File indexing completed on 2024-05-10 02:20:59
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020 #include <memory>
0021
0022
0023 #include "FWCore/Framework/interface/Frameworkfwd.h"
0024 #include "FWCore/Framework/interface/stream/EDProducer.h"
0025
0026 #include "FWCore/Framework/interface/Event.h"
0027 #include "FWCore/Framework/interface/MakerMacros.h"
0028
0029 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0030 #include "FWCore/Utilities/interface/StreamID.h"
0031
0032 #include "DataFormats/L1TCalorimeterPhase2/interface/CaloCrystalCluster.h"
0033 #include "DataFormats/L1TCalorimeterPhase2/interface/CaloTower.h"
0034 #include "DataFormats/L1TCalorimeterPhase2/interface/CaloPFCluster.h"
0035 #include "DataFormats/L1TCalorimeterPhase2/interface/Phase2L1CaloJet.h"
0036 #include "DataFormats/L1Trigger/interface/EGamma.h"
0037 #include "DataFormats/L1THGCal/interface/HGCalTower.h"
0038 #include "DataFormats/HcalDigi/interface/HcalDigiCollections.h"
0039 #include "SimDataFormats/CaloHit/interface/PCaloHitContainer.h"
0040 #include "CalibFormats/CaloTPG/interface/CaloTPGTranscoder.h"
0041 #include "CalibFormats/CaloTPG/interface/CaloTPGRecord.h"
0042 #include "L1Trigger/L1TCalorimeter/interface/CaloTools.h"
0043
0044 #include "L1Trigger/L1CaloTrigger/interface/Phase2L1CaloJetEmulator.h"
0045 #include <ap_int.h>
0046 #include <cstdio>
0047 #include <fstream>
0048 #include <iomanip>
0049 #include <iostream>
0050 #include "TF1.h"
0051
0052
0053
0054
0055
0056 class Phase2L1CaloJetEmulator : public edm::stream::EDProducer<> {
0057 public:
0058 explicit Phase2L1CaloJetEmulator(const edm::ParameterSet&);
0059 ~Phase2L1CaloJetEmulator() override;
0060
0061 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0062
0063 private:
0064 void produce(edm::Event&, const edm::EventSetup&) override;
0065 float get_jet_pt_calibration(const float& jet_pt, const float& jet_eta) const;
0066 float get_tau_pt_calibration(const float& tau_pt, const float& tau_eta) const;
0067
0068
0069 edm::EDGetTokenT<l1tp2::CaloTowerCollection> caloTowerToken_;
0070 edm::EDGetTokenT<l1t::HGCalTowerBxCollection> hgcalTowerToken_;
0071 edm::EDGetTokenT<HcalTrigPrimDigiCollection> hfToken_;
0072 edm::ESGetToken<CaloTPGTranscoder, CaloTPGRecord> decoderTag_;
0073 std::vector<edm::ParameterSet> nHits_to_nvtx_params;
0074 std::vector<edm::ParameterSet> nvtx_to_PU_sub_params;
0075 std::map<std::string, TF1> nHits_to_nvtx_funcs;
0076 std::map<std::string, TF1> hgcalEM_nvtx_to_PU_sub_funcs;
0077 std::map<std::string, TF1> hgcalHad_nvtx_to_PU_sub_funcs;
0078 std::map<std::string, TF1> hf_nvtx_to_PU_sub_funcs;
0079 std::map<std::string, std::map<std::string, TF1>> all_nvtx_to_PU_sub_funcs;
0080
0081
0082 std::vector<double> jetPtBins;
0083 std::vector<double> absEtaBinsBarrel;
0084 std::vector<double> jetCalibrationsBarrel;
0085 std::vector<double> absEtaBinsHGCal;
0086 std::vector<double> jetCalibrationsHGCal;
0087 std::vector<double> absEtaBinsHF;
0088 std::vector<double> jetCalibrationsHF;
0089
0090
0091 std::vector<double> tauPtBins;
0092 std::vector<double> tauAbsEtaBinsBarrel;
0093 std::vector<double> tauCalibrationsBarrel;
0094 std::vector<double> tauAbsEtaBinsHGCal;
0095 std::vector<double> tauCalibrationsHGCal;
0096
0097
0098 std::vector<std::vector<double>> calibrationsBarrel;
0099 std::vector<std::vector<double>> calibrationsHGCal;
0100 std::vector<std::vector<double>> calibrationsHF;
0101
0102
0103 std::vector<std::vector<double>> tauPtCalibrationsBarrel;
0104 std::vector<std::vector<double>> tauPtCalibrationsHGCal;
0105 };
0106
0107
0108
0109
0110 Phase2L1CaloJetEmulator::Phase2L1CaloJetEmulator(const edm::ParameterSet& iConfig)
0111 : caloTowerToken_(consumes<l1tp2::CaloTowerCollection>(iConfig.getParameter<edm::InputTag>("gctFullTowers"))),
0112 hgcalTowerToken_(consumes<l1t::HGCalTowerBxCollection>(iConfig.getParameter<edm::InputTag>("hgcalTowers"))),
0113 hfToken_(consumes<HcalTrigPrimDigiCollection>(iConfig.getParameter<edm::InputTag>("hcalDigis"))),
0114 decoderTag_(esConsumes<CaloTPGTranscoder, CaloTPGRecord>(edm::ESInputTag("", ""))),
0115 nHits_to_nvtx_params(iConfig.getParameter<std::vector<edm::ParameterSet>>("nHits_to_nvtx_params")),
0116 nvtx_to_PU_sub_params(iConfig.getParameter<std::vector<edm::ParameterSet>>("nvtx_to_PU_sub_params")),
0117 jetPtBins(iConfig.getParameter<std::vector<double>>("jetPtBins")),
0118 absEtaBinsBarrel(iConfig.getParameter<std::vector<double>>("absEtaBinsBarrel")),
0119 jetCalibrationsBarrel(iConfig.getParameter<std::vector<double>>("jetCalibrationsBarrel")),
0120 absEtaBinsHGCal(iConfig.getParameter<std::vector<double>>("absEtaBinsHGCal")),
0121 jetCalibrationsHGCal(iConfig.getParameter<std::vector<double>>("jetCalibrationsHGCal")),
0122 absEtaBinsHF(iConfig.getParameter<std::vector<double>>("absEtaBinsHF")),
0123 jetCalibrationsHF(iConfig.getParameter<std::vector<double>>("jetCalibrationsHF")),
0124 tauPtBins(iConfig.getParameter<std::vector<double>>("tauPtBins")),
0125 tauAbsEtaBinsBarrel(iConfig.getParameter<std::vector<double>>("tauAbsEtaBinsBarrel")),
0126 tauCalibrationsBarrel(iConfig.getParameter<std::vector<double>>("tauCalibrationsBarrel")),
0127 tauAbsEtaBinsHGCal(iConfig.getParameter<std::vector<double>>("tauAbsEtaBinsHGCal")),
0128 tauCalibrationsHGCal(iConfig.getParameter<std::vector<double>>("tauCalibrationsHGCal")) {
0129 for (uint i = 0; i < nHits_to_nvtx_params.size(); i++) {
0130 edm::ParameterSet* pset = &nHits_to_nvtx_params.at(i);
0131 std::string calo = pset->getParameter<std::string>("fit");
0132 nHits_to_nvtx_funcs[calo.c_str()] = TF1(calo.c_str(), "[0] + [1] * x");
0133 nHits_to_nvtx_funcs[calo.c_str()].SetParameter(0, pset->getParameter<std::vector<double>>("nHits_params").at(0));
0134 nHits_to_nvtx_funcs[calo.c_str()].SetParameter(1, pset->getParameter<std::vector<double>>("nHits_params").at(1));
0135 }
0136 all_nvtx_to_PU_sub_funcs["hgcalEM"] = hgcalEM_nvtx_to_PU_sub_funcs;
0137 all_nvtx_to_PU_sub_funcs["hgcalHad"] = hgcalHad_nvtx_to_PU_sub_funcs;
0138 all_nvtx_to_PU_sub_funcs["hf"] = hf_nvtx_to_PU_sub_funcs;
0139 for (uint i = 0; i < nvtx_to_PU_sub_params.size(); i++) {
0140 edm::ParameterSet* pset = &nvtx_to_PU_sub_params.at(i);
0141 std::string calo = pset->getParameter<std::string>("calo");
0142 std::string iEta = pset->getParameter<std::string>("iEta");
0143 double p1 = pset->getParameter<std::vector<double>>("nvtx_params").at(0);
0144 double p2 = pset->getParameter<std::vector<double>>("nvtx_params").at(1);
0145
0146 all_nvtx_to_PU_sub_funcs[calo.c_str()][iEta.c_str()] = TF1(calo.c_str(), "[0] + [1] * x");
0147 all_nvtx_to_PU_sub_funcs[calo.c_str()][iEta.c_str()].SetParameter(0, p1);
0148 all_nvtx_to_PU_sub_funcs[calo.c_str()][iEta.c_str()].SetParameter(1, p2);
0149 }
0150
0151
0152
0153
0154
0155
0156 int index = 0;
0157 for (unsigned int abs_eta = 0; abs_eta < absEtaBinsBarrel.size() - 1; abs_eta++) {
0158 std::vector<double> pt_bin_calibs;
0159 for (unsigned int pt = 0; pt < jetPtBins.size() - 1; pt++) {
0160 pt_bin_calibs.push_back(jetCalibrationsBarrel.at(index));
0161 index++;
0162 }
0163 calibrationsBarrel.push_back(pt_bin_calibs);
0164 }
0165
0166 index = 0;
0167 for (unsigned int abs_eta = 0; abs_eta < absEtaBinsHGCal.size() - 1; abs_eta++) {
0168 std::vector<double> pt_bin_calibs;
0169 for (unsigned int pt = 0; pt < jetPtBins.size() - 1; pt++) {
0170 pt_bin_calibs.push_back(jetCalibrationsHGCal.at(index));
0171 index++;
0172 }
0173 calibrationsHGCal.push_back(pt_bin_calibs);
0174 }
0175
0176 index = 0;
0177 for (unsigned int abs_eta = 0; abs_eta < absEtaBinsHF.size() - 1; abs_eta++) {
0178 std::vector<double> pt_bin_calibs;
0179 for (unsigned int pt = 0; pt < jetPtBins.size() - 1; pt++) {
0180 pt_bin_calibs.push_back(jetCalibrationsHF.at(index));
0181 index++;
0182 }
0183 calibrationsHF.push_back(pt_bin_calibs);
0184 }
0185
0186
0187
0188
0189
0190
0191
0192
0193
0194 index = 0;
0195 for (unsigned int abs_eta = 0; abs_eta < tauAbsEtaBinsBarrel.size() - 1; abs_eta++) {
0196 std::vector<double> pt_bin_calibs;
0197 for (unsigned int pt = 0; pt < tauPtBins.size() - 1; pt++) {
0198 pt_bin_calibs.push_back(tauCalibrationsBarrel.at(index));
0199 index++;
0200 }
0201 tauPtCalibrationsBarrel.push_back(pt_bin_calibs);
0202 }
0203
0204 index = 0;
0205 for (unsigned int abs_eta = 0; abs_eta < tauAbsEtaBinsHGCal.size() - 1; abs_eta++) {
0206 std::vector<double> pt_bin_calibs;
0207 for (unsigned int pt = 0; pt < tauPtBins.size() - 1; pt++) {
0208 pt_bin_calibs.push_back(tauCalibrationsHGCal.at(index));
0209 index++;
0210 }
0211 tauPtCalibrationsHGCal.push_back(pt_bin_calibs);
0212 }
0213
0214 produces<l1tp2::Phase2L1CaloJetCollection>("GCTJet");
0215 }
0216
0217 Phase2L1CaloJetEmulator::~Phase2L1CaloJetEmulator() {}
0218
0219
0220
0221
0222
0223
0224 void Phase2L1CaloJetEmulator::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0225 using namespace edm;
0226 std::unique_ptr<l1tp2::Phase2L1CaloJetCollection> jetCands(make_unique<l1tp2::Phase2L1CaloJetCollection>());
0227
0228
0229 edm::Handle<std::vector<l1tp2::CaloTower>> caloTowerCollection;
0230 if (!iEvent.getByToken(caloTowerToken_, caloTowerCollection))
0231 edm::LogError("Phase2L1CaloJetEmulator") << "Failed to get towers from caloTowerCollection!";
0232
0233 iEvent.getByToken(caloTowerToken_, caloTowerCollection);
0234 float GCTintTowers[nBarrelEta][nBarrelPhi];
0235 float realEta[nBarrelEta][nBarrelPhi];
0236 float realPhi[nBarrelEta][nBarrelPhi];
0237 for (const l1tp2::CaloTower& i : *caloTowerCollection) {
0238 int ieta = i.towerIEta();
0239 int iphi = i.towerIPhi();
0240 if (i.ecalTowerEt() > 1.)
0241 GCTintTowers[ieta][iphi] = i.ecalTowerEt();
0242 else
0243 GCTintTowers[ieta][iphi] = 0;
0244 realEta[ieta][iphi] = i.towerEta();
0245 realPhi[ieta][iphi] = i.towerPhi();
0246 }
0247
0248 float temporary[nBarrelEta / 2][nBarrelPhi];
0249
0250
0251 edm::Handle<l1t::HGCalTowerBxCollection> hgcalTowerCollection;
0252 if (!iEvent.getByToken(hgcalTowerToken_, hgcalTowerCollection))
0253 edm::LogError("Phase2L1CaloJetEmulator") << "Failed to get towers from hgcalTowerCollection!";
0254 l1t::HGCalTowerBxCollection hgcalTowerColl;
0255 iEvent.getByToken(hgcalTowerToken_, hgcalTowerCollection);
0256 hgcalTowerColl = (*hgcalTowerCollection.product());
0257
0258 edm::Handle<HcalTrigPrimDigiCollection> hfHandle;
0259 if (!iEvent.getByToken(hfToken_, hfHandle))
0260 edm::LogError("Phase2L1CaloJetEmulator") << "Failed to get HcalTrigPrimDigi for HF!";
0261 iEvent.getByToken(hfToken_, hfHandle);
0262
0263 int i_hgcalEM_hits_leq_threshold = 0;
0264 int i_hgcalHad_hits_leq_threshold = 0;
0265 int i_hf_hits_leq_threshold = 0;
0266 for (auto it = hgcalTowerColl.begin(0); it != hgcalTowerColl.end(0); it++) {
0267 if (it->etEm() <= 1.75 && it->etEm() >= 1.25) {
0268 i_hgcalEM_hits_leq_threshold++;
0269 }
0270 if (it->etHad() <= 1.25 && it->etHad() >= 0.75) {
0271 i_hgcalHad_hits_leq_threshold++;
0272 }
0273 }
0274 const auto& decoder = iSetup.getData(decoderTag_);
0275 for (const auto& hit : *hfHandle.product()) {
0276 double et = decoder.hcaletValue(hit.id(), hit.t0());
0277 if (abs(hit.id().ieta()) < l1t::CaloTools::kHFBegin)
0278 continue;
0279 if (abs(hit.id().ieta()) > l1t::CaloTools::kHFEnd)
0280 continue;
0281 if (et <= 15.0 && et >= 10.0)
0282 i_hf_hits_leq_threshold++;
0283 }
0284
0285 double hgcalEM_nvtx = nHits_to_nvtx_funcs["hgcalEM"].Eval(i_hgcalEM_hits_leq_threshold);
0286 if (hgcalEM_nvtx < 0)
0287 hgcalEM_nvtx = 0;
0288 double hgcalHad_nvtx = nHits_to_nvtx_funcs["hgcalHad"].Eval(i_hgcalHad_hits_leq_threshold);
0289 if (hgcalHad_nvtx < 0)
0290 hgcalHad_nvtx = 0;
0291 double hf_nvtx = nHits_to_nvtx_funcs["hf"].Eval(i_hf_hits_leq_threshold);
0292 if (hf_nvtx < 0)
0293 hf_nvtx = 0;
0294 double EstimatedNvtx = (hgcalEM_nvtx + hgcalHad_nvtx + hf_nvtx) / 3.;
0295
0296
0297 float hgcalTowers[nHgcalEta][nHgcalPhi];
0298 float hgcalEta[nHgcalEta][nHgcalPhi];
0299 float hgcalPhi[nHgcalEta][nHgcalPhi];
0300
0301 for (int iphi = 0; iphi < nHgcalPhi; iphi++) {
0302 for (int ieta = 0; ieta < nHgcalEta; ieta++) {
0303 hgcalTowers[ieta][iphi] = 0;
0304 if (ieta < nHgcalEta / 2)
0305 hgcalEta[ieta][iphi] = -3.045 + ieta * 0.087 + 0.0435;
0306 else
0307 hgcalEta[ieta][iphi] = 1.479 + (ieta - nHgcalEta / 2) * 0.087 + 0.0435;
0308 hgcalPhi[ieta][iphi] = -M_PI + (iphi * 2 * M_PI / nHgcalPhi) + (M_PI / nHgcalPhi);
0309 }
0310 }
0311
0312 for (auto it = hgcalTowerColl.begin(0); it != hgcalTowerColl.end(0); it++) {
0313 float eta = it->eta();
0314 int ieta;
0315 if (eta < 0)
0316 ieta = 19 - it->id().iEta();
0317 else
0318 ieta = 20 + it->id().iEta();
0319 if (eta > 1.479)
0320 ieta = ieta - 4;
0321 int iphi = it->id().iPhi();
0322
0323 float hgcal_etEm = it->etEm();
0324 float hgcal_etHad = it->etHad();
0325 std::string etaKey = "";
0326 if (abs(eta) <= 1.8)
0327 etaKey = "er1p4to1p8";
0328 else if (abs(eta) <= 2.1 && abs(eta) > 1.8)
0329 etaKey = "er1p8to2p1";
0330 else if (abs(eta) <= 2.4 && abs(eta) > 2.1)
0331 etaKey = "er2p1to2p4";
0332 else if (abs(eta) <= 2.7 && abs(eta) > 2.4)
0333 etaKey = "er2p4to2p7";
0334 else if (abs(eta) <= 3.1 && abs(eta) > 2.7)
0335 etaKey = "er2p7to3p1";
0336 if (!etaKey.empty()) {
0337 hgcal_etEm = it->etEm() - all_nvtx_to_PU_sub_funcs["hgcalEM"][etaKey].Eval(EstimatedNvtx);
0338 hgcal_etHad = it->etHad() - all_nvtx_to_PU_sub_funcs["hgcalHad"][etaKey].Eval(EstimatedNvtx);
0339 }
0340
0341 if (hgcal_etEm < 0)
0342 hgcal_etEm = 0;
0343 if (hgcal_etHad < 0)
0344 hgcal_etHad = 0;
0345 if ((it->etEm() + it->etHad() > 1.) && abs(eta) > 1.479)
0346 hgcalTowers[ieta][iphi] = hgcal_etEm + hgcal_etHad;
0347 }
0348
0349 float temporary_hgcal[nHgcalEta / 2][nHgcalPhi];
0350
0351
0352 float hfTowers[nHfEta][nHfPhi];
0353 float hfEta[nHfEta][nHfPhi];
0354 float hfPhi[nHfEta][nHfPhi];
0355 for (int iphi = 0; iphi < nHfPhi; iphi++) {
0356 for (int ieta = 0; ieta < nHfEta; ieta++) {
0357 hfTowers[ieta][iphi] = 0;
0358 int temp;
0359 if (ieta < nHfEta / 2)
0360 temp = ieta - l1t::CaloTools::kHFEnd;
0361 else
0362 temp = ieta - nHfEta / 2 + l1t::CaloTools::kHFBegin + 1;
0363 hfEta[ieta][iphi] = l1t::CaloTools::towerEta(temp);
0364 hfPhi[ieta][iphi] = -M_PI + (iphi * 2 * M_PI / nHfPhi) + (M_PI / nHfPhi);
0365 }
0366 }
0367
0368 for (const auto& hit : *hfHandle.product()) {
0369 double et = decoder.hcaletValue(hit.id(), hit.t0());
0370 int ieta = 0;
0371 if (abs(hit.id().ieta()) < l1t::CaloTools::kHFBegin)
0372 continue;
0373 if (abs(hit.id().ieta()) > l1t::CaloTools::kHFEnd)
0374 continue;
0375 if (hit.id().ieta() <= -(l1t::CaloTools::kHFBegin + 1)) {
0376 ieta = hit.id().ieta() + l1t::CaloTools::kHFEnd;
0377 } else if (hit.id().ieta() >= (l1t::CaloTools::kHFBegin + 1)) {
0378 ieta = nHfEta / 2 + (hit.id().ieta() - (l1t::CaloTools::kHFBegin + 1));
0379 }
0380 int iphi = 0;
0381 if (hit.id().iphi() <= nHfPhi / 2)
0382 iphi = hit.id().iphi() + (nHfPhi / 2 - 1);
0383 else if (hit.id().iphi() > nHfPhi / 2)
0384 iphi = hit.id().iphi() - (nHfPhi / 2 + 1);
0385 if (abs(hit.id().ieta()) <= 33 && abs(hit.id().ieta()) >= 29)
0386 et = et - all_nvtx_to_PU_sub_funcs["hf"]["er29to33"].Eval(EstimatedNvtx);
0387 if (abs(hit.id().ieta()) <= 37 && abs(hit.id().ieta()) >= 34)
0388 et = et - all_nvtx_to_PU_sub_funcs["hf"]["er34to37"].Eval(EstimatedNvtx);
0389 if (abs(hit.id().ieta()) <= 41 && abs(hit.id().ieta()) >= 38)
0390 et = et - all_nvtx_to_PU_sub_funcs["hf"]["er38to41"].Eval(EstimatedNvtx);
0391 if (et < 0)
0392 et = 0;
0393 if (et > 1.)
0394 hfTowers[ieta][iphi] = et;
0395 }
0396
0397 float temporary_hf[nHfEta / 2][nHfPhi];
0398
0399
0400
0401
0402
0403 vector<l1tp2::Phase2L1CaloJet> halfBarrelJets, halfHgcalJets, halfHfJets;
0404 halfBarrelJets.clear();
0405 halfHgcalJets.clear();
0406 halfHfJets.clear();
0407 vector<l1tp2::Phase2L1CaloJet> allJets;
0408 allJets.clear();
0409
0410 for (int k = 0; k < 2; k++) {
0411 halfBarrelJets.clear();
0412 halfHgcalJets.clear();
0413 halfHfJets.clear();
0414 gctobj::jetInfo jet[3 * nJets];
0415
0416
0417 for (int iphi = 0; iphi < nBarrelPhi; iphi++) {
0418 for (int ieta = 0; ieta < nBarrelEta / 2; ieta++) {
0419 if (k == 0)
0420 temporary[ieta][iphi] = GCTintTowers[ieta][iphi];
0421 else
0422 temporary[ieta][iphi] = GCTintTowers[nBarrelEta / 2 + ieta][iphi];
0423 }
0424 }
0425
0426 gctobj::GCTsupertower_t tempST[nSTEta][nSTPhi];
0427 gctobj::makeST(temporary, tempST);
0428 float TTseedThresholdBarrel = 5.;
0429
0430 for (int i = 0; i < nJets; i++) {
0431 jet[i] = gctobj::getRegion(tempST, TTseedThresholdBarrel);
0432 l1tp2::Phase2L1CaloJet tempJet;
0433 int gctjeteta = jet[i].etaCenter;
0434 int gctjetphi = jet[i].phiCenter;
0435 tempJet.setJetIEta(gctjeteta + k * nBarrelEta / 2);
0436 tempJet.setJetIPhi(gctjetphi);
0437 float jeteta = realEta[gctjeteta + k * nBarrelEta / 2][gctjetphi];
0438 float jetphi = realPhi[gctjeteta + k * nBarrelEta / 2][gctjetphi];
0439 tempJet.setJetEta(jeteta);
0440 tempJet.setJetPhi(jetphi);
0441 tempJet.setJetEt(get_jet_pt_calibration(jet[i].energy, jeteta));
0442 tempJet.setTauEt(get_tau_pt_calibration(jet[i].tauEt, jeteta));
0443 tempJet.setTowerEt(jet[i].energyMax);
0444 int gcttowereta = jet[i].etaMax;
0445 int gcttowerphi = jet[i].phiMax;
0446 tempJet.setTowerIEta(gcttowereta + k * nBarrelEta / 2);
0447 tempJet.setTowerIPhi(gcttowerphi);
0448 float towereta = realEta[gcttowereta + k * nBarrelEta / 2][gcttowerphi];
0449 float towerphi = realPhi[gcttowereta + k * nBarrelEta / 2][gcttowerphi];
0450 tempJet.setTowerEta(towereta);
0451 tempJet.setTowerPhi(towerphi);
0452 reco::Candidate::PolarLorentzVector tempJetp4;
0453 tempJetp4.SetPt(tempJet.jetEt());
0454 tempJetp4.SetEta(tempJet.jetEta());
0455 tempJetp4.SetPhi(tempJet.jetPhi());
0456 tempJetp4.SetM(0.);
0457 tempJet.setP4(tempJetp4);
0458
0459 if (jet[i].energy > 0.)
0460 halfBarrelJets.push_back(tempJet);
0461 }
0462
0463
0464 for (int iphi = 0; iphi < nHgcalPhi; iphi++) {
0465 for (int ieta = 0; ieta < nHgcalEta / 2; ieta++) {
0466 if (k == 0)
0467 temporary_hgcal[ieta][iphi] = hgcalTowers[ieta][iphi];
0468 else
0469 temporary_hgcal[ieta][iphi] = hgcalTowers[nHgcalEta / 2 + ieta][iphi];
0470 }
0471 }
0472
0473 gctobj::GCTsupertower_t tempST_hgcal[nSTEta][nSTPhi];
0474 gctobj::makeST_hgcal(temporary_hgcal, tempST_hgcal);
0475 float TTseedThresholdEndcap = 3.;
0476 for (int i = nJets; i < 2 * nJets; i++) {
0477 jet[i] = gctobj::getRegion(tempST_hgcal, TTseedThresholdEndcap);
0478 l1tp2::Phase2L1CaloJet tempJet;
0479 int hgcaljeteta = jet[i].etaCenter;
0480 int hgcaljetphi = jet[i].phiCenter;
0481 tempJet.setJetIEta(hgcaljeteta + k * nHgcalEta / 2);
0482 tempJet.setJetIPhi(hgcaljetphi);
0483 float jeteta = hgcalEta[hgcaljeteta + k * nHgcalEta / 2][hgcaljetphi];
0484 float jetphi = hgcalPhi[hgcaljeteta + k * nHgcalEta / 2][hgcaljetphi];
0485 tempJet.setJetEta(jeteta);
0486 tempJet.setJetPhi(jetphi);
0487 tempJet.setJetEt(get_jet_pt_calibration(jet[i].energy, jeteta));
0488 tempJet.setTauEt(get_tau_pt_calibration(jet[i].tauEt, jeteta));
0489 tempJet.setTowerEt(jet[i].energyMax);
0490 int hgcaltowereta = jet[i].etaMax;
0491 int hgcaltowerphi = jet[i].phiMax;
0492 tempJet.setTowerIEta(hgcaltowereta + k * nHgcalEta / 2);
0493 tempJet.setTowerIPhi(hgcaltowerphi);
0494 float towereta = hgcalEta[hgcaltowereta + k * nHgcalEta / 2][hgcaltowerphi];
0495 float towerphi = hgcalPhi[hgcaltowereta + k * nHgcalEta / 2][hgcaltowerphi];
0496 tempJet.setTowerEta(towereta);
0497 tempJet.setTowerPhi(towerphi);
0498 reco::Candidate::PolarLorentzVector tempJetp4;
0499 tempJetp4.SetPt(tempJet.jetEt());
0500 tempJetp4.SetEta(tempJet.jetEta());
0501 tempJetp4.SetPhi(tempJet.jetPhi());
0502 tempJetp4.SetM(0.);
0503 tempJet.setP4(tempJetp4);
0504
0505 if (jet[i].energy > 0.)
0506 halfHgcalJets.push_back(tempJet);
0507 }
0508
0509
0510 for (int iphi = 0; iphi < nHfPhi; iphi++) {
0511 for (int ieta = 0; ieta < nHfEta / 2; ieta++) {
0512 if (k == 0)
0513 temporary_hf[ieta][iphi] = hfTowers[ieta][iphi];
0514 else
0515 temporary_hf[ieta][iphi] = hfTowers[nHfEta / 2 + ieta][iphi];
0516 }
0517 }
0518
0519 gctobj::GCTsupertower_t tempST_hf[nSTEta][nSTPhi];
0520 gctobj::makeST_hf(temporary_hf, tempST_hf);
0521 float TTseedThresholdHF = 5.;
0522 for (int i = 2 * nJets; i < 3 * nJets; i++) {
0523 jet[i] = gctobj::getRegion(tempST_hf, TTseedThresholdHF);
0524 l1tp2::Phase2L1CaloJet tempJet;
0525 int hfjeteta = jet[i].etaCenter;
0526 int hfjetphi = jet[i].phiCenter;
0527 tempJet.setJetIEta(hfjeteta + k * nHfEta / 2);
0528 tempJet.setJetIPhi(hfjetphi);
0529 float jeteta = hfEta[hfjeteta + k * nHfEta / 2][hfjetphi];
0530 float jetphi = hfPhi[hfjeteta + k * nHfEta / 2][hfjetphi];
0531 tempJet.setJetEta(jeteta);
0532 tempJet.setJetPhi(jetphi);
0533 tempJet.setJetEt(get_jet_pt_calibration(jet[i].energy, jeteta));
0534 tempJet.setTauEt(get_tau_pt_calibration(jet[i].tauEt, jeteta));
0535 tempJet.setTowerEt(jet[i].energyMax);
0536 int hftowereta = jet[i].etaMax;
0537 int hftowerphi = jet[i].phiMax;
0538 tempJet.setTowerIEta(hftowereta + k * nHfEta / 2);
0539 tempJet.setTowerIPhi(hftowerphi);
0540 float towereta = hfEta[hftowereta + k * nHfEta / 2][hftowerphi];
0541 float towerphi = hfPhi[hftowereta + k * nHfEta / 2][hftowerphi];
0542 tempJet.setTowerEta(towereta);
0543 tempJet.setTowerPhi(towerphi);
0544 reco::Candidate::PolarLorentzVector tempJetp4;
0545 tempJetp4.SetPt(tempJet.jetEt());
0546 tempJetp4.SetEta(tempJet.jetEta());
0547 tempJetp4.SetPhi(tempJet.jetPhi());
0548 tempJetp4.SetM(0.);
0549 tempJet.setP4(tempJetp4);
0550
0551 if (jet[i].energy > 0.)
0552 halfHfJets.push_back(tempJet);
0553 }
0554
0555
0556
0557
0558
0559
0560 for (size_t i = 0; i < halfHgcalJets.size(); i++) {
0561 if (halfHgcalJets.at(i).jetIEta() >= (nHgcalEta / 2 - 2) && halfHgcalJets.at(i).jetIEta() < (nHgcalEta / 2 + 2)) {
0562 float hgcal_ieta = k * nBarrelEta + halfHgcalJets.at(i).jetIEta();
0563 for (size_t j = 0; j < halfBarrelJets.size(); j++) {
0564 float barrel_ieta = nHgcalEta / 2 + halfBarrelJets.at(j).jetIEta();
0565 if (abs(barrel_ieta - hgcal_ieta) <= 2 &&
0566 abs(halfBarrelJets.at(j).jetIPhi() - halfHgcalJets.at(i).jetIPhi()) <= 2) {
0567 float totalet = halfBarrelJets.at(j).jetEt() + halfHgcalJets.at(i).jetEt();
0568 float totalTauEt = halfBarrelJets.at(j).tauEt() + halfHgcalJets.at(i).tauEt();
0569 if (halfBarrelJets.at(j).jetEt() > halfHgcalJets.at(i).jetEt()) {
0570 halfHgcalJets.at(i).setJetEt(0.);
0571 halfHgcalJets.at(i).setTauEt(0.);
0572 halfBarrelJets.at(j).setJetEt(totalet);
0573 halfBarrelJets.at(j).setTauEt(totalTauEt);
0574 reco::Candidate::PolarLorentzVector tempJetp4;
0575 tempJetp4.SetPt(totalet);
0576 tempJetp4.SetEta(halfBarrelJets.at(j).jetEta());
0577 tempJetp4.SetPhi(halfBarrelJets.at(j).jetPhi());
0578 tempJetp4.SetM(0.);
0579 halfBarrelJets.at(j).setP4(tempJetp4);
0580 } else {
0581 halfHgcalJets.at(i).setJetEt(totalet);
0582 halfHgcalJets.at(i).setTauEt(totalTauEt);
0583 halfBarrelJets.at(j).setJetEt(0.);
0584 halfBarrelJets.at(j).setTauEt(0.);
0585 reco::Candidate::PolarLorentzVector tempJetp4;
0586 tempJetp4.SetPt(totalet);
0587 tempJetp4.SetEta(halfHgcalJets.at(i).jetEta());
0588 tempJetp4.SetPhi(halfHgcalJets.at(i).jetPhi());
0589 tempJetp4.SetM(0.);
0590 halfHgcalJets.at(i).setP4(tempJetp4);
0591 }
0592 }
0593 }
0594 } else if (halfHgcalJets.at(i).jetIEta() < 2 || halfHgcalJets.at(i).jetIEta() >= (nHgcalEta - 2)) {
0595 float hgcal_ieta = k * nBarrelEta + nHfEta / 2 + halfHgcalJets.at(i).jetIEta();
0596 for (size_t j = 0; j < halfHfJets.size(); j++) {
0597 float hf_ieta = k * nBarrelEta + k * nHgcalEta + halfHfJets.at(j).jetIEta();
0598 if (abs(hgcal_ieta - hf_ieta) < 3 && abs(halfHfJets.at(j).jetIPhi() - halfHgcalJets.at(i).jetIPhi()) < 3) {
0599 float totalet = halfHfJets.at(j).jetEt() + halfHgcalJets.at(i).jetEt();
0600 float totalTauEt = halfHfJets.at(j).tauEt() + halfHgcalJets.at(i).tauEt();
0601 if (halfHfJets.at(j).jetEt() > halfHgcalJets.at(i).jetEt()) {
0602 halfHgcalJets.at(i).setJetEt(0.);
0603 halfHgcalJets.at(i).setTauEt(0.);
0604 halfHfJets.at(j).setJetEt(totalet);
0605 halfHfJets.at(j).setTauEt(totalTauEt);
0606 reco::Candidate::PolarLorentzVector tempJetp4;
0607 tempJetp4.SetPt(totalet);
0608 tempJetp4.SetEta(halfHfJets.at(j).jetEta());
0609 tempJetp4.SetPhi(halfHfJets.at(j).jetPhi());
0610 tempJetp4.SetM(0.);
0611 halfHfJets.at(j).setP4(tempJetp4);
0612 } else {
0613 halfHgcalJets.at(i).setJetEt(totalet);
0614 halfHgcalJets.at(i).setTauEt(totalTauEt);
0615 halfHfJets.at(j).setJetEt(0.);
0616 halfHfJets.at(j).setTauEt(0.);
0617 reco::Candidate::PolarLorentzVector tempJetp4;
0618 tempJetp4.SetPt(totalet);
0619 tempJetp4.SetEta(halfHgcalJets.at(i).jetEta());
0620 tempJetp4.SetPhi(halfHgcalJets.at(i).jetPhi());
0621 tempJetp4.SetM(0.);
0622 halfHgcalJets.at(i).setP4(tempJetp4);
0623 }
0624 }
0625 }
0626 }
0627 }
0628
0629
0630
0631 vector<l1tp2::Phase2L1CaloJet> halfAllJets;
0632 halfAllJets.clear();
0633
0634 std::sort(halfBarrelJets.begin(), halfBarrelJets.end(), gctobj::compareByEt);
0635 for (size_t i = 0; i < halfBarrelJets.size(); i++) {
0636 if (halfBarrelJets.at(i).jetEt() > 0. && i < 6)
0637 halfAllJets.push_back(halfBarrelJets.at(i));
0638 }
0639
0640 std::sort(halfHgcalJets.begin(), halfHgcalJets.end(), gctobj::compareByEt);
0641 for (size_t i = 0; i < halfHgcalJets.size(); i++) {
0642 if (halfHgcalJets.at(i).jetEt() > 0. && i < 6)
0643 halfAllJets.push_back(halfHgcalJets.at(i));
0644 }
0645
0646 std::sort(halfHfJets.begin(), halfHfJets.end(), gctobj::compareByEt);
0647 for (size_t i = 0; i < halfHfJets.size(); i++) {
0648 if (halfHfJets.at(i).jetEt() > 0. && i < 6)
0649 halfAllJets.push_back(halfHfJets.at(i));
0650 }
0651
0652 std::sort(halfAllJets.begin(), halfAllJets.end(), gctobj::compareByEt);
0653 for (size_t i = 0; i < halfAllJets.size(); i++) {
0654 if (halfAllJets.at(i).jetEt() > 0. && i < 6)
0655 allJets.push_back(halfAllJets.at(i));
0656 }
0657 }
0658
0659 std::sort(allJets.begin(), allJets.end(), gctobj::compareByEt);
0660 for (size_t i = 0; i < allJets.size(); i++) {
0661 jetCands->push_back(allJets.at(i));
0662 }
0663
0664 iEvent.put(std::move(jetCands), "GCTJet");
0665 }
0666
0667
0668 float Phase2L1CaloJetEmulator::get_jet_pt_calibration(const float& jet_pt, const float& jet_eta) const {
0669 float abs_eta = std::abs(jet_eta);
0670 float tmp_jet_pt = jet_pt;
0671 if (tmp_jet_pt > 499)
0672 tmp_jet_pt = 499;
0673
0674
0675
0676 size_t eta_index = 0;
0677 size_t pt_index = 0;
0678 float calib = 1.0;
0679 if (abs_eta <= 1.5) {
0680
0681 for (unsigned int i = 1; i < absEtaBinsBarrel.size(); i++) {
0682 if (abs_eta <= absEtaBinsBarrel.at(i))
0683 break;
0684 eta_index++;
0685 }
0686
0687 for (unsigned int i = 1; i < jetPtBins.size(); i++) {
0688 if (tmp_jet_pt <= jetPtBins.at(i))
0689 break;
0690 pt_index++;
0691 }
0692 calib = calibrationsBarrel[eta_index][pt_index];
0693 }
0694 else if (abs_eta <= 3.0)
0695 {
0696
0697 for (unsigned int i = 1; i < absEtaBinsHGCal.size(); i++) {
0698 if (abs_eta <= absEtaBinsHGCal.at(i))
0699 break;
0700 eta_index++;
0701 }
0702
0703 for (unsigned int i = 1; i < jetPtBins.size(); i++) {
0704 if (tmp_jet_pt <= jetPtBins.at(i))
0705 break;
0706 pt_index++;
0707 }
0708 calib = calibrationsHGCal[eta_index][pt_index];
0709 }
0710 else
0711 {
0712
0713 for (unsigned int i = 1; i < absEtaBinsHF.size(); i++) {
0714 if (abs_eta <= absEtaBinsHF.at(i))
0715 break;
0716 eta_index++;
0717 }
0718
0719 for (unsigned int i = 1; i < jetPtBins.size(); i++) {
0720 if (tmp_jet_pt <= jetPtBins.at(i))
0721 break;
0722 pt_index++;
0723 }
0724 calib = calibrationsHF[eta_index][pt_index];
0725 }
0726
0727 return jet_pt * calib;
0728 }
0729
0730
0731 float Phase2L1CaloJetEmulator::get_tau_pt_calibration(const float& tau_pt, const float& tau_eta) const {
0732 float abs_eta = std::abs(tau_eta);
0733 float tmp_tau_pt = tau_pt;
0734 if (tmp_tau_pt > 199)
0735 tmp_tau_pt = 199;
0736
0737
0738
0739 size_t eta_index = 0;
0740 size_t pt_index = 0;
0741 float calib = 1.0;
0742 if (abs_eta <= 1.5) {
0743
0744 for (unsigned int i = 1; i < tauAbsEtaBinsBarrel.size(); i++) {
0745 if (abs_eta <= tauAbsEtaBinsBarrel.at(i))
0746 break;
0747 eta_index++;
0748 }
0749
0750 for (unsigned int i = 1; i < tauPtBins.size(); i++) {
0751 if (tmp_tau_pt <= tauPtBins.at(i))
0752 break;
0753 pt_index++;
0754 }
0755 calib = tauPtCalibrationsBarrel[eta_index][pt_index];
0756 }
0757 else if (abs_eta <= 3.0)
0758 {
0759
0760 for (unsigned int i = 1; i < tauAbsEtaBinsHGCal.size(); i++) {
0761 if (abs_eta <= tauAbsEtaBinsHGCal.at(i))
0762 break;
0763 eta_index++;
0764 }
0765
0766 for (unsigned int i = 1; i < tauPtBins.size(); i++) {
0767 if (tmp_tau_pt <= tauPtBins.at(i))
0768 break;
0769 pt_index++;
0770 }
0771 calib = tauPtCalibrationsHGCal[eta_index][pt_index];
0772 }
0773
0774 return tau_pt * calib;
0775 }
0776
0777
0778 void Phase2L1CaloJetEmulator::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0779 edm::ParameterSetDescription desc;
0780 desc.add<edm::InputTag>("gctFullTowers", edm::InputTag("l1tPhase2L1CaloEGammaEmulator", "GCTFullTowers"));
0781 desc.add<edm::InputTag>("hgcalTowers", edm::InputTag("l1tHGCalTowerProducer", "HGCalTowerProcessor"));
0782 desc.add<edm::InputTag>("hcalDigis", edm::InputTag("simHcalTriggerPrimitiveDigis"));
0783
0784 edm::ParameterSetDescription nHits_params_validator;
0785 nHits_params_validator.add<string>("fit", "type");
0786 nHits_params_validator.add<vector<double>>("nHits_params", {1., 1.});
0787 std::vector<edm::ParameterSet> nHits_params_default;
0788 edm::ParameterSet nHits_params1;
0789 nHits_params1.addParameter<string>("fit", "hgcalEM");
0790 nHits_params1.addParameter<vector<double>>("nHits_params", {157.522, 0.090});
0791 nHits_params_default.push_back(nHits_params1);
0792 edm::ParameterSet nHits_params2;
0793 nHits_params2.addParameter<string>("fit", "hgcalHad");
0794 nHits_params2.addParameter<vector<double>>("nHits_params", {159.295, 0.178});
0795 nHits_params_default.push_back(nHits_params2);
0796 edm::ParameterSet nHits_params3;
0797 nHits_params3.addParameter<string>("fit", "hf");
0798 nHits_params3.addParameter<vector<double>>("nHits_params", {165.706, 0.153});
0799 nHits_params_default.push_back(nHits_params3);
0800 desc.addVPSet("nHits_to_nvtx_params", nHits_params_validator, nHits_params_default);
0801
0802 edm::ParameterSetDescription nvtx_params_validator;
0803 nvtx_params_validator.add<string>("calo", "type");
0804 nvtx_params_validator.add<string>("iEta", "etaregion");
0805 nvtx_params_validator.add<vector<double>>("nvtx_params", {1., 1.});
0806 std::vector<edm::ParameterSet> nvtx_params_default;
0807 edm::ParameterSet nvtx_params1;
0808 nvtx_params1.addParameter<string>("calo", "hgcalEM");
0809 nvtx_params1.addParameter<string>("iEta", "er1p4to1p8");
0810 nvtx_params1.addParameter<vector<double>>("nvtx_params", {-0.011772, 0.004142});
0811 nvtx_params_default.push_back(nvtx_params1);
0812 edm::ParameterSet nvtx_params2;
0813 nvtx_params2.addParameter<string>("calo", "hgcalEM");
0814 nvtx_params2.addParameter<string>("iEta", "er1p8to2p1");
0815 nvtx_params2.addParameter<vector<double>>("nvtx_params", {-0.015488, 0.005410});
0816 nvtx_params_default.push_back(nvtx_params2);
0817 edm::ParameterSet nvtx_params3;
0818 nvtx_params3.addParameter<string>("calo", "hgcalEM");
0819 nvtx_params3.addParameter<string>("iEta", "er2p1to2p4");
0820 nvtx_params3.addParameter<vector<double>>("nvtx_params", {-0.021150, 0.006078});
0821 nvtx_params_default.push_back(nvtx_params3);
0822 edm::ParameterSet nvtx_params4;
0823 nvtx_params4.addParameter<string>("calo", "hgcalEM");
0824 nvtx_params4.addParameter<string>("iEta", "er2p4to2p7");
0825 nvtx_params4.addParameter<vector<double>>("nvtx_params", {-0.015705, 0.005339});
0826 nvtx_params_default.push_back(nvtx_params4);
0827 edm::ParameterSet nvtx_params5;
0828 nvtx_params5.addParameter<string>("calo", "hgcalEM");
0829 nvtx_params5.addParameter<string>("iEta", "er2p7to3p1");
0830 nvtx_params5.addParameter<vector<double>>("nvtx_params", {-0.018492, 0.005620});
0831 nvtx_params_default.push_back(nvtx_params5);
0832 edm::ParameterSet nvtx_params6;
0833 nvtx_params6.addParameter<string>("calo", "hgcalHad");
0834 nvtx_params6.addParameter<string>("iEta", "er1p4to1p8");
0835 nvtx_params6.addParameter<vector<double>>("nvtx_params", {0.005675, 0.000615});
0836 nvtx_params_default.push_back(nvtx_params6);
0837 edm::ParameterSet nvtx_params7;
0838 nvtx_params7.addParameter<string>("calo", "hgcalHad");
0839 nvtx_params7.addParameter<string>("iEta", "er1p8to2p1");
0840 nvtx_params7.addParameter<vector<double>>("nvtx_params", {0.004560, 0.001099});
0841 nvtx_params_default.push_back(nvtx_params7);
0842 edm::ParameterSet nvtx_params8;
0843 nvtx_params8.addParameter<string>("calo", "hgcalHad");
0844 nvtx_params8.addParameter<string>("iEta", "er2p1to2p4");
0845 nvtx_params8.addParameter<vector<double>>("nvtx_params", {0.000036, 0.001608});
0846 nvtx_params_default.push_back(nvtx_params8);
0847 edm::ParameterSet nvtx_params9;
0848 nvtx_params9.addParameter<string>("calo", "hgcalHad");
0849 nvtx_params9.addParameter<string>("iEta", "er2p4to2p7");
0850 nvtx_params9.addParameter<vector<double>>("nvtx_params", {0.000869, 0.001754});
0851 nvtx_params_default.push_back(nvtx_params9);
0852 edm::ParameterSet nvtx_params10;
0853 nvtx_params10.addParameter<string>("calo", "hgcalHad");
0854 nvtx_params10.addParameter<string>("iEta", "er2p7to3p1");
0855 nvtx_params10.addParameter<vector<double>>("nvtx_params", {-0.006574, 0.003134});
0856 nvtx_params_default.push_back(nvtx_params10);
0857 edm::ParameterSet nvtx_params11;
0858 nvtx_params11.addParameter<string>("calo", "hf");
0859 nvtx_params11.addParameter<string>("iEta", "er29to33");
0860 nvtx_params11.addParameter<vector<double>>("nvtx_params", {-0.203291, 0.044096});
0861 nvtx_params_default.push_back(nvtx_params11);
0862 edm::ParameterSet nvtx_params12;
0863 nvtx_params12.addParameter<string>("calo", "hf");
0864 nvtx_params12.addParameter<string>("iEta", "er34to37");
0865 nvtx_params12.addParameter<vector<double>>("nvtx_params", {-0.210922, 0.045628});
0866 nvtx_params_default.push_back(nvtx_params12);
0867 edm::ParameterSet nvtx_params13;
0868 nvtx_params13.addParameter<string>("calo", "hf");
0869 nvtx_params13.addParameter<string>("iEta", "er38to41");
0870 nvtx_params13.addParameter<vector<double>>("nvtx_params", {-0.229562, 0.050560});
0871 nvtx_params_default.push_back(nvtx_params13);
0872 desc.addVPSet("nvtx_to_PU_sub_params", nvtx_params_validator, nvtx_params_default);
0873
0874 desc.add<vector<double>>("jetPtBins", {0.0, 5.0, 7.5, 10.0, 12.5, 15.0, 17.5, 20.0, 22.5, 25.0, 27.5,
0875 30.0, 35.0, 40.0, 45.0, 50.0, 55.0, 60.0, 65.0, 70.0, 75.0, 80.0,
0876 85.0, 90.0, 95.0, 100.0, 110.0, 120.0, 130.0, 140.0, 150.0, 160.0, 170.0,
0877 180.0, 190.0, 200.0, 225.0, 250.0, 275.0, 300.0, 325.0, 400.0, 500.0});
0878 desc.add<vector<double>>("absEtaBinsBarrel", {0.00, 0.30, 0.60, 1.00, 1.50});
0879 desc.add<vector<double>>(
0880 "jetCalibrationsBarrel",
0881 {2.459, 2.320, 2.239, 2.166, 2.100, 2.040, 1.986, 1.937, 1.892, 1.852, 1.816, 1.768, 1.714, 1.670, 1.633, 1.603,
0882 1.578, 1.557, 1.540, 1.525, 1.513, 1.502, 1.493, 1.486, 1.479, 1.470, 1.460, 1.452, 1.445, 1.439, 1.433, 1.427,
0883 1.422, 1.417, 1.411, 1.403, 1.390, 1.377, 1.365, 1.352, 1.327, 1.284, 4.695, 3.320, 2.751, 2.361, 2.093, 1.908,
0884 1.781, 1.694, 1.633, 1.591, 1.562, 1.533, 1.511, 1.499, 1.492, 1.486, 1.482, 1.478, 1.474, 1.470, 1.467, 1.463,
0885 1.459, 1.456, 1.452, 1.447, 1.440, 1.433, 1.425, 1.418, 1.411, 1.404, 1.397, 1.390, 1.382, 1.370, 1.352, 1.334,
0886 1.316, 1.298, 1.262, 1.200, 5.100, 3.538, 2.892, 2.448, 2.143, 1.933, 1.789, 1.689, 1.620, 1.572, 1.539, 1.506,
0887 1.482, 1.469, 1.460, 1.455, 1.450, 1.446, 1.442, 1.438, 1.434, 1.431, 1.427, 1.423, 1.420, 1.414, 1.407, 1.400,
0888 1.392, 1.385, 1.378, 1.370, 1.363, 1.356, 1.348, 1.336, 1.317, 1.299, 1.281, 1.263, 1.226, 1.162, 3.850, 3.438,
0889 3.211, 3.017, 2.851, 2.708, 2.585, 2.479, 2.388, 2.310, 2.243, 2.159, 2.072, 2.006, 1.956, 1.917, 1.887, 1.863,
0890 1.844, 1.828, 1.814, 1.802, 1.791, 1.782, 1.773, 1.760, 1.744, 1.729, 1.714, 1.699, 1.685, 1.670, 1.656, 1.641,
0891 1.627, 1.602, 1.566, 1.530, 1.494, 1.458, 1.386, 1.260});
0892 desc.add<vector<double>>("absEtaBinsHGCal", {1.50, 1.90, 2.40, 3.00});
0893 desc.add<vector<double>>(
0894 "jetCalibrationsHGCal",
0895 {5.604, 4.578, 4.061, 3.647, 3.314, 3.047, 2.832, 2.660, 2.521, 2.410, 2.320, 2.216, 2.120, 2.056,
0896 2.013, 1.983, 1.961, 1.945, 1.932, 1.922, 1.913, 1.905, 1.898, 1.891, 1.884, 1.874, 1.861, 1.848,
0897 1.835, 1.822, 1.810, 1.797, 1.784, 1.771, 1.759, 1.736, 1.704, 1.673, 1.641, 1.609, 1.545, 1.434,
0898 4.385, 3.584, 3.177, 2.849, 2.584, 2.370, 2.197, 2.057, 1.944, 1.853, 1.780, 1.695, 1.616, 1.564,
0899 1.530, 1.507, 1.491, 1.480, 1.472, 1.466, 1.462, 1.459, 1.456, 1.453, 1.451, 1.447, 1.443, 1.439,
0900 1.435, 1.431, 1.427, 1.423, 1.419, 1.416, 1.412, 1.405, 1.395, 1.385, 1.376, 1.366, 1.346, 1.312,
0901 562.891, 68.647, 17.648, 5.241, 2.223, 1.490, 1.312, 1.270, 1.260, 1.259, 1.259, 1.260, 1.263, 1.265,
0902 1.267, 1.269, 1.271, 1.273, 1.275, 1.277, 1.279, 1.281, 1.283, 1.285, 1.287, 1.290, 1.295, 1.299,
0903 1.303, 1.307, 1.311, 1.315, 1.319, 1.323, 1.328, 1.335, 1.345, 1.355, 1.366, 1.376, 1.397, 1.433});
0904 desc.add<vector<double>>("absEtaBinsHF", {3.00, 3.60, 6.00});
0905 desc.add<vector<double>>(
0906 "jetCalibrationsHF",
0907 {8.169, 6.873, 6.155, 5.535, 5.001, 4.539, 4.141, 3.798, 3.501, 3.245, 3.024, 2.748, 2.463, 2.249,
0908 2.090, 1.971, 1.881, 1.814, 1.763, 1.725, 1.695, 1.673, 1.655, 1.642, 1.631, 1.618, 1.605, 1.596,
0909 1.588, 1.581, 1.575, 1.569, 1.563, 1.557, 1.551, 1.541, 1.527, 1.513, 1.498, 1.484, 1.456, 1.406,
0910 2.788, 2.534, 2.388, 2.258, 2.141, 2.037, 1.945, 1.862, 1.788, 1.722, 1.664, 1.587, 1.503, 1.436,
0911 1.382, 1.339, 1.305, 1.277, 1.255, 1.237, 1.223, 1.211, 1.201, 1.193, 1.186, 1.178, 1.170, 1.164,
0912 1.159, 1.154, 1.151, 1.147, 1.144, 1.141, 1.138, 1.133, 1.126, 1.118, 1.111, 1.104, 1.090, 1.064});
0913 desc.add<vector<double>>("tauPtBins", {0.0, 5.0, 7.5, 10.0, 12.5, 15.0, 20.0, 25.0, 30.0, 35.0,
0914 40.0, 45.0, 50.0, 55.0, 60.0, 70.0, 80.0, 100.0, 150.0, 200.0});
0915 desc.add<vector<double>>("tauAbsEtaBinsBarrel", {0.00, 0.30, 0.60, 1.00, 1.50});
0916 desc.add<vector<double>>("tauCalibrationsBarrel",
0917 {1.067, 1.067, 1.067, 1.067, 1.067, 1.067, 1.067, 1.067, 1.067, 1.067, 1.067, 1.067, 1.067,
0918 1.067, 1.067, 1.067, 1.067, 1.067, 1.067, 1.106, 1.106, 1.106, 1.106, 1.106, 1.106, 1.106,
0919 1.106, 1.106, 1.106, 1.106, 1.106, 1.106, 1.106, 1.106, 1.106, 1.106, 1.106, 1.106, 1.102,
0920 1.102, 1.102, 1.102, 1.102, 1.102, 1.102, 1.102, 1.102, 1.102, 1.102, 1.102, 1.102, 1.102,
0921 1.102, 1.102, 1.102, 1.102, 1.102, 1.139, 1.139, 1.139, 1.139, 1.139, 1.139, 1.139, 1.139});
0922 desc.add<vector<double>>("tauAbsEtaBinsHGCal", {1.50, 1.90, 2.40, 3.00});
0923 desc.add<vector<double>>(
0924 "tauCalibrationsHGCal",
0925 {1.384, 1.384, 1.384, 1.384, 1.384, 1.384, 1.384, 1.384, 1.384, 1.384, 1.384, 1.384, 1.384, 1.384, 1.384,
0926 1.384, 1.384, 1.384, 1.384, 1.473, 1.473, 1.473, 1.473, 1.473, 1.473, 1.473, 1.473, 1.473, 1.473, 1.473,
0927 1.473, 1.473, 1.473, 1.473, 1.473, 1.473, 1.473, 1.473, 1.133, 1.133, 1.133, 1.133, 1.133, 1.133, 1.133,
0928 1.133, 1.133, 1.133, 1.133, 1.133, 1.133, 1.133, 1.133, 1.133, 1.133, 1.133, 1.133});
0929
0930 descriptions.addWithDefaultLabel(desc);
0931 }
0932
0933
0934 DEFINE_FWK_MODULE(Phase2L1CaloJetEmulator);