CMS 3D CMS Logo
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
AlignmentTrackSelector.cc
Go to the documentation of this file.
1 #include "Alignment/CommonAlignmentProducer/interface/AlignmentTrackSelector.h"
2 
3 #include "FWCore/MessageLogger/interface/MessageLogger.h"
4 #include "FWCore/ParameterSet/interface/ParameterSet.h"
5 
6 #include "FWCore/Framework/interface/Event.h"
7 
8 #include "DataFormats/TrackingRecHit/interface/TrackingRecHit.h"
9 #include "DataFormats/TrackerRecHit2D/interface/SiPixelRecHit.h"
10 #include "DataFormats/TrackerRecHit2D/interface/SiStripRecHit1D.h"
11 #include "DataFormats/TrackerRecHit2D/interface/SiStripRecHit2D.h"
12 #include "DataFormats/TrackerRecHit2D/interface/ProjectedSiStripRecHit2D.h"
13 #include "DataFormats/TrackerRecHit2D/interface/SiStripMatchedRecHit2D.h"
14 #include "DataFormats/TrackerRecHit2D/interface/SiTrackerMultiRecHit.h"
15 #include "DataFormats/SiStripCluster/interface/SiStripCluster.h"
16 #include "DataFormats/TrackerRecHit2D/interface/SiStripMatchedRecHit2DCollection.h"
17 #include "DataFormats/TrackerRecHit2D/interface/SiStripRecHit2DCollection.h"
18 #include "DataFormats/Alignment/interface/AlignmentClusterFlag.h"
19 #include "DataFormats/Alignment/interface/AliClusterValueMap.h"
20 
21 #include "DataFormats/DetId/interface/DetId.h"
22 #include "DataFormats/SiStripDetId/interface/SiStripDetId.h"
23 #include "DataFormats/SiPixelDetId/interface/PixelSubdetector.h"
24 #include "DataFormats/SiStripDetId/interface/TIDDetId.h"
25 #include "DataFormats/SiStripDetId/interface/TECDetId.h"
26 #include "DataFormats/SiPixelDetId/interface/PXFDetId.h"
27 const int kBPIX = PixelSubdetector::PixelBarrel;
28 const int kFPIX = PixelSubdetector::PixelEndcap;
29 
30 // constructor ----------------------------------------------------------------
31 
32 AlignmentTrackSelector::AlignmentTrackSelector(const edm::ParameterSet & cfg) :
33  applyBasicCuts_( cfg.getParameter<bool>( "applyBasicCuts" ) ),
34  applyNHighestPt_( cfg.getParameter<bool>( "applyNHighestPt" ) ),
35  applyMultiplicityFilter_( cfg.getParameter<bool>( "applyMultiplicityFilter" ) ),
36  seedOnlyFromAbove_( cfg.getParameter<int>( "seedOnlyFrom" ) ),
37  applyIsolation_( cfg.getParameter<bool>( "applyIsolationCut" ) ),
38  chargeCheck_( cfg.getParameter<bool>( "applyChargeCheck" ) ),
39  nHighestPt_( cfg.getParameter<int>( "nHighestPt" ) ),
40  minMultiplicity_ ( cfg.getParameter<int>( "minMultiplicity" ) ),
41  maxMultiplicity_ ( cfg.getParameter<int>( "maxMultiplicity" ) ),
42  multiplicityOnInput_ ( cfg.getParameter<bool>( "multiplicityOnInput" ) ),
43  ptMin_( cfg.getParameter<double>( "ptMin" ) ),
44  ptMax_( cfg.getParameter<double>( "ptMax" ) ),
45  pMin_( cfg.getParameter<double>( "pMin" ) ),
46  pMax_( cfg.getParameter<double>( "pMax" ) ),
47  etaMin_( cfg.getParameter<double>( "etaMin" ) ),
48  etaMax_( cfg.getParameter<double>( "etaMax" ) ),
49  phiMin_( cfg.getParameter<double>( "phiMin" ) ),
50  phiMax_( cfg.getParameter<double>( "phiMax" ) ),
51  nHitMin_( cfg.getParameter<double>( "nHitMin" ) ),
52  nHitMax_( cfg.getParameter<double>( "nHitMax" ) ),
53  chi2nMax_( cfg.getParameter<double>( "chi2nMax" ) ),
54  d0Min_( cfg.getParameter<double>( "d0Min" ) ),
55  d0Max_( cfg.getParameter<double>( "d0Max" ) ),
56  dzMin_( cfg.getParameter<double>( "dzMin" ) ),
57  dzMax_( cfg.getParameter<double>( "dzMax" ) ),
58  theCharge_( cfg.getParameter<int>( "theCharge" ) ),
59  minHitChargeStrip_( cfg.getParameter<double>( "minHitChargeStrip" ) ),
60  minHitIsolation_( cfg.getParameter<double>( "minHitIsolation" ) ),
61  rphirecHitsTag_( cfg.getParameter<edm::InputTag>("rphirecHits") ),
62  matchedrecHitsTag_( cfg.getParameter<edm::InputTag>("matchedrecHits") ),
63  countStereoHitAs2D_( cfg.getParameter<bool>( "countStereoHitAs2D" ) ),
64  nHitMin2D_( cfg.getParameter<unsigned int>( "nHitMin2D" ) ),
65  // Ugly to use the same getParameter n times, but this allows const cut variables...
66  minHitsinTIB_(cfg.getParameter<edm::ParameterSet>( "minHitsPerSubDet" ).getParameter<int>( "inTIB" ) ),
67  minHitsinTOB_ (cfg.getParameter<edm::ParameterSet>( "minHitsPerSubDet" ).getParameter<int>( "inTOB" ) ),
68  minHitsinTID_ (cfg.getParameter<edm::ParameterSet>( "minHitsPerSubDet" ).getParameter<int>( "inTID" ) ),
69  minHitsinTEC_ (cfg.getParameter<edm::ParameterSet>( "minHitsPerSubDet" ).getParameter<int>( "inTEC" ) ),
70  minHitsinBPIX_ (cfg.getParameter<edm::ParameterSet>( "minHitsPerSubDet" ).getParameter<int>( "inBPIX" ) ),
71  minHitsinFPIX_ (cfg.getParameter<edm::ParameterSet>( "minHitsPerSubDet" ).getParameter<int>( "inFPIX" ) ),
72  minHitsinPIX_ (cfg.getParameter<edm::ParameterSet>( "minHitsPerSubDet" ).getParameter<int>( "inPIXEL" ) ),
73  minHitsinTIDplus_ (cfg.getParameter<edm::ParameterSet>( "minHitsPerSubDet" ).getParameter<int>( "inTIDplus" ) ),
74  minHitsinTIDminus_ (cfg.getParameter<edm::ParameterSet>( "minHitsPerSubDet" ).getParameter<int>( "inTIDminus" ) ),
75  minHitsinTECplus_ (cfg.getParameter<edm::ParameterSet>( "minHitsPerSubDet" ).getParameter<int>( "inTECplus" ) ),
76  minHitsinTECminus_ (cfg.getParameter<edm::ParameterSet>( "minHitsPerSubDet" ).getParameter<int>( "inTECminus" ) ),
77  minHitsinFPIXplus_ (cfg.getParameter<edm::ParameterSet>( "minHitsPerSubDet" ).getParameter<int>( "inFPIXplus" ) ),
78  minHitsinFPIXminus_ (cfg.getParameter<edm::ParameterSet>( "minHitsPerSubDet" ).getParameter<int>( "inFPIXminus" ) ),
79  minHitsinENDCAP_ (cfg.getParameter<edm::ParameterSet>( "minHitsPerSubDet" ).getParameter<int>( "inENDCAP" ) ),
80  minHitsinENDCAPplus_ (cfg.getParameter<edm::ParameterSet>( "minHitsPerSubDet" ).getParameter<int>( "inENDCAPplus" ) ),
81  minHitsinENDCAPminus_ (cfg.getParameter<edm::ParameterSet>( "minHitsPerSubDet" ).getParameter<int>( "inENDCAPminus" ) ),
82  maxHitDiffEndcaps_( cfg.getParameter<double>( "maxHitDiffEndcaps" ) ),
83  nLostHitMax_( cfg.getParameter<double>( "nLostHitMax" ) ),
84  clusterValueMapTag_(cfg.getParameter<edm::InputTag>("hitPrescaleMapTag")),
85  minPrescaledHits_( cfg.getParameter<int>("minPrescaledHits")),
86  applyPrescaledHitsFilter_(clusterValueMapTag_.encode().size() && minPrescaledHits_ > 0)
87 {
88 
89  //convert track quality from string to enum
90  std::vector<std::string> trkQualityStrings(cfg.getParameter<std::vector<std::string> >("trackQualities"));
91  std::string qualities;
92  if(trkQualityStrings.size()>0){
93  applyTrkQualityCheck_=true;
94  for (unsigned int i = 0; i < trkQualityStrings.size(); ++i) {
95  (qualities += trkQualityStrings[i]) += ", ";
96  trkQualities_.push_back(reco::TrackBase::qualityByName(trkQualityStrings[i]));
97  }
98  }
99  else applyTrkQualityCheck_=false;
100 
101  std::vector<std::string> trkIterStrings(cfg.getParameter<std::vector<std::string> >("iterativeTrackingSteps"));
102  if(trkIterStrings.size()>0){
103  applyIterStepCheck_=true;
104  std::string tracksteps;
105  for (unsigned int i = 0; i < trkIterStrings.size(); ++i) {
106  (tracksteps += trkIterStrings[i]) += ", ";
107  trkSteps_.push_back(reco::TrackBase::algoByName(trkIterStrings[i]));
108  }
109  }
110  else applyIterStepCheck_=false;
111 
112  if (applyBasicCuts_){
113  edm::LogInfo("AlignmentTrackSelector")
114  << "applying basic track cuts ..."
115  << "\nptmin,ptmax: " << ptMin_ << "," << ptMax_
116  << "\npmin,pmax: " << pMin_ << "," << pMax_
117  << "\netamin,etamax: " << etaMin_ << "," << etaMax_
118  << "\nphimin,phimax: " << phiMin_ << "," << phiMax_
119  << "\nnhitmin,nhitmax: " << nHitMin_ << "," << nHitMax_
120  << "\nnlosthitmax: " << nLostHitMax_
121  << "\nnhitmin2D: " << nHitMin2D_
122  << (countStereoHitAs2D_ ? "," : ", not") << " counting hits on SiStrip stereo modules as 2D"
123  << "\nchi2nmax: " << chi2nMax_;
124 
125  if (applyIsolation_)
126  edm::LogInfo("AlignmentTrackSelector")
127  << "only retain tracks isolated at least by " << minHitIsolation_ <<
128  " cm from other rechits";
129 
130  if (chargeCheck_)
131  edm::LogInfo("AlignmentTrackSelector")
132  << "only retain hits with at least " << minHitChargeStrip_ <<
133  " ADC counts of total cluster charge";
134 
135  edm::LogInfo("AlignmentTrackSelector")
136  << "Minimum number of hits in TIB/TID/TOB/TEC/BPIX/FPIX/PIXEL = "
137  << minHitsinTIB_ << "/" << minHitsinTID_ << "/" << minHitsinTOB_
138  << "/" << minHitsinTEC_ << "/" << minHitsinBPIX_ << "/" << minHitsinFPIX_ << "/" << minHitsinPIX_;
139 
140  edm::LogInfo("AlignmentTrackSelector")
141  << "Minimum number of hits in TID+/TID-/TEC+/TEC-/FPIX+/FPIX- = "
142  << minHitsinTIDplus_ << "/" << minHitsinTIDminus_
143  << "/" << minHitsinTECplus_ << "/" << minHitsinTECminus_
144  << "/" << minHitsinFPIXplus_ << "/" << minHitsinFPIXminus_;
145 
146  edm::LogInfo("AlignmentTrackSelector")
147  << "Minimum number of hits in EndCap (TID+TEC)/EndCap+/EndCap- = "
148  << minHitsinENDCAP_ << "/" << minHitsinENDCAPplus_ << "/" << minHitsinENDCAPminus_;
149 
150  edm::LogInfo("AlignmentTrackSelector")
151  << "Max value of |nHitsinENDCAPplus - nHitsinENDCAPminus| = "
152  << maxHitDiffEndcaps_;
153 
154  if (trkQualityStrings.size()) {
155  edm::LogInfo("AlignmentTrackSelector")
156  << "Select tracks with these qualities: " << qualities;
157  }
158  }
159 
160  if (applyNHighestPt_)
161  edm::LogInfo("AlignmentTrackSelector")
162  << "filter N tracks with highest Pt N=" << nHighestPt_;
163 
164  if (applyMultiplicityFilter_)
165  edm::LogInfo("AlignmentTrackSelector")
166  << "apply multiplicity filter N>= " << minMultiplicity_ << "and N<= " << maxMultiplicity_
167  << " on " << (multiplicityOnInput_ ? "input" : "output");
168 
169  if (applyPrescaledHitsFilter_) {
170  edm::LogInfo("AlignmentTrackSelector")
171  << "apply cut on number of prescaled hits N>= " << minPrescaledHits_
172  << " (prescale info from " << clusterValueMapTag_ << ")";
173 
174  }
175 
176 }
177 
178 // destructor -----------------------------------------------------------------
179 
180 AlignmentTrackSelector::~AlignmentTrackSelector()
181 {}
182 
183 
184 // do selection ---------------------------------------------------------------
185 
186 AlignmentTrackSelector::Tracks
187 AlignmentTrackSelector::select(const Tracks& tracks, const edm::Event& evt) const
188 {
189 
190  if (applyMultiplicityFilter_ && multiplicityOnInput_ &&
191  (tracks.size() < static_cast<unsigned int>(minMultiplicity_)
192  || tracks.size() > static_cast<unsigned int>(maxMultiplicity_))) {
193  return Tracks(); // empty collection
194  }
195 
196  Tracks result = tracks;
197  // apply basic track cuts (if selected)
198  if (applyBasicCuts_) result= this->basicCuts(result, evt);
199 
200  // filter N tracks with highest Pt (if selected)
201  if (applyNHighestPt_) result = this->theNHighestPtTracks(result);
202 
203  // apply minimum multiplicity requirement (if selected)
204  if (applyMultiplicityFilter_ && !multiplicityOnInput_) {
205  if (result.size() < static_cast<unsigned int>(minMultiplicity_)
206  || result.size() > static_cast<unsigned int>(maxMultiplicity_) ) {
207 
208  result.clear();
209  }
210  }
211 
212  if(applyPrescaledHitsFilter_){
213  result = this->checkPrescaledHits(result, evt);
214  }
215 
216  return result;
217 }
218 
220 bool AlignmentTrackSelector::useThisFilter()
221 {
222  return applyMultiplicityFilter_ || applyBasicCuts_ || applyNHighestPt_ || applyPrescaledHitsFilter_;
223 
224 }
225 
226 
227 // make basic cuts ------------------------------------------------------------
228 
229 AlignmentTrackSelector::Tracks
230 AlignmentTrackSelector::basicCuts(const Tracks& tracks, const edm::Event& evt) const
231 {
232  Tracks result;
233 
234  for (Tracks::const_iterator it=tracks.begin(); it != tracks.end(); ++it) {
235  const reco::Track* trackp=*it;
236  float pt=trackp->pt();
237  float p=trackp->p();
238  float eta=trackp->eta();
239  float phi=trackp->phi();
240  int nhit = trackp->numberOfValidHits();
241  int nlosthit = trackp->numberOfLostHits();
242  float chi2n = trackp->normalizedChi2();
243 
244  int q = trackp->charge();
245  bool isChargeOk = false;
246  if(theCharge_==-1 && q<0) isChargeOk = true;
247  else if (theCharge_==1 && q>0) isChargeOk = true;
248  else if (theCharge_==0) isChargeOk = true;
249 
250  float d0 = trackp->d0();
251  float dz = trackp->dz();
252 
253  // edm::LogDebug("AlignmentTrackSelector") << " pt,eta,phi,nhit: "
254  // <<pt<<","<<eta<<","<<phi<<","<<nhit;
255 
256  if (pt>ptMin_ && pt<ptMax_
257  && p>pMin_ && p<pMax_
258  && eta>etaMin_ && eta<etaMax_
259  && phi>phiMin_ && phi<phiMax_
260  && nhit>=nHitMin_ && nhit<=nHitMax_
261  && nlosthit<=nLostHitMax_
262  && chi2n<chi2nMax_
263  && isChargeOk
264  && d0>=d0Min_ && d0<=d0Max_
265  && dz>=dzMin_ && dz<=dzMax_) {
266  bool trkQualityOk=false ;
267  if (!applyTrkQualityCheck_ &&!applyIterStepCheck_)trkQualityOk=true ; // nothing required
268  else trkQualityOk = this->isOkTrkQuality(trackp);
269 
270  bool hitsCheckOk=this->detailedHitsCheck(trackp, evt);
271 
272  if (trkQualityOk && hitsCheckOk ) result.push_back(trackp);
273  }
274  }
275 
276  return result;
277 }
278 
279 //-----------------------------------------------------------------------------
280 
281 bool AlignmentTrackSelector::detailedHitsCheck(const reco::Track *trackp, const edm::Event& evt) const
282 {
283  // checking hit requirements beyond simple number of valid hits
284 
285  if (minHitsinTIB_ || minHitsinTOB_ || minHitsinTID_ || minHitsinTEC_
286  || minHitsinENDCAP_ || minHitsinENDCAPplus_ || minHitsinENDCAPminus_ || (maxHitDiffEndcaps_ < 999)
287  || minHitsinTIDplus_ || minHitsinTIDminus_
288  || minHitsinFPIXplus_ || minHitsinFPIXminus_
289  || minHitsinTECplus_ || minHitsinTECminus_
290  || minHitsinFPIX_ || minHitsinBPIX_ || minHitsinPIX_ ||nHitMin2D_ || chargeCheck_
291  || applyIsolation_ || (seedOnlyFromAbove_ == 1 || seedOnlyFromAbove_ == 2)) {
292  // any detailed hit cut is active, so have to check
293 
294  int nhitinTIB = 0, nhitinTOB = 0, nhitinTID = 0;
295  int nhitinTEC = 0, nhitinBPIX = 0, nhitinFPIX = 0, nhitinPIXEL=0;
296  int nhitinENDCAP = 0, nhitinENDCAPplus = 0, nhitinENDCAPminus = 0;
297  int nhitinTIDplus = 0, nhitinTIDminus = 0;
298  int nhitinFPIXplus = 0, nhitinFPIXminus = 0;
299  int nhitinTECplus = 0, nhitinTECminus = 0;
300  unsigned int nHit2D = 0;
301  unsigned int thishit = 0;
302 
303  for (trackingRecHit_iterator iHit = trackp->recHitsBegin(); iHit != trackp->recHitsEnd(); ++iHit) {
304  thishit++;
305  const DetId detId((*iHit)->geographicalId());
306  const int subdetId = detId.subdetId();
307 
308  // *** thishit == 1 means last hit in CTF ***
309  // (FIXME: assumption might change or not be valid for all tracking algorthms)
310  // ==> for cosmics
311  // seedOnlyFrom = 1 is TIB-TOB-TEC tracks only
312  // seedOnlyFrom = 2 is TOB-TEC tracks only
313  if (seedOnlyFromAbove_ == 1 && thishit == 1
314  && (subdetId == int(SiStripDetId::TOB) || subdetId == int(SiStripDetId::TEC))){
315  return false;
316  }
317  if (seedOnlyFromAbove_ == 2 && thishit == 1 && subdetId == int(SiStripDetId::TIB)) {
318  return false;
319  }
320 
321  if (!(*iHit)->isValid()) continue; // only real hits count as in trackp->numberOfValidHits()
322  if (detId.det() != DetId::Tracker) {
323  edm::LogError("DetectorMismatch") << "@SUB=AlignmentTrackSelector::detailedHitsCheck"
324  << "DetId.det() != DetId::Tracker (=" << DetId::Tracker
325  << "), but " << detId.det() << ".";
326  }
327  const TrackingRecHit* therechit = (*iHit).get();
328  if (chargeCheck_ && !(this->isOkCharge(therechit))) return false;
329  if (applyIsolation_ && (!this->isIsolated(therechit, evt))) return false;
330  if (SiStripDetId::TIB == subdetId) ++nhitinTIB;
331  else if (SiStripDetId::TOB == subdetId) ++nhitinTOB;
332  else if (SiStripDetId::TID == subdetId) {
333  ++nhitinTID;
334  ++nhitinENDCAP;
335  TIDDetId tidId(detId);
336  if (tidId.isZMinusSide()) {
337  ++nhitinTIDminus;
338  ++nhitinENDCAPminus;
339  }
340  else if (tidId.isZPlusSide()) {
341  ++nhitinTIDplus;
342  ++nhitinENDCAPplus;
343  }
344  }
345  else if (SiStripDetId::TEC == subdetId) {
346  ++nhitinTEC;
347  ++nhitinENDCAP;
348  TECDetId tecId(detId);
349  if (tecId.isZMinusSide()) {
350  ++nhitinTECminus;
351  ++nhitinENDCAPminus;
352  }
353  else if (tecId.isZPlusSide()) {
354  ++nhitinTECplus;
355  ++nhitinENDCAPplus;
356  }
357  }
358  else if ( kBPIX == subdetId) {++nhitinBPIX;++nhitinPIXEL;}
359  else if ( kFPIX == subdetId) {
360  ++nhitinFPIX;
361  ++nhitinPIXEL;
362  PXFDetId fpixId(detId);
363  if (fpixId.side()==1) ++nhitinFPIXminus;
364  else if (fpixId.side()==2) ++nhitinFPIXplus;
365  }
366  // Do not call isHit2D(..) if already enough 2D hits for performance reason:
367  if (nHit2D < nHitMin2D_ && this->isHit2D(**iHit)) ++nHit2D;
368  } // end loop on hits
369 
370  return (nhitinTIB >= minHitsinTIB_ && nhitinTOB >= minHitsinTOB_
371  && nhitinTID >= minHitsinTID_ && nhitinTEC >= minHitsinTEC_
372  && nhitinENDCAP >= minHitsinENDCAP_ && nhitinENDCAPplus >= minHitsinENDCAPplus_ && nhitinENDCAPminus >= minHitsinENDCAPminus_
373  && std::abs(nhitinENDCAPplus-nhitinENDCAPminus) <= maxHitDiffEndcaps_
374  && nhitinTIDplus >= minHitsinTIDplus_ && nhitinTIDminus >= minHitsinTIDminus_
375  && nhitinFPIXplus >= minHitsinFPIXplus_ && nhitinFPIXminus >= minHitsinFPIXminus_
376  && nhitinTECplus >= minHitsinTECplus_ && nhitinTECminus >= minHitsinTECminus_
377  && nhitinBPIX >= minHitsinBPIX_
378  && nhitinFPIX >= minHitsinFPIX_ && nhitinPIXEL>=minHitsinPIX_
379  && nHit2D >= nHitMin2D_);
380  } else { // no cuts set, so we are just fine and can avoid loop on hits
381  return true;
382  }
383 
384 }
385 
386 //-----------------------------------------------------------------------------
387 
388 bool AlignmentTrackSelector::isHit2D(const TrackingRecHit &hit) const
389 {
390  // we count SiStrip stereo modules as 2D if selected via countStereoHitAs2D_
391  // (since they provide theta information)
392  if (!hit.isValid() ||
393  (hit.dimension() < 2 && !countStereoHitAs2D_ && !dynamic_cast<const SiStripRecHit1D*>(&hit))){
394  return false; // real RecHit1D - but SiStripRecHit1D depends on countStereoHitAs2D_
395  } else {
396  const DetId detId(hit.geographicalId());
397  if (detId.det() == DetId::Tracker) {
398  if (detId.subdetId() == kBPIX || detId.subdetId() == kFPIX) {
399  return true; // pixel is always 2D
400  } else { // should be SiStrip now
401  const SiStripDetId stripId(detId);
402  if (stripId.stereo()) return countStereoHitAs2D_; // stereo modules
403  else if (dynamic_cast<const SiStripRecHit1D*>(&hit)
404  || dynamic_cast<const SiStripRecHit2D*>(&hit)) return false; // rphi modules hit
405  //the following two are not used any more since ages...
406  else if (dynamic_cast<const SiStripMatchedRecHit2D*>(&hit)) return true; // matched is 2D
407  else if (dynamic_cast<const ProjectedSiStripRecHit2D*>(&hit)) {
408  const ProjectedSiStripRecHit2D* pH = static_cast<const ProjectedSiStripRecHit2D*>(&hit);
409  return (countStereoHitAs2D_ && this->isHit2D(pH->originalHit())); // depends on original...
410  } else {
411  edm::LogError("UnkownType") << "@SUB=AlignmentTrackSelector::isHit2D"
412  << "Tracker hit not in pixel, neither SiStripRecHit[12]D nor "
413  << "SiStripMatchedRecHit2D nor ProjectedSiStripRecHit2D.";
414  return false;
415  }
416  }
417  } else { // not tracker??
418  edm::LogWarning("DetectorMismatch") << "@SUB=AlignmentTrackSelector::isHit2D"
419  << "Hit not in tracker with 'official' dimension >=2.";
420  return true; // dimension() >= 2 so accept that...
421  }
422  }
423  // never reached...
424 }
425 
426 //-----------------------------------------------------------------------------
427 
428 bool AlignmentTrackSelector::isOkCharge(const TrackingRecHit* hit) const
429 {
430  if (!hit || !hit->isValid()) return true;
431 
432  // check det and subdet
433  const DetId id(hit->geographicalId());
434  if (id.det() != DetId::Tracker) {
435  edm::LogWarning("DetectorMismatch") << "@SUB=isOkCharge" << "Hit not in tracker!";
436  return true;
437  }
438  if (id.subdetId() == kFPIX || id.subdetId() == kBPIX) {
439  return true; // might add some requirement...
440  }
441 
442  // We are in SiStrip now, so test normal hit:
443  const std::type_info &type = typeid(*hit);
444 
445 
446  if (type == typeid(SiStripRecHit2D)) {
447  const SiStripRecHit2D *stripHit2D = dynamic_cast<const SiStripRecHit2D*>(hit);
448  if (stripHit2D) {
449  return this->isOkChargeStripHit(*stripHit2D);
450  }
451  }
452  else if(type == typeid(SiStripRecHit1D)){
453  const SiStripRecHit1D *stripHit1D = dynamic_cast<const SiStripRecHit1D*>(hit);
454  if (stripHit1D) {
455  return this->isOkChargeStripHit(*stripHit1D);
456  }
457  }
458 
459  else if(type == typeid(SiStripMatchedRecHit2D)){ // or matched (should not occur anymore due to hit splitting since 20X)
460  const SiStripMatchedRecHit2D *matchedHit = dynamic_cast<const SiStripMatchedRecHit2D*>(hit);
461  if (matchedHit) {
462  return (this->isOkChargeStripHit(matchedHit->monoHit())
463  && this->isOkChargeStripHit(matchedHit->stereoHit()));
464  }
465  }
466  else if(type == typeid(ProjectedSiStripRecHit2D)){
467  // or projected (should not occur anymore due to hit splitting since 20X):
468  const ProjectedSiStripRecHit2D *projHit = dynamic_cast<const ProjectedSiStripRecHit2D*>(hit);
469  if (projHit) {
470  return this->isOkChargeStripHit(projHit->originalHit());
471  }
472  }
473  else{
474  edm::LogError("AlignmentTrackSelector")<< "@SUB=isOkCharge" <<"Unknown type of a valid tracker hit in Strips "
475  << " SubDet = "<<id.subdetId();
476  return false;
477  }
478 
479 
480 
481  // and now? SiTrackerMultiRecHit? Not here I guess!
482  // Should we throw instead?
483  edm::LogError("AlignmentTrackSelector")
484  << "@SUB=isOkCharge" << "Unknown valid tracker hit not in pixel, subdet " << id.subdetId()
485  << ", SiTrackerMultiRecHit " << dynamic_cast<const SiTrackerMultiRecHit*>(hit)
486  << ", BaseTrackerRecHit " << dynamic_cast<const BaseTrackerRecHit*>(hit);
487 
488  return true;
489 }
490 
491 
492 //-----------------------------------------------------------------------------
493 
494 bool AlignmentTrackSelector::isOkChargeStripHit(const SiStripRecHit2D & siStripRecHit2D) const
495 {
496  double charge = 0.;
497 
498  SiStripRecHit2D::ClusterRef cluster(siStripRecHit2D.cluster());
499  const std::vector<uint8_t> &amplitudes = cluster->amplitudes();
500 
501  for (size_t ia = 0; ia < amplitudes.size(); ++ia) {
502  charge += amplitudes[ia];
503  }
504 
505  return (charge >= minHitChargeStrip_);
506 }
507 
508 //-----------------------------------------------------------------------------
509 
510 bool AlignmentTrackSelector::isOkChargeStripHit(const SiStripRecHit1D & siStripRecHit1D) const
511 {
512  double charge = 0.;
513 
514  SiStripRecHit1D::ClusterRef cluster(siStripRecHit1D.cluster());
515  const std::vector<uint8_t> &amplitudes = cluster->amplitudes();
516 
517  for (size_t ia = 0; ia < amplitudes.size(); ++ia) {
518  charge += amplitudes[ia];
519  }
520 
521  return (charge >= minHitChargeStrip_);
522 }
523 
524 //-----------------------------------------------------------------------------
525 
526 bool AlignmentTrackSelector::isIsolated(const TrackingRecHit* therechit, const edm::Event& evt) const
527 {
528  // FIXME:
529  // adapt to changes after introduction of SiStripRecHit1D...
530  //
531  // edm::ESHandle<TrackerGeometry> tracker;
532  edm::Handle<SiStripRecHit2DCollection> rphirecHits;
533  edm::Handle<SiStripMatchedRecHit2DCollection> matchedrecHits;
534  // es.get<TrackerDigiGeometryRecord>().get(tracker);
535  evt.getByLabel( rphirecHitsTag_, rphirecHits );
536  evt.getByLabel( matchedrecHitsTag_, matchedrecHits );
537 
538  SiStripRecHit2DCollection::DataContainer::const_iterator istripSt;
539  SiStripMatchedRecHit2DCollection::DataContainer::const_iterator istripStm;
540  const SiStripRecHit2DCollection::DataContainer& stripcollSt = rphirecHits->data();
541  const SiStripMatchedRecHit2DCollection::DataContainer& stripcollStm = matchedrecHits->data();
542 
543  DetId idet = therechit->geographicalId();
544 
545  // FIXME: instead of looping the full hit collection, we should explore the features of
546  // SiStripRecHit2DCollection::rangeRphi = rphirecHits.get(idet) and loop
547  // only from rangeRphi.first until rangeRphi.second
548  for( istripSt=stripcollSt.begin(); istripSt!=stripcollSt.end(); ++istripSt ) {
549  const SiStripRecHit2D *aHit = &*(istripSt);
550  DetId mydet1 = aHit->geographicalId();
551  if (idet.rawId() != mydet1.rawId()) continue;
552  float theDistance = ( therechit->localPosition() - aHit->localPosition() ).mag();
553  // std::cout << "theDistance1 = " << theDistance << "\n";
554  if (theDistance > 0.001 && theDistance < minHitIsolation_) return false;
555  }
556 
557  // FIXME: see above
558  for( istripStm=stripcollStm.begin(); istripStm!=stripcollStm.end(); ++istripStm ) {
559  const SiStripMatchedRecHit2D *aHit = &*(istripStm);
560  DetId mydet2 = aHit->geographicalId();
561  if (idet.rawId() != mydet2.rawId()) continue;
562  float theDistance = (therechit->localPosition() - aHit->localPosition()).mag();
563  // std::cout << "theDistance1 = " << theDistance << "\n";
564  if (theDistance > 0.001 && theDistance < minHitIsolation_) return false;
565  }
566  return true;
567 }
568 
569 //-----------------------------------------------------------------------------
570 
571 AlignmentTrackSelector::Tracks
572 AlignmentTrackSelector::theNHighestPtTracks(const Tracks& tracks) const
573 {
574  Tracks sortedTracks=tracks;
575  Tracks result;
576 
577  // sort in pt
578  std::sort(sortedTracks.begin(),sortedTracks.end(),ptComparator);
579 
580  // copy theTrackMult highest pt tracks to result vector
581  int n=0;
582  for (Tracks::const_iterator it=sortedTracks.begin();
583  it!=sortedTracks.end(); ++it) {
584  if (n<nHighestPt_) { result.push_back(*it); n++; }
585  }
586 
587  return result;
588 }
589 
590 //---------
591 
592 AlignmentTrackSelector::Tracks
593 AlignmentTrackSelector::checkPrescaledHits(const Tracks& tracks, const edm::Event& evt) const
594 {
595  Tracks result;
596 
597  //take Cluster-Flag Assomap
598  edm::Handle<AliClusterValueMap> fMap;
599  evt.getByLabel( clusterValueMapTag_, fMap);
600  const AliClusterValueMap &flagMap=*fMap;
601 
602  //for each track loop on hits and count the number of taken hits
603  for (Tracks::const_iterator ittrk=tracks.begin(); ittrk != tracks.end(); ++ittrk) {
604  const reco::Track* trackp=*ittrk;
605  int ntakenhits=0;
606  // float pt=trackp->pt();
607 
608  for (trackingRecHit_iterator ith = trackp->recHitsBegin(), edh = trackp->recHitsEnd(); ith != edh; ++ith) {
609  const TrackingRecHit *hit = ith->get(); // ith is an iterator on edm::Ref to rechit
610  if(! hit->isValid())continue;
611  DetId detid = hit->geographicalId();
612  int subDet = detid.subdetId();
613  AlignmentClusterFlag flag;
614 
615  bool isPixelHit=(subDet == kFPIX || subDet == kBPIX);
616 
617  if (!isPixelHit){
618  const std::type_info &type = typeid(*hit);
619 
620  if (type == typeid(SiStripRecHit2D)) {
621  const SiStripRecHit2D* striphit=dynamic_cast<const SiStripRecHit2D*>(hit);
622  if(striphit!=0){
623  SiStripRecHit2D::ClusterRef stripclust(striphit->cluster());
624  flag = flagMap[stripclust];
625  }
626  }
627  else if(type == typeid(SiStripRecHit1D)){
628  const SiStripRecHit1D* striphit = dynamic_cast<const SiStripRecHit1D*>(hit);
629  if(striphit!=0){
630  SiStripRecHit1D::ClusterRef stripclust(striphit->cluster());
631  flag = flagMap[stripclust];
632  }
633  }
634  else{
635  edm::LogError("AlignmentTrackSelector")<<"ERROR in <AlignmentTrackSelector::checkPrescaledHits>: Dynamic cast of Strip RecHit failed!"
636  <<" Skipping this hit.";
637  continue;
638  }
639 
640 
641 
642  }//end if hit in Strips
643  else{ // test explicitely BPIX/FPIX
644  const SiPixelRecHit* pixelhit= dynamic_cast<const SiPixelRecHit*>(hit);
645  if(pixelhit!=0){
646  SiPixelRecHit::ClusterRef pixclust(pixelhit->cluster());
647  flag = flagMap[pixclust];
648  }
649  else{
650  edm::LogError("AlignmentTrackSelector")<<"ERROR in <AlignmentTrackSelector::checkPrescaledHits>: Dynamic cast of Pixel RecHit failed! ";
651  }
652  }//end else hit is in Pixel
653 
654  if(flag.isTaken())ntakenhits++;
655 
656  }//end loop on hits
657  if(ntakenhits >= minPrescaledHits_)result.push_back(trackp);
658  }//end loop on tracks
659 
660  return result;
661 }//end checkPrescaledHits
662 
663 //-----------------------------------------------------------------
664 
665 bool AlignmentTrackSelector::isOkTrkQuality(const reco::Track* track) const
666 {
667  bool qualityOk=false;
668  bool iterStepOk=false;
669 
670  //check iterative step
671  if(applyIterStepCheck_){
672  for (unsigned int i = 0; i < trkSteps_.size(); ++i) {
673  if (track->algo()==(trkSteps_[i])) {
674  iterStepOk=true;
675  }
676  }
677  }
678  else iterStepOk=true;
679 
680  //check track quality
681  if(applyTrkQualityCheck_){
682  for (unsigned int i = 0; i < trkQualities_.size(); ++i) {
683  if (track->quality(trkQualities_[i])) {
684  qualityOk=true;
685  }
686  }
687  }
688  else qualityOk=true;
689 
690  return qualityOk&&iterStepOk;
691 }//end check on track quality
SiStripRecHit1D::cluster
ClusterRef cluster() const
Definition: SiStripRecHit1D.h:36
AlignmentTrackSelector::detailedHitsCheck
bool detailedHitsCheck(const reco::Track *track, const edm::Event &evt) const
checking hit requirements beyond simple number of valid hits
Definition: AlignmentTrackSelector.cc:281
reco::TrackBase::p
double p() const
momentum vector magnitude
Definition: TrackBase.h:129
type
type
Definition: HCALResponse.h:22
edm::ParameterSet::getParameter
T getParameter(std::string const &) const
i
int i
Definition: DBlmapReader.cc:9
TrackingRecHit::dimension
virtual int dimension() const =0
edmNew::DetSetVector::const_iterator
boost::transform_iterator< IterHelp, const_IdIter > const_iterator
Definition: DetSetVectorNew.h:106
flag
long int flag
Definition: mlp_lapack.h:47
AlignmentTrackSelector::isIsolated
bool isIsolated(const TrackingRecHit *therechit, const edm::Event &evt) const
Definition: AlignmentTrackSelector.cc:526
reco::TrackBase::d0
double d0() const
dxy parameter in perigee convention (d0 = - dxy)
Definition: TrackBase.h:123
SiStripDetId::stereo
uint32_t stereo() const
Definition: SiStripDetId.h:162
edmNew::DetSetVector::DataContainer
std::vector< data_type > DataContainer
Definition: DetSetVectorNew.h:82
reco::TrackBase::normalizedChi2
double normalizedChi2() const
chi-squared divided by n.d.o.f. (or chi-squared * 1e6 if n.d.o.f. is zero)
Definition: TrackBase.h:111
mag
T mag() const
The vector magnitude. Equivalent to sqrt(vec.mag2())
Definition: Basic3DVectorLD.h:142
AlignmentTrackSelector::basicCuts
Tracks basicCuts(const Tracks &tracks, const edm::Event &evt) const
apply basic cuts on pt,eta,phi,nhit
Definition: AlignmentTrackSelector.cc:230
AlignmentTrackSelector::Tracks
std::vector< const reco::Track * > Tracks
Definition: AlignmentTrackSelector.h:22
AlignmentTrackSelector::checkPrescaledHits
Tracks checkPrescaledHits(const Tracks &tracks, const edm::Event &evt) const
Definition: AlignmentTrackSelector.cc:593
abs
#define abs(x)
Definition: mlp_lapack.h:159
AlignmentTrackSelector::trkSteps_
std::vector< reco::TrackBase::TrackAlgorithm > trkSteps_
Definition: AlignmentTrackSelector.h:90
AlignmentTrackSelector::theNHighestPtTracks
Tracks theNHighestPtTracks(const Tracks &tracks) const
filter the n highest pt tracks
Definition: AlignmentTrackSelector.cc:572
reco::TrackBase::phi
double phi() const
azimuthal angle of momentum vector
Definition: TrackBase.h:139
reco::TrackBase::numberOfLostHits
unsigned short numberOfLostHits() const
number of cases where track crossed a layer without getting a hit.
Definition: TrackBase.h:234
AlignmentTrackSelector::trkQualities_
std::vector< reco::TrackBase::TrackQuality > trkQualities_
Definition: AlignmentTrackSelector.h:88
AlignmentTrackSelector::nHitMin2D_
const unsigned int nHitMin2D_
Definition: AlignmentTrackSelector.h:75
eta
T eta() const
Definition: Basic3DVectorLD.h:172
AlignmentTrackSelector::isOkChargeStripHit
bool isOkChargeStripHit(const SiStripRecHit1D &siStripRecHit1D) const
Definition: AlignmentTrackSelector.cc:510
DeDxDiscriminatorTools::charge
double charge(const std::vector< uint8_t > &Ampls)
Definition: DeDxDiscriminatorTools.cc:43
AlignmentTrackSelector::AlignmentTrackSelector
AlignmentTrackSelector(const edm::ParameterSet &cfg)
constructor
Definition: AlignmentTrackSelector.cc:32
reco::TrackBase::algo
TrackAlgorithm algo() const
Definition: TrackBase.h:332
DetId::rawId
uint32_t rawId() const
get the raw id
Definition: DetId.h:45
kFPIX
const int kFPIX
Definition: AlignmentTrackSelector.cc:28
AlignmentTrackSelector::isHit2D
bool isHit2D(const TrackingRecHit &hit) const
Definition: AlignmentTrackSelector.cc:388
errorMatrix2Lands_multiChannel.id
tuple id
Definition: errorMatrix2Lands_multiChannel.py:151
TECDetId::isZMinusSide
bool isZMinusSide() const
Definition: TECDetId.h:86
reco::TrackBase::eta
double eta() const
pseudorapidity of momentum vector
Definition: TrackBase.h:141
AlignmentTrackSelector::select
Tracks select(const Tracks &tracks, const edm::Event &evt) const
select tracks
Definition: AlignmentTrackSelector.cc:187
reco::TrackBase::pt
double pt() const
track transverse momentum
Definition: TrackBase.h:131
query.result
tuple result
Definition: query.py:137
SiStripRecHit2D::cluster
ClusterRef cluster() const
Definition: SiStripRecHit2D.h:43
reco::TrackBase::numberOfValidHits
unsigned short numberOfValidHits() const
number of valid hits found
Definition: TrackBase.h:232
reco::Track::recHitsBegin
trackingRecHit_iterator recHitsBegin() const
Iterator to first hit on the track.
Definition: Track.h:63
AlignmentTrackSelector::matchedrecHitsTag_
const edm::InputTag matchedrecHitsTag_
Definition: AlignmentTrackSelector.h:73
TECDetId::isZPlusSide
bool isZPlusSide() const
Definition: TECDetId.h:82
DetId::subdetId
int subdetId() const
get the contents of the subdetector field (not cast into any detector&#39;s numbering enum) ...
Definition: DetId.h:39
AlignmentTrackSelector::isOkTrkQuality
bool isOkTrkQuality(const reco::Track *track) const
Definition: AlignmentTrackSelector.cc:665
edm::Event::getByLabel
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
Definition: Event.h:356
reco::TrackBase::dz
double dz() const
dz parameter (= dsz/cos(lambda)). This is the track z0 w.r.t (0,0,0) only if the refPoint is close to...
Definition: TrackBase.h:127
SiStripDetId
Detector identifier class for the strip tracker.
Definition: SiStripDetId.h:17
DetId
Definition: DetId.h:20
AlignmentTrackSelector::ptMin_
const double ptMin_
if true, cut min/maxMultiplicity on input instead of on final result
Definition: AlignmentTrackSelector.h:68
TIDDetId::isZMinusSide
bool isZMinusSide() const
Definition: TIDDetId.h:81
reco::TrackBase::qualityByName
static TrackQuality qualityByName(const std::string &name)
Definition: TrackBase.cc:46
BaseTrackerRecHit::localPosition
virtual LocalPoint localPosition() const
Definition: BaseTrackerRecHit.cc:28
SiStripMatchedRecHit2D::stereoHit
SiStripRecHit2D stereoHit() const
Definition: SiStripMatchedRecHit2D.h:18
testEve_cfg.tracks
tuple tracks
Definition: testEve_cfg.py:39
AlignmentTrackSelector::clusterValueMapTag_
const edm::InputTag clusterValueMapTag_
Definition: AlignmentTrackSelector.h:84
SiPixelRecHit::cluster
ClusterRef cluster() const
Definition: SiPixelRecHit.h:41
TrackingRecHit::isValid
bool isValid() const
Definition: TrackingRecHit.h:75
reco::TrackBase::quality
bool quality(const TrackQuality) const
Track quality.
Definition: TrackBase.h:377
SiStripMatchedRecHit2D::monoHit
SiStripRecHit2D monoHit() const
Definition: SiStripMatchedRecHit2D.h:19
AlCaHLTBitMon_ParallelJobs.p
tuple p
Definition: AlCaHLTBitMon_ParallelJobs.py:152
AlignmentTrackSelector::rphirecHitsTag_
const edm::InputTag rphirecHitsTag_
Definition: AlignmentTrackSelector.h:72
AlignmentTrackSelector::isOkCharge
bool isOkCharge(const TrackingRecHit *therechit) const
if valid, check for minimum charge (currently only in strip), if invalid give true ...
Definition: AlignmentTrackSelector.cc:428
TIDDetId::isZPlusSide
bool isZPlusSide() const
Definition: TIDDetId.h:77
reco::TrackBase::algoByName
static TrackAlgorithm algoByName(const std::string &name)
Definition: TrackBase.cc:55
PXFDetId::side
unsigned int side() const
positive or negative id
Definition: PXFDetId.h:38
reco::TrackBase::charge
int charge() const
track electric charge
Definition: TrackBase.h:113
TrackingRecHit::geographicalId
DetId geographicalId() const
Definition: TrackingRecHit.h:62
TrackingRecHit::localPosition
virtual LocalPoint localPosition() const =0
ProjectedSiStripRecHit2D::originalHit
const SiStripRecHit2D & originalHit() const
Definition: ProjectedSiStripRecHit2D.h:25
findQualityFiles.size
tuple size
Write out results.
Definition: findQualityFiles.py:442
AlignmentTrackSelector::useThisFilter
bool useThisFilter()
returns if any of the Filters is used.
Definition: AlignmentTrackSelector.cc:220
SiPixelRecHit
Our base class.
Definition: SiPixelRecHit.h:22
kBPIX
const int kBPIX
Definition: AlignmentTrackSelector.cc:27
reco::Track::recHitsEnd
trackingRecHit_iterator recHitsEnd() const
Iterator to last hit on the track.
Definition: Track.h:65
phi
Definition: DDAxes.h:10
Generated for CMSSW Reference Manual by doxygen 1.8.5