CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
GlobalTrajectoryBuilderBase.h
Go to the documentation of this file.
1 #ifndef RecoMuon_GlobalTrackingTools_GlobalTrajectoryBuilderBase_H
2 #define RecoMuon_GlobalTrackingTools_GlobalTrajectoryBuilderBase_H
3 
29 
34 class MuonServiceProxy;
35 class Trajectory;
36 class TrackTransformer;
37 class TrajectoryFitter;
40 class GlobalMuonRefitter;
41 class TrackerTopology;
42 
43 namespace edm {class ParameterSet; class Event;}
44 namespace reco {class TransientTrack;}
45 
46 // ---------------------
47 // -- Class Interface --
48 // ---------------------
49 
51 
52  public:
53 
58 
63 
64  typedef std::vector<Trajectory> TC;
65  typedef TC::const_iterator TI;
66 
67  public:
68 
71 
74 
77 
79  virtual void setEvent(const edm::Event&);
80 
81  protected:
82 
84 
88 
90  virtual std::vector<TrackCand> makeTkCandCollection(const TrackCand&) = 0;
91 
93  std::vector<TrackCand> chooseRegionalTrackerTracks(const TrackCand&,
94  const std::vector<TrackCand>&);
95 
98 
100  void checkMuonHits(const reco::Track&,
103  std::vector<int>&) const;
104 
107  const std::vector<int>&) const;
108 
111 
114  double scl_x,
115  double scl_y) const;
116 
118  const Trajectory* chooseTrajectory(const std::vector<Trajectory*>&, int) const;
119 
121  double trackProbability(const Trajectory&) const;
122 
124  void printHits(const ConstRecHitContainer&) const;
125 
127  void addTraj(TrackCand&){}
128 
131 
134  getTransientRecHits(const reco::Track&) const;
135 
138 
140  const MuonServiceProxy* service() const { return theService; }
141 
143 
146  bool barrel_a = ( a->det()->subDetector() == GeomDetEnumerators::DT ||
147  a->det()->subDetector() == GeomDetEnumerators::RPCBarrel );
148 
149  bool barrel_b = ( b->det()->subDetector() == GeomDetEnumerators::DT ||
150  b->det()->subDetector() == GeomDetEnumerators::RPCBarrel );
151 
152  if ( barrel_a && barrel_b ) return a->det()->surface().position().perp() < b->det()->surface().position().perp();
153 
154  else if ( !barrel_a && !barrel_b ) return fabs(a->globalPosition().z()) < fabs(b->globalPosition().z());
155  else if ( barrel_a && !barrel_b ) return true;
156  else if ( !barrel_a && barrel_b ) return false;
157  //shouldn;t really get here in any case (there's some sense to throw here )
158  return false;
159  }
160  };
161 
162  protected:
163 
165  float thePtCut;
166  float thePCut;
167 
168  private:
169 
176  unsigned long long theCacheId_TRH;
178 
183 
185 
188 
191 
193 };
194 #endif
MuonCandidate::CandidateContainer CandidateContainer
void printHits(const ConstRecHitContainer &) const
print all RecHits of a trajectory
GlobalMuonTrackMatcher * trackMatcher() const
ConstRecHitContainer selectTrackerHits(const ConstRecHitContainer &) const
select tracker hits; exclude some tracker hits in the global trajectory
virtual void setEvent(const edm::Event &)
pass the Event to the algo at each event
virtual std::vector< TrackCand > makeTkCandCollection(const TrackCand &)=0
make a TrackCand collection using tracker Track, Trajectory information
void addTraj(TrackCand &)
if TrackCand has only a TrackRef, attempt to add Trajectory*
std::pair< const Trajectory *, reco::TrackRef > TrackCand
GlobalTrajectoryBuilderBase(const edm::ParameterSet &, const MuonServiceProxy *)
constructor with Parameter Set and MuonServiceProxy
TransientTrackingRecHit::RecHitPointer RecHitPointer
edm::ESHandle< TransientTrackingRecHitBuilder > theMuonRecHitBuilder
MuonTransientTrackingRecHit::MuonRecHitContainer MuonRecHitContainer
MuonTransientTrackingRecHit::ConstMuonRecHitContainer ConstMuonRecHitContainer
MuonCandidate::TrajectoryContainer TrajectoryContainer
bool operator()(const TransientTrackingRecHit::ConstRecHitPointer &a, const TransientTrackingRecHit::ConstRecHitPointer &b) const
std::vector< TrackCand > chooseRegionalTrackerTracks(const TrackCand &, const std::vector< TrackCand > &)
choose tracker tracks within region of interest
MuonTransientTrackingRecHit::MuonRecHitPointer MuonRecHitPointer
std::vector< ConstRecHitPointer > RecHitContainer
void fixTEC(ConstRecHitContainer &all, double scl_x, double scl_y) const
rescale errors of outermost TEC RecHit
double trackProbability(const Trajectory &) const
calculate chi2 probability (-ln(P))
MuonTransientTrackingRecHit::ConstMuonRecHitPointer ConstMuonRecHitPointer
virtual ~GlobalTrajectoryBuilderBase()
destructor
MuonDetLayerMeasurements * theLayerMeasurements
TransientTrackingRecHit::RecHitContainer RecHitContainer
TransientTrackingRecHit::ConstRecHitContainer getTransientRecHits(const reco::Track &) const
get transient RecHits of a Track
MuonTrajectoryBuilder::TrajectoryContainer trajectories(const TrajectorySeed &)
dummy implementation, unused in this class
const MuonServiceProxy * service() const
ConstRecHitContainer selectMuonHits(const Trajectory &, const std::vector< int > &) const
select muon hits compatible with trajectory; check hits in chambers with showers
TransientTrackingRecHit::ConstRecHitContainer ConstRecHitContainer
std::vector< ConstMuonRecHitPointer > ConstMuonRecHitContainer
std::vector< ConstRecHitPointer > ConstRecHitContainer
GlobalMuonTrackMatcher * theTrackMatcher
MuonTrajectoryBuilder::CandidateContainer build(const TrackCand &, MuonTrajectoryBuilder::CandidateContainer &) const
build combined trajectory from sta Track and tracker RecHits
double b
Definition: hdecay.h:120
RectangularEtaPhiTrackingRegion defineRegionOfInterest(const reco::TrackRef &) const
define region of interest with tracker
const Trajectory * chooseTrajectory(const std::vector< Trajectory * > &, int) const
choose final trajectory
double a
Definition: hdecay.h:121
RefitDirection checkRecHitsOrdering(const ConstRecHitContainer &) const
This does nothing now.
TransientTrackingRecHit::ConstRecHitContainer ConstRecHitContainer
TransientTrackingRecHit::ConstRecHitPointer ConstRecHitPointer
MuonTrackingRegionBuilder * theRegionBuilder
std::vector< MuonRecHitPointer > MuonRecHitContainer
edm::ESHandle< TransientTrackingRecHitBuilder > theTrackerRecHitBuilder
void checkMuonHits(const reco::Track &, ConstRecHitContainer &, ConstRecHitContainer &, std::vector< int > &) const
check muon RecHits, calculate chamber occupancy and select hits to be used in the final fit ...