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
0005
0006 #include "CondFormats/ESObjects/interface/ESIntercalibConstants.h"
0007 #include "DataFormats/EcalDetId/interface/ESDetId.h"
0008 #include "CondCore/EcalPlugins/plugins/ESDrawUtils.h"
0009
0010 #include <memory>
0011 #include <sstream>
0012
0013 #include "TStyle.h"
0014 #include "TH2F.h"
0015 #include "TCanvas.h"
0016 #include "TLine.h"
0017 #include "TLatex.h"
0018
0019 namespace {
0020 enum { kESChannels = 137216 };
0021 enum { IX_MIN = 1, IY_MIN = 1, IX_MAX = 40, IY_MAX = 40 };
0022
0023
0024
0025
0026 class ESIntercalibConstantsPlot : public cond::payloadInspector::PlotImage<ESIntercalibConstants> {
0027 public:
0028 ESIntercalibConstantsPlot() : cond::payloadInspector::PlotImage<ESIntercalibConstants>("ES IntercalibConstants") {
0029 setSingleIov(true);
0030 }
0031 bool fill(const std::vector<std::tuple<cond::Time_t, cond::Hash> >& iovs) override {
0032 TH2F*** esmap = new TH2F**[2];
0033 std::string title[2][2] = {{"ES+F", "ES-F"}, {"ES+R", "ES-R"}};
0034 for (int plane = 0; plane < 2; plane++) {
0035 esmap[plane] = new TH2F*[2];
0036 for (int side = 0; side < 2; side++)
0037 esmap[plane][side] = new TH2F(
0038 Form("esmap%i%i", plane, side), title[plane][side].c_str(), IX_MAX, 0, IX_MAX, IY_MAX, 0, IY_MAX);
0039 }
0040 unsigned int run = 0;
0041 float valmin = 999.;
0042 auto iov = iovs.front();
0043 std::shared_ptr<ESIntercalibConstants> payload = fetchPayload(std::get<1>(iov));
0044 run = std::get<0>(iov);
0045 if (payload.get()) {
0046
0047 for (int id = 0; id < kESChannels; id++)
0048 if (ESDetId::validHashIndex(id)) {
0049 ESDetId myESId = ESDetId::unhashIndex(id);
0050 int side = myESId.zside();
0051 if (side < 0)
0052 side = 1;
0053 else
0054 side = 0;
0055 int plane = myESId.plane() - 1;
0056 if (side < 0 || side > 1 || plane < 0 || plane > 1) {
0057 std::cout << " channel " << id << " side " << myESId.zside() << " plane " << myESId.plane() << std::endl;
0058 return false;
0059 }
0060 ESIntercalibConstants::const_iterator IC_it = payload->find(myESId);
0061 float value = *IC_it;
0062
0063 if (myESId.strip() == 1) {
0064 esmap[plane][side]->Fill(myESId.six() - 1, myESId.siy() - 1, value);
0065 if (value < valmin)
0066 valmin = value;
0067 }
0068 }
0069 }
0070
0071 gStyle->SetOptStat(0);
0072 gStyle->SetPalette(1);
0073 TCanvas canvas("CC map", "CC map", 1680, 1320);
0074 TLatex t1;
0075 t1.SetNDC();
0076 t1.SetTextAlign(26);
0077 t1.SetTextSize(0.05);
0078 t1.DrawLatex(0.5, 0.96, Form("ES Intercalib Constants, IOV %i", run));
0079 t1.SetTextSize(0.025);
0080
0081 float xmi[2] = {0.0, 0.5};
0082 float xma[2] = {0.5, 1.0};
0083 TPad*** pad = new TPad**[2];
0084 for (int plane = 0; plane < 2; plane++) {
0085 pad[plane] = new TPad*[2];
0086 for (int side = 0; side < 2; side++) {
0087 float yma = 0.94 - (0.46 * plane);
0088 float ymi = yma - 0.44;
0089 pad[plane][side] =
0090 new TPad(Form("p_%i_%i", plane, side), Form("p_%i_%i", plane, side), xmi[side], ymi, xma[side], yma);
0091 pad[plane][side]->Draw();
0092 }
0093 }
0094
0095 int min = valmin * 10 - 1.;
0096 valmin = (float)min / 10.;
0097 for (int side = 0; side < 2; side++) {
0098 for (int plane = 0; plane < 2; plane++) {
0099 pad[plane][side]->cd();
0100 esmap[plane][side]->Draw("colz1");
0101 esmap[plane][side]->SetMinimum(valmin);
0102 DrawES(plane, side);
0103 }
0104 }
0105
0106 std::string ImageName(m_imageFileName);
0107 canvas.SaveAs(ImageName.c_str());
0108 return true;
0109 }
0110 };
0111
0112
0113
0114
0115 template <cond::payloadInspector::IOVMultiplicity nIOVs, int ntags>
0116 class ESIntercalibConstantsDiffBase : public cond::payloadInspector::PlotImage<ESIntercalibConstants, nIOVs, ntags> {
0117 public:
0118 ESIntercalibConstantsDiffBase()
0119 : cond::payloadInspector::PlotImage<ESIntercalibConstants, nIOVs, ntags>("ES IntercalibConstants difference") {}
0120 bool fill() override {
0121 TH2F*** esmap = new TH2F**[2];
0122 std::string title[2][2] = {{"ES+F", "ES-F"}, {"ES+R", "ES-R"}};
0123 for (int plane = 0; plane < 2; plane++) {
0124 esmap[plane] = new TH2F*[2];
0125 for (int side = 0; side < 2; side++)
0126 esmap[plane][side] = new TH2F(
0127 Form("esmap%i%i", plane, side), title[plane][side].c_str(), IX_MAX, 0, IX_MAX, IY_MAX, 0, IY_MAX);
0128 }
0129 unsigned int run[2];
0130 std::string l_tagname[2];
0131 float val[kESChannels], valmin = 999.;
0132 auto iovs = cond::payloadInspector::PlotBase::getTag<0>().iovs;
0133 l_tagname[0] = cond::payloadInspector::PlotBase::getTag<0>().name;
0134 auto firstiov = iovs.front();
0135 run[0] = std::get<0>(firstiov);
0136 std::tuple<cond::Time_t, cond::Hash> lastiov;
0137 if (ntags == 2) {
0138 auto tag2iovs = cond::payloadInspector::PlotBase::getTag<1>().iovs;
0139 l_tagname[1] = cond::payloadInspector::PlotBase::getTag<1>().name;
0140 lastiov = tag2iovs.front();
0141 } else {
0142 lastiov = iovs.back();
0143 l_tagname[1] = l_tagname[0];
0144 }
0145 run[1] = std::get<0>(lastiov);
0146 for (int irun = 0; irun < nIOVs; irun++) {
0147 std::shared_ptr<ESIntercalibConstants> payload;
0148 if (irun == 0) {
0149 payload = this->fetchPayload(std::get<1>(firstiov));
0150 } else {
0151 payload = this->fetchPayload(std::get<1>(lastiov));
0152 }
0153 if (payload.get()) {
0154 for (int id = 0; id < kESChannels; id++)
0155 if (ESDetId::validHashIndex(id)) {
0156 ESDetId myESId = ESDetId::unhashIndex(id);
0157
0158 ESIntercalibConstants::const_iterator IC_it = payload->find(myESId);
0159 float value = *IC_it;
0160
0161 if (irun == 0)
0162 val[id] = value;
0163 else {
0164 int side = myESId.zside();
0165 if (side < 0)
0166 side = 1;
0167 else
0168 side = 0;
0169 int plane = myESId.plane() - 1;
0170 if (side < 0 || side > 1 || plane < 0 || plane > 1) {
0171 std::cout << " channel " << id << " side " << myESId.zside() << " plane " << myESId.plane()
0172 << std::endl;
0173 return false;
0174 }
0175
0176 int diff = value - val[id];
0177 if (myESId.strip() == 1) {
0178 esmap[plane][side]->Fill(myESId.six() - 1, myESId.siy() - 1, diff);
0179 if (diff < valmin)
0180 valmin = diff;
0181 }
0182 }
0183 }
0184 }
0185 }
0186
0187 gStyle->SetOptStat(0);
0188 gStyle->SetPalette(1);
0189 TCanvas canvas("CC map", "CC map", 1680, 1320);
0190 TLatex t1;
0191 t1.SetNDC();
0192 t1.SetTextAlign(26);
0193 int len = l_tagname[0].length() + l_tagname[1].length();
0194 if (ntags == 2 && len < 58) {
0195 t1.SetTextSize(0.025);
0196 t1.DrawLatex(
0197 0.5, 0.96, Form("%s IOV %i - %s IOV %i", l_tagname[1].c_str(), run[1], l_tagname[0].c_str(), run[0]));
0198 } else {
0199 t1.SetTextSize(0.03);
0200 t1.DrawLatex(0.5, 0.96, Form("ES Intercalib Constants, IOV %i - %i", run[1], run[0]));
0201 }
0202 float xmi[2] = {0.0, 0.5};
0203 float xma[2] = {0.5, 1.0};
0204 TPad*** pad = new TPad**[2];
0205 for (int plane = 0; plane < 2; plane++) {
0206 pad[plane] = new TPad*[2];
0207 for (int side = 0; side < 2; side++) {
0208 float yma = 0.94 - (0.46 * plane);
0209 float ymi = yma - 0.44;
0210 pad[plane][side] =
0211 new TPad(Form("p_%i_%i", plane, side), Form("p_%i_%i", plane, side), xmi[side], ymi, xma[side], yma);
0212 pad[plane][side]->Draw();
0213 }
0214 }
0215
0216 int min = valmin * 10 - 1.;
0217 valmin = (float)min / 10.;
0218 for (int side = 0; side < 2; side++) {
0219 for (int plane = 0; plane < 2; plane++) {
0220 pad[plane][side]->cd();
0221 esmap[plane][side]->Draw("colz1");
0222 esmap[plane][side]->SetMinimum(valmin);
0223 DrawES(plane, side);
0224 }
0225 }
0226
0227 std::string ImageName(this->m_imageFileName);
0228 canvas.SaveAs(ImageName.c_str());
0229 return true;
0230 }
0231 };
0232 using ESIntercalibConstantsDiffOneTag = ESIntercalibConstantsDiffBase<cond::payloadInspector::SINGLE_IOV, 1>;
0233 using ESIntercalibConstantsDiffTwoTags = ESIntercalibConstantsDiffBase<cond::payloadInspector::SINGLE_IOV, 2>;
0234
0235 }
0236
0237
0238 PAYLOAD_INSPECTOR_MODULE(ESIntercalibConstants) {
0239 PAYLOAD_INSPECTOR_CLASS(ESIntercalibConstantsPlot);
0240 PAYLOAD_INSPECTOR_CLASS(ESIntercalibConstantsDiffOneTag);
0241 PAYLOAD_INSPECTOR_CLASS(ESIntercalibConstantsDiffTwoTags);
0242 }