CMS 3D CMS Logo

List of all members | Public Types | Public Member Functions | Static Public Member Functions | Private Member Functions | Private Attributes
L1TrackJetProducer Class Reference
Inheritance diagram for L1TrackJetProducer:
edm::stream::EDProducer<>

Public Types

typedef vector< L1TTTrackTypeL1TTTrackCollectionType
 
typedef TTTrack< Ref_Phase2TrackerDigi_L1TTTrackType
 
- Public Types inherited from edm::stream::EDProducer<>
using CacheTypes = CacheContexts< T... >
 
using GlobalCache = typename CacheTypes::GlobalCache
 
using HasAbility = AbilityChecker< T... >
 
using InputProcessBlockCache = typename CacheTypes::InputProcessBlockCache
 
using LuminosityBlockCache = typename CacheTypes::LuminosityBlockCache
 
using LuminosityBlockContext = LuminosityBlockContextT< LuminosityBlockCache, RunCache, GlobalCache >
 
using LuminosityBlockSummaryCache = typename CacheTypes::LuminosityBlockSummaryCache
 
using RunCache = typename CacheTypes::RunCache
 
using RunContext = RunContextT< RunCache, GlobalCache >
 
using RunSummaryCache = typename CacheTypes::RunSummaryCache
 

Public Member Functions

 L1TrackJetProducer (const ParameterSet &)
 
bool trackQualityCuts (int trk_nstub, float trk_chi2, float trk_bendchi2)
 
 ~L1TrackJetProducer () override
 
- Public Member Functions inherited from edm::stream::EDProducer<>
 EDProducer ()=default
 
 EDProducer (const EDProducer &)=delete
 
bool hasAbilityToProduceInBeginLumis () const final
 
bool hasAbilityToProduceInBeginProcessBlocks () const final
 
bool hasAbilityToProduceInBeginRuns () const final
 
bool hasAbilityToProduceInEndLumis () const final
 
bool hasAbilityToProduceInEndProcessBlocks () const final
 
bool hasAbilityToProduceInEndRuns () const final
 
const EDProduceroperator= (const EDProducer &)=delete
 

Static Public Member Functions

static void fillDescriptions (ConfigurationDescriptions &descriptions)
 

Private Member Functions

void beginStream (StreamID) override
 
void endStream () override
 
void produce (Event &, const EventSetup &) override
 

Private Attributes

const float d0CutNStubs4_
 
const float d0CutNStubs5_
 
const bool displaced_
 
const float dzPVTrk_
 
const int etaBins_
 
float etaStep_
 
const int highpTJetMinTrackMultiplicity_
 
const float highpTJetThreshold_
 
vector< Ptr< L1TTTrackType > > L1TrkPtrs_
 
const int lowpTJetMinTrackMultiplicity_
 
const float lowpTJetThreshold_
 
const double minTrkJetpT_
 
const int nDisplacedTracks_
 
const float nStubs4DisplacedBend_
 
const float nStubs4DisplacedChi2_
 
const float nStubs4PromptBend_
 
const float nStubs4PromptChi2_
 
const float nStubs5DisplacedBend_
 
const float nStubs5DisplacedChi2_
 
const float nStubs5PromptBend_
 
const float nStubs5PromptChi2_
 
const int phiBins_
 
float phiStep_
 
const edm::EDGetTokenT< std::vector< l1t::Vertex > > PVtxToken_
 
vector< int > tdtrk_
 
const EDGetTokenT< vector< TTTrack< Ref_Phase2TrackerDigi_ > > > trackToken_
 
const float trkEtaMax_
 
const int trkNPSStubMin_
 
const float trkPtMax_
 
const float trkPtMin_
 
const float trkZMax_
 
edm::ESGetToken< TrackerTopology, TrackerTopologyRcdtTopoToken_
 
const int zBins_
 
float zStep_
 

Detailed Description

Definition at line 49 of file L1TrackJetProducer.cc.

Member Typedef Documentation

◆ L1TTTrackCollectionType

Definition at line 54 of file L1TrackJetProducer.cc.

◆ L1TTTrackType

Definition at line 53 of file L1TrackJetProducer.cc.

Constructor & Destructor Documentation

◆ L1TrackJetProducer()

L1TrackJetProducer::L1TrackJetProducer ( const ParameterSet iConfig)
explicit

Definition at line 103 of file L1TrackJetProducer.cc.

