CMS 3D CMS Logo

List of all members | Classes | Public Types | Public Member Functions | Static Public Member Functions | Static Public Attributes | Private Member Functions | Private Attributes
CAHitQuadrupletGenerator Class Reference

#include <CAHitQuadrupletGenerator.h>

Inheritance diagram for CAHitQuadrupletGenerator:
HitQuadrupletGenerator OrderedHitsGenerator

Classes

class  QuantityDependsPt
 
class  QuantityDependsPtEval
 

Public Types

typedef LayerHitMapCache LayerCacheType
 
typedef OrderedHitSeeds ResultType
 

Public Member Functions

 CAHitQuadrupletGenerator (const edm::ParameterSet &cfg, edm::ConsumesCollector &&iC, bool needSeedingLayerSetsHits=true)
 
 CAHitQuadrupletGenerator (const edm::ParameterSet &cfg, edm::ConsumesCollector &iC, bool needSeedingLayerSetsHits=true)
 
void hitNtuplets (const IntermediateHitDoublets &regionDoublets, std::vector< OrderedHitSeeds > &result, const edm::EventSetup &es, const SeedingLayerSetsHits &layers)
 
virtual void hitQuadruplets (const TrackingRegion &reg, OrderedHitSeeds &quadruplets, const edm::Event &ev, const edm::EventSetup &es)
 from base class More...
 
void initEvent (const edm::Event &ev, const edm::EventSetup &es)
 
virtual ~CAHitQuadrupletGenerator ()
 
- Public Member Functions inherited from HitQuadrupletGenerator
virtual void clear () final
 
 HitQuadrupletGenerator (unsigned int size=500)
 
virtual const OrderedHitSeedsrun (const TrackingRegion &region, const edm::Event &ev, const edm::EventSetup &es) final
 
virtual ~HitQuadrupletGenerator ()
 
- Public Member Functions inherited from OrderedHitsGenerator
 OrderedHitsGenerator ()
 
virtual ~OrderedHitsGenerator ()
 

Static Public Member Functions

static void fillDescriptions (edm::ParameterSetDescription &desc)
 
static const char * fillDescriptionsLabel ()
 

Static Public Attributes

static unsigned int minLayers = 4
 

Private Member Functions

void hitQuadruplets (const TrackingRegion &reg, OrderedHitSeeds &result, std::vector< const HitDoublets * > &hitDoublets, CAGraph &g, const edm::EventSetup &es)
 

Private Attributes

const float caHardPtCut = 0.f
 
const bool caOnlyOneLastHitPerLayerFilter = false
 
const float caPhiCut = 0.1f
 
const float caThetaCut = 0.00125f
 
const float extraHitRPhitolerance
 
const bool fitFastCircle
 
const bool fitFastCircleChi2Cut
 
const QuantityDependsPt maxChi2
 
std::unique_ptr< SeedComparitortheComparitor
 
LayerCacheType theLayerCache
 
edm::EDGetTokenT< SeedingLayerSetsHitstheSeedingLayerToken
 
const bool useBendingCorrection
 

Additional Inherited Members

- Public Attributes inherited from OrderedHitsGenerator
unsigned int theMaxElement
 

Detailed Description

Definition at line 32 of file CAHitQuadrupletGenerator.h.

Member Typedef Documentation

Definition at line 34 of file CAHitQuadrupletGenerator.h.

Definition at line 37 of file CAHitQuadrupletGenerator.h.

Constructor & Destructor Documentation

CAHitQuadrupletGenerator::CAHitQuadrupletGenerator ( const edm::ParameterSet cfg,
edm::ConsumesCollector &&  iC,
bool  needSeedingLayerSetsHits = true 
)
inline

Definition at line 41 of file CAHitQuadrupletGenerator.h.

References CAHitQuadrupletGenerator_cfi::CAHitQuadrupletGenerator, and looper::cfg.

41 : CAHitQuadrupletGenerator(cfg, iC, needSeedingLayerSetsHits) {}
CAHitQuadrupletGenerator(const edm::ParameterSet &cfg, edm::ConsumesCollector &&iC, bool needSeedingLayerSetsHits=true)
CAHitQuadrupletGenerator::CAHitQuadrupletGenerator ( const edm::ParameterSet cfg,
edm::ConsumesCollector iC,
bool  needSeedingLayerSetsHits = true 
)

Definition at line 46 of file CAHitQuadrupletGenerator.cc.

References edm::ConsumesCollector::consumes(), beamerCreator::create(), edm::ParameterSet::exists(), reco::get(), edm::ParameterSet::getParameter(), AlCaHLTBitMon_QueryRunRegistry::string, theComparitor, and theSeedingLayerToken.

