===============================================================================================================================================================================================
variant2: for each run define phi-averaged A for normalization channel (Dref,16) and then, divide Rijk on it, i.e. get RRijk
eta=27
64 bool haveInput =
false;
66 bool haveOutput =
false;
70 bool cleanSimTracks =
false;
71 bool writeAllEvents =
false;
72 bool writeRecTracks =
false;
73 bool writeHitIterMasks =
false;
74 bool applyCCC =
false;
75 bool allSeeds =
false;
78 long long maxevt = -1;
80 int cutValueCCC = 1620;
83 for (
int i = 1;
i <
argc; ++
i) {
84 mArgs.push_back(
argv[
i]);
88 while (
i != mArgs.end()) {
91 if (*
i ==
"-h" || *
i ==
"-help" || *
i ==
"--help") {
94 }
else if (*
i ==
"--input") {
98 }
else if (*
i ==
"--output") {
102 }
else if (*
i ==
"--geo") {
105 }
else if (*
i ==
"--verbosity") {
108 }
else if (*
i ==
"--maxevt") {
110 maxevt = std::atoi(
i->c_str());
111 }
else if (*
i ==
"--clean-sim-tracks") {
112 cleanSimTracks =
true;
113 }
else if (*
i ==
"--write-all-events") {
114 writeAllEvents =
true;
115 }
else if (*
i ==
"--write-rec-tracks") {
116 writeRecTracks =
true;
117 }
else if (*
i ==
"--write-hit-iter-masks") {
118 writeHitIterMasks =
true;
119 }
else if (*
i ==
"--apply-ccc") {
122 cutValueCCC = std::atoi(
i->c_str());
124 }
else if (*
i ==
"--all-seeds") {
127 fprintf(
stderr,
"Error: Unknown option/argument '%s'.\n",
i->c_str());
134 if (not haveOutput
or not haveInput) {
135 fprintf(
stderr,
"Error: both input and output are required\n");
148 TTree*
t = (TTree*)
f->Get(
"trackingNtuple/tree");
150 bool hasPh2hits =
t->GetBranch(
"ph2_isLower") !=
nullptr;
154 const unsigned int nTotalLayers = lnc.nLayers();
158 std::vector<int> nhitstot(nTotalLayers, 0);
160 unsigned long long event;
161 t->SetBranchAddress(
"event", &
event);
164 std::vector<float>* sim_eta = 0;
165 std::vector<float>* sim_px = 0;
166 std::vector<float>* sim_py = 0;
167 std::vector<float>* sim_pz = 0;
168 std::vector<int>* sim_parentVtxIdx = 0;
169 std::vector<int>* sim_q = 0;
170 std::vector<int>* sim_event = 0;
171 std::vector<int>* sim_bunchCrossing = 0;
172 std::vector<int>* sim_nValid = 0;
173 t->SetBranchAddress(
"sim_eta", &sim_eta);
174 t->SetBranchAddress(
"sim_px", &sim_px);
175 t->SetBranchAddress(
"sim_py", &sim_py);
176 t->SetBranchAddress(
"sim_pz", &sim_pz);
177 t->SetBranchAddress(
"sim_parentVtxIdx", &sim_parentVtxIdx);
178 t->SetBranchAddress(
"sim_q", &sim_q);
179 t->SetBranchAddress(
"sim_event", &sim_event);
180 t->SetBranchAddress(
"sim_bunchCrossing", &sim_bunchCrossing);
181 t->SetBranchAddress(
"sim_nValid", &sim_nValid);
183 std::vector<vector<int>>* sim_trkIdx = 0;
184 t->SetBranchAddress(
"sim_trkIdx", &sim_trkIdx);
187 std::vector<float>* simvtx_x = 0;
188 std::vector<float>* simvtx_y = 0;
189 std::vector<float>* simvtx_z = 0;
190 t->SetBranchAddress(
"simvtx_x", &simvtx_x);
191 t->SetBranchAddress(
"simvtx_y", &simvtx_y);
192 t->SetBranchAddress(
"simvtx_z", &simvtx_z);
195 std::vector<short>* simhit_process = 0;
196 std::vector<int>* simhit_particle = 0;
197 std::vector<int>* simhit_simTrkIdx = 0;
198 std::vector<float>* simhit_x = 0;
199 std::vector<float>* simhit_y = 0;
200 std::vector<float>* simhit_z = 0;
201 std::vector<float>* simhit_px = 0;
202 std::vector<float>* simhit_py = 0;
203 std::vector<float>* simhit_pz = 0;
204 t->SetBranchAddress(
"simhit_process", &simhit_process);
205 t->SetBranchAddress(
"simhit_particle", &simhit_particle);
206 t->SetBranchAddress(
"simhit_simTrkIdx", &simhit_simTrkIdx);
207 t->SetBranchAddress(
"simhit_x", &simhit_x);
208 t->SetBranchAddress(
"simhit_y", &simhit_y);
209 t->SetBranchAddress(
"simhit_z", &simhit_z);
210 t->SetBranchAddress(
"simhit_px", &simhit_px);
211 t->SetBranchAddress(
"simhit_py", &simhit_py);
212 t->SetBranchAddress(
"simhit_pz", &simhit_pz);
214 std::vector<std::vector<int>>* simhit_hitIdx = 0;
215 t->SetBranchAddress(
"simhit_hitIdx", &simhit_hitIdx);
216 std::vector<std::vector<int>>* simhit_hitType = 0;
217 t->SetBranchAddress(
"simhit_hitType", &simhit_hitType);
220 std::vector<int>* trk_q = 0;
221 std::vector<unsigned int>* trk_nValid = 0;
222 std::vector<int>* trk_seedIdx = 0;
223 std::vector<unsigned long long>* trk_algoMask = 0;
224 std::vector<unsigned int>*
trk_algo = 0;
225 std::vector<unsigned int>* trk_originalAlgo = 0;
226 std::vector<float>* trk_nChi2 = 0;
227 std::vector<float>* trk_px = 0;
228 std::vector<float>* trk_py = 0;
229 std::vector<float>* trk_pz = 0;
230 std::vector<float>* trk_pt = 0;
231 std::vector<float>* trk_phi = 0;
232 std::vector<float>* trk_lambda = 0;
233 std::vector<float>* trk_refpoint_x = 0;
234 std::vector<float>* trk_refpoint_y = 0;
235 std::vector<float>* trk_refpoint_z = 0;
236 std::vector<float>* trk_dxyErr = 0;
237 std::vector<float>* trk_dzErr = 0;
238 std::vector<float>* trk_ptErr = 0;
239 std::vector<float>* trk_phiErr = 0;
240 std::vector<float>* trk_lambdaErr = 0;
241 t->SetBranchAddress(
"trk_q", &trk_q);
242 t->SetBranchAddress(
"trk_nValid", &trk_nValid);
243 t->SetBranchAddress(
"trk_seedIdx", &trk_seedIdx);
244 t->SetBranchAddress(
"trk_algoMask", &trk_algoMask);
245 t->SetBranchAddress(
"trk_algo", &
trk_algo);
246 t->SetBranchAddress(
"trk_originalAlgo", &trk_originalAlgo);
247 t->SetBranchAddress(
"trk_nChi2", &trk_nChi2);
248 t->SetBranchAddress(
"trk_px", &trk_px);
249 t->SetBranchAddress(
"trk_py", &trk_py);
250 t->SetBranchAddress(
"trk_pz", &trk_pz);
251 t->SetBranchAddress(
"trk_pt", &trk_pt);
252 t->SetBranchAddress(
"trk_phi", &trk_phi);
253 t->SetBranchAddress(
"trk_lambda", &trk_lambda);
254 t->SetBranchAddress(
"trk_refpoint_x", &trk_refpoint_x);
255 t->SetBranchAddress(
"trk_refpoint_y", &trk_refpoint_y);
256 t->SetBranchAddress(
"trk_refpoint_z", &trk_refpoint_z);
257 t->SetBranchAddress(
"trk_dxyErr", &trk_dxyErr);
258 t->SetBranchAddress(
"trk_dzErr", &trk_dzErr);
259 t->SetBranchAddress(
"trk_ptErr", &trk_ptErr);
260 t->SetBranchAddress(
"trk_phiErr", &trk_phiErr);
261 t->SetBranchAddress(
"trk_lambdaErr", &trk_lambdaErr);
263 std::vector<std::vector<int>>* trk_hitIdx = 0;
264 t->SetBranchAddress(
"trk_hitIdx", &trk_hitIdx);
265 std::vector<std::vector<int>>* trk_hitType = 0;
266 t->SetBranchAddress(
"trk_hitType", &trk_hitType);
269 std::vector<float>* see_stateTrajGlbX = 0;
270 std::vector<float>* see_stateTrajGlbY = 0;
271 std::vector<float>* see_stateTrajGlbZ = 0;
272 std::vector<float>* see_stateTrajGlbPx = 0;
273 std::vector<float>* see_stateTrajGlbPy = 0;
274 std::vector<float>* see_stateTrajGlbPz = 0;
275 std::vector<float>* see_eta = 0;
276 std::vector<float>* see_pt = 0;
277 std::vector<float>* see_stateCcov00 = 0;
278 std::vector<float>* see_stateCcov01 = 0;
279 std::vector<float>* see_stateCcov02 = 0;
280 std::vector<float>* see_stateCcov03 = 0;
281 std::vector<float>* see_stateCcov04 = 0;
282 std::vector<float>* see_stateCcov05 = 0;
283 std::vector<float>* see_stateCcov11 = 0;
284 std::vector<float>* see_stateCcov12 = 0;
285 std::vector<float>* see_stateCcov13 = 0;
286 std::vector<float>* see_stateCcov14 = 0;
287 std::vector<float>* see_stateCcov15 = 0;
288 std::vector<float>* see_stateCcov22 = 0;
289 std::vector<float>* see_stateCcov23 = 0;
290 std::vector<float>* see_stateCcov24 = 0;
291 std::vector<float>* see_stateCcov25 = 0;
292 std::vector<float>* see_stateCcov33 = 0;
293 std::vector<float>* see_stateCcov34 = 0;
294 std::vector<float>* see_stateCcov35 = 0;
295 std::vector<float>* see_stateCcov44 = 0;
296 std::vector<float>* see_stateCcov45 = 0;
297 std::vector<float>* see_stateCcov55 = 0;
298 std::vector<std::vector<float>>* see_stateCurvCov = 0;
299 std::vector<int>* see_q = 0;
300 std::vector<unsigned int>* see_algo = 0;
301 t->SetBranchAddress(
"see_stateTrajGlbX", &see_stateTrajGlbX);
302 t->SetBranchAddress(
"see_stateTrajGlbY", &see_stateTrajGlbY);
303 t->SetBranchAddress(
"see_stateTrajGlbZ", &see_stateTrajGlbZ);
304 t->SetBranchAddress(
"see_stateTrajGlbPx", &see_stateTrajGlbPx);
305 t->SetBranchAddress(
"see_stateTrajGlbPy", &see_stateTrajGlbPy);
306 t->SetBranchAddress(
"see_stateTrajGlbPz", &see_stateTrajGlbPz);
307 t->SetBranchAddress(
"see_eta", &see_eta);
308 t->SetBranchAddress(
"see_pt", &see_pt);
310 bool hasCartCov =
t->GetBranch(
"see_stateCcov00") !=
nullptr;
312 t->SetBranchAddress(
"see_stateCcov00", &see_stateCcov00);
313 t->SetBranchAddress(
"see_stateCcov01", &see_stateCcov01);
314 t->SetBranchAddress(
"see_stateCcov02", &see_stateCcov02);
315 t->SetBranchAddress(
"see_stateCcov03", &see_stateCcov03);
316 t->SetBranchAddress(
"see_stateCcov04", &see_stateCcov04);
317 t->SetBranchAddress(
"see_stateCcov05", &see_stateCcov05);
318 t->SetBranchAddress(
"see_stateCcov11", &see_stateCcov11);
319 t->SetBranchAddress(
"see_stateCcov12", &see_stateCcov12);
320 t->SetBranchAddress(
"see_stateCcov13", &see_stateCcov13);
321 t->SetBranchAddress(
"see_stateCcov14", &see_stateCcov14);
322 t->SetBranchAddress(
"see_stateCcov15", &see_stateCcov15);
323 t->SetBranchAddress(
"see_stateCcov22", &see_stateCcov22);
324 t->SetBranchAddress(
"see_stateCcov23", &see_stateCcov23);
325 t->SetBranchAddress(
"see_stateCcov24", &see_stateCcov24);
326 t->SetBranchAddress(
"see_stateCcov25", &see_stateCcov25);
327 t->SetBranchAddress(
"see_stateCcov33", &see_stateCcov33);
328 t->SetBranchAddress(
"see_stateCcov34", &see_stateCcov34);
329 t->SetBranchAddress(
"see_stateCcov35", &see_stateCcov35);
330 t->SetBranchAddress(
"see_stateCcov44", &see_stateCcov44);
331 t->SetBranchAddress(
"see_stateCcov45", &see_stateCcov45);
332 t->SetBranchAddress(
"see_stateCcov55", &see_stateCcov55);
334 t->SetBranchAddress(
"see_stateCurvCov", &see_stateCurvCov);
336 t->SetBranchAddress(
"see_q", &see_q);
337 t->SetBranchAddress(
"see_algo", &see_algo);
339 std::vector<std::vector<int>>* see_hitIdx = 0;
340 t->SetBranchAddress(
"see_hitIdx", &see_hitIdx);
341 std::vector<std::vector<int>>* see_hitType = 0;
342 t->SetBranchAddress(
"see_hitType", &see_hitType);
345 vector<unsigned short>* pix_det = 0;
346 vector<unsigned short>* pix_lay = 0;
347 vector<unsigned int>* pix_detId = 0;
348 vector<float>* pix_x = 0;
349 vector<float>* pix_y = 0;
350 vector<float>* pix_z = 0;
351 vector<float>* pix_xx = 0;
352 vector<float>* pix_xy = 0;
353 vector<float>* pix_yy = 0;
354 vector<float>* pix_yz = 0;
355 vector<float>* pix_zz = 0;
356 vector<float>* pix_zx = 0;
357 vector<int>* pix_csize_col = 0;
358 vector<int>* pix_csize_row = 0;
359 vector<uint64_t>* pix_usedMask = 0;
361 bool has910_det_lay =
t->GetBranch(
"pix_det") ==
nullptr;
362 if (has910_det_lay) {
363 t->SetBranchAddress(
"pix_subdet", &pix_det);
364 t->SetBranchAddress(
"pix_layer", &pix_lay);
366 t->SetBranchAddress(
"pix_det", &pix_det);
367 t->SetBranchAddress(
"pix_lay", &pix_lay);
369 t->SetBranchAddress(
"pix_detId", &pix_detId);
370 t->SetBranchAddress(
"pix_x", &pix_x);
371 t->SetBranchAddress(
"pix_y", &pix_y);
372 t->SetBranchAddress(
"pix_z", &pix_z);
373 t->SetBranchAddress(
"pix_xx", &pix_xx);
374 t->SetBranchAddress(
"pix_xy", &pix_xy);
375 t->SetBranchAddress(
"pix_yy", &pix_yy);
376 t->SetBranchAddress(
"pix_yz", &pix_yz);
377 t->SetBranchAddress(
"pix_zz", &pix_zz);
378 t->SetBranchAddress(
"pix_zx", &pix_zx);
379 t->SetBranchAddress(
"pix_clustSizeCol", &pix_csize_col);
380 t->SetBranchAddress(
"pix_clustSizeRow", &pix_csize_row);
381 if (writeHitIterMasks) {
382 t->SetBranchAddress(
"pix_usedMask", &pix_usedMask);
385 vector<vector<int>>* pix_simHitIdx = 0;
386 t->SetBranchAddress(
"pix_simHitIdx", &pix_simHitIdx);
387 vector<vector<float>>* pix_chargeFraction = 0;
388 t->SetBranchAddress(
"pix_chargeFraction", &pix_chargeFraction);
391 vector<short>* glu_isBarrel = 0;
392 vector<unsigned int>* glu_det = 0;
393 vector<unsigned int>* glu_lay = 0;
394 vector<unsigned int>* glu_detId = 0;
395 vector<int>* glu_monoIdx = 0;
396 vector<int>* glu_stereoIdx = 0;
397 vector<float>* glu_x = 0;
398 vector<float>* glu_y = 0;
399 vector<float>* glu_z = 0;
400 vector<float>* glu_xx = 0;
401 vector<float>* glu_xy = 0;
402 vector<float>* glu_yy = 0;
403 vector<float>* glu_yz = 0;
404 vector<float>* glu_zz = 0;
405 vector<float>* glu_zx = 0;
407 t->SetBranchAddress(
"glu_isBarrel", &glu_isBarrel);
408 if (has910_det_lay) {
409 t->SetBranchAddress(
"glu_subdet", &glu_det);
410 t->SetBranchAddress(
"glu_layer", &glu_lay);
412 t->SetBranchAddress(
"glu_det", &glu_det);
413 t->SetBranchAddress(
"glu_lay", &glu_lay);
415 t->SetBranchAddress(
"glu_detId", &glu_detId);
416 t->SetBranchAddress(
"glu_monoIdx", &glu_monoIdx);
417 t->SetBranchAddress(
"glu_stereoIdx", &glu_stereoIdx);
418 t->SetBranchAddress(
"glu_x", &glu_x);
419 t->SetBranchAddress(
"glu_y", &glu_y);
420 t->SetBranchAddress(
"glu_z", &glu_z);
421 t->SetBranchAddress(
"glu_xx", &glu_xx);
422 t->SetBranchAddress(
"glu_xy", &glu_xy);
423 t->SetBranchAddress(
"glu_yy", &glu_yy);
424 t->SetBranchAddress(
"glu_yz", &glu_yz);
425 t->SetBranchAddress(
"glu_zz", &glu_zz);
426 t->SetBranchAddress(
"glu_zx", &glu_zx);
429 vector<short>* str_isBarrel = 0;
430 vector<short>* str_isStereo = 0;
431 vector<unsigned int>* str_det = 0;
432 vector<unsigned int>* str_lay = 0;
433 vector<unsigned int>* str_detId = 0;
434 vector<unsigned int>* str_simType = 0;
435 vector<float>* str_x = 0;
436 vector<float>* str_y = 0;
437 vector<float>* str_z = 0;
438 vector<float>* str_xx = 0;
439 vector<float>* str_xy = 0;
440 vector<float>* str_yy = 0;
441 vector<float>* str_yz = 0;
442 vector<float>* str_zz = 0;
443 vector<float>* str_zx = 0;
444 vector<float>* str_chargePerCM = 0;
445 vector<int>* str_csize = 0;
446 vector<uint64_t>* str_usedMask = 0;
447 vector<vector<int>>* str_simHitIdx = 0;
448 vector<vector<float>>* str_chargeFraction = 0;
450 t->SetBranchAddress(
"str_isBarrel", &str_isBarrel);
451 t->SetBranchAddress(
"str_isStereo", &str_isStereo);
452 if (has910_det_lay) {
453 t->SetBranchAddress(
"str_subdet", &str_det);
454 t->SetBranchAddress(
"str_layer", &str_lay);
456 t->SetBranchAddress(
"str_det", &str_det);
457 t->SetBranchAddress(
"str_lay", &str_lay);
459 t->SetBranchAddress(
"str_detId", &str_detId);
460 t->SetBranchAddress(
"str_simType", &str_simType);
461 t->SetBranchAddress(
"str_x", &str_x);
462 t->SetBranchAddress(
"str_y", &str_y);
463 t->SetBranchAddress(
"str_z", &str_z);
464 t->SetBranchAddress(
"str_xx", &str_xx);
465 t->SetBranchAddress(
"str_xy", &str_xy);
466 t->SetBranchAddress(
"str_yy", &str_yy);
467 t->SetBranchAddress(
"str_yz", &str_yz);
468 t->SetBranchAddress(
"str_zz", &str_zz);
469 t->SetBranchAddress(
"str_zx", &str_zx);
470 t->SetBranchAddress(
"str_chargePerCM", &str_chargePerCM);
471 t->SetBranchAddress(
"str_clustSize", &str_csize);
472 if (writeHitIterMasks) {
473 t->SetBranchAddress(
"str_usedMask", &str_usedMask);
476 t->SetBranchAddress(
"str_simHitIdx", &str_simHitIdx);
477 t->SetBranchAddress(
"str_chargeFraction", &str_chargeFraction);
480 vector<unsigned short>* ph2_isLower = 0;
481 vector<unsigned short>* ph2_subdet = 0;
482 vector<unsigned short>* ph2_layer = 0;
483 vector<unsigned int>* ph2_detId = 0;
484 vector<unsigned short>* ph2_simType = 0;
485 vector<float>* ph2_x = 0;
486 vector<float>* ph2_y = 0;
487 vector<float>* ph2_z = 0;
488 vector<float>* ph2_xx = 0;
489 vector<float>* ph2_xy = 0;
490 vector<float>* ph2_yy = 0;
491 vector<float>* ph2_yz = 0;
492 vector<float>* ph2_zz = 0;
493 vector<float>* ph2_zx = 0;
494 vector<uint64_t>* ph2_usedMask = 0;
495 vector<vector<int>>* ph2_simHitIdx = 0;
496 if (hasPh2hits && applyCCC)
497 std::cout <<
"WARNING: applyCCC is set for Phase2 inputs: applyCCC will be ignored" << std::endl;
499 t->SetBranchAddress(
"ph2_isLower", &ph2_isLower);
500 t->SetBranchAddress(
"ph2_subdet", &ph2_subdet);
501 t->SetBranchAddress(
"ph2_layer", &ph2_layer);
502 t->SetBranchAddress(
"ph2_detId", &ph2_detId);
503 t->SetBranchAddress(
"ph2_simType", &ph2_simType);
504 t->SetBranchAddress(
"ph2_x", &ph2_x);
505 t->SetBranchAddress(
"ph2_y", &ph2_y);
506 t->SetBranchAddress(
"ph2_z", &ph2_z);
507 t->SetBranchAddress(
"ph2_xx", &ph2_xx);
508 t->SetBranchAddress(
"ph2_xy", &ph2_xy);
509 t->SetBranchAddress(
"ph2_yy", &ph2_yy);
510 t->SetBranchAddress(
"ph2_yz", &ph2_yz);
511 t->SetBranchAddress(
"ph2_zz", &ph2_zz);
512 t->SetBranchAddress(
"ph2_zx", &ph2_zx);
513 if (writeHitIterMasks) {
514 t->SetBranchAddress(
"ph2_usedMask", &ph2_usedMask);
516 t->SetBranchAddress(
"ph2_simHitIdx", &ph2_simHitIdx);
518 vector<float> ph2_chargeFraction_dummy(16, 0.
f);
527 t->SetBranchAddress(
"bsp_x", &bsp_x);
528 t->SetBranchAddress(
"bsp_y", &bsp_y);
529 t->SetBranchAddress(
"bsp_z", &bsp_z);
530 t->SetBranchAddress(
"bsp_sigmax", &bsp_sigmax);
531 t->SetBranchAddress(
"bsp_sigmay", &bsp_sigmay);
532 t->SetBranchAddress(
"bsp_sigmaz", &bsp_sigmaz);
534 long long totentries =
t->GetEntries();
535 long long savedEvents = 0;
538 int outOptions = DataFile::ES_Seeds;
540 outOptions |= DataFile::ES_CmsswTracks;
541 if (writeHitIterMasks)
542 outOptions |= DataFile::ES_HitIterMasks;
543 outOptions |= DataFile::ES_BeamSpot;
549 Event EE(0, static_cast<int>(nTotalLayers));
555 for (
long long i = 0; savedEvents < maxevt &&
i < totentries &&
i < maxevt; ++
i) {
558 cout <<
"process entry i=" <<
i <<
" out of " << totentries <<
", saved so far " << savedEvents
559 <<
", with max=" << maxevt << endl;
563 cout <<
"edm event=" <<
event << endl;
565 auto&
bs =
EE.beamSpot_;
569 bs.sigmaZ = bsp_sigmaz;
570 bs.beamWidthX = bsp_sigmax;
571 bs.beamWidthY = bsp_sigmay;
575 for (
unsigned int istr = 0; istr < str_lay->size(); ++istr) {
576 if (str_chargePerCM->at(istr) < cutValueCCC)
582 auto nSims = sim_q->size();
584 cout <<
"branches not loaded" << endl;
588 std::cout << __FILE__ <<
" " << __LINE__ <<
" nSims " << nSims <<
" nSeeds " << see_q->size() <<
" nRecT " 589 << trk_q->size() << std::endl;
592 auto bestTkIdx = [&](std::vector<int>
const& shs, std::vector<float>
const& shfs,
int rhIdx,
HitType rhType) {
602 int nshs = shs.size();
603 for (
auto const sh : shs) {
605 auto tkidx = simhit_simTrkIdx->at(sh);
610 auto hpx = simhit_px->at(sh);
611 auto hpy = simhit_py->at(sh);
612 auto hpz = simhit_pz->at(sh);
613 auto hp =
sqrt(hpx * hpx + hpy * hpy + hpz * hpz);
619 auto tpx = sim_px->at(tkidx);
620 auto tpy = sim_py->at(tkidx);
621 auto tpz = sim_pz->at(tkidx);
622 auto tp =
sqrt(tpx * tpx + tpy * tpy + tpz * tpz);
629 if (maxfrac <
hp /
tp) {
641 if (nshs == 1 && ibest >= 0) {
642 auto const& srhIdxV = simhit_hitIdx->at(shbest);
643 auto const& srhTypeV = simhit_hitType->at(shbest);
645 for (
auto itype : srhTypeV) {
647 if (
HitType(itype) == rhType && srhIdxV[ih] != rhIdx) {
654 if (ibest >= 0 &&
false) {
655 std::cout <<
" best tkIdx " << ibest <<
" rh " << rhIdx <<
" for sh " << shbest <<
" out of " << shs.size()
656 <<
" hp " << hpbest <<
" chF " << hfbest <<
" tp " << tpbest <<
" process " 657 << simhit_process->at(shbest) <<
" particle " << simhit_particle->at(shbest) << std::endl;
659 std::cout <<
" sh " << simhit_x->at(shbest) <<
", " << simhit_y->at(shbest) <<
", " << simhit_z->at(shbest)
660 <<
" rh " << str_x->at(rhIdx) <<
", " << str_y->at(rhIdx) <<
", " << str_z->at(rhIdx) << std::endl;
666 vector<Track>& simTracks_ =
EE.simTracks_;
667 vector<int> simTrackIdx_(sim_q->size(), -1);
668 vector<int> seedSimIdx(see_q->size(), -1);
669 for (
unsigned int isim = 0; isim < sim_q->size(); ++isim) {
671 auto iVtx = sim_parentVtxIdx->at(isim);
672 constexpr
float largeValF = 9999.f;
673 float sim_prodx = iVtx >= 0 ? simvtx_x->at(iVtx) : largeValF;
674 float sim_prody = iVtx >= 0 ? simvtx_y->at(iVtx) : largeValF;
675 float sim_prodz = iVtx >= 0 ? simvtx_z->at(iVtx) : largeValF;
678 vector<int>
const& trkIdxV = sim_trkIdx->at(isim);
683 const int trkIdx = trkIdxV.empty() ? -1 : trkIdxV[0];
687 std::vector<int> hitlay(nTotalLayers, 0);
688 auto const&
hits = trk_hitIdx->at(trkIdx);
689 auto const& hitTypes = trk_hitType->at(trkIdx);
691 for (
auto ihit = 0
U; ihit <
nHits; ++ihit) {
692 auto ihIdx =
hits[ihit];
693 auto const ihType =
HitType(hitTypes[ihit]);
701 lnc.convertLayerNumber(pix_det->at(ipix), pix_lay->at(ipix),
useMatched, -1, pix_z->at(ipix) > 0);
702 if (cmsswlay >= 0 && cmsswlay < static_cast<int>(nTotalLayers))
710 int cmsswlay = lnc.convertLayerNumber(
711 str_det->at(istr), str_lay->at(istr),
useMatched, str_isStereo->at(istr), str_z->at(istr) > 0);
712 if (cmsswlay >= 0 && cmsswlay < static_cast<int>(nTotalLayers))
722 lnc.convertLayerNumber(glu_det->at(iglu), glu_lay->at(iglu),
useMatched, -1, glu_z->at(iglu) > 0);
723 if (cmsswlay >= 0 && cmsswlay < static_cast<int>(nTotalLayers))
732 int cmsswlay = lnc.convertLayerNumber(
733 ph2_subdet->at(istr), ph2_layer->at(istr),
useMatched, ph2_isLower->at(istr), ph2_z->at(istr) > 0);
734 if (cmsswlay >= 0 && cmsswlay < static_cast<int>(nTotalLayers))
741 throw std::logic_error(
"Track type can not be handled");
744 for (
unsigned int i = 0;
i < nTotalLayers;
i++)
752 SVector3 mom(sim_px->at(isim), sim_py->at(isim), sim_pz->at(isim));
754 err.At(0, 0) = sim_prodx * sim_prodx;
755 err.At(1, 1) = sim_prody * sim_prody;
756 err.At(2, 2) = sim_prodz * sim_prodz;
757 err.At(3, 3) = sim_px->at(isim) * sim_px->at(isim);
758 err.At(4, 4) = sim_py->at(isim) * sim_py->at(isim);
759 err.At(5, 5) = sim_pz->at(isim) * sim_pz->at(isim);
761 state.convertFromCartesianToCCS();
765 if (sim_bunchCrossing->at(isim) == 0) {
766 if (sim_event->at(isim) == 0)
769 track.setProdType(Track::ProdType::InTimePU);
771 track.setProdType(Track::ProdType::OutOfTimePU);
774 int seedIdx = trk_seedIdx->at(trkIdx);
776 seedSimIdx[seedIdx] = simTracks_.size();
778 if (cleanSimTracks) {
782 int nRecToSimHit = 0;
785 ilay = lnc.convertLayerNumber(pix_det->at(ipix), pix_lay->at(ipix),
useMatched, -1, pix_z->at(ipix) > 0);
788 int simTkIdxNt = bestTkIdx(pix_simHitIdx->at(ipix), pix_chargeFraction->at(ipix), ipix,
HitType::Pixel);
795 ilay = lnc.convertLayerNumber(
796 ph2_subdet->at(istr), ph2_layer->at(istr),
useMatched, ph2_isLower->at(istr), ph2_z->at(istr) > 0);
801 int simTkIdxNt = bestTkIdx(ph2_simHitIdx->at(istr), ph2_chargeFraction_dummy, istr,
HitType::Phase2OT);
808 if (glu_isBarrel->at(iglu) == 0)
810 int igluMono = glu_monoIdx->at(iglu);
812 bestTkIdx(str_simHitIdx->at(igluMono), str_chargeFraction->at(igluMono), igluMono,
HitType::Strip);
819 ilay = lnc.convertLayerNumber(
820 str_det->at(istr), str_lay->at(istr),
useMatched, str_isStereo->at(istr), str_z->at(istr) > 0);
821 if (
useMatched && str_isBarrel->at(istr) == 1 && str_isStereo->at(istr))
825 int simTkIdxNt = bestTkIdx(str_simHitIdx->at(istr), str_chargeFraction->at(istr), istr,
HitType::Strip);
835 simTrackIdx_[isim] = simTracks_.size();
836 simTracks_.push_back(
track);
839 if (simTracks_.empty() and not writeAllEvents)
842 vector<Track>& seedTracks_ =
EE.seedTracks_;
843 vector<vector<int>> pixHitSeedIdx(pix_lay->size());
844 vector<vector<int>> strHitSeedIdx(hasPh2hits ? 0 : str_lay->size());
845 vector<vector<int>> gluHitSeedIdx(hasPh2hits ? 0 : glu_lay->size());
846 vector<vector<int>> ph2HitSeedIdx(hasPh2hits ? ph2_layer->size() : 0);
847 for (
unsigned int is = 0; is < see_q->size(); ++is) {
853 SVector3 pos =
SVector3(see_stateTrajGlbX->at(is), see_stateTrajGlbY->at(is), see_stateTrajGlbZ->at(is));
854 SVector3 mom =
SVector3(see_stateTrajGlbPx->at(is), see_stateTrajGlbPy->at(is), see_stateTrajGlbPz->at(is));
857 err.At(0, 0) = see_stateCcov00->at(is);
858 err.At(0, 1) = see_stateCcov01->at(is);
859 err.At(0, 2) = see_stateCcov02->at(is);
860 err.At(0, 3) = see_stateCcov03->at(is);
861 err.At(0, 4) = see_stateCcov04->at(is);
862 err.At(0, 5) = see_stateCcov05->at(is);
863 err.At(1, 1) = see_stateCcov11->at(is);
864 err.At(1, 2) = see_stateCcov12->at(is);
865 err.At(1, 3) = see_stateCcov13->at(is);
866 err.At(1, 4) = see_stateCcov14->at(is);
867 err.At(1, 5) = see_stateCcov15->at(is);
868 err.At(2, 2) = see_stateCcov22->at(is);
869 err.At(2, 3) = see_stateCcov23->at(is);
870 err.At(2, 4) = see_stateCcov24->at(is);
871 err.At(2, 5) = see_stateCcov25->at(is);
872 err.At(3, 3) = see_stateCcov33->at(is);
873 err.At(3, 4) = see_stateCcov34->at(is);
874 err.At(3, 5) = see_stateCcov35->at(is);
875 err.At(4, 4) = see_stateCcov44->at(is);
876 err.At(4, 5) = see_stateCcov45->at(is);
877 err.At(5, 5) = see_stateCcov55->at(is);
879 auto const& vCov = see_stateCurvCov->at(is);
880 assert(vCov.size() == 15);
881 auto vCovP = vCov.begin();
882 for (
int i = 0;
i < 5; ++
i)
883 for (
int j = 0;
j <=
i; ++
j)
884 err.At(
i,
j) = *(vCovP++);
888 state.convertFromCartesianToCCS();
890 state.convertFromGlbCurvilinearToCCS();
892 track.setAlgorithm(isAlgo);
893 auto const& shTypes = see_hitType->at(is);
894 auto const& shIdxs = see_hitIdx->at(is);
899 for (
unsigned int ip = 0; ip < shTypes.size(); ip++) {
900 unsigned int hidx = shIdxs[ip];
901 switch (
HitType(shTypes[ip])) {
903 pixHitSeedIdx[hidx].push_back(seedTracks_.size());
907 strHitSeedIdx[hidx].push_back(seedTracks_.size());
913 int uidx = glu_monoIdx->at(hidx);
914 strHitSeedIdx[uidx].push_back(seedTracks_.size());
915 uidx = glu_stereoIdx->at(hidx);
916 strHitSeedIdx[uidx].push_back(seedTracks_.size());
918 gluHitSeedIdx[hidx].push_back(seedTracks_.size());
923 ph2HitSeedIdx[hidx].push_back(seedTracks_.size());
929 throw std::logic_error(
"Track hit type can not be handled");
932 seedTracks_.push_back(
track);
935 if (seedTracks_.empty() and not writeAllEvents)
938 vector<Track>& cmsswTracks_ =
EE.cmsswTracks_;
939 vector<vector<int>> pixHitRecIdx(pix_lay->size());
940 vector<vector<int>> strHitRecIdx(hasPh2hits ? 0 : str_lay->size());
941 vector<vector<int>> gluHitRecIdx(hasPh2hits ? 0 : glu_lay->size());
942 vector<vector<int>> ph2HitRecIdx(hasPh2hits ? ph2_layer->size() : 0);
943 for (
unsigned int ir = 0; ir < trk_q->size(); ++ir) {
950 <<
": mask " << trk_algoMask->at(ir) <<
" origAlgo " << trk_originalAlgo->at(ir) <<
" algo " 969 float pt = trk_pt->at(ir);
970 float pz = trk_pz->at(ir);
971 float p2 =
pt *
pt + pz * pz;
972 float phi = trk_phi->at(ir);
975 float dxyErr2 = trk_dxyErr->at(ir);
977 float dzErr2 = trk_dzErr->at(ir);
979 float dzErrF2 = trk_dzErr->at(ir) * (
pt * pz /
p2);
981 err.At(0, 0) = dxyErr2 * sP * sP + dzErrF2 * cP * cP;
982 err.At(0, 1) = -dxyErr2 * cP * sP + dzErrF2 * cP * sP;
983 err.At(1, 1) = dxyErr2 * cP * cP + dzErrF2 * sP * sP;
984 err.At(0, 2) = -dzErrF2 * cP *
pt / pz;
985 err.At(1, 2) = -dzErrF2 * sP *
pt / pz;
990 SVector3 pos =
SVector3(trk_refpoint_x->at(ir), trk_refpoint_y->at(ir), trk_refpoint_z->at(ir));
993 Track track(
state, trk_nChi2->at(ir), trk_seedIdx->at(ir), 0,
nullptr);
995 auto const& hTypes = trk_hitType->at(ir);
996 auto const& hIdxs = trk_hitIdx->at(ir);
997 for (
unsigned int ip = 0; ip < hTypes.size(); ip++) {
998 unsigned int hidx = hIdxs[ip];
1002 pixHitRecIdx[hidx].push_back(cmsswTracks_.size());
1007 strHitRecIdx[hidx].push_back(cmsswTracks_.size());
1012 throw std::logic_error(
"Tracks have glued hits, but matchedHit load is not configured");
1014 gluHitRecIdx[hidx].push_back(cmsswTracks_.size());
1019 ph2HitRecIdx[hidx].push_back(cmsswTracks_.size());
1025 throw std::logic_error(
"Track hit type can not be handled");
1028 cmsswTracks_.push_back(
track);
1031 vector<vector<Hit>>& layerHits_ =
EE.layerHits_;
1032 vector<vector<uint64_t>>& layerHitMasks_ =
EE.layerHitMasks_;
1033 vector<MCHitInfo>& simHitsInfo_ =
EE.simHitsInfo_;
1035 layerHits_.resize(nTotalLayers);
1036 layerHitMasks_.resize(nTotalLayers);
1037 for (
unsigned int ipix = 0; ipix < pix_lay->size(); ++ipix) {
1039 ilay = lnc.convertLayerNumber(pix_det->at(ipix), pix_lay->at(ipix),
useMatched, -1, pix_z->at(ipix) > 0);
1043 unsigned int imoduleid = tkinfo[ilay].short_id(pix_detId->at(ipix));
1045 int simTkIdxNt = bestTkIdx(pix_simHitIdx->at(ipix), pix_chargeFraction->at(ipix), ipix,
HitType::Pixel);
1046 int simTkIdx = simTkIdxNt >= 0 ? simTrackIdx_[simTkIdxNt] : -1;
1048 SVector3 pos(pix_x->at(ipix), pix_y->at(ipix), pix_z->at(ipix));
1050 err.At(0, 0) = pix_xx->at(ipix);
1051 err.At(1, 1) = pix_yy->at(ipix);
1052 err.At(2, 2) = pix_zz->at(ipix);
1053 err.At(0, 1) = pix_xy->at(ipix);
1054 err.At(0, 2) = pix_zx->at(ipix);
1055 err.At(1, 2) = pix_yz->at(ipix);
1056 if (simTkIdx >= 0) {
1057 simTracks_[simTkIdx].addHitIdx(layerHits_[ilay].
size(), ilay, 0);
1059 for (
unsigned int is = 0; is < pixHitSeedIdx[ipix].size(); is++) {
1061 seedTracks_[pixHitSeedIdx[ipix][is]].addHitIdx(layerHits_[ilay].
size(), ilay, 0);
1063 for (
unsigned int ir = 0; ir < pixHitRecIdx[ipix].size(); ir++) {
1065 cmsswTracks_[pixHitRecIdx[ipix][ir]].addHitIdx(layerHits_[ilay].
size(), ilay, 0);
1068 hit.setupAsPixel(imoduleid, pix_csize_row->at(ipix), pix_csize_col->at(ipix));
1069 layerHits_[ilay].push_back(
hit);
1070 if (writeHitIterMasks)
1071 layerHitMasks_[ilay].push_back(pix_usedMask->at(ipix));
1072 MCHitInfo hitInfo(simTkIdx, ilay, layerHits_[ilay].
size() - 1, totHits);
1073 simHitsInfo_.push_back(hitInfo);
1078 vector<int> ph2Idx(ph2_layer->size());
1079 for (
unsigned int iph2 = 0; iph2 < ph2_layer->size(); ++iph2) {
1081 ilay = lnc.convertLayerNumber(
1082 ph2_subdet->at(iph2), ph2_layer->at(iph2),
useMatched, ph2_isLower->at(iph2), ph2_z->at(iph2) > 0);
1088 unsigned int imoduleid = tkinfo[ilay].short_id(ph2_detId->at(iph2));
1090 int simTkIdxNt = bestTkIdx(ph2_simHitIdx->at(iph2), ph2_chargeFraction_dummy, iph2,
HitType::Phase2OT);
1091 int simTkIdx = simTkIdxNt >= 0 ? simTrackIdx_[simTkIdxNt] : -1;
1093 SVector3 pos(ph2_x->at(iph2), ph2_y->at(iph2), ph2_z->at(iph2));
1095 err.At(0, 0) = ph2_xx->at(iph2);
1096 err.At(1, 1) = ph2_yy->at(iph2);
1097 err.At(2, 2) = ph2_zz->at(iph2);
1098 err.At(0, 1) = ph2_xy->at(iph2);
1099 err.At(0, 2) = ph2_zx->at(iph2);
1100 err.At(1, 2) = ph2_yz->at(iph2);
1101 if (simTkIdx >= 0) {
1102 simTracks_[simTkIdx].addHitIdx(layerHits_[ilay].
size(), ilay, 0);
1104 for (
unsigned int ir = 0; ir < ph2HitSeedIdx[iph2].size(); ir++) {
1106 seedTracks_[ph2HitSeedIdx[iph2][ir]].addHitIdx(layerHits_[ilay].
size(), ilay, 0);
1108 for (
unsigned int ir = 0; ir < ph2HitRecIdx[iph2].size(); ir++) {
1110 cmsswTracks_[ph2HitRecIdx[iph2][ir]].addHitIdx(layerHits_[ilay].
size(), ilay, 0);
1113 hit.setupAsStrip(imoduleid, 0, 1);
1114 layerHits_[ilay].push_back(
hit);
1115 if (writeHitIterMasks)
1116 layerHitMasks_[ilay].push_back(ph2_usedMask->at(iph2));
1117 MCHitInfo hitInfo(simTkIdx, ilay, layerHits_[ilay].
size() - 1, totHits);
1118 simHitsInfo_.push_back(hitInfo);
1123 for (
unsigned int iglu = 0; iglu < glu_lay->size(); ++iglu) {
1124 if (glu_isBarrel->at(iglu) == 0)
1126 int igluMono = glu_monoIdx->at(iglu);
1128 bestTkIdx(str_simHitIdx->at(igluMono), str_chargeFraction->at(igluMono), igluMono,
HitType::Strip);
1129 int simTkIdx = simTkIdxNt >= 0 ? simTrackIdx_[simTkIdxNt] : -1;
1131 int ilay = lnc.convertLayerNumber(glu_det->at(iglu), glu_lay->at(iglu),
useMatched, -1, glu_z->at(iglu) > 0);
1133 SVector3 pos(glu_x->at(iglu), glu_y->at(iglu), glu_z->at(iglu));
1135 err.At(0, 0) = glu_xx->at(iglu);
1136 err.At(1, 1) = glu_yy->at(iglu);
1137 err.At(2, 2) = glu_zz->at(iglu);
1138 err.At(0, 1) = glu_xy->at(iglu);
1139 err.At(0, 2) = glu_zx->at(iglu);
1140 err.At(1, 2) = glu_yz->at(iglu);
1141 if (simTkIdx >= 0) {
1142 simTracks_[simTkIdx].addHitIdx(layerHits_[ilay].
size(), ilay, 0);
1144 for (
unsigned int ir = 0; ir < gluHitSeedIdx[iglu].size(); ir++) {
1146 seedTracks_[gluHitSeedIdx[iglu][ir]].addHitIdx(
1147 layerHits_[ilay].
size(), ilay, 0);
1149 for (
unsigned int ir = 0; ir < gluHitRecIdx[iglu].size(); ir++) {
1151 cmsswTracks_[gluHitRecIdx[iglu][ir]].addHitIdx(
1152 layerHits_[ilay].
size(), ilay, 0);
1157 assert(
false &&
"Implement module-ids, cluster adc and spans for matched hits!");
1160 layerHits_[ilay].push_back(
hit);
1161 MCHitInfo hitInfo(simTkIdx, ilay, layerHits_[ilay].
size() - 1, totHits);
1162 simHitsInfo_.push_back(hitInfo);
1168 strIdx.resize(str_lay->size());
1169 for (
unsigned int istr = 0; istr < str_lay->size(); ++istr) {
1171 ilay = lnc.convertLayerNumber(
1172 str_det->at(istr), str_lay->at(istr),
useMatched, str_isStereo->at(istr), str_z->at(istr) > 0);
1173 if (
useMatched && str_isBarrel->at(istr) == 1 && str_isStereo->at(istr))
1178 unsigned int imoduleid = tkinfo[ilay].short_id(str_detId->at(istr));
1180 int simTkIdxNt = bestTkIdx(str_simHitIdx->at(istr), str_chargeFraction->at(istr), istr,
HitType::Strip);
1181 int simTkIdx = simTkIdxNt >= 0 ? simTrackIdx_[simTkIdxNt] : -1;
1183 bool passCCC = applyCCC ? (str_chargePerCM->at(istr) > cutValueCCC) :
true;
1186 SVector3 pos(str_x->at(istr), str_y->at(istr), str_z->at(istr));
1188 err.At(0, 0) = str_xx->at(istr);
1189 err.At(1, 1) = str_yy->at(istr);
1190 err.At(2, 2) = str_zz->at(istr);
1191 err.At(0, 1) = str_xy->at(istr);
1192 err.At(0, 2) = str_zx->at(istr);
1193 err.At(1, 2) = str_yz->at(istr);
1194 if (simTkIdx >= 0) {
1196 simTracks_[simTkIdx].addHitIdx(layerHits_[ilay].
size(), ilay, 0);
1198 simTracks_[simTkIdx].addHitIdx(-9, ilay, 0);
1200 for (
unsigned int ir = 0; ir < strHitSeedIdx[istr].size(); ir++) {
1203 seedTracks_[strHitSeedIdx[istr][ir]].addHitIdx(
1204 layerHits_[ilay].
size(), ilay, 0);
1206 seedTracks_[strHitSeedIdx[istr][ir]].addHitIdx(-9, ilay, 0);
1208 for (
unsigned int ir = 0; ir < strHitRecIdx[istr].size(); ir++) {
1211 cmsswTracks_[strHitRecIdx[istr][ir]].addHitIdx(
1212 layerHits_[ilay].
size(), ilay, 0);
1214 cmsswTracks_[strHitRecIdx[istr][ir]].addHitIdx(-9, ilay, 0);
1218 hit.setupAsStrip(imoduleid, str_chargePerCM->at(istr), str_csize->at(istr));
1219 layerHits_[ilay].push_back(
hit);
1220 if (writeHitIterMasks)
1221 layerHitMasks_[ilay].push_back(str_usedMask->at(istr));
1222 MCHitInfo hitInfo(simTkIdx, ilay, layerHits_[ilay].
size() - 1, totHits);
1223 simHitsInfo_.push_back(hitInfo);
1230 nstot += seedTracks_.size();
1231 for (
unsigned int il = 0; il < layerHits_.size(); ++il) {
1232 int nh = layerHits_[il].size();
1237 int nt = simTracks_.size();
1239 int nl = layerHits_.size();
1241 int nm = simHitsInfo_.size();
1243 int ns = seedTracks_.size();
1245 int nr = cmsswTracks_.size();
1247 printf(
"number of simTracks %i\n",
nt);
1248 printf(
"number of layerHits %i\n", nl);
1249 printf(
"number of simHitsInfo %i\n", nm);
1250 printf(
"number of seedTracks %i\n", ns);
1251 printf(
"number of recTracks %i\n",
nr);
1255 for (
int il = 0; il < nl; ++il) {
1256 int nh = layerHits_[il].size();
1257 for (
int ih = 0; ih <
nh; ++ih) {
1258 printf(
"lay=%i idx=%i mcid=%i x=(%6.3f, %6.3f, %6.3f) r=%6.3f mask=0x%lx\n",
1261 layerHits_[il][ih].mcHitID(),
1262 layerHits_[il][ih].
x(),
1263 layerHits_[il][ih].
y(),
1264 layerHits_[il][ih].
z(),
1265 sqrt(
pow(layerHits_[il][ih].
x(), 2) +
pow(layerHits_[il][ih].
y(), 2)),
1266 writeHitIterMasks ? layerHitMasks_[il][ih] : 0);
1270 for (
int i = 0;
i <
nt; ++
i) {
1273 "sim track id=%i q=%2i p=(%6.3f, %6.3f, %6.3f) x=(%6.3f, %6.3f, %6.3f) pT=%7.4f nTotal=%i nFound=%i \n",
1283 simTracks_[
i].nTotalHits(),
1284 simTracks_[
i].nFoundHits());
1285 int nh = simTracks_[
i].nTotalHits();
1286 for (
int ih = 0; ih <
nh; ++ih) {
1287 int hidx = simTracks_[
i].getHitIdx(ih);
1288 int hlay = simTracks_[
i].getHitLyr(ih);
1289 float hx = layerHits_[hlay][hidx].x();
1290 float hy = layerHits_[hlay][hidx].y();
1291 float hz = layerHits_[hlay][hidx].z();
1292 printf(
"track #%4i hit #%2i idx=%4i lay=%2i x=(% 8.3f, % 8.3f, % 8.3f) r=%8.3f\n",
1300 sqrt(hx * hx + hy * hy));
1304 for (
int i = 0;
i < ns; ++
i) {
1305 printf(
"seed id=%i label=%i algo=%i q=%2i pT=%6.3f p=(%6.3f, %6.3f, %6.3f) x=(%6.3f, %6.3f, %6.3f)\n",
1310 seedTracks_[
i].
pT(),
1311 seedTracks_[
i].
px(),
1312 seedTracks_[
i].
py(),
1313 seedTracks_[
i].pz(),
1316 seedTracks_[
i].
z());
1317 int nh = seedTracks_[
i].nTotalHits();
1318 for (
int ih = 0; ih <
nh; ++ih)
1319 printf(
"seed #%i hit #%i idx=%i\n",
i, ih, seedTracks_[
i].getHitIdx(ih));
1322 if (writeRecTracks) {
1323 for (
int i = 0;
i <
nr; ++
i) {
1324 float spt =
sqrt(
pow(cmsswTracks_[
i].
px(), 2) +
pow(cmsswTracks_[
i].
py(), 2));
1326 "rec track id=%i label=%i algo=%i chi2=%6.3f q=%2i p=(%6.3f, %6.3f, %6.3f) x=(%6.3f, %6.3f, %6.3f) " 1327 "pT=%7.4f nTotal=%i nFound=%i \n",
1331 cmsswTracks_[
i].
chi2(),
1333 cmsswTracks_[
i].
px(),
1334 cmsswTracks_[
i].
py(),
1335 cmsswTracks_[
i].pz(),
1336 cmsswTracks_[
i].
x(),
1337 cmsswTracks_[
i].
y(),
1338 cmsswTracks_[
i].
z(),
1340 cmsswTracks_[
i].nTotalHits(),
1341 cmsswTracks_[
i].nFoundHits());
1342 int nh = cmsswTracks_[
i].nTotalHits();
1343 for (
int ih = 0; ih <
nh; ++ih) {
1344 int hidx = cmsswTracks_[
i].getHitIdx(ih);
1345 int hlay = cmsswTracks_[
i].getHitLyr(ih);
1346 float hx = layerHits_[hlay][hidx].x();
1347 float hy = layerHits_[hlay][hidx].y();
1348 float hz = layerHits_[hlay][hidx].z();
1349 printf(
"track #%4i hit #%2i idx=%4i lay=%2i x=(% 8.3f, % 8.3f, % 8.3f) r=%8.3f\n",
1357 sqrt(hx * hx + hy * hy));
1364 EE.write_out(data_file);
1367 printf(
"end of event %lli\n", savedEvents);
1371 printf(
"\nSaved %lli events\n\n", savedEvents);
1373 printf(
"Average number of seeds per event %f\n",
float(nstot) /
float(savedEvents));
1374 for (
unsigned int il = 0; il < nhitstot.size(); ++il)
1375 printf(
"Average number of hits in layer %3i = %7.2f\n",
1377 float(nhitstot[il]) /
float(savedEvents));
1380 printf(
"Out of %i hits, %i failed the cut", numTotalStr, numFailCCC);
1384 printf(
"\n\n================================================================\n");
1385 printf(
"=== Max module id for %u layers\n", nTotalLayers);
1386 printf(
"================================================================\n");
1387 for (
unsigned int ii = 0;
ii < nTotalLayers; ++
ii) {
1388 printf(
"Layer%2d : %d\n",
ii, tkinfo[
ii].n_modules());
void CloseWrite(int n_written)
void next_arg_or_die(lStr_t &args, lStr_i &i)
ROOT::Math::SMatrix< float, 6, 6, ROOT::Math::MatRepSym< float, 6 > > SMatrixSym66
void read_bin_file(const std::string &fname)
bool next_arg_option(lStr_t &args, lStr_i &i)
Sin< T >::type sin(const T &t)
ROOT::Math::SVector< float, 3 > SVector3
constexpr bool useMatched
constexpr int cleanSimTrack_minSimHits
The Signals That Services Can Subscribe To This is based on ActivityRegistry and is current per Services can connect to the signals distributed by the ActivityRegistry in order to monitor the activity of the application Each possible callback has some defined which we here list in angle e< void, edm::EventID const &, edm::Timestamp const & > We also list in braces which AR_WATCH_USING_METHOD_ is used for those or
Cos< T >::type cos(const T &t)
TrackBase::TrackAlgorithm TrackAlgorithm
void printHelp(const char *av0)
void openWrite(const std::string &fname, int n_layers, int n_ev, int extra_sections=0)
constexpr int cleanSimTrack_minRecHits
std::list< std::string > lStr_t
TupleMultiplicity< TrackerTraits > const *__restrict__ uint32_t nHits
ROOT::Math::SVector< double, 3 > SVector3
ROOT::Math::SMatrix< float, 3, 3, ROOT::Math::MatRepSym< float, 3 > > SMatrixSym33