References displaced_, etaBins_, etaStep_, M_PI, phiBins_, phiStep_, trkEtaMax_, trkZMax_, zBins_, and zStep_.

104  : trkZMax_((float)iConfig.getParameter<double>("trk_zMax")),
105  trkPtMax_((float)iConfig.getParameter<double>("trk_ptMax")),
106  trkPtMin_((float)iConfig.getParameter<double>("trk_ptMin")),
107  trkEtaMax_((float)iConfig.getParameter<double>("trk_etaMax")),
108  nStubs4PromptChi2_((float)iConfig.getParameter<double>("nStubs4PromptChi2")),
109  nStubs5PromptChi2_((float)iConfig.getParameter<double>("nStubs5PromptChi2")),
110  nStubs4PromptBend_((float)iConfig.getParameter<double>("nStubs4PromptBend")),
111  nStubs5PromptBend_((float)iConfig.getParameter<double>("nStubs5PromptBend")),
112  trkNPSStubMin_((int)iConfig.getParameter<int>("trk_nPSStubMin")),
113  lowpTJetMinTrackMultiplicity_((int)iConfig.getParameter<int>("lowpTJetMinTrackMultiplicity")),
114  lowpTJetThreshold_((float)iConfig.getParameter<double>("lowpTJetThreshold")),
115  highpTJetMinTrackMultiplicity_((int)iConfig.getParameter<int>("highpTJetMinTrackMultiplicity")),
116  highpTJetThreshold_((float)iConfig.getParameter<double>("highpTJetThreshold")),
117  zBins_((int)iConfig.getParameter<int>("zBins")),
118  etaBins_((int)iConfig.getParameter<int>("etaBins")),
119  phiBins_((int)iConfig.getParameter<int>("phiBins")),
120  minTrkJetpT_(iConfig.getParameter<double>("minTrkJetpT")),
121  displaced_(iConfig.getParameter<bool>("displaced")),
122  d0CutNStubs4_((float)iConfig.getParameter<double>("d0_cutNStubs4")),
123  d0CutNStubs5_((float)iConfig.getParameter<double>("d0_cutNStubs5")),
124  nStubs4DisplacedChi2_((float)iConfig.getParameter<double>("nStubs4DisplacedChi2")),
125  nStubs5DisplacedChi2_((float)iConfig.getParameter<double>("nStubs5DisplacedChi2")),
126  nStubs4DisplacedBend_((float)iConfig.getParameter<double>("nStubs4DisplacedBend")),
127  nStubs5DisplacedBend_((float)iConfig.getParameter<double>("nStubs5DisplacedBend")),
128  nDisplacedTracks_((int)iConfig.getParameter<int>("nDisplacedTracks")),
129  dzPVTrk_((float)iConfig.getParameter<double>("MaxDzTrackPV")),
130  tTopoToken_(esConsumes<TrackerTopology, TrackerTopologyRcd>(edm::ESInputTag("", ""))),
131  trackToken_(consumes<vector<TTTrack<Ref_Phase2TrackerDigi_>>>(iConfig.getParameter<InputTag>("L1TrackInputTag"))),
132  PVtxToken_(consumes<vector<l1t::Vertex>>(iConfig.getParameter<InputTag>("L1PVertexCollection"))) {
133  zStep_ = 2.0 * trkZMax_ / (zBins_ + 1); // added +1 in denom
134  etaStep_ = 2.0 * trkEtaMax_ / etaBins_; //etaStep is the width of an etabin
135  phiStep_ = 2 * M_PI / phiBins_;
136 
137  if (displaced_)
138  produces<TkJetCollection>("L1TrackJetsExtended");
139  else
140  produces<TkJetCollection>("L1TrackJets");
141 }
const int highpTJetMinTrackMultiplicity_
T getParameter(std::string const &) const
Definition: ParameterSet.h:303
const float lowpTJetThreshold_
const float nStubs4DisplacedBend_
edm::ESGetToken< TrackerTopology, TrackerTopologyRcd > tTopoToken_
const float nStubs4DisplacedChi2_
const float nStubs5PromptChi2_
const float nStubs4PromptChi2_
const float highpTJetThreshold_
const float nStubs5DisplacedBend_
const float nStubs4PromptBend_
const float nStubs5PromptBend_
const float nStubs5DisplacedChi2_
#define M_PI
const edm::EDGetTokenT< std::vector< l1t::Vertex > > PVtxToken_
Class to store the L1 Track Trigger tracks.
Definition: TTTrack.h:29
const int lowpTJetMinTrackMultiplicity_
const EDGetTokenT< vector< TTTrack< Ref_Phase2TrackerDigi_ > > > trackToken_

