CMS 3D CMS Logo

L1TrackSelectionProducer.cc
Go to the documentation of this file.
1 // -*- C++ -*-
2 //
3 // Package: L1Trigger/L1TTrackMatch
4 // Class: L1TrackSelectionProducer
5 //
18 //
19 // Original Author: Alexx Perloff
20 // Created: Thu, 16 Dec 2021 19:02:50 GMT
21 //
22 // Updates: Claire Savard (claire.savard@colorado.edu), Nov. 2023
23 //
24 
25 // system include files
26 #include <algorithm>
27 #include <memory>
28 #include <string>
29 #include <vector>
30 
31 // Xilinx HLS includes
32 #include <ap_fixed.h>
33 #include <ap_int.h>
34 
35 // user include files
61 
62 //
63 // class declaration
64 //
65 
67 public:
69  ~L1TrackSelectionProducer() override;
70 
71  static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
72 
73 private:
74  // ----------constants, enums and typedefs ---------
75  // Relevant constants for the converted track word
77  kPtSize = TTTrack_TrackWord::TrackBitWidths::kRinvSize - 1, // Width of pt
78  kPtMagSize = 9, // Width of pt magnitude (unsigned)
79  kEtaSize = TTTrack_TrackWord::TrackBitWidths::kTanlSize, // Width of eta
80  kEtaMagSize = 3, // Width of eta magnitude (signed)
81  };
82 
84  typedef std::vector<L1Track> TTTrackCollection;
88  typedef std::unique_ptr<TTTrackRefCollection> TTTrackRefCollectionUPtr;
89 
90  // ----------member functions ----------------------
91  void printDebugInfo(const TTTrackCollectionHandle& l1TracksHandle,
92  const TTTrackRefCollectionUPtr& vTTTrackOutput,
93  const TTTrackRefCollectionUPtr& vTTTrackEmulationOutput) const;
94  void printTrackInfo(edm::LogInfo& log, const L1Track& track, bool printEmulation = false) const;
95  void produce(edm::StreamID, edm::Event&, const edm::EventSetup&) const override;
96 
97  // ----------selectors -----------------------------
98  // Based on recommendations from https://twiki.cern.ch/twiki/bin/view/CMSPublic/SWGuideGenericSelectors
101  TTTrackPtMinSelector(const edm::ParameterSet& cfg) : ptMin_(cfg.template getParameter<double>("ptMin")) {}
102  bool operator()(const L1Track& t) const { return t.momentum().perp() >= ptMin_; }
103 
104  private:
105  double ptMin_;
106  };
109  TTTrackWordPtMinSelector(const edm::ParameterSet& cfg) : ptMin_(cfg.template getParameter<double>("ptMin")) {}
110  bool operator()(const L1Track& t) const {
111  ap_uint<TrackBitWidths::kPtSize> ptEmulationBits = t.getTrackWord()(
112  TTTrack_TrackWord::TrackBitLocations::kRinvMSB - 1, TTTrack_TrackWord::TrackBitLocations::kRinvLSB);
113  ap_ufixed<TrackBitWidths::kPtSize, TrackBitWidths::kPtMagSize> ptEmulation;
114  ptEmulation.V = ptEmulationBits.range();
115  return ptEmulation.to_double() >= ptMin_;
116  }
117 
118  private:
119  double ptMin_;
120  };
124  : absEtaMax_(cfg.template getParameter<double>("absEtaMax")) {}
125  bool operator()(const L1Track& t) const { return std::abs(t.momentum().eta()) <= absEtaMax_; }
126 
127  private:
128  double absEtaMax_;
129  };
133  : absEtaMax_(cfg.template getParameter<double>("absEtaMax")) {}
134  bool operator()(const L1Track& t) const {
135  TTTrack_TrackWord::tanl_t etaEmulationBits = t.getTanlWord();
136  ap_fixed<TrackBitWidths::kEtaSize, TrackBitWidths::kEtaMagSize> etaEmulation;
137  etaEmulation.V = etaEmulationBits.range();
138  return std::abs(etaEmulation.to_double()) <= absEtaMax_;
139  }
140 
141  private:
142  double absEtaMax_;
143  };
146  TTTrackAbsZ0MaxSelector(const edm::ParameterSet& cfg) : absZ0Max_(cfg.template getParameter<double>("absZ0Max")) {}
147  bool operator()(const L1Track& t) const { return std::abs(t.z0()) <= absZ0Max_; }
148 
149  private:
150  double absZ0Max_;
151  };
155  : absZ0Max_(cfg.template getParameter<double>("absZ0Max")) {}
156  bool operator()(const L1Track& t) const {
157  double floatZ0 = t.undigitizeSignedValue(
158  t.getZ0Bits(), TTTrack_TrackWord::TrackBitWidths::kZ0Size, TTTrack_TrackWord::stepZ0, 0.0);
159  return std::abs(floatZ0) <= absZ0Max_;
160  }
161 
162  private:
163  double absZ0Max_;
164  };
168  : nStubsMin_(cfg.template getParameter<double>("nStubsMin")) {}
169  bool operator()(const L1Track& t) const { return t.getStubRefs().size() >= nStubsMin_; }
170 
171  private:
172  double nStubsMin_;
173  };
177  : nStubsMin_(cfg.template getParameter<double>("nStubsMin")) {}
178  bool operator()(const L1Track& t) const { return t.getNStubs() >= nStubsMin_; }
179 
180  private:
181  double nStubsMin_;
182  };
185  : nPSStubsMin_(nStubsMin), tTopo_(tTopo) {}
187  : nPSStubsMin_(cfg.template getParameter<double>("nPSStubsMin")), tTopo_(tTopo) {}
188  bool operator()(const L1Track& t) const {
189  int nPSStubs = 0;
190  for (const auto& stub : t.getStubRefs()) {
191  DetId detId(stub->getDetId());
192  if (detId.det() == DetId::Detector::Tracker) {
193  if ((detId.subdetId() == StripSubdetector::TOB && tTopo_.tobLayer(detId) <= 3) ||
194  (detId.subdetId() == StripSubdetector::TID && tTopo_.tidRing(detId) <= 9))
195  nPSStubs++;
196  }
197  }
198  return nPSStubs >= nPSStubsMin_;
199  }
200 
201  private:
202  double nPSStubsMin_;
204  };
208  : promptMVAMin_(cfg.template getParameter<double>("promptMVAMin")) {}
209  bool operator()(const L1Track& t) const { return t.trkMVA1() >= promptMVAMin_; }
210 
211  private:
213  };
217  : promptMVAMin_(cfg.template getParameter<double>("promptMVAMin")) {}
218  bool operator()(const L1Track& t) const {
219  return t.trkMVA1() >= promptMVAMin_;
220  } //change when mva bins in word are set
221 
222  private:
224  };
226  TTTrackBendChi2MaxSelector(double bendChi2Max) : bendChi2Max_(bendChi2Max) {}
228  : bendChi2Max_(cfg.template getParameter<double>("reducedBendChi2Max")) {}
229  bool operator()(const L1Track& t) const { return t.stubPtConsistency() < bendChi2Max_; }
230 
231  private:
232  double bendChi2Max_;
233  };
235  TTTrackWordBendChi2MaxSelector(double bendChi2Max) : bendChi2Max_(bendChi2Max) {}
237  : bendChi2Max_(cfg.template getParameter<double>("reducedBendChi2Max")) {}
238  bool operator()(const L1Track& t) const { return t.getBendChi2() < bendChi2Max_; }
239 
240  private:
241  double bendChi2Max_;
242  };
246  : reducedChi2RZMax_(cfg.template getParameter<double>("reducedChi2RZMax")) {}
247  bool operator()(const L1Track& t) const { return t.chi2ZRed() < reducedChi2RZMax_; }
248 
249  private:
251  };
255  : reducedChi2RZMax_(cfg.template getParameter<double>("reducedChi2RZMax")) {}
256  bool operator()(const L1Track& t) const { return t.getChi2RZ() < reducedChi2RZMax_; }
257 
258  private:
260  };
264  : reducedChi2RPhiMax_(cfg.template getParameter<double>("reducedChi2RPhiMax")) {}
265  bool operator()(const L1Track& t) const { return t.chi2XYRed() < reducedChi2RPhiMax_; }
266 
267  private:
269  };
273  : reducedChi2RPhiMax_(cfg.template getParameter<double>("reducedChi2RPhiMax")) {}
274  bool operator()(const L1Track& t) const { return t.getChi2RPhi() < reducedChi2RPhiMax_; }
275 
276  private:
278  };
283  : reducedChi2RZMaxNstub4_(cfg.template getParameter<double>("reducedChi2RZMaxNstub4")),
284  reducedChi2RZMaxNstub5_(cfg.template getParameter<double>("reducedChi2RZMaxNstub5")) {}
285  bool operator()(const L1Track& t) const {
286  return (((t.chi2ZRed() < reducedChi2RZMaxNstub4_) && (t.getStubRefs().size() == 4)) ||
287  ((t.chi2ZRed() < reducedChi2RZMaxNstub5_) && (t.getStubRefs().size() > 4)));
288  }
289 
290  private:
293  };
298  : reducedChi2RZMaxNstub4_(cfg.template getParameter<double>("reducedChi2RZMaxNstub4")),
299  reducedChi2RZMaxNstub5_(cfg.template getParameter<double>("reducedChi2RZMaxNstub5")) {}
300  bool operator()(const L1Track& t) const {
301  return (((t.getChi2RZ() < reducedChi2RZMaxNstub4_) && (t.getNStubs() == 4)) ||
302  ((t.getChi2RZ() < reducedChi2RZMaxNstub5_) && (t.getNStubs() > 4)));
303  }
304 
305  private:
308  };
313  : reducedChi2RPhiMaxNstub4_(cfg.template getParameter<double>("reducedChi2RPhiMaxNstub4")),
314  reducedChi2RPhiMaxNstub5_(cfg.template getParameter<double>("reducedChi2RPhiMaxNstub5")) {}
315  bool operator()(const L1Track& t) const {
316  return (((t.chi2XYRed() < reducedChi2RPhiMaxNstub4_) && (t.getStubRefs().size() == 4)) ||
317  ((t.chi2XYRed() < reducedChi2RPhiMaxNstub5_) && (t.getStubRefs().size() > 4)));
318  }
319 
320  private:
323  };
324  struct TTTrackWordChi2RPhiMaxNstubSelector { // using simulated chi2 since not implemented in track word, updates needed
328  : reducedChi2RPhiMaxNstub4_(cfg.template getParameter<double>("reducedChi2RPhiMaxNstub4")),
329  reducedChi2RPhiMaxNstub5_(cfg.template getParameter<double>("reducedChi2RPhiMaxNstub5")) {}
330  bool operator()(const L1Track& t) const {
331  return (((t.getChi2RPhi() < reducedChi2RPhiMaxNstub4_) && (t.getNStubs() == 4)) ||
332  ((t.getChi2RPhi() < reducedChi2RPhiMaxNstub5_) && (t.getNStubs() > 4)));
333  }
334 
335  private:
338  };
343  : reducedBendChi2MaxNstub4_(cfg.template getParameter<double>("reducedBendChi2MaxNstub4")),
344  reducedBendChi2MaxNstub5_(cfg.template getParameter<double>("reducedBendChi2MaxNstub5")) {}
345  bool operator()(const L1Track& t) const {
346  return (((t.stubPtConsistency() < reducedBendChi2MaxNstub4_) && (t.getStubRefs().size() == 4)) ||
347  ((t.stubPtConsistency() < reducedBendChi2MaxNstub5_) && (t.getStubRefs().size() > 4)));
348  }
349 
350  private:
353  };
358  : reducedBendChi2MaxNstub4_(cfg.template getParameter<double>("reducedBendChi2MaxNstub4")),
359  reducedBendChi2MaxNstub5_(cfg.template getParameter<double>("reducedBendChi2MaxNstub5")) {}
360  bool operator()(const L1Track& t) const {
361  return (((t.getBendChi2() < reducedBendChi2MaxNstub4_) && (t.getNStubs() == 4)) ||
362  ((t.getBendChi2() < reducedBendChi2MaxNstub5_) && (t.getNStubs() > 4)));
363  }
364 
365  private:
368  };
369 
387 
388  // ----------member data ---------------------------
398  int debug_;
399 };
400 
401 //
402 // constructors and destructor
403 //
405  : l1TracksToken_(consumes<TTTrackCollection>(iConfig.getParameter<edm::InputTag>("l1TracksInputTag"))),
407  outputCollectionName_(iConfig.getParameter<std::string>("outputCollectionName")),
408  cutSet_(iConfig.getParameter<edm::ParameterSet>("cutSet")),
409 
410  ptMin_(cutSet_.getParameter<double>("ptMin")),
411  absEtaMax_(cutSet_.getParameter<double>("absEtaMax")),
412  absZ0Max_(cutSet_.getParameter<double>("absZ0Max")),
413  promptMVAMin_(cutSet_.getParameter<double>("promptMVAMin")),
414  bendChi2Max_(cutSet_.getParameter<double>("reducedBendChi2Max")),
415  reducedChi2RZMax_(cutSet_.getParameter<double>("reducedChi2RZMax")),
416  reducedChi2RPhiMax_(cutSet_.getParameter<double>("reducedChi2RPhiMax")),
417  reducedChi2RZMaxNstub4_(cutSet_.getParameter<double>("reducedChi2RZMaxNstub4")),
418  reducedChi2RZMaxNstub5_(cutSet_.getParameter<double>("reducedChi2RZMaxNstub5")),
419  reducedChi2RPhiMaxNstub4_(cutSet_.getParameter<double>("reducedChi2RPhiMaxNstub4")),
420  reducedChi2RPhiMaxNstub5_(cutSet_.getParameter<double>("reducedChi2RPhiMaxNstub5")),
421  reducedBendChi2MaxNstub4_(cutSet_.getParameter<double>("reducedBendChi2MaxNstub4")),
422  reducedBendChi2MaxNstub5_(cutSet_.getParameter<double>("reducedBendChi2MaxNstub5")),
423  nStubsMin_(cutSet_.getParameter<int>("nStubsMin")),
424  nPSStubsMin_(cutSet_.getParameter<int>("nPSStubsMin")),
425  processSimulatedTracks_(iConfig.getParameter<bool>("processSimulatedTracks")),
426  processEmulatedTracks_(iConfig.getParameter<bool>("processEmulatedTracks")),
427  debug_(iConfig.getParameter<int>("debug")) {
428  // Confirm the the configuration makes sense
430  throw cms::Exception("You must process at least one of the track collections (simulated or emulated).");
431  }
432 
434  produces<TTTrackRefCollection>(outputCollectionName_);
435  }
437  produces<TTTrackRefCollection>(outputCollectionName_ + "Emulation");
438  }
439 }
440 
442 
443 //
444 // member functions
445 //
446 
448  const TTTrackRefCollectionUPtr& vTTTrackOutput,
449  const TTTrackRefCollectionUPtr& vTTTrackEmulationOutput) const {
450  edm::LogInfo log("L1TrackSelectionProducer");
451  log << "The original track collection (pt, eta, phi, nstub, bendchi2, chi2rz, chi2rphi, z0) values are ... \n";
452  for (const auto& track : *l1TracksHandle) {
453  printTrackInfo(log, track, debug_ >= 4);
454  }
455  log << "\t---\n\tNumber of tracks in this selection = " << l1TracksHandle->size() << "\n\n";
457  log << "The selected track collection (pt, eta, phi, nstub, bendchi2, chi2rz, chi2rphi, z0) values are ... \n";
458  for (const auto& track : *vTTTrackOutput) {
459  printTrackInfo(log, *track, debug_ >= 4);
460  }
461  log << "\t---\n\tNumber of tracks in this selection = " << vTTTrackOutput->size() << "\n\n";
462  }
464  log << "The emulation selected track collection (pt, eta, phi, nstub, bendchi2, chi2rz, chi2rphi, z0) values are "
465  "... \n";
466  for (const auto& track : *vTTTrackEmulationOutput) {
467  printTrackInfo(log, *track, debug_ >= 4);
468  }
469  log << "\t---\n\tNumber of tracks in this selection = " << vTTTrackEmulationOutput->size() << "\n\n";
470  }
472  TTTrackRefCollection inSimButNotEmu;
473  TTTrackRefCollection inEmuButNotSim;
474  std::set_difference(vTTTrackOutput->begin(),
475  vTTTrackOutput->end(),
476  vTTTrackEmulationOutput->begin(),
477  vTTTrackEmulationOutput->end(),
478  std::back_inserter(inSimButNotEmu));
479  std::set_difference(vTTTrackEmulationOutput->begin(),
480  vTTTrackEmulationOutput->end(),
481  vTTTrackOutput->begin(),
482  vTTTrackOutput->end(),
483  std::back_inserter(inEmuButNotSim));
484  log << "The set of tracks selected via cuts on the simulated values which are not in the set of tracks selected "
485  "by cutting on the emulated values ... \n";
486  for (const auto& track : inSimButNotEmu) {
487  printTrackInfo(log, *track, debug_ >= 3);
488  }
489  log << "\t---\n\tNumber of tracks in this selection = " << inSimButNotEmu.size() << "\n\n"
490  << "The set of tracks selected via cuts on the emulated values which are not in the set of tracks selected "
491  "by cutting on the simulated values ... \n";
492  for (const auto& track : inEmuButNotSim) {
493  printTrackInfo(log, *track, debug_ >= 3);
494  }
495  log << "\t---\n\tNumber of tracks in this selection = " << inEmuButNotSim.size() << "\n\n";
496  }
497 }
498 
499 void L1TrackSelectionProducer::printTrackInfo(edm::LogInfo& log, const L1Track& track, bool printEmulation) const {
500  log << "\t(" << track.momentum().perp() << ", " << track.momentum().eta() << ", " << track.momentum().phi() << ", "
501  << track.getStubRefs().size() << ", " << track.stubPtConsistency() << ", " << track.chi2ZRed() << ", "
502  << track.chi2XYRed() << ", " << track.z0() << ")\n";
503 
504  if (printEmulation) {
505  ap_uint<TrackBitWidths::kPtSize> ptEmulationBits = track.getTrackWord()(
506  TTTrack_TrackWord::TrackBitLocations::kRinvMSB - 1, TTTrack_TrackWord::TrackBitLocations::kRinvLSB);
507  ap_ufixed<TrackBitWidths::kPtSize, TrackBitWidths::kPtMagSize> ptEmulation;
508  ptEmulation.V = ptEmulationBits.range();
509  TTTrack_TrackWord::tanl_t etaEmulationBits = track.getTanlWord();
510  ap_fixed<TrackBitWidths::kEtaSize, TrackBitWidths::kEtaMagSize> etaEmulation;
511  etaEmulation.V = etaEmulationBits.range();
512  double floatTkZ0 = track.undigitizeSignedValue(
513  track.getZ0Bits(), TTTrack_TrackWord::TrackBitWidths::kZ0Size, TTTrack_TrackWord::stepZ0, 0.0);
514  double floatTkPhi = track.undigitizeSignedValue(
515  track.getPhiBits(), TTTrack_TrackWord::TrackBitWidths::kPhiSize, TTTrack_TrackWord::stepPhi0, 0.0);
516  log << "\t\t(" << ptEmulation.to_double() << ", " << etaEmulation.to_double() << ", " << floatTkPhi << ", "
517  << track.getNStubs() << ", " << track.getBendChi2() << ", " << track.getChi2RZ() << ", " << track.getChi2RPhi()
518  << ", " << floatTkZ0 << ")\n";
519  }
520 }
521 
522 // ------------ method called to produce the data ------------
524  auto vTTTrackOutput = std::make_unique<TTTrackRefCollection>();
525  auto vTTTrackEmulationOutput = std::make_unique<TTTrackRefCollection>();
526 
527  // Tracker Topology
528  const TrackerTopology& tTopo = iSetup.getData(tTopoToken_);
529 
530  TTTrackCollectionHandle l1TracksHandle;
531 
532  iEvent.getByToken(l1TracksToken_, l1TracksHandle);
533  size_t nOutputApproximate = l1TracksHandle->size();
535  vTTTrackOutput->reserve(nOutputApproximate);
536  }
538  vTTTrackEmulationOutput->reserve(nOutputApproximate);
539  }
540 
545  TTTrackNPSStubsMinSelector nPSStubsSel(nPSStubsMin_, tTopo);
554 
555  for (size_t i = 0; i < nOutputApproximate; i++) {
556  const auto& track = l1TracksHandle->at(i);
557 
558  // Select tracks based on the floating point TTTrack
559  if (processSimulatedTracks_ && kinSel(track) && nPSStubsSel(track) && chi2Sel(track) && mvaSel(track) &&
560  chi2NstubSel(track)) {
561  vTTTrackOutput->push_back(TTTrackRef(l1TracksHandle, i));
562  }
563 
564  // Select tracks based on the bitwise accurate TTTrack_TrackWord
565  if (processEmulatedTracks_ && kinSelEmu(track) && chi2SelEmu(track) && mvaSelEmu(track) && chi2NstubSelEmu(track)) {
566  vTTTrackEmulationOutput->push_back(TTTrackRef(l1TracksHandle, i));
567  }
568  }
569 
570  if (debug_ >= 2) {
571  printDebugInfo(l1TracksHandle, vTTTrackOutput, vTTTrackEmulationOutput);
572  }
573 
574  // Put the outputs into the event
576  iEvent.put(std::move(vTTTrackOutput), outputCollectionName_);
577  }
579  iEvent.put(std::move(vTTTrackEmulationOutput), outputCollectionName_ + "Emulation");
580  }
581 }
582 
583 // ------------ method fills 'descriptions' with the allowed parameters for the module ------------
585  //L1TrackSelectionProducer
587  desc.add<edm::InputTag>("l1TracksInputTag", edm::InputTag("TTTracksFromTrackletEmulation", "Level1TTTracks"));
588  desc.add<std::string>("outputCollectionName", "Level1TTTracksSelected");
589  {
590  edm::ParameterSetDescription descCutSet;
591  descCutSet.add<double>("ptMin", 2.0)->setComment("pt must be greater than this value, [GeV]");
592  descCutSet.add<double>("absEtaMax", 2.4)->setComment("absolute value of eta must be less than this value");
593  descCutSet.add<double>("absZ0Max", 15.0)->setComment("z0 must be less than this value, [cm]");
594  descCutSet.add<int>("nStubsMin", 4)->setComment("number of stubs must be greater than or equal to this value");
595  descCutSet.add<int>("nPSStubsMin", 0)
596  ->setComment("number of stubs in the PS Modules must be greater than or equal to this value");
597 
598  descCutSet.add<double>("promptMVAMin", -1.0)->setComment("MVA must be greater than this value");
599  descCutSet.add<double>("reducedBendChi2Max", 2.25)->setComment("bend chi2 must be less than this value");
600  descCutSet.add<double>("reducedChi2RZMax", 5.0)->setComment("chi2rz/dof must be less than this value");
601  descCutSet.add<double>("reducedChi2RPhiMax", 20.0)->setComment("chi2rphi/dof must be less than this value");
602  descCutSet.add<double>("reducedChi2RZMaxNstub4", 999.9)
603  ->setComment("chi2rz/dof must be less than this value in nstub==4");
604  descCutSet.add<double>("reducedChi2RZMaxNstub5", 999.9)
605  ->setComment("chi2rz/dof must be less than this value in nstub>4");
606  descCutSet.add<double>("reducedChi2RPhiMaxNstub4", 999.9)
607  ->setComment("chi2rphi/dof must be less than this value in nstub==4");
608  descCutSet.add<double>("reducedChi2RPhiMaxNstub5", 999.9)
609  ->setComment("chi2rphi/dof must be less than this value in nstub>4");
610  descCutSet.add<double>("reducedBendChi2MaxNstub4", 999.9)
611  ->setComment("bend chi2 must be less than this value in nstub==4");
612  descCutSet.add<double>("reducedBendChi2MaxNstub5", 999.9)
613  ->setComment("bend chi2 must be less than this value in nstub>4");
614 
615  desc.add<edm::ParameterSetDescription>("cutSet", descCutSet);
616  }
617  desc.add<bool>("processSimulatedTracks", true)
618  ->setComment("return selected tracks after cutting on the floating point values");
619  desc.add<bool>("processEmulatedTracks", true)
620  ->setComment("return selected tracks after cutting on the bitwise emulated values");
621  desc.add<int>("debug", 0)->setComment("Verbosity levels: 0, 1, 2, 3");
622  descriptions.addWithDefaultLabel(desc);
623 }
624 
625 //define this as a plug-in
L1TrackSelectionProducer(const edm::ParameterSet &)
void setComment(std::string const &value)
void addWithDefaultLabel(ParameterSetDescription const &psetDescription)
AndSelector< TTTrackChi2RZMaxNstubSelector, TTTrackChi2RPhiMaxNstubSelector, TTTrackBendChi2MaxNstubSelector > TTTrackChi2MaxNstubSelector
edm::ESGetToken< TrackerTopology, TrackerTopologyRcd > tTopoToken_
TTTrackBendChi2MaxNstubSelector(double reducedBendChi2MaxNstub4, double reducedBendChi2MaxNstub5)
ESGetTokenH3DDVariant esConsumes(std::string const &Record, edm::ConsumesCollector &)
Definition: DeDxTools.cc:283
unsigned int tobLayer(const DetId &id) const
edm::Handle< TTTrackCollection > TTTrackCollectionHandle
float floatZ0(z0_t z0)
Definition: datatypes.h:170
TTTrackWordBendChi2MaxNstubSelector(double reducedBendChi2MaxNstub4, double reducedBendChi2MaxNstub5)
T const & getData(const ESGetToken< T, R > &iToken) const noexcept(false)
Definition: EventSetup.h:119
AndSelector< TTTrackBendChi2MaxSelector, TTTrackChi2RZMaxSelector, TTTrackChi2RPhiMaxSelector > TTTrackBendChi2Chi2RZChi2RPhiMaxSelector
AndSelector< TTTrackWordBendChi2MaxSelector, TTTrackWordChi2RZMaxSelector, TTTrackWordChi2RPhiMaxSelector > TTTrackWordBendChi2Chi2RZChi2RPhiMaxSelector
AndSelector< TTTrackWordPtMinSelector, TTTrackWordAbsEtaMaxSelector, TTTrackWordAbsZ0MaxSelector, TTTrackWordNStubsMinSelector > TTTrackWordPtMinEtaMaxZ0MaxNStubsMinSelector
AndSelector< TTTrackPtMinSelector, TTTrackAbsEtaMaxSelector, TTTrackAbsZ0MaxSelector, TTTrackNStubsMinSelector > TTTrackPtMinEtaMaxZ0MaxNStubsMinSelector
constexpr float ptMin
TTTrackNPSStubsMinSelector(const edm::ParameterSet &cfg, const TrackerTopology &tTopo)
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
std::unique_ptr< TTTrackRefCollection > TTTrackRefCollectionUPtr
void produce(edm::StreamID, edm::Event &, const edm::EventSetup &) const override
TTTrackNPSStubsMinSelector(double nStubsMin, const TrackerTopology &tTopo)
AndSelector< TTTrackWordChi2RZMaxNstubSelector, TTTrackWordChi2RPhiMaxNstubSelector, TTTrackWordBendChi2MaxNstubSelector > TTTrackWordChi2MaxNstubSelector
int iEvent
Definition: GenABIO.cc:224
TTTrackChi2RPhiMaxNstubSelector(double reducedChi2RPhiMaxNstub4, double reducedChi2RPhiMaxNstub5)
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
edm::Ref< TTTrackCollection > TTTrackRef
def template(fileName, svg, replaceme="REPLACEME")
Definition: svgfig.py:521
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
const edm::ParameterSet cutSet_
static constexpr auto TOB
TTTrack< Ref_Phase2TrackerDigi_ > L1Track
ParameterDescriptionBase * add(U const &iLabel, T const &value)
std::vector< L1Track > TTTrackCollection
TTTrackWordChi2RPhiMaxNstubSelector(double reducedChi2RPhiMaxNstub4, double reducedChi2RPhiMaxNstub5)
Definition: DetId.h:17
TTTrackChi2RZMaxNstubSelector(double reducedChi2RZMaxNstub4, double reducedChi2RZMaxNstub5)
Class to store the L1 Track Trigger tracks.
Definition: TTTrack.h:29
void printTrackInfo(edm::LogInfo &log, const L1Track &track, bool printEmulation=false) const
static constexpr double stepPhi0
void printDebugInfo(const TTTrackCollectionHandle &l1TracksHandle, const TTTrackRefCollectionUPtr &vTTTrackOutput, const TTTrackRefCollectionUPtr &vTTTrackEmulationOutput) const
edm::RefVector< TTTrackCollection > TTTrackRefCollection
HLT enums.
ap_uint< TrackBitWidths::kTanlSize > tanl_t
unsigned int tidRing(const DetId &id) const
TTTrackWordChi2RZMaxNstubSelector(double reducedChi2RZMaxNstub4, double reducedChi2RZMaxNstub5)
static constexpr double stepZ0
static constexpr auto TID
std::vector< std::string > set_difference(std::vector< std::string > const &v1, std::vector< std::string > const &v2)
def move(src, dest)
Definition: eostools.py:511
const edm::EDGetTokenT< TTTrackCollection > l1TracksToken_