Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 13:07:30

0001 /*
0002  *  Class:DQMRivetClient 
0003  *
0004  *
0005  * 
0006  *  \author Junghwan Goh - SungKyunKwan University
0007  */
0008 
0009 #include "GeneratorInterface/RivetInterface/interface/DQMRivetClient.h"
0010 
0011 #include "DQMServices/Core/interface/DQMStore.h"
0012 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0013 #include "FWCore/ServiceRegistry/interface/Service.h"
0014 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0015 
0016 #include <TH1F.h>
0017 #include <TClass.h>
0018 #include <TString.h>
0019 #include <TPRegexp.h>
0020 
0021 #include <cmath>
0022 #include <boost/tokenizer.hpp>
0023 
0024 using namespace std;
0025 using namespace edm;
0026 
0027 typedef DQMRivetClient::MonitorElement ME;
0028 
0029 DQMRivetClient::DQMRivetClient(const ParameterSet& pset) {
0030   typedef std::vector<edm::ParameterSet> VPSet;
0031   typedef std::vector<std::string> vstring;
0032   typedef boost::escaped_list_separator<char> elsc;
0033 
0034   elsc commonEscapes("\\", " \t", "\'");
0035 
0036   // Parse Normalization commands
0037   vstring normCmds = pset.getUntrackedParameter<vstring>("normalizationToIntegral", vstring());
0038   for (vstring::const_iterator normCmd = normCmds.begin(); normCmd != normCmds.end(); ++normCmd) {
0039     if (normCmd->empty())
0040       continue;
0041     boost::tokenizer<elsc> tokens(*normCmd, commonEscapes);
0042 
0043     vector<string> args;
0044     for (boost::tokenizer<elsc>::const_iterator iToken = tokens.begin(); iToken != tokens.end(); ++iToken) {
0045       if (iToken->empty())
0046         continue;
0047       args.push_back(*iToken);
0048     }
0049 
0050     if (args.empty() or args.size() > 2) {
0051       LogInfo("DQMRivetClient") << "Wrong input to normCmds\n";
0052       continue;
0053     }
0054 
0055     NormOption opt;
0056     opt.name = args[0];
0057     opt.normHistName = args.size() == 2 ? args[1] : args[0];
0058 
0059     normOptions_.push_back(opt);
0060   }
0061 
0062   VPSet normSets = pset.getUntrackedParameter<VPSet>("normalizationToIntegralSets", VPSet());
0063   for (VPSet::const_iterator normSet = normSets.begin(); normSet != normSets.end(); ++normSet) {
0064     NormOption opt;
0065     opt.name = normSet->getUntrackedParameter<string>("name");
0066     opt.normHistName = normSet->getUntrackedParameter<string>("normalizedTo", opt.name);
0067 
0068     normOptions_.push_back(opt);
0069   }
0070 
0071   //normalize to lumi
0072   vstring lumiCmds = pset.getUntrackedParameter<vstring>("normalizationToLumi", vstring());
0073   for (vstring::const_iterator lumiCmd = lumiCmds.begin(); lumiCmd != lumiCmds.end(); ++lumiCmd) {
0074     if (lumiCmd->empty())
0075       continue;
0076     boost::tokenizer<elsc> tokens(*lumiCmd, commonEscapes);
0077 
0078     vector<string> args;
0079     for (boost::tokenizer<elsc>::const_iterator iToken = tokens.begin(); iToken != tokens.end(); ++iToken) {
0080       if (iToken->empty())
0081         continue;
0082       args.push_back(*iToken);
0083     }
0084 
0085     if (args.size() != 2) {
0086       LogInfo("DQMRivetClient") << "Wrong input to lumiCmds\n";
0087       continue;
0088     }
0089 
0090     DQMRivetClient::LumiOption opt;
0091     opt.name = args[0];
0092     opt.normHistName = args[1];
0093     opt.xsection = pset.getUntrackedParameter<double>("xsection", -1.);
0094     //opt.xsection = atof(args[2].c_str());
0095 
0096     //std::cout << opt.name << " " << opt.normHistName << " " << opt.xsection << std::endl;
0097     lumiOptions_.push_back(opt);
0098   }
0099 
0100   //multiply by a number
0101   vstring scaleCmds = pset.getUntrackedParameter<vstring>("scaleBy", vstring());
0102   for (vstring::const_iterator scaleCmd = scaleCmds.begin(); scaleCmd != scaleCmds.end(); ++scaleCmd) {
0103     if (scaleCmd->empty())
0104       continue;
0105     boost::tokenizer<elsc> tokens(*scaleCmd, commonEscapes);
0106 
0107     vector<string> args;
0108     for (boost::tokenizer<elsc>::const_iterator iToken = tokens.begin(); iToken != tokens.end(); ++iToken) {
0109       if (iToken->empty())
0110         continue;
0111       args.push_back(*iToken);
0112     }
0113 
0114     if (args.empty() or args.size() > 2) {
0115       LogInfo("DQMRivetClient") << "Wrong input to normCmds\n";
0116       continue;
0117     }
0118 
0119     ScaleFactorOption opt;
0120     opt.name = args[0];
0121     opt.scale = atof(args[1].c_str());
0122     scaleOptions_.push_back(opt);
0123   }
0124 
0125   outputFileName_ = pset.getUntrackedParameter<string>("outputFileName", "");
0126   subDirs_ = pset.getUntrackedParameter<vstring>("subDirs");
0127 }
0128 
0129 void DQMRivetClient::endRun(const edm::Run& r, const edm::EventSetup& c) {
0130   typedef vector<string> vstring;
0131 
0132   // Update 2009-09-23
0133   // Migrated all code from endJob to this function
0134   // endJob is not necessarily called in the proper sequence
0135   // and does not necessarily book histograms produced in
0136   // that step.
0137   // It more robust to do the histogram manipulation in
0138   // this endRun function
0139 
0140   theDQM = nullptr;
0141   theDQM = Service<DQMStore>().operator->();
0142 
0143   if (!theDQM) {
0144     LogInfo("DQMRivetClient") << "Cannot create DQMStore instance\n";
0145     return;
0146   }
0147 
0148   // Process wildcard in the sub-directory
0149   set<string> subDirSet;
0150 
0151   for (vstring::const_iterator iSubDir = subDirs_.begin(); iSubDir != subDirs_.end(); ++iSubDir) {
0152     string subDir = *iSubDir;
0153 
0154     if (subDir[subDir.size() - 1] == '/')
0155       subDir.erase(subDir.size() - 1);
0156 
0157     subDirSet.insert(subDir);
0158   }
0159 
0160   for (set<string>::const_iterator iSubDir = subDirSet.begin(); iSubDir != subDirSet.end(); ++iSubDir) {
0161     const string& dirName = *iSubDir;
0162     for (vector<NormOption>::const_iterator normOption = normOptions_.begin(); normOption != normOptions_.end();
0163          ++normOption) {
0164       normalizeToIntegral(dirName, normOption->name, normOption->normHistName);
0165     }
0166   }
0167 
0168   for (set<string>::const_iterator iSubDir = subDirSet.begin(); iSubDir != subDirSet.end(); ++iSubDir) {
0169     const string& dirName = *iSubDir;
0170     for (vector<LumiOption>::const_iterator lumiOption = lumiOptions_.begin(); lumiOption != lumiOptions_.end();
0171          ++lumiOption) {
0172       normalizeToLumi(dirName, lumiOption->name, lumiOption->normHistName, lumiOption->xsection);
0173     }
0174   }
0175 
0176   for (set<string>::const_iterator iSubDir = subDirSet.begin(); iSubDir != subDirSet.end(); ++iSubDir) {
0177     const string& dirName = *iSubDir;
0178     for (vector<ScaleFactorOption>::const_iterator scaleOption = scaleOptions_.begin();
0179          scaleOption != scaleOptions_.end();
0180          ++scaleOption) {
0181       scaleByFactor(dirName, scaleOption->name, scaleOption->scale);
0182     }
0183   }
0184 
0185   if (!outputFileName_.empty())
0186     theDQM->save(outputFileName_);
0187 }
0188 
0189 void DQMRivetClient::endJob() {
0190   // Update 2009-09-23
0191   // Migrated all code from here to endRun
0192 
0193   LogTrace("DQMRivetClient") << "inside of ::endJob()" << endl;
0194 }
0195 
0196 void DQMRivetClient::normalizeToIntegral(const std::string& startDir,
0197                                          const std::string& histName,
0198                                          const std::string& normHistName) {
0199   if (!theDQM->dirExists(startDir)) {
0200     LogInfo("DQMRivetClient") << "normalizeToEntries() : "
0201                               << "Cannot find sub-directory " << startDir << endl;
0202     return;
0203   }
0204 
0205   theDQM->cd();
0206 
0207   ME* element = theDQM->get(startDir + "/" + histName);
0208   ME* normME = theDQM->get(startDir + "/" + normHistName);
0209 
0210   if (!element) {
0211     LogInfo("DQMRivetClient") << "normalizeToEntries() : "
0212                               << "No such element '" << histName << "' found\n";
0213     return;
0214   }
0215 
0216   if (!normME) {
0217     LogInfo("DQMRivetClient") << "normalizeToEntries() : "
0218                               << "No such element '" << normHistName << "' found\n";
0219     return;
0220   }
0221 
0222   TH1F* hist = element->getTH1F();
0223   if (!hist) {
0224     LogInfo("DQMRivetClient") << "normalizeToEntries() : "
0225                               << "Cannot create TH1F from ME\n";
0226     return;
0227   }
0228 
0229   TH1F* normHist = normME->getTH1F();
0230   if (!normHist) {
0231     LogInfo("DQMRivetClient") << "normalizeToEntries() : "
0232                               << "Cannot create TH1F from ME\n";
0233     return;
0234   }
0235 
0236   const double entries = normHist->Integral();
0237   if (entries != 0) {
0238     hist->Scale(1. / entries, "width");
0239   } else {
0240     LogInfo("DQMRivetClient") << "normalizeToEntries() : "
0241                               << "Zero entries in histogram\n";
0242   }
0243 
0244   return;
0245 }
0246 
0247 void DQMRivetClient::normalizeToLumi(const std::string& startDir,
0248                                      const std::string& histName,
0249                                      const std::string& normHistName,
0250                                      double xsection) {
0251   normalizeToIntegral(startDir, histName, normHistName);
0252   theDQM->cd();
0253   ME* element = theDQM->get(startDir + "/" + histName);
0254   TH1F* hist = element->getTH1F();
0255   if (!hist) {
0256     LogInfo("DQMRivetClient") << "normalizeToEntries() : "
0257                               << "Cannot create TH1F from ME\n";
0258     return;
0259   }
0260   hist->Scale(xsection);
0261   return;
0262 }
0263 
0264 void DQMRivetClient::scaleByFactor(const std::string& startDir, const std::string& histName, double factor) {
0265   if (!theDQM->dirExists(startDir)) {
0266     LogInfo("DQMRivetClient") << "normalizeToEntries() : "
0267                               << "Cannot find sub-directory " << startDir << endl;
0268     return;
0269   }
0270 
0271   theDQM->cd();
0272 
0273   ME* element = theDQM->get(startDir + "/" + histName);
0274 
0275   if (!element) {
0276     LogInfo("DQMRivetClient") << "normalizeToEntries() : "
0277                               << "No such element '" << histName << "' found\n";
0278     return;
0279   }
0280 
0281   TH1F* hist = element->getTH1F();
0282   if (!hist) {
0283     LogInfo("DQMRivetClient") << "normalizeToEntries() : "
0284                               << "Cannot create TH1F from ME\n";
0285     return;
0286   }
0287   hist->Scale(factor);
0288 }
0289 
0290 #include "FWCore/Framework/interface/MakerMacros.h"
0291 
0292 DEFINE_FWK_MODULE(DQMRivetClient);