Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:07:48

0001 /*
0002  * \file L1TCompare.cc
0003  * \author P. Wittich
0004  * \brief Compare different parts of the trigger chain (e.g., RCT-GCT )
0005  *
0006  *
0007  * organized message logger
0008  *
0009  * Revision 1.12  2008/03/14 20:35:46  berryhil
0010  *
0011  *
0012  * stripped out obsolete parameter settings
0013  *
0014  * rpc tpg restored with correct dn access and dbe handling
0015  *
0016  * Revision 1.11  2008/03/12 17:24:24  berryhil
0017  *
0018  *
0019  * eliminated log files, truncated HCALTPGXana histo output
0020  *
0021  * Revision 1.10  2008/03/01 00:40:00  lat
0022  * DQM core migration.
0023  *
0024  * Revision 1.9  2008/01/22 18:56:01  muzaffar
0025  * include cleanup. Only for cc/cpp files
0026  *
0027  * Revision 1.8  2007/12/21 20:04:50  wsun
0028  * Migrated L1EtMissParticle -> L1EtMissParticleCollection.
0029  *
0030  * Revision 1.7  2007/12/21 17:41:20  berryhil
0031  *
0032  *
0033  * try/catch removal
0034  *
0035  * Revision 1.6  2007/11/19 15:08:22  lorenzo
0036  * changed top folder name
0037  *
0038  * Revision 1.5  2007/09/27 22:58:15  ratnik
0039  * QA campaign: fixes to compensate includes cleanup in  DataFormats/L1Trigger
0040  *
0041  * Revision 1.4  2007/07/19 18:05:06  berryhil
0042  *
0043  *
0044  * L1CaloRegionDetId dataformat migration for L1TCompare
0045  *
0046  * Revision 1.3  2007/06/13 11:33:39  wittich
0047  * add axis titles
0048  *
0049  * Revision 1.2  2007/06/08 08:37:43  wittich
0050  * Add ECAL TP - RCT comparisons. Lingering problems with
0051  * mismatches right now - still needs work.
0052  *
0053  * Revision 1.1  2007/06/06 14:55:51  wittich
0054  * compare within trigger subsystems
0055  *
0056  */
0057 
0058 #include "DQM/L1TMonitor/interface/L1TCompare.h"
0059 
0060 // stl
0061 #include <algorithm>
0062 
0063 using namespace l1extra;
0064 using namespace edm;
0065 
0066 const unsigned int PHIBINS = 18;
0067 const float PHIMIN = -0.5;
0068 const float PHIMAX = 17.5;
0069 
0070 // Ranks 6, 10 and 12 bits
0071 const unsigned int R6BINS = 64;
0072 const float R6MIN = -0.5;
0073 const float R6MAX = 63.5;
0074 
0075 // For GCT this should be 15 bins, -14.5 to 14.5
0076 //
0077 const unsigned int ETABINS = 22;
0078 const float ETAMIN = -0.5;
0079 const float ETAMAX = 21.5;
0080 
0081 // TPG
0082 const unsigned int TPPHIBINS = 72;
0083 const float TPPHIMIN = 0.5;
0084 const float TPPHIMAX = 72.5;
0085 
0086 const unsigned int TPETABINS = 65;
0087 const float TPETAMIN = -32.5;
0088 const float TPETAMAX = 32.5;
0089 
0090 L1TCompare::L1TCompare(const ParameterSet& ps)
0091     : rctSourceEm_token_(consumes<L1CaloEmCollection>(ps.getParameter<InputTag>("rctSource"))),
0092       rctSourceRctEmRgn_token_(consumes<L1CaloRegionCollection>(ps.getParameter<InputTag>("rctSource"))),
0093       rctSource_(ps.getParameter<InputTag>("rctSource")),
0094       gctSource_(ps.getParameter<InputTag>("gctSource")),
0095       ecalTpgSource_(ps.getParameter<edm::InputTag>("ecalTpgSource")),
0096       ecalTpgSource_token_(consumes<EcalTrigPrimDigiCollection>(ps.getParameter<edm::InputTag>("ecalTpgSource")))
0097 
0098 {
0099   // verbosity switch
0100   verbose_ = ps.getUntrackedParameter<bool>("verbose", false);
0101 
0102   if (verbose())
0103     std::cout << "L1TCompare: constructor...." << std::endl;
0104 
0105   outputFile_ = ps.getUntrackedParameter<std::string>("outputFile", "");
0106   if (!outputFile_.empty()) {
0107     std::cout << "L1T Monitoring histograms will be saved to " << outputFile_.c_str() << std::endl;
0108   }
0109 
0110   bool disable = ps.getUntrackedParameter<bool>("disableROOToutput", false);
0111   if (disable) {
0112     outputFile_ = "";
0113   }
0114 
0115   //set Token(-s)
0116   edm::InputTag gctCenJetsTag_(gctSource_.label(), "cenJets");
0117   edm::InputTag gctIsoEmCandsTag_(gctSource_.label(), "isoEm");
0118   edm::InputTag gctNonIsoEmCandsTag_(gctSource_.label(), "nonIsoEm");
0119 
0120   gctCenJetsToken_ = consumes<L1GctJetCandCollection>(gctCenJetsTag_);
0121   gctIsoEmCandsToken_ = consumes<L1GctEmCandCollection>(gctIsoEmCandsTag_);
0122   gctNonIsoEmCandsToken_ = consumes<L1GctEmCandCollection>(gctNonIsoEmCandsTag_);
0123 }
0124 
0125 L1TCompare::~L1TCompare() {}
0126 
0127 void L1TCompare::dqmBeginRun(edm::Run const&, edm::EventSetup const&) {
0128   //
0129 }
0130 
0131 void L1TCompare::bookHistograms(DQMStore::IBooker& ibooker, edm::Run const& iRun, edm::EventSetup const& iSetup) {
0132   nev_ = 0;
0133 
0134   ibooker.setCurrentFolder("L1T/Compare");
0135 
0136   // -------------------------------------------
0137   // RCT-GCT
0138   // -------------------------------------------
0139   // Isolated
0140   rctGctLeadingIsoEmRank_ =
0141       ibooker.book2D("rctGctLeadingIsoEmRank", "RCT-GCT: rank", R6BINS, R6MIN, R6MAX, R6BINS, R6MIN, R6MAX);
0142   rctGctLeadingIsoEmRank_->setAxisTitle(std::string("gct"), 1);
0143   rctGctLeadingIsoEmRank_->setAxisTitle(std::string("rct"), 2);
0144   rctGctLeadingIsoEmEta_ =
0145       ibooker.book2D("rctGctLeadingIsoEmEta", "RCT-GCT: #eta", ETABINS, ETAMIN, ETAMAX, ETABINS, ETAMIN, ETAMAX);
0146   rctGctLeadingIsoEmEta_->setAxisTitle(std::string("gct"), 1);
0147   rctGctLeadingIsoEmEta_->setAxisTitle(std::string("rct"), 2);
0148 
0149   rctGctLeadingIsoEmPhi_ =
0150       ibooker.book2D("rctGctLeadingIsoEmPhi", "RCT-GCT: #phi", PHIBINS, PHIMIN, PHIMAX, PHIBINS, PHIMIN, PHIMAX);
0151   rctGctLeadingIsoEmPhi_->setAxisTitle(std::string("gct"), 1);
0152   rctGctLeadingIsoEmPhi_->setAxisTitle(std::string("rct"), 2);
0153   // non-Isolated
0154   rctGctLeadingNonIsoEmRank_ =
0155       ibooker.book2D("rctGctLeadingNonIsoEmRank", "RCT-GCT: rank", R6BINS, R6MIN, R6MAX, R6BINS, R6MIN, R6MAX);
0156   rctGctLeadingNonIsoEmRank_->setAxisTitle(std::string("gct"), 1);
0157   rctGctLeadingNonIsoEmRank_->setAxisTitle(std::string("rct"), 2);
0158 
0159   rctGctLeadingNonIsoEmEta_ =
0160       ibooker.book2D("rctGctLeadingNonIsoEmEta", "RCT-GCT: #eta", ETABINS, ETAMIN, ETAMAX, ETABINS, ETAMIN, ETAMAX);
0161   rctGctLeadingNonIsoEmEta_->setAxisTitle(std::string("gct"), 1);
0162   rctGctLeadingNonIsoEmEta_->setAxisTitle(std::string("rct"), 2);
0163 
0164   rctGctLeadingNonIsoEmPhi_ =
0165       ibooker.book2D("rctGctLeadingNonIsoEmPhi", "RCT-GCT: #phi", PHIBINS, PHIMIN, PHIMAX, PHIBINS, PHIMIN, PHIMAX);
0166   rctGctLeadingNonIsoEmPhi_->setAxisTitle(std::string("gct"), 1);
0167   rctGctLeadingNonIsoEmPhi_->setAxisTitle(std::string("rct"), 2);
0168   // -------------------------------------------
0169   // ECAL TPG - RCT
0170   // -------------------------------------------
0171   ecalTpgRctLeadingEmRank_ =
0172       ibooker.book2D("ecalTpgRctLeadingEmRank", "ECAL TPG-RCT: rank", R6BINS, R6MIN, R6MAX, R6BINS, R6MIN, R6MAX);
0173   ecalTpgRctLeadingEmRank_->setAxisTitle(std::string("rct"), 1);
0174   ecalTpgRctLeadingEmRank_->setAxisTitle(std::string("ecal tp"), 2);
0175 
0176   ecalTpgRctLeadingEmEta_ =
0177       ibooker.book2D("ecalTpgRctLeadingEmEta", "ECAL TPG-RCT: #eta", 15, -0.5, 14.5, TPETABINS, TPETAMIN, TPETAMAX);
0178   ecalTpgRctLeadingEmEta_->setAxisTitle(std::string("rct"), 1);
0179   ecalTpgRctLeadingEmEta_->setAxisTitle(std::string("ecal tp"), 2);
0180   ecalTpgRctLeadingEmEta2_ =
0181       ibooker.book2D("ecalTpgRctLeadingEmEta2", "ECAL TPG-RCT: #eta (2)", 13, -6.5, 6.5, TPETABINS, TPETAMIN, TPETAMAX);
0182   ecalTpgRctLeadingEmEta2_->setAxisTitle(std::string("rct"), 1);
0183   ecalTpgRctLeadingEmEta2_->setAxisTitle(std::string("ecal tp"), 2);
0184   ecalTpgRctLeadingEmPhi_ = ibooker.book2D(
0185       "ecalTpgRctLeadingEmPhi", "ECAL TPG-RCT: #phi", PHIBINS, PHIMIN, PHIMAX, TPPHIBINS, TPPHIMIN, TPPHIMAX);
0186   ecalTpgRctLeadingEmPhi_->setAxisTitle(std::string("rct"), 1);
0187   ecalTpgRctLeadingEmPhi_->setAxisTitle(std::string("ecal tp"), 2);
0188   //}
0189 }
0190 
0191 void L1TCompare::analyze(const Event& e, const EventSetup& c) {
0192   ++nev_;
0193   if (verbose()) {
0194     std::cout << "L1TCompare: analyze...." << std::endl;
0195   }
0196 
0197   // L1E
0198   edm::Handle<L1EmParticleCollection> l1eIsoEm;
0199   edm::Handle<L1EmParticleCollection> l1eNonIsoEm;
0200   edm::Handle<L1JetParticleCollection> l1eCenJets;
0201   edm::Handle<L1JetParticleCollection> l1eForJets;
0202   edm::Handle<L1JetParticleCollection> l1eTauJets;
0203   //  edm::Handle < L1EtMissParticle > l1eEtMiss;
0204   edm::Handle<L1EtMissParticleCollection> l1eEtMiss;
0205   // RCT
0206   edm::Handle<L1CaloEmCollection> em;  // collection of L1CaloEmCands
0207   edm::Handle<L1CaloRegionCollection> rctEmRgn;
0208 
0209   // GCT
0210   edm::Handle<L1GctJetCandCollection> gctCenJets;
0211   edm::Handle<L1GctEmCandCollection> gctIsoEmCands;
0212   edm::Handle<L1GctEmCandCollection> gctNonIsoEmCands;
0213 
0214   e.getByToken(rctSourceEm_token_, em);
0215 
0216   if (!em.isValid()) {
0217     edm::LogInfo("DataNotFound") << "can't find L1CaloEmCollection with label " << rctSource_.label();
0218     return;
0219   }
0220 
0221   e.getByToken(rctSourceRctEmRgn_token_, rctEmRgn);
0222 
0223   if (!rctEmRgn.isValid()) {
0224     edm::LogInfo("DataNotFound") << "can't find "
0225                                  << "L1CaloRegionCollection with label " << rctSource_.label();
0226     return;
0227   }
0228 
0229   e.getByToken(gctCenJetsToken_, gctCenJets);
0230   e.getByToken(gctIsoEmCandsToken_, gctIsoEmCands);
0231   e.getByToken(gctNonIsoEmCandsToken_, gctNonIsoEmCands);
0232 
0233   if (!gctCenJets.isValid()) {
0234     std::cerr << "L1TGCT: could not find one of the classes?" << std::endl;
0235     return;
0236   }
0237   if (!gctIsoEmCands.isValid()) {
0238     std::cerr << "L1TGCT: could not find one of the classes?" << std::endl;
0239     return;
0240   }
0241   if (!gctNonIsoEmCands.isValid()) {
0242     std::cerr << "L1TGCT: could not find one of the classes?" << std::endl;
0243     return;
0244   }
0245 
0246   // GCT
0247   if (verbose()) {
0248     for (L1GctEmCandCollection::const_iterator iem = gctIsoEmCands->begin(); iem != gctIsoEmCands->end(); ++iem) {
0249       if (!iem->empty())
0250         std::cout << "GCT EM: " << iem->rank() << ", " << iem->etaIndex()
0251                   << "("
0252                   //<< int(iem->etaIndex()&0x3)*((iem->etaIndex()&0x4)?1:-1)
0253                   << "), " << iem->phiIndex() << std::endl;
0254     }
0255   }
0256   // rct phi: 0-17
0257   // rct eta: 0-21
0258 
0259   // Fill the RCT histograms
0260 
0261   // Regions
0262   RctObjectCollection rcj, rcj_iso, rcj_non_iso;
0263   for (L1CaloEmCollection::const_iterator iem = em->begin(); iem != em->end(); ++iem) {
0264     //   L1CaloRegionDetId id(false, iem->rctCrate(), iem->rctCard(),
0265     //           iem->rctRegion());
0266     L1CaloRegionDetId id(iem->rctCrate(), iem->rctCard(), iem->rctRegion());
0267 
0268     // RctObject h(id.gctEta(), id.gctPhi(), iem->rank());
0269     RctObject h(id.rctEta(), id.rctPhi(), iem->rank());
0270     if (!iem->isolated())
0271       rcj_non_iso.push_back(h);
0272     else
0273       rcj_iso.push_back(h);
0274     rcj.push_back(h);
0275   }
0276   // not so smart but ...
0277   std::sort(rcj.begin(), rcj.end(), rctObjectComp);
0278   std::sort(rcj_non_iso.begin(), rcj_non_iso.end(), rctObjectComp);
0279   std::sort(rcj_iso.begin(), rcj_iso.end(), rctObjectComp);
0280   if (verbose()) {
0281     for (RctObjectCollection::reverse_iterator ij = rcj_iso.rbegin();
0282          ij != rcj_iso.rend() && ij != rcj_iso.rbegin() + 8;
0283          ++ij) {
0284       std::cout << "RCT cj: " << ij->rank_ << ", " << ij->eta_ << ", " << ij->phi_ << std::endl;
0285     }
0286   }
0287   L1GctEmCandCollection::const_iterator lead_em = gctIsoEmCands->begin();
0288   if (!lead_em->empty()) {  // equivalent to rank == 0
0289     rctGctLeadingIsoEmEta_->Fill(lead_em->etaIndex(), rcj_iso.rbegin()->eta_);
0290     rctGctLeadingIsoEmPhi_->Fill(lead_em->phiIndex(), rcj_iso.rbegin()->phi_);
0291     rctGctLeadingIsoEmRank_->Fill(lead_em->rank(), rcj_iso.rbegin()->rank_);
0292   }
0293 
0294   // non-isolated
0295   if (verbose()) {
0296     for (L1GctEmCandCollection::const_iterator iem = gctNonIsoEmCands->begin(); iem != gctNonIsoEmCands->end(); ++iem) {
0297       if (!iem->empty())
0298         std::cout << "GCT EM non: " << iem->rank() << ", "
0299                   << iem->etaIndex()  //<< "("
0300                                       //<< int(iem->etaIndex()&0x3)*((iem->etaIndex()&0x4)?1:-1)
0301                   //<< ")"
0302                   << ", " << iem->phiIndex() << std::endl;
0303     }
0304   }
0305   if (verbose()) {
0306     for (RctObjectCollection::reverse_iterator ij = rcj_non_iso.rbegin();
0307          ij != rcj_non_iso.rend() && ij != rcj_non_iso.rbegin() + 8;
0308          ++ij) {
0309       std::cout << "RCT cj non: " << ij->rank_ << ", " << ij->eta_ << ", " << ij->phi_ << std::endl;
0310     }
0311   }
0312   lead_em = gctNonIsoEmCands->begin();
0313   if (!lead_em->empty()) {  // equivalent to rank != 0
0314     rctGctLeadingNonIsoEmEta_->Fill(lead_em->etaIndex(), rcj_non_iso.rbegin()->eta_);
0315     rctGctLeadingNonIsoEmPhi_->Fill(lead_em->phiIndex(), rcj_non_iso.rbegin()->phi_);
0316     rctGctLeadingNonIsoEmRank_->Fill(lead_em->rank(), rcj_non_iso.rbegin()->rank_);
0317   }
0318 
0319   // ECAL TPG's to RCT EM
0320   edm::Handle<EcalTrigPrimDigiCollection> eTP;
0321   e.getByToken(ecalTpgSource_token_, eTP);
0322 
0323   if (!eTP.isValid()) {
0324     edm::LogInfo("DataNotFound") << "can't find EcalTrigPrimCollection with label " << ecalTpgSource_.label();
0325     return;
0326   }
0327   RctObjectCollection ecalobs;
0328   for (EcalTrigPrimDigiCollection::const_iterator ieTP = eTP->begin(); ieTP != eTP->end(); ieTP++) {
0329     ecalobs.push_back(RctObject(ieTP->id().ieta(), ieTP->id().iphi(), ieTP->compressedEt()));
0330   }
0331   std::sort(ecalobs.begin(), ecalobs.end(), rctObjectComp);
0332   if (verbose()) {
0333     for (RctObjectCollection::reverse_iterator ij = ecalobs.rbegin();
0334          ij != ecalobs.rend() && ij != ecalobs.rbegin() + 8;
0335          ++ij) {
0336       std::cout << "ECAL cj : " << ij->rank_ << ", " << ij->eta_ << ", " << ij->phi_ << std::endl;
0337     }
0338   }
0339   // abritrary cut
0340   if (rcj.rbegin()->rank_ > 4) {
0341     ecalTpgRctLeadingEmEta_->Fill(rcj.rbegin()->eta_, ecalobs.rbegin()->eta_);
0342     int e2 = (rcj.rbegin()->eta_ & 0x7UL) * ((rcj.rbegin()->eta_ & 0x8UL) ? 1 : -1);
0343     ecalTpgRctLeadingEmEta2_->Fill(e2, ecalobs.rbegin()->eta_);
0344     ecalTpgRctLeadingEmPhi_->Fill(rcj.rbegin()->phi_, ecalobs.rbegin()->phi_);
0345     ecalTpgRctLeadingEmRank_->Fill(rcj.rbegin()->rank_, ecalobs.rbegin()->rank_);
0346   }
0347   if (verbose()) {
0348     int seta = rcj.rbegin()->eta_;
0349     seta = (seta & 0x7UL) * (seta & 0x8 ? -1 : 1);
0350     std::cout << "ZZ: " << rcj.rbegin()->eta_ << " " << rcj.rbegin()->phi_ << " " << rcj.rbegin()->rank_ << " "
0351               << (++rcj.rbegin())->rank_ << " " << ecalobs.rbegin()->eta_ << " " << ecalobs.rbegin()->phi_ << " "
0352               << ecalobs.rbegin()->rank_ << " " << (++ecalobs.rbegin())->rank_ << " " << seta << " " << std::endl;
0353   }
0354 }