ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/HbbAnalysis/src/HistosElecsBase.cc
Revision: 1.3
Committed: Tue Feb 9 14:52:25 2010 UTC (15 years, 2 months ago) by amagnan
Content type: text/plain
Branch: MAIN
Changes since 1.2: +1 -1 lines
Log Message:
export code to CMSSW_3_4_X

File Contents

# User Rev Content
1 amagnan 1.1 #include <iostream>
2     #include <fstream>
3     #include <sstream>
4     #include <cmath>
5     #include <algorithm>
6    
7     #include "UserCode/HbbAnalysis/interface/HistosElecsBase.hh"
8    
9     namespace HbbAnalysis {//namespace
10    
11     void HistosElecsBase::Initialise(TFileDirectory & aDir, std::string aName, bool aDoGenMatched){
12    
13     doGenMatched_ = aDoGenMatched;
14     CreateHistos(aName,aDir);
15    
16     p_recoOverGen = aDir.make<TH1F>("p_recoOverGen",";p_{T}^{reco}/p_{T}^{gen};N_{entries}",300,0,3);
17    
18     for (unsigned int bin(0); bin<20; bin++){//loop on eta bins
19     std::ostringstream retatitle;
20     retatitle << "p_recoOverGen_vsEta_" << bin ;
21     p_recoOverGen_vsEta[bin] = aDir.make<TH1F>(retatitle.str().c_str(),";p_{T}^{reco}/p_{T}^{gen};N_{entries}",300,0,3);
22     }//loop on eta bins
23    
24     for (unsigned int bin(0); bin<10; bin++){//loop on pt bins
25     std::ostringstream rpttitle;
26     rpttitle << "p_recoOverGen_vsPt_" << bin;
27     p_recoOverGen_vsPt[bin] = aDir.make<TH1F>(rpttitle.str().c_str(),";p_{T}^{reco}/p_{T}^{gen};N_{entries}",300,0,3);
28     }//loop on pt bins
29    
30    
31     p_nElectrons = aDir.make<TH1F>("p_nElectrons",";N_{electrons};N_{entries}",10,0,10);
32    
33     p_electronID = aDir.make<TH1F>("p_electronID",";electronID;N_{entries}",8,0,8);
34    
35     p_scSigmaEtaEta = aDir.make<TH1F>("p_scSigmaEtaEta",";#sigma_{#eta,#eta};N_{entries}",500,0,0.1);
36     p_scSigmaIEtaIEta = aDir.make<TH1F>("p_scSigmaIEtaIEta",";#sigma_{i#eta,i#eta};N_{entries}",500,0,0.1);
37     p_scE1x5 = aDir.make<TH1F>("p_scE1x5",";E_{sc}^{1#times 5} (GeV);N_{entries}",500,0,500);
38     p_scE2x5Max = aDir.make<TH1F>("p_scE2x5Max",";E_{sc}^{2#times 5} (GeV);N_{entries}",500,0,500);
39     p_scE5x5 = aDir.make<TH1F>("p_scE5x5",";E_{sc}^{5#times 5} (GeV);N_{entries}",500,0,500);
40     p_scE1x5OverscE5x5 = aDir.make<TH1F>("p_scE1x5OverscE5x5",";E_{sc}^{1#times 5}/E_{sc}^{5#times 5};N_{entries}",500,0,2);
41     p_scE2x5MaxOverscE5x5 = aDir.make<TH1F>("p_scE2x5MaxOverscE5x5",";E_{sc}^{2#times 5}/E_{sc}^{5#times 5};N_{entries}",500,0,2);
42    
43     p_eSuperClusterOverP = aDir.make<TH1F>("p_eSuperClusterOverP",";E_{sc}/p_{track};N_{entries}",100,0.5,2);
44    
45     p_HoverE = aDir.make<TH1F>("p_HoverE",";E_{had}/E_{em};N_{entries}",150,0,0.15);
46     p_deltaPhiIn = aDir.make<TH1F>("p_deltaPhiIn",";#Delta#phi_{in};N_{entries}",150,0,0.15);
47     p_deltaEtaIn = aDir.make<TH1F>("p_deltaEtaIn",";#Delta#eta_{in};N_{entries}",150,0,0.03);
48    
49     p_caloIso = aDir.make<TH1F>("p_caloIso",";caloIso;N_{entries}",400,0,20);
50     p_hcalIso = aDir.make<TH1F>("p_hcalIso",";hcalIso;N_{entries}",400,0,20);
51    
52     p_trackIso = aDir.make<TH1F>("p_trackIso",";trackIso (GeV);N_{entries}",400,0,20);
53     p_trackIsoOverEt = aDir.make<TH1F>("p_trackIsoOverEt",";trackIso/E_{T};N_{entries}",400,0,1);
54     p_ecalIso = aDir.make<TH1F>("p_ecalIso",";ecalIso (GeV);N_{entries}",400,0,20);
55     p_ecalIsoOverEt = aDir.make<TH1F>("p_ecalIsoOverEt",";ecalIso/E_{T};N_{entries}",400,0,1);
56     p_combIso = aDir.make<TH1F>("p_combIso",";combIso (GeV);N_{entries}",400,0,40);
57     p_combIsoOverEt = aDir.make<TH1F>("p_combIsoOverEt",";combIso/E_{T};N_{entries}",400,0,1);
58    
59     p_gsfTrk_pT = aDir.make<TH1F>("p_gsfTrk_pT",";p_{T}^{trk} (GeV);N_{entries}",200,0,200);
60     p_gsfTrk_IPxy = aDir.make<TH1F>("p_gsfTrk_IPxy",";IP_{xy}^{trk};N_{entries}",100,-0.1,0.1);
61     p_gsfTrk_IPz = aDir.make<TH1F>("p_gsfTrk_IPz",";IP_{z}^{trk};N_{entries}",100,-1,1);
62    
63     //efficiency curves vs cut value
64     peff_eIso[0] = aDir.make<TH1F>("peff_eIso_trkIso",";trkIso (GeV);N_{sel}/N_{tot}",50,0,5);
65     peff_eIso[1] = aDir.make<TH1F>("peff_eIso_trkIsoOverEt",";trkIso/E_{T};N_{sel}/N_{tot}",50,0,1);
66     peff_eIso[2] = aDir.make<TH1F>("peff_eIso_ecalIso",";ecalIso (GeV);N_{sel}/N_{tot}",50,0,5);
67     peff_eIso[3] = aDir.make<TH1F>("peff_eIso_ecalIsoOverEt",";ecalIso/E_{T};N_{sel}/N_{tot}",50,0,1);
68     peff_eIso[4] = aDir.make<TH1F>("peff_eIso_combIso",";combIso (GeV);N_{sel}/N_{tot}",50,0,10);
69     peff_eIso[5] = aDir.make<TH1F>("peff_eIso_combIsoOverEt",";combIso/E_{T};N_{sel}/N_{tot}",50,0,1);
70    
71     //efficiency curves vs pT muon
72     peff_eEt[0] = aDir.make<TH1F>("peff_eEt_trkIso",";E_{T} (GeV);N_{sel}/N_{tot}",18,10,100);
73     const unsigned int lNBinsPt = peff_eEt[0]->GetNbinsX();
74     const double lMinPt = peff_eEt[0]->GetXaxis()->GetXmin();
75     const double lMaxPt = peff_eEt[0]->GetXaxis()->GetXmax();
76     peff_eEt[1] = aDir.make<TH1F>("peff_eEt_trkIsoOverEt",";E_{T} (GeV);N_{sel}/N_{tot}",lNBinsPt,lMinPt,lMaxPt);
77     peff_eEt[2] = aDir.make<TH1F>("peff_eEt_ecalIso",";E_{T} (GeV);N_{sel}/N_{tot}",lNBinsPt,lMinPt,lMaxPt);
78     peff_eEt[3] = aDir.make<TH1F>("peff_eEt_ecalIsoOverEt",";E_{T} (GeV);N_{sel}/N_{tot}",lNBinsPt,lMinPt,lMaxPt);
79     peff_eEt[4] = aDir.make<TH1F>("peff_eEt_combIso",";E_{T} (GeV);N_{sel}/N_{tot}",lNBinsPt,lMinPt,lMaxPt);
80     peff_eEt[5] = aDir.make<TH1F>("peff_eEt_combIsoOverEt",";E_{T} (GeV);N_{sel}/N_{tot}",lNBinsPt,lMinPt,lMaxPt);
81    
82    
83     //efficiency curves vs eta muon
84     peff_eEta[0] = aDir.make<TH1F>("peff_eEta_trkIso",";#eta;N_{sel}/N_{tot}",42,-2.1,2.1);
85     const unsigned int lNBinsEta = peff_eEta[0]->GetNbinsX();
86     const double lMinEta = peff_eEta[0]->GetXaxis()->GetXmin();
87     const double lMaxEta = peff_eEta[0]->GetXaxis()->GetXmax();
88     peff_eEta[1] = aDir.make<TH1F>("peff_eEta_trkIsoOverEt",";#eta;N_{sel}/N_{tot}",lNBinsEta,lMinEta,lMaxEta);
89     peff_eEta[2] = aDir.make<TH1F>("peff_eEta_ecalIso",";#eta;N_{sel}/N_{tot}",lNBinsEta,lMinEta,lMaxEta);
90     peff_eEta[3] = aDir.make<TH1F>("peff_eEta_ecalIsoOverEt",";#eta;N_{sel}/N_{tot}",lNBinsEta,lMinEta,lMaxEta);
91     peff_eEta[4] = aDir.make<TH1F>("peff_eEta_combIso",";#eta;N_{sel}/N_{tot}",lNBinsEta,lMinEta,lMaxEta);
92     peff_eEta[5] = aDir.make<TH1F>("peff_eEta_combIsoOverEt",";#eta;N_{sel}/N_{tot}",lNBinsEta,lMinEta,lMaxEta);
93    
94     //efficiency calculation tool
95     for (unsigned int i(0); i<6; i++){
96     const unsigned int lNBins = peff_eIso[i]->GetNbinsX();
97     const double lMin = peff_eIso[i]->GetXaxis()->GetXmin();
98     const double lMax = peff_eIso[i]->GetXaxis()->GetXmax();
99     isoEff_[i].initialise(lNBins,lMin,lMax);
100     isoEffEt_[i].initialise(lNBinsPt,lMinPt,lMaxPt);
101     isoEffEta_[i].initialise(lNBinsEta,lMinEta,lMaxEta);
102     }
103    
104    
105     }
106    
107     void HistosElecsBase::FillEventHistograms(const std::vector<HbbAnalysis::Electron> & aCol){
108    
109     if (aCol.size() > 0) {
110     if (doGenMatched_) {
111     unsigned int nMatched = 0;
112     for (std::vector<HbbAnalysis::Electron>::const_iterator iElec = aCol.begin();
113     iElec != aCol.end();
114     iElec++)
115     {
116     if (MatchesGenElectron(*iElec)) nMatched++;
117     }
118     p_nElectrons->Fill(nMatched);
119     }
120     else p_nElectrons->Fill(aCol.size());
121    
122     }
123     }
124    
125     void HistosElecsBase::FillHistograms(const HbbAnalysis::Electron & aElec, bool isLead){//FillHistograms
126    
127     bool lIsGenMatched = !doGenMatched_ || (doGenMatched_ && MatchesGenElectron(aElec));
128     if (lIsGenMatched) {//genMatched
129     if (isLead) {//isLead
130     FillBaseHistograms(aElec.recoVars().pT,aElec.recoVars().eta,aElec.recoVars().phi,aElec.recoVars().charge);
131    
132     // if ( aElec.gsfTrack().isAvailable() && !aElec.gsfTrack().isNull() ) {
133     // p_gsfTrk_pT->Fill(aElec.gsfTrack()->pt());
134     // if ( aRecoVertices->size() >= 1 ) {
135     // const reco::Vertex & thePrimaryEventVertex = (*aRecoVertices->begin());
136     // p_gsfTrk_IPxy->Fill(aElec.gsfTrack()->dxy(thePrimaryEventVertex.position()));
137     // p_gsfTrk_IPz->Fill(aElec.gsfTrack()->dz(thePrimaryEventVertex.position()));
138     // }
139     // }
140    
141    
142     if (aElec.genVars().valid && fabs(aElec.genVars().eta)<2.1 && aElec.genVars().pT>10) {
143     int binEta = static_cast<int>((aElec.genVars().eta+2.1)/0.21);
144     int binpt;
145     if (aElec.genVars().pT < 100) binpt = static_cast<int>((aElec.genVars().pT-10)/11.25);
146     else if (aElec.genVars().pT < 120) binpt = 8;
147     else binpt = 9;
148    
149     p_recoOverGen_vsEta[binEta]->Fill(aElec.recoVars().pT/aElec.genVars().pT);
150     p_recoOverGen_vsPt[binpt]->Fill(aElec.recoVars().pT/aElec.genVars().pT);
151 amagnan 1.2 p_recoOverGen->Fill(aElec.recoVars().pT/aElec.genVars().pT);
152 amagnan 1.1 }
153    
154    
155     std::vector<std::pair<std::string,float> > lIds = aElec.idVars().electronIDs;
156     if (lIds.size() > 5) std::cout << "--- WARNING: histo will be out-of-range, number of electron IDs = " << lIds.size() << std::endl;
157     for (unsigned int iId(0); iId<lIds.size(); iId++){
158 amagnan 1.3 //const std::string & lName = lIds.at(iId).first;
159 amagnan 1.1 //std::cout << "--- ElectronID : Bin " << iId << ", id name = " << lName << std::endl;
160     p_electronID->Fill(2*iId+lIds.at(iId).second);
161     }
162    
163     p_scSigmaEtaEta->Fill(aElec.scVars().sigmaEtaEta);
164     p_scSigmaIEtaIEta->Fill(aElec.scVars().sigmaIEtaIEta);
165     p_scE1x5->Fill(aElec.scVars().e1x5);
166     p_scE2x5Max->Fill(aElec.scVars().e2x5Max);
167     p_scE5x5->Fill(aElec.scVars().e5x5);
168    
169     p_eSuperClusterOverP->Fill(aElec.scVars().eOverP);
170    
171     p_HoverE->Fill(aElec.idVars().hOverE);
172     p_deltaPhiIn->Fill(aElec.idVars().deltaPhiIn);
173     p_deltaEtaIn->Fill(aElec.idVars().deltaEtaIn);
174    
175     if (aElec.scVars().e5x5 > 0) {
176     p_scE1x5OverscE5x5->Fill(aElec.scVars().e1x5/aElec.scVars().e5x5);
177     p_scE2x5MaxOverscE5x5->Fill(aElec.scVars().e2x5Max/aElec.scVars().e5x5);
178     }
179    
180     p_caloIso->Fill(aElec.isoVars().calo);
181     p_hcalIso->Fill(aElec.isoVars().hcal);
182    
183     double lIsoVar[6] = {
184     aElec.isoVars().track,
185     aElec.isoVars().track/aElec.recoVars().pT,
186     aElec.isoVars().ecal,
187     aElec.isoVars().ecal/aElec.recoVars().pT,
188     aElec.isoVars().track+aElec.isoVars().ecal+aElec.isoVars().hcal,
189     (aElec.isoVars().track+aElec.isoVars().ecal+aElec.isoVars().hcal)/aElec.recoVars().pT
190     };
191    
192     p_trackIso->Fill(lIsoVar[0]);
193     p_trackIsoOverEt->Fill(lIsoVar[1]);
194     p_ecalIso->Fill(lIsoVar[2]);
195     p_ecalIsoOverEt->Fill(lIsoVar[3]);
196     p_combIso->Fill(lIsoVar[4]);
197     p_combIsoOverEt->Fill(lIsoVar[5]);
198    
199 amagnan 1.2 double lIsoCut[6] = {2.4,0.05,2.4,0.09,3.,0.1};
200 amagnan 1.1
201     for (unsigned int i(0); i<6; i++){//loop on iso variables
202    
203     for (unsigned int lBin(0); lBin<isoEff_[i].numberOfBins(); lBin++){
204     double lCut = isoEff_[i].xMin()+isoEff_[i].stepSize()*lBin;
205    
206     isoEff_[i].incrementTotal(lBin);
207     if (lIsoVar[i] <= lCut){
208     isoEff_[i].incrementPass(lBin);
209     }
210     }
211    
212     //eff vs pT
213     for (unsigned int lBin(0); lBin<isoEffEt_[i].numberOfBins(); lBin++){
214     double lCutMin = isoEffEt_[i].xMin()+isoEffEt_[i].stepSize()*lBin;
215     double lCutMax = isoEffEt_[i].xMin()+isoEffEt_[i].stepSize()*(lBin+1);
216    
217     if (aElec.recoVars().pT >= lCutMin && aElec.recoVars().pT < lCutMax){
218     isoEffEt_[i].incrementTotal(lBin);
219    
220     if (lIsoVar[i] <= lIsoCut[i]){
221     isoEffEt_[i].incrementPass(lBin);
222     }
223     }
224     }
225    
226     //eff vs eta
227     for (unsigned int lBin(0); lBin<isoEffEta_[i].numberOfBins(); lBin++){
228     double lCutMin = isoEffEta_[i].xMin()+isoEffEta_[i].stepSize()*lBin;
229     double lCutMax = isoEffEta_[i].xMin()+isoEffEta_[i].stepSize()*(lBin+1);
230    
231     if (aElec.recoVars().eta >= lCutMin && aElec.recoVars().eta < lCutMax){
232     isoEffEta_[i].incrementTotal(lBin);
233    
234     if (lIsoVar[i] <= lIsoCut[i]){
235     isoEffEta_[i].incrementPass(lBin);
236     }
237     }
238     }
239    
240    
241     }//loop on iso variables
242    
243    
244    
245     }//isLead
246     }//genMatched
247    
248     }//FillHistograms
249    
250    
251    
252     void HistosElecsBase::FillEffHistograms(){//FillEffHistograms
253    
254     for (unsigned int i(0); i<6; i++){//loop on iso variables
255    
256     for (unsigned int lBin(0); lBin<isoEff_[i].numberOfBins(); lBin++){
257     peff_eIso[i]->SetBinContent(lBin+1,isoEff_[i].getRatio(lBin));
258     peff_eIso[i]->SetBinError(lBin+1,isoEff_[i].getRatioError(lBin));
259     }
260     //eff vs pT
261     for (unsigned int lBin(0); lBin<isoEffEt_[i].numberOfBins(); lBin++){
262     peff_eEt[i]->SetBinContent(lBin+1,isoEffEt_[i].getRatio(lBin));
263     peff_eEt[i]->SetBinError(lBin+1,isoEffEt_[i].getRatioError(lBin));
264     }
265    
266     //eff vs eta
267     for (unsigned int lBin(0); lBin<isoEffEta_[i].numberOfBins(); lBin++){
268     peff_eEta[i]->SetBinContent(lBin+1,isoEffEta_[i].getRatio(lBin));
269     peff_eEta[i]->SetBinError(lBin+1,isoEffEta_[i].getRatioError(lBin));
270     }
271    
272     }//loop on iso variables
273    
274     }//FillEffHistograms
275    
276     bool HistosElecsBase::MatchesGenElectron(const HbbAnalysis::Electron& aElec)
277     {
278     return (aElec.genVars().valid && abs(aElec.genVars().pdgId) == 11);
279     }
280    
281     }//namespace
282    
283