CMS 3D CMS Logo

/afs/cern.ch/work/a/aaltunda/public/www/CMSSW_5_3_14/src/DQMOffline/EGamma/plugins/PhotonAnalyzer.cc

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