◆ ~L1TrackJetProducer()

L1TrackJetProducer::~L1TrackJetProducer ( )
override

Definition at line 143 of file L1TrackJetProducer.cc.

143 {}

Member Function Documentation

◆ beginStream()

void L1TrackJetProducer::beginStream ( StreamID  )
overrideprivate

Definition at line 333 of file L1TrackJetProducer.cc.

333 {}

◆ endStream()

void L1TrackJetProducer::endStream ( )
overrideprivate

Definition at line 335 of file L1TrackJetProducer.cc.

335 {}

◆ fillDescriptions()

void L1TrackJetProducer::fillDescriptions ( ConfigurationDescriptions descriptions)
static

Definition at line 357 of file L1TrackJetProducer.cc.

References edm::ConfigurationDescriptions::add(), submitPVResolutionJobs::desc, and HLT_2022v15_cff::InputTag.

357  {
358  // l1tTrackJets
360  desc.add<edm::InputTag>("L1TrackInputTag", edm::InputTag("l1tTTTracksFromTrackletEmulation", "Level1TTTracks"));
361  desc.add<edm::InputTag>("L1PVertexCollection", edm::InputTag("l1tVertexProducer", "l1vertices"));
362  desc.add<double>("MaxDzTrackPV", 1.0);
363  desc.add<double>("trk_zMax", 15.0);
364  desc.add<double>("trk_ptMax", 200.0);
365  desc.add<double>("trk_ptMin", 3.0);
366  desc.add<double>("trk_etaMax", 2.4);
367  desc.add<double>("nStubs4PromptChi2", 5.0);
368  desc.add<double>("nStubs4PromptBend", 1.7);
369  desc.add<double>("nStubs5PromptChi2", 2.75);
370  desc.add<double>("nStubs5PromptBend", 3.5);
371  desc.add<int>("trk_nPSStubMin", -1);
372  desc.add<double>("minTrkJetpT", -1.0);
373  desc.add<int>("etaBins", 24);
374  desc.add<int>("phiBins", 27);
375  desc.add<int>("zBins", 1);
376  desc.add<double>("d0_cutNStubs4", -1);
377  desc.add<double>("d0_cutNStubs5", -1);
378  desc.add<int>("lowpTJetMinTrackMultiplicity", 2);
379  desc.add<double>("lowpTJetThreshold", 50.0);
380  desc.add<int>("highpTJetMinTrackMultiplicity", 3);
381  desc.add<double>("highpTJetThreshold", 100.0);
382  desc.add<bool>("displaced", false);
383  desc.add<double>("nStubs4DisplacedChi2", 5.0);
384  desc.add<double>("nStubs4DisplacedBend", 1.7);
385  desc.add<double>("nStubs5DisplacedChi2", 2.75);
386  desc.add<double>("nStubs5DisplacedBend", 3.5);
387  desc.add<int>("nDisplacedTracks", 2);
388  descriptions.add("l1tTrackJets", desc);
389 }
void add(std::string const &label, ParameterSetDescription const &psetDescription)

◆ produce()

void L1TrackJetProducer::produce ( Event iEvent,
const EventSetup iSetup 
)
overrideprivate

Definition at line 145 of file L1TrackJetProducer.cc.

