Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-01-07 23:36:52

0001 /** \class DQMHistNormalizer
0002  *  
0003  *  Class to produce efficiency histograms by dividing nominator by denominator histograms
0004  *
0005  *  \author Christian Veelken, UC Davis
0006  */
0007 
0008 // framework & common header files
0009 #include "FWCore/Framework/interface/EDAnalyzer.h"
0010 #include "FWCore/Framework/interface/Event.h"
0011 #include "FWCore/Framework/interface/EventSetup.h"
0012 #include "FWCore/Framework/interface/MakerMacros.h"
0013 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0014 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0015 
0016 //DQM services
0017 #include "DQMServices/Core/interface/DQMStore.h"
0018 #include "FWCore/ServiceRegistry/interface/Service.h"
0019 
0020 //Regexp handling
0021 #include "classlib/utils/RegexpMatch.h"
0022 #include "classlib/utils/Regexp.h"
0023 
0024 #include <string>
0025 #include <vector>
0026 #include <map>
0027 
0028 using namespace std;
0029 
0030 class DQMHistNormalizer : public edm::EDAnalyzer {
0031 public:
0032   typedef dqm::legacy::DQMStore DQMStore;
0033   typedef dqm::legacy::MonitorElement MonitorElement;
0034 
0035   explicit DQMHistNormalizer(const edm::ParameterSet&);
0036   ~DQMHistNormalizer() override;
0037   void analyze(const edm::Event&, const edm::EventSetup&) override;
0038   void endRun(const edm::Run& r, const edm::EventSetup& c) override;
0039 
0040 private:
0041   lat::Regexp* buildRegex(const string& expr);
0042   vector<string> plotNamesToNormalize_;  //root name used by all the plots that must be normalized
0043   string reference_;
0044 };
0045 
0046 DQMHistNormalizer::DQMHistNormalizer(const edm::ParameterSet& cfg)
0047     : plotNamesToNormalize_(cfg.getParameter<std::vector<string> >("plotNamesToNormalize")),
0048       reference_(cfg.getParameter<string>("reference")) {
0049   //std::cout << "<DQMHistNormalizer::DQMHistNormalizer>:" << std::endl;
0050 }
0051 
0052 DQMHistNormalizer::~DQMHistNormalizer() {
0053   //--- nothing to be done yet
0054 }
0055 
0056 void DQMHistNormalizer::analyze(const edm::Event&, const edm::EventSetup&) {
0057   //--- nothing to be done yet
0058 }
0059 
0060 lat::Regexp* DQMHistNormalizer::buildRegex(const string& expr) {
0061   lat::Regexp* rx = nullptr;
0062   try {
0063     rx = new lat::Regexp(expr, 0, lat::Regexp::Wildcard);
0064     rx->study();
0065   } catch (lat::Error& e) {
0066     throw cms::Exception("DQMHistNormalizer")
0067         << "Invalid regular expression '" << expr.c_str() << "':" << e.explain().c_str();
0068   }
0069   return rx;
0070 }
0071 
0072 void DQMHistNormalizer::endRun(const edm::Run& r, const edm::EventSetup& c) {
0073   //std::cout << "<DQMHistNormalizer::endJob>:" << std::endl;
0074 
0075   //--- check that DQMStore service is available
0076   if (!edm::Service<DQMStore>().isAvailable()) {
0077     edm::LogError("endJob") << " Failed to access dqmStore --> histograms will NOT be plotted !!";
0078     return;
0079   }
0080 
0081   DQMStore& dqmStore = (*edm::Service<DQMStore>());
0082 
0083   vector<MonitorElement*> allOurMEs = dqmStore.getAllContents("RecoTauV/");
0084   lat::Regexp* refregex = buildRegex("*RecoTauV/*/" + reference_);
0085   vector<lat::Regexp*> toNormRegex;
0086   for (std::vector<string>::const_iterator toNorm = plotNamesToNormalize_.begin();
0087        toNorm != plotNamesToNormalize_.end();
0088        ++toNorm)
0089     toNormRegex.push_back(buildRegex("*RecoTauV/*/" + *toNorm));
0090 
0091   map<string, MonitorElement*> refsMap;
0092   vector<MonitorElement*> toNormElements;
0093 
0094   for (vector<MonitorElement*>::const_iterator element = allOurMEs.begin(); element != allOurMEs.end(); ++element) {
0095     string pathname = (*element)->getFullname();
0096     //cout << pathname << endl;
0097     //Matches reference
0098     if (refregex->match(pathname)) {
0099       //cout << "Matched to ref" << endl;
0100       string dir = pathname.substr(0, pathname.rfind('/'));
0101       if (refsMap.find(dir) != refsMap.end()) {
0102         edm::LogInfo("DQMHistNormalizer")
0103             << "DQMHistNormalizer::endRun: Warning! found multiple normalizing references for dir: " << dir << "!";
0104         edm::LogInfo("DQMHistNormalizer") << "     " << (refsMap[dir])->getFullname();
0105         edm::LogInfo("DQMHistNormalizer") << "     " << pathname;
0106       } else {
0107         refsMap[dir] = *element;
0108       }
0109     }
0110 
0111     //Matches targets
0112     for (vector<lat::Regexp*>::const_iterator reg = toNormRegex.begin(); reg != toNormRegex.end(); ++reg) {
0113       if ((*reg)->match(pathname)) {
0114         //cout << "Matched to target" << endl;
0115         toNormElements.push_back(*element);
0116         //cout << "Filled the collection" << endl;
0117       }
0118     }
0119   }
0120 
0121   delete refregex;
0122   for (vector<lat::Regexp*>::const_iterator reg = toNormRegex.begin(); reg != toNormRegex.end(); ++reg)
0123     delete *reg;
0124 
0125   for (vector<MonitorElement*>::const_iterator matchingElement = toNormElements.begin();
0126        matchingElement != toNormElements.end();
0127        ++matchingElement) {
0128     string meName = (*matchingElement)->getFullname();
0129     string dir = meName.substr(0, meName.rfind('/'));
0130 
0131     if (refsMap.find(dir) == refsMap.end()) {
0132       edm::LogInfo("DQMHistNormalizer") << "DQMHistNormalizer::endRun: Error! normalizing references for " << meName
0133                                         << " not found! Skipping...";
0134       continue;
0135     }
0136 
0137     float norm = refsMap[dir]->getTH1()->GetEntries();
0138     TH1* hist = (*matchingElement)->getTH1();
0139     if (norm != 0.) {
0140       if (!hist->GetSumw2N())
0141         hist->Sumw2();
0142       hist->Scale(1 / norm);  //use option "width" to divide the bin contents and errors by the bin width?
0143     } else {
0144       edm::LogInfo("DQMHistNormalizer") << "DQMHistNormalizer::endRun: Error! Normalization failed in "
0145                                         << hist->GetTitle() << "!";
0146     }
0147 
0148   }  //    for(vector<MonitorElement *>::const_iterator matchingElement = matchingElemts.begin(); matchingElement = matchingElemts.end(); ++matchingElement)
0149 }
0150 
0151 DEFINE_FWK_MODULE(DQMHistNormalizer);