Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-09-07 04:35:29

0001 #include "CondCore/Utilities/interface/PayloadInspectorModule.h"
0002 #include "CondCore/Utilities/interface/PayloadInspector.h"
0003 #include "CondCore/CondDB/interface/Time.h"
0004 #include "DataFormats/EcalDetId/interface/EBDetId.h"
0005 #include "DataFormats/EcalDetId/interface/EEDetId.h"
0006 #include "CondCore/EcalPlugins/plugins/EcalDrawUtils.h"
0007 
0008 // the data format of the condition to be inspected
0009 #include "CondFormats/EcalObjects/interface/EcalTPGPedestals.h"
0010 
0011 #include "TH2F.h"
0012 #include "TCanvas.h"
0013 #include "TLine.h"
0014 #include "TStyle.h"
0015 #include "TLatex.h"
0016 #include "TPave.h"
0017 #include "TPaveStats.h"
0018 #include <string>
0019 #include <fstream>
0020 
0021 namespace {
0022   enum { kEBChannels = 61200, kEEChannels = 14648, kGains = 3, kSides = 2 };
0023   enum { MIN_IETA = 1, MIN_IPHI = 1, MAX_IETA = 85, MAX_IPHI = 360 };  // barrel lower and upper bounds on eta and phi
0024   enum { IX_MIN = 1, IY_MIN = 1, IX_MAX = 100, IY_MAX = 100 };         // endcaps lower and upper bounds on x and y
0025   int gainValues[kGains] = {12, 6, 1};
0026 
0027   /**************************************************
0028      2d plot of ECAL TPGPedestals of 1 IOV
0029   **************************************************/
0030   class EcalTPGPedestalsPlot : public cond::payloadInspector::PlotImage<EcalTPGPedestals> {
0031   public:
0032     EcalTPGPedestalsPlot() : cond::payloadInspector::PlotImage<EcalTPGPedestals>("ECAL Gain Ratios - map ") {
0033       setSingleIov(true);
0034     }
0035 
0036     bool fill(const std::vector<std::tuple<cond::Time_t, cond::Hash> >& iovs) override {
0037       TH2F** barrel_m = new TH2F*[kGains];
0038       TH2F** endc_p_m = new TH2F*[kGains];
0039       TH2F** endc_m_m = new TH2F*[kGains];
0040       float mEBmin[kGains], mEEmin[kGains], mEBmax[kGains], mEEmax[kGains];
0041       for (int gainId = 0; gainId < kGains; gainId++) {
0042         barrel_m[gainId] = new TH2F(Form("EBm%i", gainId),
0043                                     Form("EB mean_x%i ", gainValues[gainId]),
0044                                     MAX_IPHI,
0045                                     0,
0046                                     MAX_IPHI,
0047                                     2 * MAX_IETA,
0048                                     -MAX_IETA,
0049                                     MAX_IETA);
0050         endc_p_m[gainId] = new TH2F(Form("EE+m%i", gainId),
0051                                     Form("EE+ mean_x%i", gainValues[gainId]),
0052                                     IX_MAX,
0053                                     IX_MIN,
0054                                     IX_MAX + 1,
0055                                     IY_MAX,
0056                                     IY_MIN,
0057                                     IY_MAX + 1);
0058         endc_m_m[gainId] = new TH2F(Form("EE-m%i", gainId),
0059                                     Form("EE- mean_x%i", gainValues[gainId]),
0060                                     IX_MAX,
0061                                     IX_MIN,
0062                                     IX_MAX + 1,
0063                                     IY_MAX,
0064                                     IY_MIN,
0065                                     IY_MAX + 1);
0066         mEBmin[gainId] = 10.;
0067         mEEmin[gainId] = 10.;
0068         mEBmax[gainId] = -10.;
0069         mEEmax[gainId] = -10.;
0070       }
0071 
0072       //      std::ofstream fout;
0073       //      fout.open("./bid.txt");
0074       auto iov = iovs.front();
0075       std::shared_ptr<EcalTPGPedestals> payload = fetchPayload(std::get<1>(iov));
0076       unsigned int run = std::get<0>(iov);
0077       if (payload.get()) {
0078         for (int sign = 0; sign < kSides; sign++) {
0079           int thesign = sign == 1 ? 1 : -1;
0080 
0081           for (int ieta = 0; ieta < MAX_IETA; ieta++) {
0082             for (int iphi = 0; iphi < MAX_IPHI; iphi++) {
0083               EBDetId id((ieta + 1) * thesign, iphi + 1);
0084               float y = -1 - ieta;
0085               if (sign == 1)
0086                 y = ieta;
0087               float val = (*payload)[id.rawId()].mean_x12;
0088               barrel_m[0]->Fill(iphi, y, val);
0089               if (val < mEBmin[0])
0090                 mEBmin[0] = val;
0091               if (val > mEBmax[0])
0092                 mEBmax[0] = val;
0093               val = (*payload)[id.rawId()].mean_x6;
0094               barrel_m[1]->Fill(iphi, y, val);
0095               if (val < mEBmin[1])
0096                 mEBmin[1] = val;
0097               if (val > mEBmax[1])
0098                 mEBmax[1] = val;
0099               val = (*payload)[id.rawId()].mean_x1;
0100               barrel_m[2]->Fill(iphi, y, val);
0101               if (val < mEBmin[2])
0102                 mEBmin[2] = val;
0103               if (val > mEBmax[2])
0104                 mEBmax[2] = val;
0105             }  // iphi
0106           }  // ieta
0107 
0108           for (int ix = 0; ix < IX_MAX; ix++) {
0109             for (int iy = 0; iy < IY_MAX; iy++) {
0110               if (!EEDetId::validDetId(ix + 1, iy + 1, thesign))
0111                 continue;
0112               EEDetId id(ix + 1, iy + 1, thesign);
0113               float val = (*payload)[id.rawId()].mean_x12;
0114               if (thesign == 1)
0115                 endc_p_m[0]->Fill(ix + 1, iy + 1, val);
0116               else
0117                 endc_m_m[0]->Fill(ix + 1, iy + 1, val);
0118               if (val < mEEmin[0])
0119                 mEEmin[0] = val;
0120               if (val > mEEmax[0])
0121                 mEEmax[0] = val;
0122               val = (*payload)[id.rawId()].mean_x6;
0123               if (thesign == 1)
0124                 endc_p_m[1]->Fill(ix + 1, iy + 1, val);
0125               else
0126                 endc_m_m[1]->Fill(ix + 1, iy + 1, val);
0127               if (val < mEEmin[1])
0128                 mEEmin[1] = val;
0129               if (val > mEEmax[1])
0130                 mEEmax[1] = val;
0131               val = (*payload)[id.rawId()].mean_x1;
0132               if (thesign == 1)
0133                 endc_p_m[2]->Fill(ix + 1, iy + 1, val);
0134               else
0135                 endc_m_m[2]->Fill(ix + 1, iy + 1, val);
0136               if (val < mEEmin[2])
0137                 mEEmin[2] = val;
0138               if (val > mEEmax[2])
0139                 mEEmax[2] = val;
0140               //          fout << " x " << ix << " y " << " val " << val << std::endl;
0141             }  // iy
0142           }  // ix
0143         }  // side
0144       }  // if payload.get()
0145       else
0146         return false;
0147       //      std::cout << " min " << rEEmin[2] << " max " << rEEmax[2] << std::endl;
0148       //      fout.close();
0149 
0150       gStyle->SetPalette(1);
0151       gStyle->SetOptStat(0);
0152       TCanvas canvas("CC map", "CC map", 1200, 900);
0153       TLatex t1;
0154       t1.SetNDC();
0155       t1.SetTextAlign(26);
0156       t1.SetTextSize(0.05);
0157       t1.DrawLatex(0.5, 0.96, Form("Ecal Gain TPGPedestals, IOV %i", run));
0158 
0159       float xmi[3] = {0.0, 0.22, 0.78};
0160       float xma[3] = {0.22, 0.78, 1.00};
0161       TPad*** pad = new TPad**[kGains];
0162       for (int gId = 0; gId < kGains; gId++) {
0163         pad[gId] = new TPad*[3];
0164         for (int obj = 0; obj < 3; obj++) {
0165           float yma = 0.94 - (0.32 * gId);
0166           float ymi = yma - 0.30;
0167           pad[gId][obj] = new TPad(Form("p_%i_%i", obj, gId), Form("p_%i_%i", obj, gId), xmi[obj], ymi, xma[obj], yma);
0168           pad[gId][obj]->Draw();
0169         }
0170       }
0171 
0172       for (int gId = 0; gId < kGains; gId++) {
0173         pad[gId][0]->cd();
0174         DrawEE(endc_m_m[gId], mEEmin[gId], mEEmax[gId]);
0175         pad[gId][1]->cd();
0176         DrawEB(barrel_m[gId], mEBmin[gId], mEBmax[gId]);
0177         pad[gId][2]->cd();
0178         DrawEE(endc_p_m[gId], mEEmin[gId], mEEmax[gId]);
0179       }
0180 
0181       std::string ImageName(m_imageFileName);
0182       canvas.SaveAs(ImageName.c_str());
0183       return true;
0184     }  // fill method
0185   };
0186 
0187   /******************************************************************
0188      2d plot of ECAL TPGPedestals difference between 2 IOVs
0189   ******************************************************************/
0190   template <cond::payloadInspector::IOVMultiplicity nIOVs, int ntags, int method>
0191   class EcalTPGPedestalsBase : public cond::payloadInspector::PlotImage<EcalTPGPedestals, nIOVs, ntags> {
0192   public:
0193     EcalTPGPedestalsBase()
0194         : cond::payloadInspector::PlotImage<EcalTPGPedestals, nIOVs, ntags>("ECAL Gain Ratios comparison") {}
0195     bool fill() override {
0196       TH2F** barrel_m = new TH2F*[kGains];
0197       TH2F** endc_p_m = new TH2F*[kGains];
0198       TH2F** endc_m_m = new TH2F*[kGains];
0199       float mEBmin[kGains], mEEmin[kGains], mEBmax[kGains], mEEmax[kGains];
0200       float mEB[kGains][kEBChannels], mEE[kGains][kEEChannels];
0201       for (int gainId = 0; gainId < kGains; gainId++) {
0202         barrel_m[gainId] = new TH2F(Form("EBm%i", gainId),
0203                                     Form("EB mean_x%i ", gainValues[gainId]),
0204                                     MAX_IPHI,
0205                                     0,
0206                                     MAX_IPHI,
0207                                     2 * MAX_IETA,
0208                                     -MAX_IETA,
0209                                     MAX_IETA);
0210         endc_p_m[gainId] = new TH2F(Form("EE+m%i", gainId),
0211                                     Form("EE+ mean_x%i", gainValues[gainId]),
0212                                     IX_MAX,
0213                                     IX_MIN,
0214                                     IX_MAX + 1,
0215                                     IY_MAX,
0216                                     IY_MIN,
0217                                     IY_MAX + 1);
0218         endc_m_m[gainId] = new TH2F(Form("EE-m%i", gainId),
0219                                     Form("EE- mean_x%i", gainValues[gainId]),
0220                                     IX_MAX,
0221                                     IX_MIN,
0222                                     IX_MAX + 1,
0223                                     IY_MAX,
0224                                     IY_MIN,
0225                                     IY_MAX + 1);
0226         mEBmin[gainId] = 10.;
0227         mEEmin[gainId] = 10.;
0228         mEBmax[gainId] = -10.;
0229         mEEmax[gainId] = -10.;
0230       }
0231 
0232       unsigned int run[2];
0233       std::string l_tagname[2];
0234       auto iovs = cond::payloadInspector::PlotBase::getTag<0>().iovs;
0235       l_tagname[0] = cond::payloadInspector::PlotBase::getTag<0>().name;
0236       auto firstiov = iovs.front();
0237       run[0] = std::get<0>(firstiov);
0238       std::tuple<cond::Time_t, cond::Hash> lastiov;
0239       if (ntags == 2) {
0240         auto tag2iovs = cond::payloadInspector::PlotBase::getTag<1>().iovs;
0241         l_tagname[1] = cond::payloadInspector::PlotBase::getTag<1>().name;
0242         lastiov = tag2iovs.front();
0243       } else {
0244         lastiov = iovs.back();
0245         l_tagname[1] = l_tagname[0];
0246       }
0247       run[1] = std::get<0>(lastiov);
0248       for (int irun = 0; irun < nIOVs; irun++) {
0249         std::shared_ptr<EcalTPGPedestals> payload;
0250         if (irun == 0) {
0251           payload = this->fetchPayload(std::get<1>(firstiov));
0252         } else {
0253           payload = this->fetchPayload(std::get<1>(lastiov));
0254         }
0255         if (payload.get()) {
0256           float dr;
0257           for (int sign = 0; sign < kSides; sign++) {
0258             int thesign = sign == 1 ? 1 : -1;
0259 
0260             for (int ieta = 0; ieta < MAX_IETA; ieta++) {
0261               for (int iphi = 0; iphi < MAX_IPHI; iphi++) {
0262                 EBDetId id((ieta + 1) * thesign, iphi + 1);
0263                 int hashindex = id.hashedIndex();
0264                 float y = -1 - ieta;
0265                 if (sign == 1)
0266                   y = ieta;
0267                 float val = (*payload)[id.rawId()].mean_x12;
0268                 if (irun == 0) {
0269                   mEB[0][hashindex] = val;
0270                 } else {
0271                   if (method == 0)
0272                     dr = val - mEB[0][hashindex];
0273                   else {
0274                     if (mEB[0][hashindex] == 0.) {
0275                       if (val == 0.)
0276                         dr = 1.;
0277                       else
0278                         dr = 9999.;
0279                     } else
0280                       dr = val / mEB[0][hashindex];
0281                   }
0282                   barrel_m[0]->Fill(iphi, y, dr);
0283                   if (dr < mEBmin[0])
0284                     mEBmin[0] = dr;
0285                   if (dr > mEBmax[0])
0286                     mEBmax[0] = dr;
0287                 }
0288                 val = (*payload)[id.rawId()].mean_x6;
0289                 if (irun == 0) {
0290                   mEB[1][hashindex] = val;
0291                 } else {
0292                   if (method == 0)
0293                     dr = val - mEB[1][hashindex];
0294                   else {
0295                     if (mEB[1][hashindex] == 0.) {
0296                       if (val == 0.)
0297                         dr = 1.;
0298                       else
0299                         dr = 9999.;
0300                     } else
0301                       dr = val / mEB[1][hashindex];
0302                   }
0303                   barrel_m[1]->Fill(iphi, y, dr);
0304                   if (dr < mEBmin[1])
0305                     mEBmin[1] = dr;
0306                   if (dr > mEBmax[1])
0307                     mEBmax[1] = dr;
0308                 }
0309                 val = (*payload)[id.rawId()].mean_x1;
0310                 if (irun == 0) {
0311                   mEB[2][hashindex] = val;
0312                 } else {
0313                   if (method == 0)
0314                     dr = val - mEB[2][hashindex];
0315                   else {
0316                     if (mEB[2][hashindex] == 0.) {
0317                       if (val == 0.)
0318                         dr = 1.;
0319                       else
0320                         dr = 9999.;
0321                     } else
0322                       dr = val / mEB[2][hashindex];
0323                   }
0324                   barrel_m[2]->Fill(iphi, y, dr);
0325                   if (dr < mEBmin[2])
0326                     mEBmin[2] = dr;
0327                   if (dr > mEBmax[2])
0328                     mEBmax[2] = dr;
0329                 }
0330               }  // iphi
0331             }  // ieta
0332 
0333             for (int ix = 0; ix < IX_MAX; ix++) {
0334               for (int iy = 0; iy < IY_MAX; iy++) {
0335                 if (!EEDetId::validDetId(ix + 1, iy + 1, thesign))
0336                   continue;
0337                 EEDetId id(ix + 1, iy + 1, thesign);
0338                 int hashindex = id.hashedIndex();
0339                 float val = (*payload)[id.rawId()].mean_x12;
0340                 if (irun == 0) {
0341                   mEE[0][hashindex] = val;
0342                 } else {
0343                   if (method == 0)
0344                     dr = val - mEE[0][hashindex];
0345                   else {
0346                     if (mEE[0][hashindex] == 0.) {
0347                       if (val == 0.)
0348                         dr = 1.;
0349                       else
0350                         dr = 9999.;
0351                     } else
0352                       dr = val / mEE[0][hashindex];
0353                   }
0354                   if (thesign == 1)
0355                     endc_p_m[0]->Fill(ix + 1, iy + 1, dr);
0356                   else
0357                     endc_m_m[0]->Fill(ix + 1, iy + 1, dr);
0358                   if (dr < mEEmin[0])
0359                     mEEmin[0] = dr;
0360                   if (dr > mEEmax[0])
0361                     mEEmax[0] = dr;
0362                 }
0363                 val = (*payload)[id.rawId()].mean_x6;
0364                 if (irun == 0) {
0365                   mEE[1][hashindex] = val;
0366                 } else {
0367                   if (method == 0)
0368                     dr = val - mEE[1][hashindex];
0369                   else {
0370                     if (mEE[1][hashindex] == 0.) {
0371                       if (val == 0.)
0372                         dr = 1.;
0373                       else
0374                         dr = 9999.;
0375                     } else
0376                       dr = val / mEE[1][hashindex];
0377                   }
0378                   if (thesign == 1)
0379                     endc_p_m[1]->Fill(ix + 1, iy + 1, dr);
0380                   else
0381                     endc_m_m[1]->Fill(ix + 1, iy + 1, dr);
0382                   if (dr < mEEmin[1])
0383                     mEEmin[1] = dr;
0384                   if (dr > mEEmax[1])
0385                     mEEmax[1] = dr;
0386                 }
0387                 val = (*payload)[id.rawId()].mean_x1;
0388                 if (irun == 0) {
0389                   mEE[2][hashindex] = val;
0390                 } else {
0391                   if (method == 0)
0392                     dr = val - mEE[2][hashindex];
0393                   else {
0394                     if (mEE[2][hashindex] == 0.) {
0395                       if (val == 0.)
0396                         dr = 1.;
0397                       else
0398                         dr = 9999.;
0399                     } else
0400                       dr = val / mEE[2][hashindex];
0401                   }
0402                   if (thesign == 1)
0403                     endc_p_m[2]->Fill(ix + 1, iy + 1, dr);
0404                   else
0405                     endc_m_m[2]->Fill(ix + 1, iy + 1, dr);
0406                   if (dr < mEEmin[2])
0407                     mEEmin[2] = dr;
0408                   if (dr > mEEmax[2])
0409                     mEEmax[2] = dr;
0410                 }
0411                 //        fout << " x " << ix << " y " << " dr " << dr << std::endl;
0412               }  // iy
0413             }  // ix
0414           }  // side
0415         }  //  if payload.get()
0416         else
0417           return false;
0418       }  // loop over IOVs
0419 
0420       gStyle->SetPalette(1);
0421       gStyle->SetOptStat(0);
0422       TCanvas canvas("CC map", "CC map", 1200, 900);
0423       TLatex t1;
0424       t1.SetNDC();
0425       t1.SetTextAlign(26);
0426       int len = l_tagname[0].length() + l_tagname[1].length();
0427       std::string dr[2] = {"-", "/"};
0428       if (ntags == 2) {
0429         if (len < 150) {
0430           t1.SetTextSize(0.03);
0431           t1.DrawLatex(
0432               0.5,
0433               0.96,
0434               Form("%s %i %s %s %i", l_tagname[1].c_str(), run[1], dr[method].c_str(), l_tagname[0].c_str(), run[0]));
0435         } else {
0436           t1.SetTextSize(0.03);
0437           t1.DrawLatex(0.5, 0.96, Form("Ecal TPGPedestals, IOV %i %s %i", run[1], dr[method].c_str(), run[0]));
0438         }
0439       } else {
0440         t1.SetTextSize(0.03);
0441         t1.DrawLatex(0.5, 0.96, Form("%s, IOV %i %s %i", l_tagname[0].c_str(), run[1], dr[method].c_str(), run[0]));
0442       }
0443 
0444       float xmi[3] = {0.0, 0.22, 0.78};
0445       float xma[3] = {0.22, 0.78, 1.00};
0446       TPad*** pad = new TPad**[kGains];
0447       for (int gId = 0; gId < kGains; gId++) {
0448         pad[gId] = new TPad*[3];
0449         for (int obj = 0; obj < 3; obj++) {
0450           float yma = 0.94 - (0.32 * gId);
0451           float ymi = yma - 0.30;
0452           pad[gId][obj] = new TPad(Form("p_%i_%i", obj, gId), Form("p_%i_%i", obj, gId), xmi[obj], ymi, xma[obj], yma);
0453           pad[gId][obj]->Draw();
0454         }
0455       }
0456 
0457       for (int gId = 0; gId < kGains; gId++) {
0458         pad[gId][0]->cd();
0459         DrawEE(endc_m_m[gId], mEEmin[gId], mEEmax[gId]);
0460         pad[gId][1]->cd();
0461         DrawEB(barrel_m[gId], mEBmin[gId], mEBmax[gId]);
0462         pad[gId][2]->cd();
0463         DrawEE(endc_p_m[gId], mEEmin[gId], mEEmax[gId]);
0464       }
0465 
0466       std::string ImageName(this->m_imageFileName);
0467       canvas.SaveAs(ImageName.c_str());
0468       return true;
0469     }  // fill method
0470   };  // class EcalTPGPedestalsDiffBase
0471   using EcalTPGPedestalsDiffOneTag = EcalTPGPedestalsBase<cond::payloadInspector::SINGLE_IOV, 1, 0>;
0472   using EcalTPGPedestalsDiffTwoTags = EcalTPGPedestalsBase<cond::payloadInspector::SINGLE_IOV, 2, 0>;
0473   using EcalTPGPedestalsRatioOneTag = EcalTPGPedestalsBase<cond::payloadInspector::SINGLE_IOV, 1, 1>;
0474   using EcalTPGPedestalsRatioTwoTags = EcalTPGPedestalsBase<cond::payloadInspector::SINGLE_IOV, 2, 1>;
0475 
0476 }  // namespace
0477 
0478 // Register the classes as boost python plugin
0479 PAYLOAD_INSPECTOR_MODULE(EcalTPGPedestals) {
0480   PAYLOAD_INSPECTOR_CLASS(EcalTPGPedestalsPlot);
0481   PAYLOAD_INSPECTOR_CLASS(EcalTPGPedestalsDiffOneTag);
0482   PAYLOAD_INSPECTOR_CLASS(EcalTPGPedestalsDiffTwoTags);
0483   PAYLOAD_INSPECTOR_CLASS(EcalTPGPedestalsRatioOneTag);
0484   PAYLOAD_INSPECTOR_CLASS(EcalTPGPedestalsRatioTwoTags);
0485 }