46  :
47 extraHitRPhitolerance(cfg.getParameter<double>("extraHitRPhitolerance")), //extra window in ThirdHitPredictionFromCircle range (divide by R to get phi)
48 maxChi2(cfg.getParameter<edm::ParameterSet>("maxChi2")),
49 fitFastCircle(cfg.getParameter<bool>("fitFastCircle")),
50 fitFastCircleChi2Cut(cfg.getParameter<bool>("fitFastCircleChi2Cut")),
51 useBendingCorrection(cfg.getParameter<bool>("useBendingCorrection")),
52 caThetaCut(cfg.getParameter<double>("CAThetaCut")),
53 caPhiCut(cfg.getParameter<double>("CAPhiCut")),
54 caHardPtCut(cfg.getParameter<double>("CAHardPtCut")),
55 caOnlyOneLastHitPerLayerFilter(cfg.getParameter<bool>("CAOnlyOneLastHitPerLayerFilter"))
56 {
57  if(needSeedingLayerSetsHits)
59 
60  if (cfg.exists("SeedComparitorPSet"))
61  {
62  edm::ParameterSet comparitorPSet =
63  cfg.getParameter<edm::ParameterSet>("SeedComparitorPSet");
64  std::string comparitorName = comparitorPSet.getParameter<std::string>("ComponentName");
65  if (comparitorName != "none")
66  {
67  theComparitor.reset(SeedComparitorFactory::get()->create(comparitorName, comparitorPSet, iC));
68  }
69  }
70 }
EDGetTokenT< ProductType > consumes(edm::InputTag const &tag)
T getParameter(std::string const &) const
def create(alignables, pedeDump, additionalData, outputFile, config)
bool exists(std::string const &parameterName) const
checks if a parameter exists
std::unique_ptr< SeedComparitor > theComparitor
const QuantityDependsPt maxChi2
edm::EDGetTokenT< SeedingLayerSetsHits > theSeedingLayerToken
T get(const Candidate &c)
Definition: component.h:55
CAHitQuadrupletGenerator::~CAHitQuadrupletGenerator ( )
virtual

Definition at line 72 of file CAHitQuadrupletGenerator.cc.

73 {
74 }

Member Function Documentation

void CAHitQuadrupletGenerator::fillDescriptions ( edm::ParameterSetDescription desc)
static

Definition at line 76 of file CAHitQuadrupletGenerator.cc.

References edm::ParameterSetDescription::add(), and AlCaHLTBitMon_QueryRunRegistry::string.

76  {
77  desc.add<double>("extraHitRPhitolerance", 0.1);
78  desc.add<bool>("fitFastCircle", false);
79  desc.add<bool>("fitFastCircleChi2Cut", false);
80  desc.add<bool>("useBendingCorrection", false);
81  desc.add<double>("CAThetaCut", 0.00125);
82  desc.add<double>("CAPhiCut", 10);
83  desc.add<double>("CAHardPtCut", 0);
84  desc.add<bool>("CAOnlyOneLastHitPerLayerFilter",false);
85  edm::ParameterSetDescription descMaxChi2;
86  descMaxChi2.add<double>("pt1", 0.2);
87  descMaxChi2.add<double>("pt2", 1.5);
88  descMaxChi2.add<double>("value1", 500);
89  descMaxChi2.add<double>("value2", 50);
90  descMaxChi2.add<bool>("enabled", true);
91  desc.add<edm::ParameterSetDescription>("maxChi2", descMaxChi2);
92 
93  edm::ParameterSetDescription descComparitor;
94  descComparitor.add<std::string>("ComponentName", "none");
95  descComparitor.setAllowAnything(); // until we have moved SeedComparitor too to EDProducers
96  desc.add<edm::ParameterSetDescription>("SeedComparitorPSet", descComparitor);
97 }
ParameterDescriptionBase * add(U const &iLabel, T const &value)
static const char* CAHitQuadrupletGenerator::fillDescriptionsLabel ( )
inlinestatic

Definition at line 47 of file CAHitQuadrupletGenerator.h.

References ev, g, LayerTriplets::layers(), and mps_fire::result.

47 { return "caHitQuadruplet"; }
void CAHitQuadrupletGenerator::hitNtuplets ( const IntermediateHitDoublets regionDoublets,
std::vector< OrderedHitSeeds > &  result,
const edm::EventSetup es,
const SeedingLayerSetsHits layers 
)

Definition at line 238 of file CAHitQuadrupletGenerator.cc.

