Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:21:25

0001 #ifndef PFTKEGALGO_REF_H
0002 #define PFTKEGALGO_REF_H
0003 
0004 #include "DataFormats/L1TParticleFlow/interface/layer1_emulator.h"
0005 #include "DataFormats/L1TParticleFlow/interface/egamma.h"
0006 #include "DataFormats/L1TParticleFlow/interface/pf.h"
0007 #include "L1Trigger/Phase2L1ParticleFlow/interface/common/inversion.h"
0008 
0009 #include "conifer.h"
0010 
0011 namespace edm {
0012   class ParameterSet;
0013   class ParameterSetDescription;
0014 }  // namespace edm
0015 
0016 namespace l1ct {
0017 
0018   struct PFTkEGAlgoEmuConfig {
0019     unsigned int nTRACK;
0020     unsigned int nTRACK_EGIN;
0021     unsigned int nEMCALO_EGIN;
0022     unsigned int nEM_EGOUT;
0023 
0024     bool filterHwQuality;
0025     bool doBremRecovery;
0026     bool writeBeforeBremRecovery;
0027     int caloHwQual;
0028     bool doEndcapHwQual;
0029     float emClusterPtMin;  // GeV
0030     float dEtaMaxBrem;
0031     float dPhiMaxBrem;
0032 
0033     std::vector<double> absEtaBoundaries;
0034     std::vector<double> dEtaValues;
0035     std::vector<double> dPhiValues;
0036     float trkQualityPtMin;  // GeV
0037     bool doCompositeTkEle;
0038     unsigned int nCompCandPerCluster;
0039     bool writeEgSta;
0040 
0041     struct IsoParameters {
0042       IsoParameters(const edm::ParameterSet &);
0043       IsoParameters(float tkQualityPtMin, float dZ, float dRMin, float dRMax)
0044           : tkQualityPtMin(Scales::makePtFromFloat(tkQualityPtMin)),
0045             dZ(Scales::makeZ0(dZ)),
0046             dRMin2(Scales::makeDR2FromFloatDR(dRMin)),
0047             dRMax2(Scales::makeDR2FromFloatDR(dRMax)) {}
0048       pt_t tkQualityPtMin;
0049       ap_int<z0_t::width + 1> dZ;
0050       int dRMin2;
0051       int dRMax2;
0052       static edm::ParameterSetDescription getParameterSetDescription();
0053     };
0054 
0055     IsoParameters tkIsoParams_tkEle;
0056     IsoParameters tkIsoParams_tkEm;
0057     IsoParameters pfIsoParams_tkEle;
0058     IsoParameters pfIsoParams_tkEm;
0059     bool doTkIso;
0060     bool doPfIso;
0061     EGIsoEleObjEmu::IsoType hwIsoTypeTkEle;
0062     EGIsoObjEmu::IsoType hwIsoTypeTkEm;
0063 
0064     struct CompIDParameters {
0065       CompIDParameters(const edm::ParameterSet &);
0066       CompIDParameters(double bdtScore_loose_wp, double bdtScore_tight_wp, const std::string &model)
0067           : bdtScore_loose_wp(bdtScore_loose_wp), bdtScore_tight_wp(bdtScore_tight_wp), conifer_model(model) {}
0068       const id_score_t bdtScore_loose_wp;  // Conifer score/4
0069       const id_score_t bdtScore_tight_wp;  // Conifer score/4
0070       const std::string conifer_model;
0071       static edm::ParameterSetDescription getParameterSetDescription();
0072     };
0073 
0074     CompIDParameters compIDparams;
0075 
0076     int debug = 0;
0077 
0078     PFTkEGAlgoEmuConfig(const edm::ParameterSet &iConfig);
0079     PFTkEGAlgoEmuConfig(unsigned int nTrack,
0080                         unsigned int nTrack_in,
0081                         unsigned int nEmCalo_in,
0082                         unsigned int nEmOut,
0083                         bool filterHwQuality,
0084                         bool doBremRecovery,
0085                         bool writeBeforeBremRecovery = false,
0086                         int caloHwQual = 4,
0087                         bool doEndcapHwQual = false,
0088                         float emClusterPtMin = 2.,
0089                         float dEtaMaxBrem = 0.02,
0090                         float dPhiMaxBrem = 0.1,
0091                         const std::vector<double> &absEtaBoundaries = {0.0, 1.5},
0092                         const std::vector<double> &dEtaValues = {0.015, 0.01},
0093                         const std::vector<double> &dPhiValues = {0.07, 0.07},
0094                         float trkQualityPtMin = 10.,
0095                         bool doCompositeTkEle = false,
0096                         unsigned int nCompCandPerCluster = 4,
0097                         bool writeEgSta = false,
0098                         const IsoParameters &tkIsoParams_tkEle = {2., 0.6, 0.03, 0.2},
0099                         const IsoParameters &tkIsoParams_tkEm = {2., 0.6, 0.07, 0.3},
0100                         const IsoParameters &pfIsoParams_tkEle = {1., 0.6, 0.03, 0.2},
0101                         const IsoParameters &pfIsoParams_tkEm = {1., 0.6, 0.07, 0.3},
0102                         bool doTkIso = true,
0103                         bool doPfIso = false,
0104                         EGIsoEleObjEmu::IsoType hwIsoTypeTkEle = EGIsoEleObjEmu::IsoType::TkIso,
0105                         EGIsoObjEmu::IsoType hwIsoTypeTkEm = EGIsoObjEmu::IsoType::TkIsoPV,
0106                         const CompIDParameters &compIDparams = {-4, 0.214844, "compositeID.json"},
0107                         int debug = 0)
0108 
0109         : nTRACK(nTrack),
0110           nTRACK_EGIN(nTrack_in),
0111           nEMCALO_EGIN(nEmCalo_in),
0112           nEM_EGOUT(nEmOut),
0113           filterHwQuality(filterHwQuality),
0114           doBremRecovery(doBremRecovery),
0115           writeBeforeBremRecovery(writeBeforeBremRecovery),
0116           caloHwQual(caloHwQual),
0117           doEndcapHwQual(doEndcapHwQual),
0118           emClusterPtMin(emClusterPtMin),
0119           dEtaMaxBrem(dEtaMaxBrem),
0120           dPhiMaxBrem(dPhiMaxBrem),
0121           absEtaBoundaries(absEtaBoundaries),
0122           dEtaValues(dEtaValues),
0123           dPhiValues(dPhiValues),
0124           trkQualityPtMin(trkQualityPtMin),
0125           doCompositeTkEle(doCompositeTkEle),
0126           nCompCandPerCluster(nCompCandPerCluster),
0127           writeEgSta(writeEgSta),
0128           tkIsoParams_tkEle(tkIsoParams_tkEle),
0129           tkIsoParams_tkEm(tkIsoParams_tkEm),
0130           pfIsoParams_tkEle(pfIsoParams_tkEle),
0131           pfIsoParams_tkEm(pfIsoParams_tkEm),
0132           doTkIso(doTkIso),
0133           doPfIso(doPfIso),
0134           hwIsoTypeTkEle(hwIsoTypeTkEle),
0135           hwIsoTypeTkEm(hwIsoTypeTkEm),
0136           compIDparams(compIDparams),
0137           debug(debug) {}
0138 
0139     static edm::ParameterSetDescription getParameterSetDescription();
0140   };
0141 
0142   class PFTkEGAlgoEmulator {
0143   public:
0144     PFTkEGAlgoEmulator(const PFTkEGAlgoEmuConfig &config);
0145 
0146     virtual ~PFTkEGAlgoEmulator() {}
0147 
0148     void toFirmware(const PFInputRegion &in, PFRegion &region, EmCaloObj calo[/*nCALO*/], TkObj track[/*nTRACK*/]) const;
0149     void toFirmware(const OutputRegion &out, EGIsoObj out_egphs[], EGIsoEleObj out_egeles[]) const;
0150     void toFirmware(const PFInputRegion &in,
0151                     const l1ct::PVObjEmu &pvin,
0152                     PFRegion &region,
0153                     TkObj track[/*nTRACK*/],
0154                     PVObj &pv) const;
0155 
0156     void run(const PFInputRegion &in, OutputRegion &out) const;
0157     void runIso(const PFInputRegion &in, const std::vector<l1ct::PVObjEmu> &pvs, OutputRegion &out) const;
0158 
0159     void setDebug(int verbose) { debug_ = verbose; }
0160 
0161     bool writeEgSta() const { return cfg.writeEgSta; }
0162 
0163     typedef ap_fixed<21, 12, AP_RND_CONV, AP_SAT> bdt_feature_t;
0164     typedef ap_fixed<12, 3, AP_RND_CONV, AP_SAT> bdt_score_t;
0165 
0166   private:
0167     void link_emCalo2emCalo(const std::vector<EmCaloObjEmu> &emcalo, std::vector<int> &emCalo2emCalo) const;
0168 
0169     void link_emCalo2tk_elliptic(const PFRegionEmu &r,
0170                                  const std::vector<EmCaloObjEmu> &emcalo,
0171                                  const std::vector<TkObjEmu> &track,
0172                                  std::vector<int> &emCalo2tk) const;
0173 
0174     void link_emCalo2tk_composite(const PFRegionEmu &r,
0175                                   const std::vector<EmCaloObjEmu> &emcalo,
0176                                   const std::vector<TkObjEmu> &track,
0177                                   std::vector<int> &emCalo2tk,
0178                                   std::vector<id_score_t> &emCaloTkBdtScore) const;
0179 
0180     struct CompositeCandidate {
0181       unsigned int cluster_idx;
0182       unsigned int track_idx;
0183       double dpt;  // For sorting
0184     };
0185 
0186     id_score_t compute_composite_score(CompositeCandidate &cand,
0187                                        const std::vector<EmCaloObjEmu> &emcalo,
0188                                        const std::vector<TkObjEmu> &track,
0189                                        const PFTkEGAlgoEmuConfig::CompIDParameters &params) const;
0190 
0191     //FIXME: still needed
0192     float deltaPhi(float phi1, float phi2) const;
0193 
0194     void sel_emCalo(unsigned int nmax_sel,
0195                     const std::vector<EmCaloObjEmu> &emcalo,
0196                     std::vector<EmCaloObjEmu> &emcalo_sel) const;
0197 
0198     void eg_algo(const PFRegionEmu &region,
0199                  const std::vector<EmCaloObjEmu> &emcalo,
0200                  const std::vector<TkObjEmu> &track,
0201                  const std::vector<int> &emCalo2emCalo,
0202                  const std::vector<int> &emCalo2tk,
0203                  const std::vector<id_score_t> &emCaloTkBdtScore,
0204                  std::vector<EGObjEmu> &egstas,
0205                  std::vector<EGIsoObjEmu> &egobjs,
0206                  std::vector<EGIsoEleObjEmu> &egeleobjs) const;
0207 
0208     void addEgObjsToPF(std::vector<EGObjEmu> &egstas,
0209                        std::vector<EGIsoObjEmu> &egobjs,
0210                        std::vector<EGIsoEleObjEmu> &egeleobjs,
0211                        const std::vector<EmCaloObjEmu> &emcalo,
0212                        const std::vector<TkObjEmu> &track,
0213                        const int calo_idx,
0214                        const unsigned int hwQual,
0215                        const pt_t ptCorr,
0216                        const int tk_idx,
0217                        const id_score_t bdtScore,
0218                        const std::vector<unsigned int> &components = {}) const;
0219 
0220     EGObjEmu &addEGStaToPF(std::vector<EGObjEmu> &egobjs,
0221                            const EmCaloObjEmu &calo,
0222                            const unsigned int hwQual,
0223                            const pt_t ptCorr,
0224                            const std::vector<unsigned int> &components) const;
0225 
0226     EGIsoObjEmu &addEGIsoToPF(std::vector<EGIsoObjEmu> &egobjs,
0227                               const EmCaloObjEmu &calo,
0228                               const unsigned int hwQual,
0229                               const pt_t ptCorr) const;
0230 
0231     EGIsoEleObjEmu &addEGIsoEleToPF(std::vector<EGIsoEleObjEmu> &egobjs,
0232                                     const EmCaloObjEmu &calo,
0233                                     const TkObjEmu &track,
0234                                     const unsigned int hwQual,
0235                                     const pt_t ptCorr,
0236                                     const id_score_t bdtScore) const;
0237 
0238     // FIXME: reimplemented from PFAlgoEmulatorBase
0239     template <typename T>
0240     void ptsort_ref(int nIn, int nOut, const std::vector<T> &in, std::vector<T> &out) const {
0241       out.resize(nOut);
0242       for (int iout = 0; iout < nOut; ++iout) {
0243         out[iout].clear();
0244       }
0245       for (int it = 0; it < nIn; ++it) {
0246         for (int iout = 0; iout < nOut; ++iout) {
0247           if (in[it].hwPt >= out[iout].hwPt) {
0248             for (int i2 = nOut - 1; i2 > iout; --i2) {
0249               out[i2] = out[i2 - 1];
0250             }
0251             out[iout] = in[it];
0252             break;
0253           }
0254         }
0255       }
0256     }
0257 
0258     template <typename T>
0259     int deltaR2(const T &charged, const EGIsoObjEmu &egphoton) const {
0260       // NOTE: we compare Tk at vertex against the calo variable...
0261       return dr2_int(charged.hwVtxEta(), charged.hwVtxPhi(), egphoton.hwEta, egphoton.hwPhi);
0262     }
0263 
0264     template <typename T>
0265     int deltaR2(const T &charged, const EGIsoEleObjEmu &egele) const {
0266       return dr2_int(charged.hwVtxEta(), charged.hwVtxPhi(), egele.hwVtxEta(), egele.hwVtxPhi());
0267     }
0268 
0269     int deltaR2(const PFNeutralObjEmu &neutral, const EGIsoObjEmu &egphoton) const {
0270       return dr2_int(neutral.hwEta, neutral.hwPhi, egphoton.hwEta, egphoton.hwPhi);
0271     }
0272 
0273     int deltaR2(const PFNeutralObjEmu &neutral, const EGIsoEleObjEmu &egele) const {
0274       // NOTE: we compare Tk at vertex against the calo variable...
0275       return dr2_int(neutral.hwEta, neutral.hwPhi, egele.hwVtxEta(), egele.hwVtxPhi());
0276     }
0277 
0278     template <typename T>
0279     ap_int<z0_t::width + 1> deltaZ0(const T &charged, const EGIsoObjEmu &egphoton, z0_t z0) const {
0280       ap_int<z0_t::width + 1> delta = charged.hwZ0 - z0;
0281       if (delta < 0)
0282         delta = -delta;
0283       return delta;
0284     }
0285 
0286     template <typename T>
0287     ap_int<z0_t::width + 1> deltaZ0(const T &charged, const EGIsoEleObjEmu &egele, z0_t z0) const {
0288       ap_int<z0_t::width + 1> delta = charged.hwZ0 - egele.hwZ0;
0289       if (delta < 0)
0290         delta = -delta;
0291       return delta;
0292     }
0293 
0294     template <typename TCH, typename TEG>
0295     void compute_sumPt(iso_t &sumPt,
0296                        iso_t &sumPtPV,
0297                        const std::vector<TCH> &objects,
0298                        unsigned int nMaxObj,
0299                        const TEG &egobj,
0300                        const PFTkEGAlgoEmuConfig::IsoParameters &params,
0301                        z0_t z0) const {
0302       for (unsigned int itk = 0; itk < std::min<unsigned>(objects.size(), nMaxObj); ++itk) {
0303         const auto &obj = objects[itk];
0304 
0305         if (obj.hwPt < params.tkQualityPtMin)
0306           continue;
0307 
0308         int dR2 = deltaR2(obj, egobj);
0309 
0310         if (dR2 > params.dRMin2 && dR2 < params.dRMax2) {
0311           sumPt += obj.hwPt;
0312           if (deltaZ0(obj, egobj, z0) < params.dZ) {
0313             sumPtPV += obj.hwPt;
0314           }
0315         }
0316       }
0317     }
0318 
0319     template <typename TEG>
0320     void compute_sumPt(iso_t &sumPt,
0321                        iso_t &sumPtPV,
0322                        const std::vector<PFNeutralObjEmu> &objects,
0323                        unsigned int nMaxObj,
0324                        const TEG &egobj,
0325                        const PFTkEGAlgoEmuConfig::IsoParameters &params,
0326                        z0_t z0) const {
0327       for (unsigned int itk = 0; itk < std::min<unsigned>(objects.size(), nMaxObj); ++itk) {
0328         const auto &obj = objects[itk];
0329 
0330         if (obj.hwPt < params.tkQualityPtMin)
0331           continue;
0332 
0333         int dR2 = deltaR2(obj, egobj);
0334 
0335         if (dR2 > params.dRMin2 && dR2 < params.dRMax2) {
0336           sumPt += obj.hwPt;
0337           // PF neutrals are not constrained by PV (since their Z0 is 0 by design)
0338           sumPtPV += obj.hwPt;
0339         }
0340       }
0341     }
0342 
0343     void compute_isolation(std::vector<EGIsoObjEmu> &egobjs,
0344                            const std::vector<TkObjEmu> &objects,
0345                            const PFTkEGAlgoEmuConfig::IsoParameters &params,
0346                            z0_t z0) const;
0347     void compute_isolation(std::vector<EGIsoEleObjEmu> &egobjs,
0348                            const std::vector<TkObjEmu> &objects,
0349                            const PFTkEGAlgoEmuConfig::IsoParameters &params,
0350                            z0_t z0) const;
0351     void compute_isolation(std::vector<EGIsoObjEmu> &egobjs,
0352                            const std::vector<PFChargedObjEmu> &charged,
0353                            const std::vector<PFNeutralObjEmu> &neutrals,
0354                            const PFTkEGAlgoEmuConfig::IsoParameters &params,
0355                            z0_t z0) const;
0356     void compute_isolation(std::vector<EGIsoEleObjEmu> &egobjs,
0357                            const std::vector<PFChargedObjEmu> &charged,
0358                            const std::vector<PFNeutralObjEmu> &neutrals,
0359                            const PFTkEGAlgoEmuConfig::IsoParameters &params,
0360                            z0_t z0) const;
0361 
0362     PFTkEGAlgoEmuConfig cfg;
0363     conifer::BDT<bdt_feature_t, ap_fixed<12, 3, AP_RND_CONV, AP_SAT>, false> *composite_bdt_;
0364     int debug_;
0365   };
0366 }  // namespace l1ct
0367 
0368 #endif