CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
SiStripPedestals_PayloadInspector.cc
Go to the documentation of this file.
1 
12 
13 // the data format of the condition to be inspected
19 
20 // needed for the tracker map
22 
23 // auxilliary functions
29 
30 #include <memory>
31 #include <sstream>
32 #include <iostream>
33 #include <iomanip>
34 #include <boost/tokenizer.hpp>
35 
36 // include ROOT
37 #include "TH2F.h"
38 #include "TLegend.h"
39 #include "TCanvas.h"
40 #include "TLine.h"
41 #include "TStyle.h"
42 #include "TLatex.h"
43 #include "TPave.h"
44 #include "TGaxis.h"
45 #include "TPaveStats.h"
46 
47 namespace {
48  using namespace cond::payloadInspector;
49 
50  class SiStripPedestalContainer : public SiStripCondObjectRepresent::SiStripDataContainer<SiStripPedestals, float> {
51  public:
52  SiStripPedestalContainer(const std::shared_ptr<SiStripPedestals>& payload,
54  const std::string& tagName)
56  payloadType_ = "SiStripPedestals";
58  }
59 
60  void storeAllValues() override {
61  std::vector<uint32_t> detid;
62  payload_->getDetIds(detid);
63 
64  for (const auto& d : detid) {
65  SiStripPedestals::Range range = payload_->getRange(d);
66  for (int it = 0; it < (range.second - range.first) * 8 / 10; ++it) {
67  // to be used to fill the histogram
68  SiStripCondData_.fillByPushBack(d, payload_->getPed(it, range));
69  }
70  }
71  }
72  };
73 
74  class SiStripPedestalCompareByPartition : public PlotImage<SiStripPedestals, MULTI_IOV, 2> {
75  public:
76  SiStripPedestalCompareByPartition()
77  : PlotImage<SiStripPedestals, MULTI_IOV, 2>("SiStrip Compare Pedestals By Partition") {}
78 
79  bool fill() override {
80  // trick to deal with the multi-ioved tag and two tag case at the same time
81  auto theIOVs = PlotBase::getTag<0>().iovs;
82  auto tagname1 = PlotBase::getTag<0>().name;
83  auto tag2iovs = PlotBase::getTag<1>().iovs;
84  auto tagname2 = PlotBase::getTag<1>().name;
85  SiStripPI::MetaData firstiov = theIOVs.front();
86  SiStripPI::MetaData lastiov = tag2iovs.front();
87 
88  std::shared_ptr<SiStripPedestals> last_payload = fetchPayload(std::get<1>(lastiov));
89  std::shared_ptr<SiStripPedestals> first_payload = fetchPayload(std::get<1>(firstiov));
90 
91  SiStripPedestalContainer* l_objContainer = new SiStripPedestalContainer(last_payload, lastiov, tagname1);
92  SiStripPedestalContainer* f_objContainer = new SiStripPedestalContainer(first_payload, firstiov, tagname2);
93 
94  l_objContainer->compare(f_objContainer);
95 
96  //l_objContainer->printAll();
97 
98  TCanvas canvas("Partition summary", "partition summary", 1400, 1000);
99  l_objContainer->fillByPartition(canvas, 300, 0.1, 300.);
100 
101  std::string fileName(m_imageFileName);
102  canvas.SaveAs(fileName.c_str());
103 
104  return true;
105  } // fill
106  };
107 
108  class SiStripPedestalDiffByPartition : public PlotImage<SiStripPedestals, MULTI_IOV, 2> {
109  public:
110  SiStripPedestalDiffByPartition()
111  : PlotImage<SiStripPedestals, MULTI_IOV, 2>("SiStrip Diff Pedestals By Partition") {}
112 
113  bool fill() override {
114  // trick to deal with the multi-ioved tag and two tag case at the same time
115  auto theIOVs = PlotBase::getTag<0>().iovs;
116  auto tagname1 = PlotBase::getTag<0>().name;
117  auto tag2iovs = PlotBase::getTag<1>().iovs;
118  auto tagname2 = PlotBase::getTag<1>().name;
119  SiStripPI::MetaData firstiov = theIOVs.front();
120  SiStripPI::MetaData lastiov = tag2iovs.front();
121 
122  std::shared_ptr<SiStripPedestals> last_payload = fetchPayload(std::get<1>(lastiov));
123  std::shared_ptr<SiStripPedestals> first_payload = fetchPayload(std::get<1>(firstiov));
124 
125  SiStripPedestalContainer* l_objContainer = new SiStripPedestalContainer(last_payload, lastiov, tagname1);
126  SiStripPedestalContainer* f_objContainer = new SiStripPedestalContainer(first_payload, firstiov, tagname2);
127 
128  l_objContainer->subtract(f_objContainer);
129 
130  //l_objContainer->printAll();
131 
132  TCanvas canvas("Partition summary", "partition summary", 1400, 1000);
133  l_objContainer->fillByPartition(canvas, 100, -30., 30.);
134 
135  std::string fileName(m_imageFileName);
136  canvas.SaveAs(fileName.c_str());
137 
138  return true;
139  } // fill
140  };
141 
142  class SiStripPedestalCorrelationByPartition : public PlotImage<SiStripPedestals> {
143  public:
144  SiStripPedestalCorrelationByPartition()
145  : PlotImage<SiStripPedestals>("SiStrip Pedestals Correlation By Partition") {
146  setSingleIov(false);
147  }
148 
149  bool fill(const std::vector<SiStripPI::MetaData>& iovs) override {
150  std::vector<SiStripPI::MetaData> sorted_iovs = iovs;
151 
152  // make absolute sure the IOVs are sortd by since
153  std::sort(begin(sorted_iovs), end(sorted_iovs), [](auto const& t1, auto const& t2) {
154  return std::get<0>(t1) < std::get<0>(t2);
155  });
156 
157  auto firstiov = sorted_iovs.front();
158  auto lastiov = sorted_iovs.back();
159 
160  std::shared_ptr<SiStripPedestals> last_payload = fetchPayload(std::get<1>(lastiov));
161  std::shared_ptr<SiStripPedestals> first_payload = fetchPayload(std::get<1>(firstiov));
162 
163  SiStripPedestalContainer* l_objContainer = new SiStripPedestalContainer(last_payload, lastiov, "");
164  SiStripPedestalContainer* f_objContainer = new SiStripPedestalContainer(first_payload, firstiov, "");
165 
166  l_objContainer->compare(f_objContainer);
167 
168  TCanvas canvas("Partition summary", "partition summary", 1200, 1200);
169  l_objContainer->fillCorrelationByPartition(canvas, 100, 0., 300.);
170 
171  std::string fileName(m_imageFileName);
172  canvas.SaveAs(fileName.c_str());
173 
174  return true;
175  } // fill
176  };
177 
178  /************************************************
179  test class
180  *************************************************/
181 
182  class SiStripPedestalsTest : public Histogram1D<SiStripPedestals, SINGLE_IOV> {
183  public:
184  SiStripPedestalsTest()
185  : Histogram1D<SiStripPedestals, SINGLE_IOV>("SiStrip Pedestals test", "SiStrip Pedestals test", 10, 0.0, 10.0),
187  edm::FileInPath("Geometry/TrackerCommonData/data/trackerParameters.xml").fullPath())} {}
188 
189  bool fill() override {
190  auto tag = PlotBase::getTag<0>();
191  for (auto const& iov : tag.iovs) {
192  std::shared_ptr<SiStripPedestals> payload = Base::fetchPayload(std::get<1>(iov));
193  if (payload.get()) {
194  fillWithValue(1.);
195 
196  std::stringstream ss;
197  ss << "Summary of strips pedestals:" << std::endl;
198 
199  //payload->printDebug(ss);
200  payload->printSummary(ss, &m_trackerTopo);
201 
202  std::vector<uint32_t> detid;
203  payload->getDetIds(detid);
204 
205  // for (const auto & d : detid) {
206  // int nstrip=0;
207  // SiStripPedestals::Range range=payload->getRange(d);
208  // for( int it=0; it < (range.second-range.first)*8/10; ++it ){
209  // auto ped = payload->getPed(it,range);
210  // nstrip++;
211  // ss << "DetId="<< d << " Strip=" << nstrip <<": "<< ped << std::endl;
212  // } // end of loop on strips
213  // } // end of loop on detIds
214 
215  std::cout << ss.str() << std::endl;
216 
217  } // payload
218  } // iovs
219  return true;
220  } // fill
221  private:
222  TrackerTopology m_trackerTopo;
223  };
224 
225  /************************************************
226  SiStrip Pedestals Profile of 1 IOV for one selected DetId
227  *************************************************/
228 
229  class SiStripPedestalPerDetId : public PlotImage<SiStripPedestals, SINGLE_IOV> {
230  public:
231  SiStripPedestalPerDetId() : PlotImage<SiStripPedestals, SINGLE_IOV>("SiStrip Pedestal values per DetId") {
232  PlotBase::addInputParam("DetIds");
233  }
234 
235  bool fill() override {
236  auto tag = PlotBase::getTag<0>();
237  auto iov = tag.iovs.front();
238  auto tagname = tag.name;
239  std::shared_ptr<SiStripPedestals> payload = fetchPayload(std::get<1>(iov));
240 
241  std::vector<uint32_t> the_detids = {};
242 
243  auto paramValues = PlotBase::inputParamValues();
244  auto ip = paramValues.find("DetIds");
245  if (ip != paramValues.end()) {
246  auto input = ip->second;
247  typedef boost::tokenizer<boost::char_separator<char>> tokenizer;
248  boost::char_separator<char> sep{","};
249  tokenizer tok{input, sep};
250  for (const auto& t : tok) {
251  the_detids.push_back(atoi(t.c_str()));
252  }
253  } else {
254  edm::LogWarning("SiStripNoisePerDetId")
255  << "\n WARNING!!!! \n The needed parameter DetIds has not been passed. Will use all Strip DetIds! \n\n";
256  the_detids.push_back(0xFFFFFFFF);
257  }
258 
259  size_t ndets = the_detids.size();
260  std::vector<std::shared_ptr<TH1F>> hpedestal;
261  std::vector<std::shared_ptr<TLegend>> legends;
262  std::vector<unsigned int> v_nAPVs;
263  std::vector<std::vector<std::shared_ptr<TLine>>> lines;
264  hpedestal.reserve(ndets);
265  legends.reserve(ndets);
266 
267  auto sides = getClosestFactors(the_detids.size());
268  edm::LogPrint("SiStripPedestalPerDetId") << "Aspect ratio: " << sides.first << ":" << sides.second << std::endl;
269 
270  if (payload.get()) {
271  //=========================
272  TCanvas canvas("ByDetId", "ByDetId", sides.second * 800, sides.first * 600);
273  canvas.Divide(sides.second, sides.first);
274  const auto detInfo =
276  for (const auto& the_detid : the_detids) {
277  edm::LogPrint("SiStripNoisePerDetId") << "DetId:" << the_detid << std::endl;
278 
279  unsigned int nAPVs = detInfo.getNumberOfApvsAndStripLength(the_detid).first;
280  if (nAPVs == 0)
281  nAPVs = 6;
282  v_nAPVs.push_back(nAPVs);
283 
284  auto histo = std::make_shared<TH1F>(
285  Form("Pedestal profile %s", std::to_string(the_detid).c_str()),
286  Form("SiStrip Pedestal profile for DetId: %s;Strip number;SiStrip Pedestal [ADC counts]",
287  std::to_string(the_detid).c_str()),
288  sistrip::STRIPS_PER_APV * nAPVs,
289  -0.5,
290  (sistrip::STRIPS_PER_APV * nAPVs) - 0.5);
291 
292  histo->SetStats(false);
293  histo->SetTitle("");
294 
295  if (the_detid != 0xFFFFFFFF) {
296  fillHisto(payload, histo, the_detid);
297  } else {
298  auto allDetIds = detInfo.getAllDetIds();
299  for (const auto& id : allDetIds) {
300  fillHisto(payload, histo, id);
301  }
302  }
303 
305  histo->GetYaxis()->SetTitleOffset(1.0);
306  hpedestal.push_back(histo);
307  } // loop on the detids
308 
309  for (size_t index = 0; index < ndets; index++) {
310  canvas.cd(index + 1);
311  canvas.cd(index + 1)->SetBottomMargin(0.11);
312  canvas.cd(index + 1)->SetTopMargin(0.06);
313  canvas.cd(index + 1)->SetLeftMargin(0.10);
314  canvas.cd(index + 1)->SetRightMargin(0.02);
315  hpedestal.at(index)->Draw();
316  hpedestal.at(index)->GetYaxis()->SetRangeUser(0, hpedestal.at(index)->GetMaximum() * 1.2);
317  canvas.cd(index)->Update();
318 
319  std::vector<int> boundaries;
320  for (size_t b = 0; b < v_nAPVs.at(index); b++) {
321  boundaries.push_back(b * sistrip::STRIPS_PER_APV);
322  }
323 
324  std::vector<std::shared_ptr<TLine>> linesVec;
325  for (const auto& bound : boundaries) {
326  auto line = std::make_shared<TLine>(hpedestal.at(index)->GetBinLowEdge(bound),
327  canvas.cd(index + 1)->GetUymin(),
328  hpedestal.at(index)->GetBinLowEdge(bound),
329  canvas.cd(index + 1)->GetUymax());
330  line->SetLineWidth(1);
331  line->SetLineStyle(9);
332  line->SetLineColor(2);
333  linesVec.push_back(line);
334  }
335  lines.push_back(linesVec);
336 
337  for (const auto& line : lines.at(index)) {
338  line->Draw("same");
339  }
340 
341  canvas.cd(index + 1);
342 
343  auto ltx = TLatex();
344  ltx.SetTextFont(62);
345  ltx.SetTextSize(0.05);
346  ltx.SetTextAlign(11);
347  ltx.DrawLatexNDC(gPad->GetLeftMargin(),
348  1 - gPad->GetTopMargin() + 0.01,
349  Form("SiStrip Pedestals profile for DetId %s", std::to_string(the_detids[index]).c_str()));
350 
351  legends.push_back(std::make_shared<TLegend>(0.45, 0.83, 0.95, 0.93));
352  legends.at(index)->SetHeader(tagname.c_str(), "C"); // option "C" allows to center the header
353  legends.at(index)->AddEntry(
354  hpedestal.at(index).get(), ("IOV: " + std::to_string(std::get<0>(iov))).c_str(), "PL");
355  legends.at(index)->SetTextSize(0.045);
356  legends.at(index)->Draw("same");
357  }
358 
359  std::string fileName(m_imageFileName);
360  canvas.SaveAs(fileName.c_str());
361  } // payload
362  return true;
363  } // fill
364 
365  private:
366  int nextPerfectSquare(int N) { return std::floor(sqrt(N)) + 1; }
367 
368  std::pair<int, int> getClosestFactors(int input) {
369  if ((input % 2 != 0) && input > 1) {
370  input += 1;
371  }
372 
373  int testNum = (int)sqrt(input);
374  while (input % testNum != 0) {
375  testNum--;
376  }
377  return std::make_pair(testNum, input / testNum);
378  }
379 
380  void fillHisto(const std::shared_ptr<SiStripPedestals> payload, std::shared_ptr<TH1F>& histo, uint32_t the_detid) {
381  int nstrip = 0;
382  SiStripPedestals::Range range = payload->getRange(the_detid);
383  for (int it = 0; it < (range.second - range.first) * 8 / 10; ++it) {
384  auto pedestal = payload->getPed(it, range);
385  nstrip++;
386  histo->AddBinContent(nstrip, pedestal);
387  } // end of loop on strips
388  }
389  };
390 
391  /************************************************
392  1d histogram of SiStripPedestals of 1 IOV
393  *************************************************/
394 
395  // inherit from one of the predefined plot class: Histogram1D
396  class SiStripPedestalsValue : public Histogram1D<SiStripPedestals, SINGLE_IOV> {
397  public:
398  SiStripPedestalsValue()
400  "SiStrip Pedestals values", "SiStrip Pedestals values", 300, 0.0, 300.0) {}
401 
402  bool fill() override {
403  auto tag = PlotBase::getTag<0>();
404  for (auto const& iov : tag.iovs) {
405  std::shared_ptr<SiStripPedestals> payload = Base::fetchPayload(std::get<1>(iov));
406  if (payload.get()) {
407  std::vector<uint32_t> detid;
408  payload->getDetIds(detid);
409 
410  for (const auto& d : detid) {
411  SiStripPedestals::Range range = payload->getRange(d);
412  for (int it = 0; it < (range.second - range.first) * 8 / 10; ++it) {
413  auto ped = payload->getPed(it, range);
414  //to be used to fill the histogram
415  fillWithValue(ped);
416  } // loop over APVs
417  } // loop over detIds
418  } // payload
419  } // iovs
420  return true;
421  } // fill
422  };
423 
424  /************************************************
425  1d histogram of SiStripPedestals of 1 IOV per Detid
426  *************************************************/
427 
428  // inherit from one of the predefined plot class: Histogram1D
429  class SiStripPedestalsValuePerDetId : public Histogram1D<SiStripPedestals, SINGLE_IOV> {
430  public:
431  SiStripPedestalsValuePerDetId()
433  "SiStrip Pedestal values per DetId", "SiStrip Pedestal values per DetId", 100, 0.0, 10.0) {
434  PlotBase::addInputParam("DetId");
435  }
436 
437  bool fill() override {
438  auto tag = PlotBase::getTag<0>();
439  for (auto const& iov : tag.iovs) {
440  std::shared_ptr<SiStripPedestals> payload = Base::fetchPayload(std::get<1>(iov));
441  unsigned int the_detid(0xFFFFFFFF);
442  auto paramValues = PlotBase::inputParamValues();
443  auto ip = paramValues.find("DetId");
444  if (ip != paramValues.end()) {
445  the_detid = std::stoul(ip->second);
446  }
447 
448  if (payload.get()) {
449  SiStripPedestals::Range range = payload->getRange(the_detid);
450  for (int it = 0; it < (range.second - range.first) * 8 / 9; ++it) {
451  auto noise = payload->getPed(it, range);
452  //to be used to fill the histogram
453  fillWithValue(noise);
454  } // loop over APVs
455  } // payload
456  } // iovs
457  return true;
458  } // fill
459  };
460 
461  /************************************************
462  templated 1d histogram of SiStripPedestals of 1 IOV
463  *************************************************/
464 
465  // inherit from one of the predefined plot class: PlotImage
466  template <SiStripPI::OpMode op_mode_>
467  class SiStripPedestalDistribution : public PlotImage<SiStripPedestals, SINGLE_IOV> {
468  public:
469  SiStripPedestalDistribution() : PlotImage<SiStripPedestals, SINGLE_IOV>("SiStrip Pedestal values") {}
470 
471  bool fill() override {
472  auto tag = PlotBase::getTag<0>();
473  auto iov = tag.iovs.front();
474  TGaxis::SetMaxDigits(3);
475  gStyle->SetOptStat("emr");
476 
477  std::shared_ptr<SiStripPedestals> payload = fetchPayload(std::get<1>(iov));
478 
479  auto mon1D = std::unique_ptr<SiStripPI::Monitor1D>(new SiStripPI::Monitor1D(
480  op_mode_,
481  "Pedestal",
482  Form("#LT Pedestal #GT per %s for IOV [%s];#LTStrip Pedestal per %s#GT [ADC counts];n. %ss",
483  opType(op_mode_).c_str(),
484  std::to_string(std::get<0>(iov)).c_str(),
485  opType(op_mode_).c_str(),
486  opType(op_mode_).c_str()),
487  300,
488  0.,
489  300.0));
490 
491  unsigned int prev_det = 0, prev_apv = 0;
492  SiStripPI::Entry epedestal;
493 
494  std::vector<uint32_t> detids;
495  payload->getDetIds(detids);
496 
497  // loop on payload
498  for (const auto& d : detids) {
499  SiStripPedestals::Range range = payload->getRange(d);
500 
501  unsigned int istrip = 0;
502 
503  for (int it = 0; it < (range.second - range.first) * 8 / 10; ++it) {
504  auto pedestal = payload->getPed(it, range);
505  bool flush = false;
506  switch (op_mode_) {
507  case (SiStripPI::APV_BASED):
508  flush = (prev_det != 0 && prev_apv != istrip / sistrip::STRIPS_PER_APV);
509  break;
511  flush = (prev_det != 0 && prev_det != d);
512  break;
513  case (SiStripPI::STRIP_BASED):
514  flush = (istrip != 0);
515  break;
516  }
517 
518  if (flush) {
519  mon1D->Fill(prev_apv, prev_det, epedestal.mean());
520  epedestal.reset();
521  }
522 
523  epedestal.add(std::min<float>(pedestal, 300.));
524  prev_apv = istrip / sistrip::STRIPS_PER_APV;
525  istrip++;
526  }
527  prev_det = d;
528  }
529 
530  //=========================
531  TCanvas canvas("Partion summary", "partition summary", 1200, 1000);
532  canvas.cd();
533  canvas.SetBottomMargin(0.11);
534  canvas.SetTopMargin(0.07);
535  canvas.SetLeftMargin(0.13);
536  canvas.SetRightMargin(0.05);
537  canvas.Modified();
538 
539  auto hist = mon1D->getHist();
541  hist.SetStats(kTRUE);
542  hist.SetFillColorAlpha(kRed, 0.35);
543  hist.Draw();
544 
545  canvas.Update();
546 
547  TPaveStats* st = (TPaveStats*)hist.GetListOfFunctions()->FindObject("stats");
548  st->SetLineColor(kRed);
549  st->SetTextColor(kRed);
550  st->SetX1NDC(.75);
551  st->SetX2NDC(.95);
552  st->SetY1NDC(.83);
553  st->SetY2NDC(.93);
554 
555  TLegend legend = TLegend(0.13, 0.83, 0.43, 0.93);
556  legend.SetHeader(Form("SiStrip Pedestal values per %s", opType(op_mode_).c_str()),
557  "C"); // option "C" allows to center the header
558  legend.AddEntry(&hist, ("IOV: " + std::to_string(std::get<0>(iov))).c_str(), "F");
559  legend.SetTextSize(0.025);
560  legend.Draw("same");
561 
562  std::string fileName(m_imageFileName);
563  canvas.SaveAs(fileName.c_str());
564 
565  return true;
566  }
567 
569  std::string types[3] = {"Strip", "APV", "Module"};
570  return types[mode];
571  }
572  };
573 
574  typedef SiStripPedestalDistribution<SiStripPI::STRIP_BASED> SiStripPedestalValuePerStrip;
575  typedef SiStripPedestalDistribution<SiStripPI::APV_BASED> SiStripPedestalValuePerAPV;
576  typedef SiStripPedestalDistribution<SiStripPI::MODULE_BASED> SiStripPedestalValuePerModule;
577 
578  /************************************************
579  template 1d histogram comparison of SiStripPedestals of 1 IOV
580  *************************************************/
581 
582  // inherit from one of the predefined plot class: PlotImage
583 
584  template <SiStripPI::OpMode op_mode_, int ntags, IOVMultiplicity nIOVs>
585  class SiStripPedestalDistributionComparisonBase : public PlotImage<SiStripPedestals, nIOVs, ntags> {
586  public:
587  SiStripPedestalDistributionComparisonBase()
588  : PlotImage<SiStripPedestals, nIOVs, ntags>("SiStrip Pedestal values comparison") {}
589 
590  bool fill() override {
591  TGaxis::SetExponentOffset(-0.1, 0.01, "y"); // X and Y offset for Y axis
592 
593  // trick to deal with the multi-ioved tag and two tag case at the same time
594  auto theIOVs = PlotBase::getTag<0>().iovs;
595  auto tagname1 = PlotBase::getTag<0>().name;
596  std::string tagname2 = "";
597  auto firstiov = theIOVs.front();
598  SiStripPI::MetaData lastiov;
599 
600  // we don't support (yet) comparison with more than 2 tags
601  assert(this->m_plotAnnotations.ntags < 3);
602 
603  if (this->m_plotAnnotations.ntags == 2) {
604  auto tag2iovs = PlotBase::getTag<1>().iovs;
605  tagname2 = PlotBase::getTag<1>().name;
606  lastiov = tag2iovs.front();
607  } else {
608  lastiov = theIOVs.back();
609  }
610 
611  std::shared_ptr<SiStripPedestals> f_payload = this->fetchPayload(std::get<1>(firstiov));
612  std::shared_ptr<SiStripPedestals> l_payload = this->fetchPayload(std::get<1>(lastiov));
613 
614  auto f_mon = std::unique_ptr<SiStripPI::Monitor1D>(new SiStripPI::Monitor1D(
615  op_mode_,
616  "f_Pedestal",
617  Form(";#LTStrip Pedestal per %s#GT [ADC counts];n. %ss", opType(op_mode_).c_str(), opType(op_mode_).c_str()),
618  300,
619  0.,
620  300.));
621 
622  auto l_mon = std::unique_ptr<SiStripPI::Monitor1D>(new SiStripPI::Monitor1D(
623  op_mode_,
624  "l_Pedestal",
625  Form(";#LTStrip Pedestal per %s#GT [ADC counts];n. %ss", opType(op_mode_).c_str(), opType(op_mode_).c_str()),
626  300,
627  0.,
628  300.));
629 
630  unsigned int prev_det = 0, prev_apv = 0;
631  SiStripPI::Entry epedestal;
632 
633  std::vector<uint32_t> f_detid;
634  f_payload->getDetIds(f_detid);
635 
636  // loop on first payload
637  for (const auto& d : f_detid) {
638  SiStripPedestals::Range range = f_payload->getRange(d);
639 
640  unsigned int istrip = 0;
641  for (int it = 0; it < (range.second - range.first) * 8 / 10; ++it) {
642  float pedestal = f_payload->getPed(it, range);
643  //to be used to fill the histogram
644 
645  bool flush = false;
646  switch (op_mode_) {
647  case (SiStripPI::APV_BASED):
648  flush = (prev_det != 0 && prev_apv != istrip / sistrip::STRIPS_PER_APV);
649  break;
651  flush = (prev_det != 0 && prev_det != d);
652  break;
653  case (SiStripPI::STRIP_BASED):
654  flush = (istrip != 0);
655  break;
656  }
657 
658  if (flush) {
659  f_mon->Fill(prev_apv, prev_det, epedestal.mean());
660  epedestal.reset();
661  }
662  epedestal.add(std::min<float>(pedestal, 300.));
663  prev_apv = istrip / sistrip::STRIPS_PER_APV;
664  istrip++;
665  }
666  prev_det = d;
667  }
668 
669  prev_det = 0, prev_apv = 0;
670  epedestal.reset();
671 
672  std::vector<uint32_t> l_detid;
673  l_payload->getDetIds(l_detid);
674 
675  // loop on first payload
676  for (const auto& d : l_detid) {
677  SiStripPedestals::Range range = l_payload->getRange(d);
678 
679  unsigned int istrip = 0;
680  for (int it = 0; it < (range.second - range.first) * 8 / 10; ++it) {
681  float pedestal = l_payload->getPed(it, range);
682 
683  bool flush = false;
684  switch (op_mode_) {
685  case (SiStripPI::APV_BASED):
686  flush = (prev_det != 0 && prev_apv != istrip / sistrip::STRIPS_PER_APV);
687  break;
689  flush = (prev_det != 0 && prev_det != d);
690  break;
691  case (SiStripPI::STRIP_BASED):
692  flush = (istrip != 0);
693  break;
694  }
695 
696  if (flush) {
697  l_mon->Fill(prev_apv, prev_det, epedestal.mean());
698  epedestal.reset();
699  }
700 
701  epedestal.add(std::min<float>(pedestal, 300.));
702  prev_apv = istrip / sistrip::STRIPS_PER_APV;
703  istrip++;
704  }
705  prev_det = d;
706  }
707 
708  auto h_first = f_mon->getHist();
709  h_first.SetStats(kFALSE);
710  auto h_last = l_mon->getHist();
711  h_last.SetStats(kFALSE);
712 
715 
716  h_first.GetYaxis()->CenterTitle(true);
717  h_last.GetYaxis()->CenterTitle(true);
718 
719  h_first.GetXaxis()->CenterTitle(true);
720  h_last.GetXaxis()->CenterTitle(true);
721 
722  h_first.SetLineWidth(2);
723  h_last.SetLineWidth(2);
724 
725  h_first.SetLineColor(kBlack);
726  h_last.SetLineColor(kBlue);
727 
728  //=========================
729  TCanvas canvas("Partion summary", "partition summary", 1200, 1000);
730  canvas.cd();
731  canvas.SetTopMargin(0.06);
732  canvas.SetBottomMargin(0.10);
733  canvas.SetLeftMargin(0.13);
734  canvas.SetRightMargin(0.05);
735  canvas.Modified();
736 
737  float theMax = (h_first.GetMaximum() > h_last.GetMaximum()) ? h_first.GetMaximum() : h_last.GetMaximum();
738 
739  h_first.SetMaximum(theMax * 1.20);
740  h_last.SetMaximum(theMax * 1.20);
741 
742  h_first.Draw();
743  h_last.SetFillColorAlpha(kBlue, 0.15);
744  h_last.Draw("same");
745 
746  TLegend legend = TLegend(0.13, 0.83, 0.95, 0.94);
747  if (this->m_plotAnnotations.ntags == 2) {
748  legend.SetHeader("#bf{Two Tags Comparison}", "C"); // option "C" allows to center the header
749  legend.AddEntry(&h_first, (tagname1 + " : " + std::to_string(std::get<0>(firstiov))).c_str(), "F");
750  legend.AddEntry(&h_last, (tagname2 + " : " + std::to_string(std::get<0>(lastiov))).c_str(), "F");
751  } else {
752  legend.SetHeader(("tag: #bf{" + tagname1 + "}").c_str(), "C"); // option "C" allows to center the header
753  legend.AddEntry(&h_first, ("IOV since: " + std::to_string(std::get<0>(firstiov))).c_str(), "F");
754  legend.AddEntry(&h_last, ("IOV since: " + std::to_string(std::get<0>(lastiov))).c_str(), "F");
755  }
756  legend.SetTextSize(0.025);
757  legend.Draw("same");
758 
759  auto ltx = TLatex();
760  ltx.SetTextFont(62);
761  ltx.SetTextSize(0.05);
762  ltx.SetTextAlign(11);
763  ltx.DrawLatexNDC(gPad->GetLeftMargin(),
764  1 - gPad->GetTopMargin() + 0.01,
765  Form("#LTSiStrip Pedestals#GT Comparison per %s", opType(op_mode_).c_str()));
766 
767  std::string fileName(this->m_imageFileName);
768  canvas.SaveAs(fileName.c_str());
769 
770  return true;
771  }
772 
774  std::string types[3] = {"Strip", "APV", "Module"};
775  return types[mode];
776  }
777  };
778 
779  template <SiStripPI::OpMode op_mode_>
780  using SiStripPedestalDistributionComparisonSingleTag =
781  SiStripPedestalDistributionComparisonBase<op_mode_, 1, MULTI_IOV>;
782 
783  template <SiStripPI::OpMode op_mode_>
784  using SiStripPedestalDistributionComparisonTwoTags =
785  SiStripPedestalDistributionComparisonBase<op_mode_, 2, SINGLE_IOV>;
786 
787  typedef SiStripPedestalDistributionComparisonSingleTag<SiStripPI::STRIP_BASED>
788  SiStripPedestalValueComparisonPerStripSingleTag;
789  typedef SiStripPedestalDistributionComparisonSingleTag<SiStripPI::APV_BASED>
790  SiStripPedestalValueComparisonPerAPVSingleTag;
791  typedef SiStripPedestalDistributionComparisonSingleTag<SiStripPI::MODULE_BASED>
792  SiStripPedestalValueComparisonPerModuleSingleTag;
793 
794  typedef SiStripPedestalDistributionComparisonTwoTags<SiStripPI::STRIP_BASED>
795  SiStripPedestalValueComparisonPerStripTwoTags;
796  typedef SiStripPedestalDistributionComparisonTwoTags<SiStripPI::APV_BASED> SiStripPedestalValueComparisonPerAPVTwoTags;
797  typedef SiStripPedestalDistributionComparisonTwoTags<SiStripPI::MODULE_BASED>
798  SiStripPedestalValueComparisonPerModuleTwoTags;
799 
800  /************************************************
801  1d histogram of fraction of Zero SiStripPedestals of 1 IOV
802  *************************************************/
803 
804  // inherit from one of the predefined plot class PlotImage
805  class SiStripZeroPedestalsFraction_TrackerMap : public PlotImage<SiStripPedestals, SINGLE_IOV> {
806  public:
807  SiStripZeroPedestalsFraction_TrackerMap()
808  : PlotImage<SiStripPedestals, SINGLE_IOV>("Tracker Map of Zero SiStripPedestals fraction per module") {}
809 
810  bool fill() override {
811  auto tag = PlotBase::getTag<0>();
812  auto iov = tag.iovs.front();
813  std::shared_ptr<SiStripPedestals> payload = fetchPayload(std::get<1>(iov));
814 
815  const auto detInfo =
817 
818  std::string titleMap =
819  "Tracker Map of Zero SiStrip Pedestals fraction per module (payload : " + std::get<1>(iov) + ")";
820 
821  std::unique_ptr<TrackerMap> tmap = std::make_unique<TrackerMap>("SiStripPedestals");
822  tmap->setTitle(titleMap);
823  tmap->setPalette(1);
824 
825  std::vector<uint32_t> detid;
826  payload->getDetIds(detid);
827 
828  std::map<uint32_t, int> zeropeds_per_detid;
829 
830  for (const auto& d : detid) {
831  int nstrips = 0;
832  SiStripPedestals::Range range = payload->getRange(d);
833  for (int it = 0; it < (range.second - range.first) * 8 / 10; ++it) {
834  nstrips++;
835  auto ped = payload->getPed(it, range);
836  if (ped == 0.) {
837  zeropeds_per_detid[d] += 1;
838  }
839  } // end of loop on strips
840  float fraction =
841  zeropeds_per_detid[d] / (sistrip::STRIPS_PER_APV * detInfo.getNumberOfApvsAndStripLength(d).first);
842  if (fraction > 0.) {
843  tmap->fill(d, fraction);
844  std::cout << "detid: " << d << " (n. APVs=" << detInfo.getNumberOfApvsAndStripLength(d).first << ") has "
845  << std::setw(4) << zeropeds_per_detid[d]
846  << " zero-pedestals strips (i.e. a fraction:" << std::setprecision(5) << fraction << ")"
847  << std::endl;
848  }
849  } // end of loop on detids
850 
851  std::string fileName(m_imageFileName);
852  tmap->save(true, 0., 0., fileName);
853 
854  return true;
855  }
856  };
857 
858  /************************************************
859  Tracker Map of SiStrip Pedestals
860  *************************************************/
861 
862  template <SiStripPI::estimator est>
863  class SiStripPedestalsTrackerMap : public PlotImage<SiStripPedestals, SINGLE_IOV> {
864  public:
865  SiStripPedestalsTrackerMap()
866  : PlotImage<SiStripPedestals, SINGLE_IOV>("Tracker Map of SiStripPedestals " + estimatorType(est) +
867  " per module") {}
868 
869  bool fill() override {
870  auto tag = PlotBase::getTag<0>();
871  auto iov = tag.iovs.front();
872  std::shared_ptr<SiStripPedestals> payload = fetchPayload(std::get<1>(iov));
873 
874  std::string titleMap =
875  "Tracker Map of SiStrip Pedestals " + estimatorType(est) + " per module (payload : " + std::get<1>(iov) + ")";
876 
877  std::unique_ptr<TrackerMap> tmap = std::make_unique<TrackerMap>("SiStripPedestals");
878  tmap->setTitle(titleMap);
879  tmap->setPalette(1);
880 
881  std::vector<uint32_t> detid;
882  payload->getDetIds(detid);
883 
884  std::map<unsigned int, float> info_per_detid;
885 
886  for (const auto& d : detid) {
887  int nstrips = 0;
888  double mean(0.), rms(0.), min(10000.), max(0.);
889  SiStripPedestals::Range range = payload->getRange(d);
890  for (int it = 0; it < (range.second - range.first) * 8 / 10; ++it) {
891  nstrips++;
892  auto ped = payload->getPed(it, range);
893  mean += ped;
894  rms += (ped * ped);
895  if (ped < min)
896  min = ped;
897  if (ped > max)
898  max = ped;
899  } // end of loop on strips
900 
901  mean /= nstrips;
902  if ((rms / nstrips - mean * mean) > 0.) {
903  rms = sqrt(rms / nstrips - mean * mean);
904  } else {
905  rms = 0.;
906  }
907 
908  switch (est) {
909  case SiStripPI::min:
910  info_per_detid[d] = min;
911  break;
912  case SiStripPI::max:
913  info_per_detid[d] = max;
914  break;
915  case SiStripPI::mean:
916  info_per_detid[d] = mean;
917  break;
918  case SiStripPI::rms:
919  info_per_detid[d] = rms;
920  break;
921  default:
922  edm::LogWarning("LogicError") << "Unknown estimator: " << est;
923  break;
924  }
925  } // end of loop on detids
926 
927  for (const auto& d : detid) {
928  tmap->fill(d, info_per_detid[d]);
929  }
930 
931  std::string fileName(m_imageFileName);
932  tmap->save(true, 0., 0., fileName);
933 
934  return true;
935  }
936  };
937 
938  typedef SiStripPedestalsTrackerMap<SiStripPI::min> SiStripPedestalsMin_TrackerMap;
939  typedef SiStripPedestalsTrackerMap<SiStripPI::max> SiStripPedestalsMax_TrackerMap;
940  typedef SiStripPedestalsTrackerMap<SiStripPI::mean> SiStripPedestalsMean_TrackerMap;
941  typedef SiStripPedestalsTrackerMap<SiStripPI::rms> SiStripPedestalsRMS_TrackerMap;
942 
943  /************************************************
944  Tracker Map of SiStrip Pedestals Summaries
945  *************************************************/
946 
947  template <SiStripPI::estimator est>
948  class SiStripPedestalsByRegion : public PlotImage<SiStripPedestals, SINGLE_IOV> {
949  public:
950  SiStripPedestalsByRegion()
951  : PlotImage<SiStripPedestals, SINGLE_IOV>("SiStrip Pedestals " + estimatorType(est) + " by Region"),
953  edm::FileInPath("Geometry/TrackerCommonData/data/trackerParameters.xml").fullPath())} {}
954 
955  bool fill() override {
956  auto tag = PlotBase::getTag<0>();
957  auto iov = tag.iovs.front();
958  std::shared_ptr<SiStripPedestals> payload = fetchPayload(std::get<1>(iov));
959 
960  SiStripDetSummary summaryPedestals{&m_trackerTopo};
961  std::vector<uint32_t> detid;
962  payload->getDetIds(detid);
963 
964  for (const auto& d : detid) {
965  int nstrips = 0;
966  double mean(0.), rms(0.), min(10000.), max(0.);
967  SiStripPedestals::Range range = payload->getRange(d);
968  for (int it = 0; it < (range.second - range.first) * 8 / 10; ++it) {
969  nstrips++;
970  auto ped = payload->getPed(it, range);
971  mean += ped;
972  rms += (ped * ped);
973  if (ped < min)
974  min = ped;
975  if (ped > max)
976  max = ped;
977  } // end of loop on strips
978 
979  mean /= nstrips;
980  if ((rms / nstrips - mean * mean) > 0.) {
981  rms = sqrt(rms / nstrips - mean * mean);
982  } else {
983  rms = 0.;
984  }
985 
986  switch (est) {
987  case SiStripPI::min:
988  summaryPedestals.add(d, min);
989  break;
990  case SiStripPI::max:
991  summaryPedestals.add(d, max);
992  break;
993  case SiStripPI::mean:
994  summaryPedestals.add(d, mean);
995  break;
996  case SiStripPI::rms:
997  summaryPedestals.add(d, rms);
998  break;
999  default:
1000  edm::LogWarning("LogicError") << "Unknown estimator: " << est;
1001  break;
1002  }
1003  } // loop on the detIds
1004 
1005  std::map<unsigned int, SiStripDetSummary::Values> map = summaryPedestals.getCounts();
1006  //=========================
1007 
1008  TCanvas canvas("Partion summary", "partition summary", 1200, 1000);
1009  canvas.cd();
1010  auto h1 = std::unique_ptr<TH1F>(new TH1F(
1011  "byRegion",
1012  Form("Average by partition of %s SiStrip Pedestals per module;;average SiStrip Pedestals %s [ADC counts]",
1013  estimatorType(est).c_str(),
1014  estimatorType(est).c_str()),
1015  map.size(),
1016  0.,
1017  map.size()));
1018  h1->SetStats(false);
1019  canvas.SetBottomMargin(0.18);
1020  canvas.SetLeftMargin(0.17);
1021  canvas.SetRightMargin(0.05);
1022  canvas.Modified();
1023 
1024  std::vector<int> boundaries;
1025  unsigned int iBin = 0;
1026 
1028  std::string currentDetector;
1029 
1030  for (const auto& element : map) {
1031  iBin++;
1032  int count = element.second.count;
1033  double mean = (element.second.mean) / count;
1034  double rms = (element.second.rms) / count - mean * mean;
1035 
1036  if (rms <= 0)
1037  rms = 0;
1038  else
1039  rms = sqrt(rms);
1040 
1041  if (currentDetector.empty())
1042  currentDetector = "TIB";
1043 
1044  switch ((element.first) / 1000) {
1045  case 1:
1046  detector = "TIB";
1047  break;
1048  case 2:
1049  detector = "TOB";
1050  break;
1051  case 3:
1052  detector = "TEC";
1053  break;
1054  case 4:
1055  detector = "TID";
1056  break;
1057  }
1058 
1059  h1->SetBinContent(iBin, mean);
1060  h1->GetXaxis()->SetBinLabel(iBin, SiStripPI::regionType(element.first).second);
1061  h1->GetXaxis()->LabelsOption("v");
1062 
1063  if (detector != currentDetector) {
1064  boundaries.push_back(iBin);
1065  currentDetector = detector;
1066  }
1067  }
1068 
1069  h1->SetMarkerStyle(20);
1070  h1->SetMarkerSize(1);
1071  h1->SetMaximum(h1->GetMaximum() * 1.1);
1072  h1->Draw("HIST");
1073  h1->Draw("Psame");
1074 
1075  canvas.Update();
1076 
1077  TLine l[boundaries.size()];
1078  unsigned int i = 0;
1079  for (const auto& line : boundaries) {
1080  l[i] = TLine(h1->GetBinLowEdge(line), canvas.GetUymin(), h1->GetBinLowEdge(line), canvas.GetUymax());
1081  l[i].SetLineWidth(1);
1082  l[i].SetLineStyle(9);
1083  l[i].SetLineColor(2);
1084  l[i].Draw("same");
1085  i++;
1086  }
1087 
1088  TLegend legend = TLegend(0.52, 0.82, 0.95, 0.9);
1089  legend.SetHeader((std::get<1>(iov)).c_str(), "C"); // option "C" allows to center the header
1090  legend.AddEntry(h1.get(), ("IOV: " + std::to_string(std::get<0>(iov))).c_str(), "PL");
1091  legend.SetTextSize(0.025);
1092  legend.Draw("same");
1093 
1094  std::string fileName(m_imageFileName);
1095  canvas.SaveAs(fileName.c_str());
1096 
1097  return true;
1098  }
1099 
1100  private:
1101  TrackerTopology m_trackerTopo;
1102  };
1103 
1104  typedef SiStripPedestalsByRegion<SiStripPI::mean> SiStripPedestalsMeanByRegion;
1105  typedef SiStripPedestalsByRegion<SiStripPI::min> SiStripPedestalsMinByRegion;
1106  typedef SiStripPedestalsByRegion<SiStripPI::max> SiStripPedestalsMaxByRegion;
1107  typedef SiStripPedestalsByRegion<SiStripPI::rms> SiStripPedestalsRMSByRegion;
1108 
1109 } // namespace
1110 
1112  PAYLOAD_INSPECTOR_CLASS(SiStripPedestalCompareByPartition);
1113  PAYLOAD_INSPECTOR_CLASS(SiStripPedestalDiffByPartition);
1114  PAYLOAD_INSPECTOR_CLASS(SiStripPedestalCorrelationByPartition);
1115  PAYLOAD_INSPECTOR_CLASS(SiStripPedestalsTest);
1116  PAYLOAD_INSPECTOR_CLASS(SiStripPedestalPerDetId);
1117  PAYLOAD_INSPECTOR_CLASS(SiStripPedestalsValue);
1118  PAYLOAD_INSPECTOR_CLASS(SiStripPedestalsValuePerDetId);
1119  PAYLOAD_INSPECTOR_CLASS(SiStripPedestalValuePerStrip);
1120  PAYLOAD_INSPECTOR_CLASS(SiStripPedestalValuePerAPV);
1121  PAYLOAD_INSPECTOR_CLASS(SiStripPedestalValuePerModule);
1122  PAYLOAD_INSPECTOR_CLASS(SiStripPedestalValueComparisonPerStripSingleTag);
1123  PAYLOAD_INSPECTOR_CLASS(SiStripPedestalValueComparisonPerAPVSingleTag);
1124  PAYLOAD_INSPECTOR_CLASS(SiStripPedestalValueComparisonPerModuleSingleTag);
1125  PAYLOAD_INSPECTOR_CLASS(SiStripPedestalValueComparisonPerStripTwoTags);
1126  PAYLOAD_INSPECTOR_CLASS(SiStripPedestalValueComparisonPerAPVTwoTags);
1127  PAYLOAD_INSPECTOR_CLASS(SiStripPedestalValueComparisonPerModuleTwoTags);
1128  PAYLOAD_INSPECTOR_CLASS(SiStripZeroPedestalsFraction_TrackerMap);
1129  PAYLOAD_INSPECTOR_CLASS(SiStripPedestalsMin_TrackerMap);
1130  PAYLOAD_INSPECTOR_CLASS(SiStripPedestalsMax_TrackerMap);
1131  PAYLOAD_INSPECTOR_CLASS(SiStripPedestalsMean_TrackerMap);
1132  PAYLOAD_INSPECTOR_CLASS(SiStripPedestalsRMS_TrackerMap);
1133  PAYLOAD_INSPECTOR_CLASS(SiStripPedestalsMeanByRegion);
1134  PAYLOAD_INSPECTOR_CLASS(SiStripPedestalsMinByRegion);
1135  PAYLOAD_INSPECTOR_CLASS(SiStripPedestalsMaxByRegion);
1136  PAYLOAD_INSPECTOR_CLASS(SiStripPedestalsRMSByRegion);
1137 }
std::pair< int, const char * > regionType(int index)
std::string to_string(const V &value)
Definition: OMSAccess.h:71
std::pair< ContainerIterator, ContainerIterator > Range
assert(be >=bs)
static std::string const input
Definition: EdmProvDump.cc:47
std::tuple< cond::Time_t, cond::Hash > MetaData
#define PAYLOAD_INSPECTOR_CLASS(CLASS_NAME)
T sqrt(T t)
Definition: SSEVec.h:19
void makeNicePlotStyle(TH1 *hist)
SiStripDetInfo read(std::string filePath)
#define PAYLOAD_INSPECTOR_MODULE(PAYLOAD_TYPENAME)
std::string estimatorType(SiStripPI::estimator e)
Log< level::Warning, true > LogPrint
d
Definition: ztail.py:151
__shared__ Hist hist
const std::map< std::string, std::string > & inputParamValues() const
void addInputParam(const std::string &paramName)
#define N
Definition: blowfish.cc:9
double b
Definition: hdecay.h:118
static const uint16_t STRIPS_PER_APV
TrackerTopology fromTrackerParametersXMLFile(const std::string &xmlFileName)
def canvas(sub, attr)
Definition: svgfig.py:482
static constexpr char const *const kDefaultFile
Log< level::Warning, false > LogWarning