Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-06-30 23:17:16

0001 #ifndef RecoLuminosity_LumiProducer_LumiCalculator_h
0002 #define RecoLuminosity_LumiProducer_LumiCalculator_h
0003 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0004 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0005 #include "FWCore/Framework/interface/Event.h"
0006 #include "FWCore/Framework/interface/Run.h"
0007 #include "FWCore/Framework/interface/LuminosityBlock.h"
0008 #include "DataFormats/Common/interface/Handle.h"
0009 #include "DataFormats/Luminosity/interface/LumiSummaryRunHeader.h"
0010 #include "DataFormats/Luminosity/interface/LumiSummary.h"
0011 #include "FWCore/Framework/interface/MakerMacros.h"
0012 #include "HLTrigger/HLTcore/interface/HLTConfigProvider.h"
0013 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0014 #include <boost/regex.hpp>
0015 #include <iostream>
0016 #include <map>
0017 #include "FWCore/Utilities/interface/EDGetToken.h"
0018 #include "FWCore/Utilities/interface/InputTag.h"
0019 struct hltPerPathInfo {
0020   hltPerPathInfo() : prescale(0) {}
0021   unsigned int prescale;
0022 };
0023 struct l1PerBitInfo {
0024   l1PerBitInfo() : prescale(0) {}
0025   unsigned int prescale;
0026 };
0027 struct MyPerLumiInfo {
0028   unsigned int lsnum;
0029   float livefraction;
0030   float intglumi;
0031   unsigned long long deadcount;
0032 };
0033 
0034 class LumiCalculator : public edm::one::EDAnalyzer<edm::one::WatchRuns, edm::one::WatchLuminosityBlocks> {
0035 public:
0036   explicit LumiCalculator(edm::ParameterSet const& pset);
0037   ~LumiCalculator() override;
0038 
0039 private:
0040   void beginJob() override;
0041   void beginRun(const edm::Run& run, const edm::EventSetup& c) override;
0042   void analyze(edm::Event const& e, edm::EventSetup const& c) override;
0043   void beginLuminosityBlock(edm::LuminosityBlock const& lumiBlock, edm::EventSetup const& c) override {}
0044   void endLuminosityBlock(edm::LuminosityBlock const& lumiBlock, edm::EventSetup const& c) override;
0045   void endRun(edm::Run const&, edm::EventSetup const&) override;
0046   void endJob() override;
0047   std::vector<std::string> splitpathstr(const std::string& strValue, const std::string separator);
0048   HLTConfigProvider hltConfig_;
0049   std::multimap<std::string, std::string> trgpathMmap_;  //key:hltpath,value:l1bit
0050   std::map<std::string, hltPerPathInfo> hltmap_;
0051   std::map<std::string, l1PerBitInfo> l1map_;  //
0052   std::vector<MyPerLumiInfo> perrunlumiinfo_;
0053 
0054 private:
0055   edm::LogInfo* log_;
0056   bool showTrgInfo_;
0057   unsigned int currentlumi_;
0058 };  //end class
0059 
0060 // -----------------------------------------------------------------
0061 
0062 LumiCalculator::LumiCalculator(edm::ParameterSet const& pset) : log_(new edm::LogInfo("LumiReport")), currentlumi_(0) {
0063   showTrgInfo_ = pset.getUntrackedParameter<bool>("showTriggerInfo", false);
0064   consumes<LumiSummary, edm::InLumi>(edm::InputTag("lumiProducer", ""));
0065   consumes<LumiSummaryRunHeader, edm::InRun>(edm::InputTag("lumiProducer", ""));
0066 }
0067 
0068 // -----------------------------------------------------------------
0069 
0070 LumiCalculator::~LumiCalculator() {
0071   delete log_;
0072   log_ = nullptr;
0073 }
0074 
0075 // -----------------------------------------------------------------
0076 
0077 void LumiCalculator::analyze(edm::Event const& e, edm::EventSetup const&) {}
0078 
0079 // -----------------------------------------------------------------
0080 
0081 void LumiCalculator::beginJob() {}
0082 
0083 // -----------------------------------------------------------------
0084 
0085 void LumiCalculator::beginRun(const edm::Run& run, const edm::EventSetup& c) {
0086   //std::cout<<"I'm in run number "<<run.run()<<std::endl;
0087   //if(!hltConfig_.init("HLT")){
0088   // throw cms::Exception("HLT process cannot be initialized");
0089   //}
0090   bool changed(true);
0091   const std::string processname("HLT");
0092   if (!hltConfig_.init(run, c, processname, changed)) {
0093     throw cms::Exception("HLT process cannot be initialized");
0094   }
0095   perrunlumiinfo_.clear();
0096   trgpathMmap_.clear();
0097   hltmap_.clear();
0098   l1map_.clear();
0099   //hltConfig_.dump("processName");
0100   //hltConfig_.dump("TableName");
0101   //hltConfig_.dump("Triggers");
0102   //hltConfig_.dump("Modules");
0103   if (showTrgInfo_) {
0104     *log_ << "======Trigger Configuration Overview======\n";
0105     *log_ << "Run " << run.run() << " Trigger Table : " << hltConfig_.tableName() << "\n";
0106   }
0107   unsigned int totaltrg = hltConfig_.size();
0108   for (unsigned int t = 0; t < totaltrg; ++t) {
0109     std::string hltname(hltConfig_.triggerName(t));
0110     std::vector<std::string> numpathmodules = hltConfig_.moduleLabels(hltname);
0111     if (showTrgInfo_) {
0112       *log_ << t << " HLT path\t" << hltname << "\n";
0113     }
0114     hltPerPathInfo hlt;
0115     hlt.prescale = 1;
0116     hltmap_.insert(std::make_pair(hltname, hlt));
0117     std::vector<std::string>::iterator hltpathBeg = numpathmodules.begin();
0118     std::vector<std::string>::iterator hltpathEnd = numpathmodules.end();
0119     unsigned int mycounter = 0;
0120     for (std::vector<std::string>::iterator numpathmodule = hltpathBeg; numpathmodule != hltpathEnd; ++numpathmodule) {
0121       if (hltConfig_.moduleType(*numpathmodule) != "HLTLevel1GTSeed") {
0122         continue;
0123       }
0124       ++mycounter;
0125 
0126       edm::ParameterSet l1GTPSet = hltConfig_.modulePSet(*numpathmodule);
0127       std::string l1pathname = l1GTPSet.getParameter<std::string>("L1SeedsLogicalExpression");
0128       //*log_<<"numpathmodule "<< *numpathmodule <<" : l1pathname "<<l1pathname<<"\n";
0129       if (mycounter > 1) {
0130         if (showTrgInfo_) {
0131           *log_ << "\tskip and erase previous seeds : multiple L1SeedsLogicalExpressions per hlt path\n";
0132         }
0133         //erase all previously calculated seeds for this path
0134         trgpathMmap_.erase(hltname);
0135         continue;
0136       }
0137       if (l1pathname.find('(') != std::string::npos) {
0138         if (showTrgInfo_) {
0139           *log_ << "  L1SeedsLogicalExpression(Complex)\t" << l1pathname << "\n";
0140           *log_ << "\tskip:contain complex logic\n";
0141         }
0142         continue;
0143       } else if (l1pathname.find("OR") != std::string::npos) {
0144         if (showTrgInfo_) {
0145           *log_ << "  L1SeedsLogicalExpression(ORed)\t" << l1pathname << "\n";
0146         }
0147         std::vector<std::string> seeds = splitpathstr(l1pathname, " OR ");
0148         if (seeds.size() > 2) {
0149           if (showTrgInfo_) {
0150             *log_ << "\tskip:contain >1 OR\n";
0151           }
0152           continue;
0153         } else {
0154           for (std::vector<std::string>::iterator i = seeds.begin(); i != seeds.end(); ++i) {
0155             if (!i->empty() && showTrgInfo_)
0156               *log_ << "\t\tseed: " << *i << "\n";
0157             if (i == seeds.begin()) {  //for now we take the first one from OR
0158               trgpathMmap_.insert(std::make_pair(hltname, *i));
0159             }
0160           }
0161         }
0162       } else if (l1pathname.find("AND") != std::string::npos) {
0163         if (showTrgInfo_) {
0164           *log_ << "  L1SeedsLogicalExpression(ANDed)\t" << l1pathname << "\n";
0165         }
0166         std::vector<std::string> seeds = splitpathstr(l1pathname, " AND ");
0167         if (seeds.size() > 2) {
0168           if (showTrgInfo_) {
0169             *log_ << "\tskip:contain >1 AND\n";
0170           }
0171           continue;
0172         } else {
0173           for (std::vector<std::string>::iterator i = seeds.begin(); i != seeds.end(); ++i) {
0174             if (!i->empty() && showTrgInfo_)
0175               *log_ << "\t\tseed: " << *i << "\n";
0176             if (i == seeds.begin()) {  //for now we take the first one
0177               trgpathMmap_.insert(std::make_pair(hltname, *i));
0178             }
0179           }
0180         }
0181       } else {
0182         if (showTrgInfo_) {
0183           *log_ << "  L1SeedsLogicalExpression(ONE)\t" << l1pathname << "\n";
0184         }
0185         if (splitpathstr(l1pathname, " NOT ").size() > 1) {
0186           if (showTrgInfo_) {
0187             *log_ << "\tskip:contain NOT\n";
0188           }
0189           continue;
0190         }
0191         trgpathMmap_.insert(std::make_pair(hltname, l1pathname));
0192       }
0193     }
0194   }
0195   if (showTrgInfo_) {
0196     *log_ << "================\n";
0197   }
0198 }
0199 
0200 // -----------------------------------------------------------------
0201 void LumiCalculator::endLuminosityBlock(edm::LuminosityBlock const& lumiBlock, edm::EventSetup const& c) {
0202   /**Integrated Luminosity per Lumi Section
0203      instantaneousLumi*93.244(sec)
0204   **/
0205   //std::cout<<"I'm in lumi block "<<lumiBlock.id()<<std::endl;
0206 
0207   edm::Handle<LumiSummary> lumiSummary;
0208   lumiBlock.getByLabel("lumiProducer", lumiSummary);
0209 
0210   edm::Handle<LumiSummaryRunHeader> lumiSummaryRH;
0211   lumiBlock.getRun().getByLabel("lumiProducer", lumiSummaryRH);
0212 
0213   MyPerLumiInfo l;
0214   l.lsnum = lumiBlock.id().luminosityBlock();
0215 
0216   //
0217   //collect lumi.
0218   //
0219   l.deadcount = lumiSummary->deadcount();
0220   l.intglumi = lumiSummary->avgInsDelLumi() * 93.244;
0221   l.livefraction = lumiSummary->liveFrac();
0222 
0223   *log_ << "====== Lumi Section " << lumiBlock.id().luminosityBlock() << " ======\n";
0224   *log_ << "\t Luminosity " << l.intglumi << "\n";
0225   *log_ << "\t Dead count " << l.deadcount << "\n";
0226   *log_ << "\t Deadtime corrected Luminosity " << l.intglumi * l.livefraction << "\n";
0227 
0228   //
0229   //print correlated hlt-l1 info, only if you ask
0230   //
0231   if (showTrgInfo_) {
0232     std::map<std::string, hltPerPathInfo>::iterator hltit;
0233     std::map<std::string, hltPerPathInfo>::iterator hltitBeg = hltmap_.begin();
0234     std::map<std::string, hltPerPathInfo>::iterator hltitEnd = hltmap_.end();
0235 
0236     typedef std::pair<std::multimap<std::string, std::string>::iterator,
0237                       std::multimap<std::string, std::string>::iterator>
0238         TRGMAPIT;
0239     unsigned int c = 0;
0240     for (hltit = hltitBeg; hltit != hltitEnd; ++hltit) {
0241       std::string hltname = hltit->first;
0242       *log_ << c << " HLT path  " << hltname << " , prescale : " << hltit->second.prescale << "\n";
0243       TRGMAPIT ppp;
0244       ppp = trgpathMmap_.equal_range(hltname);
0245       if (ppp.first == ppp.second) {
0246         *log_ << "    no L1\n";
0247       }
0248       for (std::multimap<std::string, std::string>::iterator mit = ppp.first; mit != ppp.second; ++mit) {
0249         std::string l1name = mit->second;
0250         *log_ << "    L1 name : " << l1name;
0251         LumiSummary::L1 l1result = lumiSummary->l1info(lumiSummaryRH->getL1Index(l1name));
0252         *log_ << " prescale : " << l1result.prescale << "\n";
0253         *log_ << "\n";
0254       }
0255       ++c;
0256     }
0257   }
0258   //
0259   //accumulate hlt counts. Absent for now
0260   //
0261   /**
0262      for(hltit=hltitBeg;hltit!=hltitEnd;++hltit){
0263      }
0264   **/
0265 
0266   //
0267   //accumulate l1 counts
0268   //
0269   size_t n = lumiSummary->nTriggerLine();
0270   for (size_t i = 0; i < n; ++i) {
0271     std::string l1bitname = lumiSummaryRH->getL1Name(lumiSummary->l1info(i).triggernameidx);
0272     l1PerBitInfo t;
0273     if (currentlumi_ == 0) {
0274       t.prescale = lumiSummary->l1info(i).prescale;
0275       l1map_.insert(std::make_pair(l1bitname, t));
0276     }
0277   }
0278 
0279   perrunlumiinfo_.push_back(l);
0280 
0281   ++currentlumi_;
0282 }
0283 
0284 // -----------------------------------------------------------------
0285 void LumiCalculator::endRun(edm::Run const& run, edm::EventSetup const& c) {
0286   /**Notes on calculation:
0287      
0288   1. CMS recorded Luminosity per run :
0289      sum over HLX&&HF certified LS 
0290      lumiSummary->avgInsDelLumi()*93.244*livefraction()
0291 
0292      2. Effective Luminosity per run per trigger line:
0293      For the moment, we take only the first L1 seed in case of 'OR' or 'AND'
0294      relationship between HLT and L1 seeds
0295      
0296      avgInsDelLumi()*93.244*livefraction()/(HLTprescale*L1prescale)
0297      for now HLTprescale=1
0298      
0299      3. LHC delivered:
0300      there is no point in calculating delivered when data do not contain all LS
0301   **/
0302   //std::cout<<"valid trigger lines "<<trgpathMmap_.size()<<std::endl;
0303   //std::cout<<"total lumi lines "<<perrunlumiinfo_.size()<<std::endl;
0304   std::vector<MyPerLumiInfo>::const_iterator lumiIt;
0305   std::vector<MyPerLumiInfo>::const_iterator lumiItBeg = perrunlumiinfo_.begin();
0306   std::vector<MyPerLumiInfo>::const_iterator lumiItEnd = perrunlumiinfo_.end();
0307   float recorded = 0.0;
0308 
0309   *log_ << "================ Run Summary " << run.run() << "================\n";
0310   for (lumiIt = lumiItBeg; lumiIt != lumiItEnd; ++lumiIt) {  //loop over LS
0311     recorded += lumiIt->intglumi * lumiIt->livefraction;
0312   }
0313   *log_ << "  CMS Recorded Lumi (e+27cm^-2) : " << recorded << "\n";
0314   *log_ << "  Effective Lumi (e+27cm^-2) per trigger path: "
0315         << "\n\n";
0316   std::multimap<std::string, std::string>::iterator it;
0317   std::multimap<std::string, std::string>::iterator itBeg = trgpathMmap_.begin();
0318   std::multimap<std::string, std::string>::iterator itEnd = trgpathMmap_.end();
0319   unsigned int cc = 0;
0320   for (it = itBeg; it != itEnd; ++it) {
0321     *log_ << "  " << cc << "  " << it->first << " - " << it->second << " : ";
0322     ++cc;
0323     std::map<std::string, hltPerPathInfo>::const_iterator hltIt = hltmap_.find(it->first);
0324     if (hltIt == hltmap_.end()) {
0325       std::cout << "HLT path " << it->first << " not found" << std::endl;
0326       *log_ << "\n";
0327       continue;
0328     }
0329     std::map<std::string, l1PerBitInfo>::const_iterator l1It = l1map_.find(it->second);
0330     if (l1It == l1map_.end()) {
0331       std::cout << "L1 bit " << it->second << " not found" << std::endl;
0332       *log_ << "\n";
0333       continue;
0334     }
0335     unsigned int hltprescale = hltIt->second.prescale;
0336     unsigned int l1prescale = l1It->second.prescale;
0337     if (hltprescale != 0 && l1prescale != 0) {
0338       float effectiveLumi = recorded / (hltprescale * l1prescale);
0339       *log_ << effectiveLumi << "\n";
0340     } else {
0341       *log_ << "0 prescale exception\n";
0342       continue;
0343     }
0344     *log_ << "\n";
0345   }
0346 }
0347 
0348 // -----------------------------------------------------------------
0349 void LumiCalculator::endJob() {}
0350 
0351 std::vector<std::string> LumiCalculator::splitpathstr(const std::string& strValue, const std::string separator) {
0352   std::vector<std::string> vecstrResult;
0353   boost::regex re(separator);
0354   boost::sregex_token_iterator p(strValue.begin(), strValue.end(), re, -1);
0355   boost::sregex_token_iterator end;
0356   while (p != end) {
0357     vecstrResult.push_back(*p++);
0358   }
0359   return vecstrResult;
0360 }
0361 
0362 DEFINE_FWK_MODULE(LumiCalculator);
0363 #endif