File indexing completed on 2024-09-07 04:35:27
0001 #include "CondCore/Utilities/interface/PayloadInspectorModule.h"
0002 #include "CondCore/Utilities/interface/PayloadInspector.h"
0003 #include "DataFormats/EcalDetId/interface/EBDetId.h"
0004 #include "DataFormats/EcalDetId/interface/EEDetId.h"
0005 #include "CondCore/EcalPlugins/plugins/EcalDrawUtils.h"
0006
0007
0008 #include "CondFormats/EcalObjects/interface/EcalLinearCorrections.h"
0009
0010 #include "TH2F.h"
0011 #include "TCanvas.h"
0012 #include "TLine.h"
0013 #include "TStyle.h"
0014 #include "TLatex.h"
0015
0016 #include <string>
0017
0018 namespace {
0019 constexpr int kEBChannels = 61200, kEEChannels = 14648, kSides = 2, kValues = 3, kRMS = 5;
0020 constexpr int MAX_IETA = 85, MAX_IPHI = 360;
0021 constexpr int IX_MIN = 1, IY_MIN = 1, IX_MAX = 100, IY_MAX = 100;
0022
0023
0024
0025
0026 class EcalLinearCorrectionsPlot : public cond::payloadInspector::PlotImage<EcalLinearCorrections> {
0027 public:
0028 EcalLinearCorrectionsPlot()
0029 : cond::payloadInspector::PlotImage<EcalLinearCorrections>("ECAL LinearCorrections - map ") {
0030 setSingleIov(true);
0031 }
0032
0033 bool fill(const std::vector<std::tuple<cond::Time_t, cond::Hash> >& iovs) override {
0034 TH2F** barrel = new TH2F*[kValues];
0035 TH2F** endc_p = new TH2F*[kValues];
0036 TH2F** endc_m = new TH2F*[kValues];
0037 double EBmean[kValues], EBrms[kValues], EEmean[kValues], EErms[kValues], pEBmin[kValues], pEBmax[kValues],
0038 pEEmin[kValues], pEEmax[kValues];
0039 int EBtot[kValues], EEtot[kValues];
0040 for (int valId = 0; valId < kValues; valId++) {
0041 barrel[valId] = new TH2F(
0042 Form("EBp%i", valId), Form("EBp%i", valId), MAX_IPHI, 0, MAX_IPHI, 2 * MAX_IETA, -MAX_IETA, MAX_IETA);
0043 endc_p[valId] = new TH2F(
0044 Form("EE+p%i", valId), Form("EE+p%i", valId), IX_MAX, IX_MIN, IX_MAX + 1, IY_MAX, IY_MIN, IY_MAX + 1);
0045 endc_m[valId] = new TH2F(
0046 Form("EE-p%i", valId), Form("EE-p%i", valId), IX_MAX, IX_MIN, IX_MAX + 1, IY_MAX, IY_MIN, IY_MAX + 1);
0047 EBmean[valId] = 0.;
0048 EBrms[valId] = 0.;
0049 EEmean[valId] = 0.;
0050 EErms[valId] = 0.;
0051 EBtot[valId] = 0;
0052 EEtot[valId] = 0;
0053 }
0054
0055 auto iov = iovs.front();
0056 std::shared_ptr<EcalLinearCorrections> payload = fetchPayload(std::get<1>(iov));
0057 unsigned long IOV = std::get<0>(iov);
0058 int run = 0;
0059 if (IOV < 4294967296)
0060 run = std::get<0>(iov);
0061 else {
0062 run = IOV >> 32;
0063
0064
0065 }
0066 if (payload.get()) {
0067 for (int ieta = -MAX_IETA; ieta <= MAX_IETA; ieta++) {
0068 Double_t eta = (Double_t)ieta;
0069 if (ieta == 0)
0070 continue;
0071 else if (ieta > 0)
0072 eta = eta - 0.5;
0073 else
0074 eta = eta + 0.5;
0075 for (int iphi = 1; iphi <= MAX_IPHI; iphi++) {
0076 Double_t phi = (Double_t)iphi - 0.5;
0077 EBDetId id(ieta, iphi);
0078 Double_t val = (*payload).getValueMap()[id.rawId()].p1;
0079 barrel[0]->Fill(phi, eta, val);
0080 EBmean[0] = EBmean[0] + val;
0081 EBrms[0] = EBrms[0] + val * val;
0082 EBtot[0]++;
0083 val = (*payload).getValueMap()[id.rawId()].p2;
0084 barrel[1]->Fill(phi, eta, val);
0085 EBmean[1] = EBmean[1] + val;
0086 EBrms[1] = EBrms[1] + val * val;
0087 EBtot[1]++;
0088 val = (*payload).getValueMap()[id.rawId()].p3;
0089 barrel[2]->Fill(phi, eta, val);
0090 EBmean[2] = EBmean[2] + val;
0091 EBrms[2] = EBrms[2] + val * val;
0092 EBtot[2]++;
0093 }
0094 }
0095
0096 for (int sign = 0; sign < kSides; sign++) {
0097 int thesign = sign == 1 ? 1 : -1;
0098 for (int ix = 1; ix <= IX_MAX; ix++) {
0099 for (int iy = 1; iy <= IY_MAX; iy++) {
0100 if (!EEDetId::validDetId(ix, iy, thesign))
0101 continue;
0102 EEDetId id(ix, iy, thesign);
0103 Double_t val = (*payload).getValueMap()[id.rawId()].p1;
0104 EEmean[0] = EEmean[0] + val;
0105 EErms[0] = EErms[0] + val * val;
0106 EEtot[0]++;
0107 if (thesign == 1)
0108 endc_p[0]->Fill(ix, iy, val);
0109 else
0110 endc_m[0]->Fill(ix, iy, val);
0111 val = (*payload).getValueMap()[id.rawId()].p1;
0112 EEmean[1] = EEmean[1] + val;
0113 EErms[1] = EErms[1] + val * val;
0114 EEtot[1]++;
0115 if (thesign == 1)
0116 endc_p[1]->Fill(ix, iy, val);
0117 else
0118 endc_m[1]->Fill(ix, iy, val);
0119 val = (*payload).getValueMap()[id.rawId()].p2;
0120 EEmean[2] = EEmean[2] + val;
0121 EErms[2] = EErms[2] + val * val;
0122 EEtot[2]++;
0123 if (thesign == 1)
0124 endc_p[2]->Fill(ix, iy, val);
0125 else
0126 endc_m[2]->Fill(ix, iy, val);
0127 }
0128 }
0129 }
0130 }
0131 else
0132 return false;
0133
0134 gStyle->SetPalette(1);
0135 gStyle->SetOptStat(0);
0136 TCanvas canvas("CC map", "CC map", 2800, 2600);
0137 TLatex t1;
0138 t1.SetNDC();
0139 t1.SetTextAlign(26);
0140 t1.SetTextSize(0.05);
0141
0142 if (IOV < 4294967296)
0143 t1.DrawLatex(0.5, 0.96, Form("Ecal Linear Corrections, IOV %i", run));
0144 else {
0145 time_t t = run;
0146 char buf[256];
0147 struct tm lt;
0148 localtime_r(&t, <);
0149 strftime(buf, sizeof(buf), "%F %R:%S", <);
0150 buf[sizeof(buf) - 1] = 0;
0151 t1.DrawLatex(0.5, 0.96, Form("Ecal Linear Corrections, IOV %s", buf));
0152 }
0153
0154 float xmi[3] = {0.0, 0.26, 0.74};
0155 float xma[3] = {0.26, 0.74, 1.00};
0156 TPad*** pad = new TPad**[3];
0157
0158 for (int valId = 0; valId < kValues; valId++) {
0159 pad[valId] = new TPad*[3];
0160 for (int obj = 0; obj < 3; obj++) {
0161 float yma = 0.94 - (0.32 * valId);
0162 float ymi = yma - 0.28;
0163 pad[valId][obj] =
0164 new TPad(Form("p_%i_%i", obj, valId), Form("p_%i_%i", obj, valId), xmi[obj], ymi, xma[obj], yma);
0165 pad[valId][obj]->Draw();
0166 }
0167 double vt = (double)EBtot[valId];
0168 EBmean[valId] = EBmean[valId] / vt;
0169 EBrms[valId] = (EBrms[valId] / vt) - (EBmean[valId] * EBmean[valId]);
0170 EBrms[valId] = sqrt(EBrms[valId]);
0171 if (EBrms[valId] == 0.)
0172 EBrms[valId] = 0.001;
0173 pEBmin[valId] = EBmean[valId] - kRMS * EBrms[valId];
0174 pEBmax[valId] = EBmean[valId] + kRMS * EBrms[valId];
0175
0176
0177 vt = (double)EEtot[valId];
0178 EEmean[valId] = EEmean[valId] / vt;
0179 EErms[valId] = (EErms[valId] / vt) - (EEmean[valId] * EEmean[valId]);
0180 EErms[valId] = sqrt(EErms[valId]);
0181 if (EErms[valId] == 0.)
0182 EErms[valId] = 0.001;
0183 pEEmin[valId] = EEmean[valId] - kRMS * EErms[valId];
0184 pEEmax[valId] = EEmean[valId] + kRMS * EErms[valId];
0185
0186
0187 }
0188
0189 for (int valId = 0; valId < 3; valId++) {
0190 pad[valId][0]->cd();
0191 DrawEE(endc_m[valId], pEEmin[valId], pEEmax[valId]);
0192 pad[valId][1]->cd();
0193 DrawEB(barrel[valId], pEBmin[valId], pEBmax[valId]);
0194 barrel[valId]->SetStats(false);
0195 pad[valId][2]->cd();
0196 DrawEE(endc_p[valId], pEEmin[valId], pEEmax[valId]);
0197 }
0198
0199 std::string ImageName(m_imageFileName);
0200 canvas.SaveAs(ImageName.c_str());
0201 return true;
0202 }
0203 };
0204
0205
0206
0207
0208 template <cond::payloadInspector::IOVMultiplicity nIOVs, int ntags>
0209 class EcalLinearCorrectionsDiffBase : public cond::payloadInspector::PlotImage<EcalLinearCorrections, nIOVs, ntags> {
0210 public:
0211 EcalLinearCorrectionsDiffBase()
0212 : cond::payloadInspector::PlotImage<EcalLinearCorrections, nIOVs, ntags>("ECAL LinearCorrections - map ") {}
0213
0214 bool fill() override {
0215 TH2F** barrel = new TH2F*[kValues];
0216 TH2F** endc_p = new TH2F*[kValues];
0217 TH2F** endc_m = new TH2F*[kValues];
0218 double EBmean[kValues], EBrms[kValues], EEmean[kValues], EErms[kValues], pEBmin[kValues], pEBmax[kValues],
0219 pEEmin[kValues], pEEmax[kValues];
0220 int EBtot[kValues], EEtot[kValues];
0221 for (int valId = 0; valId < kValues; valId++) {
0222 barrel[valId] = new TH2F(
0223 Form("EBp%i", valId), Form("EBp%i", valId), MAX_IPHI, 0, MAX_IPHI, 2 * MAX_IETA, -MAX_IETA, MAX_IETA);
0224 endc_p[valId] = new TH2F(
0225 Form("EE+p%i", valId), Form("EE+p%i", valId), IX_MAX, IX_MIN, IX_MAX + 1, IY_MAX, IY_MIN, IY_MAX + 1);
0226 endc_m[valId] = new TH2F(
0227 Form("EE-p%i", valId), Form("EE-p%i", valId), IX_MAX, IX_MIN, IX_MAX + 1, IY_MAX, IY_MIN, IY_MAX + 1);
0228 EBmean[valId] = 0.;
0229 EBrms[valId] = 0.;
0230 EEmean[valId] = 0.;
0231 EErms[valId] = 0.;
0232 EBtot[valId] = 0;
0233 EEtot[valId] = 0;
0234 }
0235 float vEB[kValues][kEBChannels], vEE[kValues][kEEChannels];
0236
0237 unsigned int run[2] = {0, 0};
0238 unsigned long IOV = 0;
0239 std::string l_tagname[2];
0240 auto iovs = cond::payloadInspector::PlotBase::getTag<0>().iovs;
0241 l_tagname[0] = cond::payloadInspector::PlotBase::getTag<0>().name;
0242 auto firstiov = iovs.front();
0243 IOV = std::get<0>(firstiov);
0244 if (IOV < 4294967296)
0245 run[0] = IOV;
0246 else {
0247 run[0] = IOV >> 32;
0248 }
0249 std::tuple<cond::Time_t, cond::Hash> lastiov;
0250 if (ntags == 2) {
0251 auto tag2iovs = cond::payloadInspector::PlotBase::getTag<1>().iovs;
0252 l_tagname[1] = cond::payloadInspector::PlotBase::getTag<1>().name;
0253 lastiov = tag2iovs.front();
0254 } else {
0255 lastiov = iovs.back();
0256 l_tagname[1] = l_tagname[0];
0257 }
0258 IOV = std::get<0>(lastiov);
0259 if (IOV < 4294967296)
0260 run[1] = IOV;
0261 else {
0262 run[1] = IOV >> 32;
0263 }
0264 for (int irun = 0; irun < nIOVs; irun++) {
0265 std::shared_ptr<EcalLinearCorrections> payload;
0266 if (irun == 0) {
0267 payload = this->fetchPayload(std::get<1>(firstiov));
0268 } else {
0269 payload = this->fetchPayload(std::get<1>(lastiov));
0270 }
0271 if (payload.get()) {
0272 for (int ieta = -MAX_IETA; ieta <= MAX_IETA; ieta++) {
0273 Double_t eta = (Double_t)ieta;
0274 if (ieta == 0)
0275 continue;
0276 else if (ieta > 0)
0277 eta = eta - 0.5;
0278 else
0279 eta = eta + 0.5;
0280 for (int iphi = 1; iphi <= MAX_IPHI; iphi++) {
0281 Double_t phi = (Double_t)iphi - 0.5;
0282 EBDetId id(ieta, iphi);
0283 Double_t val = (*payload).getValueMap()[id.rawId()].p1;
0284 int channel = id.hashedIndex();
0285 if (irun == 0)
0286 vEB[0][channel] = val;
0287 else {
0288 double diff = val - vEB[0][channel];
0289 barrel[0]->Fill(phi, eta, diff);
0290 EBmean[0] = EBmean[0] + diff;
0291 EBrms[0] = EBrms[0] + diff * diff;
0292 EBtot[0]++;
0293 }
0294 val = (*payload).getValueMap()[id.rawId()].p2;
0295 if (irun == 0)
0296 vEB[1][channel] = val;
0297 else {
0298 double diff = val - vEB[1][channel];
0299 barrel[1]->Fill(phi, eta, diff);
0300 EBmean[1] = EBmean[1] + diff;
0301 EBrms[1] = EBrms[1] + diff * diff;
0302 EBtot[1]++;
0303 }
0304 val = (*payload).getValueMap()[id.rawId()].p3;
0305 if (irun == 0)
0306 vEB[2][channel] = val;
0307 else {
0308 double diff = val - vEB[2][channel];
0309 barrel[2]->Fill(phi, eta, diff);
0310 EBmean[2] = EBmean[2] + diff;
0311 EBrms[2] = EBrms[2] + diff * diff;
0312 EBtot[2]++;
0313 }
0314 }
0315 }
0316
0317 for (int sign = 0; sign < kSides; sign++) {
0318 int thesign = sign == 1 ? 1 : -1;
0319 for (int ix = 1; ix <= IX_MAX; ix++) {
0320 for (int iy = 1; iy <= IY_MAX; iy++) {
0321 if (!EEDetId::validDetId(ix, iy, thesign))
0322 continue;
0323 EEDetId id(ix, iy, thesign);
0324 Double_t val = (*payload).getValueMap()[id.rawId()].p1;
0325 int channel = id.hashedIndex();
0326 if (irun == 0)
0327 vEE[0][channel] = val;
0328 else {
0329 double diff = val - vEE[0][channel];
0330 EEmean[0] = EEmean[0] + diff;
0331 EErms[0] = EErms[0] + diff * diff;
0332 EEtot[0]++;
0333 if (thesign == 1)
0334 endc_p[0]->Fill(ix, iy, diff);
0335 else
0336 endc_m[0]->Fill(ix, iy, diff);
0337 }
0338 val = (*payload).getValueMap()[id.rawId()].p1;
0339 if (irun == 0)
0340 vEE[1][channel] = val;
0341 else {
0342 double diff = val - vEE[1][channel];
0343 EEmean[1] = EEmean[1] + diff;
0344 EErms[1] = EErms[1] + diff * diff;
0345 EEtot[1]++;
0346 if (thesign == 1)
0347 endc_p[1]->Fill(ix, iy, diff);
0348 else
0349 endc_m[1]->Fill(ix, iy, diff);
0350 }
0351 val = (*payload).getValueMap()[id.rawId()].p2;
0352 if (irun == 0)
0353 vEE[2][channel] = val;
0354 else {
0355 double diff = val - vEE[2][channel];
0356 EEmean[2] = EEmean[2] + diff;
0357 EErms[2] = EErms[2] + diff * diff;
0358 EEtot[2]++;
0359 if (thesign == 1)
0360 endc_p[2]->Fill(ix, iy, diff);
0361 else
0362 endc_m[2]->Fill(ix, iy, diff);
0363 }
0364 }
0365 }
0366 }
0367 }
0368 else
0369 return false;
0370 }
0371
0372 gStyle->SetPalette(1);
0373 gStyle->SetOptStat(0);
0374 TCanvas canvas("CC map", "CC map", 2800, 2600);
0375 TLatex t1;
0376 t1.SetNDC();
0377 t1.SetTextAlign(26);
0378
0379 int len = l_tagname[0].length() + l_tagname[1].length();
0380 if (IOV < 4294967296) {
0381 t1.SetTextSize(0.05);
0382 t1.DrawLatex(0.5, 0.96, Form("Ecal Linear Corrections, IOV %i - %i", run[1], run[0]));
0383 } else {
0384 time_t t = run[0];
0385 char buf0[256], buf1[256];
0386 struct tm lt;
0387 localtime_r(&t, <);
0388 strftime(buf0, sizeof(buf0), "%F %R:%S", <);
0389 buf0[sizeof(buf0) - 1] = 0;
0390 t = run[1];
0391 localtime_r(&t, <);
0392 strftime(buf1, sizeof(buf1), "%F %R:%S", <);
0393 buf1[sizeof(buf1) - 1] = 0;
0394 if (ntags == 2) {
0395 if (len < 80) {
0396 t1.SetTextSize(0.02);
0397 t1.DrawLatex(0.5, 0.96, Form("%s %s - %s %s", l_tagname[1].c_str(), buf1, l_tagname[0].c_str(), buf0));
0398 } else {
0399 t1.SetTextSize(0.03);
0400 t1.DrawLatex(0.5, 0.96, Form("Ecal Linear Corrections, IOV %s - %s", buf1, buf0));
0401 }
0402 } else {
0403 t1.SetTextSize(0.015);
0404 t1.DrawLatex(0.5, 0.96, Form("%s, IOV %s - %s", l_tagname[0].c_str(), buf1, buf0));
0405 }
0406 }
0407
0408 float xmi[3] = {0.0, 0.26, 0.74};
0409 float xma[3] = {0.26, 0.74, 1.00};
0410 TPad*** pad = new TPad**[3];
0411
0412 for (int valId = 0; valId < kValues; valId++) {
0413 pad[valId] = new TPad*[3];
0414 for (int obj = 0; obj < 3; obj++) {
0415 float yma = 0.94 - (0.32 * valId);
0416 float ymi = yma - 0.28;
0417 pad[valId][obj] =
0418 new TPad(Form("p_%i_%i", obj, valId), Form("p_%i_%i", obj, valId), xmi[obj], ymi, xma[obj], yma);
0419 pad[valId][obj]->Draw();
0420 }
0421 double vt = (double)EBtot[valId];
0422 EBmean[valId] = EBmean[valId] / vt;
0423 EBrms[valId] = (EBrms[valId] / vt) - (EBmean[valId] * EBmean[valId]);
0424 EBrms[valId] = sqrt(EBrms[valId]);
0425 if (EBrms[valId] == 0.)
0426 EBrms[valId] = 0.001;
0427 pEBmin[valId] = EBmean[valId] - kRMS * EBrms[valId];
0428 pEBmax[valId] = EBmean[valId] + kRMS * EBrms[valId];
0429
0430
0431 vt = (double)EEtot[valId];
0432 EEmean[valId] = EEmean[valId] / vt;
0433 EErms[valId] = (EErms[valId] / vt) - (EEmean[valId] * EEmean[valId]);
0434 EErms[valId] = sqrt(EErms[valId]);
0435 if (EErms[valId] == 0.)
0436 EErms[valId] = 0.001;
0437 pEEmin[valId] = EEmean[valId] - kRMS * EErms[valId];
0438 pEEmax[valId] = EEmean[valId] + kRMS * EErms[valId];
0439
0440
0441 }
0442
0443 for (int valId = 0; valId < 3; valId++) {
0444 pad[valId][0]->cd();
0445 DrawEE(endc_m[valId], pEEmin[valId], pEEmax[valId]);
0446 pad[valId][1]->cd();
0447 DrawEB(barrel[valId], pEBmin[valId], pEBmax[valId]);
0448 barrel[valId]->SetStats(false);
0449 pad[valId][2]->cd();
0450 DrawEE(endc_p[valId], pEEmin[valId], pEEmax[valId]);
0451 }
0452
0453 std::string ImageName(this->m_imageFileName);
0454 canvas.SaveAs(ImageName.c_str());
0455 return true;
0456 }
0457 };
0458 using EcalLinearCorrectionsDiffOneTag = EcalLinearCorrectionsDiffBase<cond::payloadInspector::SINGLE_IOV, 1>;
0459 using EcalLinearCorrectionsDiffTwoTags = EcalLinearCorrectionsDiffBase<cond::payloadInspector::SINGLE_IOV, 2>;
0460
0461 }
0462
0463
0464 PAYLOAD_INSPECTOR_MODULE(EcalLinearCorrections) {
0465 PAYLOAD_INSPECTOR_CLASS(EcalLinearCorrectionsPlot);
0466 PAYLOAD_INSPECTOR_CLASS(EcalLinearCorrectionsDiffOneTag);
0467 PAYLOAD_INSPECTOR_CLASS(EcalLinearCorrectionsDiffTwoTags);
0468 }