Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-03-17 11:22:28

0001 #include "RecoTracker/MkFitCMS/standalone/buildtestMPlex.h"
0002 #include "RecoTracker/MkFitCore/standalone/ConfigStandalone.h"
0003 #include "RecoTracker/MkFitCore/src/Matrix.h"
0004 #include "RecoTracker/MkFitCore/interface/MkBuilder.h"
0005 #include "RecoTracker/MkFitCMS/interface/MkStdSeqs.h"
0006 #include "RecoTracker/MkFitCMS/standalone/MkStandaloneSeqs.h"
0007 
0008 #include "oneapi/tbb/parallel_for.h"
0009 
0010 #include <memory>
0011 
0012 namespace mkfit {
0013 
0014   inline bool sortByHitsChi2(const std::pair<Track, TrackState> &cand1, const std::pair<Track, TrackState> &cand2) {
0015     if (cand1.first.nFoundHits() == cand2.first.nFoundHits())
0016       return cand1.first.chi2() < cand2.first.chi2();
0017 
0018     return cand1.first.nFoundHits() > cand2.first.nFoundHits();
0019   }
0020 
0021   inline bool sortByPhi(const Hit &hit1, const Hit &hit2) {
0022     return std::atan2(hit1.y(), hit1.x()) < std::atan2(hit2.y(), hit2.x());
0023   }
0024 
0025   inline bool sortByEta(const Hit &hit1, const Hit &hit2) { return hit1.eta() < hit2.eta(); }
0026 
0027   inline bool sortTracksByEta(const Track &track1, const Track &track2) { return track1.momEta() < track2.momEta(); }
0028 
0029   inline bool sortTracksByPhi(const Track &track1, const Track &track2) { return track1.momPhi() < track2.momPhi(); }
0030 
0031   struct sortTracksByPhiStruct {
0032     const std::vector<std::vector<Track>> &m_track_candidates;
0033 
0034     sortTracksByPhiStruct(std::vector<std::vector<Track>> *track_candidates) : m_track_candidates(*track_candidates) {}
0035 
0036     bool operator()(const std::pair<int, int> &track1, const std::pair<int, int> &track2) {
0037       return m_track_candidates[track1.first][track1.second].posPhi() <
0038              m_track_candidates[track2.first][track2.second].posPhi();
0039     }
0040   };
0041 
0042   // within a layer with a "reasonable" geometry, ordering by Z is the same as eta
0043   inline bool sortByZ(const Hit &hit1, const Hit &hit2) { return hit1.z() < hit2.z(); }
0044 
0045   //==============================================================================
0046   // NaN and Silly track parameter check
0047   //==============================================================================
0048 
0049   namespace {
0050 
0051     int check_nan_n_silly(TrackVec &tracks, const char *prefix) {
0052       int count = 0;
0053       for (auto &t : tracks) {
0054         if (t.hasSillyValues(Const::nan_n_silly_print_bad_cands_bkfit, false, prefix)) {
0055           ++count;
0056         }
0057       }
0058       return count;
0059     }
0060 
0061     void check_nan_n_silly_candidates(Event &ev) {
0062       // MIMI -- nan_n_silly_per_layer_count is in MkBuilder, could be in MkJob.
0063       // if (Const::nan_n_silly_check_cands_every_layer)
0064       // {
0065       //   int sc = (int) ev.nan_n_silly_per_layer_count_;
0066       //   if (sc > 0)
0067       //     printf("Nan'n'Silly: Number of silly candidates over all layers = %d\n", sc);
0068       // }
0069       if (Const::nan_n_silly_check_cands_pre_bkfit) {
0070         int sc = check_nan_n_silly(ev.candidateTracks_, "Pre-bkfit silly check");
0071         if (sc > 0)
0072           printf("Nan'n'Silly: Number of silly pre-bkfit candidates = %d\n", sc);
0073       }
0074     }
0075 
0076     void check_nan_n_silly_bkfit(Event &ev) {
0077       if (Const::nan_n_silly_check_cands_post_bkfit) {
0078         int sc = check_nan_n_silly(ev.fitTracks_, "Post-bkfit silly check");
0079         if (sc > 0)
0080           printf("Nan'n'Silly: Number of silly post-bkfit candidates = %d\n", sc);
0081       }
0082     }
0083 
0084   }  // namespace
0085 
0086   //==============================================================================
0087   // runBuildTestPlexDumbCMSSW
0088   //==============================================================================
0089 
0090   void runBuildingTestPlexDumbCMSSW(Event &ev, const EventOfHits &eoh, MkBuilder &builder) {
0091     const IterationConfig &itconf = Config::ItrInfo[0];
0092 
0093     MkJob job({Config::TrkInfo, itconf, eoh, eoh.refBeamSpot()});
0094 
0095     builder.begin_event(&job, &ev, __func__);
0096 
0097     if (Config::sim_val_for_cmssw) {
0098       StdSeq::root_val_dumb_cmssw(&ev);
0099     }
0100 
0101     builder.end_event();
0102   }
0103 
0104   //==============================================================================
0105   // runBuildTestPlexBestHit
0106   //==============================================================================
0107 
0108   double runBuildingTestPlexBestHit(Event &ev, const EventOfHits &eoh, MkBuilder &builder) {
0109     const IterationConfig &itconf = Config::ItrInfo[0];
0110 
0111     const bool validation_on = (Config::sim_val || Config::quality_val);
0112 
0113     if (validation_on) {
0114       TrackVec seeds1;
0115 
0116       unsigned int algorithms[] = {4};  //only initialStep
0117 
0118       for (auto const &s : ev.seedTracks_) {
0119         //keep seeds form the first iteration for processing
0120         if (std::find(algorithms, algorithms + 1, s.algoint()) != algorithms + 1)
0121           seeds1.push_back(s);
0122       }
0123       ev.seedTracks_.swap(seeds1);  //necessary for the validation - PrepareSeeds
0124       ev.relabel_bad_seedtracks();  //necessary for the validation - PrepareSeeds
0125     }
0126 
0127     IterationMaskIfc mask_ifc;
0128 
0129     // To disable hit-masks, pass nullptr in place of &mask_ifc to MkJob ctor
0130     // and optionally comment out ev.fill_hitmask_bool_vectors() call.
0131 
0132     ev.fill_hitmask_bool_vectors(itconf.m_track_algorithm, mask_ifc.m_mask_vector);
0133 
0134     MkJob job({Config::TrkInfo, itconf, eoh, eoh.refBeamSpot(), &mask_ifc});
0135 
0136     builder.begin_event(&job, &ev, __func__);
0137 
0138     bool seeds_sorted = false;
0139     // CCCC builder.PrepareSeeds();
0140 
0141     // EventOfCandidates event_of_cands;
0142     builder.find_tracks_load_seeds_BH(ev.seedTracks_, seeds_sorted);
0143 
0144 #ifdef USE_VTUNE_PAUSE
0145     __SSC_MARK(0x111);  // use this to resume Intel SDE at the same point
0146     __itt_resume();
0147 #endif
0148 
0149     double time = dtime();
0150 
0151     builder.findTracksBestHit();
0152 
0153     time = dtime() - time;
0154 
0155 #ifdef USE_VTUNE_PAUSE
0156     __itt_pause();
0157     __SSC_MARK(0x222);  // use this to pause Intel SDE at the same point
0158 #endif
0159 
0160     // Hack, get the tracks out.
0161     ev.candidateTracks_ = builder.ref_tracks();
0162 
0163     // For best hit, the candidateTracks_ vector is the direct input to the backward fit so only need to do clean_duplicates once
0164     if (Config::quality_val || Config::sim_val || Config::cmssw_val) {
0165       //Mark tracks as duplicates; if within CMSSW, remove duplicate tracks before backward fit
0166       // CCCC if (Config::removeDuplicates) {
0167       // CCCC   StdSeq::clean_duplicates(ev.candidateTracks_);
0168       // CCCC }
0169     }
0170 
0171     job.switch_to_backward();
0172 
0173     // now do backwards fit... do we want to time this section?
0174     if (Config::backwardFit) {
0175       builder.backwardFitBH();
0176       ev.fitTracks_ = builder.ref_tracks();
0177     }
0178 
0179     if (Config::quality_val) {
0180       StdSeq::Quality qval;
0181       qval.quality_val(&ev);
0182     } else if (Config::sim_val || Config::cmssw_val) {
0183       StdSeq::root_val(&ev);
0184     }
0185 
0186     builder.end_event();
0187 
0188     // ev.print_tracks(ev.candidateTracks_, true);
0189 
0190     return time;
0191   }
0192 
0193   //==============================================================================
0194   // runBuildTestPlex Combinatorial: Standard TBB
0195   //==============================================================================
0196 
0197   double runBuildingTestPlexStandard(Event &ev, const EventOfHits &eoh, MkBuilder &builder) {
0198     const IterationConfig &itconf = Config::ItrInfo[0];
0199 
0200     const bool validation_on = (Config::sim_val || Config::quality_val);
0201 
0202     if (validation_on) {
0203       TrackVec seeds1;
0204 
0205       unsigned int algorithms[] = {4};  //only initialStep
0206 
0207       for (auto const &s : ev.seedTracks_) {
0208         //keep seeds form the first iteration for processing
0209         if (std::find(algorithms, algorithms + 1, s.algoint()) != algorithms + 1)
0210           seeds1.push_back(s);
0211       }
0212       ev.seedTracks_.swap(seeds1);  //necessary for the validation - PrepareSeeds
0213       ev.relabel_bad_seedtracks();  //necessary for the validation - PrepareSeeds
0214     }
0215 
0216     IterationMaskIfc mask_ifc;
0217 
0218     // To disable hit-masks, pass nullptr in place of &mask_ifc to MkJob ctor
0219     // and optionally comment out ev.fill_hitmask_bool_vectors() call.
0220 
0221     ev.fill_hitmask_bool_vectors(itconf.m_track_algorithm, mask_ifc.m_mask_vector);
0222 
0223     MkJob job({Config::TrkInfo, itconf, eoh, eoh.refBeamSpot(), &mask_ifc});
0224 
0225     builder.begin_event(&job, &ev, __func__);
0226 
0227     bool seeds_sorted = false;
0228     // CCCC builder.PrepareSeeds();
0229 
0230     builder.find_tracks_load_seeds(ev.seedTracks_, seeds_sorted);
0231 
0232 #ifdef USE_VTUNE_PAUSE
0233     __SSC_MARK(0x111);  // use this to resume Intel SDE at the same point
0234     __itt_resume();
0235 #endif
0236 
0237     double time = dtime();
0238 
0239     builder.findTracksStandard();
0240 
0241     time = dtime() - time;
0242 
0243 #ifdef USE_VTUNE_PAUSE
0244     __itt_pause();
0245     __SSC_MARK(0x222);  // use this to pause Intel SDE at the same point
0246 #endif
0247 
0248     check_nan_n_silly_candidates(ev);
0249 
0250     // first store candidate tracks
0251     builder.export_best_comb_cands(ev.candidateTracks_);
0252 
0253     job.switch_to_backward();
0254 
0255     // now do backwards fit... do we want to time this section?
0256     if (Config::backwardFit) {
0257       // Using the TrackVec version until we home in on THE backward fit etc.
0258       // builder.backwardFit();
0259       builder.select_best_comb_cands();
0260       builder.backwardFitBH();
0261       ev.fitTracks_ = builder.ref_tracks();
0262 
0263       check_nan_n_silly_bkfit(ev);
0264     }
0265 
0266     // CCCC StdSeq::handle_duplicates(&ev);
0267 
0268     if (Config::quality_val) {
0269       StdSeq::Quality qval;
0270       qval.quality_val(&ev);
0271     } else if (Config::sim_val || Config::cmssw_val) {
0272       StdSeq::root_val(&ev);
0273     }
0274 
0275     builder.end_event();
0276 
0277     // ev.print_tracks(ev.candidateTracks_, true);
0278 
0279     return time;
0280   }
0281 
0282   //==============================================================================
0283   // runBuildTestPlex Combinatorial: CloneEngine TBB
0284   //==============================================================================
0285 
0286   double runBuildingTestPlexCloneEngine(Event &ev, const EventOfHits &eoh, MkBuilder &builder) {
0287     const IterationConfig &itconf = Config::ItrInfo[0];
0288 
0289     const bool validation_on = (Config::sim_val || Config::quality_val);
0290 
0291     if (validation_on) {
0292       TrackVec seeds1;
0293 
0294       unsigned int algorithms[] = {4};  //only initialStep
0295 
0296       for (auto const &s : ev.seedTracks_) {
0297         //keep seeds form the first iteration for processing
0298         if (std::find(algorithms, algorithms + 1, s.algoint()) != algorithms + 1)
0299           seeds1.push_back(s);
0300       }
0301       ev.seedTracks_.swap(seeds1);  //necessary for the validation - PrepareSeeds
0302       ev.relabel_bad_seedtracks();  //necessary for the validation - PrepareSeeds
0303     }
0304 
0305     IterationMaskIfc mask_ifc;
0306 
0307     // To disable hit-masks, pass nullptr in place of &mask_ifc to MkJob ctor
0308     // and optionally comment out ev.fill_hitmask_bool_vectors() call.
0309 
0310     ev.fill_hitmask_bool_vectors(itconf.m_track_algorithm, mask_ifc.m_mask_vector);
0311 
0312     MkJob job({Config::TrkInfo, itconf, eoh, eoh.refBeamSpot(), &mask_ifc});
0313 
0314     builder.begin_event(&job, &ev, __func__);
0315 
0316     bool seeds_sorted = false;
0317     // CCCC builder.PrepareSeeds();
0318 
0319     builder.find_tracks_load_seeds(ev.seedTracks_, seeds_sorted);
0320 
0321 #ifdef USE_VTUNE_PAUSE
0322     __SSC_MARK(0x111);  // use this to resume Intel SDE at the same point
0323     __itt_resume();
0324 #endif
0325 
0326     double time = dtime();
0327 
0328     builder.findTracksCloneEngine();
0329 
0330     time = dtime() - time;
0331 
0332 #ifdef USE_VTUNE_PAUSE
0333     __itt_pause();
0334     __SSC_MARK(0x222);  // use this to pause Intel SDE at the same point
0335 #endif
0336 
0337     check_nan_n_silly_candidates(ev);
0338 
0339     // first store candidate tracks - needed for BH backward fit and root_validation
0340     builder.export_best_comb_cands(ev.candidateTracks_);
0341 
0342     job.switch_to_backward();
0343 
0344     // now do backwards fit... do we want to time this section?
0345     if (Config::backwardFit) {
0346       // a) TrackVec version:
0347       builder.select_best_comb_cands();
0348       builder.backwardFitBH();
0349       ev.fitTracks_ = builder.ref_tracks();
0350 
0351       // b) Version that runs on CombCand / TrackCand
0352       // builder.backwardFit();
0353       // builder.quality_store_tracks(ev.fitTracks_);
0354 
0355       check_nan_n_silly_bkfit(ev);
0356     }
0357 
0358     // CCCC StdSeq::handle_duplicates(&ev);
0359 
0360     // validation section
0361     if (Config::quality_val) {
0362       StdSeq::Quality qval;
0363       qval.quality_val(&ev);
0364     } else if (Config::sim_val || Config::cmssw_val) {
0365       StdSeq::root_val(&ev);
0366     }
0367 
0368     builder.end_event();
0369 
0370     // ev.print_tracks(ev.candidateTracks_, true);
0371 
0372     return time;
0373   }
0374 
0375   //==============================================================================
0376   // runBtpCe_MultiIter
0377   //
0378   // Prototype for running multiple iterations, sequentially, using the same builder.
0379   // For cmmsw seeds
0380   //
0381   // There is, in general, a mess in how tracks are processed, marked, or copied out
0382   // in various validation scenarios and export flags.
0383   //
0384   // In particular, MkBuilder::PrepareSeeds does a lot of things to whole / complete
0385   // event,seedTracks_ -- probably this would need to be split into common / and
0386   // per-iteration part.
0387   // - MkBuilder::prep_*** functions also mostly do not belong there (prep_sim is called from
0388   //   PrepareSeeds() for cmssw seeds).
0389   //
0390   // At this point we need to think about what should happen to Event before all the iterations and
0391   // after all the iterations ... from the Validation perspective.
0392   // And if we care about doing too muich work for seeds that will never get processed.
0393   //==============================================================================
0394 
0395   namespace {
0396     constexpr unsigned int algorithms[] = {4, 22, 23, 5, 24, 7, 8, 9, 10, 6};  //9 iterations
0397   }
0398 
0399   std::vector<double> runBtpCe_MultiIter(Event &ev, const EventOfHits &eoh, MkBuilder &builder, int n) {
0400     std::vector<double> timevec;
0401     if (n <= 0)
0402       return timevec;
0403     timevec.resize(n + 1, 0.0);
0404 
0405     const bool validation_on = (Config::sim_val || Config::quality_val);
0406 
0407     TrackVec seeds_used;
0408     TrackVec seeds1;
0409 
0410     if (validation_on) {
0411       for (auto const &s : ev.seedTracks_) {
0412         //keep seeds form the first n iterations for processing
0413         if (std::find(algorithms, algorithms + n, s.algoint()) != algorithms + n)
0414           seeds1.push_back(s);
0415       }
0416       ev.seedTracks_.swap(seeds1);  //necessary for the validation - PrepareSeeds
0417       ev.relabel_bad_seedtracks();  //necessary for the validation - PrepareSeeds
0418     }
0419 
0420     IterationMaskIfc mask_ifc;
0421     TrackVec seeds;
0422     TrackVec tmp_tvec;
0423 
0424     for (int it = 0; it <= n - 1; ++it) {
0425       const IterationConfig &itconf = Config::ItrInfo[it];
0426 
0427       // To disable hit-masks, pass nullptr in place of &mask_ifc to MkJob ctor
0428       // and optionally comment out ev.fill_hitmask_bool_vectors() call.
0429 
0430       ev.fill_hitmask_bool_vectors(itconf.m_track_algorithm, mask_ifc.m_mask_vector);
0431 
0432       MkJob job({Config::TrkInfo, itconf, eoh, eoh.refBeamSpot(), &mask_ifc});
0433 
0434       builder.begin_event(&job, &ev, __func__);
0435 
0436       {  // We could partition seeds once, store beg, end for each iteration in a map or vector.
0437         seeds.clear();
0438         int nc = 0;
0439         for (auto &s : ev.seedTracks_) {
0440           if (s.algoint() == itconf.m_track_algorithm) {
0441             if (itconf.m_requires_seed_hit_sorting) {
0442               s.sortHitsByLayer();
0443             }
0444             seeds.push_back(s);
0445             ++nc;
0446           } else if (nc > 0)
0447             break;
0448         }
0449       }
0450 
0451       bool do_seed_clean = bool(itconf.m_seed_cleaner);
0452 
0453       if (do_seed_clean)
0454         itconf.m_seed_cleaner(seeds, itconf, eoh.refBeamSpot());
0455 
0456       builder.seed_post_cleaning(seeds);
0457 
0458       // Add protection in case no seeds are found for iteration
0459       if (seeds.size() <= 0)
0460         continue;
0461 
0462       builder.find_tracks_load_seeds(seeds, do_seed_clean);
0463 
0464       double time = dtime();
0465 
0466       builder.findTracksCloneEngine();
0467 
0468       timevec[it] = dtime() - time;
0469       timevec[n] += timevec[it];
0470 
0471       // Print min and max size of hots vectors of CombCands.
0472       // builder.find_min_max_hots_size();
0473 
0474       if (validation_on)
0475         seeds_used.insert(seeds_used.end(), seeds.begin(), seeds.end());  //cleaned seeds need to be stored somehow
0476 
0477       // Pre backward-fit filtering.
0478       // Note -- slightly different logic than run_OneIteration as we always do nan filters for
0479       // export for validation.
0480       filter_candidates_func pre_filter;
0481       if (itconf.m_pre_bkfit_filter)
0482         pre_filter = [&](const TrackCand &tc, const MkJob &jb) -> bool {
0483           return itconf.m_pre_bkfit_filter(tc, jb) && StdSeq::qfilter_nan_n_silly<TrackCand>(tc, jb);
0484         };
0485       else
0486         pre_filter = StdSeq::qfilter_nan_n_silly<TrackCand>;
0487       // pre_filter is always at least doing nan_n_silly filter.
0488       builder.filter_comb_cands(pre_filter, true);
0489 
0490       builder.select_best_comb_cands();
0491 
0492       {
0493         builder.export_tracks(tmp_tvec);
0494         if (itconf.m_duplicate_cleaner)
0495           itconf.m_duplicate_cleaner(builder.ref_tracks_nc(), itconf);
0496         ev.candidateTracks_.reserve(ev.candidateTracks_.size() + tmp_tvec.size());
0497         for (auto &&t : tmp_tvec)
0498           ev.candidateTracks_.emplace_back(std::move(t));
0499         tmp_tvec.clear();
0500       }
0501 
0502       job.switch_to_backward();
0503 
0504       // now do backwards fit... do we want to time this section?
0505       if (Config::backwardFit) {
0506         // a) TrackVec version:
0507         // builder.backwardFitBH();
0508 
0509         // b) Version that runs on CombCand / TrackCand
0510         const bool do_backward_search = Config::backwardSearch && itconf.m_backward_search;
0511 
0512         // We copy seed-hits into Candidates ... now we have to remove them so backward fit stops
0513         // before reaching seeding region. Ideally, we wouldn't add them in the first place but
0514         // if we want to export full tracks above we need to hold on to them (alternatively, we could
0515         // have a pointer to seed track in CombCandidate and copy them from there).
0516         if (do_backward_search)
0517           builder.compactifyHitStorageForBestCand(itconf.m_backward_drop_seed_hits, itconf.m_backward_fit_min_hits);
0518 
0519         builder.backwardFit();
0520 
0521         if (do_backward_search) {
0522           builder.beginBkwSearch();
0523           builder.findTracksCloneEngine(SteeringParams::IT_BkwSearch);
0524         }
0525 
0526         // Post backward-fit filtering.
0527         // Note -- slightly different logic than run_OneIteration as we export both pre and post
0528         // backward-fit tracks.
0529         filter_candidates_func post_filter;
0530         if (itconf.m_post_bkfit_filter)
0531           post_filter = [&](const TrackCand &tc, const MkJob &jb) -> bool {
0532             return itconf.m_post_bkfit_filter(tc, jb) && StdSeq::qfilter_nan_n_silly<TrackCand>(tc, jb);
0533           };
0534         else
0535           post_filter = StdSeq::qfilter_nan_n_silly<TrackCand>;
0536         // post_filter is always at least doing nan_n_silly filter.
0537         builder.filter_comb_cands(post_filter, true);
0538 
0539         if (do_backward_search)
0540           builder.endBkwSearch();
0541 
0542         builder.select_best_comb_cands(true);  // true -> clear m_tracks as they were already filled once above
0543 
0544         if (itconf.m_duplicate_cleaner)
0545           itconf.m_duplicate_cleaner(builder.ref_tracks_nc(), itconf);
0546 
0547         builder.export_tracks(ev.fitTracks_);
0548       }
0549 
0550       builder.end_event();
0551     }
0552 
0553     // MIMI - Fake back event pointer for final processing (that should be done elsewhere)
0554     MkJob job({Config::TrkInfo, Config::ItrInfo[0], eoh, eoh.refBeamSpot()});
0555     builder.begin_event(&job, &ev, __func__);
0556 
0557     if (validation_on) {
0558       StdSeq::prep_simtracks(&ev);
0559       //swap for the cleaned seeds
0560       ev.seedTracks_.swap(seeds_used);
0561     }
0562 
0563     check_nan_n_silly_candidates(ev);
0564 
0565     if (Config::backwardFit)
0566       check_nan_n_silly_bkfit(ev);
0567 
0568     // validation section
0569     if (Config::quality_val) {
0570       StdSeq::Quality qval;
0571       qval.quality_val(&ev);
0572     } else if (Config::sim_val || Config::cmssw_val) {
0573       StdSeq::root_val(&ev);
0574     }
0575 
0576     // ev.print_tracks(ev.candidateTracks_, true);
0577 
0578     // MIMI Unfake.
0579     builder.end_event();
0580 
0581     // In CMSSW runOneIter we now release memory for comb-cands:
0582     builder.release_memory();
0583 
0584     return timevec;
0585   }
0586 
0587 }  // end namespace mkfit