GCC Code Coverage Report


Directory: ./
File: src/iguana/algorithms/clas12/PhotonGBTFilter/Algorithm.cc
Date: 2025-10-15 15:39:33
Exec Total Coverage
Lines: 194 212 91.5%
Functions: 11 11 100.0%
Branches: 184 248 74.2%

Line Branch Exec Source
1 #include "Algorithm.h"
2 #include "iguana/algorithms/TypeDefs.h"
3 #include "models/RGA_inbending_pass1.cpp"
4 #include "models/RGA_outbending_pass1.cpp"
5 #include "models/RGA_inbending_pass2.cpp"
6 #include "models/RGA_outbending_pass2.cpp"
7 #include "models/RGC_Summer2022_pass1.cpp"
8
9 namespace iguana::clas12 {
10
11 REGISTER_IGUANA_ALGORITHM(PhotonGBTFilter);
12
13 // Map for the GBT Models to use depending on pass and run number
14 std::map<std::tuple<int, int, int>, std::function<double(std::vector<float> const &)>> const PhotonGBTFilter::modelMap = {
15 {{5032, 5332, 1}, [](std::vector<float> const &data) { return ApplyCatboostModel_RGA_inbending_pass1(data); }}, // Fall2018 RGA Inbending
16 {{5032, 5332, 2}, [](std::vector<float> const &data) { return ApplyCatboostModel_RGA_inbending_pass2(data); }}, // Fall2018 RGA Inbending
17 {{5333, 5666, 1}, [](std::vector<float> const &data) { return ApplyCatboostModel_RGA_outbending_pass1(data); }}, // Fall2018 RGA Outbending
18 {{5333, 5666, 2}, [](std::vector<float> const &data) { return ApplyCatboostModel_RGA_outbending_pass2(data); }}, // Fall2018 RGA Outbending
19 2564 {{6616, 6783, 1}, [](std::vector<float> const &data) { return ApplyCatboostModel_RGA_inbending_pass1(data); }}, // Spring2019 RGA Inbending
20 {{6616, 6783, 2}, [](std::vector<float> const &data) { return ApplyCatboostModel_RGA_inbending_pass2(data); }}, // Spring2019 RGA Inbending
21 {{6156, 6603, 1}, [](std::vector<float> const &data) { return ApplyCatboostModel_RGA_inbending_pass1(data); }}, // Spring2019 RGB Inbending
22 {{6156, 6603, 2}, [](std::vector<float> const &data) { return ApplyCatboostModel_RGA_inbending_pass2(data); }}, // Spring2019 RGB Inbending
23 {{11093, 11283, 1}, [](std::vector<float> const &data) { return ApplyCatboostModel_RGA_outbending_pass1(data); }}, // Fall2019 RGB Outbending
24 {{11093, 11283, 2}, [](std::vector<float> const &data) { return ApplyCatboostModel_RGA_outbending_pass2(data); }}, // Fall2019 RGB Outbending
25 {{11284, 11300, 1}, [](std::vector<float> const &data) { return ApplyCatboostModel_RGA_inbending_pass1(data); }}, // Fall2019 RGB BAND Inbending
26 {{11284, 11300, 2}, [](std::vector<float> const &data) { return ApplyCatboostModel_RGA_inbending_pass2(data); }}, // Fall2019 RGB BAND Inbending
27 {{11323, 11571, 1}, [](std::vector<float> const &data) { return ApplyCatboostModel_RGA_inbending_pass1(data); }}, // Spring2020 RGB Inbending
28 {{11323, 11571, 2}, [](std::vector<float> const &data) { return ApplyCatboostModel_RGA_inbending_pass2(data); }}, // Spring2020 RGB Inbending
29 {{16042, 16772, 1}, [](std::vector<float> const &data) { return ApplyCatboostModel_RGC_Summer2022_pass1(data); }}, // Summer2022 RGC Inbending
30 {{16042, 16772, 2}, [](std::vector<float> const &data) { return ApplyCatboostModel_RGC_Summer2022_pass1(data); }} // Summer2022 RGC Inbending (no pass2 currently)
31 };
32
33 2 void PhotonGBTFilter::Start(hipo::banklist& banks)
34 {
35
36 2 ParseYAMLConfig();
37
38
2/4
✓ Branch 0 (4→5) taken 2 times.
✗ Branch 1 (4→43) not taken.
✓ Branch 2 (5→6) taken 2 times.
✗ Branch 3 (5→8) not taken.
2 b_particle = GetBankIndex(banks, "REC::Particle");
39
2/4
✓ Branch 0 (11→12) taken 2 times.
✗ Branch 1 (11→49) not taken.
✗ Branch 2 (12→13) not taken.
✓ Branch 3 (12→15) taken 2 times.
2 b_calorimeter = GetBankIndex(banks, "REC::Calorimeter");
40
2/4
✓ Branch 0 (18→19) taken 2 times.
✗ Branch 1 (18→55) not taken.
✓ Branch 2 (19→20) taken 2 times.
✗ Branch 3 (19→22) not taken.
2 b_config = GetBankIndex(banks, "RUN::config");
41
42
3/6
✓ Branch 0 (25→26) taken 2 times.
✗ Branch 1 (25→67) not taken.
✓ Branch 2 (26→27) taken 2 times.
✗ Branch 3 (26→61) not taken.
✓ Branch 4 (27→28) taken 2 times.
✗ Branch 5 (27→30) not taken.
4 o_pass = GetOptionScalar<int>("pass");
43
3/6
✓ Branch 0 (34→35) taken 2 times.
✗ Branch 1 (34→75) not taken.
✓ Branch 2 (35→36) taken 2 times.
✗ Branch 3 (35→69) not taken.
✓ Branch 4 (36→37) taken 2 times.
✗ Branch 5 (36→39) not taken.
4 o_threshold = GetOptionScalar<double>("threshold");
44 2 }
45
46
47
48 2000 void PhotonGBTFilter::Run(hipo::banklist& banks) const
49 {
50
51
1/2
✓ Branch 0 (3→4) taken 2000 times.
✗ Branch 1 (3→50) not taken.
2000 auto& particleBank = GetBank(banks, b_particle, "REC::Particle");
52
1/2
✓ Branch 0 (10→11) taken 2000 times.
✗ Branch 1 (10→56) not taken.
2000 auto& caloBank = GetBank(banks, b_calorimeter, "REC::Calorimeter");
53
1/2
✓ Branch 0 (17→18) taken 2000 times.
✗ Branch 1 (17→62) not taken.
2000 auto& configBank = GetBank(banks,b_config,"RUN::config");
54 2000 int runnum = configBank.getInt("run",0);
55
56 // Get CaloMap for the event
57 2000 auto calo_map = GetCaloMap(caloBank);
58
59 // dump the bank
60
2/4
✓ Branch 0 (25→26) taken 2000 times.
✗ Branch 1 (25→83) not taken.
✓ Branch 2 (28→29) taken 2000 times.
✗ Branch 3 (28→68) not taken.
4000 ShowBank(particleBank, Logger::Header("INPUT PARTICLES"));
61
62 // Loop over each photon in the particleBank to classify it
63 // Here we loop over the particleBank RowList
64 // This ensures we are only concerned with filtering photons that passed upstream filters
65
3/6
✓ Branch 0 (34→35) taken 2000 times.
✗ Branch 1 (34→83) not taken.
✓ Branch 2 (35→36) taken 2000 times.
✗ Branch 3 (35→83) not taken.
✓ Branch 4 (36→37) taken 2000 times.
✗ Branch 5 (36→74) not taken.
2000 particleBank.getMutableRowList().filter([this, &caloBank, &calo_map, runnum](auto bank, auto row) {
66 13986 auto pid = bank.getInt("pid", row);
67
2/2
✓ Branch 0 (3→4) taken 3678 times.
✓ Branch 1 (3→8) taken 10308 times.
13986 if (pid != 22) return true;
68
1/2
✓ Branch 0 (5→6) taken 3678 times.
✗ Branch 1 (5→9) not taken.
7356 return Filter(bank, caloBank, calo_map, row, runnum);
69 });
70
71 // dump the modified bank
72
2/4
✓ Branch 0 (39→40) taken 2000 times.
✗ Branch 1 (39→83) not taken.
✓ Branch 2 (42→43) taken 2000 times.
✗ Branch 3 (42→77) not taken.
6000 ShowBank(particleBank, Logger::Header("OUTPUT PARTICLES"));
73
74 2000 }
75
76 9240 bool PhotonGBTFilter::PidPurityPhotonFilter(float const E, float const Epcal, float const theta) const
77 {
78 // Apply standard pid cuts on the photon, compatible to how the models were trained
79 // 1. Minimum photon energy cut of 200 MeV
80 // 2. Photon must have deposited energy in the PCal
81 // 3. Photon must be in the Forward Detector
82
6/6
✓ Branch 0 (2→3) taken 7788 times.
✓ Branch 1 (2→6) taken 1452 times.
✓ Branch 2 (3→4) taken 7404 times.
✓ Branch 3 (3→6) taken 384 times.
✓ Branch 4 (5→6) taken 40 times.
✓ Branch 5 (5→7) taken 7364 times.
9240 if( E<0.2 || Epcal <= 0 || !(ForwardDetectorFilter(theta))) return false;
83 return true;
84 }
85
86 20669 bool PhotonGBTFilter::ForwardDetectorFilter(float const theta) const
87 {
88 // Apply forward detector cut
89
4/4
✓ Branch 0 (2→3) taken 20219 times.
✓ Branch 1 (2→4) taken 450 times.
✓ Branch 2 (3→4) taken 118 times.
✓ Branch 3 (3→5) taken 20101 times.
20669 if( theta*180/M_PI<5 || theta*180/M_PI>35 ) return false;
90 return true;
91 }
92
93 3678 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
94 {
95
96 // Set variables native to the photon we are classifying
97 3678 double gPx = particleBank.getFloat("px",row);
98 3678 double gPy = particleBank.getFloat("py",row);
99 3678 double gPz = particleBank.getFloat("pz",row);
100
101 // Set ML features intrinsic to the photon of interest
102 3678 double gE = sqrt(gPx*gPx+gPy*gPy+gPz*gPz);
103 3678 double gTheta = acos(gPz / gE);
104 3678 double gEpcal = calo_map[row].pcal_e;
105 3678 double gm2u = calo_map[row].pcal_m2u;
106 3678 double gm2v = calo_map[row].pcal_m2v;
107
108 // Apply PID purity cuts on the photon
109 // If they do not pass, then these photons are incompatible with the trained GBT model
110
2/2
✓ Branch 0 (9→11) taken 2564 times.
✓ Branch 1 (9→151) taken 1114 times.
3678 if( PidPurityPhotonFilter(gE,gEpcal,gTheta)==false ) return false;
111
112 //Define the variables "m_g" , "m_ch" , "m_nh"
113 //Should not be changed because the model was trained with this specific set of inputs
114 int const m_g = 3; // Number of neighboring gammas
115 int const m_ch = 2; // Number of neighboring charged hadrons (protons, pions, kaons)
116 int const m_nh = 2; // Number of neighboring neutral hadrons (neutrons)
117
118 double R_e = 0;
119 double dE_e = 0;
120
121 double R_gamma[m_g]; // Angular distance between calo shower centers
122 double dE_gamma[m_g]; // Energy difference
123 double Epcal_gamma[m_g]; // Energy deposited in the pcal
124 double m2u_gamma[m_g]; // Shower shape variables
125 double m2v_gamma[m_g]; // Shower shape variables
126
127 double R_ch[m_ch]; // Angular distance between calo shower centers
128 double dE_ch[m_ch]; // Energy difference
129 double Epcal_ch[m_ch]; // Energy deposited in the pcal
130 double m2u_ch[m_ch]; // Shower shape variables
131 double m2v_ch[m_ch]; // Shower shape variables
132
133 double R_nh[m_nh]; // Angular distance between calo shower centers
134 double dE_nh[m_nh]; // Energy difference
135 double Epcal_nh[m_nh]; // Energy deposited in the pcal
136 double m2u_nh[m_nh]; // Shower shape variables
137 double m2v_nh[m_nh]; // Shower shape variables
138
139 double num_photons_0_1, num_photons_0_2, num_photons_0_35;
140
141 //Initialize the arrays
142
2/2
✓ Branch 0 (11→10) taken 7692 times.
✓ Branch 1 (11→13) taken 2564 times.
10256 for (int i=0; i<m_g; ++i) {
143 7692 R_gamma[i] = 0;
144 7692 dE_gamma[i] = 0;
145 7692 Epcal_gamma[i] = 0;
146 7692 m2u_gamma[i] = 0;
147 7692 m2v_gamma[i] = 0;
148 }
149
2/2
✓ Branch 0 (13→12) taken 5128 times.
✓ Branch 1 (13→15) taken 2564 times.
7692 for (int i=0; i<m_ch; ++i){
150 5128 R_ch[i] = 0;
151 5128 dE_ch[i] = 0;
152 5128 Epcal_ch[i] = 0;
153 5128 m2u_ch[i] = 0;
154 5128 m2v_ch[i] = 0;
155 }
156
2/2
✓ Branch 0 (15→14) taken 5128 times.
✓ Branch 1 (15→16) taken 2564 times.
7692 for (int i=0; i<m_nh; ++i){
157 5128 R_nh[i] = 0;
158 5128 dE_nh[i] = 0;
159 5128 Epcal_nh[i] = 0;
160 5128 m2u_nh[i] = 0;
161 5128 m2v_nh[i] = 0;
162 }
163
164 //Set the number of photons within R<0.1, R<0.2, R<0.35
165 num_photons_0_1 = 0;
166 num_photons_0_2 = 0;
167 num_photons_0_35 = 0;
168
169
170 // Get 3-vector that points to the photon of interest's location in the calorimeter
171 2564 auto calo_POI = calo_map.at(row);
172 2564 ROOT::Math::XYZVector vPOI = GetParticleCaloVector(calo_POI);
173
174
175 // Build nearest neighbor event structure
176 // Loop over particles in the event
177 // Here we loop over particleBank.getRows(), which ignores upstream filters
178 // This is critical as the GBTs were trained on identifying nearest neighbors for the whole REC::Particle bank
179 // Only considering nearest neighbor particles that pass upstream filters would call the accuracy of the model into question
180
2/2
✓ Branch 0 (95→19) taken 22304 times.
✓ Branch 1 (95→96) taken 2564 times.
24868 for(int inner_row = 0; inner_row < particleBank.getRows(); inner_row++) {
181 // Skip over the particle if it is photon we are trying to classify
182
2/2
✓ Branch 0 (19→20) taken 2564 times.
✓ Branch 1 (19→21) taken 19740 times.
31552 if (inner_row == row) continue;
183
184 // Check if the key exists in the calo_map
185 // This skips over REC::Particle entries without a REC::Calorimeter entry
186
2/2
✓ Branch 0 (30→31) taken 8244 times.
✓ Branch 1 (30→32) taken 11496 times.
19740 if (calo_map.find(inner_row) == calo_map.end()) continue;
187 11496 auto calo_PART = calo_map.at(inner_row);
188
189 11496 auto pid = particleBank.getInt("pid",inner_row);
190
2/2
✓ Branch 0 (34→35) taken 70 times.
✓ Branch 1 (34→36) taken 11426 times.
11496 auto mass = particle::get(particle::mass, pid);
191
192 // Skip over particle if its mass was undefined
193
2/2
✓ Branch 0 (34→35) taken 70 times.
✓ Branch 1 (34→36) taken 11426 times.
11496 if (!mass.has_value()) continue;
194 11426 auto px = particleBank.getFloat("px",inner_row);
195 11426 auto py = particleBank.getFloat("py",inner_row);
196 11426 auto pz = particleBank.getFloat("pz",inner_row);
197 11426 auto p = sqrt(px*px+py*py+pz*pz);
198 11426 auto E = sqrt(p*p+mass.value()*mass.value());
199 11426 auto th = acos(pz/p);
200 // Skip over particle if it is not in the forward detector (necessary for model compatibility)
201
2/2
✓ Branch 0 (40→41) taken 172 times.
✓ Branch 1 (40→42) taken 11254 times.
11426 if (ForwardDetectorFilter(th)==false) continue;
202
203 // Get 3-vector that points to the neighboring particle's location in the calorimeter
204 11254 ROOT::Math::XYZVector vPART = GetParticleCaloVector(calo_PART);
205
206 // Get angular distance between photon of interest and particle
207 double R = ROOT::Math::VectorUtil::Angle(vPOI, vPART);
208
209 // Logic for filling nearest neighbor variables
210
2/2
✓ Branch 0 (44→45) taken 5562 times.
✓ Branch 1 (44→63) taken 5692 times.
11254 if(pid==22) {//photon
211
212 // Apply Photon Purity Cuts to ensure this neighbor can be used in classification
213
2/2
✓ Branch 0 (46→47) taken 762 times.
✓ Branch 1 (46→48) taken 4800 times.
5562 if ( PidPurityPhotonFilter(E,calo_PART.pcal_e,th)==false ) continue;
214
215
2/2
✓ Branch 0 (48→49) taken 820 times.
✓ Branch 1 (48→50) taken 3980 times.
4800 if (R < 0.1) num_photons_0_1++;
216
2/2
✓ Branch 0 (50→51) taken 504 times.
✓ Branch 1 (50→52) taken 3476 times.
4800 if (R < 0.2) num_photons_0_2++;
217
2/2
✓ Branch 0 (52→53) taken 1656 times.
✓ Branch 1 (52→54) taken 2640 times.
4800 if (R < 0.35) num_photons_0_35++;
218
219
2/2
✓ Branch 0 (62→55) taken 7200 times.
✓ Branch 1 (62→92) taken 256 times.
7456 for (int i=0; i<m_g; ++i) {
220
4/4
✓ Branch 0 (55→56) taken 1614 times.
✓ Branch 1 (55→57) taken 5586 times.
✓ Branch 2 (57→56) taken 2930 times.
✓ Branch 3 (57→61) taken 2656 times.
7200 if (R < R_gamma[i] || R_gamma[i] == 0) {
221 int j = m_g - 1;
222
2/2
✓ Branch 0 (59→58) taken 7200 times.
✓ Branch 1 (59→60) taken 4544 times.
11744 while (j > i) {
223 7200 R_gamma[j] = R_gamma[j - 1];
224 7200 dE_gamma[j] = dE_gamma[j - 1];
225 7200 Epcal_gamma[j] = Epcal_gamma[j - 1];
226 7200 m2u_gamma[j] = m2u_gamma[j - 1];
227 7200 m2v_gamma[j] = m2v_gamma[j - 1];
228 j--;
229 }
230 4544 R_gamma[i] = R;
231 4544 dE_gamma[i] = gE - E;
232 4544 Epcal_gamma[i] = calo_PART.pcal_e;
233 4544 m2u_gamma[i] = calo_PART.pcal_m2u;
234 4544 m2v_gamma[i] = calo_PART.pcal_m2v;
235 4544 break;
236 }
237 }
238 }
239
2/2
✓ Branch 0 (63→64) taken 308 times.
✓ Branch 1 (63→68) taken 5384 times.
5692 else if(pid==11){//electron
240
3/4
✓ Branch 0 (64→65) taken 306 times.
✓ Branch 1 (64→67) taken 2 times.
✗ Branch 2 (65→66) not taken.
✓ Branch 3 (65→67) taken 306 times.
308 if(R<R_e || R_e==0){
241 R_e = R;
242 308 dE_e = gE - E;
243 }
244 }
245
10/10
✓ Branch 0 (68→69) taken 5190 times.
✓ Branch 1 (68→70) taken 194 times.
✓ Branch 2 (69→70) taken 1666 times.
✓ Branch 3 (69→71) taken 3524 times.
✓ Branch 4 (71→70) taken 602 times.
✓ Branch 5 (71→72) taken 2922 times.
✓ Branch 6 (72→70) taken 16 times.
✓ Branch 7 (72→73) taken 2906 times.
✓ Branch 8 (73→70) taken 14 times.
✓ Branch 9 (73→82) taken 2892 times.
5384 else if(pid==2212||pid==-2212||pid==211||pid==-211||pid==321||pid==-321){//charged hadron
246
2/2
✓ Branch 0 (81→74) taken 2946 times.
✓ Branch 1 (81→92) taken 64 times.
3010 for (int i=0; i<m_ch; ++i) {
247
4/4
✓ Branch 0 (74→75) taken 444 times.
✓ Branch 1 (74→76) taken 2502 times.
✓ Branch 2 (76→75) taken 1984 times.
✓ Branch 3 (76→80) taken 518 times.
2946 if (R < R_ch[i] || R_ch[i] == 0) {
248 int j = m_ch - 1;
249
2/2
✓ Branch 0 (78→77) taken 2038 times.
✓ Branch 1 (78→79) taken 2428 times.
4466 while (j > i) {
250 2038 R_ch[j] = R_ch[j - 1];
251 2038 dE_ch[j] = dE_ch[j - 1];
252 2038 Epcal_ch[j] = Epcal_ch[j - 1];
253 2038 m2u_ch[j] = m2u_ch[j - 1];
254 2038 m2v_ch[j] = m2v_ch[j - 1];
255 j--;
256 }
257 2428 R_ch[i] = R;
258 2428 dE_ch[i] = gE- E;
259 2428 Epcal_ch[i] = calo_PART.pcal_e;
260 2428 m2u_ch[i] = calo_PART.pcal_m2u;
261 2428 m2v_ch[i] = calo_PART.pcal_m2v;
262 2428 break;
263 }
264 }
265 }
266
1/2
✓ Branch 0 (82→90) taken 2892 times.
✗ Branch 1 (82→91) not taken.
2892 else if(pid==2112||pid==-2112){//neutral hadron
267
2/2
✓ Branch 0 (90→83) taken 3734 times.
✓ Branch 1 (90→92) taken 264 times.
3998 for (int i=0; i<m_nh; ++i) {
268
4/4
✓ Branch 0 (83→84) taken 776 times.
✓ Branch 1 (83→85) taken 2958 times.
✓ Branch 2 (85→84) taken 1852 times.
✓ Branch 3 (85→89) taken 1106 times.
3734 if (R < R_nh[i] || R_nh[i] == 0) {
269 int j = m_nh - 1;
270
2/2
✓ Branch 0 (87→86) taken 2050 times.
✓ Branch 1 (87→88) taken 2628 times.
4678 while (j > i) {
271 2050 R_nh[j] = R_nh[j - 1];
272 2050 dE_nh[j] = dE_nh[j - 1];
273 2050 Epcal_nh[j] = Epcal_nh[j - 1];
274 2050 m2u_nh[j] = m2u_nh[j - 1];
275 2050 m2v_nh[j] = m2v_nh[j - 1];
276 j--;
277 }
278 2628 R_nh[i] = R;
279 2628 dE_nh[i] = gE - E;
280 2628 Epcal_nh[i] = calo_PART.pcal_e;
281 2628 m2u_nh[i] = calo_PART.pcal_m2u;
282 2628 m2v_nh[i] = calo_PART.pcal_m2v;
283 2628 break;
284 }
285 }
286 }
287 else{//unrecognized OR uncompatible particle type for the trained model
288 continue;
289 }
290 }
291
292 // Create and populate input_data vector for the ML model
293 std::vector<float> input_data = {
294 static_cast<float>(gE), static_cast<float>(gEpcal), static_cast<float>(gTheta),
295 2564 static_cast<float>(gm2u), static_cast<float>(gm2v), static_cast<float>(R_e),
296 2564 static_cast<float>(dE_e)
297 2564 };
298
299
300
3/4
✓ Branch 0 (98→99) taken 7692 times.
✗ Branch 1 (98→152) not taken.
✓ Branch 2 (100→98) taken 7692 times.
✓ Branch 3 (100→103) taken 2564 times.
10256 for (int i = 0; i < m_g; ++i) input_data.push_back(static_cast<float>(R_gamma[i]));
301
3/4
✓ Branch 0 (101→102) taken 7692 times.
✗ Branch 1 (101→152) not taken.
✓ Branch 2 (103→101) taken 7692 times.
✓ Branch 3 (103→106) taken 2564 times.
10256 for (int i = 0; i < m_g; ++i) input_data.push_back(static_cast<float>(dE_gamma[i]));
302
3/4
✓ Branch 0 (104→105) taken 7692 times.
✗ Branch 1 (104→152) not taken.
✓ Branch 2 (106→104) taken 7692 times.
✓ Branch 3 (106→109) taken 2564 times.
10256 for (int i = 0; i < m_g; ++i) input_data.push_back(static_cast<float>(Epcal_gamma[i]));
303
3/4
✓ Branch 0 (107→108) taken 7692 times.
✗ Branch 1 (107→152) not taken.
✓ Branch 2 (109→107) taken 7692 times.
✓ Branch 3 (109→112) taken 2564 times.
10256 for (int i = 0; i < m_g; ++i) input_data.push_back(static_cast<float>(m2u_gamma[i]));
304
3/4
✓ Branch 0 (110→111) taken 7692 times.
✗ Branch 1 (110→152) not taken.
✓ Branch 2 (112→110) taken 7692 times.
✓ Branch 3 (112→115) taken 2564 times.
10256 for (int i = 0; i < m_g; ++i) input_data.push_back(static_cast<float>(m2v_gamma[i]));
305
306
307
308
3/4
✓ Branch 0 (113→114) taken 5128 times.
✗ Branch 1 (113→152) not taken.
✓ Branch 2 (115→113) taken 5128 times.
✓ Branch 3 (115→118) taken 2564 times.
7692 for (int i = 0; i < m_ch; ++i) input_data.push_back(static_cast<float>(R_ch[i]));
309
3/4
✓ Branch 0 (116→117) taken 5128 times.
✗ Branch 1 (116→152) not taken.
✓ Branch 2 (118→116) taken 5128 times.
✓ Branch 3 (118→121) taken 2564 times.
7692 for (int i = 0; i < m_ch; ++i) input_data.push_back(static_cast<float>(dE_ch[i]));
310
3/4
✓ Branch 0 (119→120) taken 5128 times.
✗ Branch 1 (119→152) not taken.
✓ Branch 2 (121→119) taken 5128 times.
✓ Branch 3 (121→124) taken 2564 times.
7692 for (int i = 0; i < m_ch; ++i) input_data.push_back(static_cast<float>(Epcal_ch[i]));
311
3/4
✓ Branch 0 (122→123) taken 5128 times.
✗ Branch 1 (122→152) not taken.
✓ Branch 2 (124→122) taken 5128 times.
✓ Branch 3 (124→127) taken 2564 times.
7692 for (int i = 0; i < m_ch; ++i) input_data.push_back(static_cast<float>(m2u_ch[i]));
312
3/4
✓ Branch 0 (125→126) taken 5128 times.
✗ Branch 1 (125→152) not taken.
✓ Branch 2 (127→125) taken 5128 times.
✓ Branch 3 (127→130) taken 2564 times.
7692 for (int i = 0; i < m_ch; ++i) input_data.push_back(static_cast<float>(m2v_ch[i]));
313
314
315
316
3/4
✓ Branch 0 (128→129) taken 5128 times.
✗ Branch 1 (128→152) not taken.
✓ Branch 2 (130→128) taken 5128 times.
✓ Branch 3 (130→133) taken 2564 times.
7692 for (int i = 0; i < m_nh; ++i) input_data.push_back(static_cast<float>(R_nh[i]));
317
3/4
✓ Branch 0 (131→132) taken 5128 times.
✗ Branch 1 (131→152) not taken.
✓ Branch 2 (133→131) taken 5128 times.
✓ Branch 3 (133→136) taken 2564 times.
7692 for (int i = 0; i < m_nh; ++i) input_data.push_back(static_cast<float>(dE_nh[i]));
318
3/4
✓ Branch 0 (134→135) taken 5128 times.
✗ Branch 1 (134→152) not taken.
✓ Branch 2 (136→134) taken 5128 times.
✓ Branch 3 (136→139) taken 2564 times.
7692 for (int i = 0; i < m_nh; ++i) input_data.push_back(static_cast<float>(Epcal_nh[i]));
319
3/4
✓ Branch 0 (137→138) taken 5128 times.
✗ Branch 1 (137→152) not taken.
✓ Branch 2 (139→137) taken 5128 times.
✓ Branch 3 (139→142) taken 2564 times.
7692 for (int i = 0; i < m_nh; ++i) input_data.push_back(static_cast<float>(m2u_nh[i]));
320
3/4
✓ Branch 0 (140→141) taken 5128 times.
✗ Branch 1 (140→152) not taken.
✓ Branch 2 (142→140) taken 5128 times.
✓ Branch 3 (142→143) taken 2564 times.
7692 for (int i = 0; i < m_nh; ++i) input_data.push_back(static_cast<float>(m2v_nh[i]));
321
322
323
1/2
✓ Branch 0 (143→144) taken 2564 times.
✗ Branch 1 (143→152) not taken.
2564 input_data.push_back(static_cast<float>(num_photons_0_1));
324
1/2
✓ Branch 0 (144→145) taken 2564 times.
✗ Branch 1 (144→152) not taken.
2564 input_data.push_back(static_cast<float>(num_photons_0_2));
325
1/4
✓ Branch 0 (145→146) taken 2564 times.
✗ Branch 1 (145→152) not taken.
✗ Branch 2 (152→153) not taken.
✗ Branch 3 (152→155) not taken.
2564 input_data.push_back(static_cast<float>(num_photons_0_35));
326
327
1/2
✓ Branch 0 (146→147) taken 2564 times.
✗ Branch 1 (146→152) not taken.
2564 return ClassifyPhoton(input_data, runnum);
328 }
329
330 2564 bool PhotonGBTFilter::ClassifyPhoton(std::vector<float> const &input_data, int const runnum) const {
331 2564 auto modelFunction = getModelFunction(runnum);
332 double sigmoid_x = modelFunction(input_data);
333 2564 double prediction = 1 / (1 + exp(-sigmoid_x));
334
1/2
✓ Branch 0 (6→7) taken 2564 times.
✗ Branch 1 (6→8) not taken.
5128 return (prediction > o_threshold);
335 }
336
337 2000 std::map<int, PhotonGBTFilter::calo_row_data> PhotonGBTFilter::GetCaloMap(hipo::bank const& bank) const
338 {
339 std::map<int, PhotonGBTFilter::calo_row_data> calo_map;
340 // Loop over REC::Calorimeter rows
341 // Here we use bank.getRows() to purposefully ignore upstream filters
342
2/2
✓ Branch 0 (40→3) taken 12258 times.
✓ Branch 1 (40→41) taken 2000 times.
14258 for(int row = 0; row < bank.getRows(); row++){
343 12258 auto pindex = bank.getShort("pindex",row);
344 12258 auto x = bank.getFloat("x",row);
345 12258 auto y = bank.getFloat("y",row);
346 12258 auto z = bank.getFloat("z",row);
347 12258 auto m2u = bank.getFloat("m2u",row);
348 12258 auto m2v = bank.getFloat("m2v",row);
349 12258 auto layer = bank.getInt("layer",row);
350 12258 auto e = bank.getFloat("energy",row);
351
352 // Ensure an entry exists in the map for the given pindex
353
2/2
✓ Branch 0 (20→21) taken 7572 times.
✓ Branch 1 (20→23) taken 4686 times.
12258 if (calo_map.find(pindex) == calo_map.end()) {
354
1/2
✓ Branch 0 (21→22) taken 7572 times.
✗ Branch 1 (21→42) not taken.
7572 calo_map[pindex] = PhotonGBTFilter::calo_row_data();
355 }
356
357
3/4
✓ Branch 0 (23→24) taken 6292 times.
✓ Branch 1 (23→31) taken 3706 times.
✓ Branch 2 (23→35) taken 2260 times.
✗ Branch 3 (23→39) not taken.
12258 switch(layer){
358 6292 case 1: // pcal
359
1/2
✓ Branch 0 (24→25) taken 6292 times.
✗ Branch 1 (24→42) not taken.
6292 calo_map[pindex].pcal_x = x;
360
1/2
✓ Branch 0 (25→26) taken 6292 times.
✗ Branch 1 (25→42) not taken.
6292 calo_map[pindex].pcal_y = y;
361
1/2
✓ Branch 0 (26→27) taken 6292 times.
✗ Branch 1 (26→42) not taken.
6292 calo_map[pindex].pcal_z = z;
362
1/2
✓ Branch 0 (27→28) taken 6292 times.
✗ Branch 1 (27→42) not taken.
6292 calo_map[pindex].pcal_e = e;
363
1/2
✓ Branch 0 (28→29) taken 6292 times.
✗ Branch 1 (28→42) not taken.
6292 calo_map[pindex].pcal_m2u = m2u;
364
1/2
✓ Branch 0 (29→30) taken 6292 times.
✗ Branch 1 (29→42) not taken.
6292 calo_map[pindex].pcal_m2v = m2v;
365 6292 break;
366 3706 case 4: // ecin
367
1/2
✓ Branch 0 (31→32) taken 3706 times.
✗ Branch 1 (31→42) not taken.
3706 calo_map[pindex].ecin_x = x;
368
1/2
✓ Branch 0 (32→33) taken 3706 times.
✗ Branch 1 (32→42) not taken.
3706 calo_map[pindex].ecin_y = y;
369
1/2
✓ Branch 0 (33→34) taken 3706 times.
✗ Branch 1 (33→42) not taken.
3706 calo_map[pindex].ecin_z = z;
370 3706 break;
371 2260 case 7: // ecout
372
1/2
✓ Branch 0 (35→36) taken 2260 times.
✗ Branch 1 (35→42) not taken.
2260 calo_map[pindex].ecout_x = x;
373
1/2
✓ Branch 0 (36→37) taken 2260 times.
✗ Branch 1 (36→42) not taken.
2260 calo_map[pindex].ecout_y = y;
374
1/2
✓ Branch 0 (37→38) taken 2260 times.
✗ Branch 1 (37→42) not taken.
2260 calo_map[pindex].ecout_z = z;
375 2260 break;
376 }
377 }
378 2000 return calo_map;
379 }
380
381 13818 ROOT::Math::XYZVector PhotonGBTFilter::GetParticleCaloVector(PhotonGBTFilter::calo_row_data calo_row) const {
382 // Determine the 3-vector location of where the photon of interest's calo deposition is
383 // First we check the pcal coords, then ecin, then ecout
384 ROOT::Math::XYZVector v;
385
2/2
✓ Branch 0 (2→3) taken 2004 times.
✓ Branch 1 (2→6) taken 11814 times.
13818 if (calo_row.pcal_x == 0) {
386
2/2
✓ Branch 0 (3→4) taken 756 times.
✓ Branch 1 (3→5) taken 1248 times.
2004 if (calo_row.ecin_x == 0) {
387 756 v.SetXYZ(calo_row.ecout_x, calo_row.ecout_y, calo_row.ecout_z);
388 } else {
389 1248 v.SetXYZ(calo_row.ecin_x, calo_row.ecin_y, calo_row.ecin_z);
390 }
391 } else {
392 11814 v.SetXYZ(calo_row.pcal_x, calo_row.pcal_y, calo_row.pcal_z);
393 }
394 13818 return v;
395 }
396
397 2564 std::function<double(std::vector<float> const &)> PhotonGBTFilter::getModelFunction(int runnum) const {
398
399
1/2
✓ Branch 0 (8→3) taken 17948 times.
✗ Branch 1 (8→9) not taken.
17948 for (const auto &entry : modelMap) {
400
4/6
✓ Branch 0 (3→4) taken 17948 times.
✗ Branch 1 (3→7) not taken.
✓ Branch 2 (4→5) taken 2564 times.
✓ Branch 3 (4→7) taken 15384 times.
✓ Branch 4 (5→6) taken 2564 times.
✗ Branch 5 (5→7) not taken.
17948 if (runnum >= std::get<0>(entry.first) && runnum <= std::get<1>(entry.first) && o_pass == std::get<2>(entry.first)) {
401 2564 return entry.second;
402 }
403 }
404
405 // Default to RGA inbending pass1 if no match found
406 m_log->Warn("Run Number {} with pass {} has no matching PhotonGBT model...Defaulting to RGA inbending pass1...", runnum, o_pass);
407 return [](std::vector<float> const &data) { return ApplyCatboostModel_RGA_inbending_pass1(data); };
408 }
409
410 1 void PhotonGBTFilter::Stop()
411 {
412 1 }
413
414 }
415