File indexing completed on 2024-04-06 12:18:53
0001 #include "FWCore/Framework/interface/Frameworkfwd.h"
0002 #include "FWCore/Framework/interface/Event.h"
0003 #include "FWCore/Framework/interface/MakerMacros.h"
0004 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0005 #include "FWCore/ServiceRegistry/interface/Service.h"
0006 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0007 #include "FWCore/Utilities/interface/InputTag.h"
0008 #include "DQMServices/Core/interface/DQMEDHarvester.h"
0009 #include "DQMServices/Core/interface/DQMStore.h"
0010
0011 #include "TH1F.h"
0012 #include "TH2F.h"
0013 #include "TH3F.h"
0014 #include "RVersion.h"
0015
0016 #include "TEfficiency.h"
0017
0018 using namespace edm;
0019 using namespace std;
0020
0021 class HeavyFlavorHarvesting : public DQMEDHarvester {
0022 public:
0023 HeavyFlavorHarvesting(const edm::ParameterSet& pset);
0024 ~HeavyFlavorHarvesting() override;
0025
0026 private:
0027 void calculateEfficiency(const ParameterSet& pset, DQMStore::IBooker&, DQMStore::IGetter&);
0028 void calculateEfficiency1D(TH1* num, TH1* den, string name, DQMStore::IBooker&, DQMStore::IGetter&);
0029 void calculateEfficiency2D(TH2F* num, TH2F* den, string name, DQMStore::IBooker&, DQMStore::IGetter&);
0030
0031 string myDQMrootFolder;
0032 const VParameterSet efficiencies;
0033
0034 protected:
0035 void dqmEndJob(DQMStore::IBooker&, DQMStore::IGetter&) override;
0036 };
0037
0038 HeavyFlavorHarvesting::HeavyFlavorHarvesting(const edm::ParameterSet& pset)
0039 : myDQMrootFolder(pset.getUntrackedParameter<string>("MyDQMrootFolder")),
0040 efficiencies(pset.getUntrackedParameter<VParameterSet>("Efficiencies")) {}
0041
0042 void HeavyFlavorHarvesting::dqmEndJob(DQMStore::IBooker& ibooker_, DQMStore::IGetter& igetter_) {
0043 for (VParameterSet::const_iterator pset = efficiencies.begin(); pset != efficiencies.end(); pset++) {
0044 calculateEfficiency(*pset, ibooker_, igetter_);
0045 }
0046 }
0047
0048 void HeavyFlavorHarvesting::calculateEfficiency(const ParameterSet& pset,
0049 DQMStore::IBooker& ibooker_,
0050 DQMStore::IGetter& igetter_) {
0051
0052 vector<string> numDenEffMEnames = pset.getUntrackedParameter<vector<string> >("NumDenEffMEnames");
0053 if (numDenEffMEnames.size() != 3) {
0054 LogDebug("HLTriggerOfflineHeavyFlavor") << "NumDenEffMEnames must have three names" << endl;
0055 return;
0056 }
0057 string denMEname = myDQMrootFolder + "/" + numDenEffMEnames[1];
0058 string numMEname = myDQMrootFolder + "/" + numDenEffMEnames[0];
0059 MonitorElement* denME = igetter_.get(denMEname);
0060 MonitorElement* numME = igetter_.get(numMEname);
0061 if (denME == nullptr || numME == nullptr) {
0062 LogDebug("HLTriggerOfflineHeavyFlavor") << "Could not find MEs: " << denMEname << " or " << numMEname << endl;
0063 return;
0064 }
0065 TH1* den = denME->getTH1();
0066 TH1* num = numME->getTH1();
0067
0068 if (den->GetNbinsX() != num->GetNbinsX() || den->GetNbinsY() != num->GetNbinsY() ||
0069 den->GetNbinsZ() != num->GetNbinsZ()) {
0070 LogDebug("HLTriggerOfflineHeavyFlavor")
0071 << "Monitoring elements " << numMEname << " and " << denMEname << "are incompatible" << endl;
0072 return;
0073 }
0074
0075 string effName = numDenEffMEnames[2];
0076 string effDir = myDQMrootFolder;
0077 string::size_type slashPos = effName.rfind('/');
0078 if (string::npos != slashPos) {
0079 effDir += "/" + effName.substr(0, slashPos);
0080 effName.erase(0, slashPos + 1);
0081 }
0082 ibooker_.setCurrentFolder(effDir);
0083
0084 int dimensions = num->GetDimension();
0085 if (dimensions == 1) {
0086 calculateEfficiency1D(num, den, effName, ibooker_, igetter_);
0087 } else if (dimensions == 2) {
0088 calculateEfficiency2D((TH2F*)num, (TH2F*)den, effName, ibooker_, igetter_);
0089 TH1D* numX = ((TH2F*)num)->ProjectionX();
0090 TH1D* denX = ((TH2F*)den)->ProjectionX();
0091 calculateEfficiency1D(numX, denX, effName + "X", ibooker_, igetter_);
0092 delete numX;
0093 delete denX;
0094 TH1D* numY = ((TH2F*)num)->ProjectionY();
0095 TH1D* denY = ((TH2F*)den)->ProjectionY();
0096 calculateEfficiency1D(numY, denY, effName + "Y", ibooker_, igetter_);
0097 delete numY;
0098 delete denY;
0099 } else {
0100 return;
0101 }
0102 }
0103
0104 void HeavyFlavorHarvesting::calculateEfficiency1D(
0105 TH1* num, TH1* den, string effName, DQMStore::IBooker& ibooker_, DQMStore::IGetter& igetter_) {
0106 TProfile* eff;
0107 if (num->GetXaxis()->GetXbins()->GetSize() == 0) {
0108 eff = new TProfile(effName.c_str(),
0109 effName.c_str(),
0110 num->GetXaxis()->GetNbins(),
0111 num->GetXaxis()->GetXmin(),
0112 num->GetXaxis()->GetXmax());
0113 } else {
0114 eff = new TProfile(
0115 effName.c_str(), effName.c_str(), num->GetXaxis()->GetNbins(), num->GetXaxis()->GetXbins()->GetArray());
0116 }
0117 eff->SetTitle(effName.c_str());
0118 eff->SetXTitle(num->GetXaxis()->GetTitle());
0119 eff->SetYTitle("Efficiency");
0120 eff->SetOption("PE");
0121 eff->SetLineColor(2);
0122 eff->SetLineWidth(2);
0123 eff->SetMarkerStyle(20);
0124 eff->SetMarkerSize(0.8);
0125 eff->GetYaxis()->SetRangeUser(-0.001, 1.001);
0126 eff->SetStats(kFALSE);
0127 for (int i = 1; i <= num->GetNbinsX(); i++) {
0128 double e, low, high;
0129 if (int(den->GetBinContent(i)) > 0.)
0130 e = double(num->GetBinContent(i)) / double(den->GetBinContent(i));
0131 else
0132 e = 0.;
0133 low = TEfficiency::Wilson((double)den->GetBinContent(i), (double)num->GetBinContent(i), 0.683, false);
0134 high = TEfficiency::Wilson((double)den->GetBinContent(i), (double)num->GetBinContent(i), 0.683, true);
0135
0136 double err = e - low > high - e ? e - low : high - e;
0137
0138 eff->SetBinContent(i, e);
0139 eff->SetBinEntries(i, 1);
0140 eff->SetBinError(i, sqrt(e * e + err * err));
0141 }
0142 ibooker_.bookProfile(effName, eff);
0143 delete eff;
0144 }
0145
0146 void HeavyFlavorHarvesting::calculateEfficiency2D(
0147 TH2F* num, TH2F* den, string effName, DQMStore::IBooker& ibooker_, DQMStore::IGetter& igetter_) {
0148 TProfile2D* eff;
0149 if (num->GetXaxis()->GetXbins()->GetSize() == 0 && num->GetYaxis()->GetXbins()->GetSize() == 0) {
0150 eff = new TProfile2D(effName.c_str(),
0151 effName.c_str(),
0152 num->GetXaxis()->GetNbins(),
0153 num->GetXaxis()->GetXmin(),
0154 num->GetXaxis()->GetXmax(),
0155 num->GetYaxis()->GetNbins(),
0156 num->GetYaxis()->GetXmin(),
0157 num->GetYaxis()->GetXmax());
0158 } else if (num->GetXaxis()->GetXbins()->GetSize() != 0 && num->GetYaxis()->GetXbins()->GetSize() == 0) {
0159 eff = new TProfile2D(effName.c_str(),
0160 effName.c_str(),
0161 num->GetXaxis()->GetNbins(),
0162 num->GetXaxis()->GetXbins()->GetArray(),
0163 num->GetYaxis()->GetNbins(),
0164 num->GetYaxis()->GetXmin(),
0165 num->GetYaxis()->GetXmax());
0166 } else if (num->GetXaxis()->GetXbins()->GetSize() == 0 && num->GetYaxis()->GetXbins()->GetSize() != 0) {
0167 eff = new TProfile2D(effName.c_str(),
0168 effName.c_str(),
0169 num->GetXaxis()->GetNbins(),
0170 num->GetXaxis()->GetXmin(),
0171 num->GetXaxis()->GetXmax(),
0172 num->GetYaxis()->GetNbins(),
0173 num->GetYaxis()->GetXbins()->GetArray());
0174 } else {
0175 eff = new TProfile2D(effName.c_str(),
0176 effName.c_str(),
0177 num->GetXaxis()->GetNbins(),
0178 num->GetXaxis()->GetXbins()->GetArray(),
0179 num->GetYaxis()->GetNbins(),
0180 num->GetYaxis()->GetXbins()->GetArray());
0181 }
0182 eff->SetTitle(effName.c_str());
0183 eff->SetXTitle(num->GetXaxis()->GetTitle());
0184 eff->SetYTitle(num->GetYaxis()->GetTitle());
0185 eff->SetZTitle("Efficiency");
0186 eff->SetOption("colztexte");
0187 eff->GetZaxis()->SetRangeUser(-0.001, 1.001);
0188 eff->SetStats(kFALSE);
0189 for (int i = 0; i < num->GetSize(); i++) {
0190 double e, low, high;
0191 if (int(den->GetBinContent(i)) > 0.)
0192 e = double(num->GetBinContent(i)) / double(den->GetBinContent(i));
0193 else
0194 e = 0.;
0195 low = TEfficiency::Wilson((double)den->GetBinContent(i), (double)num->GetBinContent(i), 0.683, false);
0196 high = TEfficiency::Wilson((double)den->GetBinContent(i), (double)num->GetBinContent(i), 0.683, true);
0197
0198 double err = e - low > high - e ? e - low : high - e;
0199
0200 eff->SetBinContent(i, e);
0201 eff->SetBinEntries(i, 1);
0202 eff->SetBinError(i, sqrt(e * e + err * err));
0203 }
0204 ibooker_.bookProfile2D(effName, eff);
0205 delete eff;
0206 }
0207
0208 HeavyFlavorHarvesting::~HeavyFlavorHarvesting() {}
0209
0210
0211 DEFINE_FWK_MODULE(HeavyFlavorHarvesting);