Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:04:34

0001 /*
0002 Kalman Track class for 
0003 Kalman Muon Track Finder
0004 Michalis Bachtis (UCLA)
0005 Sep. 2017
0006 */
0007 
0008 #ifndef L1MuKBMTrack_H
0009 #define L1MuKBMTrack_H
0010 
0011 #include "DataFormats/Candidate/interface/LeafCandidate.h"
0012 #include "DataFormats/L1TMuon/interface/L1MuKBMTCombinedStub.h"
0013 #include "DataFormats/L1Trigger/interface/BXVector.h"
0014 
0015 class L1MuKBMTrack;
0016 typedef std::vector<L1MuKBMTrack> L1MuKBMTrackCollection;
0017 typedef BXVector<L1MuKBMTrack> L1MuKBMTrackBxCollection;
0018 
0019 class L1MuKBMTrack : public reco::LeafCandidate {
0020 public:
0021   L1MuKBMTrack();
0022   ~L1MuKBMTrack() override;
0023   L1MuKBMTrack(const L1MuKBMTCombinedStubRef&, int, int);
0024 
0025   //UnConstrained curvature at station 1
0026   int curvatureAtMuon() const;
0027   //unconstrained phi at station 1
0028   int phiAtMuon() const;
0029   //unconstrained phiB at station 1
0030   int phiBAtMuon() const;
0031   //Constrained curvature at vertex
0032   int curvatureAtVertex() const;
0033   //constrained phi at the vertex
0034   int phiAtVertex() const;
0035   //Impact parameter as calculated from the muon track
0036   int dxy() const;
0037   //Unconstrained curvature at the Muon systen
0038   int curvature() const;
0039   //Unconstrained phi at the Muon systen
0040   int positionAngle() const;
0041   //Unconstrained bending angle at the Muon systen
0042   int bendingAngle() const;
0043   //Coarse eta caluclated only using phi segments
0044   int coarseEta() const;
0045   //Approximate Chi2 metric
0046   int approxChi2() const;
0047   int trackCompatibility() const;
0048 
0049   //Approximate Chi2 metric
0050   int hitPattern() const;
0051   //step;
0052   int step() const;
0053   //sector;
0054   int sector() const;
0055   //wheel
0056   int wheel() const;
0057   //quality
0058   int quality() const;
0059 
0060   //unconstrained pt
0061   float ptUnconstrained() const;
0062 
0063   //fine eta
0064   int fineEta() const;
0065   bool hasFineEta() const;
0066 
0067   //BX
0068   int bx() const;
0069 
0070   //rank
0071   int rank() const;
0072 
0073   //Associated stubs
0074   const L1MuKBMTCombinedStubRefVector& stubs() const;
0075 
0076   //get Kalman gain
0077   const std::vector<float>& kalmanGain(unsigned int) const;
0078 
0079   //get covariance
0080   const std::vector<double>& covariance() const;
0081 
0082   //get residual
0083   int residual(uint) const;
0084 
0085   //check ogverlap
0086   bool overlapTrack(const L1MuKBMTrack&) const;
0087 
0088   bool operator==(const L1MuKBMTrack& t2) const {
0089     if (this->stubs().size() != t2.stubs().size())
0090       return false;
0091     for (unsigned int i = 0; i < this->stubs().size(); ++i) {
0092       const L1MuKBMTCombinedStubRef& s1 = this->stubs()[i];
0093       const L1MuKBMTCombinedStubRef& s2 = t2.stubs()[i];
0094       if (s1->scNum() != s2->scNum() || s1->whNum() != s2->whNum() || s1->stNum() != s2->stNum() ||
0095           s1->tag() != s2->tag())
0096         return false;
0097     }
0098     return true;
0099   }
0100 
0101   //Set coordinates general
0102   void setCoordinates(int, int, int, int);
0103 
0104   //Set coordinates at vertex
0105   void setCoordinatesAtVertex(int, int, int);
0106 
0107   //Set coordinates at muon
0108   void setCoordinatesAtMuon(int, int, int);
0109 
0110   //Set eta coarse and pattern
0111   void setCoarseEta(int);
0112 
0113   //Set phi hit pattern
0114   void setHitPattern(int);
0115 
0116   //Set chi2 like metric
0117   void setApproxChi2(int);
0118   void setTrackCompatibility(int);
0119 
0120   //Set floating point coordinates for studies
0121   void setPtEtaPhi(double, double, double);
0122   void setPtUnconstrained(float);
0123 
0124   //Add a stub
0125   void addStub(const L1MuKBMTCombinedStubRef&);
0126 
0127   //kalman gain management
0128   void setKalmanGain(
0129       unsigned int step, unsigned int K, float a1, float a2, float a3, float a4 = 0, float a5 = 0, float a6 = 0);
0130 
0131   //set covariance
0132   void setCovariance(const CovarianceMatrix&);
0133 
0134   //set fine eta
0135   void setFineEta(int);
0136 
0137   //set rank
0138   void setRank(int);
0139 
0140   //set residual
0141   void setResidual(uint, int);
0142 
0143 private:
0144   //Covariance matrix for studies
0145   std::vector<double> covariance_;
0146 
0147   L1MuKBMTCombinedStubRefVector stubs_;
0148 
0149   //vertex coordinates
0150   int curvVertex_ = 0;
0151   int phiVertex_ = 0;
0152   int dxy_ = 0;
0153 
0154   //muon coordinates
0155   int curvMuon_ = 0;
0156   int phiMuon_ = 0;
0157   int phiBMuon_ = 0;
0158 
0159   //generic coordinates
0160   int curv_ = 0;
0161   int phi_ = 0;
0162   int phiB_ = 0;
0163   //common coordinates
0164   int coarseEta_ = 0;
0165 
0166   //Approximate Chi2 metric
0167   int approxChi2_ = 0;
0168   int trackCompatibility_ = 0;
0169 
0170   //phi bitmask
0171   int hitPattern_ = 0;
0172 
0173   //propagation step
0174   int step_ = 0;
0175 
0176   //sector
0177   int sector_ = 0;
0178   //wheel
0179   int wheel_ = 0;
0180 
0181   //quality
0182   int quality_ = 0;
0183 
0184   //Fine eta
0185   int fineEta_ = 0;
0186 
0187   //has fine eta?
0188   bool hasFineEta_ = false;
0189 
0190   //BX
0191   int bx_ = 0;
0192 
0193   //rank
0194   int rank_ = 0;
0195 
0196   //Unconstrained floating point pt
0197   float ptUnconstrained_ = 0;
0198 
0199   //Kalman Gain for making LUTs
0200   std::vector<float> kalmanGain0_;
0201   std::vector<float> kalmanGain1_;
0202   std::vector<float> kalmanGain2_;
0203   std::vector<float> kalmanGain3_;
0204 
0205   std::vector<int> residuals_;
0206 };
0207 
0208 #endif