References funct::abs(), TTTrack< T >::chi2Red(), MaxZBin::clusters, filterCSVwithJSON::copy, funct::cos(), TTTrack< T >::d0(), d0CutNStubs4_, d0CutNStubs5_, displaced_, dzPVTrk_, EtaPhiBin::eta, PVValHelper::eta, PV3DBase< T, PVType, FrameType >::eta(), egammaIdentification::eta_max, egammaIdentification::eta_min, etaBins_, etaStep_, edm::EventSetup::getData(), TTTrack< T >::getStubRefs(), highpTJetMinTrackMultiplicity_, highpTJetThreshold_, MaxZBin::ht, mps_fire::i, iEvent, createfilelist::int, dqmiolumiharvest::j, reco::btau::jetEta, reco::btau::jetPhi, reco::btau::jetPt, dqmdumpme::k, L1_clustering(), L1TrkPtrs_, L2_clustering(), lowpTJetMinTrackMultiplicity_, lowpTJetThreshold_, M_PI, minTrkJetpT_, TTTrack< T >::momentum(), eostools::move(), MaxZBin::nclust, nDisplacedTracks_, EtaPhiBin::numtdtrks, EtaPhiBin::numtracks, PV3DBase< T, PVType, FrameType >::perp(), phi, EtaPhiBin::phi, phiBins_, phiStep_, EtaPhiBin::pTtot, PVtxToken_, funct::sin(), TTTrack< T >::stubPtConsistency(), tdtrk_, StripSubdetector::TID, TrackerTopology::tidRing(), StripSubdetector::TOB, TrackerTopology::tobLayer(), align::Tracker, EtaPhiBin::trackidx, trackQualityCuts(), trackToken_, trkEtaMax_, trkNPSStubMin_, trkPtMax_, trkPtMin_, trkZMax_, tTopoToken_, HLTMuonOfflineAnalyzer_cfi::z0, TTTrack< T >::z0(), MaxZBin::zbincenter, zBins_, SiStripMonitorCluster_cfi::zmax, SiStripMonitorCluster_cfi::zmin, MaxZBin::znum, and zStep_.

