Reference Manual

Go to the documentation of this file.
00001 #include <iostream>
00002 #include <iomanip>
00003 //
00004 
00005 #include "DQMOffline/EGamma/plugins/PhotonAnalyzer.h"
00006 
00007 
00019 using namespace std;
00020 
00021 
00022 PhotonAnalyzer::PhotonAnalyzer( const edm::ParameterSet& pset )
00023 {
00024 
00025     fName_                  = pset.getUntrackedParameter<string>("Name");
00026     verbosity_              = pset.getUntrackedParameter<int>("Verbosity");
00027 
00028     prescaleFactor_         = pset.getUntrackedParameter<int>("prescaleFactor",1);
00029 
00030     photonProducer_         = pset.getParameter<string>("phoProducer");
00031     photonCollection_       = pset.getParameter<string>("photonCollection");
00032 
00033     barrelRecHitProducer_   = pset.getParameter<string>("barrelRecHitProducer");
00034     barrelRecHitCollection_ = pset.getParameter<string>("barrelRecHitCollection");
00035 
00036     endcapRecHitProducer_   = pset.getParameter<string>("endcapRecHitProducer");
00037     endcapRecHitCollection_ = pset.getParameter<string>("endcapRecHitCollection");
00038 
00039     triggerEvent_           = pset.getParameter<edm::InputTag>("triggerEvent");
00040 
00041     minPhoEtCut_            = pset.getParameter<double>("minPhoEtCut");
00042     invMassEtCut_           = pset.getParameter<double>("invMassEtCut");
00043 
00044     cutStep_                = pset.getParameter<double>("cutStep");
00045     numberOfSteps_          = pset.getParameter<int>("numberOfSteps");
00046 
00047     useBinning_             = pset.getParameter<bool>("useBinning");
00048     useTriggerFiltering_    = pset.getParameter<bool>("useTriggerFiltering");
00049 
00050     standAlone_             = pset.getParameter<bool>("standAlone");
00051     outputFileName_         = pset.getParameter<string>("OutputFileName");
00052 
00053     isolationStrength_      = pset.getParameter<int>("isolationStrength");
00054 
00055     isHeavyIon_             = pset.getUntrackedParameter<bool>("isHeavyIon",false);
00056 
00057     parameters_ = pset;
00058 
00059     histo_index_photons_ = 0;
00060     histo_index_conversions_ = 0;
00061     histo_index_efficiency_ = 0;
00062     histo_index_invMass_ = 0;
00063 }
00064 
00065 
00066 
00067 PhotonAnalyzer::~PhotonAnalyzer() {}
00068 
00069 
00070 void PhotonAnalyzer::beginJob()
00071 {
00072 
00073   nEvt_=0;
00074   nEntry_=0;
00075 
00076   dbe_ = 0;
00077   dbe_ = edm::Service<DQMStore>().operator->();
00078 
00079 
00080 
00081   double eMin = parameters_.getParameter<double>("eMin");
00082   double eMax = parameters_.getParameter<double>("eMax");
00083   int    eBin = parameters_.getParameter<int>("eBin");
00084 
00085   double etMin = parameters_.getParameter<double>("etMin");
00086   double etMax = parameters_.getParameter<double>("etMax");
00087   int    etBin = parameters_.getParameter<int>("etBin");
00088 
00089   double sumMin = parameters_.getParameter<double>("sumMin");
00090   double sumMax = parameters_.getParameter<double>("sumMax");
00091   int    sumBin = parameters_.getParameter<int>("sumBin");
00092 
00093   double etaMin = parameters_.getParameter<double>("etaMin");
00094   double etaMax = parameters_.getParameter<double>("etaMax");
00095   int    etaBin = parameters_.getParameter<int>("etaBin");
00096 
00097   double phiMin = parameters_.getParameter<double>("phiMin");
00098   double phiMax = parameters_.getParameter<double>("phiMax");
00099   int    phiBin = parameters_.getParameter<int>("phiBin");
00100 
00101   double r9Min = parameters_.getParameter<double>("r9Min");
00102   double r9Max = parameters_.getParameter<double>("r9Max");
00103   int    r9Bin = parameters_.getParameter<int>("r9Bin");
00104 
00105   double hOverEMin = parameters_.getParameter<double>("hOverEMin");
00106   double hOverEMax = parameters_.getParameter<double>("hOverEMax");
00107   int    hOverEBin = parameters_.getParameter<int>("hOverEBin");
00108 
00109   double xMin = parameters_.getParameter<double>("xMin");
00110   double xMax = parameters_.getParameter<double>("xMax");
00111   int    xBin = parameters_.getParameter<int>("xBin");
00112 
00113   double yMin = parameters_.getParameter<double>("yMin");
00114   double yMax = parameters_.getParameter<double>("yMax");
00115   int    yBin = parameters_.getParameter<int>("yBin");
00116 
00117   double numberMin = parameters_.getParameter<double>("numberMin");
00118   double numberMax = parameters_.getParameter<double>("numberMax");
00119   int    numberBin = parameters_.getParameter<int>("numberBin");
00120 
00121   double zMin = parameters_.getParameter<double>("zMin");
00122   double zMax = parameters_.getParameter<double>("zMax");
00123   int    zBin = parameters_.getParameter<int>("zBin");
00124 
00125   double rMin = parameters_.getParameter<double>("rMin");
00126   double rMax = parameters_.getParameter<double>("rMax");
00127   int    rBin = parameters_.getParameter<int>("rBin");
00128 
00129   double dPhiTracksMin = parameters_.getParameter<double>("dPhiTracksMin");
00130   double dPhiTracksMax = parameters_.getParameter<double>("dPhiTracksMax");
00131   int    dPhiTracksBin = parameters_.getParameter<int>("dPhiTracksBin");
00132 
00133   double dEtaTracksMin = parameters_.getParameter<double>("dEtaTracksMin");
00134   double dEtaTracksMax = parameters_.getParameter<double>("dEtaTracksMax");
00135   int    dEtaTracksBin = parameters_.getParameter<int>("dEtaTracksBin");
00136 
00137   double sigmaIetaMin = parameters_.getParameter<double>("sigmaIetaMin");
00138   double sigmaIetaMax = parameters_.getParameter<double>("sigmaIetaMax");
00139   int    sigmaIetaBin = parameters_.getParameter<int>("sigmaIetaBin");
00140 
00141   double eOverPMin = parameters_.getParameter<double>("eOverPMin");
00142   double eOverPMax = parameters_.getParameter<double>("eOverPMax");
00143   int    eOverPBin = parameters_.getParameter<int>("eOverPBin");
00144 
00145   double chi2Min = parameters_.getParameter<double>("chi2Min");
00146   double chi2Max = parameters_.getParameter<double>("chi2Max");
00147   int    chi2Bin = parameters_.getParameter<int>("chi2Bin");
00148 
00149 
00150   int reducedEtBin = etBin/2;
00151   int reducedEtaBin = etaBin/2;
00152   int reducedSumBin = sumBin/2;
00153   int reducedR9Bin = r9Bin/2;
00154 
00155 
00156   parts_.push_back("AllEcal");
00157   parts_.push_back("Barrel");
00158   parts_.push_back("Endcaps");
00159 
00160   types_.push_back("All");
00161   types_.push_back("GoodCandidate");
00162   types_.push_back("Background");
00163 
00164 
00165 
00167 
00168   if (dbe_) {
00169 
00170     dbe_->setCurrentFolder("Egamma/PhotonAnalyzer");
00171 
00172     //int values stored in MEs to keep track of how many histograms are in each folder
00173     totalNumberOfHistos_efficiencyFolder =  dbe_->bookInt("numberOfHistogramsInEfficiencyFolder");
00174     totalNumberOfHistos_photonsFolder =     dbe_->bookInt("numberOfHistogramsInPhotonsFolder");
00175     totalNumberOfHistos_conversionsFolder = dbe_->bookInt("numberOfHistogramsInConversionsFolder");
00176     totalNumberOfHistos_invMassFolder =     dbe_->bookInt("numberOfHistogramsInInvMassFolder");
00177 
00178 
00179     //Efficiency histograms
00180 
00181     dbe_->setCurrentFolder("Egamma/PhotonAnalyzer/Efficiencies");
00182 
00183     //don't number these histograms with the "bookHisto" method, since they'll be erased in the offline client
00184     h_phoEta_Loose_ = dbe_->book1D("phoEtaLoose","Loose Photon #eta",etaBin,etaMin,etaMax);
00185     h_phoEta_Tight_ = dbe_->book1D("phoEtaTight","Tight Photon #eta",etaBin,etaMin,etaMax);
00186     h_phoEt_Loose_  = dbe_->book1D("phoEtLoose", "Loose Photon E_{T}",etBin,etMin,etMax);
00187     h_phoEt_Tight_  = dbe_->book1D("phoEtTight", "Tight Photon E_{T}",etBin,etMin,etMax);
00188 
00189 
00190     h_phoEta_preHLT_  = dbe_->book1D("phoEtaPreHLT", "Photon #eta: before HLT",etaBin,etaMin,etaMax);
00191     h_phoEta_postHLT_ = dbe_->book1D("phoEtaPostHLT","Photon #eta: after HLT",etaBin,etaMin,etaMax);
00192     h_phoEt_preHLT_   = dbe_->book1D("phoEtPreHLT",  "Photon E_{T}: before HLT",etBin,etMin,etMax);
00193     h_phoEt_postHLT_  = dbe_->book1D("phoEtPostHLT", "Photon E_{T}: after HLT",etBin,etMin,etMax);
00194 
00195     h_convEta_Loose_ = dbe_->book1D("convEtaLoose","Converted Loose Photon #eta",etaBin,etaMin,etaMax);
00196     h_convEta_Tight_ = dbe_->book1D("convEtaTight","Converted Tight Photon #eta",etaBin,etaMin,etaMax);
00197     h_convEt_Loose_  = dbe_->book1D("convEtLoose", "Converted Loose Photon E_{T}",etBin,etMin,etMax);
00198     h_convEt_Tight_  = dbe_->book1D("convEtTight", "Converted Tight Photon E_{T}",etBin,etMin,etMax);
00199 
00200     h_phoEta_Vertex_ = dbe_->book1D("phoEtaVertex","Converted Photons before valid vertex cut: #eta",etaBin,etaMin,etaMax);
00201 
00202 
00203     vector<MonitorElement*> temp1DVectorEta;
00204     vector<MonitorElement*> temp1DVectorPhi;
00205     vector<vector<MonitorElement*> > temp2DVectorPhi;
00206 
00207 
00208     for(int cut = 0; cut != numberOfSteps_; ++cut){ //looping over Et cut values
00209       for(uint type=0;type!=types_.size();++type){  //looping over isolation type
00210         currentFolder_.str("");
00211         currentFolder_ << "Egamma/PhotonAnalyzer/" << types_[type] << "Photons/Et above " << (cut+1)*cutStep_ << " GeV/Conversions";
00212         dbe_->setCurrentFolder(currentFolder_.str());
00213 
00214         temp1DVectorEta.push_back(dbe_->book1D("phoConvEtaForEfficiency","Converted Photon #eta;#eta",etaBin,etaMin,etaMax));
00215         for(uint part=0;part!=parts_.size();++part){
00216           temp1DVectorPhi.push_back(dbe_->book1D("phoConvPhiForEfficiency"+parts_[part],"Converted Photon #phi;#phi",phiBin,phiMin,phiMax));
00217         }
00218         temp2DVectorPhi.push_back(temp1DVectorPhi);
00219         temp1DVectorPhi.clear();
00220       }
00221       h_phoConvEtaForEfficiency_.push_back(temp1DVectorEta);
00222       temp1DVectorEta.clear();
00223       h_phoConvPhiForEfficiency_.push_back(temp2DVectorPhi);
00224       temp2DVectorPhi.clear();
00225     }
00226 
00227 
00228 
00229 
00230     //Invariant mass plots
00231 
00232     dbe_->setCurrentFolder("Egamma/PhotonAnalyzer/InvMass");
00233 
00234     h_invMassAllPhotons_    = bookHisto("invMassAllIsolatedPhotons","Two photon invariant mass: All isolated photons;M (GeV)",etBin,etMin,etMax);
00235     h_invMassZeroWithTracks_= bookHisto("invMassZeroWithTracks",    "Two photon invariant mass: Neither has tracks;M (GeV)",  etBin,etMin,etMax);
00236     h_invMassOneWithTracks_ = bookHisto("invMassOneWithTracks",     "Two photon invariant mass: Only one has tracks;M (GeV)", etBin,etMin,etMax);
00237     h_invMassTwoWithTracks_ = bookHisto("invMassTwoWithTracks",     "Two photon invariant mass: Both have tracks;M (GeV)",    etBin,etMin,etMax);
00238 
00239 
00241 
00242     //ENERGY VARIABLES
00243 
00244     h_phoE_  = book3DHistoVector("1D","phoE","Energy;E (GeV)",eBin,eMin,eMax);
00245     h_phoEt_ = book3DHistoVector("1D","phoEt","E_{T};E_{T} (GeV)", etBin,etMin,etMax);
00246 
00247     //NUMBER OF PHOTONS
00248 
00249     h_nPho_  = book3DHistoVector("1D","nPho","Number of Photons per Event;# #gamma",numberBin,numberMin,numberMax);
00250 
00251     //GEOMETRICAL VARIABLES
00252 
00253     //photon eta/phi
00254     h_phoEta_ = book2DHistoVector("1D","phoEta","#eta;#eta",etaBin,etaMin,etaMax) ;
00255     h_phoPhi_ = book3DHistoVector("1D","phoPhi","#phi;#phi",phiBin,phiMin,phiMax) ;
00256 
00257     //supercluster eta/phi
00258     h_scEta_ = book2DHistoVector("1D","scEta","SuperCluster #eta;#eta",etaBin,etaMin,etaMax) ;
00259     h_scPhi_ = book3DHistoVector("1D","scPhi","SuperCluster #phi;#phi",phiBin,phiMin,phiMax) ;
00260 
00261     //SHOWER SHAPE VARIABLES
00262 
00263     //r9
00264     h_r9_      = book3DHistoVector("1D","r9","R9;R9",r9Bin,r9Min, r9Max);
00265     h_r9VsEt_  = book2DHistoVector("2D","r9VsEt2D","R9 vs E_{T};E_{T} (GeV);R9",reducedEtBin,etMin,etMax,reducedR9Bin,r9Min,r9Max);
00266     p_r9VsEt_  = book2DHistoVector("Profile","r9VsEt","Avg R9 vs E_{T};E_{T} (GeV);R9",etBin,etMin,etMax,r9Bin,r9Min,r9Max);
00267     h_r9VsEta_ = book2DHistoVector("2D","r9VsEta2D","R9 vs #eta;#eta;R9",reducedEtaBin,etaMin,etaMax,reducedR9Bin,r9Min,r9Max);
00268     p_r9VsEta_ = book2DHistoVector("Profile","r9VsEta","Avg R9 vs #eta;#eta;R9",etaBin,etaMin,etaMax,r9Bin,r9Min,r9Max);
00269 
00270     //sigma ieta ieta
00271     h_phoSigmaIetaIeta_   = book3DHistoVector("1D","phoSigmaIetaIeta","#sigma_{i#etai#eta};#sigma_{i#etai#eta}",sigmaIetaBin,sigmaIetaMin,sigmaIetaMax);
00272     h_sigmaIetaIetaVsEta_ = book2DHistoVector("2D","sigmaIetaIetaVsEta2D","#sigma_{i#etai#eta} vs #eta;#eta;#sigma_{i#etai#eta}",reducedEtaBin,etaMin,etaMax,sigmaIetaBin,sigmaIetaMin,sigmaIetaMax);
00273     p_sigmaIetaIetaVsEta_ = book2DHistoVector("Profile","sigmaIetaIetaVsEta","Avg #sigma_{i#etai#eta} vs #eta;#eta;#sigma_{i#etai#eta}",etaBin,etaMin,etaMax,sigmaIetaBin,sigmaIetaMin,sigmaIetaMax);
00274 
00275     //e1x5
00276     h_e1x5VsEt_  = book2DHistoVector("2D","e1x5VsEt2D","E1x5 vs E_{T};E_{T} (GeV);E1X5 (GeV)",reducedEtBin,etMin,etMax,reducedEtBin,etMin,etMax);
00277     p_e1x5VsEt_  = book2DHistoVector("Profile","e1x5VsEt","Avg E1x5 vs E_{T};E_{T} (GeV);E1X5 (GeV)",etBin,etMin,etMax,etBin,etMin,etMax);
00278     h_e1x5VsEta_ = book2DHistoVector("2D","e1x5VsEta2D","E1x5 vs #eta;#eta;E1X5 (GeV)",reducedEtaBin,etaMin,etaMax,reducedEtBin,etMin,etMax);
00279     p_e1x5VsEta_ = book2DHistoVector("Profile","e1x5VsEta","Avg E1x5 vs #eta;#eta;E1X5 (GeV)",etaBin,etaMin,etaMax,etBin,etMin,etMax);
00280 
00281     //e2x5
00282     h_e2x5VsEt_  = book2DHistoVector("2D","e2x5VsEt2D","E2x5 vs E_{T};E_{T} (GeV);E2X5 (GeV)",reducedEtBin,etMin,etMax,reducedEtBin,etMin,etMax);
00283     p_e2x5VsEt_  = book2DHistoVector("Profile","e2x5VsEt","Avg E2x5 vs E_{T};E_{T} (GeV);E2X5 (GeV)",etBin,etMin,etMax,etBin,etMin,etMax);
00284     h_e2x5VsEta_ = book2DHistoVector("2D","e2x5VsEta2D","E2x5 vs #eta;#eta;E2X5 (GeV)",reducedEtaBin,etaMin,etaMax,reducedEtBin,etMin,etMax);
00285     p_e2x5VsEta_ = book2DHistoVector("Profile","e2x5VsEta","Avg E2x5 vs #eta;#eta;E2X5 (GeV)",etaBin,etaMin,etaMax,etBin,etMin,etMax);
00286 
00287     //r1x5
00288     h_r1x5VsEt_  = book2DHistoVector("2D","r1x5VsEt2D","R1x5 vs E_{T};E_{T} (GeV);R1X5",reducedEtBin,etMin,etMax,reducedR9Bin,r9Min,r9Max);
00289     p_r1x5VsEt_  = book2DHistoVector("Profile","r1x5VsEt","Avg R1x5 vs E_{T};E_{T} (GeV);R1X5",etBin,etMin,etMax,r9Bin,r9Min,r9Max);
00290     h_r1x5VsEta_ = book2DHistoVector("2D","r1x5VsEta2D","R1x5 vs #eta;#eta;R1X5",reducedEtaBin,etaMin,etaMax,reducedR9Bin,r9Min,r9Max);
00291     p_r1x5VsEta_ = book2DHistoVector("Profile","r1x5VsEta","Avg R1x5 vs #eta;#eta;R1X5",etaBin,etaMin,etaMax,r9Bin,r9Min,r9Max);
00292 
00293     //r2x5
00294     h_r2x5VsEt_  = book2DHistoVector("2D","r2x5VsEt2D","R2x5 vs E_{T};E_{T} (GeV);R2X5",reducedEtBin,etMin,etMax,reducedR9Bin,r9Min,r9Max);
00295     p_r2x5VsEt_  = book2DHistoVector("Profile","r2x5VsEt","Avg R2x5 vs E_{T};E_{T} (GeV);R2X5",etBin,etMin,etMax,r9Bin,r9Min,r9Max);
00296     h_r2x5VsEta_ = book2DHistoVector("2D","r2x5VsEta2D","R2x5 vs #eta;#eta;R2X5",reducedEtaBin,etaMin,etaMax,reducedR9Bin,r9Min,r9Max);
00297     p_r2x5VsEta_ = book2DHistoVector("Profile","r2x5VsEta","Avg R2x5 vs #eta;#eta;R2X5",etaBin,etaMin,etaMax,r9Bin,r9Min,r9Max);
00298 
00299     //maxEXtalOver3x3
00300     h_maxEXtalOver3x3VsEt_  = book2DHistoVector("2D","maxEXtalOver3x3VsEt2D","(Max Xtal E)/E3x3 vs E_{T};E_{T} (GeV);(Max Xtal E)/E3x3",reducedEtBin,etMin,etMax,r9Bin,r9Min,r9Max);
00301     p_maxEXtalOver3x3VsEt_  = book2DHistoVector("Profile","maxEXtalOver3x3VsEt","Avg (Max Xtal E)/E3x3 vs E_{T};E_{T} (GeV);(Max Xtal E)/E3x3",etBin,etMin,etMax,r9Bin,r9Min,r9Max);
00302     h_maxEXtalOver3x3VsEta_ = book2DHistoVector("2D","maxEXtalOver3x3VsEta2D","(Max Xtal E)/E3x3 vs #eta;#eta;(Max Xtal E)/E3x3",reducedEtaBin,etaMin,etaMax,r9Bin,r9Min,r9Max);
00303     p_maxEXtalOver3x3VsEta_ = book2DHistoVector("Profile","maxEXtalOver3x3VsEta","Avg (Max Xtal E)/E3x3 vs #eta;#eta;(Max Xtal E)/E3x3",etaBin,etaMin,etaMax,r9Bin,r9Min,r9Max);
00304 
00305 
00306     //TRACK ISOLATION VARIABLES
00307 
00308     //nTrackIsolSolid
00309     h_nTrackIsolSolid_       = book2DHistoVector("1D","nIsoTracksSolid","Number Of Tracks in the Solid Iso Cone;# tracks",numberBin,numberMin,numberMax);
00310     h_nTrackIsolSolidVsEt_   = book2DHistoVector("2D","nIsoTracksSolidVsEt2D","Number Of Tracks in the Solid Iso Cone vs E_{T};E_{T};# tracks",reducedEtBin,etMin, etMax,numberBin,numberMin,numberMax);
00311     p_nTrackIsolSolidVsEt_   = book2DHistoVector("Profile","nIsoTracksSolidVsEt","Avg Number Of Tracks in the Solid Iso Cone vs E_{T};E_{T};# tracks",etBin,etMin,etMax,numberBin,numberMin,numberMax);
00312     h_nTrackIsolSolidVsEta_  = book2DHistoVector("2D","nIsoTracksSolidVsEta2D","Number Of Tracks in the Solid Iso Cone vs #eta;#eta;# tracks",reducedEtaBin,etaMin, etaMax,numberBin,numberMin,numberMax);
00313     p_nTrackIsolSolidVsEta_  = book2DHistoVector("Profile","nIsoTracksSolidVsEta","Avg Number Of Tracks in the Solid Iso Cone vs #eta;#eta;# tracks",etaBin,etaMin, etaMax,numberBin,numberMin,numberMax);
00314 
00315     //nTrackIsolHollow
00316     h_nTrackIsolHollow_      = book2DHistoVector("1D","nIsoTracksHollow","Number Of Tracks in the Hollow Iso Cone;# tracks",numberBin,numberMin,numberMax);
00317     h_nTrackIsolHollowVsEt_  = book2DHistoVector("2D","nIsoTracksHollowVsEt2D","Number Of Tracks in the Hollow Iso Cone vs E_{T};E_{T};# tracks",reducedEtBin,etMin, etMax,numberBin,numberMin,numberMax);
00318     p_nTrackIsolHollowVsEt_  = book2DHistoVector("Profile","nIsoTracksHollowVsEt","Avg Number Of Tracks in the Hollow Iso Cone vs E_{T};E_{T};# tracks",etBin,etMin,etMax,numberBin,numberMin,numberMax);
00319     h_nTrackIsolHollowVsEta_ = book2DHistoVector("2D","nIsoTracksHollowVsEta2D","Number Of Tracks in the Hollow Iso Cone vs #eta;#eta;# tracks",reducedEtaBin,etaMin, etaMax,numberBin,numberMin,numberMax);
00320     p_nTrackIsolHollowVsEta_ = book2DHistoVector("Profile","nIsoTracksHollowVsEta","Avg Number Of Tracks in the Hollow Iso Cone vs #eta;#eta;# tracks",etaBin,etaMin, etaMax,numberBin,numberMin,numberMax);
00321 
00322     //trackPtSumSolid
00323     h_trackPtSumSolid_       = book2DHistoVector("1D","isoPtSumSolid","Track P_{T} Sum in the Solid Iso Cone;P_{T} (GeV)",sumBin,sumMin,sumMax);
00324     h_trackPtSumSolidVsEt_   = book2DHistoVector("2D","isoPtSumSolidVsEt2D","Track P_{T} Sum in the Solid Iso Cone;E_{T} (GeV);P_{T} (GeV)",reducedEtBin,etMin, etMax,reducedSumBin,sumMin,sumMax);
00325     p_trackPtSumSolidVsEt_   = book2DHistoVector("Profile","isoPtSumSolidVsEt","Avg Track P_{T} Sum in the Solid Iso Cone vs E_{T};E_{T} (GeV);P_{T} (GeV)",etBin,etMin,etMax,sumBin,sumMin,sumMax);
00326     h_trackPtSumSolidVsEta_  = book2DHistoVector("2D","isoPtSumSolidVsEta2D","Track P_{T} Sum in the Solid Iso Cone;#eta;P_{T} (GeV)",reducedEtaBin,etaMin, etaMax,reducedSumBin,sumMin,sumMax);
00327     p_trackPtSumSolidVsEta_  = book2DHistoVector("Profile","isoPtSumSolidVsEta","Avg Track P_{T} Sum in the Solid Iso Cone vs #eta;#eta;P_{T} (GeV)",etaBin,etaMin, etaMax,sumBin,sumMin,sumMax);
00328 
00329     //trackPtSumHollow
00330     h_trackPtSumHollow_      = book2DHistoVector("1D","isoPtSumHollow","Track P_{T} Sum in the Hollow Iso Cone;P_{T} (GeV)",sumBin,sumMin,sumMax);
00331     h_trackPtSumHollowVsEt_  = book2DHistoVector("2D","isoPtSumHollowVsEt2D","Track P_{T} Sum in the Hollow Iso Cone;E_{T} (GeV);P_{T} (GeV)",reducedEtBin,etMin, etMax,reducedSumBin,sumMin,sumMax);
00332     p_trackPtSumHollowVsEt_  = book2DHistoVector("Profile","isoPtSumHollowVsEt","Avg Track P_{T} Sum in the Hollow Iso Cone vs E_{T};E_{T} (GeV);P_{T} (GeV)",etBin,etMin,etMax,sumBin,sumMin,sumMax);
00333     h_trackPtSumHollowVsEta_ = book2DHistoVector("2D","isoPtSumHollowVsEta2D","Track P_{T} Sum in the Hollow Iso Cone;#eta;P_{T} (GeV)",reducedEtaBin,etaMin, etaMax,reducedSumBin,sumMin,sumMax);
00334     p_trackPtSumHollowVsEta_ = book2DHistoVector("Profile","isoPtSumHollowVsEta","Avg Track P_{T} Sum in the Hollow Iso Cone vs #eta;#eta;P_{T} (GeV)",etaBin,etaMin, etaMax,sumBin,sumMin,sumMax);
00335 
00336 
00337     //CALORIMETER ISOLATION VARIABLES
00338 
00339     //ecal sum
00340     h_ecalSum_      = book2DHistoVector("1D","ecalSum","Ecal Sum in the Iso Cone;E (GeV)",sumBin,sumMin,sumMax);
00341     h_ecalSumVsEt_  = book2DHistoVector("2D","ecalSumVsEt2D","Ecal Sum in the Iso Cone;E_{T} (GeV);E (GeV)",reducedEtBin,etMin, etMax,reducedSumBin,sumMin,sumMax);
00342     p_ecalSumVsEt_  = book3DHistoVector("Profile","ecalSumVsEt","Avg Ecal Sum in the Iso Cone vs E_{T};E_{T} (GeV);E (GeV)",etBin,etMin, etMax,sumBin,sumMin,sumMax);
00343     h_ecalSumVsEta_ = book2DHistoVector("2D","ecalSumVsEta2D","Ecal Sum in the Iso Cone;#eta;E (GeV)",reducedEtaBin,etaMin, etaMax,reducedSumBin,sumMin,sumMax);
00344     p_ecalSumVsEta_ = book2DHistoVector("Profile","ecalSumVsEta","Avg Ecal Sum in the Iso Cone vs #eta;#eta;E (GeV)",etaBin,etaMin, etaMax,sumBin,sumMin,sumMax);
00345 
00346     //hcal sum
00347     h_hcalSum_      = book2DHistoVector("1D","hcalSum","Hcal Sum in the Iso Cone;E (GeV)",sumBin,sumMin,sumMax);
00348     h_hcalSumVsEt_  = book2DHistoVector("2D","hcalSumVsEt2D","Hcal Sum in the Iso Cone;E_{T} (GeV);E (GeV)",reducedEtBin,etMin, etMax,reducedSumBin,sumMin,sumMax);
00349     p_hcalSumVsEt_  = book3DHistoVector("Profile","hcalSumVsEt","Avg Hcal Sum in the Iso Cone vs E_{T};E_{T} (GeV);E (GeV)",etBin,etMin, etMax,sumBin,sumMin,sumMax);
00350     h_hcalSumVsEta_ = book2DHistoVector("2D","hcalSumVsEta2D","Hcal Sum in the Iso Cone;#eta;E (GeV)",reducedEtaBin,etaMin, etaMax,reducedSumBin,sumMin,sumMax);
00351     p_hcalSumVsEta_ = book2DHistoVector("Profile","hcalSumVsEta","Avg Hcal Sum in the Iso Cone vs #eta;#eta;E (GeV)",etaBin,etaMin, etaMax,sumBin,sumMin,sumMax);
00352 
00353     //h over e
00354     h_hOverE_       = book3DHistoVector("1D","hOverE","H/E;H/E",hOverEBin,hOverEMin,hOverEMax);
00355     p_hOverEVsEt_   = book2DHistoVector("Profile","hOverEVsEt","Avg H/E vs Et;E_{T} (GeV);H/E",etBin,etMin,etMax,hOverEBin,hOverEMin,hOverEMax);
00356     p_hOverEVsEta_  = book2DHistoVector("Profile","hOverEVsEta","Avg H/E vs #eta;#eta;H/E",etaBin,etaMin,etaMax,hOverEBin,hOverEMin,hOverEMax);
00357     h_h1OverE_      = book3DHistoVector("1D","h1OverE","H/E for Depth 1;H/E",hOverEBin,hOverEMin,hOverEMax);
00358     h_h2OverE_      = book3DHistoVector("1D","h2OverE","H/E for Depth 2;H/E",hOverEBin,hOverEMin,hOverEMax);
00359 
00360 
00361     //OTHER VARIABLES
00362 
00363     //bad channel histograms
00364     h_phoEt_BadChannels_  = book2DHistoVector("1D","phoEtBadChannels", "Fraction Containing Bad Channels: E_{T};E_{T} (GeV)",etBin,etMin,etMax);
00365     h_phoEta_BadChannels_ = book2DHistoVector("1D","phoEtaBadChannels","Fraction Containing Bad Channels: #eta;#eta",etaBin,etaMin,etaMax);
00366     h_phoPhi_BadChannels_ = book2DHistoVector("1D","phoPhiBadChannels","Fraction Containing Bad Channels: #phi;#phi",phiBin,phiMin,phiMax);
00367 
00368 
00370 
00371     dbe_->setCurrentFolder("Egamma/PhotonAnalyzer/AllPhotons/Et Above 0 GeV/Conversions");
00372 
00373     //ENERGY VARIABLES
00374 
00375     h_phoConvE_  = book3DHistoVector("1D","phoConvE","E;E (GeV)",eBin,eMin,eMax);
00376     h_phoConvEt_ = book3DHistoVector("1D","phoConvEt","E_{T};E_{T} (GeV)",etBin,etMin,etMax);
00377 
00378     //GEOMETRICAL VARIABLES
00379 
00380     h_phoConvEta_ = book2DHistoVector("1D","phoConvEta","#eta;#eta",etaBin,etaMin,etaMax);
00381     h_phoConvPhi_ = book3DHistoVector("1D","phoConvPhi","#phi;#phi",phiBin,phiMin,phiMax);
00382 
00383     //NUMBER OF PHOTONS
00384 
00385     h_nConv_ = book3DHistoVector("1D","nConv","Number Of Conversions per Event ;# conversions",numberBin,numberMin,numberMax);
00386 
00387     //SHOWER SHAPE VARIABLES
00388 
00389     h_phoConvR9_ = book3DHistoVector("1D","phoConvR9","R9;R9",r9Bin,r9Min,r9Max);
00390 
00391     //TRACK RELATED VARIABLES
00392 
00393     h_eOverPTracks_ = book3DHistoVector("1D","eOverPTracks","E/P;E/P",eOverPBin,eOverPMin,eOverPMax);
00394     h_pOverETracks_ = book3DHistoVector("1D","pOverETracks","P/E;P/E",eOverPBin,eOverPMin,eOverPMax);
00395 
00396     h_dPhiTracksAtVtx_  = book3DHistoVector("1D","dPhiTracksAtVtx", "#Delta#phi of Tracks at Vertex;#Delta#phi",dPhiTracksBin,dPhiTracksMin,dPhiTracksMax);
00397     h_dPhiTracksAtEcal_ = book3DHistoVector("1D","dPhiTracksAtEcal", "Abs(#Delta#phi) of Tracks at Ecal;#Delta#phi",dPhiTracksBin,0.,dPhiTracksMax);
00398     h_dEtaTracksAtEcal_ = book3DHistoVector("1D","dEtaTracksAtEcal", "#Delta#eta of Tracks at Ecal;#Delta#eta",dEtaTracksBin,dEtaTracksMin,dEtaTracksMax);
00399 
00400     h_dCotTracks_      = book3DHistoVector("1D","dCotTracks","#Deltacot(#theta) of Tracks;#Deltacot(#theta)",dEtaTracksBin,dEtaTracksMin,dEtaTracksMax);
00401     p_dCotTracksVsEta_ = book2DHistoVector("Profile","dCotTracksVsEta","Avg #Deltacot(#theta) of Tracks vs #eta;#eta;#Deltacot(#theta)",etaBin,etaMin,etaMax,dEtaTracksBin,dEtaTracksMin,dEtaTracksMax);
00402 
00403     p_nHitsVsEta_ = book2DHistoVector("Profile","nHitsVsEta","Avg Number of Hits per Track vs #eta;#eta;# hits",etaBin,etaMin,etaMax,etaBin,0,16);
00404 
00405     h_tkChi2_      = book2DHistoVector("1D","tkChi2","#chi^{2} of Track Fitting;#chi^{2}",chi2Bin,chi2Min,chi2Max);
00406     p_tkChi2VsEta_ = book2DHistoVector("Profile","tkChi2VsEta","Avg #chi^{2} of Track Fitting vs #eta;#eta;#chi^{2}",etaBin,etaMin,etaMax,chi2Bin,chi2Min,chi2Max);
00407 
00408     //VERTEX RELATED VARIABLES
00409 
00410     h_convVtxRvsZ_ = book2DHistoVector("2D","convVtxRvsZ","Vertex Position;Z (cm);R (cm)",1000,zMin,zMax,rBin,rMin,rMax);
00411     h_convVtxZEndcap_    = book2DHistoVector("1D","convVtxZEndcap",   "Vertex Position: #eta > 1.5;Z (cm)",zBin,zMin,zMax);
00412     h_convVtxZ_    = book2DHistoVector("1D","convVtxZ",   "Vertex Position;Z (cm)",zBin,zMin,zMax);
00413     h_convVtxR_    = book2DHistoVector("1D","convVtxR",   "Vertex Position: #eta < 1;R (cm)",rBin,rMin,rMax);
00414     h_convVtxYvsX_ = book2DHistoVector("2D","convVtxYvsX","Vertex Position: #eta < 1;X (cm);Y (cm)",xBin,xMin,xMax,yBin,yMin,yMax);
00415 
00416 
00417 
00418     h_vertexChi2Prob_ = book2DHistoVector("1D","vertexChi2Prob","#chi^{2} Probability of Vertex Fitting;#chi^{2}",100,0.,1.0);
00419 
00420 
00421   }//end if(dbe_)
00422 
00423 
00424 }//end BeginJob
00425 
00426 
00427 
00428 void PhotonAnalyzer::analyze( const edm::Event& e, const edm::EventSetup& esup )
00429 {
00430   using namespace edm;
00431 
00432   if (nEvt_% prescaleFactor_ ) return;
00433   nEvt_++;
00434   LogInfo("PhotonAnalyzer") << "PhotonAnalyzer Analyzing event number: " << e.id() << " Global Counter " << nEvt_ <<"\n";
00435 
00436   // Get the trigger results
00437   bool validTriggerEvent=true;
00438   edm::Handle<trigger::TriggerEvent> triggerEventHandle;
00439   trigger::TriggerEvent triggerEvent;
00440   e.getByLabel(triggerEvent_,triggerEventHandle);
00441   if(!triggerEventHandle.isValid()) {
00442     edm::LogInfo("PhotonAnalyzer") << "Error! Can't get the product "<< triggerEvent_.label() << endl;
00443     validTriggerEvent=false;
00444   }
00445   if(validTriggerEvent) triggerEvent = *(triggerEventHandle.product());
00446 
00447   // Get the reconstructed photons
00448   bool validPhotons=true;
00449   Handle<reco::PhotonCollection> photonHandle;
00450   reco::PhotonCollection photonCollection;
00451   e.getByLabel(photonProducer_, photonCollection_ , photonHandle);
00452   if ( !photonHandle.isValid()) {
00453     edm::LogInfo("PhotonAnalyzer") << "Error! Can't get the product "<< photonCollection_ << endl;
00454     validPhotons=false;
00455   }
00456   if(validPhotons) photonCollection = *(photonHandle.product());
00457 
00458   // Get the PhotonId objects
00459   bool validloosePhotonID=true;
00460   Handle<edm::ValueMap<bool> > loosePhotonFlag;
00461   edm::ValueMap<bool> loosePhotonID;
00462   e.getByLabel("PhotonIDProd", "PhotonCutBasedIDLoose", loosePhotonFlag);
00463   if ( !loosePhotonFlag.isValid()) {
00464     edm::LogInfo("PhotonAnalyzer") << "Error! Can't get the product "<< "PhotonCutBasedIDLoose" << endl;
00465     validloosePhotonID=false;
00466   }
00467   if (validloosePhotonID) loosePhotonID = *(loosePhotonFlag.product());
00468 
00469   bool validtightPhotonID=true;
00470   Handle<edm::ValueMap<bool> > tightPhotonFlag;
00471   edm::ValueMap<bool> tightPhotonID;
00472   e.getByLabel("PhotonIDProd", "PhotonCutBasedIDTight", tightPhotonFlag);
00473   if ( !tightPhotonFlag.isValid()) {
00474     edm::LogInfo("PhotonAnalyzer") << "Error! Can't get the product "<< "PhotonCutBasedIDTight" << endl;
00475     validtightPhotonID=false;
00476   }
00477   if (validtightPhotonID) tightPhotonID = *(tightPhotonFlag.product());
00478 
00479 
00480 
00481   // Create array to hold #photons/event information
00482   int nPho[100][3][3];
00483 
00484   for (int cut=0; cut!=100; ++cut){
00485     for (int type=0; type!=3; ++type){
00486       for (int part=0; part!=3; ++part){
00487         nPho[cut][type][part] = 0;
00488       }
00489     }
00490   }
00491   // Create array to hold #conversions/event information
00492   int nConv[100][3][3];
00493 
00494   for (int cut=0; cut!=100; ++cut){
00495     for (int type=0; type!=3; ++type){
00496       for (int part=0; part!=3; ++part){
00497         nConv[cut][type][part] = 0;
00498       }
00499     }
00500   }
00501 
00502 
00503 
00504   //Prepare list of photon-related HLT filter names
00505 
00506   vector<int> Keys;
00507 
00508   for(uint filterIndex=0;filterIndex<triggerEvent.sizeFilters();++filterIndex){  //loop over all trigger filters in event (i.e. filters passed)
00509 
00510     string label = triggerEvent.filterTag(filterIndex).label();
00511 
00512     if(label.find( "Photon" ) != string::npos ) {  //get photon-related filters
00513 
00514       for(uint filterKeyIndex=0;filterKeyIndex<triggerEvent.filterKeys(filterIndex).size();++filterKeyIndex){  //loop over keys to objects passing this filter
00515         Keys.push_back(triggerEvent.filterKeys(filterIndex)[filterKeyIndex]);  //add keys to a vector for later reference
00516       }
00517 
00518     }
00519 
00520   }
00521 
00522   // sort Keys vector in ascending order
00523   // and erases duplicate entries from the vector
00524   sort(Keys.begin(),Keys.end());
00525   for ( uint i=0 ; i<Keys.size() ; )
00526    {
00527     if (i!=(Keys.size()-1))
00528      {
00529       if (Keys[i]==Keys[i+1]) Keys.erase(Keys.begin()+i+1) ;
00530       else ++i ;
00531      }
00532     else ++i ;
00533    }
00534 
00535   //We now have a vector of unique keys to TriggerObjects passing a photon-related filter
00536 
00537   int photonCounter = 0;
00538 
00540 
00541   for( reco::PhotonCollection::const_iterator  iPho = photonCollection.begin(); iPho != photonCollection.end(); iPho++) {
00542 
00543 
00544     //for HLT efficiency plots
00545 
00546     h_phoEta_preHLT_->Fill((*iPho).eta());
00547     h_phoEt_preHLT_->Fill( (*iPho).et());
00548 
00549 
00550     double deltaR=1000.;
00551     double deltaRMin=1000.;
00552     double deltaRMax=0.05;//sets deltaR threshold for matching photons to trigger objects
00553 
00554 
00555     for (vector<int>::const_iterator objectKey=Keys.begin();objectKey!=Keys.end();objectKey++){  //loop over keys to objects that fired photon triggers
00556 
00557       deltaR = reco::deltaR(triggerEvent.getObjects()[(*objectKey)].eta(),triggerEvent.getObjects()[(*objectKey)].phi(),(*iPho).superCluster()->eta(),(*iPho).superCluster()->phi());
00558       if(deltaR < deltaRMin) deltaRMin = deltaR;
00559 
00560     }
00561 
00562     if(deltaRMin > deltaRMax) {  //photon fails delta R cut
00563       if(useTriggerFiltering_) continue;  //throw away photons that haven't passed any photon filters
00564     }
00565 
00566     if(deltaRMin <= deltaRMax) { //photon passes delta R cut
00567       h_phoEta_postHLT_->Fill((*iPho).eta() );
00568       h_phoEt_postHLT_->Fill( (*iPho).et() );
00569     }
00570 
00571     if ((*iPho).et()  < minPhoEtCut_) continue;
00572 
00573     nEntry_++;
00574 
00575     edm::Ref<reco::PhotonCollection> photonref(photonHandle, photonCounter);
00576     photonCounter++;
00577 
00578     bool isLoosePhoton(false), isTightPhoton(false);
00579     if ( !isHeavyIon_ ) {
00580       isLoosePhoton = (loosePhotonID)[photonref];
00581       isTightPhoton = (tightPhotonID)[photonref];
00582     }
00583 
00584 
00585     //find out which part of the Ecal contains the photon
00586 
00587     bool  phoIsInBarrel=false;
00588     bool  phoIsInEndcap=false;
00589     float etaPho = (*iPho).superCluster()->eta();
00590     if ( fabs(etaPho) <  1.479 )
00591       phoIsInBarrel=true;
00592     else {
00593       phoIsInEndcap=true;
00594     }
00595 
00596     int part = 0;
00597     if ( phoIsInBarrel ) part = 1;
00598     if ( phoIsInEndcap ) part = 2;
00599 
00601     bool isIsolated=false;
00602     if ( isolationStrength_ == 0)  isIsolated = isLoosePhoton;
00603     if ( isolationStrength_ == 1)  isIsolated = isTightPhoton;
00604 
00605     int type=0;
00606     if ( isIsolated ) type=1;
00607     if ( !isIsolated ) type=2;
00608 
00609 
00610     //get rechit collection containing this photon
00611     bool validEcalRecHits=true;
00612     edm::Handle<EcalRecHitCollection>   ecalRecHitHandle;
00613     EcalRecHitCollection ecalRecHitCollection;
00614     if ( phoIsInBarrel ) {
00615       // Get handle to barrel rec hits
00616       e.getByLabel(barrelRecHitProducer_, barrelRecHitCollection_, ecalRecHitHandle);
00617       if (!ecalRecHitHandle.isValid()) {
00618         edm::LogError("PhotonAnalyzer") << "Error! Can't get the product "<<barrelRecHitProducer_;
00619         validEcalRecHits=false;
00620       }
00621     }
00622     else if ( phoIsInEndcap ) {
00623       // Get handle to endcap rec hits
00624       e.getByLabel(endcapRecHitProducer_, endcapRecHitCollection_, ecalRecHitHandle);
00625       if (!ecalRecHitHandle.isValid()) {
00626         edm::LogError("PhotonAnalyzer") << "Error! Can't get the product "<<endcapRecHitProducer_;
00627         validEcalRecHits=false;
00628       }
00629     }
00630     if (validEcalRecHits) ecalRecHitCollection = *(ecalRecHitHandle.product());
00631 
00632 
00633     //if ((*iPho).isEBEEGap()) continue;  //cut out gap photons
00634 
00635 
00636     //filling histograms to make isolation efficiencies
00637     if(isLoosePhoton){
00638       h_phoEta_Loose_->Fill((*iPho).eta());
00639       h_phoEt_Loose_->Fill( (*iPho).et() );
00640     }
00641     if(isTightPhoton){
00642       h_phoEta_Tight_->Fill((*iPho).eta());
00643       h_phoEt_Tight_->Fill( (*iPho).et() );
00644     }
00645 
00646 
00647 
00648     for (int cut = 0; cut !=numberOfSteps_; ++cut) {  //loop over different transverse energy cuts
00649       double Et =  (*iPho).et();
00650       bool passesCuts = false;
00651 
00652       //sorting the photon into the right Et-dependant folder
00653       if ( useBinning_ && Et > (cut+1)*cutStep_ && ( (Et < (cut+2)*cutStep_)  | (cut == numberOfSteps_-1) ) ){
00654         passesCuts = true;
00655       }
00656       else if ( !useBinning_ && Et > (cut+1)*cutStep_ ){
00657         passesCuts = true;
00658       }
00659 
00660       if (passesCuts){
00661 
00662         //filling isolation variable histograms
00663 
00664         //tracker isolation variables
00665 
00666         fill2DHistoVector(h_nTrackIsolSolid_, (*iPho).nTrkSolidConeDR04(), cut,type);
00667         fill2DHistoVector(h_nTrackIsolHollow_,(*iPho).nTrkHollowConeDR04(),cut,type);
00668 
00669         fill2DHistoVector(h_nTrackIsolSolidVsEta_, (*iPho).eta(),(*iPho).nTrkSolidConeDR04(), cut,type);
00670         fill2DHistoVector(p_nTrackIsolSolidVsEta_, (*iPho).eta(),(*iPho).nTrkSolidConeDR04(), cut,type);
00671         fill2DHistoVector(h_nTrackIsolHollowVsEta_,(*iPho).eta(),(*iPho).nTrkHollowConeDR04(),cut,type);
00672         fill2DHistoVector(p_nTrackIsolHollowVsEta_,(*iPho).eta(),(*iPho).nTrkHollowConeDR04(),cut,type);
00673 
00674         fill2DHistoVector(h_nTrackIsolSolidVsEt_,  (*iPho).et(), (*iPho).nTrkSolidConeDR04(), cut,type);
00675         fill2DHistoVector(p_nTrackIsolSolidVsEt_,  (*iPho).et(), (*iPho).nTrkSolidConeDR04(), cut,type);
00676         fill2DHistoVector(h_nTrackIsolHollowVsEt_, (*iPho).et(), (*iPho).nTrkHollowConeDR04(),cut,type);
00677         fill2DHistoVector(p_nTrackIsolHollowVsEt_, (*iPho).et(), (*iPho).nTrkHollowConeDR04(),cut,type);
00678 
00680 
00681         fill2DHistoVector(h_trackPtSumSolid_, (*iPho).trkSumPtSolidConeDR04(),cut,type);
00682         fill2DHistoVector(h_trackPtSumHollow_,(*iPho).trkSumPtSolidConeDR04(),cut,type);
00683 
00684         fill2DHistoVector(h_trackPtSumSolidVsEta_, (*iPho).eta(),(*iPho).trkSumPtSolidConeDR04(), cut,type);
00685         fill2DHistoVector(p_trackPtSumSolidVsEta_, (*iPho).eta(),(*iPho).trkSumPtSolidConeDR04(), cut,type);
00686         fill2DHistoVector(h_trackPtSumHollowVsEta_,(*iPho).eta(),(*iPho).trkSumPtHollowConeDR04(),cut,type);
00687         fill2DHistoVector(p_trackPtSumHollowVsEta_,(*iPho).eta(),(*iPho).trkSumPtHollowConeDR04(),cut,type);
00688 
00689         fill2DHistoVector(h_trackPtSumSolidVsEt_,  (*iPho).et(), (*iPho).trkSumPtSolidConeDR04(), cut,type);
00690         fill2DHistoVector(p_trackPtSumSolidVsEt_,  (*iPho).et(), (*iPho).trkSumPtSolidConeDR04(), cut,type);
00691         fill2DHistoVector(h_trackPtSumHollowVsEt_, (*iPho).et(), (*iPho).trkSumPtHollowConeDR04(),cut,type);
00692         fill2DHistoVector(p_trackPtSumHollowVsEt_, (*iPho).et(), (*iPho).trkSumPtHollowConeDR04(),cut,type);
00693 
00694         //calorimeter isolation variables
00695 
00696         fill2DHistoVector(h_ecalSum_,(*iPho).ecalRecHitSumEtConeDR04(),cut,type);
00697         fill2DHistoVector(h_ecalSumVsEta_,(*iPho).eta(),(*iPho).ecalRecHitSumEtConeDR04(),cut,type);
00698         fill2DHistoVector(p_ecalSumVsEta_,(*iPho).eta(),(*iPho).ecalRecHitSumEtConeDR04(),cut,type);
00699         fill2DHistoVector(h_ecalSumVsEt_, (*iPho).et(), (*iPho).ecalRecHitSumEtConeDR04(),cut,type);
00700         fill3DHistoVector(p_ecalSumVsEt_, (*iPho).et(), (*iPho).ecalRecHitSumEtConeDR04(),cut,type,part);
00701 
00703 
00704         fill2DHistoVector(h_hcalSum_,(*iPho).hcalTowerSumEtConeDR04(),cut,type);
00705         fill2DHistoVector(h_hcalSumVsEta_,(*iPho).eta(),(*iPho).hcalTowerSumEtConeDR04(),cut,type);
00706         fill2DHistoVector(p_hcalSumVsEta_,(*iPho).eta(),(*iPho).hcalTowerSumEtConeDR04(),cut,type);
00707         fill2DHistoVector(h_hcalSumVsEt_, (*iPho).et(), (*iPho).hcalTowerSumEtConeDR04(),cut,type);
00708         fill3DHistoVector(p_hcalSumVsEt_, (*iPho).et(), (*iPho).hcalTowerSumEtConeDR04(),cut,type,part);
00709 
00710         fill3DHistoVector(h_hOverE_,(*iPho).hadronicOverEm(),cut,type,part);
00711         fill2DHistoVector(p_hOverEVsEta_,(*iPho).eta(),(*iPho).hadronicOverEm(),cut,type);
00712         fill2DHistoVector(p_hOverEVsEt_, (*iPho).et(), (*iPho).hadronicOverEm(),cut,type);
00713 
00714         fill3DHistoVector(h_h1OverE_,(*iPho).hadronicDepth1OverEm(),cut,type,part);
00715         fill3DHistoVector(h_h2OverE_,(*iPho).hadronicDepth2OverEm(),cut,type,part);
00716 
00717 
00718         //filling photon histograms
00719 
00720         nPho[cut][0][0]++;
00721         nPho[cut][0][part]++;
00722         nPho[cut][type][0]++;
00723         nPho[cut][type][part]++;
00724 
00725         //energy variables
00726 
00727         fill3DHistoVector(h_phoE_, (*iPho).energy(),cut,type,part);
00728         fill3DHistoVector(h_phoEt_,(*iPho).et(),    cut,type,part);
00729 
00730         //geometrical variables
00731 
00732         fill2DHistoVector(h_phoEta_,(*iPho).eta(),cut,type);
00733         fill2DHistoVector(h_scEta_, (*iPho).superCluster()->eta(),cut,type);
00734 
00735         fill3DHistoVector(h_phoPhi_,(*iPho).phi(),cut,type,part);
00736         fill3DHistoVector(h_scPhi_, (*iPho).superCluster()->phi(),cut,type,part);
00737 
00738         //shower shape variables
00739 
00740         fill3DHistoVector(h_r9_,(*iPho).r9(),cut,type,part);
00741         fill2DHistoVector(h_r9VsEta_,(*iPho).eta(),(*iPho).r9(),cut,type);
00742         fill2DHistoVector(p_r9VsEta_,(*iPho).eta(),(*iPho).r9(),cut,type);
00743         fill2DHistoVector(h_r9VsEt_, (*iPho).et(), (*iPho).r9(),cut,type);
00744         fill2DHistoVector(p_r9VsEt_, (*iPho).et(), (*iPho).r9(),cut,type);
00745 
00746         fill2DHistoVector(h_e1x5VsEta_,(*iPho).eta(),(*iPho).e1x5(),cut,type);
00747         fill2DHistoVector(p_e1x5VsEta_,(*iPho).eta(),(*iPho).e1x5(),cut,type);
00748         fill2DHistoVector(h_e1x5VsEt_, (*iPho).et(), (*iPho).e1x5(),cut,type);
00749         fill2DHistoVector(p_e1x5VsEt_, (*iPho).et(), (*iPho).e1x5(),cut,type);
00750 
00751         fill2DHistoVector(h_e2x5VsEta_,(*iPho).eta(),(*iPho).e2x5(),cut,type);
00752         fill2DHistoVector(p_e2x5VsEta_,(*iPho).eta(),(*iPho).e2x5(),cut,type);
00753         fill2DHistoVector(h_e2x5VsEt_, (*iPho).et(), (*iPho).e2x5(),cut,type);
00754         fill2DHistoVector(p_e2x5VsEt_, (*iPho).et(), (*iPho).e2x5(),cut,type);
00755 
00756         fill2DHistoVector(h_maxEXtalOver3x3VsEta_,(*iPho).eta(),(*iPho).maxEnergyXtal()/(*iPho).e3x3(),cut,type);
00757         fill2DHistoVector(p_maxEXtalOver3x3VsEta_,(*iPho).eta(),(*iPho).maxEnergyXtal()/(*iPho).e3x3(),cut,type);
00758         fill2DHistoVector(h_maxEXtalOver3x3VsEt_, (*iPho).et(), (*iPho).maxEnergyXtal()/(*iPho).e3x3(),cut,type);
00759         fill2DHistoVector(p_maxEXtalOver3x3VsEt_, (*iPho).et(), (*iPho).maxEnergyXtal()/(*iPho).e3x3(),cut,type);
00760 
00761         fill2DHistoVector(h_r1x5VsEta_,(*iPho).eta(),(*iPho).r1x5(),cut,type);
00762         fill2DHistoVector(p_r1x5VsEta_,(*iPho).eta(),(*iPho).r1x5(),cut,type);
00763         fill2DHistoVector(h_r1x5VsEt_, (*iPho).et(), (*iPho).r1x5(),cut,type);
00764         fill2DHistoVector(p_r1x5VsEt_, (*iPho).et(), (*iPho).r1x5(),cut,type);
00765 
00766         fill2DHistoVector(h_r2x5VsEta_,(*iPho).eta(),(*iPho).r2x5(),cut,type);
00767         fill2DHistoVector(p_r2x5VsEta_,(*iPho).eta(),(*iPho).r2x5(),cut,type);
00768         fill2DHistoVector(h_r2x5VsEt_, (*iPho).et(), (*iPho).r2x5(),cut,type);
00769         fill2DHistoVector(p_r2x5VsEt_, (*iPho).et(), (*iPho).r2x5(),cut,type);
00770 
00771         fill3DHistoVector(h_phoSigmaIetaIeta_,(*iPho).sigmaIetaIeta(),cut,type,part);
00772         fill2DHistoVector(h_sigmaIetaIetaVsEta_,(*iPho).eta(),(*iPho).sigmaIetaIeta(),cut,type);
00773         fill2DHistoVector(p_sigmaIetaIetaVsEta_,(*iPho).eta(),(*iPho).sigmaIetaIeta(),cut,type);
00774 
00775 
00776 
00777         //filling histograms for photons containing a bad ECAL channel
00778 
00779         bool atLeastOneDeadChannel=false;
00780         for(reco::CaloCluster_iterator bcIt = (*iPho).superCluster()->clustersBegin();bcIt != (*iPho).superCluster()->clustersEnd(); ++bcIt) { //loop over basic clusters in SC
00781           for(vector< pair<DetId, float> >::const_iterator rhIt = (*bcIt)->hitsAndFractions().begin();rhIt != (*bcIt)->hitsAndFractions().end(); ++rhIt) { //loop over rec hits in basic cluster
00782 
00783             for(EcalRecHitCollection::const_iterator it = ecalRecHitCollection.begin(); it !=  ecalRecHitCollection.end(); ++it) { //loop over all rec hits to find the right ones
00784               if  (rhIt->first ==  (*it).id() ) { //found the matching rechit
00785                 if (  (*it).recoFlag() == 9 ) { //has a bad channel
00786                   atLeastOneDeadChannel=true;
00787                   break;
00788                 }
00789               }
00790             }
00791           }
00792         }
00793         if ( atLeastOneDeadChannel ) {
00794           fill2DHistoVector(h_phoPhi_BadChannels_,(*iPho).phi(),cut,type);
00795           fill2DHistoVector(h_phoEta_BadChannels_,(*iPho).eta(),cut,type);
00796           fill2DHistoVector(h_phoEt_BadChannels_, (*iPho).et(), cut,type);
00797         }
00798 
00799 
00800         // filling conversion-related histograms
00801         if((*iPho).hasConversionTracks()){
00802           nConv[cut][0][0]++;
00803           nConv[cut][0][part]++;
00804           nConv[cut][type][0]++;
00805           nConv[cut][type][part]++;
00806         }
00807 
00808         //loop over conversions (don't forget, we're still inside the photon loop,
00809         // i.e. these are all the conversions for this ONE photon, not for all the photons in the event)
00810         reco::ConversionRefVector conversions = (*iPho).conversions();
00811         for (unsigned int iConv=0; iConv<conversions.size(); iConv++) {
00812 
00813           reco::ConversionRef aConv=conversions[iConv];
00814 
00815           if ( aConv->nTracks() <2 ) continue;
00816 
00817           //fill histogram for denominator of vertex reconstruction efficiency plot
00818           if(cut==0) h_phoEta_Vertex_->Fill(aConv->refittedPairMomentum().eta());
00819 
00820           if ( !(aConv->conversionVertex().isValid()) ) continue;
00821 
00822           float chi2Prob = ChiSquaredProbability( aConv->conversionVertex().chi2(), aConv->conversionVertex().ndof() );
00823 
00824           if(chi2Prob<0.0005) continue;
00825 
00826           fill2DHistoVector(h_vertexChi2Prob_,chi2Prob,cut,type);
00827 
00828 
00829 
00830           fill3DHistoVector(h_phoConvE_, (*iPho).energy(),cut,type,part);
00831           fill3DHistoVector(h_phoConvEt_,(*iPho).et(),cut,type,part);
00832           fill3DHistoVector(h_phoConvR9_,(*iPho).r9(),cut,type,part);
00833 
00834           if (cut==0 && isLoosePhoton){
00835             h_convEta_Loose_->Fill((*iPho).eta());
00836             h_convEt_Loose_->Fill( (*iPho).et() );
00837           }
00838           if (cut==0 && isTightPhoton){
00839             h_convEta_Tight_->Fill((*iPho).eta());
00840             h_convEt_Tight_->Fill( (*iPho).et() );
00841           }
00842 
00843           fill2DHistoVector(h_phoConvEta_,aConv->refittedPairMomentum().eta(),cut,type);
00844           fill3DHistoVector(h_phoConvPhi_,aConv->refittedPairMomentum().phi(),cut,type,part);
00845 
00846           //we use the photon position because we'll be dividing it by a photon histogram (not a conversion histogram)
00847           fill2DHistoVector(h_phoConvEtaForEfficiency_,(*iPho).eta(),cut,type);
00848           fill3DHistoVector(h_phoConvPhiForEfficiency_,(*iPho).phi(),cut,type,part);
00849 
00850 
00851           //vertex histograms
00852           double convR= sqrt(aConv->conversionVertex().position().perp2());
00853           double scalar = aConv->conversionVertex().position().x()*aConv->refittedPairMomentum().x() + aConv->conversionVertex().position().y()*aConv->refittedPairMomentum().y();
00854           if ( scalar < 0 ) convR= -convR;
00855 
00856           fill2DHistoVector(h_convVtxRvsZ_,aConv->conversionVertex().position().z(), convR,cut,type);//trying to "see" R-Z view of tracker
00857           fill2DHistoVector(h_convVtxZ_,aConv->conversionVertex().position().z(), cut,type);
00858 
00859           if(fabs(aConv->caloCluster()[0]->eta()) > 1.5){//trying to "see" tracker endcaps
00860             fill2DHistoVector(h_convVtxZEndcap_,aConv->conversionVertex().position().z(), cut,type);
00861           }
00862           else if(fabs(aConv->caloCluster()[0]->eta()) < 1){//trying to "see" tracker barrel
00863             fill2DHistoVector(h_convVtxR_,convR,cut,type);
00864             fill2DHistoVector(h_convVtxYvsX_,aConv->conversionVertex().position().x(),aConv->conversionVertex().position().y(),cut,type);
00865           }
00866 
00867           const std::vector<edm::RefToBase<reco::Track> > tracks = aConv->tracks();
00868 
00869 
00870           for (unsigned int i=0; i<tracks.size(); i++) {
00871             fill2DHistoVector(h_tkChi2_,tracks[i]->normalizedChi2(),cut,type);
00872             fill2DHistoVector(p_tkChi2VsEta_,aConv->caloCluster()[0]->eta(),tracks[i]->normalizedChi2(),cut,type);
00873             fill2DHistoVector(p_dCotTracksVsEta_,aConv->caloCluster()[0]->eta(),aConv->pairCotThetaSeparation(),cut,type);
00874             fill2DHistoVector(p_nHitsVsEta_,aConv->caloCluster()[0]->eta(),float(tracks[i]->numberOfValidHits()),cut,type);
00875           }
00876 
00877           //calculating delta eta and delta phi of the two tracks
00878 
00879           float  DPhiTracksAtVtx = -99;
00880           float  dPhiTracksAtEcal= -99;
00881           float  dEtaTracksAtEcal= -99;
00882 
00883           float phiTk1= aConv->tracksPin()[0].phi();
00884           float phiTk2= aConv->tracksPin()[1].phi();
00885           DPhiTracksAtVtx = phiTk1-phiTk2;
00886           DPhiTracksAtVtx = phiNormalization( DPhiTracksAtVtx );
00887 
00888           if (aConv->bcMatchingWithTracks()[0].isNonnull() && aConv->bcMatchingWithTracks()[1].isNonnull() ) {
00889             float recoPhi1 = aConv->ecalImpactPosition()[0].phi();
00890             float recoPhi2 = aConv->ecalImpactPosition()[1].phi();
00891             float recoEta1 = aConv->ecalImpactPosition()[0].eta();
00892             float recoEta2 = aConv->ecalImpactPosition()[1].eta();
00893 
00894             recoPhi1 = phiNormalization(recoPhi1);
00895             recoPhi2 = phiNormalization(recoPhi2);
00896 
00897             dPhiTracksAtEcal = recoPhi1 -recoPhi2;
00898             dPhiTracksAtEcal = phiNormalization( dPhiTracksAtEcal );
00899             dEtaTracksAtEcal = recoEta1 -recoEta2;
00900 
00901           }
00902 
00903 
00904           fill3DHistoVector(h_dPhiTracksAtVtx_,DPhiTracksAtVtx,cut,type,part);
00905           fill3DHistoVector(h_dPhiTracksAtEcal_,fabs(dPhiTracksAtEcal),cut,type,part);
00906           fill3DHistoVector(h_dEtaTracksAtEcal_, dEtaTracksAtEcal,cut,type,part);
00907           fill3DHistoVector(h_eOverPTracks_,aConv->EoverPrefittedTracks(),cut,type,part);
00908           fill3DHistoVector(h_pOverETracks_,1./aConv->EoverPrefittedTracks(),cut,type,part);
00909           fill3DHistoVector(h_dCotTracks_,aConv->pairCotThetaSeparation(),cut,type,part);
00910 
00911 
00912         }//end loop over conversions
00913 
00914       }//end loop over photons passing cuts
00915     }//end loop over transverse energy cuts
00916 
00917 
00918 
00919 
00920 
00921     //make invariant mass plots
00922 
00923     if (isIsolated && iPho->et()>=invMassEtCut_){
00924 
00925       for (reco::PhotonCollection::const_iterator iPho2=iPho+1; iPho2!=photonCollection.end(); iPho2++){
00926 
00927         edm::Ref<reco::PhotonCollection> photonref2(photonHandle, photonCounter); //note: it's correct to use photonCounter and not photonCounter+1
00928                                                                                   //since it has already been incremented earlier
00929 
00930         bool  isTightPhoton2(false), isLoosePhoton2(false);
00931 
00932         if ( !isHeavyIon_ ) {
00933           isTightPhoton2 = (tightPhotonID)[photonref2];
00934           isLoosePhoton2 = (loosePhotonID)[photonref2];
00935         }
00936 
00937         bool isIsolated2=false;
00938         if ( isolationStrength_ == 0)  isIsolated2 = isLoosePhoton2;
00939         if ( isolationStrength_ == 1)  isIsolated2 = isTightPhoton2;
00940 
00941         reco::ConversionRefVector conversions = (*iPho).conversions();
00942         reco::ConversionRefVector conversions2 = (*iPho2).conversions();
00943 
00944         if(isIsolated2 && iPho2->et()>=invMassEtCut_){
00945 
00946           math::XYZTLorentzVector p12 = iPho->p4()+iPho2->p4();
00947           float gamgamMass2 = p12.Dot(p12);
00948 
00949 
00950           h_invMassAllPhotons_ -> Fill(sqrt( gamgamMass2 ));
00951 
00952           if(conversions.size()!=0 && conversions[0]->nTracks() >= 2){
00953             if(conversions2.size()!=0 && conversions2[0]->nTracks() >= 2) h_invMassTwoWithTracks_ -> Fill(sqrt( gamgamMass2 ));
00954             else h_invMassOneWithTracks_ -> Fill(sqrt( gamgamMass2 ));
00955           }
00956           else if(conversions2.size()!=0 && conversions2[0]->nTracks() >= 2) h_invMassOneWithTracks_ -> Fill(sqrt( gamgamMass2 ));
00957           else h_invMassZeroWithTracks_ -> Fill(sqrt( gamgamMass2 ));
00958         }
00959 
00960       }
00961 
00962     }
00963 
00964 
00965 
00966   }
00967 
00968 
00969   //filling number of photons/conversions per event histograms
00970   for (int cut = 0; cut !=numberOfSteps_; ++cut) {
00971     for(uint type=0;type!=types_.size();++type){
00972       for(uint part=0;part!=parts_.size();++part){
00973         h_nPho_[cut][type][part]-> Fill (float(nPho[cut][type][part]));
00974         h_nConv_[cut][type][part]-> Fill (float(nConv[cut][type][part]));
00975       }
00976     }
00977   }
00978 
00979 }//End of Analyze method
00980 
00981 void PhotonAnalyzer::endRun(const edm::Run& run, const edm::EventSetup& setup)
00982 {
00983   if(!standAlone_){
00984 
00985     dbe_->setCurrentFolder("Egamma/PhotonAnalyzer");
00986 
00987     //keep track of how many histos are in each folder
00988     totalNumberOfHistos_efficiencyFolder->Fill(histo_index_efficiency_);
00989     totalNumberOfHistos_invMassFolder->Fill(histo_index_invMass_);
00990     totalNumberOfHistos_photonsFolder->Fill(histo_index_photons_);
00991     totalNumberOfHistos_conversionsFolder->Fill(histo_index_conversions_);
00992 
00993   }
00994 
00995 }
00996 
00997 
00998 void PhotonAnalyzer::endJob()
00999 {
01000   //dbe_->showDirStructure();
01001   if(standAlone_){
01002     dbe_->setCurrentFolder("Egamma/PhotonAnalyzer");
01003 
01004     //keep track of how many histos are in each folder
01005     totalNumberOfHistos_efficiencyFolder->Fill(histo_index_efficiency_);
01006     totalNumberOfHistos_invMassFolder->Fill(histo_index_invMass_);
01007     totalNumberOfHistos_photonsFolder->Fill(histo_index_photons_);
01008     totalNumberOfHistos_conversionsFolder->Fill(histo_index_conversions_);
01009 
01010 
01011     dbe_->save(outputFileName_);
01012   }
01013 
01014 
01015 }
01016 
01018 
01019 
01020 
01021 float PhotonAnalyzer::phiNormalization(float & phi)
01022 {
01023  const float PI    = 3.1415927;
01024  const float TWOPI = 2.0*PI;
01025 
01026  if(phi >  PI) {phi = phi - TWOPI;}
01027  if(phi < -PI) {phi = phi + TWOPI;}
01028 
01029  return phi;
01030 }
01031 
01032 
01033 void  PhotonAnalyzer::fill2DHistoVector(vector<vector<MonitorElement*> >& histoVector,double x, double y, int cut, int type){
01034 
01035   histoVector[cut][0]->Fill(x,y);
01036   if(histoVector[cut].size()>1)   histoVector[cut][type]->Fill(x,y); //don't try to fill 2D histos that are only in the "AllPhotons" folder
01037 
01038 }
01039 
01040 void  PhotonAnalyzer::fill2DHistoVector(vector<vector<MonitorElement*> >& histoVector, double x, int cut, int type){
01041 
01042   histoVector[cut][0]->Fill(x);
01043   histoVector[cut][type]->Fill(x);
01044 
01045 }
01046 
01047 void  PhotonAnalyzer::fill3DHistoVector(vector<vector<vector<MonitorElement*> > >& histoVector,double x, int cut, int type, int part){
01048 
01049   histoVector[cut][0][0]->Fill(x);
01050   histoVector[cut][0][part]->Fill(x);
01051   histoVector[cut][type][0]->Fill(x);
01052   histoVector[cut][type][part]->Fill(x);
01053 
01054 }
01055 
01056 void  PhotonAnalyzer::fill3DHistoVector(vector<vector<vector<MonitorElement*> > >& histoVector,double x, double y, int cut, int type, int part){
01057 
01058   histoVector[cut][0][0]->Fill(x,y);
01059   histoVector[cut][0][part]->Fill(x,y);
01060   histoVector[cut][type][0]->Fill(x,y);
01061   histoVector[cut][type][part]->Fill(x,y);
01062 
01063 }
01064 
01065 
01066 MonitorElement* PhotonAnalyzer::bookHisto(string histoName, string title, int bin, double min, double max)
01067 {
01068 
01069   int histo_index = 0;
01070   stringstream histo_number_stream;
01071 
01072   //determining which folder we're in
01073   if(dbe_->pwd().find( "InvMass" ) != string::npos){
01074     histo_index_invMass_++;
01075     histo_index = histo_index_invMass_;
01076   }
01077   if(dbe_->pwd().find( "Efficiencies" ) != string::npos){
01078     histo_index_efficiency_++;
01079     histo_index = histo_index_efficiency_;
01080   }
01081 
01082   histo_number_stream << "h_";
01083   if(histo_index<10)   histo_number_stream << "0";
01084   histo_number_stream << histo_index;
01085 
01086   return dbe_->book1D(histo_number_stream.str()+"_"+histoName,title,bin,min,max);
01087 
01088 }
01089 
01090 
01091 vector<vector<MonitorElement*> > PhotonAnalyzer::book2DHistoVector(string histoType, string histoName, string title,
01092                                                                              int xbin, double xmin,double xmax,
01093                                                                              int ybin, double ymin, double ymax)
01094 {
01095   int histo_index = 0;
01096 
01097   vector<MonitorElement*> temp1DVector;
01098   vector<vector<MonitorElement*> > temp2DVector;
01099 
01100   //determining which folder we're in
01101   bool conversionPlot = false;
01102   if(dbe_->pwd().find( "Conversions" ) != string::npos) conversionPlot = true;
01103   bool TwoDPlot = false;
01104   if(histoName.find( "2D" ) != string::npos) TwoDPlot = true;
01105 
01106   if(conversionPlot){
01107     histo_index_conversions_++;
01108     histo_index = histo_index_conversions_;
01109   }
01110   else{
01111     histo_index_photons_++;
01112     histo_index = histo_index_photons_;
01113   }
01114 
01115   stringstream histo_number_stream;
01116   histo_number_stream << "h_";
01117   if(histo_index<10)   histo_number_stream << "0";
01118   histo_number_stream << histo_index << "_";
01119 
01120 
01121   for(int cut = 0; cut != numberOfSteps_; ++cut){ //looping over Et cut values
01122 
01123     for(uint type=0;type!=types_.size();++type){  //looping over isolation type
01124 
01125       currentFolder_.str("");
01126       currentFolder_ << "Egamma/PhotonAnalyzer/" << types_[type] << "Photons/Et above " << (cut+1)*cutStep_ << " GeV";
01127       if(conversionPlot) currentFolder_ << "/Conversions";
01128 
01129       dbe_->setCurrentFolder(currentFolder_.str());
01130 
01131       string kind;
01132       if(conversionPlot) kind = " Conversions: ";
01133       else kind = " Photons: ";
01134 
01135       if(histoType=="1D")           temp1DVector.push_back(dbe_->book1D(      histo_number_stream.str()+histoName,types_[type]+kind+title,xbin,xmin,xmax));
01136       else if(histoType=="2D"){
01137         if((TwoDPlot && type==0) || !TwoDPlot){//only book the 2D plots in the "AllPhotons" folder
01138                                     temp1DVector.push_back(dbe_->book2D(      histo_number_stream.str()+histoName,types_[type]+kind+title,xbin,xmin,xmax,ybin,ymin,ymax));
01139         }
01140       }
01141       else if(histoType=="Profile") temp1DVector.push_back(dbe_->bookProfile( histo_number_stream.str()+histoName,types_[type]+kind+title,xbin,xmin,xmax,ybin,ymin,ymax,""));
01142       else cout << "bad histoType\n";
01143     }
01144 
01145     temp2DVector.push_back(temp1DVector);
01146     temp1DVector.clear();
01147   }
01148 
01149   return temp2DVector;
01150 
01151 }
01152 
01153 
01154 vector<vector<vector<MonitorElement*> > > PhotonAnalyzer::book3DHistoVector(string histoType, string histoName, string title,
01155                                                                              int xbin, double xmin,double xmax,
01156                                                                              int ybin, double ymin, double ymax)
01157 {
01158   int histo_index = 0;
01159 
01160   vector<MonitorElement*> temp1DVector;
01161   vector<vector<MonitorElement*> > temp2DVector;
01162   vector<vector<vector<MonitorElement*> > > temp3DVector;
01163 
01164 
01165   //determining which folder we're in
01166   bool conversionPlot = false;
01167   if(dbe_->pwd().find( "Conversions" ) != string::npos) conversionPlot = true;
01168 
01169 
01170   if(conversionPlot){
01171     histo_index_conversions_++;
01172     histo_index = histo_index_conversions_;
01173   }
01174   else{
01175     histo_index_photons_++;
01176     histo_index = histo_index_photons_;
01177   }
01178 
01179   stringstream histo_number_stream;
01180   histo_number_stream << "h_";
01181   if(histo_index<10)   histo_number_stream << "0";
01182   histo_number_stream << histo_index << "_";
01183 
01184 
01185   for(int cut = 0; cut != numberOfSteps_; ++cut){     //looping over Et cut values
01186 
01187     for(uint type=0;type!=types_.size();++type){      //looping over isolation type
01188 
01189       for(uint part=0;part!=parts_.size();++part){    //looping over different parts of the ecal
01190 
01191         currentFolder_.str("");
01192         currentFolder_ << "Egamma/PhotonAnalyzer/" << types_[type] << "Photons/Et above " << (cut+1)*cutStep_ << " GeV";
01193         if(conversionPlot) currentFolder_ << "/Conversions";
01194 
01195         dbe_->setCurrentFolder(currentFolder_.str());
01196 
01197         string kind;
01198         if(conversionPlot) kind = " Conversions: ";
01199         else kind = " Photons: ";
01200 
01201         if(histoType=="1D")           temp1DVector.push_back(dbe_->book1D(      histo_number_stream.str()+histoName+parts_[part],types_[type]+kind+parts_[part]+": "+title,xbin,xmin,xmax));
01202         else if(histoType=="2D")      temp1DVector.push_back(dbe_->book2D(      histo_number_stream.str()+histoName+parts_[part],types_[type]+kind+parts_[part]+": "+title,xbin,xmin,xmax,ybin,ymin,ymax));
01203         else if(histoType=="Profile") temp1DVector.push_back(dbe_->bookProfile( histo_number_stream.str()+histoName+parts_[part],types_[type]+kind+parts_[part]+": "+title,xbin,xmin,xmax,ybin,ymin,ymax,""));
01204         else cout << "bad histoType\n";
01205 
01206 
01207       }
01208 
01209       temp2DVector.push_back(temp1DVector);
01210       temp1DVector.clear();
01211     }
01212 
01213     temp3DVector.push_back(temp2DVector);
01214     temp2DVector.clear();
01215   }
01216 
01217   return temp3DVector;
01218 }
/data/refman/pasoursint/CMSSW_4_1_8_patch13/src/DQMOffline/EGamma/plugins/PhotonAnalyzer.cc
