MkFitCore/src/MkFinder.h

0001 #ifndef RecoTracker_MkFitCore_src_MkFinder_h
0002 #define RecoTracker_MkFitCore_src_MkFinder_h
0003
0004 #include "MkBase.h"
0005 #include "RecoTracker/MkFitCore/interface/TrackerInfo.h"
0006 #include "RecoTracker/MkFitCore/interface/Track.h"
0007
0008 #include "RecoTracker/MkFitCore/interface/HitStructures.h"
0009 #include "RecoTracker/MkFitCore/interface/TrackStructures.h"
0010
0011 // Define to get printouts about track and hit chi2.
0012 // See also MkBuilder::backwardFit().
0013
0014 //#define DEBUG_BACKWARD_FIT_BH
0015 //#define DEBUG_BACKWARD_FIT
0016
0017 namespace mkfit {
0018
0019   class CandCloner;
0020   class CombCandidate;
0021   class LayerOfHits;
0022   class FindingFoos;
0023   class IterationParams;
0024   class IterationLayerConfig;
0025   class SteeringParams;
0026
0027   class Event;
0028
0029   struct UpdateIndices {
0030     int seed_idx;
0031     int cand_idx;
0032     int hit_idx;
0033
0034     UpdateIndices(int si, int ci, int hi) : seed_idx(si), cand_idx(ci), hit_idx(hi) {}
0035   };
0036
0037   class MkFinder : public MkBase {
0038     friend class MkBuilder;
0039
0040   public:
0041     static constexpr int MPlexHitIdxMax = 16;
0042
0043     using MPlexHitIdx = Matriplex::Matriplex<int, MPlexHitIdxMax, 1, NN>;
0044     using MPlexQHoT = Matriplex::Matriplex<HitOnTrack, 1, 1, NN>;
0045
0046     //----------------------------------------------------------------------------
0047
0048     MkFinder() {}
0049
0050     void setup(const PropagationConfig &pc,
0051                const IterationConfig &ic,
0052                const IterationParams &ip,
0053                const IterationLayerConfig &ilc,
0054                const SteeringParams &sp,
0055                const std::vector<bool> *ihm,
0056                const Event *ev,
0057                int region,
0058                bool infwd);
0059     void setup_bkfit(const PropagationConfig &pc, const SteeringParams &sp, const Event *ev);
0060     void release();
0061
0062     void begin_layer(const LayerOfHits &layer_of_hits);
0063     void end_layer();
0064
0065     //----------------------------------------------------------------------------
0066
0067     void inputTracksAndHitIdx(const std::vector<Track> &tracks, int beg, int end, bool inputProp);
0068
0069     void inputTracksAndHitIdx(const std::vector<Track> &tracks,
0070                               const std::vector<int> &idxs,
0071                               int beg,
0072                               int end,
0073                               bool inputProp,
0074                               int mp_offset);
0075
0076     void inputTracksAndHitIdx(const std::vector<CombCandidate> &tracks,
0077                               const std::vector<std::pair<int, int>> &idxs,
0078                               int beg,
0079                               int end,
0080                               bool inputProp);
0081
0082     void inputTracksAndHits(const std::vector<CombCandidate> &tracks,
0083                             const LayerOfHits &layer_of_hits,
0084                             const std::vector<UpdateIndices> &idxs,
0085                             int beg,
0086                             int end,
0087                             bool inputProp);
0088
0089     void inputTracksAndHitIdx(const std::vector<CombCandidate> &tracks,
0090                               const std::vector<std::pair<int, IdxChi2List>> &idxs,
0091                               int beg,
0092                               int end,
0093                               bool inputProp);
0094
0095     void outputTracksAndHitIdx(std::vector<Track> &tracks, int beg, int end, bool outputProp) const;
0096
0097     void outputTracksAndHitIdx(
0098         std::vector<Track> &tracks, const std::vector<int> &idxs, int beg, int end, bool outputProp) const;
0099
0100     void outputTrackAndHitIdx(Track &track, int itrack, bool outputProp) const {
0101       const int iO = outputProp ? iP : iC;
0102       copy_out(track, itrack, iO);
0103     }
0104
0105     void outputNonStoppedTracksAndHitIdx(
0106         std::vector<Track> &tracks, const std::vector<int> &idxs, int beg, int end, bool outputProp) const {
0107       const int iO = outputProp ? iP : iC;
0108       for (int i = beg, imp = 0; i < end; ++i, ++imp) {
0109         if (!m_Stopped[imp])
0110           copy_out(tracks[idxs[i]], imp, iO);
0111       }
0112     }
0113
0114     HitOnTrack bestHitLastHoT(int itrack) const { return m_HoTArrs[itrack][m_NHits(itrack, 0, 0) - 1]; }
0115
0116     void packModuleNormDir(const LayerOfHits &layer_of_hits, int hit_cnt, MPlexHV &norm, MPlexHV &dir, int N_proc) const;
0117
0118     //----------------------------------------------------------------------------
0119
0120     void getHitSelDynamicWindows(
0121         const float invpt, const float theta, float &min_dq, float &max_dq, float &min_dphi, float &max_dphi);
0122
0123     float getHitSelDynamicChi2Cut(const int itrk, const int ipar);
0124
0125     void selectHitIndices(const LayerOfHits &layer_of_hits, const int N_proc, bool fill_binsearch_only = false);
0126     void selectHitIndicesV2(const LayerOfHits &layer_of_hits, const int N_proc);
0127
0128     void addBestHit(const LayerOfHits &layer_of_hits, const int N_proc, const FindingFoos &fnd_foos);
0129
0130     //----------------------------------------------------------------------------
0131
0132     void findCandidates(const LayerOfHits &layer_of_hits,
0133                         std::vector<std::vector<TrackCand>> &tmp_candidates,
0134                         const int offset,
0135                         const int N_proc,
0136                         const FindingFoos &fnd_foos);
0137
0138     //----------------------------------------------------------------------------
0139
0140     void findCandidatesCloneEngine(const LayerOfHits &layer_of_hits,
0141                                    CandCloner &cloner,
0142                                    const int offset,
0143                                    const int N_proc,
0144                                    const FindingFoos &fnd_foos);
0145
0146     void updateWithLoadedHit(int N_proc, const LayerOfHits &layer_of_hits, const FindingFoos &fnd_foos);
0147
0148     void copyOutParErr(std::vector<CombCandidate> &seed_cand_vec, int N_proc, bool outputProp) const;
0149
0150     //----------------------------------------------------------------------------
0151     // Backward fit
0152
0153     void bkFitInputTracks(TrackVec &cands, int beg, int end);
0154     void bkFitOutputTracks(TrackVec &cands, int beg, int end, bool outputProp);
0155
0156     void bkFitInputTracks(EventOfCombCandidates &eocss, int beg, int end);
0157     void bkFitOutputTracks(EventOfCombCandidates &eocss, int beg, int end, bool outputProp);
0158
0159     void bkFitFitTracksBH(const EventOfHits &eventofhits,
0160                           const SteeringParams &st_par,
0161                           const int N_proc,
0162                           bool chiDebug = false);
0163
0164     void bkFitFitTracks(const EventOfHits &eventofhits,
0165                         const SteeringParams &st_par,
0166                         const int N_proc,
0167                         bool chiDebug = false);
0168
0169     void bkFitPropTracksToPCA(const int N_proc);
0170
0171     //----------------------------------------------------------------------------
0172
0173   private:
0174     void copy_in(const Track &trk, const int mslot, const int tslot) {
0175       m_Err[tslot].copyIn(mslot, trk.errors().Array());
0176       m_Par[tslot].copyIn(mslot, trk.parameters().Array());
0177
0178       m_Chg(mslot, 0, 0) = trk.charge();
0179       m_Chi2(mslot, 0, 0) = trk.chi2();
0180       m_Label(mslot, 0, 0) = trk.label();
0181
0182       m_NHits(mslot, 0, 0) = trk.nTotalHits();
0183       m_NFoundHits(mslot, 0, 0) = trk.nFoundHits();
0184
0185       m_NInsideMinusOneHits(mslot, 0, 0) = trk.nInsideMinusOneHits();
0186       m_NTailMinusOneHits(mslot, 0, 0) = trk.nTailMinusOneHits();
0187
0188       std::copy(trk.beginHitsOnTrack(), trk.endHitsOnTrack(), m_HoTArrs[mslot]);
0189     }
0190
0191     void copy_out(Track &trk, const int mslot, const int tslot) const {
0192       m_Err[tslot].copyOut(mslot, trk.errors_nc().Array());
0193       m_Par[tslot].copyOut(mslot, trk.parameters_nc().Array());
0194
0195       trk.setCharge(m_Chg(mslot, 0, 0));
0196       trk.setChi2(m_Chi2(mslot, 0, 0));
0197       trk.setLabel(m_Label(mslot, 0, 0));
0198
0199       trk.resizeHits(m_NHits(mslot, 0, 0), m_NFoundHits(mslot, 0, 0));
0200       std::copy(m_HoTArrs[mslot], &m_HoTArrs[mslot][m_NHits(mslot, 0, 0)], trk.beginHitsOnTrack_nc());
0201     }
0202
0203     void copy_in(const TrackCand &trk, const int mslot, const int tslot) {
0204       m_Err[tslot].copyIn(mslot, trk.errors().Array());
0205       m_Par[tslot].copyIn(mslot, trk.parameters().Array());
0206
0207       m_Chg(mslot, 0, 0) = trk.charge();
0208       m_Chi2(mslot, 0, 0) = trk.chi2();
0209       m_Label(mslot, 0, 0) = trk.label();
0210
0211       m_LastHitCcIndex(mslot, 0, 0) = trk.lastCcIndex();
0212       m_NFoundHits(mslot, 0, 0) = trk.nFoundHits();
0213       m_NMissingHits(mslot, 0, 0) = trk.nMissingHits();
0214       m_NOverlapHits(mslot, 0, 0) = trk.nOverlapHits();
0215
0216       m_NInsideMinusOneHits(mslot, 0, 0) = trk.nInsideMinusOneHits();
0217       m_NTailMinusOneHits(mslot, 0, 0) = trk.nTailMinusOneHits();
0218
0219       m_CombCand[mslot] = trk.combCandidate();
0220       m_TrkStatus[mslot] = trk.getStatus();
0221     }
0222
0223     void copy_out(TrackCand &trk, const int mslot, const int tslot) const {
0224       m_Err[tslot].copyOut(mslot, trk.errors_nc().Array());
0225       m_Par[tslot].copyOut(mslot, trk.parameters_nc().Array());
0226
0227       trk.setCharge(m_Chg(mslot, 0, 0));
0228       trk.setChi2(m_Chi2(mslot, 0, 0));
0229       trk.setLabel(m_Label(mslot, 0, 0));
0230
0231       trk.setLastCcIndex(m_LastHitCcIndex(mslot, 0, 0));
0232       trk.setNFoundHits(m_NFoundHits(mslot, 0, 0));
0233       trk.setNMissingHits(m_NMissingHits(mslot, 0, 0));
0234       trk.setNOverlapHits(m_NOverlapHits(mslot, 0, 0));
0235
0236       trk.setNInsideMinusOneHits(m_NInsideMinusOneHits(mslot, 0, 0));
0237       trk.setNTailMinusOneHits(m_NTailMinusOneHits(mslot, 0, 0));
0238
0239       trk.setCombCandidate(m_CombCand[mslot]);
0240       trk.setStatus(m_TrkStatus[mslot]);
0241     }
0242
0243     void add_hit(const int mslot, int index, int layer) {
0244       // Only used by BestHit.
0245       // m_NInsideMinusOneHits and m_NTailMinusOneHits are maintained here but are
0246       // not used and are not copied out (as Track does not have these members).
0247
0248       int &n_tot_hits = m_NHits(mslot, 0, 0);
0249       int &n_fnd_hits = m_NFoundHits(mslot, 0, 0);
0250
0251       if (n_tot_hits < Config::nMaxTrkHits) {
0252         m_HoTArrs[mslot][n_tot_hits++] = {index, layer};
0253         if (index >= 0) {
0254           ++n_fnd_hits;
0255           m_NInsideMinusOneHits(mslot, 0, 0) += m_NTailMinusOneHits(mslot, 0, 0);
0256           m_NTailMinusOneHits(mslot, 0, 0) = 0;
0257         } else if (index == -1) {
0258           ++m_NTailMinusOneHits(mslot, 0, 0);
0259         }
0260       } else {
0261         // printf("WARNING MkFinder::add_hit hit-on-track limit reached for label=%d\n", label_);
0262
0263         const int LH = Config::nMaxTrkHits - 1;
0264
0265         if (index >= 0) {
0266           if (m_HoTArrs[mslot][LH].index < 0)
0267             ++n_fnd_hits;
0268           m_HoTArrs[mslot][LH] = {index, layer};
0269         } else if (index == -2) {
0270           if (m_HoTArrs[mslot][LH].index >= 0)
0271             --n_fnd_hits;
0272           m_HoTArrs[mslot][LH] = {index, layer};
0273         }
0274       }
0275     }
0276
0277     int num_all_minus_one_hits(const int mslot) const {
0278       return m_NInsideMinusOneHits(mslot, 0, 0) + m_NTailMinusOneHits(mslot, 0, 0);
0279     }
0280
0281     int num_inside_minus_one_hits(const int mslot) const { return m_NInsideMinusOneHits(mslot, 0, 0); }
0282
0283     void print_par_err(int corp, int mslot) const;
0284
0285     //----------------------------------------------------------------------------
0286
0287     MPlexQF m_Chi2;
0288     MPlexQI m_Label;  // seed index in global seed vector (for MC truth match)
0289
0290     MPlexQI m_NHits;
0291     MPlexQI m_NFoundHits;
0292
0293     HitOnTrack m_HoTArrs[NN][Config::nMaxTrkHits];
0294
0295     MPlexQI m_SeedIdx;        // seed index in local thread (for bookkeeping at thread level)
0296     MPlexQI m_CandIdx;        // candidate index for the given seed (for bookkeeping of clone engine)
0297     MPlexQI m_SeedOriginIdx;  // seed index in MkBuilder seed input vector
0298
0299     MPlexQI m_Stopped;  // Flag for BestHit that a track has been stopped (and copied out already)
0300
0301     // Additions / substitutions for TrackCand copy_in/out()
0302     // One could really access the original TrackCand for all of those, especially the ones that
0303     // are STD only. This then requires access back to that TrackCand memory.
0304     // So maybe one should just have flags for CopyIn methods (or several versions). Yay, etc.
0305     MPlexQI m_NMissingHits;             // sub: m_NHits, sort of, STD only
0306     MPlexQI m_NOverlapHits;             // add: num of overlaps registered in HitOnTrack, STD only
0307     MPlexQI m_NInsideMinusOneHits;      // sub: before we copied all hit idcs and had a loop counting them only
0308     MPlexQI m_NTailMinusOneHits;        // sub: before we copied all hit idcs and had a loop counting them only
0309     MPlexQI m_LastHitCcIndex;           // add: index of last hit in m_CombCand hit tree, STD only
0310     TrackBase::Status m_TrkStatus[NN];  // STD only, status bits
0311     CombCandidate *m_CombCand[NN];
0312     // const TrackCand *m_TrkCand[NN]; // hmmh, could get all data through this guy ... but scattered
0313     // storing it in now for bkfit debug printouts
0314     TrackCand *m_TrkCand[NN];
0315
0316     // Hit indices into LayerOfHits to explore.
0317     WSR_Result m_XWsrResult[NN];  // Could also merge it with m_XHitSize. Or use smaller arrays.
0318     MPlexQI m_XHitSize;
0319     MPlexHitIdx m_XHitArr;
0320
0321     // Hit errors / parameters for hit matching, update.
0322     MPlexHS m_msErr;
0323     MPlexHV m_msPar;
0324
0325     // An idea: Do propagation to hit in FindTracksXYZZ functions.
0326     // Have some state / functions here that make this short to write.
0327     // This would simplify KalmanUtils (remove the propagate functions).
0328     // Track errors / parameters propagated to current hit.
0329     // MPlexLS    candErrAtCurrHit;
0330     // MPlexLV    candParAtCurrHit;
0331
0332     const PropagationConfig *m_prop_config = nullptr;
0333     const IterationConfig *m_iteration_config = nullptr;
0334     const IterationParams *m_iteration_params = nullptr;
0335     const IterationLayerConfig *m_iteration_layer_config = nullptr;
0336     const SteeringParams *m_steering_params = nullptr;
0337     const std::vector<bool> *m_iteration_hit_mask = nullptr;
0338     const Event *m_event = nullptr;
0339     int m_current_region = -1;
0340     bool m_in_fwd = true;
0341
0342     // Backward fit
0343     int m_CurHit[NN];
0344     const HitOnTrack *m_HoTArr[NN];
0345     int m_CurNode[NN];
0346     const HoTNode *m_HoTNodeArr[NN];
0347   };
0348
0349 }  // end namespace mkfit
0350 #endif