Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-08-22 04:57:32

0001 // -*- C++ -*-
0002 //
0003 // Package:    SiPixelPhase1ResidualsExtra
0004 // Class:      SiPixelPhase1ResidualsExtra
0005 //
0006 /**\class 
0007 
0008  Description: Create the Phsae 1 pixel DRnR plots
0009 
0010  Implementation: Introduce some computation over the PixelPhase1 residuals distributions
0011 */
0012 //
0013 // Original Author:  Alessandro Rossi
0014 //         Created:  25th May 2021
0015 //
0016 //
0017 
0018 // system includes
0019 #include <string>
0020 #include <cstdlib>
0021 #include <iostream>
0022 #include <fstream>
0023 #include <sstream>
0024 
0025 // user includes
0026 #include "DQM/SiPixelCommon/interface/SiPixelFolderOrganizer.h"
0027 #include "DQMServices/Core/interface/DQMEDHarvester.h"
0028 #include "DQMServices/Core/interface/DQMStore.h"
0029 #include "DataFormats/DetId/interface/DetId.h"
0030 #include "DataFormats/SiPixelDetId/interface/PixelBarrelNameUpgrade.h"
0031 #include "DataFormats/SiPixelDetId/interface/PixelEndcapNameUpgrade.h"
0032 #include "DataFormats/SiPixelDetId/interface/PixelSubdetector.h"
0033 #include "DataFormats/TrackerCommon/interface/PixelBarrelName.h"
0034 #include "DataFormats/TrackerCommon/interface/PixelEndcapName.h"
0035 #include "FWCore/Framework/interface/Event.h"
0036 #include "FWCore/Framework/interface/Frameworkfwd.h"
0037 #include "FWCore/Framework/interface/MakerMacros.h"
0038 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0039 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0040 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0041 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0042 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0043 #include "FWCore/ServiceRegistry/interface/Service.h"
0044 #include "Geometry/CommonDetUnit/interface/PixelGeomDetUnit.h"
0045 #include "Geometry/CommonTopologies/interface/PixelTopology.h"
0046 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0047 #include "Geometry/Records/interface/TrackerTopologyRcd.h"
0048 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0049 
0050 using namespace std;
0051 using namespace edm;
0052 
0053 class SiPixelPhase1ResidualsExtra : public DQMEDHarvester {
0054 public:
0055   explicit SiPixelPhase1ResidualsExtra(const edm::ParameterSet& conf);
0056   ~SiPixelPhase1ResidualsExtra() override;
0057 
0058   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0059 
0060 protected:
0061   // BeginRun
0062   void beginRun(edm::Run const& run, edm::EventSetup const& eSetup) override;
0063 
0064   // EndJob
0065   void dqmEndJob(DQMStore::IBooker& iBooker, DQMStore::IGetter& iGetter) override;
0066 
0067 private:
0068   const std::string topFolderName_;
0069   const std::string inputFolderName_;
0070   const int minHits_;
0071 
0072   std::map<std::string, MonitorElement*> residuals_;
0073   std::map<std::string, MonitorElement*> DRnR_;
0074 
0075   //Book Monitoring Elements
0076   void bookMEs(DQMStore::IBooker& iBooker);
0077 
0078   //Fill Monitoring Elements
0079   void fillMEs(DQMStore::IBooker& iBooker, DQMStore::IGetter& iGetter);
0080 };
0081 
0082 SiPixelPhase1ResidualsExtra::SiPixelPhase1ResidualsExtra(const edm::ParameterSet& iConfig)
0083     : DQMEDHarvester(iConfig),
0084       topFolderName_(iConfig.getParameter<std::string>("TopFolderName")),
0085       inputFolderName_(iConfig.getParameter<std::string>("InputFolderName")),
0086       minHits_(iConfig.getParameter<int>("MinHits")) {
0087   LogInfo("PixelDQM") << "SiPixelPhase1ResidualsExtra::SiPixelPhase1ResidualsExtra: Got DQM BackEnd interface" << endl;
0088 }
0089 
0090 SiPixelPhase1ResidualsExtra::~SiPixelPhase1ResidualsExtra() {
0091   LogInfo("PixelDQM") << "SiPixelPhase1ResidualsExtra::~SiPixelPhase1ResidualsExtra: Destructor" << endl;
0092 }
0093 
0094 void SiPixelPhase1ResidualsExtra::beginRun(edm::Run const& run, edm::EventSetup const& eSetup) {}
0095 
0096 void SiPixelPhase1ResidualsExtra::dqmEndJob(DQMStore::IBooker& iBooker, DQMStore::IGetter& iGetter) {
0097   bookMEs(iBooker);
0098   fillMEs(iBooker, iGetter);
0099 }
0100 
0101 //------------------------------------------------------------------
0102 // Used to book the MEs
0103 //------------------------------------------------------------------
0104 void SiPixelPhase1ResidualsExtra::bookMEs(DQMStore::IBooker& iBooker) {
0105   iBooker.cd();
0106 
0107   //New residual plots for the PXBarrel separated by inner and outer modules per layer
0108   iBooker.setCurrentFolder(topFolderName_ + "/PXBarrel");
0109 
0110   for (std::string layer : {"1", "2", "3", "4"}) {
0111     float mean_range = 20.;
0112     float rms_range = 200.;
0113     if (layer == "1") {
0114       mean_range = 50.;
0115       rms_range = 1000.;
0116     }
0117     residuals_["residual_mean_x_Inner_PXLayer_" + layer] =
0118         iBooker.book1D("residual_mean_x_Inner_PXLayer_" + layer,
0119                        "Mean of Track Residuals X Inner Modules for Layer " + layer + ";mean(x_rec-x_pred)[#mum]",
0120                        100,
0121                        -1 * mean_range,
0122                        mean_range);
0123     residuals_["residual_mean_x_Outer_PXLayer_" + layer] =
0124         iBooker.book1D("residual_mean_x_Outer_PXLayer_" + layer,
0125                        "Mean of Track Residuals X Outer Modules for Layer " + layer + ";mean(x_rec-x_pred)[#mum]",
0126                        100,
0127                        -1 * mean_range,
0128                        mean_range);
0129     residuals_["residual_mean_y_Inner_PXLayer_" + layer] =
0130         iBooker.book1D("residual_mean_y_Inner_PXLayer_" + layer,
0131                        "Mean of Track Residuals Y Inner Modules for Layer " + layer + ";mean(y_rec-y_pred)[#mum]",
0132                        100,
0133                        -1 * mean_range,
0134                        mean_range);
0135     residuals_["residual_mean_y_Outer_PXLayer_" + layer] =
0136         iBooker.book1D("residual_mean_y_Outer_PXLayer_" + layer,
0137                        "Mean of Track Residuals Y Outer Modules for Layer " + layer + ";mean(y_rec-y_pred)[#mum]",
0138                        100,
0139                        -1 * mean_range,
0140                        mean_range);
0141 
0142     residuals_["residual_rms_x_Inner_PXLayer_" + layer] =
0143         iBooker.book1D("residual_rms_x_Inner_PXLayer_" + layer,
0144                        "RMS of Track Residuals X Inner Modules for Layer " + layer + ";rms(x_rec-x_pred)[#mum]",
0145                        100,
0146                        0.,
0147                        rms_range);
0148     residuals_["residual_rms_x_Outer_PXLayer_" + layer] =
0149         iBooker.book1D("residual_rms_x_Outer_PXLayer_" + layer,
0150                        "RMS of Track Residuals X Outer Modules for Layer " + layer + ";rms(x_rec-x_pred)[#mum]",
0151                        100,
0152                        0.,
0153                        rms_range);
0154     residuals_["residual_rms_y_Inner_PXLayer_" + layer] =
0155         iBooker.book1D("residual_rms_y_Inner_PXLayer_" + layer,
0156                        "RMS of Track Residuals Y Inner Modules for Layer " + layer + ";rms(y_rec-y_pred)[#mum]",
0157                        100,
0158                        0.,
0159                        rms_range);
0160     residuals_["residual_rms_y_Outer_PXLayer_" + layer] =
0161         iBooker.book1D("residual_rms_y_Outer_PXLayer_" + layer,
0162                        "RMS of Track Residuals Y Outer Modules for Layer " + layer + ";rms(y_rec-y_pred)[#mum]",
0163                        100,
0164                        0.,
0165                        rms_range);
0166     ///Normalized Resiuduals Plots
0167     DRnR_["NormRes_mean_x_Inner_PXLayer_" + layer] = iBooker.book1D(
0168         "NormRes_mean_x_Inner_PXLayer_" + layer,
0169         "Mean of Normalized Track Residuals X Inner Modules for Layer " + layer + ";mean((x_rec-x_pred)/x_err)",
0170         100,
0171         -1 * 1,
0172         1);
0173     DRnR_["NormRes_mean_x_Outer_PXLayer_" + layer] = iBooker.book1D(
0174         "NormRes_mean_x_Outer_PXLayer_" + layer,
0175         "Mean of Normalized Track Residuals X Outer Modules for Layer " + layer + ";mean((x_rec-x_pred)/x_err)",
0176         100,
0177         -1 * 1,
0178         1);
0179     DRnR_["NormRes_mean_y_Inner_PXLayer_" + layer] = iBooker.book1D(
0180         "NormRes_mean_y_Inner_PXLayer_" + layer,
0181         "Mean of Normalized Track Residuals Y Inner Modules for Layer " + layer + ";mean((y_rec-y_pred)/y_err)",
0182         100,
0183         -1 * 1,
0184         1);
0185     DRnR_["NormRes_mean_y_Outer_PXLayer_" + layer] = iBooker.book1D(
0186         "NormRes_mean_y_Outer_PXLayer_" + layer,
0187         "Mean of Normalized Track Residuals Y Outer Modules for Layer " + layer + ";mean((y_rec-y_pred)/y_err)",
0188         100,
0189         -1 * 1,
0190         1);
0191 
0192     DRnR_["DRnR_x_Inner_PXLayer_" + layer] = iBooker.book1D(
0193         "DRnR_x_Inner_PXLayer_" + layer,
0194         "RMS of Normalized Track Residuals X Inner Modules for Layer " + layer + ";rms((x_rec-x_pred)/x_err)",
0195         100,
0196         0.,
0197         2);
0198     DRnR_["DRnR_x_Outer_PXLayer_" + layer] = iBooker.book1D(
0199         "DRnR_x_Outer_PXLayer_" + layer,
0200         "RMS of Normalized Track Residuals X Outer Modules for Layer " + layer + ";rms((x_rec-x_pred)/x_err)",
0201         100,
0202         0.,
0203         2);
0204     DRnR_["DRnR_y_Inner_PXLayer_" + layer] = iBooker.book1D(
0205         "DRnR_y_Inner_PXLayer_" + layer,
0206         "RMS of Normalized Track Residuals Y Inner Modules for Layer " + layer + ";rms((y_rec-y_pred)/y_err)",
0207         100,
0208         0.,
0209         2);
0210     DRnR_["DRnR_y_Outer_PXLayer_" + layer] = iBooker.book1D(
0211         "DRnR_y_Outer_PXLayer_" + layer,
0212         "RMS of Normalized Track Residuals Y Outer Modules for Layer " + layer + ";rms((y_rec-y_pred)/y_err)",
0213         100,
0214         0.,
0215         2);
0216   }
0217 
0218   //New residual plots for the PXForward separated by inner and outer modules
0219   iBooker.setCurrentFolder(topFolderName_ + "/PXForward");
0220 
0221   residuals_["residual_mean_x_Inner"] = iBooker.book1D(
0222       "residual_mean_x_Inner", "Mean of Track Residuals X Inner Modules;mean(x_rec-x_pred)[#mum]", 100, -20., 20.);
0223   residuals_["residual_mean_x_Outer"] = iBooker.book1D(
0224       "residual_mean_x_Outer", "Mean of Track Residuals X Outer Modules;mean(x_rec-x_pred)[#mum]", 100, -20., 20.);
0225   residuals_["residual_mean_y_Inner"] = iBooker.book1D(
0226       "residual_mean_y_Inner", "Mean of Track Residuals Y Inner Modules;mean(y_rec-y_pred)[#mum]", 100, -20., 20.);
0227   residuals_["residual_mean_y_Outer"] = iBooker.book1D(
0228       "residual_mean_y_Outer", "Mean of Track Residuals Y Outer Modules;mean(y_rec-y_pred)[#mum]", 100, -20., 20.);
0229 
0230   residuals_["residual_rms_x_Inner"] = iBooker.book1D(
0231       "residual_rms_x_Inner", "RMS of Track Residuals X Inner Modules;rms(x_rec-x_pred)[#mum]", 100, 0., 200.);
0232   residuals_["residual_rms_x_Outer"] = iBooker.book1D(
0233       "residual_rms_x_Outer", "RMS of Track Residuals X Outer Modules;rms(x_rec-x_pred)[#mum]", 100, 0., 200.);
0234   residuals_["residual_rms_y_Inner"] = iBooker.book1D(
0235       "residual_rms_y_Inner", "RMS of Track Residuals Y Inner Modules;rms(y_rec-y_pred)[#mum]", 100, 0., 200.);
0236   residuals_["residual_rms_y_Outer"] = iBooker.book1D(
0237       "residual_rms_y_Outer", "RMS of Track Residuals Y Outer Modules;rms(y_rec-y_pred)[#mum]", 100, 0., 200.);
0238   //Normalize Residuals inner/outer
0239   DRnR_["NormRes_mean_x_Inner"] =
0240       iBooker.book1D("NormRes_mean_x_Inner",
0241                      "Mean of Normalized Track Residuals X Inner Modules;mean((x_rec-x_pred)/x_err)",
0242                      100,
0243                      -1.,
0244                      1.);
0245   DRnR_["NormRes_mean_x_Outer"] =
0246       iBooker.book1D("NormRes_mean_x_Outer",
0247                      "Mean of Normalized Track Residuals X Outer Modules;mean((x_rec-x_pred)/x_err)",
0248                      100,
0249                      -1.,
0250                      1.);
0251   DRnR_["NormRes_mean_y_Inner"] =
0252       iBooker.book1D("NormRes_mean_y_Inner",
0253                      "Mean of Normalized Track Residuals Y Inner Modules;mean((y_rec-y_pred)/y_err)",
0254                      100,
0255                      -1.,
0256                      1.);
0257   DRnR_["NormRes_mean_y_Outer"] =
0258       iBooker.book1D("NormRes_mean_y_Outer",
0259                      "Mean of Normalized Track Residuals Y Outer Modules;mean((y_rec-y_pred)/y_err)",
0260                      100,
0261                      -1.,
0262                      1.);
0263 
0264   DRnR_["DRnR_x_Inner"] = iBooker.book1D(
0265       "DRnR_x_Inner", "RMS of Normalized Track Residuals X Inner Modules;rms((x_rec-x_pred)/x_err)", 100, 0., 2.);
0266   DRnR_["DRnR_x_Outer"] = iBooker.book1D(
0267       "DRnR_x_Outer", "RMS of Normalized Track Residuals X Outer Modules;rms((x_rec-x_pred)/x_err)", 100, 0., 2.);
0268   DRnR_["DRnR_y_Inner"] = iBooker.book1D(
0269       "DRnR_y_Inner", "RMS of Normalized Track Residuals Y Inner Modules;rms((y_rec-y_pred)/y_err)", 100, 0., 2.);
0270   DRnR_["DRnR_y_Outer"] = iBooker.book1D(
0271       "DRnR_y_Outer", "RMS of Normalized Track Residuals Y Outer Modules;rms((y_rec-y_pred)/y_err)", 100, 0., 2.);
0272 
0273   //New residual plots for the PXForward separated by positive and negative side
0274   iBooker.setCurrentFolder(topFolderName_ + "/PXForward");
0275 
0276   residuals_["residual_mean_x_pos"] = iBooker.book1D(
0277       "residual_mean_x_pos", "Mean of Track Residuals X pos. Side;mean(x_rec-x_pred)[#mum]", 100, -20., 20.);
0278   residuals_["residual_mean_x_neg"] = iBooker.book1D(
0279       "residual_mean_x_neg", "Mean of Track Residuals X neg. Side;mean(x_rec-x_pred)[#mum]", 100, -20., 20.);
0280   residuals_["residual_mean_y_pos"] = iBooker.book1D(
0281       "residual_mean_y_pos", "Mean of Track Residuals Y pos. Side;mean(y_rec-y_pred)[#mum]", 100, -20., 20.);
0282   residuals_["residual_mean_y_neg"] = iBooker.book1D(
0283       "residual_mean_y_neg", "Mean of Track Residuals Y neg. Side;mean(y_rec-y_pred)[#mum]", 100, -20., 20.);
0284 
0285   residuals_["residual_rms_x_pos"] =
0286       iBooker.book1D("residual_rms_x_pos", "RMS of Track Residuals X pos. Side;rms(x_rec-x_pred)[#mum]", 100, 0., 200.);
0287   residuals_["residual_rms_x_neg"] =
0288       iBooker.book1D("residual_rms_x_neg", "RMS of Track Residuals X neg. Side;rms(x_rec-x_pred)[#mum]", 100, 0., 200.);
0289   residuals_["residual_rms_y_pos"] =
0290       iBooker.book1D("residual_rms_y_pos", "RMS of Track Residuals Y pos. Side;rms(y_rec-y_pred)[#mum]", 100, 0., 200.);
0291   residuals_["residual_rms_y_neg"] =
0292       iBooker.book1D("residual_rms_y_neg", "RMS of Track Residuals Y neg. Side;rms(y_rec-y_pred)[#mum]", 100, 0., 200.);
0293   //Normalized Residuals pos/neg
0294   DRnR_["NormRes_mean_x_pos"] = iBooker.book1D(
0295       "NormRes_mean_x_pos", "Mean of Normalized Track Residuals X pos. Side;mean((x_rec-x_pred)/x_err)", 100, -1., 1.);
0296   DRnR_["NormRes_mean_x_neg"] = iBooker.book1D(
0297       "NormRes_mean_x_neg", "Mean of Normalized Track Residuals X neg. Side;mean((x_rec-x_pred)/x_err)", 100, -1., 1.);
0298   DRnR_["NormRes_mean_y_pos"] = iBooker.book1D(
0299       "NormRes_mean_y_pos", "Mean of Normalized Track Residuals Y pos. Side;mean((y_rec-y_pred)/y_err)", 100, -1., 1.);
0300   DRnR_["NormRes_mean_y_neg"] = iBooker.book1D(
0301       "NormRes_mean_y_neg", "Mean of Normalized Track Residuals Y neg. Side;mean((y_rec-y_pred)/y_err)", 100, -1., 1.);
0302 
0303   DRnR_["DRnR_x_pos"] = iBooker.book1D(
0304       "DRnR_x_pos", "RMS of Normalized Track Residuals X pos. Side;rms((x_rec-x_pred)/x_err)", 100, 0., 2.);
0305   DRnR_["DRnR_x_neg"] = iBooker.book1D(
0306       "DRnR_x_neg", "RMS of Normalized Track Residuals X neg. Side;rms((x_rec-x_pred)/x_err)", 100, 0., 2.);
0307   DRnR_["DRnR_y_pos"] = iBooker.book1D(
0308       "DRnR_y_pos", "RMS of Normalized Track Residuals Y pos. Side;rms((y_rec-y_pred)/y_err)", 100, 0., 2.);
0309   DRnR_["DRnR_y_neg"] = iBooker.book1D(
0310       "DRnR_y_neg", "RMS of Normalized Track Residuals Y neg. Side;rms((y_rec-y_pred)/y_err)", 100, 0., 2.);
0311 
0312   //Reset the iBooker
0313   iBooker.setCurrentFolder("PixelPhase1/");
0314 }
0315 
0316 //------------------------------------------------------------------
0317 // Fill the MEs
0318 //------------------------------------------------------------------
0319 void SiPixelPhase1ResidualsExtra::fillMEs(DQMStore::IBooker& iBooker, DQMStore::IGetter& iGetter) {
0320   //Fill additional residuals plots
0321   //PXBarrel
0322 
0323   //constexpr int minHits_ = 30;  //Miniminal number of hits needed for module to be filled in histograms
0324 
0325   for (std::string layer : {"1", "2", "3", "4"}) {
0326     MonitorElement* me_x =
0327         iGetter.get(inputFolderName_ + "/PXBarrel/residual_x_per_SignedModule_per_SignedLadder_PXLayer_" + layer);
0328     MonitorElement* me_y =
0329         iGetter.get(inputFolderName_ + "/PXBarrel/residual_y_per_SignedModule_per_SignedLadder_PXLayer_" + layer);
0330     MonitorElement* me2_x = iGetter.get(
0331         inputFolderName_ + "/ResidualsExtra/PXBarrel/DRnR_x_per_SignedModule_per_SignedLadder_PXLayer_" + layer);
0332     MonitorElement* me2_y = iGetter.get(
0333         inputFolderName_ + "/ResidualsExtra/PXBarrel/DRnR_y_per_SignedModule_per_SignedLadder_PXLayer_" + layer);
0334 
0335     if (me_x == nullptr || me_y == nullptr || me2_x == nullptr || me2_y == nullptr) {
0336       edm::LogWarning("SiPixelPhase1ResidualsExtra")
0337           << "Residuals plots for Pixel BPIX Layer" << layer << " not found. Skipping ResidualsExtra plots generation.";
0338       continue;
0339     }
0340 
0341     for (int i = 1; i <= me_x->getNbinsY(); i++) {
0342       if (i == (me_x->getNbinsY() / 2 + 1))
0343         continue;  //Middle bin of y axis is empty
0344 
0345       if (i <= me_x->getNbinsY() / 2) {
0346         bool iAmInner = (i % 2 == 0);  //Check whether current ladder is inner or outer ladder
0347         if (iAmInner) {
0348           for (int j : {1, 2, 3, 4, 6, 7, 8, 9}) {
0349             if (me_x->getBinEntries(j, i) < minHits_)  //Fill only if number of hits is above threshold
0350               continue;
0351             residuals_["residual_mean_x_Inner_PXLayer_" + layer]->Fill(me_x->getBinContent(j, i) * 1e4);
0352             residuals_["residual_mean_y_Inner_PXLayer_" + layer]->Fill(me_y->getBinContent(j, i) * 1e4);
0353             residuals_["residual_rms_x_Inner_PXLayer_" + layer]->Fill(me_x->getBinError(j, i) * 1e4);
0354             residuals_["residual_rms_y_Inner_PXLayer_" + layer]->Fill(me_y->getBinError(j, i) * 1e4);
0355             DRnR_["NormRes_mean_x_Inner_PXLayer_" + layer]->Fill(me2_x->getBinContent(j, i));
0356             DRnR_["NormRes_mean_y_Inner_PXLayer_" + layer]->Fill(me2_y->getBinContent(j, i));
0357             DRnR_["DRnR_x_Inner_PXLayer_" + layer]->Fill(me2_x->getBinError(j, i));
0358             DRnR_["DRnR_y_Inner_PXLayer_" + layer]->Fill(me2_y->getBinError(j, i));
0359           }
0360         } else {
0361           for (int j : {1, 2, 3, 4, 6, 7, 8, 9}) {
0362             if (me_x->getBinEntries(j, i) < minHits_)  //Fill only if number of hits is above threshold
0363               continue;
0364             residuals_["residual_mean_x_Outer_PXLayer_" + layer]->Fill(me_x->getBinContent(j, i) * 1e4);
0365             residuals_["residual_mean_y_Outer_PXLayer_" + layer]->Fill(me_y->getBinContent(j, i) * 1e4);
0366             residuals_["residual_rms_x_Outer_PXLayer_" + layer]->Fill(me_x->getBinError(j, i) * 1e4);
0367             residuals_["residual_rms_y_Outer_PXLayer_" + layer]->Fill(me_y->getBinError(j, i) * 1e4);
0368             DRnR_["NormRes_mean_x_Outer_PXLayer_" + layer]->Fill(me2_x->getBinContent(j, i));
0369             DRnR_["NormRes_mean_y_Outer_PXLayer_" + layer]->Fill(me2_y->getBinContent(j, i));
0370             DRnR_["DRnR_x_Outer_PXLayer_" + layer]->Fill(me2_x->getBinError(j, i));
0371             DRnR_["DRnR_y_Outer_PXLayer_" + layer]->Fill(me2_y->getBinError(j, i));
0372           }
0373         }
0374       } else {
0375         bool iAmInner = (i % 2 == 1);  //Check whether current ladder is inner or outer ladder
0376         if (iAmInner) {
0377           for (int j : {1, 2, 3, 4, 6, 7, 8, 9}) {
0378             if (me_x->getBinEntries(j, i) < minHits_)  //Fill only if number of hits is above threshold
0379               continue;
0380             residuals_["residual_mean_x_Inner_PXLayer_" + layer]->Fill(me_x->getBinContent(j, i) * 1e4);
0381             residuals_["residual_mean_y_Inner_PXLayer_" + layer]->Fill(me_y->getBinContent(j, i) * 1e4);
0382             residuals_["residual_rms_x_Inner_PXLayer_" + layer]->Fill(me_x->getBinError(j, i) * 1e4);
0383             residuals_["residual_rms_y_Inner_PXLayer_" + layer]->Fill(me_y->getBinError(j, i) * 1e4);
0384             DRnR_["NormRes_mean_x_Inner_PXLayer_" + layer]->Fill(me2_x->getBinContent(j, i));
0385             DRnR_["NormRes_mean_y_Inner_PXLayer_" + layer]->Fill(me2_y->getBinContent(j, i));
0386             DRnR_["DRnR_x_Inner_PXLayer_" + layer]->Fill(me2_x->getBinError(j, i));
0387             DRnR_["DRnR_y_Inner_PXLayer_" + layer]->Fill(me2_y->getBinError(j, i));
0388           }
0389         } else {
0390           for (int j : {1, 2, 3, 4, 6, 7, 8, 9}) {
0391             if (me_x->getBinEntries(j, i) < minHits_)  //Fill only if number of hits is above threshold
0392               continue;
0393             residuals_["residual_mean_x_Outer_PXLayer_" + layer]->Fill(me_x->getBinContent(j, i) * 1e4);
0394             residuals_["residual_mean_y_Outer_PXLayer_" + layer]->Fill(me_y->getBinContent(j, i) * 1e4);
0395             residuals_["residual_rms_x_Outer_PXLayer_" + layer]->Fill(me_x->getBinError(j, i) * 1e4);
0396             residuals_["residual_rms_y_Outer_PXLayer_" + layer]->Fill(me_y->getBinError(j, i) * 1e4);
0397             DRnR_["NormRes_mean_x_Outer_PXLayer_" + layer]->Fill(me2_x->getBinContent(j, i));
0398             DRnR_["NormRes_mean_y_Outer_PXLayer_" + layer]->Fill(me2_y->getBinContent(j, i));
0399             DRnR_["DRnR_x_Outer_PXLayer_" + layer]->Fill(me2_x->getBinError(j, i));
0400             DRnR_["DRnR_y_Outer_PXLayer_" + layer]->Fill(me2_y->getBinError(j, i));
0401           }
0402         }
0403       }
0404       for (int j : {1, 2, 3, 4, 6, 7, 8, 9}) {
0405         if (me_x->getBinEntries(j, i) < minHits_) {
0406           me2_x->setBinContent(j, i, 0);
0407           me2_y->setBinContent(j, i, 0);
0408           me2_x->setBinEntries(me2_x->getBin(j, i), 0);
0409           me2_y->setBinEntries(me2_y->getBin(j, i), 0);
0410         } else {
0411           me2_x->setBinContent(j, i, me2_x->getBinError(j, i));
0412           me2_y->setBinContent(j, i, me2_y->getBinError(j, i));
0413           me2_x->setBinEntries(me2_x->getBin(j, i), 1);
0414           me2_y->setBinEntries(me2_y->getBin(j, i), 1);
0415         }
0416       }
0417     }
0418   }
0419 
0420   //PXForward separating outer and inner modules as well as positive and negative side
0421   for (std::string ring : {"1", "2"}) {
0422     MonitorElement* me_x =
0423         iGetter.get(inputFolderName_ + "/PXForward/residual_x_per_PXDisk_per_SignedBladePanel_PXRing_" + ring);
0424     MonitorElement* me_y =
0425         iGetter.get(inputFolderName_ + "/PXForward/residual_y_per_PXDisk_per_SignedBladePanel_PXRing_" + ring);
0426     MonitorElement* me2_x = iGetter.get(
0427         inputFolderName_ + "/ResidualsExtra/PXForward/DRnR_x_per_PXDisk_per_SignedBladePanel_PXRing_" + ring);
0428     MonitorElement* me2_y = iGetter.get(
0429         inputFolderName_ + "/ResidualsExtra/PXForward/DRnR_y_per_PXDisk_per_SignedBladePanel_PXRing_" + ring);
0430 
0431     if (me_x == nullptr || me_y == nullptr || me2_x == nullptr || me2_y == nullptr) {
0432       edm::LogWarning("SiPixelPhase1ResidualsExtra")
0433           << "Residuals plots for Pixel FPIX Ring" << ring << " not found. Skipping ResidualsExtra plots generation.";
0434       continue;
0435     }
0436 
0437     bool posSide = false;
0438     for (int j = 1; j <= me_x->getNbinsX(); j++) {
0439       if (j == 4)
0440         continue;  //fourth x-bin in profile plots is empty
0441 
0442       if (j == 5)
0443         posSide = true;  //change to postive side
0444 
0445       for (int i = 1; i <= me_x->getNbinsY(); i++) {
0446         if (i == me_x->getNbinsY() / 2)
0447           continue;  //Middle bins of y axis is empty
0448         if (i == (me_x->getNbinsY() / 2) + 1)
0449           continue;
0450         if (me_x->getBinEntries(j, i) >= minHits_) {  //Fill only if number of hits is above threshold
0451 
0452           bool iAmInner = (i % 2 == 0);  //separate inner and outer modules
0453           if (iAmInner) {
0454             residuals_["residual_mean_x_Inner"]->Fill(me_x->getBinContent(j, i) * 1e4);
0455             residuals_["residual_mean_y_Inner"]->Fill(me_y->getBinContent(j, i) * 1e4);
0456             residuals_["residual_rms_x_Inner"]->Fill(me_x->getBinError(j, i) * 1e4);
0457             residuals_["residual_rms_y_Inner"]->Fill(me_y->getBinError(j, i) * 1e4);
0458             DRnR_["NormRes_mean_x_Inner"]->Fill(me2_x->getBinContent(j, i));
0459             DRnR_["NormRes_mean_y_Inner"]->Fill(me2_y->getBinContent(j, i));
0460             DRnR_["DRnR_x_Inner"]->Fill(me2_x->getBinError(j, i));
0461             DRnR_["DRnR_y_Inner"]->Fill(me2_y->getBinError(j, i));
0462           } else {
0463             residuals_["residual_mean_x_Outer"]->Fill(me_x->getBinContent(j, i) * 1e4);
0464             residuals_["residual_mean_y_Outer"]->Fill(me_y->getBinContent(j, i) * 1e4);
0465             residuals_["residual_rms_x_Outer"]->Fill(me_x->getBinError(j, i) * 1e4);
0466             residuals_["residual_rms_y_Outer"]->Fill(me_y->getBinError(j, i) * 1e4);
0467             DRnR_["NormRes_mean_x_Outer"]->Fill(me2_x->getBinContent(j, i));
0468             DRnR_["NormRes_mean_y_Outer"]->Fill(me2_y->getBinContent(j, i));
0469             DRnR_["DRnR_x_Outer"]->Fill(me2_x->getBinError(j, i));
0470             DRnR_["DRnR_y_Outer"]->Fill(me2_y->getBinError(j, i));
0471           }
0472 
0473           if (!posSide) {  //separate postive and negative side
0474             residuals_["residual_mean_x_neg"]->Fill(me_x->getBinContent(j, i) * 1e4);
0475             residuals_["residual_mean_y_neg"]->Fill(me_y->getBinContent(j, i) * 1e4);
0476             residuals_["residual_rms_x_neg"]->Fill(me_x->getBinError(j, i) * 1e4);
0477             residuals_["residual_rms_y_neg"]->Fill(me_y->getBinError(j, i) * 1e4);
0478             DRnR_["NormRes_mean_x_neg"]->Fill(me2_x->getBinContent(j, i));
0479             DRnR_["NormRes_mean_y_neg"]->Fill(me2_y->getBinContent(j, i));
0480             DRnR_["DRnR_x_neg"]->Fill(me2_x->getBinError(j, i));
0481             DRnR_["DRnR_y_neg"]->Fill(me2_y->getBinError(j, i));
0482           } else {
0483             residuals_["residual_mean_x_pos"]->Fill(me_x->getBinContent(j, i) * 1e4);
0484             residuals_["residual_mean_y_pos"]->Fill(me_y->getBinContent(j, i) * 1e4);
0485             residuals_["residual_rms_x_pos"]->Fill(me_x->getBinError(j, i) * 1e4);
0486             residuals_["residual_rms_y_pos"]->Fill(me_y->getBinError(j, i) * 1e4);
0487             DRnR_["NormRes_mean_x_pos"]->Fill(me2_x->getBinContent(j, i));
0488             DRnR_["NormRes_mean_y_pos"]->Fill(me2_y->getBinContent(j, i));
0489             DRnR_["DRnR_x_pos"]->Fill(me2_x->getBinError(j, i));
0490             DRnR_["DRnR_y_pos"]->Fill(me2_y->getBinError(j, i));
0491           }
0492           me2_x->setBinContent(j, i, me2_x->getBinError(j, i));
0493           me2_y->setBinContent(j, i, me2_y->getBinError(j, i));
0494           me2_x->setBinEntries(me2_x->getBin(j, i), 1);
0495           me2_y->setBinEntries(me2_y->getBin(j, i), 1);
0496         } else {
0497           me2_x->setBinContent(j, i, 0);
0498           me2_y->setBinContent(j, i, 0);
0499           me2_x->setBinEntries(me2_x->getBin(j, i), 0);
0500           me2_y->setBinEntries(me2_y->getBin(j, i), 0);
0501         }
0502       }
0503     }
0504   }
0505 }
0506 
0507 void SiPixelPhase1ResidualsExtra::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0508   edm::ParameterSetDescription desc;
0509   desc.add<std::string>("TopFolderName", "PixelPhase1/Tracks/ResidualsExtra")
0510       ->setComment("Folder in which to write output histograms");
0511   desc.add<std::string>("InputFolderName", "")->setComment("Folder from which to fetch in the input MEs");
0512   desc.add<int>("MinHits", 30)->setComment("minimum number of hits per module to fill the DRnR plots");
0513   descriptions.addWithDefaultLabel(desc);
0514 }
0515 
0516 //define this as a plug-in
0517 DEFINE_FWK_MODULE(SiPixelPhase1ResidualsExtra);