File indexing completed on 2024-10-14 22:39:29
0001
0002
0003
0004
0005
0006
0007
0008
0009 #include "CondCore/Utilities/interface/PayloadInspectorModule.h"
0010 #include "CondCore/Utilities/interface/PayloadInspector.h"
0011 #include "CondCore/CondDB/interface/Time.h"
0012
0013
0014 #include "CondFormats/SiStripObjects/interface/SiStripPedestals.h"
0015 #include "DataFormats/DetId/interface/DetId.h"
0016 #include "DataFormats/SiStripDetId/interface/StripSubdetector.h"
0017 #include "CondFormats/SiStripObjects/interface/SiStripDetSummary.h"
0018 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0019
0020
0021 #include "CommonTools/TrackerMap/interface/TrackerMap.h"
0022
0023
0024 #include "CondCore/SiStripPlugins/interface/SiStripPayloadInspectorHelper.h"
0025 #include "CalibTracker/StandaloneTrackerTopology/interface/StandaloneTrackerTopology.h"
0026 #include "CalibTracker/SiStripCommon/interface/SiStripDetInfoFileReader.h"
0027 #include "FWCore/ParameterSet/interface/FileInPath.h"
0028 #include "SiStripCondObjectRepresent.h"
0029
0030 #include <memory>
0031 #include <sstream>
0032 #include <iostream>
0033 #include <iomanip>
0034 #include <boost/tokenizer.hpp>
0035
0036
0037 #include "TH2F.h"
0038 #include "TLegend.h"
0039 #include "TCanvas.h"
0040 #include "TLine.h"
0041 #include "TStyle.h"
0042 #include "TLatex.h"
0043 #include "TPave.h"
0044 #include "TGaxis.h"
0045 #include "TPaveStats.h"
0046
0047 namespace {
0048 using namespace cond::payloadInspector;
0049
0050 class SiStripPedestalContainer : public SiStripCondObjectRepresent::SiStripDataContainer<SiStripPedestals, float> {
0051 public:
0052 SiStripPedestalContainer(const std::shared_ptr<SiStripPedestals>& payload,
0053 const SiStripPI::MetaData& metadata,
0054 const std::string& tagName)
0055 : SiStripCondObjectRepresent::SiStripDataContainer<SiStripPedestals, float>(payload, metadata, tagName) {
0056 payloadType_ = "SiStripPedestals";
0057 setGranularity(SiStripCondObjectRepresent::PERSTRIP);
0058 }
0059
0060 void storeAllValues() override {
0061 std::vector<uint32_t> detid;
0062 payload_->getDetIds(detid);
0063
0064 for (const auto& d : detid) {
0065 SiStripPedestals::Range range = payload_->getRange(d);
0066 for (int it = 0; it < (range.second - range.first) * 8 / 10; ++it) {
0067
0068 SiStripCondData_.fillByPushBack(d, payload_->getPed(it, range));
0069 }
0070 }
0071 }
0072 };
0073
0074 template <int ntags, IOVMultiplicity nIOVs>
0075 class SiStripPedestalCompareByPartitionBase : public PlotImage<SiStripPedestals, nIOVs, ntags> {
0076 public:
0077 SiStripPedestalCompareByPartitionBase()
0078 : PlotImage<SiStripPedestals, nIOVs, ntags>("SiStrip Compare Pedestals By Partition") {}
0079
0080 bool fill() override {
0081
0082 auto theIOVs = PlotBase::getTag<0>().iovs;
0083 auto tagname1 = PlotBase::getTag<0>().name;
0084 std::string tagname2{tagname1};
0085 auto firstiov = theIOVs.front();
0086 SiStripPI::MetaData lastiov;
0087
0088
0089 assert(this->m_plotAnnotations.ntags < 3);
0090
0091 if (this->m_plotAnnotations.ntags == 2) {
0092 auto tag2iovs = PlotBase::getTag<1>().iovs;
0093 tagname2 = PlotBase::getTag<1>().name;
0094 lastiov = tag2iovs.front();
0095 } else {
0096 lastiov = theIOVs.back();
0097 }
0098
0099 std::shared_ptr<SiStripPedestals> last_payload = this->fetchPayload(std::get<1>(lastiov));
0100 std::shared_ptr<SiStripPedestals> first_payload = this->fetchPayload(std::get<1>(firstiov));
0101
0102 SiStripPedestalContainer* l_objContainer = new SiStripPedestalContainer(last_payload, lastiov, tagname1);
0103 SiStripPedestalContainer* f_objContainer = new SiStripPedestalContainer(first_payload, firstiov, tagname2);
0104
0105 l_objContainer->compare(f_objContainer);
0106
0107
0108
0109 TCanvas canvas("Partition summary", "partition summary", 1400, 1000);
0110 l_objContainer->fillByPartition(canvas, 300, 0.1, 300.);
0111
0112 std::string fileName(this->m_imageFileName);
0113 canvas.SaveAs(fileName.c_str());
0114
0115 return true;
0116 }
0117 };
0118
0119 using SiStripPedestalCompareByPartitionSingleTag = SiStripPedestalCompareByPartitionBase<1, MULTI_IOV>;
0120 using SiStripPedestalCompareByPartitionTwoTags = SiStripPedestalCompareByPartitionBase<2, SINGLE_IOV>;
0121
0122 template <int ntags, IOVMultiplicity nIOVs>
0123 class SiStripPedestalDiffByPartitionBase : public PlotImage<SiStripPedestals, nIOVs, ntags> {
0124 public:
0125 SiStripPedestalDiffByPartitionBase()
0126 : PlotImage<SiStripPedestals, nIOVs, ntags>("SiStrip Diff Pedestals By Partition") {}
0127
0128 bool fill() override {
0129
0130 auto theIOVs = PlotBase::getTag<0>().iovs;
0131 auto tagname1 = PlotBase::getTag<0>().name;
0132 std::string tagname2{tagname1};
0133 auto firstiov = theIOVs.front();
0134 SiStripPI::MetaData lastiov;
0135
0136
0137 assert(this->m_plotAnnotations.ntags < 3);
0138
0139 if (this->m_plotAnnotations.ntags == 2) {
0140 auto tag2iovs = PlotBase::getTag<1>().iovs;
0141 tagname2 = PlotBase::getTag<1>().name;
0142 lastiov = tag2iovs.front();
0143 } else {
0144 lastiov = theIOVs.back();
0145 }
0146
0147 std::shared_ptr<SiStripPedestals> last_payload = this->fetchPayload(std::get<1>(lastiov));
0148 std::shared_ptr<SiStripPedestals> first_payload = this->fetchPayload(std::get<1>(firstiov));
0149
0150 SiStripPedestalContainer* l_objContainer = new SiStripPedestalContainer(last_payload, lastiov, tagname1);
0151 SiStripPedestalContainer* f_objContainer = new SiStripPedestalContainer(first_payload, firstiov, tagname2);
0152
0153 l_objContainer->subtract(f_objContainer);
0154
0155
0156
0157 TCanvas canvas("Partition summary", "partition summary", 1400, 1000);
0158 l_objContainer->fillByPartition(canvas, 100, -30., 30.);
0159
0160 std::string fileName(this->m_imageFileName);
0161 canvas.SaveAs(fileName.c_str());
0162
0163 return true;
0164 }
0165 };
0166
0167 using SiStripPedestalDiffByPartitionSingleTag = SiStripPedestalDiffByPartitionBase<1, MULTI_IOV>;
0168 using SiStripPedestalDiffByPartitionTwoTags = SiStripPedestalDiffByPartitionBase<2, SINGLE_IOV>;
0169
0170 class SiStripPedestalCorrelationByPartition : public PlotImage<SiStripPedestals> {
0171 public:
0172 SiStripPedestalCorrelationByPartition()
0173 : PlotImage<SiStripPedestals>("SiStrip Pedestals Correlation By Partition") {
0174 setSingleIov(false);
0175 }
0176
0177 bool fill(const std::vector<SiStripPI::MetaData>& iovs) override {
0178 std::vector<SiStripPI::MetaData> sorted_iovs = iovs;
0179
0180
0181 std::sort(begin(sorted_iovs), end(sorted_iovs), [](auto const& t1, auto const& t2) {
0182 return std::get<0>(t1) < std::get<0>(t2);
0183 });
0184
0185 auto firstiov = sorted_iovs.front();
0186 auto lastiov = sorted_iovs.back();
0187
0188 std::shared_ptr<SiStripPedestals> last_payload = fetchPayload(std::get<1>(lastiov));
0189 std::shared_ptr<SiStripPedestals> first_payload = fetchPayload(std::get<1>(firstiov));
0190
0191 SiStripPedestalContainer* l_objContainer = new SiStripPedestalContainer(last_payload, lastiov, "");
0192 SiStripPedestalContainer* f_objContainer = new SiStripPedestalContainer(first_payload, firstiov, "");
0193
0194 l_objContainer->compare(f_objContainer);
0195
0196 TCanvas canvas("Partition summary", "partition summary", 1200, 1200);
0197 l_objContainer->fillCorrelationByPartition(canvas, 100, 0., 300.);
0198
0199 std::string fileName(m_imageFileName);
0200 canvas.SaveAs(fileName.c_str());
0201
0202 return true;
0203 }
0204 };
0205
0206
0207
0208
0209
0210 class SiStripPedestalsTest : public Histogram1D<SiStripPedestals, SINGLE_IOV> {
0211 public:
0212 SiStripPedestalsTest()
0213 : Histogram1D<SiStripPedestals, SINGLE_IOV>("SiStrip Pedestals test", "SiStrip Pedestals test", 10, 0.0, 10.0),
0214 m_trackerTopo{StandaloneTrackerTopology::fromTrackerParametersXMLFile(
0215 edm::FileInPath("Geometry/TrackerCommonData/data/trackerParameters.xml").fullPath())} {}
0216
0217 bool fill() override {
0218 auto tag = PlotBase::getTag<0>();
0219 for (auto const& iov : tag.iovs) {
0220 std::shared_ptr<SiStripPedestals> payload = Base::fetchPayload(std::get<1>(iov));
0221 if (payload.get()) {
0222 fillWithValue(1.);
0223
0224 std::stringstream ss;
0225 ss << "Summary of strips pedestals:" << std::endl;
0226
0227
0228 payload->printSummary(ss, &m_trackerTopo);
0229
0230 std::vector<uint32_t> detid;
0231 payload->getDetIds(detid);
0232
0233
0234
0235
0236
0237
0238
0239
0240
0241
0242
0243 std::cout << ss.str() << std::endl;
0244
0245 }
0246 }
0247 return true;
0248 }
0249 private:
0250 TrackerTopology m_trackerTopo;
0251 };
0252
0253
0254
0255
0256
0257 class SiStripPedestalPerDetId : public PlotImage<SiStripPedestals, SINGLE_IOV> {
0258 public:
0259 SiStripPedestalPerDetId() : PlotImage<SiStripPedestals, SINGLE_IOV>("SiStrip Pedestal values per DetId") {
0260 PlotBase::addInputParam("DetIds");
0261 }
0262
0263 bool fill() override {
0264 auto tag = PlotBase::getTag<0>();
0265 auto iov = tag.iovs.front();
0266 auto tagname = tag.name;
0267 std::shared_ptr<SiStripPedestals> payload = fetchPayload(std::get<1>(iov));
0268
0269 std::vector<uint32_t> the_detids = {};
0270
0271 auto paramValues = PlotBase::inputParamValues();
0272 auto ip = paramValues.find("DetIds");
0273 if (ip != paramValues.end()) {
0274 auto input = ip->second;
0275 typedef boost::tokenizer<boost::char_separator<char>> tokenizer;
0276 boost::char_separator<char> sep{","};
0277 tokenizer tok{input, sep};
0278 for (const auto& t : tok) {
0279 the_detids.push_back(atoi(t.c_str()));
0280 }
0281 } else {
0282 edm::LogWarning("SiStripPedestalPerDetId")
0283 << "\n WARNING!!!! \n The needed parameter DetIds has not been passed. Will use all Strip DetIds! \n\n";
0284 the_detids.push_back(0xFFFFFFFF);
0285 }
0286
0287 size_t ndets = the_detids.size();
0288 std::vector<std::shared_ptr<TH1F>> hpedestal;
0289 std::vector<std::shared_ptr<TLegend>> legends;
0290 std::vector<unsigned int> v_nAPVs;
0291 std::vector<std::vector<std::shared_ptr<TLine>>> lines;
0292 hpedestal.reserve(ndets);
0293 legends.reserve(ndets);
0294
0295 auto sides = getClosestFactors(the_detids.size());
0296 edm::LogPrint("SiStripPedestalPerDetId") << "Aspect ratio: " << sides.first << ":" << sides.second << std::endl;
0297
0298 if (payload.get()) {
0299
0300 TCanvas canvas("ByDetId", "ByDetId", sides.second * 800, sides.first * 600);
0301 canvas.Divide(sides.second, sides.first);
0302 const auto detInfo =
0303 SiStripDetInfoFileReader::read(edm::FileInPath(SiStripDetInfoFileReader::kDefaultFile).fullPath());
0304 for (const auto& the_detid : the_detids) {
0305 edm::LogPrint("SiStripPedestalPerDetId") << "DetId:" << the_detid << std::endl;
0306
0307 unsigned int nAPVs = detInfo.getNumberOfApvsAndStripLength(the_detid).first;
0308 if (nAPVs == 0)
0309 nAPVs = 6;
0310 v_nAPVs.push_back(nAPVs);
0311
0312 auto histo = std::make_shared<TH1F>(
0313 Form("Pedestal profile %s", std::to_string(the_detid).c_str()),
0314 Form("SiStrip Pedestal profile for DetId: %s;Strip number;SiStrip Pedestal [ADC counts]",
0315 std::to_string(the_detid).c_str()),
0316 sistrip::STRIPS_PER_APV * nAPVs,
0317 -0.5,
0318 (sistrip::STRIPS_PER_APV * nAPVs) - 0.5);
0319
0320 histo->SetStats(false);
0321 histo->SetTitle("");
0322
0323 if (the_detid != 0xFFFFFFFF) {
0324 fillHisto(payload, histo, the_detid);
0325 } else {
0326 auto allDetIds = detInfo.getAllDetIds();
0327 for (const auto& id : allDetIds) {
0328 fillHisto(payload, histo, id);
0329 }
0330 }
0331
0332 SiStripPI::makeNicePlotStyle(histo.get());
0333 histo->GetYaxis()->SetTitleOffset(1.0);
0334 hpedestal.push_back(histo);
0335 }
0336
0337 for (size_t index = 0; index < ndets; index++) {
0338 canvas.cd(index + 1);
0339 canvas.cd(index + 1)->SetBottomMargin(0.11);
0340 canvas.cd(index + 1)->SetTopMargin(0.06);
0341 canvas.cd(index + 1)->SetLeftMargin(0.10);
0342 canvas.cd(index + 1)->SetRightMargin(0.02);
0343 hpedestal.at(index)->Draw();
0344 hpedestal.at(index)->GetYaxis()->SetRangeUser(0, hpedestal.at(index)->GetMaximum() * 1.2);
0345 canvas.cd(index)->Update();
0346
0347 std::vector<int> boundaries;
0348 for (size_t b = 0; b < v_nAPVs.at(index); b++) {
0349 boundaries.push_back(b * sistrip::STRIPS_PER_APV);
0350 }
0351
0352 std::vector<std::shared_ptr<TLine>> linesVec;
0353 for (const auto& bound : boundaries) {
0354 auto line = std::make_shared<TLine>(hpedestal.at(index)->GetBinLowEdge(bound),
0355 canvas.cd(index + 1)->GetUymin(),
0356 hpedestal.at(index)->GetBinLowEdge(bound),
0357 canvas.cd(index + 1)->GetUymax());
0358 line->SetLineWidth(1);
0359 line->SetLineStyle(9);
0360 line->SetLineColor(2);
0361 linesVec.push_back(line);
0362 }
0363 lines.push_back(linesVec);
0364
0365 for (const auto& line : lines.at(index)) {
0366 line->Draw("same");
0367 }
0368
0369 canvas.cd(index + 1);
0370
0371 auto ltx = TLatex();
0372 ltx.SetTextFont(62);
0373 ltx.SetTextSize(0.05);
0374 ltx.SetTextAlign(11);
0375 ltx.DrawLatexNDC(gPad->GetLeftMargin(),
0376 1 - gPad->GetTopMargin() + 0.01,
0377 Form("SiStrip Pedestals profile for DetId %s", std::to_string(the_detids[index]).c_str()));
0378
0379 legends.push_back(std::make_shared<TLegend>(0.45, 0.83, 0.95, 0.93));
0380 legends.at(index)->SetHeader(tagname.c_str(), "C");
0381 legends.at(index)->AddEntry(
0382 hpedestal.at(index).get(), ("IOV: " + std::to_string(std::get<0>(iov))).c_str(), "PL");
0383 legends.at(index)->SetTextSize(0.045);
0384 legends.at(index)->Draw("same");
0385 }
0386
0387 std::string fileName(m_imageFileName);
0388 canvas.SaveAs(fileName.c_str());
0389 }
0390 return true;
0391 }
0392
0393 private:
0394 int nextPerfectSquare(int N) { return std::floor(sqrt(N)) + 1; }
0395
0396 std::pair<int, int> getClosestFactors(int input) {
0397 if ((input % 2 != 0) && input > 1) {
0398 input += 1;
0399 }
0400
0401 int testNum = (int)sqrt(input);
0402 while (input % testNum != 0) {
0403 testNum--;
0404 }
0405 return std::make_pair(testNum, input / testNum);
0406 }
0407
0408 void fillHisto(const std::shared_ptr<SiStripPedestals> payload, std::shared_ptr<TH1F>& histo, uint32_t the_detid) {
0409 int nstrip = 0;
0410 SiStripPedestals::Range range = payload->getRange(the_detid);
0411 for (int it = 0; it < (range.second - range.first) * 8 / 10; ++it) {
0412 auto pedestal = payload->getPed(it, range);
0413 nstrip++;
0414 histo->AddBinContent(nstrip, pedestal);
0415 }
0416 }
0417 };
0418
0419
0420
0421
0422
0423
0424 class SiStripPedestalsValue : public Histogram1D<SiStripPedestals, SINGLE_IOV> {
0425 public:
0426 SiStripPedestalsValue()
0427 : Histogram1D<SiStripPedestals, SINGLE_IOV>(
0428 "SiStrip Pedestals values", "SiStrip Pedestals values", 300, 0.0, 300.0) {}
0429
0430 bool fill() override {
0431 auto tag = PlotBase::getTag<0>();
0432 for (auto const& iov : tag.iovs) {
0433 std::shared_ptr<SiStripPedestals> payload = Base::fetchPayload(std::get<1>(iov));
0434 if (payload.get()) {
0435 std::vector<uint32_t> detid;
0436 payload->getDetIds(detid);
0437
0438 for (const auto& d : detid) {
0439 SiStripPedestals::Range range = payload->getRange(d);
0440 for (int it = 0; it < (range.second - range.first) * 8 / 10; ++it) {
0441 auto ped = payload->getPed(it, range);
0442
0443 fillWithValue(ped);
0444 }
0445 }
0446 }
0447 }
0448 return true;
0449 }
0450 };
0451
0452
0453
0454
0455
0456
0457 class SiStripPedestalsValuePerDetId : public Histogram1D<SiStripPedestals, SINGLE_IOV> {
0458 public:
0459 SiStripPedestalsValuePerDetId()
0460 : Histogram1D<SiStripPedestals, SINGLE_IOV>(
0461 "SiStrip Pedestal values per DetId", "SiStrip Pedestal values per DetId", 100, 0.0, 10.0) {
0462 PlotBase::addInputParam("DetId");
0463 }
0464
0465 bool fill() override {
0466 auto tag = PlotBase::getTag<0>();
0467 for (auto const& iov : tag.iovs) {
0468 std::shared_ptr<SiStripPedestals> payload = Base::fetchPayload(std::get<1>(iov));
0469 unsigned int the_detid(0xFFFFFFFF);
0470 auto paramValues = PlotBase::inputParamValues();
0471 auto ip = paramValues.find("DetId");
0472 if (ip != paramValues.end()) {
0473 the_detid = std::stoul(ip->second);
0474 }
0475
0476 if (payload.get()) {
0477 SiStripPedestals::Range range = payload->getRange(the_detid);
0478 for (int it = 0; it < (range.second - range.first) * 8 / 9; ++it) {
0479 auto noise = payload->getPed(it, range);
0480
0481 fillWithValue(noise);
0482 }
0483 }
0484 }
0485 return true;
0486 }
0487 };
0488
0489
0490
0491
0492
0493
0494 template <SiStripPI::OpMode op_mode_>
0495 class SiStripPedestalDistribution : public PlotImage<SiStripPedestals, SINGLE_IOV> {
0496 public:
0497 SiStripPedestalDistribution() : PlotImage<SiStripPedestals, SINGLE_IOV>("SiStrip Pedestal values") {}
0498
0499 bool fill() override {
0500 auto tag = PlotBase::getTag<0>();
0501 auto iov = tag.iovs.front();
0502 TGaxis::SetMaxDigits(3);
0503 gStyle->SetOptStat("emr");
0504
0505 std::shared_ptr<SiStripPedestals> payload = fetchPayload(std::get<1>(iov));
0506
0507 auto mon1D = std::unique_ptr<SiStripPI::Monitor1D>(new SiStripPI::Monitor1D(
0508 op_mode_,
0509 "Pedestal",
0510 Form("#LT Pedestal #GT per %s for IOV [%s];#LTStrip Pedestal per %s#GT [ADC counts];n. %ss",
0511 opType(op_mode_).c_str(),
0512 std::to_string(std::get<0>(iov)).c_str(),
0513 opType(op_mode_).c_str(),
0514 opType(op_mode_).c_str()),
0515 300,
0516 0.,
0517 300.0));
0518
0519 unsigned int prev_det = 0, prev_apv = 0;
0520 SiStripPI::Entry epedestal;
0521
0522 std::vector<uint32_t> detids;
0523 payload->getDetIds(detids);
0524
0525
0526 for (const auto& d : detids) {
0527 SiStripPedestals::Range range = payload->getRange(d);
0528
0529 unsigned int istrip = 0;
0530
0531 for (int it = 0; it < (range.second - range.first) * 8 / 10; ++it) {
0532 auto pedestal = payload->getPed(it, range);
0533 bool flush = false;
0534 switch (op_mode_) {
0535 case (SiStripPI::APV_BASED):
0536 flush = (prev_det != 0 && prev_apv != istrip / sistrip::STRIPS_PER_APV);
0537 break;
0538 case (SiStripPI::MODULE_BASED):
0539 flush = (prev_det != 0 && prev_det != d);
0540 break;
0541 case (SiStripPI::STRIP_BASED):
0542 flush = (istrip != 0);
0543 break;
0544 }
0545
0546 if (flush) {
0547 mon1D->Fill(prev_apv, prev_det, epedestal.mean());
0548 epedestal.reset();
0549 }
0550
0551 epedestal.add(std::min<float>(pedestal, 300.));
0552 prev_apv = istrip / sistrip::STRIPS_PER_APV;
0553 istrip++;
0554 }
0555 prev_det = d;
0556 }
0557
0558
0559 TCanvas canvas("Partion summary", "partition summary", 1200, 1000);
0560 canvas.cd();
0561 canvas.SetBottomMargin(0.11);
0562 canvas.SetTopMargin(0.07);
0563 canvas.SetLeftMargin(0.13);
0564 canvas.SetRightMargin(0.05);
0565 canvas.Modified();
0566
0567 auto hist = mon1D->getHist();
0568 SiStripPI::makeNicePlotStyle(&hist);
0569 hist.SetStats(kTRUE);
0570 hist.SetFillColorAlpha(kRed, 0.35);
0571 hist.Draw();
0572
0573 canvas.Update();
0574
0575 TPaveStats* st = (TPaveStats*)hist.GetListOfFunctions()->FindObject("stats");
0576 st->SetLineColor(kRed);
0577 st->SetTextColor(kRed);
0578 st->SetX1NDC(.75);
0579 st->SetX2NDC(.95);
0580 st->SetY1NDC(.83);
0581 st->SetY2NDC(.93);
0582
0583 TLegend legend = TLegend(0.13, 0.83, 0.43, 0.93);
0584 legend.SetHeader(Form("SiStrip Pedestal values per %s", opType(op_mode_).c_str()),
0585 "C");
0586 legend.AddEntry(&hist, ("IOV: " + std::to_string(std::get<0>(iov))).c_str(), "F");
0587 legend.SetTextSize(0.025);
0588 legend.Draw("same");
0589
0590 std::string fileName(m_imageFileName);
0591 canvas.SaveAs(fileName.c_str());
0592
0593 return true;
0594 }
0595
0596 std::string opType(SiStripPI::OpMode mode) {
0597 std::string types[3] = {"Strip", "APV", "Module"};
0598 return types[mode];
0599 }
0600 };
0601
0602 typedef SiStripPedestalDistribution<SiStripPI::STRIP_BASED> SiStripPedestalValuePerStrip;
0603 typedef SiStripPedestalDistribution<SiStripPI::APV_BASED> SiStripPedestalValuePerAPV;
0604 typedef SiStripPedestalDistribution<SiStripPI::MODULE_BASED> SiStripPedestalValuePerModule;
0605
0606
0607
0608
0609
0610
0611
0612 template <SiStripPI::OpMode op_mode_, int ntags, IOVMultiplicity nIOVs>
0613 class SiStripPedestalDistributionComparisonBase : public PlotImage<SiStripPedestals, nIOVs, ntags> {
0614 public:
0615 SiStripPedestalDistributionComparisonBase()
0616 : PlotImage<SiStripPedestals, nIOVs, ntags>("SiStrip Pedestal values comparison") {}
0617
0618 bool fill() override {
0619 TGaxis::SetExponentOffset(-0.1, 0.01, "y");
0620
0621
0622 auto theIOVs = PlotBase::getTag<0>().iovs;
0623 auto tagname1 = PlotBase::getTag<0>().name;
0624 std::string tagname2 = "";
0625 auto firstiov = theIOVs.front();
0626 SiStripPI::MetaData lastiov;
0627
0628
0629 assert(this->m_plotAnnotations.ntags < 3);
0630
0631 if (this->m_plotAnnotations.ntags == 2) {
0632 auto tag2iovs = PlotBase::getTag<1>().iovs;
0633 tagname2 = PlotBase::getTag<1>().name;
0634 lastiov = tag2iovs.front();
0635 } else {
0636 lastiov = theIOVs.back();
0637 }
0638
0639 std::shared_ptr<SiStripPedestals> f_payload = this->fetchPayload(std::get<1>(firstiov));
0640 std::shared_ptr<SiStripPedestals> l_payload = this->fetchPayload(std::get<1>(lastiov));
0641
0642 auto f_mon = std::unique_ptr<SiStripPI::Monitor1D>(new SiStripPI::Monitor1D(
0643 op_mode_,
0644 "f_Pedestal",
0645 Form(";#LTStrip Pedestal per %s#GT [ADC counts];n. %ss", opType(op_mode_).c_str(), opType(op_mode_).c_str()),
0646 300,
0647 0.,
0648 300.));
0649
0650 auto l_mon = std::unique_ptr<SiStripPI::Monitor1D>(new SiStripPI::Monitor1D(
0651 op_mode_,
0652 "l_Pedestal",
0653 Form(";#LTStrip Pedestal per %s#GT [ADC counts];n. %ss", opType(op_mode_).c_str(), opType(op_mode_).c_str()),
0654 300,
0655 0.,
0656 300.));
0657
0658 unsigned int prev_det = 0, prev_apv = 0;
0659 SiStripPI::Entry epedestal;
0660
0661 std::vector<uint32_t> f_detid;
0662 f_payload->getDetIds(f_detid);
0663
0664
0665 for (const auto& d : f_detid) {
0666 SiStripPedestals::Range range = f_payload->getRange(d);
0667
0668 unsigned int istrip = 0;
0669 for (int it = 0; it < (range.second - range.first) * 8 / 10; ++it) {
0670 float pedestal = f_payload->getPed(it, range);
0671
0672
0673 bool flush = false;
0674 switch (op_mode_) {
0675 case (SiStripPI::APV_BASED):
0676 flush = (prev_det != 0 && prev_apv != istrip / sistrip::STRIPS_PER_APV);
0677 break;
0678 case (SiStripPI::MODULE_BASED):
0679 flush = (prev_det != 0 && prev_det != d);
0680 break;
0681 case (SiStripPI::STRIP_BASED):
0682 flush = (istrip != 0);
0683 break;
0684 }
0685
0686 if (flush) {
0687 f_mon->Fill(prev_apv, prev_det, epedestal.mean());
0688 epedestal.reset();
0689 }
0690 epedestal.add(std::min<float>(pedestal, 300.));
0691 prev_apv = istrip / sistrip::STRIPS_PER_APV;
0692 istrip++;
0693 }
0694 prev_det = d;
0695 }
0696
0697 prev_det = 0, prev_apv = 0;
0698 epedestal.reset();
0699
0700 std::vector<uint32_t> l_detid;
0701 l_payload->getDetIds(l_detid);
0702
0703
0704 for (const auto& d : l_detid) {
0705 SiStripPedestals::Range range = l_payload->getRange(d);
0706
0707 unsigned int istrip = 0;
0708 for (int it = 0; it < (range.second - range.first) * 8 / 10; ++it) {
0709 float pedestal = l_payload->getPed(it, range);
0710
0711 bool flush = false;
0712 switch (op_mode_) {
0713 case (SiStripPI::APV_BASED):
0714 flush = (prev_det != 0 && prev_apv != istrip / sistrip::STRIPS_PER_APV);
0715 break;
0716 case (SiStripPI::MODULE_BASED):
0717 flush = (prev_det != 0 && prev_det != d);
0718 break;
0719 case (SiStripPI::STRIP_BASED):
0720 flush = (istrip != 0);
0721 break;
0722 }
0723
0724 if (flush) {
0725 l_mon->Fill(prev_apv, prev_det, epedestal.mean());
0726 epedestal.reset();
0727 }
0728
0729 epedestal.add(std::min<float>(pedestal, 300.));
0730 prev_apv = istrip / sistrip::STRIPS_PER_APV;
0731 istrip++;
0732 }
0733 prev_det = d;
0734 }
0735
0736 auto h_first = f_mon->getHist();
0737 h_first.SetStats(kFALSE);
0738 auto h_last = l_mon->getHist();
0739 h_last.SetStats(kFALSE);
0740
0741 SiStripPI::makeNicePlotStyle(&h_first);
0742 SiStripPI::makeNicePlotStyle(&h_last);
0743
0744 h_first.GetYaxis()->CenterTitle(true);
0745 h_last.GetYaxis()->CenterTitle(true);
0746
0747 h_first.GetXaxis()->CenterTitle(true);
0748 h_last.GetXaxis()->CenterTitle(true);
0749
0750 h_first.SetLineWidth(2);
0751 h_last.SetLineWidth(2);
0752
0753 h_first.SetLineColor(kBlack);
0754 h_last.SetLineColor(kBlue);
0755
0756
0757 TCanvas canvas("Partion summary", "partition summary", 1200, 1000);
0758 canvas.cd();
0759 canvas.SetTopMargin(0.06);
0760 canvas.SetBottomMargin(0.10);
0761 canvas.SetLeftMargin(0.13);
0762 canvas.SetRightMargin(0.05);
0763 canvas.Modified();
0764
0765 float theMax = (h_first.GetMaximum() > h_last.GetMaximum()) ? h_first.GetMaximum() : h_last.GetMaximum();
0766
0767 h_first.SetMaximum(theMax * 1.20);
0768 h_last.SetMaximum(theMax * 1.20);
0769
0770 h_first.Draw();
0771 h_last.SetFillColorAlpha(kBlue, 0.15);
0772 h_last.Draw("same");
0773
0774 TLegend legend = TLegend(0.13, 0.83, 0.95, 0.94);
0775 if (this->m_plotAnnotations.ntags == 2) {
0776 legend.SetHeader("#bf{Two Tags Comparison}", "C");
0777 legend.AddEntry(&h_first, (tagname1 + " : " + std::to_string(std::get<0>(firstiov))).c_str(), "F");
0778 legend.AddEntry(&h_last, (tagname2 + " : " + std::to_string(std::get<0>(lastiov))).c_str(), "F");
0779 } else {
0780 legend.SetHeader(("tag: #bf{" + tagname1 + "}").c_str(), "C");
0781 legend.AddEntry(&h_first, ("IOV since: " + std::to_string(std::get<0>(firstiov))).c_str(), "F");
0782 legend.AddEntry(&h_last, ("IOV since: " + std::to_string(std::get<0>(lastiov))).c_str(), "F");
0783 }
0784 legend.SetTextSize(0.025);
0785 legend.Draw("same");
0786
0787 auto ltx = TLatex();
0788 ltx.SetTextFont(62);
0789 ltx.SetTextSize(0.05);
0790 ltx.SetTextAlign(11);
0791 ltx.DrawLatexNDC(gPad->GetLeftMargin(),
0792 1 - gPad->GetTopMargin() + 0.01,
0793 Form("#LTSiStrip Pedestals#GT Comparison per %s", opType(op_mode_).c_str()));
0794
0795 std::string fileName(this->m_imageFileName);
0796 canvas.SaveAs(fileName.c_str());
0797
0798 return true;
0799 }
0800
0801 std::string opType(SiStripPI::OpMode mode) {
0802 std::string types[3] = {"Strip", "APV", "Module"};
0803 return types[mode];
0804 }
0805 };
0806
0807 template <SiStripPI::OpMode op_mode_>
0808 using SiStripPedestalDistributionComparisonSingleTag =
0809 SiStripPedestalDistributionComparisonBase<op_mode_, 1, MULTI_IOV>;
0810
0811 template <SiStripPI::OpMode op_mode_>
0812 using SiStripPedestalDistributionComparisonTwoTags =
0813 SiStripPedestalDistributionComparisonBase<op_mode_, 2, SINGLE_IOV>;
0814
0815 typedef SiStripPedestalDistributionComparisonSingleTag<SiStripPI::STRIP_BASED>
0816 SiStripPedestalValueComparisonPerStripSingleTag;
0817 typedef SiStripPedestalDistributionComparisonSingleTag<SiStripPI::APV_BASED>
0818 SiStripPedestalValueComparisonPerAPVSingleTag;
0819 typedef SiStripPedestalDistributionComparisonSingleTag<SiStripPI::MODULE_BASED>
0820 SiStripPedestalValueComparisonPerModuleSingleTag;
0821
0822 typedef SiStripPedestalDistributionComparisonTwoTags<SiStripPI::STRIP_BASED>
0823 SiStripPedestalValueComparisonPerStripTwoTags;
0824 typedef SiStripPedestalDistributionComparisonTwoTags<SiStripPI::APV_BASED> SiStripPedestalValueComparisonPerAPVTwoTags;
0825 typedef SiStripPedestalDistributionComparisonTwoTags<SiStripPI::MODULE_BASED>
0826 SiStripPedestalValueComparisonPerModuleTwoTags;
0827
0828
0829
0830
0831
0832
0833 class SiStripZeroPedestalsFraction_TrackerMap : public PlotImage<SiStripPedestals, SINGLE_IOV> {
0834 public:
0835 SiStripZeroPedestalsFraction_TrackerMap()
0836 : PlotImage<SiStripPedestals, SINGLE_IOV>("Tracker Map of Zero SiStripPedestals fraction per module") {}
0837
0838 bool fill() override {
0839 auto tag = PlotBase::getTag<0>();
0840 auto iov = tag.iovs.front();
0841 std::shared_ptr<SiStripPedestals> payload = fetchPayload(std::get<1>(iov));
0842
0843 const auto detInfo =
0844 SiStripDetInfoFileReader::read(edm::FileInPath(SiStripDetInfoFileReader::kDefaultFile).fullPath());
0845
0846 std::string titleMap =
0847 "Tracker Map of Zero SiStrip Pedestals fraction per module (payload : " + std::get<1>(iov) + ")";
0848
0849 std::unique_ptr<TrackerMap> tmap = std::make_unique<TrackerMap>("SiStripPedestals");
0850 tmap->setTitle(titleMap);
0851 tmap->setPalette(1);
0852
0853 std::vector<uint32_t> detid;
0854 payload->getDetIds(detid);
0855
0856 std::map<uint32_t, int> zeropeds_per_detid;
0857
0858 for (const auto& d : detid) {
0859 SiStripPedestals::Range range = payload->getRange(d);
0860 for (int it = 0; it < (range.second - range.first) * 8 / 10; ++it) {
0861 auto ped = payload->getPed(it, range);
0862 if (ped == 0.) {
0863 zeropeds_per_detid[d] += 1;
0864 }
0865 }
0866 float fraction =
0867 zeropeds_per_detid[d] / (sistrip::STRIPS_PER_APV * detInfo.getNumberOfApvsAndStripLength(d).first);
0868 if (fraction > 0.) {
0869 tmap->fill(d, fraction);
0870 std::cout << "detid: " << d << " (n. APVs=" << detInfo.getNumberOfApvsAndStripLength(d).first << ") has "
0871 << std::setw(4) << zeropeds_per_detid[d]
0872 << " zero-pedestals strips (i.e. a fraction:" << std::setprecision(5) << fraction << ")"
0873 << std::endl;
0874 }
0875 }
0876
0877 std::string fileName(m_imageFileName);
0878 tmap->save(true, 0., 0., fileName);
0879
0880 return true;
0881 }
0882 };
0883
0884
0885
0886
0887
0888 template <SiStripPI::estimator est>
0889 class SiStripPedestalsTrackerMap : public PlotImage<SiStripPedestals, SINGLE_IOV> {
0890 public:
0891 SiStripPedestalsTrackerMap()
0892 : PlotImage<SiStripPedestals, SINGLE_IOV>("Tracker Map of SiStripPedestals " + estimatorType(est) +
0893 " per module") {}
0894
0895 bool fill() override {
0896 auto tag = PlotBase::getTag<0>();
0897 auto iov = tag.iovs.front();
0898 std::shared_ptr<SiStripPedestals> payload = fetchPayload(std::get<1>(iov));
0899
0900 std::string titleMap =
0901 "Tracker Map of SiStrip Pedestals " + estimatorType(est) + " per module (payload : " + std::get<1>(iov) + ")";
0902
0903 std::unique_ptr<TrackerMap> tmap = std::make_unique<TrackerMap>("SiStripPedestals");
0904 tmap->setTitle(titleMap);
0905 tmap->setPalette(1);
0906
0907 std::vector<uint32_t> detid;
0908 payload->getDetIds(detid);
0909
0910 std::map<unsigned int, float> info_per_detid;
0911
0912 for (const auto& d : detid) {
0913 int nstrips = 0;
0914 double mean(0.), rms(0.), min(10000.), max(0.);
0915 SiStripPedestals::Range range = payload->getRange(d);
0916 for (int it = 0; it < (range.second - range.first) * 8 / 10; ++it) {
0917 nstrips++;
0918 auto ped = payload->getPed(it, range);
0919 mean += ped;
0920 rms += (ped * ped);
0921 if (ped < min)
0922 min = ped;
0923 if (ped > max)
0924 max = ped;
0925 }
0926
0927 mean /= nstrips;
0928 if ((rms / nstrips - mean * mean) > 0.) {
0929 rms = sqrt(rms / nstrips - mean * mean);
0930 } else {
0931 rms = 0.;
0932 }
0933
0934 switch (est) {
0935 case SiStripPI::min:
0936 info_per_detid[d] = min;
0937 break;
0938 case SiStripPI::max:
0939 info_per_detid[d] = max;
0940 break;
0941 case SiStripPI::mean:
0942 info_per_detid[d] = mean;
0943 break;
0944 case SiStripPI::rms:
0945 info_per_detid[d] = rms;
0946 break;
0947 default:
0948 edm::LogWarning("LogicError") << "Unknown estimator: " << est;
0949 break;
0950 }
0951 }
0952
0953 for (const auto& d : detid) {
0954 tmap->fill(d, info_per_detid[d]);
0955 }
0956
0957 std::string fileName(m_imageFileName);
0958 tmap->save(true, 0., 0., fileName);
0959
0960 return true;
0961 }
0962 };
0963
0964 typedef SiStripPedestalsTrackerMap<SiStripPI::min> SiStripPedestalsMin_TrackerMap;
0965 typedef SiStripPedestalsTrackerMap<SiStripPI::max> SiStripPedestalsMax_TrackerMap;
0966 typedef SiStripPedestalsTrackerMap<SiStripPI::mean> SiStripPedestalsMean_TrackerMap;
0967 typedef SiStripPedestalsTrackerMap<SiStripPI::rms> SiStripPedestalsRMS_TrackerMap;
0968
0969
0970
0971
0972
0973 template <SiStripPI::estimator est>
0974 class SiStripPedestalsByRegion : public PlotImage<SiStripPedestals, SINGLE_IOV> {
0975 public:
0976 SiStripPedestalsByRegion()
0977 : PlotImage<SiStripPedestals, SINGLE_IOV>("SiStrip Pedestals " + estimatorType(est) + " by Region"),
0978 m_trackerTopo{StandaloneTrackerTopology::fromTrackerParametersXMLFile(
0979 edm::FileInPath("Geometry/TrackerCommonData/data/trackerParameters.xml").fullPath())} {}
0980
0981 bool fill() override {
0982 auto tag = PlotBase::getTag<0>();
0983 auto iov = tag.iovs.front();
0984 std::shared_ptr<SiStripPedestals> payload = fetchPayload(std::get<1>(iov));
0985
0986 SiStripDetSummary summaryPedestals{&m_trackerTopo};
0987 std::vector<uint32_t> detid;
0988 payload->getDetIds(detid);
0989
0990 for (const auto& d : detid) {
0991 int nstrips = 0;
0992 double mean(0.), rms(0.), min(10000.), max(0.);
0993 SiStripPedestals::Range range = payload->getRange(d);
0994 for (int it = 0; it < (range.second - range.first) * 8 / 10; ++it) {
0995 nstrips++;
0996 auto ped = payload->getPed(it, range);
0997 mean += ped;
0998 rms += (ped * ped);
0999 if (ped < min)
1000 min = ped;
1001 if (ped > max)
1002 max = ped;
1003 }
1004
1005 mean /= nstrips;
1006 if ((rms / nstrips - mean * mean) > 0.) {
1007 rms = sqrt(rms / nstrips - mean * mean);
1008 } else {
1009 rms = 0.;
1010 }
1011
1012 switch (est) {
1013 case SiStripPI::min:
1014 summaryPedestals.add(d, min);
1015 break;
1016 case SiStripPI::max:
1017 summaryPedestals.add(d, max);
1018 break;
1019 case SiStripPI::mean:
1020 summaryPedestals.add(d, mean);
1021 break;
1022 case SiStripPI::rms:
1023 summaryPedestals.add(d, rms);
1024 break;
1025 default:
1026 edm::LogWarning("LogicError") << "Unknown estimator: " << est;
1027 break;
1028 }
1029 }
1030
1031 std::map<unsigned int, SiStripDetSummary::Values> map = summaryPedestals.getCounts();
1032
1033
1034 TCanvas canvas("Partion summary", "partition summary", 1200, 1000);
1035 canvas.cd();
1036 auto h1 = std::unique_ptr<TH1F>(new TH1F(
1037 "byRegion",
1038 Form("Average by partition of %s SiStrip Pedestals per module;;average SiStrip Pedestals %s [ADC counts]",
1039 estimatorType(est).c_str(),
1040 estimatorType(est).c_str()),
1041 map.size(),
1042 0.,
1043 map.size()));
1044 h1->SetStats(false);
1045 canvas.SetBottomMargin(0.18);
1046 canvas.SetLeftMargin(0.17);
1047 canvas.SetRightMargin(0.05);
1048 canvas.Modified();
1049
1050 std::vector<int> boundaries;
1051 unsigned int iBin = 0;
1052
1053 std::string detector;
1054 std::string currentDetector;
1055
1056 for (const auto& element : map) {
1057 iBin++;
1058 int count = element.second.count;
1059 double mean = (element.second.mean) / count;
1060 double rms = (element.second.rms) / count - mean * mean;
1061
1062 if (rms <= 0)
1063 rms = 0;
1064 else
1065 rms = sqrt(rms);
1066
1067 if (currentDetector.empty())
1068 currentDetector = "TIB";
1069
1070 switch ((element.first) / 1000) {
1071 case 1:
1072 detector = "TIB";
1073 break;
1074 case 2:
1075 detector = "TOB";
1076 break;
1077 case 3:
1078 detector = "TEC";
1079 break;
1080 case 4:
1081 detector = "TID";
1082 break;
1083 }
1084
1085 h1->SetBinContent(iBin, mean);
1086 h1->GetXaxis()->SetBinLabel(iBin, SiStripPI::regionType(element.first).second);
1087 h1->GetXaxis()->LabelsOption("v");
1088
1089 if (detector != currentDetector) {
1090 boundaries.push_back(iBin);
1091 currentDetector = detector;
1092 }
1093 }
1094
1095 h1->SetMarkerStyle(20);
1096 h1->SetMarkerSize(1);
1097 h1->SetMaximum(h1->GetMaximum() * 1.1);
1098 h1->Draw("HIST");
1099 h1->Draw("Psame");
1100
1101 canvas.Update();
1102
1103 TLine l[boundaries.size()];
1104 unsigned int i = 0;
1105 for (const auto& line : boundaries) {
1106 l[i] = TLine(h1->GetBinLowEdge(line), canvas.GetUymin(), h1->GetBinLowEdge(line), canvas.GetUymax());
1107 l[i].SetLineWidth(1);
1108 l[i].SetLineStyle(9);
1109 l[i].SetLineColor(2);
1110 l[i].Draw("same");
1111 i++;
1112 }
1113
1114 TLegend legend = TLegend(0.52, 0.82, 0.95, 0.9);
1115 legend.SetHeader((std::get<1>(iov)).c_str(), "C");
1116 legend.AddEntry(h1.get(), ("IOV: " + std::to_string(std::get<0>(iov))).c_str(), "PL");
1117 legend.SetTextSize(0.025);
1118 legend.Draw("same");
1119
1120 std::string fileName(m_imageFileName);
1121 canvas.SaveAs(fileName.c_str());
1122
1123 return true;
1124 }
1125
1126 private:
1127 TrackerTopology m_trackerTopo;
1128 };
1129
1130 typedef SiStripPedestalsByRegion<SiStripPI::mean> SiStripPedestalsMeanByRegion;
1131 typedef SiStripPedestalsByRegion<SiStripPI::min> SiStripPedestalsMinByRegion;
1132 typedef SiStripPedestalsByRegion<SiStripPI::max> SiStripPedestalsMaxByRegion;
1133 typedef SiStripPedestalsByRegion<SiStripPI::rms> SiStripPedestalsRMSByRegion;
1134
1135 }
1136
1137 PAYLOAD_INSPECTOR_MODULE(SiStripPedestals) {
1138 PAYLOAD_INSPECTOR_CLASS(SiStripPedestalCompareByPartitionSingleTag);
1139 PAYLOAD_INSPECTOR_CLASS(SiStripPedestalDiffByPartitionSingleTag);
1140 PAYLOAD_INSPECTOR_CLASS(SiStripPedestalCompareByPartitionTwoTags);
1141 PAYLOAD_INSPECTOR_CLASS(SiStripPedestalDiffByPartitionTwoTags);
1142 PAYLOAD_INSPECTOR_CLASS(SiStripPedestalCorrelationByPartition);
1143 PAYLOAD_INSPECTOR_CLASS(SiStripPedestalsTest);
1144 PAYLOAD_INSPECTOR_CLASS(SiStripPedestalPerDetId);
1145 PAYLOAD_INSPECTOR_CLASS(SiStripPedestalsValue);
1146 PAYLOAD_INSPECTOR_CLASS(SiStripPedestalsValuePerDetId);
1147 PAYLOAD_INSPECTOR_CLASS(SiStripPedestalValuePerStrip);
1148 PAYLOAD_INSPECTOR_CLASS(SiStripPedestalValuePerAPV);
1149 PAYLOAD_INSPECTOR_CLASS(SiStripPedestalValuePerModule);
1150 PAYLOAD_INSPECTOR_CLASS(SiStripPedestalValueComparisonPerStripSingleTag);
1151 PAYLOAD_INSPECTOR_CLASS(SiStripPedestalValueComparisonPerAPVSingleTag);
1152 PAYLOAD_INSPECTOR_CLASS(SiStripPedestalValueComparisonPerModuleSingleTag);
1153 PAYLOAD_INSPECTOR_CLASS(SiStripPedestalValueComparisonPerStripTwoTags);
1154 PAYLOAD_INSPECTOR_CLASS(SiStripPedestalValueComparisonPerAPVTwoTags);
1155 PAYLOAD_INSPECTOR_CLASS(SiStripPedestalValueComparisonPerModuleTwoTags);
1156 PAYLOAD_INSPECTOR_CLASS(SiStripZeroPedestalsFraction_TrackerMap);
1157 PAYLOAD_INSPECTOR_CLASS(SiStripPedestalsMin_TrackerMap);
1158 PAYLOAD_INSPECTOR_CLASS(SiStripPedestalsMax_TrackerMap);
1159 PAYLOAD_INSPECTOR_CLASS(SiStripPedestalsMean_TrackerMap);
1160 PAYLOAD_INSPECTOR_CLASS(SiStripPedestalsRMS_TrackerMap);
1161 PAYLOAD_INSPECTOR_CLASS(SiStripPedestalsMeanByRegion);
1162 PAYLOAD_INSPECTOR_CLASS(SiStripPedestalsMinByRegion);
1163 PAYLOAD_INSPECTOR_CLASS(SiStripPedestalsMaxByRegion);
1164 PAYLOAD_INSPECTOR_CLASS(SiStripPedestalsRMSByRegion);
1165 }