CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
SiPixelPhase1TrackEfficiency.cc
Go to the documentation of this file.
1 // -*- C++ -*-
2 //
3 // Package: SiPixelPhase1TrackEfficiency
4 // Class: SiPixelPhase1TrackEfficiency
5 //
6 
7 // Original Author: Marcel Schneider
8 
14 
30 
38 
40 
41 namespace {
42 
43  class SiPixelPhase1TrackEfficiency final : public SiPixelPhase1Base {
44  enum { VALID, MISSING, INACTIVE, EFFICIENCY, VERTICES };
45 
46  public:
47  explicit SiPixelPhase1TrackEfficiency(const edm::ParameterSet& conf);
48  void analyze(const edm::Event&, const edm::EventSetup&) override;
49 
50  private:
54  edm::EDGetTokenT<TrajTrackAssociationCollection> trajTrackCollectionToken_;
56  bool applyVertexCut_;
57 
63  edm::ESGetToken<PixelClusterParameterEstimator, TkPixelCPERecord> pixelClusterParameterEstimatorToken_;
64 
65  const TrackerTopology* trackerTopology_;
66  const Propagator* trackerPropagator_;
67  const MeasurementEstimator* chi2MeasurementEstimator_;
68  };
69 
70  SiPixelPhase1TrackEfficiency::SiPixelPhase1TrackEfficiency(const edm::ParameterSet& iConfig)
71  : SiPixelPhase1Base(iConfig) //,
72  {
73  tracker_ = consumes<MeasurementTrackerEvent>(iConfig.getParameter<edm::InputTag>("tracker"));
74  tracksToken_ = consumes<reco::TrackCollection>(iConfig.getParameter<edm::InputTag>("tracks"));
75  vtxToken_ = consumes<reco::VertexCollection>(iConfig.getParameter<edm::InputTag>("primaryvertices"));
76  applyVertexCut_ = iConfig.getUntrackedParameter<bool>("VertexCut", true);
77  trajTrackCollectionToken_ =
78  consumes<TrajTrackAssociationCollection>(iConfig.getParameter<edm::InputTag>("trajectoryInput"));
79  clustersToken_ = consumes<edmNew::DetSetVector<SiPixelCluster>>(iConfig.getParameter<edm::InputTag>("clusters"));
80 
81  trackerTopoToken_ = esConsumes<TrackerTopology, TrackerTopologyRcd>();
82  trackerGeomToken_ = esConsumes<TrackerGeometry, TrackerDigiGeometryRecord>();
83  propagatorToken_ = esConsumes<Propagator, TrackingComponentsRecord>(edm::ESInputTag("", "PropagatorWithMaterial"));
84  chi2MeasurementEstimatorBaseToken_ =
85  esConsumes<Chi2MeasurementEstimatorBase, TrackingComponentsRecord>(edm::ESInputTag("", "Chi2"));
86  measurementTrackerToken_ = esConsumes<MeasurementTracker, CkfComponentsRecord>();
87  pixelClusterParameterEstimatorToken_ =
88  esConsumes<PixelClusterParameterEstimator, TkPixelCPERecord>(edm::ESInputTag("", "PixelCPEGeneric"));
89  }
90 
92  if (!checktrigger(iEvent, iSetup, DCS))
93  return;
94 
95  // get geometry
96  edm::ESHandle<TrackerGeometry> tracker = iSetup.getHandle(trackerGeomToken_);
97  assert(tracker.isValid());
98 
99  // get primary vertex
101  iEvent.getByToken(vtxToken_, vertices);
102 
103  // TrackerTopology for module informations
104  edm::ESHandle<TrackerTopology> trackerTopologyHandle = iSetup.getHandle(trackerTopoToken_);
105  trackerTopology_ = trackerTopologyHandle.product();
106 
107  // Tracker propagator for propagating tracks to other layers
108  edm::ESHandle<Propagator> propagatorHandle = iSetup.getHandle(propagatorToken_);
109  std::unique_ptr<Propagator> propagatorUniquePtr(propagatorHandle.product()->clone());
110  trackerPropagator_ = propagatorUniquePtr.get();
111  const_cast<Propagator*>(trackerPropagator_)->setPropagationDirection(oppositeToMomentum);
112 
113  // Measurement estimator
114  edm::ESHandle<Chi2MeasurementEstimatorBase> chi2MeasurementEstimatorHandle =
115  iSetup.getHandle(chi2MeasurementEstimatorBaseToken_);
116  chi2MeasurementEstimator_ = chi2MeasurementEstimatorHandle.product();
117 
118  //Tracker
120  iEvent.getByToken(tracker_, trackerMeas);
121 
122  edm::ESHandle<MeasurementTracker> measurementTrackerHandle = iSetup.getHandle(measurementTrackerToken_);
123 
124  //vertices
125  if (!vertices.isValid())
126  return;
127  histo[VERTICES].fill(vertices->size(), DetId(0), &iEvent);
128  if (applyVertexCut_ && vertices->empty())
129  return;
130 
131  // should be used for weird cuts
132  //const auto primaryVertex = vertices->at(0);
133 
134  // get the map
136  iEvent.getByToken(tracksToken_, tracks);
137  if (!tracks.isValid())
138  return;
139 
140  //new
141  edm::Handle<TrajTrackAssociationCollection> trajTrackCollectionHandle;
142  iEvent.getByToken(trajTrackCollectionToken_, trajTrackCollectionHandle);
143  if (!trajTrackCollectionHandle.isValid())
144  return;
145 
146  //Access Pixel Clusters
148  iEvent.getByToken(clustersToken_, siPixelClusters);
149  if (!siPixelClusters.isValid())
150  return;
151  //
152 
153  edm::ESHandle<PixelClusterParameterEstimator> cpEstimator = iSetup.getHandle(pixelClusterParameterEstimatorToken_);
154  if (!cpEstimator.isValid())
155  return;
156 
157  const PixelClusterParameterEstimator& cpe(*cpEstimator);
158  const TrackerGeometry* tkgeom = &(*tracker);
159 
161 
162  // Hp cut
163  static constexpr int TRACK_QUALITY_HIGH_PURITY_BIT = 2;
164  static constexpr int TRACK_QUALITY_HIGH_PURITY_MASK = 1 << TRACK_QUALITY_HIGH_PURITY_BIT;
165 
166  // Pt cut
167  static constexpr float TRACK_PT_CUT_VAL = 1.0f;
168 
169  // Nstrip cut
170  static constexpr int TRACK_NSTRIP_CUT_VAL = 10;
171 
172  //D0
173  static constexpr std::array<float, 4> TRACK_D0_CUT_BARREL_VAL = {{0.01f, 0.02f, 0.02f, 0.02f}};
174  static constexpr float TRACK_D0_CUT_FORWARD_VAL = 0.05f;
175 
176  //Dz
177  static constexpr float TRACK_DZ_CUT_BARREL_VAL = 0.01f;
178  static constexpr float TRACK_DZ_CUT_FORWARD_VAL = 0.5f;
179 
180  TrajectoryStateOnSurface tsosPXB2;
181  bool valid_layerFrom = false;
182 
183  const GeometricSearchTracker* gst_ = trackerMeas->geometricSearchTracker();
184  const auto* pxbLayer1_ = gst_->pixelBarrelLayers().front();
185  const LayerMeasurements* theLayerMeasurements_ = new LayerMeasurements(*measurementTrackerHandle, *trackerMeas);
186 
187  std::vector<TrajectoryMeasurement> expTrajMeasurements;
188  std::vector<std::pair<int, bool[3]>> eff_pxb1_vector;
189 
190  for (const auto& pair : *trajTrackCollectionHandle) {
191  const edm::Ref<std::vector<Trajectory>> traj = pair.key;
192  const reco::TrackRef track = pair.val;
193 
194  expTrajMeasurements.clear();
195  eff_pxb1_vector.clear();
196  //this cut is needed to be consisten with residuals calculation
197  if (applyVertexCut_ &&
198  (track->pt() < 0.75 || std::abs(track->dxy(vertices->at(0).position())) > 5 * track->dxyError()))
199  continue;
200 
201  bool isBpixtrack = false;
202  bool isFpixtrack = false;
203  int nStripHits = 0;
204  int nBpixL1Hits = 0;
205  int nBpixL2Hits = 0;
206  int nBpixL3Hits = 0;
207  int nBpixL4Hits = 0;
208  int nFpixD1Hits = 0;
209  int nFpixD2Hits = 0;
210  int nFpixD3Hits = 0;
211  bool passcuts = true;
212 
213  // first, look at the full track to see whether it is good
214  // auto const & trajParams = track.extra()->trajParams();
215 
216  auto hb = track->recHitsBegin();
217  for (unsigned int h = 0; h < track->recHitsSize(); h++) {
218  auto hit = *(hb + h);
219  if (!hit->isValid())
220  continue;
221 
222  DetId id = hit->geographicalId();
223  uint32_t subdetid = (id.subdetId());
224 
225  //Check the location of valid hit
226  if (subdetid == PixelSubdetector::PixelBarrel && hit->isValid()) {
227  isBpixtrack = true;
228  if (trackerTopology_->pxbLayer(id) == 1)
229  nBpixL1Hits++;
230  else if (trackerTopology_->pxbLayer(id) == 2)
231  nBpixL2Hits++;
232  else if (trackerTopology_->pxbLayer(id) == 3)
233  nBpixL3Hits++;
234  else if (trackerTopology_->pxbLayer(id) == 4)
235  nBpixL4Hits++;
236  } else if (subdetid == PixelSubdetector::PixelEndcap && hit->isValid()) {
237  isFpixtrack = true;
238  if (trackerTopology_->pxfDisk(id) == 1)
239  nFpixD1Hits++;
240  else if (trackerTopology_->pxfDisk(id) == 2)
241  nFpixD2Hits++;
242  else if (trackerTopology_->pxfDisk(id) == 3)
243  nFpixD3Hits++;
244  }
245 
246  // count strip hits
247  if (subdetid == StripSubdetector::TIB || subdetid == StripSubdetector::TOB ||
248  subdetid == StripSubdetector::TID || subdetid == StripSubdetector::TEC)
249  nStripHits++;
250  }
251 
252  if (!isBpixtrack && !isFpixtrack)
253  continue;
254 
255  // Hp cut
256  if (!((track->qualityMask() & TRACK_QUALITY_HIGH_PURITY_MASK) >> TRACK_QUALITY_HIGH_PURITY_BIT))
257  passcuts = false;
258 
259  // Pt cut
260  if (!(TRACK_PT_CUT_VAL < track->pt()))
261  passcuts = false;
262 
263  // Nstrip cut
264  if (!(TRACK_NSTRIP_CUT_VAL < nStripHits))
265  passcuts = false;
266 
267  // then, look at each hit
268  for (unsigned int h = 0; h < track->recHitsSize(); h++) {
269  bool passcuts_hit = true;
270  auto hit = *(hb + h);
271 
272  DetId id = hit->geographicalId();
273  uint32_t subdetid = (id.subdetId());
274  if (subdetid != PixelSubdetector::PixelBarrel && subdetid != PixelSubdetector::PixelEndcap)
275  continue;
276 
277  bool isHitValid = hit->getType() == TrackingRecHit::valid;
278  bool isHitMissing = hit->getType() == TrackingRecHit::missing;
279  bool isHitInactive = hit->getType() == TrackingRecHit::inactive;
280 
281  //D0
282  if (subdetid == PixelSubdetector::PixelBarrel) {
283  if (!((std::abs(track->dxy(vertices->at(0).position())) * -1.0) <
284  TRACK_D0_CUT_BARREL_VAL[trackerTopology_->pxbLayer(id) - 1]))
285  passcuts_hit = false;
286  } else if (subdetid == PixelSubdetector::PixelEndcap) {
287  if (!((std::abs(track->dxy(vertices->at(0).position())) * -1.0) < TRACK_D0_CUT_FORWARD_VAL))
288  passcuts_hit = false;
289  }
290 
291  //Dz
292  if (subdetid == PixelSubdetector::PixelBarrel) {
293  if (!(std::abs(track->dz(vertices->at(0).position())) < TRACK_DZ_CUT_BARREL_VAL))
294  passcuts_hit = false;
295  } else if (subdetid == PixelSubdetector::PixelEndcap) {
296  if (!(std::abs(track->dz(vertices->at(0).position())) < TRACK_DZ_CUT_FORWARD_VAL))
297  passcuts_hit = false;
298  }
299 
300  // Pixhit cut
301  if (subdetid == PixelSubdetector::PixelBarrel) {
302  if (trackerTopology_->pxbLayer(id) == 1) {
303  if (!((nBpixL2Hits > 0 && nBpixL3Hits > 0 && nBpixL4Hits > 0) ||
304  (nBpixL2Hits > 0 && nBpixL3Hits > 0 && nFpixD1Hits > 0) ||
305  (nBpixL2Hits > 0 && nFpixD1Hits > 0 && nFpixD2Hits > 0) ||
306  (nFpixD1Hits > 0 && nFpixD2Hits > 0 && nFpixD3Hits > 0)))
307  passcuts_hit = false;
308  } else if (trackerTopology_->pxbLayer(id) == 2) {
309  if (!((nBpixL1Hits > 0 && nBpixL3Hits > 0 && nBpixL4Hits > 0) ||
310  (nBpixL1Hits > 0 && nBpixL3Hits > 0 && nFpixD1Hits > 0) ||
311  (nBpixL1Hits > 0 && nFpixD1Hits > 0 && nFpixD2Hits > 0)))
312  passcuts_hit = false;
313  } else if (trackerTopology_->pxbLayer(id) == 3) {
314  if (!((nBpixL1Hits > 0 && nBpixL2Hits > 0 && nBpixL4Hits > 0) ||
315  (nBpixL1Hits > 0 && nBpixL2Hits > 0 && nFpixD1Hits > 0)))
316  passcuts_hit = false;
317  } else if (trackerTopology_->pxbLayer(id) == 4)
318  if (!((nBpixL1Hits > 0 && nBpixL2Hits > 0 && nBpixL3Hits > 0)))
319  passcuts_hit = false;
320  } else if (subdetid == PixelSubdetector::PixelEndcap) {
321  if (trackerTopology_->pxfDisk(id) == 1) {
322  if (!((nBpixL1Hits > 0 && nBpixL2Hits > 0 && nBpixL3Hits > 0) ||
323  (nBpixL1Hits > 0 && nBpixL2Hits > 0 && nFpixD2Hits > 0) ||
324  (nBpixL1Hits > 0 && nFpixD2Hits > 0 && nFpixD3Hits > 0)))
325  passcuts_hit = false;
326  } else if (trackerTopology_->pxfDisk(id) == 2) {
327  if (!((nBpixL1Hits > 0 && nBpixL2Hits > 0 && nFpixD1Hits > 0) ||
328  (nBpixL1Hits > 0 && nFpixD1Hits > 0 && nFpixD3Hits > 0)))
329  passcuts_hit = false;
330  } else if (trackerTopology_->pxfDisk(id) == 3) {
331  if (!((nBpixL1Hits > 0 && nFpixD1Hits > 0 && nFpixD2Hits > 0)))
332  passcuts_hit = false;
333  }
334  }
335  /*
336  //Fiducial Cut - will work on it later
337  const SiPixelRecHit* pixhit = dynamic_cast<const SiPixelRecHit*>(hit);
338  const PixelGeomDetUnit* geomdetunit = dynamic_cast<const PixelGeomDetUnit*> ( tracker->idToDet(id) );
339  const PixelTopology& topol = geomdetunit->specificTopology();
340 
341  LocalPoint lp;
342  if (pixhit) {
343  lp = pixhit->localPosition();
344  }
345 
346  MeasurementPoint mp = topol.measurementPosition(lp);
347  int row = (int) mp.x() % 80;
348  int col = (int) mp.y() % 52;
349 
350  int centerrow = 40;
351  int centercol = 26;
352 
353  if (!((col < (centercol + 10)) && (col > (centercol - 10)) && (row < (centerrow + 10)) && (row > (centerrow -10 )))) passcuts_hit = false;
354  */
355 
356  if (passcuts_hit && passcuts) {
357  if (!(subdetid == PixelSubdetector::PixelBarrel && trackerTopology_->pxbLayer(id) == 1)) {
358  if (isHitValid) {
359  histo[VALID].fill(id, &iEvent);
360  histo[EFFICIENCY].fill(1, id, &iEvent);
361  }
362  if (isHitMissing) {
363  histo[MISSING].fill(id, &iEvent);
364  histo[EFFICIENCY].fill(0, id, &iEvent);
365  }
366  if (isHitInactive) {
367  histo[INACTIVE].fill(id, &iEvent);
368  }
369  }
370  }
371  }
372 
374  valid_layerFrom = false;
375 
376  //propagation only from PXB2 and PXD1, more cuts later
377  for (const auto& tm : traj->measurements()) {
378  if (tm.recHit().get() && tm.recHitR().isValid()) {
379  DetId where = tm.recHitR().geographicalId();
380  int source_det = where.subdetId();
381 
382  if (source_det == PixelSubdetector::SubDetector::PixelBarrel) {
383  int source_layer = trackerTopology_->pxbLayer(where);
384  if (source_layer == 2) {
385  if (tm.updatedState().isValid()) {
386  tsosPXB2 = tm.updatedState();
387  valid_layerFrom = true;
388  }
389  }
390  }
391 
392  if (source_det == PixelSubdetector::SubDetector::PixelEndcap) {
393  int source_layer = trackerTopology_->pxfDisk(where);
394  if (source_layer == 1) {
395  if (tm.updatedState().isValid()) {
396  tsosPXB2 = tm.updatedState();
397  valid_layerFrom = true;
398  }
399  }
400  }
401  }
402  } //uodated tsosPXB2 here
403 
404  if (!valid_layerFrom)
405  continue;
406  if (!tsosPXB2.isValid())
407  continue;
408 
409  //propagation A: Calculate the efficiency by the distance to the closest cluster
410  expTrajMeasurements =
411  theLayerMeasurements_->measurements(*pxbLayer1_, tsosPXB2, *trackerPropagator_, *chi2MeasurementEstimator_);
412  auto compDets = pxbLayer1_->compatibleDets(tsosPXB2, *trackerPropagator_, *chi2MeasurementEstimator_);
413  std::pair<int, bool[3]> eff_map;
414 
415  //Fiducial Cut, only calculate the efficiency of the central pixels
416  for (uint p = 0; p < expTrajMeasurements.size(); p++) {
417  bool valid = false;
418  bool missing = false;
419  bool passcuts_hit = true;
420  TrajectoryMeasurement pxb1TM(expTrajMeasurements[p]);
421  const auto& pxb1Hit = pxb1TM.recHit();
422  bool inactive = (pxb1Hit->getType() == TrackingRecHit::inactive);
423  int detidHit = pxb1Hit->geographicalId();
424  if (detidHit == 0)
425  continue;
426 
427  const SiPixelRecHit* pixhit = dynamic_cast<const SiPixelRecHit*>(pxb1Hit->hit());
428  const PixelGeomDetUnit* geomdetunit = dynamic_cast<const PixelGeomDetUnit*>(tracker->idToDet(detidHit));
429  const PixelTopology& topol = geomdetunit->specificTopology();
430 
431  if (!pixhit)
432  continue;
433 
434  LocalPoint lp = pixhit->localPosition();
435  MeasurementPoint mp = topol.measurementPosition(lp);
436  const int nRows = topol.rowsperroc();
437  const int nColumns = topol.colsperroc();
438  int row = (int)mp.x() % nRows;
439  int col = (int)mp.y() % nColumns;
440 
441  int centerrow = nRows / 2;
442  int centercol = nColumns / 2;
443 
444  if (!((col < (centercol + 10)) && (col > (centercol - 10)) && (row < (centerrow + 10)) &&
445  (row > (centerrow - 10))))
446  passcuts_hit = false;
447 
448  //Access the distance to the closest cluster
449  for (const auto& detAndState : compDets) {
450  const auto& pXb1_lpos = detAndState.second.localPosition();
451  if (pxb1Hit->geographicalId().rawId() != detAndState.first->geographicalId().rawId())
452  continue;
453  int detid = detAndState.first->geographicalId().rawId();
454 
456  iter_cl != siPixelClusters->end();
457  iter_cl++) {
458  DetId detId(iter_cl->id());
459  float minD[2], minDist = 10000;
460  minD[0] = minD[1] = 10000.;
461  if (detId.rawId() != detAndState.first->geographicalId().rawId())
462  continue;
463 
464  const PixelGeomDetUnit* pixdet = (const PixelGeomDetUnit*)tkgeom->idToDetUnit(detId);
465  edmNew::DetSet<SiPixelCluster>::const_iterator itCluster = iter_cl->begin();
466  if (passcuts_hit) {
467  for (; itCluster != iter_cl->end(); ++itCluster) {
468  LocalPoint lp(itCluster->x(), itCluster->y(), 0.);
469  PixelClusterParameterEstimator::ReturnType params = cpe.getParameters(*itCluster, *pixdet);
470  lp = std::get<0>(params);
471 
472  float Xdist = abs(lp.x() - pXb1_lpos.x());
473  float Ydist = abs(lp.y() - pXb1_lpos.y());
474  float dist = sqrt(Xdist * Xdist + Ydist * Ydist);
475  if (dist < minDist) {
476  minDist = dist;
477  minD[0] = Xdist;
478  minD[1] = Ydist;
479  }
480  }
481 
482  if ((minD[0] < 0.02) && (minD[1] < 0.02)) {
483  valid = true;
484  missing = false;
485  inactive = false;
486  } else if (inactive) {
487  valid = false;
488  missing = false;
489  } else {
490  missing = true;
491  valid = false;
492  }
493  }
494  }
495 
496  //cuts: exactly the same as for other hits but assuming PXB1
497 
498  //D0
499  if (!((std::abs(track->dxy(vertices->at(0).position())) * -1.0) <
500  TRACK_D0_CUT_BARREL_VAL[trackerTopology_->pxbLayer(detid) - 1]))
501  passcuts_hit = false;
502  //Dz
503  if (!(std::abs(track->dz(vertices->at(0).position())) < TRACK_DZ_CUT_BARREL_VAL))
504  passcuts_hit = false;
505  // Pixhit cut
506  if (!((nBpixL2Hits > 0 && nBpixL3Hits > 0 && nBpixL4Hits > 0) ||
507  (nBpixL2Hits > 0 && nBpixL3Hits > 0 && nFpixD1Hits > 0) ||
508  (nBpixL2Hits > 0 && nFpixD1Hits > 0 && nFpixD2Hits > 0) ||
509  (nFpixD1Hits > 0 && nFpixD2Hits > 0 && nFpixD3Hits > 0)))
510  passcuts_hit = false;
511  bool found_det = false;
512 
513  if (passcuts && passcuts_hit) {
514  for (unsigned int i_eff = 0; i_eff < eff_pxb1_vector.size(); i_eff++) {
515  //in case found hit in the same det, take only the valid hit
516  if (eff_pxb1_vector[i_eff].first == detid) {
517  found_det = true;
518  if (eff_pxb1_vector[i_eff].second[0] == false && valid == true) {
519  eff_pxb1_vector[i_eff].second[0] = valid;
520  eff_pxb1_vector[i_eff].second[1] = missing;
521  eff_pxb1_vector[i_eff].second[2] = inactive;
522  }
523  }
524  }
525  if (!found_det) {
526  eff_map.first = detid;
527  eff_map.second[0] = valid;
528  eff_map.second[1] = missing;
529  eff_map.second[2] = inactive;
530  eff_pxb1_vector.push_back(eff_map);
531  }
532  }
533  }
534  }
535 
536  if (eff_pxb1_vector.size() == 1) {
537  //eff map is filled -> decide what to do for double hits, ie eff_pxb1_vector.size>1 ... if 1 just use MISSING and VALID as usual
538 
539  if (eff_pxb1_vector[0].second[0]) {
540  histo[VALID].fill(eff_pxb1_vector[0].first, &iEvent);
541  histo[EFFICIENCY].fill(1, eff_pxb1_vector[0].first, &iEvent);
542  }
543  if (eff_pxb1_vector[0].second[1]) {
544  histo[MISSING].fill(eff_pxb1_vector[0].first, &iEvent);
545  histo[EFFICIENCY].fill(0, eff_pxb1_vector[0].first, &iEvent);
546  }
547 
548  if (eff_pxb1_vector[0].second[2]) {
549  histo[INACTIVE].fill(eff_pxb1_vector[0].first, &iEvent);
550  }
551  }
552 
553  } //trajTrackHandle
554 
555  histo[VALID].executePerEventHarvesting(&iEvent);
556  histo[MISSING].executePerEventHarvesting(&iEvent);
557  histo[INACTIVE].executePerEventHarvesting(&iEvent);
558  }
559 
560 } // namespace
561 
562 DEFINE_FWK_MODULE(SiPixelPhase1TrackEfficiency);
static constexpr auto TEC
T getParameter(std::string const &) const
Definition: ParameterSet.h:307
static const char tracker_[]
const TrackerGeomDet * idToDetUnit(DetId) const override
Return the pointer to the GeomDetUnit corresponding to a given DetId.
virtual int rowsperroc() const =0
virtual Propagator * clone() const =0
data_type const * const_iterator
Definition: DetSetNew.h:31
T x() const
Definition: PV2DBase.h:43
assert(be >=bs)
example_stream void analyze(const edm::Event &, const edm::EventSetup &) override
key_type key() const
Accessor for product key.
Definition: Ref.h:250
missing
Definition: combine.py:5
T y() const
Definition: PV2DBase.h:44
U second(std::pair< T, U > const &p)
T x() const
Definition: PV3DBase.h:59
T y() const
Definition: PV3DBase.h:60
int iEvent
Definition: GenABIO.cc:224
T const * product() const
Definition: ESHandle.h:86
virtual int colsperroc() const =0
T sqrt(T t)
Definition: SSEVec.h:19
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
constexpr int subdetId() const
get the contents of the subdetector field (not cast into any detector&#39;s numbering enum) ...
Definition: DetId.h:48
virtual MeasurementPoint measurementPosition(const LocalPoint &) const =0
static constexpr auto TOB
ESHandle< T > getHandle(const ESGetToken< T, R > &iToken) const
Definition: EventSetup.h:130
bool isValid() const
Definition: ESHandle.h:44
std::vector< BarrelDetLayer const * > const & pixelBarrelLayers() const
boost::transform_iterator< IterHelp, const_IdIter > const_iterator
Definition: DetId.h:17
static constexpr auto TIB
void analyze(edm::Event const &e, edm::EventSetup const &) override=0
bool isValid() const
Definition: HandleBase.h:70
const GeometricSearchTracker * geometricSearchTracker() const
iterator end()
Definition: DetSetNew.h:53
virtual const PixelTopology & specificTopology() const
Returns a reference to the pixel proxy topology.
LocalPoint localPosition() const override
col
Definition: cuy.py:1009
std::vector< TrajectoryMeasurement > measurements(const DetLayer &layer, const TrajectoryStateOnSurface &startingState, const Propagator &prop, const MeasurementEstimator &est) const
The Signals That Services Can Subscribe To This is based on ActivityRegistry h
Helper function to determine trigger accepts.
Definition: Activities.doc:4
std::tuple< LocalPoint, LocalError, SiPixelRecHitQuality::QualWordType > ReturnType
static constexpr auto TID
Our base class.
Definition: SiPixelRecHit.h:23