CMS 3D CMS Logo

HLTPixelClusterShapeFilter.cc
Go to the documentation of this file.
2 
6 
8 
9 //
10 // class declaration
11 //
12 
14 public:
16  ~HLTPixelClusterShapeFilter() override;
17  static void fillDescriptions(edm::ConfigurationDescriptions & descriptions);
18 
19 private:
20 
22  edm::InputTag inputTag_; // input tag identifying product containing pixel clusters
23  double minZ_; // beginning z-vertex position
24  double maxZ_; // end z-vertex position
25  double zStep_; // size of steps in z-vertex test
26 
27  std::vector<double> clusterPars_; //pixel cluster polynomial pars for vertex compatibility cut
28  int nhitsTrunc_; //maximum pixel clusters to apply compatibility check
29  double clusterTrunc_; //maximum vertex compatibility value for event rejection
30 
31  struct VertexHit
32  {
33  float z;
34  float r;
35  float w;
36  };
37 
38  bool hltFilter(edm::Event&, const edm::EventSetup&, trigger::TriggerFilterObjectWithRefs & filterproduct) const override;
39  int getContainedHits(const std::vector<VertexHit> &hits, double z0, double &chi) const;
40 
41 };
42 
49 
59 
60 //
61 // constructors and destructor
62 //
63 
65  inputTag_ (config.getParameter<edm::InputTag>("inputTag")),
66  minZ_ (config.getParameter<double>("minZ")),
67  maxZ_ (config.getParameter<double>("maxZ")),
68  zStep_ (config.getParameter<double>("zStep")),
69  clusterPars_ (config.getParameter< std::vector<double> >("clusterPars")),
70  nhitsTrunc_ (config.getParameter<int>("nhitsTrunc")),
71  clusterTrunc_ (config.getParameter<double>("clusterTrunc"))
72 {
73  inputToken_ = consumes<SiPixelRecHitCollection>(inputTag_);
74  LogDebug("") << "Using the " << inputTag_ << " input collection";
75 }
76 
78 
79 void
83  desc.add<edm::InputTag>("inputTag",edm::InputTag("hltSiPixelRecHits"));
84  desc.add<double>("minZ",-20.0);
85  desc.add<double>("maxZ",20.05);
86  desc.add<double>("zStep",0.2);
87  std::vector<double> temp; temp.push_back(0.0); temp.push_back(0.0045);
88  desc.add<std::vector<double> >("clusterPars",temp);
89  desc.add<int>("nhitsTrunc",150.);
90  desc.add<double>("clusterTrunc",2.0);
91  descriptions.add("hltPixelClusterShapeFilter",desc);
92 }
93 
94 //
95 // member functions
96 //
97 
98 // ------------ method called to produce the data ------------
100 {
101  // All HLT filters must create and fill an HLT filter object,
102  // recording any reconstructed physics objects satisfying (or not)
103  // this HLT filter, and place it in the Event.
104 
105  // The filter object
106  if (saveTags()) filterproduct.addCollectionTag(inputTag_);
107  bool accept = true;
108 
109  // get hold of products from Event
111  event.getByToken(inputToken_, hRecHits);
112 
113  // get tracker geometry
114  if (hRecHits.isValid()) {
115  edm::ESHandle<TrackerGeometry> trackerHandle;
116  iSetup.get<TrackerDigiGeometryRecord>().get(trackerHandle);
117  const TrackerGeometry *tgeo = trackerHandle.product();
118  const SiPixelRecHitCollection *hits = hRecHits.product();
119 
120  // loop over pixel rechits
121  int nPxlHits=0;
122  std::vector<VertexHit> vhits;
123  for(auto const & hit : hits->data()) {
124  if (!hit.isValid())
125  continue;
126  ++nPxlHits;
127  DetId id(hit.geographicalId());
128  if(id.subdetId() != int(PixelSubdetector::PixelBarrel))
129  continue;
130  const PixelGeomDetUnit *pgdu = static_cast<const PixelGeomDetUnit*>(tgeo->idToDet(id));
131  if (true) {
132  const PixelTopology *pixTopo = &(pgdu->specificTopology());
133  std::vector<SiPixelCluster::Pixel> pixels(hit.cluster()->pixels());
134  bool pixelOnEdge = false;
135  for(std::vector<SiPixelCluster::Pixel>::const_iterator pixel = pixels.begin();
136  pixel != pixels.end(); ++pixel) {
137  int pixelX = pixel->x;
138  int pixelY = pixel->y;
139  if(pixTopo->isItEdgePixelInX(pixelX) || pixTopo->isItEdgePixelInY(pixelY)) {
140  pixelOnEdge = true;
141  break;
142  }
143  }
144  if (pixelOnEdge)
145  continue;
146  }
147 
148  LocalPoint lpos = LocalPoint(hit.localPosition().x(),
149  hit.localPosition().y(),
150  hit.localPosition().z());
151  GlobalPoint gpos = pgdu->toGlobal(lpos);
152  VertexHit vh;
153  vh.z = gpos.z();
154  vh.r = gpos.perp();
155  vh.w = hit.cluster()->sizeY();
156  vhits.push_back(vh);
157  }
158 
159  // estimate z-position from cluster lengths
160  double zest = 0.0;
161  int nhits = 0, nhits_max = 0;
162  double chi = 0, chi_max = 1e+9;
163  for(double z0 = minZ_; z0 <= maxZ_; z0 += zStep_) {
164  nhits = getContainedHits(vhits, z0, chi);
165  if(nhits == 0)
166  continue;
167  if(nhits > nhits_max) {
168  chi_max = 1e+9;
169  nhits_max = nhits;
170  }
171  if(nhits >= nhits_max && chi < chi_max) {
172  chi_max = chi;
173  zest = z0;
174  }
175  }
176 
177  chi = 0;
178  int nbest=0, nminus=0, nplus=0;
179  nbest = getContainedHits(vhits,zest,chi);
180  nminus = getContainedHits(vhits,zest-10.,chi);
181  nplus = getContainedHits(vhits,zest+10.,chi);
182 
183  double clusVtxQual=0.0;
184  if ((nminus+nplus)> 0)
185  clusVtxQual = (2.0*nbest)/(nminus+nplus); // A/B
186  else if (nbest>0)
187  clusVtxQual = 1000.0; // A/0 (set to arbitrarily large number)
188  else
189  clusVtxQual = 0; // 0/0 (already the default)
190 
191  // construct polynomial cut on cluster vertex quality vs. npixelhits
192  double polyCut=0;
193  for(unsigned int i=0; i < clusterPars_.size(); i++) {
194  polyCut += clusterPars_[i]*std::pow((double)nPxlHits,(int)i);
195  }
196  if(nPxlHits < nhitsTrunc_)
197  polyCut=0; // don't use cut below nhitsTrunc_ pixel hits
198  if(polyCut > clusterTrunc_ && clusterTrunc_ > 0)
199  polyCut=clusterTrunc_; // no cut above clusterTrunc_
200 
201  if (clusVtxQual < polyCut) accept = false;
202  }
203 
204  // return with final filter decision
205  return accept;
206 }
207 
208 int HLTPixelClusterShapeFilter::getContainedHits(const std::vector<VertexHit> &hits, double z0, double &chi) const
209 {
210  // Calculate number of hits contained in v-shaped window in cluster y-width vs. z-position.
211  int n = 0;
212  chi = 0.;
213 
214  for(auto hit : hits) {
215  double p = 2 * fabs(hit.z - z0)/hit.r + 0.5; // FIXME
216  if(fabs(p - hit.w) <= 1.) {
217  chi += fabs(p - hit.w);
218  n++;
219  }
220  }
221  return n;
222 }
223 
224 
225 // define as a framework module
#define LogDebug(id)
Point3DBase< Scalar, LocalTag > LocalPoint
Definition: Definitions.h:32
virtual bool isItEdgePixelInX(int ixbin) const =0
~HLTPixelClusterShapeFilter() override
GlobalPoint toGlobal(const Local2DPoint &lp) const
Conversion to the global R.F. from the R.F. of the GeomDet.
Definition: GeomDet.h:54
virtual bool isItEdgePixelInY(int iybin) const =0
int getContainedHits(const std::vector< VertexHit > &hits, double z0, double &chi) const
Definition: config.py:1
bool hltFilter(edm::Event &, const edm::EventSetup &, trigger::TriggerFilterObjectWithRefs &filterproduct) const override
bool accept(const edm::Event &event, const edm::TriggerResults &triggerTable, const std::string &triggerPath)
Definition: TopDQMHelpers.h:30
edm::EDGetTokenT< SiPixelRecHitCollection > inputToken_
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
data_type const * data(size_t cell) const
ParameterDescriptionBase * add(U const &iLabel, T const &value)
bool isValid() const
Definition: HandleBase.h:74
Definition: DetId.h:18
static void makeHLTFilterDescription(edm::ParameterSetDescription &desc)
Definition: HLTFilter.cc:29
T const * product() const
Definition: Handle.h:74
virtual const PixelTopology & specificTopology() const
Returns a reference to the pixel proxy topology.
void addCollectionTag(const edm::InputTag &collectionTag)
collectionTags
void add(std::string const &label, ParameterSetDescription const &psetDescription)
bool saveTags() const
Definition: HLTFilter.h:45
HLT enums.
T get() const
Definition: EventSetup.h:71
const TrackerGeomDet * idToDet(DetId) const override
HLTPixelClusterShapeFilter(const edm::ParameterSet &)
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
T const * product() const
Definition: ESHandle.h:86
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:40
Definition: event.py:1