File indexing completed on 2024-09-07 04:35:31
0001 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0002
0003 #include "CondCore/Utilities/interface/PayloadInspectorModule.h"
0004 #include "CondCore/Utilities/interface/PayloadInspector.h"
0005 #include "CondCore/CondDB/interface/Time.h"
0006
0007
0008 #include "CondFormats/JetMETObjects/interface/JetResolution.h"
0009 #include "CondFormats/JetMETObjects/interface/JetResolutionObject.h"
0010 #include "CondFormats/JetMETObjects/interface/Utilities.h"
0011
0012 #include <memory>
0013 #include <sstream>
0014 #include <fstream>
0015 #include <iostream>
0016
0017
0018 #include "TH2F.h"
0019 #include "TLegend.h"
0020 #include "TCanvas.h"
0021 #include "TLine.h"
0022 #include "TStyle.h"
0023 #include "TLatex.h"
0024 #include "TPave.h"
0025 #include "TPaveStats.h"
0026
0027 #define MIN_ETA -5.05
0028 #define MAX_ETA 5.05
0029 #define NBIN_ETA 51
0030 #define MIN_PT -5.
0031 #define MAX_PT 3005.
0032 #define NBIN_PT 301
0033
0034 namespace JME {
0035
0036 using namespace cond::payloadInspector;
0037
0038 enum index { NORM = 0, DOWN = 1, UP = 2 };
0039
0040
0041
0042
0043
0044
0045
0046
0047
0048 class JetResolutionVsEta : public cond::payloadInspector::Histogram1D<JetResolutionObject, SINGLE_IOV> {
0049 public:
0050 JetResolutionVsEta()
0051 : cond::payloadInspector::Histogram1D<JetResolutionObject, SINGLE_IOV>(
0052 "Jet Resolution", "#eta", NBIN_ETA, MIN_ETA, MAX_ETA, "Resolution") {
0053 cond::payloadInspector::PlotBase::addInputParam("Jet_Pt");
0054 cond::payloadInspector::PlotBase::addInputParam("Jet_Rho");
0055 }
0056
0057 bool fill() override {
0058 double par_Pt = 100.;
0059 double par_Eta = 1.;
0060 double par_Rho = 20.;
0061
0062
0063 auto paramValues = cond::payloadInspector::PlotBase::inputParamValues();
0064 auto ip = paramValues.find("Jet_Pt");
0065 if (ip != paramValues.end()) {
0066 par_Pt = std::stod(ip->second);
0067 edm::LogPrint("JER_PI") << "Jet Pt: " << par_Pt;
0068 }
0069 ip = paramValues.find("Jet_Rho");
0070 if (ip != paramValues.end()) {
0071 par_Rho = std::stod(ip->second);
0072 edm::LogPrint("JER_PI") << "Jet Rho: " << par_Rho;
0073 }
0074
0075 auto tag = PlotBase::getTag<0>();
0076 for (auto const& iov : tag.iovs) {
0077 std::shared_ptr<JetResolutionObject> payload = Base::fetchPayload(std::get<1>(iov));
0078 if (payload.get()) {
0079 if (!payload->getRecords().empty() &&
0080 payload->getDefinition().getFormulaString().compare("") == 0)
0081 return false;
0082
0083 for (size_t idx = 0; idx <= NBIN_ETA; idx++) {
0084 par_Eta = (idx + 0.5) * (MAX_ETA - MIN_ETA) / NBIN_ETA + MIN_ETA;
0085
0086 JetParameters j_param;
0087 j_param.setJetPt(par_Pt);
0088 j_param.setRho(par_Rho);
0089 j_param.setJetEta(par_Eta);
0090
0091 if (payload->getRecord(j_param) == nullptr) {
0092 continue;
0093 }
0094
0095 const JetResolutionObject::Record record = *(payload->getRecord(j_param));
0096 float res = payload->evaluateFormula(record, j_param);
0097 fillWithBinAndValue(idx, res);
0098 }
0099 } else
0100 return false;
0101 }
0102 return true;
0103 }
0104 };
0105
0106 class JetResolutionVsPt : public cond::payloadInspector::Histogram1D<JetResolutionObject, SINGLE_IOV> {
0107 public:
0108 JetResolutionVsPt()
0109 : cond::payloadInspector::Histogram1D<JetResolutionObject, SINGLE_IOV>(
0110 "Jet Energy Resolution", "p_T", NBIN_PT, MIN_PT, MAX_PT, "Resolution") {
0111 cond::payloadInspector::PlotBase::addInputParam("Jet_Eta");
0112 cond::payloadInspector::PlotBase::addInputParam("Jet_Rho");
0113 }
0114
0115 bool fill() override {
0116 double par_Pt = 100.;
0117 double par_Eta = 1.;
0118 double par_Rho = 20.;
0119
0120 auto paramValues = cond::payloadInspector::PlotBase::inputParamValues();
0121 auto ip = paramValues.find("Jet_Eta");
0122 if (ip != paramValues.end()) {
0123 par_Eta = std::stod(ip->second);
0124 edm::LogPrint("JER_PI") << "Jet Eta: " << par_Eta;
0125 }
0126 ip = paramValues.find("Jet_Rho");
0127 if (ip != paramValues.end()) {
0128 par_Rho = std::stod(ip->second);
0129 edm::LogPrint("JER_PI") << "Jet Rho: " << par_Rho;
0130 }
0131
0132 auto tag = PlotBase::getTag<0>();
0133 for (auto const& iov : tag.iovs) {
0134 std::shared_ptr<JetResolutionObject> payload = Base::fetchPayload(std::get<1>(iov));
0135 if (payload.get()) {
0136 if (!payload->getRecords().empty() &&
0137 payload->getDefinition().getFormulaString().compare("") == 0)
0138 return false;
0139
0140 for (size_t idx = 0; idx <= NBIN_PT; idx++) {
0141 par_Pt = (idx + 0.5) * (MAX_PT - MIN_PT) / NBIN_PT + MIN_PT;
0142
0143 JetParameters j_param;
0144 j_param.setJetEta(par_Eta);
0145 j_param.setRho(par_Rho);
0146 j_param.setJetPt(par_Pt);
0147
0148 if (payload->getRecord(j_param) == nullptr) {
0149 continue;
0150 }
0151
0152 const JetResolutionObject::Record record = *(payload->getRecord(j_param));
0153 float res = payload->evaluateFormula(record, j_param);
0154 fillWithBinAndValue(idx + 1, res);
0155 }
0156 } else
0157 return false;
0158 }
0159 return true;
0160 }
0161 };
0162
0163 class JetResolutionSummary : public cond::payloadInspector::PlotImage<JetResolutionObject, SINGLE_IOV> {
0164 public:
0165 JetResolutionSummary()
0166 : cond::payloadInspector::PlotImage<JetResolutionObject, SINGLE_IOV>("Jet Resolution Summary") {
0167 cond::payloadInspector::PlotBase::addInputParam("Jet_Pt");
0168 cond::payloadInspector::PlotBase::addInputParam("Jet_Eta");
0169 cond::payloadInspector::PlotBase::addInputParam("Jet_Rho");
0170 }
0171
0172 bool fill() override {
0173 double par_Pt = 100.;
0174 double par_Eta = 1.;
0175 double par_Rho = 20.;
0176
0177 auto paramValues = cond::payloadInspector::PlotBase::inputParamValues();
0178 auto ip = paramValues.find("Jet_Pt");
0179 if (ip != paramValues.end()) {
0180 par_Pt = std::stod(ip->second);
0181 }
0182 ip = paramValues.find("Jet_Eta");
0183 if (ip != paramValues.end()) {
0184 par_Eta = std::stod(ip->second);
0185 }
0186 ip = paramValues.find("Jet_Rho");
0187 if (ip != paramValues.end()) {
0188 par_Rho = std::stod(ip->second);
0189 }
0190
0191 TH1D* resol_eta = new TH1D("Jet Resolution vs #eta", "", NBIN_ETA, MIN_ETA, MAX_ETA);
0192 TH1D* resol_pt = new TH1D("Jet Resolution vs p_T", "", NBIN_PT, MIN_PT, MAX_PT);
0193 TLegend* leg_eta = new TLegend(0.26, 0.73, 0.935, 0.90);
0194 TLegend* leg_pt = new TLegend(0.26, 0.73, 0.935, 0.90);
0195
0196 leg_eta->SetBorderSize(0);
0197 leg_eta->SetLineStyle(0);
0198 leg_eta->SetFillStyle(0);
0199
0200 leg_eta->SetTextFont(42);
0201 leg_pt->SetBorderSize(0);
0202 leg_pt->SetLineStyle(0);
0203 leg_pt->SetFillStyle(0);
0204
0205 auto tag = PlotBase::getTag<0>();
0206 auto iov = tag.iovs.front();
0207 std::shared_ptr<JetResolutionObject> payload = fetchPayload(std::get<1>(iov));
0208 unsigned int run = std::get<0>(iov);
0209 std::string tagname = tag.name;
0210 std::stringstream ss_tagname(tag.name);
0211 std::string stmp;
0212
0213 std::string tag_ver;
0214 std::string tag_res;
0215 std::string tag_jet;
0216
0217 getline(ss_tagname, stmp, '_');
0218 getline(ss_tagname, stmp, '_');
0219 tag_ver = stmp;
0220 getline(ss_tagname, stmp, '_');
0221 tag_ver += '_' + stmp;
0222 getline(ss_tagname, stmp, '_');
0223 tag_ver += '_' + stmp;
0224 getline(ss_tagname, stmp, '_');
0225 tag_ver += '_' + stmp;
0226 getline(ss_tagname, stmp, '_');
0227 tag_res = stmp;
0228 getline(ss_tagname, stmp, '_');
0229 tag_jet = stmp;
0230
0231 if (payload.get()) {
0232 if (!payload->getRecords().empty() &&
0233 payload->getDefinition().getFormulaString().compare("") == 0)
0234 return false;
0235
0236 for (size_t idx = 0; idx <= NBIN_ETA; idx++) {
0237 double x_axis = (idx + 0.5) * (MAX_ETA - MIN_ETA) / NBIN_ETA + MIN_ETA;
0238
0239 JetParameters j_param;
0240 j_param.setJetPt(par_Pt);
0241 j_param.setRho(par_Rho);
0242 j_param.setJetEta(x_axis);
0243
0244 if (payload->getRecord(j_param) == nullptr) {
0245 continue;
0246 }
0247
0248 const JetResolutionObject::Record record = *(payload->getRecord(j_param));
0249 float res = payload->evaluateFormula(record, j_param);
0250 resol_eta->SetBinContent(idx + 1, res);
0251 }
0252
0253 for (size_t idx = 0; idx <= NBIN_PT; idx++) {
0254 double x_axis = (idx + 0.5) * (MAX_PT - MIN_PT) / NBIN_PT + MIN_PT;
0255
0256 JetParameters j_param;
0257 j_param.setJetEta(par_Eta);
0258 j_param.setRho(par_Rho);
0259 j_param.setJetPt(x_axis);
0260
0261 if (payload->getRecord(j_param) == nullptr) {
0262 continue;
0263 }
0264
0265 const JetResolutionObject::Record record = *(payload->getRecord(j_param));
0266 float res = payload->evaluateFormula(record, j_param);
0267 resol_pt->SetBinContent(idx + 1, res);
0268 }
0269
0270 gStyle->SetOptStat(0);
0271 gStyle->SetLabelFont(42, "XYZ");
0272 gStyle->SetLabelSize(0.05, "XYZ");
0273 gStyle->SetFrameLineWidth(3);
0274
0275 std::string title = Form("Summary Run %i", run);
0276 TCanvas canvas("Jet Resolution Summary", title.c_str(), 800, 1200);
0277 canvas.Divide(1, 2);
0278
0279 canvas.cd(1);
0280 resol_eta->SetTitle("Jet Resolution vs. #eta");
0281 resol_eta->SetXTitle("#eta");
0282 resol_eta->SetYTitle("Resolution");
0283 resol_eta->SetLineWidth(3);
0284 resol_eta->SetMaximum(resol_eta->GetMaximum() * 1.25);
0285 resol_eta->Draw("");
0286
0287 leg_eta->AddEntry(resol_eta, (tag_ver + '_' + tag_jet).c_str(), "l");
0288 leg_eta->AddEntry((TObject*)nullptr, Form("JetPt=%.2f; JetRho=%.2f", par_Pt, par_Rho), "");
0289 leg_eta->Draw();
0290
0291 canvas.cd(2);
0292 resol_pt->SetTitle("Jet Resolution vs. p_{T}");
0293 resol_pt->SetXTitle("p_{T} [GeV]");
0294 resol_pt->SetYTitle("Resolution");
0295 resol_pt->SetLineWidth(3);
0296 resol_pt->Draw("][");
0297
0298 leg_pt->AddEntry(resol_pt, (tag_ver + '_' + tag_jet).c_str(), "l");
0299 leg_pt->AddEntry((TObject*)nullptr, Form("JetEta=%.2f; JetRho=%.2f", par_Eta, par_Rho), "");
0300 leg_pt->Draw();
0301
0302 canvas.SaveAs(m_imageFileName.c_str());
0303
0304 return true;
0305 } else
0306 return false;
0307 }
0308
0309 };
0310
0311 class JetResolutionCompare : public cond::payloadInspector::PlotImage<JetResolutionObject, SINGLE_IOV, 2> {
0312 public:
0313 JetResolutionCompare()
0314 : cond::payloadInspector::PlotImage<JetResolutionObject, SINGLE_IOV, 2>("Jet Resolution Compare") {
0315 cond::payloadInspector::PlotBase::addInputParam("Jet_Pt");
0316 cond::payloadInspector::PlotBase::addInputParam("Jet_Eta");
0317 cond::payloadInspector::PlotBase::addInputParam("Jet_Rho");
0318 }
0319
0320 bool fill() override {
0321 double par_Pt = 100.;
0322 double par_Eta = 1.;
0323 double par_Rho = 20.;
0324
0325 auto paramValues = cond::payloadInspector::PlotBase::inputParamValues();
0326 auto ip = paramValues.find("Jet_Pt");
0327 if (ip != paramValues.end()) {
0328 par_Pt = std::stod(ip->second);
0329 }
0330 ip = paramValues.find("Jet_Eta");
0331 if (ip != paramValues.end()) {
0332 par_Eta = std::stod(ip->second);
0333 }
0334 ip = paramValues.find("Jet_Rho");
0335 if (ip != paramValues.end()) {
0336 par_Rho = std::stod(ip->second);
0337 }
0338
0339 TH1D* resol_eta1 = new TH1D("Jet Resolution vs #eta one", "", NBIN_ETA, MIN_ETA, MAX_ETA);
0340 TH1D* resol_eta2 = new TH1D("Jet Resolution vs #eta two", "", NBIN_ETA, MIN_ETA, MAX_ETA);
0341 TH1D* resol_pt1 = new TH1D("Jet Resolution vs p_T one", "", NBIN_PT, MIN_PT, MAX_PT);
0342 TH1D* resol_pt2 = new TH1D("Jet Resolution vs p_T two", "", NBIN_PT, MIN_PT, MAX_PT);
0343 TLegend* leg_eta = new TLegend(0.26, 0.73, 0.935, 0.90);
0344 TLegend* leg_pt = new TLegend(0.26, 0.73, 0.935, 0.90);
0345
0346 leg_eta->SetBorderSize(0);
0347 leg_eta->SetLineStyle(0);
0348 leg_eta->SetFillStyle(0);
0349
0350 leg_eta->SetTextFont(42);
0351 leg_pt->SetBorderSize(0);
0352 leg_pt->SetLineStyle(0);
0353 leg_pt->SetFillStyle(0);
0354
0355 auto tag1 = PlotBase::getTag<0>();
0356 auto iov1 = tag1.iovs.front();
0357 std::shared_ptr<JetResolutionObject> payload1 = fetchPayload(std::get<1>(iov1));
0358 auto tag2 = PlotBase::getTag<1>();
0359 auto iov2 = tag2.iovs.front();
0360 std::shared_ptr<JetResolutionObject> payload2 = fetchPayload(std::get<1>(iov2));
0361
0362 std::stringstream ss_tagname1(tag1.name);
0363 std::stringstream ss_tagname2(tag2.name);
0364 std::string stmp;
0365
0366 std::string tag_ver1;
0367 std::string tag_ver2;
0368
0369 getline(ss_tagname1, stmp, '_');
0370 getline(ss_tagname1, stmp);
0371 tag_ver1 = stmp;
0372
0373 getline(ss_tagname2, stmp, '_');
0374 getline(ss_tagname2, stmp);
0375 tag_ver2 = stmp;
0376
0377 if (payload1.get() && payload2.get()) {
0378 if (!payload1->getRecords().empty() &&
0379 payload1->getDefinition().getFormulaString().compare("") == 0)
0380 return false;
0381
0382 for (size_t idx = 0; idx <= NBIN_ETA; idx++) {
0383 double x_axis = (idx + 0.5) * (MAX_ETA - MIN_ETA) / NBIN_ETA + MIN_ETA;
0384
0385 JetParameters j_param;
0386 j_param.setJetPt(par_Pt);
0387 j_param.setRho(par_Rho);
0388 j_param.setJetEta(x_axis);
0389
0390 if (payload1->getRecord(j_param) == nullptr || payload2->getRecord(j_param) == nullptr) {
0391 continue;
0392 }
0393
0394 const JetResolutionObject::Record record1 = *(payload1->getRecord(j_param));
0395 const JetResolutionObject::Record record2 = *(payload2->getRecord(j_param));
0396 float res1 = payload1->evaluateFormula(record1, j_param);
0397 resol_eta1->SetBinContent(idx + 1, res1);
0398 float res2 = payload2->evaluateFormula(record2, j_param);
0399 resol_eta2->SetBinContent(idx + 1, res2);
0400 }
0401
0402 for (size_t idx = 0; idx <= NBIN_PT; idx++) {
0403 double x_axis = (idx + 0.5) * (MAX_PT - MIN_PT) / NBIN_PT + MIN_PT;
0404
0405 JetParameters j_param;
0406 j_param.setJetEta(par_Eta);
0407 j_param.setRho(par_Rho);
0408 j_param.setJetPt(x_axis);
0409
0410 if (payload1->getRecord(j_param) == nullptr || payload2->getRecord(j_param) == nullptr) {
0411 continue;
0412 }
0413
0414 const JetResolutionObject::Record record1 = *(payload1->getRecord(j_param));
0415 const JetResolutionObject::Record record2 = *(payload2->getRecord(j_param));
0416 float res1 = payload1->evaluateFormula(record1, j_param);
0417 resol_pt1->SetBinContent(idx + 1, res1);
0418 float res2 = payload2->evaluateFormula(record2, j_param);
0419 resol_pt2->SetBinContent(idx + 1, res2);
0420 }
0421
0422 gStyle->SetOptStat(0);
0423 gStyle->SetLabelFont(42, "XYZ");
0424 gStyle->SetLabelSize(0.05, "XYZ");
0425 gStyle->SetFrameLineWidth(3);
0426
0427 std::string title = Form("Comparison between %s and %s", tag1.name.c_str(), tag2.name.c_str());
0428 TCanvas canvas("Jet Resolution Comparison", title.c_str(), 800, 1200);
0429 canvas.Divide(1, 2);
0430
0431 canvas.cd(1);
0432 resol_eta1->SetTitle("Jet Resolution Comparison vs. #eta");
0433 resol_eta1->SetXTitle("#eta");
0434 resol_eta1->SetYTitle("Resolution");
0435 resol_eta1->SetLineWidth(3);
0436 resol_eta1->SetMaximum(resol_eta1->GetMaximum() * 1.25);
0437 resol_eta1->Draw("][");
0438
0439 resol_eta2->SetLineColor(2);
0440 resol_eta2->SetLineWidth(3);
0441 resol_eta2->SetLineStyle(2);
0442 resol_eta2->Draw("][same");
0443
0444 leg_eta->AddEntry(resol_eta1, tag_ver1.c_str(), "l");
0445 leg_eta->AddEntry(resol_eta2, tag_ver2.c_str(), "l");
0446 leg_eta->AddEntry((TObject*)nullptr, Form("JetPt=%.2f; JetRho=%.2f", par_Pt, par_Rho), "");
0447 leg_eta->Draw();
0448
0449 canvas.cd(2);
0450 resol_pt1->SetTitle("Jet Resolution Comparison vs. p_{T}");
0451 resol_pt1->SetXTitle("p_{T} [GeV]");
0452 resol_pt1->SetYTitle("Resolution");
0453 resol_pt1->SetLineWidth(3);
0454 resol_pt1->Draw("][");
0455
0456 resol_pt2->SetLineColor(2);
0457 resol_pt2->SetLineWidth(3);
0458 resol_pt2->SetLineStyle(2);
0459 resol_pt2->Draw("][same");
0460
0461 leg_pt->AddEntry(resol_pt1, tag_ver1.c_str(), "l");
0462 leg_pt->AddEntry(resol_pt2, tag_ver2.c_str(), "l");
0463 leg_pt->AddEntry((TObject*)nullptr, Form("JetEta=%.2f; JetRho=%.2f", par_Eta, par_Rho), "");
0464 leg_pt->Draw();
0465
0466 canvas.SaveAs(m_imageFileName.c_str());
0467
0468 return true;
0469 } else
0470 return false;
0471 }
0472
0473 };
0474
0475 template <index ii>
0476 class JetScaleFactorVsEta : public cond::payloadInspector::Histogram1D<JetResolutionObject, SINGLE_IOV> {
0477 public:
0478 JetScaleFactorVsEta()
0479 : cond::payloadInspector::Histogram1D<JetResolutionObject, SINGLE_IOV>(
0480 "Jet Energy Scale Factor", "#eta", NBIN_ETA, MIN_ETA, MAX_ETA, "Scale Factor") {
0481 cond::payloadInspector::PlotBase::addInputParam("Jet_Pt");
0482 cond::payloadInspector::PlotBase::addInputParam("Jet_Eta");
0483 }
0484
0485 bool fill() override {
0486 double par_Pt = 100.;
0487 double par_Eta = 1.;
0488
0489 auto paramValues = cond::payloadInspector::PlotBase::inputParamValues();
0490 auto ip = paramValues.find("Jet_Pt");
0491 if (ip != paramValues.end()) {
0492 par_Pt = std::stod(ip->second);
0493 edm::LogPrint("JER_PI") << "Jet Pt: " << par_Pt;
0494 }
0495 ip = paramValues.find("Jet_Eta");
0496 if (ip != paramValues.end()) {
0497 par_Eta = std::stod(ip->second);
0498 edm::LogPrint("JER_PI") << "Jet Eta: " << par_Eta;
0499 }
0500
0501 auto tag = PlotBase::getTag<0>();
0502 for (auto const& iov : tag.iovs) {
0503 std::shared_ptr<JetResolutionObject> payload = Base::fetchPayload(std::get<1>(iov));
0504 if (payload.get()) {
0505 if (!payload->getRecords().empty() &&
0506 payload->getDefinition().getFormulaString().compare("") != 0)
0507 return false;
0508
0509 for (size_t idx = 0; idx <= NBIN_ETA; idx++) {
0510 par_Eta = (idx + 0.5) * (MAX_ETA - MIN_ETA) / NBIN_ETA + MIN_ETA;
0511
0512 JetParameters j_param;
0513 j_param.setJetPt(par_Pt);
0514 j_param.setJetEta(par_Eta);
0515
0516 if (payload->getRecord(j_param) == nullptr) {
0517 continue;
0518 }
0519
0520 const JetResolutionObject::Record record = *(payload->getRecord(j_param));
0521 if (record.getParametersValues().size() == 3) {
0522 double sf = record.getParametersValues()[ii];
0523 fillWithBinAndValue(idx + 1, sf);
0524 }
0525 }
0526 } else
0527 return false;
0528 }
0529 return true;
0530 }
0531 };
0532
0533 typedef JetScaleFactorVsEta<NORM> JetScaleFactorVsEtaNORM;
0534 typedef JetScaleFactorVsEta<DOWN> JetScaleFactorVsEtaDOWN;
0535 typedef JetScaleFactorVsEta<UP> JetScaleFactorVsEtaUP;
0536
0537 template <index ii>
0538 class JetScaleFactorVsPt : public cond::payloadInspector::Histogram1D<JetResolutionObject, SINGLE_IOV> {
0539 public:
0540 JetScaleFactorVsPt()
0541 : cond::payloadInspector::Histogram1D<JetResolutionObject, SINGLE_IOV>(
0542 "Jet Energy Scale Factor", "p_T", NBIN_PT, MIN_PT, MAX_PT, "Scale Factor") {
0543 cond::payloadInspector::PlotBase::addInputParam("Jet_Pt");
0544 cond::payloadInspector::PlotBase::addInputParam("Jet_Eta");
0545 }
0546
0547 bool fill() override {
0548 double par_Pt = 100.;
0549 double par_Eta = 1.;
0550
0551 auto paramValues = cond::payloadInspector::PlotBase::inputParamValues();
0552 auto ip = paramValues.find("Jet_Pt");
0553 if (ip != paramValues.end()) {
0554 par_Pt = std::stod(ip->second);
0555 edm::LogPrint("JER_PI") << "Jet Pt: " << par_Pt;
0556 }
0557 ip = paramValues.find("Jet_Eta");
0558 if (ip != paramValues.end()) {
0559 par_Eta = std::stod(ip->second);
0560 edm::LogPrint("JER_PI") << "Jet Eta: " << par_Eta;
0561 }
0562
0563 auto tag = PlotBase::getTag<0>();
0564 for (auto const& iov : tag.iovs) {
0565 std::shared_ptr<JetResolutionObject> payload = Base::fetchPayload(std::get<1>(iov));
0566 if (payload.get()) {
0567 if (!payload->getRecords().empty() &&
0568 payload->getDefinition().getFormulaString().compare("") != 0)
0569 return false;
0570
0571 for (size_t idx = 0; idx <= NBIN_PT; idx++) {
0572 par_Pt = (idx + 0.5) * (MAX_PT - MIN_PT) / NBIN_PT + MIN_PT;
0573
0574 JetParameters j_param;
0575 j_param.setJetEta(par_Eta);
0576 j_param.setJetPt(par_Pt);
0577
0578 if (payload->getRecord(j_param) == nullptr) {
0579 continue;
0580 }
0581
0582 const JetResolutionObject::Record record = *(payload->getRecord(j_param));
0583 if (record.getParametersValues().size() == 3) {
0584 double sf = record.getParametersValues()[ii];
0585 fillWithBinAndValue(idx + 1, sf);
0586 }
0587 }
0588 } else
0589 return false;
0590 }
0591 return true;
0592 }
0593 };
0594
0595 typedef JetScaleFactorVsPt<NORM> JetScaleFactorVsPtNORM;
0596 typedef JetScaleFactorVsPt<DOWN> JetScaleFactorVsPtDOWN;
0597 typedef JetScaleFactorVsPt<UP> JetScaleFactorVsPtUP;
0598
0599 class JetScaleFactorSummary : public cond::payloadInspector::PlotImage<JetResolutionObject, SINGLE_IOV> {
0600 public:
0601 JetScaleFactorSummary()
0602 : cond::payloadInspector::PlotImage<JetResolutionObject, SINGLE_IOV>("Jet ScaleFactor Summary") {
0603 cond::payloadInspector::PlotBase::addInputParam("Jet_Pt");
0604 cond::payloadInspector::PlotBase::addInputParam("Jet_Eta");
0605 }
0606
0607 bool fill() override {
0608 double par_Pt = 100.;
0609 double par_Eta = 1.;
0610
0611 auto paramValues = cond::payloadInspector::PlotBase::inputParamValues();
0612 auto ip = paramValues.find("Jet_Pt");
0613 if (ip != paramValues.end()) {
0614 par_Pt = std::stod(ip->second);
0615 }
0616 ip = paramValues.find("Jet_Eta");
0617 if (ip != paramValues.end()) {
0618 par_Eta = std::stod(ip->second);
0619 }
0620
0621 TH1D* sf_eta_norm = new TH1D("Jet SF vs #eta NORM", "", NBIN_ETA, MIN_ETA, MAX_ETA);
0622 TH1D* sf_eta_down = new TH1D("Jet SF vs #eta DOWN", "", NBIN_ETA, MIN_ETA, MAX_ETA);
0623 TH1D* sf_eta_up = new TH1D("Jet SF vs #eta UP", "", NBIN_ETA, MIN_ETA, MAX_ETA);
0624 TH1D* sf_pt_norm = new TH1D("Jet SF vs p_T NORM", "", NBIN_PT, MIN_PT, MAX_PT);
0625 TH1D* sf_pt_down = new TH1D("Jet SF vs p_T DOWN", "", NBIN_PT, MIN_PT, MAX_PT);
0626 TH1D* sf_pt_up = new TH1D("Jet SF vs p_T UP", "", NBIN_PT, MIN_PT, MAX_PT);
0627
0628 TLegend* leg_eta = new TLegend(0.26, 0.73, 0.935, 0.90);
0629 TLegend* leg_pt = new TLegend(0.26, 0.73, 0.935, 0.90);
0630
0631 leg_eta->SetBorderSize(0);
0632 leg_eta->SetLineStyle(0);
0633 leg_eta->SetFillStyle(0);
0634
0635 leg_eta->SetTextFont(42);
0636 leg_pt->SetBorderSize(0);
0637 leg_pt->SetLineStyle(0);
0638 leg_pt->SetFillStyle(0);
0639
0640 auto tag = PlotBase::getTag<0>();
0641 auto iov = tag.iovs.front();
0642 std::shared_ptr<JetResolutionObject> payload = fetchPayload(std::get<1>(iov));
0643 unsigned int run = std::get<0>(iov);
0644 std::string tagname = tag.name;
0645 std::stringstream ss_tagname(tag.name);
0646 std::string stmp;
0647
0648 std::string tag_ver;
0649 std::string tag_res;
0650 std::string tag_jet;
0651
0652 getline(ss_tagname, stmp, '_');
0653 getline(ss_tagname, stmp, '_');
0654 tag_ver = stmp;
0655 getline(ss_tagname, stmp, '_');
0656 tag_ver += '_' + stmp;
0657 getline(ss_tagname, stmp, '_');
0658 tag_ver += '_' + stmp;
0659 getline(ss_tagname, stmp, '_');
0660 tag_ver += '_' + stmp;
0661 getline(ss_tagname, stmp, '_');
0662 tag_res = stmp;
0663 getline(ss_tagname, stmp, '_');
0664 tag_jet = stmp;
0665
0666 bool is_2bin = false;
0667
0668 if (payload.get()) {
0669 if (!payload->getRecords().empty() &&
0670 payload->getDefinition().getFormulaString().compare("") != 0)
0671 return false;
0672
0673 is_2bin = false;
0674 for (size_t idx = 0; idx <= NBIN_ETA; idx++) {
0675 double x_axis = (idx + 0.5) * (MAX_ETA - MIN_ETA) / NBIN_ETA + MIN_ETA;
0676
0677 JetParameters j_param;
0678 j_param.setJetPt(par_Pt);
0679 j_param.setJetEta(x_axis);
0680
0681 if (payload->getRecord(j_param) == nullptr) {
0682 continue;
0683 }
0684
0685 const JetResolutionObject::Record record = *(payload->getRecord(j_param));
0686 if (record.getBinsRange().size() > 1)
0687 is_2bin = true;
0688
0689 if (record.getParametersValues().size() == 3) {
0690 sf_eta_norm->SetBinContent(idx + 1, record.getParametersValues()[0]);
0691 sf_eta_down->SetBinContent(idx + 1, record.getParametersValues()[1]);
0692 sf_eta_up->SetBinContent(idx + 1, record.getParametersValues()[2]);
0693 }
0694 }
0695
0696 for (size_t idx = 0; idx <= NBIN_PT; idx++) {
0697 double x_axis = (idx + 0.5) * (MAX_PT - MIN_PT) / NBIN_PT + MIN_PT;
0698
0699 JetParameters j_param;
0700 j_param.setJetEta(par_Eta);
0701 j_param.setJetPt(x_axis);
0702
0703 if (payload->getRecord(j_param) == nullptr) {
0704 continue;
0705 }
0706
0707 const JetResolutionObject::Record record = *(payload->getRecord(j_param));
0708 if (record.getParametersValues().size() == 3) {
0709 sf_pt_norm->SetBinContent(idx + 1, record.getParametersValues()[0]);
0710 sf_pt_down->SetBinContent(idx + 1, record.getParametersValues()[1]);
0711 sf_pt_up->SetBinContent(idx + 1, record.getParametersValues()[2]);
0712 }
0713 }
0714
0715 gStyle->SetOptStat(0);
0716 gStyle->SetLabelFont(42, "XYZ");
0717 gStyle->SetLabelSize(0.05, "XYZ");
0718 gStyle->SetFrameLineWidth(3);
0719
0720 std::string title = Form("Summary Run %i", run);
0721 TCanvas canvas("Jet ScaleFactor Summary", title.c_str(), 800, 1200);
0722 canvas.Divide(1, 2);
0723
0724 canvas.cd(1);
0725 sf_eta_up->SetTitle("ScaleFactor vs. #eta");
0726 sf_eta_up->SetXTitle("#eta");
0727 sf_eta_up->SetYTitle("Scale Factor");
0728 sf_eta_up->SetLineStyle(7);
0729 sf_eta_up->SetLineWidth(3);
0730 sf_eta_up->SetFillColorAlpha(kGray, 0.5);
0731 sf_eta_up->SetMaximum(sf_eta_up->GetMaximum() * 1.25);
0732 sf_eta_up->SetMinimum(0.);
0733 sf_eta_up->Draw("][");
0734
0735 sf_eta_down->SetLineStyle(7);
0736 sf_eta_down->SetLineWidth(3);
0737 sf_eta_down->SetFillColorAlpha(kWhite, 1);
0738 sf_eta_down->Draw("][ same");
0739
0740 sf_eta_norm->SetLineStyle(1);
0741 sf_eta_norm->SetLineWidth(5);
0742 sf_eta_norm->SetFillColor(0);
0743 sf_eta_norm->Draw("][ same");
0744 sf_eta_norm->Draw("axis same");
0745
0746 leg_eta->AddEntry(sf_eta_norm, (tag_ver + '_' + tag_jet).c_str(), "l");
0747 leg_eta->AddEntry((TObject*)nullptr, Form("JetPt=%.2f", par_Pt), "");
0748 leg_eta->Draw();
0749
0750 if (is_2bin == true) {
0751 canvas.cd(2);
0752 sf_pt_up->SetTitle("ScaleFactor vs. p_{T}");
0753 sf_pt_up->SetXTitle("p_{T} [GeV]");
0754 sf_pt_up->SetYTitle("Scale Factor");
0755 sf_pt_up->SetLineStyle(7);
0756 sf_pt_up->SetLineWidth(3);
0757 sf_pt_up->SetFillColorAlpha(kGray, 0.5);
0758 sf_pt_up->SetMaximum(sf_pt_up->GetMaximum() * 1.25);
0759 sf_pt_up->SetMinimum(0.);
0760 sf_pt_up->Draw("][");
0761
0762 sf_pt_down->SetLineStyle(7);
0763 sf_pt_down->SetLineWidth(3);
0764 sf_pt_down->SetFillColorAlpha(kWhite, 1);
0765 sf_pt_down->Draw("][ same");
0766
0767 sf_pt_norm->SetLineStyle(1);
0768 sf_pt_norm->SetLineWidth(5);
0769 sf_pt_norm->SetFillColor(0);
0770 sf_pt_norm->Draw("][ same");
0771 sf_pt_norm->Draw("axis same");
0772
0773 leg_pt->AddEntry(sf_pt_norm, (tag_ver + '_' + tag_jet).c_str(), "l");
0774 leg_pt->AddEntry((TObject*)nullptr, Form("JetEta=%.2f", par_Eta), "");
0775 leg_pt->Draw();
0776 }
0777
0778 canvas.SaveAs(m_imageFileName.c_str());
0779
0780 return true;
0781 } else
0782 return false;
0783 }
0784
0785 };
0786
0787 class JetScaleFactorCompare : public cond::payloadInspector::PlotImage<JetResolutionObject, SINGLE_IOV, 2> {
0788 public:
0789 JetScaleFactorCompare()
0790 : cond::payloadInspector::PlotImage<JetResolutionObject, SINGLE_IOV, 2>("Jet ScaleFactor Summary") {
0791 cond::payloadInspector::PlotBase::addInputParam("Jet_Pt");
0792 cond::payloadInspector::PlotBase::addInputParam("Jet_Eta");
0793 }
0794
0795 bool fill() override {
0796 double par_Pt = 100.;
0797 double par_Eta = 1.;
0798
0799 auto paramValues = cond::payloadInspector::PlotBase::inputParamValues();
0800 auto ip = paramValues.find("Jet_Pt");
0801 if (ip != paramValues.end()) {
0802 par_Pt = std::stod(ip->second);
0803 }
0804 ip = paramValues.find("Jet_Eta");
0805 if (ip != paramValues.end()) {
0806 par_Eta = std::stod(ip->second);
0807 }
0808
0809 TH1D* sf_eta_one = new TH1D("Jet SF vs #eta one", "", NBIN_ETA, MIN_ETA, MAX_ETA);
0810 TH1D* sf_eta_two = new TH1D("Jet SF vs #eta two", "", NBIN_ETA, MIN_ETA, MAX_ETA);
0811 TH1D* sf_pt_one = new TH1D("Jet SF vs p_T one", "", NBIN_PT, MIN_PT, MAX_PT);
0812 TH1D* sf_pt_two = new TH1D("Jet SF vs p_T two", "", NBIN_PT, MIN_PT, MAX_PT);
0813
0814 TLegend* leg_eta = new TLegend(0.26, 0.73, 0.935, 0.90);
0815 TLegend* leg_pt = new TLegend(0.26, 0.73, 0.935, 0.90);
0816
0817 leg_eta->SetBorderSize(0);
0818 leg_eta->SetLineStyle(0);
0819 leg_eta->SetFillStyle(0);
0820
0821 leg_eta->SetTextFont(42);
0822 leg_pt->SetBorderSize(0);
0823 leg_pt->SetLineStyle(0);
0824 leg_pt->SetFillStyle(0);
0825
0826 auto tag1 = PlotBase::getTag<0>();
0827 auto iov1 = tag1.iovs.front();
0828 std::shared_ptr<JetResolutionObject> payload1 = fetchPayload(std::get<1>(iov1));
0829
0830 auto tag2 = PlotBase::getTag<1>();
0831 auto iov2 = tag2.iovs.front();
0832 std::shared_ptr<JetResolutionObject> payload2 = fetchPayload(std::get<1>(iov2));
0833
0834 std::stringstream ss_tagname1(tag1.name);
0835 std::stringstream ss_tagname2(tag2.name);
0836 std::string stmp;
0837
0838 std::string tag_ver1;
0839 std::string tag_ver2;
0840
0841 getline(ss_tagname1, stmp, '_');
0842 getline(ss_tagname1, stmp);
0843 tag_ver1 = stmp;
0844
0845 getline(ss_tagname2, stmp, '_');
0846 getline(ss_tagname2, stmp);
0847 tag_ver2 = stmp;
0848
0849 bool is_2bin = false;
0850
0851 if (payload1.get() && payload2.get()) {
0852 if (!payload1->getRecords().empty() &&
0853 payload1->getDefinition().getFormulaString().compare("") != 0)
0854 return false;
0855
0856 is_2bin = false;
0857 for (size_t idx = 0; idx <= NBIN_ETA; idx++) {
0858 double x_axis = (idx + 0.5) * (MAX_ETA - MIN_ETA) / NBIN_ETA + MIN_ETA;
0859
0860 JetParameters j_param;
0861 j_param.setJetPt(par_Pt);
0862 j_param.setJetEta(x_axis);
0863
0864 const JetResolutionObject::Record* rcd_p1 = payload1->getRecord(j_param);
0865 const JetResolutionObject::Record* rcd_p2 = payload2->getRecord(j_param);
0866 if (rcd_p1 == nullptr || rcd_p2 == nullptr) {
0867 continue;
0868 }
0869
0870 const JetResolutionObject::Record record1 = *(rcd_p1);
0871 const JetResolutionObject::Record record2 = *(rcd_p2);
0872
0873 if (record1.getBinsRange().size() > 1 && record2.getBinsRange().size() > 1)
0874 is_2bin = true;
0875
0876 if (record1.getParametersValues().size() == 3) {
0877 sf_eta_one->SetBinContent(idx + 1, record1.getParametersValues()[0]);
0878 }
0879 if (record2.getParametersValues().size() == 3) {
0880 sf_eta_two->SetBinContent(idx + 1, record2.getParametersValues()[0]);
0881 }
0882 }
0883
0884 for (size_t idx = 0; idx <= NBIN_PT; idx++) {
0885 double x_axis = (idx + 0.5) * (MAX_PT - MIN_PT) / NBIN_PT + MIN_PT;
0886
0887 JetParameters j_param;
0888 j_param.setJetEta(par_Eta);
0889 j_param.setJetPt(x_axis);
0890
0891 const JetResolutionObject::Record* rcd_p1 = payload1->getRecord(j_param);
0892 const JetResolutionObject::Record* rcd_p2 = payload2->getRecord(j_param);
0893 if (rcd_p1 == nullptr || rcd_p2 == nullptr) {
0894 continue;
0895 }
0896
0897 const JetResolutionObject::Record record1 = *(rcd_p1);
0898 const JetResolutionObject::Record record2 = *(rcd_p2);
0899
0900 if (record1.getParametersValues().size() == 3) {
0901 sf_pt_one->SetBinContent(idx + 1, record1.getParametersValues()[0]);
0902 }
0903 if (record2.getParametersValues().size() == 3) {
0904 sf_pt_two->SetBinContent(idx + 1, record2.getParametersValues()[0]);
0905 }
0906 }
0907
0908 gStyle->SetOptStat(0);
0909 gStyle->SetLabelFont(42, "XYZ");
0910 gStyle->SetLabelSize(0.05, "XYZ");
0911 gStyle->SetFrameLineWidth(3);
0912
0913 std::string title = Form("Comparison");
0914 TCanvas canvas("Jet ScaleFactor", title.c_str(), 800, 1200);
0915 canvas.Divide(1, 2);
0916
0917 canvas.cd(1);
0918 sf_eta_one->SetTitle("Jet ScaleFactor Comparison vs. #eta");
0919 sf_eta_one->SetXTitle("#eta");
0920 sf_eta_one->SetYTitle("Scale Factor");
0921 sf_eta_one->SetLineWidth(3);
0922 sf_eta_one->SetMaximum(sf_eta_one->GetMaximum() * 1.25);
0923 sf_eta_one->SetMinimum(0.);
0924 sf_eta_one->Draw("][");
0925
0926 sf_eta_two->SetLineColor(2);
0927 sf_eta_two->SetLineWidth(3);
0928 sf_eta_two->SetLineStyle(2);
0929 sf_eta_two->Draw("][ same");
0930
0931 leg_eta->AddEntry(sf_eta_one, tag_ver1.c_str(), "l");
0932 leg_eta->AddEntry(sf_eta_two, tag_ver2.c_str(), "l");
0933 leg_eta->AddEntry((TObject*)nullptr, Form("JetPt=%.2f", par_Pt), "");
0934 leg_eta->Draw();
0935
0936 if (is_2bin == true) {
0937 canvas.cd(2);
0938 sf_pt_one->SetTitle("Jet ScaleFactor Comparison vs. p_{T}");
0939 sf_pt_one->SetXTitle("p_{T} [GeV]");
0940 sf_pt_one->SetYTitle("Scale Factor");
0941 sf_pt_one->SetLineWidth(3);
0942 sf_pt_one->SetMaximum(sf_pt_one->GetMaximum() * 1.25);
0943 sf_pt_one->SetMinimum(0.);
0944 sf_pt_one->Draw("][");
0945
0946 sf_pt_two->SetLineColor(2);
0947 sf_pt_two->SetLineWidth(3);
0948 sf_pt_two->SetLineStyle(2);
0949 sf_pt_two->Draw("][ same");
0950
0951 leg_pt->AddEntry(sf_pt_one, tag_ver1.c_str(), "l");
0952 leg_pt->AddEntry(sf_pt_two, tag_ver2.c_str(), "l");
0953 leg_pt->AddEntry((TObject*)nullptr, Form("JetEta=%.2f", par_Eta), "");
0954 leg_pt->Draw();
0955 }
0956
0957 canvas.SaveAs(m_imageFileName.c_str());
0958
0959 return true;
0960 } else
0961 return false;
0962 }
0963
0964 };
0965
0966
0967 PAYLOAD_INSPECTOR_MODULE(JetResolutionObject) {
0968 PAYLOAD_INSPECTOR_CLASS(JetResolutionVsEta);
0969 PAYLOAD_INSPECTOR_CLASS(JetResolutionVsPt);
0970 PAYLOAD_INSPECTOR_CLASS(JetScaleFactorVsEtaNORM);
0971 PAYLOAD_INSPECTOR_CLASS(JetScaleFactorVsEtaDOWN);
0972 PAYLOAD_INSPECTOR_CLASS(JetScaleFactorVsEtaUP);
0973 PAYLOAD_INSPECTOR_CLASS(JetScaleFactorVsPtNORM);
0974 PAYLOAD_INSPECTOR_CLASS(JetScaleFactorVsPtDOWN);
0975 PAYLOAD_INSPECTOR_CLASS(JetScaleFactorVsPtUP);
0976
0977 PAYLOAD_INSPECTOR_CLASS(JetResolutionSummary);
0978 PAYLOAD_INSPECTOR_CLASS(JetScaleFactorSummary);
0979 PAYLOAD_INSPECTOR_CLASS(JetResolutionCompare);
0980 PAYLOAD_INSPECTOR_CLASS(JetScaleFactorCompare);
0981 }
0982
0983 }