Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-03-17 11:28:14

0001 // system include files
0002 #include <memory>
0003 
0004 // user include files
0005 #include "FWCore/Framework/interface/Event.h"
0006 #include "FWCore/Framework/interface/Frameworkfwd.h"
0007 #include "FWCore/Framework/interface/MakerMacros.h"
0008 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0009 #include "FWCore/ServiceRegistry/interface/Service.h"
0010 #include "FWCore/Utilities/interface/InputTag.h"
0011 #include "TGraphAsymmErrors.h"
0012 
0013 /// Geometry
0014 #include "Geometry/CommonDetUnit/interface/GeomDet.h"
0015 #include "Geometry/Records/interface/MuonGeometryRecord.h"
0016 
0017 #include "Geometry/CommonTopologies/interface/StripTopology.h"
0018 #include "Geometry/GEMGeometry/interface/ME0EtaPartition.h"
0019 #include "Geometry/GEMGeometry/interface/ME0EtaPartitionSpecs.h"
0020 #include "Geometry/GEMGeometry/interface/ME0Geometry.h"
0021 
0022 #include "DQMServices/Core/interface/DQMStore.h"
0023 
0024 /// Log messages
0025 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0026 #include "FWCore/ServiceRegistry/interface/Service.h"
0027 
0028 #include "Validation/MuonME0Validation/plugins/MuonME0SegHarvestor.h"
0029 
0030 MuonME0SegHarvestor::MuonME0SegHarvestor(const edm::ParameterSet &ps) {
0031   dbe_path_ = std::string("MuonME0RecHitsV/ME0SegmentsTask/");
0032 }
0033 
0034 MuonME0SegHarvestor::~MuonME0SegHarvestor() {}
0035 
0036 TProfile *MuonME0SegHarvestor::ComputeEff(TH1F *num, TH1F *denum, std::string nameHist) {
0037   std::string name = "eff_" + nameHist;
0038   std::string title = "Segment Efficiency" + std::string(num->GetTitle());
0039   TProfile *efficHist = new TProfile(name.c_str(),
0040                                      title.c_str(),
0041                                      denum->GetXaxis()->GetNbins(),
0042                                      denum->GetXaxis()->GetXmin(),
0043                                      denum->GetXaxis()->GetXmax());
0044 
0045   for (int i = 1; i <= denum->GetNbinsX(); i++) {
0046     double nNum = num->GetBinContent(i);
0047     double nDenum = denum->GetBinContent(i);
0048     if (nDenum == 0 || nNum == 0) {
0049       continue;
0050     }
0051     if (nNum > nDenum) {
0052       double temp = nDenum;
0053       nDenum = nNum;
0054       nNum = temp;
0055       edm::LogWarning("MuonME0SegHarvestor")
0056           << "Alert! specific bin's num is bigger than denum " << i << " " << nNum << " " << nDenum;
0057     }
0058     const double effVal = nNum / nDenum;
0059     efficHist->SetBinContent(i, effVal);
0060     efficHist->SetBinEntries(i, 1);
0061     efficHist->SetBinError(i, 0);
0062     const double errLo = TEfficiency::ClopperPearson((int)nDenum, (int)nNum, 0.683, false);
0063     const double errUp = TEfficiency::ClopperPearson((int)nDenum, (int)nNum, 0.683, true);
0064     const double errVal = (effVal - errLo > errUp - effVal) ? effVal - errLo : errUp - effVal;
0065     efficHist->SetBinError(i, sqrt(effVal * effVal + errVal * errVal));
0066   }
0067   return efficHist;
0068 }
0069 
0070 void MuonME0SegHarvestor::ProcessBooking(
0071     DQMStore::IBooker &ibooker, DQMStore::IGetter &ig, std::string nameHist, TH1F *num, TH1F *den) {
0072   if (num != nullptr && den != nullptr) {
0073     TProfile *profile = ComputeEff(num, den, nameHist);
0074 
0075     TString x_axis_title = TString(num->GetXaxis()->GetTitle());
0076     TString title = TString::Format("Segment Efficiency;%s;Eff.", x_axis_title.Data());
0077 
0078     profile->SetTitle(title.Data());
0079     ibooker.bookProfile(profile->GetName(), profile);
0080 
0081     delete profile;
0082 
0083   } else {
0084     edm::LogWarning("MuonME0SegHarvestor") << "Can not find histograms";
0085     if (num == nullptr)
0086       edm::LogWarning("MuonME0SegHarvestor") << "num not found";
0087     if (den == nullptr)
0088       edm::LogWarning("MuonME0SegHarvestor") << "den not found";
0089   }
0090   return;
0091 }
0092 
0093 void MuonME0SegHarvestor::dqmEndJob(DQMStore::IBooker &ibooker, DQMStore::IGetter &ig) {
0094   ig.setCurrentFolder(dbe_path_);
0095 
0096   TString eta_label_den = TString(dbe_path_) + "me0_simsegment_eta";
0097   TString eta_label_num = TString(dbe_path_) + "me0_matchedsimsegment_eta";
0098   TString pt_label_den = TString(dbe_path_) + "me0_simsegment_pt";
0099   TString pt_label_num = TString(dbe_path_) + "me0_matchedsimsegment_pt";
0100   TString phi_label_den = TString(dbe_path_) + "me0_simsegment_phi";
0101   TString phi_label_num = TString(dbe_path_) + "me0_matchedsimsegment_phi";
0102 
0103   if (ig.get(eta_label_num.Data()) != nullptr && ig.get(eta_label_den.Data()) != nullptr) {
0104     TH1F *num_vs_eta = (TH1F *)ig.get(eta_label_num.Data())->getTH1F()->Clone();
0105     num_vs_eta->Sumw2();
0106     TH1F *den_vs_eta = (TH1F *)ig.get(eta_label_den.Data())->getTH1F()->Clone();
0107     den_vs_eta->Sumw2();
0108 
0109     ProcessBooking(ibooker, ig, "me0segment_eff_vs_eta", num_vs_eta, den_vs_eta);
0110 
0111     delete num_vs_eta;
0112     delete den_vs_eta;
0113 
0114   } else
0115     edm::LogWarning("MuonME0SegHarvestor") << "Can not find histograms: " << eta_label_num << " or " << eta_label_den;
0116 
0117   if (ig.get(pt_label_num.Data()) != nullptr && ig.get(pt_label_den.Data()) != nullptr) {
0118     TH1F *num_vs_pt = (TH1F *)ig.get(pt_label_num.Data())->getTH1F()->Clone();
0119     num_vs_pt->Sumw2();
0120     TH1F *den_vs_pt = (TH1F *)ig.get(pt_label_den.Data())->getTH1F()->Clone();
0121     den_vs_pt->Sumw2();
0122 
0123     ProcessBooking(ibooker, ig, "me0segment_eff_vs_pt", num_vs_pt, den_vs_pt);
0124 
0125     delete num_vs_pt;
0126     delete den_vs_pt;
0127 
0128   } else
0129     edm::LogWarning("MuonME0SegHarvestor") << "Can not find histograms: " << pt_label_num << " or " << pt_label_den;
0130 
0131   if (ig.get(phi_label_num.Data()) != nullptr && ig.get(phi_label_den.Data()) != nullptr) {
0132     TH1F *num_vs_phi = (TH1F *)ig.get(phi_label_num.Data())->getTH1F()->Clone();
0133     num_vs_phi->Sumw2();
0134     TH1F *den_vs_phi = (TH1F *)ig.get(phi_label_den.Data())->getTH1F()->Clone();
0135     den_vs_phi->Sumw2();
0136 
0137     ProcessBooking(ibooker, ig, "me0segment_eff_vs_phi", num_vs_phi, den_vs_phi);
0138 
0139     delete num_vs_phi;
0140     delete den_vs_phi;
0141 
0142   } else
0143     edm::LogWarning("MuonME0SegHarvestor") << "Can not find histograms: " << phi_label_num << " or " << phi_label_den;
0144 }
0145 
0146 // define this as a plug-in
0147 DEFINE_FWK_MODULE(MuonME0SegHarvestor);