Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-07-31 22:58:34

0001 
0002 /**\class L1TrackerHTMissEmulatorProducer L1TrackerHTMissEmulatorProducer.cc
0003  L1Trigger/L1TTrackMatch/plugins/L1TrackerHTMissEmulatorProducer.cc
0004  Description: Takes L1TTkJets and performs a integer emulation of Track-based missing HT, outputting a collection of EtSum 
0005 */
0006 
0007 // Original Author:  Hardik Routray
0008 //         Created:  Mon, 11 Oct 2021
0009 // Update: George Karathanasis, CU Boulder
0010 //         2/4/2024
0011 
0012 // system include files
0013 #include <memory>
0014 #include <numeric>
0015 
0016 // user include files
0017 #include "FWCore/Framework/interface/Frameworkfwd.h"
0018 #include "FWCore/Framework/interface/stream/EDProducer.h"
0019 #include "FWCore/Framework/interface/Event.h"
0020 #include "FWCore/Framework/interface/MakerMacros.h"
0021 #include "FWCore/Framework/interface/EventSetup.h"
0022 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0023 #include "DataFormats/Math/interface/LorentzVector.h"
0024 #include "DataFormats/L1TCorrelator/interface/TkHTMiss.h"
0025 #include "DataFormats/L1TCorrelator/interface/TkHTMissFwd.h"
0026 #include "DataFormats/L1Trigger/interface/EtSum.h"
0027 #include "DataFormats/L1Trigger/interface/TkJetWord.h"
0028 #include "L1Trigger/L1TTrackMatch/interface/L1TkHTMissEmulatorProducer.h"
0029 
0030 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0031 
0032 using namespace l1t;
0033 
0034 class L1TkHTMissEmulatorProducer : public edm::stream::EDProducer<> {
0035 public:
0036   explicit L1TkHTMissEmulatorProducer(const edm::ParameterSet&);
0037   ~L1TkHTMissEmulatorProducer() override = default;
0038 
0039 private:
0040   virtual void beginJob();
0041   void produce(edm::Event&, const edm::EventSetup&) override;
0042   virtual void endJob();
0043 
0044   // ----------member data ---------------------------
0045 
0046   bool debug_ = false;
0047   bool displaced_;
0048 
0049   int cosLUTbins;
0050   int aSteps = 8;
0051 
0052   l1tmhtemu::pt_t jetMinPt_;
0053   l1tmhtemu::eta_t jetMaxEta_;
0054   l1tmhtemu::ntracks_t minNtracksHighPt_;
0055   l1tmhtemu::ntracks_t minNtracksLowPt_;
0056 
0057   std::vector<l1tmhtemu::phi_t> cosLUT_;
0058   std::vector<l1tmhtemu::MHTphi_t> atanLUT_;
0059   std::vector<l1tmhtemu::Et_t> magNormalisationLUT_;
0060 
0061   std::string L1MHTCollectionName_;
0062 
0063   const edm::EDGetTokenT<TkJetWordCollection> jetToken_;
0064 };
0065 
0066 L1TkHTMissEmulatorProducer::L1TkHTMissEmulatorProducer(const edm::ParameterSet& iConfig)
0067     : jetToken_(consumes<TkJetWordCollection>(iConfig.getParameter<edm::InputTag>("L1TkJetEmulationInputTag"))) {
0068   debug_ = iConfig.getParameter<bool>("debug");
0069   displaced_ = iConfig.getParameter<bool>("displaced");
0070 
0071   jetMinPt_ = l1tmhtemu::digitizeSignedValue<l1tmhtemu::pt_t>(
0072       (float)iConfig.getParameter<double>("jet_minPt"), l1tmhtemu::kInternalPtWidth, l1tmhtemu::kStepPt);
0073   jetMaxEta_ = l1tmhtemu::digitizeSignedValue<l1tmhtemu::eta_t>(
0074       (float)iConfig.getParameter<double>("jet_maxEta"), l1tmhtemu::kInternalEtaWidth, l1tmhtemu::kStepEta);
0075   minNtracksHighPt_ = (l1tmhtemu::ntracks_t)iConfig.getParameter<int>("jet_minNtracksHighPt");
0076   minNtracksLowPt_ = (l1tmhtemu::ntracks_t)iConfig.getParameter<int>("jet_minNtracksLowPt");
0077 
0078   cosLUTbins = floor(l1tmhtemu::kMaxCosLUTPhi / l1tmhtemu::kStepPhi);
0079   cosLUT_ = l1tmhtemu::generateCosLUT(cosLUTbins);
0080 
0081   atanLUT_ = l1tmhtemu::generateaTanLUT(aSteps);
0082   magNormalisationLUT_ = l1tmhtemu::generatemagNormalisationLUT(aSteps);
0083 
0084   // Name of output ED Product
0085   L1MHTCollectionName_ = (std::string)iConfig.getParameter<std::string>("L1MHTCollectionName");
0086 
0087   produces<std::vector<EtSum>>(L1MHTCollectionName_);
0088 
0089   if (debug_) {
0090     edm::LogVerbatim("L1TrackerHTMissEmulatorProducer")
0091         << "-------------------------------------------------------------------------\n"
0092         << "====BITWIDTHS====\n"
0093         << "pt: " << l1t::TkJetWord::TkJetBitWidths::kPtSize << " eta: " << l1t::TkJetWord::TkJetBitWidths::kGlbEtaSize
0094         << " phi:" << l1t::TkJetWord::TkJetBitWidths::kGlbPhiSize << "\n"
0095         << "====CUT AP_INTS====\n"
0096         << "minpt: " << jetMinPt_ << " maxeta: " << jetMaxEta_ << " minNtracksHighPt: " << minNtracksHighPt_
0097         << " minNtracksLowPt: " << minNtracksLowPt_ << "\n"
0098         << "====CUT AP_INTS TO FLOATS====\n"
0099         << "minpt: " << (float)jetMinPt_ * l1tmhtemu::kStepPt << " maxeta: " << (float)jetMaxEta_ * l1tmhtemu::kStepEta
0100         << " minNtracksHighPt: " << (int)minNtracksHighPt_ << " minNtracksLowPt: " << (int)minNtracksLowPt_ << "\n"
0101         << "-------------------------------------------------------------------------\n";
0102   }
0103 }
0104 
0105 void L1TkHTMissEmulatorProducer::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0106   using namespace edm;
0107   std::unique_ptr<std::vector<l1t::EtSum>> MHTCollection(new std::vector<l1t::EtSum>(0));
0108 
0109   // L1 track-trigger jets
0110   edm::Handle<TkJetWordCollection> L1TkJetsHandle;
0111   iEvent.getByToken(jetToken_, L1TkJetsHandle);
0112   std::vector<TkJetWord>::const_iterator jetIter;
0113 
0114   if (!L1TkJetsHandle.isValid() && !displaced_) {
0115     LogError("TkHTMissEmulatorProducer") << "\nWarning: TkJetCollection not found in the event. Exit\n";
0116     return;
0117   }
0118 
0119   if (!L1TkJetsHandle.isValid() && displaced_) {
0120     LogError("TkHTMissEmulatorProducer") << "\nWarning: TkJetExtendedCollection not found in the event. Exit\n";
0121     return;
0122   }
0123 
0124   // floats used for debugging
0125   float sumPx_ = 0;
0126   float sumPy_ = 0;
0127   float HT_ = 0;
0128 
0129   l1tmhtemu::Et_t sumPx = 0;
0130   l1tmhtemu::Et_t sumPy = 0;
0131   l1tmhtemu::MHT_t HT = 0;
0132 
0133   // loop over jets
0134   int jetn = 0;
0135 
0136   for (jetIter = L1TkJetsHandle->begin(); jetIter != L1TkJetsHandle->end(); ++jetIter) {
0137     // floats used for debugging
0138     float tmp_jet_px_ = jetIter->pt() * cos(jetIter->glbphi());
0139     float tmp_jet_py_ = jetIter->pt() * sin(jetIter->glbphi());
0140     //float tmp_jet_et_ = jetIter->pt();  // FIXME Get Et from the emulated jets
0141     float tmp_jet_pt_ = jetIter->pt();
0142 
0143     bool tmp_jet_isDisplaced_ = jetIter->dispflag();
0144 
0145     l1tmhtemu::pt_t tmp_jet_pt =
0146         l1tmhtemu::digitizeSignedValue<l1tmhtemu::pt_t>(jetIter->pt(), l1tmhtemu::kInternalPtWidth, l1tmhtemu::kStepPt);
0147     l1tmhtemu::eta_t tmp_jet_eta = l1tmhtemu::digitizeSignedValue<l1tmhtemu::eta_t>(
0148         jetIter->glbeta(), l1tmhtemu::kInternalEtaWidth, l1tmhtemu::kStepEta);
0149     l1tmhtemu::phi_t tmp_jet_phi = l1tmhtemu::digitizeSignedValue<l1tmhtemu::phi_t>(
0150         jetIter->glbphi(), l1tmhtemu::kInternalPhiWidth, l1tmhtemu::kStepPhi);
0151     l1tmhtemu::ntracks_t tmp_jet_nt = l1tmhtemu::ntracks_t(jetIter->nt());
0152 
0153     l1tmhtemu::phi_t tmp_jet_cos_phi = l1tmhtemu::phi_t(-999);
0154     l1tmhtemu::phi_t tmp_jet_sin_phi = l1tmhtemu::phi_t(-999);
0155 
0156     if (tmp_jet_phi >= 0) {
0157       tmp_jet_cos_phi = cosLUT_[tmp_jet_phi];
0158 
0159       if (cosLUTbins / 2 + 1 - tmp_jet_phi >= 0)
0160         tmp_jet_sin_phi = cosLUT_[cosLUTbins / 2 + 1 - tmp_jet_phi];
0161       else
0162         tmp_jet_sin_phi = cosLUT_[-1 * (cosLUTbins / 2 + 1 - tmp_jet_phi)];
0163 
0164     } else {
0165       tmp_jet_cos_phi = cosLUT_[-1 * tmp_jet_phi];
0166 
0167       if (cosLUTbins / 2 + 1 - (-1 * tmp_jet_phi) >= 0)
0168         tmp_jet_sin_phi = -1 * cosLUT_[cosLUTbins / 2 + 1 - (-1 * tmp_jet_phi)];
0169       else
0170         tmp_jet_sin_phi = -1 * cosLUT_[-1 * (cosLUTbins / 2 + 1 - (-1 * tmp_jet_phi))];
0171     }
0172 
0173     l1tmhtemu::Et_t tmp_jet_px = tmp_jet_pt * tmp_jet_cos_phi;
0174     l1tmhtemu::Et_t tmp_jet_py = tmp_jet_pt * tmp_jet_sin_phi;
0175 
0176     jetn++;
0177 
0178     if (debug_) {
0179       edm::LogVerbatim("L1TrackerHTMissEmulatorProducer")
0180           << "****JET EMULATION" << jetn << "****\n"
0181           << "FLOATS ORIGINAL\n"
0182           << "PT: " << jetIter->pt() << "| ETA: " << jetIter->glbeta() << "| PHI: " << jetIter->glbphi()
0183           << "| NTRACKS: " << jetIter->nt() << "| COS(PHI): " << cos(jetIter->glbphi())
0184           << "| SIN(PHI): " << sin(jetIter->glbphi()) << "| Px: " << jetIter->pt() * cos(jetIter->glbphi())
0185           << "| Py: " << jetIter->pt() * sin(jetIter->glbphi()) << "\n"
0186           << "AP_INTS RAW\n"
0187           << "PT: " << jetIter->ptWord() << "| ETA: " << jetIter->glbEtaWord() << "| PHI: " << jetIter->glbPhiWord()
0188           << "| NTRACKS: " << jetIter->ntWord() << "\n"
0189           << "AP_INTS NEW\n"
0190           << "PT: " << tmp_jet_pt << "| ETA: " << tmp_jet_eta << "| PHI: " << tmp_jet_phi << "| NTRACKS: " << tmp_jet_nt
0191           << "| COS(PHI): " << tmp_jet_cos_phi << "| SIN(PHI): " << tmp_jet_sin_phi << "| Px: " << tmp_jet_px
0192           << "| Py: " << tmp_jet_py << "\n"
0193           << "AP_INTS NEW TO FLOATS\n"
0194           << "PT: " << (float)tmp_jet_pt * l1tmhtemu::kStepPt << "| ETA: " << (float)tmp_jet_eta * l1tmhtemu::kStepEta
0195           << "| PHI: " << (float)tmp_jet_phi * l1tmhtemu::kStepPhi << "| NTRACKS: " << (int)tmp_jet_nt
0196           << "| COS(PHI): " << (float)tmp_jet_cos_phi * l1tmhtemu::kStepPhi
0197           << "| SIN(PHI): " << (float)tmp_jet_sin_phi * l1tmhtemu::kStepPhi
0198           << "| Px: " << (float)tmp_jet_px * l1tmhtemu::kStepPt * l1tmhtemu::kStepPhi
0199           << "| Py: " << (float)tmp_jet_py * l1tmhtemu::kStepPt * l1tmhtemu::kStepPhi << "\n"
0200           << "-------------------------------------------------------------------------\n";
0201     }
0202 
0203     if (tmp_jet_pt < jetMinPt_)
0204       continue;
0205     if (tmp_jet_eta > jetMaxEta_ or tmp_jet_eta < -1 * jetMaxEta_)
0206       continue;
0207     if (tmp_jet_nt < minNtracksLowPt_ && tmp_jet_pt > 200)
0208       continue;
0209     if (tmp_jet_nt < minNtracksHighPt_ && tmp_jet_pt > 400)
0210       continue;
0211     if (displaced_ && !tmp_jet_isDisplaced_)
0212       continue;
0213 
0214     if (debug_) {
0215       sumPx_ += tmp_jet_px_;
0216       sumPy_ += tmp_jet_py_;
0217       HT_ += tmp_jet_pt_;
0218     }
0219 
0220     sumPx += tmp_jet_pt * tmp_jet_cos_phi;
0221     sumPy += tmp_jet_pt * tmp_jet_sin_phi;
0222     HT += tmp_jet_pt;
0223 
0224   }  // end jet loop
0225 
0226   // define missing HT
0227 
0228   // Perform cordic sqrt, take x,y and converts to polar coordinate r,phi where
0229   // r=sqrt(x**2+y**2) and phi = atan(y/x)
0230   l1tmhtemu::EtMiss EtMiss = l1tmhtemu::cordicSqrt(sumPx, sumPy, aSteps, atanLUT_, magNormalisationLUT_);
0231   math::XYZTLorentzVector missingEt(-sumPx, -sumPy, 0, EtMiss.Et);
0232 
0233   l1tmhtemu::MHTphi_t phi = 0;
0234 
0235   if ((sumPx < 0) && (sumPy < 0))
0236     phi = EtMiss.Phi - l1tmhtemu::kMHTPhiBins / 2;
0237   else if ((sumPx >= 0) && (sumPy >= 0))
0238     phi = (EtMiss.Phi) + l1tmhtemu::kMHTPhiBins / 2;
0239   else if ((sumPx >= 0) && (sumPy < 0))
0240     phi = EtMiss.Phi - l1tmhtemu::kMHTPhiBins / 2;
0241   else if ((sumPx < 0) && (sumPy >= 0))
0242     phi = EtMiss.Phi - 3 * l1tmhtemu::kMHTPhiBins / 2;
0243 
0244   if (debug_) {
0245     edm::LogVerbatim("L1TrackerHTMissEmulatorProducer")
0246         << "-------------------------------------------------------------------------\n"
0247         << "====MHT FLOATS====\n"
0248         << "sumPx: " << sumPx_ << "| sumPy: " << sumPy_ << "| ET: " << sqrt(sumPx_ * sumPx_ + sumPy_ * sumPy_)
0249         << "| HT: " << HT_ << "| PHI: " << atan2(sumPy_, sumPx_) << "\n"
0250         << "====MHT AP_INTS====\n"
0251         << "sumPx: " << sumPx << "| sumPy: " << sumPy << "| ET: " << EtMiss.Et << "| HT: " << HT << "| PHI: " << phi
0252         << "\n"
0253         << "====MHT AP_INTS TO FLOATS====\n"
0254         << "sumPx: " << (float)sumPx * l1tmhtemu::kStepPt * l1tmhtemu::kStepPhi
0255         << "| sumPy: " << (float)sumPy * l1tmhtemu::kStepPt * l1tmhtemu::kStepPhi << "| ET: " << EtMiss.Et.to_double()
0256         << "| HT: " << (float)HT * l1tmhtemu::kStepPt << "| PHI: " << (float)phi * l1tmhtemu::kStepMHTPhi - M_PI << "\n"
0257         << "-------------------------------------------------------------------------\n";
0258   }
0259   //rescale HT to correct output range
0260   HT = HT / (int)(1 / l1tmhtemu::kStepPt);
0261 
0262   EtSum L1HTSum(missingEt, EtSum::EtSumType::kMissingHt, (int)HT.range(), 0, (int)phi, (int)jetn);
0263 
0264   MHTCollection->push_back(L1HTSum);
0265   iEvent.put(std::move(MHTCollection), L1MHTCollectionName_);
0266 
0267 }  //end producer
0268 
0269 void L1TkHTMissEmulatorProducer::beginJob() {}
0270 
0271 void L1TkHTMissEmulatorProducer::endJob() {}
0272 
0273 DEFINE_FWK_MODULE(L1TkHTMissEmulatorProducer);