File indexing completed on 2022-04-14 02:41:01
0001 #include "FWCore/Framework/interface/stream/EDProducer.h"
0002 #include "FWCore/Framework/interface/Event.h"
0003 #include "FWCore/Framework/interface/MakerMacros.h"
0004
0005 #include "DataFormats/L1THGCal/interface/HGCalMulticluster.h"
0006 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0007 #include "DataFormats/L1Trigger/interface/EGamma.h"
0008 #include "L1Trigger/L1CaloTrigger/interface/L1EGammaEECalibrator.h"
0009 #include "DataFormats/Math/interface/deltaPhi.h"
0010
0011 namespace l1tp2 {
0012
0013 bool compare_cluster_pt(const l1t::HGCalMulticluster *cl1, const l1t::HGCalMulticluster *cl2) {
0014 return cl1->pt() > cl2->pt();
0015 }
0016 };
0017
0018 int etaBin(const l1t::HGCalMulticluster *cl) {
0019 static float constexpr eta_min = 1.;
0020 static float constexpr eta_max = 4.;
0021 static unsigned constexpr n_eta_bins = 150;
0022 int eta_bin = floor((std::abs(cl->eta()) - eta_min) / ((eta_max - eta_min) / n_eta_bins));
0023 if (cl->eta() < 0)
0024 return -1 * eta_bin;
0025 return eta_bin;
0026 }
0027
0028 int get_phi_bin(const l1t::HGCalMulticluster *cl) {
0029 static constexpr float phi_min = -M_PI;
0030 static constexpr float phi_max = M_PI;
0031 static constexpr unsigned n_phi_bins = 63;
0032 return floor(std::abs(reco::deltaPhi(cl->phi(), phi_min)) / ((phi_max - phi_min) / n_phi_bins));
0033 }
0034
0035 pair<int, int> get_eta_phi_bin(const l1t::HGCalMulticluster *cl) { return std::make_pair(etaBin(cl), get_phi_bin(cl)); }
0036
0037 class L1EGammaEEProducer : public edm::stream::EDProducer<> {
0038 public:
0039 explicit L1EGammaEEProducer(const edm::ParameterSet &);
0040
0041 private:
0042 void produce(edm::Event &, const edm::EventSetup &) override;
0043
0044 edm::EDGetToken multiclusters_token_;
0045 L1EGammaEECalibrator calibrator_;
0046 };
0047
0048 L1EGammaEEProducer::L1EGammaEEProducer(const edm::ParameterSet &iConfig)
0049 : multiclusters_token_(
0050 consumes<l1t::HGCalMulticlusterBxCollection>(iConfig.getParameter<edm::InputTag>("Multiclusters"))),
0051 calibrator_(iConfig.getParameter<edm::ParameterSet>("calibrationConfig")) {
0052 produces<BXVector<l1t::EGamma>>("L1EGammaCollectionBXVWithCuts");
0053 }
0054
0055 void L1EGammaEEProducer::produce(edm::Event &iEvent, const edm::EventSetup &iSetup) {
0056 float minEt_ = 0;
0057
0058 std::unique_ptr<BXVector<l1t::EGamma>> l1EgammaBxCollection(new l1t::EGammaBxCollection);
0059
0060
0061 edm::Handle<l1t::HGCalMulticlusterBxCollection> multiclusters_h;
0062 iEvent.getByToken(multiclusters_token_, multiclusters_h);
0063 const l1t::HGCalMulticlusterBxCollection &multiclusters = *multiclusters_h;
0064
0065 std::vector<const l1t::HGCalMulticluster *> selected_multiclusters;
0066 std::map<std::pair<int, int>, std::vector<const l1t::HGCalMulticluster *>> etaphi_bins;
0067
0068
0069 for (auto cl3d = multiclusters.begin(0); cl3d != multiclusters.end(0); cl3d++) {
0070 if (cl3d->hwQual()) {
0071 if (cl3d->et() > minEt_) {
0072 int hw_quality = 1;
0073 if (std::abs(cl3d->eta()) >= 1.52) {
0074 hw_quality = 2;
0075 }
0076
0077 float calib_factor = calibrator_.calibrationFactor(cl3d->pt(), cl3d->eta());
0078 l1t::EGamma eg =
0079 l1t::EGamma(reco::Candidate::PolarLorentzVector(cl3d->pt() / calib_factor, cl3d->eta(), cl3d->phi(), 0.));
0080 eg.setHwQual(hw_quality);
0081 eg.setHwIso(1);
0082 eg.setIsoEt(-1);
0083 l1EgammaBxCollection->push_back(0, eg);
0084 if (hw_quality == 2) {
0085
0086 l1t::EGamma eg_emint = l1t::EGamma(reco::Candidate::PolarLorentzVector(
0087 cl3d->iPt(l1t::HGCalMulticluster::EnergyInterpretation::EM), cl3d->eta(), cl3d->phi(), 0.));
0088 eg_emint.setHwQual(4);
0089 eg_emint.setHwIso(1);
0090 eg_emint.setIsoEt(-1);
0091 l1EgammaBxCollection->push_back(0, eg_emint);
0092
0093 selected_multiclusters.push_back(&(*cl3d));
0094 auto eta_phi_bin = get_eta_phi_bin(&(*cl3d));
0095 auto bucket = etaphi_bins.find(eta_phi_bin);
0096 if (bucket == etaphi_bins.end()) {
0097 std::vector<const l1t::HGCalMulticluster *> vec;
0098 vec.push_back(&(*cl3d));
0099 etaphi_bins[eta_phi_bin] = vec;
0100 } else {
0101 bucket->second.push_back(&(*cl3d));
0102 }
0103 }
0104 }
0105 }
0106 }
0107
0108 std::sort(selected_multiclusters.begin(), selected_multiclusters.end(), l1tp2::compare_cluster_pt);
0109 std::set<const l1t::HGCalMulticluster *> used_clusters;
0110 for (const auto &cl3d : selected_multiclusters) {
0111 if (used_clusters.find(cl3d) == used_clusters.end()) {
0112 float pt = cl3d->pt();
0113
0114 if (cl3d->hOverE() != -1)
0115 pt = cl3d->pt() / (1 + cl3d->hOverE());
0116 reco::Candidate::PolarLorentzVector mom(pt, cl3d->eta(), cl3d->phi(), 0.);
0117 reco::Candidate::PolarLorentzVector mom_eint(
0118 cl3d->iPt(l1t::HGCalMulticluster::EnergyInterpretation::EM), cl3d->eta(), cl3d->phi(), 0.);
0119
0120
0121 used_clusters.insert(cl3d);
0122 auto eta_phi_bin = get_eta_phi_bin(cl3d);
0123
0124 for (int eta_bin : {eta_phi_bin.first - 1, eta_phi_bin.first, eta_phi_bin.first + 1}) {
0125 for (int phi_bin : {eta_phi_bin.second - 1, eta_phi_bin.second, eta_phi_bin.second + 1}) {
0126 auto bucket = etaphi_bins.find(std::make_pair(eta_bin, phi_bin));
0127 if (bucket != etaphi_bins.end()) {
0128
0129 for (const auto &other_cl_ptr : bucket->second) {
0130 if (used_clusters.find(other_cl_ptr) == used_clusters.end()) {
0131 if (std::abs(other_cl_ptr->eta() - cl3d->eta()) < 0.02) {
0132 if (std::abs(reco::deltaPhi(other_cl_ptr->phi(), cl3d->phi())) < 0.1) {
0133 float pt_other = other_cl_ptr->pt();
0134 if (other_cl_ptr->hOverE() != -1)
0135 pt_other = other_cl_ptr->pt() / (1 + other_cl_ptr->hOverE());
0136 mom += reco::Candidate::PolarLorentzVector(pt_other, other_cl_ptr->eta(), other_cl_ptr->phi(), 0.);
0137 mom_eint += reco::Candidate::PolarLorentzVector(
0138 other_cl_ptr->iPt(l1t::HGCalMulticluster::EnergyInterpretation::EM),
0139 other_cl_ptr->eta(),
0140 other_cl_ptr->phi(),
0141 0.);
0142 used_clusters.insert(other_cl_ptr);
0143 }
0144 }
0145 }
0146 }
0147 }
0148 }
0149 }
0150 float calib_factor = calibrator_.calibrationFactor(mom.pt(), mom.eta());
0151 l1t::EGamma eg =
0152 l1t::EGamma(reco::Candidate::PolarLorentzVector(mom.pt() / calib_factor, mom.eta(), mom.phi(), 0.));
0153 eg.setHwQual(3);
0154 eg.setHwIso(1);
0155 l1EgammaBxCollection->push_back(0, eg);
0156
0157 l1t::EGamma eg_emint_brec =
0158 l1t::EGamma(reco::Candidate::PolarLorentzVector(mom_eint.pt(), mom_eint.eta(), mom_eint.phi(), 0.));
0159 eg_emint_brec.setHwQual(5);
0160 eg_emint_brec.setHwIso(1);
0161 l1EgammaBxCollection->push_back(0, eg_emint_brec);
0162 }
0163 }
0164
0165 iEvent.put(std::move(l1EgammaBxCollection), "L1EGammaCollectionBXVWithCuts");
0166 }
0167
0168 DEFINE_FWK_MODULE(L1EGammaEEProducer);