File indexing completed on 2023-05-05 02:47:53
0001
0002 #include <memory>
0003 #include <algorithm>
0004
0005
0006 #include "FWCore/Framework/interface/Frameworkfwd.h"
0007 #include "FWCore/Framework/interface/stream/EDProducer.h"
0008 #include "FWCore/Framework/interface/Event.h"
0009 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0010 #include "FWCore/ParameterSet/interface/FileInPath.h"
0011 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0012
0013 #include "DataFormats/HcalDetId/interface/HcalTrigTowerDetId.h"
0014 #include "DataFormats/HcalDigi/interface/HcalDigiCollections.h"
0015 #include "CalibFormats/CaloTPG/interface/CaloTPGTranscoder.h"
0016 #include "CalibFormats/CaloTPG/interface/CaloTPGRecord.h"
0017 #include "L1Trigger/L1TCalorimeter/interface/CaloTools.h"
0018
0019 #include "DataFormats/L1THGCal/interface/HGCalTower.h"
0020
0021 #include "DataFormats/L1TCalorimeterPhase2/interface/CaloTower.h"
0022
0023 #include "DataFormats/Math/interface/deltaPhi.h"
0024
0025 #include "L1Trigger/Phase2L1ParticleFlow/interface/corrector.h"
0026 #include "L1Trigger/Phase2L1ParticleFlow/interface/ParametricResolution.h"
0027 #include "L1Trigger/Phase2L1ParticleFlow/interface/CaloClusterer.h"
0028
0029
0030 class L1TPFCaloProducer : public edm::stream::EDProducer<> {
0031 public:
0032 explicit L1TPFCaloProducer(const edm::ParameterSet &);
0033
0034 private:
0035 bool ecalOnly_, debug_;
0036 std::vector<edm::EDGetTokenT<reco::CandidateView>> ecalCands_;
0037 std::vector<edm::EDGetTokenT<reco::CandidateView>> hcalCands_;
0038
0039 std::vector<edm::EDGetTokenT<HcalTrigPrimDigiCollection>> hcalDigis_;
0040 edm::ESGetToken<CaloTPGTranscoder, CaloTPGRecord> decoderTag_;
0041 bool hcalDigisBarrel_, hcalDigisHF_;
0042 std::vector<edm::EDGetTokenT<l1tp2::CaloTowerCollection>> phase2barrelTowers_;
0043 std::vector<edm::EDGetTokenT<l1t::HGCalTowerBxCollection>> hcalHGCTowers_;
0044 bool hcalHGCTowersHadOnly_;
0045
0046 l1tpf::corrector emCorrector_;
0047 l1tpf::corrector hcCorrector_;
0048 l1tpf::corrector hadCorrector_;
0049
0050 l1tpf_calo::SingleCaloClusterer ecalClusterer_, hcalClusterer_;
0051 std::unique_ptr<l1tpf_calo::SimpleCaloLinkerBase> caloLinker_;
0052
0053 l1tpf::ParametricResolution resol_;
0054
0055 void produce(edm::Event &, const edm::EventSetup &) override;
0056
0057 void readHcalDigis_(edm::Event &event, const edm::EventSetup &);
0058 void readPhase2BarrelCaloTowers_(edm::Event &event, const edm::EventSetup &);
0059 void readHcalHGCTowers_(edm::Event &event, const edm::EventSetup &);
0060 struct SimpleHGCTC {
0061 float et, eta, phi;
0062 SimpleHGCTC(float aet, float aeta, float aphi) : et(aet), eta(aeta), phi(aphi) {}
0063 };
0064 };
0065
0066
0067
0068
0069 L1TPFCaloProducer::L1TPFCaloProducer(const edm::ParameterSet &iConfig)
0070 : ecalOnly_(iConfig.existsAs<bool>("ecalOnly") ? iConfig.getParameter<bool>("ecalOnly") : false),
0071 debug_(iConfig.getUntrackedParameter<int>("debug", 0)),
0072 decoderTag_(esConsumes<CaloTPGTranscoder, CaloTPGRecord>(edm::ESInputTag("", ""))),
0073 emCorrector_(iConfig.getParameter<std::string>("emCorrector"), -1, debug_),
0074 hcCorrector_(iConfig.getParameter<std::string>("hcCorrector"), -1, debug_),
0075 hadCorrector_(iConfig.getParameter<std::string>("hadCorrector"),
0076 iConfig.getParameter<double>("hadCorrectorEmfMax"),
0077 debug_),
0078 ecalClusterer_(iConfig.getParameter<edm::ParameterSet>("ecalClusterer")),
0079 hcalClusterer_(iConfig.getParameter<edm::ParameterSet>("hcalClusterer")),
0080 caloLinker_(l1tpf_calo::makeCaloLinker(
0081 iConfig.getParameter<edm::ParameterSet>("linker"), ecalClusterer_, hcalClusterer_)),
0082 resol_(iConfig.getParameter<edm::ParameterSet>("resol")) {
0083 produces<l1t::PFClusterCollection>("ecalCells");
0084
0085 produces<l1t::PFClusterCollection>("emCalibrated");
0086 produces<l1t::PFClusterCollection>("emUncalibrated");
0087
0088 for (auto &tag : iConfig.getParameter<std::vector<edm::InputTag>>("ecalCandidates")) {
0089 ecalCands_.push_back(consumes<reco::CandidateView>(tag));
0090 }
0091
0092 if (ecalOnly_)
0093 return;
0094
0095 produces<l1t::PFClusterCollection>("hcalCells");
0096
0097 produces<l1t::PFClusterCollection>("hcalUnclustered");
0098 produces<l1t::PFClusterCollection>("hcalUncalibrated");
0099 produces<l1t::PFClusterCollection>("hcalCalibrated");
0100
0101 produces<l1t::PFClusterCollection>("uncalibrated");
0102 produces<l1t::PFClusterCollection>("calibrated");
0103
0104 for (auto &tag : iConfig.getParameter<std::vector<edm::InputTag>>("hcalCandidates")) {
0105 hcalCands_.push_back(consumes<reco::CandidateView>(tag));
0106 }
0107
0108 for (auto &tag : iConfig.getParameter<std::vector<edm::InputTag>>("hcalDigis")) {
0109 hcalDigis_.push_back(consumes<HcalTrigPrimDigiCollection>(tag));
0110 }
0111 if (!hcalDigis_.empty()) {
0112 hcalDigisBarrel_ = iConfig.getParameter<bool>("hcalDigisBarrel");
0113 hcalDigisHF_ = iConfig.getParameter<bool>("hcalDigisHF");
0114 }
0115
0116 for (auto &tag : iConfig.getParameter<std::vector<edm::InputTag>>("phase2barrelCaloTowers")) {
0117 phase2barrelTowers_.push_back(consumes<l1tp2::CaloTowerCollection>(tag));
0118 }
0119
0120 for (auto &tag : iConfig.getParameter<std::vector<edm::InputTag>>("hcalHGCTowers")) {
0121 hcalHGCTowers_.push_back(consumes<l1t::HGCalTowerBxCollection>(tag));
0122 }
0123 if (!hcalHGCTowers_.empty())
0124 hcalHGCTowersHadOnly_ = iConfig.getParameter<bool>("hcalHGCTowersHadOnly");
0125 }
0126
0127
0128 void L1TPFCaloProducer::produce(edm::Event &iEvent, const edm::EventSetup &iSetup) {
0129
0130 edm::Handle<reco::CandidateView> ecals;
0131 for (const auto &token : ecalCands_) {
0132 iEvent.getByToken(token, ecals);
0133 for (const reco::Candidate &it : *ecals) {
0134 if (debug_)
0135 edm::LogWarning("L1TPFCaloProducer")
0136 << "adding ECal input pt " << it.pt() << ", eta " << it.eta() << ", phi " << it.phi() << "\n";
0137 ecalClusterer_.add(it);
0138 }
0139 }
0140
0141
0142 if (!ecalOnly_) {
0143 edm::Handle<reco::CandidateView> hcals;
0144 for (const auto &token : hcalCands_) {
0145 iEvent.getByToken(token, hcals);
0146 for (const reco::Candidate &it : *hcals) {
0147 if (debug_)
0148 edm::LogWarning("L1TPFCaloProducer")
0149 << "adding HCal cand input pt " << it.pt() << ", eta " << it.eta() << ", phi " << it.phi() << "\n";
0150 hcalClusterer_.add(it);
0151 }
0152 }
0153 if (!hcalDigis_.empty()) {
0154 readHcalDigis_(iEvent, iSetup);
0155 }
0156 if (!phase2barrelTowers_.empty()) {
0157 readPhase2BarrelCaloTowers_(iEvent, iSetup);
0158 }
0159 if (!hcalHGCTowers_.empty()) {
0160 readHcalHGCTowers_(iEvent, iSetup);
0161 }
0162 }
0163
0164
0165 ecalClusterer_.run();
0166
0167 auto ecalCellsH = iEvent.put(ecalClusterer_.fetchCells(), "ecalCells");
0168
0169 iEvent.put(ecalClusterer_.fetch(ecalCellsH), "emUncalibrated");
0170
0171 if (emCorrector_.valid()) {
0172 ecalClusterer_.correct(
0173 [&](const l1tpf_calo::Cluster &c) -> float { return emCorrector_.correctedPt(0., c.et, std::abs(c.eta)); });
0174 }
0175
0176 std::unique_ptr<l1t::PFClusterCollection> corrEcal = ecalClusterer_.fetch(ecalCellsH);
0177
0178 if (debug_) {
0179 for (const l1t::PFCluster &it : *corrEcal) {
0180 edm::LogWarning("L1TPFCaloProducer")
0181 << "corrected ECal cluster pt " << it.pt() << ", eta " << it.eta() << ", phi " << it.phi() << "\n";
0182 }
0183 }
0184
0185 auto ecalClustH = iEvent.put(std::move(corrEcal), "emCalibrated");
0186
0187 if (ecalOnly_) {
0188 ecalClusterer_.clear();
0189 return;
0190 }
0191
0192 hcalClusterer_.run();
0193
0194 auto hcalCellsH = iEvent.put(hcalClusterer_.fetchCells(), "hcalCells");
0195
0196
0197 iEvent.put(hcalClusterer_.fetchCells(true), "hcalUnclustered");
0198
0199 iEvent.put(hcalClusterer_.fetch(hcalCellsH), "hcalUncalibrated");
0200
0201 if (hcCorrector_.valid()) {
0202 hcalClusterer_.correct(
0203 [&](const l1tpf_calo::Cluster &c) -> float { return hcCorrector_.correctedPt(c.et, 0., std::abs(c.eta)); });
0204 }
0205
0206 auto hcalClustH = iEvent.put(hcalClusterer_.fetch(hcalCellsH), "hcalCalibrated");
0207
0208
0209 caloLinker_->run();
0210
0211 iEvent.put(caloLinker_->fetch(ecalClustH, hcalClustH), "uncalibrated");
0212
0213 if (hadCorrector_.valid()) {
0214 caloLinker_->correct([&](const l1tpf_calo::CombinedCluster &c) -> float {
0215 if (debug_)
0216 edm::LogWarning("L1TPFCaloProducer") << "raw linked cluster pt " << c.et << ", eta " << c.eta << ", phi "
0217 << c.phi << ", emPt " << c.ecal_et << "\n";
0218 return hadCorrector_.correctedPt(c.et, c.ecal_et, std::abs(c.eta));
0219 });
0220 }
0221
0222 std::unique_ptr<l1t::PFClusterCollection> clusters = caloLinker_->fetch(ecalClustH, hcalClustH);
0223 for (l1t::PFCluster &c : *clusters) {
0224 c.setPtError(resol_(c.pt(), std::abs(c.eta())));
0225 if (debug_)
0226 edm::LogWarning("L1TPFCaloProducer") << "calibrated linked cluster pt " << c.pt() << ", eta " << c.eta()
0227 << ", phi " << c.phi() << ", emPt " << c.emEt() << "\n";
0228 }
0229 iEvent.put(std::move(clusters), "calibrated");
0230
0231 ecalClusterer_.clear();
0232 hcalClusterer_.clear();
0233 caloLinker_->clear();
0234 }
0235
0236 void L1TPFCaloProducer::readHcalDigis_(edm::Event &iEvent, const edm::EventSetup &iSetup) {
0237 const auto &decoder = iSetup.getData(decoderTag_);
0238 edm::Handle<HcalTrigPrimDigiCollection> hcalTPs;
0239 for (const auto &token : hcalDigis_) {
0240 iEvent.getByToken(token, hcalTPs);
0241 for (const auto &itr : *hcalTPs) {
0242 HcalTrigTowerDetId id = itr.id();
0243 double et = decoder.hcaletValue(itr.id(), itr.t0());
0244 if (et <= 0)
0245 continue;
0246 float towerEta = l1t::CaloTools::towerEta(id.ieta());
0247 float towerPhi = l1t::CaloTools::towerPhi(id.ieta(), id.iphi());
0248 if (!hcalDigisBarrel_ && std::abs(towerEta) < 2)
0249 continue;
0250 if (!hcalDigisHF_ && std::abs(towerEta) > 2)
0251 continue;
0252 if (debug_)
0253 edm::LogWarning("L1TPFCaloProducer")
0254 << "adding HCal digi input pt " << et << ", eta " << towerEta << ", phi " << towerPhi << "\n";
0255 hcalClusterer_.add(et, towerEta, towerPhi);
0256 }
0257 }
0258 }
0259
0260 void L1TPFCaloProducer::readPhase2BarrelCaloTowers_(edm::Event &event, const edm::EventSetup &) {
0261 edm::Handle<l1tp2::CaloTowerCollection> towers;
0262 for (const auto &token : phase2barrelTowers_) {
0263 event.getByToken(token, towers);
0264 for (const auto &t : *towers) {
0265
0266 if ((int)t.towerIEta() == -1016 && (int)t.towerIPhi() == -962)
0267 continue;
0268 if (debug_ && (t.hcalTowerEt() > 0 || t.ecalTowerEt() > 0)) {
0269 edm::LogWarning("L1TPFCaloProducer")
0270 << "adding phase2 L1 CaloTower eta " << t.towerEta() << " phi " << t.towerPhi() << " ieta "
0271 << t.towerIEta() << " iphi " << t.towerIPhi() << " ecal " << t.ecalTowerEt() << " hcal "
0272 << t.hcalTowerEt() << "\n";
0273 }
0274 hcalClusterer_.add(t.hcalTowerEt(), t.towerEta(), t.towerPhi());
0275 ecalClusterer_.add(t.ecalTowerEt(), t.towerEta(), t.towerPhi());
0276 }
0277 }
0278 }
0279
0280 void L1TPFCaloProducer::readHcalHGCTowers_(edm::Event &iEvent, const edm::EventSetup &iSetup) {
0281 edm::Handle<l1t::HGCalTowerBxCollection> hgcTowers;
0282
0283 for (const auto &token : hcalHGCTowers_) {
0284 iEvent.getByToken(token, hgcTowers);
0285 for (auto it = hgcTowers->begin(0), ed = hgcTowers->end(0); it != ed; ++it) {
0286 if (debug_)
0287 edm::LogWarning("L1TPFCaloProducer")
0288 << "adding HGC Tower hadEt " << it->etHad() << ", emEt " << it->etEm() << ", pt " << it->pt() << ", eta "
0289 << it->eta() << ", phi " << it->phi() << "\n";
0290 hcalClusterer_.add(it->etHad(), it->eta(), it->phi());
0291 if (!hcalHGCTowersHadOnly_)
0292 ecalClusterer_.add(it->etEm(), it->eta(), it->phi());
0293 }
0294 }
0295 }
0296
0297
0298 #include "FWCore/Framework/interface/MakerMacros.h"
0299 DEFINE_FWK_MODULE(L1TPFCaloProducer);