References funct::abs(), EnergyCorrector::c, caHardPtCut, caOnlyOneLastHitPerLayerFilter, caPhiCut, caThetaCut, HiEvtPlane_cfi::chi2, FastCircleFit::chi2(), RZLine::chi2(), ThirdHitPredictionFromCircle::curvature(), PixelRecoUtilities::curvature(), GlobalErrorBase< T, ErrorWeightType >::czz(), relativeConstraints::error, CAHitQuadrupletGenerator::QuantityDependsPt::evaluator(), extraHitRPhitolerance, fitFastCircle, fitFastCircleChi2Cut, runEdmFileComparison::found, g, edm::EventSetup::get(), mps_fire::i, hcalTTPDigis_cfi::id, SeedingLayerSetsHits::SeedingLayer::index(), diffTreeTool::index, SurfaceOrientation::inner, IntermediateHitDoublets::LayerPairHitDoublets::innerLayerIndex(), PixelRecoUtilities::inversePt(), gedGsfElectrons_cfi::isBarrel, edm::isNotFinite(), TrackerTopology::layer(), geometryCSVtoXML::line, hpstanc_transforms::max, maxChi2, TrackingRegion::origin(), SurfaceOrientation::outer, IntermediateHitDoublets::LayerPairHitDoublets::outerLayerIndex(), PixelSubdetector::PixelBarrel, point, edm::ESHandle< T >::product(), EnergyCorrector::pt, GlobalErrorBase< T, ErrorWeightType >::rerr(), TrackerTopology::side(), funct::sqr(), mathSSE::sqrt(), GeomDetEnumerators::subDetId, theComparitor, useBendingCorrection, CAHitQuadrupletGenerator::QuantityDependsPtEval::value(), x, PV3DBase< T, PVType, FrameType >::x(), y, PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

