File indexing completed on 2024-04-06 12:33:14
0001
0002
0003
0004
0005
0006
0007
0008
0009 #include "FWCore/Framework/interface/LuminosityBlock.h"
0010 #include "FWCore/Framework/interface/Run.h"
0011 #include "FWCore/Framework/interface/Event.h"
0012 #include "FWCore/Framework/interface/MakerMacros.h"
0013 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0014 #include "FWCore/ServiceRegistry/interface/Service.h"
0015
0016 #include "DQMServices/Core/interface/DQMEDHarvester.h"
0017 #include "DQMServices/Core/interface/DQMStore.h"
0018
0019
0020
0021
0022
0023 class OffsetDQMPostProcessor : public DQMEDHarvester {
0024 public:
0025 explicit OffsetDQMPostProcessor(const edm::ParameterSet&);
0026 ~OffsetDQMPostProcessor() override;
0027
0028 private:
0029 void dqmEndJob(DQMStore::IBooker&, DQMStore::IGetter&) override;
0030
0031 std::string offsetPlotBaseName;
0032 std::string offsetDir;
0033 float offsetR;
0034 std::vector<std::string> pftypes;
0035 std::vector<std::string> offsetVariableTypes;
0036 int muHigh;
0037 int npvHigh;
0038
0039 bool debug = false;
0040 };
0041
0042
0043
0044
0045 OffsetDQMPostProcessor::OffsetDQMPostProcessor(const edm::ParameterSet& iConfig) {
0046 offsetPlotBaseName = iConfig.getParameter<std::string>("offsetPlotBaseName");
0047 offsetDir = iConfig.getParameter<std::string>("offsetDir");
0048 offsetVariableTypes = iConfig.getParameter<std::vector<std::string> >("offsetVariableTypes");
0049 offsetR = iConfig.getUntrackedParameter<double>("offsetR");
0050 pftypes = iConfig.getParameter<std::vector<std::string> >("pftypes");
0051 muHigh = iConfig.getUntrackedParameter<int>("muHigh");
0052 npvHigh = iConfig.getUntrackedParameter<int>("npvHigh");
0053 };
0054
0055 OffsetDQMPostProcessor::~OffsetDQMPostProcessor() {}
0056
0057
0058 void OffsetDQMPostProcessor::dqmEndJob(DQMStore::IBooker& ibook_, DQMStore::IGetter& iget_) {
0059 iget_.setCurrentFolder(offsetDir);
0060
0061 std::string stitle;
0062 std::vector<MonitorElement*> vME;
0063 std::vector<std::string> MEStrings = iget_.getMEs();
0064 std::for_each(MEStrings.begin(), MEStrings.end(), [&](auto& s) { s.insert(0, offsetDir); });
0065
0066
0067 MonitorElement* mtmp;
0068 TProfile* hproftmp;
0069 TH1F* htmp;
0070 TH1F* hscaled;
0071
0072
0073
0074
0075 for (std::vector<std::string>::const_iterator i = offsetVariableTypes.begin(); i != offsetVariableTypes.end(); ++i) {
0076
0077
0078
0079 stitle = offsetDir + (*i);
0080 std::vector<std::string>::const_iterator it = std::find(MEStrings.begin(), MEStrings.end(), stitle);
0081 if (it == MEStrings.end())
0082 continue;
0083 mtmp = iget_.get(stitle);
0084 float avg = mtmp->getMean();
0085 int iavg = int(avg + 0.5);
0086
0087 if (avg < 1.)
0088 avg = 1.;
0089
0090 if (iavg < 0)
0091 iavg = 0;
0092 if (*i == "npv" && iavg >= npvHigh)
0093 iavg = npvHigh - 1;
0094 else if (*i == "mu" && iavg >= muHigh)
0095 iavg = muHigh - 1;
0096
0097
0098
0099
0100 stitle = (*i) + "_mean";
0101 MonitorElement* MEmean = ibook_.bookFloat(stitle);
0102 MEmean->Fill(avg);
0103 vME.push_back(MEmean);
0104
0105
0106
0107
0108 for (std::vector<std::string>::const_iterator j = pftypes.begin(); j != pftypes.end(); ++j) {
0109
0110 std::string str_base = *i + std::to_string(iavg);
0111 if ((*i) == "npv")
0112 stitle = offsetDir + "npvPlots/" + str_base + "/" + offsetPlotBaseName + "_" + str_base + "_" + (*j);
0113 else if ((*i) == "mu")
0114 stitle = offsetDir + "muPlots/" + str_base + "/" + offsetPlotBaseName + "_" + str_base + "_" + (*j);
0115 else
0116 return;
0117
0118
0119 mtmp = iget_.get(stitle);
0120 hproftmp = (TProfile*)mtmp->getTProfile();
0121 htmp = (TH1F*)hproftmp->ProjectionX();
0122 TAxis* xaxis = (TAxis*)htmp->GetXaxis();
0123 stitle = offsetPlotBaseName + "_" + str_base + "_" + *j;
0124 hscaled = new TH1F(stitle.c_str(), stitle.c_str(), xaxis->GetNbins(), xaxis->GetXbins()->GetArray());
0125
0126 htmp->Scale(pow(offsetR, 2) / 2. / float(avg));
0127 for (int ibin = 1; ibin <= hscaled->GetNbinsX(); ibin++) {
0128 hscaled->SetBinContent(ibin, htmp->GetBinContent(ibin) / htmp->GetBinWidth(ibin));
0129 hscaled->SetBinError(ibin, htmp->GetBinError(ibin) / htmp->GetBinWidth(ibin));
0130 }
0131
0132
0133 stitle = offsetPlotBaseName + "_" + *i + "_" + *j;
0134 mtmp = ibook_.book1D(stitle.c_str(), hscaled);
0135 vME.push_back(mtmp);
0136 }
0137 }
0138
0139
0140
0141
0142 if (debug) {
0143 for (std::vector<MonitorElement*>::const_iterator i = vME.begin(); i != vME.end(); ++i)
0144 (*i)->getTH1F()->Print();
0145 }
0146 }
0147
0148 #include "FWCore/Framework/interface/MakerMacros.h"
0149 DEFINE_FWK_MODULE(OffsetDQMPostProcessor);