Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-07-01 02:23:05

0001 #ifndef VertexReco_Vertex_h
0002 #define VertexReco_Vertex_h
0003 /** \class reco::Vertex
0004  *  
0005  * A reconstructed Vertex providing position, error, chi2, ndof 
0006  * and reconstrudted tracks.
0007  * The vertex can be valid, fake, or invalid.
0008  * A valid vertex is one which has been obtained from a vertex fit of tracks, 
0009  * and all data is meaningful
0010  * A fake vertex is a vertex which was not made out of a proper fit with
0011  * tracks, but still has a position and error (chi2 and ndof are null). 
0012  * For a primary vertex, it could simply be the beam line.
0013  * A fake vertex is considered valid.
0014  * An invalid vertex has no meaningful data.
0015  *
0016  * \author Luca Lista, INFN
0017  *
0018  *
0019  */
0020 #include <Rtypes.h>
0021 #include "DataFormats/VertexReco/interface/VertexFwd.h"
0022 #include "DataFormats/Math/interface/Error.h"
0023 #include "DataFormats/Math/interface/Point3D.h"
0024 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0025 #include "DataFormats/TrackReco/interface/Track.h"
0026 #include "DataFormats/Common/interface/RefToBase.h"
0027 #include <Math/GenVector/PxPyPzE4D.h>
0028 #include <Math/GenVector/PxPyPzM4D.h>
0029 #include "DataFormats/Math/interface/LorentzVector.h"
0030 
0031 namespace reco {
0032 
0033   class Track;
0034 
0035   class Vertex {
0036   public:
0037     /// The iteratator for the vector<TrackRef>
0038     typedef std::vector<TrackBaseRef>::const_iterator trackRef_iterator;
0039     /// point in the space
0040     typedef math::XYZPoint Point;
0041     /// error matrix dimension
0042     enum { dimension = 3, dimension4D = 4 };
0043     /// covariance error matrix (3x3)
0044     typedef math::Error<dimension>::type Error;
0045     /// covariance error matrix (3x3)
0046     typedef math::Error<dimension>::type CovarianceMatrix;
0047     /// covariance error matrix (4x4)
0048     typedef math::Error<dimension4D>::type Error4D;
0049     /// covariance error matrix (4x4)
0050     typedef math::Error<dimension4D>::type CovarianceMatrix4D;
0051     /// matix size
0052     enum { size = dimension * (dimension + 1) / 2, size4D = (dimension4D) * (dimension4D + 1) / 2 };
0053     /// index type
0054     typedef unsigned int index;
0055     /// default constructor - The vertex will not be valid. Position, error,
0056     /// chi2, ndof will have random entries, and the vectors of tracks will be empty
0057     /// Use the isValid method to check that your vertex is valid.
0058     Vertex() : chi2_(0.0), ndof_(0), position_(0., 0., 0.), time_(0.) {
0059       validity_ = false;
0060       for (int i = 0; i < size4D; ++i)
0061         covariance_[i] = 0.;
0062     }
0063     /// Constructor for a fake vertex.
0064     Vertex(const Point &, const Error &);
0065     /// Constructor for a fake vertex. 4D
0066     Vertex(const Point &, const Error4D &, double);
0067     /// constructor for a valid vertex, with all data
0068     Vertex(const Point &, const Error &, double chi2, double ndof, size_t size);
0069     /// constructor for a valid vertex, with all data 4D
0070     Vertex(const Point &, const Error4D &, double time, double chi2, double ndof, size_t size);
0071     /// Tells whether the vertex is valid.
0072     bool isValid() const { return validity_; }
0073     /// Tells whether a Vertex is fake, i.e. not a vertex made out of a proper
0074     /// fit with tracks.
0075     /// For a primary vertex, it could simply be the beam line.
0076     bool isFake() const { return (chi2_ == 0 && ndof_ == 0 && tracks_.empty()); }
0077     /// reserve space for the tracks
0078     void reserve(int size, bool refitAsWell = false) {
0079       tracks_.reserve(size);
0080       if (refitAsWell)
0081         refittedTracks_.reserve(size);
0082       weights_.reserve(size);
0083     }
0084     /// add a reference to a Track
0085     template <typename Ref>
0086     void add(Ref const &r, float w = 1.0) {
0087       tracks_.emplace_back(r);
0088       weights_.emplace_back(w * 255.f);
0089     }
0090     /// add the original a Track(reference) and the smoothed Track
0091     void add(const TrackBaseRef &r, const Track &refTrack, float w = 1.0);
0092     void removeTracks();
0093 
0094     ///returns the weight with which a Track has contributed to the vertex-fit.
0095     template <typename TREF>
0096     float trackWeight(const TREF &r) const {
0097       int i = 0;
0098       for (auto const &t : tracks_) {
0099         if ((r.id() == t.id()) & (t.key() == r.key()))
0100           return weights_[i] / 255.f;
0101         ++i;
0102       }
0103       return 0;
0104     }
0105     // track collections
0106     auto const &tracks() const { return tracks_; }
0107     /// first iterator over tracks
0108     trackRef_iterator tracks_begin() const { return tracks_.begin(); }
0109     /// last iterator over tracks
0110     trackRef_iterator tracks_end() const { return tracks_.end(); }
0111     /// number of tracks
0112     size_t tracksSize() const { return tracks_.size(); }
0113     /// python friendly track getting
0114     const TrackBaseRef &trackRefAt(size_t idx) const { return tracks_[idx]; }
0115     /// chi-squares
0116     double chi2() const { return chi2_; }
0117     /** Number of degrees of freedom
0118      *  Meant to be Double32_t for soft-assignment fitters: 
0119      *  tracks may contribute to the vertex with fractional weights.
0120      *  The ndof is then = to the sum of the track weights.
0121      *  see e.g. CMS NOTE-2006/032, CMS NOTE-2004/002
0122      */
0123     double ndof() const { return ndof_; }
0124     /// chi-squared divided by n.d.o.f.
0125     double normalizedChi2() const { return ndof_ != 0 ? chi2_ / ndof_ : chi2_ * 1e6; }
0126     /// position
0127     const Point &position() const { return position_; }
0128     /// x coordinate
0129     double x() const { return position_.X(); }
0130     /// y coordinate
0131     double y() const { return position_.Y(); }
0132     /// z coordinate
0133     double z() const { return position_.Z(); }
0134     /// t coordinate
0135     double t() const { return time_; }
0136     /// error on x
0137     double xError() const { return sqrt(covariance(0, 0)); }
0138     /// error on y
0139     double yError() const { return sqrt(covariance(1, 1)); }
0140     /// error on z
0141     double zError() const { return sqrt(covariance(2, 2)); }
0142     /// error on t
0143     double tError() const { return sqrt(covariance(3, 3)); }
0144     /// (i, j)-th element of error matrix, i, j = 0, ... 2
0145     // Note that:
0146     //   double error( int i, int j ) const
0147     // is OBSOLETE, use covariance(i, j)
0148     double covariance(int i, int j) const { return covariance_[idx(i, j)]; }
0149     /// return SMatrix
0150     CovarianceMatrix covariance() const {
0151       Error m;
0152       fill(m);
0153       return m;
0154     }
0155     /// return SMatrix 4D
0156     CovarianceMatrix4D covariance4D() const {
0157       Error4D m;
0158       fill(m);
0159       return m;
0160     }
0161 
0162     /// return SMatrix
0163     Error error() const {
0164       Error m;
0165       fill(m);
0166       return m;
0167     }
0168     /// return SMatrix
0169     Error4D error4D() const {
0170       Error4D m;
0171       fill(m);
0172       return m;
0173     }
0174 
0175     /// fill SMatrix
0176     void fill(CovarianceMatrix &v) const;
0177     /// 4D version
0178     void fill(CovarianceMatrix4D &v) const;
0179 
0180     /// Checks whether refitted tracks are stored.
0181     bool hasRefittedTracks() const { return !refittedTracks_.empty(); }
0182 
0183     /// Returns the original track which corresponds to a particular refitted Track
0184     /// Throws an exception if now refitted tracks are stored ot the track is not found in the list
0185     TrackBaseRef originalTrack(const Track &refTrack) const;
0186 
0187     /// Returns the refitted track which corresponds to a particular original Track
0188     /// Throws an exception if now refitted tracks are stored ot the track is not found in the list
0189     Track refittedTrack(const TrackBaseRef &track) const;
0190 
0191     /// Returns the refitted track which corresponds to a particular original Track
0192     /// Throws an exception if now refitted tracks are stored ot the track is not found in the list
0193     Track refittedTrack(const TrackRef &track) const;
0194 
0195     /// Returns the container of refitted tracks
0196     const std::vector<Track> &refittedTracks() const { return refittedTracks_; }
0197 
0198     /// Returns the four momentum of the sum of the tracks, assuming the given mass for the decay products
0199     math::XYZTLorentzVectorD p4(float mass = 0.13957018, float minWeight = 0.5) const;
0200 
0201     /// Returns the number of tracks in the vertex with weight above minWeight
0202     unsigned int nTracks(float minWeight = 0.5) const;
0203 
0204     class TrackEqual {
0205     public:
0206       TrackEqual(const Track &t) : track_(t) {}
0207       bool operator()(const Track &t) const { return t.pt() == track_.pt(); }
0208 
0209     private:
0210       const Track &track_;
0211     };
0212 
0213   private:
0214     /// chi-sqared
0215     float chi2_;
0216     /// number of degrees of freedom
0217     float ndof_;
0218     /// position
0219     Point position_;
0220     /// covariance matrix (4x4) as vector
0221     float covariance_[size4D];
0222     /// reference to tracks
0223     std::vector<TrackBaseRef> tracks_;
0224     /// The vector of refitted tracks
0225     std::vector<Track> refittedTracks_;
0226     std::vector<uint8_t> weights_;
0227     /// tells wether the vertex is really valid.
0228     bool validity_;
0229     double time_;
0230 
0231     /// position index
0232     index idx(index i, index j) const {
0233       int a = (i <= j ? i : j), b = (i <= j ? j : i);
0234       return b * (b + 1) / 2 + a;
0235     }
0236   };
0237 
0238 }  // namespace reco
0239 
0240 #endif