File indexing completed on 2024-07-02 00:53:27
0001 #include "CondCore/Utilities/interface/PayloadInspectorModule.h"
0002 #include "CondCore/Utilities/interface/PayloadInspector.h"
0003 #include "CondCore/CondDB/interface/Time.h"
0004 #include "DataFormats/HcalDetId/interface/HcalDetId.h"
0005 #include "Geometry/CaloTopology/interface/HcalTopology.h"
0006 #include "CondCore/HcalPlugins/interface/HcalObjRepresent.h"
0007
0008
0009 #include "CondFormats/HcalObjects/interface/HcalRespCorrs.h"
0010
0011 #include "TLegend.h"
0012 #include "TH2F.h"
0013 #include "TCanvas.h"
0014 #include "TLine.h"
0015 #include "TStyle.h"
0016 #include "TLatex.h"
0017 #include "TPave.h"
0018 #include "TPaveStats.h"
0019 #include <string>
0020 #include <fstream>
0021
0022 namespace {
0023
0024 using namespace cond::payloadInspector;
0025
0026 class HcalRespCorrContainer : public HcalObjRepresent::HcalDataContainer<HcalRespCorrs, HcalRespCorr> {
0027 public:
0028 HcalRespCorrContainer(std::shared_ptr<HcalRespCorrs> payload, unsigned int run)
0029 : HcalObjRepresent::HcalDataContainer<HcalRespCorrs, HcalRespCorr>(payload, run) {}
0030 float getValue(const HcalRespCorr* rCor) override { return rCor->getValue(); }
0031 };
0032
0033
0034
0035
0036 class HcalRespCorrsPlotAll : public cond::payloadInspector::PlotImage<HcalRespCorrs> {
0037 public:
0038 HcalRespCorrsPlotAll() : cond::payloadInspector::PlotImage<HcalRespCorrs>("HCAL RespCorr Ratios - map ") {
0039 setSingleIov(true);
0040 }
0041
0042 bool fill(const std::vector<std::tuple<cond::Time_t, cond::Hash>>& iovs) override {
0043 auto iov = iovs.front();
0044 std::shared_ptr<HcalRespCorrs> payload = fetchPayload(std::get<1>(iov));
0045 if (payload.get()) {
0046 HcalRespCorrContainer* objContainer = new HcalRespCorrContainer(payload, std::get<0>(iov));
0047 std::string ImageName(m_imageFileName);
0048 objContainer->getCanvasAll()->SaveAs(ImageName.c_str());
0049 return true;
0050 } else
0051 return false;
0052 }
0053 };
0054
0055
0056
0057
0058 class HcalRespCorrsRatioAll : public cond::payloadInspector::PlotImage<HcalRespCorrs> {
0059 public:
0060 HcalRespCorrsRatioAll() : cond::payloadInspector::PlotImage<HcalRespCorrs>("HCAL RespCorr Ratios difference") {
0061 setSingleIov(false);
0062 }
0063
0064 bool fill(const std::vector<std::tuple<cond::Time_t, cond::Hash>>& iovs) override {
0065 auto iov1 = iovs.front();
0066 auto iov2 = iovs.back();
0067
0068 std::shared_ptr<HcalRespCorrs> payload1 = fetchPayload(std::get<1>(iov1));
0069 std::shared_ptr<HcalRespCorrs> payload2 = fetchPayload(std::get<1>(iov2));
0070
0071 if (payload1.get() && payload2.get()) {
0072 HcalRespCorrContainer* objContainer1 = new HcalRespCorrContainer(payload1, std::get<0>(iov1));
0073 HcalRespCorrContainer* objContainer2 = new HcalRespCorrContainer(payload2, std::get<0>(iov2));
0074 objContainer2->Divide(objContainer1);
0075 std::string ImageName(m_imageFileName);
0076 objContainer2->getCanvasAll()->SaveAs(ImageName.c_str());
0077 return true;
0078 } else
0079 return false;
0080 }
0081 };
0082
0083
0084
0085
0086 template <int ntags, IOVMultiplicity nIOVs>
0087 class HcalRespCorrsComparatorBase : public cond::payloadInspector::PlotImage<HcalRespCorrs, nIOVs, ntags> {
0088 public:
0089 HcalRespCorrsComparatorBase()
0090 : cond::payloadInspector::PlotImage<HcalRespCorrs, nIOVs, ntags>("HcalRespCorrs Comparison") {}
0091
0092 using MetaData = std::tuple<cond::Time_t, cond::Hash>;
0093
0094 bool fill() override {
0095
0096 auto theIOVs = PlotBase::getTag<0>().iovs;
0097 auto tagname1 = PlotBase::getTag<0>().name;
0098 std::string tagname2 = "";
0099 auto firstiov = theIOVs.front();
0100 MetaData lastiov;
0101
0102
0103 assert(this->m_plotAnnotations.ntags < 3);
0104
0105 if (this->m_plotAnnotations.ntags == 2) {
0106 auto tag2iovs = PlotBase::getTag<1>().iovs;
0107 tagname2 = PlotBase::getTag<1>().name;
0108 lastiov = tag2iovs.front();
0109 } else {
0110 lastiov = theIOVs.back();
0111 }
0112
0113 std::shared_ptr<HcalRespCorrs> last_payload = this->fetchPayload(std::get<1>(lastiov));
0114 std::shared_ptr<HcalRespCorrs> first_payload = this->fetchPayload(std::get<1>(firstiov));
0115
0116 std::string lastIOVsince = std::to_string(std::get<0>(lastiov));
0117 std::string firstIOVsince = std::to_string(std::get<0>(firstiov));
0118
0119 HcalRespCorrContainer* last_objContainer = new HcalRespCorrContainer(last_payload, std::get<0>(lastiov));
0120 HcalRespCorrContainer* first_objContainer = new HcalRespCorrContainer(first_payload, std::get<0>(firstiov));
0121
0122 const auto& lastItems = last_objContainer->getAllItems();
0123 const auto& firstItems = first_objContainer->getAllItems();
0124
0125 assert(lastItems.size() == firstItems.size());
0126
0127 TCanvas canvas("Partition summary", "partition summary", 1400, 1000);
0128 canvas.Divide(3, 2);
0129
0130 std::map<std::string, std::shared_ptr<TH1F>> first_plots;
0131 std::map<std::string, std::shared_ptr<TH1F>> last_plots;
0132
0133
0134 std::vector<std::string> parts(lastItems.size());
0135
0136
0137 std::transform(lastItems.begin(),
0138 lastItems.end(),
0139 parts.begin(),
0140 [](const std::pair<std::string, std::vector<HcalRespCorr>>& pair) {
0141 return pair.first;
0142 });
0143
0144 auto legend = TLegend(0.13, 0.85, 0.98, 0.95);
0145 legend.SetTextSize(0.032);
0146
0147 unsigned int count{0};
0148 for (const auto& part : parts) {
0149 count++;
0150 first_plots[part] =
0151 std::make_shared<TH1F>(Form("f_%s_%s", part.c_str(), firstIOVsince.c_str()),
0152 Form("Response corrections [%s];correction factor;entries", part.c_str()),
0153 100,
0154 0.,
0155 3.);
0156
0157
0158 auto it = std::find_if(
0159 firstItems.begin(),
0160 firstItems.end(),
0161 [&part](const std::pair<std::string, std::vector<HcalRespCorr>>& pair) { return pair.first == part; });
0162
0163
0164 if (it != firstItems.end()) {
0165 const std::vector<HcalRespCorr>& result = it->second;
0166 if (DEBUG) {
0167 std::cout << "Found vector<HcalRespCorr> for key: " << part << std::endl;
0168 }
0169 for (auto& item : result) {
0170 HcalDetId detId = HcalDetId(item.rawId());
0171 if (DEBUG) {
0172 int iphi = detId.iphi();
0173 int ieta = detId.ieta();
0174 int depth = detId.depth();
0175
0176 std::cout << detId << " iphi" << iphi << " ieta: " << ieta << " depth:" << depth << std::endl;
0177 }
0178 if (detId != HcalDetId())
0179 first_plots[part]->Fill(first_objContainer->getValue(&item));
0180 }
0181
0182
0183 } else {
0184 std::cout << "Key not found: " << part << std::endl;
0185 }
0186
0187 last_plots[part] =
0188 std::make_shared<TH1F>(Form("l_%s_%s", part.c_str(), lastIOVsince.c_str()),
0189 Form("Response Corrections [%s];correction factor;entries", part.c_str()),
0190 100,
0191 0,
0192 3.);
0193
0194
0195 auto it2 = std::find_if(
0196 lastItems.begin(), lastItems.end(), [&part](const std::pair<std::string, std::vector<HcalRespCorr>>& pair) {
0197 return pair.first == part;
0198 });
0199
0200
0201 if (it2 != lastItems.end()) {
0202 const std::vector<HcalRespCorr>& result = it2->second;
0203 if (DEBUG) {
0204 std::cout << "Found vector<HcalRespCorr> for key: " << part << std::endl;
0205 }
0206 for (auto& item : result) {
0207 HcalDetId detId = HcalDetId(item.rawId());
0208 if (DEBUG) {
0209 int iphi = detId.iphi();
0210 int ieta = detId.ieta();
0211 int depth = detId.depth();
0212
0213 std::cout << detId << " iphi" << iphi << " ieta: " << ieta << " depth:" << depth << std::endl;
0214 }
0215 if (detId != HcalDetId())
0216 last_plots[part]->Fill(last_objContainer->getValue(&item));
0217 }
0218
0219
0220 } else {
0221 std::cout << "Key not found: " << part << std::endl;
0222 }
0223
0224 canvas.cd(count);
0225
0226 canvas.cd(count)->SetTopMargin(0.05);
0227 canvas.cd(count)->SetLeftMargin(0.13);
0228 canvas.cd(count)->SetRightMargin(0.02);
0229
0230 if (count == 1) {
0231 legend.AddEntry(first_plots[part].get(), (std::get<1>(firstiov)).c_str(), "L");
0232 legend.AddEntry(last_plots[part].get(), (std::get<1>(lastiov)).c_str(), "L");
0233 }
0234
0235 const auto& extrema = getExtrema(first_plots[part].get(), last_plots[part].get());
0236 first_plots[part]->SetMaximum(1.1 * extrema.second);
0237
0238 beautifyPlot(first_plots[part], kBlue);
0239 first_plots[part]->Draw();
0240 beautifyPlot(last_plots[part], kRed);
0241 last_plots[part]->Draw("same");
0242
0243
0244 addStatisticsPaveText(first_plots[part].get(), "Resp Corr", 0.80);
0245
0246
0247 addStatisticsPaveText(last_plots[part].get(), "Resp Corr", 0.70);
0248
0249 legend.Draw("same");
0250 }
0251
0252 std::string fileName(this->m_imageFileName);
0253 canvas.SaveAs(fileName.c_str());
0254
0255 return true;
0256 }
0257
0258 private:
0259 bool DEBUG{false};
0260
0261 void beautifyPlot(std::shared_ptr<TH1F> hist, int kColor) {
0262 hist->SetStats(kFALSE);
0263 hist->SetLineWidth(2);
0264 hist->SetLineColor(kColor);
0265 hist->GetXaxis()->CenterTitle(true);
0266 hist->GetYaxis()->CenterTitle(true);
0267 hist->GetXaxis()->SetTitleFont(42);
0268 hist->GetYaxis()->SetTitleFont(42);
0269 hist->GetXaxis()->SetTitleSize(0.05);
0270 hist->GetYaxis()->SetTitleSize(0.05);
0271 hist->GetXaxis()->SetTitleOffset(0.9);
0272 hist->GetYaxis()->SetTitleOffset(1.5);
0273 hist->GetXaxis()->SetLabelFont(42);
0274 hist->GetYaxis()->SetLabelFont(42);
0275 hist->GetYaxis()->SetLabelSize(.05);
0276 hist->GetXaxis()->SetLabelSize(.05);
0277 }
0278
0279 std::pair<float, float> getExtrema(TH1* h1, TH1* h2) {
0280 float theMax(-9999.);
0281 float theMin(9999.);
0282 theMax = h1->GetMaximum() > h2->GetMaximum() ? h1->GetMaximum() : h2->GetMaximum();
0283 theMin = h1->GetMinimum() < h2->GetMaximum() ? h1->GetMinimum() : h2->GetMinimum();
0284
0285 float add_min = theMin > 0. ? -0.05 : 0.05;
0286 float add_max = theMax > 0. ? 0.05 : -0.05;
0287
0288 auto result = std::make_pair(theMin * (1 + add_min), theMax * (1 + add_max));
0289 return result;
0290 }
0291
0292
0293 void addStatisticsPaveText(TH1* hist, const std::string& label, double yPosition) {
0294
0295 double mean = hist->GetMean();
0296 double rms = hist->GetRMS();
0297
0298
0299 std::ostringstream meanStream;
0300 meanStream << std::setprecision(4) << std::scientific << mean;
0301 std::string meanStr = meanStream.str();
0302
0303 std::ostringstream rmsStream;
0304 rmsStream << std::setprecision(4) << std::scientific << rms;
0305 std::string rmsStr = rmsStream.str();
0306
0307
0308 double x1 = 0.55;
0309 double x2 = 0.95;
0310 double y1 = yPosition - 0.1;
0311 double y2 = yPosition;
0312
0313
0314 TPaveText* statsPave = new TPaveText(x1, y1, x2, y2, "NDC");
0315 statsPave->SetFillColor(0);
0316 statsPave->SetBorderSize(1);
0317 statsPave->SetLineColor(hist->GetLineColor());
0318 statsPave->SetTextColor(hist->GetLineColor());
0319 statsPave->SetTextAlign(12);
0320
0321
0322 statsPave->AddText((label + " Mean: " + meanStr).c_str());
0323 statsPave->AddText((label + " RMS: " + rmsStr).c_str());
0324
0325
0326 statsPave->Draw();
0327 }
0328 };
0329
0330 using HcalRespCorrsComparatorSingleTag = HcalRespCorrsComparatorBase<1, MULTI_IOV>;
0331 using HcalRespCorrsComparatorTwoTags = HcalRespCorrsComparatorBase<2, SINGLE_IOV>;
0332
0333
0334
0335
0336 template <int ntags, IOVMultiplicity nIOVs>
0337 class HcalRespCorrsCorrelationBase : public cond::payloadInspector::PlotImage<HcalRespCorrs, nIOVs, ntags> {
0338 public:
0339 HcalRespCorrsCorrelationBase()
0340 : cond::payloadInspector::PlotImage<HcalRespCorrs, nIOVs, ntags>("HcalRespCorrs Comparison") {}
0341
0342 using MetaData = std::tuple<cond::Time_t, cond::Hash>;
0343
0344 bool fill() override {
0345
0346 auto theIOVs = PlotBase::getTag<0>().iovs;
0347 auto tagname1 = PlotBase::getTag<0>().name;
0348 std::string tagname2 = "";
0349 auto firstiov = theIOVs.front();
0350 MetaData lastiov;
0351
0352
0353 assert(this->m_plotAnnotations.ntags < 3);
0354
0355 if (this->m_plotAnnotations.ntags == 2) {
0356 auto tag2iovs = PlotBase::getTag<1>().iovs;
0357 tagname2 = PlotBase::getTag<1>().name;
0358 lastiov = tag2iovs.front();
0359 } else {
0360 lastiov = theIOVs.back();
0361 }
0362
0363 std::shared_ptr<HcalRespCorrs> last_payload = this->fetchPayload(std::get<1>(lastiov));
0364 std::shared_ptr<HcalRespCorrs> first_payload = this->fetchPayload(std::get<1>(firstiov));
0365
0366 std::string lastIOVsince = std::to_string(std::get<0>(lastiov));
0367 std::string firstIOVsince = std::to_string(std::get<0>(firstiov));
0368
0369 HcalRespCorrContainer* last_objContainer = new HcalRespCorrContainer(last_payload, std::get<0>(lastiov));
0370 HcalRespCorrContainer* first_objContainer = new HcalRespCorrContainer(first_payload, std::get<0>(firstiov));
0371
0372 const auto& lastItems = last_objContainer->getAllItems();
0373 const auto& firstItems = first_objContainer->getAllItems();
0374
0375 assert(lastItems.size() == firstItems.size());
0376
0377 TCanvas canvas("Partition summary", "partition summary", 1400, 1000);
0378 canvas.Divide(3, 2);
0379
0380 std::map<std::string, std::shared_ptr<TH2F>> plots;
0381
0382
0383 std::vector<std::string> parts(lastItems.size());
0384
0385
0386 std::transform(lastItems.begin(),
0387 lastItems.end(),
0388 parts.begin(),
0389 [](const std::pair<std::string, std::vector<HcalRespCorr>>& pair) {
0390 return pair.first;
0391 });
0392
0393 auto legend = TLegend(0.13, 0.85, 0.98, 0.95);
0394 legend.SetTextSize(0.032);
0395
0396 unsigned int count{0};
0397 for (const auto& part : parts) {
0398 count++;
0399 plots[part] =
0400 std::make_shared<TH2F>(Form("%s_%s_vs_%s", part.c_str(), firstIOVsince.c_str(), lastIOVsince.c_str()),
0401 Form("%s;correction factor (#bf{#color[2]{%s}}, IOV #bf{#color[2]{%s}});correction "
0402 "factor (#bf{#color[4]{%s}}, IOV #bf{#color[4]{%s}})",
0403 part.c_str(),
0404 tagname1.c_str(),
0405 firstIOVsince.c_str(),
0406 tagname2.c_str(),
0407 lastIOVsince.c_str()),
0408 100,
0409 0.,
0410 3.,
0411 100,
0412 0.,
0413 3.);
0414
0415
0416 auto it = std::find_if(
0417 firstItems.begin(),
0418 firstItems.end(),
0419 [&part](const std::pair<std::string, std::vector<HcalRespCorr>>& pair) { return pair.first == part; });
0420
0421
0422 auto it2 = std::find_if(
0423 lastItems.begin(), lastItems.end(), [&part](const std::pair<std::string, std::vector<HcalRespCorr>>& pair) {
0424 return pair.first == part;
0425 });
0426
0427
0428 if (it != firstItems.end() && it2 != lastItems.end()) {
0429 const std::vector<HcalRespCorr>& result = it->second;
0430 const std::vector<HcalRespCorr>& result2 = it2->second;
0431
0432 for (auto& item : result) {
0433 HcalDetId detId = HcalDetId(item.rawId());
0434 if (detId == HcalDetId()) {
0435 continue;
0436 }
0437 for (auto& item2 : result2) {
0438 HcalDetId detId2 = HcalDetId(item2.rawId());
0439 if (detId == detId2) {
0440 plots[part]->Fill(first_objContainer->getValue(&item), last_objContainer->getValue(&item2));
0441 }
0442 }
0443 }
0444 }
0445
0446 canvas.cd(count);
0447
0448 canvas.cd(count)->SetTopMargin(0.05);
0449 canvas.cd(count)->SetLeftMargin(0.13);
0450 canvas.cd(count)->SetRightMargin(0.02);
0451
0452 beautifyPlot(plots[part]);
0453 plots[part]->Draw();
0454 }
0455
0456 std::string fileName(this->m_imageFileName);
0457 canvas.SaveAs(fileName.c_str());
0458
0459 return true;
0460 }
0461
0462 void beautifyPlot(std::shared_ptr<TH2F> hist) {
0463 hist->SetStats(kFALSE);
0464 hist->GetXaxis()->CenterTitle(true);
0465 hist->GetYaxis()->CenterTitle(true);
0466 hist->GetXaxis()->SetTitleFont(42);
0467 hist->GetYaxis()->SetTitleFont(42);
0468 hist->GetXaxis()->SetTitleSize(0.035);
0469 hist->GetYaxis()->SetTitleSize(0.035);
0470 hist->GetXaxis()->SetTitleOffset(1.55);
0471 hist->GetYaxis()->SetTitleOffset(1.55);
0472 hist->GetXaxis()->SetLabelFont(42);
0473 hist->GetYaxis()->SetLabelFont(42);
0474 hist->GetYaxis()->SetLabelSize(.05);
0475 hist->GetXaxis()->SetLabelSize(.05);
0476 }
0477 };
0478
0479 using HcalRespCorrsCorrelationSingleTag = HcalRespCorrsCorrelationBase<1, MULTI_IOV>;
0480 using HcalRespCorrsCorrelationTwoTags = HcalRespCorrsCorrelationBase<2, SINGLE_IOV>;
0481
0482
0483
0484
0485 class HcalRespCorrsEtaPlotAll : public cond::payloadInspector::PlotImage<HcalRespCorrs> {
0486 public:
0487 HcalRespCorrsEtaPlotAll() : cond::payloadInspector::PlotImage<HcalRespCorrs>("HCAL RespCorr Ratios - map ") {
0488 setSingleIov(true);
0489 }
0490
0491 bool fill(const std::vector<std::tuple<cond::Time_t, cond::Hash>>& iovs) override {
0492 auto iov = iovs.front();
0493 std::shared_ptr<HcalRespCorrs> payload = fetchPayload(std::get<1>(iov));
0494 if (payload.get()) {
0495 HcalRespCorrContainer* objContainer = new HcalRespCorrContainer(payload, std::get<0>(iov));
0496 std::string ImageName(m_imageFileName);
0497 objContainer->getCanvasAll("EtaProfile")->SaveAs(ImageName.c_str());
0498 return true;
0499 } else
0500 return false;
0501 }
0502 };
0503
0504
0505
0506
0507 class HcalRespCorrsEtaRatioAll : public cond::payloadInspector::PlotImage<HcalRespCorrs> {
0508 public:
0509 HcalRespCorrsEtaRatioAll() : cond::payloadInspector::PlotImage<HcalRespCorrs>("HCAL RespCorr Ratios difference") {
0510 setSingleIov(false);
0511 }
0512
0513 bool fill(const std::vector<std::tuple<cond::Time_t, cond::Hash>>& iovs) override {
0514 auto iov1 = iovs.front();
0515 auto iov2 = iovs.back();
0516
0517 std::shared_ptr<HcalRespCorrs> payload1 = fetchPayload(std::get<1>(iov1));
0518 std::shared_ptr<HcalRespCorrs> payload2 = fetchPayload(std::get<1>(iov2));
0519
0520 if (payload1.get() && payload2.get()) {
0521 HcalRespCorrContainer* objContainer1 = new HcalRespCorrContainer(payload1, std::get<0>(iov1));
0522 HcalRespCorrContainer* objContainer2 = new HcalRespCorrContainer(payload2, std::get<0>(iov2));
0523 objContainer2->Divide(objContainer1);
0524 std::string ImageName(m_imageFileName);
0525 objContainer2->getCanvasAll("EtaProfile")->SaveAs(ImageName.c_str());
0526 return true;
0527 } else
0528 return false;
0529 }
0530 };
0531
0532
0533
0534 class HcalRespCorrsPhiPlotAll : public cond::payloadInspector::PlotImage<HcalRespCorrs> {
0535 public:
0536 HcalRespCorrsPhiPlotAll() : cond::payloadInspector::PlotImage<HcalRespCorrs>("HCAL RespCorr Ratios - map ") {
0537 setSingleIov(true);
0538 }
0539
0540 bool fill(const std::vector<std::tuple<cond::Time_t, cond::Hash>>& iovs) override {
0541 auto iov = iovs.front();
0542 std::shared_ptr<HcalRespCorrs> payload = fetchPayload(std::get<1>(iov));
0543 if (payload.get()) {
0544 HcalRespCorrContainer* objContainer = new HcalRespCorrContainer(payload, std::get<0>(iov));
0545 std::string ImageName(m_imageFileName);
0546 objContainer->getCanvasAll("PhiProfile")->SaveAs(ImageName.c_str());
0547 return true;
0548 } else
0549 return false;
0550 }
0551 };
0552
0553
0554
0555
0556 class HcalRespCorrsPhiRatioAll : public cond::payloadInspector::PlotImage<HcalRespCorrs> {
0557 public:
0558 HcalRespCorrsPhiRatioAll() : cond::payloadInspector::PlotImage<HcalRespCorrs>("HCAL RespCorr Ratios difference") {
0559 setSingleIov(false);
0560 }
0561
0562 bool fill(const std::vector<std::tuple<cond::Time_t, cond::Hash>>& iovs) override {
0563 auto iov1 = iovs.front();
0564 auto iov2 = iovs.back();
0565
0566 std::shared_ptr<HcalRespCorrs> payload1 = fetchPayload(std::get<1>(iov1));
0567 std::shared_ptr<HcalRespCorrs> payload2 = fetchPayload(std::get<1>(iov2));
0568
0569 if (payload1.get() && payload2.get()) {
0570 HcalRespCorrContainer* objContainer1 = new HcalRespCorrContainer(payload1, std::get<0>(iov1));
0571 HcalRespCorrContainer* objContainer2 = new HcalRespCorrContainer(payload2, std::get<0>(iov2));
0572 objContainer2->Divide(objContainer1);
0573 std::string ImageName(m_imageFileName);
0574 objContainer2->getCanvasAll("PhiProfile")->SaveAs(ImageName.c_str());
0575 return true;
0576 } else
0577 return false;
0578 }
0579 };
0580
0581
0582
0583 class HcalRespCorrsPlotHBHO : public cond::payloadInspector::PlotImage<HcalRespCorrs> {
0584 public:
0585 HcalRespCorrsPlotHBHO() : cond::payloadInspector::PlotImage<HcalRespCorrs>("HCAL RespCorr Ratios - map ") {
0586 setSingleIov(true);
0587 }
0588
0589 bool fill(const std::vector<std::tuple<cond::Time_t, cond::Hash>>& iovs) override {
0590 auto iov = iovs.front();
0591 std::shared_ptr<HcalRespCorrs> payload = fetchPayload(std::get<1>(iov));
0592 if (payload.get()) {
0593 HcalRespCorrContainer* objContainer = new HcalRespCorrContainer(payload, std::get<0>(iov));
0594 std::string ImageName(m_imageFileName);
0595 objContainer->getCanvasHBHO()->SaveAs(ImageName.c_str());
0596 return true;
0597 } else
0598 return false;
0599 }
0600 };
0601
0602
0603
0604
0605 class HcalRespCorrsRatioHBHO : public cond::payloadInspector::PlotImage<HcalRespCorrs> {
0606 public:
0607 HcalRespCorrsRatioHBHO() : cond::payloadInspector::PlotImage<HcalRespCorrs>("HCAL RespCorr Ratios difference") {
0608 setSingleIov(false);
0609 }
0610
0611 bool fill(const std::vector<std::tuple<cond::Time_t, cond::Hash>>& iovs) override {
0612 auto iov1 = iovs.front();
0613 auto iov2 = iovs.back();
0614
0615 std::shared_ptr<HcalRespCorrs> payload1 = fetchPayload(std::get<1>(iov1));
0616 std::shared_ptr<HcalRespCorrs> payload2 = fetchPayload(std::get<1>(iov2));
0617
0618 if (payload1.get() && payload2.get()) {
0619 HcalRespCorrContainer* objContainer1 = new HcalRespCorrContainer(payload1, std::get<0>(iov1));
0620 HcalRespCorrContainer* objContainer2 = new HcalRespCorrContainer(payload2, std::get<0>(iov2));
0621 objContainer2->Divide(objContainer1);
0622 std::string ImageName(m_imageFileName);
0623 objContainer2->getCanvasHBHO()->SaveAs(ImageName.c_str());
0624 return true;
0625 } else
0626 return false;
0627 }
0628 };
0629
0630
0631
0632 class HcalRespCorrsPlotHE : public cond::payloadInspector::PlotImage<HcalRespCorrs> {
0633 public:
0634 HcalRespCorrsPlotHE() : cond::payloadInspector::PlotImage<HcalRespCorrs>("HCAL RespCorr Ratios - map ") {
0635 setSingleIov(true);
0636 }
0637
0638 bool fill(const std::vector<std::tuple<cond::Time_t, cond::Hash>>& iovs) override {
0639 auto iov = iovs.front();
0640 std::shared_ptr<HcalRespCorrs> payload = fetchPayload(std::get<1>(iov));
0641 if (payload.get()) {
0642 HcalRespCorrContainer* objContainer = new HcalRespCorrContainer(payload, std::get<0>(iov));
0643 std::string ImageName(m_imageFileName);
0644 objContainer->getCanvasHE()->SaveAs(ImageName.c_str());
0645 return true;
0646 } else
0647 return false;
0648 }
0649 };
0650
0651
0652
0653
0654 class HcalRespCorrsRatioHE : public cond::payloadInspector::PlotImage<HcalRespCorrs> {
0655 public:
0656 HcalRespCorrsRatioHE() : cond::payloadInspector::PlotImage<HcalRespCorrs>("HCAL RespCorr Ratios difference") {
0657 setSingleIov(false);
0658 }
0659
0660 bool fill(const std::vector<std::tuple<cond::Time_t, cond::Hash>>& iovs) override {
0661 auto iov1 = iovs.front();
0662 auto iov2 = iovs.back();
0663
0664 std::shared_ptr<HcalRespCorrs> payload1 = fetchPayload(std::get<1>(iov1));
0665 std::shared_ptr<HcalRespCorrs> payload2 = fetchPayload(std::get<1>(iov2));
0666
0667 if (payload1.get() && payload2.get()) {
0668 HcalRespCorrContainer* objContainer1 = new HcalRespCorrContainer(payload1, std::get<0>(iov1));
0669 HcalRespCorrContainer* objContainer2 = new HcalRespCorrContainer(payload2, std::get<0>(iov2));
0670 objContainer2->Divide(objContainer1);
0671 std::string ImageName(m_imageFileName);
0672 objContainer2->getCanvasHE()->SaveAs(ImageName.c_str());
0673 return true;
0674 } else
0675 return false;
0676 }
0677 };
0678
0679
0680
0681 class HcalRespCorrsPlotHF : public cond::payloadInspector::PlotImage<HcalRespCorrs> {
0682 public:
0683 HcalRespCorrsPlotHF() : cond::payloadInspector::PlotImage<HcalRespCorrs>("HCAL RespCorr Ratios - map ") {
0684 setSingleIov(true);
0685 }
0686
0687 bool fill(const std::vector<std::tuple<cond::Time_t, cond::Hash>>& iovs) override {
0688 auto iov = iovs.front();
0689 std::shared_ptr<HcalRespCorrs> payload = fetchPayload(std::get<1>(iov));
0690 if (payload.get()) {
0691 HcalRespCorrContainer* objContainer = new HcalRespCorrContainer(payload, std::get<0>(iov));
0692 std::string ImageName(m_imageFileName);
0693 objContainer->getCanvasHF()->SaveAs(ImageName.c_str());
0694 return true;
0695 } else
0696 return false;
0697 }
0698 };
0699
0700
0701
0702
0703 class HcalRespCorrsRatioHF : public cond::payloadInspector::PlotImage<HcalRespCorrs> {
0704 public:
0705 HcalRespCorrsRatioHF() : cond::payloadInspector::PlotImage<HcalRespCorrs>("HCAL RespCorr Ratios difference") {
0706 setSingleIov(false);
0707 }
0708
0709 bool fill(const std::vector<std::tuple<cond::Time_t, cond::Hash>>& iovs) override {
0710 auto iov1 = iovs.front();
0711 auto iov2 = iovs.back();
0712
0713 std::shared_ptr<HcalRespCorrs> payload1 = fetchPayload(std::get<1>(iov1));
0714 std::shared_ptr<HcalRespCorrs> payload2 = fetchPayload(std::get<1>(iov2));
0715
0716 if (payload1.get() && payload2.get()) {
0717 HcalRespCorrContainer* objContainer1 = new HcalRespCorrContainer(payload1, std::get<0>(iov1));
0718 HcalRespCorrContainer* objContainer2 = new HcalRespCorrContainer(payload2, std::get<0>(iov2));
0719 objContainer2->Divide(objContainer1);
0720 std::string ImageName(m_imageFileName);
0721 objContainer2->getCanvasHF()->SaveAs(ImageName.c_str());
0722 return true;
0723 } else
0724 return false;
0725 }
0726 };
0727 }
0728
0729
0730 PAYLOAD_INSPECTOR_MODULE(HcalRespCorrs) {
0731 PAYLOAD_INSPECTOR_CLASS(HcalRespCorrsPlotAll);
0732 PAYLOAD_INSPECTOR_CLASS(HcalRespCorrsRatioAll);
0733 PAYLOAD_INSPECTOR_CLASS(HcalRespCorrsComparatorSingleTag);
0734 PAYLOAD_INSPECTOR_CLASS(HcalRespCorrsComparatorTwoTags);
0735 PAYLOAD_INSPECTOR_CLASS(HcalRespCorrsCorrelationSingleTag);
0736 PAYLOAD_INSPECTOR_CLASS(HcalRespCorrsCorrelationTwoTags);
0737 PAYLOAD_INSPECTOR_CLASS(HcalRespCorrsEtaPlotAll);
0738 PAYLOAD_INSPECTOR_CLASS(HcalRespCorrsEtaRatioAll);
0739 PAYLOAD_INSPECTOR_CLASS(HcalRespCorrsPhiPlotAll);
0740 PAYLOAD_INSPECTOR_CLASS(HcalRespCorrsPhiRatioAll);
0741 PAYLOAD_INSPECTOR_CLASS(HcalRespCorrsPlotHBHO);
0742 PAYLOAD_INSPECTOR_CLASS(HcalRespCorrsRatioHBHO);
0743 PAYLOAD_INSPECTOR_CLASS(HcalRespCorrsPlotHE);
0744 PAYLOAD_INSPECTOR_CLASS(HcalRespCorrsRatioHE);
0745 PAYLOAD_INSPECTOR_CLASS(HcalRespCorrsPlotHF);
0746 PAYLOAD_INSPECTOR_CLASS(HcalRespCorrsRatioHF);
0747 }