145  {
146  unique_ptr<TkJetCollection> L1L1TrackJetProducer(new TkJetCollection);
147 
148  // Read inputs
149  const TrackerTopology &tTopo = iSetup.getData(tTopoToken_);
150 
152  iEvent.getByToken(trackToken_, TTTrackHandle);
153 
155  iEvent.getByToken(PVtxToken_, PVtx);
156  float PVz = (PVtx->at(0)).z0();
157 
158  L1TrkPtrs_.clear();
159  tdtrk_.clear();
160 
161  for (unsigned int this_l1track = 0; this_l1track < TTTrackHandle->size(); this_l1track++) {
162  edm::Ptr<L1TTTrackType> trkPtr(TTTrackHandle, this_l1track);
163  float trk_pt = trkPtr->momentum().perp();
164  int trk_nstubs = (int)trkPtr->getStubRefs().size();
165  float trk_chi2dof = trkPtr->chi2Red();
166  float trk_d0 = trkPtr->d0();
167  float trk_bendchi2 = trkPtr->stubPtConsistency();
168  float trk_z0 = trkPtr->z0();
169 
170  int trk_nPS = 0;
171  for (int istub = 0; istub < trk_nstubs; istub++) { // loop over the stubs
172  DetId detId(trkPtr->getStubRefs().at(istub)->getDetId());
173  if (detId.det() == DetId::Detector::Tracker) {
174  if ((detId.subdetId() == StripSubdetector::TOB && tTopo.tobLayer(detId) <= 3) ||
175  (detId.subdetId() == StripSubdetector::TID && tTopo.tidRing(detId) <= 9))
176  trk_nPS++;
177  }
178  }
179  // select tracks
180  if (trk_nPS < trkNPSStubMin_)
181  continue;
182  if (!trackQualityCuts(trk_nstubs, trk_chi2dof, trk_bendchi2))
183  continue;
184  if (std::abs(PVz - trk_z0) > dzPVTrk_ && dzPVTrk_ > 0)
185  continue;
186  if (std::abs(trk_z0) > trkZMax_)
187  continue;
188  if (std::abs(trkPtr->momentum().eta()) > trkEtaMax_)
189  continue;
190  if (trk_pt < trkPtMin_)
191  continue;
192  L1TrkPtrs_.push_back(trkPtr);
193 
194  if ((std::abs(trk_d0) > d0CutNStubs5_ && trk_nstubs >= 5 && d0CutNStubs5_ >= 0) ||
195  (trk_nstubs == 4 && std::abs(trk_d0) > d0CutNStubs4_ && d0CutNStubs4_ >= 0))
196  tdtrk_.push_back(1); //displaced track
197  else
198  tdtrk_.push_back(0); // not displaced track
199  }
200 
201  if (!L1TrkPtrs_.empty()) {
202  MaxZBin mzb;
203  mzb.znum = -1;
204  mzb.nclust = -1;
205  mzb.ht = -1;
206  EtaPhiBin epbins_default[phiBins_][etaBins_]; // create grid of phiBins
207 
208  float phi = -1.0 * M_PI;
209  for (int i = 0; i < phiBins_; ++i) {
210  float eta = -1.0 * trkEtaMax_;
211  for (int j = 0; j < etaBins_; ++j) {
212  epbins_default[i][j].phi = (phi + (phi + phiStep_)) / 2.0; // phimin=phi; phimax= phimin+step
213  epbins_default[i][j].eta = (eta + (eta + etaStep_)) / 2.0; // phimin=phi; phimax phimin+step
214  eta += etaStep_;
215  } // for each etabin
216  phi += phiStep_;
217  } // for each phibin (finished creating epbins)
218 
219  std::vector<float> zmins, zmaxs;
220  for (int zbin = 0; zbin < zBins_; zbin++) {
221  zmins.push_back(-1.0 * trkZMax_ + zStep_ * zbin);
222  zmaxs.push_back(-1.0 * trkZMax_ + zStep_ * zbin + zStep_ * 2);
223  }
224 
225  // create vectors that hold data
226  std::vector<std::vector<EtaPhiBin>> L1clusters;
227  L1clusters.reserve(phiBins_);
228  std::vector<EtaPhiBin> L2clusters;
229 
230  for (unsigned int zbin = 0; zbin < zmins.size(); ++zbin) {
231  // initialize matrices
232  float zmin = zmins[zbin];
233  float zmax = zmaxs[zbin];
234  EtaPhiBin epbins[phiBins_][etaBins_];
235  std::copy(&epbins_default[0][0], &epbins_default[0][0] + phiBins_ * etaBins_, &epbins[0][0]);
236 
237  //clear containers
238  L1clusters.clear();
239  L2clusters.clear();
240 
241  // fill grid
242  for (unsigned int k = 0; k < L1TrkPtrs_.size(); ++k) {
243  float trkZ = L1TrkPtrs_[k]->z0();
244  if (zmax < trkZ)
245  continue;
246  if (zmin > trkZ)
247  continue;
248  if (zbin == 0 && zmin == trkZ)
249  continue;
250  float trkpt = L1TrkPtrs_[k]->momentum().perp();
251  float trketa = L1TrkPtrs_[k]->momentum().eta();
252  float trkphi = L1TrkPtrs_[k]->momentum().phi();
253  for (int i = 0; i < phiBins_; ++i) {
254  for (int j = 0; j < etaBins_; ++j) {
255  float eta_min = epbins[i][j].eta - etaStep_ / 2.0; //eta min
256  float eta_max = epbins[i][j].eta + etaStep_ / 2.0; //eta max
257  float phi_min = epbins[i][j].phi - phiStep_ / 2.0; //phi min
258  float phi_max = epbins[i][j].phi + phiStep_ / 2.0; //phi max
259 
260  if ((trketa < eta_min) || (trketa > eta_max) || (trkphi < phi_min) || (trkphi > phi_max))
261  continue;
262 
263  if (trkpt < trkPtMax_)
264  epbins[i][j].pTtot += trkpt;
265  else
266  epbins[i][j].pTtot += trkPtMax_;
267  epbins[i][j].numtdtrks += tdtrk_[k];
268  epbins[i][j].trackidx.push_back(k);
269  ++epbins[i][j].numtracks;
270  } // for each phibin
271  } // for each phibin
272  } //end loop over tracks
273 
274  // first layer clustering - in eta for all phi bins
275  for (int phibin = 0; phibin < phiBins_; ++phibin) {
276  L1clusters.push_back(L1_clustering(epbins[phibin], etaBins_, etaStep_));
277  }
278 
279  //second layer clustering in phi for all eta clusters
280  L2clusters = L2_clustering(L1clusters, phiBins_, phiStep_, etaStep_);
281 
282  // sum all cluster pTs in this zbin to find max
283  float sum_pt = 0;
284  for (unsigned int k = 0; k < L2clusters.size(); ++k) {
285  if (L2clusters[k].pTtot > lowpTJetThreshold_ && L2clusters[k].numtracks < lowpTJetMinTrackMultiplicity_)
286  continue;
287  if (L2clusters[k].pTtot > highpTJetThreshold_ && L2clusters[k].numtracks < highpTJetMinTrackMultiplicity_)
288  continue;
289  if (L2clusters[k].pTtot > minTrkJetpT_)
290  sum_pt += L2clusters[k].pTtot;
291  }
292 
293  if (sum_pt < mzb.ht)
294  continue;
295  // if ht is larger than previous max, this is the new vertex zbin
296  mzb.ht = sum_pt;
297  mzb.znum = zbin;
298  mzb.clusters = L2clusters;
299  mzb.nclust = L2clusters.size();
300  mzb.zbincenter = (zmin + zmax) / 2.0;
301  } //zbin loop
302 
303  vector<Ptr<L1TTTrackType>> L1TrackAssocJet;
304  for (unsigned int j = 0; j < mzb.clusters.size(); ++j) {
305  float jetEta = mzb.clusters[j].eta;
306  float jetPhi = mzb.clusters[j].phi;
307  float jetPt = mzb.clusters[j].pTtot;
308  float jetPx = jetPt * cos(jetPhi);
309  float jetPy = jetPt * sin(jetPhi);
310  float jetPz = jetPt * sinh(jetEta);
311  float jetP = jetPt * cosh(jetEta);
312  int totalDisptrk = mzb.clusters[j].numtdtrks;
313  bool isDispJet = (totalDisptrk > nDisplacedTracks_ || totalDisptrk == nDisplacedTracks_);
314 
315  math::XYZTLorentzVector jetP4(jetPx, jetPy, jetPz, jetP);
316  L1TrackAssocJet.clear();
317  for (unsigned int itrk = 0; itrk < mzb.clusters[j].trackidx.size(); itrk++)
318  L1TrackAssocJet.push_back(L1TrkPtrs_[mzb.clusters[j].trackidx[itrk]]);
319 
320  TkJet trkJet(jetP4, L1TrackAssocJet, mzb.zbincenter, mzb.clusters[j].numtracks, 0, totalDisptrk, 0, isDispJet);
321 
322  if (!L1TrackAssocJet.empty())
323  L1L1TrackJetProducer->push_back(trkJet);
324  }
325 
326  if (displaced_)
327  iEvent.put(std::move(L1L1TrackJetProducer), "L1TrackJetsExtended");
328  else
329  iEvent.put(std::move(L1L1TrackJetProducer), "L1TrackJets");
330  }
331 }
bool trackQualityCuts(int trk_nstub, float trk_chi2, float trk_bendchi2)
const int highpTJetMinTrackMultiplicity_
unsigned int tobLayer(const DetId &id) const
std::vector< EtaPhiBin > L2_clustering(std::vector< std::vector< EtaPhiBin >> &L1clusters, int phiBins_, float phiStep_, float etaStep_)
Definition: L1Clustering.h:136
const float lowpTJetThreshold_
T const & getData(const ESGetToken< T, R > &iToken) const noexcept(false)
Definition: EventSetup.h:119
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
float phi
Definition: L1Clustering.h:17
edm::ESGetToken< TrackerTopology, TrackerTopologyRcd > tTopoToken_
std::vector< unsigned int > trackidx
Definition: L1Clustering.h:19
float ht
Definition: L1Clustering.h:28
XYZTLorentzVectorD XYZTLorentzVector
Lorentz vector with cylindrical internal representation using pseudorapidity.
Definition: LorentzVector.h:29
int nclust
Definition: L1Clustering.h:25
int iEvent
Definition: GenABIO.cc:224
const float highpTJetThreshold_
float zbincenter
Definition: L1Clustering.h:26
vector< Ptr< L1TTTrackType > > L1TrkPtrs_
std::vector< EtaPhiBin > L1_clustering(EtaPhiBin *phislice, int etaBins_, float etaStep_)
Definition: L1Clustering.h:31
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
float pTtot
Definition: L1Clustering.h:11
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
int numtdtrks
Definition: L1Clustering.h:14
static constexpr auto TOB
#define M_PI
const edm::EDGetTokenT< std::vector< l1t::Vertex > > PVtxToken_
Definition: DetId.h:17
float eta
Definition: L1Clustering.h:18
int numtracks
Definition: L1Clustering.h:12
const int lowpTJetMinTrackMultiplicity_
const EDGetTokenT< vector< TTTrack< Ref_Phase2TrackerDigi_ > > > trackToken_
unsigned int tidRing(const DetId &id) const
std::vector< EtaPhiBin > clusters
Definition: L1Clustering.h:27
std::vector< TkJet > TkJetCollection
Definition: TkJetFwd.h:20
static constexpr auto TID
def move(src, dest)
Definition: eostools.py:511
int znum
Definition: L1Clustering.h:24