241  {
242  CAGraph g;
243 
244  std::vector<const HitDoublets *> hitDoublets;
245 
246  auto layerPairEqual = [](const IntermediateHitDoublets::LayerPairHitDoublets& pair,
249  return pair.innerLayerIndex() == inner && pair.outerLayerIndex() == outer;
250  };
252  es.get<TrackerTopologyRcd>().get(tTopoHand);
253  const TrackerTopology *tTopo=tTopoHand.product();
254 
255  const int numberOfHitsInNtuplet = 4;
256  std::vector<CACell::CAntuplet> foundQuadruplets;
257 
258 
259 
260  int index =0;
261  for(const auto& regionLayerPairs: regionDoublets) {
262 
263  const TrackingRegion& region = regionLayerPairs.region();
264  hitDoublets.clear();
265  foundQuadruplets.clear();
266  if (index == 0){
267  createGraphStructure(layers, g);
268  }
269  else{
270  clearGraphStructure(layers, g);
271  }
272 
273  fillGraph(layers, g, hitDoublets,
274  [&](const SeedingLayerSetsHits::SeedingLayer& inner,
276  std::vector<const HitDoublets *>& hitDoublets) {
277  using namespace std::placeholders;
278  auto found = std::find_if(regionLayerPairs.begin(), regionLayerPairs.end(),
279  std::bind(layerPairEqual, _1, inner.index(), outer.index()));
280  if(found != regionLayerPairs.end()) {
281  hitDoublets.emplace_back(&(found->doublets()));
282  return true;
283  }
284  return false;
285  });
286 
287  CellularAutomaton ca(g);
288 
289  ca.createAndConnectCells(hitDoublets, region, caThetaCut,
291 
292  ca.evolve(numberOfHitsInNtuplet);
293 
294  ca.findNtuplets(foundQuadruplets, numberOfHitsInNtuplet);
295 
296  auto & allCells = ca.getAllCells();
297 
298  const QuantityDependsPtEval maxChi2Eval = maxChi2.evaluator(es);
299 
300  // re-used thoughout
301  std::array<float, 4> bc_r;
302  std::array<float, 4> bc_z;
303  std::array<float, 4> bc_errZ2;
304  std::array<GlobalPoint, 4> gps;
305  std::array<GlobalError, 4> ges;
306  std::array<bool, 4> barrels;
307  bool hasAlreadyPushedACandidate = false;
308  float selectedChi2 = std::numeric_limits<float>::max();
309  unsigned int previousfourthLayerId = 0;
310  std::array<unsigned int, 2> previousCellIds ={{0,0}};
311  unsigned int previousSideId = 0;
312  int previousSubDetId = 0;
313 
314  unsigned int numberOfFoundQuadruplets = foundQuadruplets.size();
315 
316  // Loop over quadruplets
317  for (unsigned int quadId = 0; quadId < numberOfFoundQuadruplets; ++quadId)
318  {
319 
320  auto isBarrel = [](const unsigned id) -> bool
321  {
322  return id == PixelSubdetector::PixelBarrel;
323  };
324  for(unsigned int i = 0; i< 3; ++i)
325  {
326  auto const& ahit = allCells[foundQuadruplets[quadId][i]].getInnerHit();
327  gps[i] = ahit->globalPosition();
328  ges[i] = ahit->globalPositionError();
329  barrels[i] = isBarrel(ahit->geographicalId().subdetId());
330  }
331 
332  auto const& ahit = allCells[foundQuadruplets[quadId][2]].getOuterHit();
333  gps[3] = ahit->globalPosition();
334  ges[3] = ahit->globalPositionError();
335  barrels[3] = isBarrel(ahit->geographicalId().subdetId());
337  {
338  const auto fourthLayerId = tTopo->layer(ahit->geographicalId());
339  const auto sideId = tTopo->side(ahit->geographicalId());
340  const auto subDetId = ahit->geographicalId().subdetId();
341  const auto isTheSameTriplet = (quadId != 0) && (foundQuadruplets[quadId][0] == previousCellIds[0]) && (foundQuadruplets[quadId][1] == previousCellIds[1]);
342  const auto isTheSameFourthLayer = (quadId != 0) && (fourthLayerId == previousfourthLayerId) && (subDetId == previousSubDetId) && (sideId == previousSideId);
343 
344  previousCellIds = {{foundQuadruplets[quadId][0], foundQuadruplets[quadId][1]}};
345  previousfourthLayerId = fourthLayerId;
346 
347 
348  if(!(isTheSameTriplet && isTheSameFourthLayer ))
349  {
350  selectedChi2 = std::numeric_limits<float>::max();
351  hasAlreadyPushedACandidate = false;
352  }
353 
354  }
355  // TODO:
356  // - 'line' is not used for anything
357  // - if we decide to always do the circle fit for 4 hits, we don't
358  // need ThirdHitPredictionFromCircle for the curvature; then we
359  // could remove extraHitRPhitolerance configuration parameter
360  PixelRecoLineRZ line(gps[0], gps[2]);
361  ThirdHitPredictionFromCircle predictionRPhi(gps[0], gps[2], extraHitRPhitolerance);
362  const float curvature = predictionRPhi.curvature(ThirdHitPredictionFromCircle::Vector2D(gps[1].x(), gps[1].y()));
363  const float abscurv = std::abs(curvature);
364  const float thisMaxChi2 = maxChi2Eval.value(abscurv);
365  if (theComparitor)
366  {
367  SeedingHitSet tmpTriplet(allCells[foundQuadruplets[quadId][0]].getInnerHit(),
368  allCells[foundQuadruplets[quadId][2]].getInnerHit(),
369  allCells[foundQuadruplets[quadId][2]].getOuterHit());
370 
371 
372  if (!theComparitor->compatible(tmpTriplet) )
373  {
374  continue;
375  }
376  }
377 
378  float chi2 = std::numeric_limits<float>::quiet_NaN();
379  // TODO: Do we have any use case to not use bending correction?
381  {
382  // Following PixelFitterByConformalMappingAndLine
383  const float simpleCot = ( gps.back().z() - gps.front().z() ) / (gps.back().perp() - gps.front().perp() );
384  const float pt = 1.f / PixelRecoUtilities::inversePt(abscurv, es);
385  for (int i=0; i < 4; ++i)
386  {
387  const GlobalPoint & point = gps[i];
388  const GlobalError & error = ges[i];
389  bc_r[i] = sqrt( sqr(point.x() - region.origin().x()) + sqr(point.y() - region.origin().y()) );
390  bc_r[i] += pixelrecoutilities::LongitudinalBendingCorrection(pt, es)(bc_r[i]);
391  bc_z[i] = point.z() - region.origin().z();
392  bc_errZ2[i] = (barrels[i]) ? error.czz() : error.rerr(point)*sqr(simpleCot);
393  }
394  RZLine rzLine(bc_r, bc_z, bc_errZ2, RZLine::ErrZ2_tag());
395  chi2 = rzLine.chi2();
396  }
397  else
398  {
399  RZLine rzLine(gps, ges, barrels);
400  chi2 = rzLine.chi2();
401  }
402  if (edm::isNotFinite(chi2) || chi2 > thisMaxChi2)
403  {
404  continue;
405  }
406  // TODO: Do we have any use case to not use circle fit? Maybe
407  // HLT where low-pT inefficiency is not a problem?
408  if (fitFastCircle)
409  {
410  FastCircleFit c(gps, ges);
411  chi2 += c.chi2();
412  if (edm::isNotFinite(chi2)) continue;
413  if (fitFastCircleChi2Cut && chi2 > thisMaxChi2) continue;
414  }
416  {
417  if (chi2 < selectedChi2)
418  {
419  selectedChi2 = chi2;
420  if(hasAlreadyPushedACandidate)
421  {
422  result[index].pop_back();
423  }
424  result[index].emplace_back(allCells[foundQuadruplets[quadId][0]].getInnerHit(),
425  allCells[foundQuadruplets[quadId][1]].getInnerHit(),
426  allCells[foundQuadruplets[quadId][2]].getInnerHit(),
427  allCells[foundQuadruplets[quadId][2]].getOuterHit());
428  hasAlreadyPushedACandidate = true;
429  }
430  }
431  else
432  {
433  result[index].emplace_back(allCells[foundQuadruplets[quadId][0]].getInnerHit(),
434  allCells[foundQuadruplets[quadId][1]].getInnerHit(),
435  allCells[foundQuadruplets[quadId][2]].getInnerHit(),
436  allCells[foundQuadruplets[quadId][2]].getOuterHit());
437  }
438  }
439  index++;
440  }
441 
442 }
GlobalPoint const & origin() const
unsigned int side(const DetId &id) const
T y() const
Definition: PV3DBase.h:63
std::unique_ptr< SeedComparitor > theComparitor
The Signals That Services Can Subscribe To This is based on ActivityRegistry and is current per Services can connect to the signals distributed by the ActivityRegistry in order to monitor the activity of the application Each possible callback has some defined which we here list in angle e g
Definition: Activities.doc:4
T inversePt(T curvature, const edm::EventSetup &iSetup)
unsigned int subDetId[21]
bool isNotFinite(T x)
Definition: isFinite.h:10
T curvature(T InversePt, const edm::EventSetup &iSetup)
SeedingLayerSetsHits::LayerIndex innerLayerIndex() const
T sqrt(T t)
Definition: SSEVec.h:18
T z() const
Definition: PV3DBase.h:64
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
Definition: RZLine.h:12
const QuantityDependsPt maxChi2
T rerr(const GlobalPoint &aPoint) const
const T & get() const
Definition: EventSetup.h:56
unsigned int layer(const DetId &id) const
SeedingLayerSetsHits::LayerIndex outerLayerIndex() const
Square< F >::type sqr(const F &f)
Definition: Square.h:13
QuantityDependsPtEval evaluator(const edm::EventSetup &es) const
T x() const
Definition: PV3DBase.h:62
T const * product() const
Definition: ESHandle.h:86
*vegas h *****************************************************used in the default bin number in original ***version of VEGAS is ***a higher bin number might help to derive a more precise ***grade subtle point
Definition: invegas.h:5
void CAHitQuadrupletGenerator::hitQuadruplets ( const TrackingRegion reg,
OrderedHitSeeds quadruplets,
const edm::Event ev,
const edm::EventSetup es 
)
virtual

from base class

Implements HitQuadrupletGenerator.

Definition at line 198 of file CAHitQuadrupletGenerator.cc.

References LayerHitMapCache::clear(), HitPairGeneratorFromLayerPair::doublets(), MillePedeFileConverter_cfg::e, g, edm::Event::getByToken(), SurfaceOrientation::inner, LayerTriplets::layers(), SeedingLayerSetsHits::numberOfLayersInSet(), SurfaceOrientation::outer, theComparitor, theLayerCache, and theSeedingLayerToken.

201 {
203  ev.getByToken(theSeedingLayerToken, hlayers);
204  const SeedingLayerSetsHits& layers = *hlayers;
205  if (layers.numberOfLayersInSet() != 4)
206  throw cms::Exception("Configuration")
207  << "CAHitQuadrupletsGenerator expects SeedingLayerSetsHits::numberOfLayersInSet() to be 4, got "
208  << layers.numberOfLayersInSet();
209 
210  CAGraph g;
211 
212  std::vector<HitDoublets> hitDoublets;
213 
214 
215  HitPairGeneratorFromLayerPair thePairGenerator(0, 1, &theLayerCache);
216 
217  createGraphStructure(layers, g);
218  fillGraph(layers, g, hitDoublets,
221  std::vector<HitDoublets>& hitDoublets) {
222  hitDoublets.emplace_back(thePairGenerator.doublets(region, ev, es, inner, outer));
223  return true;
224  });
225 
226  if (theComparitor)
227  theComparitor->init(ev, es);
228 
229  std::vector<const HitDoublets *> hitDoubletsPtr;
230  hitDoubletsPtr.reserve(hitDoublets.size());
231  for(const auto& e: hitDoublets)
232  hitDoubletsPtr.emplace_back(&e);
233 
234  hitQuadruplets(region, result, hitDoubletsPtr, g, es);
236 }
unsigned short numberOfLayersInSet() const
Get number of layers in each SeedingLayerSets.
std::vector< LayerSetAndLayers > layers(const SeedingLayerSetsHits &sets)
Definition: LayerTriplets.cc:4
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:460
std::unique_ptr< SeedComparitor > theComparitor
The Signals That Services Can Subscribe To This is based on ActivityRegistry and is current per Services can connect to the signals distributed by the ActivityRegistry in order to monitor the activity of the application Each possible callback has some defined which we here list in angle e g
Definition: Activities.doc:4
virtual void hitQuadruplets(const TrackingRegion &reg, OrderedHitSeeds &quadruplets, const edm::Event &ev, const edm::EventSetup &es)
from base class
edm::EDGetTokenT< SeedingLayerSetsHits > theSeedingLayerToken
void CAHitQuadrupletGenerator::hitQuadruplets ( const TrackingRegion reg,
OrderedHitSeeds result,
std::vector< const HitDoublets * > &  hitDoublets,
CAGraph g,
const edm::EventSetup es 
)
private

