CMS 3D CMS Logo

StandAloneMuonFilter.h
Go to the documentation of this file.
1 #ifndef RecoMuon_StandAloneTrackFinder_StandAloneMuonFilter_H
2 #define RecoMuon_StandAloneTrackFinder_StandAloneMuonFilter_H
3 
18 
21 
22 class Propagator;
23 class DetLayer;
25 class Trajectory;
28 class MuonServiceProxy;
29 
30 namespace edm {class ParameterSet; class EventSetup; class Event;}
31 
33 
34  public:
37 
39  virtual ~StandAloneMuonFilter();
40 
41  // Operations
42 
44  void refit(const TrajectoryStateOnSurface& initialState, const DetLayer*, Trajectory& trajectory);
45 
47  FreeTrajectoryState lastUpdatedFTS() const {return *theLastUpdatedTSOS.freeTrajectoryState();}
48 
50  FreeTrajectoryState lastButOneUpdatedFTS() const {return *theLastButOneUpdatedTSOS.freeTrajectoryState();}
51 
53  TrajectoryStateOnSurface lastCompatibleTSOS() const {return theLastCompatibleTSOS;}
54 
56  TrajectoryStateOnSurface lastUpdatedTSOS() const {return theLastUpdatedTSOS;}
57 
59  TrajectoryStateOnSurface lastButOneUpdatedTSOS() const {return theLastButOneUpdatedTSOS;}
60 
61  void reset();
62 
64  virtual void setEvent(const edm::Event& event);
65 
66  int getTotalChamberUsed() const {return totalChambers;}
67  int getDTChamberUsed() const {return dtChambers;}
68  int getCSCChamberUsed() const {return cscChambers;}
69  int getRPCChamberUsed() const {return rpcChambers;}
70  int getGEMChamberUsed() const {return gemChambers;}
71  int getME0ChamberUsed() const {return me0Chambers;}
72 
73  int getTotalCompatibleChambers() const {return totalCompatibleChambers;}
74  int getDTCompatibleChambers() const {return dtCompatibleChambers;}
75  int getCSCCompatibleChambers() const {return cscCompatibleChambers;}
76  int getRPCCompatibleChambers() const {return rpcCompatibleChambers;}
77  int getGEMCompatibleChambers() const {return gemCompatibleChambers;}
78  int getME0CompatibleChambers() const {return me0CompatibleChambers;}
79 
80 
81  inline bool goodState() const {return totalChambers >= 2 &&
82  ((dtChambers + cscChambers + gemChambers + me0Chambers) >0 ||
83  onlyRPC());}
84 
85  inline bool isCompatibilitySatisfied() const {return totalCompatibleChambers >= 2 &&
86  ((dtCompatibleChambers + cscCompatibleChambers + gemCompatibleChambers + me0CompatibleChambers) >0 ||
87  onlyRPC());}
88 
90  std::vector<const DetLayer*> layers() const {return theDetLayers;}
91 
93  const DetLayer* lastDetLayer() const {return theDetLayers.back();}
94 
96  PropagationDirection propagationDirection() const;
97 
99  NavigationDirection fitDirection() const {return theFitDirection;}
100 
102  bool onlyRPC() const {return theRPCLoneliness;}
103 
105  MeasurementEstimator *estimator() const {return theEstimator;}
106 
108  const Propagator *propagator() const;
109 
111  MuonTrajectoryUpdator *updator() const {return theMuonUpdator;}
112 
113  void createDefaultTrajectory(const Trajectory &, Trajectory &);
114 
115 
116 protected:
117 
118 private:
119 
121  void setLastCompatibleTSOS(TrajectoryStateOnSurface tsos) { theLastCompatibleTSOS = tsos;}
122 
124  void setLastUpdatedTSOS(TrajectoryStateOnSurface tsos) { theLastUpdatedTSOS = tsos;}
125 
127  void setLastButOneUpdatedTSOS(TrajectoryStateOnSurface tsos) { theLastButOneUpdatedTSOS = tsos;}
128 
130  void incrementChamberCounters(const DetLayer *layer);
131  void incrementCompatibleChamberCounters(const DetLayer *layer);
132 
134  std::vector<const DetLayer*> compatibleLayers(const DetLayer *initialLayer,
135  const FreeTrajectoryState& fts,
136  PropagationDirection propDir);
137 
138  bool update(const DetLayer * layer,
139  const TrajectoryMeasurement * meas,
140  Trajectory & trajectory);
141 
142  std::vector<TrajectoryMeasurement>
143  findBestMeasurements(const DetLayer * layer, const TrajectoryStateOnSurface & tsos);
144 
151 
154 
157 
162 
166  MuonBestMeasurementFinder *bestMeasurementFinder() const {return theBestMeasurementFinder;}
167 
169  double theMaxChi2;
172  double theNSigma;
173 
176 
178  std::vector<const DetLayer*> theDetLayers;
179 
182 
186 
189 
196 
203 
206 };
207 #endif
208 
std::vector< const DetLayer * > layers() const
return the layer used for the refit
MeasurementEstimator * estimator() const
access at the estimator
bool onlyRPC() const
True if there are only the RPC measurements.
NavigationDirection fitDirection() const
Return the fit direction.
double theMaxChi2
The max allowed chi2 to accept a rechit in the fit.
MuonTrajectoryUpdator * theMuonUpdator
the muon updator (it doesn&#39;t inhert from an updator, but it has one!)
bool isCompatibilitySatisfied() const
void setLastCompatibleTSOS(TrajectoryStateOnSurface tsos)
Set the last compatible TSOS.
int getRPCCompatibleChambers() const
PropagationDirection
TrajectoryStateOnSurface lastUpdatedTSOS() const
the Trajectory state on the last surface of the fitting
TrajectoryStateOnSurface theLastCompatibleTSOS
the trajectory state on the last compatible surface
bool theRPCLoneliness
True if there are only the RPC measurements.
std::string theMuonUpdatorName
its name
void setLastUpdatedTSOS(TrajectoryStateOnSurface tsos)
Set the last TSOS.
std::string thePropagatorName
the propagator name
MuonDetLayerMeasurements * theMeasurementExtractor
The Measurement extractor.
int getCSCCompatibleChambers() const
int getDTCompatibleChambers() const
TrajectoryStateOnSurface lastButOneUpdatedTSOS() const
the Trajectory state on the last surface of the fitting
void setLastButOneUpdatedTSOS(TrajectoryStateOnSurface tsos)
Set the last but one TSOS.
int getTotalChamberUsed() const
MuonTrajectoryUpdator * updator() const
access at the muon updator
std::vector< const DetLayer * > theDetLayers
the det layer used in the reconstruction
NavigationDirection theFitDirection
the propagation direction
int getGEMCompatibleChambers() const
MeasurementEstimator * theEstimator
The Estimator.
const DetLayer * lastDetLayer() const
return the last det layer
TrajectoryStateOnSurface theLastButOneUpdatedTSOS
the trajectory state on the last but one available surface
int getME0CompatibleChambers() const
HLT enums.
#define update(a, b)
TrajectoryStateOnSurface theLastUpdatedTSOS
the trajectory state on the last available surface
int getTotalCompatibleChambers() const
MuonBestMeasurementFinder * bestMeasurementFinder() const
Access to the best measurement finder.
FreeTrajectoryState lastUpdatedFTS() const
the last free trajectory state
void reset(double vett[256])
Definition: TPedValues.cc:11
TrajectoryStateOnSurface lastCompatibleTSOS() const
the Trajectory state on the last compatible surface
Definition: event.py:1
MuonBestMeasurementFinder * theBestMeasurementFinder
The best measurement finder: search for the best measurement among the TMs available.
const MuonServiceProxy * theService
FreeTrajectoryState lastButOneUpdatedFTS() const
the last but one free trajectory state