Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-04-12 07:44:04

0001 /*!
0002   \file SiStripApvGains_PayloadInspector
0003   \Payload Inspector Plugin for SiStrip Gain
0004   \author M. Musich
0005   \version $Revision: 1.0 $
0006   \date $Date: 2017/07/02 17:59:56 $
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 // the data format of the condition to be inspected
0014 #include "CondFormats/SiStripObjects/interface/SiStripApvGain.h"
0015 #include "DataFormats/DetId/interface/DetId.h"
0016 #include "DataFormats/SiStripDetId/interface/StripSubdetector.h"
0017 #include "CondFormats/SiStripObjects/interface/SiStripDetSummary.h"
0018 
0019 // needed for the tracker map
0020 #include "CommonTools/TrackerMap/interface/TrackerMap.h"
0021 
0022 // auxilliary functions
0023 #include "CondCore/SiStripPlugins/interface/SiStripPayloadInspectorHelper.h"
0024 #include "CalibTracker/StandaloneTrackerTopology/interface/StandaloneTrackerTopology.h"
0025 #include "SiStripCondObjectRepresent.h"
0026 
0027 #include <memory>
0028 #include <sstream>
0029 #include <iostream>
0030 
0031 // include ROOT
0032 #include "TProfile.h"
0033 #include "TH2F.h"
0034 #include "THStack.h"
0035 #include "TLegend.h"
0036 #include "TCanvas.h"
0037 #include "TLine.h"
0038 #include "TStyle.h"
0039 #include "TLatex.h"
0040 #include "TPave.h"
0041 #include "TPaveStats.h"
0042 
0043 namespace {
0044 
0045   using namespace cond::payloadInspector;
0046 
0047   class SiStripApvGainContainer : public SiStripCondObjectRepresent::SiStripDataContainer<SiStripApvGain, float> {
0048   public:
0049     SiStripApvGainContainer(const std::shared_ptr<SiStripApvGain>& payload,
0050                             const SiStripPI::MetaData& metadata,
0051                             const std::string& tagName)
0052         : SiStripCondObjectRepresent::SiStripDataContainer<SiStripApvGain, float>(payload, metadata, tagName) {
0053       payloadType_ = "SiStripApvGain";
0054       setGranularity(SiStripCondObjectRepresent::PERAPV);
0055     }
0056 
0057     void storeAllValues() override {
0058       std::vector<uint32_t> detid;
0059       payload_->getDetIds(detid);
0060 
0061       for (const auto& d : detid) {
0062         SiStripApvGain::Range range = payload_->getRange(d);
0063         for (int it = 0; it < range.second - range.first; it++) {
0064           // to be used to fill the histogram
0065           SiStripCondData_.fillByPushBack(d, payload_->getApvGain(it, range));
0066         }
0067       }
0068     }
0069   };
0070 
0071   /************************************************
0072     testing the machinery
0073   ************************************************/
0074   class SiStripApvGainTest : public Histogram1D<SiStripApvGain, SINGLE_IOV> {
0075   public:
0076     SiStripApvGainTest()
0077         : Histogram1D<SiStripApvGain, SINGLE_IOV>("SiStrip ApvGain values", "SiStrip ApvGain values", 1, 0.0, 1.) {}
0078 
0079     bool fill() override {
0080       auto tag = PlotBase::getTag<0>();
0081       auto tagname = tag.name;
0082       for (auto const& iov : tag.iovs) {
0083         std::shared_ptr<SiStripApvGain> payload = Base::fetchPayload(std::get<1>(iov));
0084         if (payload.get()) {
0085           SiStripApvGainContainer* objContainer = new SiStripApvGainContainer(payload, iov, tagname);
0086           objContainer->printAll();
0087 
0088         }  // payload
0089       }    // iovs
0090       return true;
0091     }  // fill
0092   };
0093 
0094   class SiStripApvGainByPartition : public PlotImage<SiStripApvGain, SINGLE_IOV> {
0095   public:
0096     SiStripApvGainByPartition() : PlotImage<SiStripApvGain, SINGLE_IOV>("SiStrip ApvGains By Partition") {}
0097 
0098     bool fill() override {
0099       auto tag = PlotBase::getTag<0>();
0100       auto iov = tag.iovs.front();
0101       auto tagname = tag.name;
0102       std::shared_ptr<SiStripApvGain> payload = fetchPayload(std::get<1>(iov));
0103       if (payload.get()) {
0104         SiStripApvGainContainer* objContainer = new SiStripApvGainContainer(payload, iov, tagname);
0105         //objContainer->printAll();
0106 
0107         TCanvas canvas("Partition summary", "partition summary", 1400, 1000);
0108         objContainer->fillByPartition(canvas, 100, 0., 2.);
0109 
0110         std::string fileName(m_imageFileName);
0111         canvas.SaveAs(fileName.c_str());
0112       }  // payload
0113       return true;
0114     }  // fill
0115   };
0116 
0117   class SiStripApvGainCompareByPartition : public PlotImage<SiStripApvGain, MULTI_IOV, 2> {
0118   public:
0119     SiStripApvGainCompareByPartition()
0120         : PlotImage<SiStripApvGain, MULTI_IOV, 2>("SiStrip Compare ApvGains By Partition") {}
0121 
0122     bool fill() override {
0123       // trick to deal with the multi-ioved tag and two tag case at the same time
0124       auto theIOVs = PlotBase::getTag<0>().iovs;
0125       auto tagname1 = PlotBase::getTag<0>().name;
0126       auto tag2iovs = PlotBase::getTag<1>().iovs;
0127       auto tagname2 = PlotBase::getTag<1>().name;
0128       SiStripPI::MetaData firstiov = theIOVs.front();
0129       SiStripPI::MetaData lastiov = tag2iovs.front();
0130 
0131       std::shared_ptr<SiStripApvGain> last_payload = fetchPayload(std::get<1>(lastiov));
0132       std::shared_ptr<SiStripApvGain> first_payload = fetchPayload(std::get<1>(firstiov));
0133 
0134       SiStripApvGainContainer* l_objContainer = new SiStripApvGainContainer(last_payload, lastiov, tagname1);
0135       SiStripApvGainContainer* f_objContainer = new SiStripApvGainContainer(first_payload, firstiov, tagname2);
0136 
0137       l_objContainer->compare(f_objContainer);
0138 
0139       //l_objContainer->printAll();
0140 
0141       TCanvas canvas("Partition summary", "partition summary", 1400, 1000);
0142       l_objContainer->fillByPartition(canvas, 100, 0.5, 1.5);
0143 
0144       std::string fileName(m_imageFileName);
0145       canvas.SaveAs(fileName.c_str());
0146 
0147       return true;
0148     }  // fill
0149   };
0150 
0151   class SiStripApvGainRatioByPartition : public PlotImage<SiStripApvGain, MULTI_IOV, 2> {
0152   public:
0153     SiStripApvGainRatioByPartition() : PlotImage<SiStripApvGain, MULTI_IOV, 2>("SiStrip Ratio ApvGains By Partition") {}
0154 
0155     bool fill() override {
0156       // trick to deal with the multi-ioved tag and two tag case at the same time
0157       auto theIOVs = PlotBase::getTag<0>().iovs;
0158       auto tagname1 = PlotBase::getTag<0>().name;
0159       auto tag2iovs = PlotBase::getTag<1>().iovs;
0160       auto tagname2 = PlotBase::getTag<1>().name;
0161       SiStripPI::MetaData firstiov = theIOVs.front();
0162       SiStripPI::MetaData lastiov = tag2iovs.front();
0163 
0164       std::shared_ptr<SiStripApvGain> last_payload = fetchPayload(std::get<1>(lastiov));
0165       std::shared_ptr<SiStripApvGain> first_payload = fetchPayload(std::get<1>(firstiov));
0166 
0167       SiStripApvGainContainer* l_objContainer = new SiStripApvGainContainer(last_payload, lastiov, tagname1);
0168       SiStripApvGainContainer* f_objContainer = new SiStripApvGainContainer(first_payload, firstiov, tagname2);
0169 
0170       l_objContainer->divide(f_objContainer);
0171 
0172       //l_objContainer->printAll();
0173 
0174       TCanvas canvas("Partition summary", "partition summary", 1400, 1000);
0175       l_objContainer->fillByPartition(canvas, 200, 0.5, 1.5);
0176       //for (int i = 1; i <= 4; i++)
0177       //  canvas.cd(i)->SetLogy();
0178 
0179       std::string fileName(m_imageFileName);
0180       canvas.SaveAs(fileName.c_str());
0181 
0182       return true;
0183     }  // fill
0184   };
0185 
0186   class SiStripApvGainDiffByPartition : public PlotImage<SiStripApvGain, MULTI_IOV, 2> {
0187   public:
0188     SiStripApvGainDiffByPartition() : PlotImage<SiStripApvGain, MULTI_IOV, 2>("SiStrip Diff ApvGains By Partition") {}
0189 
0190     bool fill() override {
0191       // trick to deal with the multi-ioved tag and two tag case at the same time
0192       auto theIOVs = PlotBase::getTag<0>().iovs;
0193       auto tagname1 = PlotBase::getTag<0>().name;
0194       auto tag2iovs = PlotBase::getTag<1>().iovs;
0195       auto tagname2 = PlotBase::getTag<1>().name;
0196       SiStripPI::MetaData firstiov = theIOVs.front();
0197       SiStripPI::MetaData lastiov = tag2iovs.front();
0198 
0199       std::shared_ptr<SiStripApvGain> last_payload = fetchPayload(std::get<1>(lastiov));
0200       std::shared_ptr<SiStripApvGain> first_payload = fetchPayload(std::get<1>(firstiov));
0201 
0202       SiStripApvGainContainer* l_objContainer = new SiStripApvGainContainer(last_payload, lastiov, tagname1);
0203       SiStripApvGainContainer* f_objContainer = new SiStripApvGainContainer(first_payload, firstiov, tagname2);
0204 
0205       l_objContainer->subtract(f_objContainer);
0206 
0207       //l_objContainer->printAll();
0208 
0209       TCanvas canvas("Partition summary", "partition summary", 1400, 1000);
0210       l_objContainer->fillByPartition(canvas, 100, -0.1, 0.1);
0211 
0212       std::string fileName(m_imageFileName);
0213       canvas.SaveAs(fileName.c_str());
0214 
0215       return true;
0216     }  // fill
0217   };
0218 
0219   /************************************************
0220     1d histogram of SiStripApvGains of 1 IOV 
0221   *************************************************/
0222 
0223   // inherit from one of the predefined plot class: Histogram1D
0224   class SiStripApvGainsValue : public Histogram1D<SiStripApvGain, SINGLE_IOV> {
0225   public:
0226     SiStripApvGainsValue()
0227         : Histogram1D<SiStripApvGain, SINGLE_IOV>("SiStripApv Gains values", "SiStripApv Gains values", 200, 0.0, 2.0) {
0228     }
0229 
0230     bool fill() override {
0231       auto tag = PlotBase::getTag<0>();
0232       for (auto const& iov : tag.iovs) {
0233         std::shared_ptr<SiStripApvGain> payload = Base::fetchPayload(std::get<1>(iov));
0234         if (payload.get()) {
0235           std::vector<uint32_t> detid;
0236           payload->getDetIds(detid);
0237 
0238           for (const auto& d : detid) {
0239             SiStripApvGain::Range range = payload->getRange(d);
0240             for (int it = 0; it < range.second - range.first; it++) {
0241               // to be used to fill the histogram
0242               fillWithValue(payload->getApvGain(it, range));
0243 
0244             }  // loop over APVs
0245           }    // loop over detIds
0246         }      // payload
0247       }        // iovs
0248       return true;
0249     }  // fill
0250   };
0251 
0252   /************************************************
0253     1d histogram of means of SiStripApvGains
0254     for Tracker Barrel of 1 IOV 
0255   *************************************************/
0256 
0257   // inherit from one of the predefined plot class: Histogram1D
0258   class SiStripApvBarrelGainsByLayer : public Histogram1D<SiStripApvGain, SINGLE_IOV> {
0259   public:
0260     SiStripApvBarrelGainsByLayer()
0261         : Histogram1D<SiStripApvGain, SINGLE_IOV>("SiStripApv Gains averages by Barrel layer",
0262                                                   "Barrel layer (0-3: TIB), (4-9: TOB)",
0263                                                   10,
0264                                                   0,
0265                                                   10,
0266                                                   "average SiStripApv Gain") {}
0267 
0268     bool fill() override {
0269       auto tag = PlotBase::getTag<0>();
0270       for (auto const& iov : tag.iovs) {
0271         std::shared_ptr<SiStripApvGain> payload = Base::fetchPayload(std::get<1>(iov));
0272         if (payload.get()) {
0273           TrackerTopology tTopo = StandaloneTrackerTopology::fromTrackerParametersXMLFile(
0274               edm::FileInPath("Geometry/TrackerCommonData/data/trackerParameters.xml").fullPath());
0275 
0276           std::vector<uint32_t> detid;
0277           payload->getDetIds(detid);
0278 
0279           std::map<int, std::pair<float, float>> sumOfGainsByLayer;
0280 
0281           for (const auto& d : detid) {
0282             int subid = DetId(d).subdetId();
0283             int layer(-1);
0284             if (subid != StripSubdetector::TIB && subid != StripSubdetector::TOB)
0285               continue;
0286             if (subid == StripSubdetector::TIB) {
0287               layer = tTopo.tibLayer(d);
0288             } else if (subid == StripSubdetector::TOB) {
0289               // layers of TOB start at 5th bin
0290               layer = tTopo.tobLayer(d);
0291               layer += 4;
0292             }
0293 
0294             SiStripApvGain::Range range = payload->getRange(d);
0295             for (int it = 0; it < range.second - range.first; it++) {
0296               sumOfGainsByLayer[layer].first += payload->getApvGain(it, range);
0297               sumOfGainsByLayer[layer].second += 1.;
0298             }  // loop over APVs
0299           }    // loop over detIds
0300 
0301           // loop on the map to fill the plot
0302           for (auto& data : sumOfGainsByLayer) {
0303             fillWithBinAndValue(data.first - 1, (data.second.first / data.second.second));
0304           }
0305 
0306         }  // payload
0307       }    // iovs
0308       return true;
0309     }  // fill
0310   };
0311 
0312   /************************************************
0313     2d histogram of absolute (i.e. not average)
0314     SiStripApvGains for Tracker Barrel of 1 IOV
0315   *************************************************/
0316 
0317   class SiStripApvAbsoluteBarrelGainsByLayer : public Histogram2D<SiStripApvGain, SINGLE_IOV> {
0318   public:
0319     SiStripApvAbsoluteBarrelGainsByLayer()
0320         : Histogram2D<SiStripApvGain, SINGLE_IOV>("SiStripApv Gains by Barrel layer",
0321                                                   "Barrel layer (0-3: TIB), (4-9: TOB)",
0322                                                   10,
0323                                                   0,
0324                                                   10,
0325                                                   "SiStripApv Gain",
0326                                                   200,
0327                                                   0.0,
0328                                                   2.0) {}
0329     bool fill() override {
0330       auto tag = PlotBase::getTag<0>();
0331       for (auto const& iov : tag.iovs) {
0332         std::shared_ptr<SiStripApvGain> payload = Base::fetchPayload(std::get<1>(iov));
0333         if (payload.get()) {
0334           TrackerTopology tTopo = StandaloneTrackerTopology::fromTrackerParametersXMLFile(
0335               edm::FileInPath("Geometry/TrackerCommonData/data/trackerParameters.xml").fullPath());
0336 
0337           std::vector<uint32_t> detid;
0338           payload->getDetIds(detid);
0339           for (const auto& d : detid) {
0340             int subid = DetId(d).subdetId();
0341             if (subid != 3 && subid != 5)
0342               continue;
0343 
0344             SiStripApvGain::Range range = payload->getRange(d);
0345             for (int it = 0; it < range.second - range.first; it++) {
0346               float gain = payload->getApvGain(it, range);
0347               fillWithValue(static_cast<float>((subid == 5) ? tTopo.tobLayer(d) + 3 : tTopo.tibLayer(d) - 1),
0348                             (gain > 2.0) ? 2.0 : gain);
0349             }
0350           }  //loop over detIds
0351         }    // loop over payloads
0352       }      // loop over iovs
0353       return true;
0354     }  // fill
0355   };
0356 
0357   /************************************************
0358     1d histogram of means of SiStripApvGains
0359     for Tracker Endcaps (minus side) of 1 IOV 
0360   *************************************************/
0361 
0362   // inherit from one of the predefined plot class: Histogram1D
0363   class SiStripApvEndcapMinusGainsByDisk : public Histogram1D<SiStripApvGain, SINGLE_IOV> {
0364   public:
0365     SiStripApvEndcapMinusGainsByDisk()
0366         : Histogram1D<SiStripApvGain, SINGLE_IOV>("SiStripApv Gains averages by Endcap (minus) disk",
0367                                                   "Endcap (minus) disk (0-2: TID), (3-11: TEC)",
0368                                                   12,
0369                                                   0,
0370                                                   12,
0371                                                   "average SiStripApv Gain") {}
0372 
0373     bool fill() override {
0374       auto tag = PlotBase::getTag<0>();
0375       for (auto const& iov : tag.iovs) {
0376         std::shared_ptr<SiStripApvGain> payload = Base::fetchPayload(std::get<1>(iov));
0377         if (payload.get()) {
0378           TrackerTopology tTopo = StandaloneTrackerTopology::fromTrackerParametersXMLFile(
0379               edm::FileInPath("Geometry/TrackerCommonData/data/trackerParameters.xml").fullPath());
0380 
0381           std::vector<uint32_t> detid;
0382           payload->getDetIds(detid);
0383 
0384           std::map<int, std::pair<float, float>> sumOfGainsByDisk;
0385 
0386           for (const auto& d : detid) {
0387             int disk = -1;
0388             int side = -1;
0389             int subid = DetId(d).subdetId();
0390             if (subid != StripSubdetector::TID && subid != StripSubdetector::TEC)
0391               continue;
0392 
0393             if (subid == StripSubdetector::TID) {
0394               side = tTopo.tidSide(d);
0395               disk = tTopo.tidWheel(d);
0396             } else {
0397               side = tTopo.tecSide(d);
0398               disk = tTopo.tecWheel(d);
0399               // disks of TEC start at 4th bin
0400               disk += 3;
0401             }
0402 
0403             // only negative side
0404             if (side != 1)
0405               continue;
0406 
0407             SiStripApvGain::Range range = payload->getRange(d);
0408             for (int it = 0; it < range.second - range.first; it++) {
0409               sumOfGainsByDisk[disk].first += payload->getApvGain(it, range);
0410               sumOfGainsByDisk[disk].second += 1.;
0411             }  // loop over APVs
0412           }    // loop over detIds
0413 
0414           // loop on the map to fill the plot
0415           for (auto& data : sumOfGainsByDisk) {
0416             fillWithBinAndValue(data.first - 1, (data.second.first / data.second.second));
0417           }
0418 
0419         }  // payload
0420       }    // iovs
0421       return true;
0422     }  // fill
0423   };
0424 
0425   /************************************************
0426     1d histogram of means of SiStripApvGains
0427     for Tracker Endcaps (plus side) of 1 IOV 
0428   *************************************************/
0429 
0430   // inherit from one of the predefined plot class: Histogram1D
0431   class SiStripApvEndcapPlusGainsByDisk : public Histogram1D<SiStripApvGain, SINGLE_IOV> {
0432   public:
0433     SiStripApvEndcapPlusGainsByDisk()
0434         : Histogram1D<SiStripApvGain, SINGLE_IOV>("SiStripApv Gains averages by Endcap (plus) disk",
0435                                                   "Endcap (plus) disk (0-2: TID), (3-11: TEC)",
0436                                                   12,
0437                                                   0,
0438                                                   12,
0439                                                   "average SiStripApv Gain") {}
0440 
0441     bool fill() override {
0442       auto tag = PlotBase::getTag<0>();
0443       for (auto const& iov : tag.iovs) {
0444         std::shared_ptr<SiStripApvGain> payload = Base::fetchPayload(std::get<1>(iov));
0445         if (payload.get()) {
0446           TrackerTopology tTopo = StandaloneTrackerTopology::fromTrackerParametersXMLFile(
0447               edm::FileInPath("Geometry/TrackerCommonData/data/trackerParameters.xml").fullPath());
0448 
0449           std::vector<uint32_t> detid;
0450           payload->getDetIds(detid);
0451 
0452           std::map<int, std::pair<float, float>> sumOfGainsByDisk;
0453 
0454           for (const auto& d : detid) {
0455             int disk = -1;
0456             int side = -1;
0457             int subid = DetId(d).subdetId();
0458             if (subid != StripSubdetector::TID && subid != StripSubdetector::TEC)
0459               continue;
0460 
0461             if (subid == StripSubdetector::TID) {
0462               side = tTopo.tidSide(d);
0463               disk = tTopo.tidWheel(d);
0464               ;
0465             } else {
0466               side = tTopo.tecSide(d);
0467               disk = tTopo.tecWheel(d);
0468               // disks of TEC start at 4th bin
0469               disk += 3;
0470             }
0471 
0472             // only positive side
0473             if (side != 2)
0474               continue;
0475 
0476             SiStripApvGain::Range range = payload->getRange(d);
0477             for (int it = 0; it < range.second - range.first; it++) {
0478               sumOfGainsByDisk[disk].first += payload->getApvGain(it, range);
0479               sumOfGainsByDisk[disk].second += 1.;
0480             }  // loop over APVs
0481           }    // loop over detIds
0482 
0483           // loop on the map to fill the plot
0484           for (auto& data : sumOfGainsByDisk) {
0485             fillWithBinAndValue(data.first - 1, (data.second.first / data.second.second));
0486           }
0487 
0488         }  // payload
0489       }    // iovs
0490       return true;
0491     }  // fill
0492   };
0493 
0494   /************************************************
0495     2D histogram of absolute (i.e. not average)
0496     SiStripApv Gains on the Endcap- for 1 IOV
0497    ************************************************/
0498   class SiStripApvAbsoluteEndcapMinusGainsByDisk : public Histogram2D<SiStripApvGain, SINGLE_IOV> {
0499   public:
0500     SiStripApvAbsoluteEndcapMinusGainsByDisk()
0501         : Histogram2D<SiStripApvGain, SINGLE_IOV>("SiStripApv Gains averages by Endcap (minus) disk",
0502                                                   "Endcap (minus) disk (0-2: TID), (3-11: TEC)",
0503                                                   12,
0504                                                   0,
0505                                                   12,
0506                                                   "SiStripApv Gain",
0507                                                   200,
0508                                                   0.0,
0509                                                   2.0) {}
0510 
0511     bool fill() override {
0512       auto tag = PlotBase::getTag<0>();
0513       for (auto const& iov : tag.iovs) {
0514         std::shared_ptr<SiStripApvGain> payload = Base::fetchPayload(std::get<1>(iov));
0515         if (payload.get()) {
0516           TrackerTopology tTopo = StandaloneTrackerTopology::fromTrackerParametersXMLFile(
0517               edm::FileInPath("Geometry/TrackerCommonData/data/trackerParameters.xml").fullPath());
0518 
0519           std::vector<uint32_t> detid;
0520           payload->getDetIds(detid);
0521 
0522           for (const auto& d : detid) {
0523             int subid = DetId(d).subdetId(), side = -1, disk = -1;
0524 
0525             switch (subid) {
0526               case 4:
0527                 side = tTopo.tidSide(d);
0528                 disk = tTopo.tidWheel(d);
0529                 break;
0530               case 6:
0531                 side = tTopo.tecSide(d);
0532                 disk = tTopo.tecWheel(d) + 4;
0533                 break;
0534               default:
0535                 continue;
0536             }
0537 
0538             if (side != 1)
0539               continue;
0540             SiStripApvGain::Range range = payload->getRange(d);
0541             for (int it = 0; it < range.second - range.first; it++) {
0542               float gain = payload->getApvGain(it, range);
0543               fillWithValue((float)disk - 1, (gain > 2.0) ? 2.0 : gain);
0544             }  // apvs
0545           }    // detids
0546         }
0547       }  // iovs
0548       return true;
0549     }  // fill
0550   };
0551 
0552   /************************************************
0553     2D histogram of absolute (i.e. not average)
0554     SiStripApv Gains on the Endcap+ for 1 IOV
0555    ************************************************/
0556   class SiStripApvAbsoluteEndcapPlusGainsByDisk : public Histogram2D<SiStripApvGain, SINGLE_IOV> {
0557   public:
0558     SiStripApvAbsoluteEndcapPlusGainsByDisk()
0559         : Histogram2D<SiStripApvGain, SINGLE_IOV>("SiStripApv Gains averages by Endcap (plus) disk",
0560                                                   "Endcap (plus) disk (0-2: TID), (3-11: TEC)",
0561                                                   12,
0562                                                   0,
0563                                                   12,
0564                                                   "SiStripApv Gain",
0565                                                   200,
0566                                                   0.0,
0567                                                   2.0) {}
0568     bool fill() override {
0569       auto tag = PlotBase::getTag<0>();
0570       for (auto const& iov : tag.iovs) {
0571         std::shared_ptr<SiStripApvGain> payload = Base::fetchPayload(std::get<1>(iov));
0572         if (payload.get()) {
0573           TrackerTopology tTopo = StandaloneTrackerTopology::fromTrackerParametersXMLFile(
0574               edm::FileInPath("Geometry/TrackerCommonData/data/trackerParameters.xml").fullPath());
0575 
0576           std::vector<uint32_t> detid;
0577           payload->getDetIds(detid);
0578 
0579           for (const auto& d : detid) {
0580             int subid = DetId(d).subdetId(), side = -1, disk = -1;
0581 
0582             switch (subid) {
0583               case 4:
0584                 side = tTopo.tidSide(d);
0585                 disk = tTopo.tidWheel(d);
0586                 break;
0587               case 6:
0588                 side = tTopo.tecSide(d);
0589                 disk = tTopo.tecWheel(d) + 4;
0590                 break;
0591               default:
0592                 continue;
0593             }
0594 
0595             if (side != 2)
0596               continue;
0597             SiStripApvGain::Range range = payload->getRange(d);
0598             for (int it = 0; it < range.second - range.first; it++) {
0599               float gain = payload->getApvGain(it, range);
0600               fillWithValue((float)disk - 1, (gain > 2.0) ? 2.0 : gain);
0601             }  //apvs
0602           }    //detids
0603         }
0604       }  // iovs
0605       return true;
0606     }  // fill
0607   };
0608 
0609   /************************************************
0610     TrackerMap of SiStripApvGains (average gain per detid)
0611   *************************************************/
0612   class SiStripApvGainsAverageTrackerMap : public PlotImage<SiStripApvGain, SINGLE_IOV> {
0613   public:
0614     SiStripApvGainsAverageTrackerMap() : PlotImage<SiStripApvGain, SINGLE_IOV>("Tracker Map of average SiStripGains") {}
0615 
0616     bool fill() override {
0617       auto tag = PlotBase::getTag<0>();
0618       auto iov = tag.iovs.front();
0619       std::shared_ptr<SiStripApvGain> payload = fetchPayload(std::get<1>(iov));
0620 
0621       std::string titleMap = "SiStrip APV Gain average per module (payload : " + std::get<1>(iov) + ")";
0622 
0623       std::unique_ptr<TrackerMap> tmap = std::make_unique<TrackerMap>("SiStripApvGains");
0624       tmap->setTitle(titleMap);
0625       tmap->setPalette(1);
0626 
0627       std::vector<uint32_t> detid;
0628       payload->getDetIds(detid);
0629 
0630       std::map<uint32_t, float> store;
0631 
0632       for (const auto& d : detid) {
0633         SiStripApvGain::Range range = payload->getRange(d);
0634         float sumOfGains = 0;
0635         float nAPVsPerModule = 0.;
0636         for (int it = 0; it < range.second - range.first; it++) {
0637           nAPVsPerModule += 1;
0638           sumOfGains += payload->getApvGain(it, range);
0639         }  // loop over APVs
0640         // fill the tracker map taking the average gain on a single DetId
0641         store[d] = (sumOfGains / nAPVsPerModule);
0642         tmap->fill(d, (sumOfGains / nAPVsPerModule));
0643       }  // loop over detIds
0644 
0645       //=========================
0646       // saturate at 2 std deviations
0647       auto range = SiStripPI::getTheRange(store, 2);
0648 
0649       std::string fileName(m_imageFileName);
0650       tmap->save(true, range.first, range.second, fileName);
0651 
0652       return true;
0653     }
0654   };
0655 
0656   /************************************************
0657     TrackerMap of SiStripApvGains (module with default)
0658   *************************************************/
0659   class SiStripApvGainsDefaultTrackerMap : public PlotImage<SiStripApvGain, SINGLE_IOV> {
0660   public:
0661     SiStripApvGainsDefaultTrackerMap()
0662         : PlotImage<SiStripApvGain, SINGLE_IOV>("Tracker Map of SiStripGains to default") {}
0663 
0664     bool fill() override {
0665       auto tag = PlotBase::getTag<0>();
0666       auto iov = tag.iovs.front();
0667       std::shared_ptr<SiStripApvGain> payload = fetchPayload(std::get<1>(iov));
0668 
0669       std::unique_ptr<TrackerMap> tmap = std::make_unique<TrackerMap>("SiStripApvGains");
0670 
0671       tmap->setPalette(1);
0672 
0673       std::vector<uint32_t> detid;
0674       payload->getDetIds(detid);
0675 
0676       /*
0677     the defaul G1 value comes from the ratio of DefaultTickHeight/GainNormalizationFactor
0678     as defined in the default of the O2O producer: OnlineDB/SiStripESSources/src/SiStripCondObjBuilderFromDb.cc
0679       */
0680 
0681       float G1default = 690. / 640.;
0682       float G2default = 1.;
0683 
0684       int totalG1DefaultAPVs = 0;
0685       int totalG2DefaultAPVs = 0;
0686 
0687       for (const auto& d : detid) {
0688         SiStripApvGain::Range range = payload->getRange(d);
0689         float sumOfGains = 0;
0690         float nAPVsPerModule = 0.;
0691         int countDefaults = 0;
0692         for (int it = 0; it < range.second - range.first; it++) {
0693           nAPVsPerModule += 1;
0694           sumOfGains += payload->getApvGain(it, range);
0695           if ((payload->getApvGain(it, range)) == G1default || (payload->getApvGain(it, range)) == G2default)
0696             countDefaults++;
0697         }  // loop over APVs
0698         // fill the tracker map taking the average gain on a single DetId
0699         if (countDefaults > 0.) {
0700           tmap->fill(d, countDefaults);
0701 
0702           if (std::fmod((sumOfGains / countDefaults), G1default) == 0.) {
0703             totalG1DefaultAPVs += countDefaults;
0704           } else if (std::fmod((sumOfGains / countDefaults), G2default) == 0.) {
0705             totalG2DefaultAPVs += countDefaults;
0706           }
0707         }
0708       }  // loop over detIds
0709 
0710       //=========================
0711 
0712       std::string gainType = totalG1DefaultAPVs == 0 ? "G2 value (=1)" : "G1 value (=690./640.)";
0713 
0714       std::string titleMap = "# of APVs/module w/ default " + gainType + " (payload : " + std::get<1>(iov) + ")";
0715       tmap->setTitle(titleMap);
0716 
0717       std::string fileName(m_imageFileName);
0718       tmap->save(true, 0, 0, fileName);
0719 
0720       return true;
0721     }
0722   };
0723 
0724   /************************************************
0725     TrackerMap of SiStripApvGains (ratio with previous gain per detid)
0726   *************************************************/
0727 
0728   template <int ntags, IOVMultiplicity nIOVs>
0729   class SiStripApvGainsRatioTrackerMapBase : public PlotImage<SiStripApvGain, nIOVs, ntags> {
0730   public:
0731     SiStripApvGainsRatioTrackerMapBase()
0732         : PlotImage<SiStripApvGain, nIOVs, ntags>("Tracker Map of ratio of SiStripGains with previous IOV") {
0733       PlotBase::addInputParam("nsigma");
0734     }
0735 
0736     bool fill() override {
0737       // determine n. sigmas
0738       unsigned int nsigma(1);
0739 
0740       auto paramValues = PlotBase::inputParamValues();
0741       auto ip = paramValues.find("nsigma");
0742       if (ip != paramValues.end()) {
0743         nsigma = std::stoul(ip->second);
0744       }
0745 
0746       // trick to deal with the multi-ioved tag and two tag case at the same time
0747       auto theIOVs = PlotBase::getTag<0>().iovs;
0748       auto tagname1 = PlotBase::getTag<0>().name;
0749       std::string tagname2 = "";
0750       auto firstiov = theIOVs.front();
0751       SiStripPI::MetaData lastiov;
0752 
0753       // we don't support (yet) comparison with more than 2 tags
0754       assert(this->m_plotAnnotations.ntags < 3);
0755 
0756       if (this->m_plotAnnotations.ntags == 2) {
0757         auto tag2iovs = PlotBase::getTag<1>().iovs;
0758         tagname2 = PlotBase::getTag<1>().name;
0759         lastiov = tag2iovs.front();
0760       } else {
0761         lastiov = theIOVs.back();
0762       }
0763 
0764       std::shared_ptr<SiStripApvGain> last_payload = this->fetchPayload(std::get<1>(lastiov));
0765       std::shared_ptr<SiStripApvGain> first_payload = this->fetchPayload(std::get<1>(firstiov));
0766 
0767       std::string titleMap = "SiStrip APV Gain ratio per module average (IOV: ";
0768 
0769       titleMap += std::to_string(std::get<0>(firstiov));
0770       titleMap += "/ IOV:";
0771       titleMap += std::to_string(std::get<0>(lastiov));
0772       titleMap += ")";
0773 
0774       titleMap += +" " + std::to_string(nsigma) + " std. dev. saturation";
0775 
0776       std::unique_ptr<TrackerMap> tmap = std::make_unique<TrackerMap>("SiStripApvGains");
0777       tmap->setTitle(titleMap);
0778       tmap->setPalette(1);
0779 
0780       std::map<uint32_t, float> lastmap, firstmap;
0781 
0782       std::vector<uint32_t> detid;
0783       last_payload->getDetIds(detid);
0784 
0785       // cache the last IOV
0786       for (const auto& d : detid) {
0787         SiStripApvGain::Range range = last_payload->getRange(d);
0788         float Gain = 0;
0789         float nAPV = 0;
0790         for (int it = 0; it < range.second - range.first; it++) {
0791           nAPV += 1;
0792           Gain += last_payload->getApvGain(it, range);
0793         }  // loop over APVs
0794         lastmap[d] = (Gain / nAPV);
0795       }  // loop over detIds
0796 
0797       detid.clear();
0798 
0799       first_payload->getDetIds(detid);
0800 
0801       // cache the first IOV
0802       for (const auto& d : detid) {
0803         SiStripApvGain::Range range = first_payload->getRange(d);
0804         float Gain = 0;
0805         float nAPV = 0;
0806         for (int it = 0; it < range.second - range.first; it++) {
0807           nAPV += 1;
0808           Gain += first_payload->getApvGain(it, range);
0809         }  // loop over APVs
0810         firstmap[d] = (Gain / nAPV);
0811       }  // loop over detIds
0812 
0813       std::map<uint32_t, float> cachedRatio;
0814       for (const auto& d : detid) {
0815         float ratio = firstmap[d] / lastmap[d];
0816         tmap->fill(d, ratio);
0817         cachedRatio[d] = ratio;
0818       }
0819 
0820       //=========================
0821       auto range = SiStripPI::getTheRange(cachedRatio, nsigma);
0822 
0823       std::string fileName(this->m_imageFileName);
0824       tmap->save(true, range.first, range.second, fileName);
0825 
0826       return true;
0827     }
0828   };
0829 
0830   using SiStripApvGainsAvgDeviationRatioWithPreviousIOVTrackerMap = SiStripApvGainsRatioTrackerMapBase<1, MULTI_IOV>;
0831   using SiStripApvGainsAvgDeviationRatioTrackerMapTwoTags = SiStripApvGainsRatioTrackerMapBase<2, SINGLE_IOV>;
0832 
0833   /************************************************
0834    TrackerMap of SiStripApvGains (ratio for largest deviation with previous gain per detid)
0835   *************************************************/
0836 
0837   template <int ntags, IOVMultiplicity nIOVs>
0838   class SiStripApvGainsRatioMaxDeviationTrackerMapBase : public PlotImage<SiStripApvGain, nIOVs, ntags> {
0839   public:
0840     SiStripApvGainsRatioMaxDeviationTrackerMapBase()
0841         : PlotImage<SiStripApvGain, nIOVs, ntags>(
0842               "Tracker Map of ratio (for largest deviation) of SiStripGains with previous IOV") {
0843       PlotBase::addInputParam("nsigma");
0844     }
0845 
0846     bool fill() override {
0847       unsigned int nsigma(1);
0848       auto paramValues = PlotBase::inputParamValues();
0849       auto ip = paramValues.find("nsigma");
0850       if (ip != paramValues.end()) {
0851         nsigma = std::stoul(ip->second);
0852         edm::LogPrint("SiStripApvGain_PayloadInspector")
0853             << "using custom z-axis saturation: " << nsigma << " sigmas" << std::endl;
0854       } else {
0855         edm::LogPrint("SiStripApvGain_PayloadInspector")
0856             << "using default saturation: " << nsigma << " sigmas" << std::endl;
0857       }
0858 
0859       // trick to deal with the multi-ioved tag and two tag case at the same time
0860       auto theIOVs = PlotBase::getTag<0>().iovs;
0861       auto tagname1 = PlotBase::getTag<0>().name;
0862       std::string tagname2 = "";
0863       auto firstiov = theIOVs.front();
0864       SiStripPI::MetaData lastiov;
0865 
0866       // we don't support (yet) comparison with more than 2 tags
0867       assert(this->m_plotAnnotations.ntags < 3);
0868 
0869       if (this->m_plotAnnotations.ntags == 2) {
0870         auto tag2iovs = PlotBase::getTag<1>().iovs;
0871         tagname2 = PlotBase::getTag<1>().name;
0872         lastiov = tag2iovs.front();
0873       } else {
0874         lastiov = theIOVs.back();
0875       }
0876 
0877       std::shared_ptr<SiStripApvGain> last_payload = this->fetchPayload(std::get<1>(lastiov));
0878       std::shared_ptr<SiStripApvGain> first_payload = this->fetchPayload(std::get<1>(firstiov));
0879 
0880       std::string titleMap = "SiStrip APV Gain ratio for largest deviation per module (IOV: ";
0881 
0882       titleMap += std::to_string(std::get<0>(firstiov));
0883       titleMap += "/ IOV:";
0884       titleMap += std::to_string(std::get<0>(lastiov));
0885       titleMap += ") ";
0886 
0887       titleMap += +" - " + std::to_string(nsigma) + " std. dev. saturation";
0888 
0889       std::unique_ptr<TrackerMap> tmap = std::make_unique<TrackerMap>("SiStripApvGains");
0890       tmap->setTitle(titleMap);
0891       tmap->setPalette(1);
0892 
0893       std::map<std::pair<uint32_t, int>, float> lastmap, firstmap;
0894 
0895       std::vector<uint32_t> detid;
0896       last_payload->getDetIds(detid);
0897 
0898       // cache the last IOV
0899       for (const auto& d : detid) {
0900         SiStripApvGain::Range range = last_payload->getRange(d);
0901         float nAPV = 0;
0902         for (int it = 0; it < range.second - range.first; it++) {
0903           nAPV += 1;
0904           float Gain = last_payload->getApvGain(it, range);
0905           std::pair<uint32_t, int> index = std::make_pair(d, nAPV);
0906           lastmap[index] = Gain;
0907         }  // loop over APVs
0908       }    // loop over detIds
0909 
0910       detid.clear();
0911 
0912       first_payload->getDetIds(detid);
0913 
0914       // cache the first IOV
0915       for (const auto& d : detid) {
0916         SiStripApvGain::Range range = first_payload->getRange(d);
0917         float nAPV = 0;
0918         for (int it = 0; it < range.second - range.first; it++) {
0919           nAPV += 1;
0920           float Gain = first_payload->getApvGain(it, range);
0921           std::pair<uint32_t, int> index = std::make_pair(d, nAPV);
0922           firstmap[index] = Gain;
0923         }  // loop over APVs
0924       }    // loop over detIds
0925 
0926       // find the largest deviation
0927       std::map<uint32_t, float> cachedRatio;
0928 
0929       for (const auto& item : firstmap) {
0930         // packed index (detid,APV)
0931         auto index = item.first;
0932         auto mod = item.first.first;
0933 
0934         float ratio = firstmap[index] / lastmap[index];
0935         // if we have already cached something
0936         if (cachedRatio[mod]) {
0937           // if the discrepancy with 1 of ratio is larger than the cached value
0938           if (std::abs(ratio - 1.) > std::abs(cachedRatio[mod] - 1.)) {
0939             cachedRatio[mod] = ratio;
0940           }
0941         } else {
0942           cachedRatio[mod] = ratio;
0943         }
0944       }
0945 
0946       for (const auto& element : cachedRatio) {
0947         tmap->fill(element.first, element.second);
0948       }
0949 
0950       // get the range of the TrackerMap (saturate at +/-n std deviations)
0951       auto range = SiStripPI::getTheRange(cachedRatio, nsigma);
0952 
0953       //=========================
0954 
0955       std::string fileName(this->m_imageFileName);
0956       tmap->save(true, range.first, range.second, fileName);
0957 
0958       return true;
0959     }
0960   };
0961 
0962   using SiStripApvGainsMaxDeviationRatioWithPreviousIOVTrackerMap =
0963       SiStripApvGainsRatioMaxDeviationTrackerMapBase<1, MULTI_IOV>;
0964 
0965   using SiStripApvGainsMaxDeviationRatioTrackerMapTwoTags =
0966       SiStripApvGainsRatioMaxDeviationTrackerMapBase<2, SINGLE_IOV>;
0967 
0968   /************************************************
0969     TrackerMap of SiStripApvGains (maximum gain per detid)
0970   *************************************************/
0971   class SiStripApvGainsMaximumTrackerMap : public PlotImage<SiStripApvGain, SINGLE_IOV> {
0972   public:
0973     SiStripApvGainsMaximumTrackerMap()
0974         : PlotImage<SiStripApvGain, SINGLE_IOV>("Tracker Map of SiStripAPVGains (maximum per DetId)") {}
0975 
0976     bool fill() override {
0977       auto tag = PlotBase::getTag<0>();
0978       auto iov = tag.iovs.front();
0979       std::shared_ptr<SiStripApvGain> payload = fetchPayload(std::get<1>(iov));
0980 
0981       std::string titleMap = "SiStrip APV Gain maximum per module (payload : " + std::get<1>(iov) + ")";
0982 
0983       std::unique_ptr<TrackerMap> tmap = std::make_unique<TrackerMap>("SiStripApvGains");
0984       tmap->setTitle(titleMap);
0985       tmap->setPalette(1);
0986 
0987       std::vector<uint32_t> detid;
0988       payload->getDetIds(detid);
0989 
0990       for (const auto& d : detid) {
0991         SiStripApvGain::Range range = payload->getRange(d);
0992         float theMaxGain = 0;
0993         for (int it = 0; it < range.second - range.first; it++) {
0994           float currentGain = payload->getApvGain(it, range);
0995           if (currentGain > theMaxGain) {
0996             theMaxGain = currentGain;
0997           }
0998         }  // loop over APVs
0999         // fill the tracker map taking the maximum gain on a single DetId
1000         tmap->fill(d, theMaxGain);
1001       }  // loop over detIds
1002 
1003       //=========================
1004 
1005       std::pair<float, float> extrema = tmap->getAutomaticRange();
1006 
1007       std::string fileName(m_imageFileName);
1008 
1009       // protect against uniform values (gains are defined positive)
1010       if (extrema.first != extrema.second) {
1011         tmap->save(true, 0, 0, fileName);
1012       } else {
1013         tmap->save(true, extrema.first * 0.95, extrema.first * 1.05, fileName);
1014       }
1015 
1016       return true;
1017     }
1018   };
1019 
1020   /************************************************
1021     TrackerMap of SiStripApvGains (minimum gain per detid)
1022   *************************************************/
1023   class SiStripApvGainsMinimumTrackerMap : public PlotImage<SiStripApvGain, SINGLE_IOV> {
1024   public:
1025     SiStripApvGainsMinimumTrackerMap()
1026         : PlotImage<SiStripApvGain, SINGLE_IOV>("Tracker Map of SiStripAPVGains (minimum per DetId)") {}
1027 
1028     bool fill() override {
1029       auto tag = PlotBase::getTag<0>();
1030       auto iov = tag.iovs.front();
1031 
1032       std::shared_ptr<SiStripApvGain> payload = fetchPayload(std::get<1>(iov));
1033 
1034       std::string titleMap = "SiStrip APV Gain minumum per module (payload : " + std::get<1>(iov) + ")";
1035 
1036       std::unique_ptr<TrackerMap> tmap = std::make_unique<TrackerMap>("SiStripApvGains");
1037       tmap->setTitle(titleMap);
1038       tmap->setPalette(1);
1039 
1040       std::vector<uint32_t> detid;
1041       payload->getDetIds(detid);
1042 
1043       for (const auto& d : detid) {
1044         SiStripApvGain::Range range = payload->getRange(d);
1045         float theMinGain = 999.;
1046         for (int it = 0; it < range.second - range.first; it++) {
1047           float currentGain = payload->getApvGain(it, range);
1048           if (currentGain < theMinGain) {
1049             theMinGain = currentGain;
1050           }
1051         }  // loop over APVs
1052         // fill the tracker map taking the minimum gain on a single DetId
1053         tmap->fill(d, theMinGain);
1054       }  // loop over detIds
1055 
1056       //=========================
1057 
1058       std::pair<float, float> extrema = tmap->getAutomaticRange();
1059 
1060       std::string fileName(m_imageFileName);
1061 
1062       // protect against uniform values (gains are defined positive)
1063       if (extrema.first != extrema.second) {
1064         tmap->save(true, 0, 0, fileName);
1065       } else {
1066         tmap->save(true, extrema.first * 0.95, extrema.first * 1.05, fileName);
1067       }
1068 
1069       return true;
1070     }
1071   };
1072 
1073   /************************************************
1074     time history histogram of SiStripApvGains 
1075   *************************************************/
1076 
1077   class SiStripApvGainByRunMeans : public HistoryPlot<SiStripApvGain, float> {
1078   public:
1079     SiStripApvGainByRunMeans()
1080         : HistoryPlot<SiStripApvGain, float>("SiStripApv Gains average", "average Strip APV gain value") {}
1081     ~SiStripApvGainByRunMeans() override = default;
1082 
1083     float getFromPayload(SiStripApvGain& payload) override {
1084       std::vector<uint32_t> detid;
1085       payload.getDetIds(detid);
1086 
1087       float nAPVs = 0;
1088       float sumOfGains = 0;
1089 
1090       for (const auto& d : detid) {
1091         SiStripApvGain::Range range = payload.getRange(d);
1092         for (int it = 0; it < range.second - range.first; it++) {
1093           nAPVs += 1;
1094           sumOfGains += payload.getApvGain(it, range);
1095         }  // loop over APVs
1096       }    // loop over detIds
1097 
1098       return sumOfGains / nAPVs;
1099     }  // payload
1100   };
1101 
1102   /************************************************
1103     time history of SiStripApvGains properties
1104   *************************************************/
1105 
1106   template <SiStripPI::estimator est>
1107   class SiStripApvGainProperties : public HistoryPlot<SiStripApvGain, float> {
1108   public:
1109     SiStripApvGainProperties()
1110         : HistoryPlot<SiStripApvGain, float>("SiStripApv Gains " + estimatorType(est),
1111                                              estimatorType(est) + " Strip APV gain value") {}
1112     ~SiStripApvGainProperties() override = default;
1113 
1114     float getFromPayload(SiStripApvGain& payload) override {
1115       std::vector<uint32_t> detid;
1116       payload.getDetIds(detid);
1117 
1118       float nAPVs = 0;
1119       float sumOfGains = 0;
1120       float meanOfGains = 0;
1121       float rmsOfGains = 0;
1122       float min(0.), max(0.);
1123 
1124       for (const auto& d : detid) {
1125         SiStripApvGain::Range range = payload.getRange(d);
1126         for (int it = 0; it < range.second - range.first; it++) {
1127           nAPVs += 1;
1128           float gain = payload.getApvGain(it, range);
1129           if (gain < min)
1130             min = gain;
1131           if (gain > max)
1132             max = gain;
1133           sumOfGains += gain;
1134           rmsOfGains += (gain * gain);
1135         }  // loop over APVs
1136       }    // loop over detIds
1137 
1138       meanOfGains = sumOfGains / nAPVs;
1139 
1140       switch (est) {
1141         case SiStripPI::min:
1142           return min;
1143           break;
1144         case SiStripPI::max:
1145           return max;
1146           break;
1147         case SiStripPI::mean:
1148           return meanOfGains;
1149           break;
1150         case SiStripPI::rms:
1151           if ((rmsOfGains / nAPVs - meanOfGains * meanOfGains) > 0.) {
1152             return sqrt(rmsOfGains / nAPVs - meanOfGains * meanOfGains);
1153           } else {
1154             return 0.;
1155           }
1156           break;
1157         default:
1158           edm::LogWarning("LogicError") << "Unknown estimator: " << est;
1159           break;
1160       }
1161       return 0.;
1162     }  // payload
1163   };
1164 
1165   typedef SiStripApvGainProperties<SiStripPI::min> SiStripApvGainMin_History;
1166   typedef SiStripApvGainProperties<SiStripPI::max> SiStripApvGainMax_History;
1167   typedef SiStripApvGainProperties<SiStripPI::mean> SiStripApvGainMean_History;
1168   typedef SiStripApvGainProperties<SiStripPI::rms> SiStripApvGainRMS_History;
1169 
1170   /************************************************
1171     time history histogram of TIB SiStripApvGains 
1172   *************************************************/
1173 
1174   class SiStripApvTIBGainByRunMeans : public HistoryPlot<SiStripApvGain, float> {
1175   public:
1176     SiStripApvTIBGainByRunMeans()
1177         : HistoryPlot<SiStripApvGain, float>("SiStripApv Gains average",
1178                                              "average Tracker Inner Barrel APV gain value") {}
1179     ~SiStripApvTIBGainByRunMeans() override = default;
1180 
1181     float getFromPayload(SiStripApvGain& payload) override {
1182       std::vector<uint32_t> detid;
1183       payload.getDetIds(detid);
1184 
1185       float nAPVs = 0;
1186       float sumOfGains = 0;
1187 
1188       for (const auto& d : detid) {
1189         int subid = DetId(d).subdetId();
1190         if (subid != StripSubdetector::TIB)
1191           continue;
1192 
1193         SiStripApvGain::Range range = payload.getRange(d);
1194         for (int it = 0; it < range.second - range.first; it++) {
1195           nAPVs += 1;
1196           sumOfGains += payload.getApvGain(it, range);
1197         }  // loop over APVs
1198       }    // loop over detIds
1199 
1200       return sumOfGains / nAPVs;
1201 
1202     }  // payload
1203   };
1204 
1205   /************************************************
1206     time history histogram of TOB SiStripApvGains 
1207   *************************************************/
1208 
1209   class SiStripApvTOBGainByRunMeans : public HistoryPlot<SiStripApvGain, float> {
1210   public:
1211     SiStripApvTOBGainByRunMeans()
1212         : HistoryPlot<SiStripApvGain, float>("SiStripApv Gains average", "average Tracker Outer Barrel gain value") {}
1213     ~SiStripApvTOBGainByRunMeans() override = default;
1214 
1215     float getFromPayload(SiStripApvGain& payload) override {
1216       std::vector<uint32_t> detid;
1217       payload.getDetIds(detid);
1218 
1219       float nAPVs = 0;
1220       float sumOfGains = 0;
1221 
1222       for (const auto& d : detid) {
1223         int subid = DetId(d).subdetId();
1224         if (subid != StripSubdetector::TOB)
1225           continue;
1226 
1227         SiStripApvGain::Range range = payload.getRange(d);
1228         for (int it = 0; it < range.second - range.first; it++) {
1229           nAPVs += 1;
1230           sumOfGains += payload.getApvGain(it, range);
1231         }  // loop over APVs
1232       }    // loop over detIds
1233 
1234       return sumOfGains / nAPVs;
1235 
1236     }  // payload
1237   };
1238 
1239   /************************************************
1240     time history histogram of TID SiStripApvGains 
1241   *************************************************/
1242 
1243   class SiStripApvTIDGainByRunMeans : public HistoryPlot<SiStripApvGain, float> {
1244   public:
1245     SiStripApvTIDGainByRunMeans()
1246         : HistoryPlot<SiStripApvGain, float>("SiStripApv Gains average", "average Tracker Inner Disks APV gain value") {
1247     }
1248     ~SiStripApvTIDGainByRunMeans() override = default;
1249 
1250     float getFromPayload(SiStripApvGain& payload) override {
1251       std::vector<uint32_t> detid;
1252       payload.getDetIds(detid);
1253 
1254       float nAPVs = 0;
1255       float sumOfGains = 0;
1256       for (const auto& d : detid) {
1257         int subid = DetId(d).subdetId();
1258         if (subid != StripSubdetector::TID)
1259           continue;
1260 
1261         SiStripApvGain::Range range = payload.getRange(d);
1262         for (int it = 0; it < range.second - range.first; it++) {
1263           nAPVs += 1;
1264           sumOfGains += payload.getApvGain(it, range);
1265         }  // loop over APVs
1266       }    // loop over detIds
1267 
1268       return sumOfGains / nAPVs;
1269 
1270     }  // payload
1271   };
1272 
1273   /************************************************
1274     time history histogram of TEC SiStripApvGains 
1275   *************************************************/
1276 
1277   class SiStripApvTECGainByRunMeans : public HistoryPlot<SiStripApvGain, float> {
1278   public:
1279     SiStripApvTECGainByRunMeans()
1280         : HistoryPlot<SiStripApvGain, float>("SiStripApv Gains average in TEC",
1281                                              "average Tracker Endcaps APV gain value") {}
1282     ~SiStripApvTECGainByRunMeans() override = default;
1283 
1284     float getFromPayload(SiStripApvGain& payload) override {
1285       std::vector<uint32_t> detid;
1286       payload.getDetIds(detid);
1287 
1288       float nAPVs = 0;
1289       float sumOfGains = 0;
1290 
1291       for (const auto& d : detid) {
1292         int subid = DetId(d).subdetId();
1293         if (subid != StripSubdetector::TEC)
1294           continue;
1295 
1296         SiStripApvGain::Range range = payload.getRange(d);
1297         for (int it = 0; it < range.second - range.first; it++) {
1298           nAPVs += 1;
1299           sumOfGains += payload.getApvGain(it, range);
1300         }  // loop over APVs
1301       }    // loop over detIds
1302 
1303       return sumOfGains / nAPVs;
1304 
1305     }  // payload
1306   };
1307 
1308   /************************************************
1309     test class
1310   *************************************************/
1311 
1312   class SiStripApvGainsTest : public Histogram1D<SiStripApvGain, SINGLE_IOV> {
1313   public:
1314     SiStripApvGainsTest()
1315         : Histogram1D<SiStripApvGain, SINGLE_IOV>("SiStripApv Gains test", "SiStripApv Gains test", 10, 0.0, 10.0),
1316           m_trackerTopo{StandaloneTrackerTopology::fromTrackerParametersXMLFile(
1317               edm::FileInPath("Geometry/TrackerCommonData/data/trackerParameters.xml").fullPath())} {}
1318 
1319     bool fill() override {
1320       auto tag = PlotBase::getTag<0>();
1321       for (auto const& iov : tag.iovs) {
1322         std::shared_ptr<SiStripApvGain> payload = Base::fetchPayload(std::get<1>(iov));
1323         if (payload.get()) {
1324           std::vector<uint32_t> detid;
1325           payload->getDetIds(detid);
1326 
1327           SiStripDetSummary summaryGain{&m_trackerTopo};
1328 
1329           for (const auto& d : detid) {
1330             SiStripApvGain::Range range = payload->getRange(d);
1331             for (int it = 0; it < range.second - range.first; ++it) {
1332               summaryGain.add(d, payload->getApvGain(it, range));
1333               fillWithValue(payload->getApvGain(it, range));
1334             }
1335           }
1336           std::map<unsigned int, SiStripDetSummary::Values> map = summaryGain.getCounts();
1337 
1338           //SiStripPI::printSummary(map);
1339 
1340           std::stringstream ss;
1341           ss << "Summary of gain values:" << std::endl;
1342           summaryGain.print(ss, true);
1343           std::cout << ss.str() << std::endl;
1344 
1345         }  // payload
1346       }    // iovs
1347       return true;
1348     }  // fill
1349   private:
1350     TrackerTopology m_trackerTopo;
1351   };
1352 
1353   /************************************************
1354     Compare Gains from 2 IOVs, 2 pads canvas, firsr for ratio, second for scatter plot
1355   *************************************************/
1356 
1357   template <int ntags, IOVMultiplicity nIOVs>
1358   class SiStripApvGainsComparatorBase : public PlotImage<SiStripApvGain, nIOVs, ntags> {
1359   public:
1360     SiStripApvGainsComparatorBase() : PlotImage<SiStripApvGain, nIOVs, ntags>("SiStripGains Comparison") {}
1361 
1362     bool fill() override {
1363       // trick to deal with the multi-ioved tag and two tag case at the same time
1364       auto theIOVs = PlotBase::getTag<0>().iovs;
1365       auto tagname1 = PlotBase::getTag<0>().name;
1366       std::string tagname2 = "";
1367       auto firstiov = theIOVs.front();
1368       SiStripPI::MetaData lastiov;
1369 
1370       // we don't support (yet) comparison with more than 2 tags
1371       assert(this->m_plotAnnotations.ntags < 3);
1372 
1373       if (this->m_plotAnnotations.ntags == 2) {
1374         auto tag2iovs = PlotBase::getTag<1>().iovs;
1375         tagname2 = PlotBase::getTag<1>().name;
1376         lastiov = tag2iovs.front();
1377       } else {
1378         lastiov = theIOVs.back();
1379       }
1380 
1381       std::shared_ptr<SiStripApvGain> last_payload = this->fetchPayload(std::get<1>(lastiov));
1382       std::shared_ptr<SiStripApvGain> first_payload = this->fetchPayload(std::get<1>(firstiov));
1383 
1384       std::string lastIOVsince = std::to_string(std::get<0>(lastiov));
1385       std::string firstIOVsince = std::to_string(std::get<0>(firstiov));
1386 
1387       std::vector<uint32_t> detid;
1388       last_payload->getDetIds(detid);
1389 
1390       std::map<std::pair<uint32_t, int>, float> lastmap, firstmap;
1391 
1392       // loop on the last payload
1393       for (const auto& d : detid) {
1394         SiStripApvGain::Range range = last_payload->getRange(d);
1395         float Gain = 0;
1396         float nAPV = 0;
1397         for (int it = 0; it < range.second - range.first; ++it) {
1398           nAPV += 1;
1399           Gain = last_payload->getApvGain(it, range);
1400           std::pair<uint32_t, int> index = std::make_pair(d, nAPV);
1401           lastmap[index] = Gain;
1402         }  // end loop on APVs
1403       }    // end loop on detids
1404 
1405       detid.clear();
1406       first_payload->getDetIds(detid);
1407 
1408       // loop on the first payload
1409       for (const auto& d : detid) {
1410         SiStripApvGain::Range range = first_payload->getRange(d);
1411         float Gain = 0;
1412         float nAPV = 0;
1413         for (int it = 0; it < range.second - range.first; ++it) {
1414           nAPV += 1;
1415           Gain = first_payload->getApvGain(it, range);
1416           std::pair<uint32_t, int> index = std::make_pair(d, nAPV);
1417           firstmap[index] = Gain;
1418         }  // end loop on APVs
1419       }    // end loop on detids
1420 
1421       TCanvas canvas("Payload comparison", "payload comparison", 1400, 1000);
1422       canvas.Divide(2, 1);
1423 
1424       std::map<std::string, std::shared_ptr<TH1F>> ratios;
1425       std::map<std::string, std::shared_ptr<TH2F>> scatters;
1426       std::map<std::string, int> colormap;
1427       std::map<std::string, int> markermap;
1428       colormap["TIB"] = kRed;
1429       markermap["TIB"] = kFullCircle;
1430       colormap["TOB"] = kGreen;
1431       markermap["TOB"] = kFullTriangleUp;
1432       colormap["TID"] = kBlack;
1433       markermap["TID"] = kFullSquare;
1434       colormap["TEC"] = kBlue;
1435       markermap["TEC"] = kFullTriangleDown;
1436 
1437       std::vector<std::string> parts = {"TEC", "TOB", "TIB", "TID"};
1438 
1439       for (const auto& part : parts) {
1440         ratios[part] = std::make_shared<TH1F>(
1441             Form("hRatio_%s", part.c_str()),
1442             Form("Gains ratio IOV: %s/ IOV: %s ;Previous Gain (%s) / New Gain (%s);Number of APV",
1443                  firstIOVsince.c_str(),
1444                  lastIOVsince.c_str(),
1445                  firstIOVsince.c_str(),
1446                  lastIOVsince.c_str()),
1447             200,
1448             0.,
1449             2.);
1450         scatters[part] =
1451             std::make_shared<TH2F>(Form("hScatter_%s", part.c_str()),
1452                                    Form("new Gain (%s) vs previous Gain (%s);Previous Gain (%s);New Gain (%s)",
1453                                         lastIOVsince.c_str(),
1454                                         firstIOVsince.c_str(),
1455                                         firstIOVsince.c_str(),
1456                                         lastIOVsince.c_str()),
1457                                    100,
1458                                    0.5,
1459                                    1.8,
1460                                    100,
1461                                    0.5,
1462                                    1.8);
1463       }
1464 
1465       // now loop on the cached maps
1466       for (const auto& item : firstmap) {
1467         // packed index (detid,APV)
1468         auto index = item.first;
1469         auto mod = item.first.first;
1470 
1471         int subid = DetId(mod).subdetId();
1472         float ratio = firstmap[index] / lastmap[index];
1473 
1474         if (subid == StripSubdetector::TIB) {
1475           ratios["TIB"]->Fill(ratio);
1476           scatters["TIB"]->Fill(firstmap[index], lastmap[index]);
1477         }
1478 
1479         if (subid == StripSubdetector::TOB) {
1480           ratios["TOB"]->Fill(ratio);
1481           scatters["TOB"]->Fill(firstmap[index], lastmap[index]);
1482         }
1483 
1484         if (subid == StripSubdetector::TID) {
1485           ratios["TID"]->Fill(ratio);
1486           scatters["TID"]->Fill(firstmap[index], lastmap[index]);
1487         }
1488 
1489         if (subid == StripSubdetector::TEC) {
1490           ratios["TEC"]->Fill(ratio);
1491           scatters["TEC"]->Fill(firstmap[index], lastmap[index]);
1492         }
1493       }
1494 
1495       auto legend = TLegend(0.60, 0.8, 0.92, 0.95);
1496       legend.SetTextSize(0.05);
1497       canvas.cd(1)->SetLogy();
1498       canvas.cd(1)->SetTopMargin(0.05);
1499       canvas.cd(1)->SetLeftMargin(0.13);
1500       canvas.cd(1)->SetRightMargin(0.08);
1501 
1502       for (const auto& part : parts) {
1503         SiStripPI::makeNicePlotStyle(ratios[part].get());
1504         ratios[part]->SetMinimum(1.);
1505         ratios[part]->SetStats(false);
1506         ratios[part]->SetLineWidth(2);
1507         ratios[part]->SetLineColor(colormap[part]);
1508         if (part == "TEC")
1509           ratios[part]->Draw();
1510         else
1511           ratios[part]->Draw("same");
1512         legend.AddEntry(ratios[part].get(), part.c_str(), "L");
1513       }
1514 
1515       legend.Draw("same");
1516       SiStripPI::drawStatBox(ratios, colormap, parts);
1517 
1518       auto legend2 = TLegend(0.60, 0.8, 0.92, 0.95);
1519       legend2.SetTextSize(0.05);
1520       canvas.cd(2);
1521       canvas.cd(2)->SetTopMargin(0.05);
1522       canvas.cd(2)->SetLeftMargin(0.13);
1523       canvas.cd(2)->SetRightMargin(0.08);
1524 
1525       for (const auto& part : parts) {
1526         SiStripPI::makeNicePlotStyle(scatters[part].get());
1527         scatters[part]->SetStats(false);
1528         scatters[part]->SetMarkerColor(colormap[part]);
1529         scatters[part]->SetMarkerStyle(markermap[part]);
1530         scatters[part]->SetMarkerSize(0.5);
1531 
1532         auto temp = (TH2F*)(scatters[part]->Clone());
1533         temp->SetMarkerSize(1.3);
1534 
1535         if (part == "TEC")
1536           scatters[part]->Draw("P");
1537         else
1538           scatters[part]->Draw("Psame");
1539 
1540         legend2.AddEntry(temp, part.c_str(), "P");
1541       }
1542 
1543       TLine diagonal(0.5, 0.5, 1.8, 1.8);
1544       diagonal.SetLineWidth(3);
1545       diagonal.SetLineStyle(2);
1546       diagonal.Draw("same");
1547 
1548       legend2.Draw("same");
1549 
1550       std::string fileName(this->m_imageFileName);
1551       canvas.SaveAs(fileName.c_str());
1552 
1553       return true;
1554     }
1555   };
1556 
1557   using SiStripApvGainsComparatorSingleTag = SiStripApvGainsComparatorBase<1, MULTI_IOV>;
1558   using SiStripApvGainsComparatorTwoTags = SiStripApvGainsComparatorBase<2, SINGLE_IOV>;
1559 
1560   /************************************************
1561     Plot stack of gain by region
1562   *************************************************/
1563 
1564   class SiStripApvGainsTHStack : public PlotImage<SiStripApvGain, SINGLE_IOV> {
1565   public:
1566     SiStripApvGainsTHStack()
1567         : PlotImage<SiStripApvGain, SINGLE_IOV>("Stack of SiStrip APV gains values"),
1568           m_trackerTopo{StandaloneTrackerTopology::fromTrackerParametersXMLFile(
1569               edm::FileInPath("Geometry/TrackerCommonData/data/trackerParameters.xml").fullPath())} {}
1570 
1571     bool fill() override {
1572       //TH1F::SetDefaultSumw2(true);
1573       auto tag = PlotBase::getTag<0>();
1574       auto iov = tag.iovs.front();
1575 
1576       std::shared_ptr<SiStripApvGain> payload = this->fetchPayload(std::get<1>(iov));
1577       std::string IOVsince = std::to_string(std::get<0>(iov));
1578 
1579       std::map<partition, std::shared_ptr<TH1F>> h_gains;
1580 
1581       //std::vector< SiStripPI::TrackerRegion > regions;
1582       std::vector<partition> regions;
1583 
1584       std::vector<uint32_t> detid;
1585       payload->getDetIds(detid);
1586 
1587       // fill the vector of regions
1588       for (const auto& d : detid) {
1589         //auto region = this->getTheRegion(d);
1590         auto region = this->getThePartition(d);
1591         if (std::find(regions.begin(), regions.end(), region) == regions.end()) {
1592           regions.push_back(region);
1593         }
1594       }
1595 
1596       LogDebug("SiStripApvGainsTHStack") << "regions.size()=" << regions.size() << std::endl;
1597 
1598       for (const auto& r : regions) {
1599         //auto part = std::string(SiStripPI::regionType(r).second);
1600 
1601         auto part = std::string(this->partitionName(r));
1602 
1603         h_gains[r] = std::make_shared<TH1F>(Form("hGains_%s", part.c_str()),
1604                                             Form("Gains values for IOV: %s ;Gain;Number of APV", IOVsince.c_str()),
1605                                             100,
1606                                             0.5,
1607                                             1.5);
1608       }
1609 
1610       // loop on the payload
1611       for (const auto& d : detid) {
1612         //auto region = this->getTheRegion(d);
1613         auto region = this->getThePartition(d);
1614         SiStripApvGain::Range range = payload->getRange(d);
1615         for (int it = 0; it < range.second - range.first; ++it) {
1616           float gain = payload->getApvGain(it, range);
1617           h_gains[region]->Fill(gain);
1618         }  // end loop on APVs
1619       }    // end loop on detids
1620 
1621       TCanvas canvas("Payload breakout", "payload breakout", 1200, 800);
1622       canvas.Divide(2, 1);
1623 
1624       std::array<int, 6> colors = {{kRed, kBlue, kGreen, kCyan, 8, kMagenta}};
1625 
1626       THStack* hs = new THStack("hs", Form("Gains values for IOV: %s;Gain;Number of APV", IOVsince.c_str()));
1627       int colorCounter = 0;
1628       for (const auto& r : regions) {
1629         hs->Add(h_gains[r].get());
1630         SiStripPI::makeNicePlotStyle(h_gains[r].get());
1631         h_gains[r]->SetFillColor(colors[colorCounter]);
1632         //h_gains[r]->SetLineColor(colorCounter);
1633         h_gains[r]->SetLineWidth(2);
1634         colorCounter++;
1635       }
1636 
1637       TLegend legend = TLegend(0.60, 0.65, 0.95, 0.93);
1638       legend.SetTextSize(0.05);
1639       legend.SetHeader("Gain break-out", "C");  // option "C" allows to center the header
1640       for (const auto& r : regions) {
1641         auto part = std::string(this->partitionName(r));
1642         legend.AddEntry(h_gains[r].get(), part.c_str(), "F");
1643       }
1644 
1645       canvas.cd(1)->SetLogy();
1646       canvas.cd(1)->SetTopMargin(0.07);
1647       canvas.cd(1)->SetBottomMargin(0.10);
1648       canvas.cd(1)->SetLeftMargin(0.15);
1649       canvas.cd(1)->SetRightMargin(0.05);
1650       //      hs->Draw("NOSTACKB");
1651 
1652       int count(0);
1653       auto stack = hs->GetHists();
1654       double maximum = hs->GetMaximum("nostack");  //SiStripPI::getMaximum(stack);
1655 
1656       TLegend legend2 = TLegend(0.70, 0.65, 0.95, 0.93);
1657       legend2.SetTextSize(0.05);
1658       legend2.SetHeader("Partition", "C");  // option "C" allows to center the header
1659 
1660       for (const auto&& elem : *stack) {
1661         auto clone = (TH1F*)(elem->Clone(Form("hclone_%s", elem->GetName())));
1662         SiStripPI::makeNicePlotStyle(clone);
1663         clone->SetFillColor(0);
1664         clone->SetMarkerStyle(20);
1665         clone->SetLineColor(colors[count]);
1666         clone->SetMarkerColor(colors[count]);
1667         clone->SetMaximum(maximum * 10);
1668         TString candName = clone->GetName();
1669         legend2.AddEntry(clone, candName.ReplaceAll("hclone_hGains_", ""), "L");
1670         if (count == 0) {
1671           clone->Draw("HIST");
1672         } else {
1673           clone->Draw("HISTsame");
1674         }
1675         count++;
1676       }
1677 
1678       legend2.Draw("same");
1679 
1680       canvas.cd(2);  //->SetLogy();
1681       canvas.cd(2)->SetTopMargin(0.07);
1682       canvas.cd(2)->SetBottomMargin(0.10);
1683       canvas.cd(2)->SetLeftMargin(0.12);
1684       canvas.cd(2)->SetRightMargin(0.05);
1685       hs->Draw();
1686       // all graphics manipulations *after* drawing the stack!
1687       hs->GetYaxis()->SetMaxDigits(2);
1688       SiStripPI::makeNiceStyle<THStack>(hs);
1689       legend.Draw("same");
1690 
1691       std::string fileName(this->m_imageFileName);
1692       canvas.SaveAs(fileName.c_str());
1693 
1694       return true;
1695     }
1696 
1697   private:
1698     TrackerTopology m_trackerTopo;
1699     enum partition { TIB = 30, TIDP = 41, TIDM = 42, TOB = 50, TECP = 61, TECM = 62, END_OF_PARTS };
1700 
1701     const char* partitionName(partition part) {
1702       std::map<partition, const char*> mapping = {{partition::TIB, "TIB"},
1703                                                   {partition::TIDP, "TIPp"},
1704                                                   {partition::TIDM, "TIDm"},
1705                                                   {partition::TOB, "TOB"},
1706                                                   {partition::TECP, "TECp"},
1707                                                   {partition::TECM, "TECm"}};
1708 
1709       if (mapping.find(part) == mapping.end()) {
1710         throw cms::Exception("Invalid Partition passed");
1711       } else {
1712         return mapping[part];
1713       }
1714     }
1715 
1716     partition getThePartition(DetId detid) {
1717       int detNum = 0;
1718       int side = 0;
1719       switch (detid.subdetId()) {
1720         case StripSubdetector::TIB:
1721           detNum = 30;
1722           break;
1723         case StripSubdetector::TOB:
1724           detNum = 50;
1725           break;
1726         case StripSubdetector::TEC:
1727           // is this module in TEC+ or TEC-?
1728           side = m_trackerTopo.tecSide(detid);
1729           detNum = 60;
1730           break;
1731         case StripSubdetector::TID:
1732           // is this module in TID+ or TID-?
1733           side = m_trackerTopo.tidSide(detid);
1734           detNum = 40;
1735           break;
1736       }
1737 
1738       detNum += side;
1739       return static_cast<partition>(detNum);
1740     }
1741 
1742     SiStripPI::TrackerRegion getTheRegion(DetId detid) {
1743       int layer = 0;
1744       int stereo = 0;
1745       int detNum = 0;
1746 
1747       switch (detid.subdetId()) {
1748         case StripSubdetector::TIB:
1749           layer = m_trackerTopo.tibLayer(detid);
1750           stereo = m_trackerTopo.tibStereo(detid);
1751           detNum = 1000;
1752           break;
1753         case StripSubdetector::TOB:
1754           layer = m_trackerTopo.tobLayer(detid);
1755           stereo = m_trackerTopo.tobStereo(detid);
1756           detNum = 2000;
1757           break;
1758         case StripSubdetector::TEC:
1759           // is this module in TEC+ or TEC-?
1760           layer = m_trackerTopo.tecWheel(detid);
1761           stereo = m_trackerTopo.tecStereo(detid);
1762           detNum = 3000;
1763           break;
1764         case StripSubdetector::TID:
1765           // is this module in TID+ or TID-?
1766           layer = m_trackerTopo.tidWheel(detid);
1767           stereo = m_trackerTopo.tidStereo(detid);
1768           detNum = 4000;
1769           break;
1770       }
1771 
1772       detNum += layer * 10 + stereo;
1773       return static_cast<SiStripPI::TrackerRegion>(detNum);
1774     }
1775   };
1776 
1777   //*******************************************//
1778   // Compare Gains from 2 IOVs
1779   //******************************************//
1780 
1781   template <int ntags, IOVMultiplicity nIOVs>
1782   class SiStripApvGainsValuesComparatorBase : public PlotImage<SiStripApvGain, nIOVs, ntags> {
1783   public:
1784     SiStripApvGainsValuesComparatorBase()
1785         : PlotImage<SiStripApvGain, nIOVs, ntags>("Comparison of SiStrip APV gains values") {}
1786 
1787     bool fill() override {
1788       TH1F::SetDefaultSumw2(true);
1789 
1790       // trick to deal with the multi-ioved tag and two tag case at the same time
1791       auto theIOVs = PlotBase::getTag<0>().iovs;
1792       auto tagname1 = PlotBase::getTag<0>().name;
1793       std::string tagname2 = "";
1794       auto firstiov = theIOVs.front();
1795       SiStripPI::MetaData lastiov;
1796 
1797       // we don't support (yet) comparison with more than 2 tags
1798       assert(this->m_plotAnnotations.ntags < 3);
1799 
1800       if (this->m_plotAnnotations.ntags == 2) {
1801         auto tag2iovs = PlotBase::getTag<1>().iovs;
1802         tagname2 = PlotBase::getTag<1>().name;
1803         lastiov = tag2iovs.front();
1804       } else {
1805         lastiov = theIOVs.back();
1806       }
1807 
1808       std::shared_ptr<SiStripApvGain> last_payload = this->fetchPayload(std::get<1>(lastiov));
1809       std::shared_ptr<SiStripApvGain> first_payload = this->fetchPayload(std::get<1>(firstiov));
1810 
1811       std::string lastIOVsince = std::to_string(std::get<0>(lastiov));
1812       std::string firstIOVsince = std::to_string(std::get<0>(firstiov));
1813 
1814       std::vector<uint32_t> detid;
1815       last_payload->getDetIds(detid);
1816 
1817       std::map<std::pair<uint32_t, int>, float> lastmap, firstmap;
1818 
1819       // loop on the last payload
1820       for (const auto& d : detid) {
1821         SiStripApvGain::Range range = last_payload->getRange(d);
1822         float nAPV = 0;
1823         for (int it = 0; it < range.second - range.first; ++it) {
1824           nAPV += 1;
1825           auto index = std::make_pair(d, nAPV);
1826           lastmap[index] = last_payload->getApvGain(it, range);
1827         }  // end loop on APVs
1828       }    // end loop on detids
1829 
1830       detid.clear();
1831       first_payload->getDetIds(detid);
1832 
1833       // loop on the first payload
1834       for (const auto& d : detid) {
1835         SiStripApvGain::Range range = first_payload->getRange(d);
1836         float nAPV = 0;
1837         for (int it = 0; it < range.second - range.first; ++it) {
1838           nAPV += 1;
1839           auto index = std::make_pair(d, nAPV);
1840           firstmap[index] = last_payload->getApvGain(it, range);
1841         }  // end loop on APVs
1842       }    // end loop on detids
1843 
1844       TCanvas canvas("Payload comparison", "payload comparison", 1000, 1000);
1845       canvas.cd();
1846 
1847       TPad pad1("pad1", "pad1", 0, 0.3, 1, 1.0);
1848       pad1.SetBottomMargin(0.02);  // Upper and lower plot are joined
1849       pad1.SetTopMargin(0.07);
1850       pad1.SetRightMargin(0.05);
1851       pad1.SetLeftMargin(0.15);
1852       pad1.Draw();  // Draw the upper pad: pad1
1853       pad1.cd();    // pad1 becomes the current pad
1854 
1855       auto h_firstGains =
1856           std::make_shared<TH1F>("hFirstGains", "SiStrip APV gains values; APV Gains;n. APVs", 200, 0.2, 1.8);
1857       auto h_lastGains =
1858           std::make_shared<TH1F>("hLastGains", "SiStrip APV gains values; APV Gains;n. APVs", 200, 0.2, 1.8);
1859 
1860       for (const auto& item : firstmap) {
1861         h_firstGains->Fill(item.second);
1862       }
1863 
1864       for (const auto& item : lastmap) {
1865         h_lastGains->Fill(item.second);
1866       }
1867 
1868       SiStripPI::makeNicePlotStyle(h_lastGains.get());
1869       SiStripPI::makeNicePlotStyle(h_firstGains.get());
1870 
1871       TH1F* hratio = (TH1F*)h_firstGains->Clone("hratio");
1872 
1873       h_firstGains->SetLineColor(kRed);
1874       h_lastGains->SetLineColor(kBlue);
1875 
1876       h_firstGains->SetMarkerColor(kRed);
1877       h_lastGains->SetMarkerColor(kBlue);
1878 
1879       h_firstGains->SetMarkerSize(1.);
1880       h_lastGains->SetMarkerSize(1.);
1881 
1882       h_firstGains->SetLineWidth(1);
1883       h_lastGains->SetLineWidth(1);
1884 
1885       h_firstGains->SetMarkerStyle(20);
1886       h_lastGains->SetMarkerStyle(21);
1887 
1888       h_firstGains->GetXaxis()->SetLabelOffset(2.);
1889       h_lastGains->GetXaxis()->SetLabelOffset(2.);
1890 
1891       h_firstGains->Draw("HIST");
1892       h_lastGains->Draw("HISTsame");
1893 
1894       TLegend legend = TLegend(0.70, 0.7, 0.95, 0.9);
1895       legend.SetHeader("Gain Comparison", "C");  // option "C" allows to center the header
1896       legend.AddEntry(h_firstGains.get(), ("IOV: " + std::to_string(std::get<0>(firstiov))).c_str(), "PL");
1897       legend.AddEntry(h_lastGains.get(), ("IOV: " + std::to_string(std::get<0>(lastiov))).c_str(), "PL");
1898       legend.Draw("same");
1899 
1900       // lower plot will be in pad
1901       canvas.cd();  // Go back to the main canvas before defining pad2
1902       TPad pad2("pad2", "pad2", 0, 0.005, 1, 0.3);
1903       pad2.SetTopMargin(0.01);
1904       pad2.SetBottomMargin(0.2);
1905       pad2.SetRightMargin(0.05);
1906       pad2.SetLeftMargin(0.15);
1907       pad2.SetGridy();  // horizontal grid
1908       pad2.Draw();
1909       pad2.cd();  // pad2 becomes the current pad
1910 
1911       // Define the ratio plot
1912       hratio->SetLineColor(kBlack);
1913       hratio->SetMarkerColor(kBlack);
1914       hratio->SetTitle("");
1915       hratio->SetMinimum(0.55);  // Define Y ..
1916       hratio->SetMaximum(1.55);  // .. range
1917       hratio->SetStats(false);   // No statistics on lower plot
1918       hratio->Divide(h_lastGains.get());
1919       hratio->SetMarkerStyle(20);
1920       hratio->Draw("ep");  // Draw the ratio plot
1921 
1922       // Y axis ratio plot settings
1923       hratio->GetYaxis()->SetTitle(
1924           ("ratio " + std::to_string(std::get<0>(firstiov)) + " / " + std::to_string(std::get<0>(lastiov))).c_str());
1925 
1926       hratio->GetYaxis()->SetNdivisions(505);
1927 
1928       SiStripPI::makeNicePlotStyle(hratio);
1929 
1930       hratio->GetYaxis()->SetTitleSize(25);
1931       hratio->GetXaxis()->SetLabelSize(25);
1932 
1933       hratio->GetYaxis()->SetTitleFont(43);
1934       hratio->GetYaxis()->SetTitleOffset(2.5);
1935       hratio->GetYaxis()->SetLabelFont(43);  // Absolute font size in pixel (precision 3)
1936       hratio->GetYaxis()->SetLabelSize(25);
1937 
1938       // X axis ratio plot settings
1939       hratio->GetXaxis()->SetTitleSize(30);
1940       hratio->GetXaxis()->SetTitleFont(43);
1941       hratio->GetXaxis()->SetTitle("SiStrip APV Gains");
1942       hratio->GetXaxis()->SetLabelFont(43);  // Absolute font size in pixel (precision 3)
1943       hratio->GetXaxis()->SetTitleOffset(3.);
1944 
1945       std::string fileName(this->m_imageFileName);
1946       canvas.SaveAs(fileName.c_str());
1947 
1948       return true;
1949     }
1950   };
1951 
1952   using SiStripApvGainsValuesComparatorSingleTag = SiStripApvGainsValuesComparatorBase<1, MULTI_IOV>;
1953   using SiStripApvGainsValuesComparatorTwoTags = SiStripApvGainsValuesComparatorBase<2, SINGLE_IOV>;
1954 
1955   //*******************************************//
1956   // Compare Gains ratio from 2 IOVs, region by region
1957   //******************************************//
1958 
1959   template <int ntags, IOVMultiplicity nIOVs>
1960   class SiStripApvGainsRatioComparatorByRegionBase : public PlotImage<SiStripApvGain, nIOVs, ntags> {
1961   public:
1962     SiStripApvGainsRatioComparatorByRegionBase()
1963         : PlotImage<SiStripApvGain, nIOVs, ntags>("Module by Module Comparison of SiStrip APV gains"),
1964           m_trackerTopo{StandaloneTrackerTopology::fromTrackerParametersXMLFile(
1965               edm::FileInPath("Geometry/TrackerCommonData/data/trackerParameters.xml").fullPath())} {}
1966 
1967     bool fill() override {
1968       //gStyle->SetPalette(5);
1969       SiStripPI::setPaletteStyle(SiStripPI::GRAY);
1970 
1971       // trick to deal with the multi-ioved tag and two tag case at the same time
1972       auto theIOVs = PlotBase::getTag<0>().iovs;
1973       auto tagname1 = PlotBase::getTag<0>().name;
1974       std::string tagname2 = "";
1975       auto firstiov = theIOVs.front();
1976       SiStripPI::MetaData lastiov;
1977 
1978       // we don't support (yet) comparison with more than 2 tags
1979       assert(this->m_plotAnnotations.ntags < 3);
1980 
1981       if (this->m_plotAnnotations.ntags == 2) {
1982         auto tag2iovs = PlotBase::getTag<1>().iovs;
1983         tagname2 = PlotBase::getTag<1>().name;
1984         lastiov = tag2iovs.front();
1985       } else {
1986         lastiov = theIOVs.back();
1987       }
1988 
1989       std::shared_ptr<SiStripApvGain> last_payload = this->fetchPayload(std::get<1>(lastiov));
1990       std::shared_ptr<SiStripApvGain> first_payload = this->fetchPayload(std::get<1>(firstiov));
1991 
1992       std::string lastIOVsince = std::to_string(std::get<0>(lastiov));
1993       std::string firstIOVsince = std::to_string(std::get<0>(firstiov));
1994 
1995       std::vector<uint32_t> detid;
1996       last_payload->getDetIds(detid);
1997 
1998       std::map<std::pair<uint32_t, int>, float> lastmap, firstmap;
1999 
2000       // loop on the last payload
2001       for (const auto& d : detid) {
2002         SiStripApvGain::Range range = last_payload->getRange(d);
2003         float Gain = 0;
2004         float nAPV = 0;
2005         for (int it = 0; it < range.second - range.first; ++it) {
2006           nAPV += 1;
2007           Gain = last_payload->getApvGain(it, range);
2008           std::pair<uint32_t, int> index = std::make_pair(d, nAPV);
2009           lastmap[index] = Gain;
2010         }  // end loop on APVs
2011       }    // end loop on detids
2012 
2013       detid.clear();
2014       first_payload->getDetIds(detid);
2015 
2016       // loop on the first payload
2017       for (const auto& d : detid) {
2018         SiStripApvGain::Range range = first_payload->getRange(d);
2019         float Gain = 0;
2020         float nAPV = 0;
2021         for (int it = 0; it < range.second - range.first; ++it) {
2022           nAPV += 1;
2023           Gain = first_payload->getApvGain(it, range);
2024           std::pair<uint32_t, int> index = std::make_pair(d, nAPV);
2025           firstmap[index] = Gain;
2026         }  // end loop on APVs
2027       }    // end loop on detids
2028 
2029       TCanvas canvas("Payload comparison by Tracker Region", "payload comparison by Tracker Region", 1800, 800);
2030       canvas.Divide(2, 1);
2031 
2032       auto h2first = std::make_unique<TH2F>(
2033           "byRegion1", "SiStrip APV Gain values by region;; average SiStrip Gain", 38, 1., 39., 100., 0., 2.);
2034       auto h2last = std::make_unique<TH2F>(
2035           "byRegion2", "SiStrip APV Gain values by region;; average SiStrip Gain", 38, 1., 39., 100., 0., 2.);
2036 
2037       auto h2ratio = std::make_unique<TH2F>("byRegionRatio",
2038                                             Form("SiStrip APV Gains ratio by region;; Gains ratio IOV: %s/ IOV %s",
2039                                                  lastIOVsince.c_str(),
2040                                                  firstIOVsince.c_str()),
2041                                             38,
2042                                             1.,
2043                                             39.,
2044                                             100.,
2045                                             0.85,
2046                                             1.15);
2047 
2048       h2first->SetStats(false);
2049       h2last->SetStats(false);
2050       h2ratio->SetStats(false);
2051 
2052       canvas.cd(1)->SetBottomMargin(0.18);
2053       canvas.cd(1)->SetLeftMargin(0.12);
2054       canvas.cd(1)->SetRightMargin(0.05);
2055       canvas.Modified();
2056 
2057       std::vector<int> boundaries;
2058       std::string detector;
2059       std::string currentDetector;
2060 
2061       for (const auto& element : lastmap) {
2062         auto region = this->getTheRegion(element.first.first);
2063         auto bin = SiStripPI::regionType(region).first;
2064         auto label = SiStripPI::regionType(region).second;
2065 
2066         h2last->Fill(bin, element.second);
2067         h2last->GetXaxis()->SetBinLabel(bin, label);
2068         h2ratio->Fill(bin, element.second / firstmap[element.first]);
2069         h2ratio->GetXaxis()->SetBinLabel(bin, label);
2070       }
2071 
2072       for (const auto& element : firstmap) {
2073         auto region = this->getTheRegion(element.first.first);
2074         auto bin = SiStripPI::regionType(region).first;
2075         auto label = SiStripPI::regionType(region).second;
2076 
2077         h2first->Fill(bin, element.second);
2078         h2first->GetXaxis()->SetBinLabel(bin, label);
2079       }
2080 
2081       h2first->GetXaxis()->LabelsOption("v");
2082       h2last->GetXaxis()->LabelsOption("v");
2083       h2ratio->GetXaxis()->LabelsOption("v");
2084 
2085       h2last->SetLineColor(kBlue);
2086       h2first->SetLineColor(kRed);
2087       h2first->SetFillColor(kRed);
2088 
2089       h2first->SetMarkerStyle(20);
2090       h2last->SetMarkerStyle(21);
2091 
2092       h2first->SetMarkerColor(kRed);
2093       h2last->SetMarkerColor(kBlue);
2094 
2095       canvas.cd(1);
2096       h2first->Draw("BOX");
2097       h2last->Draw("BOXsame");
2098 
2099       TLegend legend = TLegend(0.70, 0.8, 0.95, 0.9);
2100       legend.SetHeader("Gain Comparison", "C");  // option "C" allows to center the header
2101       legend.AddEntry(h2first.get(), ("IOV: " + std::to_string(std::get<0>(firstiov))).c_str(), "F");
2102       legend.AddEntry(h2last.get(), ("IOV: " + std::to_string(std::get<0>(lastiov))).c_str(), "F");
2103       legend.Draw("same");
2104 
2105       canvas.cd(2);
2106       canvas.cd(2)->SetBottomMargin(0.18);
2107       canvas.cd(2)->SetLeftMargin(0.12);
2108       canvas.cd(2)->SetRightMargin(0.12);
2109 
2110       h2ratio->Draw("COLZ");
2111       auto hpfx_tmp = (TProfile*)(h2ratio->ProfileX("_pfx", 1, -1, "o"));
2112       hpfx_tmp->SetStats(kFALSE);
2113       hpfx_tmp->SetMarkerColor(kRed);
2114       hpfx_tmp->SetLineColor(kRed);
2115       hpfx_tmp->SetMarkerSize(1.2);
2116       hpfx_tmp->SetMarkerStyle(20);
2117       hpfx_tmp->Draw("same");
2118 
2119       std::string fileName(this->m_imageFileName);
2120       canvas.SaveAs(fileName.c_str());
2121 
2122       delete hpfx_tmp;
2123       return true;
2124     }
2125 
2126   private:
2127     TrackerTopology m_trackerTopo;
2128 
2129     SiStripPI::TrackerRegion getTheRegion(DetId detid) {
2130       int layer = 0;
2131       int stereo = 0;
2132       int detNum = 0;
2133 
2134       switch (detid.subdetId()) {
2135         case StripSubdetector::TIB:
2136           layer = m_trackerTopo.tibLayer(detid);
2137           stereo = m_trackerTopo.tibStereo(detid);
2138           detNum = 1000;
2139           break;
2140         case StripSubdetector::TOB:
2141           layer = m_trackerTopo.tobLayer(detid);
2142           stereo = m_trackerTopo.tobStereo(detid);
2143           detNum = 2000;
2144           break;
2145         case StripSubdetector::TEC:
2146           // is this module in TEC+ or TEC-?
2147           layer = m_trackerTopo.tecWheel(detid);
2148           stereo = m_trackerTopo.tecStereo(detid);
2149           detNum = 3000;
2150           break;
2151         case StripSubdetector::TID:
2152           // is this module in TID+ or TID-?
2153           layer = m_trackerTopo.tidWheel(detid);
2154           stereo = m_trackerTopo.tidStereo(detid);
2155           detNum = 4000;
2156           break;
2157       }
2158 
2159       detNum += layer * 10 + stereo;
2160       return static_cast<SiStripPI::TrackerRegion>(detNum);
2161     }
2162   };
2163 
2164   using SiStripApvGainsRatioComparatorByRegionSingleTag = SiStripApvGainsRatioComparatorByRegionBase<1, MULTI_IOV>;
2165   using SiStripApvGainsRatioComparatorByRegionTwoTags = SiStripApvGainsRatioComparatorByRegionBase<2, SINGLE_IOV>;
2166 
2167   /************************************************
2168     Compare Gains for each tracker region
2169   *************************************************/
2170 
2171   template <int ntags, IOVMultiplicity nIOVs>
2172   class SiStripApvGainsComparatorByRegionBase : public PlotImage<SiStripApvGain, nIOVs, ntags> {
2173   public:
2174     SiStripApvGainsComparatorByRegionBase()
2175         : PlotImage<SiStripApvGain, nIOVs, ntags>("SiStripGains Comparison By Region"),
2176           m_trackerTopo{StandaloneTrackerTopology::fromTrackerParametersXMLFile(
2177               edm::FileInPath("Geometry/TrackerCommonData/data/trackerParameters.xml").fullPath())} {}
2178 
2179     bool fill() override {
2180       // trick to deal with the multi-ioved tag and two tag case at the same time
2181       auto theIOVs = PlotBase::getTag<0>().iovs;
2182       auto tagname1 = PlotBase::getTag<0>().name;
2183       std::string tagname2 = "";
2184       auto firstiov = theIOVs.front();
2185       SiStripPI::MetaData lastiov;
2186 
2187       // we don't support (yet) comparison with more than 2 tags
2188       assert(this->m_plotAnnotations.ntags < 3);
2189 
2190       if (this->m_plotAnnotations.ntags == 2) {
2191         auto tag2iovs = PlotBase::getTag<1>().iovs;
2192         tagname2 = PlotBase::getTag<1>().name;
2193         lastiov = tag2iovs.front();
2194       } else {
2195         lastiov = theIOVs.back();
2196       }
2197 
2198       std::shared_ptr<SiStripApvGain> last_payload = this->fetchPayload(std::get<1>(lastiov));
2199       std::shared_ptr<SiStripApvGain> first_payload = this->fetchPayload(std::get<1>(firstiov));
2200 
2201       std::vector<uint32_t> detid;
2202       last_payload->getDetIds(detid);
2203 
2204       SiStripDetSummary summaryLastGain{&m_trackerTopo};
2205 
2206       for (const auto& d : detid) {
2207         SiStripApvGain::Range range = last_payload->getRange(d);
2208         for (int it = 0; it < range.second - range.first; ++it) {
2209           summaryLastGain.add(d, last_payload->getApvGain(it, range));
2210         }
2211       }
2212 
2213       SiStripDetSummary summaryFirstGain{&m_trackerTopo};
2214 
2215       for (const auto& d : detid) {
2216         SiStripApvGain::Range range = first_payload->getRange(d);
2217         for (int it = 0; it < range.second - range.first; ++it) {
2218           summaryFirstGain.add(d, first_payload->getApvGain(it, range));
2219         }
2220       }
2221 
2222       std::map<unsigned int, SiStripDetSummary::Values> firstmap = summaryFirstGain.getCounts();
2223       std::map<unsigned int, SiStripDetSummary::Values> lastmap = summaryLastGain.getCounts();
2224       //=========================
2225 
2226       TCanvas canvas("Region summary", "region summary", 1200, 1000);
2227       canvas.cd();
2228 
2229       auto hfirst = std::make_unique<TH1F>("byRegion1",
2230                                            "SiStrip APV Gain average by region;; average SiStrip Gain",
2231                                            firstmap.size(),
2232                                            0.,
2233                                            firstmap.size());
2234       auto hlast = std::make_unique<TH1F>(
2235           "byRegion2", "SiStrip APV Gain average by region;; average SiStrip Gain", lastmap.size(), 0., lastmap.size());
2236 
2237       hfirst->SetStats(false);
2238       hlast->SetStats(false);
2239 
2240       canvas.SetBottomMargin(0.18);
2241       canvas.SetLeftMargin(0.12);
2242       canvas.SetRightMargin(0.05);
2243       canvas.Modified();
2244 
2245       std::vector<int> boundaries;
2246       unsigned int iBin = 0;
2247 
2248       std::string detector;
2249       std::string currentDetector;
2250 
2251       for (const auto& element : lastmap) {
2252         iBin++;
2253         int count = element.second.count;
2254         double mean = (element.second.mean) / count;
2255 
2256         if (currentDetector.empty())
2257           currentDetector = "TIB";
2258 
2259         switch ((element.first) / 1000) {
2260           case 1:
2261             detector = "TIB";
2262             break;
2263           case 2:
2264             detector = "TOB";
2265             break;
2266           case 3:
2267             detector = "TEC";
2268             break;
2269           case 4:
2270             detector = "TID";
2271             break;
2272         }
2273 
2274         hlast->SetBinContent(iBin, mean);
2275         hlast->SetBinError(iBin, mean / 10000.);
2276         hlast->GetXaxis()->SetBinLabel(iBin, SiStripPI::regionType(element.first).second);
2277         hlast->GetXaxis()->LabelsOption("v");
2278 
2279         if (detector != currentDetector) {
2280           boundaries.push_back(iBin);
2281           currentDetector = detector;
2282         }
2283       }
2284 
2285       // reset the count
2286       iBin = 0;
2287 
2288       for (const auto& element : firstmap) {
2289         iBin++;
2290         int count = element.second.count;
2291         double mean = (element.second.mean) / count;
2292 
2293         hfirst->SetBinContent(iBin, mean);
2294         hfirst->SetBinError(iBin, mean / 10000.);
2295         hfirst->GetXaxis()->SetBinLabel(iBin, SiStripPI::regionType(element.first).second);
2296         hfirst->GetXaxis()->LabelsOption("v");
2297       }
2298 
2299       auto extrema = SiStripPI::getExtrema(hfirst.get(), hlast.get());
2300       hlast->GetYaxis()->SetRangeUser(extrema.first, extrema.second);
2301 
2302       hlast->SetMarkerStyle(20);
2303       hlast->SetMarkerSize(1);
2304       hlast->Draw("E1");
2305       hlast->Draw("Psame");
2306 
2307       hfirst->SetMarkerStyle(18);
2308       hfirst->SetMarkerSize(1);
2309       hfirst->SetLineColor(kBlue);
2310       hfirst->SetMarkerColor(kBlue);
2311       hfirst->Draw("E1same");
2312       hfirst->Draw("Psame");
2313 
2314       canvas.Update();
2315       canvas.cd();
2316 
2317       TLine l[boundaries.size()];
2318       unsigned int i = 0;
2319       for (const auto& line : boundaries) {
2320         l[i] = TLine(
2321             hfirst->GetBinLowEdge(line), canvas.cd()->GetUymin(), hfirst->GetBinLowEdge(line), canvas.cd()->GetUymax());
2322         l[i].SetLineWidth(1);
2323         l[i].SetLineStyle(9);
2324         l[i].SetLineColor(2);
2325         l[i].Draw("same");
2326         i++;
2327       }
2328 
2329       TLegend legend = TLegend(0.70, 0.8, 0.95, 0.9);
2330       legend.SetHeader("Gain Comparison", "C");  // option "C" allows to center the header
2331       legend.AddEntry(hfirst.get(), ("IOV: " + std::to_string(std::get<0>(firstiov))).c_str(), "PL");
2332       legend.AddEntry(hlast.get(), ("IOV: " + std::to_string(std::get<0>(lastiov))).c_str(), "PL");
2333       legend.Draw("same");
2334 
2335       std::string fileName(this->m_imageFileName);
2336       canvas.SaveAs(fileName.c_str());
2337 
2338       return true;
2339     }
2340 
2341   private:
2342     TrackerTopology m_trackerTopo;
2343   };
2344 
2345   using SiStripApvGainsComparatorByRegionSingleTag = SiStripApvGainsComparatorByRegionBase<1, MULTI_IOV>;
2346   using SiStripApvGainsComparatorByRegionTwoTags = SiStripApvGainsComparatorByRegionBase<2, SINGLE_IOV>;
2347 
2348   /************************************************
2349     Plot gain averages by region 
2350   *************************************************/
2351 
2352   class SiStripApvGainsByRegion : public PlotImage<SiStripApvGain, SINGLE_IOV> {
2353   public:
2354     SiStripApvGainsByRegion()
2355         : PlotImage<SiStripApvGain, SINGLE_IOV>("SiStripGains By Region"),
2356           m_trackerTopo{StandaloneTrackerTopology::fromTrackerParametersXMLFile(
2357               edm::FileInPath("Geometry/TrackerCommonData/data/trackerParameters.xml").fullPath())} {}
2358 
2359     bool fill() override {
2360       auto tag = PlotBase::getTag<0>();
2361       auto iov = tag.iovs.front();
2362       std::shared_ptr<SiStripApvGain> payload = fetchPayload(std::get<1>(iov));
2363 
2364       std::vector<uint32_t> detid;
2365       payload->getDetIds(detid);
2366 
2367       SiStripDetSummary summaryGain{&m_trackerTopo};
2368 
2369       for (const auto& d : detid) {
2370         SiStripApvGain::Range range = payload->getRange(d);
2371         for (int it = 0; it < range.second - range.first; ++it) {
2372           summaryGain.add(d, payload->getApvGain(it, range));
2373         }
2374       }
2375 
2376       std::map<unsigned int, SiStripDetSummary::Values> map = summaryGain.getCounts();
2377       //=========================
2378 
2379       TCanvas canvas("Region summary", "region summary", 1200, 1000);
2380       canvas.cd();
2381       auto h1 = std::make_unique<TH1F>(
2382           "byRegion", "SiStrip Gain average by region;; average SiStrip Gain", map.size(), 0., map.size());
2383       h1->SetStats(false);
2384       canvas.SetBottomMargin(0.18);
2385       canvas.SetLeftMargin(0.12);
2386       canvas.SetRightMargin(0.05);
2387       canvas.Modified();
2388 
2389       std::vector<int> boundaries;
2390       unsigned int iBin = 0;
2391 
2392       std::string detector;
2393       std::string currentDetector;
2394 
2395       for (const auto& element : map) {
2396         iBin++;
2397         int count = element.second.count;
2398         double mean = (element.second.mean) / count;
2399 
2400         if (currentDetector.empty())
2401           currentDetector = "TIB";
2402 
2403         switch ((element.first) / 1000) {
2404           case 1:
2405             detector = "TIB";
2406             break;
2407           case 2:
2408             detector = "TOB";
2409             break;
2410           case 3:
2411             detector = "TEC";
2412             break;
2413           case 4:
2414             detector = "TID";
2415             break;
2416         }
2417 
2418         h1->SetBinContent(iBin, mean);
2419         h1->GetXaxis()->SetBinLabel(iBin, SiStripPI::regionType(element.first).second);
2420         h1->GetXaxis()->LabelsOption("v");
2421 
2422         if (detector != currentDetector) {
2423           boundaries.push_back(iBin);
2424           currentDetector = detector;
2425         }
2426       }
2427 
2428       h1->SetMarkerStyle(20);
2429       h1->SetMarkerSize(1);
2430       h1->Draw("HIST");
2431       h1->Draw("Psame");
2432 
2433       canvas.Update();
2434 
2435       TLine l[boundaries.size()];
2436       unsigned int i = 0;
2437       for (const auto& line : boundaries) {
2438         l[i] = TLine(h1->GetBinLowEdge(line), canvas.GetUymin(), h1->GetBinLowEdge(line), canvas.GetUymax());
2439         l[i].SetLineWidth(1);
2440         l[i].SetLineStyle(9);
2441         l[i].SetLineColor(2);
2442         l[i].Draw("same");
2443         i++;
2444       }
2445 
2446       TLegend legend = TLegend(0.52, 0.82, 0.95, 0.9);
2447       legend.SetHeader((std::get<1>(iov)).c_str(), "C");  // option "C" allows to center the header
2448       legend.AddEntry(h1.get(), ("IOV: " + std::to_string(std::get<0>(iov))).c_str(), "PL");
2449       legend.SetTextSize(0.025);
2450       legend.Draw("same");
2451 
2452       std::string fileName(m_imageFileName);
2453       canvas.SaveAs(fileName.c_str());
2454 
2455       return true;
2456     }
2457 
2458   private:
2459     TrackerTopology m_trackerTopo;
2460   };
2461 
2462 }  // namespace
2463 
2464 // Register the classes as boost python plugin
2465 PAYLOAD_INSPECTOR_MODULE(SiStripApvGain) {
2466   PAYLOAD_INSPECTOR_CLASS(SiStripApvGainsValue);
2467   PAYLOAD_INSPECTOR_CLASS(SiStripApvGainTest);
2468   PAYLOAD_INSPECTOR_CLASS(SiStripApvGainByPartition);
2469   PAYLOAD_INSPECTOR_CLASS(SiStripApvGainCompareByPartition);
2470   PAYLOAD_INSPECTOR_CLASS(SiStripApvGainRatioByPartition);
2471   PAYLOAD_INSPECTOR_CLASS(SiStripApvGainDiffByPartition);
2472   PAYLOAD_INSPECTOR_CLASS(SiStripApvGainsTest);
2473   PAYLOAD_INSPECTOR_CLASS(SiStripApvGainsByRegion);
2474   PAYLOAD_INSPECTOR_CLASS(SiStripApvGainsComparatorSingleTag);
2475   PAYLOAD_INSPECTOR_CLASS(SiStripApvGainsComparatorTwoTags);
2476   PAYLOAD_INSPECTOR_CLASS(SiStripApvGainsValuesComparatorSingleTag);
2477   PAYLOAD_INSPECTOR_CLASS(SiStripApvGainsValuesComparatorTwoTags);
2478   PAYLOAD_INSPECTOR_CLASS(SiStripApvGainsComparatorByRegionSingleTag);
2479   PAYLOAD_INSPECTOR_CLASS(SiStripApvGainsComparatorByRegionTwoTags);
2480   PAYLOAD_INSPECTOR_CLASS(SiStripApvGainsRatioComparatorByRegionSingleTag);
2481   PAYLOAD_INSPECTOR_CLASS(SiStripApvGainsRatioComparatorByRegionTwoTags);
2482   PAYLOAD_INSPECTOR_CLASS(SiStripApvGainsTHStack);
2483   PAYLOAD_INSPECTOR_CLASS(SiStripApvBarrelGainsByLayer);
2484   PAYLOAD_INSPECTOR_CLASS(SiStripApvAbsoluteBarrelGainsByLayer);
2485   PAYLOAD_INSPECTOR_CLASS(SiStripApvEndcapMinusGainsByDisk);
2486   PAYLOAD_INSPECTOR_CLASS(SiStripApvEndcapPlusGainsByDisk);
2487   PAYLOAD_INSPECTOR_CLASS(SiStripApvAbsoluteEndcapMinusGainsByDisk);
2488   PAYLOAD_INSPECTOR_CLASS(SiStripApvAbsoluteEndcapPlusGainsByDisk);
2489   PAYLOAD_INSPECTOR_CLASS(SiStripApvGainsAverageTrackerMap);
2490   PAYLOAD_INSPECTOR_CLASS(SiStripApvGainsDefaultTrackerMap);
2491   PAYLOAD_INSPECTOR_CLASS(SiStripApvGainsMaximumTrackerMap);
2492   PAYLOAD_INSPECTOR_CLASS(SiStripApvGainsMinimumTrackerMap);
2493   PAYLOAD_INSPECTOR_CLASS(SiStripApvGainsAvgDeviationRatioWithPreviousIOVTrackerMap);
2494   PAYLOAD_INSPECTOR_CLASS(SiStripApvGainsAvgDeviationRatioTrackerMapTwoTags);
2495   PAYLOAD_INSPECTOR_CLASS(SiStripApvGainsMaxDeviationRatioWithPreviousIOVTrackerMap);
2496   PAYLOAD_INSPECTOR_CLASS(SiStripApvGainsMaxDeviationRatioTrackerMapTwoTags);
2497   PAYLOAD_INSPECTOR_CLASS(SiStripApvGainByRunMeans);
2498   PAYLOAD_INSPECTOR_CLASS(SiStripApvGainMin_History);
2499   PAYLOAD_INSPECTOR_CLASS(SiStripApvGainMax_History);
2500   PAYLOAD_INSPECTOR_CLASS(SiStripApvGainMean_History);
2501   PAYLOAD_INSPECTOR_CLASS(SiStripApvGainRMS_History);
2502   PAYLOAD_INSPECTOR_CLASS(SiStripApvTIBGainByRunMeans);
2503   PAYLOAD_INSPECTOR_CLASS(SiStripApvTIDGainByRunMeans);
2504   PAYLOAD_INSPECTOR_CLASS(SiStripApvTOBGainByRunMeans);
2505   PAYLOAD_INSPECTOR_CLASS(SiStripApvTECGainByRunMeans);
2506 }