EcalPlugins/plugins/EcalPulseCovariances_PayloadInspector.cc

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 // the data format of the condition to be inspected
0008 #include "CondFormats/EcalObjects/interface/EcalPulseCovariances.h"
0009
0010 #include "TH2F.h"
0011 #include "TCanvas.h"
0012 #include "TStyle.h"
0013 #include "TLine.h"
0014 #include "TLatex.h"
0015
0016 #include <string>
0017
0018 namespace {
0019   enum { kEBChannels = 61200, kEEChannels = 14648, kSides = 2, kRMS = 5, TEMPLATESAMPLES = 12 };
0020   enum { MIN_IETA = 1, MIN_IPHI = 1, MAX_IETA = 85, MAX_IPHI = 360 };  // barrel lower and upper bounds on eta and phi
0021   enum { IX_MIN = 1, IY_MIN = 1, IX_MAX = 100, IY_MAX = 100 };         // endcaps lower and upper bounds on x and y
0022
0023   /*****************************************************
0024      2d plot of ECAL PulseCovariances of 1 IOV
0025   *****************************************************/
0026   class EcalPulseCovariancesPlot : public cond::payloadInspector::PlotImage<EcalPulseCovariances> {
0027   public:
0028     EcalPulseCovariancesPlot()
0029         : cond::payloadInspector::PlotImage<EcalPulseCovariances>("ECAL PulseCovariances - 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*[TEMPLATESAMPLES];
0035       TH2F** endc_p = new TH2F*[TEMPLATESAMPLES];
0036       TH2F** endc_m = new TH2F*[TEMPLATESAMPLES];
0037       //      long double EBmean[TEMPLATESAMPLES], EBrms[TEMPLATESAMPLES], EEmean[TEMPLATESAMPLES], EErms[TEMPLATESAMPLES];
0038       //      int EBtot[TEMPLATESAMPLES], EEtot[TEMPLATESAMPLES];
0039       double pEBmin[TEMPLATESAMPLES], pEBmax[TEMPLATESAMPLES], pEEmin[TEMPLATESAMPLES], pEEmax[TEMPLATESAMPLES];
0040       for (int s = 0; s < TEMPLATESAMPLES; ++s) {
0041         barrel[s] = new TH2F(
0042             Form("EBs%i", s), Form("sample %i EB", s), MAX_IPHI, 0, MAX_IPHI, 2 * MAX_IETA, -MAX_IETA, MAX_IETA);
0043         endc_p[s] = new TH2F(
0044             Form("EE+s%i", s), Form("sample %i EE+", s), IX_MAX, IX_MIN, IX_MAX + 1, IY_MAX, IY_MIN, IY_MAX + 1);
0045         endc_m[s] = new TH2F(
0046             Form("EE-s%i", s), Form("sample %i EE-", s), IX_MAX, IX_MIN, IX_MAX + 1, IY_MAX, IY_MIN, IY_MAX + 1);
0047         /*  EBmean[s] = 0.;
0048     EBrms[s] = 0.;
0049     EEmean[s] = 0.;
0050     EErms[s] = 0.;
0051     EBtot[s] = 0;
0052     EEtot[s] = 0;*/
0053         pEBmin[s] = 10.;
0054         pEBmax[s] = -10.;
0055         pEEmin[s] = 10.;
0056         pEEmax[s] = -10.;
0057       }
0058
0059       auto iov = iovs.front();
0060       std::shared_ptr<EcalPulseCovariances> payload = fetchPayload(std::get<1>(iov));
0061       unsigned int run = std::get<0>(iov);
0062       if (payload.get()) {
0063         //  int chan = 0;
0064         for (int ieta = -MAX_IETA; ieta <= MAX_IETA; ieta++) {
0065           Double_t eta = (Double_t)ieta;
0066           if (ieta == 0)
0067             continue;
0068           else if (ieta > 0.)
0069             eta = eta - 0.5;  //   0.5 to 84.5
0070           else
0071             eta = eta + 0.5;  //  -84.5 to -0.5
0072           for (int iphi = 1; iphi <= MAX_IPHI; iphi++) {
0073             Double_t phi = (Double_t)iphi - 0.5;
0074             EBDetId id(ieta, iphi);
0075             /*
0076           chan++;
0077         if(chan < 20) std::cout << " channel " << chan << " hash " << id.hashedIndex() << " raw " << id.rawId() << std::endl;
0078         int rawId = id.rawId();
0079         EcalPulseCovariance pulse = (*payload)[rawId];
0080         if(chan < 20) {
0081           for(int i = 0; i < TEMPLATESAMPLES; ++i) {
0082         std::cout << std::setw(2) << i;
0083         for(int j = 0; j < TEMPLATESAMPLES; ++j) {
0084           double val = pulse.covval[i][j];
0085           std::cout << " " << std::setw(10) << val ;
0086         }
0087         std::cout << std::endl;
0088           }
0089         }
0090         */
0091             for (int s = 0; s < TEMPLATESAMPLES; ++s) {
0092               double val = (*payload)[id.rawId()].covval[0][s];
0093               barrel[s]->Fill(phi, eta, val);
0094               //          EBmean[s] = EBmean[s] + val;
0095               //          EBrms[s] = EBrms[s] + val * val;
0096               //          EBtot[s]++;
0097               if (val < pEBmin[s])
0098                 pEBmin[s] = val;
0099               if (val > pEBmax[s])
0100                 pEBmax[s] = val;
0101             }
0102           }
0103         }
0104
0105         for (int sign = 0; sign < kSides; sign++) {
0106           int thesign = sign == 1 ? 1 : -1;
0107           for (int ix = 1; ix <= IX_MAX; ix++) {
0108             for (int iy = 1; iy <= IY_MAX; iy++) {
0109               if (!EEDetId::validDetId(ix, iy, thesign))
0110                 continue;
0111               EEDetId id(ix, iy, thesign);
0112               for (int s = 0; s < TEMPLATESAMPLES; ++s) {
0113                 double val = (*payload)[id.rawId()].covval[0][s];
0114                 //      EEmean[s] = EEmean[s] + val;
0115                 //      EErms[s] = EErms[s] + val * val;
0116                 //      EEtot[s]++;
0117                 if (val < pEEmin[s])
0118                   pEEmin[s] = val;
0119                 if (val > pEEmax[s])
0120                   pEEmax[s] = val;
0121                 if (thesign == 1)
0122                   endc_p[s]->Fill(ix, iy, val);
0123                 else
0124                   endc_m[s]->Fill(ix, iy, val);
0125               }
0126             }  // iy
0127           }    // ix
0128         }      // side
0129       }        // payload
0130       /*
0131       for(int s = 0; s < TEMPLATESAMPLES; ++s) {
0132     std::cout << "EB sample " << s << " mean " << EBmean[s] << " rms " << EBrms[s] << " entries " << EBtot[s]
0133           << " EE mean " << EEmean[s] << " rms " << EErms[s] << " entries " << EEtot[s]  << std::endl;
0134     double vt =(double)EBtot[s];
0135     EBmean[s] = EBmean[s] / vt;
0136     EBrms[s] = EBrms[s] / vt - (EBmean[s] * EBmean[s]);
0137     if(EBrms[s] > 0) EBrms[s] = sqrt(EBrms[s]);
0138     else EBrms[s] = 1.e-06;
0139     pEBmin[s] = EBmean[s] - kRMS * EBrms[s];
0140     pEBmax[s] = EBmean[s] + kRMS * EBrms[s];
0141     std::cout << "EB sample " << s << " mean " << EBmean[s] << " rms " << EBrms[s] << " entries " << EBtot[s] << " min " << pEBmin[s] << " max " << pEBmax[s] << std::endl;
0142     //  if(pEBmin[s] <= 0.) pEBmin[s] = 1.e-06;
0143     vt =(double)EEtot[s];
0144     EEmean[s] = EEmean[s] / vt;
0145     EErms[s] = EErms[s] / vt - (EEmean[s] * EEmean[s]);
0146     if(EErms[s] > 0) EErms[s] = sqrt(EErms[s]);
0147     else EErms[s] = 1.e-06;
0148     pEEmin[s] = EEmean[s] - kRMS * EErms[s];
0149     pEEmax[s] = EEmean[s] + kRMS * EErms[s];
0150     std::cout << "EE sample " << s  << " mean " << EEmean[s] << " rms " << EErms[s] << " entries " << EEtot[s] << " min " << pEEmin[s] << " max " << pEEmax[s] << std::endl;
0151     //  if(pEEmin[s] <= 0.) pEEmin[s] = 1.e-06;
0152       }
0153       */
0154       //      for(int s = 0; s < TEMPLATESAMPLES; ++s)
0155       //    std::cout << " sample " << s << " EB min " << pEBmin[s] << " max " << pEBmax[s] << " EE min " << pEEmin[s] << " max " << pEEmax[s] << std::endl;
0156       gStyle->SetPalette(1);
0157       gStyle->SetOptStat(0);
0158       TCanvas canvas("CC map", "CC map", 1600, 2800);
0159       TLatex t1;
0160       t1.SetNDC();
0161       t1.SetTextAlign(26);
0162       t1.SetTextSize(0.05);
0163       t1.DrawLatex(0.5, 0.96, Form("Ecal PulseCovariances, IOV %i", run));
0164
0165       float xmi[3] = {0.0, 0.24, 0.76};
0166       float xma[3] = {0.24, 0.76, 1.00};
0167       TPad*** pad = new TPad**[6];
0168       for (int s = 0; s < 6; s++) {
0169         pad[s] = new TPad*[3];
0170         for (int obj = 0; obj < 3; obj++) {
0171           float yma = 0.94 - (0.16 * s);
0172           float ymi = yma - 0.14;
0173           pad[s][obj] = new TPad(Form("p_%i_%i", obj, s), Form("p_%i_%i", obj, s), xmi[obj], ymi, xma[obj], yma);
0174           pad[s][obj]->Draw();
0175         }
0176       }
0177
0178       int ipad = 0;
0179       for (int s = 0; s < 7; ++s) {  // plot only the measured ones, not the extrapolated
0180         // for(int s = 7; s<12; ++s) {  // plot only the extrapolated ones
0181         if (s == 2)
0182           continue;  // do not plot the maximum sample, which is 0 by default
0183         pad[ipad][0]->cd();
0184         if (pEBmin[s] == pEBmax[s]) {  // same values everywhere!..
0185           pEBmin[s] = pEBmin[s] - 1.e-06;
0186           pEBmax[s] = pEBmax[s] + 1.e-06;
0187         }
0188         if (pEEmin[s] == pEEmax[s]) {
0189           pEEmin[s] = pEEmin[s] - 1.e-06;
0190           pEEmax[s] = pEEmax[s] + 1.e-06;
0191         }
0192         DrawEE(endc_m[s], pEEmin[s], pEEmax[s]);
0193         //  pad[ipad][0]->SetLogz(1);
0194         pad[ipad][1]->cd();
0195         DrawEB(barrel[s], pEBmin[s], pEBmax[s]);
0196         //  pad[ipad][1]->SetLogz(1);
0197         pad[ipad][2]->cd();
0198         DrawEE(endc_p[s], pEEmin[s], pEEmax[s]);
0199         //  pad[ipad][2]->SetLogz(1);
0200         ipad++;
0201       }
0202
0203       std::string ImageName(m_imageFileName);
0204       canvas.SaveAs(ImageName.c_str());
0205       return true;
0206     }  // fill method
0207   };
0208
0209   /*****************************************************
0210      2d plot of ECAL PulseCovariances matrix of 1 IOV
0211   *****************************************************/
0212   class EcalPulseCovariancesMatrix : public cond::payloadInspector::PlotImage<EcalPulseCovariances> {
0213   public:
0214     EcalPulseCovariancesMatrix()
0215         : cond::payloadInspector::PlotImage<EcalPulseCovariances>("ECAL PulseCovariances - matrix ") {
0216       setSingleIov(true);
0217     }
0218
0219     bool fill(const std::vector<std::tuple<cond::Time_t, cond::Hash> >& iovs) override {
0220       double EBmean[TEMPLATESAMPLES][TEMPLATESAMPLES], EBrms[TEMPLATESAMPLES][TEMPLATESAMPLES],
0221           EEmean[TEMPLATESAMPLES][TEMPLATESAMPLES], EErms[TEMPLATESAMPLES][TEMPLATESAMPLES];
0222       for (int i = 0; i < TEMPLATESAMPLES; i++) {
0223         for (int j = 0; j < TEMPLATESAMPLES; j++) {
0224           EBmean[i][j] = 0.;
0225           EBrms[i][j] = 0.;
0226           EEmean[i][j] = 0.;
0227           EErms[i][j] = 0.;
0228         }
0229       }
0230
0231       auto iov = iovs.front();
0232       std::shared_ptr<EcalPulseCovariances> payload = fetchPayload(std::get<1>(iov));
0233       unsigned int run = std::get<0>(iov);
0234       if (payload.get()) {
0235         for (int ieta = -MAX_IETA; ieta <= MAX_IETA; ieta++) {
0236           for (int iphi = 1; iphi <= MAX_IPHI; iphi++) {
0237             EBDetId id(ieta, iphi);
0238             for (int i = 0; i < TEMPLATESAMPLES; ++i) {
0239               for (int j = 0; j < TEMPLATESAMPLES; ++j) {
0240                 double val = (*payload)[id.rawId()].covval[i][j];
0241                 EBmean[i][j] = EBmean[i][j] + val;
0242                 EBrms[i][j] = EBrms[i][j] + val * val;
0243               }
0244             }
0245           }
0246         }
0247
0248         for (int sign = 0; sign < kSides; sign++) {
0249           int thesign = sign == 1 ? 1 : -1;
0250           for (int ix = 1; ix <= IX_MAX; ix++) {
0251             for (int iy = 1; iy <= IY_MAX; iy++) {
0252               if (!EEDetId::validDetId(ix, iy, thesign))
0253                 continue;
0254               EEDetId id(ix, iy, thesign);
0255               for (int i = 0; i < TEMPLATESAMPLES; i++) {
0256                 for (int j = 0; j < TEMPLATESAMPLES; j++) {
0257                   double val = (*payload)[id.rawId()].covval[i][j];
0258                   EEmean[i][j] = EEmean[i][j] + val;
0259                   EErms[i][j] = EErms[i][j] + val * val;
0260                 }
0261               }
0262             }  // iy
0263           }    // ix
0264         }      // side
0265       }        // payload
0266
0267       TH2F* barrel_m =
0268           new TH2F("EBm", "EB mean", TEMPLATESAMPLES, 0, TEMPLATESAMPLES, TEMPLATESAMPLES, 0., TEMPLATESAMPLES);
0269       TH2F* barrel_r =
0270           new TH2F("EBr", "EB rms", TEMPLATESAMPLES, 0, TEMPLATESAMPLES, TEMPLATESAMPLES, 0., TEMPLATESAMPLES);
0271       TH2F* endcap_m =
0272           new TH2F("EEm", "EE mean", TEMPLATESAMPLES, 0, TEMPLATESAMPLES, TEMPLATESAMPLES, 0., TEMPLATESAMPLES);
0273       TH2F* endcap_r =
0274           new TH2F("EEr", "EE rms", TEMPLATESAMPLES, 0, TEMPLATESAMPLES, TEMPLATESAMPLES, 0., TEMPLATESAMPLES);
0275       for (int i = 0; i < TEMPLATESAMPLES; i++) {
0276         //  std::cout << "EB sample " << i << std::endl
0277         //  std::cout << "EB sample " << i;
0278         for (int j = 0; j < TEMPLATESAMPLES; j++) {
0279           double vt = (double)kEBChannels;
0280           EBmean[i][j] = EBmean[i][j] / vt;
0281           barrel_m->Fill(i, j, EBmean[i][j]);
0282           EBrms[i][j] = EBrms[i][j] / vt - (EBmean[i][j] * EBmean[i][j]);
0283           if (EBrms[i][j] >= 0)
0284             EBrms[i][j] = sqrt(EBrms[i][j]);
0285           else
0286             EBrms[i][j] = 0.;
0287           barrel_r->Fill(i, j, EBrms[i][j]);
0288           //      std::cout << "EB sample " << j << " mean " << EBmean[i][j] << " rms " << EBrms[i][j] << std::endl;
0289           //      std::cout << " " << std::setw(10) << EBrms[i][j];
0290           vt = (double)kEEChannels;
0291           EEmean[i][j] = EEmean[i][j] / vt;
0292           endcap_m->Fill(i, j, EEmean[i][j]);
0293           EErms[i][j] = EErms[i][j] / vt - (EEmean[i][j] * EEmean[i][j]);
0294           if (EErms[i][j] >= 0)
0295             EErms[i][j] = sqrt(EErms[i][j]);
0296           else
0297             EErms[i][j] = 0.;
0298           endcap_r->Fill(i, j, EErms[i][j]);
0299           //      std::cout << "EE sample " << j  << " mean " << EEmean[i][j] << " rms " << EErms[i][j] << std::endl;
0300         }
0301         //  std::cout << std::endl;
0302       }
0303
0304       gStyle->SetPalette(1);
0305       gStyle->SetOptStat(0);
0306       TCanvas canvas("CC map", "CC map", 1440, 1500);
0307       TLatex t1;
0308       t1.SetNDC();
0309       t1.SetTextAlign(26);
0310       t1.SetTextSize(0.05);
0311       t1.DrawLatex(0.5, 0.96, Form("Ecal PulseCovariances, IOV %i", run));
0312
0313       float xmi[2] = {0.0, 0.5};
0314       float xma[2] = {0.5, 1.0};
0315       TPad*** pad = new TPad**[2];
0316       for (int s = 0; s < 2; s++) {
0317         pad[s] = new TPad*[2];
0318         for (int obj = 0; obj < 2; obj++) {
0319           float yma = 0.94 - (0.47 * s);
0320           float ymi = yma - 0.45;
0321           pad[s][obj] = new TPad(Form("p_%i_%i", obj, s), Form("p_%i_%i", obj, s), xmi[obj], ymi, xma[obj], yma);
0322           pad[s][obj]->Draw();
0323         }
0324       }
0325
0326       pad[0][0]->cd();
0327       //      barrel_m->Draw("COLZ1");
0328       grid(barrel_m);
0329       pad[0][1]->cd();
0330       //      endcap_m->Draw("COLZ1");
0331       grid(endcap_m);
0332       pad[1][0]->cd();
0333       //      barrel_r->Draw("COLZ1");
0334       grid(barrel_r);
0335       pad[1][1]->cd();
0336       //      endcap_r->Draw("COLZ1");
0337       grid(endcap_r);
0338
0339       std::string ImageName(m_imageFileName);
0340       canvas.SaveAs(ImageName.c_str());
0341       return true;
0342     }  // fill method
0343
0344     void grid(TH2F* matrix) {
0345       matrix->Draw("COLZ1");
0346       TLine* lh = new TLine(1., 0., 1., 1.);
0347       lh->SetLineWidth(2);
0348       TLine* lv = new TLine(1., 0., 1., 1.);
0349       lv->SetLineWidth(2);
0350       //      double max = (double)TEMPLATESAMPLES;
0351       for (int i = 1; i < TEMPLATESAMPLES; i++) {
0352         //  double x = (double)i;
0353         lv = new TLine(i, 0., i, TEMPLATESAMPLES);
0354         lv->Draw();
0355         lh = new TLine(0., i, TEMPLATESAMPLES, i);
0356         lh->Draw();
0357       }
0358     }
0359   };
0360
0361 }  // namespace
0362
0363 // Register the classes as boost python plugin
0364 PAYLOAD_INSPECTOR_MODULE(EcalPulseCovariances) {
0365   PAYLOAD_INSPECTOR_CLASS(EcalPulseCovariancesPlot);
0366   PAYLOAD_INSPECTOR_CLASS(EcalPulseCovariancesMatrix);
0367 }