Definition at line 444 of file CAHitQuadrupletGenerator.cc.

References funct::abs(), EnergyCorrector::c, caHardPtCut, caOnlyOneLastHitPerLayerFilter, caPhiCut, caThetaCut, HiEvtPlane_cfi::chi2, FastCircleFit::chi2(), RZLine::chi2(), CellularAutomaton::createAndConnectCells(), ThirdHitPredictionFromCircle::curvature(), PixelRecoUtilities::curvature(), GlobalErrorBase< T, ErrorWeightType >::czz(), relativeConstraints::error, CAHitQuadrupletGenerator::QuantityDependsPt::evaluator(), CellularAutomaton::evolve(), extraHitRPhitolerance, CellularAutomaton::findNtuplets(), fitFastCircle, fitFastCircleChi2Cut, edm::EventSetup::get(), CellularAutomaton::getAllCells(), mps_fire::i, hcalTTPDigis_cfi::id, PixelRecoUtilities::inversePt(), gedGsfElectrons_cfi::isBarrel, edm::isNotFinite(), TrackerTopology::layer(), geometryCSVtoXML::line, hpstanc_transforms::max, maxChi2, TrackingRegion::origin(), PixelSubdetector::PixelBarrel, point, edm::ESHandle< T >::product(), EnergyCorrector::pt, GlobalErrorBase< T, ErrorWeightType >::rerr(), TrackerTopology::side(), funct::sqr(), mathSSE::sqrt(), GeomDetEnumerators::subDetId, theComparitor, useBendingCorrection, CAHitQuadrupletGenerator::QuantityDependsPtEval::value(), x, PV3DBase< T, PVType, FrameType >::x(), y, PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

