| Line | Branch | Exec | Source |
|---|---|---|---|
| 1 | #include "Algorithm.h" | ||
| 2 | |||
| 3 | namespace iguana::clas12 { | ||
| 4 | |||
| 5 | REGISTER_IGUANA_ALGORITHM(PhotonGBTFilter); | ||
| 6 | |||
| 7 | 2 | void PhotonGBTFilter::Start(hipo::banklist& banks) | |
| 8 | { | ||
| 9 | 2 | ParseYAMLConfig(); | |
| 10 | |||
| 11 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
2 | b_particle = GetBankIndex(banks, "REC::Particle"); |
| 12 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
2 | b_calorimeter = GetBankIndex(banks, "REC::Calorimeter"); |
| 13 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
2 | b_config = GetBankIndex(banks, "RUN::config"); |
| 14 | |||
| 15 |
2/4✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
|
2 | o_pass = GetOptionScalar<int>("pass",{}); |
| 16 |
2/4✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
|
2 | o_threshold = GetOptionScalar<double>("threshold",{}); |
| 17 | 2 | } | |
| 18 | |||
| 19 | |||
| 20 | |||
| 21 | 2000 | void PhotonGBTFilter::Run(hipo::banklist& banks) const | |
| 22 | { | ||
| 23 | |||
| 24 |
1/2✓ Branch 0 taken 2000 times.
✗ Branch 1 not taken.
|
2000 | auto& particleBank = GetBank(banks, b_particle, "REC::Particle"); |
| 25 |
1/2✓ Branch 0 taken 2000 times.
✗ Branch 1 not taken.
|
2000 | auto& caloBank = GetBank(banks, b_calorimeter, "REC::Calorimeter"); |
| 26 |
1/2✓ Branch 0 taken 2000 times.
✗ Branch 1 not taken.
|
2000 | auto& configBank = GetBank(banks,b_config,"RUN::config"); |
| 27 | 2000 | int runnum = configBank.getInt("run",0); | |
| 28 | |||
| 29 | // Get CaloMap for the event | ||
| 30 | 2000 | auto calo_map = GetCaloMap(caloBank); | |
| 31 | |||
| 32 | // dump the bank | ||
| 33 |
2/4✓ Branch 0 taken 2000 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2000 times.
✗ Branch 3 not taken.
|
2000 | ShowBank(particleBank, Logger::Header("INPUT PARTICLES")); |
| 34 | |||
| 35 | // Loop over each photon in the particleBank to classify it | ||
| 36 | // Here we loop over the particleBank RowList | ||
| 37 | // This ensures we are only concerned with filtering photons that passed upstream filters | ||
| 38 |
4/8✓ Branch 0 taken 2000 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2000 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 2000 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2000 times.
✗ Branch 7 not taken.
|
4000 | particleBank.getMutableRowList().filter([this, &caloBank, &calo_map, runnum](auto bank, auto row) { |
| 39 | 14286 | auto pid = bank.getInt("pid", row); | |
| 40 |
2/2✓ Branch 0 taken 3730 times.
✓ Branch 1 taken 10556 times.
|
14286 | if (pid != 22) return true; |
| 41 |
1/2✓ Branch 0 taken 3730 times.
✗ Branch 1 not taken.
|
7460 | return Filter(bank, caloBank, calo_map, row, runnum); |
| 42 | }); | ||
| 43 | |||
| 44 | // dump the modified bank | ||
| 45 |
2/4✓ Branch 0 taken 2000 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2000 times.
✗ Branch 3 not taken.
|
4000 | ShowBank(particleBank, Logger::Header("OUTPUT PARTICLES")); |
| 46 | |||
| 47 | 2000 | } | |
| 48 | |||
| 49 | 9782 | bool PhotonGBTFilter::PidPurityPhotonFilter(float const E, float const Epcal, float const theta) const | |
| 50 | { | ||
| 51 | // Apply standard pid cuts on the photon, compatible to how the models were trained | ||
| 52 | // 1. Minimum photon energy cut of 200 MeV | ||
| 53 | // 2. Photon must have deposited energy in the PCal | ||
| 54 | // 3. Photon must be in the Forward Detector | ||
| 55 |
6/6✓ Branch 0 taken 8332 times.
✓ Branch 1 taken 1450 times.
✓ Branch 2 taken 7992 times.
✓ Branch 3 taken 340 times.
✓ Branch 4 taken 22 times.
✓ Branch 5 taken 7970 times.
|
9782 | if( E<0.2 || Epcal <= 0 || !(ForwardDetectorFilter(theta))) return false; |
| 56 | return true; | ||
| 57 | } | ||
| 58 | |||
| 59 | 21899 | bool PhotonGBTFilter::ForwardDetectorFilter(float const theta) const | |
| 60 | { | ||
| 61 | // Apply forward detector cut | ||
| 62 |
4/4✓ Branch 0 taken 21492 times.
✓ Branch 1 taken 407 times.
✓ Branch 2 taken 150 times.
✓ Branch 3 taken 21342 times.
|
21899 | if( theta*180/M_PI<5 || theta*180/M_PI>35 ) return false; |
| 63 | return true; | ||
| 64 | } | ||
| 65 | |||
| 66 | 3730 | bool PhotonGBTFilter::Filter(hipo::bank const &particleBank, hipo::bank const &caloBank, std::map<int, PhotonGBTFilter::calo_row_data> calo_map, int const row, int const runnum) const | |
| 67 | { | ||
| 68 | |||
| 69 | // Set variables native to the photon we are classifying | ||
| 70 | 3730 | double gPx = particleBank.getFloat("px",row); | |
| 71 | 3730 | double gPy = particleBank.getFloat("py",row); | |
| 72 | 3730 | double gPz = particleBank.getFloat("pz",row); | |
| 73 | |||
| 74 | // Set ML features intrinsic to the photon of interest | ||
| 75 | 3730 | double gE = sqrt(gPx*gPx+gPy*gPy+gPz*gPz); | |
| 76 | 3730 | double gTheta = acos(gPz / gE); | |
| 77 | 3730 | double gEpcal = calo_map[row].pcal_e; | |
| 78 | 3730 | double gm2u = calo_map[row].pcal_m2u; | |
| 79 | 3730 | double gm2v = calo_map[row].pcal_m2v; | |
| 80 | |||
| 81 | // Apply PID purity cuts on the photon | ||
| 82 | // If they do not pass, then these photons are incompatible with the trained GBT model | ||
| 83 |
2/2✓ Branch 0 taken 2654 times.
✓ Branch 1 taken 1076 times.
|
3730 | if( PidPurityPhotonFilter(gE,gEpcal,gTheta)==false ) return false; |
| 84 | |||
| 85 | //Define the variables "m_g" , "m_ch" , "m_nh" | ||
| 86 | //Should not be changed because the model was trained with this specific set of inputs | ||
| 87 | int const m_g = 3; // Number of neighboring gammas | ||
| 88 | int const m_ch = 2; // Number of neighboring charged hadrons (protons, pions, kaons) | ||
| 89 | int const m_nh = 2; // Number of neighboring neutral hadrons (neutrons) | ||
| 90 | |||
| 91 | double R_e = 0; | ||
| 92 | double dE_e = 0; | ||
| 93 | |||
| 94 | double R_gamma[m_g]; // Angular distance between calo shower centers | ||
| 95 | double dE_gamma[m_g]; // Energy difference | ||
| 96 | double Epcal_gamma[m_g]; // Energy deposited in the pcal | ||
| 97 | double m2u_gamma[m_g]; // Shower shape variables | ||
| 98 | double m2v_gamma[m_g]; // Shower shape variables | ||
| 99 | |||
| 100 | double R_ch[m_ch]; // Angular distance between calo shower centers | ||
| 101 | double dE_ch[m_ch]; // Energy difference | ||
| 102 | double Epcal_ch[m_ch]; // Energy deposited in the pcal | ||
| 103 | double m2u_ch[m_ch]; // Shower shape variables | ||
| 104 | double m2v_ch[m_ch]; // Shower shape variables | ||
| 105 | |||
| 106 | double R_nh[m_nh]; // Angular distance between calo shower centers | ||
| 107 | double dE_nh[m_nh]; // Energy difference | ||
| 108 | double Epcal_nh[m_nh]; // Energy deposited in the pcal | ||
| 109 | double m2u_nh[m_nh]; // Shower shape variables | ||
| 110 | double m2v_nh[m_nh]; // Shower shape variables | ||
| 111 | |||
| 112 | double num_photons_0_1, num_photons_0_2, num_photons_0_35; | ||
| 113 | |||
| 114 | //Initialize the arrays | ||
| 115 |
2/2✓ Branch 0 taken 7962 times.
✓ Branch 1 taken 2654 times.
|
10616 | for (int i=0; i<m_g; ++i) { |
| 116 | 7962 | R_gamma[i] = 0; | |
| 117 | 7962 | dE_gamma[i] = 0; | |
| 118 | 7962 | Epcal_gamma[i] = 0; | |
| 119 | 7962 | m2u_gamma[i] = 0; | |
| 120 | 7962 | m2v_gamma[i] = 0; | |
| 121 | } | ||
| 122 |
2/2✓ Branch 0 taken 5308 times.
✓ Branch 1 taken 2654 times.
|
7962 | for (int i=0; i<m_ch; ++i){ |
| 123 | 5308 | R_ch[i] = 0; | |
| 124 | 5308 | dE_ch[i] = 0; | |
| 125 | 5308 | Epcal_ch[i] = 0; | |
| 126 | 5308 | m2u_ch[i] = 0; | |
| 127 | 5308 | m2v_ch[i] = 0; | |
| 128 | } | ||
| 129 |
2/2✓ Branch 0 taken 5308 times.
✓ Branch 1 taken 2654 times.
|
7962 | for (int i=0; i<m_nh; ++i){ |
| 130 | 5308 | R_nh[i] = 0; | |
| 131 | 5308 | dE_nh[i] = 0; | |
| 132 | 5308 | Epcal_nh[i] = 0; | |
| 133 | 5308 | m2u_nh[i] = 0; | |
| 134 | 5308 | m2v_nh[i] = 0; | |
| 135 | } | ||
| 136 | |||
| 137 | //Set the number of photons within R<0.1, R<0.2, R<0.35 | ||
| 138 | num_photons_0_1 = 0; | ||
| 139 | num_photons_0_2 = 0; | ||
| 140 | num_photons_0_35 = 0; | ||
| 141 | |||
| 142 | |||
| 143 | // Get 3-vector that points to the photon of interest's location in the calorimeter | ||
| 144 | 2654 | auto calo_POI = calo_map.at(row); | |
| 145 | 2654 | ROOT::Math::XYZVector vPOI = GetParticleCaloVector(calo_POI); | |
| 146 | |||
| 147 | |||
| 148 | // Build nearest neighbor event structure | ||
| 149 | // Loop over particles in the event | ||
| 150 | // Here we loop over particleBank.getRows(), which ignores upstream filters | ||
| 151 | // This is critical as the GBTs were trained on identifying nearest neighbors for the whole REC::Particle bank | ||
| 152 | // Only considering nearest neighbor particles that pass upstream filters would call the accuracy of the model into question | ||
| 153 |
2/2✓ Branch 0 taken 23966 times.
✓ Branch 1 taken 2654 times.
|
26620 | for(int inner_row = 0; inner_row < particleBank.getRows(); inner_row++) { |
| 154 | // Skip over the particle if it is photon we are trying to classify | ||
| 155 |
2/2✓ Branch 0 taken 2654 times.
✓ Branch 1 taken 21312 times.
|
34172 | if (inner_row == row) continue; |
| 156 | |||
| 157 | // Check if the key exists in the calo_map | ||
| 158 | // This skips over REC::Particle entries without a REC::Calorimeter entry | ||
| 159 |
2/2✓ Branch 0 taken 9160 times.
✓ Branch 1 taken 12152 times.
|
21312 | if (calo_map.find(inner_row) == calo_map.end()) continue; |
| 160 | 12152 | auto calo_PART = calo_map.at(inner_row); | |
| 161 | |||
| 162 | 12152 | auto pid = particleBank.getInt("pid",inner_row); | |
| 163 | 12152 | auto mass = GetMass(pid); | |
| 164 | |||
| 165 | // Skip over particle if its mass was undefined | ||
| 166 |
2/2✓ Branch 0 taken 110 times.
✓ Branch 1 taken 12042 times.
|
12152 | if (mass == -1.0) continue; |
| 167 | 12042 | auto px = particleBank.getFloat("px",inner_row); | |
| 168 | 12042 | auto py = particleBank.getFloat("py",inner_row); | |
| 169 | 12042 | auto pz = particleBank.getFloat("pz",inner_row); | |
| 170 | 12042 | auto p = sqrt(px*px+py*py+pz*pz); | |
| 171 | 12042 | auto E = sqrt(p*p+mass*mass); | |
| 172 | 12042 | auto th = acos(pz/p); | |
| 173 | // Skip over particle if it is not in the forward detector (necessary for model compatibility) | ||
| 174 |
2/2✓ Branch 0 taken 200 times.
✓ Branch 1 taken 11842 times.
|
12042 | if (ForwardDetectorFilter(th)==false) continue; |
| 175 | |||
| 176 | // Get 3-vector that points to the neighboring particle's location in the calorimeter | ||
| 177 | 11842 | ROOT::Math::XYZVector vPART = GetParticleCaloVector(calo_PART); | |
| 178 | |||
| 179 | // Get angular distance between photon of interest and particle | ||
| 180 | double R = ROOT::Math::VectorUtil::Angle(vPOI, vPART); | ||
| 181 | |||
| 182 | // Logic for filling nearest neighbor variables | ||
| 183 |
2/2✓ Branch 0 taken 6052 times.
✓ Branch 1 taken 5790 times.
|
11842 | if(pid==22) {//photon |
| 184 | |||
| 185 | // Apply Photon Purity Cuts to ensure this neighbor can be used in classification | ||
| 186 |
2/2✓ Branch 0 taken 736 times.
✓ Branch 1 taken 5316 times.
|
6052 | if ( PidPurityPhotonFilter(E,calo_PART.pcal_e,th)==false ) continue; |
| 187 | |||
| 188 |
2/2✓ Branch 0 taken 836 times.
✓ Branch 1 taken 4480 times.
|
5316 | if (R < 0.1) num_photons_0_1++; |
| 189 |
2/2✓ Branch 0 taken 536 times.
✓ Branch 1 taken 3944 times.
|
5316 | if (R < 0.2) num_photons_0_2++; |
| 190 |
2/2✓ Branch 0 taken 1812 times.
✓ Branch 1 taken 2968 times.
|
5316 | if (R < 0.35) num_photons_0_35++; |
| 191 | |||
| 192 |
2/2✓ Branch 0 taken 8216 times.
✓ Branch 1 taken 380 times.
|
8596 | for (int i=0; i<m_g; ++i) { |
| 193 |
4/4✓ Branch 0 taken 1698 times.
✓ Branch 1 taken 6518 times.
✓ Branch 2 taken 3238 times.
✓ Branch 3 taken 3280 times.
|
8216 | if (R < R_gamma[i] || R_gamma[i] == 0) { |
| 194 | int j = m_g - 1; | ||
| 195 |
2/2✓ Branch 0 taken 7732 times.
✓ Branch 1 taken 4936 times.
|
12668 | while (j > i) { |
| 196 | 7732 | R_gamma[j] = R_gamma[j - 1]; | |
| 197 | 7732 | dE_gamma[j] = dE_gamma[j - 1]; | |
| 198 | 7732 | Epcal_gamma[j] = Epcal_gamma[j - 1]; | |
| 199 | 7732 | m2u_gamma[j] = m2u_gamma[j - 1]; | |
| 200 | 7732 | m2v_gamma[j] = m2v_gamma[j - 1]; | |
| 201 | j--; | ||
| 202 | } | ||
| 203 | 4936 | R_gamma[i] = R; | |
| 204 | 4936 | dE_gamma[i] = gE - E; | |
| 205 | 4936 | Epcal_gamma[i] = calo_PART.pcal_e; | |
| 206 | 4936 | m2u_gamma[i] = calo_PART.pcal_m2u; | |
| 207 | 4936 | m2v_gamma[i] = calo_PART.pcal_m2v; | |
| 208 | 4936 | break; | |
| 209 | } | ||
| 210 | } | ||
| 211 | } | ||
| 212 |
2/2✓ Branch 0 taken 278 times.
✓ Branch 1 taken 5512 times.
|
5790 | else if(pid==11){//electron |
| 213 |
4/4✓ Branch 0 taken 276 times.
✓ Branch 1 taken 2 times.
✓ Branch 2 taken 4 times.
✓ Branch 3 taken 272 times.
|
278 | if(R<R_e || R_e==0){ |
| 214 | R_e = R; | ||
| 215 | 274 | dE_e = gE - E; | |
| 216 | } | ||
| 217 | } | ||
| 218 |
6/6✓ Branch 0 taken 5308 times.
✓ Branch 1 taken 204 times.
✓ Branch 2 taken 2190 times.
✓ Branch 3 taken 3118 times.
✓ Branch 4 taken 28 times.
✓ Branch 5 taken 3090 times.
|
5512 | else if(pid==2212||pid==-2212||pid==211||pid==-211||pid==321||pid==-321){//charged hadron |
| 219 |
2/2✓ Branch 0 taken 2746 times.
✓ Branch 1 taken 22 times.
|
2768 | for (int i=0; i<m_ch; ++i) { |
| 220 |
4/4✓ Branch 0 taken 382 times.
✓ Branch 1 taken 2364 times.
✓ Branch 2 taken 2018 times.
✓ Branch 3 taken 346 times.
|
2746 | if (R < R_ch[i] || R_ch[i] == 0) { |
| 221 | int j = m_ch - 1; | ||
| 222 |
2/2✓ Branch 0 taken 2098 times.
✓ Branch 1 taken 2400 times.
|
4498 | while (j > i) { |
| 223 | 2098 | R_ch[j] = R_ch[j - 1]; | |
| 224 | 2098 | dE_ch[j] = dE_ch[j - 1]; | |
| 225 | 2098 | Epcal_ch[j] = Epcal_ch[j - 1]; | |
| 226 | 2098 | m2u_ch[j] = m2u_ch[j - 1]; | |
| 227 | 2098 | m2v_ch[j] = m2v_ch[j - 1]; | |
| 228 | j--; | ||
| 229 | } | ||
| 230 | 2400 | R_ch[i] = R; | |
| 231 | 2400 | dE_ch[i] = gE- E; | |
| 232 | 2400 | Epcal_ch[i] = calo_PART.pcal_e; | |
| 233 | 2400 | m2u_ch[i] = calo_PART.pcal_m2u; | |
| 234 | 2400 | m2v_ch[i] = calo_PART.pcal_m2v; | |
| 235 | 2400 | break; | |
| 236 | } | ||
| 237 | } | ||
| 238 | } | ||
| 239 |
1/2✓ Branch 0 taken 3090 times.
✗ Branch 1 not taken.
|
3090 | else if(pid==2112||pid==-2112){//neutral hadron |
| 240 |
2/2✓ Branch 0 taken 4064 times.
✓ Branch 1 taken 316 times.
|
4380 | for (int i=0; i<m_nh; ++i) { |
| 241 |
4/4✓ Branch 0 taken 814 times.
✓ Branch 1 taken 3250 times.
✓ Branch 2 taken 1960 times.
✓ Branch 3 taken 1290 times.
|
4064 | if (R < R_nh[i] || R_nh[i] == 0) { |
| 242 | int j = m_nh - 1; | ||
| 243 |
2/2✓ Branch 0 taken 2116 times.
✓ Branch 1 taken 2774 times.
|
4890 | while (j > i) { |
| 244 | 2116 | R_nh[j] = R_nh[j - 1]; | |
| 245 | 2116 | dE_nh[j] = dE_nh[j - 1]; | |
| 246 | 2116 | Epcal_nh[j] = Epcal_nh[j - 1]; | |
| 247 | 2116 | m2u_nh[j] = m2u_nh[j - 1]; | |
| 248 | 2116 | m2v_nh[j] = m2v_nh[j - 1]; | |
| 249 | j--; | ||
| 250 | } | ||
| 251 | 2774 | R_nh[i] = R; | |
| 252 | 2774 | dE_nh[i] = gE - E; | |
| 253 | 2774 | Epcal_nh[i] = calo_PART.pcal_e; | |
| 254 | 2774 | m2u_nh[i] = calo_PART.pcal_m2u; | |
| 255 | 2774 | m2v_nh[i] = calo_PART.pcal_m2v; | |
| 256 | 2774 | break; | |
| 257 | } | ||
| 258 | } | ||
| 259 | } | ||
| 260 | else{//unrecognized OR uncompatible particle type for the trained model | ||
| 261 | ✗ | continue; | |
| 262 | } | ||
| 263 | } | ||
| 264 | |||
| 265 | // Create and populate input_data vector for the ML model | ||
| 266 | std::vector<float> input_data = { | ||
| 267 | static_cast<float>(gE), static_cast<float>(gEpcal), static_cast<float>(gTheta), | ||
| 268 | 2654 | static_cast<float>(gm2u), static_cast<float>(gm2v), static_cast<float>(R_e), | |
| 269 | 2654 | static_cast<float>(dE_e) | |
| 270 | 2654 | }; | |
| 271 | |||
| 272 | |||
| 273 |
3/4✓ Branch 0 taken 7962 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 7962 times.
✓ Branch 3 taken 2654 times.
|
10616 | for (int i = 0; i < m_g; ++i) input_data.push_back(static_cast<float>(R_gamma[i])); |
| 274 |
3/4✓ Branch 0 taken 7962 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 7962 times.
✓ Branch 3 taken 2654 times.
|
10616 | for (int i = 0; i < m_g; ++i) input_data.push_back(static_cast<float>(dE_gamma[i])); |
| 275 |
3/4✓ Branch 0 taken 7962 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 7962 times.
✓ Branch 3 taken 2654 times.
|
10616 | for (int i = 0; i < m_g; ++i) input_data.push_back(static_cast<float>(Epcal_gamma[i])); |
| 276 |
3/4✓ Branch 0 taken 7962 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 7962 times.
✓ Branch 3 taken 2654 times.
|
10616 | for (int i = 0; i < m_g; ++i) input_data.push_back(static_cast<float>(m2u_gamma[i])); |
| 277 |
3/4✓ Branch 0 taken 7962 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 7962 times.
✓ Branch 3 taken 2654 times.
|
10616 | for (int i = 0; i < m_g; ++i) input_data.push_back(static_cast<float>(m2v_gamma[i])); |
| 278 | |||
| 279 | |||
| 280 | |||
| 281 |
3/4✓ Branch 0 taken 5308 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 5308 times.
✓ Branch 3 taken 2654 times.
|
7962 | for (int i = 0; i < m_ch; ++i) input_data.push_back(static_cast<float>(R_ch[i])); |
| 282 |
3/4✓ Branch 0 taken 5308 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 5308 times.
✓ Branch 3 taken 2654 times.
|
7962 | for (int i = 0; i < m_ch; ++i) input_data.push_back(static_cast<float>(dE_ch[i])); |
| 283 |
3/4✓ Branch 0 taken 5308 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 5308 times.
✓ Branch 3 taken 2654 times.
|
7962 | for (int i = 0; i < m_ch; ++i) input_data.push_back(static_cast<float>(Epcal_ch[i])); |
| 284 |
3/4✓ Branch 0 taken 5308 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 5308 times.
✓ Branch 3 taken 2654 times.
|
7962 | for (int i = 0; i < m_ch; ++i) input_data.push_back(static_cast<float>(m2u_ch[i])); |
| 285 |
3/4✓ Branch 0 taken 5308 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 5308 times.
✓ Branch 3 taken 2654 times.
|
7962 | for (int i = 0; i < m_ch; ++i) input_data.push_back(static_cast<float>(m2v_ch[i])); |
| 286 | |||
| 287 | |||
| 288 | |||
| 289 |
3/4✓ Branch 0 taken 5308 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 5308 times.
✓ Branch 3 taken 2654 times.
|
7962 | for (int i = 0; i < m_nh; ++i) input_data.push_back(static_cast<float>(R_nh[i])); |
| 290 |
3/4✓ Branch 0 taken 5308 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 5308 times.
✓ Branch 3 taken 2654 times.
|
7962 | for (int i = 0; i < m_nh; ++i) input_data.push_back(static_cast<float>(dE_nh[i])); |
| 291 |
3/4✓ Branch 0 taken 5308 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 5308 times.
✓ Branch 3 taken 2654 times.
|
7962 | for (int i = 0; i < m_nh; ++i) input_data.push_back(static_cast<float>(Epcal_nh[i])); |
| 292 |
3/4✓ Branch 0 taken 5308 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 5308 times.
✓ Branch 3 taken 2654 times.
|
7962 | for (int i = 0; i < m_nh; ++i) input_data.push_back(static_cast<float>(m2u_nh[i])); |
| 293 |
3/4✓ Branch 0 taken 5308 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 5308 times.
✓ Branch 3 taken 2654 times.
|
7962 | for (int i = 0; i < m_nh; ++i) input_data.push_back(static_cast<float>(m2v_nh[i])); |
| 294 | |||
| 295 | |||
| 296 |
1/2✓ Branch 0 taken 2654 times.
✗ Branch 1 not taken.
|
2654 | input_data.push_back(static_cast<float>(num_photons_0_1)); |
| 297 |
1/2✓ Branch 0 taken 2654 times.
✗ Branch 1 not taken.
|
2654 | input_data.push_back(static_cast<float>(num_photons_0_2)); |
| 298 |
1/4✓ Branch 0 taken 2654 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
2654 | input_data.push_back(static_cast<float>(num_photons_0_35)); |
| 299 | |||
| 300 |
1/2✓ Branch 0 taken 2654 times.
✗ Branch 1 not taken.
|
2654 | return ClassifyPhoton(input_data, runnum); |
| 301 | } | ||
| 302 | |||
| 303 | 2654 | bool PhotonGBTFilter::ClassifyPhoton(std::vector<float> const &input_data, int const runnum) const { | |
| 304 | 2654 | auto modelFunction = getModelFunction(runnum); | |
| 305 | double sigmoid_x = modelFunction(input_data); | ||
| 306 | 2654 | double prediction = 1 / (1 + exp(-sigmoid_x)); | |
| 307 |
1/2✓ Branch 0 taken 2654 times.
✗ Branch 1 not taken.
|
5308 | return (prediction > o_threshold); |
| 308 | } | ||
| 309 | |||
| 310 | 2000 | std::map<int, PhotonGBTFilter::calo_row_data> PhotonGBTFilter::GetCaloMap(hipo::bank const& bank) const | |
| 311 | { | ||
| 312 | std::map<int, PhotonGBTFilter::calo_row_data> calo_map; | ||
| 313 | // Loop over REC::Calorimeter rows | ||
| 314 | // Here we use bank.getRows() to purposefully ignore upstream filters | ||
| 315 |
2/2✓ Branch 0 taken 12212 times.
✓ Branch 1 taken 2000 times.
|
14212 | for(int row = 0; row < bank.getRows(); row++){ |
| 316 | 12212 | auto pindex = bank.getShort("pindex",row); | |
| 317 | 12212 | auto x = bank.getFloat("x",row); | |
| 318 | 12212 | auto y = bank.getFloat("y",row); | |
| 319 | 12212 | auto z = bank.getFloat("z",row); | |
| 320 | 12212 | auto m2u = bank.getFloat("m2u",row); | |
| 321 | 12212 | auto m2v = bank.getFloat("m2v",row); | |
| 322 | 12212 | auto layer = bank.getInt("layer",row); | |
| 323 | 12212 | auto e = bank.getFloat("energy",row); | |
| 324 | |||
| 325 | // Ensure an entry exists in the map for the given pindex | ||
| 326 |
2/2✓ Branch 0 taken 7598 times.
✓ Branch 1 taken 4614 times.
|
12212 | if (calo_map.find(pindex) == calo_map.end()) { |
| 327 |
1/2✓ Branch 0 taken 7598 times.
✗ Branch 1 not taken.
|
7598 | calo_map[pindex] = PhotonGBTFilter::calo_row_data(); |
| 328 | } | ||
| 329 | |||
| 330 |
3/4✓ Branch 0 taken 6204 times.
✓ Branch 1 taken 3546 times.
✓ Branch 2 taken 2462 times.
✗ Branch 3 not taken.
|
12212 | switch(layer){ |
| 331 | 6204 | case 1: // pcal | |
| 332 |
1/2✓ Branch 0 taken 6204 times.
✗ Branch 1 not taken.
|
6204 | calo_map[pindex].pcal_x = x; |
| 333 |
1/2✓ Branch 0 taken 6204 times.
✗ Branch 1 not taken.
|
6204 | calo_map[pindex].pcal_y = y; |
| 334 |
1/2✓ Branch 0 taken 6204 times.
✗ Branch 1 not taken.
|
6204 | calo_map[pindex].pcal_z = z; |
| 335 |
1/2✓ Branch 0 taken 6204 times.
✗ Branch 1 not taken.
|
6204 | calo_map[pindex].pcal_e = e; |
| 336 |
1/2✓ Branch 0 taken 6204 times.
✗ Branch 1 not taken.
|
6204 | calo_map[pindex].pcal_m2u = m2u; |
| 337 |
1/2✓ Branch 0 taken 6204 times.
✗ Branch 1 not taken.
|
6204 | calo_map[pindex].pcal_m2v = m2v; |
| 338 | 6204 | break; | |
| 339 | 3546 | case 4: // ecin | |
| 340 |
1/2✓ Branch 0 taken 3546 times.
✗ Branch 1 not taken.
|
3546 | calo_map[pindex].ecin_x = x; |
| 341 |
1/2✓ Branch 0 taken 3546 times.
✗ Branch 1 not taken.
|
3546 | calo_map[pindex].ecin_y = y; |
| 342 |
1/2✓ Branch 0 taken 3546 times.
✗ Branch 1 not taken.
|
3546 | calo_map[pindex].ecin_z = z; |
| 343 | 3546 | break; | |
| 344 | 2462 | case 7: // ecout | |
| 345 |
1/2✓ Branch 0 taken 2462 times.
✗ Branch 1 not taken.
|
2462 | calo_map[pindex].ecout_x = x; |
| 346 |
1/2✓ Branch 0 taken 2462 times.
✗ Branch 1 not taken.
|
2462 | calo_map[pindex].ecout_y = y; |
| 347 |
1/2✓ Branch 0 taken 2462 times.
✗ Branch 1 not taken.
|
2462 | calo_map[pindex].ecout_z = z; |
| 348 | 2462 | break; | |
| 349 | } | ||
| 350 | } | ||
| 351 | 2000 | return calo_map; | |
| 352 | } | ||
| 353 | |||
| 354 |
2/2✓ Branch 0 taken 2152 times.
✓ Branch 1 taken 12344 times.
|
14496 | ROOT::Math::XYZVector PhotonGBTFilter::GetParticleCaloVector(PhotonGBTFilter::calo_row_data calo_row) const { |
| 355 | // Determine the 3-vector location of where the photon of interest's calo deposition is | ||
| 356 | // First we check the pcal coords, then ecin, then ecout | ||
| 357 | ROOT::Math::XYZVector v; | ||
| 358 |
2/2✓ Branch 0 taken 2152 times.
✓ Branch 1 taken 12344 times.
|
14496 | if (calo_row.pcal_x == 0) { |
| 359 |
2/2✓ Branch 0 taken 842 times.
✓ Branch 1 taken 1310 times.
|
2152 | if (calo_row.ecin_x == 0) { |
| 360 | 842 | v.SetXYZ(calo_row.ecout_x, calo_row.ecout_y, calo_row.ecout_z); | |
| 361 | } else { | ||
| 362 | 1310 | v.SetXYZ(calo_row.ecin_x, calo_row.ecin_y, calo_row.ecin_z); | |
| 363 | } | ||
| 364 | } else { | ||
| 365 | 12344 | v.SetXYZ(calo_row.pcal_x, calo_row.pcal_y, calo_row.pcal_z); | |
| 366 | } | ||
| 367 | 14496 | return v; | |
| 368 | } | ||
| 369 | |||
| 370 | 12152 | double PhotonGBTFilter::GetMass(int pid) const { | |
| 371 |
2/2✓ Branch 0 taken 12042 times.
✓ Branch 1 taken 110 times.
|
12152 | auto it = particle::mass.find(static_cast<particle::PDG>(pid)); |
| 372 |
2/2✓ Branch 0 taken 12042 times.
✓ Branch 1 taken 110 times.
|
12152 | if (it != particle::mass.end()) { |
| 373 | 12042 | return it->second; | |
| 374 | } else { | ||
| 375 | return -1.0; // Default mass if pid not found | ||
| 376 | } | ||
| 377 | } | ||
| 378 | |||
| 379 | 2654 | std::function<double(std::vector<float> const &)> PhotonGBTFilter::getModelFunction(int runnum) const { | |
| 380 | |||
| 381 |
1/2✓ Branch 0 taken 18578 times.
✗ Branch 1 not taken.
|
18578 | for (const auto &entry : modelMap) { |
| 382 |
4/6✓ Branch 0 taken 18578 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2654 times.
✓ Branch 3 taken 15924 times.
✓ Branch 4 taken 2654 times.
✗ Branch 5 not taken.
|
18578 | if (runnum >= std::get<0>(entry.first) && runnum <= std::get<1>(entry.first) && o_pass == std::get<2>(entry.first)) { |
| 383 | 2654 | return entry.second; | |
| 384 | } | ||
| 385 | } | ||
| 386 | |||
| 387 | // Default to RGA inbending pass1 if no match found | ||
| 388 | ✗ | m_log->Warn("Run Number {} with pass {} has no matching PhotonGBT model...Defaulting to RGA inbending pass1...", runnum, o_pass); | |
| 389 | ✗ | return [](std::vector<float> const &data) { return ApplyCatboostModel_RGA_inbending_pass1(data); }; | |
| 390 | } | ||
| 391 | |||
| 392 | 1 | void PhotonGBTFilter::Stop() | |
| 393 | { | ||
| 394 | 1 | } | |
| 395 | |||
| 396 | } | ||
| 397 |