Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-06-29 02:25:58

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 // class declaration
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     //virtual void beginRun(edm::Run const&, edm::EventSetup const&) override;
0042     //virtual void endRun(edm::Run const&, edm::EventSetup const&) override;
0043     //virtual void beginLuminosityBlock(edm::LuminosityBlock const&, edm::EventSetup const&) override;
0044     //virtual void endLuminosityBlock(edm::LuminosityBlock const&, edm::EventSetup const&) override;
0045 
0046     // ----------member data ---------------------------
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   // constants, enums and typedefs
0145   //
0146 
0147   //
0148   // static data member definitions
0149   //
0150 
0151   //
0152   // constructors and destructor
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     //now do what ever initialization is needed
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     // register what you consume and keep token for later access:
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     // do anything here that needs to be done at desctruction time
0301     // (e.g. close files, deallocate resources etc.)
0302   }
0303 
0304   //
0305   // member functions
0306   //
0307 
0308   // ------------ method called for each event  ------------
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     // get TPs ?
0336     // get regions ?
0337     // get RCT clusters ?
0338 
0339     //check mpbx and dmxbx
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     // get towers
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     // get cluster
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     // get EG
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     // get tau
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     // get jet
0451     int njetmp = 0;
0452     std::vector<l1t::Jet> thejets_poseta;
0453     std::vector<l1t::Jet> thejets_negeta;
0454 
0455     if (m_doMPJets) {
0456       Handle<BXVector<l1t::Jet> > mpjets;
0457       iEvent.getByToken(m_mpJetToken, mpjets);
0458 
0459       //Handle<BXVector<l1t::Jet>> jets;
0460       //iEvent.getByToken(m_jetToken,jets);
0461       //if (mpjets->size(0) == jets->size(0)) std::cout<<"******notequal"<<std::endl;
0462       for (int ibx = mpjets->getFirstBX(); ibx <= mpjets->getLastBX(); ++ibx) {
0463         if (!m_allBx && ibx != m_mpBx)
0464           continue;
0465 
0466         for (auto itr = mpjets->begin(ibx); itr != mpjets->end(ibx); ++itr) {
0467           njetmp += 1;
0468           hbx_.at(MPJet)->Fill(ibx);
0469           het_.at(MPJet)->Fill(itr->hwPt());
0470           heta_.at(MPJet)->Fill(itr->hwEta());
0471           hphi_.at(MPJet)->Fill(itr->hwPhi());
0472           hetaphi_.at(MPJet)->Fill(itr->hwEta(), itr->hwPhi(), itr->hwPt());
0473 
0474           text << "MP Jet : "
0475                << " BX=" << ibx << " ipt=" << itr->hwPt() << " ieta=" << itr->hwEta() << " iphi=" << itr->hwPhi()
0476                << std::endl;
0477 
0478           if (m_doEvtDisp)
0479             hEvtMPJet->Fill(itr->hwEta(), itr->hwPhi(), itr->hwPt());
0480 
0481           itr->hwEta() > 0 ? thejets_poseta.push_back(*itr) : thejets_negeta.push_back(*itr);
0482         }
0483       }
0484     }
0485 
0486     if (!thejets_poseta.empty()) {
0487       for (unsigned int i = 0; i < thejets_poseta.size() - 1; i++) {
0488         for (unsigned int j = i + 1; j < thejets_poseta.size(); j++) {
0489           hsortMP_->Fill(thejets_poseta.at(i).hwPt() - thejets_poseta.at(j).hwPt());
0490         }
0491       }
0492     }
0493 
0494     if (!thejets_negeta.empty()) {
0495       for (unsigned int i = 0; i < thejets_negeta.size() - 1; i++) {
0496         for (unsigned int j = i + 1; j < thejets_negeta.size(); j++) {
0497           hsortMP_->Fill(thejets_negeta.at(i).hwPt() - thejets_negeta.at(j).hwPt());
0498         }
0499       }
0500     }
0501 
0502     //std::cout<<"njetmp "<<njetmp<<std::endl;
0503 
0504     // get sums
0505     if (m_doMPSums) {
0506       Handle<BXVector<l1t::EtSum> > mpsums;
0507       iEvent.getByToken(m_mpSumToken, mpsums);
0508 
0509       for (int ibx = mpsums->getFirstBX(); ibx <= mpsums->getLastBX(); ++ibx) {
0510         if (!m_allBx && ibx != m_mpBx)
0511           continue;
0512 
0513         for (auto itr = mpsums->begin(ibx); itr != mpsums->end(ibx); ++itr) {
0514           switch (itr->getType()) {
0515             case l1t::EtSum::EtSumType::kTotalEt:
0516               het_.at(MPSumET)->Fill(itr->hwPt());
0517               break;
0518             case l1t::EtSum::EtSumType::kTotalEtHF:
0519               het_.at(MPSumETHF)->Fill(itr->hwPt());
0520               break;
0521             case l1t::EtSum::EtSumType::kTotalEtEm:
0522               het_.at(MPSumETEm)->Fill(itr->hwPt());
0523               break;
0524             case l1t::EtSum::EtSumType::kTotalEtx:
0525               het_.at(MPSumMETx)->Fill(itr->hwPt());
0526               break;
0527             case l1t::EtSum::EtSumType::kTotalEtxHF:
0528               het_.at(MPSumMETxHF)->Fill(itr->hwPt());
0529               break;
0530             case l1t::EtSum::EtSumType::kTotalEty:
0531               het_.at(MPSumMETy)->Fill(itr->hwPt());
0532               break;
0533             case l1t::EtSum::EtSumType::kTotalEtyHF:
0534               het_.at(MPSumMETyHF)->Fill(itr->hwPt());
0535               break;
0536             case l1t::EtSum::EtSumType::kTotalHt:
0537               het_.at(MPSumHT)->Fill(itr->hwPt());
0538               break;
0539             case l1t::EtSum::EtSumType::kTotalHtHF:
0540               het_.at(MPSumHTHF)->Fill(itr->hwPt());
0541               break;
0542             case l1t::EtSum::EtSumType::kTotalHtx:
0543               het_.at(MPSumMHTx)->Fill(itr->hwPt());
0544               break;
0545             case l1t::EtSum::EtSumType::kTotalHtxHF:
0546               het_.at(MPSumMHTxHF)->Fill(itr->hwPt());
0547               break;
0548             case l1t::EtSum::EtSumType::kTotalHty:
0549               het_.at(MPSumMHTy)->Fill(itr->hwPt());
0550               break;
0551             case l1t::EtSum::EtSumType::kTotalHtyHF:
0552               het_.at(MPSumMHTyHF)->Fill(itr->hwPt());
0553               break;
0554             case l1t::EtSum::EtSumType::kMinBiasHFP0:
0555               het_.at(MPMinBiasHFP0)->Fill(itr->hwPt());
0556               break;
0557             case l1t::EtSum::EtSumType::kMinBiasHFM0:
0558               het_.at(MPMinBiasHFM0)->Fill(itr->hwPt());
0559               break;
0560             case l1t::EtSum::EtSumType::kMinBiasHFP1:
0561               het_.at(MPMinBiasHFP1)->Fill(itr->hwPt());
0562               break;
0563             case l1t::EtSum::EtSumType::kMinBiasHFM1:
0564               het_.at(MPMinBiasHFM1)->Fill(itr->hwPt());
0565               break;
0566             case l1t::EtSum::EtSumType::kTowerCount:
0567               het_.at(MPSumHITowCount)->Fill(itr->hwPt());
0568               break;
0569 
0570             default:
0571               std::cout << "wrong type of MP sum" << std::endl;
0572           }
0573           text << "MP Sum : "
0574                << " type=" << itr->getType() << " BX=" << ibx << " ipt=" << itr->hwPt() << " ieta=" << itr->hwEta()
0575                << " iphi=" << itr->hwPhi() << std::endl;
0576         }
0577       }
0578     }
0579 
0580     // get EG
0581     if (m_doEGs) {
0582       Handle<BXVector<l1t::EGamma> > egs;
0583       iEvent.getByToken(m_egToken, egs);
0584 
0585       for (int ibx = egs->getFirstBX(); ibx <= egs->getLastBX(); ++ibx) {
0586         if (!m_allBx && ibx != m_dmxBx)
0587           continue;
0588 
0589         for (auto itr = egs->begin(ibx); itr != egs->end(ibx); ++itr) {
0590           hbx_.at(EG)->Fill(ibx);
0591           het_.at(EG)->Fill(itr->hwPt());
0592           heta_.at(EG)->Fill(itr->hwEta());
0593           hphi_.at(EG)->Fill(itr->hwPhi());
0594           hetaphi_.at(EG)->Fill(itr->hwEta(), itr->hwPhi(), itr->hwPt());
0595           hiso_.at(EG)->Fill(itr->hwIso());
0596 
0597           text << "EG : "
0598                << " BX=" << ibx << " ipt=" << itr->hwPt() << " ieta=" << itr->hwEta() << " iphi=" << itr->hwPhi()
0599                << std::endl;
0600 
0601           if (m_doEvtDisp)
0602             hEvtDemuxEG->Fill(itr->hwEta(), itr->hwPhi(), itr->hwPt());
0603         }
0604       }
0605     }
0606 
0607     // get tau
0608     if (m_doTaus) {
0609       Handle<BXVector<l1t::Tau> > taus;
0610       iEvent.getByToken(m_tauToken, taus);
0611 
0612       for (int ibx = taus->getFirstBX(); ibx <= taus->getLastBX(); ++ibx) {
0613         if (!m_allBx && ibx != m_dmxBx)
0614           continue;
0615 
0616         for (auto itr = taus->begin(ibx); itr != taus->end(ibx); ++itr) {
0617           hbx_.at(Tau)->Fill(ibx);
0618           het_.at(Tau)->Fill(itr->hwPt());
0619           heta_.at(Tau)->Fill(itr->hwEta());
0620           hphi_.at(Tau)->Fill(itr->hwPhi());
0621           hetaphi_.at(Tau)->Fill(itr->hwEta(), itr->hwPhi(), itr->hwPt());
0622 
0623           text << "Tau : "
0624                << " BX=" << ibx << " ipt=" << itr->hwPt() << " ieta=" << itr->hwEta() << " iphi=" << itr->hwPhi()
0625                << std::endl;
0626 
0627           if (m_doEvtDisp)
0628             hEvtDemuxTau->Fill(itr->hwEta(), itr->hwPhi(), itr->hwPt());
0629         }
0630       }
0631     }
0632 
0633     // get jet
0634     int njetdem = 0;
0635     std::vector<l1t::Jet> thejets;
0636 
0637     if (m_doJets) {
0638       Handle<BXVector<l1t::Jet> > jets;
0639       iEvent.getByToken(m_jetToken, jets);
0640 
0641       for (int ibx = jets->getFirstBX(); ibx <= jets->getLastBX(); ++ibx) {
0642         if (!m_allBx && ibx != m_dmxBx)
0643           continue;
0644 
0645         for (auto itr = jets->begin(ibx); itr != jets->end(ibx); ++itr) {
0646           njetdem += 1;
0647           hbx_.at(Jet)->Fill(ibx);
0648           het_.at(Jet)->Fill(itr->hwPt());
0649           heta_.at(Jet)->Fill(itr->hwEta());
0650           hphi_.at(Jet)->Fill(itr->hwPhi());
0651           hetaphi_.at(Jet)->Fill(itr->hwEta(), itr->hwPhi(), itr->hwPt());
0652 
0653           text << "Jet : "
0654                << " BX=" << ibx << " ipt=" << itr->hwPt() << " ieta=" << itr->hwEta() << " iphi=" << itr->hwPhi()
0655                << std::endl;
0656 
0657           if (m_doEvtDisp)
0658             hEvtDemuxJet->Fill(itr->hwEta(), itr->hwPhi(), itr->hwPt());
0659 
0660           thejets.push_back(*itr);
0661         }
0662       }
0663     }
0664 
0665     if (!thejets.empty()) {
0666       for (unsigned int i = 0; i < thejets.size() - 1; i++) {
0667         for (unsigned int j = i + 1; j < thejets.size(); j++) {
0668           hsort_->Fill(thejets.at(i).hwPt() - thejets.at(j).hwPt());
0669         }
0670       }
0671     }
0672 
0673     // get sums
0674     if (m_doSums) {
0675       Handle<BXVector<l1t::EtSum> > sums;
0676       iEvent.getByToken(m_sumToken, sums);
0677 
0678       for (int ibx = sums->getFirstBX(); ibx <= sums->getLastBX(); ++ibx) {
0679         if (!m_allBx && ibx != m_dmxBx)
0680           continue;
0681 
0682         for (auto itr = sums->begin(ibx); itr != sums->end(ibx); ++itr) {
0683           switch (itr->getType()) {
0684             case l1t::EtSum::EtSumType::kTotalEt:
0685               het_.at(SumET)->Fill(itr->hwPt());
0686               break;
0687             case l1t::EtSum::EtSumType::kTotalEtHF:
0688               break;
0689             case l1t::EtSum::EtSumType::kTotalEtEm:
0690               het_.at(SumETEm)->Fill(itr->hwPt());
0691               break;
0692             case l1t::EtSum::EtSumType::kTotalHt:
0693               het_.at(SumHT)->Fill(itr->hwPt());
0694               break;
0695             case l1t::EtSum::EtSumType::kTotalHtHF:
0696               break;
0697             case l1t::EtSum::EtSumType::kMissingEt:
0698               het_.at(SumMET)->Fill(itr->hwPt());
0699               hphi_.at(SumMET)->Fill(itr->hwPhi());
0700               break;
0701             case l1t::EtSum::EtSumType::kMissingEtHF:
0702               het_.at(SumMETHF)->Fill(itr->hwPt());
0703               hphi_.at(SumMETHF)->Fill(itr->hwPhi());
0704               break;
0705             case l1t::EtSum::EtSumType::kMissingHt:
0706               het_.at(SumMHT)->Fill(itr->hwPt());
0707               hphi_.at(SumMHT)->Fill(itr->hwPhi());
0708               break;
0709             case l1t::EtSum::EtSumType::kMissingHtHF:
0710               het_.at(SumMHTHF)->Fill(itr->hwPt());
0711               hphi_.at(SumMHTHF)->Fill(itr->hwPhi());
0712               break;
0713             case l1t::EtSum::EtSumType::kMinBiasHFP0:
0714               het_.at(MinBiasHFP0)->Fill(itr->hwPt());
0715               break;
0716             case l1t::EtSum::EtSumType::kMinBiasHFM0:
0717               het_.at(MinBiasHFM0)->Fill(itr->hwPt());
0718               break;
0719             case l1t::EtSum::EtSumType::kMinBiasHFP1:
0720               het_.at(MinBiasHFP1)->Fill(itr->hwPt());
0721               break;
0722             case l1t::EtSum::EtSumType::kMinBiasHFM1:
0723               het_.at(MinBiasHFM1)->Fill(itr->hwPt());
0724               break;
0725             case l1t::EtSum::EtSumType::kTowerCount:
0726               het_.at(SumHITowCount)->Fill(itr->hwPt());
0727               break;
0728             case l1t::EtSum::EtSumType::kCentrality:
0729               het_.at(SumCentrality)->Fill(itr->hwPt());
0730               break;
0731             case l1t::EtSum::EtSumType::kAsymEt:
0732               het_.at(SumAsymEt)->Fill(itr->hwPt());
0733               break;
0734             case l1t::EtSum::EtSumType::kAsymHt:
0735               het_.at(SumAsymHt)->Fill(itr->hwPt());
0736               break;
0737             case l1t::EtSum::EtSumType::kAsymEtHF:
0738               het_.at(SumAsymEtHF)->Fill(itr->hwPt());
0739               break;
0740             case l1t::EtSum::EtSumType::kAsymHtHF:
0741               het_.at(SumAsymHtHF)->Fill(itr->hwPt());
0742               break;
0743 
0744             default:
0745               std::cout << "wrong type of demux sum: " << itr->getType() << std::endl;
0746           }
0747 
0748           text << "Sum : "
0749                << " type=" << itr->getType() << " BX=" << ibx << " ipt=" << itr->hwPt() << " ieta=" << itr->hwEta()
0750                << " iphi=" << itr->hwPhi() << std::endl;
0751         }
0752       }
0753     }
0754 
0755     if (doText_)
0756       edm::LogVerbatim("L1TCaloEvents") << text.str();
0757   }
0758 
0759   // ------------ method called once each job just before starting event loop  ------------
0760   void L1TStage2CaloAnalyzer::beginJob() {
0761     edm::Service<TFileService> fs;
0762 
0763     auto itr = types_.cbegin();
0764     auto str = typeStr_.cbegin();
0765 
0766     for (; itr != types_.end(); ++itr, ++str) {
0767       dirs_.insert(std::pair<ObjectType, TFileDirectory>(*itr, fs->mkdir(*str)));
0768 
0769       if (*itr == MPSumMETx || *itr == MPSumMETxHF || *itr == MPSumMETy || *itr == MPSumMETyHF) {
0770         het_.insert(
0771             std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("et", "", 2200000, -2199999500, 2200000500)));
0772       } else if (*itr == MPSumMHTx || *itr == MPSumMHTxHF || *itr == MPSumMHTy || *itr == MPSumMHTyHF) {
0773         het_.insert(
0774             std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("et", "", 2200000, -2199999500, 2200000500)));
0775       } else if (*itr == SumET || *itr == SumETEm || *itr == MPSumETEm || *itr == MPSumET || *itr == MPSumETHF ||
0776                  *itr == SumHT || *itr == MPSumHT || *itr == MPSumHTHF) {
0777         het_.insert(std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("et", "", 100000, -0.5, 99999.5)));
0778       } else if (*itr == MPMinBiasHFP0 || *itr == MPMinBiasHFM0 || *itr == MPMinBiasHFP1 || *itr == MPMinBiasHFM1 ||
0779                  *itr == MinBiasHFP0 || *itr == MinBiasHFM0 || *itr == MinBiasHFP1 || *itr == MinBiasHFM1) {
0780         het_.insert(std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("et", "", 16, -0.5, 15.5)));
0781       } else if (*itr == MPSumHITowCount || *itr == SumHITowCount) {
0782         het_.insert(std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("et", "", 5904, -0.5, 5903.5)));
0783       } else if (*itr == EG || *itr == Tau || *itr == MPEG || *itr == MPTau) {
0784         het_.insert(std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("et", "", 1500, -0.5, 1499.5)));
0785       } else {
0786         het_.insert(std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("et", "", 100000, -0.5, 99999.5)));
0787       }
0788 
0789       hbx_.insert(std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("bx", "", 11, -5.5, 5.5)));
0790 
0791       if (*itr == EG || *itr == Jet || *itr == Tau) {
0792         heta_.insert(std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("eta", "", 227, -113.5, 113.5)));
0793         hphi_.insert(std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("phi", "", 144, -0.5, 143.5)));
0794         hetaphi_.insert(std::pair<ObjectType, TH2F*>(
0795             *itr, dirs_.at(*itr).make<TH2F>("etaphi", "", 227, -113.5, 113.5, 144, -0.5, 143.5)));
0796         if (*itr == EG)
0797           hiso_.insert(std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("iso", "", 4, -0.5, 3.5)));
0798       } else if (*itr == Tower || *itr == Cluster || *itr == MPEG || *itr == MPJet || *itr == MPTau) {
0799         heta_.insert(std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("eta", "", 83, -41.5, 41.5)));
0800         hphi_.insert(std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("phi", "", 72, 0.5, 72.5)));
0801         hetaphi_.insert(
0802             std::pair<ObjectType, TH2F*>(*itr, dirs_.at(*itr).make<TH2F>("etaphi", "", 83, -41.5, 41.5, 72, .5, 72.5)));
0803         if (*itr == MPEG)
0804           hiso_.insert(std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("iso", "", 4, -0.5, 3.5)));
0805       } else if (*itr == SumMET || *itr == SumMETHF || *itr == SumMHT || *itr == SumMHTHF) {
0806         hphi_.insert(std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("phi", "", 1008, -0.5, 1007.5)));
0807       }
0808 
0809       if (*itr == Tower) {
0810         hem_.insert(std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("em", "", 101, -0.5, 100.5)));
0811         hhad_.insert(std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("had", "", 101, -0.5, 100.5)));
0812         hratio_.insert(std::pair<ObjectType, TH1F*>(*itr, dirs_.at(*itr).make<TH1F>("ratio", "", 11, -0.5, 10.5)));
0813       }
0814     }
0815 
0816     if (m_doEvtDisp) {
0817       evtDispDir_ = fs->mkdir("Events");
0818     }
0819 
0820     hsort_ = fs->make<TH1F>("sort", "", 201, -100.5, 100.5);
0821     hsortMP_ = fs->make<TH1F>("sortMP", "", 201, -100.5, 100.5);
0822   }
0823 
0824   // ------------ method called once each job just after ending the event loop  ------------
0825   void L1TStage2CaloAnalyzer::endJob() {}
0826 
0827   // ------------ method called when starting to processes a run  ------------
0828   /*
0829     void 
0830     L1TStage2CaloAnalyzer::beginRun(edm::Run const&, edm::EventSetup const&)
0831     {
0832     }
0833   */
0834 
0835   // ------------ method called when ending the processing of a run  ------------
0836   /*
0837     void 
0838     L1TStage2CaloAnalyzer::endRun(edm::Run const&, edm::EventSetup const&)
0839     {
0840     }
0841   */
0842 
0843   // ------------ method called when starting to processes a luminosity block  ------------
0844   /*
0845     void 
0846     L1TStage2CaloAnalyzer::beginLuminosityBlock(edm::LuminosityBlock const&, edm::EventSetup const&)
0847     {
0848     }
0849   */
0850 
0851   // ------------ method called when ending the processing of a luminosity block  ------------
0852   /*
0853     void 
0854     L1TStage2CaloAnalyzer::endLuminosityBlock(edm::LuminosityBlock const&, edm::EventSetup const&)
0855     {
0856     }
0857   */
0858 
0859   // ------------ method fills 'descriptions' with the allowed parameters for the module  ------------
0860   void L1TStage2CaloAnalyzer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0861     //The following says we do not know what parameters are allowed so do no validation
0862     // Please change this to state exactly what you do use, even if it is no parameters
0863     edm::ParameterSetDescription desc;
0864     desc.setUnknown();
0865     descriptions.addDefault(desc);
0866   }
0867 
0868 }  // namespace l1t
0869 
0870 using namespace l1t;
0871 
0872 //define this as a plug-in
0873 DEFINE_FWK_MODULE(L1TStage2CaloAnalyzer);