File indexing completed on 2024-09-07 04:35:28
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/EcalPulseShapes.h"
0009
0010 #include "TH2F.h"
0011 #include "TProfile.h"
0012 #include "TCanvas.h"
0013 #include "TStyle.h"
0014 #include "TLine.h"
0015 #include "TLatex.h"
0016
0017 #include <string>
0018
0019 namespace {
0020 constexpr int kSides = 2, kRMS = 3, TEMPLATESAMPLES = 12;
0021 constexpr int MAX_IETA = 85, MAX_IPHI = 360;
0022 constexpr int IX_MIN = 1, IY_MIN = 1, IX_MAX = 100, IY_MAX = 100;
0023
0024
0025
0026
0027 class EcalPulseShapesPlot : public cond::payloadInspector::PlotImage<EcalPulseShapes> {
0028 public:
0029 EcalPulseShapesPlot() : cond::payloadInspector::PlotImage<EcalPulseShapes>("ECAL PulseShapes - 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 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
0054
0055
0056
0057 }
0058
0059 auto iov = iovs.front();
0060 std::shared_ptr<EcalPulseShapes> payload = fetchPayload(std::get<1>(iov));
0061 unsigned int run = std::get<0>(iov);
0062 if (payload.get()) {
0063
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;
0070 else
0071 eta = eta + 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 for (int s = 0; s < TEMPLATESAMPLES; s++) {
0076 double val = (*payload)[id.rawId()].pdfval[s];
0077 barrel[s]->Fill(phi, eta, val);
0078 EBmean[s] = EBmean[s] + val;
0079 EBrms[s] = EBrms[s] + val * val;
0080 EBtot[s]++;
0081
0082
0083 }
0084 }
0085 }
0086
0087 for (int sign = 0; sign < kSides; sign++) {
0088 int thesign = sign == 1 ? 1 : -1;
0089 for (int ix = 1; ix <= IX_MAX; ix++) {
0090 for (int iy = 1; iy <= IY_MAX; iy++) {
0091 if (!EEDetId::validDetId(ix, iy, thesign))
0092 continue;
0093 EEDetId id(ix, iy, thesign);
0094 for (int s = 0; s < TEMPLATESAMPLES; s++) {
0095 double val = (*payload)[id.rawId()].pdfval[s];
0096 EEmean[s] = EEmean[s] + val;
0097 EErms[s] = EErms[s] + val * val;
0098 EEtot[s]++;
0099
0100
0101 if (thesign == 1)
0102 endc_p[s]->Fill(ix, iy, val);
0103 else
0104 endc_m[s]->Fill(ix, iy, val);
0105 }
0106 }
0107 }
0108 }
0109 }
0110
0111 for (int s = 0; s < TEMPLATESAMPLES; s++) {
0112 std::cout << "EB sample " << s << " mean " << EBmean[s] << " rms " << EBrms[s] << " entries " << EBtot[s]
0113 << " EE mean " << EEmean[s] << " rms " << EErms[s] << " entries " << EEtot[s] << std::endl;
0114 double vt = (double)EBtot[s];
0115 EBmean[s] = EBmean[s] / vt;
0116 EBrms[s] = EBrms[s] / vt - (EBmean[s] * EBmean[s]);
0117 if (EBrms[s] > 0)
0118 EBrms[s] = sqrt(EBrms[s]);
0119 else
0120 EBrms[s] = 1.e-06;
0121 pEBmin[s] = EBmean[s] - kRMS * EBrms[s];
0122 pEBmax[s] = EBmean[s] + kRMS * EBrms[s];
0123 std::cout << "EB sample " << s << " mean " << EBmean[s] << " rms " << EBrms[s] << " entries " << EBtot[s]
0124 << " min " << pEBmin[s] << " max " << pEBmax[s] << std::endl;
0125
0126 vt = (double)EEtot[s];
0127 EEmean[s] = EEmean[s] / vt;
0128 EErms[s] = EErms[s] / vt - (EEmean[s] * EEmean[s]);
0129 if (EErms[s] > 0)
0130 EErms[s] = sqrt(EErms[s]);
0131 else
0132 EErms[s] = 1.e-06;
0133 pEEmin[s] = EEmean[s] - kRMS * EErms[s];
0134 pEEmax[s] = EEmean[s] + kRMS * EErms[s];
0135 std::cout << "EE sample " << s << " mean " << EEmean[s] << " rms " << EErms[s] << " entries " << EEtot[s]
0136 << " min " << pEEmin[s] << " max " << pEEmax[s] << std::endl;
0137
0138 }
0139
0140
0141
0142 gStyle->SetPalette(1);
0143 gStyle->SetOptStat(0);
0144 TCanvas canvas("CC map", "CC map", 1600, 5600);
0145 TLatex t1;
0146 t1.SetNDC();
0147 t1.SetTextAlign(26);
0148 t1.SetTextSize(0.05);
0149 t1.DrawLatex(0.5, 0.98, Form("Ecal PulseShapes, IOV %i", run));
0150
0151 float xmi[3] = {0.0, 0.24, 0.76};
0152 float xma[3] = {0.24, 0.76, 1.00};
0153 TPad*** pad = new TPad**[6];
0154 for (int s = 0; s < TEMPLATESAMPLES; s++) {
0155 pad[s] = new TPad*[3];
0156 for (int obj = 0; obj < 3; obj++) {
0157 float yma = 0.96 - (0.08 * s);
0158 float ymi = yma - 0.07;
0159 pad[s][obj] = new TPad(Form("p_%i_%i", obj, s), Form("p_%i_%i", obj, s), xmi[obj], ymi, xma[obj], yma);
0160 pad[s][obj]->Draw();
0161 }
0162 }
0163
0164 int ipad = 0;
0165 for (int s = 0; s < TEMPLATESAMPLES; s++) {
0166 pad[ipad][0]->cd();
0167 if (pEBmin[s] == pEBmax[s]) {
0168 pEBmin[s] = pEBmin[s] - 1.e-06;
0169 pEBmax[s] = pEBmax[s] + 1.e-06;
0170 }
0171 if (pEEmin[s] == pEEmax[s]) {
0172 pEEmin[s] = pEEmin[s] - 1.e-06;
0173 pEEmax[s] = pEEmax[s] + 1.e-06;
0174 }
0175 DrawEE(endc_m[s], pEEmin[s], pEEmax[s]);
0176 pad[ipad][1]->cd();
0177 DrawEB(barrel[s], pEBmin[s], pEBmax[s]);
0178 pad[ipad][2]->cd();
0179 DrawEE(endc_p[s], pEEmin[s], pEEmax[s]);
0180 ipad++;
0181 }
0182
0183 std::string ImageName(m_imageFileName);
0184 canvas.SaveAs(ImageName.c_str());
0185
0186 return true;
0187 }
0188 };
0189
0190
0191
0192
0193 class EcalPulseShapesProfile : public cond::payloadInspector::PlotImage<EcalPulseShapes> {
0194 public:
0195 EcalPulseShapesProfile() : cond::payloadInspector::PlotImage<EcalPulseShapes>("ECAL PulseShapes - Profile ") {
0196 setSingleIov(true);
0197 }
0198
0199 bool fill(const std::vector<std::tuple<cond::Time_t, cond::Hash> >& iovs) override {
0200 TProfile* barrel = new TProfile("EB", "EB profile", TEMPLATESAMPLES, 0, TEMPLATESAMPLES);
0201 TProfile* endcap = new TProfile("EE", "EE profile", TEMPLATESAMPLES, 0, TEMPLATESAMPLES);
0202
0203 auto iov = iovs.front();
0204 std::shared_ptr<EcalPulseShapes> payload = fetchPayload(std::get<1>(iov));
0205 unsigned int run = std::get<0>(iov);
0206 if (payload.get()) {
0207 for (int ieta = -MAX_IETA; ieta <= MAX_IETA; ieta++) {
0208 for (int iphi = 1; iphi <= MAX_IPHI; iphi++) {
0209 EBDetId id(ieta, iphi);
0210 EcalPulseShape pulse = (*payload)[id.rawId()];
0211 for (int s = 0; s < TEMPLATESAMPLES; s++) {
0212 double val = pulse.pdfval[s];
0213 barrel->Fill(s, val);
0214 }
0215 }
0216 }
0217
0218 for (int sign = 0; sign < kSides; sign++) {
0219 int thesign = sign == 1 ? 1 : -1;
0220 for (int ix = 1; ix <= IX_MAX; ix++) {
0221 for (int iy = 1; iy <= IY_MAX; iy++) {
0222 if (!EEDetId::validDetId(ix, iy, thesign))
0223 continue;
0224 EEDetId id(ix, iy, thesign);
0225 EcalPulseShape pulse = (*payload)[id.rawId()];
0226 for (int s = 0; s < TEMPLATESAMPLES; s++) {
0227 double val = pulse.pdfval[s];
0228 endcap->Fill(s, val);
0229 }
0230 }
0231 }
0232 }
0233 }
0234
0235 gStyle->SetPalette(1);
0236 gStyle->SetOptStat(0);
0237 TCanvas canvas("CC map", "CC map", 500, 1000);
0238 TLatex t1;
0239 t1.SetNDC();
0240 t1.SetTextAlign(26);
0241 t1.SetTextSize(0.05);
0242 t1.DrawLatex(0.5, 0.98, Form("Ecal PulseShapes, IOV %i", run));
0243
0244 TPad** pad = new TPad*[2];
0245 for (int s = 0; s < 2; s++) {
0246 float yma = 0.96 - (0.47 * s);
0247 float ymi = yma - 0.45;
0248 pad[s] = new TPad(Form("p_%i", s), Form("p_%i", s), 0.0, ymi, 1.0, yma);
0249 pad[s]->Draw();
0250 }
0251
0252 pad[0]->cd();
0253 barrel->Draw("");
0254 pad[1]->cd();
0255 endcap->Draw("");
0256
0257 std::string ImageName(m_imageFileName);
0258 canvas.SaveAs(ImageName.c_str());
0259 return true;
0260 }
0261 };
0262
0263 }
0264
0265
0266 PAYLOAD_INSPECTOR_MODULE(EcalPulseShapes) {
0267 PAYLOAD_INSPECTOR_CLASS(EcalPulseShapesPlot);
0268 PAYLOAD_INSPECTOR_CLASS(EcalPulseShapesProfile);
0269 }