◆ trackQualityCuts()

bool L1TrackJetProducer::trackQualityCuts ( int  trk_nstub,
float  trk_chi2,
float  trk_bendchi2 
)

Definition at line 337 of file L1TrackJetProducer.cc.

References displaced_, nStubs4DisplacedBend_, nStubs4DisplacedChi2_, nStubs4PromptBend_, nStubs4PromptChi2_, nStubs5DisplacedBend_, nStubs5DisplacedChi2_, nStubs5PromptBend_, and nStubs5PromptChi2_.

Referenced by produce().

337  {
338  bool PassQuality = false;
339  if (!displaced_) {
340  if (trk_nstub == 4 && trk_bendchi2 < nStubs4PromptBend_ &&
341  trk_chi2 < nStubs4PromptChi2_) // 4 stubs are the lowest track quality and have different cuts
342  PassQuality = true;
343  if (trk_nstub > 4 && trk_bendchi2 < nStubs5PromptBend_ &&
344  trk_chi2 < nStubs5PromptChi2_) // above 4 stubs diffent selection imposed (genrally looser)
345  PassQuality = true;
346  } else {
347  if (trk_nstub == 4 && trk_bendchi2 < nStubs4DisplacedBend_ &&
348  trk_chi2 < nStubs4DisplacedChi2_) // 4 stubs are the lowest track quality and have different cuts
349  PassQuality = true;
350  if (trk_nstub > 4 && trk_bendchi2 < nStubs5DisplacedBend_ &&
351  trk_chi2 < nStubs5DisplacedChi2_) // above 4 stubs diffent selection imposed (genrally looser)
352  PassQuality = true;
353  }
354  return PassQuality;
355 }
const float nStubs4DisplacedBend_
const float nStubs4DisplacedChi2_
const float nStubs5PromptChi2_
const float nStubs4PromptChi2_
const float nStubs5DisplacedBend_
const float nStubs4PromptBend_
const float nStubs5PromptBend_
const float nStubs5DisplacedChi2_

