File indexing completed on 2024-09-11 04:32:52
0001 #include "FWCore/Framework/interface/Frameworkfwd.h"
0002 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0003 #include "FWCore/Framework/interface/Event.h"
0004 #include "FWCore/Framework/interface/MakerMacros.h"
0005 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0006 #include "FWCore/ServiceRegistry/interface/Service.h"
0007 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0008 #include "FWCore/Utilities/interface/InputTag.h"
0009 #include "DQMServices/Core/interface/DQMStore.h"
0010
0011 #include "DQMOffline/Trigger/plugins/GenericTnPFitter.h"
0012
0013 #include <TString.h>
0014 #include <TPRegexp.h>
0015
0016 using namespace edm;
0017 using namespace dqmTnP;
0018 using namespace std;
0019
0020 using vstring = std::vector<std::string>;
0021
0022 class DQMGenericTnPClient : public edm::one::EDAnalyzer<edm::one::SharedResources, edm::one::WatchRuns> {
0023 public:
0024 typedef dqm::legacy::MonitorElement MonitorElement;
0025 typedef dqm::legacy::DQMStore DQMStore;
0026
0027 DQMGenericTnPClient(const edm::ParameterSet& pset);
0028 ~DQMGenericTnPClient() override;
0029 void analyze(const edm::Event& event, const edm::EventSetup& eventSetup) override {}
0030 void beginRun(edm::Run const&, edm::EventSetup const&) override {}
0031 void endRun(const edm::Run& run, const edm::EventSetup& setup) override;
0032 void calculateEfficiency(const std::string& dirName, const ParameterSet& pset);
0033 void findAllSubdirectories(const std::string& dir, std::set<std::string>* myList, TString pattern);
0034
0035 private:
0036 DQMStore* dqmStore;
0037 TFile* plots;
0038 vstring subDirs;
0039 std::string myDQMrootFolder;
0040 bool verbose;
0041 const VParameterSet efficiencies;
0042 GaussianPlusLinearFitter* GPLfitter;
0043 VoigtianPlusExponentialFitter* VPEfitter;
0044 };
0045
0046 DQMGenericTnPClient::DQMGenericTnPClient(const edm::ParameterSet& pset)
0047 : subDirs(pset.getUntrackedParameter<vstring>("subDirs", vstring())),
0048 myDQMrootFolder(pset.getUntrackedParameter<std::string>("MyDQMrootFolder", "")),
0049 verbose(pset.getUntrackedParameter<bool>("Verbose", false)),
0050 efficiencies(pset.getUntrackedParameter<VParameterSet>("Efficiencies")) {
0051 usesResource("DQMStore");
0052 TString savePlotsInRootFileName = pset.getUntrackedParameter<string>("SavePlotsInRootFileName", "");
0053 plots = savePlotsInRootFileName != "" ? new TFile(savePlotsInRootFileName, "recreate") : nullptr;
0054 GPLfitter = new GaussianPlusLinearFitter(verbose);
0055 VPEfitter = new VoigtianPlusExponentialFitter(verbose);
0056 }
0057
0058 void DQMGenericTnPClient::endRun(const edm::Run& run, const edm::EventSetup& setup) {
0059 TPRegexp metacharacters("[\\^\\$\\.\\*\\+\\?\\|\\(\\)\\{\\}\\[\\]]");
0060
0061 dqmStore = Service<DQMStore>().operator->();
0062 if (!dqmStore) {
0063 LogError("DQMGenericTnPClient") << "Could not find DQMStore service\n";
0064 return;
0065 }
0066 dqmStore->setCurrentFolder(myDQMrootFolder);
0067
0068 set<std::string> subDirSet;
0069
0070 if (!myDQMrootFolder.empty())
0071 subDirSet.insert(myDQMrootFolder);
0072 else {
0073 for (auto iSubDir = subDirs.begin(); iSubDir != subDirs.end(); ++iSubDir) {
0074 std::string subDir = *iSubDir;
0075 if (subDir[subDir.size() - 1] == '/')
0076 subDir.erase(subDir.size() - 1);
0077 if (TString(subDir).Contains(metacharacters)) {
0078 const string::size_type shiftPos = subDir.rfind('/');
0079 const string searchPath = subDir.substr(0, shiftPos);
0080 const string pattern = subDir.substr(shiftPos + 1, subDir.length());
0081 findAllSubdirectories(searchPath, &subDirSet, pattern);
0082 } else {
0083 subDirSet.insert(subDir);
0084 }
0085 }
0086 }
0087
0088 for (auto const& dirName : subDirSet) {
0089 for (auto const& efficiencie : efficiencies) {
0090 calculateEfficiency(dirName, efficiencie);
0091 }
0092 }
0093 }
0094
0095 void DQMGenericTnPClient::calculateEfficiency(const std::string& dirName, const ParameterSet& pset) {
0096
0097 string allMEname = dirName + "/" + pset.getUntrackedParameter<string>("DenominatorMEname");
0098 string passMEname = dirName + "/" + pset.getUntrackedParameter<string>("NumeratorMEname");
0099 MonitorElement* allME = dqmStore->get(allMEname);
0100 MonitorElement* passME = dqmStore->get(passMEname);
0101 if (allME == nullptr || passME == nullptr) {
0102 LogDebug("DQMGenericTnPClient") << "Could not find MEs: " << allMEname << " or " << passMEname << endl;
0103 return;
0104 }
0105 TH1* all = allME->getTH1();
0106 TH1* pass = passME->getTH1();
0107
0108 string fitFunction = pset.getUntrackedParameter<string>("FitFunction");
0109 AbstractFitter* fitter = nullptr;
0110 if (fitFunction == "GaussianPlusLinear") {
0111 GPLfitter->setup(pset.getUntrackedParameter<double>("ExpectedMean"),
0112 pset.getUntrackedParameter<double>("FitRangeLow"),
0113 pset.getUntrackedParameter<double>("FitRangeHigh"),
0114 pset.getUntrackedParameter<double>("ExpectedSigma"));
0115 fitter = GPLfitter;
0116 } else if (fitFunction == "VoigtianPlusExponential") {
0117 VPEfitter->setup(pset.getUntrackedParameter<double>("ExpectedMean"),
0118 pset.getUntrackedParameter<double>("FitRangeLow"),
0119 pset.getUntrackedParameter<double>("FitRangeHigh"),
0120 pset.getUntrackedParameter<double>("ExpectedSigma"),
0121 pset.getUntrackedParameter<double>("FixedWidth"));
0122 fitter = VPEfitter;
0123 } else {
0124 LogError("DQMGenericTnPClient") << "Fit function: " << fitFunction << " does not exist" << endl;
0125 return;
0126 }
0127
0128 int dimensions = all->GetDimension();
0129 int massDimension = pset.getUntrackedParameter<int>("MassDimension");
0130 if (massDimension > dimensions) {
0131 LogError("DQMGenericTnPClient") << "Monitoring Element " << allMEname << " has smaller dimension than "
0132 << massDimension << endl;
0133 return;
0134 }
0135
0136 string effName = pset.getUntrackedParameter<string>("EfficiencyMEname");
0137 string effDir = dirName;
0138 string::size_type slashPos = effName.rfind('/');
0139 if (string::npos != slashPos) {
0140 effDir += "/" + effName.substr(0, slashPos);
0141 effName.erase(0, slashPos + 1);
0142 }
0143 dqmStore->setCurrentFolder(effDir);
0144 TString prefix(effDir.c_str());
0145 prefix.ReplaceAll('/', '_');
0146
0147 if (dimensions == 2) {
0148 TProfile* eff = nullptr;
0149 TProfile* effChi2 = nullptr;
0150 TString error = fitter->calculateEfficiency(
0151 (TH2*)pass, (TH2*)all, massDimension, eff, effChi2, plots ? prefix + effName.c_str() : "");
0152 if (error != "") {
0153 LogError("DQMGenericTnPClient") << error << endl;
0154 return;
0155 }
0156 dqmStore->bookProfile(effName, eff);
0157 dqmStore->bookProfile(effName + "Chi2", effChi2);
0158 delete eff;
0159 delete effChi2;
0160 } else if (dimensions == 3) {
0161 TProfile2D* eff = nullptr;
0162 TProfile2D* effChi2 = nullptr;
0163 TString error = fitter->calculateEfficiency(
0164 (TH3*)pass, (TH3*)all, massDimension, eff, effChi2, plots ? prefix + effName.c_str() : "");
0165 if (error != "") {
0166 LogError("DQMGenericTnPClient") << error << endl;
0167 return;
0168 }
0169 dqmStore->bookProfile2D(effName, eff);
0170 dqmStore->bookProfile2D(effName + "Chi2", effChi2);
0171 delete eff;
0172 delete effChi2;
0173 }
0174 }
0175
0176 DQMGenericTnPClient::~DQMGenericTnPClient() {
0177 delete GPLfitter;
0178 if (plots) {
0179 plots->Close();
0180 }
0181 }
0182
0183 void DQMGenericTnPClient::findAllSubdirectories(const std::string& dir,
0184 std::set<std::string>* myList,
0185 TString pattern = "") {
0186 if (!dqmStore->dirExists(dir)) {
0187 LogError("DQMGenericTnPClient") << " DQMGenericTnPClient::findAllSubdirectories ==> Missing folder " << dir
0188 << " !!!";
0189 return;
0190 }
0191 TPRegexp nonPerlWildcard("\\w\\*|^\\*");
0192 if (pattern != "") {
0193 if (pattern.Contains(nonPerlWildcard))
0194 pattern.ReplaceAll("*", ".*");
0195 TPRegexp regexp(pattern);
0196 dqmStore->cd(dir);
0197 vector<string> foundDirs = dqmStore->getSubdirs();
0198 for (auto iDir = foundDirs.begin(); iDir != foundDirs.end(); ++iDir) {
0199 TString dirName = iDir->substr(iDir->rfind('/') + 1, iDir->length());
0200 if (dirName.Contains(regexp))
0201 findAllSubdirectories(*iDir, myList);
0202 }
0203 } else if (dqmStore->dirExists(dir)) {
0204 myList->insert(dir);
0205 dqmStore->cd(dir);
0206 findAllSubdirectories(dir, myList, "*");
0207 } else {
0208 LogInfo("DQMGenericClient") << "Trying to find sub-directories of " << dir << " failed because " << dir
0209 << " does not exist";
0210 }
0211 return;
0212 }
0213
0214
0215 DEFINE_FWK_MODULE(DQMGenericTnPClient);