Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-06-10 01:53:57

0001 #include "L1Trigger/Phase2L1ParticleFlow/interface/pf/pfalgo2hgc_ref.h"
0002 #include "L1Trigger/Phase2L1ParticleFlow/interface/dbgPrintf.h"
0003 
0004 #include <cmath>
0005 #include <cstdio>
0006 #include <algorithm>
0007 #include <memory>
0008 
0009 #ifdef CMSSW_GIT_HASH
0010 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0011 
0012 l1ct::PFAlgo2HGCEmulator::PFAlgo2HGCEmulator(const edm::ParameterSet& iConfig)
0013     : PFAlgoEmulatorBase(iConfig.getParameter<uint32_t>("nTrack"),
0014                          iConfig.getParameter<uint32_t>("nCalo"),
0015                          iConfig.getParameter<uint32_t>("nMu"),
0016                          iConfig.getParameter<uint32_t>("nSelCalo"),
0017                          l1ct::Scales::makeDR2FromFloatDR(iConfig.getParameter<double>("trackMuDR")),
0018                          l1ct::Scales::makeDR2FromFloatDR(iConfig.getParameter<double>("trackCaloDR")),
0019                          l1ct::Scales::makePtFromFloat(iConfig.getParameter<double>("maxInvisiblePt")),
0020                          l1ct::Scales::makePtFromFloat(iConfig.getParameter<double>("tightTrackMaxInvisiblePt"))) {
0021   debug_ = iConfig.getUntrackedParameter<bool>("debug", false);
0022   loadPtErrBins(iConfig);
0023 }
0024 
0025 #endif
0026 
0027 void l1ct::PFAlgo2HGCEmulator::toFirmware(const PFInputRegion& in,
0028                                           PFRegion& region,
0029                                           HadCaloObj calo[/*nCALO*/],
0030                                           TkObj track[/*nTRACK*/],
0031                                           MuObj mu[/*nMU*/]) const {
0032   region = in.region;
0033   l1ct::toFirmware(in.track, nTRACK_, track);
0034   l1ct::toFirmware(in.hadcalo, nCALO_, calo);
0035   l1ct::toFirmware(in.muon, nMU_, mu);
0036 }
0037 
0038 void l1ct::PFAlgo2HGCEmulator::toFirmware(const OutputRegion& out,
0039                                           PFChargedObj outch[/*nTRACK*/],
0040                                           PFNeutralObj outne[/*nSELCALO*/],
0041                                           PFChargedObj outmu[/*nMU*/]) const {
0042   l1ct::toFirmware(out.pfcharged, nTRACK_, outch);
0043   l1ct::toFirmware(out.pfneutral, nSELCALO_, outne);
0044   l1ct::toFirmware(out.pfmuon, nMU_, outmu);
0045 }
0046 
0047 void l1ct::PFAlgo2HGCEmulator::run(const PFInputRegion& in, OutputRegion& out) const {
0048   unsigned int nTRACK = std::min<unsigned>(nTRACK_, in.track.size());
0049   unsigned int nCALO = std::min<unsigned>(nCALO_, in.hadcalo.size());
0050   unsigned int nSELCALO = std::min<unsigned>(nSELCALO_, in.hadcalo.size());
0051   unsigned int nMU = std::min<unsigned>(nMU_, in.muon.size());
0052 
0053   if (debug_) {
0054     dbgPrintf("FW\nFW  \t region eta %+5.2f [ %+5.2f , %+5.2f ], phi %+5.2f [ %+5.2f , %+5.2f ]   packed %s\n",
0055               in.region.floatEtaCenter(),
0056               in.region.floatEtaMinExtra(),
0057               in.region.floatEtaMaxExtra(),
0058               in.region.floatPhiCenter(),
0059               in.region.floatPhiCenter() - in.region.floatPhiHalfWidthExtra(),
0060               in.region.floatPhiCenter() + in.region.floatPhiHalfWidthExtra(),
0061               in.region.pack().to_string(16).c_str());
0062 
0063     dbgPrintf("FW  \t N(track) %3lu   N(calo) %3lu   N(mu) %3lu\n", in.track.size(), in.hadcalo.size(), in.muon.size());
0064 
0065     for (unsigned int i = 0; i < nTRACK; ++i) {
0066       if (in.track[i].hwPt == 0)
0067         continue;
0068       dbgPrintf(
0069           "FW  \t track %3d: pt %8.2f [ %8d ]  calo eta %+5.2f [ %+5d ]  calo phi %+5.2f [ %+5d ]  vtx eta %+5.2f  "
0070           "vtx phi %+5.2f  charge %+2d  qual %d  fid %d  glb eta %+5.2f phi %+5.2f  packed %s\n",
0071           i,
0072           in.track[i].floatPt(),
0073           in.track[i].intPt(),
0074           in.track[i].floatEta(),
0075           in.track[i].intEta(),
0076           in.track[i].floatPhi(),
0077           in.track[i].intPhi(),
0078           in.track[i].floatVtxEta(),
0079           in.track[i].floatVtxPhi(),
0080           in.track[i].intCharge(),
0081           int(in.track[i].hwQuality),
0082           int(in.region.isFiducial(in.track[i].hwEta, in.track[i].hwPhi)),
0083           in.region.floatGlbEta(in.track[i].hwVtxEta()),
0084           in.region.floatGlbPhi(in.track[i].hwVtxPhi()),
0085           in.track[i].pack().to_string(16).c_str());
0086     }
0087     for (unsigned int i = 0; i < nCALO; ++i) {
0088       if (in.hadcalo[i].hwPt == 0)
0089         continue;
0090       dbgPrintf(
0091           "FW  \t calo  %3d: pt %8.2f [ %8d ]  calo eta %+5.2f [ %+5d ]  calo phi %+5.2f [ %+5d ]  calo emPt %8.2f [ "
0092           "%6d ]   emID %2d  fid %d  glb eta %+5.2f phi %+5.2f  packed %s \n",
0093           i,
0094           in.hadcalo[i].floatPt(),
0095           in.hadcalo[i].intPt(),
0096           in.hadcalo[i].floatEta(),
0097           in.hadcalo[i].intEta(),
0098           in.hadcalo[i].floatPhi(),
0099           in.hadcalo[i].intPhi(),
0100           in.hadcalo[i].floatEmPt(),
0101           in.hadcalo[i].intEmPt(),
0102           in.hadcalo[i].hwEmID.to_int(),
0103           int(in.region.isFiducial(in.hadcalo[i].hwEta, in.hadcalo[i].hwPhi)),
0104           in.region.floatGlbEtaOf(in.hadcalo[i]),
0105           in.region.floatGlbPhiOf(in.hadcalo[i]),
0106           in.hadcalo[i].pack().to_string(16).c_str());
0107     }
0108     for (unsigned int i = 0; i < nMU; ++i) {
0109       if (in.muon[i].hwPt == 0)
0110         continue;
0111       dbgPrintf(
0112           "FW  \t muon  %3d: pt %8.2f [ %8d ]  calo eta %+5.2f [ %+5d ]  calo phi %+5.2f [ %+5d ]  "
0113           "vtx eta %+5.2f  vtx phi %+5.2f  charge %+2d  qual %2d  glb eta %+5.2f phi %+5.2f  packed %s \n",
0114           i,
0115           in.muon[i].floatPt(),
0116           in.muon[i].intPt(),
0117           in.muon[i].floatEta(),
0118           in.muon[i].intEta(),
0119           in.muon[i].floatPhi(),
0120           in.muon[i].intPhi(),
0121           in.muon[i].floatVtxEta(),
0122           in.muon[i].floatVtxPhi(),
0123           in.muon[i].intCharge(),
0124           int(in.muon[i].hwQuality),
0125           in.region.floatGlbEta(in.muon[i].hwVtxEta()),
0126           in.region.floatGlbPhi(in.muon[i].hwVtxPhi()),
0127           in.muon[i].pack().to_string(16).c_str());
0128     }
0129   }
0130 
0131   ////////////////////////////////////////////////////
0132   // TK-MU Linking
0133   std::vector<int> iMu;
0134   pfalgo_mu_ref(in, out, iMu);
0135 
0136   ////////////////////////////////////////////////////
0137   // TK-HAD Linking
0138 
0139   // initialize sum track pt
0140   std::vector<pt_t> calo_sumtk(nCALO), calo_subpt(nCALO);
0141   std::vector<pt2_t> calo_sumtkErr2(nCALO);
0142   for (unsigned int ic = 0; ic < nCALO; ++ic) {
0143     calo_sumtk[ic] = 0;
0144     calo_sumtkErr2[ic] = 0;
0145   }
0146 
0147   // initialize good track bit
0148   std::vector<bool> track_good(nTRACK, false);
0149   std::vector<bool> isEle(nTRACK, false);
0150   for (unsigned int it = 0; it < nTRACK; ++it) {
0151     if (!in.track[it].isPFLoose())
0152       continue;
0153     pt_t ptInv = in.track[it].isPFTight() ? tk_MAXINVPT_TIGHT_ : tk_MAXINVPT_LOOSE_;
0154     track_good[it] = (in.track[it].hwPt < ptInv) || (iMu[it] != -1);
0155     isEle[it] = false;
0156   }
0157 
0158   // initialize output
0159   out.pfcharged.resize(nTRACK);
0160   out.pfneutral.resize(nSELCALO);
0161   for (unsigned int ipf = 0; ipf < nTRACK; ++ipf)
0162     out.pfcharged[ipf].clear();
0163   for (unsigned int ipf = 0; ipf < nSELCALO; ++ipf)
0164     out.pfneutral[ipf].clear();
0165 
0166   // for each track, find the closest calo
0167   std::vector<int> tk2calo(nTRACK, -1);
0168   for (unsigned int it = 0; it < nTRACK; ++it) {
0169     if (in.track[it].hwPt > 0 && in.track[it].isPFLoose() && iMu[it] == -1) {
0170       pt_t tkCaloPtErr = ptErr_ref(in.region, in.track[it]);
0171       int ibest = best_match_with_pt_ref(dR2MAX_TK_CALO_, in.hadcalo, in.track[it], tkCaloPtErr);
0172       if (ibest != -1) {
0173         if (debug_)
0174           dbgPrintf("FW  \t track  %3d pt %8.2f caloPtErr %6.2f matched to calo %3d pt %8.2f\n",
0175                     it,
0176                     in.track[it].floatPt(),
0177                     Scales::floatPt(tkCaloPtErr),
0178                     ibest,
0179                     in.hadcalo[ibest].floatPt());
0180         track_good[it] = true;
0181         isEle[it] = in.hadcalo[ibest].hwIsEM();
0182         calo_sumtk[ibest] += in.track[it].hwPt;
0183         calo_sumtkErr2[ibest] += tkCaloPtErr * tkCaloPtErr;
0184       }
0185       tk2calo[it] = ibest;  // for emulator info
0186     }
0187   }
0188 
0189   for (unsigned int ic = 0; ic < nCALO; ++ic) {
0190     if (calo_sumtk[ic] > 0) {
0191       pt_t ptdiff = in.hadcalo[ic].hwPt - calo_sumtk[ic];
0192       pt2_t sigmamult =
0193           calo_sumtkErr2[ic];  //  + (calo_sumtkErr2[ic] >> 1)); // this multiplies by 1.5 = sqrt(1.5)^2 ~ (1.2)^2
0194       if (debug_ && (in.hadcalo[ic].hwPt > 0)) {
0195         dbgPrintf(
0196             "FW  \t calo  %3d pt %8.2f [ %7d ] eta %+5.2f [ %+5d ] has a sum track pt %8.2f, difference %7.2f +- %.2f "
0197             "\n",
0198             ic,
0199             in.hadcalo[ic].floatPt(),
0200             in.hadcalo[ic].intPt(),
0201             in.hadcalo[ic].floatEta(),
0202             in.hadcalo[ic].intEta(),
0203             Scales::floatPt(calo_sumtk[ic]),
0204             Scales::floatPt(ptdiff),
0205             std::sqrt(Scales::floatPt(calo_sumtkErr2[ic])));
0206       }
0207       if (ptdiff > 0 && ptdiff * ptdiff > sigmamult) {
0208         calo_subpt[ic] = ptdiff;
0209       } else {
0210         calo_subpt[ic] = 0;
0211       }
0212     } else {
0213       calo_subpt[ic] = in.hadcalo[ic].hwPt;
0214     }
0215     if (debug_ && (in.hadcalo[ic].hwPt > 0))
0216       dbgPrintf(
0217           "FW  \t calo'  %3d pt %8.2f ---> %8.2f \n", ic, in.hadcalo[ic].floatPt(), Scales::floatPt(calo_subpt[ic]));
0218   }
0219 
0220   // copy out charged hadrons
0221   for (unsigned int it = 0; it < nTRACK; ++it) {
0222     if (in.track[it].hwPt > 0 && track_good[it]) {
0223       fillPFCand(in.track[it], out.pfcharged[it], /*isMu=*/(iMu[it] != -1), isEle[it]);
0224       // extra emulator information
0225       if (tk2calo[it] != -1)
0226         out.pfcharged[it].srcCluster = in.hadcalo[tk2calo[it]].src;
0227       if (iMu[it] != -1)
0228         out.pfcharged[it].srcMu = in.muon[iMu[it]].src;
0229     }
0230   }
0231 
0232   // copy out neutral hadrons with sorting and cropping
0233   std::vector<PFNeutralObjEmu> outne_all(nCALO);
0234   for (unsigned int ipf = 0; ipf < nCALO; ++ipf)
0235     outne_all[ipf].clear();
0236   for (unsigned int ic = 0; ic < nCALO; ++ic) {
0237     if (calo_subpt[ic] > 0) {
0238       fillPFCand(in.hadcalo[ic], outne_all[ic], in.hadcalo[ic].hwIsEM());
0239       outne_all[ic].hwPt = calo_subpt[ic];
0240       outne_all[ic].hwEmPt = in.hadcalo[ic].hwIsEM() ? calo_subpt[ic] : pt_t(0);  // FIXME
0241     }
0242   }
0243 
0244   if (nCALO_ == nSELCALO_) {
0245     std::swap(outne_all, out.pfneutral);
0246   } else {
0247     ptsort_ref(nCALO, nSELCALO, outne_all, out.pfneutral);
0248   }
0249 
0250   if (debug_) {
0251     for (unsigned int i = 0; i < nTRACK; ++i) {
0252       if (out.pfcharged[i].hwPt == 0)
0253         continue;
0254       dbgPrintf(
0255           "FW  \t outch %3d: pt %8.2f [ %8d ]  calo eta %+5.2f [ %+5d ]  calo phi %+5.2f [ %+5d ]  pid %d  packed %s\n",
0256           i,
0257           out.pfcharged[i].floatPt(),
0258           out.pfcharged[i].intPt(),
0259           out.pfcharged[i].floatEta(),
0260           out.pfcharged[i].intEta(),
0261           out.pfcharged[i].floatPhi(),
0262           out.pfcharged[i].intPhi(),
0263           out.pfcharged[i].intId(),
0264           out.pfcharged[i].pack().to_string(16).c_str());
0265     }
0266     for (unsigned int i = 0; i < nSELCALO; ++i) {
0267       if (out.pfneutral[i].hwPt == 0)
0268         continue;
0269       dbgPrintf(
0270           "FW  \t outne %3d: pt %8.2f [ %8d ]  calo eta %+5.2f [ %+5d ]  calo phi %+5.2f [ %+5d ]  pid %d  packed %s\n",
0271           i,
0272           out.pfneutral[i].floatPt(),
0273           out.pfneutral[i].intPt(),
0274           out.pfneutral[i].floatEta(),
0275           out.pfneutral[i].intEta(),
0276           out.pfneutral[i].floatPhi(),
0277           out.pfneutral[i].intPhi(),
0278           out.pfneutral[i].intId(),
0279           out.pfneutral[i].pack().to_string(16).c_str());
0280     }
0281   }
0282 }