TrackBase.h

Go to the documentation of this file.

#ifndef TrackReco_TrackBase_h
#define TrackReco_TrackBase_h

#include "DataFormats/Math/interface/Vector.h"
#include "DataFormats/Math/interface/Error.h"
#include "DataFormats/Math/interface/Vector3D.h"
#include "DataFormats/Math/interface/Point3D.h"
#include "DataFormats/Math/interface/Error.h"
#include "DataFormats/TrackReco/interface/HitPattern.h"
#include "DataFormats/BeamSpot/interface/BeamSpot.h"

namespace reco {

  class TrackBase {
  public:
    enum { dimension = 5 };
    enum { covarianceSize = dimension * ( dimension + 1 ) / 2 };
    typedef math::Vector<dimension>::type ParameterVector;
    typedef math::Error<dimension>::type CovarianceMatrix;
    typedef math::XYZVector Vector;
    typedef math::XYZPoint Point;
    enum { i_qoverp = 0 , i_lambda, i_phi, i_dxy, i_dsz }; 
    typedef unsigned int index;
    enum TrackAlgorithm { undefAlgorithm=0, ctf=1, rs=2, cosmics=3, iter0=4, 
              iter1=5, iter2=6, iter3=7, iter4=8, iter5=9, iter6=10, iter7=11, iter8=12, iter9=13,iter10=14,
              outInEcalSeededConv=15, inOutEcalSeededConv=16, 
              nuclInter=17,
              standAloneMuon=18,globalMuon=19,cosmicStandAloneMuon=20,cosmicGlobalMuon=21,
                          iter1LargeD0=22,iter2LargeD0=23,iter3LargeD0=24,iter4LargeD0=25,iter5LargeD0=26,
              bTagGhostTracks=27,
              beamhalo=28, 
              gsf=29,
              algoSize=30 };
    static const std::string algoNames[];

    enum TrackQuality { undefQuality=-1, loose=0, tight=1, highPurity=2, confirmed=3, goodIterative=4, looseSetWithPV=5, highPuritySetWithPV=6, qualitySize=7 };
    static const std::string qualityNames[];

    TrackBase();
    TrackBase( double chi2, double ndof, const Point & referencePoint,
           const Vector & momentum, int charge, const CovarianceMatrix &,
           TrackAlgorithm=undefAlgorithm, TrackQuality quality=undefQuality,signed char nloops=0);
    ~TrackBase();
    double chi2() const { return chi2_; }
    double ndof() const { return ndof_; }
    double normalizedChi2() const { return ndof_ != 0 ? chi2_ / ndof_ : chi2_ * 1e6; }
    int charge() const { return charge_; }
    double qoverp() const { return charge() / p(); }
    double theta() const { return momentum_.theta(); }
    double lambda() const { return M_PI/2 - momentum_.theta(); }
    double dxy() const { return ( - vx() * py() + vy() * px() ) / pt(); }
    double d0() const { return - dxy(); }
    double dsz() const { return vz()*pt()/p() - (vx()*px()+vy()*py())/pt() * pz()/p(); }
    double dz() const { return vz() - (vx()*px()+vy()*py())/pt() * (pz()/pt()); }
    double p() const { return momentum_.R(); }
    double pt() const { return sqrt( momentum_.Perp2() ); }
    double px() const { return momentum_.x(); }
    double py() const { return momentum_.y(); }
    double pz() const { return momentum_.z(); }
    double phi() const { return momentum_.Phi(); }
    double eta() const { return momentum_.Eta(); }
    double vx() const { return vertex_.x(); }
    double vy() const { return vertex_.y(); }
    double vz() const { return vertex_.z(); }

    const Vector & momentum() const { return momentum_; }

    const Point & referencePoint() const { return vertex_; }

    const Point & vertex() const { return vertex_; }

    double dxy(const Point& myBeamSpot) const { 
      return ( - (vx()-myBeamSpot.x()) * py() + (vy()-myBeamSpot.y()) * px() ) / pt(); 
    }

    double dxy(const BeamSpot& theBeamSpot) const { 
      return dxy(theBeamSpot.position(vz()));
    }

    double dsz(const Point& myBeamSpot) const { 
      return (vz()-myBeamSpot.z())*pt()/p() - ((vx()-myBeamSpot.x())*px()+(vy()-myBeamSpot.y())*py())/pt() *pz()/p(); 
    }
    double dz(const Point& myBeamSpot) const { 
      return (vz()-myBeamSpot.z()) - ((vx()-myBeamSpot.x())*px()+(vy()-myBeamSpot.y())*py())/pt() * pz()/pt(); 
    }

