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