Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 11:56:34

0001 /**
0002  * @package   Alignment/MillePedeAlignmentAlgorithm
0003  * @file      MillePedeDQMModule.cc
0004  *
0005  * @author    Max Stark (max.stark@cern.ch)
0006  * @date      Feb 19, 2016
0007  */
0008 
0009 /*** header-file ***/
0010 #include "Alignment/MillePedeAlignmentAlgorithm/plugins/MillePedeDQMModule.h"
0011 
0012 /*** ROOT objects ***/
0013 #include "TH1F.h"
0014 
0015 /*** Core framework functionality ***/
0016 #include "FWCore/Framework/interface/EventSetup.h"
0017 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0018 
0019 /*** Alignment ***/
0020 #include "Alignment/MillePedeAlignmentAlgorithm/interface/PedeLabelerBase.h"
0021 #include "Alignment/MillePedeAlignmentAlgorithm/interface/PedeLabelerPluginFactory.h"
0022 #include "Alignment/TrackerAlignment/interface/AlignableTracker.h"
0023 
0024 /*** Necessary Framework infrastructure ***/
0025 #include "FWCore/Framework/interface/ProcessBlock.h"
0026 #include "DataFormats/Alignment/interface/AlignmentToken.h"
0027 
0028 MillePedeDQMModule ::MillePedeDQMModule(const edm::ParameterSet& config)
0029     : tTopoToken_(esConsumes<edm::Transition::BeginRun>()),
0030       gDetToken_(esConsumes<edm::Transition::BeginRun>()),
0031       ptpToken_(esConsumes<edm::Transition::BeginRun>()),
0032       aliThrToken_(esConsumes<edm::Transition::BeginRun>()),
0033       geomToken_(esConsumes<edm::Transition::BeginRun>()),
0034       outputFolder_(config.getParameter<std::string>("outputFolder")),
0035       mpReaderConfig_(config.getParameter<edm::ParameterSet>("MillePedeFileReader")),
0036       isHG_(mpReaderConfig_.getParameter<bool>("isHG")) {
0037   consumes<AlignmentToken, edm::InProcess>(config.getParameter<edm::InputTag>("alignmentTokenSrc"));
0038 }
0039 
0040 MillePedeDQMModule ::~MillePedeDQMModule() {}
0041 
0042 //=============================================================================
0043 //===   INTERFACE IMPLEMENTATION                                            ===
0044 //=============================================================================
0045 
0046 void MillePedeDQMModule ::bookHistograms(DQMStore::IBooker& booker) {
0047   edm::LogInfo("MillePedeDQMModule") << "Booking histograms";
0048 
0049   booker.cd();
0050   if (!isHG_) {
0051     if (outputFolder_.find("HG") != std::string::npos) {
0052       throw cms::Exception("LogicError")
0053           << "MillePedeDQMModule is configured as Low Granularity but the outputfolder is for High Granularity";
0054     }
0055 
0056     booker.setCurrentFolder(outputFolder_);
0057     h_xPos = booker.book1D("Xpos", "Alignment fit #DeltaX;;#mum", 36, 0., 36.);
0058     h_xRot = booker.book1D("Xrot", "Alignment fit #Delta#theta_{X};;#murad", 36, 0., 36.);
0059     h_yPos = booker.book1D("Ypos", "Alignment fit #DeltaY;;#mum", 36, 0., 36.);
0060     h_yRot = booker.book1D("Yrot", "Alignment fit #Delta#theta_{Y};;#murad", 36, 0., 36.);
0061     h_zPos = booker.book1D("Zpos", "Alignment fit #DeltaZ;;#mum", 36, 0., 36.);
0062     h_zRot = booker.book1D("Zrot", "Alignment fit #Delta#theta_{Z};;#murad", 36, 0., 36.);
0063     statusResults = booker.book2D("statusResults", "Status of SiPixelAli PCL workflow;;", 6, 0., 6., 1, 0., 1.);
0064   } else {
0065     if (outputFolder_.find("HG") == std::string::npos) {
0066       throw cms::Exception("LogicError")
0067           << "MillePedeDQMModule is configured as High Granularity but the outputfolder is for Low Granularity";
0068     }
0069 
0070     booker.setCurrentFolder(outputFolder_);
0071 
0072     layerVec = {{"Layer1", pixelTopologyMap_->getPXBLadders(1)},
0073                 {"Layer2", pixelTopologyMap_->getPXBLadders(2)},
0074                 {"Layer3", pixelTopologyMap_->getPXBLadders(3)},
0075                 {"Layer4", pixelTopologyMap_->getPXBLadders(4)},
0076                 {"Disk-3", pixelTopologyMap_->getPXFBlades(-3) * 2},
0077                 {"Disk-2", pixelTopologyMap_->getPXFBlades(-2) * 2},
0078                 {"Disk-1", pixelTopologyMap_->getPXFBlades(-1) * 2},
0079                 {"Disk1", pixelTopologyMap_->getPXFBlades(1) * 2},
0080                 {"Disk2", pixelTopologyMap_->getPXFBlades(2) * 2},
0081                 {"Disk3", pixelTopologyMap_->getPXFBlades(3) * 2}};
0082 
0083     for (const auto& layer : layerVec) {
0084       h_xPos_HG[layer.first] = booker.book1D("Xpos_HG_" + layer.first,
0085                                              "Alignment fit #DeltaX for " + layer.first + ";;#mum",
0086                                              layer.second + 5,
0087                                              0.,
0088                                              layer.second + 5);
0089       h_xRot_HG[layer.first] = booker.book1D("Xrot_HG_" + layer.first,
0090                                              "Alignment fit #Delta#theta_{X} for " + layer.first + ";;#murad",
0091                                              layer.second + 5,
0092                                              0.,
0093                                              layer.second + 5);
0094       h_yPos_HG[layer.first] = booker.book1D("Ypos_HG_" + layer.first,
0095                                              "Alignment fit #DeltaY for " + layer.first + ";;#mum",
0096                                              layer.second + 5,
0097                                              0.,
0098                                              layer.second + 5);
0099       h_yRot_HG[layer.first] = booker.book1D("Yrot_HG_" + layer.first,
0100                                              "Alignment fit #Delta#theta_{Y} for " + layer.first + ";;#murad",
0101                                              layer.second + 5,
0102                                              0.,
0103                                              layer.second + 5);
0104       h_zPos_HG[layer.first] = booker.book1D("Zpos_HG_" + layer.first,
0105                                              "Alignment fit #DeltaZ for " + layer.first + ";;#mum",
0106                                              layer.second + 5,
0107                                              0.,
0108                                              layer.second + 5);
0109       h_zRot_HG[layer.first] = booker.book1D("Zrot_HG_" + layer.first,
0110                                              "Alignment fit #Delta#theta_{Z} for " + layer.first + ";;#murad",
0111                                              layer.second + 5,
0112                                              0.,
0113                                              layer.second + 5);
0114     }
0115 
0116     statusResults =
0117         booker.book2D("statusResults", "Fraction threshold check for SiPixelAliHG PCL;;", 6, 0., 6., 10, 0., 10.);
0118   }
0119 
0120   binariesAvalaible = booker.bookInt("BinariesFound");
0121   exitCode = booker.bookString("PedeExitCode", "");
0122   isVetoed = booker.bookString("IsVetoed", "");
0123 
0124   booker.cd();
0125 }
0126 
0127 void MillePedeDQMModule ::dqmEndJob(DQMStore::IBooker& booker, DQMStore::IGetter&) {
0128   bookHistograms(booker);
0129   if (mpReader_) {
0130     mpReader_->read();
0131   } else {
0132     throw cms::Exception("LogicError") << "@SUB=MillePedeDQMModule::dqmEndJob\n"
0133                                        << "Try to read MillePede results before initializing MillePedeFileReader";
0134   }
0135   if (!isHG_) {
0136     fillExpertHistos();
0137     fillStatusHisto(statusResults);
0138   } else {
0139     fillExpertHistos_HG();
0140     fillStatusHistoHG(statusResults);
0141   }
0142   binariesAvalaible->Fill(mpReader_->binariesAmount());
0143   auto theResults = mpReader_->getResults();
0144   std::string exitCodeStr = theResults.getExitMessage();
0145 
0146   std::string vetoStr{};
0147   if (mpReader_->storeAlignments()) {
0148     vetoStr = "DB Updated!"; /* easy peasy, fait accompli an alignment is there */
0149   } else {
0150     if (theResults.isHighGranularity()) { /* HG case */
0151       if (theResults.getDBVetoed() && theResults.getDBUpdated()) {
0152         vetoStr = "DB Update Vetoed"; /* this can happen in the HG PCL case */
0153       } else {
0154         vetoStr = "N/A";
0155       }
0156     } else { /* LG case */
0157       if (theResults.exceedsCutoffs()) {
0158         vetoStr = "DB Update Vetoed"; /* this can happen in the LG PCL case */
0159       } else {
0160         vetoStr = "N/A";
0161       }  // if the alignment exceeds the cutoffs
0162     }    // LG case
0163   }      // if the alignment was not stored
0164 
0165   exitCode->Fill(exitCodeStr);
0166   isVetoed->Fill(vetoStr);
0167 }
0168 
0169 //=============================================================================
0170 //===   PRIVATE METHOD IMPLEMENTATION                                       ===
0171 //=============================================================================
0172 
0173 void MillePedeDQMModule ::beginRun(const edm::Run&, const edm::EventSetup& setup) {
0174   if (!setupChanged(setup))
0175     return;
0176 
0177   const TrackerTopology* const tTopo = &setup.getData(tTopoToken_);
0178   const GeometricDet* geometricDet = &setup.getData(gDetToken_);
0179   const PTrackerParameters* ptp = &setup.getData(ptpToken_);
0180   const TrackerGeometry* geom = &setup.getData(geomToken_);
0181 
0182   pixelTopologyMap_ = std::make_shared<PixelTopologyMap>(geom, tTopo);
0183 
0184   // take the thresholds from DB
0185   const auto& thresholds_ = &setup.getData(aliThrToken_);
0186 
0187   auto myThresholds = std::make_shared<AlignPCLThresholdsHG>();
0188   myThresholds->setAlignPCLThresholds(thresholds_->getNrecords(), thresholds_->getThreshold_Map());
0189   myThresholds->setFloatMap(thresholds_->getFloatMap());
0190 
0191   TrackerGeomBuilderFromGeometricDet builder;
0192 
0193   const auto trackerGeometry = builder.build(geometricDet, *ptp, tTopo);
0194   tracker_ = std::make_unique<AlignableTracker>(trackerGeometry, tTopo);
0195 
0196   const std::string labelerPlugin{"PedeLabeler"};
0197   edm::ParameterSet labelerConfig{};
0198   labelerConfig.addUntrackedParameter("plugin", labelerPlugin);
0199   labelerConfig.addUntrackedParameter("RunRangeSelection", edm::VParameterSet{});
0200 
0201   std::shared_ptr<PedeLabelerBase> pedeLabeler{PedeLabelerPluginFactory::get()->create(
0202       labelerPlugin, PedeLabelerBase::TopLevelAlignables(tracker_.get(), nullptr, nullptr), labelerConfig)};
0203 
0204   mpReader_ = std::make_unique<MillePedeFileReader>(
0205       mpReaderConfig_, pedeLabeler, std::shared_ptr<const AlignPCLThresholdsHG>(myThresholds), pixelTopologyMap_);
0206 }
0207 
0208 void MillePedeDQMModule ::fillStatusHisto(MonitorElement* statusHisto) {
0209   TH2F* histo_status = statusHisto->getTH2F();
0210   auto theResults = mpReader_->getResults();
0211   theResults.print();
0212   histo_status->SetBinContent(1, 1, theResults.getDBUpdated());
0213   histo_status->GetXaxis()->SetBinLabel(1, "DB updated");
0214   histo_status->SetBinContent(2, 1, theResults.exceedsCutoffs());
0215   histo_status->GetXaxis()->SetBinLabel(2, "significant movement");
0216   histo_status->SetBinContent(3, 1, theResults.getDBVetoed());
0217   histo_status->GetXaxis()->SetBinLabel(3, "DB update vetoed");
0218   histo_status->SetBinContent(4, 1, !theResults.exceedsThresholds());
0219   histo_status->GetXaxis()->SetBinLabel(4, "within max movement");
0220   histo_status->SetBinContent(5, 1, !theResults.exceedsMaxError());
0221   histo_status->GetXaxis()->SetBinLabel(5, "within max error");
0222   histo_status->SetBinContent(6, 1, !theResults.belowSignificance());
0223   histo_status->GetXaxis()->SetBinLabel(6, "above significance");
0224 }
0225 
0226 void MillePedeDQMModule ::fillStatusHistoHG(MonitorElement* statusHisto) {
0227   TH2F* histo_status = statusHisto->getTH2F();
0228   auto& theResults = mpReader_->getResultsHG();
0229   histo_status->GetXaxis()->SetBinLabel(1, "#DeltaX");
0230   histo_status->GetXaxis()->SetBinLabel(2, "#DeltaY");
0231   histo_status->GetXaxis()->SetBinLabel(3, "#DeltaZ");
0232   histo_status->GetXaxis()->SetBinLabel(4, "#Delta#theta_{X}");
0233   histo_status->GetXaxis()->SetBinLabel(5, "#Delta#theta_{Y}");
0234   histo_status->GetXaxis()->SetBinLabel(6, "#Delta#theta_{Z}");
0235 
0236   int i = 0;
0237   for (const auto& result : theResults) {
0238     histo_status->GetYaxis()->SetBinLabel(i + 1, result.first.data());
0239     for (std::size_t j = 0; j < result.second.size(); ++j) {
0240       histo_status->SetBinContent(j + 1, i + 1, result.second[j]);
0241     }
0242     i++;
0243   }
0244 }
0245 
0246 void MillePedeDQMModule ::fillExpertHistos() {
0247   std::array<double, SIZE_INDEX> Xcut_, sigXcut_, maxMoveXcut_, maxErrorXcut_;
0248   std::array<double, SIZE_INDEX> tXcut_, sigtXcut_, maxMovetXcut_, maxErrortXcut_;
0249 
0250   std::array<double, SIZE_INDEX> Ycut_, sigYcut_, maxMoveYcut_, maxErrorYcut_;
0251   std::array<double, SIZE_INDEX> tYcut_, sigtYcut_, maxMovetYcut_, maxErrortYcut_;
0252 
0253   std::array<double, SIZE_INDEX> Zcut_, sigZcut_, maxMoveZcut_, maxErrorZcut_;
0254   std::array<double, SIZE_INDEX> tZcut_, sigtZcut_, maxMovetZcut_, maxErrortZcut_;
0255 
0256   auto myMap = mpReader_->getThresholdMap();
0257 
0258   std::vector<std::string> alignablesList;
0259   for (auto it = myMap.begin(); it != myMap.end(); ++it) {
0260     alignablesList.push_back(it->first);
0261   }
0262 
0263   for (auto& alignable : alignablesList) {
0264     int detIndex = getIndexFromString(alignable);
0265 
0266     Xcut_[detIndex] = myMap[alignable].getXcut();
0267     sigXcut_[detIndex] = myMap[alignable].getSigXcut();
0268     maxMoveXcut_[detIndex] = myMap[alignable].getMaxMoveXcut();
0269     maxErrorXcut_[detIndex] = myMap[alignable].getErrorXcut();
0270 
0271     Ycut_[detIndex] = myMap[alignable].getYcut();
0272     sigYcut_[detIndex] = myMap[alignable].getSigYcut();
0273     maxMoveYcut_[detIndex] = myMap[alignable].getMaxMoveYcut();
0274     maxErrorYcut_[detIndex] = myMap[alignable].getErrorYcut();
0275 
0276     Zcut_[detIndex] = myMap[alignable].getZcut();
0277     sigZcut_[detIndex] = myMap[alignable].getSigZcut();
0278     maxMoveZcut_[detIndex] = myMap[alignable].getMaxMoveZcut();
0279     maxErrorZcut_[detIndex] = myMap[alignable].getErrorZcut();
0280 
0281     tXcut_[detIndex] = myMap[alignable].getThetaXcut();
0282     sigtXcut_[detIndex] = myMap[alignable].getSigThetaXcut();
0283     maxMovetXcut_[detIndex] = myMap[alignable].getMaxMoveThetaXcut();
0284     maxErrortXcut_[detIndex] = myMap[alignable].getErrorThetaXcut();
0285 
0286     tYcut_[detIndex] = myMap[alignable].getThetaYcut();
0287     sigtYcut_[detIndex] = myMap[alignable].getSigThetaYcut();
0288     maxMovetYcut_[detIndex] = myMap[alignable].getMaxMoveThetaYcut();
0289     maxErrortYcut_[detIndex] = myMap[alignable].getErrorThetaYcut();
0290 
0291     tZcut_[detIndex] = myMap[alignable].getThetaZcut();
0292     sigtZcut_[detIndex] = myMap[alignable].getSigThetaZcut();
0293     maxMovetZcut_[detIndex] = myMap[alignable].getMaxMoveThetaZcut();
0294     maxErrortZcut_[detIndex] = myMap[alignable].getErrorThetaZcut();
0295   }
0296 
0297   fillExpertHisto(h_xPos, Xcut_, sigXcut_, maxMoveXcut_, maxErrorXcut_, mpReader_->getXobs(), mpReader_->getXobsErr());
0298   fillExpertHisto(
0299       h_xRot, tXcut_, sigtXcut_, maxMovetXcut_, maxErrortXcut_, mpReader_->getTXobs(), mpReader_->getTXobsErr());
0300 
0301   fillExpertHisto(h_yPos, Ycut_, sigYcut_, maxMoveYcut_, maxErrorYcut_, mpReader_->getYobs(), mpReader_->getYobsErr());
0302   fillExpertHisto(
0303       h_yRot, tYcut_, sigtYcut_, maxMovetYcut_, maxErrortYcut_, mpReader_->getTYobs(), mpReader_->getTYobsErr());
0304 
0305   fillExpertHisto(h_zPos, Zcut_, sigZcut_, maxMoveZcut_, maxErrorZcut_, mpReader_->getZobs(), mpReader_->getZobsErr());
0306   fillExpertHisto(
0307       h_zRot, tZcut_, sigtZcut_, maxMovetZcut_, maxErrortZcut_, mpReader_->getTZobs(), mpReader_->getTZobsErr());
0308 }
0309 
0310 void MillePedeDQMModule ::fillExpertHisto(MonitorElement* histo,
0311                                           const std::array<double, SIZE_INDEX>& cut,
0312                                           const std::array<double, SIZE_INDEX>& sigCut,
0313                                           const std::array<double, SIZE_INDEX>& maxMoveCut,
0314                                           const std::array<double, SIZE_INDEX>& maxErrorCut,
0315                                           const std::array<double, SIZE_LG_STRUCTS>& obs,
0316                                           const std::array<double, SIZE_LG_STRUCTS>& obsErr) {
0317   TH1F* histo_0 = histo->getTH1F();
0318 
0319   double max_ = *std::max_element(maxMoveCut.begin(), maxMoveCut.end());
0320 
0321   histo_0->SetMinimum(-(max_));
0322   histo_0->SetMaximum(max_);
0323 
0324   //  Schematics of the bin contents
0325   //
0326   //  XX XX XX XX XX XX    OO OO OO OO    II II II II
0327   // |--|--|--|--|--|--|--|--|--|--|--|--|--|--|--|--|--|
0328   // | 1| 2| 3| 4| 5| 6| 7| 8| 9|10|11|12|13|14|15|16|17|  ...
0329   //
0330   // |-----------------|  |-----------|  |-----------|
0331   // |observed movement|  |thresholds1|  |thresholds2|
0332 
0333   for (size_t i = 0; i < obs.size(); ++i) {
0334     // fist obs.size() bins for observed movements
0335     histo_0->SetBinContent(i + 1, obs[i]);
0336     histo_0->SetBinError(i + 1, obsErr[i]);
0337 
0338     // then at bin 8,8+5,8+10,... for cutoffs
0339     // 5 bins is the space allocated for the 4 other thresholds + 1 empty separation bin
0340     histo_0->SetBinContent(8 + i * 5, cut[i]);
0341 
0342     // then at bin 9,9+5,9+10,... for significances
0343     histo_0->SetBinContent(9 + i * 5, sigCut[i]);
0344 
0345     // then at bin 10,10+5,10+10,... for maximum movements
0346     histo_0->SetBinContent(10 + i * 5, maxMoveCut[i]);
0347 
0348     // then at bin 11,11+5,11+10,... for maximum errors
0349     histo_0->SetBinContent(11 + i * 5, maxErrorCut[i]);
0350   }
0351 }
0352 
0353 void MillePedeDQMModule ::fillExpertHistos_HG() {
0354   std::array<double, SIZE_INDEX> Xcut_, sigXcut_, maxMoveXcut_, maxErrorXcut_;
0355   std::array<double, SIZE_INDEX> tXcut_, sigtXcut_, maxMovetXcut_, maxErrortXcut_;
0356 
0357   std::array<double, SIZE_INDEX> Ycut_, sigYcut_, maxMoveYcut_, maxErrorYcut_;
0358   std::array<double, SIZE_INDEX> tYcut_, sigtYcut_, maxMovetYcut_, maxErrortYcut_;
0359 
0360   std::array<double, SIZE_INDEX> Zcut_, sigZcut_, maxMoveZcut_, maxErrorZcut_;
0361   std::array<double, SIZE_INDEX> tZcut_, sigtZcut_, maxMovetZcut_, maxErrortZcut_;
0362 
0363   auto myMap = mpReader_->getThresholdMap();
0364 
0365   std::vector<std::string> alignablesList;
0366   for (auto it = myMap.begin(); it != myMap.end(); ++it) {
0367     alignablesList.push_back(it->first);
0368   }
0369 
0370   for (auto& alignable : alignablesList) {
0371     int detIndex = getIndexFromString(alignable);
0372 
0373     Xcut_[detIndex] = myMap[alignable].getXcut();
0374     sigXcut_[detIndex] = myMap[alignable].getSigXcut();
0375     maxMoveXcut_[detIndex] = myMap[alignable].getMaxMoveXcut();
0376     maxErrorXcut_[detIndex] = myMap[alignable].getErrorXcut();
0377 
0378     Ycut_[detIndex] = myMap[alignable].getYcut();
0379     sigYcut_[detIndex] = myMap[alignable].getSigYcut();
0380     maxMoveYcut_[detIndex] = myMap[alignable].getMaxMoveYcut();
0381     maxErrorYcut_[detIndex] = myMap[alignable].getErrorYcut();
0382 
0383     Zcut_[detIndex] = myMap[alignable].getZcut();
0384     sigZcut_[detIndex] = myMap[alignable].getSigZcut();
0385     maxMoveZcut_[detIndex] = myMap[alignable].getMaxMoveZcut();
0386     maxErrorZcut_[detIndex] = myMap[alignable].getErrorZcut();
0387 
0388     tXcut_[detIndex] = myMap[alignable].getThetaXcut();
0389     sigtXcut_[detIndex] = myMap[alignable].getSigThetaXcut();
0390     maxMovetXcut_[detIndex] = myMap[alignable].getMaxMoveThetaXcut();
0391     maxErrortXcut_[detIndex] = myMap[alignable].getErrorThetaXcut();
0392 
0393     tYcut_[detIndex] = myMap[alignable].getThetaYcut();
0394     sigtYcut_[detIndex] = myMap[alignable].getSigThetaYcut();
0395     maxMovetYcut_[detIndex] = myMap[alignable].getMaxMoveThetaYcut();
0396     maxErrortYcut_[detIndex] = myMap[alignable].getErrorThetaYcut();
0397 
0398     tZcut_[detIndex] = myMap[alignable].getThetaZcut();
0399     sigtZcut_[detIndex] = myMap[alignable].getSigThetaZcut();
0400     maxMovetZcut_[detIndex] = myMap[alignable].getMaxMoveThetaZcut();
0401     maxErrortZcut_[detIndex] = myMap[alignable].getErrorThetaZcut();
0402   }
0403 
0404   fillExpertHisto_HG(
0405       h_xPos_HG, Xcut_, sigXcut_, maxMoveXcut_, maxErrorXcut_, mpReader_->getXobs_HG(), mpReader_->getXobsErr_HG());
0406   fillExpertHisto_HG(h_xRot_HG,
0407                      tXcut_,
0408                      sigtXcut_,
0409                      maxMovetXcut_,
0410                      maxErrortXcut_,
0411                      mpReader_->getTXobs_HG(),
0412                      mpReader_->getTXobsErr_HG());
0413 
0414   fillExpertHisto_HG(
0415       h_yPos_HG, Ycut_, sigYcut_, maxMoveYcut_, maxErrorYcut_, mpReader_->getYobs_HG(), mpReader_->getYobsErr_HG());
0416   fillExpertHisto_HG(h_yRot_HG,
0417                      tYcut_,
0418                      sigtYcut_,
0419                      maxMovetYcut_,
0420                      maxErrortYcut_,
0421                      mpReader_->getTYobs_HG(),
0422                      mpReader_->getTYobsErr_HG());
0423 
0424   fillExpertHisto_HG(
0425       h_zPos_HG, Zcut_, sigZcut_, maxMoveZcut_, maxErrorZcut_, mpReader_->getZobs_HG(), mpReader_->getZobsErr_HG());
0426   fillExpertHisto_HG(h_zRot_HG,
0427                      tZcut_,
0428                      sigtZcut_,
0429                      maxMovetZcut_,
0430                      maxErrortZcut_,
0431                      mpReader_->getTZobs_HG(),
0432                      mpReader_->getTZobsErr_HG());
0433 }
0434 
0435 void MillePedeDQMModule ::fillExpertHisto_HG(std::map<std::string, MonitorElement*>& histo_map,
0436                                              const std::array<double, SIZE_INDEX>& cut,
0437                                              const std::array<double, SIZE_INDEX>& sigCut,
0438                                              const std::array<double, SIZE_INDEX>& maxMoveCut,
0439                                              const std::array<double, SIZE_INDEX>& maxErrorCut,
0440                                              const std::array<double, SIZE_HG_STRUCTS>& obs,
0441                                              const std::array<double, SIZE_HG_STRUCTS>& obsErr) {
0442   int currentStart = 0;
0443   int bin = 0;
0444   double max_ = 0;
0445 
0446   for (const auto& layer : layerVec) {
0447     TH1F* histo_0 = histo_map[layer.first]->getTH1F();
0448 
0449     max_ = -1;
0450     for (int i = currentStart; i < (currentStart + layer.second); ++i) {
0451       // first obs.size() bins for observed movements
0452       bin = i - currentStart + 1;
0453 
0454       // fill observed values
0455       histo_0->SetBinContent(bin, obs[i]);
0456       histo_0->SetBinError(bin, obsErr[i]);
0457 
0458       if (std::abs(obs[i]) > max_) {
0459         max_ = std::abs(obs[i]);
0460       }
0461     }
0462 
0463     // five extra bins at the end, one empty, one with threshold, one with sigCut, one with maxMoveCut, one with MaxErrorCut
0464     histo_0->SetBinContent(bin + 1, 0);
0465     histo_0->SetBinError(bin + 1, 0);
0466 
0467     int detIndex;
0468     if (layer.first.find("Disk") != std::string::npos) {
0469       // 7 is the detId for panels, see getIndexFromString
0470       detIndex = 7;
0471       histo_0->GetXaxis()->SetTitle("Panel");
0472     } else {
0473       // 6 is the detId for ladders, see getIndexFromString
0474       detIndex = 6;
0475       histo_0->GetXaxis()->SetTitle("Ladder");
0476     }
0477 
0478     histo_0->SetBinContent(bin + 2, cut[detIndex]);
0479     histo_0->SetBinError(bin + 2, 0);
0480     histo_0->SetBinContent(bin + 3, sigCut[detIndex]);
0481     histo_0->SetBinError(bin + 3, 0);
0482     histo_0->SetBinContent(bin + 4, maxMoveCut[detIndex]);
0483     histo_0->SetBinError(bin + 4, 0);
0484     histo_0->SetBinContent(bin + 5, maxErrorCut[detIndex]);
0485     histo_0->SetBinError(bin + 5, 0);
0486 
0487     // always scale so the cutoff is visible
0488     max_ = std::max(cut[detIndex] * 1.2, max_);
0489 
0490     histo_0->SetMinimum(-(max_) * 1.2);
0491     histo_0->SetMaximum(max_ * 1.2);
0492 
0493     currentStart += layer.second;
0494   }
0495 }
0496 
0497 bool MillePedeDQMModule ::setupChanged(const edm::EventSetup& setup) {
0498   bool changed{false};
0499 
0500   if (watchIdealGeometryRcd_.check(setup))
0501     changed = true;
0502   if (watchTrackerTopologyRcd_.check(setup))
0503     changed = true;
0504   if (watchPTrackerParametersRcd_.check(setup))
0505     changed = true;
0506 
0507   return changed;
0508 }
0509 
0510 int MillePedeDQMModule ::getIndexFromString(const std::string& alignableId) {
0511   if (alignableId == "TPBHalfBarrelXminus") {
0512     return 3;
0513   } else if (alignableId == "TPBHalfBarrelXplus") {
0514     return 2;
0515   } else if (alignableId == "TPEHalfCylinderXminusZminus") {
0516     return 1;
0517   } else if (alignableId == "TPEHalfCylinderXplusZminus") {
0518     return 0;
0519   } else if (alignableId == "TPEHalfCylinderXminusZplus") {
0520     return 5;
0521   } else if (alignableId == "TPEHalfCylinderXplusZplus") {
0522     return 4;
0523   } else if (alignableId.rfind("TPBLadder", 0) == 0) {
0524     return 6;
0525   } else if (alignableId.rfind("TPEPanel", 0) == 0) {
0526     return 7;
0527   } else {
0528     throw cms::Exception("LogicError") << "@SUB=MillePedeDQMModule::getIndexFromString\n"
0529                                        << "Retrieving conversion for not supported Alignable partition" << alignableId;
0530   }
0531 }