Member Data Documentation

◆ d0CutNStubs4_

const float L1TrackJetProducer::d0CutNStubs4_
private

Definition at line 89 of file L1TrackJetProducer.cc.

Referenced by produce().

◆ d0CutNStubs5_

const float L1TrackJetProducer::d0CutNStubs5_
private

Definition at line 90 of file L1TrackJetProducer.cc.

Referenced by produce().

◆ displaced_

const bool L1TrackJetProducer::displaced_
private

Definition at line 88 of file L1TrackJetProducer.cc.

Referenced by L1TrackJetProducer(), produce(), and trackQualityCuts().

◆ dzPVTrk_

const float L1TrackJetProducer::dzPVTrk_
private

Definition at line 96 of file L1TrackJetProducer.cc.

Referenced by produce().

◆ etaBins_

const int L1TrackJetProducer::etaBins_
private

Definition at line 82 of file L1TrackJetProducer.cc.

Referenced by L1TrackJetProducer(), and produce().

◆ etaStep_

float L1TrackJetProducer::etaStep_
private

Definition at line 86 of file L1TrackJetProducer.cc.

Referenced by L1TrackJetProducer(), and produce().

◆ highpTJetMinTrackMultiplicity_

const int L1TrackJetProducer::highpTJetMinTrackMultiplicity_
private

Definition at line 79 of file L1TrackJetProducer.cc.

Referenced by produce().

◆ highpTJetThreshold_

const float L1TrackJetProducer::highpTJetThreshold_
private

Definition at line 80 of file L1TrackJetProducer.cc.

Referenced by produce().

◆ L1TrkPtrs_

vector<Ptr<L1TTTrackType> > L1TrackJetProducer::L1TrkPtrs_
private

Definition at line 66 of file L1TrackJetProducer.cc.

Referenced by produce().

◆ lowpTJetMinTrackMultiplicity_

const int L1TrackJetProducer::lowpTJetMinTrackMultiplicity_
private

Definition at line 77 of file L1TrackJetProducer.cc.

Referenced by produce().

◆ lowpTJetThreshold_

const float L1TrackJetProducer::lowpTJetThreshold_
private

Definition at line 78 of file L1TrackJetProducer.cc.

Referenced by produce().

◆ minTrkJetpT_

const double L1TrackJetProducer::minTrkJetpT_
private

Definition at line 84 of file L1TrackJetProducer.cc.

Referenced by produce().

◆ nDisplacedTracks_

const int L1TrackJetProducer::nDisplacedTracks_
private

