CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
TrajSeedMatcher.cc
Go to the documentation of this file.
2 
7 
11 
15 
17 
20 
22 
24  : cacheIDMagField_(0),
25  paramMagFieldLabel_(pset.getParameter<std::string>("paramMagField")),
26  navSchoolLabel_(pset.getParameter<std::string>("navSchool")),
27  detLayerGeomLabel_(pset.getParameter<std::string>("detLayerGeom")),
28  useRecoVertex_(pset.getParameter<bool>("useRecoVertex")),
29  enableHitSkipping_(pset.getParameter<bool>("enableHitSkipping")),
30  requireExactMatchCount_(pset.getParameter<bool>("requireExactMatchCount")),
31  useParamMagFieldIfDefined_(pset.getParameter<bool>("useParamMagFieldIfDefined")),
32  minNrHits_(pset.getParameter<std::vector<unsigned int> >("minNrHits")),
33  minNrHitsValidLayerBins_(pset.getParameter<std::vector<int> >("minNrHitsValidLayerBins")) {
34  const auto cutsPSets = pset.getParameter<std::vector<edm::ParameterSet> >("matchingCuts");
35  for (const auto& cutPSet : cutsPSets) {
36  int version = cutPSet.getParameter<int>("version");
37  switch (version) {
38  case 1:
39  matchingCuts_.emplace_back(std::make_unique<MatchingCutsV1>(cutPSet));
40  break;
41  case 2:
42  matchingCuts_.emplace_back(std::make_unique<MatchingCutsV2>(cutPSet));
43  break;
44  default:
45  throw cms::Exception("InvalidConfig") << " Error TrajSeedMatcher::TrajSeedMatcher pixel match cuts version "
46  << version << " not recognised" << std::endl;
47  }
48  }
49 
50  if (minNrHitsValidLayerBins_.size() + 1 != minNrHits_.size()) {
51  throw cms::Exception("InvalidConfig")
52  << " TrajSeedMatcher::TrajSeedMatcher minNrHitsValidLayerBins should be 1 less than minNrHits when its "
53  << minNrHitsValidLayerBins_.size() << " vs " << minNrHits_.size();
54  }
55 }
56 
59  desc.add<bool>("useRecoVertex", false);
60  desc.add<bool>("enableHitSkipping", false);
61  desc.add<bool>("requireExactMatchCount", true);
62  desc.add<bool>("useParamMagFieldIfDefined", true);
63  desc.add<std::string>("paramMagField", "ParabolicMf");
64  desc.add<std::string>("navSchool", "SimpleNavigationSchool");
65  desc.add<std::string>("detLayerGeom", "hltESPGlobalDetLayerGeometry");
66  desc.add<std::vector<int> >("minNrHitsValidLayerBins", {4});
67  desc.add<std::vector<unsigned int> >("minNrHits", {2, 3});
68 
70  auto cutDescCases = 1 >> (edm::ParameterDescription<double>("dPhiMax", 0.04, true) and
71  edm::ParameterDescription<double>("dRZMax", 0.09, true) and
72  edm::ParameterDescription<double>("dRZMaxLowEtThres", 20., true) and
73  edm::ParameterDescription<std::vector<double> >(
74  "dRZMaxLowEtEtaBins", std::vector<double>{1., 1.5}, true) and
75  edm::ParameterDescription<std::vector<double> >(
76  "dRZMaxLowEt", std::vector<double>{0.09, 0.15, 0.09}, true)) or
77  2 >> (edm::ParameterDescription<std::vector<double> >("dPhiMaxHighEt", {0.003}, true) and
78  edm::ParameterDescription<std::vector<double> >("dPhiMaxHighEtThres", {0.0}, true) and
79  edm::ParameterDescription<std::vector<double> >("dPhiMaxLowEtGrad", {0.0}, true) and
80  edm::ParameterDescription<std::vector<double> >("dRZMaxHighEt", {0.005}, true) and
81  edm::ParameterDescription<std::vector<double> >("dRZMaxHighEtThres", {30}, true) and
82  edm::ParameterDescription<std::vector<double> >("dRZMaxLowEtGrad", {-0.002}, true) and
83  edm::ParameterDescription<std::vector<double> >("etaBins", {}, true));
84  cutsDesc.ifValue(edm::ParameterDescription<int>("version", 1, true), std::move(cutDescCases));
85 
87  defaults.addParameter<double>("dPhiMax", 0.04);
88  defaults.addParameter<double>("dRZMax", 0.09);
89  defaults.addParameter<double>("dRZMaxLowEtThres", 0.09);
90  defaults.addParameter<std::vector<double> >("dRZMaxLowEtEtaBins", std::vector<double>{1., 1.5});
91  defaults.addParameter<std::vector<double> >("dRZMaxLowEt", std::vector<double>{0.09, 0.09, 0.09});
92  defaults.addParameter<int>("version", 1);
93  desc.addVPSet("matchingCuts", cutsDesc, std::vector<edm::ParameterSet>{defaults, defaults, defaults});
94  return desc;
95 }
96 
99  iSetup.get<IdealMagneticFieldRecord>().get(magField_);
102  cacheIDMagField_ = iSetup.get<IdealMagneticFieldRecord>().cacheIdentifier();
103  forwardPropagator_ = std::make_unique<PropagatorWithMaterial>(alongMomentum, kElectronMass_, &*(magField_));
104  backwardPropagator_ = std::make_unique<PropagatorWithMaterial>(oppositeToMomentum, kElectronMass_, &*(magField_));
105  }
108 }
109 
110 std::vector<TrajSeedMatcher::SeedWithInfo> TrajSeedMatcher::compatibleSeeds(const TrajectorySeedCollection& seeds,
111  const GlobalPoint& candPos,
112  const GlobalPoint& vprim,
113  const float energy) {
115  throw cms::Exception("LogicError") << __FUNCTION__
116  << " can not make pixel seeds as event setup has not properly been called";
117  }
118 
119  clearCache();
120 
121  std::vector<SeedWithInfo> matchedSeeds;
122 
123  //these are super expensive functions
124  TrajectoryStateOnSurface scTrajStateOnSurfNeg = makeTrajStateOnSurface(candPos, vprim, energy, -1);
125  TrajectoryStateOnSurface scTrajStateOnSurfPos = makeTrajStateOnSurface(candPos, vprim, energy, 1);
126 
127  for (const auto& seed : seeds) {
128  std::vector<SCHitMatch> matchedHitsNeg = processSeed(seed, candPos, vprim, energy, scTrajStateOnSurfNeg);
129  std::vector<SCHitMatch> matchedHitsPos = processSeed(seed, candPos, vprim, energy, scTrajStateOnSurfPos);
130 
131  int nrValidLayersPos = 0;
132  int nrValidLayersNeg = 0;
133  if (matchedHitsNeg.size() >= 2) {
134  nrValidLayersNeg = getNrValidLayersAlongTraj(matchedHitsNeg[0], matchedHitsNeg[1], candPos, vprim, energy, -1);
135  }
136  if (matchedHitsPos.size() >= 2) {
137  nrValidLayersPos = getNrValidLayersAlongTraj(matchedHitsPos[0], matchedHitsPos[1], candPos, vprim, energy, +1);
138  }
139 
140  int nrValidLayers = std::max(nrValidLayersNeg, nrValidLayersPos);
141  size_t nrHitsRequired = getNrHitsRequired(nrValidLayers);
142  bool matchCountPasses;
144  // If the input seed collection is not cross-cleaned, an exact match is necessary to
145  // prevent redundant seeds.
146  matchCountPasses = matchedHitsNeg.size() == nrHitsRequired || matchedHitsPos.size() == nrHitsRequired;
147  } else {
148  matchCountPasses = matchedHitsNeg.size() >= nrHitsRequired || matchedHitsPos.size() >= nrHitsRequired;
149  }
150  if (matchCountPasses) {
151  matchedSeeds.push_back({seed, matchedHitsPos, matchedHitsNeg, nrValidLayers});
152  }
153  }
154  return matchedSeeds;
155 }
156 
157 std::vector<TrajSeedMatcher::SCHitMatch> TrajSeedMatcher::processSeed(const TrajectorySeed& seed,
158  const GlobalPoint& candPos,
159  const GlobalPoint& vprim,
160  const float energy,
161  const TrajectoryStateOnSurface& initialTrajState) {
162  //next try passing these variables in once...
163  const float candEta = candPos.eta();
164  const float candEt = energy * std::sin(candPos.theta());
165  const int charge = initialTrajState.charge();
166 
167  std::vector<SCHitMatch> matchedHits;
168  FreeTrajectoryState firstHitFreeTraj;
169  GlobalPoint prevHitPos;
171  for (size_t hitNr = 0; hitNr < matchingCuts_.size() && hitNr < seed.nHits(); hitNr++) {
172  bool matched = false;
173  if (matchedHits.empty()) { // First hit is special
174  SCHitMatch firstHit = matchFirstHit(seed, initialTrajState, vprim, *backwardPropagator_);
175  firstHit.setExtra(candEt, candEta, candPos.phi(), charge, 1);
176  if (passesMatchSel(firstHit, 0)) {
177  matched = true;
178  matchedHits.push_back(firstHit);
179  //now we can figure out the z vertex
180  double zVertex = useRecoVertex_ ? vprim.z() : getZVtxFromExtrapolation(vprim, firstHit.hitPos(), candPos);
181  vertex = GlobalPoint(vprim.x(), vprim.y(), zVertex);
182  firstHitFreeTraj =
184  prevHitPos = firstHit.hitPos();
185  }
186  } else { // All subsequent hits are handled the same
187  SCHitMatch hit = match2ndToNthHit(seed, firstHitFreeTraj, hitNr, prevHitPos, vertex, *forwardPropagator_);
188  hit.setExtra(candEt, candEta, candPos.phi(), charge, 1);
189  if (passesMatchSel(hit, matchedHits.size())) {
190  matched = true;
191  matchedHits.push_back(hit);
192  prevHitPos = hit.hitPos();
193  }
194  }
195  if (!matched and !enableHitSkipping_)
196  break;
197  }
198  return matchedHits;
199 }
200 
201 // compute the z vertex from the candidate position and the found pixel hit
203  const GlobalPoint& hitPos,
204  const GlobalPoint& candPos) {
205  auto sq = [](float x) { return x * x; };
206  auto calRDiff = [sq](const GlobalPoint& p1, const GlobalPoint& p2) {
207  return std::sqrt(sq(p2.x() - p1.x()) + sq(p2.y() - p1.y()));
208  };
209  const double r1Diff = calRDiff(primeVtxPos, hitPos);
210  const double r2Diff = calRDiff(hitPos, candPos);
211  return hitPos.z() - r1Diff * (candPos.z() - hitPos.z()) / r2Diff;
212 }
213 
214 bool TrajSeedMatcher::passTrajPreSel(const GlobalPoint& hitPos, const GlobalPoint& candPos) const {
215  float dt = hitPos.x() * candPos.x() + hitPos.y() * candPos.y();
216  if (dt < 0)
217  return false;
218  if (dt < kPhiCut_ * (candPos.perp() * hitPos.perp()))
219  return false;
220  return true;
221 }
222 
224  const TrajectoryStateOnSurface& initialState,
226  auto& trajStateFromVtxCache =
228 
229  auto key = hit.det()->gdetIndex();
230  auto res = trajStateFromVtxCache.find(key);
231  if (res != trajStateFromVtxCache.end())
232  return res->second;
233  else { //doesnt exist, need to make it
234  //FIXME: check for efficiency
235  auto val = trajStateFromVtxCache.emplace(key, propagator.propagate(initialState, hit.det()->surface()));
236  return val.first->second;
237  }
238 }
239 
241  const FreeTrajectoryState& initialState,
242  const GlobalPoint& point,
244  auto& trajStateFromPointCache =
246 
247  auto key = std::make_pair(hit.det()->gdetIndex(), point);
248  auto res = trajStateFromPointCache.find(key);
249  if (res != trajStateFromPointCache.end())
250  return res->second;
251  else { //doesnt exist, need to make it
252  //FIXME: check for efficiency
253  auto val = trajStateFromPointCache.emplace(key, propagator.propagate(initialState, hit.det()->surface()));
254  return val.first->second;
255  }
256 }
257 
259  const GlobalPoint& vtx,
260  const float energy,
261  const int charge) const {
262  FreeTrajectoryState freeTS = FTSFromVertexToPointFactory::get(getMagField(pos), pos, vtx, energy, charge);
264  return TrajectoryStateOnSurface(freeTS, *bpb(freeTS.position(), freeTS.momentum()));
265 }
266 
268  const TrajectoryStateOnSurface& initialState,
269  const GlobalPoint& vtxPos,
271  const TrajectorySeed::range& hits = seed.recHits();
272  auto hitIt = hits.first;
273 
274  if (hitIt->isValid()) {
275  const TrajectoryStateOnSurface& trajStateFromVtx = getTrajStateFromVtx(*hitIt, initialState, propagator);
276  if (trajStateFromVtx.isValid())
277  return SCHitMatch(vtxPos, trajStateFromVtx, *hitIt);
278  }
279  return SCHitMatch();
280 }
281 
283  const FreeTrajectoryState& initialState,
284  const size_t hitNr,
285  const GlobalPoint& prevHitPos,
286  const GlobalPoint& vtxPos,
288  const TrajectorySeed::range& hits = seed.recHits();
289  auto hitIt = hits.first + hitNr;
290 
291  if (hitIt->isValid()) {
292  const TrajectoryStateOnSurface& trajState = getTrajStateFromPoint(*hitIt, initialState, prevHitPos, propagator);
293 
294  if (trajState.isValid()) {
295  return SCHitMatch(vtxPos, trajState, *hitIt);
296  }
297  }
298  return SCHitMatch();
299 }
300 
306 }
307 
308 bool TrajSeedMatcher::passesMatchSel(const TrajSeedMatcher::SCHitMatch& hit, const size_t hitNr) const {
309  if (hitNr < matchingCuts_.size()) {
310  return (*matchingCuts_[hitNr])(hit);
311  } else {
312  throw cms::Exception("LogicError") << " Error, attempting to apply selection to hit " << hitNr
313  << " but only cuts for " << matchingCuts_.size() << " defined";
314  }
315 }
316 
318  const SCHitMatch& hit2,
319  const GlobalPoint& candPos,
320  const GlobalPoint& vprim,
321  const float energy,
322  const int charge) {
323  double zVertex = useRecoVertex_ ? vprim.z() : getZVtxFromExtrapolation(vprim, hit1.hitPos(), candPos);
324  GlobalPoint vertex(vprim.x(), vprim.y(), zVertex);
325 
326  FreeTrajectoryState firstHitFreeTraj =
328  const TrajectoryStateOnSurface& secondHitTraj =
329  getTrajStateFromPoint(*hit2.hit(), firstHitFreeTraj, hit1.hitPos(), *forwardPropagator_);
330  return getNrValidLayersAlongTraj(hit2.hit()->geographicalId(), secondHitTraj);
331 }
332 
333 int TrajSeedMatcher::getNrValidLayersAlongTraj(const DetId& hitId, const TrajectoryStateOnSurface& hitTrajState) const {
334  const DetLayer* detLayer = detLayerGeom_->idToLayer(hitId);
335  if (detLayer == nullptr)
336  return 0;
337 
338  const FreeTrajectoryState& hitFreeState = *hitTrajState.freeState();
339  const std::vector<const DetLayer*> inLayers =
340  navSchool_->compatibleLayers(*detLayer, hitFreeState, oppositeToMomentum);
341  const std::vector<const DetLayer*> outLayers = navSchool_->compatibleLayers(*detLayer, hitFreeState, alongMomentum);
342 
343  int nrValidLayers = 1; //because our current hit is also valid and wont be included in the count otherwise
344  int nrPixInLayers = 0;
345  int nrPixOutLayers = 0;
346  for (auto layer : inLayers) {
347  if (GeomDetEnumerators::isTrackerPixel(layer->subDetector())) {
348  nrPixInLayers++;
349  if (layerHasValidHits(*layer, hitTrajState, *backwardPropagator_))
350  nrValidLayers++;
351  }
352  }
353  for (auto layer : outLayers) {
354  if (GeomDetEnumerators::isTrackerPixel(layer->subDetector())) {
355  nrPixOutLayers++;
356  if (layerHasValidHits(*layer, hitTrajState, *forwardPropagator_))
357  nrValidLayers++;
358  }
359  }
360  return nrValidLayers;
361 }
362 
364  const TrajectoryStateOnSurface& hitSurState,
365  const Propagator& propToLayerFromState) const {
366  //FIXME: do not hardcode with werid magic numbers stolen from ancient tracking code
367  //its taken from https://cmssdt.cern.ch/dxr/CMSSW/source/RecoTracker/TrackProducer/interface/TrackProducerBase.icc#165
368  //which inspires this code
369  Chi2MeasurementEstimator estimator(30., -3.0, 0.5, 2.0, 0.5, 1.e12); // same as defauts....
370 
371  const std::vector<GeometricSearchDet::DetWithState>& detWithState =
372  layer.compatibleDets(hitSurState, propToLayerFromState, estimator);
373  if (detWithState.empty())
374  return false;
375  else {
376  DetId id = detWithState.front().first->geographicalId();
378  if (measDet.isActive())
379  return true;
380  else
381  return false;
382  }
383 }
384 
385 size_t TrajSeedMatcher::getNrHitsRequired(const int nrValidLayers) const {
386  for (size_t binNr = 0; binNr < minNrHitsValidLayerBins_.size(); binNr++) {
387  if (nrValidLayers < minNrHitsValidLayerBins_[binNr])
388  return minNrHits_[binNr];
389  }
390  return minNrHits_.back();
391 }
392 
394  const TrajectoryStateOnSurface& trajState,
395  const TrackingRecHit& hit)
396  : detId_(hit.geographicalId()),
397  hitPos_(hit.globalPosition()),
398  hit_(&hit),
399  et_(0),
400  eta_(0),
401  phi_(0),
402  charge_(0),
403  nrClus_(0) {
404  EleRelPointPair pointPair(hitPos_, trajState.globalParameters().position(), vtxPos);
405  dRZ_ = detId_.subdetId() == PixelSubdetector::PixelBarrel ? pointPair.dZ() : pointPair.dPerp();
406  dPhi_ = pointPair.dPhi();
407 }
408 
410  const std::vector<SCHitMatch>& posCharge,
411  const std::vector<SCHitMatch>& negCharge,
412  int nrValidLayers)
413  : seed_(seed), nrValidLayers_(nrValidLayers) {
414  size_t nrHitsMax = std::max(posCharge.size(), negCharge.size());
415  for (size_t hitNr = 0; hitNr < nrHitsMax; hitNr++) {
416  DetId detIdPos = hitNr < posCharge.size() ? posCharge[hitNr].detId() : DetId(0);
417  float dRZPos = hitNr < posCharge.size() ? posCharge[hitNr].dRZ() : std::numeric_limits<float>::max();
418  float dPhiPos = hitNr < posCharge.size() ? posCharge[hitNr].dPhi() : std::numeric_limits<float>::max();
419 
420  DetId detIdNeg = hitNr < negCharge.size() ? negCharge[hitNr].detId() : DetId(0);
421  float dRZNeg = hitNr < negCharge.size() ? negCharge[hitNr].dRZ() : std::numeric_limits<float>::max();
422  float dPhiNeg = hitNr < negCharge.size() ? negCharge[hitNr].dPhi() : std::numeric_limits<float>::max();
423 
424  if (detIdPos != detIdNeg && (detIdPos.rawId() != 0 && detIdNeg.rawId() != 0)) {
425  cms::Exception("LogicError")
426  << " error in " << __FILE__ << ", " << __LINE__
427  << " hits to be combined have different detIDs, this should not be possible and nothing good will come of it";
428  }
429  DetId detId = detIdPos.rawId() != 0 ? detIdPos : detIdNeg;
430  matchInfo_.push_back(MatchInfo(detId, dRZPos, dRZNeg, dPhiPos, dPhiNeg));
431  }
432 }
433 
435  : dPhiMax_(pset.getParameter<double>("dPhiMax")),
436  dRZMax_(pset.getParameter<double>("dRZMax")),
437  dRZMaxLowEtThres_(pset.getParameter<double>("dRZMaxLowEtThres")),
438  dRZMaxLowEtEtaBins_(pset.getParameter<std::vector<double> >("dRZMaxLowEtEtaBins")),
439  dRZMaxLowEt_(pset.getParameter<std::vector<double> >("dRZMaxLowEt")) {
440  if (dRZMaxLowEtEtaBins_.size() + 1 != dRZMaxLowEt_.size()) {
441  throw cms::Exception("InvalidConfig") << " dRZMaxLowEtEtaBins should be 1 less than dRZMaxLowEt when its "
442  << dRZMaxLowEtEtaBins_.size() << " vs " << dRZMaxLowEt_.size();
443  }
444 }
445 
447  if (dPhiMax_ >= 0 && std::abs(scHitMatch.dPhi()) > dPhiMax_)
448  return false;
449 
450  const float dRZMax = getDRZCutValue(scHitMatch.et(), scHitMatch.eta());
451  if (dRZMax_ >= 0 && std::abs(scHitMatch.dRZ()) > dRZMax)
452  return false;
453 
454  return true;
455 }
456 
457 float TrajSeedMatcher::MatchingCutsV1::getDRZCutValue(const float scEt, const float scEta) const {
458  if (scEt >= dRZMaxLowEtThres_)
459  return dRZMax_;
460  else {
461  const float absEta = std::abs(scEta);
462  for (size_t etaNr = 0; etaNr < dRZMaxLowEtEtaBins_.size(); etaNr++) {
463  if (absEta < dRZMaxLowEtEtaBins_[etaNr])
464  return dRZMaxLowEt_[etaNr];
465  }
466  return dRZMaxLowEt_.back();
467  }
468 }
469 
471  : dPhiHighEt_(pset.getParameter<std::vector<double> >("dPhiMaxHighEt")),
472  dPhiHighEtThres_(pset.getParameter<std::vector<double> >("dPhiMaxHighEtThres")),
473  dPhiLowEtGrad_(pset.getParameter<std::vector<double> >("dPhiMaxLowEtGrad")),
474  dRZHighEt_(pset.getParameter<std::vector<double> >("dRZMaxHighEt")),
475  dRZHighEtThres_(pset.getParameter<std::vector<double> >("dRZMaxHighEtThres")),
476  dRZLowEtGrad_(pset.getParameter<std::vector<double> >("dRZMaxLowEtGrad")),
477  etaBins_(pset.getParameter<std::vector<double> >("etaBins")) {
478  auto binSizeCheck = [](size_t sizeEtaBins, const std::vector<double>& vec, const std::string& name) {
479  if (vec.size() != sizeEtaBins + 1) {
480  throw cms::Exception("InvalidConfig")
481  << " when constructing TrajSeedMatcher::MatchingCutsV2 " << name << " has " << vec.size()
482  << " bins, it should be equal to #bins of etaBins+1" << sizeEtaBins + 1;
483  }
484  };
485  binSizeCheck(etaBins_.size(), dPhiHighEt_, "dPhiMaxHighEt");
486  binSizeCheck(etaBins_.size(), dPhiHighEtThres_, "dPhiMaxHighEtThres");
487  binSizeCheck(etaBins_.size(), dPhiLowEtGrad_, "dPhiMaxLowEtGrad");
488  binSizeCheck(etaBins_.size(), dRZHighEt_, "dRZMaxHighEt");
489  binSizeCheck(etaBins_.size(), dRZHighEtThres_, "dRZMaxHighEtThres");
490  binSizeCheck(etaBins_.size(), dRZLowEtGrad_, "dRZMaxLowEtGrad");
491 }
492 
494  size_t binNr = getBinNr(scHitMatch.eta());
495  float dPhiMax = getCutValue(scHitMatch.et(), dPhiHighEt_[binNr], dPhiHighEtThres_[binNr], dPhiLowEtGrad_[binNr]);
496  if (dPhiMax >= 0 && std::abs(scHitMatch.dPhi()) > dPhiMax)
497  return false;
498  float dRZMax = getCutValue(scHitMatch.et(), dRZHighEt_[binNr], dRZHighEtThres_[binNr], dRZLowEtGrad_[binNr]);
499  if (dRZMax >= 0 && std::abs(scHitMatch.dRZ()) > dRZMax)
500  return false;
501 
502  return true;
503 }
504 
505 //eta bins is exactly 1 smaller than the vectors which will be accessed by this bin nr
507  const float absEta = std::abs(eta);
508  for (size_t etaNr = 0; etaNr < etaBins_.size(); etaNr++) {
509  if (absEta < etaBins_[etaNr])
510  return etaNr;
511  }
512  return etaBins_.size();
513 }
MatchingCutsV1(const edm::ParameterSet &pset)
const GlobalPoint & hitPos() const
T getParameter(std::string const &) const
unsigned long long cacheIdentifier() const
MeasurementDetWithData idToDet(const DetId &id) const
Previous MeasurementDetSystem interface.
void doEventSetup(const edm::EventSetup &iSetup)
float dt
Definition: AMPTWrapper.h:136
bool useParamMagFieldIfDefined_
IntGlobalPointPairUnorderedMap< TrajectoryStateOnSurface > trajStateFromPointNegChargeCache_
std::pair< const_iterator, const_iterator > range
std::vector< TrajSeedMatcher::SeedWithInfo > compatibleSeeds(const TrajectorySeedCollection &seeds, const GlobalPoint &candPos, const GlobalPoint &vprim, const float energy)
static FreeTrajectoryState get(MagneticField const &magField, GlobalPoint const &xmeas, GlobalPoint const &xvert, float momentum, TrackCharge charge)
ParameterDescriptionBase * addVPSet(U const &iLabel, ParameterSetDescription const &validator, std::vector< ParameterSet > const &defaults)
T perp() const
Definition: PV3DBase.h:69
MatchingCutsV2(const edm::ParameterSet &pset)
DetId detId(size_t hitNr) const
std::vector< double > dRZLowEtGrad_
const TrajectoryStateOnSurface & getTrajStateFromVtx(const TrackingRecHit &hit, const TrajectoryStateOnSurface &initialState, const PropagatorWithMaterial &propagator)
std::vector< double > dPhiLowEtGrad_
int gdetIndex() const
Definition: GeomDet.h:87
const TrackingRecHit * hit() const
const MagneticField & getMagField(const GlobalPoint &point) const
TrajectoryStateOnSurface makeTrajStateOnSurface(const GlobalPoint &pos, const GlobalPoint &vtx, const float energy, const int charge) const
edm::ESHandle< MagneticField > magFieldParam_
IntGlobalPointPairUnorderedMap< TrajectoryStateOnSurface > trajStateFromPointPosChargeCache_
bool passesMatchSel(const SCHitMatch &hit, const size_t hitNr) const
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
edm::ESHandle< DetLayerGeometry > detLayerGeom_
Geom::Phi< T > phi() const
Definition: PV3DBase.h:66
float dRZPos(size_t hitNr) const
bool passTrajPreSel(const GlobalPoint &hitPos, const GlobalPoint &candPos) const
Global3DPoint GlobalPoint
Definition: GlobalPoint.h:10
constexpr uint32_t rawId() const
get the raw id
Definition: DetId.h:57
std::unique_ptr< PropagatorWithMaterial > forwardPropagator_
T y() const
Definition: PV3DBase.h:60
const std::vector< int > minNrHitsValidLayerBins_
TrackCharge charge() const
const std::vector< double > dRZMaxLowEt_
static float getZVtxFromExtrapolation(const GlobalPoint &primeVtxPos, const GlobalPoint &hitPos, const GlobalPoint &candPos)
const Plane & surface() const
The nominal surface of the GeomDet.
Definition: GeomDet.h:37
virtual std::vector< DetWithState > compatibleDets(const TrajectoryStateOnSurface &startingState, const Propagator &prop, const MeasurementEstimator &est) const
Definition: Electron.h:6
const TrajectoryStateOnSurface & getTrajStateFromPoint(const TrackingRecHit &hit, const FreeTrajectoryState &initialState, const GlobalPoint &point, const PropagatorWithMaterial &propagator)
Geom::Theta< T > theta() const
Definition: PV3DBase.h:72
edm::ESHandle< NavigationSchool > navSchool_
std::vector< double > etaBins_
unsigned long long cacheIDMagField_
static edm::ParameterSetDescription makePSetDescription()
std::vector< TrajectorySeed > TrajectorySeedCollection
edm::ESHandle< MagneticField > magField_
const GeomDet * det() const
T sqrt(T t)
Definition: SSEVec.h:19
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< void, edm::EventID const &, edm::Timestamp const & > We also list in braces which AR_WATCH_USING_METHOD_ is used for those or
Definition: Activities.doc:12
std::vector< MatchInfo > matchInfo_
FreeTrajectoryState const * freeState(bool withErrors=true) const
void addParameter(std::string const &name, T const &value)
Definition: ParameterSet.h:124
T z() const
Definition: PV3DBase.h:61
std::vector< const DetLayer * > compatibleLayers(const DetLayer &detLayer, Args &&...args) const
Returns all layers compatible.
constexpr int subdetId() const
get the contents of the subdetector field (not cast into any detector&#39;s numbering enum) ...
Definition: DetId.h:48
std::vector< std::unique_ptr< MatchingCuts > > matchingCuts_
std::string paramMagFieldLabel_
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
std::unordered_map< int, TrajectoryStateOnSurface > trajStateFromVtxPosChargeCache_
void setExtra(float et, float eta, float phi, int charge, int nrClus)
float dPhiPos(size_t hitNr) const
bool operator()(const SCHitMatch &scHitMatch) const override
bool operator()(const SCHitMatch &scHitMatch) const override
ParameterDescriptionBase * add(U const &iLabel, T const &value)
double p2[4]
Definition: TauolaWrapper.h:90
std::vector< double > dPhiHighEt_
size_t getNrHitsRequired(const int nrValidLayers) const
GlobalVector momentum() const
size_t getBinNr(float eta) const
std::string navSchoolLabel_
float getDRZCutValue(const float scEt, const float scEta) const
Definition: DetId.h:17
GlobalPoint position() const
const std::vector< double > dRZMaxLowEtEtaBins_
const GlobalTrajectoryParameters & globalParameters() const
edm::Handle< MeasurementTrackerEvent > measTkEvt_
static float kElectronMass_
range recHits() const
std::vector< double > dPhiHighEtThres_
TrajectoryStateOnSurface propagate(STA const &state, SUR const &surface) const
Definition: Propagator.h:50
TrajSeedMatcher::SCHitMatch match2ndToNthHit(const TrajectorySeed &seed, const FreeTrajectoryState &trajState, const size_t hitNr, const GlobalPoint &prevHitPos, const GlobalPoint &vtxPos, const PropagatorWithMaterial &propagator)
T eta() const
Definition: PV3DBase.h:73
std::vector< SCHitMatch > processSeed(const TrajectorySeed &seed, const GlobalPoint &candPos, const GlobalPoint &vprim, const float energy, const TrajectoryStateOnSurface &initialTrajState)
unsigned int nHits() const
SeedWithInfo(const TrajectorySeed &seed, const std::vector< SCHitMatch > &posCharge, const std::vector< SCHitMatch > &negCharge, int nrValidLayers)
double p1[4]
Definition: TauolaWrapper.h:89
T get() const
Definition: EventSetup.h:73
TrajSeedMatcher(const edm::ParameterSet &pset)
std::string detLayerGeomLabel_
float dRZNeg(size_t hitNr) const
DetId geographicalId() const
std::unordered_map< int, TrajectoryStateOnSurface > trajStateFromVtxNegChargeCache_
bool isValid() const
Definition: ESHandle.h:44
virtual const DetLayer * idToLayer(const DetId &detId) const
int getNrValidLayersAlongTraj(const SCHitMatch &hit1, const SCHitMatch &hit2, const GlobalPoint &candPos, const GlobalPoint &vprim, const float energy, const int charge)
T x() const
Definition: PV3DBase.h:59
std::unique_ptr< PropagatorWithMaterial > backwardPropagator_
bool layerHasValidHits(const DetLayer &layer, const TrajectoryStateOnSurface &hitSurState, const Propagator &propToLayerFromState) const
TrajSeedMatcher::SCHitMatch matchFirstHit(const TrajectorySeed &seed, const TrajectoryStateOnSurface &trajState, const GlobalPoint &vtxPos, const PropagatorWithMaterial &propagator)
def move(src, dest)
Definition: eostools.py:511
static float kPhiCut_
bool isTrackerPixel(GeomDetEnumerators::SubDetector m)
#define constexpr
*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
std::vector< double > dRZHighEtThres_
std::vector< double > dRZHighEt_
float dPhiNeg(size_t hitNr) const
float getCutValue(float et, float highEt, float highEtThres, float lowEtGrad) const
const std::vector< unsigned int > minNrHits_