From bceb401be55f10a6260b26273f89a5a5a5658b1c Mon Sep 17 00:00:00 2001 From: Sunanda Date: Fri, 24 Jul 2020 19:46:09 +0200 Subject: [PATCH 1/3] Modify the scripts for study radiation damage using muons --- .../HcalCalibAlgos/macros/AnalyzeLepTree.C | 102 +++++++++++++----- .../HcalCalibAlgos/macros/HBHEMuonHighEta.C | 96 ++++++++++++++--- .../macros/HBHEMuonOfflineAnalyzer.C | 34 +++--- 3 files changed, 174 insertions(+), 58 deletions(-) diff --git a/Calibration/HcalCalibAlgos/macros/AnalyzeLepTree.C b/Calibration/HcalCalibAlgos/macros/AnalyzeLepTree.C index d50d81acc0ade..327a1aca015ae 100644 --- a/Calibration/HcalCalibAlgos/macros/AnalyzeLepTree.C +++ b/Calibration/HcalCalibAlgos/macros/AnalyzeLepTree.C @@ -1,5 +1,5 @@ ////////////////////////////////////////////////////////////////////////////// -// This class analyzes the "Lep Tree" created by HBHEOfflineAnalyzer +// This class analyzes the "Lep Tree" created by HBHEMuonOfflineAnalyzer // It has two constructors using either a pointer to the tree chain or // the file name where the tree resides // There are 2 additional arguments: @@ -284,6 +284,8 @@ void AnalyzeLepTree::Init(TChain *tree) { fChain->SetBranchAddress("t_depth", &t_depth, &b_t_depth); Notify(); + t_ediff = 0; + bookHisto(); } @@ -342,6 +344,7 @@ void AnalyzeLepTree::Loop(Long64_t nmax, bool debug) { const double ethr = 0.00001; // Threshold of 10 keV Long64_t nbytes = 0, nb = 0; + int32_t n15(0), n16(0); for (Long64_t jentry=0; jentry 0 && pbin >= 0 && vbin >= 0) { if (kdepth_ == 0) { for (unsigned int k=0; ksize(); ++k) { + if (eta == 15) ++n15; + else if (eta == 16) ++n16; int depth = (*t_depth)[k]; unsigned int id = packID(zside,eta,phi,depth+1,vbin,pbin); double ene = (*t_ene)[k]; @@ -412,29 +417,50 @@ void AnalyzeLepTree::Loop(Long64_t nmax, bool debug) { } } } else if (kdepth_ == 1) { - double ene(0), enec(0), actl(0), charge(0); - unsigned int id = packID(zside,eta,phi,1,vbin,pbin); + double ene[2], enec[2], actl[2], charge[2]; + for (unsigned int k=0; k<2; ++k) { + ene[k] = enec[k] = actl[k] = charge[k] = 0; + } for (unsigned int k=0; ksize(); ++k) { if ((*t_ene)[k] > 0 && (*t_actln)[k] > 0) { - ene += (*t_ene)[k]; - enec += (*t_enec)[k]; - charge += (*t_charge)[k]; - actl += (*t_actln)[k]; + int dep = (*t_depth)[k]; + int depth = ((eta != 16) ? 0 : ((dep > 1) ? 1 : 0)); + ene[depth] += (*t_ene)[k]; + enec[depth] += (*t_enec)[k]; + charge[depth] += (*t_charge)[k]; + actl[depth] += (*t_actln)[k]; } } - if (ene > ethr && actl > 0 && charge > 0 && t_ediff < cutEdiff_) { - std::map::iterator it1 = h_Energy_.find(id); - if (it1 != h_Energy_.end()) (it1->second)->Fill(ene); - std::map::iterator it2 = h_Ecorr_.find(id); - if (it2 != h_Ecorr_.end()) (it2->second)->Fill(ene/actl); - std::map::iterator it3 = h_EnergyC_.find(id); - if (it3 != h_EnergyC_.end()) (it3->second)->Fill(enec); - std::map::iterator it4 = h_EcorrC_.find(id); - if (it4 != h_EcorrC_.end()) (it4->second)->Fill(enec/actl); - std::map::iterator it5 = h_Charge_.find(id); - if (it5 != h_Charge_.end()) (it5->second)->Fill(charge); - std::map::iterator it6 = h_Chcorr_.find(id); - if (it6 != h_Chcorr_.end()) (it6->second)->Fill(charge/actl); + int nDepth = (eta == 16) ? 2 : 1; + for (int k=0; k ethr && actl[k] > 0 && charge[k] > 0 && + t_ediff < cutEdiff_) { + if (eta == 15) ++n15; + else if (eta == 16) ++n16; + int depth = k + 1; + unsigned int id = packID(zside,eta,phi,depth,vbin,pbin); + std::map::iterator it1 = h_Energy_.find(id); + if (it1 != h_Energy_.end()) (it1->second)->Fill(ene[k]); + std::map::iterator it2 = h_Ecorr_.find(id); + if (it2 != h_Ecorr_.end()) (it2->second)->Fill(ene[k]/actl[k]); + std::map::iterator it3 = h_EnergyC_.find(id); + if (it3 != h_EnergyC_.end()) (it3->second)->Fill(enec[k]); + std::map::iterator it4 = h_EcorrC_.find(id); + if (it4 != h_EcorrC_.end()) (it4->second)->Fill(enec[k]/actl[k]); + std::map::iterator it5 = h_Charge_.find(id); + if (it5 != h_Charge_.end()) (it5->second)->Fill(charge[k]); + std::map::iterator it6 = h_Chcorr_.find(id); + if (it6 != h_Chcorr_.end()) (it6->second)->Fill(charge[k]/actl[k]); + if (((eta == 15) || (eta == 16)) && debug_) + std::cout << zside << ":" << eta << ":" << phi << ":" + << t_iphi << ":" << depth << ":" << vbin << ":" + << pbin << " ID " << std::hex << id << std::dec + << " Flags " << (it1 != h_Energy_.end()) << ":" + << (it2 != h_Ecorr_.end()) << ":" + << (it3 != h_Charge_.end()) << ":" + << (it4 != h_Chcorr_.end()) << " E " << ene + << " C " << charge << " L " << actl << std::endl; + } } } else { double ene[3], enec[3], actl[3], charge[3]; @@ -454,7 +480,10 @@ void AnalyzeLepTree::Loop(Long64_t nmax, bool debug) { for (int k=0; k ethr && actl[k] > 0 && charge[k] > 0 && t_ediff < cutEdiff_) { - unsigned int id = packID(zside,eta,phi,k+1,vbin,pbin); + if (eta == 15) ++n15; + else if (eta == 16) ++n16; + int depth = k + 1; + unsigned int id = packID(zside,eta,phi,depth,vbin,pbin); std::map::iterator it1 = h_Energy_.find(id); if (it1 != h_Energy_.end()) (it1->second)->Fill(ene[k]); std::map::iterator it2 = h_Ecorr_.find(id); @@ -467,12 +496,23 @@ void AnalyzeLepTree::Loop(Long64_t nmax, bool debug) { if (it5 != h_Charge_.end()) (it5->second)->Fill(charge[k]); std::map::iterator it6 = h_Chcorr_.find(id); if (it6 != h_Chcorr_.end()) (it6->second)->Fill(charge[k]/actl[k]); + if (((eta == 15) || (eta == 16)) && debug_) + std::cout << zside << ":" << eta << ":" << phi << ":" + << t_iphi << ":" << depth << ":" << vbin << ":" + << pbin << " ID " << std::hex << id << std::dec + << " Flags " << (it1 != h_Energy_.end()) << ":" + << (it2 != h_Ecorr_.end()) << ":" + << (it3 != h_Charge_.end()) << ":" + << (it4 != h_Chcorr_.end()) << " E " << ene[k] + << " C " << charge[k] << " L " << actl[k] << std::endl; } } } } } } + std::cout << "Number of events with eta15: " << n15 << " and eta16: " + << n16 << std::endl; } bool AnalyzeLepTree::fillChain(TChain *chain, const char* inputFileList) { @@ -535,7 +575,8 @@ void AnalyzeLepTree::bookHisto() { int eta = (ieta>0) ? ieta : -ieta; if (eta != 0) { int ndepth = ((kdepth_ == 0) ? nDepthBins(eta, 63, modeLHC_) : - (kdepth_ == 1) ? 1 : nDepthBins(eta, 63, 0)); + ((kdepth_ != 1) ? nDepthBins(eta, 63, 0) : + (eta == 16) ? 2 : 1)); std::cout << "Eta " << ieta << " with " << nPhiBins(eta) << " phi bins " << ndepth << " maximum depths and " << nPBins(eta) << " p bins" << std::endl; @@ -638,11 +679,16 @@ void AnalyzeLepTree::bookHisto() { sprintf (phis, "All i#phi"); } int ndepth = ((kdepth_ == 0) ? nDepthBins(eta, phi0, modeLHC_) : - (kdepth_ == 1) ? 1 : nDepthBins(eta, phi0, 0)); + ((kdepth_ != 1) ? nDepthBins(eta, phi0, 0) : + (eta == 16) ? 2 : 1)); for (int depth=0; depth 66 || eta < 0))) { } else { - if (dep == 3) depth = 2; - else if (dep == 4) depth = 3; + if (dep > 2) depth = 3; } } else if (ieta < 26) { if (modeLHC_ == 0 || (modeLHC_ == 3 && (phi < 63 || phi > 66 || eta < 0))) { diff --git a/Calibration/HcalCalibAlgos/macros/HBHEMuonHighEta.C b/Calibration/HcalCalibAlgos/macros/HBHEMuonHighEta.C index 911c2ab251d8d..37583bf95f6b3 100644 --- a/Calibration/HcalCalibAlgos/macros/HBHEMuonHighEta.C +++ b/Calibration/HcalCalibAlgos/macros/HBHEMuonHighEta.C @@ -1,9 +1,19 @@ -////////////////////////////////////////////////////////// -// This class has been automatically generated on -// Fri Oct 25 18:21:22 2019 by ROOT version 6.14/09 -// from TTree HBHEMuonHighEta/HBHEMuonHighEta -// found on file: muonHighEta.root -////////////////////////////////////////////////////////// +/////////////////////////////////////////////////////////////////////////////// +// +// HBHEMuonHighEta h1(infile, outfile,t mode, debug) +// h1.Loop() +// +// infile const char* Name of the input file +// outfile const char* Name of the output file +// (dyll_PU20_25_output_10.root) +// mode int Geometry file used 0:(defined by maxDHB/HE); +// 1 (Run 1; valid till 2016); 2 (Run 2; 2018); +// 3 (Run 3; post LS2); 4 (2017 Plan 1); +// 5 (Run 4; post LS3); (default: 0) +// debug bool Debug flag (default: false) +// +/////////////////////////////////////////////////////////////////////////////// + #include #include #include @@ -25,7 +35,7 @@ class HBHEMuonHighEta { public : HBHEMuonHighEta(const char *infile, const char *outfile, - const int mode=1, const bool debug=false); + const int mode=0, const bool debug=false); virtual ~HBHEMuonHighEta(); virtual Int_t Cut(Long64_t entry); virtual Int_t GetEntry(Long64_t entry); @@ -38,9 +48,12 @@ public : private: void BookHistograms(const char* fname); void Close(); + int nDepthBins(int eta, int phi); + TTree *fChain; //!pointer to the analyzed TTree or TChain + Int_t fCurrent; //!current Tree number in a TChain - TTree *fChain; //!pointer to the analyzed TTree or TChain - Int_t fCurrent; //!current Tree number in a TChain + static const int maxDepthHB_ = 7; + static const int maxDepthHE_ = 4; // Fixed size dimensions of array or collections stored in the TTree if any. // Declaration of leaf types @@ -54,7 +67,7 @@ private: std::vector *IsolationR03; std::vector *ecal_3into3; std::vector *hcal_3into3; - std::vector *ho_3into3; + std::vector *tracker_3into3; std::vector *emaxNearP; UInt_t Run_No; UInt_t Event_No; @@ -157,7 +170,7 @@ private: TBranch *b_IsolationR03; TBranch *b_ecal_3into3; TBranch *b_hcal_3into3; - TBranch *b_ho_3into3; + TBranch *b_tracker_3into3; TBranch *b_emaxNearP; TBranch *b_Run_No; TBranch *b_Event_No; @@ -312,7 +325,7 @@ void HBHEMuonHighEta::Init(TTree *tree) { IsolationR03 = 0; ecal_3into3 = 0; hcal_3into3 = 0; - ho_3into3 = 0; + tracker_3into3 = 0; emaxNearP = 0; matchedId = 0; hcal_cellHot = 0; @@ -417,7 +430,7 @@ void HBHEMuonHighEta::Init(TTree *tree) { fChain->SetBranchAddress("IsolationR03", &IsolationR03, &b_IsolationR03); fChain->SetBranchAddress("ecal_3into3", &ecal_3into3, &b_ecal_3into3); fChain->SetBranchAddress("hcal_3into3", &hcal_3into3, &b_hcal_3into3); - fChain->SetBranchAddress("ho_3into3", &ho_3into3, &b_ho_3into3); + fChain->SetBranchAddress("tracker_3into3", &tracker_3into3, &b_tracker_3into3); fChain->SetBranchAddress("emaxNearP", &emaxNearP, &b_emaxNearP); fChain->SetBranchAddress("Run_No", &Run_No, &b_Run_No); fChain->SetBranchAddress("Event_No", &Event_No, &b_Event_No); @@ -585,3 +598,60 @@ void HBHEMuonHighEta::Close() { output_file->Close(); if (debug_) std::cout << "now doing return" << std::endl; } + + +int HBHEMuonHighEta::nDepthBins(int eta, int phi) { + // Run 1 scenario + int nDepthR1[29]={1,1,1,1,1,1,1,1,1,1,1,1,1,1,2,3,1,2,2,2,2,2,2,2,2,2,3,3,2}; + // Run 2 scenario from 2018 + int nDepthR2[29]={1,1,1,1,1,1,1,1,1,1,1,1,1,1,2,4,3,5,6,6,6,6,6,6,6,7,7,7,3}; + // Run 3 scenario + int nDepthR3[29]={4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,3,5,6,6,6,6,6,6,6,7,7,7,3}; + // Run 4 scenario + int nDepthR4[29]={4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,7,7,7,7,7,7,7,7,7,7,7,7,7}; + // for 2021 scenario multi depth segmentation + // int nDepth[29]={3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,5,5,5,5,5,5,5,5,5,5,5,5,5}; + // modeLHC_ = 0 --> nbin defined maxDepthHB/HE + // = 1 --> corresponds to Run 1 (valid till 2016) + // = 2 --> corresponds to Run 2 (2018 geometry) + // = 3 --> corresponds to Run 3 (post LS2) + // = 4 --> corresponds to 2017 (Plan 1) + // = 5 --> corresponds to Run 4 (post LS3) + int nbin(0); + if (modeLHC_ == 0) { + if (eta<=15) { + nbin = maxDepthHB_; + } else if (eta == 16) { + nbin = 4; + } else { + nbin = maxDepthHE_; + } + } else if (modeLHC_ == 1) { + nbin = nDepthR1[eta-1]; + } else if (modeLHC_ == 2) { + nbin = nDepthR2[eta-1]; + } else if (modeLHC_ == 3) { + nbin = nDepthR3[eta-1]; + } else if (modeLHC_ == 4) { + if (phi > 0) { + if (eta >= 16 && phi >= 63 && phi <= 66) { + nbin = nDepthR2[eta-1]; + } else { + nbin = nDepthR1[eta-1]; + } + } else { + if (eta >= 16) { + nbin = (nDepthR2[eta-1] > nDepthR1[eta-1]) ? nDepthR2[eta-1] : nDepthR1[eta-1]; + } else { + nbin = nDepthR1[eta-1]; + } + } + } else { + if (eta > 0 && eta < 30) { + nbin = nDepthR4[eta-1]; + } else { + nbin = nDepthR4[28]; + } + } + return nbin; +} diff --git a/Calibration/HcalCalibAlgos/macros/HBHEMuonOfflineAnalyzer.C b/Calibration/HcalCalibAlgos/macros/HBHEMuonOfflineAnalyzer.C index 25d28b50099c1..7cf5515401f40 100644 --- a/Calibration/HcalCalibAlgos/macros/HBHEMuonOfflineAnalyzer.C +++ b/Calibration/HcalCalibAlgos/macros/HBHEMuonOfflineAnalyzer.C @@ -384,19 +384,19 @@ public : int& zside, int& ieta, int& iphi, int& depth); private: - static const int maxDep=7; - static const int maxEta=29; - static const int maxPhi=72; + static const int maxDep_=7; + static const int maxEta_=29; + static const int maxPhi_=72; //3x16x72x2 + 5x4x72x2 + 5x9x36x2 - static const int maxHist=20000;//13032; + static const int maxHist_=20000;//13032; static const int nCut_=1; static const unsigned int nmax_=10; const bool debug_; int modeLHC_, maxDepthHB_, maxDepthHE_, maxDepth_; int runLo_, runHi_, etaMin_, etaMax_; bool cFactor_, useCorrect_, mergeDepth_; - int nHist, nDepths[maxEta], nDepthsPhi[maxEta]; - int indxEta[maxEta][maxDep][maxPhi]; + int nHist, nDepths[maxEta_], nDepthsPhi[maxEta_]; + int indxEta[maxEta_][maxDep_][maxPhi_]; TFile *output_file; std::map corrFac_[nmax_]; std::vector runlow_; @@ -415,19 +415,19 @@ private: TH1D *h_LongImpactParameter[3], *h_LongImpactParameterBin1[3], *h_LongImpactParameterBin2[3]; TH1D *h_TransImpactParameter[3], *h_TransImpactParameterBin1[3], *h_TransImpactParameterBin2[3]; - TH1D *h_Hot_MuonEnergy_hcal_ClosestCell[3][maxHist] , *h_Hot_MuonEnergy_hcal_HotCell[3][maxHist] , *h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[3][maxHist], *h_HotCell_MuonEnergy_phi[3][maxHist], *h_active_length_Fill[3][maxHist], *h_p_muon_ineta[3][maxHist], *h_charge_signal[3][maxHist], *h_charge_bg[3][maxHist]; + TH1D *h_Hot_MuonEnergy_hcal_ClosestCell[3][maxHist_] , *h_Hot_MuonEnergy_hcal_HotCell[3][maxHist_] , *h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[3][maxHist_], *h_HotCell_MuonEnergy_phi[3][maxHist_], *h_active_length_Fill[3][maxHist_], *h_p_muon_ineta[3][maxHist_], *h_charge_signal[3][maxHist_], *h_charge_bg[3][maxHist_]; TH2D *h_2D_Bin1[3], *h_2D_Bin2[3]; TH1D *h_ecal_energy[3], *h_hcal_energy[3], *h_3x3_ecal[3], *h_1x1_hcal[3]; - TH1D *h_MuonHittingEcal[3], *h_HotCell[3], *h_MuonEnergy_hcal[3][maxHist]; - TH1D *h_Hot_MuonEnergy_hcal[3][maxHist]; + TH1D *h_MuonHittingEcal[3], *h_HotCell[3], *h_MuonEnergy_hcal[3][maxHist_]; + TH1D *h_Hot_MuonEnergy_hcal[3][maxHist_]; TH2D *hcal_ietaVsEnergy[3]; TProfile *h_EtaX_hcal[3], *h_PhiY_hcal[3], *h_EtaX_ecal[3], *h_PhiY_ecal[3]; TProfile *h_Eta_ecal[3], *h_Phi_ecal[3]; - TProfile *h_MuonEnergy_eta[3][maxDep], *h_MuonEnergy_phi[3][maxDep], *h_MuonEnergy_muon_eta[3][maxDep]; - TProfile *h_Hot_MuonEnergy_eta[3][maxDep], *h_Hot_MuonEnergy_phi[3][maxDep], *h_Hot_MuonEnergy_muon_eta[3][maxDep]; - TProfile *h_IsoHot_MuonEnergy_eta[3][maxDep], *h_IsoHot_MuonEnergy_phi[3][maxDep], *h_IsoHot_MuonEnergy_muon_eta[3][maxDep]; - TProfile *h_IsoWithoutHot_MuonEnergy_eta[3][maxDep], *h_IsoWithoutHot_MuonEnergy_phi[3][maxDep], *h_IsoWithoutHot_MuonEnergy_muon_eta[3][maxDep]; - TProfile *h_HotWithoutIso_MuonEnergy_eta[3][maxDep], *h_HotWithoutIso_MuonEnergy_phi[3][maxDep], *h_HotWithoutIso_MuonEnergy_muon_eta[3][maxDep]; + TProfile *h_MuonEnergy_eta[3][maxDep_], *h_MuonEnergy_phi[3][maxDep_], *h_MuonEnergy_muon_eta[3][maxDep_]; + TProfile *h_Hot_MuonEnergy_eta[3][maxDep_], *h_Hot_MuonEnergy_phi[3][maxDep_], *h_Hot_MuonEnergy_muon_eta[3][maxDep_]; + TProfile *h_IsoHot_MuonEnergy_eta[3][maxDep_], *h_IsoHot_MuonEnergy_phi[3][maxDep_], *h_IsoHot_MuonEnergy_muon_eta[3][maxDep_]; + TProfile *h_IsoWithoutHot_MuonEnergy_eta[3][maxDep_], *h_IsoWithoutHot_MuonEnergy_phi[3][maxDep_], *h_IsoWithoutHot_MuonEnergy_muon_eta[3][maxDep_]; + TProfile *h_HotWithoutIso_MuonEnergy_eta[3][maxDep_], *h_HotWithoutIso_MuonEnergy_phi[3][maxDep_], *h_HotWithoutIso_MuonEnergy_muon_eta[3][maxDep_]; }; @@ -1176,8 +1176,8 @@ void HBHEMuonOfflineAnalyzer::bookHistograms(const char* fname) { } } } - if (nHist >= maxHist) { - std::cout << "Problem here " << nHist << ":" << maxHist << std::endl; + if (nHist >= maxHist_) { + std::cout << "Problem here " << nHist << ":" << maxHist_ << std::endl; } // TDirectory *d_output_file[nCut_][29]; //output_file->cd(); @@ -1714,7 +1714,7 @@ int HBHEMuonOfflineAnalyzer::nDepthBins(int eta, int phi) { int nDepthR3[29]={4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,3,5,6,6,6,6,6,6,6,7,7,7,3}; // Run 4 scenario int nDepthR4[29]={4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,7,7,7,7,7,7,7,7,7,7,7,7,7}; - // for 2021 scenario multi depth segmentation + // for a test scenario with multi depth segmentation considered during Run 1 // int nDepth[29]={3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,5,5,5,5,5,5,5,5,5,5,5,5,5}; // modeLHC_ = 0 --> nbin defined maxDepthHB/HE // = 1 --> corresponds to Run 1 (valid till 2016) From c107188783f98e58c97deb0eac0deba5498e8505 Mon Sep 17 00:00:00 2001 From: Sunanda Date: Fri, 24 Jul 2020 19:58:04 +0200 Subject: [PATCH 2/3] Code check --- .../HcalCalibAlgos/macros/AnalyzeLepTree.C | 1455 ++++----- .../HcalCalibAlgos/macros/HBHEMuonHighEta.C | 535 ++-- .../macros/HBHEMuonOfflineAnalyzer.C | 2640 +++++++++-------- 3 files changed, 2396 insertions(+), 2234 deletions(-) diff --git a/Calibration/HcalCalibAlgos/macros/AnalyzeLepTree.C b/Calibration/HcalCalibAlgos/macros/AnalyzeLepTree.C index 327a1aca015ae..c217a18f7aafd 100644 --- a/Calibration/HcalCalibAlgos/macros/AnalyzeLepTree.C +++ b/Calibration/HcalCalibAlgos/macros/AnalyzeLepTree.C @@ -5,13 +5,13 @@ // There are 2 additional arguments: // mode (a packed integer) with bits specifying some action // Bit 0 : 0 will produce set of histograms of energy deposit, -// active length, active length corrected energy; +// active length, active length corrected energy; // 1 will produce plots of p, nvx and scatter plots for // each ieta // 1 : 0 ignores the information on number of vertex // 1 groups ranges of # of vertex 0:15:20:25:30:100 // 2 : 0 ignores the information on track momentum -// 1 separate plots for certain momentum range +// 1 separate plots for certain momentum range // (the range depends on ieta) // 3-4: 0 separate plot for each depth // 1 sums up all depths @@ -39,12 +39,12 @@ // a1.writeMeanError(outfileMean); // // outfile (char*) Name of the file where histograms to be written -// outfileMean (char*) Name of the file where means with errors to be +// outfileMean (char*) Name of the file where means with errors to be // written when it is run with mode bit 0 set to 1 // drawStatBox (bool) If Statbox to be drawn or not // type (int) Each bit says what plots to be drawn // If bit 0 of "mode" is 1 -// (0: momentum for each ieta; +// (0: momentum for each ieta; // 1: number of vertex for each ieta; // 2: 2D plot for nvx vs p for each ieta; // 3: number of vertex for each ieta & depth; @@ -56,7 +56,7 @@ // value and "phi" as of (type/256)&127 is 0 or // the specified value. Bit 0 set plots the energy // distribution; 1 set plots active length corrected -// energy; 2 set plots charge; 3 set plots active +// energy; 2 set plots charge; 3 set plots active // length corrected charge distributions // save (bool) If plots to be saved as pdf file or not /////////////////////////////////////////////////////////////////////////////// @@ -81,83 +81,77 @@ #include class AnalyzeLepTree { -public : - AnalyzeLepTree(TChain *tree, int mode=0, int modeLHC=3); - AnalyzeLepTree(const char* fname, int mode=0, int modeLHC=3); +public: + AnalyzeLepTree(TChain* tree, int mode = 0, int modeLHC = 3); + AnalyzeLepTree(const char* fname, int mode = 0, int modeLHC = 3); virtual ~AnalyzeLepTree(); - virtual Int_t Cut(Long64_t entry); - virtual Int_t GetEntry(Long64_t entry); - virtual Long64_t LoadTree(Long64_t entry); - virtual void Init(TChain *tree); - virtual void Loop(Long64_t nmax=-1, bool debug=false); - virtual Bool_t Notify(); - virtual void Show(Long64_t entry = -1); - void writeHisto(const char* outfile); - void writeMeanError(const char* outfile); - std::vector plotHisto(bool drawStatBox, int type, bool save=false); + virtual Int_t Cut(Long64_t entry); + virtual Int_t GetEntry(Long64_t entry); + virtual Long64_t LoadTree(Long64_t entry); + virtual void Init(TChain* tree); + virtual void Loop(Long64_t nmax = -1, bool debug = false); + virtual Bool_t Notify(); + virtual void Show(Long64_t entry = -1); + void writeHisto(const char* outfile); + void writeMeanError(const char* outfile); + std::vector plotHisto(bool drawStatBox, int type, bool save = false); + private: - bool fillChain(TChain* chain, const char* fname); - void bookHisto(); - void plotHisto(std::map hists, - std::vector& cvs, bool save); - void plotHisto(std::map hists, int phi0, - int depth0, std::vector& cvs, - bool save); - TCanvas* plotHisto(TH1D* hist); - void plot2DHisto(std::map hists, - std::vector& cvs, bool save); - int getCollapsedDepth(int eta, int phi, int depth); - int getRBX(int eta); - int getPBin(int eta); - int getVxBin(); - int getDepthBin(int depth); - int getPhiBin(int eta); - void makeVxBins(int modeLHC); - int nDepthBins(int eta, int phi, int modeLHC); - int nPhiBins(int eta); - int nPBins(int eta); - int nVxBins(); - unsigned int packID(int zside, int eta, int phi, int depth, - int nvx, int ipbin); - void unpackID(unsigned int id, int& zside, int& eta, - int& phi, int& depth, int& nvx, int& ipbin); - void getBins(int type, int eta, int phi, int depth, - int& nbins, double& xmax); + bool fillChain(TChain* chain, const char* fname); + void bookHisto(); + void plotHisto(std::map hists, std::vector& cvs, bool save); + void plotHisto(std::map hists, int phi0, int depth0, std::vector& cvs, bool save); + TCanvas* plotHisto(TH1D* hist); + void plot2DHisto(std::map hists, std::vector& cvs, bool save); + int getCollapsedDepth(int eta, int phi, int depth); + int getRBX(int eta); + int getPBin(int eta); + int getVxBin(); + int getDepthBin(int depth); + int getPhiBin(int eta); + void makeVxBins(int modeLHC); + int nDepthBins(int eta, int phi, int modeLHC); + int nPhiBins(int eta); + int nPBins(int eta); + int nVxBins(); + unsigned int packID(int zside, int eta, int phi, int depth, int nvx, int ipbin); + void unpackID(unsigned int id, int& zside, int& eta, int& phi, int& depth, int& nvx, int& ipbin); + void getBins(int type, int eta, int phi, int depth, int& nbins, double& xmax); private: - TChain *fChain; //!pointer to the analyzed TTree or TChain - Int_t fCurrent; //!current Tree number in a TChain - + TChain* fChain; //!pointer to the analyzed TTree or TChain + Int_t fCurrent; //!current Tree number in a TChain + // Declaration of leaf types - Int_t t_ieta; - Int_t t_iphi; - Int_t t_nvtx; - Double_t t_p; - Double_t t_ediff; - std::vector *t_ene; - std::vector *t_enec; - std::vector *t_charge; - std::vector *t_actln; - std::vector *t_depth; - + Int_t t_ieta; + Int_t t_iphi; + Int_t t_nvtx; + Double_t t_p; + Double_t t_ediff; + std::vector* t_ene; + std::vector* t_enec; + std::vector* t_charge; + std::vector* t_actln; + std::vector* t_depth; + // List of branches - TBranch *b_t_ieta; //! - TBranch *b_t_iphi; //! - TBranch *b_t_nvtx; //! - TBranch *b_t_p; //! - TBranch *b_t_ediff; //! - TBranch *b_t_ene; //! - TBranch *b_t_enec; //! - TBranch *b_t_charge; //! - TBranch *b_t_actln; //! - TBranch *b_t_depth; //! + TBranch* b_t_ieta; //! + TBranch* b_t_iphi; //! + TBranch* b_t_nvtx; //! + TBranch* b_t_p; //! + TBranch* b_t_ediff; //! + TBranch* b_t_ene; //! + TBranch* b_t_enec; //! + TBranch* b_t_charge; //! + TBranch* b_t_actln; //! + TBranch* b_t_depth; //! - static const int etamax_=26, npbin_=9, nvbin_=6; - static const bool debug_=false; - int mode_, modeLHC_, exRBX_, kphi_, kdepth_; - std::vector npvbin_, iprange_; - std::vector prange_[5]; - double cutEdiff_; + static const int etamax_ = 26, npbin_ = 9, nvbin_ = 6; + static const bool debug_ = false; + int mode_, modeLHC_, exRBX_, kphi_, kdepth_; + std::vector npvbin_, iprange_; + std::vector prange_[5]; + double cutEdiff_; std::map h_p_, h_nv_; std::map h_pnv_; std::map h_p2_, h_nv2_; @@ -166,43 +160,44 @@ private: std::map h_ediff_, h_ediff_nvtx_; }; -AnalyzeLepTree::AnalyzeLepTree(TChain *tree, int mode1, - int mode2) : mode_(mode1), modeLHC_(mode2) { - std::cout << "Proceed with a tree chain with " << tree->GetEntries() - << " entries" << std::endl; +AnalyzeLepTree::AnalyzeLepTree(TChain* tree, int mode1, int mode2) : mode_(mode1), modeLHC_(mode2) { + std::cout << "Proceed with a tree chain with " << tree->GetEntries() << " entries" << std::endl; Init(tree); } -AnalyzeLepTree::AnalyzeLepTree(const char* fname, int mode1, - int mode2) : mode_(mode1), modeLHC_(mode2) { - TChain *chain = new TChain("Lep_Tree"); - if (!fillChain(chain,fname)) { +AnalyzeLepTree::AnalyzeLepTree(const char* fname, int mode1, int mode2) : mode_(mode1), modeLHC_(mode2) { + TChain* chain = new TChain("Lep_Tree"); + if (!fillChain(chain, fname)) { std::cout << "*****No valid tree chain can be obtained*****" << std::endl; } else { - std::cout << "Proceed with a tree chain with " << chain->GetEntries() - << " entries" << std::endl; + std::cout << "Proceed with a tree chain with " << chain->GetEntries() << " entries" << std::endl; Init(chain); } } AnalyzeLepTree::~AnalyzeLepTree() { - if (!fChain) return; + if (!fChain) + return; delete fChain->GetCurrentFile(); } Int_t AnalyzeLepTree::GetEntry(Long64_t entry) { // Read contents of entry. - if (!fChain) return 0; + if (!fChain) + return 0; return fChain->GetEntry(entry); } Long64_t AnalyzeLepTree::LoadTree(Long64_t entry) { // Set the environment to read one entry - if (!fChain) return -5; + if (!fChain) + return -5; Long64_t centry = fChain->LoadTree(entry); - if (centry < 0) return centry; - if (!fChain->InheritsFrom(TChain::Class())) return centry; - TChain *chain = (TChain*)fChain; + if (centry < 0) + return centry; + if (!fChain->InheritsFrom(TChain::Class())) + return centry; + TChain* chain = (TChain*)fChain; if (chain->GetTreeNumber() != fCurrent) { fCurrent = chain->GetTreeNumber(); Notify(); @@ -210,7 +205,7 @@ Long64_t AnalyzeLepTree::LoadTree(Long64_t entry) { return centry; } -void AnalyzeLepTree::Init(TChain *tree) { +void AnalyzeLepTree::Init(TChain* tree) { // The Init() function is called when the selector needs to initialize // a new tree or chain. Typically here the branch addresses and branch // pointers of the tree will be set. @@ -220,71 +215,85 @@ void AnalyzeLepTree::Init(TChain *tree) { // (once per file to be processed). makeVxBins(modeLHC_); - exRBX_ = (mode_/128)%32; - kphi_ = (mode_/32)%4; - kdepth_ = (mode_/8)%4; - if ((mode_%2) == 0) std::cout << "Produce set of histograms of energy, " - << " active length, active length corrected " - << "energy for 3 types" << std::endl; - else std::cout << "Produce plots of p, nvx and scatter plots " - << "for each ieta" << std::endl; - if (((mode_/2)%2) == 0) { - std::cout << "Ignore the information on number of vertex iformation" - << std::endl; + exRBX_ = (mode_ / 128) % 32; + kphi_ = (mode_ / 32) % 4; + kdepth_ = (mode_ / 8) % 4; + if ((mode_ % 2) == 0) + std::cout << "Produce set of histograms of energy, " + << " active length, active length corrected " + << "energy for 3 types" << std::endl; + else + std::cout << "Produce plots of p, nvx and scatter plots " + << "for each ieta" << std::endl; + if (((mode_ / 2) % 2) == 0) { + std::cout << "Ignore the information on number of vertex iformation" << std::endl; } else { std::cout << "Group ranges of # of vertex "; - for (unsigned int k=0; kSetMakeClass(1); - fChain->SetBranchAddress("t_ieta", &t_ieta, &b_t_ieta); - fChain->SetBranchAddress("t_iphi", &t_iphi, &b_t_iphi); - fChain->SetBranchAddress("t_nvtx", &t_nvtx, &b_t_nvtx); - fChain->SetBranchAddress("t_p", &t_p, &b_t_p); - fChain->SetBranchAddress("t_ediff", &t_ediff, &b_t_ediff); - fChain->SetBranchAddress("t_ene", &t_ene, &b_t_ene); - fChain->SetBranchAddress("t_enec", &t_enec, &b_t_enec); + fChain->SetBranchAddress("t_ieta", &t_ieta, &b_t_ieta); + fChain->SetBranchAddress("t_iphi", &t_iphi, &b_t_iphi); + fChain->SetBranchAddress("t_nvtx", &t_nvtx, &b_t_nvtx); + fChain->SetBranchAddress("t_p", &t_p, &b_t_p); + fChain->SetBranchAddress("t_ediff", &t_ediff, &b_t_ediff); + fChain->SetBranchAddress("t_ene", &t_ene, &b_t_ene); + fChain->SetBranchAddress("t_enec", &t_enec, &b_t_enec); fChain->SetBranchAddress("t_charge", &t_charge, &b_t_charge); - fChain->SetBranchAddress("t_actln", &t_actln, &b_t_actln); - fChain->SetBranchAddress("t_depth", &t_depth, &b_t_depth); + fChain->SetBranchAddress("t_actln", &t_actln, &b_t_actln); + fChain->SetBranchAddress("t_depth", &t_depth, &b_t_depth); Notify(); - t_ediff = 0; + t_ediff = 0; bookHisto(); } @@ -301,11 +310,12 @@ Bool_t AnalyzeLepTree::Notify() { void AnalyzeLepTree::Show(Long64_t entry) { // Print contents of entry. // If entry is not specified, print current entry - if (!fChain) return; + if (!fChain) + return; fChain->Show(entry); } -Int_t AnalyzeLepTree::Cut(Long64_t ) { +Int_t AnalyzeLepTree::Cut(Long64_t) { // This function may be called from Loop. // returns 1 if entry is accepted. // returns -1 otherwise. @@ -336,224 +346,244 @@ void AnalyzeLepTree::Loop(Long64_t nmax, bool debug) { // METHOD2: replace line // fChain->GetEntry(jentry); //read all branches //by b_branchname->GetEntry(ientry); //read only this branch - if (fChain == 0) return; + if (fChain == 0) + return; Long64_t nentries = fChain->GetEntriesFast(); std::cout << "Number of entries: " << nentries << ":" << nmax << std::endl; - if (nmax > 0 && nmax < nentries) nentries = nmax; - const double ethr = 0.00001; // Threshold of 10 keV + if (nmax > 0 && nmax < nentries) + nentries = nmax; + const double ethr = 0.00001; // Threshold of 10 keV Long64_t nbytes = 0, nb = 0; int32_t n15(0), n16(0); - for (Long64_t jentry=0; jentryGetEntry(jentry); nbytes += nb; + if (ientry < 0) + break; + nb = fChain->GetEntry(jentry); + nbytes += nb; // if (Cut(ientry) < 0) continue; int zside = (t_ieta > 0) ? 1 : -1; - int eta = (t_ieta > 0) ? t_ieta : -t_ieta; - int phi = getPhiBin(eta); - int pbin = getPBin(eta); - int vbin = getVxBin(); - if ((mode_/1)%2 == 1) { - unsigned int id0 = packID(zside,eta,1,1,1,1); - std::map::iterator it1 = h_p_.find(id0); - if (it1 != h_p_.end()) (it1->second)->Fill(t_p); - std::map::iterator it2 = h_nv_.find(id0); - if (it2 != h_nv_.end()) (it2->second)->Fill(t_nvtx); - std::map::iterator it3 = h_pnv_.find(id0); - if (it3 != h_pnv_.end()) (it3->second)->Fill(t_nvtx,t_p); - unsigned int id1 = packID(zside,eta,1,1,1,pbin); - std::map::iterator it4 = h_p2_.find(id1); - if (it4 != h_p2_.end()) (it4->second)->Fill(t_p); - unsigned int id2 = packID(zside,eta,1,1,vbin,1); - std::map::iterator it5 = h_nv2_.find(id2); - if (it5 != h_nv2_.end()) (it5->second)->Fill(t_nvtx); - std::map::iterator it6 = h_ediff_.find(id0); - if (it6 != h_ediff_.end()) (it6->second)->Fill(t_ediff); - std::map::iterator it7 = h_ediff_nvtx_.find(id2); - if (it7 != h_ediff_nvtx_.end()) (it7->second)->Fill(t_ediff); + int eta = (t_ieta > 0) ? t_ieta : -t_ieta; + int phi = getPhiBin(eta); + int pbin = getPBin(eta); + int vbin = getVxBin(); + if ((mode_ / 1) % 2 == 1) { + unsigned int id0 = packID(zside, eta, 1, 1, 1, 1); + std::map::iterator it1 = h_p_.find(id0); + if (it1 != h_p_.end()) + (it1->second)->Fill(t_p); + std::map::iterator it2 = h_nv_.find(id0); + if (it2 != h_nv_.end()) + (it2->second)->Fill(t_nvtx); + std::map::iterator it3 = h_pnv_.find(id0); + if (it3 != h_pnv_.end()) + (it3->second)->Fill(t_nvtx, t_p); + unsigned int id1 = packID(zside, eta, 1, 1, 1, pbin); + std::map::iterator it4 = h_p2_.find(id1); + if (it4 != h_p2_.end()) + (it4->second)->Fill(t_p); + unsigned int id2 = packID(zside, eta, 1, 1, vbin, 1); + std::map::iterator it5 = h_nv2_.find(id2); + if (it5 != h_nv2_.end()) + (it5->second)->Fill(t_nvtx); + std::map::iterator it6 = h_ediff_.find(id0); + if (it6 != h_ediff_.end()) + (it6->second)->Fill(t_ediff); + std::map::iterator it7 = h_ediff_nvtx_.find(id2); + if (it7 != h_ediff_nvtx_.end()) + (it7->second)->Fill(t_ediff); } else { if (phi > 0 && pbin >= 0 && vbin >= 0) { - if (kdepth_ == 0) { - for (unsigned int k=0; ksize(); ++k) { - if (eta == 15) ++n15; - else if (eta == 16) ++n16; - int depth = (*t_depth)[k]; - unsigned int id = packID(zside,eta,phi,depth+1,vbin,pbin); - double ene = (*t_ene)[k]; - double enec = (*t_enec)[k]; - double charge = (*t_charge)[k]; - double actl = (*t_actln)[k]; - if (ene > ethr && actl > 0 && charge > 0 && t_ediff < cutEdiff_) { - std::map::iterator it1 = h_Energy_.find(id); - if (it1 != h_Energy_.end()) (it1->second)->Fill(ene); - std::map::iterator it2 = h_Ecorr_.find(id); - if (it2 != h_Ecorr_.end()) (it2->second)->Fill(ene/actl); - std::map::iterator it3 = h_EnergyC_.find(id); - if (it3 != h_EnergyC_.end()) (it3->second)->Fill(enec); - std::map::iterator it4 = h_EcorrC_.find(id); - if (it4 != h_EcorrC_.end()) (it4->second)->Fill(enec/actl); - std::map::iterator it5 = h_Charge_.find(id); - if (it5 != h_Charge_.end()) (it5->second)->Fill(charge); - std::map::iterator it6 = h_Chcorr_.find(id); - if (it6 != h_Chcorr_.end()) (it6->second)->Fill(charge/actl); - if (debug_) { -// if ((eta>20 && (t_iphi > 35)) || (t_iphi > 71)) std::cout << zside << ":" << eta << ":" << phi << ":" << t_iphi << ":" << depth+1 << ":" << vbin << ":" << pbin << " ID " << std::hex << id << std::dec << " Flags " << (it1 != h_Energy_.end()) << ":" << (it2 != h_Ecorr_.end()) << ":" << (it3 != h_EnergyC_.end()) << ":" << (it4 != h_EcorrC_.end()) << ":" << (it5 != h_Charge_.end()) << ":" << (it6 != h_Chcorr_.end()) << " E " << ene << " C " << charge << " L " << actl << std::endl; - if ((it1 == h_Energy_.end()) || (it2 == h_Ecorr_.end()) || - (it3 == h_EnergyC_.end()) || (it4 == h_EcorrC_.end()) || - (it5 == h_Charge_.end()) || (it6 == h_Chcorr_.end())) - std::cout << zside << ":" << eta << ":" << phi << ":" - << t_iphi << ":" << depth+1 << ":" << vbin << ":" - << pbin << " ID " << std::hex << id << std::dec - << " Flags " << (it1 != h_Energy_.end()) << ":" - << (it2 != h_Ecorr_.end()) << ":" - << (it3 != h_Charge_.end()) << ":" - << (it4 != h_Chcorr_.end()) << " E " << ene - << " C " << charge << " L " << actl << std::endl; - } - } - } - } else if (kdepth_ == 1) { - double ene[2], enec[2], actl[2], charge[2]; - for (unsigned int k=0; k<2; ++k) { - ene[k] = enec[k] = actl[k] = charge[k] = 0; - } - for (unsigned int k=0; ksize(); ++k) { - if ((*t_ene)[k] > 0 && (*t_actln)[k] > 0) { - int dep = (*t_depth)[k]; - int depth = ((eta != 16) ? 0 : ((dep > 1) ? 1 : 0)); - ene[depth] += (*t_ene)[k]; - enec[depth] += (*t_enec)[k]; - charge[depth] += (*t_charge)[k]; - actl[depth] += (*t_actln)[k]; - } - } - int nDepth = (eta == 16) ? 2 : 1; - for (int k=0; k ethr && actl[k] > 0 && charge[k] > 0 && - t_ediff < cutEdiff_) { - if (eta == 15) ++n15; - else if (eta == 16) ++n16; - int depth = k + 1; - unsigned int id = packID(zside,eta,phi,depth,vbin,pbin); - std::map::iterator it1 = h_Energy_.find(id); - if (it1 != h_Energy_.end()) (it1->second)->Fill(ene[k]); - std::map::iterator it2 = h_Ecorr_.find(id); - if (it2 != h_Ecorr_.end()) (it2->second)->Fill(ene[k]/actl[k]); - std::map::iterator it3 = h_EnergyC_.find(id); - if (it3 != h_EnergyC_.end()) (it3->second)->Fill(enec[k]); - std::map::iterator it4 = h_EcorrC_.find(id); - if (it4 != h_EcorrC_.end()) (it4->second)->Fill(enec[k]/actl[k]); - std::map::iterator it5 = h_Charge_.find(id); - if (it5 != h_Charge_.end()) (it5->second)->Fill(charge[k]); - std::map::iterator it6 = h_Chcorr_.find(id); - if (it6 != h_Chcorr_.end()) (it6->second)->Fill(charge[k]/actl[k]); - if (((eta == 15) || (eta == 16)) && debug_) - std::cout << zside << ":" << eta << ":" << phi << ":" - << t_iphi << ":" << depth << ":" << vbin << ":" - << pbin << " ID " << std::hex << id << std::dec - << " Flags " << (it1 != h_Energy_.end()) << ":" - << (it2 != h_Ecorr_.end()) << ":" - << (it3 != h_Charge_.end()) << ":" - << (it4 != h_Chcorr_.end()) << " E " << ene - << " C " << charge << " L " << actl << std::endl; - } - } - } else { - double ene[3], enec[3], actl[3], charge[3]; - for (unsigned int k=0; k<3; ++k) { - ene[k] = enec[k] = actl[k] = charge[k] = 0; - } - for (unsigned int k=0; ksize(); ++k) { - if ((*t_ene)[k] > 0 && (*t_actln)[k] > 0) { - int dep = (*t_depth)[k]; - int depth = getCollapsedDepth(zside*eta, phi, dep) - 1; - ene[depth] += (*t_ene)[k]; - enec[depth] += (*t_enec)[k]; - charge[depth] += (*t_charge)[k]; - actl[depth] += (*t_actln)[k]; - } - } - for (int k=0; k ethr && actl[k] > 0 && charge[k] > 0 && - t_ediff < cutEdiff_) { - if (eta == 15) ++n15; - else if (eta == 16) ++n16; - int depth = k + 1; - unsigned int id = packID(zside,eta,phi,depth,vbin,pbin); - std::map::iterator it1 = h_Energy_.find(id); - if (it1 != h_Energy_.end()) (it1->second)->Fill(ene[k]); - std::map::iterator it2 = h_Ecorr_.find(id); - if (it2 != h_Ecorr_.end()) (it2->second)->Fill(ene[k]/actl[k]); - std::map::iterator it3 = h_EnergyC_.find(id); - if (it3 != h_EnergyC_.end()) (it3->second)->Fill(enec[k]); - std::map::iterator it4 = h_EcorrC_.find(id); - if (it4 != h_EcorrC_.end()) (it4->second)->Fill(enec[k]/actl[k]); - std::map::iterator it5 = h_Charge_.find(id); - if (it5 != h_Charge_.end()) (it5->second)->Fill(charge[k]); - std::map::iterator it6 = h_Chcorr_.find(id); - if (it6 != h_Chcorr_.end()) (it6->second)->Fill(charge[k]/actl[k]); - if (((eta == 15) || (eta == 16)) && debug_) - std::cout << zside << ":" << eta << ":" << phi << ":" - << t_iphi << ":" << depth << ":" << vbin << ":" - << pbin << " ID " << std::hex << id << std::dec - << " Flags " << (it1 != h_Energy_.end()) << ":" - << (it2 != h_Ecorr_.end()) << ":" - << (it3 != h_Charge_.end()) << ":" - << (it4 != h_Chcorr_.end()) << " E " << ene[k] - << " C " << charge[k] << " L " << actl[k] << std::endl; - } - } - } + if (kdepth_ == 0) { + for (unsigned int k = 0; k < t_depth->size(); ++k) { + if (eta == 15) + ++n15; + else if (eta == 16) + ++n16; + int depth = (*t_depth)[k]; + unsigned int id = packID(zside, eta, phi, depth + 1, vbin, pbin); + double ene = (*t_ene)[k]; + double enec = (*t_enec)[k]; + double charge = (*t_charge)[k]; + double actl = (*t_actln)[k]; + if (ene > ethr && actl > 0 && charge > 0 && t_ediff < cutEdiff_) { + std::map::iterator it1 = h_Energy_.find(id); + if (it1 != h_Energy_.end()) + (it1->second)->Fill(ene); + std::map::iterator it2 = h_Ecorr_.find(id); + if (it2 != h_Ecorr_.end()) + (it2->second)->Fill(ene / actl); + std::map::iterator it3 = h_EnergyC_.find(id); + if (it3 != h_EnergyC_.end()) + (it3->second)->Fill(enec); + std::map::iterator it4 = h_EcorrC_.find(id); + if (it4 != h_EcorrC_.end()) + (it4->second)->Fill(enec / actl); + std::map::iterator it5 = h_Charge_.find(id); + if (it5 != h_Charge_.end()) + (it5->second)->Fill(charge); + std::map::iterator it6 = h_Chcorr_.find(id); + if (it6 != h_Chcorr_.end()) + (it6->second)->Fill(charge / actl); + if (debug_) { + // if ((eta>20 && (t_iphi > 35)) || (t_iphi > 71)) std::cout << zside << ":" << eta << ":" << phi << ":" << t_iphi << ":" << depth+1 << ":" << vbin << ":" << pbin << " ID " << std::hex << id << std::dec << " Flags " << (it1 != h_Energy_.end()) << ":" << (it2 != h_Ecorr_.end()) << ":" << (it3 != h_EnergyC_.end()) << ":" << (it4 != h_EcorrC_.end()) << ":" << (it5 != h_Charge_.end()) << ":" << (it6 != h_Chcorr_.end()) << " E " << ene << " C " << charge << " L " << actl << std::endl; + if ((it1 == h_Energy_.end()) || (it2 == h_Ecorr_.end()) || (it3 == h_EnergyC_.end()) || + (it4 == h_EcorrC_.end()) || (it5 == h_Charge_.end()) || (it6 == h_Chcorr_.end())) + std::cout << zside << ":" << eta << ":" << phi << ":" << t_iphi << ":" << depth + 1 << ":" << vbin + << ":" << pbin << " ID " << std::hex << id << std::dec << " Flags " + << (it1 != h_Energy_.end()) << ":" << (it2 != h_Ecorr_.end()) << ":" + << (it3 != h_Charge_.end()) << ":" << (it4 != h_Chcorr_.end()) << " E " << ene << " C " + << charge << " L " << actl << std::endl; + } + } + } + } else if (kdepth_ == 1) { + double ene[2], enec[2], actl[2], charge[2]; + for (unsigned int k = 0; k < 2; ++k) { + ene[k] = enec[k] = actl[k] = charge[k] = 0; + } + for (unsigned int k = 0; k < t_depth->size(); ++k) { + if ((*t_ene)[k] > 0 && (*t_actln)[k] > 0) { + int dep = (*t_depth)[k]; + int depth = ((eta != 16) ? 0 : ((dep > 1) ? 1 : 0)); + ene[depth] += (*t_ene)[k]; + enec[depth] += (*t_enec)[k]; + charge[depth] += (*t_charge)[k]; + actl[depth] += (*t_actln)[k]; + } + } + int nDepth = (eta == 16) ? 2 : 1; + for (int k = 0; k < nDepth; ++k) { + if (ene[k] > ethr && actl[k] > 0 && charge[k] > 0 && t_ediff < cutEdiff_) { + if (eta == 15) + ++n15; + else if (eta == 16) + ++n16; + int depth = k + 1; + unsigned int id = packID(zside, eta, phi, depth, vbin, pbin); + std::map::iterator it1 = h_Energy_.find(id); + if (it1 != h_Energy_.end()) + (it1->second)->Fill(ene[k]); + std::map::iterator it2 = h_Ecorr_.find(id); + if (it2 != h_Ecorr_.end()) + (it2->second)->Fill(ene[k] / actl[k]); + std::map::iterator it3 = h_EnergyC_.find(id); + if (it3 != h_EnergyC_.end()) + (it3->second)->Fill(enec[k]); + std::map::iterator it4 = h_EcorrC_.find(id); + if (it4 != h_EcorrC_.end()) + (it4->second)->Fill(enec[k] / actl[k]); + std::map::iterator it5 = h_Charge_.find(id); + if (it5 != h_Charge_.end()) + (it5->second)->Fill(charge[k]); + std::map::iterator it6 = h_Chcorr_.find(id); + if (it6 != h_Chcorr_.end()) + (it6->second)->Fill(charge[k] / actl[k]); + if (((eta == 15) || (eta == 16)) && debug_) + std::cout << zside << ":" << eta << ":" << phi << ":" << t_iphi << ":" << depth << ":" << vbin << ":" + << pbin << " ID " << std::hex << id << std::dec << " Flags " << (it1 != h_Energy_.end()) + << ":" << (it2 != h_Ecorr_.end()) << ":" << (it3 != h_Charge_.end()) << ":" + << (it4 != h_Chcorr_.end()) << " E " << ene << " C " << charge << " L " << actl << std::endl; + } + } + } else { + double ene[3], enec[3], actl[3], charge[3]; + for (unsigned int k = 0; k < 3; ++k) { + ene[k] = enec[k] = actl[k] = charge[k] = 0; + } + for (unsigned int k = 0; k < t_depth->size(); ++k) { + if ((*t_ene)[k] > 0 && (*t_actln)[k] > 0) { + int dep = (*t_depth)[k]; + int depth = getCollapsedDepth(zside * eta, phi, dep) - 1; + ene[depth] += (*t_ene)[k]; + enec[depth] += (*t_enec)[k]; + charge[depth] += (*t_charge)[k]; + actl[depth] += (*t_actln)[k]; + } + } + for (int k = 0; k < nDepthBins(eta, phi, 0); ++k) { + if (ene[k] > ethr && actl[k] > 0 && charge[k] > 0 && t_ediff < cutEdiff_) { + if (eta == 15) + ++n15; + else if (eta == 16) + ++n16; + int depth = k + 1; + unsigned int id = packID(zside, eta, phi, depth, vbin, pbin); + std::map::iterator it1 = h_Energy_.find(id); + if (it1 != h_Energy_.end()) + (it1->second)->Fill(ene[k]); + std::map::iterator it2 = h_Ecorr_.find(id); + if (it2 != h_Ecorr_.end()) + (it2->second)->Fill(ene[k] / actl[k]); + std::map::iterator it3 = h_EnergyC_.find(id); + if (it3 != h_EnergyC_.end()) + (it3->second)->Fill(enec[k]); + std::map::iterator it4 = h_EcorrC_.find(id); + if (it4 != h_EcorrC_.end()) + (it4->second)->Fill(enec[k] / actl[k]); + std::map::iterator it5 = h_Charge_.find(id); + if (it5 != h_Charge_.end()) + (it5->second)->Fill(charge[k]); + std::map::iterator it6 = h_Chcorr_.find(id); + if (it6 != h_Chcorr_.end()) + (it6->second)->Fill(charge[k] / actl[k]); + if (((eta == 15) || (eta == 16)) && debug_) + std::cout << zside << ":" << eta << ":" << phi << ":" << t_iphi << ":" << depth << ":" << vbin << ":" + << pbin << " ID " << std::hex << id << std::dec << " Flags " << (it1 != h_Energy_.end()) + << ":" << (it2 != h_Ecorr_.end()) << ":" << (it3 != h_Charge_.end()) << ":" + << (it4 != h_Chcorr_.end()) << " E " << ene[k] << " C " << charge[k] << " L " << actl[k] + << std::endl; + } + } + } } } } - std::cout << "Number of events with eta15: " << n15 << " and eta16: " - << n16 << std::endl; + std::cout << "Number of events with eta15: " << n15 << " and eta16: " << n16 << std::endl; } -bool AnalyzeLepTree::fillChain(TChain *chain, const char* inputFileList) { - +bool AnalyzeLepTree::fillChain(TChain* chain, const char* inputFileList) { int kount(0); std::string fname(inputFileList); - if (fname.substr(fname.size()-5,5) == ".root") { + if (fname.substr(fname.size() - 5, 5) == ".root") { chain->Add(fname.c_str()); } else { ifstream infile(inputFileList); if (!infile.is_open()) { - std::cout << "** ERROR: Can't open '" << inputFileList << "' for input" - << std::endl; + std::cout << "** ERROR: Can't open '" << inputFileList << "' for input" << std::endl; return false; } while (1) { infile >> fname; - if (!infile.good()) break; + if (!infile.good()) + break; chain->Add(fname.c_str()); ++kount; } infile.close(); } - std::cout << "Adds " << kount << " files in the chain from " << fname - << std::endl; + std::cout << "Adds " << kount << " files in the chain from " << fname << std::endl; return true; } void AnalyzeLepTree::bookHisto() { - - for (int i=0; i<5; ++i) prange_[i].clear(); - int ipbrng[30]= {0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,2,2,2,2,2,2,3,3,3,4,4,4,4}; - for (int i=0; i<30; ++i) iprange_.push_back(ipbrng[i]); - double prange0[npbin_] = {0,30,45,55,75,100,125,150,500}; - double prange1[npbin_] = {0,50,75,100,125,150,200,300,500}; - double prange2[npbin_] = {0,60,75,100,125,150,200,300,500}; - double prange3[npbin_] = {0,100,125,150,175,200,300,400,500}; - double prange4[npbin_] = {0,125,150,175,200,250,300,400,500}; - double prangeX[npbin_] = {125,150,200,250,300,350,400,500}; - for (int i=0; i0) ? ieta : -ieta; - if (eta != 0) { - int ndepth = ((kdepth_ == 0) ? nDepthBins(eta, 63, modeLHC_) : - ((kdepth_ != 1) ? nDepthBins(eta, 63, 0) : - (eta == 16) ? 2 : 1)); - std::cout << "Eta " << ieta << " with " << nPhiBins(eta) - << " phi bins " << ndepth << " maximum depths and " - << nPBins(eta) << " p bins" << std::endl; - } + std::cout << "Eta range " << -etamax_ << ":" << etamax_ << " # of vtx bins " << nVxBins() << std::endl; + if ((mode_ / 1) % 2 == 0) { + for (int ieta = -etamax_; ieta <= etamax_; ++ieta) { + int eta = (ieta > 0) ? ieta : -ieta; + if (eta != 0) { + int ndepth = ((kdepth_ == 0) ? nDepthBins(eta, 63, modeLHC_) + : ((kdepth_ != 1) ? nDepthBins(eta, 63, 0) : (eta == 16) ? 2 : 1)); + std::cout << "Eta " << ieta << " with " << nPhiBins(eta) << " phi bins " << ndepth << " maximum depths and " + << nPBins(eta) << " p bins" << std::endl; + } } } } char name[100], title[200]; unsigned int book1(0), book2(0); - if ((mode_/1)%2 == 1) { - h_p_.clear(); h_nv_.clear(); h_pnv_.clear(); h_nv2_.clear(); h_p2_.clear(); - h_ediff_.clear(); h_ediff_nvtx_.clear(); - for (int ieta=-etamax_; ieta<=etamax_; ++ieta) { + if ((mode_ / 1) % 2 == 1) { + h_p_.clear(); + h_nv_.clear(); + h_pnv_.clear(); + h_nv2_.clear(); + h_p2_.clear(); + h_ediff_.clear(); + h_ediff_nvtx_.clear(); + for (int ieta = -etamax_; ieta <= etamax_; ++ieta) { if (ieta != 0) { - int zside = (ieta>0) ? 1 : -1; - int eta = (ieta>0) ? ieta : -ieta; - unsigned int id0 = packID(zside,eta,1,1,1,1); - sprintf(name, "peta%d", ieta); - sprintf(title,"Momentum for i#eta = %d (GeV)", ieta); - h_p_[id0] = new TH1D(name, title, 100, 0.0, 500.0); - ++book1; - sprintf(name, "Ediff_eta%d", ieta); - sprintf(title,"Energy difference for i#eta = %d (GeV)", ieta); + int zside = (ieta > 0) ? 1 : -1; + int eta = (ieta > 0) ? ieta : -ieta; + unsigned int id0 = packID(zside, eta, 1, 1, 1, 1); + sprintf(name, "peta%d", ieta); + sprintf(title, "Momentum for i#eta = %d (GeV)", ieta); + h_p_[id0] = new TH1D(name, title, 100, 0.0, 500.0); + ++book1; + sprintf(name, "Ediff_eta%d", ieta); + sprintf(title, "Energy difference for i#eta = %d (GeV)", ieta); h_ediff_[id0] = new TH1D(name, title, 1000, -500.0, 500.0); ++book1; - sprintf(name, "nveta%d", ieta); - sprintf(title,"Number of Vertex for i#eta = %d", ieta); - h_nv_[id0] = new TH1D(name, title, 100, 0.0, 100.0); - ++book1; - sprintf(name, "pnveta%d", ieta); - sprintf(title,"i#eta = %d", ieta); - TH2D* h2 = new TH2D(name, title, 100, 0.0, 100.0, 100, 0.0, 500.0); - ++book2; - h2->GetXaxis()->SetTitle("Number of Vertex"); - h2->GetYaxis()->SetTitle("Momentum (GeV)"); - h_pnv_[id0] = h2; - ++book1; - char etas[10]; - sprintf (etas, "i#eta=%d", ieta); - char name[100], title[200]; - for (int pbin=0; pbin= 0 && eta < (int)(iprange_.size())) ? - iprange_[eta]-1 : iprange_[0]; - sprintf (ps, "p=%d:%d", (int)prange_[np][pbin], - (int)prange_[np][pbin+1]); - }; - unsigned int id = packID(zside,eta,1,1,1,pbin); - sprintf (name,"pvx%d111%d",ieta,pbin); - sprintf (title,"Momentum for %s %s", etas, ps); - h_p2_[id] = new TH1D(name,title,100,0.0,500.0); - h_p2_[id]->Sumw2(); - ++book1; - } - for (int vbin=0; vbinSumw2(); - ++book1; - sprintf (name,"Ediff_nvx%d11%d1",ieta,vbin); - sprintf (title,"Energy difference for %s %s", etas, vtx); - h_ediff_nvtx_[id] = new TH1D(name,title,1000,-500.0,500.0); + sprintf(name, "nveta%d", ieta); + sprintf(title, "Number of Vertex for i#eta = %d", ieta); + h_nv_[id0] = new TH1D(name, title, 100, 0.0, 100.0); + ++book1; + sprintf(name, "pnveta%d", ieta); + sprintf(title, "i#eta = %d", ieta); + TH2D* h2 = new TH2D(name, title, 100, 0.0, 100.0, 100, 0.0, 500.0); + ++book2; + h2->GetXaxis()->SetTitle("Number of Vertex"); + h2->GetYaxis()->SetTitle("Momentum (GeV)"); + h_pnv_[id0] = h2; + ++book1; + char etas[10]; + sprintf(etas, "i#eta=%d", ieta); + char name[100], title[200]; + for (int pbin = 0; pbin < nPBins(eta); ++pbin) { + char ps[20]; + if ((mode_ / 4) % 2 == 1) { + int np = (eta >= 0 && eta < (int)(iprange_.size())) ? iprange_[eta] - 1 : iprange_[0]; + sprintf(ps, "p=%d:%d", (int)prange_[np][pbin], (int)prange_[np][pbin + 1]); + }; + unsigned int id = packID(zside, eta, 1, 1, 1, pbin); + sprintf(name, "pvx%d111%d", ieta, pbin); + sprintf(title, "Momentum for %s %s", etas, ps); + h_p2_[id] = new TH1D(name, title, 100, 0.0, 500.0); + h_p2_[id]->Sumw2(); + ++book1; + } + for (int vbin = 0; vbin < nVxBins(); ++vbin) { + char vtx[12]; + if ((mode_ / 2) % 2 == 1) { + sprintf(vtx, "N_{vtx}=%d:%d", npvbin_[vbin], npvbin_[vbin + 1]); + } else { + sprintf(vtx, "all N_{vtx}"); + } + unsigned int id = packID(zside, eta, 1, 1, vbin, 1); + sprintf(name, "nvx%d11%d1", ieta, vbin); + sprintf(title, "Number of vertex for %s %s", etas, vtx); + h_nv2_[id] = new TH1D(name, title, 100, 0.0, 100.0); + h_nv2_[id]->Sumw2(); + ++book1; + sprintf(name, "Ediff_nvx%d11%d1", ieta, vbin); + sprintf(title, "Energy difference for %s %s", etas, vtx); + h_ediff_nvtx_[id] = new TH1D(name, title, 1000, -500.0, 500.0); h_ediff_nvtx_[id]->Sumw2(); ++book1; - } + } } } } else { - h_Energy_.clear(); h_Ecorr_.clear(); h_Charge_.clear(); h_Chcorr_.clear(); - h_EnergyC_.clear(); h_EcorrC_.clear(); - for (int ieta=-etamax_; ieta<=etamax_; ++ieta) { + h_Energy_.clear(); + h_Ecorr_.clear(); + h_Charge_.clear(); + h_Chcorr_.clear(); + h_EnergyC_.clear(); + h_EcorrC_.clear(); + for (int ieta = -etamax_; ieta <= etamax_; ++ieta) { if (ieta != 0) { - int zside = (ieta>0) ? 1 : -1; - int eta = (ieta>0) ? ieta : -ieta; - char etas[20]; - sprintf (etas, "i#eta=%d", ieta); - for (int iphi=0; iphi= 0 && eta < (int)(iprange_.size())) ? - iprange_[eta]-1 : iprange_[0]; - sprintf (ps, "p=%d:%d", (int)prange_[np][pbin], (int)prange_[np][pbin+1]); - } else { - sprintf (ps, "all p"); - }; - for (int vbin=0; vbin 0) ? 1 : -1; + int eta = (ieta > 0) ? ieta : -ieta; + char etas[20]; + sprintf(etas, "i#eta=%d", ieta); + for (int iphi = 0; iphi < nPhiBins(eta); ++iphi) { + char phis[20]; + int phi(1), phi0(63); + if (kphi_ == 1) { + phi0 = phi = (eta <= 20) ? iphi + 1 : 2 * iphi + 1; + sprintf(phis, "i#phi=%d", phi); + } else if (kphi_ == 2) { + phi0 = 4 * iphi + 1; + phi = iphi + 1; + sprintf(phis, "RBX=%d", iphi + 1); + } else if (kphi_ == 3) { + sprintf(phis, "All except RBX %d", exRBX_); + } else { + sprintf(phis, "All i#phi"); + } + int ndepth = ((kdepth_ == 0) ? nDepthBins(eta, phi0, modeLHC_) + : ((kdepth_ != 1) ? nDepthBins(eta, phi0, 0) : (eta == 16) ? 2 : 1)); + for (int depth = 0; depth < ndepth; ++depth) { + char deps[20]; + if (kdepth_ == 1) { + if (depth == 0) + sprintf(deps, "all depths"); + else + sprintf(deps, "all endcap depths"); + } else { + sprintf(deps, "Depth=%d", depth + 1); + } + for (int pbin = 0; pbin < nPBins(eta); ++pbin) { + char ps[30]; + if ((mode_ / 4) % 2 == 1) { + int np = (eta >= 0 && eta < (int)(iprange_.size())) ? iprange_[eta] - 1 : iprange_[0]; + sprintf(ps, "p=%d:%d", (int)prange_[np][pbin], (int)prange_[np][pbin + 1]); + } else { + sprintf(ps, "all p"); + }; + for (int vbin = 0; vbin < nVxBins(); ++vbin) { + int nbin(4000); + double xmax(10.0); + char vtx[20]; + if ((mode_ / 2) % 2 == 1) { + sprintf(vtx, "N_{vtx}=%d:%d", npvbin_[vbin], npvbin_[vbin + 1]); + } else { + sprintf(vtx, "all N_{vtx}"); + } + unsigned int id = packID(zside, eta, phi, depth + 1, vbin, pbin); + char name[100], title[200]; + sprintf(name, "EdepE%dF%dD%dV%dP%d", ieta, phi, depth, vbin, pbin); + sprintf(title, "Deposited energy for %s %s %s %s %s (GeV)", etas, phis, deps, ps, vtx); + getBins(0, ieta, phi0, depth + 1, nbin, xmax); + h_Energy_[id] = new TH1D(name, title, nbin, 0.0, xmax); + ++book1; + sprintf(name, "EcorE%dF%dD%dV%dP%d", ieta, phi, depth, vbin, pbin); + sprintf(title, "Active length corrected energy for %s %s %s %s %s (GeV/cm)", etas, phis, deps, ps, vtx); + getBins(1, ieta, phi0, depth + 1, nbin, xmax); + h_Ecorr_[id] = new TH1D(name, title, nbin, 0.0, xmax); + ++book1; + sprintf(name, "EdepCE%dF%dD%dV%dP%d", ieta, phi, depth, vbin, pbin); + sprintf( + title, "Response Corrected deposited energy for %s %s %s %s %s (GeV)", etas, phis, deps, ps, vtx); + getBins(2, ieta, phi0, depth + 1, nbin, xmax); + h_EnergyC_[id] = new TH1D(name, title, nbin, 0.0, xmax); + ++book1; + sprintf(name, "EcorCE%dF%dD%dV%dP%d", ieta, phi, depth, vbin, pbin); + sprintf(title, + "Response and active length corrected energy for %s %s %s %s %s (GeV/cm)", + etas, + phis, + deps, + ps, + vtx); + getBins(3, ieta, phi0, depth + 1, nbin, xmax); + h_EcorrC_[id] = new TH1D(name, title, nbin, 0.0, xmax); + ++book1; + sprintf(name, "ChrgE%dF%dD%dV%dP%d", ieta, phi, depth, vbin, pbin); + sprintf(title, "Measured charge for %s %s %s %s %s (cm)", etas, phis, deps, ps, vtx); + getBins(4, ieta, phi0, depth + 1, nbin, xmax); + h_Charge_[id] = new TH1D(name, title, nbin, 0.0, xmax); + ++book1; + sprintf(name, "ChcorE%dF%dD%dV%dP%d", ieta, phi, depth, vbin, pbin); + sprintf(title, "Active length corrected charge for %s %s %s %s %s (cm)", etas, phis, deps, ps, vtx); + getBins(5, ieta, phi0, depth + 1, nbin, xmax); + h_Chcorr_[id] = new TH1D(name, title, nbin, 0.0, xmax); + ++book1; + } + } + } + } } } } @@ -750,65 +791,69 @@ void AnalyzeLepTree::bookHisto() { } void AnalyzeLepTree::writeHisto(const char* outfile) { - - TFile* output_file = TFile::Open(outfile,"RECREATE"); + TFile* output_file = TFile::Open(outfile, "RECREATE"); output_file->cd(); - if ((mode_/1)%2 == 1) { - for (std::map::const_iterator itr = h_p_.begin(); - itr != h_p_.end(); ++itr) (itr->second)->Write(); - for (std::map::const_iterator itr = h_nv_.begin(); - itr != h_nv_.end(); ++itr) (itr->second)->Write(); - for (std::map::const_iterator itr = h_pnv_.begin(); - itr != h_pnv_.end(); ++itr) (itr->second)->Write(); - for (std::map::const_iterator itr = h_p2_.begin(); - itr != h_p2_.end(); ++itr) (itr->second)->Write(); - for (std::map::const_iterator itr = h_nv2_.begin(); - itr != h_nv2_.end(); ++itr) (itr->second)->Write(); - for (std::map::const_iterator itr = h_ediff_.begin(); - itr != h_ediff_.end(); ++itr) (itr->second)->Write(); - for (std::map::const_iterator itr = h_ediff_nvtx_.begin(); - itr != h_ediff_nvtx_.end(); ++itr) (itr->second)->Write(); + if ((mode_ / 1) % 2 == 1) { + for (std::map::const_iterator itr = h_p_.begin(); itr != h_p_.end(); ++itr) + (itr->second)->Write(); + for (std::map::const_iterator itr = h_nv_.begin(); itr != h_nv_.end(); ++itr) + (itr->second)->Write(); + for (std::map::const_iterator itr = h_pnv_.begin(); itr != h_pnv_.end(); ++itr) + (itr->second)->Write(); + for (std::map::const_iterator itr = h_p2_.begin(); itr != h_p2_.end(); ++itr) + (itr->second)->Write(); + for (std::map::const_iterator itr = h_nv2_.begin(); itr != h_nv2_.end(); ++itr) + (itr->second)->Write(); + for (std::map::const_iterator itr = h_ediff_.begin(); itr != h_ediff_.end(); ++itr) + (itr->second)->Write(); + for (std::map::const_iterator itr = h_ediff_nvtx_.begin(); itr != h_ediff_nvtx_.end(); ++itr) + (itr->second)->Write(); } else { - for (int ieta=-etamax_; ieta<=etamax_; ++ieta) { + for (int ieta = -etamax_; ieta <= etamax_; ++ieta) { if (ieta != 0) { - char dirname[50]; - sprintf (dirname, "DieMuonEta%d", ieta); - TDirectory* d_output = output_file->mkdir(dirname); - d_output->cd(); - int zside = (ieta>0) ? 1 : -1; - int eta = (ieta>0) ? ieta : -ieta; - for (int iphi=0; iphi::const_iterator itr; - itr = h_Energy_.find(id); - if (itr != h_Energy_.end()) (itr->second)->Write(); - itr = h_Ecorr_.find(id); - if (itr != h_Ecorr_.end()) (itr->second)->Write(); - itr = h_EnergyC_.find(id); - if (itr != h_EnergyC_.end()) (itr->second)->Write(); - itr = h_EcorrC_.find(id); - if (itr != h_EcorrC_.end()) (itr->second)->Write(); - itr = h_Charge_.find(id); - if (itr != h_Charge_.end()) (itr->second)->Write(); - itr = h_Chcorr_.find(id); - if (itr != h_Chcorr_.end()) (itr->second)->Write(); - } - } - } - } + char dirname[50]; + sprintf(dirname, "DieMuonEta%d", ieta); + TDirectory* d_output = output_file->mkdir(dirname); + d_output->cd(); + int zside = (ieta > 0) ? 1 : -1; + int eta = (ieta > 0) ? ieta : -ieta; + for (int iphi = 0; iphi < nPhiBins(eta); ++iphi) { + int phi(1), phi0(1); + if (kphi_ == 1) { + phi0 = phi = (eta <= 20) ? iphi + 1 : 2 * iphi + 1; + } else if (kphi_ == 2) { + phi0 = 4 * iphi + 1; + phi = iphi + 1; + }; + int ndepth = ((kdepth_ == 0) ? nDepthBins(eta, phi0, modeLHC_) + : ((kdepth_ != 1) ? nDepthBins(eta, phi0, 0) : (eta == 16) ? 2 : 1)); + for (int depth = 0; depth < ndepth; ++depth) { + for (int pbin = 0; pbin < nPBins(eta); ++pbin) { + for (int vbin = 0; vbin < nVxBins(); ++vbin) { + unsigned int id = packID(zside, eta, phi, depth + 1, vbin, pbin); + std::map::const_iterator itr; + itr = h_Energy_.find(id); + if (itr != h_Energy_.end()) + (itr->second)->Write(); + itr = h_Ecorr_.find(id); + if (itr != h_Ecorr_.end()) + (itr->second)->Write(); + itr = h_EnergyC_.find(id); + if (itr != h_EnergyC_.end()) + (itr->second)->Write(); + itr = h_EcorrC_.find(id); + if (itr != h_EcorrC_.end()) + (itr->second)->Write(); + itr = h_Charge_.find(id); + if (itr != h_Charge_.end()) + (itr->second)->Write(); + itr = h_Chcorr_.find(id); + if (itr != h_Chcorr_.end()) + (itr->second)->Write(); + } + } + } + } } output_file->cd(); } @@ -816,118 +861,123 @@ void AnalyzeLepTree::writeHisto(const char* outfile) { } void AnalyzeLepTree::writeMeanError(const char* outfile) { - - if ((mode_/1)%2 == 1) { + if ((mode_ / 1) % 2 == 1) { std::ofstream fOutput(outfile); - for (int vbin=0; vbin0) ? 1 : -1; - int eta = (ieta>0) ? ieta : -ieta; - unsigned int id = packID(zside,eta,1,1,vbin,1); - char name[100]; - sprintf (name,"nvx%d11%d1",ieta,vbin); - std::map::iterator itr = h_nv2_.find(id); - if (itr != h_nv2_.end()) { - double mean = (itr->second)->GetMean(); - double error= (itr->second)->GetMeanError(); - char vtx[12]; - if ((mode_/2)%2 == 1) { - sprintf(vtx, "Nvtx=%3d:%3d", npvbin_[vbin], npvbin_[vbin+1]); - } else { - sprintf(vtx, "all Nvtx"); - } - fOutput << "Eta " << std::setw(3) << ieta << " " << vtx << " " - << mean << " +- " << error << std::endl; - } - } + for (int vbin = 0; vbin < nVxBins(); ++vbin) { + for (int ieta = -etamax_; ieta <= etamax_; ++ieta) { + if (ieta != 0) { + int zside = (ieta > 0) ? 1 : -1; + int eta = (ieta > 0) ? ieta : -ieta; + unsigned int id = packID(zside, eta, 1, 1, vbin, 1); + char name[100]; + sprintf(name, "nvx%d11%d1", ieta, vbin); + std::map::iterator itr = h_nv2_.find(id); + if (itr != h_nv2_.end()) { + double mean = (itr->second)->GetMean(); + double error = (itr->second)->GetMeanError(); + char vtx[12]; + if ((mode_ / 2) % 2 == 1) { + sprintf(vtx, "Nvtx=%3d:%3d", npvbin_[vbin], npvbin_[vbin + 1]); + } else { + sprintf(vtx, "all Nvtx"); + } + fOutput << "Eta " << std::setw(3) << ieta << " " << vtx << " " << mean << " +- " << error << std::endl; + } + } } } fOutput.close(); } } -std::vector AnalyzeLepTree::plotHisto(bool drawStatBox, int type, - bool save) { - +std::vector AnalyzeLepTree::plotHisto(bool drawStatBox, int type, bool save) { std::vector cvs; - gStyle->SetCanvasBorderMode(0); gStyle->SetCanvasColor(kWhite); - gStyle->SetPadColor(kWhite); gStyle->SetFillColor(kWhite); + gStyle->SetCanvasBorderMode(0); + gStyle->SetCanvasColor(kWhite); + gStyle->SetPadColor(kWhite); + gStyle->SetFillColor(kWhite); gStyle->SetOptTitle(0); if (drawStatBox) { - gStyle->SetOptStat(111110); gStyle->SetOptFit(1); + gStyle->SetOptStat(111110); + gStyle->SetOptFit(1); } else { - gStyle->SetOptStat(0); gStyle->SetOptFit(0); + gStyle->SetOptStat(0); + gStyle->SetOptFit(0); } - if ((mode_/1)%2 == 1) { - if (type%2 > 0) plotHisto(h_p_, cvs, save); - if ((type/2)%2 > 0) plotHisto(h_nv_, cvs, save); - if ((type/4)%2 > 0) plot2DHisto(h_pnv_, cvs, save); - if ((type/8)%2 > 0) plotHisto(h_nv2_, cvs, save); - if ((type/16)%2 > 0) plotHisto(h_p2_, cvs, save); - if ((type/32)%2 > 0) plotHisto(h_ediff_, cvs, save); - if ((type/32)%2 > 0) plotHisto(h_ediff_nvtx_, cvs, save); + if ((mode_ / 1) % 2 == 1) { + if (type % 2 > 0) + plotHisto(h_p_, cvs, save); + if ((type / 2) % 2 > 0) + plotHisto(h_nv_, cvs, save); + if ((type / 4) % 2 > 0) + plot2DHisto(h_pnv_, cvs, save); + if ((type / 8) % 2 > 0) + plotHisto(h_nv2_, cvs, save); + if ((type / 16) % 2 > 0) + plotHisto(h_p2_, cvs, save); + if ((type / 32) % 2 > 0) + plotHisto(h_ediff_, cvs, save); + if ((type / 32) % 2 > 0) + plotHisto(h_ediff_nvtx_, cvs, save); } else { - int depth0 = (type/16)&15; - int phi0 = (type/256)&127; - bool doEn = ((type/1)%2 > 0); - bool doEnL = ((type/2)%2 > 0); - bool doChg = ((type/4)%2 > 0); - bool doChL = ((type/8)%2 > 0); - if (doEn) plotHisto(h_Energy_, phi0, depth0, cvs, save); - if (doEn) plotHisto(h_EnergyC_, phi0, depth0, cvs, save); - if (doEnL) plotHisto(h_Ecorr_, phi0, depth0, cvs, save); - if (doEnL) plotHisto(h_EcorrC_, phi0, depth0, cvs, save); - if (doChg) plotHisto(h_Charge_, phi0, depth0, cvs, save); - if (doChL) plotHisto(h_Chcorr_, phi0, depth0, cvs, save); + int depth0 = (type / 16) & 15; + int phi0 = (type / 256) & 127; + bool doEn = ((type / 1) % 2 > 0); + bool doEnL = ((type / 2) % 2 > 0); + bool doChg = ((type / 4) % 2 > 0); + bool doChL = ((type / 8) % 2 > 0); + if (doEn) + plotHisto(h_Energy_, phi0, depth0, cvs, save); + if (doEn) + plotHisto(h_EnergyC_, phi0, depth0, cvs, save); + if (doEnL) + plotHisto(h_Ecorr_, phi0, depth0, cvs, save); + if (doEnL) + plotHisto(h_EcorrC_, phi0, depth0, cvs, save); + if (doChg) + plotHisto(h_Charge_, phi0, depth0, cvs, save); + if (doChL) + plotHisto(h_Chcorr_, phi0, depth0, cvs, save); } return cvs; } -void AnalyzeLepTree::plotHisto(std::map hists, - std::vector& cvs, bool save) { - - for (std::map::const_iterator itr = hists.begin(); - itr != hists.end(); ++itr) { +void AnalyzeLepTree::plotHisto(std::map hists, std::vector& cvs, bool save) { + for (std::map::const_iterator itr = hists.begin(); itr != hists.end(); ++itr) { TH1D* hist = (itr->second); if (hist != 0) { TCanvas* pad = plotHisto(hist); cvs.push_back(pad); if (save) { - char name[100]; - sprintf (name, "c_%s.pdf", pad->GetName()); - pad->Print(name); - } + char name[100]; + sprintf(name, "c_%s.pdf", pad->GetName()); + pad->Print(name); + } } } } -void AnalyzeLepTree::plotHisto(std::map hists, int phi0, - int depth0, std::vector& cvs, - bool save) { - - for (std::map::const_iterator itr=hists.begin(); - itr != hists.end(); ++itr) { +void AnalyzeLepTree::plotHisto( + std::map hists, int phi0, int depth0, std::vector& cvs, bool save) { + for (std::map::const_iterator itr = hists.begin(); itr != hists.end(); ++itr) { int zside, eta, phi, depth, pbin, vbin; - unpackID(itr->first,zside,eta,phi,depth,vbin,pbin); + unpackID(itr->first, zside, eta, phi, depth, vbin, pbin); TH1D* hist = itr->second; - if (((depth == depth0) || (depth0 == 0)) && - ((phi == phi0) || (phi0 == 0)) && (hist != 0)) { + if (((depth == depth0) || (depth0 == 0)) && ((phi == phi0) || (phi0 == 0)) && (hist != 0)) { TCanvas* pad = plotHisto(hist); cvs.push_back(pad); if (save) { - char name[100]; - sprintf (name, "c_%s.pdf", pad->GetName()); - pad->Print(name); - } + char name[100]; + sprintf(name, "c_%s.pdf", pad->GetName()); + pad->Print(name); + } } } } TCanvas* AnalyzeLepTree::plotHisto(TH1D* hist) { - - TCanvas *pad = new TCanvas(hist->GetName(), hist->GetName(), 700, 500); + TCanvas* pad = new TCanvas(hist->GetName(), hist->GetName(), 700, 500); pad->SetRightMargin(0.10); pad->SetTopMargin(0.10); hist->GetXaxis()->SetTitleSize(0.04); @@ -943,22 +993,22 @@ TCanvas* AnalyzeLepTree::plotHisto(TH1D* hist) { pad->Update(); TPaveStats* st1 = (TPaveStats*)hist->GetListOfFunctions()->FindObject("stats"); if (st1 != NULL) { - st1->SetY1NDC(0.70); st1->SetY2NDC(0.90); - st1->SetX1NDC(0.65); st1->SetX2NDC(0.90); + st1->SetY1NDC(0.70); + st1->SetY2NDC(0.90); + st1->SetX1NDC(0.65); + st1->SetX2NDC(0.90); } pad->Modified(); pad->Update(); return pad; } -void AnalyzeLepTree::plot2DHisto(std::map hists, - std::vector& cvs, bool save) { +void AnalyzeLepTree::plot2DHisto(std::map hists, std::vector& cvs, bool save) { char name[100]; - for (std::map::const_iterator itr = hists.begin(); - itr != hists.end(); ++itr) { + for (std::map::const_iterator itr = hists.begin(); itr != hists.end(); ++itr) { TH2D* hist = (itr->second); if (hist != 0) { - TCanvas *pad = new TCanvas(hist->GetName(), hist->GetName(), 700, 700); + TCanvas* pad = new TCanvas(hist->GetName(), hist->GetName(), 700, 700); pad->SetRightMargin(0.10); pad->SetTopMargin(0.10); hist->GetXaxis()->SetTitleSize(0.04); @@ -970,53 +1020,65 @@ void AnalyzeLepTree::plot2DHisto(std::map hists, pad->Update(); TPaveStats* st1 = (TPaveStats*)hist->GetListOfFunctions()->FindObject("stats"); if (st1 != NULL) { - st1->SetY1NDC(0.65); st1->SetY2NDC(0.90); - st1->SetX1NDC(0.65); st1->SetX2NDC(0.90); + st1->SetY1NDC(0.65); + st1->SetY2NDC(0.90); + st1->SetX1NDC(0.65); + st1->SetX2NDC(0.90); } pad->Modified(); pad->Update(); cvs.push_back(pad); if (save) { - sprintf (name, "c_%s.pdf", pad->GetName()); - pad->Print(name); + sprintf(name, "c_%s.pdf", pad->GetName()); + pad->Print(name); } } - } + } } int AnalyzeLepTree::getCollapsedDepth(int eta, int phi, int dep) { - int depth = dep+1; - int ieta = (eta > 0) ? eta : -eta; + int depth = dep + 1; + int ieta = (eta > 0) ? eta : -eta; if (ieta <= 14 || ieta == 17) { depth = 1; } else if (ieta == 15) { if (modeLHC_ > 3) { - if (dep > 3) depth = 2; - else depth = 1; + if (dep > 3) + depth = 2; + else + depth = 1; } } else if (ieta == 16) { if (modeLHC_ == 0 || (modeLHC_ == 3 && (phi < 63 || phi > 66 || eta < 0))) { } else { - if (dep > 2) depth = 3; + if (dep > 2) + depth = 3; } } else if (ieta < 26) { if (modeLHC_ == 0 || (modeLHC_ == 3 && (phi < 63 || phi > 66 || eta < 0))) { } else { - if (dep < 3) depth = 1; - else depth = 2; + if (dep < 3) + depth = 1; + else + depth = 2; } } else if (ieta == 26) { if (modeLHC_ == 0 || (modeLHC_ == 3 && (phi < 63 || phi > 66 || eta < 0))) { } else { - if (dep < 4) depth = 1; - else depth = 2; + if (dep < 4) + depth = 1; + else + depth = 2; } } else { if (modeLHC_ == 0 || (modeLHC_ == 3 && (phi < 63 || phi > 66 || eta < 0))) { } else { - if (dep < 3) depth = 1; - else if (dep == 3) depth = 2; - else depth = 3; + if (dep < 3) + depth = 1; + else if (dep == 3) + depth = 2; + else + depth = 3; } } return depth; @@ -1024,17 +1086,19 @@ int AnalyzeLepTree::getCollapsedDepth(int eta, int phi, int dep) { int AnalyzeLepTree::getRBX(int eta) { int rbx(1); - int phi = (eta > 20) ? (2*t_iphi + 1) : (t_iphi + 1); - if (phi > 2 && phi < 71) rbx = (phi+5)/4; + int phi = (eta > 20) ? (2 * t_iphi + 1) : (t_iphi + 1); + if (phi > 2 && phi < 71) + rbx = (phi + 5) / 4; return rbx; } int AnalyzeLepTree::getPBin(int eta) { int bin(0); - if ((mode_/4)%2 == 1) { + if ((mode_ / 4) % 2 == 1) { int np = (eta >= 0 && eta < (int)(iprange_.size())) ? iprange_[eta] : iprange_[0]; - for (unsigned int k=1; k 20) ? (2*t_iphi + 1) : (t_iphi + 1); + bin = (eta > 20) ? (2 * t_iphi + 1) : (t_iphi + 1); } else if (kphi_ == 2) { bin = getRBX(eta); } else if (kphi_ == 3) { - if (exRBX_ == getRBX(eta)) bin = -1; + if (exRBX_ == getRBX(eta)) + bin = -1; } return bin; } @@ -1075,52 +1141,56 @@ void AnalyzeLepTree::makeVxBins(int modeLHC) { int npvbin2[nvbin_] = {0, 25, 30, 35, 40, 100}; int npvbin3[nvbin_] = {0, 30, 40, 50, 70, 200}; npvbin_.clear(); - for (int i=0; i corresponds to Run 1 (valid till 2016) // = 1 --> corresponds to Run 2 (2018 geometry) // = 2 --> corresponds to Run 3 (post LS2) // = 3 --> corresponds to 2017 (Plan 1) // = 4 --> corresponds to Run4 (post LS3) - int nbin(0); + int nbin(0); if (modeLHC == 0) { - nbin = nDepthR1[eta-1]; + nbin = nDepthR1[eta - 1]; } else if (modeLHC == 1) { - nbin = nDepthR2[eta-1]; + nbin = nDepthR2[eta - 1]; } else if (modeLHC == 2) { - nbin = nDepthR3[eta-1]; + nbin = nDepthR3[eta - 1]; } else if (modeLHC == 3) { if (phi > 0) { if (eta >= 16 && phi >= 63 && phi <= 66) { - nbin = nDepthR2[eta-1]; + nbin = nDepthR2[eta - 1]; } else { - nbin = nDepthR1[eta-1]; + nbin = nDepthR1[eta - 1]; } } else { if (eta >= 16) { - nbin = (nDepthR2[eta-1] > nDepthR1[eta-1]) ? nDepthR2[eta-1] : nDepthR1[eta-1]; + nbin = (nDepthR2[eta - 1] > nDepthR1[eta - 1]) ? nDepthR2[eta - 1] : nDepthR1[eta - 1]; } else { - nbin = nDepthR1[eta-1]; + nbin = nDepthR1[eta - 1]; } } } else { if (eta > 0 && eta < 30) { - nbin = nDepthR4[eta-1]; + nbin = nDepthR4[eta - 1]; } else { nbin = nDepthR4[28]; } @@ -1130,55 +1200,54 @@ int AnalyzeLepTree::nDepthBins(int eta, int phi, int modeLHC) { int AnalyzeLepTree::nPhiBins(int eta) { int nphi = (eta <= 20) ? 72 : 36; - if (modeLHC_ == 4 && eta > 16) nphi = 360; - if (kphi_ == 0) nphi = 1; - else if (kphi_ == 2) nphi = 18; - else if (kphi_ == 3) nphi = 1; + if (modeLHC_ == 4 && eta > 16) + nphi = 360; + if (kphi_ == 0) + nphi = 1; + else if (kphi_ == 2) + nphi = 18; + else if (kphi_ == 3) + nphi = 1; return nphi; } int AnalyzeLepTree::nPBins(int eta) { int bin(1); - if ((mode_/4)%2 == 1) { - int np = (eta >= 0 && eta < (int)(iprange_.size())) ? iprange_[eta]-1 : iprange_[0]; - bin = (int)(prange_[np].size())-1; + if ((mode_ / 4) % 2 == 1) { + int np = (eta >= 0 && eta < (int)(iprange_.size())) ? iprange_[eta] - 1 : iprange_[0]; + bin = (int)(prange_[np].size()) - 1; } return bin; } int AnalyzeLepTree::nVxBins() { - int nvx = ((mode_/2)%2 == 1) ? ((int)(npvbin_.size())-1) : 1; + int nvx = ((mode_ / 2) % 2 == 1) ? ((int)(npvbin_.size()) - 1) : 1; return nvx; } -unsigned int AnalyzeLepTree::packID(int zside, int eta, int phi, int depth, - int nvx, int ipbin) { +unsigned int AnalyzeLepTree::packID(int zside, int eta, int phi, int depth, int nvx, int ipbin) { unsigned int id = (zside > 0) ? 1 : 0; - id |= (((nvx&7)<<19) | ((ipbin&7)<<16) | ((depth&7)<<13) | ((eta&31)<<8) | - ((phi&127)<<1)); + id |= (((nvx & 7) << 19) | ((ipbin & 7) << 16) | ((depth & 7) << 13) | ((eta & 31) << 8) | ((phi & 127) << 1)); return id; } -void AnalyzeLepTree::unpackID(unsigned int id, int& zside, int& eta, int& phi, - int& depth, int& nvx, int& ipbin) { - zside = (id%2 == 0) ? -1 : 1; - phi = (id >> 1) & 127; - eta = (id >> 8) & 31; +void AnalyzeLepTree::unpackID(unsigned int id, int& zside, int& eta, int& phi, int& depth, int& nvx, int& ipbin) { + zside = (id % 2 == 0) ? -1 : 1; + phi = (id >> 1) & 127; + eta = (id >> 8) & 31; depth = (id >> 13) & 7; ipbin = (id >> 16) & 7; - nvx = (id >> 19) & 7; + nvx = (id >> 19) & 7; } -void AnalyzeLepTree::getBins(int type, int ieta, int phi, int depth, int& nbin, - double& xmax) { - int eta = (ieta >= 0) ? ieta : -ieta; - bool barrel = (eta < 16) || ((eta == 16) && (depth <= 2)); - bool rbx17 = (phi >= 63) && (phi <= 66) && (ieta >= 16) && (!barrel); +void AnalyzeLepTree::getBins(int type, int ieta, int phi, int depth, int& nbin, double& xmax) { + int eta = (ieta >= 0) ? ieta : -ieta; + bool barrel = (eta < 16) || ((eta == 16) && (depth <= 2)); + bool rbx17 = (phi >= 63) && (phi <= 66) && (ieta >= 16) && (!barrel); nbin = 5000; xmax = 10.0; if (type >= 4) { - if ((modeLHC_ == 0) || (((modeLHC_ == 1) || (modeLHC_ == 3)) && barrel) || - ((modeLHC_ == 3) && (!rbx17))) { + if ((modeLHC_ == 0) || (((modeLHC_ == 1) || (modeLHC_ == 3)) && barrel) || ((modeLHC_ == 3) && (!rbx17))) { // HPD Channels xmax = 50.0; } else { diff --git a/Calibration/HcalCalibAlgos/macros/HBHEMuonHighEta.C b/Calibration/HcalCalibAlgos/macros/HBHEMuonHighEta.C index 37583bf95f6b3..98551fdc150d8 100644 --- a/Calibration/HcalCalibAlgos/macros/HBHEMuonHighEta.C +++ b/Calibration/HcalCalibAlgos/macros/HBHEMuonHighEta.C @@ -4,7 +4,7 @@ // h1.Loop() // // infile const char* Name of the input file -// outfile const char* Name of the output file +// outfile const char* Name of the output file // (dyll_PU20_25_output_10.root) // mode int Geometry file used 0:(defined by maxDHB/HE); // 1 (Run 1; valid till 2016); 2 (Run 2; 2018); @@ -32,272 +32,272 @@ #include class HBHEMuonHighEta { - -public : - HBHEMuonHighEta(const char *infile, const char *outfile, - const int mode=0, const bool debug=false); +public: + HBHEMuonHighEta(const char *infile, const char *outfile, const int mode = 0, const bool debug = false); virtual ~HBHEMuonHighEta(); - virtual Int_t Cut(Long64_t entry); - virtual Int_t GetEntry(Long64_t entry); + virtual Int_t Cut(Long64_t entry); + virtual Int_t GetEntry(Long64_t entry); virtual Long64_t LoadTree(Long64_t entry); - virtual void Init(TTree *tree); - virtual void Loop(); - virtual Bool_t Notify(); - virtual void Show(Long64_t entry = -1); + virtual void Init(TTree *tree); + virtual void Loop(); + virtual Bool_t Notify(); + virtual void Show(Long64_t entry = -1); private: - void BookHistograms(const char* fname); - void Close(); - int nDepthBins(int eta, int phi); - TTree *fChain; //!pointer to the analyzed TTree or TChain - Int_t fCurrent; //!current Tree number in a TChain + void BookHistograms(const char *fname); + void Close(); + int nDepthBins(int eta, int phi); + TTree *fChain; //!pointer to the analyzed TTree or TChain + Int_t fCurrent; //!current Tree number in a TChain - static const int maxDepthHB_ = 7; - static const int maxDepthHE_ = 4; + static const int maxDepthHB_ = 7; + static const int maxDepthHE_ = 4; // Fixed size dimensions of array or collections stored in the TTree if any. // Declaration of leaf types - std::vector *pt_of_muon; - std::vector *eta_of_muon; - std::vector *phi_of_muon; - std::vector *energy_of_muon; - std::vector *p_of_muon; - std::vector *MediumMuon; - std::vector *IsolationR04; - std::vector *IsolationR03; - std::vector *ecal_3into3; - std::vector *hcal_3into3; - std::vector *tracker_3into3; - std::vector *emaxNearP; - UInt_t Run_No; - UInt_t Event_No; - UInt_t GoodVertex; - std::vector *matchedId; - std::vector *hcal_cellHot; - std::vector *ecal_3x3; - std::vector *hcal_1x1; + std::vector *pt_of_muon; + std::vector *eta_of_muon; + std::vector *phi_of_muon; + std::vector *energy_of_muon; + std::vector *p_of_muon; + std::vector *MediumMuon; + std::vector *IsolationR04; + std::vector *IsolationR03; + std::vector *ecal_3into3; + std::vector *hcal_3into3; + std::vector *tracker_3into3; + std::vector *emaxNearP; + UInt_t Run_No; + UInt_t Event_No; + UInt_t GoodVertex; + std::vector *matchedId; + std::vector *hcal_cellHot; + std::vector *ecal_3x3; + std::vector *hcal_1x1; std::vector *ecal_detID; std::vector *hcal_detID; std::vector *ehcal_detID; - std::vector *hcal_ieta; - std::vector *hcal_iphi; - std::vector *hcal_edepth1; - std::vector *hcal_activeL1; - std::vector *hcal_edepthHot1; - std::vector *hcal_activeHotL1; - std::vector *hcal_cdepthHot1; - std::vector *hcal_cdepthHotBG1; - std::vector *hcal_edepthCorrect1; - std::vector *hcal_edepthHotCorrect1; - std::vector *hcal_depthMatch1; - std::vector *hcal_depthMatchHot1; - std::vector *hcal_edepth2; - std::vector *hcal_activeL2; - std::vector *hcal_edepthHot2; - std::vector *hcal_activeHotL2; - std::vector *hcal_cdepthHot2; - std::vector *hcal_cdepthHotBG2; - std::vector *hcal_edepthCorrect2; - std::vector *hcal_edepthHotCorrect2; - std::vector *hcal_depthMatch2; - std::vector *hcal_depthMatchHot2; - std::vector *hcal_edepth3; - std::vector *hcal_activeL3; - std::vector *hcal_edepthHot3; - std::vector *hcal_activeHotL3; - std::vector *hcal_cdepthHot3; - std::vector *hcal_cdepthHotBG3; - std::vector *hcal_edepthCorrect3; - std::vector *hcal_edepthHotCorrect3; - std::vector *hcal_depthMatch3; - std::vector *hcal_depthMatchHot3; - std::vector *hcal_edepth4; - std::vector *hcal_activeL4; - std::vector *hcal_edepthHot4; - std::vector *hcal_activeHotL4; - std::vector *hcal_cdepthHot4; - std::vector *hcal_cdepthHotBG4; - std::vector *hcal_edepthCorrect4; - std::vector *hcal_edepthHotCorrect4; - std::vector *hcal_depthMatch4; - std::vector *hcal_depthMatchHot4; - std::vector *hcal_edepth5; - std::vector *hcal_activeL5; - std::vector *hcal_edepthHot5; - std::vector *hcal_activeHotL5; - std::vector *hcal_cdepthHot5; - std::vector *hcal_cdepthHotBG5; - std::vector *hcal_edepthCorrect5; - std::vector *hcal_edepthHotCorrect5; - std::vector *hcal_depthMatch5; - std::vector *hcal_depthMatchHot5; - std::vector *hcal_edepth6; - std::vector *hcal_activeL6; - std::vector *hcal_edepthHot6; - std::vector *hcal_activeHotL6; - std::vector *hcal_cdepthHot6; - std::vector *hcal_cdepthHotBG6; - std::vector *hcal_edepthCorrect6; - std::vector *hcal_edepthHotCorrect6; - std::vector *hcal_depthMatch6; - std::vector *hcal_depthMatchHot6; - std::vector *hcal_edepth7; - std::vector *hcal_activeL7; - std::vector *hcal_edepthHot7; - std::vector *hcal_activeHotL7; - std::vector *hcal_cdepthHot7; - std::vector *hcal_cdepthHotBG7; - std::vector *hcal_edepthCorrect7; - std::vector *hcal_edepthHotCorrect7; - std::vector *hcal_depthMatch7; - std::vector *hcal_depthMatchHot7; - std::vector *activeLength; - std::vector *activeLengthHot; - std::vector *trackDz; - std::vector *trackLayerCrossed; - std::vector *trackOuterHit; - std::vector *trackMissedInnerHits; - std::vector *trackMissedOuterHits; + std::vector *hcal_ieta; + std::vector *hcal_iphi; + std::vector *hcal_edepth1; + std::vector *hcal_activeL1; + std::vector *hcal_edepthHot1; + std::vector *hcal_activeHotL1; + std::vector *hcal_cdepthHot1; + std::vector *hcal_cdepthHotBG1; + std::vector *hcal_edepthCorrect1; + std::vector *hcal_edepthHotCorrect1; + std::vector *hcal_depthMatch1; + std::vector *hcal_depthMatchHot1; + std::vector *hcal_edepth2; + std::vector *hcal_activeL2; + std::vector *hcal_edepthHot2; + std::vector *hcal_activeHotL2; + std::vector *hcal_cdepthHot2; + std::vector *hcal_cdepthHotBG2; + std::vector *hcal_edepthCorrect2; + std::vector *hcal_edepthHotCorrect2; + std::vector *hcal_depthMatch2; + std::vector *hcal_depthMatchHot2; + std::vector *hcal_edepth3; + std::vector *hcal_activeL3; + std::vector *hcal_edepthHot3; + std::vector *hcal_activeHotL3; + std::vector *hcal_cdepthHot3; + std::vector *hcal_cdepthHotBG3; + std::vector *hcal_edepthCorrect3; + std::vector *hcal_edepthHotCorrect3; + std::vector *hcal_depthMatch3; + std::vector *hcal_depthMatchHot3; + std::vector *hcal_edepth4; + std::vector *hcal_activeL4; + std::vector *hcal_edepthHot4; + std::vector *hcal_activeHotL4; + std::vector *hcal_cdepthHot4; + std::vector *hcal_cdepthHotBG4; + std::vector *hcal_edepthCorrect4; + std::vector *hcal_edepthHotCorrect4; + std::vector *hcal_depthMatch4; + std::vector *hcal_depthMatchHot4; + std::vector *hcal_edepth5; + std::vector *hcal_activeL5; + std::vector *hcal_edepthHot5; + std::vector *hcal_activeHotL5; + std::vector *hcal_cdepthHot5; + std::vector *hcal_cdepthHotBG5; + std::vector *hcal_edepthCorrect5; + std::vector *hcal_edepthHotCorrect5; + std::vector *hcal_depthMatch5; + std::vector *hcal_depthMatchHot5; + std::vector *hcal_edepth6; + std::vector *hcal_activeL6; + std::vector *hcal_edepthHot6; + std::vector *hcal_activeHotL6; + std::vector *hcal_cdepthHot6; + std::vector *hcal_cdepthHotBG6; + std::vector *hcal_edepthCorrect6; + std::vector *hcal_edepthHotCorrect6; + std::vector *hcal_depthMatch6; + std::vector *hcal_depthMatchHot6; + std::vector *hcal_edepth7; + std::vector *hcal_activeL7; + std::vector *hcal_edepthHot7; + std::vector *hcal_activeHotL7; + std::vector *hcal_cdepthHot7; + std::vector *hcal_cdepthHotBG7; + std::vector *hcal_edepthCorrect7; + std::vector *hcal_edepthHotCorrect7; + std::vector *hcal_depthMatch7; + std::vector *hcal_depthMatchHot7; + std::vector *activeLength; + std::vector *activeLengthHot; + std::vector *trackDz; + std::vector *trackLayerCrossed; + std::vector *trackOuterHit; + std::vector *trackMissedInnerHits; + std::vector *trackMissedOuterHits; // List of branches - TBranch *b_pt_of_muon; - TBranch *b_eta_of_muon; - TBranch *b_phi_of_muon; - TBranch *b_energy_of_muon; - TBranch *b_p_of_muon; - TBranch *b_MediumMuon; - TBranch *b_IsolationR04; - TBranch *b_IsolationR03; - TBranch *b_ecal_3into3; - TBranch *b_hcal_3into3; - TBranch *b_tracker_3into3; - TBranch *b_emaxNearP; - TBranch *b_Run_No; - TBranch *b_Event_No; - TBranch *b_GoodVertex; - TBranch *b_matchedId; - TBranch *b_hcal_cellHot; - TBranch *b_ecal_3x3; - TBranch *b_hcal_1x1; - TBranch *b_ecal_detID; - TBranch *b_hcal_detID; - TBranch *b_ehcal_detID; - TBranch *b_hcal_ieta; - TBranch *b_hcal_iphi; - TBranch *b_hcal_edepth1; - TBranch *b_hcal_activeL1; - TBranch *b_hcal_edepthHot1; - TBranch *b_hcal_activeHotL1; - TBranch *b_hcal_cdepthHot1; - TBranch *b_hcal_cdepthHotBG1; - TBranch *b_hcal_edepthCorrect1; - TBranch *b_hcal_edepthHotCorrect1; - TBranch *b_hcal_depthMatch1; - TBranch *b_hcal_depthMatchHot1; - TBranch *b_hcal_edepth2; - TBranch *b_hcal_activeL2; - TBranch *b_hcal_edepthHot2; - TBranch *b_hcal_activeHotL2; - TBranch *b_hcal_cdepthHot2; - TBranch *b_hcal_cdepthHotBG2; - TBranch *b_hcal_edepthCorrect2; - TBranch *b_hcal_edepthHotCorrect2; - TBranch *b_hcal_depthMatch2; - TBranch *b_hcal_depthMatchHot2; - TBranch *b_hcal_edepth3; - TBranch *b_hcal_activeL3; - TBranch *b_hcal_edepthHot3; - TBranch *b_hcal_activeHotL3; - TBranch *b_hcal_cdepthHot3; - TBranch *b_hcal_cdepthHotBG3; - TBranch *b_hcal_edepthCorrect3; - TBranch *b_hcal_edepthHotCorrect3; - TBranch *b_hcal_depthMatch3; - TBranch *b_hcal_depthMatchHot3; - TBranch *b_hcal_edepth4; - TBranch *b_hcal_activeL4; - TBranch *b_hcal_edepthHot4; - TBranch *b_hcal_activeHotL4; - TBranch *b_hcal_cdepthHot4; - TBranch *b_hcal_cdepthHotBG4; - TBranch *b_hcal_edepthCorrect4; - TBranch *b_hcal_edepthHotCorrect4; - TBranch *b_hcal_depthMatch4; - TBranch *b_hcal_depthMatchHot4; - TBranch *b_hcal_edepth5; - TBranch *b_hcal_activeL5; - TBranch *b_hcal_edepthHot5; - TBranch *b_hcal_activeHotL5; - TBranch *b_hcal_cdepthHot5; - TBranch *b_hcal_cdepthHotBG5; - TBranch *b_hcal_edepthCorrect5; - TBranch *b_hcal_edepthHotCorrect5; - TBranch *b_hcal_depthMatch5; - TBranch *b_hcal_depthMatchHot5; - TBranch *b_hcal_edepth6; - TBranch *b_hcal_activeL6; - TBranch *b_hcal_edepthHot6; - TBranch *b_hcal_activeHotL6; - TBranch *b_hcal_cdepthHot6; - TBranch *b_hcal_cdepthHotBG6; - TBranch *b_hcal_edepthCorrect6; - TBranch *b_hcal_edepthHotCorrect6; - TBranch *b_hcal_depthMatch6; - TBranch *b_hcal_depthMatchHot6; - TBranch *b_hcal_edepth7; - TBranch *b_hcal_activeL7; - TBranch *b_hcal_edepthHot7; - TBranch *b_hcal_activeHotL7; - TBranch *b_hcal_cdepthHot7; - TBranch *b_hcal_cdepthHotBG7; - TBranch *b_hcal_edepthCorrect7; - TBranch *b_hcal_edepthHotCorrect7; - TBranch *b_hcal_depthMatch7; - TBranch *b_hcal_depthMatchHot7; - TBranch *b_activeLength; - TBranch *b_activeLengthHot; - TBranch *b_trackDz; - TBranch *b_trackLayerCrossed; - TBranch *b_trackOuterHit; - TBranch *b_trackMissedInnerHits; - TBranch *b_trackMissedOuterHits; + TBranch *b_pt_of_muon; + TBranch *b_eta_of_muon; + TBranch *b_phi_of_muon; + TBranch *b_energy_of_muon; + TBranch *b_p_of_muon; + TBranch *b_MediumMuon; + TBranch *b_IsolationR04; + TBranch *b_IsolationR03; + TBranch *b_ecal_3into3; + TBranch *b_hcal_3into3; + TBranch *b_tracker_3into3; + TBranch *b_emaxNearP; + TBranch *b_Run_No; + TBranch *b_Event_No; + TBranch *b_GoodVertex; + TBranch *b_matchedId; + TBranch *b_hcal_cellHot; + TBranch *b_ecal_3x3; + TBranch *b_hcal_1x1; + TBranch *b_ecal_detID; + TBranch *b_hcal_detID; + TBranch *b_ehcal_detID; + TBranch *b_hcal_ieta; + TBranch *b_hcal_iphi; + TBranch *b_hcal_edepth1; + TBranch *b_hcal_activeL1; + TBranch *b_hcal_edepthHot1; + TBranch *b_hcal_activeHotL1; + TBranch *b_hcal_cdepthHot1; + TBranch *b_hcal_cdepthHotBG1; + TBranch *b_hcal_edepthCorrect1; + TBranch *b_hcal_edepthHotCorrect1; + TBranch *b_hcal_depthMatch1; + TBranch *b_hcal_depthMatchHot1; + TBranch *b_hcal_edepth2; + TBranch *b_hcal_activeL2; + TBranch *b_hcal_edepthHot2; + TBranch *b_hcal_activeHotL2; + TBranch *b_hcal_cdepthHot2; + TBranch *b_hcal_cdepthHotBG2; + TBranch *b_hcal_edepthCorrect2; + TBranch *b_hcal_edepthHotCorrect2; + TBranch *b_hcal_depthMatch2; + TBranch *b_hcal_depthMatchHot2; + TBranch *b_hcal_edepth3; + TBranch *b_hcal_activeL3; + TBranch *b_hcal_edepthHot3; + TBranch *b_hcal_activeHotL3; + TBranch *b_hcal_cdepthHot3; + TBranch *b_hcal_cdepthHotBG3; + TBranch *b_hcal_edepthCorrect3; + TBranch *b_hcal_edepthHotCorrect3; + TBranch *b_hcal_depthMatch3; + TBranch *b_hcal_depthMatchHot3; + TBranch *b_hcal_edepth4; + TBranch *b_hcal_activeL4; + TBranch *b_hcal_edepthHot4; + TBranch *b_hcal_activeHotL4; + TBranch *b_hcal_cdepthHot4; + TBranch *b_hcal_cdepthHotBG4; + TBranch *b_hcal_edepthCorrect4; + TBranch *b_hcal_edepthHotCorrect4; + TBranch *b_hcal_depthMatch4; + TBranch *b_hcal_depthMatchHot4; + TBranch *b_hcal_edepth5; + TBranch *b_hcal_activeL5; + TBranch *b_hcal_edepthHot5; + TBranch *b_hcal_activeHotL5; + TBranch *b_hcal_cdepthHot5; + TBranch *b_hcal_cdepthHotBG5; + TBranch *b_hcal_edepthCorrect5; + TBranch *b_hcal_edepthHotCorrect5; + TBranch *b_hcal_depthMatch5; + TBranch *b_hcal_depthMatchHot5; + TBranch *b_hcal_edepth6; + TBranch *b_hcal_activeL6; + TBranch *b_hcal_edepthHot6; + TBranch *b_hcal_activeHotL6; + TBranch *b_hcal_cdepthHot6; + TBranch *b_hcal_cdepthHotBG6; + TBranch *b_hcal_edepthCorrect6; + TBranch *b_hcal_edepthHotCorrect6; + TBranch *b_hcal_depthMatch6; + TBranch *b_hcal_depthMatchHot6; + TBranch *b_hcal_edepth7; + TBranch *b_hcal_activeL7; + TBranch *b_hcal_edepthHot7; + TBranch *b_hcal_activeHotL7; + TBranch *b_hcal_cdepthHot7; + TBranch *b_hcal_cdepthHotBG7; + TBranch *b_hcal_edepthCorrect7; + TBranch *b_hcal_edepthHotCorrect7; + TBranch *b_hcal_depthMatch7; + TBranch *b_hcal_depthMatchHot7; + TBranch *b_activeLength; + TBranch *b_activeLengthHot; + TBranch *b_trackDz; + TBranch *b_trackLayerCrossed; + TBranch *b_trackOuterHit; + TBranch *b_trackMissedInnerHits; + TBranch *b_trackMissedOuterHits; - int modeLHC_; - bool debug_; - TFile *output_file; + int modeLHC_; + bool debug_; + TFile *output_file; }; -HBHEMuonHighEta::HBHEMuonHighEta(const char *infile, const char *outfile, - const int mode, const bool debug) { +HBHEMuonHighEta::HBHEMuonHighEta(const char *infile, const char *outfile, const int mode, const bool debug) { modeLHC_ = mode; debug_ = debug; - TFile *file = new TFile(infile); - TDirectory *dir = (TDirectory*)(file->FindObjectAny("hcalHBHEMuonHighEta")); - TTree *tree = (TTree*)(dir->FindObjectAny("HBHEMuonHighEta")); - std::cout << "Attaches tree HBHEMuonHighEta at " << tree << " in file " - << infile << std::endl; - + TFile *file = new TFile(infile); + TDirectory *dir = (TDirectory *)(file->FindObjectAny("hcalHBHEMuonHighEta")); + TTree *tree = (TTree *)(dir->FindObjectAny("HBHEMuonHighEta")); + std::cout << "Attaches tree HBHEMuonHighEta at " << tree << " in file " << infile << std::endl; + BookHistograms(outfile); Init(tree); } HBHEMuonHighEta::~HBHEMuonHighEta() { Close(); - if (!fChain) return; + if (!fChain) + return; delete fChain->GetCurrentFile(); } Int_t HBHEMuonHighEta::GetEntry(Long64_t entry) { // Read contents of entry. - if (!fChain) return 0; + if (!fChain) + return 0; return fChain->GetEntry(entry); } Long64_t HBHEMuonHighEta::LoadTree(Long64_t entry) { // Set the environment to read one entry - if (!fChain) return -5; + if (!fChain) + return -5; Long64_t centry = fChain->LoadTree(entry); - if (centry < 0) return centry; + if (centry < 0) + return centry; if (fChain->GetTreeNumber() != fCurrent) { fCurrent = fChain->GetTreeNumber(); Notify(); @@ -413,13 +413,14 @@ void HBHEMuonHighEta::Init(TTree *tree) { trackOuterHit = 0; trackMissedInnerHits = 0; trackMissedOuterHits = 0; - + // Set branch addresses and branch pointers - if (!tree) return; + if (!tree) + return; fChain = tree; fCurrent = -1; fChain->SetMakeClass(1); - + fChain->SetBranchAddress("pt_of_muon", &pt_of_muon, &b_pt_of_muon); fChain->SetBranchAddress("eta_of_muon", &eta_of_muon, &b_eta_of_muon); fChain->SetBranchAddress("phi_of_muon", &phi_of_muon, &b_phi_of_muon); @@ -530,18 +531,19 @@ Bool_t HBHEMuonHighEta::Notify() { // is started when using PROOF. It is normally not necessary to make changes // to the generated code, but the routine can be extended by the // user if needed. The return value is currently not used. - - return kTRUE; + + return kTRUE; } void HBHEMuonHighEta::Show(Long64_t entry) { // Print contents of entry. // If entry is not specified, print current entry - if (!fChain) return; + if (!fChain) + return; fChain->Show(entry); } -Int_t HBHEMuonHighEta::Cut(Long64_t ) { +Int_t HBHEMuonHighEta::Cut(Long64_t) { // This function may be called from Loop. // returns 1 if entry is accepted. // returns -1 otherwise. @@ -557,7 +559,7 @@ void HBHEMuonHighEta::Loop() { // root> t.Show(16); // Read and show values of entry 16 // root> t.Loop(); // Loop on all entries // - + // This is the loop skeleton where: // jentry is the global entry number in the chain // ientry is the entry number in the current Tree @@ -572,43 +574,44 @@ void HBHEMuonHighEta::Loop() { // METHOD2: replace line // fChain->GetEntry(jentry); //read all branches //by b_branchname->GetEntry(ientry); //read only this branch - if (fChain == 0) return; + if (fChain == 0) + return; Long64_t nentries = fChain->GetEntriesFast(); Long64_t nbytes = 0, nb = 0; - for (Long64_t jentry=0; jentryGetEntry(jentry); nbytes += nb; + if (ientry < 0) + break; + nb = fChain->GetEntry(jentry); + nbytes += nb; // if (Cut(ientry) < 0) continue; - } + } } -void HBHEMuonHighEta::BookHistograms(const char* fname) { - - output_file = TFile::Open(fname,"RECREATE"); -} +void HBHEMuonHighEta::BookHistograms(const char *fname) { output_file = TFile::Open(fname, "RECREATE"); } void HBHEMuonHighEta::Close() { output_file->cd(); - if (debug_) std::cout << "file yet to be Written" << std::endl; + if (debug_) + std::cout << "file yet to be Written" << std::endl; output_file->Write(); std::cout << "output file Written" << std::endl; output_file->Close(); - if (debug_) std::cout << "now doing return" << std::endl; + if (debug_) + std::cout << "now doing return" << std::endl; } - int HBHEMuonHighEta::nDepthBins(int eta, int phi) { // Run 1 scenario - int nDepthR1[29]={1,1,1,1,1,1,1,1,1,1,1,1,1,1,2,3,1,2,2,2,2,2,2,2,2,2,3,3,2}; + int nDepthR1[29] = {1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 3, 1, 2, 2, 2, 2, 2, 2, 2, 2, 2, 3, 3, 2}; // Run 2 scenario from 2018 - int nDepthR2[29]={1,1,1,1,1,1,1,1,1,1,1,1,1,1,2,4,3,5,6,6,6,6,6,6,6,7,7,7,3}; + int nDepthR2[29] = {1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 5, 6, 6, 6, 6, 6, 6, 6, 7, 7, 7, 3}; // Run 3 scenario - int nDepthR3[29]={4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,3,5,6,6,6,6,6,6,6,7,7,7,3}; + int nDepthR3[29] = {4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 3, 5, 6, 6, 6, 6, 6, 6, 6, 7, 7, 7, 3}; // Run 4 scenario - int nDepthR4[29]={4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,7,7,7,7,7,7,7,7,7,7,7,7,7}; + int nDepthR4[29] = {4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7}; // for 2021 scenario multi depth segmentation // int nDepth[29]={3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,5,5,5,5,5,5,5,5,5,5,5,5,5}; // modeLHC_ = 0 --> nbin defined maxDepthHB/HE @@ -617,9 +620,9 @@ int HBHEMuonHighEta::nDepthBins(int eta, int phi) { // = 3 --> corresponds to Run 3 (post LS2) // = 4 --> corresponds to 2017 (Plan 1) // = 5 --> corresponds to Run 4 (post LS3) - int nbin(0); + int nbin(0); if (modeLHC_ == 0) { - if (eta<=15) { + if (eta <= 15) { nbin = maxDepthHB_; } else if (eta == 16) { nbin = 4; @@ -627,28 +630,28 @@ int HBHEMuonHighEta::nDepthBins(int eta, int phi) { nbin = maxDepthHE_; } } else if (modeLHC_ == 1) { - nbin = nDepthR1[eta-1]; + nbin = nDepthR1[eta - 1]; } else if (modeLHC_ == 2) { - nbin = nDepthR2[eta-1]; + nbin = nDepthR2[eta - 1]; } else if (modeLHC_ == 3) { - nbin = nDepthR3[eta-1]; + nbin = nDepthR3[eta - 1]; } else if (modeLHC_ == 4) { if (phi > 0) { if (eta >= 16 && phi >= 63 && phi <= 66) { - nbin = nDepthR2[eta-1]; + nbin = nDepthR2[eta - 1]; } else { - nbin = nDepthR1[eta-1]; + nbin = nDepthR1[eta - 1]; } } else { if (eta >= 16) { - nbin = (nDepthR2[eta-1] > nDepthR1[eta-1]) ? nDepthR2[eta-1] : nDepthR1[eta-1]; + nbin = (nDepthR2[eta - 1] > nDepthR1[eta - 1]) ? nDepthR2[eta - 1] : nDepthR1[eta - 1]; } else { - nbin = nDepthR1[eta-1]; + nbin = nDepthR1[eta - 1]; } } } else { if (eta > 0 && eta < 30) { - nbin = nDepthR4[eta-1]; + nbin = nDepthR4[eta - 1]; } else { nbin = nDepthR4[28]; } diff --git a/Calibration/HcalCalibAlgos/macros/HBHEMuonOfflineAnalyzer.C b/Calibration/HcalCalibAlgos/macros/HBHEMuonOfflineAnalyzer.C index 7cf5515401f40..5f46b5b9e68a4 100644 --- a/Calibration/HcalCalibAlgos/macros/HBHEMuonOfflineAnalyzer.C +++ b/Calibration/HcalCalibAlgos/macros/HBHEMuonOfflineAnalyzer.C @@ -8,12 +8,12 @@ // // tree TTree* Pointer to the tree chain // infile const char* Name of the input file -// outfile const char* Name of the output file +// outfile const char* Name of the output file // (dyll_PU20_25_output_10.root) // rcorFile consr char* name of the text file having the correction factors // as a function of run numbers to be used for raddam // correction (default="", no corr.) -// flag int Flag of 2 digits ("to"): to where "o" decides if +// flag int Flag of 2 digits ("to"): to where "o" decides if // corrected (1) or default (0) energy to be used; // "t" decides if all depths to be merged (1) or not // (0) (default is 0) @@ -51,434 +51,461 @@ #include class HBHEMuonOfflineAnalyzer { - -public : - TChain *fChain; //!pointer to the analyzed TTree/TChain - Int_t fCurrent; //!current Tree number in a TChain +public: + TChain *fChain; //!pointer to the analyzed TTree/TChain + Int_t fCurrent; //!current Tree number in a TChain // Fixed size dimensions of array or collections stored in the TTree if any. // Declaration of leaf types - UInt_t Event_No; - UInt_t Run_No; - UInt_t LumiNumber; - UInt_t BXNumber; - UInt_t GoodVertex; - std::vector *PF_Muon; - std::vector *Global_Muon; - std::vector *Tracker_muon; - std::vector *MuonIsTight; - std::vector *MuonIsMedium; - std::vector *pt_of_muon; - std::vector *eta_of_muon; - std::vector *phi_of_muon; - std::vector *energy_of_muon; - std::vector *p_of_muon; - std::vector *muon_trkKink; - std::vector *muon_chi2LocalPosition; - std::vector *muon_segComp; - std::vector *TrackerLayer; - std::vector *NumPixelLayers; - std::vector *InnerTrackPixelHits; - std::vector *innerTrack; - std::vector *chiTracker; - std::vector *DxyTracker; - std::vector *DzTracker; - std::vector *innerTrackpt; - std::vector *innerTracketa; - std::vector *innerTrackphi; - std::vector *tight_validFraction; - std::vector *OuterTrack; - std::vector *OuterTrackPt; - std::vector *OuterTrackEta; - std::vector *OuterTrackPhi; - std::vector *OuterTrackChi; - std::vector *OuterTrackHits; - std::vector *OuterTrackRHits; - std::vector *GlobalTrack; - std::vector *GlobalTrckPt; - std::vector *GlobalTrckEta; - std::vector *GlobalTrckPhi; - std::vector *Global_Muon_Hits; - std::vector *MatchedStations; - std::vector *GlobTrack_Chi; - std::vector *Tight_LongitudinalImpactparameter; - std::vector *Tight_TransImpactparameter; - std::vector *IsolationR04; - std::vector *IsolationR03; - std::vector *ecal_3into3; - std::vector *hcal_3into3; - std::vector *tracker_3into3; - std::vector *matchedId; - std::vector *hcal_cellHot; - std::vector *ecal_3x3; - std::vector *hcal_1x1; + UInt_t Event_No; + UInt_t Run_No; + UInt_t LumiNumber; + UInt_t BXNumber; + UInt_t GoodVertex; + std::vector *PF_Muon; + std::vector *Global_Muon; + std::vector *Tracker_muon; + std::vector *MuonIsTight; + std::vector *MuonIsMedium; + std::vector *pt_of_muon; + std::vector *eta_of_muon; + std::vector *phi_of_muon; + std::vector *energy_of_muon; + std::vector *p_of_muon; + std::vector *muon_trkKink; + std::vector *muon_chi2LocalPosition; + std::vector *muon_segComp; + std::vector *TrackerLayer; + std::vector *NumPixelLayers; + std::vector *InnerTrackPixelHits; + std::vector *innerTrack; + std::vector *chiTracker; + std::vector *DxyTracker; + std::vector *DzTracker; + std::vector *innerTrackpt; + std::vector *innerTracketa; + std::vector *innerTrackphi; + std::vector *tight_validFraction; + std::vector *OuterTrack; + std::vector *OuterTrackPt; + std::vector *OuterTrackEta; + std::vector *OuterTrackPhi; + std::vector *OuterTrackChi; + std::vector *OuterTrackHits; + std::vector *OuterTrackRHits; + std::vector *GlobalTrack; + std::vector *GlobalTrckPt; + std::vector *GlobalTrckEta; + std::vector *GlobalTrckPhi; + std::vector *Global_Muon_Hits; + std::vector *MatchedStations; + std::vector *GlobTrack_Chi; + std::vector *Tight_LongitudinalImpactparameter; + std::vector *Tight_TransImpactparameter; + std::vector *IsolationR04; + std::vector *IsolationR03; + std::vector *ecal_3into3; + std::vector *hcal_3into3; + std::vector *tracker_3into3; + std::vector *matchedId; + std::vector *hcal_cellHot; + std::vector *ecal_3x3; + std::vector *hcal_1x1; std::vector *ecal_detID; std::vector *hcal_detID; std::vector *ehcal_detID; - std::vector *hcal_ieta; - std::vector *hcal_iphi; - std::vector *hcal_edepth1; - std::vector *hcal_activeL1; - std::vector *hcal_edepthHot1; - std::vector *hcal_activeHotL1; - std::vector *hcal_edepthCorrect1; - std::vector *hcal_edepthHotCorrect1; - std::vector *hcal_cdepthHot1; - std::vector *hcal_cdepthHotBG1; - std::vector *hcal_depthMatch1; - std::vector *hcal_depthMatchHot1; - std::vector *hcal_edepth2; - std::vector *hcal_activeL2; - std::vector *hcal_edepthHot2; - std::vector *hcal_activeHotL2; - std::vector *hcal_edepthCorrect2; - std::vector *hcal_edepthHotCorrect2; - std::vector *hcal_cdepthHot2; - std::vector *hcal_cdepthHotBG2; - std::vector *hcal_depthMatch2; - std::vector *hcal_depthMatchHot2; - std::vector *hcal_edepth3; - std::vector *hcal_activeL3; - std::vector *hcal_edepthHot3; - std::vector *hcal_activeHotL3; - std::vector *hcal_edepthCorrect3; - std::vector *hcal_edepthHotCorrect3; - std::vector *hcal_cdepthHot3; - std::vector *hcal_cdepthHotBG3; - std::vector *hcal_depthMatch3; - std::vector *hcal_depthMatchHot3; - std::vector *hcal_edepth4; - std::vector *hcal_activeL4; - std::vector *hcal_edepthHot4; - std::vector *hcal_activeHotL4; - std::vector *hcal_edepthCorrect4; - std::vector *hcal_edepthHotCorrect4; - std::vector *hcal_cdepthHot4; - std::vector *hcal_cdepthHotBG4; - std::vector *hcal_depthMatch4; - std::vector *hcal_depthMatchHot4; - std::vector *hcal_edepth5; - std::vector *hcal_activeL5; - std::vector *hcal_edepthHot5; - std::vector *hcal_activeHotL5; - std::vector *hcal_edepthCorrect5; - std::vector *hcal_edepthHotCorrect5; - std::vector *hcal_cdepthHot5; - std::vector *hcal_cdepthHotBG5; - std::vector *hcal_depthMatch5; - std::vector *hcal_depthMatchHot5; - std::vector *hcal_edepth6; - std::vector *hcal_activeL6; - std::vector *hcal_edepthHot6; - std::vector *hcal_activeHotL6; - std::vector *hcal_edepthCorrect6; - std::vector *hcal_edepthHotCorrect6; - std::vector *hcal_cdepthHot6; - std::vector *hcal_cdepthHotBG6; - std::vector *hcal_depthMatch6; - std::vector *hcal_depthMatchHot6; - std::vector *hcal_edepth7; - std::vector *hcal_activeL7; - std::vector *hcal_edepthHot7; - std::vector *hcal_activeHotL7; - std::vector *hcal_edepthCorrect7; - std::vector *hcal_edepthHotCorrect7; - std::vector *hcal_cdepthHot7; - std::vector *hcal_cdepthHotBG7; - std::vector *hcal_depthMatch7; - std::vector *hcal_depthMatchHot7; - std::vector *activeLength; - std::vector *activeLengthHot; - std::vector *hltresults; - std::vector *all_triggers; - + std::vector *hcal_ieta; + std::vector *hcal_iphi; + std::vector *hcal_edepth1; + std::vector *hcal_activeL1; + std::vector *hcal_edepthHot1; + std::vector *hcal_activeHotL1; + std::vector *hcal_edepthCorrect1; + std::vector *hcal_edepthHotCorrect1; + std::vector *hcal_cdepthHot1; + std::vector *hcal_cdepthHotBG1; + std::vector *hcal_depthMatch1; + std::vector *hcal_depthMatchHot1; + std::vector *hcal_edepth2; + std::vector *hcal_activeL2; + std::vector *hcal_edepthHot2; + std::vector *hcal_activeHotL2; + std::vector *hcal_edepthCorrect2; + std::vector *hcal_edepthHotCorrect2; + std::vector *hcal_cdepthHot2; + std::vector *hcal_cdepthHotBG2; + std::vector *hcal_depthMatch2; + std::vector *hcal_depthMatchHot2; + std::vector *hcal_edepth3; + std::vector *hcal_activeL3; + std::vector *hcal_edepthHot3; + std::vector *hcal_activeHotL3; + std::vector *hcal_edepthCorrect3; + std::vector *hcal_edepthHotCorrect3; + std::vector *hcal_cdepthHot3; + std::vector *hcal_cdepthHotBG3; + std::vector *hcal_depthMatch3; + std::vector *hcal_depthMatchHot3; + std::vector *hcal_edepth4; + std::vector *hcal_activeL4; + std::vector *hcal_edepthHot4; + std::vector *hcal_activeHotL4; + std::vector *hcal_edepthCorrect4; + std::vector *hcal_edepthHotCorrect4; + std::vector *hcal_cdepthHot4; + std::vector *hcal_cdepthHotBG4; + std::vector *hcal_depthMatch4; + std::vector *hcal_depthMatchHot4; + std::vector *hcal_edepth5; + std::vector *hcal_activeL5; + std::vector *hcal_edepthHot5; + std::vector *hcal_activeHotL5; + std::vector *hcal_edepthCorrect5; + std::vector *hcal_edepthHotCorrect5; + std::vector *hcal_cdepthHot5; + std::vector *hcal_cdepthHotBG5; + std::vector *hcal_depthMatch5; + std::vector *hcal_depthMatchHot5; + std::vector *hcal_edepth6; + std::vector *hcal_activeL6; + std::vector *hcal_edepthHot6; + std::vector *hcal_activeHotL6; + std::vector *hcal_edepthCorrect6; + std::vector *hcal_edepthHotCorrect6; + std::vector *hcal_cdepthHot6; + std::vector *hcal_cdepthHotBG6; + std::vector *hcal_depthMatch6; + std::vector *hcal_depthMatchHot6; + std::vector *hcal_edepth7; + std::vector *hcal_activeL7; + std::vector *hcal_edepthHot7; + std::vector *hcal_activeHotL7; + std::vector *hcal_edepthCorrect7; + std::vector *hcal_edepthHotCorrect7; + std::vector *hcal_cdepthHot7; + std::vector *hcal_cdepthHotBG7; + std::vector *hcal_depthMatch7; + std::vector *hcal_depthMatchHot7; + std::vector *activeLength; + std::vector *activeLengthHot; + std::vector *hltresults; + std::vector *all_triggers; + // List of branches - TBranch *b_Event_No; //! - TBranch *b_Run_No; //! - TBranch *b_LumiNumber; //! - TBranch *b_BXNumber; //! - TBranch *b_GoodVertex; //! - TBranch *b_PF_Muon; //! - TBranch *b_Global_Muon; //! - TBranch *b_Tracker_muon; //! - TBranch *b_MuonIsTight; //! - TBranch *b_MuonIsMedium; //! - TBranch *b_pt_of_muon; //! - TBranch *b_eta_of_muon; //! - TBranch *b_phi_of_muon; //! - TBranch *b_energy_of_muon; //! - TBranch *b_p_of_muon; //! - TBranch *b_muon_trkKink; //! - TBranch *b_muon_chi2LocalPosition; //! - TBranch *b_muon_segComp; //! - TBranch *b_TrackerLayer; //! - TBranch *b_NumPixelLayers; //! - TBranch *b_InnerTrackPixelHits; //! - TBranch *b_innerTrack; //! - TBranch *b_chiTracker; //! - TBranch *b_DxyTracker; //! - TBranch *b_DzTracker; //! - TBranch *b_innerTrackpt; //! - TBranch *b_innerTracketa; //! - TBranch *b_innerTrackphi; //! - TBranch *b_tight_validFraction; //! - TBranch *b_OuterTrack; //! - TBranch *b_OuterTrackPt; //! - TBranch *b_OuterTrackEta; //! - TBranch *b_OuterTrackPhi; //! - TBranch *b_OuterTrackChi; //! - TBranch *b_OuterTrackHits; //! - TBranch *b_OuterTrackRHits; //! - TBranch *b_GlobalTrack; //! - TBranch *b_GlobalTrckPt; //! - TBranch *b_GlobalTrckEta; //! - TBranch *b_GlobalTrckPhi; //! - TBranch *b_Global_Muon_Hits; //! - TBranch *b_MatchedStations; //! - TBranch *b_GlobTrack_Chi; //! - TBranch *b_Tight_LongitudinalImpactparameter; //! - TBranch *b_Tight_TransImpactparameter; //! - TBranch *b_IsolationR04; //! - TBranch *b_IsolationR03; //! - TBranch *b_ecal_3into3; //! - TBranch *b_hcal_3into3; //! - TBranch *b_tracker_3into3; //! - TBranch *b_matchedId; //! - TBranch *b_hcal_cellHot; //! - TBranch *b_ecal_3x3; //! - TBranch *b_hcal_1x1; //! - TBranch *b_hcal_detID; //! - TBranch *b_ecal_detID; //! - TBranch *b_ehcal_detID; //! - TBranch *b_hcal_ieta; //! - TBranch *b_hcal_iphi; //! - TBranch *b_hcal_edepth1; //! - TBranch *b_hcal_activeL1; //! - TBranch *b_hcal_edepthHot1; //! - TBranch *b_hcal_activeHotL1; //! - TBranch *b_hcal_edepthCorrect1; //! - TBranch *b_hcal_edepthHotCorrect1; //! - TBranch *b_hcal_cdepthHot1; //! - TBranch *b_hcal_cdepthHotBG1; //! - TBranch *b_hcal_depthMatch1; //! - TBranch *b_hcal_depthMatchHot1; //! - TBranch *b_hcal_edepth2; //! - TBranch *b_hcal_activeL2; //! - TBranch *b_hcal_edepthHot2; //! - TBranch *b_hcal_activeHotL2; //! - TBranch *b_hcal_edepthCorrect2; //! - TBranch *b_hcal_edepthHotCorrect2; //! - TBranch *b_hcal_cdepthHot2; //! - TBranch *b_hcal_cdepthHotBG2; //! - TBranch *b_hcal_depthMatch2; //! - TBranch *b_hcal_depthMatchHot2; //! - TBranch *b_hcal_edepth3; //! - TBranch *b_hcal_activeL3; //! - TBranch *b_hcal_edepthHot3; //! - TBranch *b_hcal_activeHotL3; //! - TBranch *b_hcal_edepthCorrect3; //! - TBranch *b_hcal_edepthHotCorrect3; //! - TBranch *b_hcal_cdepthHot3; //! - TBranch *b_hcal_cdepthHotBG3; //! - TBranch *b_hcal_depthMatch3; //! - TBranch *b_hcal_depthMatchHot3; //! - TBranch *b_hcal_edepth4; //! - TBranch *b_hcal_activeL4; //! - TBranch *b_hcal_edepthHot4; //! - TBranch *b_hcal_activeHotL4; //! - TBranch *b_hcal_edepthCorrect4; //! - TBranch *b_hcal_edepthHotCorrect4; //! - TBranch *b_hcal_cdepthHot4; //! - TBranch *b_hcal_cdepthHotBG4; //! - TBranch *b_hcal_depthMatch4; //! - TBranch *b_hcal_depthMatchHot4; //! - TBranch *b_hcal_edepth5; //! - TBranch *b_hcal_activeL5; //! - TBranch *b_hcal_edepthHot5; //! - TBranch *b_hcal_activeHotL5; //! - TBranch *b_hcal_edepthCorrect5; //! - TBranch *b_hcal_edepthHotCorrect5; //! - TBranch *b_hcal_cdepthHot5; //! - TBranch *b_hcal_cdepthHotBG5; //! - TBranch *b_hcal_depthMatch5; //! - TBranch *b_hcal_depthMatchHot5; //! - TBranch *b_hcal_edepth6; //! - TBranch *b_hcal_activeL6; //! - TBranch *b_hcal_edepthHot6; //! - TBranch *b_hcal_activeHotL6; //! - TBranch *b_hcal_edepthCorrect6; //! - TBranch *b_hcal_edepthHotCorrect6; //! - TBranch *b_hcal_cdepthHot6; //! - TBranch *b_hcal_cdepthHotBG6; //! - TBranch *b_hcal_depthMatch6; //! - TBranch *b_hcal_depthMatchHot6; //! - TBranch *b_hcal_edepth7; //! - TBranch *b_hcal_activeL7; //! - TBranch *b_hcal_edepthHot7; //! - TBranch *b_hcal_activeHotL7; //! - TBranch *b_hcal_edepthCorrect7; //! - TBranch *b_hcal_edepthHotCorrect7; //! - TBranch *b_hcal_cdepthHot7; //! - TBranch *b_hcal_cdepthHotBG7; //! - TBranch *b_hcal_depthMatch7; //! - TBranch *b_hcal_depthMatchHot7; //! - TBranch *b_activeLength; //! - TBranch *b_activeLengthHot; //! - TBranch *b_hltresults; //! - TBranch *b_all_triggers; //! - - - HBHEMuonOfflineAnalyzer(TChain *tree=0, - const char* outfile="dyll_PU20_25_output_10.root", - const char* rcorFileName="", int flag=0, int mode=2, - int maxDHB=4, int maxDHE=7, int runLo=1, - int runHi=99999999, int etaMin=1, int etaMax=29, - bool debug=false); - HBHEMuonOfflineAnalyzer(const char* infile, - const char* outfile="dyll_PU20_25_output_10.root", - const char* rcorFileName="", int flag=0, int mode=2, - int maxDHB=4, int maxDHE=7, int runLo=1, - int runHi=99999999, int etaMin=1, int etaMax=29, - bool debug=false); + TBranch *b_Event_No; //! + TBranch *b_Run_No; //! + TBranch *b_LumiNumber; //! + TBranch *b_BXNumber; //! + TBranch *b_GoodVertex; //! + TBranch *b_PF_Muon; //! + TBranch *b_Global_Muon; //! + TBranch *b_Tracker_muon; //! + TBranch *b_MuonIsTight; //! + TBranch *b_MuonIsMedium; //! + TBranch *b_pt_of_muon; //! + TBranch *b_eta_of_muon; //! + TBranch *b_phi_of_muon; //! + TBranch *b_energy_of_muon; //! + TBranch *b_p_of_muon; //! + TBranch *b_muon_trkKink; //! + TBranch *b_muon_chi2LocalPosition; //! + TBranch *b_muon_segComp; //! + TBranch *b_TrackerLayer; //! + TBranch *b_NumPixelLayers; //! + TBranch *b_InnerTrackPixelHits; //! + TBranch *b_innerTrack; //! + TBranch *b_chiTracker; //! + TBranch *b_DxyTracker; //! + TBranch *b_DzTracker; //! + TBranch *b_innerTrackpt; //! + TBranch *b_innerTracketa; //! + TBranch *b_innerTrackphi; //! + TBranch *b_tight_validFraction; //! + TBranch *b_OuterTrack; //! + TBranch *b_OuterTrackPt; //! + TBranch *b_OuterTrackEta; //! + TBranch *b_OuterTrackPhi; //! + TBranch *b_OuterTrackChi; //! + TBranch *b_OuterTrackHits; //! + TBranch *b_OuterTrackRHits; //! + TBranch *b_GlobalTrack; //! + TBranch *b_GlobalTrckPt; //! + TBranch *b_GlobalTrckEta; //! + TBranch *b_GlobalTrckPhi; //! + TBranch *b_Global_Muon_Hits; //! + TBranch *b_MatchedStations; //! + TBranch *b_GlobTrack_Chi; //! + TBranch *b_Tight_LongitudinalImpactparameter; //! + TBranch *b_Tight_TransImpactparameter; //! + TBranch *b_IsolationR04; //! + TBranch *b_IsolationR03; //! + TBranch *b_ecal_3into3; //! + TBranch *b_hcal_3into3; //! + TBranch *b_tracker_3into3; //! + TBranch *b_matchedId; //! + TBranch *b_hcal_cellHot; //! + TBranch *b_ecal_3x3; //! + TBranch *b_hcal_1x1; //! + TBranch *b_hcal_detID; //! + TBranch *b_ecal_detID; //! + TBranch *b_ehcal_detID; //! + TBranch *b_hcal_ieta; //! + TBranch *b_hcal_iphi; //! + TBranch *b_hcal_edepth1; //! + TBranch *b_hcal_activeL1; //! + TBranch *b_hcal_edepthHot1; //! + TBranch *b_hcal_activeHotL1; //! + TBranch *b_hcal_edepthCorrect1; //! + TBranch *b_hcal_edepthHotCorrect1; //! + TBranch *b_hcal_cdepthHot1; //! + TBranch *b_hcal_cdepthHotBG1; //! + TBranch *b_hcal_depthMatch1; //! + TBranch *b_hcal_depthMatchHot1; //! + TBranch *b_hcal_edepth2; //! + TBranch *b_hcal_activeL2; //! + TBranch *b_hcal_edepthHot2; //! + TBranch *b_hcal_activeHotL2; //! + TBranch *b_hcal_edepthCorrect2; //! + TBranch *b_hcal_edepthHotCorrect2; //! + TBranch *b_hcal_cdepthHot2; //! + TBranch *b_hcal_cdepthHotBG2; //! + TBranch *b_hcal_depthMatch2; //! + TBranch *b_hcal_depthMatchHot2; //! + TBranch *b_hcal_edepth3; //! + TBranch *b_hcal_activeL3; //! + TBranch *b_hcal_edepthHot3; //! + TBranch *b_hcal_activeHotL3; //! + TBranch *b_hcal_edepthCorrect3; //! + TBranch *b_hcal_edepthHotCorrect3; //! + TBranch *b_hcal_cdepthHot3; //! + TBranch *b_hcal_cdepthHotBG3; //! + TBranch *b_hcal_depthMatch3; //! + TBranch *b_hcal_depthMatchHot3; //! + TBranch *b_hcal_edepth4; //! + TBranch *b_hcal_activeL4; //! + TBranch *b_hcal_edepthHot4; //! + TBranch *b_hcal_activeHotL4; //! + TBranch *b_hcal_edepthCorrect4; //! + TBranch *b_hcal_edepthHotCorrect4; //! + TBranch *b_hcal_cdepthHot4; //! + TBranch *b_hcal_cdepthHotBG4; //! + TBranch *b_hcal_depthMatch4; //! + TBranch *b_hcal_depthMatchHot4; //! + TBranch *b_hcal_edepth5; //! + TBranch *b_hcal_activeL5; //! + TBranch *b_hcal_edepthHot5; //! + TBranch *b_hcal_activeHotL5; //! + TBranch *b_hcal_edepthCorrect5; //! + TBranch *b_hcal_edepthHotCorrect5; //! + TBranch *b_hcal_cdepthHot5; //! + TBranch *b_hcal_cdepthHotBG5; //! + TBranch *b_hcal_depthMatch5; //! + TBranch *b_hcal_depthMatchHot5; //! + TBranch *b_hcal_edepth6; //! + TBranch *b_hcal_activeL6; //! + TBranch *b_hcal_edepthHot6; //! + TBranch *b_hcal_activeHotL6; //! + TBranch *b_hcal_edepthCorrect6; //! + TBranch *b_hcal_edepthHotCorrect6; //! + TBranch *b_hcal_cdepthHot6; //! + TBranch *b_hcal_cdepthHotBG6; //! + TBranch *b_hcal_depthMatch6; //! + TBranch *b_hcal_depthMatchHot6; //! + TBranch *b_hcal_edepth7; //! + TBranch *b_hcal_activeL7; //! + TBranch *b_hcal_edepthHot7; //! + TBranch *b_hcal_activeHotL7; //! + TBranch *b_hcal_edepthCorrect7; //! + TBranch *b_hcal_edepthHotCorrect7; //! + TBranch *b_hcal_cdepthHot7; //! + TBranch *b_hcal_cdepthHotBG7; //! + TBranch *b_hcal_depthMatch7; //! + TBranch *b_hcal_depthMatchHot7; //! + TBranch *b_activeLength; //! + TBranch *b_activeLengthHot; //! + TBranch *b_hltresults; //! + TBranch *b_all_triggers; //! + + HBHEMuonOfflineAnalyzer(TChain *tree = 0, + const char *outfile = "dyll_PU20_25_output_10.root", + const char *rcorFileName = "", + int flag = 0, + int mode = 2, + int maxDHB = 4, + int maxDHE = 7, + int runLo = 1, + int runHi = 99999999, + int etaMin = 1, + int etaMax = 29, + bool debug = false); + HBHEMuonOfflineAnalyzer(const char *infile, + const char *outfile = "dyll_PU20_25_output_10.root", + const char *rcorFileName = "", + int flag = 0, + int mode = 2, + int maxDHB = 4, + int maxDHE = 7, + int runLo = 1, + int runHi = 99999999, + int etaMin = 1, + int etaMax = 29, + bool debug = false); // mode of LHC is kept 1 for 2017 scenario as no change in depth segmentation // mode of LHC is 0 for 2021 virtual ~HBHEMuonOfflineAnalyzer(); - virtual Int_t Cut(Long64_t entry); - virtual Int_t GetEntry(Long64_t entry); - virtual Long64_t LoadTree(Long64_t entry); - virtual void Init(TChain *tree, const char* rcorFileName, - int flag, int mode, int maxDHB, int maxDHE, - int runLo, int runHi, int etaMin, int etaMax); - virtual void Loop(); - virtual Bool_t Notify(); - virtual void Show(Long64_t entry = -1); - - bool fillChain(TChain *chain, const char* inputFileList); - bool readCorr(const char* rcorFileName); - void bookHistograms(const char* ); - bool getEnergy(unsigned int ml, int dep, double& enb, - double& enu, double& enh, double& enc, - double& chgS, double& chgB, double& actL); - void writeHistograms(); - bool looseMuon(unsigned int ml); - bool tightMuon(unsigned int ml); - bool softMuon(unsigned int ml); - bool mediumMuon2016(unsigned int ml); - void etaPhiHcal(unsigned int detId, int &eta, - int &phi, int &depth); - void etaPhiEcal(unsigned int detId, int& type, int& zside, - int& etaX, int& phiY, int& plane, - int& strip); - void calculateP(double pt, double eta, double& pM); - void close(); - int nDepthBins(int ieta, int iphi); - int nPhiBins(int ieta); - float getCorr(int run, unsigned int id); - std::vector splitString(const std::string&); - unsigned int getDetIdHBHE(int ieta, int iphi, int depth); - unsigned int getDetId(int subdet, int ieta, int iphi, int depth); - unsigned int correctDetId(const unsigned int& detId); - void unpackDetId(unsigned int detId, int& subdet, - int& zside, int& ieta, int& iphi, - int& depth); + virtual Int_t Cut(Long64_t entry); + virtual Int_t GetEntry(Long64_t entry); + virtual Long64_t LoadTree(Long64_t entry); + virtual void Init(TChain *tree, + const char *rcorFileName, + int flag, + int mode, + int maxDHB, + int maxDHE, + int runLo, + int runHi, + int etaMin, + int etaMax); + virtual void Loop(); + virtual Bool_t Notify(); + virtual void Show(Long64_t entry = -1); + + bool fillChain(TChain *chain, const char *inputFileList); + bool readCorr(const char *rcorFileName); + void bookHistograms(const char *); + bool getEnergy(unsigned int ml, + int dep, + double &enb, + double &enu, + double &enh, + double &enc, + double &chgS, + double &chgB, + double &actL); + void writeHistograms(); + bool looseMuon(unsigned int ml); + bool tightMuon(unsigned int ml); + bool softMuon(unsigned int ml); + bool mediumMuon2016(unsigned int ml); + void etaPhiHcal(unsigned int detId, int &eta, int &phi, int &depth); + void etaPhiEcal(unsigned int detId, int &type, int &zside, int &etaX, int &phiY, int &plane, int &strip); + void calculateP(double pt, double eta, double &pM); + void close(); + int nDepthBins(int ieta, int iphi); + int nPhiBins(int ieta); + float getCorr(int run, unsigned int id); + std::vector splitString(const std::string &); + unsigned int getDetIdHBHE(int ieta, int iphi, int depth); + unsigned int getDetId(int subdet, int ieta, int iphi, int depth); + unsigned int correctDetId(const unsigned int &detId); + void unpackDetId(unsigned int detId, int &subdet, int &zside, int &ieta, int &iphi, int &depth); + private: - static const int maxDep_=7; - static const int maxEta_=29; - static const int maxPhi_=72; + static const int maxDep_ = 7; + static const int maxEta_ = 29; + static const int maxPhi_ = 72; //3x16x72x2 + 5x4x72x2 + 5x9x36x2 - static const int maxHist_=20000;//13032; - static const int nCut_=1; - static const unsigned int nmax_=10; - const bool debug_; - int modeLHC_, maxDepthHB_, maxDepthHE_, maxDepth_; - int runLo_, runHi_, etaMin_, etaMax_; - bool cFactor_, useCorrect_, mergeDepth_; - int nHist, nDepths[maxEta_], nDepthsPhi[maxEta_]; - int indxEta[maxEta_][maxDep_][maxPhi_]; - TFile *output_file; - std::map corrFac_[nmax_]; - std::vector runlow_; - - TTree *outtree_; - int t_ieta, t_iphi, t_nvtx; - double t_p, t_ediff; - std::vector t_ene, t_enec, t_actln, t_charge; - std::vector t_depth; - - TH1D *h_evtype, *h_Pt_Muon[3], *h_Eta_Muon[3], *h_Phi_Muon[3], *h_P_Muon[3]; - TH1D *h_PF_Muon[3], *h_GlobTrack_Chi[3], *h_Global_Muon_Hits[3]; - TH1D *h_MatchedStations[3], *h_Tight_TransImpactparameter[3]; - TH1D *h_Tight_LongitudinalImpactparameter[3], *h_InnerTrackPixelHits[3]; - TH1D *h_TrackerLayer[3], *h_IsolationR04[3] , *h_Global_Muon[3]; - TH1D *h_LongImpactParameter[3], *h_LongImpactParameterBin1[3], *h_LongImpactParameterBin2[3]; - - TH1D *h_TransImpactParameter[3], *h_TransImpactParameterBin1[3], *h_TransImpactParameterBin2[3]; - TH1D *h_Hot_MuonEnergy_hcal_ClosestCell[3][maxHist_] , *h_Hot_MuonEnergy_hcal_HotCell[3][maxHist_] , *h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[3][maxHist_], *h_HotCell_MuonEnergy_phi[3][maxHist_], *h_active_length_Fill[3][maxHist_], *h_p_muon_ineta[3][maxHist_], *h_charge_signal[3][maxHist_], *h_charge_bg[3][maxHist_]; - TH2D *h_2D_Bin1[3], *h_2D_Bin2[3]; - TH1D *h_ecal_energy[3], *h_hcal_energy[3], *h_3x3_ecal[3], *h_1x1_hcal[3]; - TH1D *h_MuonHittingEcal[3], *h_HotCell[3], *h_MuonEnergy_hcal[3][maxHist_]; - TH1D *h_Hot_MuonEnergy_hcal[3][maxHist_]; - TH2D *hcal_ietaVsEnergy[3]; + static const int maxHist_ = 20000; //13032; + static const int nCut_ = 1; + static const unsigned int nmax_ = 10; + const bool debug_; + int modeLHC_, maxDepthHB_, maxDepthHE_, maxDepth_; + int runLo_, runHi_, etaMin_, etaMax_; + bool cFactor_, useCorrect_, mergeDepth_; + int nHist, nDepths[maxEta_], nDepthsPhi[maxEta_]; + int indxEta[maxEta_][maxDep_][maxPhi_]; + TFile *output_file; + std::map corrFac_[nmax_]; + std::vector runlow_; + + TTree *outtree_; + int t_ieta, t_iphi, t_nvtx; + double t_p, t_ediff; + std::vector t_ene, t_enec, t_actln, t_charge; + std::vector t_depth; + + TH1D *h_evtype, *h_Pt_Muon[3], *h_Eta_Muon[3], *h_Phi_Muon[3], *h_P_Muon[3]; + TH1D *h_PF_Muon[3], *h_GlobTrack_Chi[3], *h_Global_Muon_Hits[3]; + TH1D *h_MatchedStations[3], *h_Tight_TransImpactparameter[3]; + TH1D *h_Tight_LongitudinalImpactparameter[3], *h_InnerTrackPixelHits[3]; + TH1D *h_TrackerLayer[3], *h_IsolationR04[3], *h_Global_Muon[3]; + TH1D *h_LongImpactParameter[3], *h_LongImpactParameterBin1[3], *h_LongImpactParameterBin2[3]; + + TH1D *h_TransImpactParameter[3], *h_TransImpactParameterBin1[3], *h_TransImpactParameterBin2[3]; + TH1D *h_Hot_MuonEnergy_hcal_ClosestCell[3][maxHist_], *h_Hot_MuonEnergy_hcal_HotCell[3][maxHist_], + *h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[3][maxHist_], *h_HotCell_MuonEnergy_phi[3][maxHist_], + *h_active_length_Fill[3][maxHist_], *h_p_muon_ineta[3][maxHist_], *h_charge_signal[3][maxHist_], + *h_charge_bg[3][maxHist_]; + TH2D *h_2D_Bin1[3], *h_2D_Bin2[3]; + TH1D *h_ecal_energy[3], *h_hcal_energy[3], *h_3x3_ecal[3], *h_1x1_hcal[3]; + TH1D *h_MuonHittingEcal[3], *h_HotCell[3], *h_MuonEnergy_hcal[3][maxHist_]; + TH1D *h_Hot_MuonEnergy_hcal[3][maxHist_]; + TH2D *hcal_ietaVsEnergy[3]; TProfile *h_EtaX_hcal[3], *h_PhiY_hcal[3], *h_EtaX_ecal[3], *h_PhiY_ecal[3]; TProfile *h_Eta_ecal[3], *h_Phi_ecal[3]; TProfile *h_MuonEnergy_eta[3][maxDep_], *h_MuonEnergy_phi[3][maxDep_], *h_MuonEnergy_muon_eta[3][maxDep_]; - TProfile *h_Hot_MuonEnergy_eta[3][maxDep_], *h_Hot_MuonEnergy_phi[3][maxDep_], *h_Hot_MuonEnergy_muon_eta[3][maxDep_]; - TProfile *h_IsoHot_MuonEnergy_eta[3][maxDep_], *h_IsoHot_MuonEnergy_phi[3][maxDep_], *h_IsoHot_MuonEnergy_muon_eta[3][maxDep_]; - TProfile *h_IsoWithoutHot_MuonEnergy_eta[3][maxDep_], *h_IsoWithoutHot_MuonEnergy_phi[3][maxDep_], *h_IsoWithoutHot_MuonEnergy_muon_eta[3][maxDep_]; - TProfile *h_HotWithoutIso_MuonEnergy_eta[3][maxDep_], *h_HotWithoutIso_MuonEnergy_phi[3][maxDep_], *h_HotWithoutIso_MuonEnergy_muon_eta[3][maxDep_]; - + TProfile *h_Hot_MuonEnergy_eta[3][maxDep_], *h_Hot_MuonEnergy_phi[3][maxDep_], *h_Hot_MuonEnergy_muon_eta[3][maxDep_]; + TProfile *h_IsoHot_MuonEnergy_eta[3][maxDep_], *h_IsoHot_MuonEnergy_phi[3][maxDep_], + *h_IsoHot_MuonEnergy_muon_eta[3][maxDep_]; + TProfile *h_IsoWithoutHot_MuonEnergy_eta[3][maxDep_], *h_IsoWithoutHot_MuonEnergy_phi[3][maxDep_], + *h_IsoWithoutHot_MuonEnergy_muon_eta[3][maxDep_]; + TProfile *h_HotWithoutIso_MuonEnergy_eta[3][maxDep_], *h_HotWithoutIso_MuonEnergy_phi[3][maxDep_], + *h_HotWithoutIso_MuonEnergy_muon_eta[3][maxDep_]; }; - HBHEMuonOfflineAnalyzer::HBHEMuonOfflineAnalyzer(TChain *tree, - const char* outFileName, - const char* rcorFileName, - int flag, int mode, - int maxDHB, int maxDHE, - int runLo, int runHi, - int etaMin, int etaMax, - bool deb) : debug_(deb), - cFactor_(false) { - - Init(tree, rcorFileName, flag, mode, maxDHB, maxDHE, runLo, runHi, - etaMin, etaMax); - + const char *outFileName, + const char *rcorFileName, + int flag, + int mode, + int maxDHB, + int maxDHE, + int runLo, + int runHi, + int etaMin, + int etaMax, + bool deb) + : debug_(deb), cFactor_(false) { + Init(tree, rcorFileName, flag, mode, maxDHB, maxDHE, runLo, runHi, etaMin, etaMax); + //Now book histograms bookHistograms(outFileName); } -HBHEMuonOfflineAnalyzer::HBHEMuonOfflineAnalyzer(const char* infile, - const char* outFileName, - const char* rcorFileName, - int flag, int mode, - int maxDHB, int maxDHE, - int runLo, int runHi, - int etaMin, int etaMax, - bool deb) : debug_(deb), - cFactor_(false) { - +HBHEMuonOfflineAnalyzer::HBHEMuonOfflineAnalyzer(const char *infile, + const char *outFileName, + const char *rcorFileName, + int flag, + int mode, + int maxDHB, + int maxDHE, + int runLo, + int runHi, + int etaMin, + int etaMax, + bool deb) + : debug_(deb), cFactor_(false) { TChain *chain = new TChain("hcalHBHEMuon/TREE"); - if (!fillChain(chain,infile)) { + if (!fillChain(chain, infile)) { std::cout << "*****No valid tree chain can be obtained*****" << std::endl; } else { - std::cout << "Proceed with a tree chain with " << chain->GetEntries() - << " entries" << std::endl; - Init(chain, rcorFileName, flag, mode, maxDHB, maxDHE, runLo, runHi, - etaMin, etaMax); - + std::cout << "Proceed with a tree chain with " << chain->GetEntries() << " entries" << std::endl; + Init(chain, rcorFileName, flag, mode, maxDHB, maxDHE, runLo, runHi, etaMin, etaMax); + //Now book histograms bookHistograms(outFileName); } } HBHEMuonOfflineAnalyzer::~HBHEMuonOfflineAnalyzer() { - if (!fChain) return; + if (!fChain) + return; delete fChain->GetCurrentFile(); } -Int_t HBHEMuonOfflineAnalyzer::Cut(Long64_t ) { +Int_t HBHEMuonOfflineAnalyzer::Cut(Long64_t) { // This function may be called from Loop. // returns 1 if entry is accepted. // returns -1 otherwise. @@ -487,15 +514,18 @@ Int_t HBHEMuonOfflineAnalyzer::Cut(Long64_t ) { Int_t HBHEMuonOfflineAnalyzer::GetEntry(Long64_t entry) { // Read contents of entry. - if (!fChain) return 0; + if (!fChain) + return 0; return fChain->GetEntry(entry); } Long64_t HBHEMuonOfflineAnalyzer::LoadTree(Long64_t entry) { // Set the environment to read one entry - if (!fChain) return -5; + if (!fChain) + return -5; Long64_t centry = fChain->LoadTree(entry); - if (centry < 0) return centry; + if (centry < 0) + return centry; if (fChain->GetTreeNumber() != fCurrent) { fCurrent = fChain->GetTreeNumber(); Notify(); @@ -503,27 +533,35 @@ Long64_t HBHEMuonOfflineAnalyzer::LoadTree(Long64_t entry) { return centry; } -void HBHEMuonOfflineAnalyzer::Init(TChain *tree, const char* rcorFileName, - int flag, int mode, int maxDHB, int maxDHE, - int runLo, int runHi, int etaMin, - int etaMax) { - - modeLHC_ = mode; +void HBHEMuonOfflineAnalyzer::Init(TChain *tree, + const char *rcorFileName, + int flag, + int mode, + int maxDHB, + int maxDHE, + int runLo, + int runHi, + int etaMin, + int etaMax) { + modeLHC_ = mode; maxDepthHB_ = maxDHB; maxDepthHE_ = maxDHE; - maxDepth_ = (maxDepthHB_ > maxDepthHE_) ? maxDepthHB_ : maxDepthHE_; - runLo_ = runLo; - runHi_ = runHi; - etaMin_ = (etaMin > 0) ? etaMin : 1; - etaMax_ = (etaMax <= 29) ? etaMax : 29; + maxDepth_ = (maxDepthHB_ > maxDepthHE_) ? maxDepthHB_ : maxDepthHE_; + runLo_ = runLo; + runHi_ = runHi; + etaMin_ = (etaMin > 0) ? etaMin : 1; + etaMax_ = (etaMax <= 29) ? etaMax : 29; if (etaMax_ <= etaMin_) { - if (etaMax_ == 29) etaMin_ = etaMax_ - 1; - else etaMax_ = etaMin_ + 1; + if (etaMax_ == 29) + etaMin_ = etaMax_ - 1; + else + etaMax_ = etaMin_ + 1; } - useCorrect_ = ((flag%10) > 0); - mergeDepth_ = (((flag/10)%10) > 0); - if (std::string(rcorFileName) != "") cFactor_ = readCorr(rcorFileName); - + useCorrect_ = ((flag % 10) > 0); + mergeDepth_ = (((flag / 10) % 10) > 0); + if (std::string(rcorFileName) != "") + cFactor_ = readCorr(rcorFileName); + // The Init() function is called when the selector needs to initialize // a new tree or chain. Typically here the branch addresses and branch // pointers of the tree will be set. @@ -531,135 +569,136 @@ void HBHEMuonOfflineAnalyzer::Init(TChain *tree, const char* rcorFileName, // code, but the routine can be extended by the user if needed. // Init() will be called many times when running on PROOF // (once per file to be processed). - + // Set object pointer - PF_Muon = 0; - Global_Muon = 0; - Tracker_muon = 0; - pt_of_muon = 0; - eta_of_muon = 0; - phi_of_muon = 0; - energy_of_muon = 0; - p_of_muon = 0; - muon_trkKink = 0; - muon_chi2LocalPosition = 0; - muon_segComp = 0; - TrackerLayer = 0; - NumPixelLayers = 0; - InnerTrackPixelHits = 0; - innerTrack = 0; - chiTracker = 0; - DxyTracker = 0; - DzTracker = 0; - innerTrackpt = 0; - innerTracketa = 0; - innerTrackphi = 0; - tight_validFraction = 0; - OuterTrack = 0; - OuterTrackPt = 0; - OuterTrackEta = 0; - OuterTrackPhi = 0; - OuterTrackHits = 0; - OuterTrackRHits = 0; - OuterTrackChi = 0; - GlobalTrack = 0; - GlobalTrckPt = 0; - GlobalTrckEta = 0; - GlobalTrckPhi = 0; - Global_Muon_Hits = 0; - MatchedStations = 0; - GlobTrack_Chi = 0; + PF_Muon = 0; + Global_Muon = 0; + Tracker_muon = 0; + pt_of_muon = 0; + eta_of_muon = 0; + phi_of_muon = 0; + energy_of_muon = 0; + p_of_muon = 0; + muon_trkKink = 0; + muon_chi2LocalPosition = 0; + muon_segComp = 0; + TrackerLayer = 0; + NumPixelLayers = 0; + InnerTrackPixelHits = 0; + innerTrack = 0; + chiTracker = 0; + DxyTracker = 0; + DzTracker = 0; + innerTrackpt = 0; + innerTracketa = 0; + innerTrackphi = 0; + tight_validFraction = 0; + OuterTrack = 0; + OuterTrackPt = 0; + OuterTrackEta = 0; + OuterTrackPhi = 0; + OuterTrackHits = 0; + OuterTrackRHits = 0; + OuterTrackChi = 0; + GlobalTrack = 0; + GlobalTrckPt = 0; + GlobalTrckEta = 0; + GlobalTrckPhi = 0; + Global_Muon_Hits = 0; + MatchedStations = 0; + GlobTrack_Chi = 0; Tight_LongitudinalImpactparameter = 0; - Tight_TransImpactparameter = 0; - IsolationR04 = 0; - IsolationR03 = 0; - ecal_3into3 = 0; - hcal_3into3 = 0; - tracker_3into3 = 0; - matchedId = 0; - hcal_cellHot = 0; - ecal_3x3 = 0; - hcal_1x1 = 0; - ecal_detID = 0; - hcal_detID = 0; - ehcal_detID = 0; - hcal_edepth1 = 0; - hcal_activeL1 = 0; - hcal_edepthHot1 = 0; - hcal_activeHotL1 = 0; - hcal_edepthCorrect1 = 0; - hcal_edepthHotCorrect1 = 0; - hcal_cdepthHot1 = 0; - hcal_cdepthHotBG1 = 0; - hcal_depthMatch1 = 0; - hcal_depthMatchHot1 = 0; - hcal_edepth2 = 0; - hcal_activeL2 = 0; - hcal_edepthHot2 = 0; - hcal_activeHotL2 = 0; - hcal_edepthCorrect2 = 0; - hcal_edepthHotCorrect2 = 0; - hcal_cdepthHot2 = 0; - hcal_cdepthHotBG2 = 0; - hcal_depthMatch2 = 0; - hcal_depthMatchHot2 = 0; - hcal_edepth3 = 0; - hcal_activeL3 = 0; - hcal_edepthHot3 = 0; - hcal_activeHotL3 = 0; - hcal_edepthCorrect3 = 0; - hcal_edepthHotCorrect3 = 0; - hcal_cdepthHot3 = 0; - hcal_cdepthHotBG3 = 0; - hcal_depthMatch3 = 0; - hcal_depthMatchHot3 = 0; - hcal_edepth4 = 0; - hcal_activeL4 = 0; - hcal_edepthHot4 = 0; - hcal_activeHotL4 = 0; - hcal_edepthCorrect4 = 0; - hcal_edepthHotCorrect4 = 0; - hcal_cdepthHot4 = 0; - hcal_cdepthHotBG4 = 0; - hcal_depthMatch4 = 0; - hcal_depthMatchHot4 = 0; - hcal_edepth5 = 0; - hcal_activeL5 = 0; - hcal_edepthHot5 = 0; - hcal_activeHotL5 = 0; - hcal_edepthCorrect5 = 0; - hcal_edepthHotCorrect5 = 0; - hcal_cdepthHot5 = 0; - hcal_cdepthHotBG5 = 0; - hcal_depthMatch5 = 0; - hcal_depthMatchHot5 = 0; - hcal_edepth6 = 0; - hcal_activeL6 = 0; - hcal_edepthHot6 = 0; - hcal_activeHotL6 = 0; - hcal_edepthCorrect6 = 0; - hcal_edepthHotCorrect6 = 0; - hcal_cdepthHot6 = 0; - hcal_cdepthHotBG6 = 0; - hcal_depthMatch6 = 0; - hcal_depthMatchHot6 = 0; - hcal_edepth7 = 0; - hcal_activeL7 = 0; - hcal_edepthHot7 = 0; - hcal_activeHotL7 = 0; - hcal_edepthCorrect7 = 0; - hcal_edepthHotCorrect7 = 0; - hcal_cdepthHot7 = 0; - hcal_cdepthHotBG7 = 0; - hcal_depthMatch7 = 0; - hcal_depthMatchHot7 = 0; - activeLength = 0; - activeLengthHot = 0; - hltresults = 0; - all_triggers = 0; + Tight_TransImpactparameter = 0; + IsolationR04 = 0; + IsolationR03 = 0; + ecal_3into3 = 0; + hcal_3into3 = 0; + tracker_3into3 = 0; + matchedId = 0; + hcal_cellHot = 0; + ecal_3x3 = 0; + hcal_1x1 = 0; + ecal_detID = 0; + hcal_detID = 0; + ehcal_detID = 0; + hcal_edepth1 = 0; + hcal_activeL1 = 0; + hcal_edepthHot1 = 0; + hcal_activeHotL1 = 0; + hcal_edepthCorrect1 = 0; + hcal_edepthHotCorrect1 = 0; + hcal_cdepthHot1 = 0; + hcal_cdepthHotBG1 = 0; + hcal_depthMatch1 = 0; + hcal_depthMatchHot1 = 0; + hcal_edepth2 = 0; + hcal_activeL2 = 0; + hcal_edepthHot2 = 0; + hcal_activeHotL2 = 0; + hcal_edepthCorrect2 = 0; + hcal_edepthHotCorrect2 = 0; + hcal_cdepthHot2 = 0; + hcal_cdepthHotBG2 = 0; + hcal_depthMatch2 = 0; + hcal_depthMatchHot2 = 0; + hcal_edepth3 = 0; + hcal_activeL3 = 0; + hcal_edepthHot3 = 0; + hcal_activeHotL3 = 0; + hcal_edepthCorrect3 = 0; + hcal_edepthHotCorrect3 = 0; + hcal_cdepthHot3 = 0; + hcal_cdepthHotBG3 = 0; + hcal_depthMatch3 = 0; + hcal_depthMatchHot3 = 0; + hcal_edepth4 = 0; + hcal_activeL4 = 0; + hcal_edepthHot4 = 0; + hcal_activeHotL4 = 0; + hcal_edepthCorrect4 = 0; + hcal_edepthHotCorrect4 = 0; + hcal_cdepthHot4 = 0; + hcal_cdepthHotBG4 = 0; + hcal_depthMatch4 = 0; + hcal_depthMatchHot4 = 0; + hcal_edepth5 = 0; + hcal_activeL5 = 0; + hcal_edepthHot5 = 0; + hcal_activeHotL5 = 0; + hcal_edepthCorrect5 = 0; + hcal_edepthHotCorrect5 = 0; + hcal_cdepthHot5 = 0; + hcal_cdepthHotBG5 = 0; + hcal_depthMatch5 = 0; + hcal_depthMatchHot5 = 0; + hcal_edepth6 = 0; + hcal_activeL6 = 0; + hcal_edepthHot6 = 0; + hcal_activeHotL6 = 0; + hcal_edepthCorrect6 = 0; + hcal_edepthHotCorrect6 = 0; + hcal_cdepthHot6 = 0; + hcal_cdepthHotBG6 = 0; + hcal_depthMatch6 = 0; + hcal_depthMatchHot6 = 0; + hcal_edepth7 = 0; + hcal_activeL7 = 0; + hcal_edepthHot7 = 0; + hcal_activeHotL7 = 0; + hcal_edepthCorrect7 = 0; + hcal_edepthHotCorrect7 = 0; + hcal_cdepthHot7 = 0; + hcal_cdepthHotBG7 = 0; + hcal_depthMatch7 = 0; + hcal_depthMatchHot7 = 0; + activeLength = 0; + activeLengthHot = 0; + hltresults = 0; + all_triggers = 0; // Set branch addresses and branch pointers - if (!tree) return; + if (!tree) + return; fChain = tree; fCurrent = -1; fChain->SetMakeClass(1); @@ -707,7 +746,8 @@ void HBHEMuonOfflineAnalyzer::Init(TChain *tree, const char* rcorFileName, fChain->SetBranchAddress("Global_Muon_Hits", &Global_Muon_Hits, &b_Global_Muon_Hits); fChain->SetBranchAddress("MatchedStations", &MatchedStations, &b_MatchedStations); fChain->SetBranchAddress("GlobTrack_Chi", &GlobTrack_Chi, &b_GlobTrack_Chi); - fChain->SetBranchAddress("Tight_LongitudinalImpactparameter", &Tight_LongitudinalImpactparameter, &b_Tight_LongitudinalImpactparameter); + fChain->SetBranchAddress( + "Tight_LongitudinalImpactparameter", &Tight_LongitudinalImpactparameter, &b_Tight_LongitudinalImpactparameter); fChain->SetBranchAddress("Tight_TransImpactparameter", &Tight_TransImpactparameter, &b_Tight_TransImpactparameter); fChain->SetBranchAddress("IsolationR04", &IsolationR04, &b_IsolationR04); fChain->SetBranchAddress("IsolationR03", &IsolationR03, &b_IsolationR03); @@ -797,244 +837,263 @@ void HBHEMuonOfflineAnalyzer::Init(TChain *tree, const char* rcorFileName, fChain->SetBranchAddress("activeLengthHot", &activeLengthHot, &b_activeLengthHot); fChain->SetBranchAddress("hltresults", &hltresults, &b_hltresults); fChain->SetBranchAddress("all_triggers", &all_triggers, &b_all_triggers); - + Notify(); } void HBHEMuonOfflineAnalyzer::Loop() { - //declarations - if (fChain == 0) return; - + if (fChain == 0) + return; + Long64_t nentries = fChain->GetEntriesFast(); - - if (debug_) std::cout << "nevent = " << nentries << std::endl; - + + if (debug_) + std::cout << "nevent = " << nentries << std::endl; + Long64_t nbytes = 0, nb = 0; - - for (Long64_t jentry=0; jentryGetEntry(jentry); nbytes += nb; - if ((int)(Run_No) < runLo_ || (int)(Run_No) > runHi_) continue; - if (debug_) std::cout << "Run " << Run_No << " Event " << Event_No << " Muons " << pt_of_muon->size() << std::endl; + if (ientry < 0) + break; + nb = fChain->GetEntry(jentry); + nbytes += nb; + if ((int)(Run_No) < runLo_ || (int)(Run_No) > runHi_) + continue; + if (debug_) + std::cout << "Run " << Run_No << " Event " << Event_No << " Muons " << pt_of_muon->size() << std::endl; bool loose(false), soft(false), tight(false), pcut(false), ptcut(false); - for (unsigned int ml = 0; ml< pt_of_muon->size(); ml++) { - - t_ene.clear(); t_enec.clear(); t_charge.clear(); t_actln.clear(); + for (unsigned int ml = 0; ml < pt_of_muon->size(); ml++) { + t_ene.clear(); + t_enec.clear(); + t_charge.clear(); + t_actln.clear(); t_depth.clear(); - - if(debug_) std::cout << "ecal_det_id " << ecal_detID->at(ml) << std::endl; - + + if (debug_) + std::cout << "ecal_det_id " << ecal_detID->at(ml) << std::endl; + int typeEcal, etaXEcal, phiYEcal, zsideEcal, planeEcal, stripEcal; - etaPhiEcal(ecal_detID->at(ml),typeEcal,zsideEcal,etaXEcal,phiYEcal,planeEcal,stripEcal); - double etaEcal = (etaXEcal-0.5)*zsideEcal; - double phiEcal = phiYEcal-0.5; - - if (debug_) std::cout << "hcal_det_id " << std::hex << hcal_detID->at(ml) - << std::dec; - - int etaHcal, phiHcal, depthHcal; - etaPhiHcal(hcal_detID->at(ml),etaHcal,phiHcal,depthHcal); - - int eta = (etaHcal > 0) ? etaHcal-1 : -(1+etaHcal); - double etaXHcal = (etaHcal > 0) ? etaHcal-0.5 : etaHcal+0.5; - int nDepth = nDepthBins(eta+1,phiHcal); - int nPhi = nPhiBins(eta+1); - int PHI = (nPhi > 36) ? (phiHcal-1) : (phiHcal-1)/2; - double phiYHcal = (phiHcal-0.5); - t_ieta = etaHcal; - t_iphi = PHI; - t_p = p_of_muon->at(ml); - t_ediff = hcal_3into3->at(ml) - hcal_1x1->at(ml); - t_nvtx = GoodVertex; - if (p_of_muon->at(ml) > 20) pcut = true; - if (pt_of_muon->at(ml) > 20) ptcut = true; - if (looseMuon(ml)) loose = true; - if (softMuon(ml)) soft = true; - if (tightMuon(ml)) tight = true; - - if (debug_) - std::cout << " etaHcal " << etaHcal << ":" << etaXHcal << " phiHcal " - << phiHcal << ":" << phiYHcal << ":" << PHI << " Depth " - << nDepth << " Muon Pt " << pt_of_muon->at(ml) << " Isol " - << IsolationR04->at(ml) << std::endl; - - for (int cut=0; cut= etaMin_) && ((eta+1) <= etaMax_)) { - // h_P_Muon[cut]->Fill(p_of_muon->at(ml)); - h_Pt_Muon[cut]->Fill(pt_of_muon->at(ml)); - h_Eta_Muon[cut]->Fill(eta_of_muon->at(ml)); - h_Phi_Muon[cut]->Fill(phi_of_muon->at(ml)); - h_PF_Muon[cut]->Fill(PF_Muon->at(ml)); - h_GlobTrack_Chi[cut]->Fill(GlobTrack_Chi->at(ml)); - h_Global_Muon_Hits[cut]->Fill(Global_Muon_Hits->at(ml)); - h_MatchedStations[cut]->Fill(MatchedStations->at(ml)); - h_Tight_TransImpactparameter[cut]->Fill(Tight_TransImpactparameter->at(ml)); - h_Tight_LongitudinalImpactparameter[cut]->Fill(Tight_LongitudinalImpactparameter->at(ml)); - h_InnerTrackPixelHits[cut]->Fill(InnerTrackPixelHits->at(ml)); - h_TrackerLayer[cut]->Fill(TrackerLayer->at(ml)); - h_IsolationR04[cut]->Fill(IsolationR04->at(ml)); - h_Global_Muon[cut]->Fill(Global_Muon->at(ml)); - - h_TransImpactParameter[cut]->Fill(Tight_TransImpactparameter->at(ml)); - h_LongImpactParameter[cut]->Fill(Tight_LongitudinalImpactparameter->at(ml)); - - //in Phi Bins - if(((phi_of_muon->at(ml)) >= -1.5) || ((phi_of_muon->at(ml)) <= 0.5)) { - h_TransImpactParameterBin1[cut]->Fill(Tight_TransImpactparameter->at(ml)); - h_LongImpactParameterBin1[cut]->Fill(Tight_LongitudinalImpactparameter->at(ml)); - h_2D_Bin1[cut]->Fill(Tight_TransImpactparameter->at(ml),Tight_LongitudinalImpactparameter->at(ml)); - } - - if((phi_of_muon->at(ml) > 0.5) || (phi_of_muon->at(ml) < -1.5)) { - h_TransImpactParameterBin2[cut]->Fill(Tight_TransImpactparameter->at(ml)); - h_LongImpactParameterBin2[cut]->Fill(Tight_LongitudinalImpactparameter->at(ml)); - h_2D_Bin2[cut]->Fill(Tight_TransImpactparameter->at(ml),Tight_LongitudinalImpactparameter->at(ml)); - } - - h_ecal_energy[cut]->Fill(ecal_3into3->at(ml)); - h_3x3_ecal[cut]->Fill(ecal_3x3->at(ml)); - h_Eta_ecal[cut]->Fill(eta_of_muon->at(ml),ecal_3x3->at(ml)); - h_Phi_ecal[cut]->Fill(phi_of_muon->at(ml),ecal_3x3->at(ml)); - h_MuonHittingEcal[cut]->Fill(typeEcal); - if (typeEcal == 1) { - h_EtaX_ecal[cut]->Fill(etaEcal,ecal_3x3->at(ml)); - h_PhiY_ecal[cut]->Fill(phiEcal,ecal_3x3->at(ml)); - } - - h_hcal_energy[cut]->Fill(hcal_3into3->at(ml)); - h_1x1_hcal[cut]->Fill(hcal_1x1->at(ml)); - h_EtaX_hcal[cut]->Fill(etaXHcal,hcal_1x1->at(ml)); - h_PhiY_hcal[cut]->Fill(phiYHcal,hcal_1x1->at(ml)); - h_HotCell[cut]->Fill(hcal_cellHot->at(ml)); - if (mergeDepth_) { - double en1(0), en2(0), energyFill(0), chargeS(0), chargeBG(0); - double enh(0), enc(0); - for (int dep=0; dep 0) ? indxEta[eta][0][PHI] : 1+indxEta[eta][0][PHI]; - if (debug_)// || eta==15 || eta==17) - std::cout << "Matched Id " << matchedId->at(ml) << " Hot " - << hcal_cellHot->at(ml) << " eta " << etaHcal << ":" - << eta << " phi " << phiHcal << ":" << PHI - << " Index " << ind << " E " << en1 << ":" << en2 - << ":" << enh << ":" << enc << " L " << energyFill - << " Charge " << chargeS << ":" << chargeBG <at(ml))) continue; - if (hcal_cellHot->at(ml)==1) { - if (energyFill > 0) { - h_Hot_MuonEnergy_hcal_HotCell[cut][ind]->Fill(en2); - h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[cut][ind]->Fill(en2/energyFill); - h_active_length_Fill[cut][ind]->Fill(energyFill); - h_p_muon_ineta[cut][ind]->Fill(p_of_muon->at(ml)); - h_charge_signal[cut][ind]->Fill(chargeS); - h_charge_bg[cut][ind]->Fill(chargeBG); - - t_ene.push_back(enh); - t_enec.push_back(enc); - t_charge.push_back(chargeS); - t_actln.push_back(energyFill); - t_depth.push_back(0); - - outtree_->Fill(); - } - } - } else { - bool fillTree(false); - for (int dep=0; dep 0) ? indxEta[eta][dep][PHI] : 1+indxEta[eta][dep][PHI]; - if (debug_)// || eta==15 || eta==17) - std::cout << "Matched Id " << matchedId->at(ml) << " Hot " - << hcal_cellHot->at(ml) << " eta " << etaHcal << ":" - << eta << " phi " << phiHcal << ":" << PHI - << " depth " << dep << " Index " << ind << " E " - << en1 << ":" << en2 << ":" << enh << ":" << enc - << " L " << energyFill << " Charge " << chargeS - << ":" << chargeBG << std::endl; - if (!(matchedId->at(ml))) continue; - if (ok1) { - if (debug_) std::cout<<"enter ok1"<at(ml)==1) { - if (energyFill > 0) { - h_Hot_MuonEnergy_hcal_HotCell[cut][ind]->Fill(en2); - h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[cut][ind]->Fill(en2/energyFill); - h_active_length_Fill[cut][ind]->Fill(energyFill); - h_p_muon_ineta[cut][ind]->Fill(p_of_muon->at(ml)); - h_charge_signal[cut][ind]->Fill(chargeS); - h_charge_bg[cut][ind]->Fill(chargeBG); - t_ene.push_back(enh); - t_enec.push_back(enc); - t_charge.push_back(chargeS); - t_actln.push_back(energyFill); - // added depth vector AmanKalsi - t_depth.push_back(dep); - fillTree = true; - } else { - t_ene.push_back(-999.0); - t_enec.push_back(-999.0); - t_charge.push_back(-999.0); - t_actln.push_back(-999.0); - t_depth.push_back(-999.0); - } - if(debug_) std::cout<<"enter hot cell"<at(ml)!=1) { - } - } - - if(debug_) std::cout<<"ETA \t"<Fill(); - } - } + etaPhiEcal(ecal_detID->at(ml), typeEcal, zsideEcal, etaXEcal, phiYEcal, planeEcal, stripEcal); + double etaEcal = (etaXEcal - 0.5) * zsideEcal; + double phiEcal = phiYEcal - 0.5; + + if (debug_) + std::cout << "hcal_det_id " << std::hex << hcal_detID->at(ml) << std::dec; + + int etaHcal, phiHcal, depthHcal; + etaPhiHcal(hcal_detID->at(ml), etaHcal, phiHcal, depthHcal); + + int eta = (etaHcal > 0) ? etaHcal - 1 : -(1 + etaHcal); + double etaXHcal = (etaHcal > 0) ? etaHcal - 0.5 : etaHcal + 0.5; + int nDepth = nDepthBins(eta + 1, phiHcal); + int nPhi = nPhiBins(eta + 1); + int PHI = (nPhi > 36) ? (phiHcal - 1) : (phiHcal - 1) / 2; + double phiYHcal = (phiHcal - 0.5); + t_ieta = etaHcal; + t_iphi = PHI; + t_p = p_of_muon->at(ml); + t_ediff = hcal_3into3->at(ml) - hcal_1x1->at(ml); + t_nvtx = GoodVertex; + if (p_of_muon->at(ml) > 20) + pcut = true; + if (pt_of_muon->at(ml) > 20) + ptcut = true; + if (looseMuon(ml)) + loose = true; + if (softMuon(ml)) + soft = true; + if (tightMuon(ml)) + tight = true; + + if (debug_) + std::cout << " etaHcal " << etaHcal << ":" << etaXHcal << " phiHcal " << phiHcal << ":" << phiYHcal << ":" + << PHI << " Depth " << nDepth << " Muon Pt " << pt_of_muon->at(ml) << " Isol " << IsolationR04->at(ml) + << std::endl; + + for (int cut = 0; cut < nCut_; ++cut) { + bool select(false); + if (cut == 0) + select = tightMuon(ml); + else if (cut == 1) + select = softMuon(ml); + else + select = looseMuon(ml); + + if (select && ((eta + 1) >= etaMin_) && ((eta + 1) <= etaMax_)) { + // h_P_Muon[cut]->Fill(p_of_muon->at(ml)); + h_Pt_Muon[cut]->Fill(pt_of_muon->at(ml)); + h_Eta_Muon[cut]->Fill(eta_of_muon->at(ml)); + h_Phi_Muon[cut]->Fill(phi_of_muon->at(ml)); + h_PF_Muon[cut]->Fill(PF_Muon->at(ml)); + h_GlobTrack_Chi[cut]->Fill(GlobTrack_Chi->at(ml)); + h_Global_Muon_Hits[cut]->Fill(Global_Muon_Hits->at(ml)); + h_MatchedStations[cut]->Fill(MatchedStations->at(ml)); + h_Tight_TransImpactparameter[cut]->Fill(Tight_TransImpactparameter->at(ml)); + h_Tight_LongitudinalImpactparameter[cut]->Fill(Tight_LongitudinalImpactparameter->at(ml)); + h_InnerTrackPixelHits[cut]->Fill(InnerTrackPixelHits->at(ml)); + h_TrackerLayer[cut]->Fill(TrackerLayer->at(ml)); + h_IsolationR04[cut]->Fill(IsolationR04->at(ml)); + h_Global_Muon[cut]->Fill(Global_Muon->at(ml)); + + h_TransImpactParameter[cut]->Fill(Tight_TransImpactparameter->at(ml)); + h_LongImpactParameter[cut]->Fill(Tight_LongitudinalImpactparameter->at(ml)); + + //in Phi Bins + if (((phi_of_muon->at(ml)) >= -1.5) || ((phi_of_muon->at(ml)) <= 0.5)) { + h_TransImpactParameterBin1[cut]->Fill(Tight_TransImpactparameter->at(ml)); + h_LongImpactParameterBin1[cut]->Fill(Tight_LongitudinalImpactparameter->at(ml)); + h_2D_Bin1[cut]->Fill(Tight_TransImpactparameter->at(ml), Tight_LongitudinalImpactparameter->at(ml)); + } + + if ((phi_of_muon->at(ml) > 0.5) || (phi_of_muon->at(ml) < -1.5)) { + h_TransImpactParameterBin2[cut]->Fill(Tight_TransImpactparameter->at(ml)); + h_LongImpactParameterBin2[cut]->Fill(Tight_LongitudinalImpactparameter->at(ml)); + h_2D_Bin2[cut]->Fill(Tight_TransImpactparameter->at(ml), Tight_LongitudinalImpactparameter->at(ml)); + } + + h_ecal_energy[cut]->Fill(ecal_3into3->at(ml)); + h_3x3_ecal[cut]->Fill(ecal_3x3->at(ml)); + h_Eta_ecal[cut]->Fill(eta_of_muon->at(ml), ecal_3x3->at(ml)); + h_Phi_ecal[cut]->Fill(phi_of_muon->at(ml), ecal_3x3->at(ml)); + h_MuonHittingEcal[cut]->Fill(typeEcal); + if (typeEcal == 1) { + h_EtaX_ecal[cut]->Fill(etaEcal, ecal_3x3->at(ml)); + h_PhiY_ecal[cut]->Fill(phiEcal, ecal_3x3->at(ml)); + } + + h_hcal_energy[cut]->Fill(hcal_3into3->at(ml)); + h_1x1_hcal[cut]->Fill(hcal_1x1->at(ml)); + h_EtaX_hcal[cut]->Fill(etaXHcal, hcal_1x1->at(ml)); + h_PhiY_hcal[cut]->Fill(phiYHcal, hcal_1x1->at(ml)); + h_HotCell[cut]->Fill(hcal_cellHot->at(ml)); + if (mergeDepth_) { + double en1(0), en2(0), energyFill(0), chargeS(0), chargeBG(0); + double enh(0), enc(0); + for (int dep = 0; dep < nDepth; ++dep) { + double enb(0), enu(0), eh0(0), ec0(0), chgS(0), chgB(0), actL(0); + getEnergy(ml, dep, enb, enu, eh0, ec0, chgS, chgB, actL); + en1 += ((useCorrect_) ? enu : enb); + en2 += ((useCorrect_) ? ec0 : eh0); + enh += (eh0); + enc += (ec0); + energyFill += (actL); + chargeS += (chgS); + chargeBG += (chgB); + } + int ind = (etaHcal > 0) ? indxEta[eta][0][PHI] : 1 + indxEta[eta][0][PHI]; + if (debug_) // || eta==15 || eta==17) + std::cout << "Matched Id " << matchedId->at(ml) << " Hot " << hcal_cellHot->at(ml) << " eta " << etaHcal + << ":" << eta << " phi " << phiHcal << ":" << PHI << " Index " << ind << " E " << en1 << ":" + << en2 << ":" << enh << ":" << enc << " L " << energyFill << " Charge " << chargeS << ":" + << chargeBG << std::endl; + if (!(matchedId->at(ml))) + continue; + if (hcal_cellHot->at(ml) == 1) { + if (energyFill > 0) { + h_Hot_MuonEnergy_hcal_HotCell[cut][ind]->Fill(en2); + h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[cut][ind]->Fill(en2 / energyFill); + h_active_length_Fill[cut][ind]->Fill(energyFill); + h_p_muon_ineta[cut][ind]->Fill(p_of_muon->at(ml)); + h_charge_signal[cut][ind]->Fill(chargeS); + h_charge_bg[cut][ind]->Fill(chargeBG); + + t_ene.push_back(enh); + t_enec.push_back(enc); + t_charge.push_back(chargeS); + t_actln.push_back(energyFill); + t_depth.push_back(0); + + outtree_->Fill(); + } + } + } else { + bool fillTree(false); + for (int dep = 0; dep < nDepth; ++dep) { + if (debug_) + std::cout << "dep:" << dep << std::endl; + + double energyFill(0), chargeS(-9999), chargeBG(-9999); + double enh(-9999), enc(-9999), enb(0), enu(0); + bool ok1 = getEnergy(ml, dep, enb, enu, enh, enc, chargeS, chargeBG, energyFill); + double en1 = ((useCorrect_) ? enu : enb); + double en2 = ((useCorrect_) ? enc : enh); + if (debug_) + std::cout << "Hello in " << dep + 1 << " " << en1 << ":" << en2 << ":" << energyFill << std::endl; + + bool ok2 = ok1; + + if (debug_) + std::cout << "Before Index " << ok1 << ":" << ok2 << std::endl; + + int ind = (etaHcal > 0) ? indxEta[eta][dep][PHI] : 1 + indxEta[eta][dep][PHI]; + if (debug_) // || eta==15 || eta==17) + std::cout << "Matched Id " << matchedId->at(ml) << " Hot " << hcal_cellHot->at(ml) << " eta " << etaHcal + << ":" << eta << " phi " << phiHcal << ":" << PHI << " depth " << dep << " Index " << ind + << " E " << en1 << ":" << en2 << ":" << enh << ":" << enc << " L " << energyFill << " Charge " + << chargeS << ":" << chargeBG << std::endl; + if (!(matchedId->at(ml))) + continue; + if (ok1) { + if (debug_) + std::cout << "enter ok1" << std::endl; + + if (hcal_cellHot->at(ml) == 1) { + if (energyFill > 0) { + h_Hot_MuonEnergy_hcal_HotCell[cut][ind]->Fill(en2); + h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[cut][ind]->Fill(en2 / energyFill); + h_active_length_Fill[cut][ind]->Fill(energyFill); + h_p_muon_ineta[cut][ind]->Fill(p_of_muon->at(ml)); + h_charge_signal[cut][ind]->Fill(chargeS); + h_charge_bg[cut][ind]->Fill(chargeBG); + t_ene.push_back(enh); + t_enec.push_back(enc); + t_charge.push_back(chargeS); + t_actln.push_back(energyFill); + // added depth vector AmanKalsi + t_depth.push_back(dep); + fillTree = true; + } else { + t_ene.push_back(-999.0); + t_enec.push_back(-999.0); + t_charge.push_back(-999.0); + t_actln.push_back(-999.0); + t_depth.push_back(-999.0); + } + if (debug_) + std::cout << "enter hot cell" << std::endl; + } + } + + if (ok2) { + if (debug_) + std::cout << "enter ok2" << std::endl; + if (hcal_cellHot->at(ml) != 1) { + } + } + + if (debug_) + std::cout << "ETA \t" << eta << "DEPTH \t" << dep << std::endl; + } + if (fillTree) + outtree_->Fill(); + } + } } } int evtype(0); - if (pcut) evtype += 1; - if (ptcut) evtype += 2; - if (loose) evtype += 4; - if (soft) evtype += 8; - if (tight) evtype += 16; + if (pcut) + evtype += 1; + if (ptcut) + evtype += 2; + if (loose) + evtype += 4; + if (soft) + evtype += 8; + if (tight) + evtype += 16; h_evtype->Fill(evtype); } close(); @@ -1046,133 +1105,126 @@ Bool_t HBHEMuonOfflineAnalyzer::Notify() { // is started when using PROOF. It is normally not necessary to make changes // to the generated code, but the routine can be extended by the // user if needed. The return value is currently not used. - + return kTRUE; } void HBHEMuonOfflineAnalyzer::Show(Long64_t entry) { // Print contents of entry. // If entry is not specified, print current entry - if (!fChain) return; + if (!fChain) + return; fChain->Show(entry); } -bool HBHEMuonOfflineAnalyzer::fillChain(TChain *chain, - const char* inputFileList) { - +bool HBHEMuonOfflineAnalyzer::fillChain(TChain *chain, const char *inputFileList) { std::string fname(inputFileList); - if (fname.substr(fname.size()-5,5) == ".root") { + if (fname.substr(fname.size() - 5, 5) == ".root") { chain->Add(fname.c_str()); } else { std::ifstream infile(inputFileList); if (!infile.is_open()) { - std::cout << "** ERROR: Can't open '" << inputFileList << "' for input" - << std::endl; + std::cout << "** ERROR: Can't open '" << inputFileList << "' for input" << std::endl; return false; } while (1) { infile >> fname; - if (!infile.good()) break; + if (!infile.good()) + break; chain->Add(fname.c_str()); } infile.close(); } - std::cout << "No. of Entries in this tree : " << chain->GetEntries() - << std::endl; + std::cout << "No. of Entries in this tree : " << chain->GetEntries() << std::endl; return true; } -bool HBHEMuonOfflineAnalyzer::readCorr(const char* infile) { - +bool HBHEMuonOfflineAnalyzer::readCorr(const char *infile) { std::ifstream fInput(infile); unsigned int ncorr(0), all(0), good(0); if (!fInput.good()) { std::cout << "Cannot open file " << infile << std::endl; } else { - char buffer [1024]; + char buffer[1024]; while (fInput.getline(buffer, 1024)) { ++all; std::string bufferString(buffer); - if (bufferString.substr(0,5) == "#IOVs") { - std::vector items = splitString(bufferString.substr(6)); - ncorr = items.size() - 1; - for (unsigned int n=0; n items = splitString(bufferString.substr(6)); + ncorr = items.size() - 1; + for (unsigned int n = 0; n < ncorr; ++n) { + int run = std::atoi(items[n].c_str()); + runlow_.push_back(run); + } + std::cout << ncorr << ":" << runlow_.size() << " Run ranges" << std::endl; + for (unsigned int n = 0; n < runlow_.size(); ++n) + std::cout << " [" << n << "] " << runlow_[n]; + std::cout << std::endl; + } else if (buffer[0] == '#') { + continue; //ignore other comments } else { - std::vector items = splitString(bufferString); - if (items.size () != ncorr+3) { - std::cout << "Ignore line: " << buffer << std::endl; - } else { - ++good; - int ieta = std::atoi (items[0].c_str()); - int iphi = std::atoi (items[1].c_str()); - int depth = std::atoi (items[2].c_str()); - unsigned int id = getDetIdHBHE(ieta,iphi,depth); - for (unsigned int n=0; n items = splitString(bufferString); + if (items.size() != ncorr + 3) { + std::cout << "Ignore line: " << buffer << std::endl; + } else { + ++good; + int ieta = std::atoi(items[0].c_str()); + int iphi = std::atoi(items[1].c_str()); + int depth = std::atoi(items[2].c_str()); + unsigned int id = getDetIdHBHE(ieta, iphi, depth); + for (unsigned int n = 0; n < ncorr; ++n) { + float corrf = std::atof(items[n + 3].c_str()); + if (n < nmax_) + corrFac_[n][id] = corrf; + } + if (debug_) { + std::cout << "ID " << std::hex << id << std::dec << ":" << id << " (eta " << ieta << " phi " << iphi + << " depth " << depth << ")"; + for (unsigned int n = 0; n < ncorr; ++n) + std::cout << " " << corrFac_[n][id]; + std::cout << std::endl; + } + } } } fInput.close(); - std::cout << "Reads total of " << all << " and " << good << " good records" - << std::endl; + std::cout << "Reads total of " << all << " and " << good << " good records" << std::endl; } return (good > 0); } -void HBHEMuonOfflineAnalyzer::bookHistograms(const char* fname) { - +void HBHEMuonOfflineAnalyzer::bookHistograms(const char *fname) { std::cout << "BookHistograms" << std::endl; - output_file = TFile::Open(fname,"RECREATE"); + output_file = TFile::Open(fname, "RECREATE"); output_file->cd(); - outtree_ = new TTree("Lep_Tree","Lep_Tree"); - outtree_->Branch("t_ieta", &t_ieta); - outtree_->Branch("t_iphi", &t_iphi); - outtree_->Branch("t_nvtx", &t_nvtx); - outtree_->Branch("t_p", &t_p); - outtree_->Branch("t_ediff", &t_ediff); - outtree_->Branch("t_ene", &t_ene); - outtree_->Branch("t_enec", &t_enec); - outtree_->Branch("t_charge", &t_charge); - outtree_->Branch("t_actln", &t_actln); - outtree_->Branch("t_depth", &t_depth); - - std::string type[]={"tight"};//,"soft","loose"}; + outtree_ = new TTree("Lep_Tree", "Lep_Tree"); + outtree_->Branch("t_ieta", &t_ieta); + outtree_->Branch("t_iphi", &t_iphi); + outtree_->Branch("t_nvtx", &t_nvtx); + outtree_->Branch("t_p", &t_p); + outtree_->Branch("t_ediff", &t_ediff); + outtree_->Branch("t_ene", &t_ene); + outtree_->Branch("t_enec", &t_enec); + outtree_->Branch("t_charge", &t_charge); + outtree_->Branch("t_actln", &t_actln); + outtree_->Branch("t_depth", &t_depth); + + std::string type[] = {"tight"}; //,"soft","loose"}; char name[128], title[500]; nHist = 0; - for (int eta=etaMin_; eta<=etaMax_; ++eta) { - - int nDepth = nDepthBins(eta,-1); - int nPhi = nPhiBins(eta); + for (int eta = etaMin_; eta <= etaMax_; ++eta) { + int nDepth = nDepthBins(eta, -1); + int nPhi = nPhiBins(eta); //std::cout<<"problem 2"<cd(); - h_evtype = new TH1D("EvType", "Event Type", 100,0,100); - for (int i=0; i 0.5 and #phi< -1.5 ", type[i].c_str()); - h_TransImpactParameterBin2[i] = new TH1D(name, title,100,0,0.5); + h_evtype = new TH1D("EvType", "Event Type", 100, 0, 100); + for (int i = 0; i < nCut_; ++i) { + sprintf(name, "h_Pt_Muon_%s", type[i].c_str()); + sprintf(title, "p_{T} of %s muons (GeV)", type[i].c_str()); + h_Pt_Muon[i] = new TH1D(name, title, 100, 0, 200); + + sprintf(name, "h_Eta_Muon_%s", type[i].c_str()); + sprintf(title, "#eta of %s muons", type[i].c_str()); + h_Eta_Muon[i] = new TH1D(name, title, 50, -2.5, 2.5); + + sprintf(name, "h_Phi_Muon_%s", type[i].c_str()); + sprintf(title, "#phi of %s muons", type[i].c_str()); + h_Phi_Muon[i] = new TH1D(name, title, 100, -3.1415926, 3.1415926); + + sprintf(name, "h_P_Muon_%s", type[i].c_str()); + sprintf(title, "p of %s muons (GeV)", type[i].c_str()); + h_P_Muon[i] = new TH1D(name, title, 100, 0, 200); + + sprintf(name, "h_PF_Muon_%s", type[i].c_str()); + sprintf(title, "PF %s muons (GeV)", type[i].c_str()); + h_PF_Muon[i] = new TH1D(name, title, 2, 0, 2); + + sprintf(name, "h_Global_Muon_Chi2_%s", type[i].c_str()); + sprintf(title, "Chi2 Global %s muons (GeV)", type[i].c_str()); + h_GlobTrack_Chi[i] = new TH1D(name, title, 15, 0, 15); + + sprintf(name, "h_Global_Muon_Hits_%s", type[i].c_str()); + sprintf(title, "Global Hits %s muons (GeV)", type[i].c_str()); + h_Global_Muon_Hits[i] = new TH1D(name, title, 10, 0, 10); + + sprintf(name, "h_Matched_Stations_%s", type[i].c_str()); + sprintf(title, "Matched Stations %s muons (GeV)", type[i].c_str()); + h_MatchedStations[i] = new TH1D(name, title, 10, 0, 10); + + sprintf(name, "h_Transverse_ImpactParameter_%s", type[i].c_str()); + sprintf(title, "Transverse_ImpactParameter of %s muons (GeV)", type[i].c_str()); + h_Tight_TransImpactparameter[i] = new TH1D(name, title, 50, 0, 10); + + sprintf(name, "h_Longitudinal_ImpactParameter_%s", type[i].c_str()); + sprintf(title, "Longitudinal_ImpactParameter of %s muons (GeV)", type[i].c_str()); + h_Tight_LongitudinalImpactparameter[i] = new TH1D(name, title, 20, 0, 10); + + sprintf(name, "h_InnerTrack_PixelHits_%s", type[i].c_str()); + sprintf(title, "InnerTrack_PixelHits of %s muons (GeV)", type[i].c_str()); + h_InnerTrackPixelHits[i] = new TH1D(name, title, 20, 0, 20); + + sprintf(name, "h_TrackLayers_%s", type[i].c_str()); + sprintf(title, "No. of Tracker Layers of %s muons (GeV)", type[i].c_str()); + h_TrackerLayer[i] = new TH1D(name, title, 20, 0, 20); + ; + + sprintf(name, "h_IsolationR04_%s", type[i].c_str()); + sprintf(title, "IsolationR04 %s muons (GeV)", type[i].c_str()); + h_IsolationR04[i] = new TH1D(name, title, 45, 0, 5); + ; + + sprintf(name, "h_Global_Muon_%s", type[i].c_str()); + sprintf(title, "Global %s muons (GeV)", type[i].c_str()); + h_Global_Muon[i] = new TH1D(name, title, 2, 0, 2); + + sprintf(name, "h_TransImpactParameter_%s", type[i].c_str()); + sprintf(title, "TransImpactParameter of %s muons (GeV)", type[i].c_str()); + h_TransImpactParameter[i] = new TH1D(name, title, 100, 0, 0.5); + + sprintf(name, "h_TransImpactParameterBin1_%s", type[i].c_str()); + sprintf(title, "TransImpactParameter of %s muons (GeV) in -1.5 <= #phi <= 0.5", type[i].c_str()); + h_TransImpactParameterBin1[i] = new TH1D(name, title, 100, 0, 0.5); + + sprintf(name, "h_TransImpactParameterBin2_%s", type[i].c_str()); + sprintf(title, "TransImpactParameter of %s muons (GeV) in #phi> 0.5 and #phi< -1.5 ", type[i].c_str()); + h_TransImpactParameterBin2[i] = new TH1D(name, title, 100, 0, 0.5); // - sprintf (name, "h_LongImpactParameter_%s", type[i].c_str()); - sprintf (title, "LongImpactParameter of %s muons (GeV)", type[i].c_str()); - h_LongImpactParameter[i] = new TH1D(name, title,100,0,30); - - sprintf (name, "h_LongImpactParameterBin1_%s", type[i].c_str()); - sprintf (title, "LongImpactParameter of %s muons (GeV) in -1.5 <= #phi <= 0.5", type[i].c_str()); - h_LongImpactParameterBin1[i] = new TH1D(name, title,100,0,30); - - sprintf (name, "h_LongImpactParameterBin2_%s", type[i].c_str()); - sprintf (title, "LongImpactParameter of %s muons (GeV) in #phi> 0.5 and #phi< -1.5 ", type[i].c_str()); - h_LongImpactParameterBin2[i] = new TH1D(name, title,100,0,30); - - sprintf (name, "h_2D_Bin1_%s", type[i].c_str()); - sprintf (title, "Trans/Long ImpactParameter of %s muons (GeV) in -1.5 <= #phi< 0.5 ", type[i].c_str()); - h_2D_Bin1[i] = new TH2D(name, title, 100,0,0.5,100,0,30); - - sprintf (name, "h_2D_Bin2_%s", type[i].c_str()); - sprintf (title, "Trans/Long ImpactParameter of %s muons (GeV) in #phi> 0.5 and #phi< -1.5 ", type[i].c_str()); - h_2D_Bin2[i] = new TH2D(name, title, 100,0,0.5,100,0,30); - - sprintf (name, "h_ecal_energy_%s", type[i].c_str()); - sprintf (title, "ECAL energy for %s muons", type[i].c_str()); - h_ecal_energy[i] = new TH1D(name, title,1000,-10.0,90.0); - - sprintf (name, "h_hcal_energy_%s", type[i].c_str()); - sprintf (title, "HCAL energy for %s muons", type[i].c_str()); - h_hcal_energy[i] = new TH1D(name, title,500,-10.0,90.0); - - sprintf (name, "h_3x3_ecal_%s", type[i].c_str()); - sprintf (title, "ECAL energy in 3x3 for %s muons", type[i].c_str()); - h_3x3_ecal[i] = new TH1D(name, title,1000,-10.0,90.0); - - sprintf (name, "h_1x1_hcal_%s", type[i].c_str()); - sprintf (title, "HCAL energy in 1x1 for %s muons", type[i].c_str()); - h_1x1_hcal[i] = new TH1D(name, title,500,-10.0,90.0); - - sprintf (name, "h_EtaX_hcal_%s", type[i].c_str()); - sprintf (title, "HCAL energy as a function of i#eta for %s muons", type[i].c_str()); - h_EtaX_hcal[i] = new TProfile(name, title,60,-30.0,30.0); - - sprintf (name, "h_PhiY_hcal_%s", type[i].c_str()); - sprintf (title, "HCAL energy as a function of i#phi for %s muons", type[i].c_str()); - h_PhiY_hcal[i] = new TProfile(name, title,72,0,72); - - sprintf (name, "h_EtaX_ecal_%s", type[i].c_str()); - sprintf (title, "EB energy as a function of i#eta for %s muons", type[i].c_str()); - h_EtaX_ecal[i] = new TProfile(name, title,170,-85.0,85.0); - - sprintf (name, "h_PhiY_ecal_%s", type[i].c_str()); - sprintf (title, "EB energy as a function of i#phi for %s muons", type[i].c_str()); - h_PhiY_ecal[i] = new TProfile(name, title,360,0,360); - - sprintf (name, "h_Eta_ecal_%s", type[i].c_str()); - sprintf (title, "ECAL energy as a function of #eta for %s muons", type[i].c_str()); - h_Eta_ecal[i] = new TProfile(name, title,100,-2.5,2.5); - - sprintf (name, "h_Phi_ecal_%s", type[i].c_str()); - sprintf (title, "ECAL energy as a function of #phi for %s muons", type[i].c_str()); - h_Phi_ecal[i] = new TProfile(name, title,100,-3.1415926,3.1415926); - - sprintf (name, "h_MuonHittingEcal_%s", type[i].c_str()); - sprintf (title, "%s muons hitting ECAL", type[i].c_str()); - h_MuonHittingEcal[i] = new TH1D(name, title,100,0,5.0); - - sprintf (name, "h_HotCell_%s", type[i].c_str()); - sprintf (title, "Hot cell for %s muons", type[i].c_str()); - h_HotCell[i] = new TH1D(name, title,100,0,2); + sprintf(name, "h_LongImpactParameter_%s", type[i].c_str()); + sprintf(title, "LongImpactParameter of %s muons (GeV)", type[i].c_str()); + h_LongImpactParameter[i] = new TH1D(name, title, 100, 0, 30); + + sprintf(name, "h_LongImpactParameterBin1_%s", type[i].c_str()); + sprintf(title, "LongImpactParameter of %s muons (GeV) in -1.5 <= #phi <= 0.5", type[i].c_str()); + h_LongImpactParameterBin1[i] = new TH1D(name, title, 100, 0, 30); + + sprintf(name, "h_LongImpactParameterBin2_%s", type[i].c_str()); + sprintf(title, "LongImpactParameter of %s muons (GeV) in #phi> 0.5 and #phi< -1.5 ", type[i].c_str()); + h_LongImpactParameterBin2[i] = new TH1D(name, title, 100, 0, 30); + + sprintf(name, "h_2D_Bin1_%s", type[i].c_str()); + sprintf(title, "Trans/Long ImpactParameter of %s muons (GeV) in -1.5 <= #phi< 0.5 ", type[i].c_str()); + h_2D_Bin1[i] = new TH2D(name, title, 100, 0, 0.5, 100, 0, 30); + + sprintf(name, "h_2D_Bin2_%s", type[i].c_str()); + sprintf(title, "Trans/Long ImpactParameter of %s muons (GeV) in #phi> 0.5 and #phi< -1.5 ", type[i].c_str()); + h_2D_Bin2[i] = new TH2D(name, title, 100, 0, 0.5, 100, 0, 30); + + sprintf(name, "h_ecal_energy_%s", type[i].c_str()); + sprintf(title, "ECAL energy for %s muons", type[i].c_str()); + h_ecal_energy[i] = new TH1D(name, title, 1000, -10.0, 90.0); + + sprintf(name, "h_hcal_energy_%s", type[i].c_str()); + sprintf(title, "HCAL energy for %s muons", type[i].c_str()); + h_hcal_energy[i] = new TH1D(name, title, 500, -10.0, 90.0); + + sprintf(name, "h_3x3_ecal_%s", type[i].c_str()); + sprintf(title, "ECAL energy in 3x3 for %s muons", type[i].c_str()); + h_3x3_ecal[i] = new TH1D(name, title, 1000, -10.0, 90.0); + + sprintf(name, "h_1x1_hcal_%s", type[i].c_str()); + sprintf(title, "HCAL energy in 1x1 for %s muons", type[i].c_str()); + h_1x1_hcal[i] = new TH1D(name, title, 500, -10.0, 90.0); + + sprintf(name, "h_EtaX_hcal_%s", type[i].c_str()); + sprintf(title, "HCAL energy as a function of i#eta for %s muons", type[i].c_str()); + h_EtaX_hcal[i] = new TProfile(name, title, 60, -30.0, 30.0); + + sprintf(name, "h_PhiY_hcal_%s", type[i].c_str()); + sprintf(title, "HCAL energy as a function of i#phi for %s muons", type[i].c_str()); + h_PhiY_hcal[i] = new TProfile(name, title, 72, 0, 72); + + sprintf(name, "h_EtaX_ecal_%s", type[i].c_str()); + sprintf(title, "EB energy as a function of i#eta for %s muons", type[i].c_str()); + h_EtaX_ecal[i] = new TProfile(name, title, 170, -85.0, 85.0); + + sprintf(name, "h_PhiY_ecal_%s", type[i].c_str()); + sprintf(title, "EB energy as a function of i#phi for %s muons", type[i].c_str()); + h_PhiY_ecal[i] = new TProfile(name, title, 360, 0, 360); + + sprintf(name, "h_Eta_ecal_%s", type[i].c_str()); + sprintf(title, "ECAL energy as a function of #eta for %s muons", type[i].c_str()); + h_Eta_ecal[i] = new TProfile(name, title, 100, -2.5, 2.5); + + sprintf(name, "h_Phi_ecal_%s", type[i].c_str()); + sprintf(title, "ECAL energy as a function of #phi for %s muons", type[i].c_str()); + h_Phi_ecal[i] = new TProfile(name, title, 100, -3.1415926, 3.1415926); + + sprintf(name, "h_MuonHittingEcal_%s", type[i].c_str()); + sprintf(title, "%s muons hitting ECAL", type[i].c_str()); + h_MuonHittingEcal[i] = new TH1D(name, title, 100, 0, 5.0); + + sprintf(name, "h_HotCell_%s", type[i].c_str()); + sprintf(title, "Hot cell for %s muons", type[i].c_str()); + h_HotCell[i] = new TH1D(name, title, 100, 0, 2); // output_file->cd(); - for (int eta=etaMin_; eta<=etaMax_; ++eta) { - int nDepth = nDepthBins(eta,-1); - int nPhi = nPhiBins(eta); + for (int eta = etaMin_; eta <= etaMax_; ++eta) { + int nDepth = nDepthBins(eta, -1); + int nPhi = nPhiBins(eta); //sprintf(name, "Dir_muon_type_%s_ieta%d",type[i].c_str(), eta); //d_output_file[i][eta]= output_file->mkdir(name); //output_file->cd(name); //d_output_file[i][eta]->cd(); - for (int depth=0; depthSumw2(); - - sprintf (name, "h_Hot_MuonEnergy_hc_%d_%d_%d_%s_HotCell_ByActiveLength", eta, (depth+1), PHI0, type[i].c_str()); - sprintf (title, "HCAL energy in hot tower (i#eta=%d, depth=%d, i#phi = %d) for extrapolated %s muons (Hot Cell) divided by Active Length", eta, (depth+1), PHI0, type[i].c_str()); - h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[i][ih] = new TH1D(name, title,4000,0.0,10.0); - h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[i][ih]->Sumw2(); - - sprintf (name, "h_active_length_Fill_%d_%d_%d_%s", eta, (depth+1), PHI0, type[i].c_str()); - sprintf (title, "active_length%d_%d_%d_%s", eta, (depth+1), PHI0, type[i].c_str()); - h_active_length_Fill[i][ih] = new TH1D(name, title,20,0,20); - h_active_length_Fill[i][ih]->Sumw2(); - - sprintf (name, "h_p_muon_in_%d_%d_%d_%s", eta, (depth+1), PHI0, type[i].c_str()); - sprintf (title, "p_muon_in%d_%d_%d_%s", eta, (depth+1), PHI0, type[i].c_str()); - h_p_muon_ineta[i][ih] = new TH1D(name, title, 500,0,500); - h_p_muon_ineta[i][ih]->Sumw2(); - - sprintf (name, "h_charge_signal_in_%d_%d_%d_%s", eta, (depth+1), PHI0, type[i].c_str()); - sprintf (name, "charge_signal_in_%d_%d_%d_%s", eta, (depth+1), PHI0, type[i].c_str()); - h_charge_signal[i][ih] = new TH1D(name, title, 500,0,500); - h_charge_signal[i][ih]->Sumw2(); - - sprintf (name, "h_charge_bg_in_%d_%d_%d_%s", eta, (depth+1), PHI0, type[i].c_str()); - sprintf (name, "charge_bg_in_%d_%d_%d_%s", eta, (depth+1), PHI0, type[i].c_str()); - h_charge_bg[i][ih] = new TH1D(name, title, 500,0,500); - h_charge_bg[i][ih]->Sumw2(); - - ih++; - - sprintf (name, "h_Hot_MuonEnergy_hc_%d_%d_%d_%s_HotCell", -eta, (depth+1), PHI0, type[i].c_str()); - sprintf (title, "HCAL energy in hot tower (i#eta=%d, depth=%d, i#phi = %d) for extrapolated %s muons (Hot Cell)", -eta, (depth+1), PHI0, type[i].c_str()); - h_Hot_MuonEnergy_hcal_HotCell[i][ih] = new TH1D(name, title,4000,0.0,10.0); - h_Hot_MuonEnergy_hcal_HotCell[i][ih] ->Sumw2(); - - sprintf (name, "h_Hot_MuonEnergy_hc_%d_%d_%d_%s_HotCell_ByActiveLength", -eta, (depth+1),PHI0, type[i].c_str()); - sprintf (title, "HCAL energy in hot tower (i#eta=%d, depth=%d, i#phi=%d) for extrapolated %s muons (Hot Cell) divided by Active Length", -eta, (depth+1), PHI0, type[i].c_str()); - h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[i][ih] = new TH1D(name, title,4000,0.0,10.0); - h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[i][ih]->Sumw2(); - - sprintf (name, "h_active_length_Fill_%d_%d_%d_%s",-eta, (depth+1), PHI0, type[i].c_str()); - sprintf (title, "active_length%d_%d_%d_%s",-eta, (depth+1), PHI0, type[i].c_str()); - h_active_length_Fill[i][ih] = new TH1D(name, title,20,0,20); - h_active_length_Fill[i][ih]->Sumw2(); - - sprintf (name, "h_p_muon_in_%d_%d_%d_%s",-eta, (depth+1), PHI0, type[i].c_str()); - sprintf (title, "p_muon_in%d_%d_%d_%s",-eta, (depth+1), PHI0, type[i].c_str()); - h_p_muon_ineta[i][ih] = new TH1D(name, title, 500,0,500); - h_p_muon_ineta[i][ih]->Sumw2(); - - sprintf (name, "h_charge_signal_in_%d_%d_%d_%s",-eta, (depth+1), PHI0, type[i].c_str()); - sprintf (name, "charge_signal_in_%d_%d_%d_%s",-eta, (depth+1), PHI0, type[i].c_str()); - h_charge_signal[i][ih] = new TH1D(name, title, 500,0,500); - h_charge_signal[i][ih]->Sumw2(); - - sprintf (name, "h_charge_bg_in_%d_%d_%d_%s",-eta, (depth+1), PHI0, type[i].c_str()); - sprintf (name, "charge_bg_in_%d_%d_%d_%s",-eta, (depth+1), PHI0, type[i].c_str()); - h_charge_bg[i][ih] = new TH1D(name, title, 500,0,500); - h_charge_bg[i][ih]->Sumw2(); - - } + for (int depth = 0; depth < nDepth; ++depth) { + for (int PHI = 0; PHI < nPhi; ++PHI) { + int PHI0 = (nPhi == 72) ? PHI + 1 : 2 * PHI + 1; + int ih = indxEta[eta - 1][depth][PHI]; + if (debug_) + std::cout << "eta:" << eta << " depth:" << depth << " PHI:" << PHI << ":" << PHI0 << " ih:" << ih + << std::endl; + + sprintf(name, "h_Hot_MuonEnergy_hc_%d_%d_%d_%s_HotCell", eta, (depth + 1), PHI0, type[i].c_str()); + sprintf(title, + "HCAL energy in hot tower (i#eta=%d, depth=%d, i#phi = %d) for extrapolated %s muons (Hot Cell)", + eta, + (depth + 1), + PHI0, + type[i].c_str()); + h_Hot_MuonEnergy_hcal_HotCell[i][ih] = new TH1D(name, title, 4000, 0.0, 10.0); + h_Hot_MuonEnergy_hcal_HotCell[i][ih]->Sumw2(); + + sprintf( + name, "h_Hot_MuonEnergy_hc_%d_%d_%d_%s_HotCell_ByActiveLength", eta, (depth + 1), PHI0, type[i].c_str()); + sprintf(title, + "HCAL energy in hot tower (i#eta=%d, depth=%d, i#phi = %d) for extrapolated %s muons (Hot Cell) " + "divided by Active Length", + eta, + (depth + 1), + PHI0, + type[i].c_str()); + h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[i][ih] = new TH1D(name, title, 4000, 0.0, 10.0); + h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[i][ih]->Sumw2(); + + sprintf(name, "h_active_length_Fill_%d_%d_%d_%s", eta, (depth + 1), PHI0, type[i].c_str()); + sprintf(title, "active_length%d_%d_%d_%s", eta, (depth + 1), PHI0, type[i].c_str()); + h_active_length_Fill[i][ih] = new TH1D(name, title, 20, 0, 20); + h_active_length_Fill[i][ih]->Sumw2(); + + sprintf(name, "h_p_muon_in_%d_%d_%d_%s", eta, (depth + 1), PHI0, type[i].c_str()); + sprintf(title, "p_muon_in%d_%d_%d_%s", eta, (depth + 1), PHI0, type[i].c_str()); + h_p_muon_ineta[i][ih] = new TH1D(name, title, 500, 0, 500); + h_p_muon_ineta[i][ih]->Sumw2(); + + sprintf(name, "h_charge_signal_in_%d_%d_%d_%s", eta, (depth + 1), PHI0, type[i].c_str()); + sprintf(name, "charge_signal_in_%d_%d_%d_%s", eta, (depth + 1), PHI0, type[i].c_str()); + h_charge_signal[i][ih] = new TH1D(name, title, 500, 0, 500); + h_charge_signal[i][ih]->Sumw2(); + + sprintf(name, "h_charge_bg_in_%d_%d_%d_%s", eta, (depth + 1), PHI0, type[i].c_str()); + sprintf(name, "charge_bg_in_%d_%d_%d_%s", eta, (depth + 1), PHI0, type[i].c_str()); + h_charge_bg[i][ih] = new TH1D(name, title, 500, 0, 500); + h_charge_bg[i][ih]->Sumw2(); + + ih++; + + sprintf(name, "h_Hot_MuonEnergy_hc_%d_%d_%d_%s_HotCell", -eta, (depth + 1), PHI0, type[i].c_str()); + sprintf(title, + "HCAL energy in hot tower (i#eta=%d, depth=%d, i#phi = %d) for extrapolated %s muons (Hot Cell)", + -eta, + (depth + 1), + PHI0, + type[i].c_str()); + h_Hot_MuonEnergy_hcal_HotCell[i][ih] = new TH1D(name, title, 4000, 0.0, 10.0); + h_Hot_MuonEnergy_hcal_HotCell[i][ih]->Sumw2(); + + sprintf( + name, "h_Hot_MuonEnergy_hc_%d_%d_%d_%s_HotCell_ByActiveLength", -eta, (depth + 1), PHI0, type[i].c_str()); + sprintf(title, + "HCAL energy in hot tower (i#eta=%d, depth=%d, i#phi=%d) for extrapolated %s muons (Hot Cell) " + "divided by Active Length", + -eta, + (depth + 1), + PHI0, + type[i].c_str()); + h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[i][ih] = new TH1D(name, title, 4000, 0.0, 10.0); + h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[i][ih]->Sumw2(); + + sprintf(name, "h_active_length_Fill_%d_%d_%d_%s", -eta, (depth + 1), PHI0, type[i].c_str()); + sprintf(title, "active_length%d_%d_%d_%s", -eta, (depth + 1), PHI0, type[i].c_str()); + h_active_length_Fill[i][ih] = new TH1D(name, title, 20, 0, 20); + h_active_length_Fill[i][ih]->Sumw2(); + + sprintf(name, "h_p_muon_in_%d_%d_%d_%s", -eta, (depth + 1), PHI0, type[i].c_str()); + sprintf(title, "p_muon_in%d_%d_%d_%s", -eta, (depth + 1), PHI0, type[i].c_str()); + h_p_muon_ineta[i][ih] = new TH1D(name, title, 500, 0, 500); + h_p_muon_ineta[i][ih]->Sumw2(); + + sprintf(name, "h_charge_signal_in_%d_%d_%d_%s", -eta, (depth + 1), PHI0, type[i].c_str()); + sprintf(name, "charge_signal_in_%d_%d_%d_%s", -eta, (depth + 1), PHI0, type[i].c_str()); + h_charge_signal[i][ih] = new TH1D(name, title, 500, 0, 500); + h_charge_signal[i][ih]->Sumw2(); + + sprintf(name, "h_charge_bg_in_%d_%d_%d_%s", -eta, (depth + 1), PHI0, type[i].c_str()); + sprintf(name, "charge_bg_in_%d_%d_%d_%s", -eta, (depth + 1), PHI0, type[i].c_str()); + h_charge_bg[i][ih] = new TH1D(name, title, 500, 0, 500); + h_charge_bg[i][ih]->Sumw2(); + } } //output_file->cd(); - } + } } //output_file->cd(); } -bool HBHEMuonOfflineAnalyzer::getEnergy(unsigned int ml, int dep, double& enb, - double& enu, double& enh, double& enc, - double& chgS, double& chgB, - double& actL) { +bool HBHEMuonOfflineAnalyzer::getEnergy(unsigned int ml, + int dep, + double &enb, + double &enu, + double &enh, + double &enc, + double &chgS, + double &chgB, + double &actL) { double cfac(1.0); - bool flag(true); + bool flag(true); if (cFactor_) { - int ieta = hcal_ieta->at(ml); - int iphi = hcal_iphi->at(ml); - unsigned int detId = getDetIdHBHE(ieta, iphi, dep+1); - cfac = getCorr(Run_No, detId); + int ieta = hcal_ieta->at(ml); + int iphi = hcal_iphi->at(ml); + unsigned int detId = getDetIdHBHE(ieta, iphi, dep + 1); + cfac = getCorr(Run_No, detId); } if (dep == 0) { - enb = cfac * hcal_edepth1->at(ml); - enu = cfac * hcal_edepthCorrect1->at(ml); - enh = cfac * hcal_edepthHot1->at(ml); - enc = cfac * hcal_edepthHotCorrect1->at(ml); + enb = cfac * hcal_edepth1->at(ml); + enu = cfac * hcal_edepthCorrect1->at(ml); + enh = cfac * hcal_edepthHot1->at(ml); + enc = cfac * hcal_edepthHotCorrect1->at(ml); chgS = hcal_cdepthHot1->at(ml); actL = hcal_activeHotL1->at(ml); chgB = hcal_cdepthHotBG1->at(ml); } else if (dep == 1) { - enb = cfac * hcal_edepth2->at(ml); - enu = cfac * hcal_edepthCorrect2->at(ml); - enh = cfac * hcal_edepthHot2->at(ml); - enc = cfac * hcal_edepthHotCorrect2->at(ml); + enb = cfac * hcal_edepth2->at(ml); + enu = cfac * hcal_edepthCorrect2->at(ml); + enh = cfac * hcal_edepthHot2->at(ml); + enc = cfac * hcal_edepthHotCorrect2->at(ml); chgS = hcal_cdepthHot2->at(ml); actL = hcal_activeHotL2->at(ml); chgB = hcal_cdepthHotBG2->at(ml); } else if (dep == 2) { - enb = cfac * hcal_edepth3->at(ml); - enu = cfac * hcal_edepthCorrect3->at(ml); - enh = cfac * hcal_edepthHot3->at(ml); - enc = cfac * hcal_edepthHotCorrect3->at(ml); + enb = cfac * hcal_edepth3->at(ml); + enu = cfac * hcal_edepthCorrect3->at(ml); + enh = cfac * hcal_edepthHot3->at(ml); + enc = cfac * hcal_edepthHotCorrect3->at(ml); chgS = hcal_cdepthHot3->at(ml); actL = hcal_activeHotL3->at(ml); chgB = hcal_cdepthHotBG3->at(ml); } else if (dep == 3) { - enb = cfac * hcal_edepth4->at(ml); - enu = cfac * hcal_edepthCorrect4->at(ml); - enh = cfac * hcal_edepthHot4->at(ml); - enc = cfac * hcal_edepthHotCorrect4->at(ml); + enb = cfac * hcal_edepth4->at(ml); + enu = cfac * hcal_edepthCorrect4->at(ml); + enh = cfac * hcal_edepthHot4->at(ml); + enc = cfac * hcal_edepthHotCorrect4->at(ml); chgS = hcal_cdepthHot4->at(ml); actL = hcal_activeHotL4->at(ml); chgB = hcal_cdepthHotBG4->at(ml); } else if (dep == 4) { if (hcal_edepthCorrect5->size() > ml) { - enb = cfac * hcal_edepth5->at(ml); - enu = cfac * hcal_edepthCorrect5->at(ml); - enh = cfac * hcal_edepthHot5->at(ml); - enc = cfac * hcal_edepthHotCorrect5->at(ml); + enb = cfac * hcal_edepth5->at(ml); + enu = cfac * hcal_edepthCorrect5->at(ml); + enh = cfac * hcal_edepthHot5->at(ml); + enc = cfac * hcal_edepthHotCorrect5->at(ml); chgS = hcal_cdepthHot5->at(ml); actL = hcal_activeHotL5->at(ml); chgB = hcal_cdepthHotBG5->at(ml); } else { - enb = enu = enh = enc = chgS = actL = chgB = 0; + enb = enu = enh = enc = chgS = actL = chgB = 0; flag = false; } } else if (dep == 5) { if (hcal_edepthCorrect6->size() > ml) { - enb = cfac * hcal_edepth6->at(ml); - enu = cfac * hcal_edepthCorrect6->at(ml); - enh = cfac * hcal_edepthHot6->at(ml); - enc = cfac * hcal_edepthHotCorrect6->at(ml); + enb = cfac * hcal_edepth6->at(ml); + enu = cfac * hcal_edepthCorrect6->at(ml); + enh = cfac * hcal_edepthHot6->at(ml); + enc = cfac * hcal_edepthHotCorrect6->at(ml); chgS = hcal_cdepthHot6->at(ml); actL = hcal_activeHotL6->at(ml); chgB = hcal_cdepthHotBG6->at(ml); } else { - enb = enu = enh = enc = chgS = actL = chgB = 0; + enb = enu = enh = enc = chgS = actL = chgB = 0; flag = false; } } else if (dep == 6) { if (hcal_edepthCorrect7->size() > ml) { - enb = cfac * hcal_edepth7->at(ml); - enu = cfac * hcal_edepthCorrect7->at(ml); - enh = cfac * hcal_edepthHot7->at(ml); - enc = cfac * hcal_edepthHotCorrect7->at(ml); + enb = cfac * hcal_edepth7->at(ml); + enu = cfac * hcal_edepthCorrect7->at(ml); + enh = cfac * hcal_edepthHot7->at(ml); + enc = cfac * hcal_edepthHotCorrect7->at(ml); chgS = hcal_cdepthHot7->at(ml); actL = hcal_activeHotL7->at(ml); chgB = hcal_cdepthHotBG7->at(ml); } else { - enb = enu = enh = enc = chgS = actL = chgB = 0; + enb = enu = enh = enc = chgS = actL = chgB = 0; flag = false; } } else { - enb = enu = enh = enc = chgS = actL = chgB = 0; + enb = enu = enh = enc = chgS = actL = chgB = 0; flag = false; } return flag; @@ -1496,98 +1580,93 @@ bool HBHEMuonOfflineAnalyzer::getEnergy(unsigned int ml, int dep, double& enb, bool HBHEMuonOfflineAnalyzer::looseMuon(unsigned int ml) { if (pt_of_muon->at(ml) > 20.) { if (mediumMuon2016(ml)) { - if (IsolationR04->at(ml) < 0.25) { - return true; + if (IsolationR04->at(ml) < 0.25) { + return true; } - } - } - return false; -} + } + } + return false; +} bool HBHEMuonOfflineAnalyzer::softMuon(unsigned int ml) { if (pt_of_muon->at(ml) > 20.) { if (mediumMuon2016(ml)) { - if (IsolationR03->at(ml) < 0.10) { - return true; + if (IsolationR03->at(ml) < 0.10) { + return true; } - } - } - return false; + } + } + return false; } bool HBHEMuonOfflineAnalyzer::tightMuon(unsigned int ml) { if (pt_of_muon->at(ml) > 20.) { if (mediumMuon2016(ml)) { - if (IsolationR04->at(ml) < 0.15) { - return true; + if (IsolationR04->at(ml) < 0.15) { + return true; } - } + } } return false; } bool HBHEMuonOfflineAnalyzer::mediumMuon2016(unsigned int ml) { - bool medium16 = (((PF_Muon->at(ml)) && - (Global_Muon->at(ml) || Tracker_muon->at(ml))) && - (tight_validFraction->at(ml) > 0.49)); - if (!medium16) return medium16; - - bool goodGlob = (Global_Muon->at(ml) && - GlobTrack_Chi->at(ml) < 3 && - muon_chi2LocalPosition->at(ml) < 12 && - muon_trkKink->at(ml) < 20); - medium16 = muon_segComp->at(ml) > (goodGlob ? 0.303 : 0.451); + bool medium16 = + (((PF_Muon->at(ml)) && (Global_Muon->at(ml) || Tracker_muon->at(ml))) && (tight_validFraction->at(ml) > 0.49)); + if (!medium16) + return medium16; + + bool goodGlob = (Global_Muon->at(ml) && GlobTrack_Chi->at(ml) < 3 && muon_chi2LocalPosition->at(ml) < 12 && + muon_trkKink->at(ml) < 20); + medium16 = muon_segComp->at(ml) > (goodGlob ? 0.303 : 0.451); return medium16; } -void HBHEMuonOfflineAnalyzer::etaPhiHcal(unsigned int detId, int &eta, - int &phi, int &depth) { +void HBHEMuonOfflineAnalyzer::etaPhiHcal(unsigned int detId, int &eta, int &phi, int &depth) { int zside, etaAbs; - if ((detId&0x1000000)==0) { - zside = (detId&0x2000)?(1):(-1); - etaAbs = (detId>>7)&0x3F; - phi = detId&0x7F; - depth = (detId>>14)&0x1F; + if ((detId & 0x1000000) == 0) { + zside = (detId & 0x2000) ? (1) : (-1); + etaAbs = (detId >> 7) & 0x3F; + phi = detId & 0x7F; + depth = (detId >> 14) & 0x1F; } else { - zside = (detId&0x80000)?(1):(-1); - etaAbs = (detId>>10)&0x1FF; - phi = detId&0x3FF; - depth = (detId>>20)&0xF; + zside = (detId & 0x80000) ? (1) : (-1); + etaAbs = (detId >> 10) & 0x1FF; + phi = detId & 0x3FF; + depth = (detId >> 20) & 0xF; } - eta = etaAbs*zside; + eta = etaAbs * zside; } -void HBHEMuonOfflineAnalyzer::etaPhiEcal(unsigned int detId, int& type, - int& zside, int& etaX, int& phiY, - int& plane, int& strip) { - - type = ((detId>>25)&0x7); - // std::cout<<"type"<> 25) & 0x7); + // std::cout<<"type"<>9)&0x7F; - phiY = detId&0x1FF; - } else if (type==2) { - zside = (detId&0x4000)?(1):(-1); - etaX = (detId>>7)&0x7F; - phiY = (detId&0x7F); - } else if (type==3) { - zside = (detId&0x80000)?(1):(-1); - etaX = (detId>>6)&0x3F; + zside = (detId & 0x10000) ? (1) : (-1); + etaX = (detId >> 9) & 0x7F; + phiY = detId & 0x1FF; + } else if (type == 2) { + zside = (detId & 0x4000) ? (1) : (-1); + etaX = (detId >> 7) & 0x7F; + phiY = (detId & 0x7F); + } else if (type == 3) { + zside = (detId & 0x80000) ? (1) : (-1); + etaX = (detId >> 6) & 0x3F; /** get the sensor iy */ - phiY = (detId>>12)&0x3F; + phiY = (detId >> 12) & 0x3F; /** get the strip */ - plane = ((detId>>18)&0x1)+1; - strip = detId&0x3F; + plane = ((detId >> 18) & 0x1) + 1; + strip = detId & 0x3F; } else { zside = etaX = phiY = 0; } } -void HBHEMuonOfflineAnalyzer::calculateP(double pt, double eta, double& pM) { - pM = (pt*cos(2*(1/atan(exp(eta))))); +void HBHEMuonOfflineAnalyzer::calculateP(double pt, double eta, double &pM) { + pM = (pt * cos(2 * (1 / atan(exp(eta))))); } void HBHEMuonOfflineAnalyzer::close() { @@ -1601,24 +1680,23 @@ void HBHEMuonOfflineAnalyzer::close() { } void HBHEMuonOfflineAnalyzer::writeHistograms() { - //output_file->cd(); - std::string type[]={"tight"};//,"soft","loose"}; + std::string type[] = {"tight"}; //,"soft","loose"}; char name[128]; - - std::cout<<"WriteHistograms"<Write(); //output_file->cd(); - for (int i=0; iWrite(); h_Eta_Muon[i]->Write(); h_Phi_Muon[i]->Write(); - + h_P_Muon[i]->Write(); h_PF_Muon[i]->Write(); h_GlobTrack_Chi[i]->Write(); h_Global_Muon_Hits[i]->Write(); - + h_MatchedStations[i]->Write(); - h_Tight_TransImpactparameter[i]->Write(); + h_Tight_TransImpactparameter[i]->Write(); h_Tight_LongitudinalImpactparameter[i]->Write(); h_InnerTrackPixelHits[i]->Write(); h_TrackerLayer[i]->Write(); h_IsolationR04[i]->Write(); - + h_Global_Muon[i]->Write(); - h_TransImpactParameter[i]->Write();; + h_TransImpactParameter[i]->Write(); + ; h_TransImpactParameterBin1[i]->Write(); h_TransImpactParameterBin2[i]->Write(); // h_LongImpactParameter[i]->Write(); h_LongImpactParameterBin1[i]->Write(); h_LongImpactParameterBin2[i]->Write(); - + h_ecal_energy[i]->Write(); - h_hcal_energy[i]->Write();; + h_hcal_energy[i]->Write(); + ; h_3x3_ecal[i]->Write(); - h_1x1_hcal[i]->Write();; + h_1x1_hcal[i]->Write(); + ; h_EtaX_hcal[i]->Write(); - h_PhiY_hcal[i]->Write();; - - h_EtaX_ecal[i]->Write();; - h_PhiY_ecal[i]->Write();; - h_Eta_ecal[i]->Write();; - h_Phi_ecal[i]->Write();; - h_MuonHittingEcal[i]->Write();; - h_HotCell[i]->Write();; + h_PhiY_hcal[i]->Write(); + ; + + h_EtaX_ecal[i]->Write(); + ; + h_PhiY_ecal[i]->Write(); + ; + h_Eta_ecal[i]->Write(); + ; + h_Phi_ecal[i]->Write(); + ; + h_MuonHittingEcal[i]->Write(); + ; + h_HotCell[i]->Write(); + ; output_file->cd(); - for (int eta=etaMin_; eta<=etaMax_; ++eta) { - int nDepth = nDepthBins(eta,-1); - int nPhi = nPhiBins(eta); - sprintf(name, "Dir_muon_type_%s_ieta%d",type[i].c_str(), eta); - d_output_file[i][eta-1]= output_file->mkdir(name); + for (int eta = etaMin_; eta <= etaMax_; ++eta) { + int nDepth = nDepthBins(eta, -1); + int nPhi = nPhiBins(eta); + sprintf(name, "Dir_muon_type_%s_ieta%d", type[i].c_str(), eta); + d_output_file[i][eta - 1] = output_file->mkdir(name); //output_file->cd(name); - d_output_file[i][eta-1]->cd(); - for (int depth=0; depthWrite(); - - h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[i][ih]->Write(); - h_active_length_Fill[i][ih]->Write(); - h_p_muon_ineta[i][ih]->Write(); - h_charge_signal[i][ih]->Write(); - h_charge_bg[i][ih]->Write(); - ih++; - h_Hot_MuonEnergy_hcal_HotCell[i][ih]->Write(); - h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[i][ih]->Write(); - h_active_length_Fill[i][ih]->Write(); - h_p_muon_ineta[i][ih]->Write(); - h_charge_signal[i][ih]->Write(); - h_charge_bg[i][ih]->Write(); - } + d_output_file[i][eta - 1]->cd(); + for (int depth = 0; depth < nDepth; ++depth) { + for (int PHI = 0; PHI < nPhi; ++PHI) { + // std::cout<<"eta:"<Write(); + + h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[i][ih]->Write(); + h_active_length_Fill[i][ih]->Write(); + h_p_muon_ineta[i][ih]->Write(); + h_charge_signal[i][ih]->Write(); + h_charge_bg[i][ih]->Write(); + ih++; + h_Hot_MuonEnergy_hcal_HotCell[i][ih]->Write(); + h_Hot_MuonEnergy_hcal_HotCell_VsActiveLength[i][ih]->Write(); + h_active_length_Fill[i][ih]->Write(); + h_p_muon_ineta[i][ih]->Write(); + h_charge_signal[i][ih]->Write(); + h_charge_bg[i][ih]->Write(); + } } output_file->cd(); - } + } } output_file->cd(); } int HBHEMuonOfflineAnalyzer::nDepthBins(int eta, int phi) { // Run 1 scenario - int nDepthR1[29]={1,1,1,1,1,1,1,1,1,1,1,1,1,1,2,3,1,2,2,2,2,2,2,2,2,2,3,3,2}; + int nDepthR1[29] = {1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 3, 1, 2, 2, 2, 2, 2, 2, 2, 2, 2, 3, 3, 2}; // Run 2 scenario from 2018 - int nDepthR2[29]={1,1,1,1,1,1,1,1,1,1,1,1,1,1,2,4,3,5,6,6,6,6,6,6,6,7,7,7,3}; + int nDepthR2[29] = {1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 5, 6, 6, 6, 6, 6, 6, 6, 7, 7, 7, 3}; // Run 3 scenario - int nDepthR3[29]={4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,3,5,6,6,6,6,6,6,6,7,7,7,3}; + int nDepthR3[29] = {4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 3, 5, 6, 6, 6, 6, 6, 6, 6, 7, 7, 7, 3}; // Run 4 scenario - int nDepthR4[29]={4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,7,7,7,7,7,7,7,7,7,7,7,7,7}; + int nDepthR4[29] = {4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7}; // for a test scenario with multi depth segmentation considered during Run 1 // int nDepth[29]={3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,5,5,5,5,5,5,5,5,5,5,5,5,5}; // modeLHC_ = 0 --> nbin defined maxDepthHB/HE @@ -1722,9 +1808,9 @@ int HBHEMuonOfflineAnalyzer::nDepthBins(int eta, int phi) { // = 3 --> corresponds to Run 3 (post LS2) // = 4 --> corresponds to 2017 (Plan 1) // = 5 --> corresponds to Run 4 (post LS3) - int nbin(0); + int nbin(0); if (modeLHC_ == 0) { - if (eta<=15) { + if (eta <= 15) { nbin = maxDepthHB_; } else if (eta == 16) { nbin = 4; @@ -1732,28 +1818,28 @@ int HBHEMuonOfflineAnalyzer::nDepthBins(int eta, int phi) { nbin = maxDepthHE_; } } else if (modeLHC_ == 1) { - nbin = nDepthR1[eta-1]; + nbin = nDepthR1[eta - 1]; } else if (modeLHC_ == 2) { - nbin = nDepthR2[eta-1]; + nbin = nDepthR2[eta - 1]; } else if (modeLHC_ == 3) { - nbin = nDepthR3[eta-1]; + nbin = nDepthR3[eta - 1]; } else if (modeLHC_ == 4) { if (phi > 0) { if (eta >= 16 && phi >= 63 && phi <= 66) { - nbin = nDepthR2[eta-1]; + nbin = nDepthR2[eta - 1]; } else { - nbin = nDepthR1[eta-1]; + nbin = nDepthR1[eta - 1]; } } else { if (eta >= 16) { - nbin = (nDepthR2[eta-1] > nDepthR1[eta-1]) ? nDepthR2[eta-1] : nDepthR1[eta-1]; + nbin = (nDepthR2[eta - 1] > nDepthR1[eta - 1]) ? nDepthR2[eta - 1] : nDepthR1[eta - 1]; } else { - nbin = nDepthR1[eta-1]; + nbin = nDepthR1[eta - 1]; } } } else { if (eta > 0 && eta < 30) { - nbin = nDepthR4[eta-1]; + nbin = nDepthR4[eta - 1]; } else { nbin = nDepthR4[28]; } @@ -1763,17 +1849,19 @@ int HBHEMuonOfflineAnalyzer::nDepthBins(int eta, int phi) { int HBHEMuonOfflineAnalyzer::nPhiBins(int eta) { int nphi = (eta <= 20) ? 72 : 36; - if (modeLHC_ == 5 && eta > 16) nphi = 360; + if (modeLHC_ == 5 && eta > 16) + nphi = 360; return nphi; } float HBHEMuonOfflineAnalyzer::getCorr(int run, unsigned int id) { float cfac(1.0); int ip(-1); - for (unsigned int k=0; k= runlow_[i]) { - ip = (int)(i); break; + ip = (int)(i); + break; } } if (debug_) { @@ -1781,92 +1869,94 @@ float HBHEMuonOfflineAnalyzer::getCorr(int run, unsigned int id) { } unsigned idx = correctDetId(id); if (ip >= 0) { - std::map::iterator itr = corrFac_[ip].find(idx); - if (itr != corrFac_[ip].end()) cfac = itr->second; + std::map::iterator itr = corrFac_[ip].find(idx); + if (itr != corrFac_[ip].end()) + cfac = itr->second; } if (debug_) { - int subdet, zside, ieta, iphi, depth; + int subdet, zside, ieta, iphi, depth; unpackDetId(idx, subdet, zside, ieta, iphi, depth); - std::cout << "ID " << std::hex << id << std::dec << " (Sub " << subdet - << " eta " << zside*ieta << " phi " << iphi << " depth " << depth - << ") Factor " << cfac << std::endl; + std::cout << "ID " << std::hex << id << std::dec << " (Sub " << subdet << " eta " << zside * ieta << " phi " << iphi + << " depth " << depth << ") Factor " << cfac << std::endl; } return cfac; } -std::vector HBHEMuonOfflineAnalyzer::splitString(const std::string& fLine) { - std::vector result; +std::vector HBHEMuonOfflineAnalyzer::splitString(const std::string &fLine) { + std::vector result; int start = 0; bool empty = true; - for (unsigned i = 0; i <= fLine.size (); i++) { - if (fLine [i] == ' ' || i == fLine.size ()) { + for (unsigned i = 0; i <= fLine.size(); i++) { + if (fLine[i] == ' ' || i == fLine.size()) { if (!empty) { - std::string item (fLine, start, i-start); - result.push_back (item); - empty = true; + std::string item(fLine, start, i - start); + result.push_back(item); + empty = true; } - start = i+1; + start = i + 1; } else { - if (empty) empty = false; + if (empty) + empty = false; } } return result; } -unsigned int HBHEMuonOfflineAnalyzer::getDetIdHBHE(int ieta, int iphi, - int depth) { - int eta = std::abs(ieta); +unsigned int HBHEMuonOfflineAnalyzer::getDetIdHBHE(int ieta, int iphi, int depth) { + int eta = std::abs(ieta); int subdet(0); - if (eta > 16) subdet = 2; - else if (eta == 16 && depth > 2) subdet = 2; - else subdet = 1; - return getDetId(subdet,ieta,iphi,depth); + if (eta > 16) + subdet = 2; + else if (eta == 16 && depth > 2) + subdet = 2; + else + subdet = 1; + return getDetId(subdet, ieta, iphi, depth); } -unsigned int HBHEMuonOfflineAnalyzer::getDetId(int subdet, int ieta, int iphi, - int depth) { - // All numbers used here are described as masks/offsets in +unsigned int HBHEMuonOfflineAnalyzer::getDetId(int subdet, int ieta, int iphi, int depth) { + // All numbers used here are described as masks/offsets in // DataFormats/HcalDetId/interface/HcalDetId.h - unsigned int id = ((4<<28)|((subdet&0x7)<<25)); - id |= ((0x1000000) | ((depth&0xF)<<20) | - ((ieta>0)?(0x80000|(ieta<<10)):((-ieta)<<10)) | - (iphi&0x3FF)); + unsigned int id = ((4 << 28) | ((subdet & 0x7) << 25)); + id |= ((0x1000000) | ((depth & 0xF) << 20) | ((ieta > 0) ? (0x80000 | (ieta << 10)) : ((-ieta) << 10)) | + (iphi & 0x3FF)); return id; } -unsigned int HBHEMuonOfflineAnalyzer::correctDetId(const unsigned int & detId) { +unsigned int HBHEMuonOfflineAnalyzer::correctDetId(const unsigned int &detId) { int subdet, ieta, zside, depth, iphi; unpackDetId(detId, subdet, zside, ieta, iphi, depth); if (subdet == 0) { - if (ieta > 16) subdet = 2; - else if (ieta == 16 && depth > 2) subdet = 2; - else subdet = 1; + if (ieta > 16) + subdet = 2; + else if (ieta == 16 && depth > 2) + subdet = 2; + else + subdet = 1; } - unsigned int id = getDetId(subdet,ieta*zside,iphi,depth); + unsigned int id = getDetId(subdet, ieta * zside, iphi, depth); if ((id != detId) && debug_) { - std::cout << "Correct Id " << std::hex << detId << " to " << id << std::dec - << "(Sub " << subdet << " eta " << ieta*zside << " phi " << iphi - << " depth " << depth << ")" << std::endl; + std::cout << "Correct Id " << std::hex << detId << " to " << id << std::dec << "(Sub " << subdet << " eta " + << ieta * zside << " phi " << iphi << " depth " << depth << ")" << std::endl; } return id; } -void HBHEMuonOfflineAnalyzer::unpackDetId(unsigned int detId, int& subdet, - int& zside, int& ieta, int& iphi, - int& depth) { +void HBHEMuonOfflineAnalyzer::unpackDetId( + unsigned int detId, int &subdet, int &zside, int &ieta, int &iphi, int &depth) { // The maskings are defined in DataFormats/DetId/interface/DetId.h // and in DataFormats/HcalDetId/interface/HcalDetId.h // The macro does not invoke the classes there and use them subdet = ((detId >> 25) & (0x7)); - if ((detId&0x1000000) == 0) { - ieta = ((detId >> 7) & 0x3F); - zside = (detId&0x2000)?(1):(-1); - depth = ((detId >> 14) & 0x1F); - iphi = (detId & 0x3F); + if ((detId & 0x1000000) == 0) { + ieta = ((detId >> 7) & 0x3F); + zside = (detId & 0x2000) ? (1) : (-1); + depth = ((detId >> 14) & 0x1F); + iphi = (detId & 0x3F); } else { - ieta = ((detId >> 10) & 0x1FF); - zside = (detId&0x80000)?(1):(-1); - depth = ((detId >> 20) & 0xF); - iphi = (detId & 0x3FF); + ieta = ((detId >> 10) & 0x1FF); + zside = (detId & 0x80000) ? (1) : (-1); + depth = ((detId >> 20) & 0xF); + iphi = (detId & 0x3FF); } } From 721a15c3a94a22bffb78ca74e7285e530fed3c20 Mon Sep 17 00:00:00 2001 From: Sunanda Date: Mon, 27 Jul 2020 18:02:21 +0200 Subject: [PATCH 3/3] Change default for 2018 runs --- Calibration/HcalCalibAlgos/macros/AnalyzeLepTree.C | 14 +++++++------- 1 file changed, 7 insertions(+), 7 deletions(-) diff --git a/Calibration/HcalCalibAlgos/macros/AnalyzeLepTree.C b/Calibration/HcalCalibAlgos/macros/AnalyzeLepTree.C index c217a18f7aafd..cb0f0db73c1f4 100644 --- a/Calibration/HcalCalibAlgos/macros/AnalyzeLepTree.C +++ b/Calibration/HcalCalibAlgos/macros/AnalyzeLepTree.C @@ -23,12 +23,12 @@ // 7-11: RBX # to be excluded (maximum 5 bits needed for RBX) // 12: 0 varying ranges of p depending on |ieta| // 1 constant p ranges -// 13-15: 0 no cut on ediff; 1-4 cuts at 5, 10, 15, 20 GeV +// 13-15: 0 no cut on ediff; 1-14 cuts at 5, 10, 15, 20 GeV // modeLHC (integer) specifies the detector condition // 0 Run1 detector (till 2016) -// 1 Plan36 detector -// 2 Phase1 detector -// 3 Plan1 detector +// 1 Plan36 detector (2018) +// 2 Phase1 detector (Run3) +// 3 Plan1 detector (2017) // 4 Phase2 detector // // AnalyzeLepTree a1(tree, mode, modeLHC); @@ -82,8 +82,8 @@ class AnalyzeLepTree { public: - AnalyzeLepTree(TChain* tree, int mode = 0, int modeLHC = 3); - AnalyzeLepTree(const char* fname, int mode = 0, int modeLHC = 3); + AnalyzeLepTree(TChain* tree, int mode = 0, int modeLHC = 1); + AnalyzeLepTree(const char* fname, int mode = 0, int modeLHC = 1); virtual ~AnalyzeLepTree(); virtual Int_t Cut(Long64_t entry); virtual Int_t GetEntry(Long64_t entry); @@ -262,7 +262,7 @@ void AnalyzeLepTree::Init(TChain* tree) { else if (modeLHC_ == 3) std::cout << "This is Plan1 detector (2017)\n"; else - std::cout << "This is Phase2 detector (after 2024)\n"; + std::cout << "This is Phase2 detector (after 2026)\n"; static const double cuts[8] = {200, 5, 10, 15, 20, 25, 30, 40}; int cutE = (mode_ / 4096) % 8; cutEdiff_ = cuts[cutE];