CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
SiPixelPhase1RecHits.cc
Go to the documentation of this file.
1 // -*- C++ -*-
2 //
3 // Package: SiPixelPhase1RecHits
4 // Class: SiPixelPhase1RecHits
5 //
6 
7 // Original Author: Marcel Schneider
8 
21 
23  SiPixelPhase1Base(iConfig)
24 {
25  srcToken_ = consumes<reco::TrackCollection>(iConfig.getParameter<edm::InputTag>("src"));
26 
27  offlinePrimaryVerticesToken_ = consumes<reco::VertexCollection>(std::string("offlinePrimaryVertices"));
28 
29  onlyValid_=iConfig.getParameter<bool>("onlyValidHits");
30 
31  applyVertexCut_=iConfig.getUntrackedParameter<bool>("VertexCut",true);
32 
33 }
34 
36  if( !checktrigger(iEvent,iSetup,DCS) ) return;
37 
39  iSetup.get<TrackerDigiGeometryRecord>().get(tracker);
40  assert(tracker.isValid());
41 
43  iEvent.getByToken( srcToken_, tracks);
44  if (!tracks.isValid()) return;
45 
47  iEvent.getByToken(offlinePrimaryVerticesToken_, vertices);
48 
49  if (applyVertexCut_ && (!vertices.isValid() || vertices->empty())) return;
50 
51 
52  for (auto const & track : *tracks) {
53 
54  if (applyVertexCut_ && (track.pt() < 0.75 || std::abs( track.dxy(vertices->at(0).position()) ) > 5*track.dxyError())) continue;
55 
56  bool isBpixtrack = false, isFpixtrack = false;
57 
58  auto const & trajParams = track.extra()->trajParams();
59  auto hb = track.recHitsBegin();
60  for(unsigned int h=0;h<track.recHitsSize();h++){
61 
62  auto hit = *(hb+h);
63  if(!hit->isValid()) continue;
64 
65  DetId id = hit->geographicalId();
66  uint32_t subdetid = (id.subdetId());
67 
68  if (subdetid == PixelSubdetector::PixelBarrel) isBpixtrack = true;
69  if (subdetid == PixelSubdetector::PixelEndcap) isFpixtrack = true;
70  }
71 
72  if (!isBpixtrack && !isFpixtrack) continue;
73 
74  // then, look at each hit
75  for(unsigned int h=0;h<track.recHitsSize();h++){
76  auto rechit = *(hb+h);
77 
78  if(!rechit->isValid()) continue;
79 
80  //continue if not a Pixel recHit
81  DetId id = rechit->geographicalId();
82  uint32_t subdetid = (id.subdetId());
83 
84  if ( subdetid != PixelSubdetector::PixelBarrel
85  && subdetid != PixelSubdetector::PixelEndcap) continue;
86 
87  bool isHitValid = rechit->getType()==TrackingRecHit::valid;
88  if (onlyValid_ && !isHitValid) continue; //useful to run on cosmics where the TrackEfficiency plugin is not used
89 
90  const SiPixelRecHit* prechit = dynamic_cast<const SiPixelRecHit*>(rechit);//to be used to get the associated cluster and the cluster probability
91 
92  int sizeX=0, sizeY=0;
93 
94  if (isHitValid){
95  SiPixelRecHit::ClusterRef const& clust = prechit->cluster();
96  sizeX = (*clust).sizeX();
97  sizeY = (*clust).sizeY();
98  }
99 
100  const PixelGeomDetUnit* geomdetunit = dynamic_cast<const PixelGeomDetUnit*> ( tracker->idToDet(id) );
101  const PixelTopology& topol = geomdetunit->specificTopology();
102 
103  LocalPoint lp = trajParams[h].position();
104  MeasurementPoint mp = topol.measurementPosition(lp);
105 
106  int row = (int) mp.x();
107  int col = (int) mp.y();
108 
109  float rechit_x = lp.x();
110  float rechit_y = lp.y();
111 
112  LocalError lerr = rechit->localPositionError();
113  float lerr_x = sqrt(lerr.xx());
114  float lerr_y = sqrt(lerr.yy());
115 
116  histo[NRECHITS].fill(id, &iEvent, col, row); //in general a inclusive counter of missing/valid/inactive hits
117 
118  if (isHitValid){
119  histo[CLUST_X].fill(sizeX, id, &iEvent, col, row);
120  histo[CLUST_Y].fill(sizeY, id, &iEvent, col, row);
121  }
122 
123  histo[ERROR_X].fill(lerr_x, id, &iEvent);
124  histo[ERROR_Y].fill(lerr_y, id, &iEvent);
125 
126  histo[POS].fill(rechit_x, rechit_y, id, &iEvent);
127 
128  if (isHitValid){
129  double clusterProbability= prechit->clusterProbability(0);
130  if (clusterProbability > 0)
131  histo[CLUSTER_PROB].fill(log10(clusterProbability), id, &iEvent);
132  }
133  }
134  }
135 
136  histo[NRECHITS].executePerEventHarvesting(&iEvent);
137 }
138 
140 
T getParameter(std::string const &) const
T getUntrackedParameter(std::string const &, T const &) const
float xx() const
Definition: LocalError.h:24
float clusterProbability(unsigned int flags=0) const
T y() const
Definition: PV2DBase.h:46
FWCore Framework interface EventSetupRecordImplementation h
Helper function to determine trigger accepts.
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:519
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:17
SiPixelPhase1RecHits(const edm::ParameterSet &conf)
T y() const
Definition: PV3DBase.h:63
bool checktrigger(const edm::Event &iEvent, const edm::EventSetup &iSetup, const unsigned trgidx) const
int iEvent
Definition: GenABIO.cc:230
float yy() const
Definition: LocalError.h:26
T sqrt(T t)
Definition: SSEVec.h:18
virtual MeasurementPoint measurementPosition(const LocalPoint &) const =0
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
bool isValid() const
Definition: HandleBase.h:74
Definition: DetId.h:18
ClusterRef cluster() const
Definition: SiPixelRecHit.h:49
const T & get() const
Definition: EventSetup.h:58
virtual const PixelTopology & specificTopology() const
Returns a reference to the pixel proxy topology.
std::vector< HistogramManager > histo
col
Definition: cuy.py:1008
const TrackerGeomDet * idToDet(DetId) const override
bool isValid() const
Definition: ESHandle.h:47
T x() const
Definition: PV2DBase.h:45
T x() const
Definition: PV3DBase.h:62
edm::EDGetTokenT< reco::VertexCollection > offlinePrimaryVerticesToken_
void analyze(const edm::Event &, const edm::EventSetup &) override
edm::EDGetTokenT< reco::TrackCollection > srcToken_
Our base class.
Definition: SiPixelRecHit.h:23