CMS 3D CMS Logo

HitPixelLayersTPSelector.h
Go to the documentation of this file.
1 #ifndef HitPixelLayersTrackSelection_h
2 #define HitPixelLayersTrackSelection_h
3 
6 
13 
22 {
23 
24  public:
25  // input collection type
27 
28 
29  // output collection type
31 
32  // iterator over result collection type.
34 
35  // constructor from parameter set configurability
37  tripletSeedOnly_(iConfig.getParameter<bool>("tripletSeedOnly")),
38  ptMin_(iConfig.getParameter<double>("ptMin")),
39  minRapidity_(iConfig.getParameter<double>("minRapidity")),
40  maxRapidity_(iConfig.getParameter<double>("maxRapidity")),
41  tip_(iConfig.getParameter<double>("tip")),
42  lip_(iConfig.getParameter<double>("lip")),
43  minHit_(iConfig.getParameter<int>("minHit")),
44  signalOnly_(iConfig.getParameter<bool>("signalOnly")),
45  chargedOnly_(iConfig.getParameter<bool>("chargedOnly")),
46  primaryOnly_(iConfig.getParameter<bool>("primaryOnly")),
47  tpStatusBased_(iConfig.getParameter<bool>("tpStatusBased")),
48  pdgId_(iConfig.getParameter< std::vector<int> >("pdgId"))
49  {};
50 
51  // select object from a collection and
52  // possibly event content
53  void select( const edm::Handle<collection> & TPCH, const edm::Event & iEvent, const edm::EventSetup & iSetup)
54  {
55  selected_.clear();
56  //Retrieve tracker topology from geometry
58  iSetup.get<TrackerTopologyRcd>().get(tTopoHand);
59  const TrackerTopology *tTopo=tTopoHand.product();
60 
61 
62  const collection & tpc = *(TPCH.product());
63 
64  for (TrackingParticleCollection::size_type i=0; i<tpc.size(); i++)
65  {
66  TrackingParticleRef tpr(TPCH, i);
67 
68  // quickly reject if it is from pile-up
69  if (signalOnly_ && !(tpr->eventId().bunchCrossing()==0 && tpr->eventId().event()==0) ) continue;
70  if (chargedOnly_ && tpr->charge()==0) continue; //select only if charge!=0
71  if (tpStatusBased_ && primaryOnly_ && tpr->status()!=1 ) continue; // TP status based sel primary
72  if ((!tpStatusBased_) && primaryOnly_ && tpr->parentVertex()->nSourceTracks()!=0 ) continue; // vertex based sel for primary
73 
74  // loop over specified PID values
75  bool testId = false;
76  unsigned int idSize = pdgId_.size();
77  if (idSize==0) testId = true;
78  else for (unsigned int it=0;it!=idSize;++it){
79  if (tpr->pdgId()==pdgId_[it]) testId = true;
80  }
81 
82  // selection criteria
83  if ( tpr->numberOfTrackerLayers() >= minHit_ &&
84  sqrt(tpr->momentum().perp2()) >= ptMin_ &&
85  tpr->momentum().eta() >= minRapidity_ && tpr->momentum().eta() <= maxRapidity_ &&
86  sqrt(tpr->vertex().perp2()) <= tip_ &&
87  fabs(tpr->vertex().z()) <= lip_ &&
88  testId)
89  {
90  if (tripletSeedOnly_ && !goodHitPattern(pixelHitPattern(tpr,tTopo)) ) continue; //findable triplet seed
91  selected_.push_back(tpr);
92  }
93 
94  }
95  }
96 
97  // return pixel layer hit pattern
98  std::vector<bool> pixelHitPattern( const TrackingParticleRef& simTrack, const TrackerTopology *tTopo )
99  {
100  std::vector<bool> hitpattern(5,false); // PXB 0,1,2 PXF 0,1
101  // This currently will always return false, since we can no loger use the sim hits to check for triplets. This would need to be fixed if we want to enable this feature, but it's not being used at the moment, since tripletSeedOnly is always set to False - Matt Nguyen, 24/7/2013
102 
103  return hitpattern;
104  }
105 
106  // test whether hit pattern would give a pixel triplet seed
107  bool goodHitPattern( const std::vector<bool>& hitpattern )
108  {
109  if( (hitpattern[0] && hitpattern[1] && hitpattern[2]) ||
110  (hitpattern[0] && hitpattern[1] && hitpattern[3]) ||
111  (hitpattern[0] && hitpattern[3] && hitpattern[4]) )
112  return true;
113  else
114  return false;
115  }
116 
117  // iterators over selected objects: collection begin
118  const_iterator begin() const
119  {
120  return selected_.begin();
121  }
122 
123  // iterators over selected objects: collection end
124  const_iterator end() const
125  {
126  return selected_.end();
127  }
128 
129  // true if no object has been selected
130  size_t size() const
131  {
132  return selected_.size();
133  }
134 
135  //private:
136 
137  container selected_;
139  double ptMin_;
140  double minRapidity_;
141  double maxRapidity_;
142  double tip_;
143  double lip_;
144  int minHit_;
149  std::vector<int> pdgId_;
150 
151 
152 };
153 
154 
155 #endif
std::vector< TrackingParticle > TrackingParticleCollection
std::vector< bool > pixelHitPattern(const TrackingParticleRef &simTrack, const TrackerTopology *tTopo)
const_iterator end() const
Termination of iteration.
Definition: RefVector.h:253
uint16_t size_type
const_iterator begin() const
Initialize an iterator over the RefVector.
Definition: RefVector.h:248
int iEvent
Definition: GenABIO.cc:230
TrackingParticleCollection collection
void select(const edm::Handle< collection > &TPCH, const edm::Event &iEvent, const edm::EventSetup &iSetup)
TrackingParticleRefVector container
simTrack
per collection params
T sqrt(T t)
Definition: SSEVec.h:18
container::const_iterator const_iterator
T const * product() const
Definition: Handle.h:81
void clear()
Clear the vector.
Definition: RefVector.h:147
bool goodHitPattern(const std::vector< bool > &hitpattern)
HitPixelLayersTPSelector(const edm::ParameterSet &iConfig, edm::ConsumesCollector &&iC)
T get() const
Definition: EventSetup.h:63
void push_back(value_type const &ref)
Add a Ref<C, T> to the RefVector.
Definition: RefVector.h:69
size_type size() const
Size of the RefVector.
Definition: RefVector.h:107
T const * product() const
Definition: ESHandle.h:86
const_iterator begin() const