Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-07-18 08:35:34

0001 // Original Author: Samuel Edwin Leigh, Tyler Wu,
0002 //                  Rutgers, the State University of New Jersey
0003 //
0004 // Rewritting/Improvements:      George Karathanasis,
0005 //                          georgios.karathanasis@cern.ch, CU Boulder
0006 //
0007 //         Created:  Wed, 01 Aug 2018 14:01:41 GMT
0008 //         Latest update: Nov 2022 (by GK)
0009 //
0010 // Track jets are clustered in a two-layer process, first by clustering in phi,
0011 // then by clustering in eta. The code proceeds as following: putting all tracks// in a grid of eta vs phi space, and then cluster them. Finally we merge the cl
0012 // usters when needed. The code is improved to use the same module between emula
0013 // tion and simulation was also improved, with bug fixes and being faster.
0014 // Introduction to object (p10-13):
0015 // https://indico.cern.ch/event/791517/contributions/3341650/attachments/1818736/2973771/TrackBasedAlgos_L1TMadrid_MacDonald.pdf
0016 // New and improved version: https://indico.cern.ch/event/1203796/contributions/5073056/attachments/2519806/4333006/trackjet_emu.pdf
0017 
0018 // L1T include files
0019 #include "DataFormats/L1TCorrelator/interface/TkJet.h"
0020 #include "DataFormats/L1TCorrelator/interface/TkJetFwd.h"
0021 #include "DataFormats/L1TrackTrigger/interface/TTTypes.h"
0022 #include "DataFormats/L1TrackTrigger/interface/TTTrack.h"
0023 #include "DataFormats/L1TrackTrigger/interface/TTTrack_TrackWord.h"
0024 #include "DataFormats/L1Trigger/interface/TkJetWord.h"
0025 #include "DataFormats/L1Trigger/interface/VertexWord.h"
0026 
0027 // system include files
0028 #include "DataFormats/Common/interface/Handle.h"
0029 #include "FWCore/Framework/interface/Frameworkfwd.h"
0030 #include "FWCore/Framework/interface/stream/EDProducer.h"
0031 #include "FWCore/Framework/interface/Event.h"
0032 #include "FWCore/Framework/interface/EventSetup.h"
0033 #include "FWCore/Framework/interface/MakerMacros.h"
0034 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0035 #include "FWCore/Utilities/interface/StreamID.h"
0036 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0037 #include "DataFormats/Common/interface/Ref.h"
0038 #include "DataFormats/Common/interface/RefVector.h"
0039 
0040 //own headers
0041 #include "L1TrackJetClustering.h"
0042 
0043 //general
0044 #include <ap_int.h>
0045 
0046 using namespace std;
0047 using namespace edm;
0048 using namespace l1t;
0049 using namespace l1ttrackjet;
0050 
0051 class L1TrackJetEmulatorProducer : public stream::EDProducer<> {
0052 public:
0053   explicit L1TrackJetEmulatorProducer(const ParameterSet &);
0054   ~L1TrackJetEmulatorProducer() override = default;
0055   typedef TTTrack<Ref_Phase2TrackerDigi_> L1TTTrackType;
0056   typedef vector<L1TTTrackType> L1TTTrackCollectionType;
0057   typedef edm::RefVector<L1TTTrackCollectionType> L1TTTrackRefCollectionType;
0058   static void fillDescriptions(ConfigurationDescriptions &descriptions);
0059 
0060 private:
0061   void produce(Event &, const EventSetup &) override;
0062 
0063   // ----------member data ---------------------------
0064 
0065   std::vector<edm::Ptr<L1TTTrackType>> L1TrkPtrs_;
0066   vector<int> tdtrk_;
0067   const float trkZMax_;
0068   const float trkPtMax_;
0069   const float trkPtMin_;
0070   const float trkEtaMax_;
0071   const float nStubs4PromptChi2_;
0072   const float nStubs5PromptChi2_;
0073   const float nStubs4PromptBend_;
0074   const float nStubs5PromptBend_;
0075   const int trkNPSStubMin_;
0076   const int lowpTJetMinTrackMultiplicity_;
0077   const float lowpTJetThreshold_;
0078   const int highpTJetMinTrackMultiplicity_;
0079   const float highpTJetThreshold_;
0080   const int zBins_;
0081   const int etaBins_;
0082   const int phiBins_;
0083   const double minTrkJetpT_;
0084   const bool displaced_;
0085   const float d0CutNStubs4_;
0086   const float d0CutNStubs5_;
0087   const float nStubs4DisplacedChi2_;
0088   const float nStubs5DisplacedChi2_;
0089   const float nStubs4DisplacedBend_;
0090   const float nStubs5DisplacedBend_;
0091   const int nDisplacedTracks_;
0092   const float dzPVTrk_;
0093 
0094   float PVz;
0095   float zStep_;
0096   glbeta_intern etaStep_;
0097   glbphi_intern phiStep_;
0098 
0099   TTTrack_TrackWord trackword;
0100 
0101   edm::ESGetToken<TrackerTopology, TrackerTopologyRcd> tTopoToken_;
0102   const EDGetTokenT<L1TTTrackRefCollectionType> trackToken_;
0103   const EDGetTokenT<l1t::VertexWordCollection> PVtxToken_;
0104 };
0105 
0106 //constructor
0107 L1TrackJetEmulatorProducer::L1TrackJetEmulatorProducer(const ParameterSet &iConfig)
0108     : trkZMax_(iConfig.getParameter<double>("trk_zMax")),
0109       trkPtMax_(iConfig.getParameter<double>("trk_ptMax")),
0110       trkPtMin_(iConfig.getParameter<double>("trk_ptMin")),
0111       trkEtaMax_(iConfig.getParameter<double>("trk_etaMax")),
0112       nStubs4PromptChi2_(iConfig.getParameter<double>("nStubs4PromptChi2")),
0113       nStubs5PromptChi2_(iConfig.getParameter<double>("nStubs5PromptChi2")),
0114       nStubs4PromptBend_(iConfig.getParameter<double>("nStubs4PromptBend")),
0115       nStubs5PromptBend_(iConfig.getParameter<double>("nStubs5PromptBend")),
0116       trkNPSStubMin_(iConfig.getParameter<int>("trk_nPSStubMin")),
0117       lowpTJetMinTrackMultiplicity_(iConfig.getParameter<int>("lowpTJetMinTrackMultiplicity")),
0118       lowpTJetThreshold_(iConfig.getParameter<double>("lowpTJetThreshold")),
0119       highpTJetMinTrackMultiplicity_(iConfig.getParameter<int>("highpTJetMinTrackMultiplicity")),
0120       highpTJetThreshold_(iConfig.getParameter<double>("highpTJetThreshold")),
0121       zBins_(iConfig.getParameter<int>("zBins")),
0122       etaBins_(iConfig.getParameter<int>("etaBins")),
0123       phiBins_(iConfig.getParameter<int>("phiBins")),
0124       minTrkJetpT_(iConfig.getParameter<double>("minTrkJetpT")),
0125       displaced_(iConfig.getParameter<bool>("displaced")),
0126       d0CutNStubs4_(iConfig.getParameter<double>("d0_cutNStubs4")),
0127       d0CutNStubs5_(iConfig.getParameter<double>("d0_cutNStubs5")),
0128       nStubs4DisplacedChi2_(iConfig.getParameter<double>("nStubs4DisplacedChi2")),
0129       nStubs5DisplacedChi2_(iConfig.getParameter<double>("nStubs5DisplacedChi2")),
0130       nStubs4DisplacedBend_(iConfig.getParameter<double>("nStubs4DisplacedBend")),
0131       nStubs5DisplacedBend_(iConfig.getParameter<double>("nStubs5DisplacedBend")),
0132       nDisplacedTracks_(iConfig.getParameter<int>("nDisplacedTracks")),
0133       dzPVTrk_(iConfig.getParameter<double>("MaxDzTrackPV")),
0134       tTopoToken_(esConsumes<TrackerTopology, TrackerTopologyRcd>(edm::ESInputTag("", ""))),
0135       trackToken_(consumes<L1TTTrackRefCollectionType>(iConfig.getParameter<InputTag>("L1TrackInputTag"))),
0136       PVtxToken_(consumes<l1t::VertexWordCollection>(iConfig.getParameter<InputTag>("L1PVertexInputTag"))) {
0137   zStep_ = 2.0 * trkZMax_ / (zBins_ + 1);                 // added +1 in denom
0138   etaStep_ = glbeta_intern(2.0 * trkEtaMax_ / etaBins_);  //etaStep is the width of an etabin
0139   phiStep_ = DoubleToBit(2.0 * (M_PI) / phiBins_,
0140                          TTTrack_TrackWord::TrackBitWidths::kPhiSize + kExtraGlobalPhiBit,
0141                          TTTrack_TrackWord::stepPhi0);  ///phiStep is the width of a phibin
0142 
0143   if (displaced_)
0144     produces<l1t::TkJetWordCollection>("L1TrackJetsExtended");
0145   else
0146     produces<l1t::TkJetWordCollection>("L1TrackJets");
0147 }
0148 
0149 void L1TrackJetEmulatorProducer::produce(Event &iEvent, const EventSetup &iSetup) {
0150   unique_ptr<l1t::TkJetWordCollection> L1TrackJetContainer(new l1t::TkJetWordCollection);
0151   // Read inputs
0152   const TrackerTopology &tTopo = iSetup.getData(tTopoToken_);
0153 
0154   edm::Handle<L1TTTrackRefCollectionType> TTTrackHandle;
0155   iEvent.getByToken(trackToken_, TTTrackHandle);
0156 
0157   edm::Handle<l1t::VertexWordCollection> PVtx;
0158   iEvent.getByToken(PVtxToken_, PVtx);
0159   float PVz = (PVtx->at(0)).z0();
0160 
0161   L1TrkPtrs_.clear();
0162   tdtrk_.clear();
0163   // track selection
0164   for (unsigned int this_l1track = 0; this_l1track < TTTrackHandle->size(); this_l1track++) {
0165     edm::Ptr<L1TTTrackType> trkPtr(TTTrackHandle, this_l1track);
0166     float trk_pt = trkPtr->momentum().perp();
0167     int trk_nstubs = (int)trkPtr->getStubRefs().size();
0168     float trk_chi2dof = trkPtr->chi2Red();
0169     float trk_bendchi2 = trkPtr->stubPtConsistency();
0170     int trk_nPS = 0;
0171     for (int istub = 0; istub < trk_nstubs; istub++) {
0172       DetId detId(trkPtr->getStubRefs().at(istub)->getDetId());
0173       if (detId.det() == DetId::Detector::Tracker) {
0174         if ((detId.subdetId() == StripSubdetector::TOB && tTopo.tobLayer(detId) <= 3) ||
0175             (detId.subdetId() == StripSubdetector::TID && tTopo.tidRing(detId) <= 9))
0176           trk_nPS++;
0177       }
0178     }
0179     // selection tracks - supposed to happen on seperate module (kept for legacy/debug reasons)
0180     if (trk_nPS < trkNPSStubMin_)
0181       continue;
0182     if (!TrackQualitySelection(trk_nstubs,
0183                                trk_chi2dof,
0184                                trk_bendchi2,
0185                                nStubs4PromptBend_,
0186                                nStubs5PromptBend_,
0187                                nStubs4PromptChi2_,
0188                                nStubs5PromptChi2_,
0189                                nStubs4DisplacedBend_,
0190                                nStubs5DisplacedBend_,
0191                                nStubs4DisplacedChi2_,
0192                                nStubs5DisplacedChi2_,
0193                                displaced_))
0194       continue;
0195     if (std::abs(PVz - trkPtr->z0()) > dzPVTrk_ && dzPVTrk_ > 0)
0196       continue;
0197     if (std::abs(trkPtr->z0()) > trkZMax_)
0198       continue;
0199     if (std::abs(trkPtr->momentum().eta()) > trkEtaMax_)
0200       continue;
0201     if (trk_pt < trkPtMin_)
0202       continue;
0203     L1TrkPtrs_.push_back(trkPtr);
0204   }
0205 
0206   // if no tracks pass selection return empty containers
0207   if (L1TrkPtrs_.empty()) {
0208     if (displaced_)
0209       iEvent.put(std::move(L1TrackJetContainer), "L1TrackJetsExtended");
0210     else
0211       iEvent.put(std::move(L1TrackJetContainer), "L1TrackJets");
0212     return;
0213   }
0214 
0215   TrackJetEmulationMaxZBin mzb;
0216   mzb.znum = 0;
0217   mzb.nclust = 0;
0218   mzb.ht = 0;
0219 
0220   TrackJetEmulationEtaPhiBin epbins_default[phiBins_][etaBins_];  // create grid of phiBins
0221   glbphi_intern phi = DoubleToBit(
0222       -1.0 * M_PI, TTTrack_TrackWord::TrackBitWidths::kPhiSize + kExtraGlobalPhiBit, TTTrack_TrackWord::stepPhi0);
0223   for (int i = 0; i < phiBins_; ++i) {
0224     glbeta_intern eta = -1 * trkEtaMax_;
0225     for (int j = 0; j < etaBins_; ++j) {
0226       epbins_default[i][j].phi = (phi + (phi + phiStep_)) / 2;  // phi coord of bin center
0227       epbins_default[i][j].eta = (eta + (eta + etaStep_)) / 2;  // eta coord of bin center
0228       eta += etaStep_;
0229     }  // for each etabin
0230     phi += phiStep_;
0231   }  // for each phibin (finished creating epbins)
0232 
0233   // bins for z if multibin option is selected
0234   std::vector<z0_intern> zmins, zmaxs;
0235   for (int zbin = 0; zbin < zBins_; zbin++) {
0236     zmins.push_back(DoubleToBit(
0237         -1.0 * trkZMax_ + zStep_ * zbin, TTTrack_TrackWord::TrackBitWidths::kZ0Size, TTTrack_TrackWord::stepZ0));
0238     zmaxs.push_back(DoubleToBit(-1.0 * trkZMax_ + zStep_ * zbin + 2 * zStep_,
0239                                 TTTrack_TrackWord::TrackBitWidths::kZ0Size,
0240                                 TTTrack_TrackWord::stepZ0));
0241   }
0242 
0243   // create vectors that hold clusters
0244   std::vector<std::vector<TrackJetEmulationEtaPhiBin>> L1clusters;
0245   L1clusters.reserve(phiBins_);
0246   std::vector<TrackJetEmulationEtaPhiBin> L2clusters;
0247 
0248   // default (empty) grid
0249   for (int i = 0; i < phiBins_; ++i) {
0250     for (int j = 0; j < etaBins_; ++j) {
0251       epbins_default[i][j].pTtot = 0;
0252       epbins_default[i][j].used = false;
0253       epbins_default[i][j].ntracks = 0;
0254       epbins_default[i][j].nxtracks = 0;
0255       epbins_default[i][j].trackidx.clear();
0256     }
0257   }
0258 
0259   // logic: loop over z bins find tracks in this bin and arrange them in a 2D eta-phi matrix
0260   for (unsigned int zbin = 0; zbin < zmins.size(); ++zbin) {
0261     // initialize matrices for every z bin
0262     z0_intern zmin = zmins[zbin];
0263     z0_intern zmax = zmaxs[zbin];
0264     TrackJetEmulationEtaPhiBin epbins[phiBins_][etaBins_];
0265     std::copy(&epbins_default[0][0], &epbins_default[0][0] + phiBins_ * etaBins_, &epbins[0][0]);
0266 
0267     //clear containers
0268     L1clusters.clear();
0269     L2clusters.clear();
0270 
0271     // fill grid
0272     for (unsigned int k = 0; k < L1TrkPtrs_.size(); ++k) {
0273       //// conversions
0274       //-z0
0275       z0_intern trkZ = L1TrkPtrs_[k]->getZ0Word();
0276 
0277       if (zmax < trkZ)
0278         continue;
0279       if (zmin > trkZ)
0280         continue;
0281       if (zbin == 0 && zmin == trkZ)
0282         continue;
0283 
0284       //-pt
0285       ap_uint<TTTrack_TrackWord::TrackBitWidths::kRinvSize - 1> ptEmulationBits = L1TrkPtrs_[k]->getRinvWord();
0286       pt_intern trkpt;
0287       trkpt.V = ptEmulationBits.range();
0288 
0289       //-eta
0290       TTTrack_TrackWord::tanl_t etaEmulationBits = L1TrkPtrs_[k]->getTanlWord();
0291       glbeta_intern trketa;
0292       trketa.V = etaEmulationBits.range();
0293 
0294       //-phi
0295       int Sector = L1TrkPtrs_[k]->phiSector();
0296       double sector_phi_value = 0;
0297       if (Sector < 5) {
0298         sector_phi_value = 2.0 * M_PI * Sector / 9.0;
0299       } else {
0300         sector_phi_value = (-1.0 * M_PI + M_PI / 9.0 + (Sector - 5) * 2.0 * M_PI / 9.0);
0301       }
0302 
0303       glbphi_intern trkphiSector = DoubleToBit(sector_phi_value,
0304                                                TTTrack_TrackWord::TrackBitWidths::kPhiSize + kExtraGlobalPhiBit,
0305                                                TTTrack_TrackWord::stepPhi0);
0306       glbphi_intern local_phi = 0;
0307       local_phi.V = L1TrkPtrs_[k]->getPhiWord();
0308       glbphi_intern local_phi2 =
0309           DoubleToBit(BitToDouble(local_phi, TTTrack_TrackWord::TrackBitWidths::kPhiSize, TTTrack_TrackWord::stepPhi0),
0310                       TTTrack_TrackWord::TrackBitWidths::kPhiSize + kExtraGlobalPhiBit,
0311                       TTTrack_TrackWord::stepPhi0);
0312       glbphi_intern trkphi = local_phi2 + trkphiSector;
0313 
0314       //-d0
0315       d0_intern abs_trkD0 = L1TrkPtrs_[k]->getD0Word();
0316 
0317       //-nstub
0318       int trk_nstubs = (int)L1TrkPtrs_[k]->getStubRefs().size();
0319 
0320       // now fill the 2D grid with tracks
0321       for (int i = 0; i < phiBins_; ++i) {
0322         for (int j = 0; j < etaBins_; ++j) {
0323           glbeta_intern eta_min = epbins[i][j].eta - etaStep_ / 2;  //eta min
0324           glbeta_intern eta_max = epbins[i][j].eta + etaStep_ / 2;  //eta max
0325           glbphi_intern phi_min = epbins[i][j].phi - phiStep_ / 2;  //phi min
0326           glbphi_intern phi_max = epbins[i][j].phi + phiStep_ / 2;  //phi max
0327           if ((trketa < eta_min) && j != 0)
0328             continue;
0329           if ((trketa > eta_max) && j != etaBins_ - 1)
0330             continue;
0331           if ((trkphi < phi_min) && i != 0)
0332             continue;
0333           if ((trkphi > phi_max) && i != phiBins_ - 1)
0334             continue;
0335 
0336           if (trkpt < pt_intern(trkPtMax_))
0337             epbins[i][j].pTtot += trkpt;
0338           else
0339             epbins[i][j].pTtot += pt_intern(trkPtMax_);
0340           if ((abs_trkD0 >
0341                    DoubleToBit(d0CutNStubs5_, TTTrack_TrackWord::TrackBitWidths::kD0Size, TTTrack_TrackWord::stepD0) &&
0342                trk_nstubs >= 5 && d0CutNStubs5_ >= 0) ||
0343               (abs_trkD0 >
0344                    DoubleToBit(d0CutNStubs4_, TTTrack_TrackWord::TrackBitWidths::kD0Size, TTTrack_TrackWord::stepD0) &&
0345                trk_nstubs == 4 && d0CutNStubs4_ >= 0))
0346             epbins[i][j].nxtracks += 1;
0347 
0348           epbins[i][j].trackidx.push_back(k);
0349           ++epbins[i][j].ntracks;
0350         }  // for each etabin
0351       }    // for each phibin
0352     }      //end loop over tracks
0353 
0354     // first layer clustering - in eta using grid
0355     for (int phibin = 0; phibin < phiBins_; ++phibin) {
0356       L1clusters.push_back(L1_clustering<TrackJetEmulationEtaPhiBin, pt_intern, glbeta_intern, glbphi_intern>(
0357           epbins[phibin], etaBins_, etaStep_));
0358     }
0359 
0360     //second layer clustering in phi for all eta clusters
0361     L2clusters = L2_clustering<TrackJetEmulationEtaPhiBin, pt_intern, glbeta_intern, glbphi_intern>(
0362         L1clusters, phiBins_, phiStep_, etaStep_);
0363 
0364     // sum all cluster pTs in this zbin to find max
0365     pt_intern sum_pt = 0;
0366     for (unsigned int k = 0; k < L2clusters.size(); ++k) {
0367       if (L2clusters[k].pTtot > pt_intern(highpTJetThreshold_) && L2clusters[k].ntracks < lowpTJetMinTrackMultiplicity_)
0368         continue;
0369       if (L2clusters[k].pTtot > pt_intern(highpTJetThreshold_) &&
0370           L2clusters[k].ntracks < highpTJetMinTrackMultiplicity_)
0371         continue;
0372 
0373       if (L2clusters[k].pTtot > pt_intern(minTrkJetpT_))
0374         sum_pt += L2clusters[k].pTtot;
0375     }
0376     if (sum_pt < mzb.ht)
0377       continue;
0378     // if ht is larger than previous max, this is the new vertex zbin
0379     mzb.ht = sum_pt;
0380     mzb.znum = zbin;
0381     mzb.clusters = L2clusters;
0382     mzb.nclust = L2clusters.size();
0383     mzb.zbincenter = (zmin + zmax) / 2.0;
0384   }  //zbin loop
0385 
0386   vector<edm::Ptr<L1TTTrackType>> L1TrackAssocJet;
0387   for (unsigned int j = 0; j < mzb.clusters.size(); ++j) {
0388     if (mzb.clusters[j].pTtot < pt_intern(trkPtMin_))
0389       continue;
0390 
0391     l1t::TkJetWord::glbeta_t jetEta = DoubleToBit(double(mzb.clusters[j].eta),
0392                                                   TkJetWord::TkJetBitWidths::kGlbEtaSize,
0393                                                   TkJetWord::MAX_ETA / (1 << TkJetWord::TkJetBitWidths::kGlbEtaSize));
0394     l1t::TkJetWord::glbphi_t jetPhi = DoubleToBit(
0395         BitToDouble(mzb.clusters[j].phi, TTTrack_TrackWord::TrackBitWidths::kPhiSize + 4, TTTrack_TrackWord::stepPhi0),
0396         TkJetWord::TkJetBitWidths::kGlbPhiSize,
0397         (2. * std::abs(M_PI)) / (1 << TkJetWord::TkJetBitWidths::kGlbPhiSize));
0398     l1t::TkJetWord::z0_t jetZ0 = 0;
0399     l1t::TkJetWord::pt_t jetPt = mzb.clusters[j].pTtot;
0400     l1t::TkJetWord::nt_t total_ntracks = mzb.clusters[j].ntracks;
0401     l1t::TkJetWord::nx_t total_disptracks = mzb.clusters[j].nxtracks;
0402     l1t::TkJetWord::dispflag_t dispflag = 0;
0403     l1t::TkJetWord::tkjetunassigned_t unassigned = 0;
0404 
0405     if (total_disptracks > nDisplacedTracks_ || total_disptracks == nDisplacedTracks_)
0406       dispflag = 1;
0407     L1TrackAssocJet.clear();
0408     for (unsigned int itrk = 0; itrk < mzb.clusters[j].trackidx.size(); itrk++)
0409       L1TrackAssocJet.push_back(L1TrkPtrs_[mzb.clusters[j].trackidx[itrk]]);
0410 
0411     l1t::TkJetWord trkJet(jetPt, jetEta, jetPhi, jetZ0, total_ntracks, total_disptracks, dispflag, unassigned);
0412 
0413     L1TrackJetContainer->push_back(trkJet);
0414   }
0415 
0416   std::sort(L1TrackJetContainer->begin(), L1TrackJetContainer->end(), [](auto &a, auto &b) { return a.pt() > b.pt(); });
0417   if (displaced_)
0418     iEvent.put(std::move(L1TrackJetContainer), "L1TrackJetsExtended");
0419   else
0420     iEvent.put(std::move(L1TrackJetContainer), "L1TrackJets");
0421 }
0422 
0423 void L1TrackJetEmulatorProducer::fillDescriptions(ConfigurationDescriptions &descriptions) {
0424   //The following says we do not know what parameters are allowed so do no validation
0425   // Please change this to state exactly what you do use, even if it is no parameters
0426   ParameterSetDescription desc;
0427   desc.add<edm::InputTag>("L1TrackInputTag", edm::InputTag("l1tTTTracksFromTrackletEmulation", "Level1TTTracks"));
0428   desc.add<edm::InputTag>("L1PVertexInputTag", edm::InputTag("l1tVertexFinderEmulator", "L1VerticesEmulation"));
0429   desc.add<double>("MaxDzTrackPV", 1.0);
0430   desc.add<double>("trk_zMax", 15.0);
0431   desc.add<double>("trk_ptMax", 200.0);
0432   desc.add<double>("trk_ptMin", 3.0);
0433   desc.add<double>("trk_etaMax", 2.4);
0434   desc.add<double>("nStubs4PromptChi2", 5.0);
0435   desc.add<double>("nStubs4PromptBend", 1.7);
0436   desc.add<double>("nStubs5PromptChi2", 2.75);
0437   desc.add<double>("nStubs5PromptBend", 3.5);
0438   desc.add<int>("trk_nPSStubMin", -1);
0439   desc.add<double>("minTrkJetpT", -1.0);
0440   desc.add<int>("etaBins", 24);
0441   desc.add<int>("phiBins", 27);
0442   desc.add<int>("zBins", 1);
0443   desc.add<double>("d0_cutNStubs4", -1);
0444   desc.add<double>("d0_cutNStubs5", -1);
0445   desc.add<int>("lowpTJetMinTrackMultiplicity", 2);
0446   desc.add<double>("lowpTJetThreshold", 50.0);
0447   desc.add<int>("highpTJetMinTrackMultiplicity", 3);
0448   desc.add<double>("highpTJetThreshold", 100.0);
0449   desc.add<bool>("displaced", false);
0450   desc.add<double>("nStubs4DisplacedChi2", 5.0);
0451   desc.add<double>("nStubs4DisplacedBend", 1.7);
0452   desc.add<double>("nStubs5DisplacedChi2", 2.75);
0453   desc.add<double>("nStubs5DisplacedBend", 3.5);
0454   desc.add<int>("nDisplacedTracks", 2);
0455   descriptions.add("l1tTrackJetsEmulator", desc);
0456 }
0457 
0458 //define this as a plug-in
0459 DEFINE_FWK_MODULE(L1TrackJetEmulatorProducer);