sPhenix code displayed by LXR |
|
|
|
|||
File indexing completed on 2025-08-06 08:13:58
0001 #pragma once 0002 #include "include_tracking/least_square2.cc" 0003 0004 class truth 0005 { 0006 public: 0007 double m_truthpx; 0008 double m_truthpy; 0009 double m_truthpz; 0010 double m_truthpt; 0011 double m_trutheta; 0012 double m_truththeta; 0013 double m_truthphi; 0014 double m_truthz_out; 0015 0016 double m_xvtx; 0017 double m_yvtx; 0018 double m_zvtx; 0019 0020 int m_truthpid; 0021 int m_status; 0022 0023 TF1 *truth_xy = nullptr; 0024 TF1 *truth_rz = nullptr; 0025 0026 double dca_mean[3]; 0027 0028 double center[2]; 0029 double m_rad; 0030 0031 bool is_intt = false; 0032 bool is_charged = false; 0033 bool have_cluster = false; 0034 bool have_track = false; 0035 0036 void set_truth(double M_truthpx, double M_truthpy, double M_truthpz, double M_truthpt, int M_status, double M_trutheta, int M_truthpid, double M_truththeta, double M_truthphi, double M_xvtx, double M_yvtx, double M_zvtx) 0037 { 0038 m_truthpx = M_truthpx; 0039 m_truthpy = M_truthpy; 0040 m_truthpz = M_truthpz; 0041 m_truthpt = M_truthpt; 0042 m_status = M_status; 0043 m_trutheta = M_trutheta; 0044 m_truthpid = M_truthpid; 0045 m_truthphi = M_truthphi; 0046 m_truththeta = M_truththeta; 0047 m_xvtx = M_xvtx; 0048 m_yvtx = M_yvtx; 0049 m_zvtx = M_zvtx; 0050 } 0051 0052 double getzout() 0053 { 0054 double y = m_truthpy / fabs(m_truthpy) * 10; 0055 // cout << "truth_rz->GetParameter(0) : " << truth_rz->GetParameter(0) << endl; 0056 // cout<<"y : "<<y<<endl; 0057 double z_outer = (y - truth_rz->GetParameter(0)) / truth_rz->GetParameter(1); 0058 // cout << "z_outer : " << z_outer << endl; 0059 0060 return z_outer; 0061 } 0062 0063 double getpointr(int mode) 0064 { 0065 double r = 0; 0066 if (mode == 1) // truth 0067 { 0068 r = m_truthpy / fabs(m_truthpy) * m_truthpt; 0069 } 0070 else if (mode == 2) // dca_mean 0071 { 0072 r = dca_mean[1] / fabs(dca_mean[1]) * sqrt(dca_mean[0] * dca_mean[0] + dca_mean[1] * dca_mean[1]); 0073 } 0074 else if (mode == 3) // vertex 0075 { 0076 r = /*m_yvtx / fabs(m_yvtx) **/ sqrt(m_xvtx * m_xvtx + m_yvtx * m_yvtx); 0077 } 0078 0079 return r; 0080 } 0081 0082 void calc_line() 0083 { 0084 if (truth_xy != nullptr) 0085 { 0086 delete truth_xy; 0087 truth_xy = nullptr; 0088 } 0089 0090 truth_xy = new TF1("truth_xy", "pol1", -1000, 1000); 0091 0092 truth_xy->SetParameter(1, m_truthpy / m_truthpx); 0093 truth_xy->SetParameter(0, -m_yvtx / m_truthpx * m_xvtx + m_yvtx); 0094 0095 if (truth_rz != nullptr) 0096 { 0097 delete truth_rz; 0098 truth_rz = nullptr; 0099 } 0100 0101 truth_rz = new TF1("truth_rz", "pol1", -1000, 1000); 0102 0103 // from the truth vertex 0104 double slope = getpointr(1) / m_truthpz; 0105 double intercept = -slope * m_zvtx + getpointr(3); 0106 // cout << slope << " " << m_zvtx << " " << getpointr(3) << endl; 0107 0108 truth_rz->SetParameter(1, slope); 0109 truth_rz->SetParameter(0, intercept); 0110 0111 m_truthz_out = getzout(); 0112 } 0113 0114 void calc_center() 0115 { 0116 m_rad = m_truthpt / (0.3 * 1.4) * 100; 0117 double sign = m_truthpid / fabs(m_truthpid); 0118 center[0] = -sign * (m_xvtx + m_rad * (-sin(m_truthphi))); 0119 center[1] = -sign * (m_yvtx + m_rad * cos(m_truthphi)); 0120 } 0121 0122 ///////////////////////////////////////// 0123 0124 // vector<double> get_crossline(double x1, double y1, double r1) 0125 // { 0126 // vector<double> cross_co; 0127 // // ax + by + c = 0 0128 // cross_co.push_back(2 * (cx - x1)); // a 0129 // cross_co.push_back(2 * (cy - y1)); // b 0130 // cross_co.push_back((x1 + cx) * (x1 - cx) + (y1 + cy) * (y1 - cy) + (rad + r1) * (rad - r1)); // c 0131 0132 // return cross_co; 0133 // } 0134 0135 // vector<double> get_crosspoint(const vector<double> &V) // ax + by + c = 0 0136 // { 0137 // vector<double> cross; 0138 // double a = V[0]; 0139 // double b = V[1]; 0140 // double c = V[2]; 0141 0142 // double d = fabs(a * cx + b * cy + c) / sqrt(a * a + b * b); 0143 // cout << "d ; " << d << " " << rad << endl; 0144 0145 // double theta = atan2(b, a); 0146 0147 // if (d > rad) 0148 // { 0149 // cout << "d > rad" << endl; 0150 // } 0151 // else if (d == rad) 0152 // { 0153 // cout << "d == rad" << endl; 0154 0155 // if (a * cx + b * cy + c > 0) 0156 // theta += M_PI; 0157 // cross.push_back(rad * cos(theta) + cx); 0158 // cross.push_back(rad * sin(theta) + cy); 0159 // } 0160 // else 0161 // { 0162 // cout << "d < rad" << endl; 0163 // double alpha, beta, phi; 0164 // phi = acos(d / rad); 0165 // alpha = theta - phi; 0166 // beta = theta + phi; 0167 // if (a * cx + b * cy + c > 0) 0168 // { 0169 // alpha += M_PI; 0170 // beta += M_PI; 0171 // } 0172 // // double temp_cross[2][2]; 0173 // vector<vector<double>> temp_cross = {{rad * cos(alpha) + cx, rad * cos(beta) + cx}, {rad * sin(alpha) + cy, rad * sin(beta) + cy}}; 0174 // // for (unsigned int j = 0; j < temp_cross.at(0).size(); j++) 0175 // // { 0176 // // cout << "temp_cross : "; 0177 // // for (unsigned int i = 0; i < temp_cross.size(); i++) 0178 // // cout << temp_cross.at(j).at(i) << " "; 0179 // // cout << endl; 0180 // // } 0181 // cross = get_closestpoint(temp_cross); 0182 // } 0183 // // cout << "cross size : " << cross.size() << endl; 0184 // // for (int i = 0; i < cross.size(); i++) 0185 // // cout << cross[i] << endl; 0186 // return cross; 0187 // } 0188 0189 // vector<double> get_closestpoint(const vector<vector<double>> VV) 0190 // { 0191 // vector<double> closest; 0192 // double pre_phi = 999; 0193 // double phi_p1 = atan2(p1.y() - dca_mean[1], p1.x() - dca_mean[0]); 0194 // for (unsigned int i = 0; i < VV.at(0).size(); i++) 0195 // { 0196 // cout << "dca_mean : " << dca_mean[0] << " " << dca_mean[1] << endl; 0197 // cout << "VV : " << VV[0][i] << " " << VV[1][i] << endl; 0198 // cout << "pVV : " << VV[0][i] - dca_mean[0] << " " << VV[1][i] - dca_mean[1] << endl; 0199 // cout << "p1 : " << p1.x() - dca_mean[0] << " " << p1.y() - dca_mean[1] << endl; 0200 // // cout << VV.at(i).at(0) << endl; 0201 // double dot = (VV[0][i] - dca_mean[0]) * (p1.x() - dca_mean[0]) + (VV[1][i] - dca_mean[1]) * (p1.y() - dca_mean[1]); 0202 // double temp_phi = atan2(VV[1][i] - dca_mean[1], VV[0][i] - dca_mean[0]); 0203 // cout << "dot : " << dot << endl; 0204 // if (dot > 0) 0205 // { 0206 // if (fabs(temp_phi - phi_p1) < fabs(pre_phi - phi_p1)) 0207 // { 0208 // closest.erase(closest.begin(), closest.end()); 0209 // cout << "VV size : " << VV.size() << " " << VV.at(0).size() << endl; 0210 // for (int j = 0; j < 2; j++) 0211 // closest.push_back(VV[j][i]); 0212 // pre_phi = temp_phi; 0213 // } 0214 // } 0215 // } 0216 // // closest = temp_closest; 0217 // // for (unsigned int i = 0; i < VV.at(0).size(); i++) 0218 // // { 0219 // // for (int j = 0; j < 2; j++) 0220 // // closest.push_back(VV[j][i]); 0221 // // } 0222 // // for (unsigned int i = 0; i < closest.size(); i++) 0223 // // cout << closest[i] << endl; 0224 0225 // return closest; 0226 // } 0227 0228 // vector<double> get_crossframe(int Detect) 0229 // { 0230 // cout << endl; 0231 // cout << "calc frame" << endl; 0232 // double size; 0233 // if (Detect == 0) 0234 // size = 15; 0235 // else if (Detect == 1) 0236 // size = 120; 0237 // vector<vector<double>> crossframe(2, vector<double>(0)); 0238 // for (int j = 0; j < 2; j++) 0239 // { 0240 // for (int i = 0; i < 2; i++) 0241 // { 0242 // double a_ = (j == 0) ? 1 : 0; 0243 // double b_ = (j == 0) ? 0 : 1; 0244 // double c_ = (2 * i - 1) * size; 0245 // cout << "kesu : " << a_ << " " << b_ << " " << c_ << endl; 0246 // vector<double> co_ = {a_, b_, c_}; 0247 // vector<double> temp_crossframe = get_crosspoint(co_); 0248 // if (temp_crossframe.size() == 2) 0249 // { 0250 // cout << "temp_crossframe : " << temp_crossframe[0] << " " << temp_crossframe[1] << endl; 0251 // crossframe.at(0).push_back(temp_crossframe[0]); 0252 // crossframe.at(1).push_back(temp_crossframe[1]); 0253 // } 0254 // } 0255 // } 0256 // cout << "crossframe" << endl; 0257 // cout << "size : " << crossframe.size() << " " << crossframe.at(0).size() << endl; 0258 // for (int j = 0; j < crossframe.at(0).size(); j++) 0259 // { 0260 // for (int i = 0; i < crossframe.size(); i++) 0261 // { 0262 // cout << crossframe[i][j] << ", "; 0263 // } 0264 // cout << endl; 0265 // } 0266 // vector<double> result; 0267 // if (crossframe.at(0).size() != 0) 0268 // result = get_closestpoint(crossframe); 0269 // cout << "closest frame" << endl; 0270 // for (int i = 0; i < result.size(); i++) 0271 // { 0272 // cout << result[i] << ", "; 0273 // } 0274 // cout << endl; 0275 // cout << "calc end" << endl; 0276 0277 // return result; 0278 // } 0279 };
| [ Source navigation ] | [ Diff markup ] | [ Identifier search ] | [ general search ] |
|
This page was automatically generated by the 2.3.7 LXR engine. The LXR team |
|