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