CMS 3D CMS Logo

CSCSegAlgoDF.cc
Go to the documentation of this file.
1 
11 #include "CSCSegAlgoDF.h"
12 #include "CSCSegFit.h"
13 
15 
17 
20 
22 
24 #include "CSCSegAlgoShowering.h"
25 
26 #include <algorithm>
27 #include <cmath>
28 #include <iostream>
29 #include <string>
30 
31 /* Constructor
32  *
33  */
35  : CSCSegmentAlgorithm(ps), myName("CSCSegAlgoDF"), sfit_(nullptr) {
36  debug = ps.getUntrackedParameter<bool>("CSCSegmentDebug");
37  minLayersApart = ps.getParameter<int>("minLayersApart");
38  minHitsPerSegment = ps.getParameter<int>("minHitsPerSegment");
39  dRPhiFineMax = ps.getParameter<double>("dRPhiFineMax");
40  dPhiFineMax = ps.getParameter<double>("dPhiFineMax");
41  tanThetaMax = ps.getParameter<double>("tanThetaMax");
42  tanPhiMax = ps.getParameter<double>("tanPhiMax");
43  chi2Max = ps.getParameter<double>("chi2Max");
44  preClustering = ps.getUntrackedParameter<bool>("preClustering");
45  minHitsForPreClustering = ps.getParameter<int>("minHitsForPreClustering");
46  nHitsPerClusterIsShower = ps.getParameter<int>("nHitsPerClusterIsShower");
47  Pruning = ps.getUntrackedParameter<bool>("Pruning");
48  maxRatioResidual = ps.getParameter<double>("maxRatioResidualPrune");
49 
52 }
53 
54 /* Destructor
55  *
56  */
58  delete preCluster_;
59  delete showering_;
60 }
61 
62 /* run
63  *
64  */
65 std::vector<CSCSegment> CSCSegAlgoDF::run(const CSCChamber* aChamber, const ChamberHitContainer& rechits) {
66  // Store chamber info in temp memory
67  theChamber = aChamber;
68 
69  int nHits = rechits.size();
70 
71  // Segments prior to pruning
72  std::vector<CSCSegment> segments_temp;
73 
75  // This is where the segment origin is in the chamber on avg.
76  std::vector<CSCSegment> testSegments;
77  std::vector<ChamberHitContainer> clusteredHits = preCluster_->clusterHits(theChamber, rechits);
78  // loop over the found clusters:
79  for (std::vector<ChamberHitContainer>::iterator subrechits = clusteredHits.begin();
80  subrechits != clusteredHits.end();
81  ++subrechits) {
82  // build the subset of segments:
83  std::vector<CSCSegment> segs = buildSegments((*subrechits));
84  // add the found subset of segments to the collection of all segments in this chamber:
85  segments_temp.insert(segments_temp.end(), segs.begin(), segs.end());
86  }
87  } else {
88  std::vector<CSCSegment> segs = buildSegments(rechits);
89  // add the found subset of segments to the collection of all segments in this chamber:
90  segments_temp.insert(segments_temp.end(), segs.begin(), segs.end());
91  }
92 
93  return segments_temp;
94 }
95 
96 /* This builds segments by first creating proto-segments from at least 3 hits.
97  * We intend to try all possible pairs of hits to start segment building. 'All possible'
98  * means each hit lies on different layers in the chamber. Once a hit has been assigned
99  * to a segment, we don't consider it again, THAT IS, FOR THE FIRST PASS ONLY !
100  * In fact, this is one of the possible flaw with the SK algorithms as it sometimes manages
101  * to build segments with the wrong starting points. In the DF algorithm, the endpoints
102  * are tested as the best starting points in a 2nd loop.
103  *
104  * Also, a maximum of 5 segments is allowed in the chamber (and then it just gives up.)
105  *
106  * @@ There are 7 return's in this function!
107  *
108  */
109 std::vector<CSCSegment> CSCSegAlgoDF::buildSegments(const ChamberHitContainer& _rechits) {
110  ChamberHitContainer rechits = _rechits;
111  // Clear buffer for segment vector
112  std::vector<CSCSegment> segmentInChamber;
113  segmentInChamber.clear();
114 
115  unsigned nHitInChamber = rechits.size();
116 
117  // std::cout << "[CSCSegAlgoDF::buildSegments] address of chamber = " << theChamber << std::endl;
118  // std::cout << "[CSCSegAlgoDF::buildSegments] starting in " << theChamber->id()
119  // << " with " << nHitInChamber << " rechits" << std::endl;
120 
121  // Return #1 - OK, there aren't enough rechits to build a segment
122  if (nHitInChamber < 3)
123  return segmentInChamber;
124 
125  LayerIndex layerIndex(nHitInChamber);
126 
127  size_t nLayers = 0;
128  size_t old_layer = 0;
129  for (size_t i = 0; i < nHitInChamber; ++i) {
130  size_t this_layer = rechits[i]->cscDetId().layer();
131  // std::cout << "[CSCSegAlgoDF::buildSegments] this_layer = " << this_layer << std::endl;
132  layerIndex[i] = this_layer;
133  // std::cout << "[CSCSegAlgoDF::buildSegments] layerIndex[" << i << "] = " << layerIndex[i] << std::endl;
134  if (this_layer != old_layer) {
135  old_layer = this_layer;
136  ++nLayers;
137  }
138  }
139 
140  // std::cout << "[CSCSegAlgoDF::buildSegments] layers with rechits = " << nLayers << std::endl;
141 
142  // Return #2 - OK, there aren't enough hit layers to build a segment
143  if (nLayers < 3)
144  return segmentInChamber;
145 
146  double z1 = theChamber->layer(1)->position().z();
147  double z6 = theChamber->layer(6)->position().z();
148 
149  if (z1 > 0.) {
150  if (z1 > z6) {
151  reverse(layerIndex.begin(), layerIndex.end());
152  reverse(rechits.begin(), rechits.end());
153  }
154  } else if (z1 < 0.) {
155  if (z1 < z6) {
156  reverse(layerIndex.begin(), layerIndex.end());
157  reverse(rechits.begin(), rechits.end());
158  }
159  }
160 
161  // std::cout << "[CSCSegAlgoDF::buildSegments] rechits have been ordered" << std::endl;
162 
163  // Showering muon
164  if (preClustering && int(nHitInChamber) > nHitsPerClusterIsShower && nLayers > 2) {
165  // std::cout << "[CSCSegAlgoDF::buildSegments] showering block" << std::endl;
166 
168 
169  // Return #3 - OK, this is now 'effectve' rechits
170  // Make sure have at least 3 hits...
171  if (segShower.nRecHits() < 3)
172  return segmentInChamber;
173 
174  segmentInChamber.push_back(segShower);
175  if (debug)
176  dumpSegment(segShower);
177 
178  // Return #4 - OK, only get one try at building a segment from shower
179  return segmentInChamber;
180  }
181 
182  // Initialize flags that a given hit has been allocated to a segment
183  BoolContainer used_ini(rechits.size(), false);
184  usedHits = used_ini;
185 
187  ChamberHitContainerCIt ie = rechits.end();
188 
189  // std::cout << "[CSCSegAlgoDF::buildSegments] entering rechit loop" << std::endl;
190 
191  // Now Loop over hits within the chamber to find 1st seed for segment building
192  for (ChamberHitContainerCIt i1 = ib; i1 < ie; ++i1) {
193  if (usedHits[i1 - ib])
194  continue;
195 
196  const CSCRecHit2D* h1 = *i1;
197  int layer1 = layerIndex[i1 - ib];
198  const CSCLayer* l1 = theChamber->layer(layer1);
199  GlobalPoint gp1 = l1->toGlobal(h1->localPosition());
200  LocalPoint lp1 = theChamber->toLocal(gp1);
201 
202  // Loop over hits backward to find 2nd seed for segment building
203  for (ChamberHitContainerCIt i2 = ie - 1; i2 > ib; --i2) {
204  if (usedHits[i2 - ib])
205  continue; // Hit has been used already
206 
207  int layer2 = layerIndex[i2 - ib];
208  if ((layer2 - layer1) < minLayersApart)
209  continue;
210 
211  const CSCRecHit2D* h2 = *i2;
212  const CSCLayer* l2 = theChamber->layer(layer2);
213  GlobalPoint gp2 = l2->toGlobal(h2->localPosition());
214  LocalPoint lp2 = theChamber->toLocal(gp2);
215 
216  // Clear proto segment so it can be (re)-filled
217  protoSegment.clear();
218 
219  // We want hit wrt chamber (and local z will be != 0)
220  float dz = gp2.z() - gp1.z();
221  float slope_u = (lp2.x() - lp1.x()) / dz;
222  float slope_v = (lp2.y() - lp1.y()) / dz;
223 
224  // Test if entrance angle is roughly pointing towards IP
225  if (fabs(slope_v) > tanThetaMax)
226  continue;
227  if (fabs(slope_u) > tanPhiMax)
228  continue;
229 
230  protoSegment.push_back(h1);
231  protoSegment.push_back(h2);
232 
233  // std::cout << "[CSCSegAlgoDF::buildSegments] about to fit 2 hits on layers "
234  // << layer1 << " and " << layer2 << std::endl;
235 
236  // protoSegment has just 2 hits - but need to create a CSCSegFit to hold it in case
237  // we fail to add any more hits
239 
240  // Try adding hits to proto segment
241  tryAddingHitsToSegment(rechits, i1, i2, layerIndex);
242 
243  // Check no. of hits on segment to see if segment is large enough
244  bool segok = true;
245  unsigned iadd = 0;
246 
247  if (protoSegment.size() < minHitsPerSegment + iadd)
248  segok = false;
249 
250  if (Pruning && segok)
252 
253  // Check if segment satisfies chi2 requirement
254  if (sfit_->chi2() > chi2Max)
255  segok = false;
256 
257  if (segok) {
258  // Create an actual CSCSegment - retrieve all info from the current fit
260  // std::cout << "[CSCSegAlgoDF::buildSegments] about to delete sfit= = " << sfit_ << std::endl;
261  delete sfit_;
262  sfit_ = nullptr; // avoid possibility of attempting a second delete later
263 
264  segmentInChamber.push_back(temp);
265  if (debug)
266  dumpSegment(temp);
267 
268  // Return #5 - OK, fewer than 3 rechits not on this segment left in chamber
269  if (nHitInChamber - protoSegment.size() < 3)
270  return segmentInChamber;
271  // Return $6 - already have more than 4 segments in this chamber
272  if (segmentInChamber.size() > 4)
273  return segmentInChamber;
274 
275  // Flag used hits
277  }
278  }
279  }
280  // Return #7
281  return segmentInChamber;
282 }
283 
284 /* Method tryAddingHitsToSegment
285  *
286  * Look at left over hits and try to add them to proto segment by looking how far they
287  * are from the segment in terms of the hit error matrix (so how many sigmas away).
288  *
289  */
293  const LayerIndex& layerIndex) {
294  /* Iterate over the layers with hits in the chamber
295  * Skip the layers containing the segment endpoints on first pass, but then
296  * try hits on layer containing the segment starting points on 2nd pass
297  * if segment has >2 hits. Once a hit is added to a layer, don't replace it
298  * until 2nd iteration
299  */
300 
301  // std::cout << "[CSCSegAlgoDF::tryAddingHitsToSegment] entering"
302  // << " with rechits.size() = " << rechits.size() << std::endl;
303 
305  ChamberHitContainerCIt ie = rechits.end();
306  closeHits.clear();
307 
308  // int counter1 = 0;
309  // int counter2 = 0;
310 
311  for (ChamberHitContainerCIt i = ib; i != ie; ++i) {
312  // std::cout << "counter1 = " << ++counter1 << std::endl;
313  if (i == i1 || i == i2)
314  continue;
315  if (usedHits[i - ib])
316  continue; // Don't use hits already part of a segment.
317 
318  // std::cout << "counter2 = " << ++counter2 << std::endl;
319  const CSCRecHit2D* h = *i;
320  int layer = layerIndex[i - ib];
321  int layer1 = layerIndex[i1 - ib];
322  int layer2 = layerIndex[i2 - ib];
323 
324  // std::cout << "[CSCSegAlgoDF::tryAddingHitsToSegment] layer, layer1, layer2 = "
325  // << layer << ", " << layer1 << ", " << layer2 << std::endl;
326 
327  // Low multiplicity case
328  // only adds hit to protoSegment if no hit on layer already; otherwise adds to closeHits
329  if (rechits.size() < 9) {
330  // std::cout << "low mult" << std::endl;
331  if (isHitNearSegment(h)) {
332  if (!hasHitOnLayer(layer)) {
333  addHit(h, layer);
334  } else {
335  closeHits.push_back(h);
336  }
337  }
338 
339  // High multiplicity case
340  // only adds hit to protoSegment if no hit on layer already AND then refits; otherwise adds to closeHits
341  } else {
342  // std::cout << "high mult" << std::endl;
343  if (isHitNearSegment(h)) {
344  // std::cout << "near seg" << std::endl;
345  if (!hasHitOnLayer(layer)) {
346  // std::cout << "no hit on layer" << std::endl;
347  if (addHit(h, layer)) {
348  // std::cout << "update fit" << std::endl;
350  }
351  // Don't change the starting points at this stage !!!
352  } else {
353  // std::cout << "already hit on layer" << std::endl;
354  closeHits.push_back(h);
355  if (layer != layer1 && layer != layer2)
357  }
358  }
359  }
360  }
361 
362  if (int(protoSegment.size()) < 3)
363  return;
364  // std::cout << "final fit" << std::endl;
366 
367  // 2nd pass to remove biases
368  // This time, also consider changing the endpoints
369  for (ChamberHitContainerCIt i = closeHits.begin(); i != closeHits.end(); ++i) {
370  // std::cout << "2nd pass" << std::endl;
371  const CSCRecHit2D* h = *i;
372  int layer = (*i)->cscDetId().layer();
374  }
375 }
376 
377 /* isHitNearSegment
378  *
379  * Compare rechit with expected position from protosegment
380  */
382  const CSCLayer* layer = theChamber->layer(hit->cscDetId().layer());
383 
384  // hit phi position in global coordinates
385  GlobalPoint Hgp = layer->toGlobal(hit->localPosition());
386  float Hphi = Hgp.phi();
387  if (Hphi < 0.)
388  Hphi += 2. * M_PI;
389  LocalPoint Hlp = theChamber->toLocal(Hgp);
390  float z = Hlp.z();
391 
392  float LocalX = sfit_->xfit(z);
393  float LocalY = sfit_->yfit(z);
394 
395  LocalPoint Slp(LocalX, LocalY, z);
396  GlobalPoint Sgp = theChamber->toGlobal(Slp);
397  float Sphi = Sgp.phi();
398  if (Sphi < 0.)
399  Sphi += 2. * M_PI;
400  float R = sqrt(Sgp.x() * Sgp.x() + Sgp.y() * Sgp.y());
401 
402  float deltaPhi = Sphi - Hphi;
403  if (deltaPhi > 2. * M_PI)
404  deltaPhi -= 2. * M_PI;
405  if (deltaPhi < -2. * M_PI)
406  deltaPhi += 2. * M_PI;
407  if (deltaPhi < 0.)
408  deltaPhi = -deltaPhi;
409 
410  float RdeltaPhi = R * deltaPhi;
411 
412  if (RdeltaPhi < dRPhiFineMax && deltaPhi < dPhiFineMax)
413  return true;
414 
415  return false;
416 }
417 
418 /* Method addHit
419  *
420  * Test if can add hit to proto segment. If so, try to add it.
421  *
422  */
423 bool CSCSegAlgoDF::addHit(const CSCRecHit2D* aHit, int layer) {
424  // std::cout << "[CSCSegAlgoDF::addHit] on layer " << layer << " to protoSegment.size() = "
425  // << protoSegment.size() << std::endl;
426 
427  // Return true if hit was added successfully and then parameters are updated.
428  // Return false if there is already a hit on the same layer, or insert failed.
429 
430  if (protoSegment.size() > 5)
431  return false; //@@ can only have 6 hits at most
432 
433  // Test that we are not trying to add the same hit again
434  for (ChamberHitContainer::const_iterator it = protoSegment.begin(); it != protoSegment.end(); ++it)
435  if (aHit == (*it))
436  return false;
437 
438  protoSegment.push_back(aHit);
439 
440  return true;
441 }
442 
443 /* Method updateParameters
444  *
445  * Perform a simple Least Square Fit on proto segment to determine slope and intercept
446  *
447  */
449  // Delete existing CSCSegFit, create a new one and make the fit
450  // Uses internal variables - theChamber & protoSegment
451 
452  // std::cout << "[CSCSegAlgoDF::updateParameters] about to delete sfit_ = " << sfit_ << std::endl;
453  delete sfit_;
454  // std::cout << "[CSCSegAlgoDF::updateParameters] protoSegment.size() = "
455  // << protoSegment.size() << std::endl;
456  // std::cout << "[CSCSegAlgoDF::updateParameters] theChamber = " << theChamber << std::endl;
458  // std::cout << "[CSCSegAlgoDF::updateParameters] new sfit_ = " << sfit_ << std::endl;
459  sfit_->fit();
460 }
461 
462 /* hasHitOnLayer
463  *
464  * Just make sure hit to be added to layer comes from different layer than those in proto segment
465  */
467  // std::cout << "[CSCSegAlgoDF::hasHitOnLayer] on layer " << layer << std::endl;
468 
469  // Is there already a hit on this layer?
470  for (ChamberHitContainerCIt it = protoSegment.begin(); it != protoSegment.end(); it++)
471  if ((*it)->cscDetId().layer() == layer)
472  return true;
473 
474  return false;
475 }
476 
477 /* Method compareProtoSegment
478  *
479  * For hit coming from the same layer of an existing hit within the proto segment
480  * test if achieve better chi^2 by using this hit than the other
481  *
482  */
484  // std::cout << "[CSCSegAlgoDF::compareProtoSegment] for hit on layer " << layer
485  // << " with protoSegment.size() = " << protoSegment.size() << std::endl;
486 
487  // Try adding the hit to existing segment, and remove old one existing in same layer
488  ChamberHitContainer::iterator it;
489  for (it = protoSegment.begin(); it != protoSegment.end();) {
490  if ((*it)->cscDetId().layer() == layer) {
491  it = protoSegment.erase(it);
492  } else {
493  ++it;
494  }
495  }
496 
497  // std::cout << "[CSCSegAlgoDF::compareProtoSegment] about to add hit on layer " << layer
498  // << " with protoSegment.size() = " << protoSegment.size() << std::endl;
499 
500  bool ok = addHit(h, layer);
501 
502  CSCSegFit* newfit = nullptr;
503  if (ok) {
504  newfit = new CSCSegFit(theChamber, protoSegment);
505  // std::cout << "[CSCSegAlgoDF::compareProtoSegment] newfit = " << newfit << std::endl;
506  newfit->fit();
507  }
508  if (!ok || (newfit->chi2() > sfit_->chi2())) {
509  // std::cout << "[CSCSegAlgoDF::compareProtoSegment] about to delete newfit = " << newfit << std::endl;
510  delete newfit; // failed to add a hit or new fit is worse
511  } else {
512  // std::cout << "[CSCSegAlgoDF::compareProtoSegment] about to delete sfit_ = " << sfit_ << std::endl;
513  delete sfit_; // new fit is better
514  sfit_ = newfit;
515  // std::cout << "[CSCSegAlgoDF::compareProtoSegment] reset sfit_ = " << sfit_ << std::endl;
516  }
517 }
518 
519 /* Method flagHitsAsUsed
520  *
521  * Flag hits which have entered segment building so we don't reuse them.
522  * Also flag does which were very close to segment to reduce combinatorics
523  */
524 void CSCSegAlgoDF::flagHitsAsUsed(const ChamberHitContainer& rechitsInChamber) {
525  // Flag hits on segment as used
526  ChamberHitContainerCIt ib = rechitsInChamber.begin();
528 
529  for (hi = protoSegment.begin(); hi != protoSegment.end(); ++hi) {
530  for (iu = ib; iu != rechitsInChamber.end(); ++iu) {
531  if (*hi == *iu)
532  usedHits[iu - ib] = true;
533  }
534  }
535  // Don't reject hits marked as "nearby" for now.
536  // So this is bypassed at all times for now !!!
537  // Perhaps add back to speed up algorithm some more
538  if (!closeHits.empty())
539  return;
540  for (hi = closeHits.begin(); hi != closeHits.end(); ++hi) {
541  for (iu = ib; iu != rechitsInChamber.end(); ++iu) {
542  if (*hi == *iu)
543  usedHits[iu - ib] = true;
544  }
545  }
546 }
547 
548 // Try to clean up segments by quickly looking at residuals
550  // Only prune if have at least 5 hits
551  if (protoSegment.size() < 5)
552  return;
553 
554  // Now Study residuals
555 
556  float maxResidual = 0.;
557  float sumResidual = 0.;
558  int nHits = 0;
559  int badIndex = -1;
560  int j = 0;
561 
562  ChamberHitContainer::const_iterator ih;
563 
564  for (ih = protoSegment.begin(); ih != protoSegment.end(); ++ih) {
565  const CSCRecHit2D& hit = (**ih);
566  const CSCLayer* layer = theChamber->layer(hit.cscDetId().layer());
567  GlobalPoint gp = layer->toGlobal(hit.localPosition());
569 
570  float residual = sfit_->Rdev(lp.x(), lp.y(), lp.z());
571 
572  sumResidual += residual;
573  nHits++;
574  if (residual > maxResidual) {
575  maxResidual = residual;
576  badIndex = j;
577  }
578  j++;
579  }
580 
581  float corrAvgResidual = (sumResidual - maxResidual) / (nHits - 1);
582 
583  // Keep all hits
584  if (maxResidual / corrAvgResidual < maxRatioResidual)
585  return;
586 
587  // Drop worse hit and recompute segment properties + fill
588 
589  ChamberHitContainer newProtoSegment;
590 
591  j = 0;
592  for (ih = protoSegment.begin(); ih != protoSegment.end(); ++ih) {
593  if (j != badIndex)
594  newProtoSegment.push_back(*ih);
595  j++;
596  }
597 
598  protoSegment.clear();
599 
600  for (ih = newProtoSegment.begin(); ih != newProtoSegment.end(); ++ih) {
601  protoSegment.push_back(*ih);
602  }
603 
604  // Compute segment parameters
606 }
607 
608 void CSCSegAlgoDF::dumpSegment(const CSCSegment& seg) const {
609  edm::LogVerbatim("CSCSegment") << "CSCSegment in " << theChamber->id() << "\nlocal position = " << seg.localPosition()
610  << "\nerror = " << seg.localPositionError()
611  << "\nlocal direction = " << seg.localDirection()
612  << "\nerror =" << seg.localDirectionError() << "\ncovariance matrix"
613  << seg.parametersError() << "chi2/ndf = " << seg.chi2() << "/"
614  << seg.degreesOfFreedom() << "\n#rechits = " << seg.specificRecHits().size()
615  << "\ntime = " << seg.time();
616 }
void tryAddingHitsToSegment(const ChamberHitContainer &rechitsInChamber, const ChamberHitContainerCIt i1, const ChamberHitContainerCIt i2, const LayerIndex &layerIndex)
Utility functions.
Log< level::Info, true > LogVerbatim
void fit(void)
Definition: CSCSegFit.cc:13
ChamberHitContainer protoSegment
Definition: CSCSegAlgoDF.h:118
std::vector< const CSCRecHit2D * > ChamberHitContainer
Definition: CSCSegAlgoDF.h:53
std::vector< const CSCRecHit2D * >::const_iterator ChamberHitContainerCIt
Definition: CSCSegAlgoDF.h:54
T getParameter(std::string const &) const
Definition: ParameterSet.h:307
const CSCLayer * layer(CSCDetId id) const
Return the layer corresponding to the given id.
Definition: CSCChamber.cc:30
double chi2(void) const
Definition: CSCSegFit.h:82
bool isHitNearSegment(const CSCRecHit2D *h) const
float time() const
Definition: CSCSegment.cc:144
LocalPoint localPosition() const override
Definition: CSCSegment.h:39
LocalPoint toLocal(const GlobalPoint &gp) const
Conversion to the R.F. of the GeomDet.
Definition: GeomDet.h:58
T z() const
Definition: PV3DBase.h:61
CSCSegAlgoPreClustering * preCluster_
Definition: CSCSegAlgoDF.h:139
CSCSegment showerSeg(const CSCChamber *aChamber, const ChamberHitContainer &rechits)
Geom::Phi< T > phi() const
Definition: PV3DBase.h:66
LocalPoint intercept() const
Definition: CSCSegFit.h:84
ChamberHitContainer closeHits
Definition: CSCSegAlgoDF.h:116
bool addHit(const CSCRecHit2D *hit, int layer)
float yfit(float z) const
Definition: CSCSegFit.cc:432
bool hasHitOnLayer(int layer) const
CSCSegAlgoDF(const edm::ParameterSet &ps)
Constructor.
Definition: CSCSegAlgoDF.cc:34
CSCSegAlgoShowering * showering_
Definition: CSCSegAlgoDF.h:140
double chi2() const override
Chi2 of the segment fit.
Definition: CSCSegment.h:58
bool preClustering
Definition: CSCSegAlgoDF.h:123
constexpr std::array< uint8_t, layerIndexSize > layer
std::vector< CSCSegment > buildSegments(const ChamberHitContainer &rechits)
T getUntrackedParameter(std::string const &, T const &) const
LocalVector localDirection() const override
Local direction.
Definition: CSCSegment.h:42
void dumpSegment(const CSCSegment &seg) const
T x() const
Definition: PV3DBase.h:59
T y() const
Definition: PV3DBase.h:60
CSCSegFit * sfit_
Definition: CSCSegAlgoDF.h:141
void updateParameters(void)
std::deque< bool > BoolContainer
Definition: CSCSegAlgoDF.h:55
Definition: EPCuts.h:4
T sqrt(T t)
Definition: SSEVec.h:19
void compareProtoSegment(const CSCRecHit2D *h, int layer)
int nRecHits() const
Definition: CSCSegment.h:68
std::vector< int > LayerIndex
Typedefs.
Definition: CSCSegAlgoDF.h:52
int degreesOfFreedom() const override
Degrees of freedom of the segment fit.
Definition: CSCSegment.h:62
float maxRatioResidual
Definition: CSCSegAlgoDF.h:137
int minHitsForPreClustering
Definition: CSCSegAlgoDF.h:124
GlobalPoint toGlobal(const Local2DPoint &lp) const
Conversion to the global R.F. from the R.F. of the GeomDet.
Definition: GeomDet.h:49
#define M_PI
LocalPoint localPosition() const override
Definition: CSCRecHit2D.h:56
LocalError localDirectionError() const override
Error on the local direction.
Definition: CSCSegment.cc:52
void flagHitsAsUsed(const ChamberHitContainer &rechitsInChamber)
float tanThetaMax
Definition: CSCSegAlgoDF.h:135
void pruneFromResidual(void)
AlgebraicSymMatrix covarianceMatrix(void)
Definition: CSCSegFit.cc:352
double dPhiFineMax
Definition: CSCSegAlgoDF.h:133
AlgebraicSymMatrix parametersError() const override
Covariance matrix of parameters()
Definition: CSCSegment.h:49
caConstants::TupleMultiplicity const CAHitNtupletGeneratorKernelsGPU::HitToTuple const cms::cuda::AtomicPairCounter GPUCACell const *__restrict__ uint32_t const *__restrict__ gpuPixelDoublets::CellNeighborsVector const gpuPixelDoublets::CellTracksVector const GPUCACell::OuterHitOfCell const int32_t nHits
std::vector< std::vector< const CSCRecHit2D * > > clusterHits(const CSCChamber *aChamber, const ChamberHitContainer &rechits)
clusterize
float xfit(float z) const
Definition: CSCSegFit.cc:426
const Surface::PositionType & position() const
The position (origin of the R.F.)
Definition: GeomDet.h:43
std::vector< CSCSegment > run(const CSCChamber *aChamber, const ChamberHitContainer &rechits) override
Definition: CSCSegAlgoDF.cc:65
double dRPhiFineMax
Definition: CSCSegAlgoDF.h:132
~CSCSegAlgoDF() override
Destructor.
Definition: CSCSegAlgoDF.cc:57
const std::vector< CSCRecHit2D > & specificRecHits() const
Definition: CSCSegment.h:66
LocalVector localdir() const
Definition: CSCSegFit.h:85
const CSCChamber * theChamber
Definition: CSCSegAlgoDF.h:113
CSCSetOfHits hits(void) const
Definition: CSCSegFit.h:79
LocalError localPositionError() const override
Definition: CSCSegment.cc:48
int minHitsPerSegment
Definition: CSCSegAlgoDF.h:130
The Signals That Services Can Subscribe To This is based on ActivityRegistry h
Helper function to determine trigger accepts.
Definition: Activities.doc:4
int nHitsPerClusterIsShower
Definition: CSCSegAlgoDF.h:128
float Rdev(float x, float y, float z) const
Definition: CSCSegFit.cc:438
BoolContainer usedHits
Definition: CSCSegAlgoDF.h:114
ib
Definition: cuy.py:661
CSCDetId id() const
Get the (concrete) DetId.
Definition: CSCChamber.h:34