CMS 3D CMS Logo

EcalTPGLinearizationConst_PayloadInspector.cc
Go to the documentation of this file.
7 
8 // the data format of the condition to be inspected
10 
11 #include "TH2F.h"
12 #include "TCanvas.h"
13 #include "TLine.h"
14 #include "TStyle.h"
15 #include "TLatex.h"
16 #include "TPave.h"
17 #include "TPaveStats.h"
18 #include <string>
19 #include <fstream>
20 
21 namespace {
22  enum {kEBChannels = 61200, kEEChannels = 14648, kGains = 3, kSides = 2};
23  enum {MIN_IETA = 1, MIN_IPHI = 1, MAX_IETA = 85, MAX_IPHI = 360}; // barrel lower and upper bounds on eta and phi
24  enum {IX_MIN = 1, IY_MIN = 1, IX_MAX = 100, IY_MAX = 100}; // endcaps lower and upper bounds on x and y
25  int gainValues[kGains] = {12, 6, 1};
26 
27  /**************************************************
28  2d plot of ECAL TPGLinearizationConst of 1 IOV
29  **************************************************/
30  class EcalTPGLinearizationConstPlot : public cond::payloadInspector::PlotImage<EcalTPGLinearizationConst> {
31  public:
32  EcalTPGLinearizationConstPlot() : cond::payloadInspector::PlotImage<EcalTPGLinearizationConst>("ECAL Gain Ratios - map ") {
33  setSingleIov( true );
34  }
35 
36  bool fill( const std::vector<std::tuple<cond::Time_t,cond::Hash> >& iovs ) override{
37  TH2F** barrel_m = new TH2F*[kGains];
38  TH2F** endc_p_m = new TH2F*[kGains];
39  TH2F** endc_m_m = new TH2F*[kGains];
40  TH2F** barrel_r = new TH2F*[kGains];
41  TH2F** endc_p_r = new TH2F*[kGains];
42  TH2F** endc_m_r = new TH2F*[kGains];
43  float mEBmin[kGains], mEEmin[kGains], mEBmax[kGains], mEEmax[kGains], rEBmin[kGains], rEEmin[kGains], rEBmax[kGains], rEEmax[kGains];
44  for (int gainId = 0; gainId < kGains; gainId++) {
45  barrel_m[gainId] = new TH2F(Form("EBm%i", gainId), Form("EB mult_x%i ", gainValues[gainId]), MAX_IPHI, 0, MAX_IPHI, 2 * MAX_IETA, -MAX_IETA, MAX_IETA);
46  endc_p_m[gainId] = new TH2F(Form("EE+m%i", gainId), Form("EE+ mult_x%i", gainValues[gainId]), IX_MAX, IX_MIN, IX_MAX + 1, IY_MAX, IY_MIN, IY_MAX + 1);
47  endc_m_m[gainId] = new TH2F(Form("EE-m%i", gainId), Form("EE- mult_x%i", gainValues[gainId]), IX_MAX, IX_MIN, IX_MAX + 1, IY_MAX, IY_MIN, IY_MAX + 1);
48  barrel_r[gainId] = new TH2F(Form("EBr%i", gainId), Form("EB shift_x%i", gainValues[gainId]), MAX_IPHI, 0, MAX_IPHI, 2 * MAX_IETA, -MAX_IETA, MAX_IETA);
49  endc_p_r[gainId] = new TH2F(Form("EE+r%i",gainId), Form("EE+ shift_x%i", gainValues[gainId]), IX_MAX, IX_MIN, IX_MAX + 1, IY_MAX, IY_MIN, IY_MAX + 1);
50  endc_m_r[gainId] = new TH2F(Form("EE-r%i",gainId), Form("EE- shift_x%i", gainValues[gainId]), IX_MAX, IX_MIN, IX_MAX + 1, IY_MAX, IY_MIN, IY_MAX + 1);
51  mEBmin[gainId] = 10.;
52  mEEmin[gainId] = 10.;
53  mEBmax[gainId] = -10.;
54  mEEmax[gainId] = -10.;
55  rEBmin[gainId] = 10.;
56  rEEmin[gainId] = 10.;
57  rEBmax[gainId] = -10.;
58  rEEmax[gainId] = -10.;
59  }
60 
61  // std::ofstream fout;
62  // fout.open("./bid.txt");
63  auto iov = iovs.front();
64  std::shared_ptr<EcalTPGLinearizationConst> payload = fetchPayload( std::get<1>(iov) );
65  unsigned int run = std::get<0>(iov);
66  if( payload.get() ){
67  for (int sign=0; sign < kSides; sign++) {
68  int thesign = sign==1 ? 1:-1;
69 
70  for (int ieta = 0; ieta < MAX_IETA; ieta++) {
71  for (int iphi = 0; iphi < MAX_IPHI; iphi++) {
72  EBDetId id((ieta+1)*thesign, iphi+1);
73  float y = -1 - ieta;
74  if(sign == 1) y = ieta;
75  float val = (*payload)[id.rawId()].mult_x12;
76  barrel_m[0]->Fill(iphi, y, val);
77  if(val < mEBmin[0]) mEBmin[0] = val;
78  if(val > mEBmax[0]) mEBmax[0] = val;
79  val = (*payload)[id.rawId()].shift_x12;
80  barrel_r[0]->Fill(iphi, y, val);
81  if(val < rEBmin[0]) rEBmin[0] = val;
82  if(val > rEBmax[0]) rEBmax[0] = val;
83  val = (*payload)[id.rawId()].mult_x6;
84  barrel_m[1]->Fill(iphi, y, val);
85  if(val < mEBmin[1]) mEBmin[1] = val;
86  if(val > mEBmax[1]) mEBmax[1] = val;
87  val = (*payload)[id.rawId()].shift_x6;
88  barrel_r[1]->Fill(iphi, y, val);
89  if(val < rEBmin[1]) rEBmin[1] = val;
90  if(val > rEBmax[1]) rEBmax[1] = val;
91  val = (*payload)[id.rawId()].mult_x1;
92  barrel_m[2]->Fill(iphi, y, val);
93  if(val < mEBmin[2]) mEBmin[2] = val;
94  if(val > mEBmax[2]) mEBmax[2] = val;
95  val = (*payload)[id.rawId()].shift_x1;
96  barrel_r[2]->Fill(iphi, y, val);
97  if(val < rEBmin[2]) rEBmin[2] = val;
98  if(val > rEBmax[2]) rEBmax[2] = val;
99  } // iphi
100  } // ieta
101 
102  for (int ix = 0; ix < IX_MAX; ix++) {
103  for (int iy = 0; iy < IY_MAX; iy++) {
104  if (! EEDetId::validDetId(ix+1,iy+1,thesign)) continue;
105  EEDetId id(ix+1,iy+1,thesign);
106  float val = (*payload)[id.rawId()].mult_x12;
107  if (thesign==1) endc_p_m[0]->Fill(ix + 1, iy + 1, val);
108  else endc_m_m[0]->Fill(ix + 1, iy + 1, val);
109  if(val < mEEmin[0]) mEEmin[0] = val;
110  if(val > mEEmax[0]) mEEmax[0] = val;
111  val = (*payload)[id.rawId()].shift_x12;
112  if (thesign==1) endc_p_r[0]->Fill(ix + 1, iy + 1, val);
113  else endc_m_r[0]->Fill(ix + 1, iy + 1, val);
114  if(val < rEEmin[0]) rEEmin[0] = val;
115  if(val > rEEmax[0]) rEEmax[0] = val;
116  val = (*payload)[id.rawId()].mult_x6;
117  if (thesign==1) endc_p_m[1]->Fill(ix + 1, iy + 1, val);
118  else endc_m_m[1]->Fill(ix + 1, iy + 1, val);
119  if(val < mEEmin[1]) mEEmin[1] = val;
120  if(val > mEEmax[1]) mEEmax[1] = val;
121  val = (*payload)[id.rawId()].shift_x6;
122  if (thesign==1) endc_p_r[1]->Fill(ix + 1, iy + 1, val);
123  else endc_m_r[1]->Fill(ix + 1, iy + 1, val);
124  if(val < rEEmin[1]) rEEmin[1] = val;
125  if(val > rEEmax[1]) rEEmax[1] = val;
126  val = (*payload)[id.rawId()].mult_x1;
127  if (thesign==1) endc_p_m[2]->Fill(ix + 1, iy + 1, val);
128  else endc_m_m[2]->Fill(ix + 1, iy + 1, val);
129  if(val < mEEmin[2]) mEEmin[2] = val;
130  if(val > mEEmax[2]) mEEmax[2] = val;
131  val = (*payload)[id.rawId()].shift_x1;
132  if (thesign==1) endc_p_r[2]->Fill(ix + 1, iy + 1, val);
133  else endc_m_r[2]->Fill(ix + 1, iy + 1, val);
134  if(val < rEEmin[2]) rEEmin[2] = val;
135  if(val > rEEmax[2]) rEEmax[2] = val;
136  // fout << " x " << ix << " y " << " val " << val << std::endl;
137  } // iy
138  } // ix
139  } // side
140  } // if payload.get()
141  else return false;
142  // std::cout << " min " << rEEmin[2] << " max " << rEEmax[2] << std::endl;
143  // fout.close();
144 
145  gStyle->SetPalette(1);
146  gStyle->SetOptStat(0);
147  TCanvas canvas("CC map","CC map",1200,1800);
148  TLatex t1;
149  t1.SetNDC();
150  t1.SetTextAlign(26);
151  t1.SetTextSize(0.05);
152  t1.DrawLatex(0.5, 0.96, Form("Ecal Gain TPGLinearizationConst, IOV %i", run));
153 
154  float xmi[3] = {0.0 , 0.22, 0.78};
155  float xma[3] = {0.22, 0.78, 1.00};
156  TPad*** pad = new TPad**[6];
157  for (int gId = 0; gId < 6; gId++) {
158  pad[gId] = new TPad*[3];
159  for (int obj = 0; obj < 3; obj++) {
160  float yma = 0.94 - (0.16 * gId);
161  float ymi = yma - 0.14;
162  pad[gId][obj] = new TPad(Form("p_%i_%i", obj, gId),Form("p_%i_%i", obj, gId),
163  xmi[obj], ymi, xma[obj], yma);
164  pad[gId][obj]->Draw();
165  }
166  }
167 
168  for (int gId = 0; gId < kGains; gId++) {
169  pad[gId][0]->cd();
170  DrawEE(endc_m_m[gId], mEEmin[gId], mEEmax[gId]);
171  pad[gId + 3][0]->cd();
172  DrawEE(endc_m_r[gId], rEEmin[gId], rEEmax[gId]);
173  pad[gId][1]->cd();
174  DrawEB(barrel_m[gId], mEBmin[gId], mEBmax[gId]);
175  pad[gId + 3][1]->cd();
176  DrawEB(barrel_r[gId], rEBmin[gId], rEBmax[gId]);
177  pad[gId][2]->cd();
178  DrawEE(endc_p_m[gId], mEEmin[gId], mEEmax[gId]);
179  pad[gId + 3][2]->cd();
180  DrawEE(endc_p_r[gId], rEEmin[gId], rEEmax[gId]);
181  }
182 
183  std::string ImageName(m_imageFileName);
184  canvas.SaveAs(ImageName.c_str());
185  return true;
186  }// fill method
187  };
188 
189  /******************************************************************
190  2d plot of ECAL TPGLinearizationConst difference between 2 IOVs
191  ******************************************************************/
192  class EcalTPGLinearizationConstDiff : public cond::payloadInspector::PlotImage<EcalTPGLinearizationConst> {
193 
194  public:
195  EcalTPGLinearizationConstDiff() : cond::payloadInspector::PlotImage<EcalTPGLinearizationConst>("ECAL Gain Ratios difference") {
196  setSingleIov(false);
197  }
198 
199  bool fill( const std::vector<std::tuple<cond::Time_t,cond::Hash> >& iovs ) override{
200  TH2F** barrel_m = new TH2F*[kGains];
201  TH2F** endc_p_m = new TH2F*[kGains];
202  TH2F** endc_m_m = new TH2F*[kGains];
203  TH2F** barrel_r = new TH2F*[kGains];
204  TH2F** endc_p_r = new TH2F*[kGains];
205  TH2F** endc_m_r = new TH2F*[kGains];
206  float mEBmin[kGains], mEEmin[kGains], mEBmax[kGains], mEEmax[kGains], rEBmin[kGains], rEEmin[kGains], rEBmax[kGains], rEEmax[kGains];
207  float mEB[kGains][kEBChannels], mEE[kGains][kEEChannels], rEB[kGains][kEBChannels], rEE[kGains][kEEChannels];
208  for (int gainId = 0; gainId < kGains; gainId++) {
209  barrel_m[gainId] = new TH2F(Form("EBm%i", gainId), Form("EB mult_x%i ", gainValues[gainId]), MAX_IPHI, 0, MAX_IPHI, 2 * MAX_IETA, -MAX_IETA, MAX_IETA);
210  endc_p_m[gainId] = new TH2F(Form("EE+m%i", gainId), Form("EE+ mult_x%i", gainValues[gainId]), IX_MAX, IX_MIN, IX_MAX + 1, IY_MAX, IY_MIN, IY_MAX + 1);
211  endc_m_m[gainId] = new TH2F(Form("EE-m%i", gainId), Form("EE- mult_x%i", gainValues[gainId]), IX_MAX, IX_MIN, IX_MAX + 1, IY_MAX, IY_MIN, IY_MAX + 1);
212  barrel_r[gainId] = new TH2F(Form("EBr%i", gainId), Form("EB shift_x%i", gainValues[gainId]), MAX_IPHI, 0, MAX_IPHI, 2 * MAX_IETA, -MAX_IETA, MAX_IETA);
213  endc_p_r[gainId] = new TH2F(Form("EE+r%i",gainId), Form("EE+ shift_x%i", gainValues[gainId]), IX_MAX, IX_MIN, IX_MAX + 1, IY_MAX, IY_MIN, IY_MAX + 1);
214  endc_m_r[gainId] = new TH2F(Form("EE-r%i",gainId), Form("EE- shift_x%i", gainValues[gainId]), IX_MAX, IX_MIN, IX_MAX + 1, IY_MAX, IY_MIN, IY_MAX + 1);
215  mEBmin[gainId] = 10.;
216  mEEmin[gainId] = 10.;
217  mEBmax[gainId] = -10.;
218  mEEmax[gainId] = -10.;
219  rEBmin[gainId] = 10.;
220  rEEmin[gainId] = 10.;
221  rEBmax[gainId] = -10.;
222  rEEmax[gainId] = -10.;
223  }
224 
225  unsigned int run[2], irun = 0;
226  //float gEB[3][kEBChannels], gEE[3][kEEChannels];
227  for ( auto const & iov: iovs) {
228  std::shared_ptr<EcalTPGLinearizationConst> payload = fetchPayload( std::get<1>(iov) );
229  run[irun] = std::get<0>(iov);
230  if( payload.get() ){
231  for (int sign=0; sign < kSides; sign++) {
232  int thesign = sign==1 ? 1:-1;
233 
234  for (int ieta = 0; ieta < MAX_IETA; ieta++) {
235  for (int iphi = 0; iphi < MAX_IPHI; iphi++) {
236  EBDetId id((ieta+1)*thesign, iphi+1);
237  int hashindex = id.hashedIndex();
238  float y = -1 - ieta;
239  if(sign == 1) y = ieta;
240  float val = (*payload)[id.rawId()].mult_x12;
241  if(irun == 0) {
242  mEB[0][hashindex] = val;
243  }
244  else {
245  float diff = val - mEB[0][hashindex];
246  barrel_m[0]->Fill(iphi, y, diff);
247  if(diff < mEBmin[0]) mEBmin[0] = diff;
248  if(diff > mEBmax[0]) mEBmax[0] = diff;
249  }
250  val = (*payload)[id.rawId()].shift_x12;
251  if(irun == 0) {
252  rEB[0][hashindex] = val;
253  }
254  else {
255  float diff = val - rEB[0][hashindex];
256  barrel_r[0]->Fill(iphi, y, diff);
257  if(diff < rEBmin[0]) rEBmin[0] = diff;
258  if(diff > rEBmax[0]) rEBmax[0] = diff;
259  }
260  val = (*payload)[id.rawId()].mult_x6;
261  if(irun == 0) {
262  mEB[1][hashindex] = val;
263  }
264  else {
265  float diff = val - mEB[1][hashindex];
266  barrel_m[1]->Fill(iphi, y, diff);
267  if(diff < mEBmin[1]) mEBmin[1] = diff;
268  if(diff > mEBmax[1]) mEBmax[1] = diff;
269  }
270  val = (*payload)[id.rawId()].shift_x6;
271  if(irun == 0) {
272  rEB[1][hashindex] = val;
273  }
274  else {
275  float diff = val - rEB[1][hashindex];
276  barrel_r[1]->Fill(iphi, y, diff);
277  if(diff < rEBmin[1]) rEBmin[1] = diff;
278  if(diff > rEBmax[1]) rEBmax[1] = diff;
279  }
280  val = (*payload)[id.rawId()].mult_x1;
281  if(irun == 0) {
282  mEB[2][hashindex] = val;
283  }
284  else {
285  float diff = val - mEB[2][hashindex];
286  barrel_m[2]->Fill(iphi, y, diff);
287  if(diff < mEBmin[2]) mEBmin[2] = diff;
288  if(diff > mEBmax[2]) mEBmax[2] = diff;
289  }
290  val = (*payload)[id.rawId()].shift_x1;
291  if(irun == 0) {
292  rEB[2][hashindex] = val;
293  }
294  else {
295  float diff = val - rEB[2][hashindex];
296  barrel_r[2]->Fill(iphi, y, diff);
297  if(diff < rEBmin[2]) rEBmin[2] = diff;
298  if(diff > rEBmax[2]) rEBmax[2] = diff;
299  }
300  } // iphi
301  } // ieta
302 
303  for (int ix = 0; ix < IX_MAX; ix++) {
304  for (int iy = 0; iy < IY_MAX; iy++) {
305  if (! EEDetId::validDetId(ix+1,iy+1,thesign)) continue;
306  EEDetId id(ix+1,iy+1,thesign);
307  int hashindex = id.hashedIndex();
308  float val = (*payload)[id.rawId()].mult_x12;
309  if(irun == 0) {
310  mEE[0][hashindex] = val;
311  }
312  else {
313  float diff = val - mEE[0][hashindex];
314  if (thesign==1) endc_p_m[0]->Fill(ix + 1, iy + 1, diff);
315  else endc_m_m[0]->Fill(ix + 1, iy + 1, diff);
316  if(diff < mEEmin[0]) mEEmin[0] = diff;
317  if(diff > mEEmax[0]) mEEmax[0] = diff;
318  }
319  val = (*payload)[id.rawId()].shift_x12;
320  if(irun == 0) {
321  rEE[0][hashindex] = val;
322  }
323  else {
324  float diff = val - rEE[0][hashindex];
325  if (thesign==1) endc_p_r[0]->Fill(ix + 1, iy + 1, diff);
326  else endc_m_r[0]->Fill(ix + 1, iy + 1, diff);
327  if(diff < rEEmin[0]) rEEmin[0] = diff;
328  if(diff > rEEmax[0]) rEEmax[0] = diff;
329  }
330  val = (*payload)[id.rawId()].mult_x6;
331  if(irun == 0) {
332  mEE[1][hashindex] = val;
333  }
334  else {
335  float diff = val - mEE[1][hashindex];
336  if (thesign==1) endc_p_m[1]->Fill(ix + 1, iy + 1, diff);
337  else endc_m_m[1]->Fill(ix + 1, iy + 1, diff);
338  if(diff < mEEmin[1]) mEEmin[1] = diff;
339  if(diff > mEEmax[1]) mEEmax[1] = diff;
340  }
341  val = (*payload)[id.rawId()].shift_x6;
342  if(irun == 0) {
343  rEE[1][hashindex] = val;
344  }
345  else {
346  float diff = val - rEE[1][hashindex];
347  if (thesign==1) endc_p_r[1]->Fill(ix + 1, iy + 1, diff);
348  else endc_m_r[1]->Fill(ix + 1, iy + 1, diff);
349  if(diff < rEEmin[1]) rEEmin[1] = diff;
350  if(diff > rEEmax[1]) rEEmax[1] = diff;
351  }
352  val = (*payload)[id.rawId()].mult_x1;
353  if(irun == 0) {
354  mEE[2][hashindex] = val;
355  }
356  else {
357  float diff = val - mEE[2][hashindex];
358  if (thesign==1) endc_p_m[2]->Fill(ix + 1, iy + 1, diff);
359  else endc_m_m[2]->Fill(ix + 1, iy + 1, diff);
360  if(diff < mEEmin[2]) mEEmin[2] = diff;
361  if(diff > mEEmax[2]) mEEmax[2] = diff;
362  }
363  val = (*payload)[id.rawId()].shift_x1;
364  if(irun == 0) {
365  rEE[2][hashindex] = val;
366  }
367  else {
368  float diff = val - rEE[2][hashindex];
369  if (thesign==1) endc_p_r[2]->Fill(ix + 1, iy + 1, diff);
370  else endc_m_r[2]->Fill(ix + 1, iy + 1, diff);
371  if(diff < rEEmin[2]) rEEmin[2] = diff;
372  if(diff > rEEmax[2]) rEEmax[2] = diff;
373  }
374  // fout << " x " << ix << " y " << " diff " << diff << std::endl;
375  } // iy
376  } // ix
377  } // side
378  } // if payload.get()
379  else return false;
380  irun++;
381  } // loop over IOVs
382 
383  gStyle->SetPalette(1);
384  gStyle->SetOptStat(0);
385  TCanvas canvas("CC map","CC map",1200,1800);
386  TLatex t1;
387  t1.SetNDC();
388  t1.SetTextAlign(26);
389  t1.SetTextSize(0.05);
390  t1.DrawLatex(0.5, 0.96, Form("Ecal TPGLinearizationConst, IOV %i - %i", run[1], run[0]));
391 
392  float xmi[3] = {0.0 , 0.22, 0.78};
393  float xma[3] = {0.22, 0.78, 1.00};
394  TPad*** pad = new TPad**[6];
395  for (int gId = 0; gId < 6; gId++) {
396  pad[gId] = new TPad*[3];
397  for (int obj = 0; obj < 3; obj++) {
398  float yma = 0.94 - (0.16 * gId);
399  float ymi = yma - 0.14;
400  pad[gId][obj] = new TPad(Form("p_%i_%i", obj, gId),Form("p_%i_%i", obj, gId),
401  xmi[obj], ymi, xma[obj], yma);
402  pad[gId][obj]->Draw();
403  }
404  }
405 
406  for (int gId = 0; gId < kGains; gId++) {
407  pad[gId][0]->cd();
408  DrawEE(endc_m_m[gId], mEEmin[gId], mEEmax[gId]);
409  pad[gId + 3][0]->cd();
410  DrawEE(endc_m_r[gId], rEEmin[gId], rEEmax[gId]);
411  pad[gId][1]->cd();
412  DrawEB(barrel_m[gId], mEBmin[gId], mEBmax[gId]);
413  pad[gId + 3][1]->cd();
414  DrawEB(barrel_r[gId], rEBmin[gId], rEBmax[gId]);
415  pad[gId][2]->cd();
416  DrawEE(endc_p_m[gId], mEEmin[gId], mEEmax[gId]);
417  pad[gId + 3][2]->cd();
418  DrawEE(endc_p_r[gId], rEEmin[gId], rEEmax[gId]);
419  }
420 
421  std::string ImageName(m_imageFileName);
422  canvas.SaveAs(ImageName.c_str());
423  return true;
424  }// fill method
425  };
426 
427 } // close namespace
428 
429  // Register the classes as boost python plugin
431  PAYLOAD_INSPECTOR_CLASS(EcalTPGLinearizationConstPlot);
432  PAYLOAD_INSPECTOR_CLASS(EcalTPGLinearizationConstDiff);
433 }
std::shared_ptr< PayloadType > fetchPayload(const cond::Hash &payloadHash)
static const int kSides
const Int_t gainValues[kGains]
void DrawEE(TH2F *endc, float min, float max)
Definition: EcalDrawUtils.h:29
void DrawEB(TH2F *ebmap, float min, float max)
Definition: EcalDrawUtils.h:4
#define PAYLOAD_INSPECTOR_CLASS(CLASS_NAME)
virtual bool fill(const std::vector< std::tuple< cond::Time_t, cond::Hash > > &iovs)=0
#define PAYLOAD_INSPECTOR_MODULE(PAYLOAD_TYPENAME)
static bool validDetId(int crystal_ix, int crystal_iy, int iz)
Definition: EEDetId.h:248
constexpr int gainId(sample_type sample)
get the gainId (2 bits)
Definition: plugin.cc:24
def canvas(sub, attr)
Definition: svgfig.py:482
const Int_t kGains