Definition at line 95 of file L1TrackJetProducer.cc.

Referenced by produce().

◆ nStubs4DisplacedBend_

const float L1TrackJetProducer::nStubs4DisplacedBend_
private

Definition at line 93 of file L1TrackJetProducer.cc.

Referenced by trackQualityCuts().

◆ nStubs4DisplacedChi2_

const float L1TrackJetProducer::nStubs4DisplacedChi2_
private

Definition at line 91 of file L1TrackJetProducer.cc.

Referenced by trackQualityCuts().

◆ nStubs4PromptBend_

const float L1TrackJetProducer::nStubs4PromptBend_
private

Definition at line 74 of file L1TrackJetProducer.cc.

Referenced by trackQualityCuts().

◆ nStubs4PromptChi2_

const float L1TrackJetProducer::nStubs4PromptChi2_
private

Definition at line 72 of file L1TrackJetProducer.cc.

Referenced by trackQualityCuts().

◆ nStubs5DisplacedBend_

const float L1TrackJetProducer::nStubs5DisplacedBend_
private

Definition at line 94 of file L1TrackJetProducer.cc.

Referenced by trackQualityCuts().

◆ nStubs5DisplacedChi2_

const float L1TrackJetProducer::nStubs5DisplacedChi2_
private

Definition at line 92 of file L1TrackJetProducer.cc.

Referenced by trackQualityCuts().

◆ nStubs5PromptBend_

const float L1TrackJetProducer::nStubs5PromptBend_
private

Definition at line 75 of file L1TrackJetProducer.cc.

Referenced by trackQualityCuts().

◆ nStubs5PromptChi2_

const float L1TrackJetProducer::nStubs5PromptChi2_
private

Definition at line 73 of file L1TrackJetProducer.cc.

Referenced by trackQualityCuts().

◆ phiBins_

const int L1TrackJetProducer::phiBins_
private

Definition at line 83 of file L1TrackJetProducer.cc.

Referenced by L1TrackJetProducer(), and produce().

◆ phiStep_

float L1TrackJetProducer::phiStep_
private

Definition at line 87 of file L1TrackJetProducer.cc.

Referenced by L1TrackJetProducer(), and produce().

◆ PVtxToken_

const edm::EDGetTokenT<std::vector<l1t::Vertex> > L1TrackJetProducer::PVtxToken_
private

Definition at line 100 of file L1TrackJetProducer.cc.

Referenced by produce().

◆ tdtrk_

vector<int> L1TrackJetProducer::tdtrk_
private

Definition at line 67 of file L1TrackJetProducer.cc.

Referenced by produce().

◆ trackToken_

const EDGetTokenT<vector<TTTrack<Ref_Phase2TrackerDigi_> > > L1TrackJetProducer::trackToken_
private

Definition at line 99 of file L1TrackJetProducer.cc.

Referenced by produce().

◆ trkEtaMax_

const float L1TrackJetProducer::trkEtaMax_
private

Definition at line 71 of file L1TrackJetProducer.cc.

Referenced by L1TrackJetProducer(), and produce().

◆ trkNPSStubMin_

const int L1TrackJetProducer::trkNPSStubMin_
private

Definition at line 76 of file L1TrackJetProducer.cc.

Referenced by produce().

◆ trkPtMax_

const float L1TrackJetProducer::trkPtMax_
private

Definition at line 69 of file L1TrackJetProducer.cc.

Referenced by produce().

◆ trkPtMin_

const float L1TrackJetProducer::trkPtMin_
private

Definition at line 70 of file L1TrackJetProducer.cc.

Referenced by produce().

◆ trkZMax_

const float L1TrackJetProducer::trkZMax_
private

Definition at line 68 of file L1TrackJetProducer.cc.

Referenced by L1TrackJetProducer(), and produce().

◆ tTopoToken_

edm::ESGetToken<TrackerTopology, TrackerTopologyRcd> L1TrackJetProducer::tTopoToken_
private

Definition at line 98 of file L1TrackJetProducer.cc.

Referenced by produce().

◆ zBins_

const int L1TrackJetProducer::zBins_
private

Definition at line 81 of file L1TrackJetProducer.cc.

Referenced by L1TrackJetProducer(), and produce().

◆ zStep_

float L1TrackJetProducer::zStep_
private

Definition at line 85 of file L1TrackJetProducer.cc.

Referenced by L1TrackJetProducer(), and produce().