Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-09-02 03:45:33

0001 // -*- C++ -*-
0002 //
0003 // Package:    SiPixelPhase1EfficiencyExtras
0004 // Class:      SiPixelPhase1EfficiencyExtras
0005 //
0006 /**\class 
0007 
0008  Description: Create the Phase 1 extra efficiency trend plots
0009 
0010  Implementation:
0011      <Notes on implementation>
0012 */
0013 //
0014 // Original Author:  Jack Sisson, Julie Hogan
0015 //         Created:  7 July, 2021
0016 //
0017 //
0018 
0019 // Framework
0020 #include "FWCore/ServiceRegistry/interface/Service.h"
0021 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0022 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0023 #include "FWCore/Framework/interface/Frameworkfwd.h"
0024 #include "FWCore/Framework/interface/MakerMacros.h"
0025 #include "FWCore/Framework/interface/Event.h"
0026 #include "FWCore/Framework/interface/EventSetup.h"
0027 
0028 // DQM Framework
0029 #include "DQMServices/Core/interface/DQMStore.h"
0030 #include "DQMServices/Core/interface/DQMEDHarvester.h"
0031 
0032 using namespace std;
0033 using namespace edm;
0034 
0035 class SiPixelPhase1EfficiencyExtras : public DQMEDHarvester {
0036 public:
0037   explicit SiPixelPhase1EfficiencyExtras(const edm::ParameterSet& conf);
0038   ~SiPixelPhase1EfficiencyExtras() override;
0039 
0040 protected:
0041   void beginRun(edm::Run const& run, edm::EventSetup const& eSetup) override;
0042 
0043   void dqmEndJob(DQMStore::IBooker& iBooker, DQMStore::IGetter& iGetter) override;
0044 
0045   const std::string effFolderName_;
0046   const std::string vtxFolderName_;
0047   const std::string instLumiFolderName_;
0048 };
0049 
0050 SiPixelPhase1EfficiencyExtras::SiPixelPhase1EfficiencyExtras(const edm::ParameterSet& iConfig)
0051     : effFolderName_(iConfig.getParameter<std::string>("EffFolderName")),
0052       vtxFolderName_(iConfig.getParameter<std::string>("VtxFolderName")),
0053       instLumiFolderName_(iConfig.getParameter<std::string>("InstLumiFolderName")) {
0054   LogInfo("PixelDQM") << "SiPixelPhase1EfficiencyExtras::SiPixelPhase1EfficiencyExtras: Hello!" << endl;
0055 }
0056 
0057 SiPixelPhase1EfficiencyExtras::~SiPixelPhase1EfficiencyExtras() {
0058   LogInfo("PixelDQM") << "SiPixelPhase1EfficiencyExtras::~SiPixelPhase1EfficiencyExtras: Destructor" << endl;
0059 }
0060 
0061 void SiPixelPhase1EfficiencyExtras::beginRun(edm::Run const& run, edm::EventSetup const& eSetup) {}
0062 
0063 //------------------------------------------------------------------
0064 // Method called for every event
0065 //------------------------------------------------------------------
0066 void SiPixelPhase1EfficiencyExtras::dqmEndJob(DQMStore::IBooker& iBooker, DQMStore::IGetter& iGetter) {
0067   iBooker.setCurrentFolder(effFolderName_);
0068 
0069   //Get the existing histos
0070   MonitorElement* vtx_v_lumi = iGetter.get(vtxFolderName_ + "/NumberOfGoodPVtxVsLS_GenTk");
0071 
0072   MonitorElement* scalLumi_v_lumi = iGetter.get(instLumiFolderName_ + "/lumiVsLS");
0073 
0074   MonitorElement* eff_v_lumi_forward =
0075       iGetter.get(effFolderName_ + "/hitefficiency_per_Lumisection_per_PXDisk_PXForward");
0076 
0077   MonitorElement* eff_v_lumi_barrel =
0078       iGetter.get(effFolderName_ + "/hitefficiency_per_Lumisection_per_PXLayer_PXBarrel");
0079 
0080   //set up some booleans that will tell us which graphs to create
0081   bool createNvtx = true;
0082   bool createInstLumi = true;
0083 
0084   //check which of the MEs exist and respond appropriately
0085   if (!eff_v_lumi_forward) {
0086     edm::LogWarning("SiPixelPhase1EfficiencyExtras")
0087         << "no hitefficiency_per_Lumisection_per_PXDisk_PXForward ME is available in " << effFolderName_ << std::endl;
0088     return;
0089   }
0090   if (!eff_v_lumi_barrel) {
0091     edm::LogWarning("SiPixelPhase1EfficiencyExtras")
0092         << "no hitefficiency_per_Lumisection_per_PXLayer_PXBarrel ME is available in " << effFolderName_ << std::endl;
0093     return;
0094   }
0095   if (!vtx_v_lumi) {
0096     edm::LogWarning("SiPixelPhase1EfficiencyExtras")
0097         << "no NumberOfGoodPVtxVsLS_GenTK ME is available in " << vtxFolderName_ << std::endl;
0098     createNvtx = false;
0099   }
0100   if (!scalLumi_v_lumi) {
0101     edm::LogWarning("SiPixelPhase1EfficiencyExtras")
0102         << "no lumiVsLS ME is available in " << instLumiFolderName_ << std::endl;
0103     createInstLumi = false;
0104   }
0105 
0106   //If the existing MEs are empty, set the boolean to skip booking
0107   if (vtx_v_lumi && vtx_v_lumi->getEntries() == 0)
0108     createNvtx = false;
0109   if (scalLumi_v_lumi && scalLumi_v_lumi->getEntries() == 0)
0110     createInstLumi = false;
0111 
0112   double eff = 0.0;
0113 
0114   //Will pass if nvtx ME exists and is not empty
0115   if (createNvtx) {
0116     //Book new histos
0117     MonitorElement* eff_v_vtx_barrel =
0118         iBooker.book2D("hitefficiency_per_meanNvtx_per_PXLayer_PXBarrel",
0119                        "hitefficiency_per_meanNvtx_per_PXLayer_PXBarrel; meanNvtx; PXLayer",
0120                        500,
0121                        0,
0122                        100,
0123                        3,
0124                        .5,
0125                        3.5);
0126 
0127     MonitorElement* eff_v_vtx_forward =
0128         iBooker.book2D("hitefficiency_per_meanNvtx_per_PXDisk_PXForward",
0129                        "hitefficiency_per_meanNvtx_per_PXDisk_PXForward; meanNvtx; PXDisk",
0130                        500,
0131                        0,
0132                        100,
0133                        7,
0134                        -3.5,
0135                        3.5);
0136 
0137     //initialize variables
0138     int numLumiNvtx = int(vtx_v_lumi->getNbinsX());
0139     double nvtx = 0.0;
0140     int binNumVtx = 0;
0141 
0142     //For loop to loop through lumisections
0143     for (int iLumi = 1; iLumi < numLumiNvtx - 1; iLumi++) {
0144       //get the meanNvtx for each lumi
0145       nvtx = vtx_v_lumi->getBinContent(iLumi);
0146 
0147       //Filter out useless iterations
0148       if (nvtx != 0) {
0149         //Grab the bin number for the nvtx
0150         binNumVtx = eff_v_vtx_barrel->getTH2F()->FindBin(nvtx);
0151 
0152         //loop through the layers
0153         for (int iLayer = 1; iLayer < 8; iLayer++) {
0154           //get the eff at the lumisection and layer
0155           eff = eff_v_lumi_forward->getBinContent(iLumi - 1, iLayer);
0156 
0157           //set the efficiency in the new histo
0158           eff_v_vtx_forward->setBinContent(binNumVtx, iLayer, eff);
0159         }
0160 
0161         //loop through the layers
0162         for (int iLayer = 1; iLayer < 5; iLayer++) {
0163           //get the efficiency for each lumi at each layer
0164           eff = eff_v_lumi_barrel->getBinContent(iLumi - 1, iLayer);
0165 
0166           //set the efficiency
0167           eff_v_vtx_barrel->setBinContent(binNumVtx, iLayer, eff);
0168         }
0169       }
0170     }
0171   }
0172   // Will pass if InstLumi ME exists and is not empty
0173   if (createInstLumi) {
0174     //Get the max value of inst lumi for plot
0175     int yMax2 = scalLumi_v_lumi->getTProfile()->GetMaximum();
0176     yMax2 = yMax2 + yMax2 * .1;
0177 
0178     //Book new histos
0179     MonitorElement* eff_v_scalLumi_barrel =
0180         iBooker.book2D("hitefficiency_per_scalLumi_per_PXLayer_PXBarrel",
0181                        "hitefficiency_per_scalLumi_per_PXLayer_PXBarrel; scal inst lumi E30; PXLayer",
0182                        500,
0183                        0,
0184                        yMax2,
0185                        3,
0186                        .5,
0187                        3.5);
0188 
0189     MonitorElement* eff_v_scalLumi_forward =
0190         iBooker.book2D("hitefficiency_per_scalLumi_per_PXDisk_PXForward",
0191                        "hitefficiency_per_scalLumi_per_PXDisk_PXForward; scal inst lumi E30; PXDisk",
0192                        500,
0193                        0,
0194                        yMax2,
0195                        7,
0196                        -3.5,
0197                        3.5);
0198 
0199     //initialize variables
0200     int numLumiScal = int(scalLumi_v_lumi->getNbinsX());
0201     double scalLumi = 0.0;
0202     int binNumScal = 0;
0203 
0204     //For loop to loop through lumisections
0205     for (int iLumi = 1; iLumi < numLumiScal - 1; iLumi++) {
0206       //get the inst lumi for each lumi
0207       scalLumi = scalLumi_v_lumi->getBinContent(iLumi);
0208 
0209       //Filter out useless iterations
0210       if (scalLumi != 0) {
0211         //Grab the bin number for the inst lumi
0212         binNumScal = eff_v_scalLumi_barrel->getTH2F()->FindBin(scalLumi);
0213 
0214         //loop through the layers
0215         for (int iLayer = 1; iLayer < 8; iLayer++) {
0216           //get the eff at the lumisection and layer
0217           eff = eff_v_lumi_forward->getBinContent(iLumi - 1, iLayer);
0218 
0219           //set the efficiency in the new histo
0220           eff_v_scalLumi_forward->setBinContent(binNumScal, iLayer, eff);
0221         }
0222 
0223         //loop through the layers
0224         for (int iLayer = 1; iLayer < 5; iLayer++) {
0225           //get the eff at the lumisection and layer
0226           eff = eff_v_lumi_barrel->getBinContent(iLumi - 1, iLayer);
0227 
0228           //set the efficiency in the new histo
0229           eff_v_scalLumi_barrel->setBinContent(binNumScal, iLayer, eff);
0230         }
0231       }
0232     }
0233   } else
0234     return;
0235 }
0236 
0237 //define this as a plug-in
0238 DEFINE_FWK_MODULE(SiPixelPhase1EfficiencyExtras);