File indexing completed on 2023-05-26 01:15:26
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 (const auto& record : payload->getRecords()) {
0084
0085 if (!record.getVariablesRange().empty() && payload->getDefinition().getVariableName(0) == "JetPt" &&
0086 record.getVariablesRange()[0].is_inside(par_Pt)) {
0087 if (record.getBinsRange().size() > 1 && payload->getDefinition().getBinName(1) == "Rho" &&
0088 record.getBinsRange()[1].is_inside(par_Rho)) {
0089 if (!record.getBinsRange().empty() && payload->getDefinition().getBinName(0) == "JetEta") {
0090 reco::FormulaEvaluator f(payload->getDefinition().getFormulaString());
0091
0092 for (size_t idx = 0; idx <= NBIN_ETA; idx++) {
0093 par_Eta = (idx + 0.5) * (MAX_ETA - MIN_ETA) / NBIN_ETA + MIN_ETA;
0094 if (record.getBinsRange()[0].is_inside(par_Eta)) {
0095 std::vector<double> var = {par_Pt};
0096 std::vector<double> param;
0097 for (size_t i = 0; i < record.getParametersValues().size(); i++) {
0098 double par = record.getParametersValues()[i];
0099 param.push_back(par);
0100 }
0101 float res = f.evaluate(var, param);
0102 fillWithBinAndValue(idx + 1, res);
0103 }
0104 }
0105 }
0106 }
0107 }
0108 }
0109 return true;
0110 }
0111 }
0112 return false;
0113 }
0114 };
0115
0116 class JetResolutionVsPt : public cond::payloadInspector::Histogram1D<JetResolutionObject, SINGLE_IOV> {
0117 public:
0118 JetResolutionVsPt()
0119 : cond::payloadInspector::Histogram1D<JetResolutionObject, SINGLE_IOV>(
0120 "Jet Energy Resolution", "p_T", NBIN_PT, MIN_PT, MAX_PT, "Resolution") {
0121 cond::payloadInspector::PlotBase::addInputParam("Jet_Eta");
0122 cond::payloadInspector::PlotBase::addInputParam("Jet_Rho");
0123 }
0124
0125 bool fill() override {
0126 double par_Pt = 100.;
0127 double par_Eta = 1.;
0128 double par_Rho = 20.;
0129
0130 auto paramValues = cond::payloadInspector::PlotBase::inputParamValues();
0131 auto ip = paramValues.find("Jet_Eta");
0132 if (ip != paramValues.end()) {
0133 par_Eta = std::stod(ip->second);
0134 edm::LogPrint("JER_PI") << "Jet Eta: " << par_Eta;
0135 }
0136 ip = paramValues.find("Jet_Rho");
0137 if (ip != paramValues.end()) {
0138 par_Rho = std::stod(ip->second);
0139 edm::LogPrint("JER_PI") << "Jet Rho: " << par_Rho;
0140 }
0141
0142 auto tag = PlotBase::getTag<0>();
0143 for (auto const& iov : tag.iovs) {
0144 std::shared_ptr<JetResolutionObject> payload = Base::fetchPayload(std::get<1>(iov));
0145 if (payload.get()) {
0146 if (!payload->getRecords().empty() &&
0147 payload->getDefinition().getFormulaString().compare("") == 0)
0148 return false;
0149
0150 for (const auto& record : payload->getRecords()) {
0151
0152 if (!record.getBinsRange().empty() && payload->getDefinition().getBinName(0) == "JetEta" &&
0153 record.getBinsRange()[0].is_inside(par_Eta)) {
0154 if (record.getBinsRange().size() > 1 && payload->getDefinition().getBinName(1) == "Rho" &&
0155 record.getBinsRange()[1].is_inside(par_Rho)) {
0156 if (!record.getVariablesRange().empty() && payload->getDefinition().getVariableName(0) == "JetPt") {
0157 reco::FormulaEvaluator f(payload->getDefinition().getFormulaString());
0158
0159 for (size_t idx = 0; idx <= NBIN_PT; idx++) {
0160 par_Pt = (idx + 0.5) * (MAX_PT - MIN_PT) / NBIN_PT + MIN_PT;
0161 if (record.getVariablesRange()[0].is_inside(par_Pt)) {
0162 std::vector<double> var = {par_Pt};
0163 std::vector<double> param;
0164 for (size_t i = 0; i < record.getParametersValues().size(); i++) {
0165 double par = record.getParametersValues()[i];
0166 param.push_back(par);
0167 }
0168 float res = f.evaluate(var, param);
0169 fillWithBinAndValue(idx + 1, res);
0170 }
0171 }
0172 }
0173 }
0174 }
0175 }
0176 return true;
0177 }
0178 }
0179 return false;
0180 }
0181 };
0182
0183 class JetResolutionSummary : public cond::payloadInspector::PlotImage<JetResolutionObject, SINGLE_IOV> {
0184 public:
0185 JetResolutionSummary()
0186 : cond::payloadInspector::PlotImage<JetResolutionObject, SINGLE_IOV>("Jet Resolution Summary") {
0187 cond::payloadInspector::PlotBase::addInputParam("Jet_Pt");
0188 cond::payloadInspector::PlotBase::addInputParam("Jet_Eta");
0189 cond::payloadInspector::PlotBase::addInputParam("Jet_Rho");
0190 }
0191
0192 bool fill() override {
0193 double par_Pt = 100.;
0194 double par_Eta = 1.;
0195 double par_Rho = 20.;
0196
0197 auto paramValues = cond::payloadInspector::PlotBase::inputParamValues();
0198 auto ip = paramValues.find("Jet_Pt");
0199 if (ip != paramValues.end()) {
0200 par_Pt = std::stod(ip->second);
0201 }
0202 ip = paramValues.find("Jet_Eta");
0203 if (ip != paramValues.end()) {
0204 par_Eta = std::stod(ip->second);
0205 }
0206 ip = paramValues.find("Jet_Rho");
0207 if (ip != paramValues.end()) {
0208 par_Rho = std::stod(ip->second);
0209 }
0210
0211 TH1D* resol_eta = new TH1D("Jet Resolution vs #eta", "", NBIN_ETA, MIN_ETA, MAX_ETA);
0212 TH1D* resol_pt = new TH1D("Jet Resolution vs p_T", "", NBIN_PT, MIN_PT, MAX_PT);
0213 TLegend* leg_eta = new TLegend(0.26, 0.73, 0.935, 0.90);
0214 TLegend* leg_pt = new TLegend(0.26, 0.73, 0.935, 0.90);
0215
0216 leg_eta->SetBorderSize(0);
0217 leg_eta->SetLineStyle(0);
0218 leg_eta->SetFillStyle(0);
0219
0220 leg_eta->SetTextFont(42);
0221 leg_pt->SetBorderSize(0);
0222 leg_pt->SetLineStyle(0);
0223 leg_pt->SetFillStyle(0);
0224
0225 auto tag = PlotBase::getTag<0>();
0226 auto iov = tag.iovs.front();
0227 std::shared_ptr<JetResolutionObject> payload = fetchPayload(std::get<1>(iov));
0228 unsigned int run = std::get<0>(iov);
0229 std::string tagname = tag.name;
0230 std::stringstream ss_tagname(tag.name);
0231 std::string stmp;
0232
0233 std::string tag_ver;
0234 std::string tag_res;
0235 std::string tag_jet;
0236
0237 getline(ss_tagname, stmp, '_');
0238 getline(ss_tagname, stmp, '_');
0239 tag_ver = stmp;
0240 getline(ss_tagname, stmp, '_');
0241 tag_ver += '_' + stmp;
0242 getline(ss_tagname, stmp, '_');
0243 tag_ver += '_' + stmp;
0244 getline(ss_tagname, stmp, '_');
0245 tag_ver += '_' + stmp;
0246 getline(ss_tagname, stmp, '_');
0247 tag_res = stmp;
0248 getline(ss_tagname, stmp, '_');
0249 tag_jet = stmp;
0250
0251 if (payload.get()) {
0252 if (!payload->getRecords().empty() &&
0253 payload->getDefinition().getFormulaString().compare("") == 0)
0254 return false;
0255
0256 for (const auto& record : payload->getRecords()) {
0257
0258 if (!record.getVariablesRange().empty() && payload->getDefinition().getVariableName(0) == "JetPt" &&
0259 record.getVariablesRange()[0].is_inside(par_Pt)) {
0260 if (record.getBinsRange().size() > 1 && payload->getDefinition().getBinName(1) == "Rho" &&
0261 record.getBinsRange()[1].is_inside(par_Rho)) {
0262 if (!record.getBinsRange().empty() && payload->getDefinition().getBinName(0) == "JetEta") {
0263 reco::FormulaEvaluator f(payload->getDefinition().getFormulaString());
0264
0265 for (size_t idx = 0; idx <= NBIN_ETA; idx++) {
0266 double x_axis = (idx + 0.5) * (MAX_ETA - MIN_ETA) / NBIN_ETA + MIN_ETA;
0267 if (record.getBinsRange()[0].is_inside(x_axis)) {
0268 std::vector<double> var = {par_Pt};
0269 std::vector<double> param;
0270 for (size_t i = 0; i < record.getParametersValues().size(); i++) {
0271 double par = record.getParametersValues()[i];
0272 param.push_back(par);
0273 }
0274 float res = f.evaluate(var, param);
0275 resol_eta->SetBinContent(idx + 1, res);
0276 }
0277 }
0278 }
0279 }
0280 }
0281
0282 if (!record.getBinsRange().empty() && payload->getDefinition().getBinName(0) == "JetEta" &&
0283 record.getBinsRange()[0].is_inside(par_Eta)) {
0284 if (record.getBinsRange().size() > 1 && payload->getDefinition().getBinName(1) == "Rho" &&
0285 record.getBinsRange()[1].is_inside(par_Rho)) {
0286 if (!record.getVariablesRange().empty() && payload->getDefinition().getVariableName(0) == "JetPt") {
0287 reco::FormulaEvaluator f(payload->getDefinition().getFormulaString());
0288
0289 for (size_t idx = 0; idx <= NBIN_PT + 2; idx++) {
0290 double x_axis = (idx + 0.5) * (MAX_PT - MIN_PT) / NBIN_PT + MIN_PT;
0291 if (record.getVariablesRange()[0].is_inside(x_axis)) {
0292 std::vector<double> var = {x_axis};
0293 std::vector<double> param;
0294 for (size_t i = 0; i < record.getParametersValues().size(); i++) {
0295 double par = record.getParametersValues()[i];
0296 param.push_back(par);
0297 }
0298 float res = f.evaluate(var, param);
0299 resol_pt->SetBinContent(idx + 1, res);
0300 }
0301 }
0302 }
0303 }
0304 }
0305 }
0306
0307 gStyle->SetOptStat(0);
0308 gStyle->SetLabelFont(42, "XYZ");
0309 gStyle->SetLabelSize(0.05, "XYZ");
0310 gStyle->SetFrameLineWidth(3);
0311
0312 std::string title = Form("Summary Run %i", run);
0313 TCanvas canvas("Jet Resolution Summary", title.c_str(), 800, 1200);
0314 canvas.Divide(1, 2);
0315
0316 canvas.cd(1);
0317 resol_eta->SetTitle(tag_res.c_str());
0318 resol_eta->SetXTitle("#eta");
0319 resol_eta->SetYTitle("Resolution");
0320 resol_eta->SetLineWidth(3);
0321 resol_eta->SetMaximum(resol_eta->GetMaximum() * 1.25);
0322 resol_eta->Draw("");
0323
0324 leg_eta->AddEntry(resol_eta, (tag_ver + '_' + tag_jet).c_str(), "l");
0325 leg_eta->AddEntry((TObject*)nullptr, Form("JetPt=%.2f; JetRho=%.2f", par_Pt, par_Rho), "");
0326 leg_eta->Draw();
0327
0328 canvas.cd(2);
0329 resol_pt->SetXTitle("p_{T} [GeV]");
0330 resol_pt->SetYTitle("Resolution");
0331 resol_pt->SetLineWidth(3);
0332 resol_pt->Draw("][");
0333
0334 leg_pt->AddEntry(resol_pt, (tag_ver + '_' + tag_jet).c_str(), "l");
0335 leg_pt->AddEntry((TObject*)nullptr, Form("JetEta=%.2f; JetRho=%.2f", par_Eta, par_Rho), "");
0336 leg_pt->Draw();
0337
0338 canvas.SaveAs(m_imageFileName.c_str());
0339
0340 return true;
0341 } else
0342 return false;
0343 }
0344
0345 };
0346
0347 template <index ii>
0348 class JetScaleFactorVsEta : public cond::payloadInspector::Histogram1D<JetResolutionObject, SINGLE_IOV> {
0349 public:
0350 JetScaleFactorVsEta()
0351 : cond::payloadInspector::Histogram1D<JetResolutionObject, SINGLE_IOV>(
0352 "Jet Energy Scale Factor", "#eta", NBIN_ETA, MIN_ETA, MAX_ETA, "Scale Factor") {
0353 cond::payloadInspector::PlotBase::addInputParam("Jet_Pt");
0354 cond::payloadInspector::PlotBase::addInputParam("Jet_Eta");
0355 }
0356
0357 bool fill() override {
0358 double par_Pt = 100.;
0359 double par_Eta = 1.;
0360
0361 auto paramValues = cond::payloadInspector::PlotBase::inputParamValues();
0362 auto ip = paramValues.find("Jet_Pt");
0363 if (ip != paramValues.end()) {
0364 par_Pt = std::stod(ip->second);
0365 edm::LogPrint("JER_PI") << "Jet Pt: " << par_Pt;
0366 }
0367 ip = paramValues.find("Jet_Eta");
0368 if (ip != paramValues.end()) {
0369 par_Eta = std::stod(ip->second);
0370 edm::LogPrint("JER_PI") << "Jet Eta: " << par_Eta;
0371 }
0372
0373 auto tag = PlotBase::getTag<0>();
0374 for (auto const& iov : tag.iovs) {
0375 std::shared_ptr<JetResolutionObject> payload = Base::fetchPayload(std::get<1>(iov));
0376 if (payload.get()) {
0377 if (!payload->getRecords().empty() &&
0378 payload->getDefinition().getFormulaString().compare("") != 0)
0379 return false;
0380
0381 for (const auto& record : payload->getRecords()) {
0382 if (!record.getBinsRange().empty() && payload->getDefinition().getBinName(0) == "JetEta" &&
0383 record.getParametersValues().size() == 3) {
0384
0385 if (record.getBinsRange().size() > 1 && payload->getDefinition().getBinName(1) == "JetPt" &&
0386 !record.getBinsRange()[1].is_inside(par_Pt))
0387 continue;
0388
0389 for (size_t it = 0; it <= NBIN_ETA; it++) {
0390 par_Eta = (it + 0.5) * (MAX_ETA - MIN_ETA) / NBIN_ETA + MIN_ETA;
0391 if (record.getBinsRange()[0].is_inside(par_Eta)) {
0392 double sf = 0.;
0393 sf = record.getParametersValues()[ii];
0394 fillWithBinAndValue(it, sf);
0395 }
0396 }
0397 }
0398 }
0399 return true;
0400 } else
0401 return false;
0402 }
0403 return false;
0404 }
0405 };
0406
0407 typedef JetScaleFactorVsEta<NORM> JetScaleFactorVsEtaNORM;
0408 typedef JetScaleFactorVsEta<DOWN> JetScaleFactorVsEtaDOWN;
0409 typedef JetScaleFactorVsEta<UP> JetScaleFactorVsEtaUP;
0410
0411 template <index ii>
0412 class JetScaleFactorVsPt : public cond::payloadInspector::Histogram1D<JetResolutionObject, SINGLE_IOV> {
0413 public:
0414 JetScaleFactorVsPt()
0415 : cond::payloadInspector::Histogram1D<JetResolutionObject, SINGLE_IOV>(
0416 "Jet Energy Scale Factor", "p_T", NBIN_PT, MIN_PT, MAX_PT, "Scale Factor") {
0417 cond::payloadInspector::PlotBase::addInputParam("Jet_Pt");
0418 cond::payloadInspector::PlotBase::addInputParam("Jet_Eta");
0419 }
0420
0421 bool fill() override {
0422 double par_Pt = 100.;
0423 double par_Eta = 1.;
0424
0425 auto paramValues = cond::payloadInspector::PlotBase::inputParamValues();
0426 auto ip = paramValues.find("Jet_Pt");
0427 if (ip != paramValues.end()) {
0428 par_Pt = std::stod(ip->second);
0429 edm::LogPrint("JER_PI") << "Jet Pt: " << par_Pt;
0430 }
0431 ip = paramValues.find("Jet_Eta");
0432 if (ip != paramValues.end()) {
0433 par_Eta = std::stod(ip->second);
0434 edm::LogPrint("JER_PI") << "Jet Eta: " << par_Eta;
0435 }
0436
0437 auto tag = PlotBase::getTag<0>();
0438 for (auto const& iov : tag.iovs) {
0439 std::shared_ptr<JetResolutionObject> payload = Base::fetchPayload(std::get<1>(iov));
0440 if (payload.get()) {
0441 if (!payload->getRecords().empty() &&
0442 payload->getDefinition().getFormulaString().compare("") != 0)
0443 return false;
0444
0445 for (const auto& record : payload->getRecords()) {
0446 if (record.getBinsRange().size() > 1 && payload->getDefinition().getBinName(0) == "JetEta" &&
0447 record.getBinsRange()[0].is_inside(par_Eta) &&
0448 payload->getDefinition().getBinName(1) == "JetPt" &&
0449 record.getParametersValues().size() == 3) {
0450
0451 for (size_t it = 0; it <= NBIN_PT; it++) {
0452 par_Pt = (it + 0.5) * (MAX_PT - MIN_PT) / NBIN_PT + MIN_PT;
0453 if (record.getBinsRange()[1].is_inside(par_Pt)) {
0454 double sf = 0.;
0455 sf = record.getParametersValues()[ii];
0456 fillWithBinAndValue(it, sf);
0457 }
0458 }
0459 }
0460 }
0461 return true;
0462 } else
0463 return false;
0464 }
0465 return false;
0466 }
0467 };
0468
0469 typedef JetScaleFactorVsPt<NORM> JetScaleFactorVsPtNORM;
0470 typedef JetScaleFactorVsPt<DOWN> JetScaleFactorVsPtDOWN;
0471 typedef JetScaleFactorVsPt<UP> JetScaleFactorVsPtUP;
0472
0473 class JetScaleFactorSummary : public cond::payloadInspector::PlotImage<JetResolutionObject, SINGLE_IOV> {
0474 public:
0475 JetScaleFactorSummary()
0476 : cond::payloadInspector::PlotImage<JetResolutionObject, SINGLE_IOV>("Jet ScaleFactor Summary") {
0477 cond::payloadInspector::PlotBase::addInputParam("Jet_Pt");
0478 cond::payloadInspector::PlotBase::addInputParam("Jet_Eta");
0479 }
0480
0481 bool fill() override {
0482 double par_Pt = 100.;
0483 double par_Eta = 1.;
0484
0485 auto paramValues = cond::payloadInspector::PlotBase::inputParamValues();
0486 auto ip = paramValues.find("Jet_Pt");
0487 if (ip != paramValues.end()) {
0488 par_Pt = std::stod(ip->second);
0489 }
0490 ip = paramValues.find("Jet_Eta");
0491 if (ip != paramValues.end()) {
0492 par_Eta = std::stod(ip->second);
0493 }
0494
0495 TH1D* sf_eta_norm = new TH1D("Jet SF vs #eta NORM", "", NBIN_ETA, MIN_ETA, MAX_ETA);
0496 TH1D* sf_eta_down = new TH1D("Jet SF vs #eta DOWN", "", NBIN_ETA, MIN_ETA, MAX_ETA);
0497 TH1D* sf_eta_up = new TH1D("Jet SF vs #eta UP", "", NBIN_ETA, MIN_ETA, MAX_ETA);
0498 TH1D* sf_pt_norm = new TH1D("Jet SF vs p_T NORM", "", NBIN_PT, MIN_PT, MAX_PT);
0499 TH1D* sf_pt_down = new TH1D("Jet SF vs p_T DOWN", "", NBIN_PT, MIN_PT, MAX_PT);
0500 TH1D* sf_pt_up = new TH1D("Jet SF vs p_T UP", "", NBIN_PT, MIN_PT, MAX_PT);
0501
0502 TLegend* leg_eta = new TLegend(0.26, 0.73, 0.935, 0.90);
0503 TLegend* leg_pt = new TLegend(0.26, 0.73, 0.935, 0.90);
0504
0505 leg_eta->SetBorderSize(0);
0506 leg_eta->SetLineStyle(0);
0507 leg_eta->SetFillStyle(0);
0508
0509 leg_eta->SetTextFont(42);
0510 leg_pt->SetBorderSize(0);
0511 leg_pt->SetLineStyle(0);
0512 leg_pt->SetFillStyle(0);
0513
0514 auto tag = PlotBase::getTag<0>();
0515 auto iov = tag.iovs.front();
0516 std::shared_ptr<JetResolutionObject> payload = fetchPayload(std::get<1>(iov));
0517 unsigned int run = std::get<0>(iov);
0518 std::string tagname = tag.name;
0519 std::stringstream ss_tagname(tag.name);
0520 std::string stmp;
0521
0522 std::string tag_ver;
0523 std::string tag_res;
0524 std::string tag_jet;
0525
0526 getline(ss_tagname, stmp, '_');
0527 getline(ss_tagname, stmp, '_');
0528 tag_ver = stmp;
0529 getline(ss_tagname, stmp, '_');
0530 tag_ver += '_' + stmp;
0531 getline(ss_tagname, stmp, '_');
0532 tag_ver += '_' + stmp;
0533 getline(ss_tagname, stmp, '_');
0534 tag_ver += '_' + stmp;
0535 getline(ss_tagname, stmp, '_');
0536 tag_res = stmp;
0537 getline(ss_tagname, stmp, '_');
0538 tag_jet = stmp;
0539
0540 bool is_2bin = false;
0541
0542 if (payload.get()) {
0543 if (!payload->getRecords().empty() &&
0544 payload->getDefinition().getFormulaString().compare("") != 0)
0545 return false;
0546
0547 is_2bin = false;
0548 for (const auto& record : payload->getRecords()) {
0549 if (record.getBinsRange().size() > 1)
0550 is_2bin = true;
0551
0552 if (!record.getBinsRange().empty() && payload->getDefinition().getBinName(0) == "JetEta" &&
0553 record.getParametersValues().size() == 3) {
0554
0555 for (size_t it = 0; it <= NBIN_ETA; it++) {
0556 double x_axis = (it + 0.5) * (MAX_ETA - MIN_ETA) / NBIN_ETA + MIN_ETA;
0557 if (((is_2bin == false) || (is_2bin == true && record.getBinsRange()[1].is_inside(par_Pt))) &&
0558 record.getBinsRange()[0].is_inside(x_axis)) {
0559 sf_eta_norm->SetBinContent(it + 1, record.getParametersValues()[0]);
0560 sf_eta_down->SetBinContent(it + 1, record.getParametersValues()[1]);
0561 sf_eta_up->SetBinContent(it + 1, record.getParametersValues()[2]);
0562 }
0563 }
0564 }
0565
0566 if (record.getBinsRange().size() > 1 && payload->getDefinition().getBinName(0) == "JetEta" &&
0567 record.getBinsRange()[0].is_inside(par_Eta) &&
0568 payload->getDefinition().getBinName(1) == "JetPt" &&
0569 record.getParametersValues().size() == 3) {
0570
0571 is_2bin = true;
0572
0573 for (size_t it = 0; it <= NBIN_PT; it++) {
0574 double x_axis = (it + 0.5) * (MAX_PT - MIN_PT) / NBIN_PT + MIN_PT;
0575 if (record.getBinsRange()[1].is_inside(x_axis)) {
0576 sf_pt_norm->SetBinContent(it + 1, record.getParametersValues()[0]);
0577 sf_pt_down->SetBinContent(it + 1, record.getParametersValues()[1]);
0578 sf_pt_up->SetBinContent(it + 1, record.getParametersValues()[2]);
0579 }
0580 }
0581 }
0582 }
0583
0584 gStyle->SetOptStat(0);
0585 gStyle->SetLabelFont(42, "XYZ");
0586 gStyle->SetLabelSize(0.05, "XYZ");
0587 gStyle->SetFrameLineWidth(3);
0588
0589 std::string title = Form("Summary Run %i", run);
0590 TCanvas canvas("Jet ScaleFactor Summary", title.c_str(), 800, 1200);
0591 canvas.Divide(1, 2);
0592
0593 canvas.cd(1);
0594 sf_eta_up->SetTitle("ScaleFactor vs. #eta");
0595 sf_eta_up->SetXTitle("#eta");
0596 sf_eta_up->SetYTitle("Scale Factor");
0597 sf_eta_up->SetLineStyle(7);
0598 sf_eta_up->SetLineWidth(3);
0599 sf_eta_up->SetFillColorAlpha(kGray, 0.5);
0600 sf_eta_up->SetMaximum(sf_eta_up->GetMaximum() * 1.25);
0601 sf_eta_up->SetMinimum(0.);
0602 sf_eta_up->Draw("][");
0603
0604 sf_eta_down->SetLineStyle(7);
0605 sf_eta_down->SetLineWidth(3);
0606 sf_eta_down->SetFillColorAlpha(kWhite, 1);
0607 sf_eta_down->Draw("][ same");
0608
0609 sf_eta_norm->SetLineStyle(1);
0610 sf_eta_norm->SetLineWidth(5);
0611 sf_eta_norm->SetFillColor(0);
0612 sf_eta_norm->Draw("][ same");
0613 sf_eta_norm->Draw("axis same");
0614
0615 leg_eta->AddEntry(sf_eta_norm, (tag_ver + '_' + tag_jet).c_str(), "l");
0616 leg_eta->AddEntry((TObject*)nullptr, Form("JetPt=%.2f", par_Pt), "");
0617 leg_eta->Draw();
0618
0619 if (is_2bin == true) {
0620 canvas.cd(2);
0621 sf_pt_up->SetTitle("ScaleFactor vs. p_{T}");
0622 sf_pt_up->SetXTitle("p_{T} [GeV]");
0623 sf_pt_up->SetYTitle("Scale Factor");
0624 sf_pt_up->SetLineStyle(7);
0625 sf_pt_up->SetLineWidth(3);
0626 sf_pt_up->SetFillColorAlpha(kGray, 0.5);
0627 sf_pt_up->SetMaximum(sf_pt_up->GetMaximum() * 1.25);
0628 sf_pt_up->SetMinimum(0.);
0629 sf_pt_up->Draw("][");
0630
0631 sf_pt_down->SetLineStyle(7);
0632 sf_pt_down->SetLineWidth(3);
0633 sf_pt_down->SetFillColorAlpha(kWhite, 1);
0634 sf_pt_down->Draw("][ same");
0635
0636 sf_pt_norm->SetLineStyle(1);
0637 sf_pt_norm->SetLineWidth(5);
0638 sf_pt_norm->SetFillColor(0);
0639 sf_pt_norm->Draw("][ same");
0640 sf_pt_norm->Draw("axis same");
0641
0642 leg_pt->AddEntry(sf_pt_norm, (tag_ver + '_' + tag_jet).c_str(), "l");
0643 leg_pt->AddEntry((TObject*)nullptr, Form("JetEta=%.2f", par_Eta), "");
0644 leg_pt->Draw();
0645 }
0646
0647 canvas.SaveAs(m_imageFileName.c_str());
0648
0649 return true;
0650 } else
0651 return false;
0652 }
0653
0654 };
0655
0656
0657 PAYLOAD_INSPECTOR_MODULE(JetResolutionObject) {
0658 PAYLOAD_INSPECTOR_CLASS(JetResolutionVsEta);
0659 PAYLOAD_INSPECTOR_CLASS(JetResolutionVsPt);
0660 PAYLOAD_INSPECTOR_CLASS(JetScaleFactorVsEtaNORM);
0661 PAYLOAD_INSPECTOR_CLASS(JetScaleFactorVsEtaDOWN);
0662 PAYLOAD_INSPECTOR_CLASS(JetScaleFactorVsEtaUP);
0663 PAYLOAD_INSPECTOR_CLASS(JetScaleFactorVsPtNORM);
0664 PAYLOAD_INSPECTOR_CLASS(JetScaleFactorVsPtDOWN);
0665 PAYLOAD_INSPECTOR_CLASS(JetScaleFactorVsPtUP);
0666 PAYLOAD_INSPECTOR_CLASS(JetResolutionSummary);
0667 PAYLOAD_INSPECTOR_CLASS(JetScaleFactorSummary);
0668 }
0669
0670 }