447  {
448  //Retrieve tracker topology from geometry
450  es.get<TrackerTopologyRcd>().get(tTopoHand);
451  const TrackerTopology *tTopo=tTopoHand.product();
452 
453  const int numberOfHitsInNtuplet = 4;
454  std::vector<CACell::CAntuplet> foundQuadruplets;
455 
456  CellularAutomaton ca(g);
457 
458  ca.createAndConnectCells(hitDoublets, region, caThetaCut,
460 
461  ca.evolve(numberOfHitsInNtuplet);
462 
463  ca.findNtuplets(foundQuadruplets, numberOfHitsInNtuplet);
464 
465  auto & allCells = ca.getAllCells();
466 
467  const QuantityDependsPtEval maxChi2Eval = maxChi2.evaluator(es);
468 
469  // re-used thoughout
470  std::array<float, 4> bc_r;
471  std::array<float, 4> bc_z;
472  std::array<float, 4> bc_errZ2;
473  std::array<GlobalPoint, 4> gps;
474  std::array<GlobalError, 4> ges;
475  std::array<bool, 4> barrels;
476  bool hasAlreadyPushedACandidate = false;
477  float selectedChi2 = std::numeric_limits<float>::max();
478  unsigned int previousfourthLayerId = 0;
479  std::array<unsigned int, 2> previousCellIds ={{0,0}};
480  unsigned int previousSideId = 0;
481  int previousSubDetId = 0;
482 
483  unsigned int numberOfFoundQuadruplets = foundQuadruplets.size();
484 
485  // Loop over quadruplets
486  for (unsigned int quadId = 0; quadId < numberOfFoundQuadruplets; ++quadId)
487  {
488 
489  auto isBarrel = [](const unsigned id) -> bool
490  {
491  return id == PixelSubdetector::PixelBarrel;
492  };
493  for(unsigned int i = 0; i< 3; ++i)
494  {
495  auto const& ahit = allCells[foundQuadruplets[quadId][i]].getInnerHit();
496  gps[i] = ahit->globalPosition();
497  ges[i] = ahit->globalPositionError();
498  barrels[i] = isBarrel(ahit->geographicalId().subdetId());
499  }
500 
501  auto const& ahit = allCells[foundQuadruplets[quadId][2]].getOuterHit();
502  gps[3] = ahit->globalPosition();
503  ges[3] = ahit->globalPositionError();
504  barrels[3] = isBarrel(ahit->geographicalId().subdetId());
505 
507  {
508  const auto fourthLayerId = tTopo->layer(ahit->geographicalId());
509  const auto sideId = tTopo->side(ahit->geographicalId());
510  const auto subDetId = ahit->geographicalId().subdetId();
511  const auto isTheSameTriplet = (quadId != 0) && (foundQuadruplets[quadId][0] == previousCellIds[0]) && (foundQuadruplets[quadId][1] == previousCellIds[1]);
512  const auto isTheSameFourthLayer = (quadId != 0) && (fourthLayerId == previousfourthLayerId) && (subDetId == previousSubDetId) && (sideId == previousSideId);
513 
514  previousCellIds = {{foundQuadruplets[quadId][0], foundQuadruplets[quadId][1]}};
515  previousfourthLayerId = fourthLayerId;
516 
517 
518  if(!(isTheSameTriplet && isTheSameFourthLayer ))
519  {
520  selectedChi2 = std::numeric_limits<float>::max();
521  hasAlreadyPushedACandidate = false;
522  }
523 
524  }
525 
526 
527  // TODO:
528  // - 'line' is not used for anything
529  // - if we decide to always do the circle fit for 4 hits, we don't
530  // need ThirdHitPredictionFromCircle for the curvature; then we
531  // could remove extraHitRPhitolerance configuration parameter
532  PixelRecoLineRZ line(gps[0], gps[2]);
533  ThirdHitPredictionFromCircle predictionRPhi(gps[0], gps[2], extraHitRPhitolerance);
534  const float curvature = predictionRPhi.curvature(ThirdHitPredictionFromCircle::Vector2D(gps[1].x(), gps[1].y()));
535  const float abscurv = std::abs(curvature);
536  const float thisMaxChi2 = maxChi2Eval.value(abscurv);
537 
538  if (theComparitor)
539  {
540  SeedingHitSet tmpTriplet(allCells[foundQuadruplets[quadId][0]].getInnerHit(),
541  allCells[foundQuadruplets[quadId][2]].getInnerHit(),
542  allCells[foundQuadruplets[quadId][2]].getOuterHit());
543 
544 
545  if (!theComparitor->compatible(tmpTriplet) )
546  {
547  continue;
548  }
549  }
550 
551  float chi2 = std::numeric_limits<float>::quiet_NaN();
552  // TODO: Do we have any use case to not use bending correction?
554  {
555  // Following PixelFitterByConformalMappingAndLine
556  const float simpleCot = ( gps.back().z() - gps.front().z() ) / (gps.back().perp() - gps.front().perp() );
557  const float pt = 1.f / PixelRecoUtilities::inversePt(abscurv, es);
558  for (int i=0; i < 4; ++i)
559  {
560  const GlobalPoint & point = gps[i];
561  const GlobalError & error = ges[i];
562  bc_r[i] = sqrt( sqr(point.x() - region.origin().x()) + sqr(point.y() - region.origin().y()) );
563  bc_r[i] += pixelrecoutilities::LongitudinalBendingCorrection(pt, es)(bc_r[i]);
564  bc_z[i] = point.z() - region.origin().z();
565  bc_errZ2[i] = (barrels[i]) ? error.czz() : error.rerr(point)*sqr(simpleCot);
566  }
567  RZLine rzLine(bc_r, bc_z, bc_errZ2, RZLine::ErrZ2_tag());
568  chi2 = rzLine.chi2();
569  }
570  else
571  {
572  RZLine rzLine(gps, ges, barrels);
573  chi2 = rzLine.chi2();
574  }
575  if (edm::isNotFinite(chi2) || chi2 > thisMaxChi2)
576  {
577  continue;
578 
579  }
580  // TODO: Do we have any use case to not use circle fit? Maybe
581  // HLT where low-pT inefficiency is not a problem?
582  if (fitFastCircle)
583  {
584  FastCircleFit c(gps, ges);
585  chi2 += c.chi2();
586  if (edm::isNotFinite(chi2))
587  continue;
588  if (fitFastCircleChi2Cut && chi2 > thisMaxChi2)
589  continue;
590  }
591 
593  {
594  if (chi2 < selectedChi2)
595  {
596  selectedChi2 = chi2;
597 
598  if(hasAlreadyPushedACandidate)
599  {
600  result.pop_back();
601 
602  }
603 
604  result.emplace_back(allCells[foundQuadruplets[quadId][0]].getInnerHit(),
605  allCells[foundQuadruplets[quadId][1]].getInnerHit(),
606  allCells[foundQuadruplets[quadId][2]].getInnerHit(),
607  allCells[foundQuadruplets[quadId][2]].getOuterHit());
608  hasAlreadyPushedACandidate = true;
609 
610  }
611  }
612  else
613  {
614  result.emplace_back(allCells[foundQuadruplets[quadId][0]].getInnerHit(),
615  allCells[foundQuadruplets[quadId][1]].getInnerHit(),
616  allCells[foundQuadruplets[quadId][2]].getInnerHit(),
617  allCells[foundQuadruplets[quadId][2]].getOuterHit());
618  }
619 
620  }
621 
622 }
unsigned int side(const DetId &id) const
T y() const
Definition: PV3DBase.h:63
std::unique_ptr< SeedComparitor > theComparitor
T inversePt(T curvature, const edm::EventSetup &iSetup)
unsigned int subDetId[21]
bool isNotFinite(T x)
Definition: isFinite.h:10
T curvature(T InversePt, const edm::EventSetup &iSetup)
T sqrt(T t)
Definition: SSEVec.h:18
T z() const
Definition: PV3DBase.h:64
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
Definition: RZLine.h:12
const QuantityDependsPt maxChi2
T rerr(const GlobalPoint &aPoint) const
const T & get() const
Definition: EventSetup.h:56
unsigned int layer(const DetId &id) const
Square< F >::type sqr(const F &f)
Definition: Square.h:13
QuantityDependsPtEval evaluator(const edm::EventSetup &es) const
T x() const
Definition: PV3DBase.h:62
T const * product() const
Definition: ESHandle.h:86
*vegas h *****************************************************used in the default bin number in original ***version of VEGAS is ***a higher bin number might help to derive a more precise ***grade subtle point
Definition: invegas.h:5
void CAHitQuadrupletGenerator::initEvent ( const edm::Event ev,
const edm::EventSetup es 
)

Member Data Documentation

const float CAHitQuadrupletGenerator::caHardPtCut = 0.f
private

Definition at line 148 of file CAHitQuadrupletGenerator.h.

Referenced by hitNtuplets(), and hitQuadruplets().

const bool CAHitQuadrupletGenerator::caOnlyOneLastHitPerLayerFilter = false
private

Definition at line 149 of file CAHitQuadrupletGenerator.h.

Referenced by hitNtuplets(), and hitQuadruplets().

const float CAHitQuadrupletGenerator::caPhiCut = 0.1f
private

Definition at line 147 of file CAHitQuadrupletGenerator.h.

Referenced by hitNtuplets(), and hitQuadruplets().

const float CAHitQuadrupletGenerator::caThetaCut = 0.00125f
private

Definition at line 146 of file CAHitQuadrupletGenerator.h.

Referenced by hitNtuplets(), and hitQuadruplets().

const float CAHitQuadrupletGenerator::extraHitRPhitolerance
private

Definition at line 139 of file CAHitQuadrupletGenerator.h.

Referenced by hitNtuplets(), and hitQuadruplets().

const bool CAHitQuadrupletGenerator::fitFastCircle
private

Definition at line 142 of file CAHitQuadrupletGenerator.h.

Referenced by hitNtuplets(), and hitQuadruplets().

const bool CAHitQuadrupletGenerator::fitFastCircleChi2Cut
private

Definition at line 143 of file CAHitQuadrupletGenerator.h.

Referenced by hitNtuplets(), and hitQuadruplets().

const QuantityDependsPt CAHitQuadrupletGenerator::maxChi2
private

Definition at line 141 of file CAHitQuadrupletGenerator.h.

Referenced by hitNtuplets(), and hitQuadruplets().

unsigned int CAHitQuadrupletGenerator::minLayers = 4
static

Definition at line 36 of file CAHitQuadrupletGenerator.h.

std::unique_ptr<SeedComparitor> CAHitQuadrupletGenerator::theComparitor
private
LayerCacheType CAHitQuadrupletGenerator::theLayerCache
private

Definition at line 71 of file CAHitQuadrupletGenerator.h.

Referenced by hitQuadruplets().

edm::EDGetTokenT<SeedingLayerSetsHits> CAHitQuadrupletGenerator::theSeedingLayerToken
private

Definition at line 69 of file CAHitQuadrupletGenerator.h.

Referenced by CAHitQuadrupletGenerator(), and hitQuadruplets().

const bool CAHitQuadrupletGenerator::useBendingCorrection
private

Definition at line 144 of file CAHitQuadrupletGenerator.h.

Referenced by hitNtuplets(), and hitQuadruplets().