CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_5_3_0/src/RecoVertex/AdaptiveVertexFit/interface/AdaptiveVertexFitter.h

Go to the documentation of this file.
00001 #ifndef AdaptiveVertexFitter_H
00002 #define AdaptiveVertexFitter_H
00003 
00004 #include "RecoVertex/VertexPrimitives/interface/VertexFitter.h"
00005 #include "RecoVertex/LinearizationPointFinders/interface/DefaultLinearizationPointFinder.h"
00006 #include "RecoVertex/VertexTools/interface/GeometricAnnealing.h"
00007 #include "RecoVertex/KalmanVertexFit/interface/KalmanVertexUpdator.h"
00008 #include "RecoVertex/VertexTools/interface/DummyVertexSmoother.h"
00009 // #include "RecoVertex/AdaptiveVertexFit/interface/KalmanVertexSmoother.h"
00010 #include "RecoVertex/VertexPrimitives/interface/VertexState.h"
00011 #include "RecoVertex/VertexPrimitives/interface/VertexTrackCompatibilityEstimator.h"
00012 #include "RecoVertex/VertexTools/interface/LinearizedTrackStateFactory.h"
00013 #include "RecoVertex/KalmanVertexFit/interface/KalmanVertexTrackCompatibilityEstimator.h"
00014 #include "FWCore/ParameterSet/interface/ParameterSet.h"
00015 
00029 class AdaptiveVertexFitter : public VertexFitter<5> {
00030 
00031 public:
00032 
00033   typedef ReferenceCountingPointer<VertexTrack<5> > RefCountedVertexTrack;
00034   typedef ReferenceCountingPointer<LinearizedTrackState<5> > RefCountedLinearizedTrackState;
00035 
00041   AdaptiveVertexFitter(
00042       const AnnealingSchedule & ann = GeometricAnnealing(),
00043       const LinearizationPointFinder & linP =
00044              DefaultLinearizationPointFinder(),
00045       const VertexUpdator<5> & updator = KalmanVertexUpdator<5>(),
00046       const VertexTrackCompatibilityEstimator<5> & estor =
00047              KalmanVertexTrackCompatibilityEstimator<5>(),
00048       const VertexSmoother<5> & smoother = DummyVertexSmoother<5>(),
00049       const AbstractLTSFactory<5> & ltsf = LinearizedTrackStateFactory() );
00050 
00051   AdaptiveVertexFitter( const AdaptiveVertexFitter & original );
00052 
00053   virtual ~AdaptiveVertexFitter();
00054 
00061   virtual CachingVertex<5> vertex( const std::vector<reco::TransientTrack> & ) const;
00062 
00071   virtual CachingVertex<5> vertex(const std::vector<RefCountedVertexTrack> & ) const;
00072 
00076   virtual CachingVertex<5> vertex(const std::vector<RefCountedVertexTrack> &,
00077      const reco::BeamSpot & spot ) const;
00078 
00082   virtual CachingVertex<5> vertex( const std::vector<reco::TransientTrack> &,
00083                                 const GlobalPoint& linPoint ) const;
00084 
00090   virtual CachingVertex<5> vertex( const std::vector<reco::TransientTrack> &,
00091                                 const GlobalPoint & priorPos,
00092                                 const GlobalError & priorError ) const;
00093 
00098   virtual CachingVertex<5> vertex(const std::vector<reco::TransientTrack> & tracks,
00099                 const reco::BeamSpot& beamSpot) const;
00100 
00101 
00106   virtual CachingVertex<5> vertex( const std::vector<RefCountedVertexTrack> &,
00107                                 const GlobalPoint & priorPos,
00108                                 const GlobalError & priorError ) const;
00109 
00110   AdaptiveVertexFitter * clone() const;
00111 
00118   void setWeightThreshold ( float w );
00119 
00131   void setParameters( double maxshift=0.0001, double maxlpshift=0.1, 
00132                       unsigned maxstep=30, double weightthreshold=.001 );
00133 
00139   void setParameters ( const edm::ParameterSet & );
00140 
00141   void gsfIntermediarySmoothing(bool sm) { gsfIntermediarySmoothing_ = sm; }
00142 
00143   bool gsfIntermediarySmoothing() const { return gsfIntermediarySmoothing_;}
00144 
00145 private:
00156   std::vector<RefCountedVertexTrack> reLinearizeTracks(
00157                 const std::vector<RefCountedVertexTrack> & tracks,
00158                 const CachingVertex<5> & vertex ) const;
00159 
00160 
00165   std::vector<RefCountedVertexTrack> reWeightTracks(
00166                         const std::vector<RefCountedLinearizedTrackState> &,
00167                         const CachingVertex<5> & seed ) const;
00168 
00174   std::vector<RefCountedVertexTrack> reWeightTracks(
00175                         const std::vector<RefCountedVertexTrack> &,
00176                         const CachingVertex<5> & seed) const;
00177 
00178 
00184   std::vector<RefCountedVertexTrack> weightTracks(
00185                         const std::vector<RefCountedLinearizedTrackState> &,
00186                         const VertexState & seed ) const;
00187 
00191   std::vector<RefCountedVertexTrack> linearizeTracks(
00192                         const std::vector<reco::TransientTrack> &,
00193                         const VertexState & ) const;
00197   CachingVertex<5> fit( const std::vector<RefCountedVertexTrack> & tracks,
00198                      const VertexState & priorSeed,
00199                      bool withPrior) const;
00200 
00201   double getWeight ( float chi2 ) const;
00202 private:
00203   double theMaxShift;
00204   double theMaxLPShift;
00205   int theMaxStep;
00206   double theWeightThreshold;
00207   mutable int theNr;
00208 
00209   LinearizationPointFinder * theLinP;
00210   VertexUpdator<5> * theUpdator;
00211   VertexSmoother<5> * theSmoother;
00212   AnnealingSchedule * theAssProbComputer;
00213   VertexTrackCompatibilityEstimator<5> * theComp;
00214   const AbstractLTSFactory<5> * theLinTrkFactory;
00215   bool gsfIntermediarySmoothing_;
00216   mutable int mctr_;
00217 };
00218 
00219 #endif