CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_5_3_10_patch2/src/DQMOffline/PFTau/src/MatchCandidateBenchmark.cc

Go to the documentation of this file.
00001 #include "DQMOffline/PFTau/interface/MatchCandidateBenchmark.h"
00002 
00003 #include "DataFormats/Candidate/interface/Candidate.h"
00004 
00005 
00006 // #include "DQMServices/Core/interface/MonitorElement.h"
00007 // #include "DQMServices/Core/interface/DQMStore.h"
00008 
00009 #include <TROOT.h>
00010 #include <TFile.h>
00011 #include <TH1.h>
00012 #include <TH2.h>
00013 
00014 
00015 using namespace std;
00016 
00017 MatchCandidateBenchmark::MatchCandidateBenchmark(Mode mode)  : Benchmark(mode) {
00018   delta_et_Over_et_VS_et_ = 0; 
00019   delta_et_VS_et_         = 0; 
00020   delta_eta_VS_et_        = 0; 
00021   delta_phi_VS_et_        = 0;
00022 
00023   histogramBooked_ = false;
00024 }
00025 
00026 MatchCandidateBenchmark::~MatchCandidateBenchmark() {}
00027 
00028 
00029 void MatchCandidateBenchmark::setup() {
00030   if (!histogramBooked_) {
00031     PhaseSpace ptPS;
00032     PhaseSpace dptOvptPS;
00033     PhaseSpace dptPS;
00034     PhaseSpace detaPS;
00035     PhaseSpace dphiPS;
00036     switch(mode_) {
00037     case VALIDATION:
00038       ptPS = PhaseSpace(100,0,1000);
00039       dptOvptPS = PhaseSpace( 200, -1, 1);
00040       dphiPS = PhaseSpace( 200, -1, 1);
00041       detaPS = PhaseSpace( 200, -1, 1);
00042       dptPS = PhaseSpace( 100, -100, 100);
00043       break;
00044     case DQMOFFLINE:
00045     default:
00046       ptPS = PhaseSpace(50,0,100);
00047       dptOvptPS = PhaseSpace( 50, -1, 1);
00048       dphiPS = PhaseSpace( 50, -1, 1);
00049       detaPS = PhaseSpace( 50, -1, 1);
00050       dptPS = PhaseSpace( 50, -50, 50);
00051       break;
00052     }
00053     float ptBins[11] = {0, 1, 2, 5, 10, 20, 50, 100, 200, 400, 1000};
00054     
00055     delta_et_Over_et_VS_et_ = book2D("delta_et_Over_et_VS_et_", 
00056                                      ";E_{T, true} (GeV);#DeltaE_{T}/E_{T}",
00057                                      10, ptBins, 
00058                                      dptOvptPS.n, dptOvptPS.m, dptOvptPS.M );
00059     
00060     
00061     delta_et_VS_et_ = book2D("delta_et_VS_et_", 
00062                              ";E_{T, true} (GeV);#DeltaE_{T}",
00063                              10, ptBins,
00064                              dptPS.n, dptPS.m, dptPS.M );
00065     
00066     delta_eta_VS_et_ = book2D("delta_eta_VS_et_", 
00067                               ";#E_{T, true} (GeV);#Delta#eta",
00068                               10, ptBins,
00069                               detaPS.n, detaPS.m, detaPS.M );
00070     
00071     delta_phi_VS_et_ = book2D("delta_phi_VS_et_", 
00072                               ";E_{T, true} (GeV);#Delta#phi",
00073                               10, ptBins,
00074                               dphiPS.n, dphiPS.m, dphiPS.M );
00075     histogramBooked_ = true;
00076   } 
00077 }
00078 void MatchCandidateBenchmark::setup(const edm::ParameterSet& parameterSet) {
00079 
00080   if (!histogramBooked_) {
00081     
00082     edm::ParameterSet dptPS = parameterSet.getParameter<edm::ParameterSet>("DeltaPtHistoParameter");
00083     edm::ParameterSet dptOvptPS = parameterSet.getParameter<edm::ParameterSet>("DeltaPtOvPtHistoParameter");
00084     edm::ParameterSet detaPS = parameterSet.getParameter<edm::ParameterSet>("DeltaEtaHistoParameter");
00085     edm::ParameterSet dphiPS = parameterSet.getParameter<edm::ParameterSet>("DeltaPhiHistoParameter");
00086     
00087     std::vector<double> ptBinsPS = parameterSet.getParameter< std::vector<double> >( "VariablePtBins" );
00088     float* ptBins = new float[ptBinsPS.size()];
00089     for (size_t i = 0; i < ptBinsPS.size(); i++) {
00090       ptBins[i] = ptBinsPS[i];
00091     }
00092     
00093     if (dptOvptPS.getParameter<bool>("switchOn")) {
00094       delta_et_Over_et_VS_et_ = book2D("delta_et_Over_et_VS_et_", 
00095                                        ";E_{T, true} (GeV);#DeltaE_{T}/E_{T}",
00096                                        ptBinsPS.size()-1, ptBins, 
00097                                        dptOvptPS.getParameter<int32_t>("nBin"), 
00098                                        dptOvptPS.getParameter<double>("xMin"), 
00099                                        dptOvptPS.getParameter<double>("xMax"));
00100     }
00101     
00102     if (dptPS.getParameter<bool>("switchOn")) {
00103       delta_et_VS_et_ = book2D("delta_et_VS_et_", 
00104                                ";E_{T, true} (GeV);#DeltaE_{T}",
00105                                ptBinsPS.size()-1, ptBins,
00106                                dptPS.getParameter<int32_t>("nBin"), 
00107                                dptPS.getParameter<double>("xMin"), 
00108                                dptPS.getParameter<double>("xMax"));
00109     }
00110     
00111     if (detaPS.getParameter<bool>("switchOn")) {
00112       delta_eta_VS_et_ = book2D("delta_eta_VS_et_", 
00113                                 ";#E_{T, true} (GeV);#Delta#eta",
00114                                 ptBinsPS.size()-1, ptBins,
00115                                 detaPS.getParameter<int32_t>("nBin"), 
00116                                 detaPS.getParameter<double>("xMin"), 
00117                                 detaPS.getParameter<double>("xMax"));
00118     }
00119     
00120     if (dphiPS.getParameter<bool>("switchOn")) {
00121       delta_phi_VS_et_ = book2D("delta_phi_VS_et_", 
00122                                 ";E_{T, true} (GeV);#Delta#phi",
00123                                 ptBinsPS.size()-1, ptBins,
00124                                 dphiPS.getParameter<int32_t>("nBin"), 
00125                                 dphiPS.getParameter<double>("xMin"),
00126                                 dphiPS.getParameter<double>("xMax"));
00127     }
00128     histogramBooked_ = true;
00129     delete ptBins;
00130   }
00131 }
00132 
00133 void MatchCandidateBenchmark::fillOne(const reco::Candidate& cand,
00134                                       const reco::Candidate& matchedCand) {
00135   
00136   if( !isInRange(cand.pt(), cand.eta(), cand.phi() ) ) return;
00137 
00138   if (histogramBooked_) {
00139     if (delta_et_Over_et_VS_et_) delta_et_Over_et_VS_et_->Fill( matchedCand.pt(), (cand.pt() - matchedCand.pt())/matchedCand.pt() );
00140     if (delta_et_VS_et_) delta_et_VS_et_->Fill( matchedCand.pt(), cand.pt() - matchedCand.pt() );
00141     if (delta_eta_VS_et_) delta_eta_VS_et_->Fill( matchedCand.pt(), cand.eta() - matchedCand.eta() );
00142     if (delta_phi_VS_et_) delta_phi_VS_et_->Fill( matchedCand.pt(), cand.phi() - matchedCand.phi() );
00143   }  
00144 
00145 }