CMS 3D CMS Logo

EcalLaserAPDPNRatios_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 "TStyle.h"
14 #include "TLine.h"
15 #include "TLatex.h"
16 
17 #include <memory>
18 #include <sstream>
19 
20 
21 
22 namespace {
23  enum {kEBChannels = 61200, kEEChannels = 14648};
24  enum {MIN_IETA = 1, MIN_IPHI = 1, MAX_IETA = 85, MAX_IPHI = 360, EBhistEtaMax = 171}; // barrel lower and upper bounds on eta and phi
25  enum {IX_MIN = 1, IY_MIN = 1, IX_MAX = 100, IY_MAX = 100, EEhistXMax = 220}; // endcaps lower and upper bounds on x and y
26 
27  /*******************************************************
28 
29  2d histogram of ECAL barrel APDPNRatios of 1 IOV
30 
31  *******************************************************/
32 
33 
34  // inherit from one of the predefined plot class: Histogram2D
35  class EcalLaserAPDPNRatiosEBMap : public cond::payloadInspector::Histogram2D<EcalLaserAPDPNRatios> {
36 
37  public:
38  EcalLaserAPDPNRatiosEBMap() : cond::payloadInspector::Histogram2D<EcalLaserAPDPNRatios>( "ECAL Barrel APDPNRatios - map ",
39  "iphi", MAX_IPHI, MIN_IPHI, MAX_IPHI + 1,
40  "ieta", EBhistEtaMax, -MAX_IETA, MAX_IETA + 1) {
41  Base::setSingleIov( true );
42  }
43 
44  // Histogram2D::fill (virtual) needs be overridden - the implementation should use fillWithValue
45  bool fill( const std::vector<std::tuple<cond::Time_t,cond::Hash> >& iovs ) override{
46 
47  for (auto const & iov: iovs) {
48  std::shared_ptr<EcalLaserAPDPNRatios> payload = Base::fetchPayload( std::get<1>(iov) );
49  if( payload.get() ){
50  // set to -1 for ieta 0 (no crystal)
51  for(int iphi = 1; iphi < 361; iphi++) fillWithValue(iphi, 0, -1);
52 
53  for(int cellid = EBDetId::MIN_HASH; cellid < EBDetId::kSizeForDenseIndexing; ++cellid) {
54  uint32_t rawid = EBDetId::unhashIndex(cellid);
55  float p2 = (payload->getLaserMap())[rawid].p2;
56  // fill the Histogram2D here
57  fillWithValue((EBDetId(rawid)).iphi(), (EBDetId(rawid)).ieta(), p2);
58  } // loop over cellid
59  } // if payload.get()
60  } // loop over IOV's (1 in this case)
61 
62  return true;
63  }// fill method
64  };
65 
66  class EcalLaserAPDPNRatiosEEMap : public cond::payloadInspector::Histogram2D<EcalLaserAPDPNRatios> {
67 
68  private:
69  int EEhistSplit = 20;
70 
71  public:
72  EcalLaserAPDPNRatiosEEMap() : cond::payloadInspector::Histogram2D<EcalLaserAPDPNRatios>( "ECAL Endcap APDPNRatios - map ",
73  "ix", EEhistXMax, IX_MIN, EEhistXMax + 1,
74  "iy", IY_MAX, IY_MIN, IY_MAX + 1) {
75  Base::setSingleIov( true );
76  }
77 
78  bool fill( const std::vector<std::tuple<cond::Time_t,cond::Hash> >& iovs ) override{
79 
80  for (auto const & iov: iovs) {
81  std::shared_ptr<EcalLaserAPDPNRatios> payload = Base::fetchPayload( std::get<1>(iov) );
82  if( payload.get() ){
83 
84  // set to -1 everywhwere
85  for(int ix = IX_MIN; ix < EEhistXMax + 1; ix++)
86  for(int iy = IY_MAX; iy < IY_MAX + 1; iy++)
87  fillWithValue(ix, iy, -1);
88 
89  for(int cellid = 0; cellid < EEDetId::kSizeForDenseIndexing; ++cellid) {
90  if(!EEDetId::validHashIndex(cellid)) continue;
91  uint32_t rawid = EEDetId::unhashIndex(cellid);
92  float p2 = (payload->getLaserMap())[rawid].p2;
93  EEDetId myEEId(rawid);
94  if(myEEId.zside() == -1)
95  fillWithValue(myEEId.ix(), myEEId.iy(), p2);
96  else
97  fillWithValue(myEEId.ix() + IX_MAX + EEhistSplit, myEEId.iy(), p2);
98  } // loop over cellid
99  } // payload
100  } // loop over IOV's (1 in this case)
101  return true;
102  }// fill method
103  };
104 
105  /*************************************************
106  2d plot of ECAL IntercalibConstants of 1 IOV
107  *************************************************/
108  class EcalLaserAPDPNRatiosPlot : public cond::payloadInspector::PlotImage<EcalLaserAPDPNRatios> {
109 
110  public:
111  EcalLaserAPDPNRatiosPlot() : cond::payloadInspector::PlotImage<EcalLaserAPDPNRatios>("ECAL Laser APDPNRatios - map ") {
112  setSingleIov( true );
113  }
114 
115  bool fill( const std::vector<std::tuple<cond::Time_t,cond::Hash> >& iovs ) override{
116  TH2F** barrel = new TH2F*[3];
117  TH2F** endc_p = new TH2F*[3];
118  TH2F** endc_m = new TH2F*[3];
119  float pEBmin[3], pEEmin[3];
120 
121  for (int i = 0; i < 3; i++) {
122  barrel[i] = new TH2F(Form("EBp%i", i),Form("EB p%i", i+1),MAX_IPHI, 0, MAX_IPHI, 2 * MAX_IETA, -MAX_IETA, MAX_IETA);
123  endc_p[i] = new TH2F(Form("EE+p%i",i),Form("EE+ p%i",i+1),IX_MAX, IX_MIN, IX_MAX + 1, IY_MAX, IY_MIN, IY_MAX + 1);
124  endc_m[i] = new TH2F(Form("EE-p%i",i),Form("EE- p%i",i+1),IX_MAX, IX_MIN, IX_MAX + 1, IY_MAX, IY_MIN, IY_MAX + 1);
125  pEBmin[i] = 10.;
126  pEEmin[i] = 10.;
127  }
128 
129  auto iov = iovs.front();
130  std::shared_ptr<EcalLaserAPDPNRatios> payload = fetchPayload( std::get<1>(iov) );
131  unsigned long IOV = std::get<0>(iov);
132  int run = 0;
133  if(IOV < 4294967296) run = std::get<0>(iov);
134  else // time type IOV
135  run = IOV >> 32;
136  if( payload.get() ){
137  // looping over the EB channels, via the dense-index, mapped into EBDetId's
138  for(int cellid = 0; cellid < EBDetId::kSizeForDenseIndexing; ++cellid) { // loop on EB cells
139  uint32_t rawid = EBDetId::unhashIndex(cellid);
140  Double_t phi = (Double_t)(EBDetId(rawid)).iphi() - 0.5;
141  Double_t eta = (Double_t)(EBDetId(rawid)).ieta();
142  if(eta > 0.) eta = eta - 0.5; // 0.5 to 84.5
143  else eta = eta + 0.5; // -84.5 to -0.5
144  float p1 = (payload->getLaserMap())[rawid].p1;
145  if(p1 < pEBmin[0]) pEBmin[0] = p1;
146  float p2 = (payload->getLaserMap())[rawid].p2;
147  if(p2 < pEBmin[1]) pEBmin[1] = p2;
148  float p3 = (payload->getLaserMap())[rawid].p3;
149  if(p3 < pEBmin[2]) pEBmin[2] = p3;
150  barrel[0]->Fill(phi, eta, p1);
151  barrel[1]->Fill(phi, eta, p2);
152  barrel[2]->Fill(phi, eta, p3);
153  } // loop over cellid
154 
155  // looping over the EE channels
156  for(int cellid = 0; cellid < EEDetId::kSizeForDenseIndexing; ++cellid) {
157  if(!EEDetId::validHashIndex(cellid)) continue;
158  uint32_t rawid = EEDetId::unhashIndex(cellid);
159  EEDetId myEEId(rawid);
160  float p1 = (payload->getLaserMap())[rawid].p1;
161  if(p1 < pEEmin[0]) pEEmin[0] = p1;
162  float p2 = (payload->getLaserMap())[rawid].p2;
163  if(p2 < pEEmin[1]) pEEmin[1] = p2;
164  float p3 = (payload->getLaserMap())[rawid].p3;
165  if(p3 < pEEmin[2]) pEEmin[2] = p3;
166  if (myEEId.zside() == 1) {
167  endc_p[0]->Fill(myEEId.ix(), myEEId.iy(), p1);
168  endc_p[1]->Fill(myEEId.ix(), myEEId.iy(), p2);
169  endc_p[2]->Fill(myEEId.ix(), myEEId.iy(), p3);
170  }
171  else {
172  endc_m[0]->Fill(myEEId.ix(), myEEId.iy(), p1);
173  endc_m[1]->Fill(myEEId.ix(), myEEId.iy(), p2);
174  endc_m[2]->Fill(myEEId.ix(), myEEId.iy(), p3);
175  }
176  } // validDetId
177  } // if payload.get()
178  else return false;
179 
180  gStyle->SetPalette(1);
181  gStyle->SetOptStat(0);
182  TCanvas canvas("CC map","CC map",2800,2600);
183  TLatex t1;
184  t1.SetNDC();
185  t1.SetTextAlign(26);
186  t1.SetTextSize(0.05);
187  if(IOV < 4294967296)
188  t1.DrawLatex(0.5, 0.96, Form("Ecal Laser APD/PN, IOV %i", run));
189  else { // time type IOV
190  time_t t = run;
191  char buf[256];
192  struct tm lt;
193  localtime_r(&t, &lt);
194  strftime(buf, sizeof(buf), "%F %R:%S", &lt);
195  buf[sizeof(buf)-1] = 0;
196  t1.DrawLatex(0.5, 0.96, Form("Ecal Laser APD/PN, IOV %s", buf));
197  }
198 
199  float xmi[3] = {0.0 , 0.26, 0.74};
200  float xma[3] = {0.26, 0.74, 1.00};
201  TPad*** pad = new TPad**[3];
202  for (int i = 0; i < 3; i++) {
203  pad[i] = new TPad*[3];
204  for (int obj = 0; obj < 3; obj++) {
205  float yma = 0.94 - (0.32 * i);
206  float ymi = yma - 0.28;
207  pad[i][obj] = new TPad(Form("p_%i_%i", obj, i),Form("p_%i_%i", obj, i),
208  xmi[obj], ymi, xma[obj], yma);
209  pad[i][obj]->Draw();
210  }
211  }
212 
213  for (int i = 0; i < 3; i++) {
214  // compute histo limits with some rounding
215  // std::cout << " before " << pEBmin[i];
216  float xmin = pEBmin[i] * 10.;
217  int min = (int)xmin;
218  pEBmin[i] = (float)min /10.;
219  // std::cout << " after " << pEBmin[i] << std::endl << " before " << pEEmin[i];
220  xmin = pEEmin[i] * 10.;
221  min = (int)xmin;
222  pEEmin[i] = (float)min /10.;
223  // std::cout << " after " << pEEmin[i] << std::endl;
224  pad[i][0]->cd();
225  DrawEE(endc_m[i], pEEmin[i], 1.1);
226  pad[i][1]->cd();
227  DrawEB(barrel[i], pEBmin[i], 1.1);
228  pad[i][2]->cd();
229  DrawEE(endc_p[i], pEEmin[i], 1.1);
230  }
231 
232  std::string ImageName(m_imageFileName);
233  canvas.SaveAs(ImageName.c_str());
234  return true;
235  }// fill method
236  };
237 
238  /*****************************************************************
239  2d plot of ECAL IntercalibConstants difference between 2 IOVs
240  ******************************************************************/
241  class EcalLaserAPDPNRatiosDiff : public cond::payloadInspector::PlotImage<EcalLaserAPDPNRatios> {
242 
243  public:
244  EcalLaserAPDPNRatiosDiff() : cond::payloadInspector::PlotImage<EcalLaserAPDPNRatios>("ECAL Laser APDPNRatios difference") {
245  setSingleIov(false);
246  }
247 
248  bool fill( const std::vector<std::tuple<cond::Time_t,cond::Hash> >& iovs ) override{
249  TH2F** barrel = new TH2F*[3];
250  TH2F** endc_p = new TH2F*[3];
251  TH2F** endc_m = new TH2F*[3];
252  float pEBmin[3], pEEmin[3], pEBmax[3], pEEmax[3];
253 
254 
255  for (int i = 0; i < 3; i++) {
256  barrel[i] = new TH2F(Form("EBp%i", i),Form("EB p%i", i+1),MAX_IPHI, 0, MAX_IPHI, 2 * MAX_IETA, -MAX_IETA, MAX_IETA);
257  endc_p[i] = new TH2F(Form("EE+p%i",i),Form("EE+ p%i",i+1),IX_MAX, IX_MIN, IX_MAX + 1, IY_MAX, IY_MIN, IY_MAX + 1);
258  endc_m[i] = new TH2F(Form("EE-p%i",i),Form("EE- p%i",i+1),IX_MAX, IX_MIN, IX_MAX + 1, IY_MAX, IY_MIN, IY_MAX + 1);
259  pEBmin[i] = 10.;
260  pEEmin[i] = 10.;
261  pEBmax[i] = -10.;
262  pEEmax[i] = -10.;
263  }
264  unsigned int run[2] = {0, 0}, irun = 0;
265  unsigned long IOV = 0;
266  float pEB[3][kEBChannels], pEE[3][kEEChannels];
267  for ( auto const & iov: iovs) {
268  std::shared_ptr<EcalLaserAPDPNRatios> payload = fetchPayload( std::get<1>(iov) );
269  IOV = std::get<0>(iov);
270  if(IOV < 4294967296) run[irun] = std::get<0>(iov);
271  else // time type IOV
272  run[irun] = IOV >> 32;
273  if( payload.get() ){
274  // looping over the EB channels, via the dense-index, mapped into EBDetId's
275  for(int cellid = 0; cellid < EBDetId::kSizeForDenseIndexing; ++cellid) { // loop on EB cells
276  uint32_t rawid = EBDetId::unhashIndex(cellid);
277  if(irun == 0) {
278  pEB[0][cellid] = (payload->getLaserMap())[rawid].p1;
279  pEB[1][cellid] = (payload->getLaserMap())[rawid].p2;
280  pEB[2][cellid] = (payload->getLaserMap())[rawid].p3;
281  }
282  else {
283  Double_t phi = (Double_t)(EBDetId(rawid)).iphi() - 0.5;
284  Double_t eta = (Double_t)(EBDetId(rawid)).ieta();
285  if(eta > 0.) eta = eta - 0.5; // 0.5 to 84.5
286  else eta = eta + 0.5; // -84.5 to -0.5
287  double diff = (payload->getLaserMap())[rawid].p1 - pEB[0][cellid];
288  if(diff < pEBmin[0]) pEBmin[0] = diff;
289  if(diff > pEBmax[0]) pEBmax[0] = diff;
290  barrel[0]->Fill(phi, eta, diff);
291  diff = (payload->getLaserMap())[rawid].p2 - pEB[1][cellid];
292  if(diff < pEBmin[1]) pEBmin[1] = diff;
293  if(diff > pEBmax[1]) pEBmax[1] = diff;
294  barrel[1]->Fill(phi, eta, diff);
295  diff = (payload->getLaserMap())[rawid].p3 - pEB[2][cellid];
296  if(diff < pEBmin[2]) pEBmin[2] = diff;
297  if(diff > pEBmax[2]) pEBmax[2] = diff;
298  barrel[2]->Fill(phi, eta, diff);
299  }
300  } // loop over cellid
301 
302  // looping over the EE channels
303  for(int cellid = 0; cellid < EEDetId::kSizeForDenseIndexing; ++cellid) {
304  if(!EEDetId::validHashIndex(cellid)) continue;
305  uint32_t rawid = EEDetId::unhashIndex(cellid);
306  EEDetId myEEId(rawid);
307  if(irun == 0) {
308  pEE[0][cellid] = (payload->getLaserMap())[rawid].p1;
309  pEE[1][cellid] = (payload->getLaserMap())[rawid].p2;
310  pEE[2][cellid] = (payload->getLaserMap())[rawid].p3;
311  }
312  else {
313  double diff1 = (payload->getLaserMap())[rawid].p1 - pEE[0][cellid];
314  if(diff1 < pEEmin[0]) pEEmin[0] = diff1;
315  if(diff1 > pEEmax[0]) pEEmax[0] = diff1;
316  double diff2 = (payload->getLaserMap())[rawid].p2 - pEE[1][cellid];
317  if(diff2 < pEEmin[1]) pEEmin[1] = diff2;
318  if(diff2 > pEEmax[1]) pEEmax[1] = diff2;
319  double diff3 = (payload->getLaserMap())[rawid].p3 - pEE[2][cellid];
320  if(diff3 < pEEmin[2]) pEEmin[2] = diff3;
321  if(diff3 > pEEmax[2]) pEEmax[2] = diff3;
322  if (myEEId.zside() == 1) {
323  endc_p[0]->Fill(myEEId.ix(), myEEId.iy(), diff1);
324  endc_p[1]->Fill(myEEId.ix(), myEEId.iy(), diff2);
325  endc_p[2]->Fill(myEEId.ix(), myEEId.iy(), diff3);
326  }
327  else {
328  endc_m[0]->Fill(myEEId.ix(), myEEId.iy(), diff1);
329  endc_m[1]->Fill(myEEId.ix(), myEEId.iy(), diff2);
330  endc_m[2]->Fill(myEEId.ix(), myEEId.iy(), diff3);
331  }
332  }
333  } // loop over cellid
334  } // if payload.get()
335  else return false;
336  irun++;
337  } // loop over IOVs
338 
339  gStyle->SetPalette(1);
340  gStyle->SetOptStat(0);
341  TCanvas canvas("CC map","CC map",2800,2600);
342  TLatex t1;
343  t1.SetNDC();
344  t1.SetTextAlign(26);
345  t1.SetTextSize(0.05);
346  if(IOV < 4294967296) {
347  t1.SetTextSize(0.05);
348  t1.DrawLatex(0.5, 0.96, Form("Ecal Laser APD/PN, IOV %i - %i", run[1], run[0]));
349  }
350  else { // time type IOV
351  time_t t = run[0];
352  char buf0[256], buf1[256];
353  struct tm lt;
354  localtime_r(&t, &lt);
355  strftime(buf0, sizeof(buf0), "%F %R:%S", &lt);
356  buf0[sizeof(buf0)-1] = 0;
357  t = run[1];
358  localtime_r(&t, &lt);
359  strftime(buf1, sizeof(buf1), "%F %R:%S", &lt);
360  buf1[sizeof(buf1)-1] = 0;
361  t1.SetTextSize(0.015);
362  t1.DrawLatex(0.5, 0.96, Form("Ecal Laser APD/PN, IOV %s - %s", buf1, buf0));
363  }
364 
365  float xmi[3] = {0.0 , 0.24, 0.76};
366  float xma[3] = {0.24, 0.76, 1.00};
367  TPad*** pad = new TPad**[3];
368  for (int i = 0; i < 3; i++) {
369  pad[i] = new TPad*[3];
370  for (int obj = 0; obj < 3; obj++) {
371  float yma = 0.94 - (0.32 * i);
372  float ymi = yma - 0.28;
373  pad[i][obj] = new TPad(Form("p_%i_%i", obj, i),Form("p_%i_%i", obj, i),
374  xmi[obj], ymi, xma[obj], yma);
375  pad[i][obj]->Draw();
376  }
377  }
378 
379  for (int i = 0; i < 3; i++) {
380  // compute histo limits with some rounding
381  // std::cout << " before min " << pEBmin[i] << " max " << pEBmax[i];
382  float xmin = (pEBmin[i] - 0.009) * 100.;
383  int min = (int)xmin;
384  pEBmin[i] = (float)min /100.;
385  float xmax = (pEBmax[i] + 0.009) * 100.;
386  int max = (int)xmax;
387  pEBmax[i] = (float)max /100.;
388  // std::cout << " after min " << pEBmin[i] << " max " << pEBmax[i] << std::endl << " before min " << pEEmin[i] << " max " << pEEmax[i];
389  xmin = (pEEmin[i] + 0.009) * 100.;
390  min = (int)xmin;
391  pEEmin[i] = (float)min /100.;
392  xmax = (pEEmax[i] + 0.009) * 100.;
393  max = (int)xmax;
394  pEEmax[i] = (float)max /100.;
395  // std::cout << " after min " << pEEmin[i] << " max " << pEEmax[i]<< std::endl;
396  pad[i][0]->cd();
397  DrawEE(endc_m[i], pEEmin[i], pEEmax[i]);
398  endc_m[i]->GetZaxis()->SetLabelSize(0.02);
399  pad[i][1]->cd();
400  DrawEB(barrel[i], pEBmin[i], pEBmax[i]);
401  pad[i][2]->cd();
402  DrawEE(endc_p[i], pEEmin[i], pEEmax[i]);
403  endc_p[i]->GetZaxis()->SetLabelSize(0.02);
404  }
405 
406  std::string ImageName(m_imageFileName);
407  canvas.SaveAs(ImageName.c_str());
408  return true;
409  }// fill method
410  };
411 
412 } // close namespace
413 
414 // Register the classes as boost python plugin
416  PAYLOAD_INSPECTOR_CLASS(EcalLaserAPDPNRatiosEBMap);
417  PAYLOAD_INSPECTOR_CLASS(EcalLaserAPDPNRatiosEEMap);
418  PAYLOAD_INSPECTOR_CLASS(EcalLaserAPDPNRatiosPlot);
419  PAYLOAD_INSPECTOR_CLASS(EcalLaserAPDPNRatiosDiff);
420 }
void fillWithValue(float xvalue, float yvalue, float weight=1)
static EEDetId unhashIndex(int hi)
Definition: EEDetId.cc:65
static bool validHashIndex(int i)
Definition: EEDetId.h:239
void DrawEE(TH2F *endc, float min, float max)
Definition: EcalDrawUtils.h:29
void DrawEB(TH2F *ebmap, float min, float max)
Definition: EcalDrawUtils.h:4
bool fill(const std::vector< std::tuple< cond::Time_t, cond::Hash > > &iovs) override
#define PAYLOAD_INSPECTOR_CLASS(CLASS_NAME)
static const int MIN_HASH
Definition: EBDetId.h:149
T min(T a, T b)
Definition: MathUtil.h:58
#define PAYLOAD_INSPECTOR_MODULE(PAYLOAD_TYPENAME)
double p2[4]
Definition: TauolaWrapper.h:90
static EBDetId unhashIndex(int hi)
get a DetId from a compact index for arrays
Definition: EBDetId.h:110
Definition: plugin.cc:24
double p1[4]
Definition: TauolaWrapper.h:89
def canvas(sub, attr)
Definition: svgfig.py:482
std::shared_ptr< PayloadType > fetchPayload(const cond::Hash &payloadHash)
double p3[4]
Definition: TauolaWrapper.h:91