    ParameterVector parameters() const { 
      return ParameterVector(qoverp(),lambda(),phi(),dxy(),dsz());
    }
    CovarianceMatrix covariance() const { CovarianceMatrix m; fill( m ); return m; }

    double parameter(int i) const { return parameters()[i]; }
    double covariance( int i, int j ) const { return covariance_[ covIndex( i, j ) ]; }
    double error( int i ) const { return sqrt( covariance_[ covIndex( i, i ) ] ); }

    double qoverpError() const { return error( i_qoverp ); }
    double ptError() const { 
      return (charge()!=0) ?  sqrt( 
            pt()*pt()*p()*p()/charge()/charge() * covariance(i_qoverp,i_qoverp)
          + 2*pt()*p()/charge()*pz() * covariance(i_qoverp,i_lambda)
          + pz()*pz() * covariance(i_lambda,i_lambda) ) : 1.e6;
    }
    double thetaError() const { return error( i_lambda ); }
    double lambdaError() const { return error( i_lambda ); }
    double etaError() const { return error( i_lambda ) * p()/pt(); }
    double phiError() const { return error( i_phi ); }
    double dxyError() const { return error( i_dxy ); }
    double d0Error() const { return error( i_dxy ); }
    double dszError() const { return error( i_dsz ); }
    double dzError() const { return error( i_dsz ) * p()/pt(); }

    CovarianceMatrix & fill( CovarianceMatrix & v ) const;
    static index covIndex( index i, index j );
   
    const HitPattern & hitPattern() const { return hitPattern_; }
    const HitPattern & trackerExpectedHitsInner() const { return trackerExpectedHitsInner_; }
    const HitPattern & trackerExpectedHitsOuter() const { return trackerExpectedHitsOuter_; }



    unsigned short numberOfValidHits() const { return hitPattern_.numberOfValidHits(); }
    unsigned short numberOfLostHits() const { return hitPattern_.numberOfLostHits(); }
    double validFraction() const {
      int valid = hitPattern_.numberOfValidTrackerHits();
      int lost  = hitPattern_.numberOfLostTrackerHits();
      int lostIn = trackerExpectedHitsInner_.numberOfLostTrackerHits();
      int lostOut = trackerExpectedHitsOuter_.numberOfLostTrackerHits();
      if ((valid+lost+lostIn+lostOut)==0) return -1;
      return valid/(1.0*(valid+lost+lostIn+lostOut));
    }
    template<typename C>
    void setHitPattern( const C & c ) { hitPattern_.set( c.begin(), c.end() ); }
    template<typename C>
    void setTrackerExpectedHitsInner( const C & c ) { trackerExpectedHitsInner_.set( c.begin(), c.end() ); }
    template<typename C>
    void setTrackerExpectedHitsOuter( const C & c ) { trackerExpectedHitsOuter_.set( c.begin(), c.end() ); }
    
    template<typename I>
    void setHitPattern( const I & begin, const I & end ) { hitPattern_.set( begin, end ); }
    template<typename I>
    void setTrackerExpectedHitsInner( const I & begin, const I & end ) { trackerExpectedHitsInner_.set( begin, end ); }
    template<typename I>
    void setTrackerExpectedHitsOuter( const I & begin, const I & end ) { trackerExpectedHitsOuter_.set( begin, end ); }

    void setHitPattern( const TrackingRecHit & hit, size_t i ) { hitPattern_.set( hit, i ); }
    void appendHitPattern( const TrackingRecHit & hit) { hitPattern_.appendHit( hit); }
    void setTrackerExpectedHitsInner( const TrackingRecHit & hit, size_t i ) { trackerExpectedHitsInner_.set( hit, i ); }
    void setTrackerExpectedHitsOuter( const TrackingRecHit & hit, size_t i ) { trackerExpectedHitsOuter_.set( hit, i ); }

    void setHitPattern( const HitPattern& hitP) {hitPattern_ = hitP;}
    void setTrackerExpectedHitsInner( const HitPattern& hitP ) {trackerExpectedHitsInner_ = hitP;}
    void setTrackerExpectedHitsOuter( const HitPattern& hitP ) {trackerExpectedHitsOuter_ = hitP;}


