CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
SiStripElectronSeedGenerator.h
Go to the documentation of this file.
1 #ifndef SiStripElectronSeedGenerator_H
2 #define SiStripElectronSeedGenerator_H
3 
19 
24 
30 
36 
42 
46 
51 
56 
59 
61 class KFUpdator;
62 class MeasurementTracker;
64 class NavigationSchool;
65 
67 {
68 public:
69 
70  struct Tokens {
73  };
74 
79 
81  const Tokens&);
82 
84 
85  void setupES(const edm::EventSetup& setup);
86  void run(edm::Event&, const edm::EventSetup& setup,
89 
90 private:
91  double normalPhi(double phi) const {
92  while (phi > 2.* M_PI) { phi -= 2.*M_PI; }
93  while (phi < 0) { phi += 2.*M_PI; }
94  return phi;
95  }
96 
97  double phiDiff(double phi1, double phi2){
98  double result = normalPhi(phi1) - normalPhi(phi2);
99  if(result > M_PI) result -= 2*M_PI;
100  if(result < -M_PI) result += 2*M_PI;
101  return result;
102  }
103 
104  double unwrapPhi(double phi) const {
105  while (phi > M_PI) { phi -= 2.*M_PI; }
106  while (phi < -M_PI) { phi += 2.*M_PI; }
107  return phi;
108  }
109 
111  const MeasurementTrackerEvent &trackerData,
113 
114  int whichSubdetector(std::vector<const SiStripMatchedRecHit2D*>::const_iterator hit);
115 
116  bool preselection(GlobalPoint position,GlobalPoint superCluster,double phiVsRSlope, int hitLayer);
117  //hitLayer: 1 = TIB, 2 = TID, 3 = TEC, 4 = Mono
118 
119  bool checkHitsAndTSOS(std::vector<const SiStripMatchedRecHit2D*>::const_iterator hit1,
120  std::vector<const SiStripMatchedRecHit2D*>::const_iterator hit2,
121  double scr,double scz,double pT,double scEta);
122 
123  bool altCheckHitsAndTSOS(std::vector<const SiStripMatchedRecHit2D*>::const_iterator hit1,
124  std::vector<const SiStripRecHit2D*>::const_iterator hit2,
125  double scr,double scz,double pT,double scEta);
126 
127  const SiStripMatchedRecHit2D* matchedHitConverter(ConstRecHitPointer crhp);
128  const SiStripRecHit2D* backupHitConverter(ConstRecHitPointer crhp);
129 
130  std::vector<bool> useDetLayer(double scEta);
131 
137 
141 
146 
149 
150  // member vectors to hold the good hits found between hit selection and combinatorics
151  std::vector<const SiStripMatchedRecHit2D*> layer1Hits_;
152  std::vector<const SiStripMatchedRecHit2D*> layer2Hits_;
153  std::vector<const SiStripRecHit2D*> backupLayer2Hits_;
154 
156 
157  unsigned long long cacheIDMagField_;
158  unsigned long long cacheIDCkfComp_;
159  unsigned long long cacheIDTrkGeom_;
160 
176  double tidEtaUsage_;
181 
182 };
183 
184 #endif // SiStripElectronSeedGenerator_H
185 
186 
bool checkHitsAndTSOS(std::vector< const SiStripMatchedRecHit2D * >::const_iterator hit1, std::vector< const SiStripMatchedRecHit2D * >::const_iterator hit2, double scr, double scz, double pT, double scEta)
std::vector< bool > useDetLayer(double scEta)
const MeasurementTracker * theMeasurementTracker
edm::ESHandle< MagneticField > theMagField
TransientTrackingRecHit::ConstRecHitPointer ConstRecHitPointer
Chi2MeasurementEstimator * theEstimator
const SiStripMatchedRecHit2D * matchedHitConverter(ConstRecHitPointer crhp)
edm::EDGetTokenT< reco::BeamSpot > beamSpotTag_
std::vector< ConstRecHitPointer > RecHitContainer
static int position[TOTALCHAMBERS][3]
Definition: ReadPGInfo.cc:509
edm::EDGetTokenT< MeasurementTrackerEvent > theMeasurementTrackerEventTag
edm::EDGetTokenT< MeasurementTrackerEvent > token_mte
TransientTrackingRecHit::RecHitPointer RecHitPointer
tuple result
Definition: query.py:137
edm::OwnVector< TrackingRecHit > PRecHitContainer
SiStripElectronSeedGenerator(const edm::ParameterSet &, const Tokens &)
std::vector< ElectronSeed > ElectronSeedCollection
collection of ElectronSeed objects
edm::ESHandle< TrackerGeometry > trackerGeometryHandle
edm::Handle< reco::BeamSpot > theBeamSpot
int whichSubdetector(std::vector< const SiStripMatchedRecHit2D * >::const_iterator hit)
void findSeedsFromCluster(edm::Ref< reco::SuperClusterCollection >, edm::Handle< reco::BeamSpot >, const MeasurementTrackerEvent &trackerData, reco::ElectronSeedCollection &)
bool preselection(GlobalPoint position, GlobalPoint superCluster, double phiVsRSlope, int hitLayer)
#define M_PI
Definition: BFit3D.cc:3
const SiStripRecHit2D * backupHitConverter(ConstRecHitPointer crhp)
bool altCheckHitsAndTSOS(std::vector< const SiStripMatchedRecHit2D * >::const_iterator hit1, std::vector< const SiStripRecHit2D * >::const_iterator hit2, double scr, double scz, double pT, double scEta)
void setupES(const edm::EventSetup &setup)
const SiStripRecHitMatcher * theMatcher_
std::vector< const SiStripMatchedRecHit2D * > layer2Hits_
double phiDiff(double phi1, double phi2)
edm::ESHandle< MeasurementTracker > measurementTrackerHandle
std::vector< const SiStripRecHit2D * > backupLayer2Hits_
void run(edm::Event &, const edm::EventSetup &setup, const edm::Handle< reco::SuperClusterCollection > &, reco::ElectronSeedCollection &)
std::vector< const SiStripMatchedRecHit2D * > layer1Hits_
TransientTrackingRecHit::RecHitContainer RecHitContainer
void setup(std::vector< TH2F > &depth, std::string name, std::string units="")
edm::EDGetTokenT< reco::BeamSpot > token_bs
Definition: DDAxes.h:10