File indexing completed on 2023-01-16 23:37:03
0001 #include "FWCore/Framework/interface/Event.h"
0002 #include "FWCore/Framework/interface/MakerMacros.h"
0003 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0004
0005 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0006
0007 #include "FWCore/ServiceRegistry/interface/Service.h"
0008 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0009
0010 #include "CondFormats/L1TObjects/interface/CaloParams.h"
0011 #include "CondFormats/DataRecord/interface/L1TCaloParamsRcd.h"
0012
0013 #include "DataFormats/L1TCalorimeter/interface/CaloTower.h"
0014 #include "DataFormats/L1TCalorimeter/interface/CaloCluster.h"
0015 #include "DataFormats/L1Trigger/interface/EGamma.h"
0016 #include "DataFormats/L1Trigger/interface/Tau.h"
0017 #include "DataFormats/L1Trigger/interface/Jet.h"
0018 #include "DataFormats/L1Trigger/interface/EtSum.h"
0019
0020 #include "TH1F.h"
0021 #include "TH2F.h"
0022
0023
0024
0025
0026
0027 namespace l1t {
0028
0029 class L1TStage2CaloAnalyzer : public edm::one::EDAnalyzer<edm::one::SharedResources> {
0030 public:
0031 explicit L1TStage2CaloAnalyzer(const edm::ParameterSet&);
0032 ~L1TStage2CaloAnalyzer() override;
0033
0034 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0035
0036 private:
0037 void beginJob() override;
0038 void analyze(const edm::Event&, const edm::EventSetup&) override;
0039 void endJob() override;
0040
0041
0042
0043
0044
0045
0046
0047 edm::EDGetToken m_towerToken;
0048 edm::EDGetToken m_clusterToken;
0049 edm::EDGetToken m_mpEGToken;
0050 edm::EDGetToken m_mpTauToken;
0051 edm::EDGetToken m_mpJetToken;
0052 edm::EDGetToken m_mpSumToken;
0053 edm::EDGetToken m_egToken;
0054 edm::EDGetToken m_tauToken;
0055 edm::EDGetToken m_jetToken;
0056 edm::EDGetToken m_sumToken;
0057
0058 bool m_doTowers;
0059 bool m_doClusters;
0060 bool m_doMPEGs;
0061 bool m_doMPTaus;
0062 bool m_doMPJets;
0063 bool m_doMPSums;
0064 bool m_doEGs;
0065 bool m_doTaus;
0066 bool m_doJets;
0067 bool m_doSums;
0068
0069 bool doText_;
0070 bool doHistos_;
0071
0072 enum ObjectType {
0073 Tower = 1,
0074 Cluster = 2,
0075 EG = 3,
0076 Tau = 4,
0077 Jet = 5,
0078 SumET = 6,
0079 SumETEm = 7,
0080 SumHT = 8,
0081 SumMET = 9,
0082 SumMETHF = 10,
0083 SumMHT = 11,
0084 SumMHTHF = 12,
0085 MPEG = 13,
0086 MPTau = 14,
0087 MPJet = 15,
0088 MPSum = 16,
0089 MPSumET = 17,
0090 MPSumETHF = 18,
0091 MPSumMETx = 19,
0092 MPSumMETxHF = 20,
0093 MPSumMETy = 21,
0094 MPSumMETyHF = 22,
0095 MPSumHT = 23,
0096 MPSumHTHF = 24,
0097 MPSumMHTx = 25,
0098 MPSumMHTxHF = 26,
0099 MPSumMHTy = 27,
0100 MPSumMHTyHF = 28,
0101 MPMinBiasHFP1 = 29,
0102 MPMinBiasHFM1 = 30,
0103 MPMinBiasHFP0 = 31,
0104 MPMinBiasHFM0 = 32,
0105 MinBiasHFP1 = 33,
0106 MinBiasHFM1 = 34,
0107 MinBiasHFP0 = 35,
0108 MinBiasHFM0 = 36,
0109 MPSumETEm = 37,
0110 MPSumHITowCount = 38,
0111 SumHITowCount = 39,
0112 SumCentrality = 40,
0113 SumAsymEt = 41,
0114 SumAsymEtHF = 42,
0115 SumAsymHt = 43,
0116 SumAsymHtHF = 44
0117 };
0118
0119 std::vector<ObjectType> types_;
0120 std::vector<std::string> typeStr_;
0121
0122 std::map<ObjectType, TFileDirectory> dirs_;
0123 std::map<ObjectType, TH1F*> het_;
0124 std::map<ObjectType, TH1F*> heta_;
0125 std::map<ObjectType, TH1F*> hphi_;
0126 std::map<ObjectType, TH1F*> hbx_;
0127 std::map<ObjectType, TH1F*> hem_;
0128 std::map<ObjectType, TH1F*> hhad_;
0129 std::map<ObjectType, TH1F*> hratio_;
0130 std::map<ObjectType, TH2F*> hetaphi_;
0131 std::map<ObjectType, TH1F*> hiso_;
0132
0133 TFileDirectory evtDispDir_;
0134
0135 TH1F *hsortMP_, *hsort_;
0136
0137 int m_mpBx = 0;
0138 int m_dmxBx = 0;
0139 bool m_allBx = false;
0140 bool m_doEvtDisp = false;
0141 };
0142
0143
0144
0145
0146
0147
0148
0149
0150
0151
0152
0153
0154 L1TStage2CaloAnalyzer::L1TStage2CaloAnalyzer(const edm::ParameterSet& iConfig)
0155 : doText_(iConfig.getUntrackedParameter<bool>("doText", true)),
0156 doHistos_(iConfig.getUntrackedParameter<bool>("doHistos", true)) {
0157 usesResource(TFileService::kSharedResource);
0158
0159
0160 m_mpBx = iConfig.getParameter<int>("mpBx");
0161 m_dmxBx = iConfig.getParameter<int>("dmxBx");
0162 m_allBx = iConfig.getParameter<bool>("allBx");
0163 m_doEvtDisp = iConfig.getParameter<bool>("doEvtDisp");
0164
0165
0166 edm::InputTag nullTag("None");
0167
0168 edm::InputTag towerTag = iConfig.getParameter<edm::InputTag>("towerToken");
0169 m_towerToken = consumes<l1t::CaloTowerBxCollection>(towerTag);
0170 m_doTowers = !(towerTag == nullTag);
0171
0172 edm::InputTag clusterTag = iConfig.getParameter<edm::InputTag>("clusterToken");
0173 m_clusterToken = consumes<l1t::CaloClusterBxCollection>(clusterTag);
0174 m_doClusters = !(clusterTag == nullTag);
0175
0176 edm::InputTag mpEGTag = iConfig.getParameter<edm::InputTag>("mpEGToken");
0177 m_mpEGToken = consumes<l1t::EGammaBxCollection>(mpEGTag);
0178 m_doMPEGs = !(mpEGTag == nullTag);
0179
0180 edm::InputTag mpTauTag = iConfig.getParameter<edm::InputTag>("mpTauToken");
0181 m_mpTauToken = consumes<l1t::TauBxCollection>(mpTauTag);
0182 m_doMPTaus = !(mpTauTag == nullTag);
0183
0184 edm::InputTag mpJetTag = iConfig.getParameter<edm::InputTag>("mpJetToken");
0185 m_mpJetToken = consumes<l1t::JetBxCollection>(mpJetTag);
0186 m_doMPJets = !(mpJetTag == nullTag);
0187
0188 edm::InputTag mpSumTag = iConfig.getParameter<edm::InputTag>("mpEtSumToken");
0189 m_mpSumToken = consumes<l1t::EtSumBxCollection>(mpSumTag);
0190 m_doMPSums = !(mpSumTag == nullTag);
0191
0192 edm::InputTag egTag = iConfig.getParameter<edm::InputTag>("egToken");
0193 m_egToken = consumes<l1t::EGammaBxCollection>(egTag);
0194 m_doEGs = !(egTag == nullTag);
0195
0196 edm::InputTag tauTag = iConfig.getParameter<edm::InputTag>("tauToken");
0197 m_tauToken = consumes<l1t::TauBxCollection>(tauTag);
0198 m_doTaus = !(tauTag == nullTag);
0199
0200 edm::InputTag jetTag = iConfig.getParameter<edm::InputTag>("jetToken");
0201 m_jetToken = consumes<l1t::JetBxCollection>(jetTag);
0202 m_doJets = !(jetTag == nullTag);
0203
0204 edm::InputTag sumTag = iConfig.getParameter<edm::InputTag>("etSumToken");
0205 m_sumToken = consumes<l1t::EtSumBxCollection>(sumTag);
0206 m_doSums = !(sumTag == nullTag);
0207
0208 std::cout << "Processing " << sumTag.label() << std::endl;
0209
0210 types_.push_back(Tower);
0211 types_.push_back(Cluster);
0212 types_.push_back(MPEG);
0213 types_.push_back(MPTau);
0214 types_.push_back(MPJet);
0215 types_.push_back(MPSumET);
0216 types_.push_back(MPSumETHF);
0217 types_.push_back(MPSumETEm);
0218 types_.push_back(MPSumMETx);
0219 types_.push_back(MPSumMETxHF);
0220 types_.push_back(MPSumMETy);
0221 types_.push_back(MPSumMETyHF);
0222 types_.push_back(MPSumHT);
0223 types_.push_back(MPSumHTHF);
0224 types_.push_back(MPSumMHTx);
0225 types_.push_back(MPSumMHTxHF);
0226 types_.push_back(MPSumMHTy);
0227 types_.push_back(MPSumMHTyHF);
0228 types_.push_back(EG);
0229 types_.push_back(Tau);
0230 types_.push_back(Jet);
0231 types_.push_back(SumET);
0232 types_.push_back(SumETEm);
0233 types_.push_back(SumHT);
0234 types_.push_back(SumMET);
0235 types_.push_back(SumMETHF);
0236 types_.push_back(SumMHT);
0237 types_.push_back(SumMHTHF);
0238 types_.push_back(MPMinBiasHFP0);
0239 types_.push_back(MPMinBiasHFM0);
0240 types_.push_back(MPMinBiasHFP1);
0241 types_.push_back(MPMinBiasHFM1);
0242 types_.push_back(MinBiasHFP0);
0243 types_.push_back(MinBiasHFM0);
0244 types_.push_back(MinBiasHFP1);
0245 types_.push_back(MinBiasHFM1);
0246 types_.push_back(MPSumHITowCount);
0247 types_.push_back(SumHITowCount);
0248 types_.push_back(SumCentrality);
0249 types_.push_back(SumAsymEt);
0250 types_.push_back(SumAsymEtHF);
0251 types_.push_back(SumAsymHt);
0252 types_.push_back(SumAsymHtHF);
0253
0254 typeStr_.push_back("tower");
0255 typeStr_.push_back("cluster");
0256 typeStr_.push_back("mpeg");
0257 typeStr_.push_back("mptau");
0258 typeStr_.push_back("mpjet");
0259 typeStr_.push_back("mpsumet");
0260 typeStr_.push_back("mpsumethf");
0261 typeStr_.push_back("mpsumetem");
0262 typeStr_.push_back("mpsummetx");
0263 typeStr_.push_back("mpsummetxhf");
0264 typeStr_.push_back("mpsummety");
0265 typeStr_.push_back("mpsummetyhf");
0266 typeStr_.push_back("mpsumht");
0267 typeStr_.push_back("mpsumhthf");
0268 typeStr_.push_back("mpsummhtx");
0269 typeStr_.push_back("mpsummhtxhf");
0270 typeStr_.push_back("mpsummhty");
0271 typeStr_.push_back("mpsummhtyhf");
0272 typeStr_.push_back("eg");
0273 typeStr_.push_back("tau");
0274 typeStr_.push_back("jet");
0275 typeStr_.push_back("sumet");
0276 typeStr_.push_back("sumetem");
0277 typeStr_.push_back("sumht");
0278 typeStr_.push_back("summet");
0279 typeStr_.push_back("summethf");
0280 typeStr_.push_back("summht");
0281 typeStr_.push_back("summhthf");
0282 typeStr_.push_back("mpminbiashfp0");
0283 typeStr_.push_back("mpminbiashfm0");
0284 typeStr_.push_back("mpminbiashfp1");
0285 typeStr_.push_back("mpminbiashfm1");
0286 typeStr_.push_back("minbiashfp0");
0287 typeStr_.push_back("minbiashfm0");
0288 typeStr_.push_back("minbiashfp1");
0289 typeStr_.push_back("minbiashfm1");
0290 typeStr_.push_back("mpsumhitowercount");
0291 typeStr_.push_back("sumhitowercount");
0292 typeStr_.push_back("sumcentrality");
0293 typeStr_.push_back("sumasymet");
0294 typeStr_.push_back("sumasymethf");
0295 typeStr_.push_back("sumasymht");
0296 typeStr_.push_back("sumasymhthf");
0297 }
0298
0299 L1TStage2CaloAnalyzer::~L1TStage2CaloAnalyzer() {
0300
0301
0302 }
0303
0304
0305
0306
0307
0308
0309 void L1TStage2CaloAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0310 using namespace edm;
0311
0312 std::stringstream text;
0313
0314 TH2F* hEvtTow = new TH2F();
0315 TH2F* hEvtMPEG = new TH2F();
0316 TH2F* hEvtMPTau = new TH2F();
0317 TH2F* hEvtMPJet = new TH2F();
0318 TH2F* hEvtDemuxEG = new TH2F();
0319 TH2F* hEvtDemuxTau = new TH2F();
0320 TH2F* hEvtDemuxJet = new TH2F();
0321
0322 if (m_doEvtDisp) {
0323 std::stringstream ss;
0324 ss << iEvent.run() << "-" << iEvent.id().event();
0325 TFileDirectory dir = evtDispDir_.mkdir(ss.str());
0326 hEvtTow = dir.make<TH2F>("Tower", "", 83, -41.5, 41.5, 72, .5, 72.5);
0327 hEvtMPEG = dir.make<TH2F>("MPEG", "", 83, -41.5, 41.5, 72, .5, 72.5);
0328 hEvtMPTau = dir.make<TH2F>("MPTau", "", 83, -41.5, 41.5, 72, .5, 72.5);
0329 hEvtMPJet = dir.make<TH2F>("MPJet", "", 83, -41.5, 41.5, 72, .5, 72.5);
0330 hEvtDemuxEG = dir.make<TH2F>("DemuxEG", "", 227, -113.5, 113.5, 144, -0.5, 143.5);
0331 hEvtDemuxTau = dir.make<TH2F>("DemuxTau", "", 227, -113.5, 113.5, 144, -0.5, 143.5);
0332 hEvtDemuxJet = dir.make<TH2F>("DemuxJet", "", 227, -113.5, 113.5, 144, -0.5, 143.5);
0333 }
0334
0335
0336
0337
0338
0339
0340 if (m_mpBx < -2 || m_mpBx > 2 || m_dmxBx < -2 || m_dmxBx > 2)
0341 edm::LogError("L1T")
0342 << "Selected MP Bx or Demux Bx to fill histograms is outside of range -2,2. Histos will be empty!";
0343
0344
0345 if (m_doTowers) {
0346 Handle<BXVector<l1t::CaloTower> > towers;
0347 iEvent.getByToken(m_towerToken, towers);
0348
0349 for (int ibx = towers->getFirstBX(); ibx <= towers->getLastBX(); ++ibx) {
0350 if (!m_allBx && ibx != m_mpBx)
0351 continue;
0352
0353 for (auto itr = towers->begin(ibx); itr != towers->end(ibx); ++itr) {
0354 if (itr->hwPt() <= 0)
0355 continue;
0356
0357 hbx_.at(Tower)->Fill(ibx);
0358 het_.at(Tower)->Fill(itr->hwPt());
0359 heta_.at(Tower)->Fill(itr->hwEta());
0360 hphi_.at(Tower)->Fill(itr->hwPhi());
0361 hem_.at(Tower)->Fill(itr->hwEtEm());
0362 hhad_.at(Tower)->Fill(itr->hwEtHad());
0363 hratio_.at(Tower)->Fill(itr->hwEtRatio());
0364 hetaphi_.at(Tower)->Fill(itr->hwEta(), itr->hwPhi(), itr->hwPt());
0365
0366 text << "Tower : "
0367 << " BX=" << ibx << " ipt=" << itr->hwPt() << " ieta=" << itr->hwEta() << " iphi=" << itr->hwPhi()
0368 << " iem=" << itr->hwEtEm() << " ihad=" << itr->hwEtHad() << " iratio=" << itr->hwEtRatio() << std::endl;
0369
0370 if (m_doEvtDisp)
0371 hEvtTow->Fill(itr->hwEta(), itr->hwPhi(), itr->hwPt());
0372 }
0373 }
0374 }
0375
0376 if (m_doClusters) {
0377 Handle<BXVector<l1t::CaloCluster> > clusters;
0378 iEvent.getByToken(m_clusterToken, clusters);
0379
0380 for (int ibx = clusters->getFirstBX(); ibx <= clusters->getLastBX(); ++ibx) {
0381 if (!m_allBx && ibx != m_mpBx)
0382 continue;
0383
0384 for (auto itr = clusters->begin(ibx); itr != clusters->end(ibx); ++itr) {
0385 hbx_.at(Cluster)->Fill(ibx);
0386 het_.at(Cluster)->Fill(itr->hwPt());
0387 heta_.at(Cluster)->Fill(itr->hwEta());
0388 hphi_.at(Cluster)->Fill(itr->hwPhi());
0389 hetaphi_.at(Cluster)->Fill(itr->hwEta(), itr->hwPhi(), itr->hwPt());
0390 text << "Cluster : "
0391 << " BX=" << ibx << " ipt=" << itr->hwPt() << " ieta=" << itr->hwEta() << " iphi=" << itr->hwPhi()
0392 << std::endl;
0393 }
0394 }
0395 }
0396
0397
0398 if (m_doMPEGs) {
0399 Handle<BXVector<l1t::EGamma> > mpegs;
0400 iEvent.getByToken(m_mpEGToken, mpegs);
0401
0402 for (int ibx = mpegs->getFirstBX(); ibx <= mpegs->getLastBX(); ++ibx) {
0403 if (!m_allBx && ibx != m_mpBx)
0404 continue;
0405
0406 for (auto itr = mpegs->begin(ibx); itr != mpegs->end(ibx); ++itr) {
0407 hbx_.at(MPEG)->Fill(ibx);
0408 het_.at(MPEG)->Fill(itr->hwPt());
0409 heta_.at(MPEG)->Fill(itr->hwEta());
0410 hphi_.at(MPEG)->Fill(itr->hwPhi());
0411 hetaphi_.at(MPEG)->Fill(itr->hwEta(), itr->hwPhi(), itr->hwPt());
0412 hiso_.at(MPEG)->Fill(itr->hwIso());
0413
0414 text << "MP EG : "
0415 << " BX=" << ibx << " ipt=" << itr->hwPt() << " ieta=" << itr->hwEta() << " iphi=" << itr->hwPhi()
0416 << std::endl;
0417
0418 if (m_doEvtDisp)
0419 hEvtMPEG->Fill(itr->hwEta(), itr->hwPhi(), itr->hwPt());
0420 }
0421 }
0422 }
0423
0424
0425 if (m_doMPTaus) {
0426 Handle<BXVector<l1t::Tau> > mptaus;
0427 iEvent.getByToken(m_mpTauToken, mptaus);
0428
0429 for (int ibx = mptaus->getFirstBX(); ibx <= mptaus->getLastBX(); ++ibx) {
0430 if (!m_allBx && ibx != m_mpBx)
0431 continue;
0432
0433 for (auto itr = mptaus->begin(ibx); itr != mptaus->end(ibx); ++itr) {
0434 hbx_.at(MPTau)->Fill(ibx);
0435 het_.at(MPTau)->Fill(itr->hwPt());
0436 heta_.at(MPTau)->Fill(itr->hwEta());
0437 hphi_.at(MPTau)->Fill(itr->hwPhi());
0438 hetaphi_.at(MPTau)->Fill(itr->hwEta(), itr->hwPhi(), itr->hwPt());
0439
0440 text << "MP Tau : "
0441 << " BX=" << ibx << " ipt=" << itr->hwPt() << " ieta=" << itr->hwEta() << " iphi=" << itr->hwPhi()
0442 << std::endl;
0443
0444 if (m_doEvtDisp)
0445 hEvtMPTau->Fill(itr->hwEta(), itr->hwPhi(), itr->hwPt());
0446 }
0447 }
0448 }
0449
0450
0451 std::vector<l1t::Jet> thejets_poseta;
0452 std::vector<l1t::Jet> thejets_negeta;
0453
0454 if (m_doMPJets) {
0455 Handle<BXVector<l1t::Jet> > mpjets;
0456 iEvent.getByToken(m_mpJetToken, mpjets);
0457
0458
0459
0460
0461 for (int ibx = mpjets->getFirstBX(); ibx <= mpjets->getLastBX(); ++ibx) {
0462 if (!m_allBx && ibx != m_mpBx)
0463 continue;
0464
0465 for (auto itr = mpjets->begin(ibx); itr != mpjets->end(ibx); ++itr) {
0466 hbx_.at(MPJet)->Fill(ibx);
0467 het_.at(MPJet)->Fill(itr->hwPt());
0468 heta_.at(MPJet)->Fill(itr->hwEta());
0469 hphi_.at(MPJet)->Fill(itr->hwPhi());
0470 hetaphi_.at(MPJet)->Fill(itr->hwEta(), itr->hwPhi(), itr->hwPt());
0471
0472 text << "MP Jet : "
0473 << " BX=" << ibx << " ipt=" << itr->hwPt() << " ieta=" << itr->hwEta() << " iphi=" << itr->hwPhi()
0474 << std::endl;
0475
0476 if (m_doEvtDisp)
0477 hEvtMPJet->Fill(itr->hwEta(), itr->hwPhi(), itr->hwPt());
0478
0479 itr->hwEta() > 0 ? thejets_poseta.push_back(*itr) : thejets_negeta.push_back(*itr);
0480 }
0481 }
0482 }
0483
0484 if (!thejets_poseta.empty()) {
0485 for (unsigned int i = 0; i < thejets_poseta.size() - 1; i++) {
0486 for (unsigned int j = i + 1; j < thejets_poseta.size(); j++) {
0487 hsortMP_->Fill(thejets_poseta.at(i).hwPt() - thejets_poseta.at(j).hwPt());
0488 }
0489 }
0490 }
0491
0492 if (!thejets_negeta.empty()) {
0493 for (unsigned int i = 0; i < thejets_negeta.size() - 1; i++) {
0494 for (unsigned int j = i + 1; j < thejets_negeta.size(); j++) {
0495 hsortMP_->Fill(thejets_negeta.at(i).hwPt() - thejets_negeta.at(j).hwPt());
0496 }
0497 }
0498 }
0499
0500
0501 if (m_doMPSums) {
0502 Handle<BXVector<l1t::EtSum> > mpsums;
0503 iEvent.getByToken(m_mpSumToken, mpsums);
0504
0505 for (int ibx = mpsums->getFirstBX(); ibx <= mpsums->getLastBX(); ++ibx) {
0506 if (!m_allBx && ibx != m_mpBx)
0507 continue;
0508
0509 for (auto itr = mpsums->begin(ibx); itr != mpsums->end(ibx); ++itr) {
0510 switch (itr->getType()) {
0511 case l1t::EtSum::EtSumType::kTotalEt:
0512 het_.at(MPSumET)->Fill(itr->hwPt());
0513 break;
0514 case l1t::EtSum::EtSumType::kTotalEtHF:
0515 het_.at(MPSumETHF)->Fill(itr->hwPt());
0516 break;
0517 case l1t::EtSum::EtSumType::kTotalEtEm:
0518 het_.at(MPSumETEm)->Fill(itr->hwPt());
0519 break;
0520 case l1t::EtSum::EtSumType::kTotalEtx:
0521 het_.at(MPSumMETx)->Fill(itr->hwPt());
0522 break;
0523 case l1t::EtSum::EtSumType::kTotalEtxHF:
0524 het_.at(MPSumMETxHF)->Fill(itr->hwPt());
0525 break;
0526 case l1t::EtSum::EtSumType::kTotalEty:
0527 het_.at(MPSumMETy)->Fill(itr->hwPt());
0528 break;
0529 case l1t::EtSum::EtSumType::kTotalEtyHF:
0530 het_.at(MPSumMETyHF)->Fill(itr->hwPt());
0531 break;
0532 case l1t::EtSum::EtSumType::kTotalHt:
0533 het_.at(MPSumHT)->Fill(itr->hwPt());
0534 break;
0535 case l1t::EtSum::EtSumType::kTotalHtHF:
0536 het_.at(MPSumHTHF)->Fill(itr->hwPt());
0537 break;
0538 case l1t::EtSum::EtSumType::kTotalHtx:
0539 het_.at(MPSumMHTx)->Fill(itr->hwPt());
0540 break;
0541 case l1t::EtSum::EtSumType::kTotalHtxHF:
0542 het_.at(MPSumMHTxHF)->Fill(itr->hwPt());
0543 break;
0544 case l1t::EtSum::EtSumType::kTotalHty:
0545 het_.at(MPSumMHTy)->Fill(itr->hwPt());
0546 break;
0547 case l1t::EtSum::EtSumType::kTotalHtyHF:
0548 het_.at(MPSumMHTyHF)->Fill(itr->hwPt());
0549 break;
0550 case l1t::EtSum::EtSumType::kMinBiasHFP0:
0551 het_.at(MPMinBiasHFP0)->Fill(itr->hwPt());
0552 break;
0553 case l1t::EtSum::EtSumType::kMinBiasHFM0:
0554 het_.at(MPMinBiasHFM0)->Fill(itr->hwPt());
0555 break;
0556 case l1t::EtSum::EtSumType::kMinBiasHFP1:
0557 het_.at(MPMinBiasHFP1)->Fill(itr->hwPt());
0558 break;
0559 case l1t::EtSum::EtSumType::kMinBiasHFM1:
0560 het_.at(MPMinBiasHFM1)->Fill(itr->hwPt());
0561 break;
0562 case l1t::EtSum::EtSumType::kTowerCount:
0563 het_.at(MPSumHITowCount)->Fill(itr->hwPt());
0564 break;
0565
0566 default:
0567 std::cout << "wrong type of MP sum" << std::endl;
0568 }
0569 text << "MP Sum : "
0570 << " type=" << itr->getType() << " BX=" << ibx << " ipt=" << itr->hwPt() << " ieta=" << itr->hwEta()
0571 << " iphi=" << itr->hwPhi() << std::endl;
0572 }
0573 }
0574 }
0575
0576
0577 if (m_doEGs) {
0578 Handle<BXVector<l1t::EGamma> > egs;
0579 iEvent.getByToken(m_egToken, egs);
0580
0581 for (int ibx = egs->getFirstBX(); ibx <= egs->getLastBX(); ++ibx) {
0582 if (!m_allBx && ibx != m_dmxBx)
0583 continue;
0584
0585 for (auto itr = egs->begin(ibx); itr != egs->end(ibx); ++itr) {
0586 hbx_.at(EG)->Fill(ibx);
0587 het_.at(EG)->Fill(itr->hwPt());
0588 heta_.at(EG)->Fill(itr->hwEta());
0589 hphi_.at(EG)->Fill(itr->hwPhi());
0590 hetaphi_.at(EG)->Fill(itr->hwEta(), itr->hwPhi(), itr->hwPt());
0591 hiso_.at(EG)->Fill(itr->hwIso());
0592
0593 text << "EG : "
0594 << " BX=" << ibx << " ipt=" << itr->hwPt() << " ieta=" << itr->hwEta() << " iphi=" << itr->hwPhi()
0595 << std::endl;
0596
0597 if (m_doEvtDisp)
0598 hEvtDemuxEG->Fill(itr->hwEta(), itr->hwPhi(), itr->hwPt());
0599 }
0600 }
0601 }
0602
0603
0604 if (m_doTaus) {
0605 Handle<BXVector<l1t::Tau> > taus;
0606 iEvent.getByToken(m_tauToken, taus);
0607
0608 for (int ibx = taus->getFirstBX(); ibx <= taus->getLastBX(); ++ibx) {
0609 if (!m_allBx && ibx != m_dmxBx)
0610 continue;
0611
0612 for (auto itr = taus->begin(ibx); itr != taus->end(ibx); ++itr) {
0613 hbx_.at(Tau)->Fill(ibx);
0614 het_.at(Tau)->Fill(itr->hwPt());
0615 heta_.at(Tau)->Fill(itr->hwEta());
0616 hphi_.at(Tau)->Fill(itr->hwPhi());
0617 hetaphi_.at(Tau)->Fill(itr->hwEta(), itr->hwPhi(), itr->hwPt());
0618
0619 text << "Tau : "
0620 << " BX=" << ibx << " ipt=" << itr->hwPt() << " ieta=" << itr->hwEta() << " iphi=" << itr->hwPhi()
0621 << std::endl;
0622
0623 if (m_doEvtDisp)
0624 hEvtDemuxTau->Fill(itr->hwEta(), itr->hwPhi(), itr->hwPt());
0625 }
0626 }
0627 }
0628
0629
0630 std::vector<l1t::Jet> thejets;
0631
0632 if (m_doJets) {
0633 Handle<BXVector<l1t::Jet> > jets;
0634 iEvent.getByToken(m_jetToken, jets);
0635
0636 for (int ibx = jets->getFirstBX(); ibx <= jets->getLastBX(); ++ibx) {
0637 if (!m_allBx && ibx != m_dmxBx)
0638 continue;
0639
0640 for (auto itr = jets->begin(ibx); itr != jets->end(ibx); ++itr) {
0641 hbx_.at(Jet)->Fill(ibx);
0642 het_.at(Jet)->Fill(itr->hwPt());
0643 heta_.at(Jet)->Fill(itr->hwEta());
0644 hphi_.at(Jet)->Fill(itr->hwPhi());
0645 hetaphi_.at(Jet)->Fill(itr->hwEta(), itr->hwPhi(), itr->hwPt());
0646
0647 text << "Jet : "
0648 << " BX=" << ibx << " ipt=" << itr->hwPt() << " ieta=" << itr->hwEta() << " iphi=" << itr->hwPhi()
0649 << std::endl;
0650
0651 if (m_doEvtDisp)
0652 hEvtDemuxJet->Fill(itr->hwEta(), itr->hwPhi(), itr->hwPt());
0653
0654 thejets.push_back(*itr);
0655 }
0656 }
0657 }
0658
0659 if (!thejets.empty()) {
0660 for (unsigned int i = 0; i < thejets.size() - 1; i++) {
0661 for (unsigned int j = i + 1; j < thejets.size(); j++) {
0662 hsort_->Fill(thejets.at(i).hwPt() - thejets.at(j).hwPt());
0663 }
0664 }
0665 }
0666
0667
0668 if (m_doSums) {
0669 Handle<BXVector<l1t::EtSum> > sums;
0670 iEvent.getByToken(m_sumToken, sums);
0671
0672 for (int ibx = sums->getFirstBX(); ibx <= sums->getLastBX(); ++ibx) {
0673 if (!m_allBx && ibx != m_dmxBx)
0674 continue;
0675
0676 for (auto itr = sums->begin(ibx); itr != sums->end(ibx); ++itr) {
0677 switch (itr->getType()) {
0678 case l1t::EtSum::EtSumType::kTotalEt:
0679 het_.at(SumET)->Fill(itr->hwPt());
0680 break;
0681 case l1t::EtSum::EtSumType::kTotalEtHF:
0682 break;
0683 case l1t::EtSum::EtSumType::kTotalEtEm:
0684 het_.at(SumETEm)->Fill(itr->hwPt());
0685 break;
0686 case l1t::EtSum::EtSumType::kTotalHt:
0687 het_.at(SumHT)->Fill(itr->hwPt());
0688 break;
0689 case l1t::EtSum::EtSumType::kTotalHtHF:
0690 break;
0691 case l1t::EtSum::EtSumType::kMissingEt:
0692 het_.at(SumMET)->Fill(itr->hwPt());
0693 hphi_.at(SumMET)->Fill(itr->hwPhi());
0694 break;
0695 case l1t::EtSum::EtSumType::kMissingEtHF:
0696 het_.at(SumMETHF)->Fill(itr->hwPt());
0697 hphi_.at(SumMETHF)->Fill(itr->hwPhi());
0698 break;
0699 case l1t::EtSum::EtSumType::kMissingHt:
0700 het_.at(SumMHT)->Fill(itr->hwPt());
0701 hphi_.at(SumMHT)->Fill(itr->hwPhi());
0702 break;
0703 case l1t::EtSum::EtSumType::kMissingHtHF:
0704 het_.at(SumMHTHF)->Fill(itr->hwPt());
0705 hphi_.at(SumMHTHF)->Fill(itr->hwPhi());
0706 break;
0707 case l1t::EtSum::EtSumType::kMinBiasHFP0:
0708 het_.at(MinBiasHFP0)->Fill(itr->hwPt());
0709 break;
0710 case l1t::EtSum::EtSumType::kMinBiasHFM0:
0711 het_.at(MinBiasHFM0)->Fill(itr->hwPt());
0712 break;
0713 case l1t::EtSum::EtSumType::kMinBiasHFP1:
0714 het_.at(MinBiasHFP1)->Fill(itr->hwPt());
0715 break;
0716 case l1t::EtSum::EtSumType::kMinBiasHFM1:
0717 het_.at(MinBiasHFM1)->Fill(itr->hwPt());
0718 break;
0719 case l1t::EtSum::EtSumType::kTowerCount:
0720 het_.at(SumHITowCount)->Fill(itr->hwPt());
0721 break;
0722 case l1t::EtSum::EtSumType::kCentrality:
0723 het_.at(SumCentrality)->Fill(itr->hwPt());
0724 break;
0725 case l1t::EtSum::EtSumType::kAsymEt:
0726 het_.at(SumAsymEt)->Fill(itr->hwPt());
0727 break;
0728 case l1t::EtSum::EtSumType::kAsymHt:
0729 het_.at(SumAsymHt)->Fill(itr->hwPt());
0730 break;
0731 case l1t::EtSum::EtSumType::kAsymEtHF:
0732 het_.at(SumAsymEtHF)->Fill(itr->hwPt());
0733 break;
0734 case l1t::EtSum::EtSumType::kAsymHtHF:
0735 het_.at(SumAsymHtHF)->Fill(itr->hwPt());
0736 break;
0737
0738 default:
0739 std::cout << "wrong type of demux sum: " << itr->getType() << std::endl;
0740 }
0741
0742 text << "Sum : "
0743 << " type=" << itr->getType() << " BX=" << ibx << " ipt=" << itr->hwPt() << " ieta=" << itr->hwEta()
0744 << " iphi=" << itr->hwPhi() << std::endl;
0745 }
0746 }
0747 }
0748
0749 if (doText_)
0750 edm::LogVerbatim("L1TCaloEvents") << text.str();
0751 }
0752
0753
0754 void L1TStage2CaloAnalyzer::beginJob() {
0755 edm::Service<TFileService> fs;
0756
0757 auto itr = types_.cbegin();
0758 auto str = typeStr_.cbegin();
0759
0760 for (; itr != types_.end(); ++itr, ++str) {
0761 dirs_.insert(std::pair<ObjectType, TFileDirectory>(*itr, fs->mkdir(*str)));
0762
0763 if (*itr == MPSumMETx || *itr == MPSumMETxHF || *itr == MPSumMETy || *itr == MPSumMETyHF) {
0764 het_.insert(
0765 std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("et", "", 2200000, -2199999500, 2200000500)));
0766 } else if (*itr == MPSumMHTx || *itr == MPSumMHTxHF || *itr == MPSumMHTy || *itr == MPSumMHTyHF) {
0767 het_.insert(
0768 std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("et", "", 2200000, -2199999500, 2200000500)));
0769 } else if (*itr == SumET || *itr == SumETEm || *itr == MPSumETEm || *itr == MPSumET || *itr == MPSumETHF ||
0770 *itr == SumHT || *itr == MPSumHT || *itr == MPSumHTHF) {
0771 het_.insert(std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("et", "", 100000, -0.5, 99999.5)));
0772 } else if (*itr == MPMinBiasHFP0 || *itr == MPMinBiasHFM0 || *itr == MPMinBiasHFP1 || *itr == MPMinBiasHFM1 ||
0773 *itr == MinBiasHFP0 || *itr == MinBiasHFM0 || *itr == MinBiasHFP1 || *itr == MinBiasHFM1) {
0774 het_.insert(std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("et", "", 16, -0.5, 15.5)));
0775 } else if (*itr == MPSumHITowCount || *itr == SumHITowCount) {
0776 het_.insert(std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("et", "", 5904, -0.5, 5903.5)));
0777 } else if (*itr == EG || *itr == Tau || *itr == MPEG || *itr == MPTau) {
0778 het_.insert(std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("et", "", 1500, -0.5, 1499.5)));
0779 } else {
0780 het_.insert(std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("et", "", 100000, -0.5, 99999.5)));
0781 }
0782
0783 hbx_.insert(std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("bx", "", 11, -5.5, 5.5)));
0784
0785 if (*itr == EG || *itr == Jet || *itr == Tau) {
0786 heta_.insert(std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("eta", "", 227, -113.5, 113.5)));
0787 hphi_.insert(std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("phi", "", 144, -0.5, 143.5)));
0788 hetaphi_.insert(std::pair<ObjectType, TH2F*>(
0789 *itr, dirs_.at(*itr).make<TH2F>("etaphi", "", 227, -113.5, 113.5, 144, -0.5, 143.5)));
0790 if (*itr == EG)
0791 hiso_.insert(std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("iso", "", 4, -0.5, 3.5)));
0792 } else if (*itr == Tower || *itr == Cluster || *itr == MPEG || *itr == MPJet || *itr == MPTau) {
0793 heta_.insert(std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("eta", "", 83, -41.5, 41.5)));
0794 hphi_.insert(std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("phi", "", 72, 0.5, 72.5)));
0795 hetaphi_.insert(
0796 std::pair<ObjectType, TH2F*>(*itr, dirs_.at(*itr).make<TH2F>("etaphi", "", 83, -41.5, 41.5, 72, .5, 72.5)));
0797 if (*itr == MPEG)
0798 hiso_.insert(std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("iso", "", 4, -0.5, 3.5)));
0799 } else if (*itr == SumMET || *itr == SumMETHF || *itr == SumMHT || *itr == SumMHTHF) {
0800 hphi_.insert(std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("phi", "", 1008, -0.5, 1007.5)));
0801 }
0802
0803 if (*itr == Tower) {
0804 hem_.insert(std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("em", "", 101, -0.5, 100.5)));
0805 hhad_.insert(std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("had", "", 101, -0.5, 100.5)));
0806 hratio_.insert(std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("ratio", "", 11, -0.5, 10.5)));
0807 }
0808 }
0809
0810 if (m_doEvtDisp) {
0811 evtDispDir_ = fs->mkdir("Events");
0812 }
0813
0814 hsort_ = fs->make<TH1F>("sort", "", 201, -100.5, 100.5);
0815 hsortMP_ = fs->make<TH1F>("sortMP", "", 201, -100.5, 100.5);
0816 }
0817
0818
0819 void L1TStage2CaloAnalyzer::endJob() {}
0820
0821
0822
0823
0824
0825
0826
0827
0828
0829
0830
0831
0832
0833
0834
0835
0836
0837
0838
0839
0840
0841
0842
0843
0844
0845
0846
0847
0848
0849
0850
0851
0852
0853
0854 void L1TStage2CaloAnalyzer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0855
0856
0857 edm::ParameterSetDescription desc;
0858 desc.setUnknown();
0859 descriptions.addDefault(desc);
0860 }
0861
0862 }
0863
0864 using namespace l1t;
0865
0866
0867 DEFINE_FWK_MODULE(L1TStage2CaloAnalyzer);