Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:29:15

0001 #ifndef RecoVertex_PrimaryVertexProducer_AdaptiveChisquarePrimaryVertexFitter_h
0002 #define RecoVertex_PrimaryVertexProducer_AdaptiveChisquarePrimaryVertexFitter_h
0003 
0004 /**\class AdaptiveChisquarePrimaryVertexFitter
0005 
0006   Description: simultaneous chisquared fit of primary vertices 
0007 
0008 */
0009 #include <vector>
0010 
0011 #include "RecoVertex/VertexPrimitives/interface/TransientVertex.h"
0012 #include "TrackingTools/TransientTrack/interface/TransientTrack.h"
0013 #include "RecoVertex/PrimaryVertexProducer/interface/PrimaryVertexFitterBase.h"
0014 
0015 class AdaptiveChisquarePrimaryVertexFitter : public PrimaryVertexFitterBase {
0016 public:
0017   AdaptiveChisquarePrimaryVertexFitter(double chicutoff = 2.5,
0018                                        double zcutoff = 1.0,
0019                                        double mintrkweight = 0.4,
0020                                        bool multivertexfit = false);
0021   ~AdaptiveChisquarePrimaryVertexFitter() override = default;
0022 
0023   std::vector<TransientVertex> fit(const std::vector<reco::TransientTrack> &,
0024                                    const std::vector<TransientVertex> &,
0025                                    const reco::BeamSpot &,
0026                                    const bool) override;
0027 
0028   using Error3 = ROOT::Math::SMatrix<double, 3>;
0029 
0030 protected:
0031   void verify() {  // DEBUG only
0032     unsigned int nt = trackinfo_.size();
0033     unsigned int nv = xv_.size();
0034     assert((yv_.size() == nv) && "yv size");
0035     assert((zv_.size() == nv) && "zv size");
0036     assert((tkfirstv_.size() == (nv + 1)) && "tkfirstv size");
0037     assert((tkmap_.size() == tkweight_.size()) && "tkmapsize <> tkweightssize");
0038     for (unsigned int k = 0; k < nv; k++) {
0039       assert((tkfirstv_[k] < tkweight_.size()) && "tkfirst[k]");
0040       assert((tkfirstv_[k + 1] <= tkweight_.size()) && "tkfirst[k+1]");
0041       assert((tkfirstv_[k] <= tkfirstv_[k + 1]) && "tkfirst[k/k+1]");
0042       for (unsigned int j = tkfirstv_[k]; j < tkfirstv_[k + 1]; k++) {
0043         assert((j < tkmap_.size()) && "illegal tkfirst entry");
0044         unsigned int i = tkmap_[j];
0045         assert((i < nt) && "illegal tkmap entry");
0046         assert((tkweight_[i] >= 0) && "negative tkweight or nan");
0047         assert((tkweight_[i] <= 1) && "tkweight > 1 or nan");
0048       }
0049     }
0050   };
0051 
0052   struct TrackInfo {
0053     double S11, S22, S12;  // inverse of the covariance (sub-)matrix
0054     Error3 C;              // H^T S H
0055     double g[3];
0056     double H1[3], H2[3];
0057     double b1, b2;
0058     double zpca, dzError;
0059   };
0060 
0061   std::vector<TransientVertex> vertices(const std::vector<reco::TransientTrack> &,
0062                                         const std::vector<TransientVertex> &,
0063                                         const reco::BeamSpot &,
0064                                         const bool);
0065   TransientVertex refit(const TransientVertex &, const reco::BeamSpot &, const bool);
0066   double track_in_vertex_chsq(const TrackInfo &, const double, const double, const double);
0067   void fill_trackinfo(const std::vector<reco::TransientTrack> &, const reco::BeamSpot &);
0068   void fill_weights(const reco::BeamSpot &, const double beta = 1.);
0069   TransientVertex get_TransientVertex(const unsigned int,
0070                                       const std::vector<std::pair<unsigned int, float>> &,
0071                                       const std::vector<reco::TransientTrack> &,
0072                                       const float,
0073                                       const reco::BeamSpot &);
0074   Error3 get_inverse_beam_covariance(const reco::BeamSpot &);
0075   double update(const reco::BeamSpot &, float beam_weight, const bool fill_covariances = false);
0076   void make_vtx_trk_map(const double);
0077   bool clean();
0078   void remove_vertex(unsigned int);
0079 
0080   // track information
0081   std::vector<TrackInfo> trackinfo_;
0082 
0083   // vertex lists:
0084   std::vector<double> xv_;
0085   std::vector<double> yv_;
0086   std::vector<double> zv_;
0087   std::vector<Error3> covv_;
0088   // track-vertex-mapping and weights after a coarse z-cut:
0089   std::vector<unsigned int> tkfirstv_;  // parallel to the vertex list
0090   std::vector<unsigned int> tkmap_;     // parallel to tkweight
0091   std::vector<double> tkweight_;        // parallel to tkmap
0092   // configuration
0093   double chi_cutoff_;
0094   double z_cutoff_;
0095   double min_trackweight_;
0096   double multivertexfit_;
0097 };
0098 #endif