    void setAlgorithm(const TrackAlgorithm a, bool set=true) { if (set) algorithm_=a; else algorithm_=undefAlgorithm;}
    TrackAlgorithm algo() const ;
    std::string algoName() const;
    static std::string algoName(TrackAlgorithm);
    static TrackAlgorithm algoByName(const std::string &name);

    bool quality(const TrackQuality ) const;
    void setQuality(const TrackQuality, bool set=true); 
    static std::string qualityName(TrackQuality);
    static TrackQuality qualityByName(const std::string &name);


    int qualityMask() const { return quality_; }
    void setQualityMask(int qualMask) {quality_ = qualMask;}

   
    void setNLoops(signed char value) { nLoops_=value;}

    bool isLooper() const { return (nLoops_>0);}
    signed char nLoops() const {return nLoops_;}
   
  private:
    HitPattern hitPattern_;
    HitPattern trackerExpectedHitsInner_; 
    HitPattern trackerExpectedHitsOuter_;

    float covariance_[ covarianceSize ];
    float chi2_;
    Point vertex_;
    Vector momentum_;

    float ndof_;

    char charge_;
    uint8_t algorithm_;
    uint8_t quality_;

    signed char nLoops_; 

  };

  inline TrackBase::index TrackBase::covIndex( index i, index j )  {
    int a = ( i <= j ? i : j ), b = ( i <= j ? j : i );
    return b * ( b + 1 ) / 2 + a;
  }

  inline TrackBase::TrackAlgorithm TrackBase::algo() const {
    return (TrackAlgorithm) algorithm_;
  }

  inline std::string TrackBase::algoName() const{
    // I'd like to do:
    // return TrackBase::algoName(algorithm_);
    // but I cannot define a const static function. Why???

    switch(algorithm_)
      {
      case undefAlgorithm: return "undefAlgorithm";
      case ctf: return "ctf";
      case rs: return "rs";
      case cosmics: return "cosmics";
      case beamhalo: return "beamhalo";
      case iter0: return "iter0";
      case iter1: return "iter1";
      case iter2: return "iter2";
      case iter3: return "iter3";
      case iter4: return "iter4";
      case iter5: return "iter5";
      case iter6: return "iter6";
      case iter7: return "iter7";
      case iter8: return "iter8";
      case iter9: return "iter9";
      case iter10: return "iter10";
      case outInEcalSeededConv: return "outInEcalSeededConv";
      case inOutEcalSeededConv: return "inOutEcalSeededConv";
      case nuclInter: return "nuclInter";
      case standAloneMuon: return "standAloneMuon";
      case globalMuon: return "globalMuon";
      case cosmicStandAloneMuon: return "cosmicStandAloneMuon";
      case cosmicGlobalMuon: return "cosmicGlobalMuon";
      case iter1LargeD0: return "iter1LargeD0";
      case iter2LargeD0: return "iter2LargeD0";
      case iter3LargeD0: return "iter3LargeD0";
      case iter4LargeD0: return "iter4LargeD0";
      case iter5LargeD0: return "iter5LargeD0";
      case bTagGhostTracks: return "bTagGhostTracks";
      case gsf: return "gsf";
      }
    return "undefAlgorithm";
  }
      
  inline bool TrackBase::quality( const TrackBase::TrackQuality q) const{
    switch(q){
    case undefQuality: return (quality_==0);
    case goodIterative: 
      {
        return ( ((quality_ & (1<<TrackBase::confirmed))>>TrackBase::confirmed) ||
                 ((quality_ & (1<<TrackBase::highPurity))>>TrackBase::highPurity));
      }
     default:
      {
    return (quality_ & (1<<q))>>q;
      }
    }
    return false;
  }

  inline void TrackBase::setQuality( const TrackBase::TrackQuality q, bool set){
    switch(q){
    case undefQuality: quality_=0;
    default: 
      {
    if (set)//regular OR if setting value to true
      quality_= quality_ | (1<<q) ;
    else // doing "half-XOR" if unsetting value
      quality_= quality_ & (~(1<<q));
      }
    }
  }

  inline std::string TrackBase::qualityName(TrackQuality q) {
    if(int(q) < int(qualitySize) && int(q)>=0) return qualityNames[int(q)];
    return "undefQuality";
  }

  inline std::string TrackBase::algoName(TrackAlgorithm a){
    if(int(a) < int(algoSize) && int(a)>0) return algoNames[int(a)];
    return "undefAlgorithm";
  }



}

#endif