Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-10-25 10:07:00

0001 #include "Validation/RecoEgamma/plugins/EgammaObjects.h"
0002 
0003 #include "FWCore/Utilities/interface/Exception.h"
0004 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0005 
0006 #include "DataFormats/Common/interface/Handle.h"
0007 #include "DataFormats/EgammaReco/interface/BasicCluster.h"
0008 #include "DataFormats/EgammaReco/interface/SuperCluster.h"
0009 #include "DataFormats/EgammaReco/interface/ClusterShape.h"
0010 #include "DataFormats/EgammaCandidates/interface/Photon.h"
0011 #include "DataFormats/EgammaCandidates/interface/PhotonFwd.h"
0012 #include "DataFormats/EgammaCandidates/interface/GsfElectron.h"
0013 #include "DataFormats/EgammaCandidates/interface/GsfElectronFwd.h"
0014 
0015 #include "SimDataFormats/GeneratorProducts/interface/HepMCProduct.h"
0016 
0017 #include "TF1.h"
0018 
0019 EgammaObjects::EgammaObjects(const edm::ParameterSet& ps) {
0020   particleID = ps.getParameter<int>("particleID");
0021   EtCut = ps.getParameter<int>("EtCut");
0022 
0023   if (particleID == 22)
0024     particleString = "Photon";
0025   else if (particleID == 11)
0026     particleString = "Electron";
0027   else
0028     throw(std::runtime_error(
0029         "\n\nEgammaObjects: Only 11 = Photon and 22 = Electron are acceptable parictleIDs! Exiting...\n\n"));
0030 
0031   loadCMSSWObjects(ps);
0032   loadHistoParameters(ps);
0033 
0034   rootFile_ = TFile::Open(ps.getParameter<std::string>("outputFile").c_str(), "RECREATE");
0035 
0036   hist_EtaEfficiency_ = nullptr;
0037   hist_EtaNumRecoOverNumTrue_ = nullptr;
0038   hist_deltaEtaVsEt_ = nullptr;
0039   hist_deltaEtaVsE_ = nullptr;
0040   hist_deltaEtaVsEta_ = nullptr;
0041   hist_deltaEtaVsPhi_ = nullptr;
0042   hist_resolutionEtaVsEt_ = nullptr;
0043   hist_resolutionEtaVsE_ = nullptr;
0044   hist_resolutionEtaVsEta_ = nullptr;
0045   hist_resolutionEtaVsPhi_ = nullptr;
0046 
0047   hist_Phi_ = nullptr;
0048   hist_PhiOverTruth_ = nullptr;
0049   hist_PhiEfficiency_ = nullptr;
0050   hist_PhiNumRecoOverNumTrue_ = nullptr;
0051   hist_deltaPhiVsEt_ = nullptr;
0052   hist_deltaPhiVsE_ = nullptr;
0053   hist_deltaPhiVsEta_ = nullptr;
0054   hist_deltaPhiVsPhi_ = nullptr;
0055   hist_resolutionPhiVsEt_ = nullptr;
0056   hist_resolutionPhiVsE_ = nullptr;
0057   hist_resolutionPhiVsEta_ = nullptr;
0058   hist_resolutionPhiVsPhi_ = nullptr;
0059 
0060   hist_All_recoMass_ = nullptr;
0061   hist_BarrelOnly_recoMass_ = nullptr;
0062   hist_EndcapOnly_recoMass_ = nullptr;
0063   hist_Mixed_recoMass_ = nullptr;
0064 
0065   hist_recoMass_withBackgroud_NoEtCut_ = nullptr;
0066   hist_recoMass_withBackgroud_5EtCut_ = nullptr;
0067   hist_recoMass_withBackgroud_10EtCut_ = nullptr;
0068   hist_recoMass_withBackgroud_20EtCut_ = nullptr;
0069 
0070   _TEMP_scatterPlot_EtOverTruthVsEt_ = nullptr;
0071   _TEMP_scatterPlot_EtOverTruthVsE_ = nullptr;
0072   _TEMP_scatterPlot_EtOverTruthVsEta_ = nullptr;
0073   _TEMP_scatterPlot_EtOverTruthVsPhi_ = nullptr;
0074 
0075   _TEMP_scatterPlot_EOverTruthVsEt_ = nullptr;
0076   _TEMP_scatterPlot_EOverTruthVsE_ = nullptr;
0077   _TEMP_scatterPlot_EOverTruthVsEta_ = nullptr;
0078   _TEMP_scatterPlot_EOverTruthVsPhi_ = nullptr;
0079 
0080   _TEMP_scatterPlot_deltaEtaVsEt_ = nullptr;
0081   _TEMP_scatterPlot_deltaEtaVsE_ = nullptr;
0082   _TEMP_scatterPlot_deltaEtaVsEta_ = nullptr;
0083   _TEMP_scatterPlot_deltaEtaVsPhi_ = nullptr;
0084 
0085   _TEMP_scatterPlot_deltaPhiVsEt_ = nullptr;
0086   _TEMP_scatterPlot_deltaPhiVsE_ = nullptr;
0087   _TEMP_scatterPlot_deltaPhiVsEta_ = nullptr;
0088   _TEMP_scatterPlot_deltaPhiVsPhi_ = nullptr;
0089 }
0090 
0091 void EgammaObjects::loadCMSSWObjects(const edm::ParameterSet& ps) {
0092   MCTruthCollectionT_ = consumes<edm::HepMCProduct>(ps.getParameter<edm::InputTag>("MCTruthCollection"));
0093   RecoCollectionT_ = consumes<reco::GsfElectronCollection>(ps.getParameter<edm::InputTag>("RecoCollection"));
0094 }
0095 
0096 void EgammaObjects::loadHistoParameters(const edm::ParameterSet& ps) {
0097   hist_min_Et_ = ps.getParameter<double>("hist_min_Et");
0098   hist_max_Et_ = ps.getParameter<double>("hist_max_Et");
0099   hist_bins_Et_ = ps.getParameter<int>("hist_bins_Et");
0100 
0101   hist_min_E_ = ps.getParameter<double>("hist_min_E");
0102   hist_max_E_ = ps.getParameter<double>("hist_max_E");
0103   hist_bins_E_ = ps.getParameter<int>("hist_bins_E");
0104 
0105   hist_min_Eta_ = ps.getParameter<double>("hist_min_Eta");
0106   hist_max_Eta_ = ps.getParameter<double>("hist_max_Eta");
0107   hist_bins_Eta_ = ps.getParameter<int>("hist_bins_Eta");
0108 
0109   hist_min_Phi_ = ps.getParameter<double>("hist_min_Phi");
0110   hist_max_Phi_ = ps.getParameter<double>("hist_max_Phi");
0111   hist_bins_Phi_ = ps.getParameter<int>("hist_bins_Phi");
0112 
0113   hist_min_EtOverTruth_ = ps.getParameter<double>("hist_min_EtOverTruth");
0114   hist_max_EtOverTruth_ = ps.getParameter<double>("hist_max_EtOverTruth");
0115   hist_bins_EtOverTruth_ = ps.getParameter<int>("hist_bins_EtOverTruth");
0116 
0117   hist_min_EOverTruth_ = ps.getParameter<double>("hist_min_EOverTruth");
0118   hist_max_EOverTruth_ = ps.getParameter<double>("hist_max_EOverTruth");
0119   hist_bins_EOverTruth_ = ps.getParameter<int>("hist_bins_EOverTruth");
0120 
0121   hist_min_EtaOverTruth_ = ps.getParameter<double>("hist_min_EtaOverTruth");
0122   hist_max_EtaOverTruth_ = ps.getParameter<double>("hist_max_EtaOverTruth");
0123   hist_bins_EtaOverTruth_ = ps.getParameter<int>("hist_bins_EtaOverTruth");
0124 
0125   hist_min_PhiOverTruth_ = ps.getParameter<double>("hist_min_PhiOverTruth");
0126   hist_max_PhiOverTruth_ = ps.getParameter<double>("hist_max_PhiOverTruth");
0127   hist_bins_PhiOverTruth_ = ps.getParameter<int>("hist_bins_PhiOverTruth");
0128 
0129   hist_min_deltaEta_ = ps.getParameter<double>("hist_min_deltaEta");
0130   hist_max_deltaEta_ = ps.getParameter<double>("hist_max_deltaEta");
0131   hist_bins_deltaEta_ = ps.getParameter<int>("hist_bins_deltaEta");
0132 
0133   hist_min_deltaPhi_ = ps.getParameter<double>("hist_min_deltaPhi");
0134   hist_max_deltaPhi_ = ps.getParameter<double>("hist_max_deltaPhi");
0135   hist_bins_deltaPhi_ = ps.getParameter<int>("hist_bins_deltaPhi");
0136 
0137   hist_min_recoMass_ = ps.getParameter<double>("hist_min_recoMass");
0138   hist_max_recoMass_ = ps.getParameter<double>("hist_max_recoMass");
0139   hist_bins_recoMass_ = ps.getParameter<int>("hist_bins_recoMass");
0140 }
0141 
0142 EgammaObjects::~EgammaObjects() { delete rootFile_; }
0143 
0144 void EgammaObjects::beginJob() {
0145   TH1::SetDefaultSumw2(true);
0146 
0147   createBookedHistoObjects();
0148   createTempHistoObjects();
0149 }
0150 
0151 void EgammaObjects::createBookedHistoObjects() {
0152   hist_Et_ =
0153       new TH1D("hist_Et_", ("Et Distribution of " + particleString).c_str(), hist_bins_Et_, hist_min_Et_, hist_max_Et_);
0154   hist_EtOverTruth_ = new TH1D("hist_EtOverTruth_",
0155                                ("Reco Et over True Et of " + particleString).c_str(),
0156                                hist_bins_EtOverTruth_,
0157                                hist_min_EtOverTruth_,
0158                                hist_max_EtOverTruth_);
0159   hist_EtEfficiency_ = new TH1D("hist_EtEfficiency_",
0160                                 ("# of True " + particleString + " Reconstructed over # of True " + particleString +
0161                                  " VS Et of " + particleString)
0162                                     .c_str(),
0163                                 hist_bins_Et_,
0164                                 hist_min_Et_,
0165                                 hist_max_Et_);
0166   hist_EtNumRecoOverNumTrue_ = new TH1D(
0167       "hist_EtNumRecoOverNumTrue_",
0168       ("# of Reco " + particleString + " over # of True " + particleString + " VS Et of " + particleString).c_str(),
0169       hist_bins_Et_,
0170       hist_min_Et_,
0171       hist_max_Et_);
0172 
0173   hist_E_ =
0174       new TH1D("hist_E_", ("E Distribution of " + particleString).c_str(), hist_bins_E_, hist_min_E_, hist_max_E_);
0175   hist_EOverTruth_ = new TH1D("hist_EOverTruth_",
0176                               ("Reco E over True E of " + particleString).c_str(),
0177                               hist_bins_EOverTruth_,
0178                               hist_min_EOverTruth_,
0179                               hist_max_EOverTruth_);
0180   hist_EEfficiency_ = new TH1D(
0181       "hist_EEfficiency_",
0182       ("# of True " + particleString + " Reconstructed over # of True " + particleString + " VS E of " + particleString)
0183           .c_str(),
0184       hist_bins_E_,
0185       hist_min_E_,
0186       hist_max_E_);
0187   hist_ENumRecoOverNumTrue_ = new TH1D(
0188       "hist_ENumRecoOverNumTrue_",
0189       ("# of Reco " + particleString + " over # of True " + particleString + " VS E of " + particleString).c_str(),
0190       hist_bins_E_,
0191       hist_min_E_,
0192       hist_max_E_);
0193 
0194   hist_Eta_ = new TH1D(
0195       "hist_Eta_", ("Eta Distribution of " + particleString).c_str(), hist_bins_Eta_, hist_min_Eta_, hist_max_Eta_);
0196   hist_EtaOverTruth_ = new TH1D("hist_EtaOverTruth_",
0197                                 ("Reco Eta over True Eta of " + particleString).c_str(),
0198                                 hist_bins_EtaOverTruth_,
0199                                 hist_min_EtaOverTruth_,
0200                                 hist_max_EtaOverTruth_);
0201   hist_EtaEfficiency_ = new TH1D("hist_EtaEfficiency_",
0202                                  ("# of True " + particleString + " Reconstructed over # of True " + particleString +
0203                                   " VS Eta of " + particleString)
0204                                      .c_str(),
0205                                  hist_bins_Eta_,
0206                                  hist_min_Eta_,
0207                                  hist_max_Eta_);
0208   hist_EtaNumRecoOverNumTrue_ = new TH1D(
0209       "hist_EtaNumRecoOverNumTrue_",
0210       ("# of Reco " + particleString + " over # of True " + particleString + " VS Eta of " + particleString).c_str(),
0211       hist_bins_Eta_,
0212       hist_min_Eta_,
0213       hist_max_Eta_);
0214 
0215   hist_Phi_ = new TH1D(
0216       "hist_Phi_", ("Phi Distribution of " + particleString).c_str(), hist_bins_Phi_, hist_min_Phi_, hist_max_Phi_);
0217   hist_PhiOverTruth_ = new TH1D("hist_PhiOverTruth_",
0218                                 ("Reco Phi over True Phi of " + particleString).c_str(),
0219                                 hist_bins_PhiOverTruth_,
0220                                 hist_min_PhiOverTruth_,
0221                                 hist_max_PhiOverTruth_);
0222   hist_PhiEfficiency_ = new TH1D("hist_PhiEfficiency_",
0223                                  ("# of True " + particleString + " Reconstructed over # of True " + particleString +
0224                                   " VS Phi of " + particleString)
0225                                      .c_str(),
0226                                  hist_bins_Phi_,
0227                                  hist_min_Phi_,
0228                                  hist_max_Phi_);
0229   hist_PhiNumRecoOverNumTrue_ = new TH1D(
0230       "hist_PhiNumRecoOverNumTrue_",
0231       ("# of Reco " + particleString + " over # of True " + particleString + " VS Phi of " + particleString).c_str(),
0232       hist_bins_Phi_,
0233       hist_min_Phi_,
0234       hist_max_Phi_);
0235 
0236   std::string recoParticleName;
0237 
0238   if (particleID == 22)
0239     recoParticleName = "Higgs";
0240   else if (particleID == 11)
0241     recoParticleName = "Z";
0242 
0243   hist_All_recoMass_ = new TH1D("hist_All_recoMass_",
0244                                 (recoParticleName + " Mass from " + particleString + " in All Regions").c_str(),
0245                                 hist_bins_recoMass_,
0246                                 hist_min_recoMass_,
0247                                 hist_max_recoMass_);
0248   hist_BarrelOnly_recoMass_ = new TH1D("hist_BarrelOnly_recoMass_",
0249                                        (recoParticleName + " Mass from " + particleString + " in Barrel").c_str(),
0250                                        hist_bins_recoMass_,
0251                                        hist_min_recoMass_,
0252                                        hist_max_recoMass_);
0253   hist_EndcapOnly_recoMass_ = new TH1D("hist_EndcapOnly_recoMass_",
0254                                        (recoParticleName + " Mass from " + particleString + " in EndCap").c_str(),
0255                                        hist_bins_recoMass_,
0256                                        hist_min_recoMass_,
0257                                        hist_max_recoMass_);
0258   hist_Mixed_recoMass_ = new TH1D("hist_Mixed_recoMass_",
0259                                   (recoParticleName + " Mass from " + particleString + " in Split Detectors").c_str(),
0260                                   hist_bins_recoMass_,
0261                                   hist_min_recoMass_,
0262                                   hist_max_recoMass_);
0263 
0264   hist_recoMass_withBackgroud_NoEtCut_ =
0265       new TH1D("hist_recoMass_withBackgroud_NoEtCut_",
0266                (recoParticleName + " Mass from " + particleString + " with Background, No Et Cut").c_str(),
0267                hist_bins_recoMass_,
0268                hist_min_recoMass_,
0269                hist_max_recoMass_);
0270   hist_recoMass_withBackgroud_5EtCut_ =
0271       new TH1D("hist_recoMass_withBackgroud_5EtCut_",
0272                (recoParticleName + " Mass from " + particleString + " with Background, 5 Et Cut").c_str(),
0273                hist_bins_recoMass_,
0274                hist_min_recoMass_,
0275                hist_max_recoMass_);
0276   hist_recoMass_withBackgroud_10EtCut_ =
0277       new TH1D("hist_recoMass_withBackgroud_10EtCut_",
0278                (recoParticleName + " Mass from " + particleString + " with Background, 10 Et Cut").c_str(),
0279                hist_bins_recoMass_,
0280                hist_min_recoMass_,
0281                hist_max_recoMass_);
0282   hist_recoMass_withBackgroud_20EtCut_ =
0283       new TH1D("hist_recoMass_withBackgroud_20EtCut_",
0284                (recoParticleName + " Mass from " + particleString + " with Background, 20 Et Cut").c_str(),
0285                hist_bins_recoMass_,
0286                hist_min_recoMass_,
0287                hist_max_recoMass_);
0288 }
0289 
0290 void EgammaObjects::createTempHistoObjects() {
0291   _TEMP_scatterPlot_EtOverTruthVsEt_ = new TH2D("_TEMP_scatterPlot_EtOverTruthVsEt_",
0292                                                 "_TEMP_scatterPlot_EtOverTruthVsEt_",
0293                                                 hist_bins_Et_,
0294                                                 hist_min_Et_,
0295                                                 hist_max_Et_,
0296                                                 hist_bins_EtOverTruth_,
0297                                                 hist_min_EtOverTruth_,
0298                                                 hist_max_EtOverTruth_);
0299   _TEMP_scatterPlot_EtOverTruthVsE_ = new TH2D("_TEMP_scatterPlot_EtOverTruthVsE_",
0300                                                "_TEMP_scatterPlot_EtOverTruthVsE_",
0301                                                hist_bins_E_,
0302                                                hist_min_E_,
0303                                                hist_max_E_,
0304                                                hist_bins_EtOverTruth_,
0305                                                hist_min_EtOverTruth_,
0306                                                hist_max_EtOverTruth_);
0307   _TEMP_scatterPlot_EtOverTruthVsEta_ = new TH2D("_TEMP_scatterPlot_EtOverTruthVsEta_",
0308                                                  "_TEMP_scatterPlot_EtOverTruthVsEta_",
0309                                                  hist_bins_Eta_,
0310                                                  hist_min_Eta_,
0311                                                  hist_max_Eta_,
0312                                                  hist_bins_EtOverTruth_,
0313                                                  hist_min_EtOverTruth_,
0314                                                  hist_max_EtOverTruth_);
0315   _TEMP_scatterPlot_EtOverTruthVsPhi_ = new TH2D("_TEMP_scatterPlot_EtOverTruthVsPhi_",
0316                                                  "_TEMP_scatterPlot_EtOverTruthVsPhi_",
0317                                                  hist_bins_Phi_,
0318                                                  hist_min_Phi_,
0319                                                  hist_max_Phi_,
0320                                                  hist_bins_EtOverTruth_,
0321                                                  hist_min_EtOverTruth_,
0322                                                  hist_max_EtOverTruth_);
0323 
0324   _TEMP_scatterPlot_EOverTruthVsEt_ = new TH2D("_TEMP_scatterPlot_EOverTruthVsEt_",
0325                                                "_TEMP_scatterPlot_EOverTruthVsEt_",
0326                                                hist_bins_Et_,
0327                                                hist_min_Et_,
0328                                                hist_max_Et_,
0329                                                hist_bins_EOverTruth_,
0330                                                hist_min_EOverTruth_,
0331                                                hist_max_EOverTruth_);
0332   _TEMP_scatterPlot_EOverTruthVsE_ = new TH2D("_TEMP_scatterPlot_EOverTruthVsE_",
0333                                               "_TEMP_scatterPlot_EOverTruthVsE_",
0334                                               hist_bins_E_,
0335                                               hist_min_E_,
0336                                               hist_max_E_,
0337                                               hist_bins_EOverTruth_,
0338                                               hist_min_EOverTruth_,
0339                                               hist_max_EOverTruth_);
0340   _TEMP_scatterPlot_EOverTruthVsEta_ = new TH2D("_TEMP_scatterPlot_EOverTruthVsEta_",
0341                                                 "_TEMP_scatterPlot_EOverTruthVsEta_",
0342                                                 hist_bins_Eta_,
0343                                                 hist_min_Eta_,
0344                                                 hist_max_Eta_,
0345                                                 hist_bins_EOverTruth_,
0346                                                 hist_min_EOverTruth_,
0347                                                 hist_max_EOverTruth_);
0348   _TEMP_scatterPlot_EOverTruthVsPhi_ = new TH2D("_TEMP_scatterPlot_EOverTruthVsPhi_",
0349                                                 "_TEMP_scatterPlot_EOverTruthVsPhi_",
0350                                                 hist_bins_Phi_,
0351                                                 hist_min_Phi_,
0352                                                 hist_max_Phi_,
0353                                                 hist_bins_EOverTruth_,
0354                                                 hist_min_EOverTruth_,
0355                                                 hist_max_EOverTruth_);
0356 
0357   _TEMP_scatterPlot_deltaEtaVsEt_ = new TH2D("_TEMP_scatterPlot_deltaEtaVsEt_",
0358                                              "_TEMP_scatterPlot_deltaEtaVsEt_",
0359                                              hist_bins_Et_,
0360                                              hist_min_Et_,
0361                                              hist_max_Et_,
0362                                              hist_bins_deltaEta_,
0363                                              hist_min_deltaEta_,
0364                                              hist_max_deltaEta_);
0365   _TEMP_scatterPlot_deltaEtaVsE_ = new TH2D("_TEMP_scatterPlot_deltaEtaVsE_",
0366                                             "_TEMP_scatterPlot_deltaEtaVsE_",
0367                                             hist_bins_E_,
0368                                             hist_min_E_,
0369                                             hist_max_E_,
0370                                             hist_bins_deltaEta_,
0371                                             hist_min_deltaEta_,
0372                                             hist_max_deltaEta_);
0373   _TEMP_scatterPlot_deltaEtaVsEta_ = new TH2D("_TEMP_scatterPlot_deltaEtaVsEta_",
0374                                               "_TEMP_scatterPlot_deltaEtaVsEta_",
0375                                               hist_bins_Eta_,
0376                                               hist_min_Eta_,
0377                                               hist_max_Eta_,
0378                                               hist_bins_deltaEta_,
0379                                               hist_min_deltaEta_,
0380                                               hist_max_deltaEta_);
0381   _TEMP_scatterPlot_deltaEtaVsPhi_ = new TH2D("_TEMP_scatterPlot_deltaEtaVsPhi_",
0382                                               "_TEMP_scatterPlot_deltaEtaVsPhi_",
0383                                               hist_bins_Phi_,
0384                                               hist_min_Phi_,
0385                                               hist_max_Phi_,
0386                                               hist_bins_deltaEta_,
0387                                               hist_min_deltaEta_,
0388                                               hist_max_deltaEta_);
0389 
0390   _TEMP_scatterPlot_deltaPhiVsEt_ = new TH2D("_TEMP_scatterPlot_deltaPhiVsEt_",
0391                                              "_TEMP_scatterPlot_deltaPhiVsEt_",
0392                                              hist_bins_Et_,
0393                                              hist_min_Et_,
0394                                              hist_max_Et_,
0395                                              hist_bins_deltaPhi_,
0396                                              hist_min_deltaPhi_,
0397                                              hist_max_deltaPhi_);
0398   _TEMP_scatterPlot_deltaPhiVsE_ = new TH2D("_TEMP_scatterPlot_deltaPhiVsE_",
0399                                             "_TEMP_scatterPlot_deltaPhiVsE_",
0400                                             hist_bins_E_,
0401                                             hist_min_E_,
0402                                             hist_max_E_,
0403                                             hist_bins_deltaPhi_,
0404                                             hist_min_deltaPhi_,
0405                                             hist_max_deltaPhi_);
0406   _TEMP_scatterPlot_deltaPhiVsEta_ = new TH2D("_TEMP_scatterPlot_deltaPhiVsEta_",
0407                                               "_TEMP_scatterPlot_deltaPhiVsEta_",
0408                                               hist_bins_Eta_,
0409                                               hist_min_Eta_,
0410                                               hist_max_Eta_,
0411                                               hist_bins_deltaPhi_,
0412                                               hist_min_deltaPhi_,
0413                                               hist_max_deltaPhi_);
0414   _TEMP_scatterPlot_deltaPhiVsPhi_ = new TH2D("_TEMP_scatterPlot_deltaPhiVsPhi_",
0415                                               "_TEMP_scatterPlot_deltaPhiVsPhi_",
0416                                               hist_bins_Phi_,
0417                                               hist_min_Phi_,
0418                                               hist_max_Phi_,
0419                                               hist_bins_deltaPhi_,
0420                                               hist_min_deltaPhi_,
0421                                               hist_max_deltaPhi_);
0422 }
0423 
0424 void EgammaObjects::analyze(const edm::Event& evt, const edm::EventSetup& es) {
0425   if (particleID == 22)
0426     analyzePhotons(evt, es);
0427   else if (particleID == 11)
0428     analyzeElectrons(evt, es);
0429 }
0430 
0431 void EgammaObjects::analyzePhotons(const edm::Event& evt, const edm::EventSetup& es) {
0432   edm::Handle<reco::PhotonCollection> pPhotons;
0433   evt.getByToken(RecoCollectionT_, pPhotons);
0434   if (!pPhotons.isValid()) {
0435     Labels l;
0436     labelsForToken(RecoCollectionT_, l);
0437     edm::LogError("EgammaObjects") << "Error! can't get collection with label " << l.module;
0438   }
0439 
0440   const reco::PhotonCollection* photons = pPhotons.product();
0441   std::vector<reco::Photon> photonsMCMatched;
0442 
0443   for (reco::PhotonCollection::const_iterator aClus = photons->begin(); aClus != photons->end(); aClus++) {
0444     if (aClus->et() >= EtCut) {
0445       hist_Et_->Fill(aClus->et());
0446       hist_E_->Fill(aClus->energy());
0447       hist_Eta_->Fill(aClus->eta());
0448       hist_Phi_->Fill(aClus->phi());
0449     }
0450   }
0451 
0452   for (int firstPhoton = 0, numPhotons = photons->size(); firstPhoton < numPhotons - 1; firstPhoton++)
0453     for (int secondPhoton = firstPhoton + 1; secondPhoton < numPhotons; secondPhoton++) {
0454       reco::Photon pOne = (*photons)[firstPhoton];
0455       reco::Photon pTwo = (*photons)[secondPhoton];
0456 
0457       double recoMass = findRecoMass(pOne, pTwo);
0458 
0459       hist_recoMass_withBackgroud_NoEtCut_->Fill(recoMass);
0460 
0461       if (pOne.et() >= 5 && pTwo.et() >= 5)
0462         hist_recoMass_withBackgroud_5EtCut_->Fill(recoMass);
0463 
0464       if (pOne.et() >= 10 && pTwo.et() >= 10)
0465         hist_recoMass_withBackgroud_10EtCut_->Fill(recoMass);
0466 
0467       if (pOne.et() >= 20 && pTwo.et() >= 20)
0468         hist_recoMass_withBackgroud_20EtCut_->Fill(recoMass);
0469     }
0470 
0471   edm::Handle<edm::HepMCProduct> pMCTruth;
0472   evt.getByToken(MCTruthCollectionT_, pMCTruth);
0473   if (!pMCTruth.isValid()) {
0474     Labels l;
0475     labelsForToken(MCTruthCollectionT_, l);
0476     edm::LogError("EgammaObjects") << "Error! can't get collection with label " << l.module;
0477   }
0478 
0479   const HepMC::GenEvent* genEvent = pMCTruth->GetEvent();
0480 
0481   for (HepMC::GenEvent::particle_const_iterator currentParticle = genEvent->particles_begin();
0482        currentParticle != genEvent->particles_end();
0483        currentParticle++) {
0484     if (abs((*currentParticle)->pdg_id()) == 22 && (*currentParticle)->status() == 1 &&
0485         (*currentParticle)->momentum().e() /
0486                 cosh(ecalEta((*currentParticle)->momentum().eta(),
0487                              (*currentParticle)->production_vertex()->position().z() / 10.,
0488                              (*currentParticle)->production_vertex()->position().perp() / 10.)) >=
0489             EtCut) {
0490       HepMC::FourVector vtx = (*currentParticle)->production_vertex()->position();
0491       double phiTrue = (*currentParticle)->momentum().phi();
0492       double etaTrue = ecalEta((*currentParticle)->momentum().eta(), vtx.z() / 10., vtx.perp() / 10.);
0493       double eTrue = (*currentParticle)->momentum().e();
0494       double etTrue = (*currentParticle)->momentum().e() / cosh(etaTrue);
0495 
0496       double etaCurrent, etaFound = -999;
0497       double phiCurrent, phiFound = -999;
0498       double etCurrent, etFound = -999;
0499       double eCurrent, eFound = -999;
0500 
0501       reco::Photon bestMatchPhoton;
0502 
0503       double closestParticleDistance = 999;
0504 
0505       for (reco::PhotonCollection::const_iterator aClus = photons->begin(); aClus != photons->end(); aClus++) {
0506         if (aClus->et() > EtCut) {
0507           etaCurrent = aClus->eta();
0508           phiCurrent = aClus->phi();
0509           etCurrent = aClus->et();
0510           eCurrent = aClus->energy();
0511 
0512           double deltaPhi = phiCurrent - phiTrue;
0513           if (deltaPhi > Geom::pi())
0514             deltaPhi -= 2. * Geom::pi();
0515           if (deltaPhi < -Geom::pi())
0516             deltaPhi += 2. * Geom::pi();
0517           double deltaR = std::sqrt(std::pow(etaCurrent - etaTrue, 2) + std::pow(deltaPhi, 2));
0518 
0519           if (deltaR < closestParticleDistance) {
0520             etFound = etCurrent;
0521             eFound = eCurrent;
0522             etaFound = etaCurrent;
0523             phiFound = phiCurrent;
0524             closestParticleDistance = deltaR;
0525             bestMatchPhoton = *aClus;
0526           }
0527         }
0528       }
0529 
0530       if (closestParticleDistance < 0.05 && etFound / etTrue > .5 && etFound / etTrue < 1.5) {
0531         hist_EtOverTruth_->Fill(etFound / etTrue);
0532         hist_EOverTruth_->Fill(eFound / eTrue);
0533         hist_EtaOverTruth_->Fill(etaFound / etaTrue);
0534         hist_PhiOverTruth_->Fill(phiFound / phiTrue);
0535 
0536         hist_EtEfficiency_->Fill(etTrue);
0537         hist_EEfficiency_->Fill(eTrue);
0538         hist_EtaEfficiency_->Fill(etaTrue);
0539         hist_PhiEfficiency_->Fill(phiTrue);
0540 
0541         double deltaPhi = phiFound - phiTrue;
0542         if (deltaPhi > Geom::pi())
0543           deltaPhi -= 2. * Geom::pi();
0544         if (deltaPhi < -Geom::pi())
0545           deltaPhi += 2. * Geom::pi();
0546 
0547         _TEMP_scatterPlot_EtOverTruthVsEt_->Fill(etFound, etFound / etTrue);
0548         _TEMP_scatterPlot_EtOverTruthVsE_->Fill(eFound, etFound / etTrue);
0549         _TEMP_scatterPlot_EtOverTruthVsEta_->Fill(etaFound, etFound / etTrue);
0550         _TEMP_scatterPlot_EtOverTruthVsPhi_->Fill(phiFound, etFound / etTrue);
0551 
0552         _TEMP_scatterPlot_EOverTruthVsEt_->Fill(etFound, eFound / eTrue);
0553         _TEMP_scatterPlot_EOverTruthVsE_->Fill(eFound, eFound / eTrue);
0554         _TEMP_scatterPlot_EOverTruthVsEta_->Fill(etaFound, eFound / eTrue);
0555         _TEMP_scatterPlot_EOverTruthVsPhi_->Fill(phiFound, eFound / eTrue);
0556 
0557         _TEMP_scatterPlot_deltaEtaVsEt_->Fill(etFound, etaFound - etaTrue);
0558         _TEMP_scatterPlot_deltaEtaVsE_->Fill(eFound, etaFound - etaTrue);
0559         _TEMP_scatterPlot_deltaEtaVsEta_->Fill(etaFound, etaFound - etaTrue);
0560         _TEMP_scatterPlot_deltaEtaVsPhi_->Fill(phiFound, etaFound - etaTrue);
0561 
0562         _TEMP_scatterPlot_deltaPhiVsEt_->Fill(etFound, deltaPhi);
0563         _TEMP_scatterPlot_deltaPhiVsE_->Fill(eFound, deltaPhi);
0564         _TEMP_scatterPlot_deltaPhiVsEta_->Fill(etaFound, deltaPhi);
0565         _TEMP_scatterPlot_deltaPhiVsPhi_->Fill(phiFound, deltaPhi);
0566 
0567         photonsMCMatched.push_back(bestMatchPhoton);
0568       }
0569 
0570       hist_EtNumRecoOverNumTrue_->Fill(etTrue);
0571       hist_ENumRecoOverNumTrue_->Fill(eTrue);
0572       hist_EtaNumRecoOverNumTrue_->Fill(etaTrue);
0573       hist_PhiNumRecoOverNumTrue_->Fill(phiTrue);
0574     }
0575   }
0576 
0577   if (photonsMCMatched.size() == 2) {
0578     reco::Photon pOne = photonsMCMatched[0];
0579     reco::Photon pTwo = photonsMCMatched[1];
0580 
0581     double recoMass = findRecoMass(pOne, pTwo);
0582 
0583     hist_All_recoMass_->Fill(recoMass);
0584 
0585     if (pOne.superCluster()->seed()->algo() == 1 && pTwo.superCluster()->seed()->algo() == 1)
0586       hist_BarrelOnly_recoMass_->Fill(recoMass);
0587     else if (pOne.superCluster()->seed()->algo() == 0 && pTwo.superCluster()->seed()->algo() == 0)
0588       hist_EndcapOnly_recoMass_->Fill(recoMass);
0589     else
0590       hist_Mixed_recoMass_->Fill(recoMass);
0591   }
0592 }
0593 
0594 void EgammaObjects::analyzeElectrons(const edm::Event& evt, const edm::EventSetup& es) {
0595   edm::Handle<reco::GsfElectronCollection> pElectrons;
0596   evt.getByToken(RecoCollectionT_, pElectrons);
0597   if (!pElectrons.isValid()) {
0598     Labels l;
0599     labelsForToken(RecoCollectionT_, l);
0600     edm::LogError("DOEPlotsProducerElectrons") << "Error! can't get collection with label " << l.module;
0601   }
0602 
0603   const reco::GsfElectronCollection* electrons = pElectrons.product();
0604   std::vector<reco::GsfElectron> electronsMCMatched;
0605 
0606   for (reco::GsfElectronCollection::const_iterator aClus = electrons->begin(); aClus != electrons->end(); aClus++) {
0607     if (aClus->et() >= EtCut) {
0608       hist_Et_->Fill(aClus->et());
0609       hist_E_->Fill(aClus->energy());
0610       hist_Eta_->Fill(aClus->eta());
0611       hist_Phi_->Fill(aClus->phi());
0612     }
0613   }
0614 
0615   for (int firstElectron = 0, numElectrons = electrons->size(); firstElectron < numElectrons - 1; firstElectron++)
0616     for (int secondElectron = firstElectron + 1; secondElectron < numElectrons; secondElectron++) {
0617       reco::GsfElectron eOne = (*electrons)[firstElectron];
0618       reco::GsfElectron eTwo = (*electrons)[secondElectron];
0619 
0620       double recoMass = findRecoMass(eOne, eTwo);
0621 
0622       hist_recoMass_withBackgroud_NoEtCut_->Fill(recoMass);
0623 
0624       if (eOne.et() >= 5 && eTwo.et() >= 5)
0625         hist_recoMass_withBackgroud_5EtCut_->Fill(recoMass);
0626 
0627       if (eOne.et() >= 10 && eTwo.et() >= 10)
0628         hist_recoMass_withBackgroud_10EtCut_->Fill(recoMass);
0629 
0630       if (eOne.et() >= 20 && eTwo.et() >= 20)
0631         hist_recoMass_withBackgroud_20EtCut_->Fill(recoMass);
0632     }
0633 
0634   edm::Handle<edm::HepMCProduct> pMCTruth;
0635   evt.getByToken(MCTruthCollectionT_, pMCTruth);
0636   if (!pMCTruth.isValid()) {
0637     Labels l;
0638     labelsForToken(MCTruthCollectionT_, l);
0639     edm::LogError("DOEPlotsProducerElectrons") << "Error! can't get collection with label " << l.module;
0640   }
0641 
0642   const HepMC::GenEvent* genEvent = pMCTruth->GetEvent();
0643   for (HepMC::GenEvent::particle_const_iterator currentParticle = genEvent->particles_begin();
0644        currentParticle != genEvent->particles_end();
0645        currentParticle++) {
0646     if (abs((*currentParticle)->pdg_id()) == 11 && (*currentParticle)->status() == 1 &&
0647         (*currentParticle)->momentum().e() / cosh((*currentParticle)->momentum().eta()) >= EtCut) {
0648       double phiTrue = (*currentParticle)->momentum().phi();
0649       double etaTrue = (*currentParticle)->momentum().eta();
0650       double eTrue = (*currentParticle)->momentum().e();
0651       double etTrue = (*currentParticle)->momentum().e() / cosh(etaTrue);
0652 
0653       double etaCurrent, etaFound = -999;
0654       double phiCurrent, phiFound = -999;
0655       double etCurrent, etFound = -999;
0656       double eCurrent, eFound = -999;
0657 
0658       reco::GsfElectron bestMatchElectron;
0659 
0660       double closestParticleDistance = 999;
0661 
0662       for (reco::GsfElectronCollection::const_iterator aClus = electrons->begin(); aClus != electrons->end(); aClus++) {
0663         if (aClus->et() > EtCut) {
0664           etaCurrent = aClus->eta();
0665           phiCurrent = aClus->phi();
0666           etCurrent = aClus->et();
0667           eCurrent = aClus->energy();
0668 
0669           double deltaPhi = phiCurrent - phiTrue;
0670           if (deltaPhi > Geom::pi())
0671             deltaPhi -= 2. * Geom::pi();
0672           if (deltaPhi < -Geom::pi())
0673             deltaPhi += 2. * Geom::pi();
0674           double deltaR = std::sqrt(std::pow(etaCurrent - etaTrue, 2) + std::pow(deltaPhi, 2));
0675 
0676           if (deltaR < closestParticleDistance) {
0677             etFound = etCurrent;
0678             eFound = eCurrent;
0679             etaFound = etaCurrent;
0680             phiFound = phiCurrent;
0681             closestParticleDistance = deltaR;
0682             bestMatchElectron = *aClus;
0683           }
0684         }
0685       }
0686 
0687       if (closestParticleDistance < 0.05 && etFound / etTrue > .5 && etFound / etTrue < 1.5) {
0688         hist_EtOverTruth_->Fill(etFound / etTrue);
0689         hist_EOverTruth_->Fill(eFound / eTrue);
0690         hist_EtaOverTruth_->Fill(etaFound / etaTrue);
0691         hist_PhiOverTruth_->Fill(phiFound / phiTrue);
0692 
0693         hist_EtEfficiency_->Fill(etTrue);
0694         hist_EEfficiency_->Fill(eTrue);
0695         hist_EtaEfficiency_->Fill(etaTrue);
0696         hist_PhiEfficiency_->Fill(phiTrue);
0697 
0698         double deltaPhi = phiFound - phiTrue;
0699         if (deltaPhi > Geom::pi())
0700           deltaPhi -= 2. * Geom::pi();
0701         if (deltaPhi < -Geom::pi())
0702           deltaPhi += 2. * Geom::pi();
0703 
0704         _TEMP_scatterPlot_EtOverTruthVsEt_->Fill(etFound, etFound / etTrue);
0705         _TEMP_scatterPlot_EtOverTruthVsE_->Fill(eFound, etFound / etTrue);
0706         _TEMP_scatterPlot_EtOverTruthVsEta_->Fill(etaFound, etFound / etTrue);
0707         _TEMP_scatterPlot_EtOverTruthVsPhi_->Fill(phiFound, etFound / etTrue);
0708 
0709         _TEMP_scatterPlot_EOverTruthVsEt_->Fill(etFound, eFound / eTrue);
0710         _TEMP_scatterPlot_EOverTruthVsE_->Fill(eFound, eFound / eTrue);
0711         _TEMP_scatterPlot_EOverTruthVsEta_->Fill(etaFound, eFound / eTrue);
0712         _TEMP_scatterPlot_EOverTruthVsPhi_->Fill(phiFound, eFound / eTrue);
0713 
0714         _TEMP_scatterPlot_deltaEtaVsEt_->Fill(etFound, etaFound - etaTrue);
0715         _TEMP_scatterPlot_deltaEtaVsE_->Fill(eFound, etaFound - etaTrue);
0716         _TEMP_scatterPlot_deltaEtaVsEta_->Fill(etaFound, etaFound - etaTrue);
0717         _TEMP_scatterPlot_deltaEtaVsPhi_->Fill(phiFound, etaFound - etaTrue);
0718 
0719         _TEMP_scatterPlot_deltaPhiVsEt_->Fill(etFound, deltaPhi);
0720         _TEMP_scatterPlot_deltaPhiVsE_->Fill(eFound, deltaPhi);
0721         _TEMP_scatterPlot_deltaPhiVsEta_->Fill(etaFound, deltaPhi);
0722         _TEMP_scatterPlot_deltaPhiVsPhi_->Fill(phiFound, deltaPhi);
0723 
0724         electronsMCMatched.push_back(bestMatchElectron);
0725       }
0726 
0727       hist_EtNumRecoOverNumTrue_->Fill(etTrue);
0728       hist_ENumRecoOverNumTrue_->Fill(eTrue);
0729       hist_EtaNumRecoOverNumTrue_->Fill(etaTrue);
0730       hist_PhiNumRecoOverNumTrue_->Fill(phiTrue);
0731     }
0732   }
0733 
0734   if (electronsMCMatched.size() == 2) {
0735     reco::GsfElectron eOne = electronsMCMatched[0];
0736     reco::GsfElectron eTwo = electronsMCMatched[1];
0737 
0738     double recoMass = findRecoMass(eOne, eTwo);
0739 
0740     hist_All_recoMass_->Fill(recoMass);
0741 
0742     if (eOne.superCluster()->seed()->algo() == 1 && eTwo.superCluster()->seed()->algo() == 1)
0743       hist_BarrelOnly_recoMass_->Fill(recoMass);
0744     else if (eOne.superCluster()->seed()->algo() == 0 && eTwo.superCluster()->seed()->algo() == 0)
0745       hist_EndcapOnly_recoMass_->Fill(recoMass);
0746     else
0747       hist_Mixed_recoMass_->Fill(recoMass);
0748   }
0749 }
0750 
0751 double EgammaObjects::findRecoMass(const reco::Photon& pOne, const reco::Photon& pTwo) {
0752   double cosTheta = (cos(pOne.superCluster()->phi() - pTwo.superCluster()->phi()) +
0753                      sinh(pOne.superCluster()->eta()) * sinh(pTwo.superCluster()->eta())) /
0754                     (cosh(pOne.superCluster()->eta()) * cosh(pTwo.superCluster()->eta()));
0755 
0756   double recoMass = sqrt(2 * (pOne.superCluster())->energy() * (pTwo.superCluster())->energy() * (1 - cosTheta));
0757 
0758   return recoMass;
0759 }
0760 
0761 double EgammaObjects::findRecoMass(const reco::GsfElectron& eOne, const reco::GsfElectron& eTwo) {
0762   double cosTheta = (cos(eOne.caloPosition().phi() - eTwo.caloPosition().phi()) +
0763                      sinh(eOne.caloPosition().eta()) * sinh(eTwo.caloPosition().eta())) /
0764                     (cosh(eOne.caloPosition().eta()) * cosh(eTwo.caloPosition().eta()));
0765 
0766   double recoMass = sqrt(2 * eOne.caloEnergy() * eTwo.caloEnergy() * (1 - cosTheta));
0767 
0768   return recoMass;
0769 }
0770 
0771 float EgammaObjects::ecalEta(float EtaParticle, float Zvertex, float plane_Radius) {
0772   const float R_ECAL = 136.5;
0773   const float Z_Endcap = 328.0;
0774   const float etaBarrelEndcap = 1.479;
0775 
0776   if (EtaParticle != 0.) {
0777     float Theta = 0.0;
0778     float ZEcal = (R_ECAL - plane_Radius) * sinh(EtaParticle) + Zvertex;
0779 
0780     if (ZEcal != 0.0)
0781       Theta = atan(R_ECAL / ZEcal);
0782     if (Theta < 0.0)
0783       Theta = Theta + Geom::pi();
0784 
0785     float ETA = -log(tan(0.5 * Theta));
0786 
0787     if (std::abs(ETA) > etaBarrelEndcap) {
0788       float Zend = Z_Endcap;
0789       if (EtaParticle < 0.0)
0790         Zend = -Zend;
0791       float Zlen = Zend - Zvertex;
0792       float RR = Zlen / sinh(EtaParticle);
0793       Theta = atan((RR + plane_Radius) / Zend);
0794       if (Theta < 0.0)
0795         Theta = Theta + Geom::pi();
0796       ETA = -log(tan(0.5 * Theta));
0797     }
0798 
0799     return ETA;
0800   } else {
0801     edm::LogWarning("") << "[EgammaObjects::ecalEta] Warning: Eta equals to zero, not correcting";
0802     return EtaParticle;
0803   }
0804 }
0805 
0806 void EgammaObjects::endJob() {
0807   rootFile_->cd();
0808   rootFile_->mkdir(particleString.c_str());
0809 
0810   getDeltaResHistosViaSlicing();
0811   getEfficiencyHistosViaDividing();
0812   fitHistos();
0813 
0814   applyLabels();
0815   setDrawOptions();
0816   saveHistos();
0817   rootFile_->Close();
0818 }
0819 
0820 void EgammaObjects::getDeltaResHistosViaSlicing() {
0821   _TEMP_scatterPlot_EtOverTruthVsEt_->FitSlicesY(nullptr, 1, hist_bins_Et_, 10, "QRG3");
0822   _TEMP_scatterPlot_EtOverTruthVsE_->FitSlicesY(nullptr, 1, hist_bins_E_, 10, "QRG3");
0823   _TEMP_scatterPlot_EtOverTruthVsEta_->FitSlicesY(nullptr, 1, hist_bins_Eta_, 10, "QRG2");
0824   _TEMP_scatterPlot_EtOverTruthVsPhi_->FitSlicesY(nullptr, 1, hist_bins_Phi_, 10, "QRG2");
0825 
0826   _TEMP_scatterPlot_EOverTruthVsEt_->FitSlicesY(nullptr, 1, hist_bins_Et_, 10, "QRG3");
0827   _TEMP_scatterPlot_EOverTruthVsE_->FitSlicesY(nullptr, 1, hist_bins_E_, 10, "QRG3");
0828   _TEMP_scatterPlot_EOverTruthVsEta_->FitSlicesY(nullptr, 1, hist_bins_Eta_, 10, "QRG2");
0829   _TEMP_scatterPlot_EOverTruthVsPhi_->FitSlicesY(nullptr, 1, hist_bins_Phi_, 10, "QRG2");
0830 
0831   _TEMP_scatterPlot_deltaEtaVsEt_->FitSlicesY(nullptr, 1, hist_bins_Et_, 10, "QRG3");
0832   _TEMP_scatterPlot_deltaEtaVsE_->FitSlicesY(nullptr, 1, hist_bins_E_, 10, "QRG3");
0833   _TEMP_scatterPlot_deltaEtaVsEta_->FitSlicesY(nullptr, 1, hist_bins_Eta_, 10, "QRG2");
0834   _TEMP_scatterPlot_deltaEtaVsPhi_->FitSlicesY(nullptr, 1, hist_bins_Phi_, 10, "QRG2");
0835 
0836   _TEMP_scatterPlot_deltaPhiVsEt_->FitSlicesY(nullptr, 1, hist_bins_Et_, 10, "QRG3");
0837   _TEMP_scatterPlot_deltaPhiVsE_->FitSlicesY(nullptr, 1, hist_bins_E_, 10, "QRG3");
0838   _TEMP_scatterPlot_deltaPhiVsEta_->FitSlicesY(nullptr, 1, hist_bins_Eta_, 10, "QRG2");
0839   _TEMP_scatterPlot_deltaPhiVsPhi_->FitSlicesY(nullptr, 1, hist_bins_Phi_, 10, "QRG2");
0840 
0841   hist_EtOverTruthVsEt_ = (TH1D*)gDirectory->Get("_TEMP_scatterPlot_EtOverTruthVsEt__1");
0842   hist_EtOverTruthVsE_ = (TH1D*)gDirectory->Get("_TEMP_scatterPlot_EtOverTruthVsE__1");
0843   hist_EtOverTruthVsEta_ = (TH1D*)gDirectory->Get("_TEMP_scatterPlot_EtOverTruthVsEta__1");
0844   hist_EtOverTruthVsPhi_ = (TH1D*)gDirectory->Get("_TEMP_scatterPlot_EtOverTruthVsPhi__1");
0845 
0846   hist_EOverTruthVsEt_ = (TH1D*)gDirectory->Get("_TEMP_scatterPlot_EOverTruthVsEt__1");
0847   hist_EOverTruthVsE_ = (TH1D*)gDirectory->Get("_TEMP_scatterPlot_EOverTruthVsE__1");
0848   hist_EOverTruthVsEta_ = (TH1D*)gDirectory->Get("_TEMP_scatterPlot_EOverTruthVsEta__1");
0849   hist_EOverTruthVsPhi_ = (TH1D*)gDirectory->Get("_TEMP_scatterPlot_EOverTruthVsPhi__1");
0850 
0851   hist_deltaEtaVsEt_ = (TH1D*)gDirectory->Get("_TEMP_scatterPlot_deltaEtaVsEt__1");
0852   hist_deltaEtaVsE_ = (TH1D*)gDirectory->Get("_TEMP_scatterPlot_deltaEtaVsE__1");
0853   hist_deltaEtaVsEta_ = (TH1D*)gDirectory->Get("_TEMP_scatterPlot_deltaEtaVsEta__1");
0854   hist_deltaEtaVsPhi_ = (TH1D*)gDirectory->Get("_TEMP_scatterPlot_deltaEtaVsPhi__1");
0855 
0856   hist_deltaPhiVsEt_ = (TH1D*)gDirectory->Get("_TEMP_scatterPlot_deltaPhiVsEt__1");
0857   hist_deltaPhiVsE_ = (TH1D*)gDirectory->Get("_TEMP_scatterPlot_deltaPhiVsE__1");
0858   hist_deltaPhiVsEta_ = (TH1D*)gDirectory->Get("_TEMP_scatterPlot_deltaPhiVsEta__1");
0859   hist_deltaPhiVsPhi_ = (TH1D*)gDirectory->Get("_TEMP_scatterPlot_deltaPhiVsPhi__1");
0860 
0861   hist_EtOverTruthVsEt_->SetNameTitle("hist_EtOverTruthVsEt_",
0862                                       ("Reco Et over True Et VS Et of " + particleString).c_str());
0863   hist_EtOverTruthVsE_->SetNameTitle("hist_EtOverTruthVsE_",
0864                                      ("Reco Et over True Et VS E of " + particleString).c_str());
0865   hist_EtOverTruthVsEta_->SetNameTitle("hist_EtOverTruthVsEta_",
0866                                        ("Reco Et over True Et VS Eta of " + particleString).c_str());
0867   hist_EtOverTruthVsPhi_->SetNameTitle("hist_EtOverTruthVsPhi_",
0868                                        ("Reco Et over True Et VS Phi of " + particleString).c_str());
0869 
0870   hist_EOverTruthVsEt_->SetNameTitle("hist_EOverTruthVsEt_", ("Reco E over True E VS Et of " + particleString).c_str());
0871   hist_EOverTruthVsE_->SetNameTitle("hist_EOverTruthVsE_", ("Reco E over True E VS E of " + particleString).c_str());
0872   hist_EOverTruthVsEta_->SetNameTitle("hist_EOverTruthVsEta_",
0873                                       ("Reco E over True E VS Eta of " + particleString).c_str());
0874   hist_EOverTruthVsPhi_->SetNameTitle("hist_EOverTruthVsPhi_",
0875                                       ("Reco E over True E VS Phi of " + particleString).c_str());
0876 
0877   hist_deltaEtaVsEt_->SetNameTitle("hist_deltaEtaVsEt_", ("delta Eta VS Et of " + particleString).c_str());
0878   hist_deltaEtaVsE_->SetNameTitle("hist_deltaEtaVsE_", ("delta Eta VS E of " + particleString).c_str());
0879   hist_deltaEtaVsEta_->SetNameTitle("hist_deltaEtaVsEta_", ("delta Eta VS Eta of " + particleString).c_str());
0880   hist_deltaEtaVsPhi_->SetNameTitle("hist_deltaEtaVsPhi_", ("delta Eta VS Phi of " + particleString).c_str());
0881 
0882   hist_deltaPhiVsEt_->SetNameTitle("hist_deltaPhiVsEt_", ("delta Phi VS Et of " + particleString).c_str());
0883   hist_deltaPhiVsE_->SetNameTitle("hist_deltaPhiVsE_", ("delta Phi VS E of " + particleString).c_str());
0884   hist_deltaPhiVsEta_->SetNameTitle("hist_deltaPhiVsEta_", ("delta Phi VS Eta of " + particleString).c_str());
0885   hist_deltaPhiVsPhi_->SetNameTitle("hist_deltaPhiVsPhi_", ("delta Phi VS Phi of " + particleString).c_str());
0886 
0887   hist_resolutionEtVsEt_ = (TH1D*)gDirectory->Get("_TEMP_scatterPlot_EtOverTruthVsEt__2");
0888   hist_resolutionEtVsE_ = (TH1D*)gDirectory->Get("_TEMP_scatterPlot_EtOverTruthVsE__2");
0889   hist_resolutionEtVsEta_ = (TH1D*)gDirectory->Get("_TEMP_scatterPlot_EtOverTruthVsEta__2");
0890   hist_resolutionEtVsPhi_ = (TH1D*)gDirectory->Get("_TEMP_scatterPlot_EtOverTruthVsPhi__2");
0891 
0892   hist_resolutionEVsEt_ = (TH1D*)gDirectory->Get("_TEMP_scatterPlot_EOverTruthVsEt__2");
0893   hist_resolutionEVsE_ = (TH1D*)gDirectory->Get("_TEMP_scatterPlot_EOverTruthVsE__2");
0894   hist_resolutionEVsEta_ = (TH1D*)gDirectory->Get("_TEMP_scatterPlot_EOverTruthVsEta__2");
0895   hist_resolutionEVsPhi_ = (TH1D*)gDirectory->Get("_TEMP_scatterPlot_EOverTruthVsPhi__2");
0896 
0897   hist_resolutionEtaVsEt_ = (TH1D*)gDirectory->Get("_TEMP_scatterPlot_deltaEtaVsEt__2");
0898   hist_resolutionEtaVsE_ = (TH1D*)gDirectory->Get("_TEMP_scatterPlot_deltaEtaVsE__2");
0899   hist_resolutionEtaVsEta_ = (TH1D*)gDirectory->Get("_TEMP_scatterPlot_deltaEtaVsEta__2");
0900   hist_resolutionEtaVsPhi_ = (TH1D*)gDirectory->Get("_TEMP_scatterPlot_deltaEtaVsPhi__2");
0901 
0902   hist_resolutionPhiVsEt_ = (TH1D*)gDirectory->Get("_TEMP_scatterPlot_deltaPhiVsEt__2");
0903   hist_resolutionPhiVsE_ = (TH1D*)gDirectory->Get("_TEMP_scatterPlot_deltaPhiVsE__2");
0904   hist_resolutionPhiVsEta_ = (TH1D*)gDirectory->Get("_TEMP_scatterPlot_deltaPhiVsEta__2");
0905   hist_resolutionPhiVsPhi_ = (TH1D*)gDirectory->Get("_TEMP_scatterPlot_deltaPhiVsPhi__2");
0906 
0907   hist_resolutionEtVsEt_->SetNameTitle("hist_resolutionEtVsEt_",
0908                                        ("#sigma of Reco Et over True Et VS Et of " + particleString).c_str());
0909   hist_resolutionEtVsE_->SetNameTitle("hist_resolutionEtVsE_",
0910                                       ("#sigma of Reco Et over True Et VS E of " + particleString).c_str());
0911   hist_resolutionEtVsEta_->SetNameTitle("hist_resolutionEtVsEta_",
0912                                         ("#sigma of Reco Et over True Et VS Eta of " + particleString).c_str());
0913   hist_resolutionEtVsPhi_->SetNameTitle("hist_resolutionEtVsPhi_",
0914                                         ("#sigma of Reco Et over True Et VS Phi of " + particleString).c_str());
0915 
0916   hist_resolutionEVsEt_->SetNameTitle("hist_resolutionEVsEt_",
0917                                       ("#sigma of Reco E over True E VS Et of " + particleString).c_str());
0918   hist_resolutionEVsE_->SetNameTitle("hist_resolutionEVsE_",
0919                                      ("#sigma of Reco E over True E VS E of " + particleString).c_str());
0920   hist_resolutionEVsEta_->SetNameTitle("hist_resolutionEVsEta_",
0921                                        ("#sigma of Reco E over True E VS Eta of " + particleString).c_str());
0922   hist_resolutionEVsPhi_->SetNameTitle("hist_resolutionEVsPhi_",
0923                                        ("#sigma of Reco E over True E VS Phi of " + particleString).c_str());
0924 
0925   hist_resolutionEtaVsEt_->SetNameTitle("hist_resolutionEtaVsEt_",
0926                                         ("#sigma of delta Eta VS Et of " + particleString).c_str());
0927   hist_resolutionEtaVsE_->SetNameTitle("hist_resolutionEtaVsE_",
0928                                        ("#sigma of delta Eta VS E of " + particleString).c_str());
0929   hist_resolutionEtaVsEta_->SetNameTitle("hist_resolutionEtaVsEta_",
0930                                          ("#sigma of delta Eta VS Eta of " + particleString).c_str());
0931   hist_resolutionEtaVsPhi_->SetNameTitle("hist_resolutionEtaVsPhi_",
0932                                          ("#sigma of delta Eta VS Phi of " + particleString).c_str());
0933 
0934   hist_resolutionPhiVsEt_->SetNameTitle("hist_resolutionPhiVsEt_",
0935                                         ("#sigma of delta Phi VS Et of " + particleString).c_str());
0936   hist_resolutionPhiVsE_->SetNameTitle("hist_resolutionPhiVsE_",
0937                                        ("#sigma of delta Phi VS E of " + particleString).c_str());
0938   hist_resolutionPhiVsEta_->SetNameTitle("hist_resolutionPhiVsEta_",
0939                                          ("#sigma of delta Phi VS Eta of " + particleString).c_str());
0940   hist_resolutionPhiVsPhi_->SetNameTitle("hist_resolutionPhiVsPhi_",
0941                                          ("#sigma of delta Phi VS Phi of " + particleString).c_str());
0942 }
0943 
0944 void EgammaObjects::getEfficiencyHistosViaDividing() {
0945   hist_EtEfficiency_->Divide(hist_EtEfficiency_, hist_EtNumRecoOverNumTrue_, 1, 1);
0946   hist_EEfficiency_->Divide(hist_EEfficiency_, hist_ENumRecoOverNumTrue_, 1, 1);
0947   hist_EtaEfficiency_->Divide(hist_EtaEfficiency_, hist_EtaNumRecoOverNumTrue_, 1, 1);
0948   hist_PhiEfficiency_->Divide(hist_PhiEfficiency_, hist_PhiNumRecoOverNumTrue_, 1, 1);
0949 
0950   hist_EtNumRecoOverNumTrue_->Divide(hist_Et_, hist_EtNumRecoOverNumTrue_, 1, 1);
0951   hist_ENumRecoOverNumTrue_->Divide(hist_E_, hist_ENumRecoOverNumTrue_, 1, 1);
0952   hist_EtaNumRecoOverNumTrue_->Divide(hist_Eta_, hist_EtaNumRecoOverNumTrue_, 1, 1);
0953   hist_PhiNumRecoOverNumTrue_->Divide(hist_Phi_, hist_PhiNumRecoOverNumTrue_, 1, 1);
0954 }
0955 
0956 void EgammaObjects::fitHistos() {
0957   //Use our own copy for thread safety
0958   TF1 gaus("mygaus", "gaus");
0959   hist_EtOverTruth_->Fit(&gaus, "QEM");
0960   //  hist_EtNumRecoOverNumTrue_->Fit("pol1","QEM");
0961 
0962   hist_EOverTruth_->Fit(&gaus, "QEM");
0963   //  hist_ENumRecoOverNumTrue_->Fit("pol1","QEM");
0964 
0965   hist_EtaOverTruth_->Fit(&gaus, "QEM");
0966   //  hist_EtaNumRecoOverNumTrue_->Fit("pol1","QEM");
0967 
0968   hist_PhiOverTruth_->Fit(&gaus, "QEM");
0969   //  hist_PhiNumRecoOverNumTrue_->Fit("pol1","QEM");
0970 
0971   /*
0972   hist_EtOverTruthVsEt_->Fit("pol1","QEM");
0973   hist_EtOverTruthVsEta_->Fit("pol1","QEM");
0974   hist_EtOverTruthVsPhi_->Fit("pol1","QEM");
0975   hist_resolutionEtVsEt_->Fit("pol1","QEM");
0976   hist_resolutionEtVsEta_->Fit("pol1","QEM");
0977   hist_resolutionEtVsPhi_->Fit("pol1","QEM");
0978 
0979   hist_EOverTruthVsEt_->Fit("pol1","QEM");
0980   hist_EOverTruthVsEta_->Fit("pol1","QEM");
0981   hist_EOverTruthVsPhi_->Fit("pol1","QEM");
0982   hist_resolutionEVsEt_->Fit("pol1","QEM");
0983   hist_resolutionEVsEta_->Fit("pol1","QEM");
0984   hist_resolutionEVsPhi_->Fit("pol1","QEM");
0985 
0986   hist_deltaEtaVsEt_->Fit("pol1","QEM");
0987   hist_deltaEtaVsEta_->Fit("pol1","QEM");
0988   hist_deltaEtaVsPhi_->Fit("pol1","QEM");
0989   hist_resolutionEtaVsEt_->Fit("pol1","QEM");
0990   hist_resolutionEtaVsEta_->Fit("pol1","QEM");
0991   hist_resolutionEtaVsPhi_->Fit("pol1","QEM");
0992 
0993   hist_deltaPhiVsEt_->Fit("pol1","QEM");
0994   hist_deltaPhiVsEta_->Fit("pol1","QEM");
0995   hist_deltaPhiVsPhi_->Fit("pol1","QEM");
0996   hist_resolutionPhiVsEt_->Fit("pol1","QEM");
0997   hist_resolutionPhiVsEta_->Fit("pol1","QEM");
0998   hist_resolutionPhiVsPhi_->Fit("pol1","QEM");
0999   */
1000 }
1001 
1002 void EgammaObjects::applyLabels() {
1003   hist_Et_->GetXaxis()->SetTitle("Et (GeV)");
1004   hist_Et_->GetYaxis()->SetTitle("# per Et Bin");
1005   hist_EtOverTruth_->GetXaxis()->SetTitle("Reco Et/True Et");
1006   hist_EtOverTruth_->GetYaxis()->SetTitle("# per Ratio Bin");
1007   hist_EtEfficiency_->GetXaxis()->SetTitle("Et (GeV)");
1008   hist_EtEfficiency_->GetYaxis()->SetTitle("# True Reconstructed/# True per Et Bin");
1009   hist_EtNumRecoOverNumTrue_->GetXaxis()->SetTitle("Et (GeV)");
1010   hist_EtNumRecoOverNumTrue_->GetYaxis()->SetTitle("# Reco/# True per Et Bin");
1011   hist_EtOverTruthVsEt_->GetXaxis()->SetTitle("Et (GeV)");
1012   hist_EtOverTruthVsEt_->GetYaxis()->SetTitle("Reco Et/True Et per Et Bin");
1013   hist_EtOverTruthVsE_->GetXaxis()->SetTitle("E (GeV)");
1014   hist_EtOverTruthVsE_->GetYaxis()->SetTitle("Reco Et/True Et per E Bin");
1015   hist_EtOverTruthVsEta_->GetXaxis()->SetTitle("#eta (Radians)");
1016   hist_EtOverTruthVsEta_->GetYaxis()->SetTitle("Reco Et/True Et per Eta Bin");
1017   hist_EtOverTruthVsPhi_->GetXaxis()->SetTitle("#phi (Radians)");
1018   hist_EtOverTruthVsPhi_->GetYaxis()->SetTitle("Reco Et/True Et per Phi Bin");
1019   hist_resolutionEtVsEt_->GetXaxis()->SetTitle("Et (GeV)");
1020   hist_resolutionEtVsEt_->GetYaxis()->SetTitle("#sigma of Reco Et/True Et per Et Bin");
1021   hist_resolutionEtVsE_->GetXaxis()->SetTitle("E (GeV)");
1022   hist_resolutionEtVsE_->GetYaxis()->SetTitle("#sigma of Reco Et/True Et per E Bin");
1023   hist_resolutionEtVsEta_->GetXaxis()->SetTitle("#eta (Radians)");
1024   hist_resolutionEtVsEta_->GetYaxis()->SetTitle("#sigma of Reco Et/True Et per Eta Bin");
1025   hist_resolutionEtVsPhi_->GetXaxis()->SetTitle("#phi (Radians)");
1026   hist_resolutionEtVsPhi_->GetYaxis()->SetTitle("#sigma of Reco Et/True Et per Phi Bin");
1027 
1028   hist_E_->GetXaxis()->SetTitle("E (GeV)");
1029   hist_E_->GetYaxis()->SetTitle("# per E Bin");
1030   hist_EOverTruth_->GetXaxis()->SetTitle("Reco E/True E");
1031   hist_EOverTruth_->GetYaxis()->SetTitle("# per Ratio Bin");
1032   hist_EEfficiency_->GetXaxis()->SetTitle("E (GeV)");
1033   hist_EEfficiency_->GetYaxis()->SetTitle("# True Reconstructed/# True per E Bin");
1034   hist_ENumRecoOverNumTrue_->GetXaxis()->SetTitle("E (GeV)");
1035   hist_ENumRecoOverNumTrue_->GetYaxis()->SetTitle("# Reco/# True per E Bin");
1036   hist_EOverTruthVsEt_->GetXaxis()->SetTitle("Et (GeV)");
1037   hist_EOverTruthVsEt_->GetYaxis()->SetTitle("Reco E/True E per Et Bin");
1038   hist_EOverTruthVsE_->GetXaxis()->SetTitle("E (GeV)");
1039   hist_EOverTruthVsE_->GetYaxis()->SetTitle("Reco E/True E per E Bin");
1040   hist_EOverTruthVsEta_->GetXaxis()->SetTitle("#eta (Radians)");
1041   hist_EOverTruthVsEta_->GetYaxis()->SetTitle("Reco E/True E per Eta Bin");
1042   hist_EOverTruthVsPhi_->GetXaxis()->SetTitle("#phi (Radians)");
1043   hist_EOverTruthVsPhi_->GetYaxis()->SetTitle("Reco E/True E per Phi Bin");
1044   hist_resolutionEVsEt_->GetXaxis()->SetTitle("Et (GeV)");
1045   hist_resolutionEVsEt_->GetYaxis()->SetTitle("#sigma of Reco E/True E per Et Bin");
1046   hist_resolutionEVsE_->GetXaxis()->SetTitle("E (GeV)");
1047   hist_resolutionEVsE_->GetYaxis()->SetTitle("#sigma of Reco E/True E per E Bin");
1048   hist_resolutionEVsEta_->GetXaxis()->SetTitle("#eta (Radians)");
1049   hist_resolutionEVsEta_->GetYaxis()->SetTitle("#sigma of Reco E/True E per Eta Bin");
1050   hist_resolutionEVsPhi_->GetXaxis()->SetTitle("#phi (Radians)");
1051   hist_resolutionEVsPhi_->GetYaxis()->SetTitle("#sigma of Reco E/True E per Phi Bin");
1052 
1053   hist_Eta_->GetXaxis()->SetTitle("#eta (Radians)");
1054   hist_Eta_->GetYaxis()->SetTitle("# per Eta Bin");
1055   hist_EtaOverTruth_->GetXaxis()->SetTitle("Reco Eta/True Eta");
1056   hist_EtaOverTruth_->GetYaxis()->SetTitle("# per Ratio Bin");
1057   hist_EtaEfficiency_->GetXaxis()->SetTitle("#eta (Radians)");
1058   hist_EtaEfficiency_->GetYaxis()->SetTitle("# True Reconstructed/# True per Eta Bin");
1059   hist_EtaNumRecoOverNumTrue_->GetXaxis()->SetTitle("#eta (Radians)");
1060   hist_EtaNumRecoOverNumTrue_->GetYaxis()->SetTitle("# Reco/# True per Eta Bin");
1061   hist_deltaEtaVsEt_->GetXaxis()->SetTitle("Et (GeV)");
1062   hist_deltaEtaVsEt_->GetYaxis()->SetTitle("Reco Eta - True Eta per Et Bin");
1063   hist_deltaEtaVsE_->GetXaxis()->SetTitle("E (GeV)");
1064   hist_deltaEtaVsE_->GetYaxis()->SetTitle("Reco Eta - True Eta per E Bin");
1065   hist_deltaEtaVsEta_->GetXaxis()->SetTitle("#eta (Radians)");
1066   hist_deltaEtaVsEta_->GetYaxis()->SetTitle("Reco Eta - True Eta per Eta Bin");
1067   hist_deltaEtaVsPhi_->GetXaxis()->SetTitle("#phi (Radians)");
1068   hist_deltaEtaVsPhi_->GetYaxis()->SetTitle("Reco Eta - True Eta per Phi Bin");
1069   hist_resolutionEtaVsEt_->GetXaxis()->SetTitle("Et (GeV)");
1070   hist_resolutionEtaVsEt_->GetYaxis()->SetTitle("#sigma of Reco Eta - True Eta per Et Bin");
1071   hist_resolutionEtaVsE_->GetXaxis()->SetTitle("E (GeV)");
1072   hist_resolutionEtaVsE_->GetYaxis()->SetTitle("#sigma of Reco Eta - True Eta per E Bin");
1073   hist_resolutionEtaVsEta_->GetXaxis()->SetTitle("#eta (Radians)");
1074   hist_resolutionEtaVsEta_->GetYaxis()->SetTitle("#sigma of Reco Eta - True Eta per Eta Bin");
1075   hist_resolutionEtaVsPhi_->GetXaxis()->SetTitle("#phi (Radians)");
1076   hist_resolutionEtaVsPhi_->GetYaxis()->SetTitle("#sigma of Reco Eta - True Eta per Phi Bin");
1077 
1078   hist_Phi_->GetXaxis()->SetTitle("#phi (Radians)");
1079   hist_Phi_->GetYaxis()->SetTitle("# per Phi Bin");
1080   hist_PhiOverTruth_->GetXaxis()->SetTitle("Reco Phi/True Phi");
1081   hist_PhiOverTruth_->GetYaxis()->SetTitle("# per Ratio Bin");
1082   hist_PhiEfficiency_->GetXaxis()->SetTitle("#phi (Radians)");
1083   hist_PhiEfficiency_->GetYaxis()->SetTitle("# True Reconstructed/# True per Phi Bin");
1084   hist_PhiNumRecoOverNumTrue_->GetXaxis()->SetTitle("#Phi (Radians)");
1085   hist_PhiNumRecoOverNumTrue_->GetYaxis()->SetTitle("# Reco/# True per Phi Bin");
1086   hist_deltaPhiVsEt_->GetXaxis()->SetTitle("Et (GeV)");
1087   hist_deltaPhiVsEt_->GetYaxis()->SetTitle("Reco Phi - True Phi per Et Bin");
1088   hist_deltaPhiVsE_->GetXaxis()->SetTitle("E (GeV)");
1089   hist_deltaPhiVsE_->GetYaxis()->SetTitle("Reco Phi - True Phi per E Bin");
1090   hist_deltaPhiVsEta_->GetXaxis()->SetTitle("#eta (Radians)");
1091   hist_deltaPhiVsEta_->GetYaxis()->SetTitle("Reco Phi - True Phi per Eta Bin");
1092   hist_deltaPhiVsPhi_->GetXaxis()->SetTitle("#phi (Radians)");
1093   hist_deltaPhiVsPhi_->GetYaxis()->SetTitle("Reco Phi - True Phi per Phi Bin");
1094   hist_resolutionPhiVsEt_->GetXaxis()->SetTitle("Et (GeV)");
1095   hist_resolutionPhiVsEt_->GetYaxis()->SetTitle("#sigma of Reco Phi - True Phi per Et Bin");
1096   hist_resolutionPhiVsE_->GetXaxis()->SetTitle("E (GeV)");
1097   hist_resolutionPhiVsE_->GetYaxis()->SetTitle("#sigma of Reco Phi - True Phi per E Bin");
1098   hist_resolutionPhiVsEta_->GetXaxis()->SetTitle("#eta (Radians)");
1099   hist_resolutionPhiVsEta_->GetYaxis()->SetTitle("#sigma of Reco Phi - True Phi per Eta Bin");
1100   hist_resolutionPhiVsPhi_->GetXaxis()->SetTitle("#phi (Radians)");
1101   hist_resolutionPhiVsPhi_->GetYaxis()->SetTitle("#sigma of Reco Phi - True Phi per Phi Bin");
1102 
1103   std::string recoParticleName;
1104 
1105   if (particleID == 22)
1106     recoParticleName = "Reco Higgs";
1107   else if (particleID == 11)
1108     recoParticleName = "Reco Z";
1109 
1110   hist_All_recoMass_->GetXaxis()->SetTitle((recoParticleName + " Mass (GeV)").c_str());
1111   hist_All_recoMass_->GetYaxis()->SetTitle("# of Reco Masses per Mass Bin");
1112   hist_BarrelOnly_recoMass_->GetXaxis()->SetTitle((recoParticleName + " Mass (GeV)").c_str());
1113   hist_BarrelOnly_recoMass_->GetYaxis()->SetTitle("# of Reco Masses per Mass Bin");
1114   hist_EndcapOnly_recoMass_->GetXaxis()->SetTitle((recoParticleName + " Mass (GeV)").c_str());
1115   hist_EndcapOnly_recoMass_->GetYaxis()->SetTitle("# of Reco Masses per Mass Bin");
1116   hist_Mixed_recoMass_->GetXaxis()->SetTitle((recoParticleName + " Mass (GeV)").c_str());
1117   hist_Mixed_recoMass_->GetYaxis()->SetTitle("# of Reco Masses per Mass Bin");
1118   hist_recoMass_withBackgroud_NoEtCut_->GetXaxis()->SetTitle((recoParticleName + " Mass (GeV)").c_str());
1119   hist_recoMass_withBackgroud_NoEtCut_->GetYaxis()->SetTitle("# of Reco Masses per Mass Bin");
1120   hist_recoMass_withBackgroud_5EtCut_->GetXaxis()->SetTitle((recoParticleName + " Mass (GeV)").c_str());
1121   hist_recoMass_withBackgroud_5EtCut_->GetYaxis()->SetTitle("# of Reco Masses per Mass Bin");
1122   hist_recoMass_withBackgroud_10EtCut_->GetXaxis()->SetTitle((recoParticleName + " Mass (GeV)").c_str());
1123   hist_recoMass_withBackgroud_10EtCut_->GetYaxis()->SetTitle("# of Reco Masses per Mass Bin");
1124   hist_recoMass_withBackgroud_20EtCut_->GetXaxis()->SetTitle((recoParticleName + " Mass (GeV)").c_str());
1125   hist_recoMass_withBackgroud_20EtCut_->GetYaxis()->SetTitle("# of Reco Masses per Mass Bin");
1126 }
1127 
1128 void EgammaObjects::setDrawOptions() {
1129   hist_Et_->SetOption("e");
1130   hist_EtOverTruth_->SetOption("e");
1131   hist_EtEfficiency_->SetOption("e");
1132   hist_EtNumRecoOverNumTrue_->SetOption("e");
1133   hist_EtOverTruthVsEt_->SetOption("e");
1134   hist_EtOverTruthVsE_->SetOption("e");
1135   hist_EtOverTruthVsEta_->SetOption("e");
1136   hist_EtOverTruthVsPhi_->SetOption("e");
1137   hist_resolutionEtVsEt_->SetOption("e");
1138   hist_resolutionEtVsE_->SetOption("e");
1139   hist_resolutionEtVsEta_->SetOption("e");
1140   hist_resolutionEtVsPhi_->SetOption("e");
1141 
1142   hist_E_->SetOption("e");
1143   hist_EOverTruth_->SetOption("e");
1144   hist_EEfficiency_->SetOption("e");
1145   hist_ENumRecoOverNumTrue_->SetOption("e");
1146   hist_EOverTruthVsEt_->SetOption("e");
1147   hist_EOverTruthVsE_->SetOption("e");
1148   hist_EOverTruthVsEta_->SetOption("e");
1149   hist_EOverTruthVsPhi_->SetOption("e");
1150   hist_resolutionEVsEt_->SetOption("e");
1151   hist_resolutionEVsE_->SetOption("e");
1152   hist_resolutionEVsEta_->SetOption("e");
1153   hist_resolutionEVsPhi_->SetOption("e");
1154 
1155   hist_Eta_->SetOption("e");
1156   hist_EtaOverTruth_->SetOption("e");
1157   hist_EtaEfficiency_->SetOption("e");
1158   hist_EtaNumRecoOverNumTrue_->SetOption("e");
1159   hist_deltaEtaVsEt_->SetOption("e");
1160   hist_deltaEtaVsE_->SetOption("e");
1161   hist_deltaEtaVsEta_->SetOption("e");
1162   hist_deltaEtaVsPhi_->SetOption("e");
1163   hist_resolutionEtaVsEt_->SetOption("e");
1164   hist_resolutionEtaVsE_->SetOption("e");
1165   hist_resolutionEtaVsEta_->SetOption("e");
1166   hist_resolutionEtaVsPhi_->SetOption("e");
1167 
1168   hist_Phi_->SetOption("e");
1169   hist_PhiOverTruth_->SetOption("e");
1170   hist_PhiEfficiency_->SetOption("e");
1171   hist_PhiNumRecoOverNumTrue_->SetOption("e");
1172   hist_deltaPhiVsEt_->SetOption("e");
1173   hist_deltaPhiVsE_->SetOption("e");
1174   hist_deltaPhiVsEta_->SetOption("e");
1175   hist_deltaPhiVsPhi_->SetOption("e");
1176   hist_resolutionPhiVsEt_->SetOption("e");
1177   hist_resolutionPhiVsE_->SetOption("e");
1178   hist_resolutionPhiVsEta_->SetOption("e");
1179   hist_resolutionPhiVsPhi_->SetOption("e");
1180 
1181   hist_All_recoMass_->SetOption("e");
1182   hist_BarrelOnly_recoMass_->SetOption("e");
1183   hist_EndcapOnly_recoMass_->SetOption("e");
1184   hist_Mixed_recoMass_->SetOption("e");
1185   hist_recoMass_withBackgroud_NoEtCut_->SetOption("e");
1186   hist_recoMass_withBackgroud_5EtCut_->SetOption("e");
1187   hist_recoMass_withBackgroud_10EtCut_->SetOption("e");
1188   hist_recoMass_withBackgroud_20EtCut_->SetOption("e");
1189 }
1190 
1191 void EgammaObjects::saveHistos() {
1192   rootFile_->cd();
1193   rootFile_->GetDirectory(particleString.c_str())->mkdir("ET");
1194   rootFile_->cd(("/" + particleString + "/ET").c_str());
1195 
1196   hist_Et_->Write();
1197   hist_EtOverTruth_->Write();
1198   hist_EtEfficiency_->Write();
1199   hist_EtNumRecoOverNumTrue_->Write();
1200   hist_EtOverTruthVsEt_->Write();
1201   hist_EtOverTruthVsE_->Write();
1202   hist_EtOverTruthVsEta_->Write();
1203   hist_EtOverTruthVsPhi_->Write();
1204   hist_resolutionEtVsEt_->Write();
1205   hist_resolutionEtVsE_->Write();
1206   hist_resolutionEtVsEta_->Write();
1207   hist_resolutionEtVsPhi_->Write();
1208 
1209   rootFile_->cd();
1210   rootFile_->GetDirectory(particleString.c_str())->mkdir("E");
1211   rootFile_->cd(("/" + particleString + "/E").c_str());
1212 
1213   hist_E_->Write();
1214   hist_EOverTruth_->Write();
1215   hist_EEfficiency_->Write();
1216   hist_ENumRecoOverNumTrue_->Write();
1217   hist_EOverTruthVsEt_->Write();
1218   hist_EOverTruthVsE_->Write();
1219   hist_EOverTruthVsEta_->Write();
1220   hist_EOverTruthVsPhi_->Write();
1221   hist_resolutionEVsEt_->Write();
1222   hist_resolutionEVsE_->Write();
1223   hist_resolutionEVsEta_->Write();
1224   hist_resolutionEVsPhi_->Write();
1225 
1226   rootFile_->cd();
1227   rootFile_->GetDirectory(particleString.c_str())->mkdir("Eta");
1228   rootFile_->cd(("/" + particleString + "/Eta").c_str());
1229 
1230   hist_Eta_->Write();
1231   hist_EtaOverTruth_->Write();
1232   hist_EtaEfficiency_->Write();
1233   hist_EtaNumRecoOverNumTrue_->Write();
1234   hist_deltaEtaVsEt_->Write();
1235   hist_deltaEtaVsE_->Write();
1236   hist_deltaEtaVsEta_->Write();
1237   hist_deltaEtaVsPhi_->Write();
1238   hist_resolutionEtaVsEt_->Write();
1239   hist_resolutionEtaVsE_->Write();
1240   hist_resolutionEtaVsEta_->Write();
1241   hist_resolutionEtaVsPhi_->Write();
1242 
1243   rootFile_->cd();
1244   rootFile_->GetDirectory(particleString.c_str())->mkdir("Phi");
1245   rootFile_->cd(("/" + particleString + "/Phi").c_str());
1246 
1247   hist_Phi_->Write();
1248   hist_PhiOverTruth_->Write();
1249   hist_PhiEfficiency_->Write();
1250   hist_PhiNumRecoOverNumTrue_->Write();
1251   hist_deltaPhiVsEt_->Write();
1252   hist_deltaPhiVsE_->Write();
1253   hist_deltaPhiVsEta_->Write();
1254   hist_deltaPhiVsPhi_->Write();
1255   hist_resolutionPhiVsEt_->Write();
1256   hist_resolutionPhiVsE_->Write();
1257   hist_resolutionPhiVsEta_->Write();
1258   hist_resolutionPhiVsPhi_->Write();
1259 
1260   std::string recoParticleName;
1261 
1262   if (particleID == 22)
1263     recoParticleName = "HiggsRecoMass";
1264   else if (particleID == 11)
1265     recoParticleName = "ZRecoMass";
1266 
1267   rootFile_->cd();
1268   rootFile_->GetDirectory(particleString.c_str())->mkdir(recoParticleName.c_str());
1269   rootFile_->cd(("/" + particleString + "/" + recoParticleName).c_str());
1270 
1271   hist_All_recoMass_->Write();
1272   hist_BarrelOnly_recoMass_->Write();
1273   hist_EndcapOnly_recoMass_->Write();
1274   hist_Mixed_recoMass_->Write();
1275   hist_recoMass_withBackgroud_NoEtCut_->Write();
1276   hist_recoMass_withBackgroud_5EtCut_->Write();
1277   hist_recoMass_withBackgroud_10EtCut_->Write();
1278   hist_recoMass_withBackgroud_20EtCut_->Write();
1279 
1280   rootFile_->cd();
1281   rootFile_->GetDirectory(particleString.c_str())->mkdir("_TempScatterPlots");
1282   rootFile_->cd(("/" + particleString + "/_TempScatterPlots").c_str());
1283 
1284   _TEMP_scatterPlot_EtOverTruthVsEt_->Write();
1285   _TEMP_scatterPlot_EtOverTruthVsE_->Write();
1286   _TEMP_scatterPlot_EtOverTruthVsEta_->Write();
1287   _TEMP_scatterPlot_EtOverTruthVsPhi_->Write();
1288 
1289   _TEMP_scatterPlot_EOverTruthVsEt_->Write();
1290   _TEMP_scatterPlot_EOverTruthVsE_->Write();
1291   _TEMP_scatterPlot_EOverTruthVsEta_->Write();
1292   _TEMP_scatterPlot_EOverTruthVsPhi_->Write();
1293 
1294   _TEMP_scatterPlot_deltaEtaVsEt_->Write();
1295   _TEMP_scatterPlot_deltaEtaVsE_->Write();
1296   _TEMP_scatterPlot_deltaEtaVsEta_->Write();
1297   _TEMP_scatterPlot_deltaEtaVsPhi_->Write();
1298 
1299   _TEMP_scatterPlot_deltaPhiVsEt_->Write();
1300   _TEMP_scatterPlot_deltaPhiVsE_->Write();
1301   _TEMP_scatterPlot_deltaPhiVsEta_->Write();
1302   _TEMP_scatterPlot_deltaPhiVsPhi_->Write();
1303 
1304   rootFile_->cd();
1305 }