ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/Hammad/macros/EfficiencyPlot.h
Revision: 1.1
Committed: Mon Aug 31 11:53:54 2009 UTC (15 years, 8 months ago) by ghammad
Content type: text/plain
Branch: MAIN
CVS Tags: HEAD
Error occurred while calculating annotation data.
Log Message:
personal code : first version

File Contents

# Content
1 #include "/home/gregory/Analysis/TtSemiLeptonic/tdrstyle.C"
2 #include "map.h"
3
4 void EfficiencyPlot(){
5
6 //using namespace std;
7
8 //
9 // load framework libraries
10 //
11 //gSystem->Load( "libFWCoreFWLite" );
12 //AutoLibraryLoader::enable();
13 //
14 // Set P-TDR plots style
15 //
16
17 setTDRStyle();
18
19 //TFile *InputFile = TFile::Open("../rootfiles/multijets_EvtSelect_NonIso_4jetsPt30_JetComb_MediumBtagging_obs.root");
20 //TFile *InputFile = TFile::Open("../rootfiles/multijets_EvtSelect_NonIso_VetoIso_4jetsPt30_JetComb_MediumBtagging_obs.root");
21 TFile *InputFile = TFile::Open("../rootfiles/allprocesses_EvtSelect_NonIso_VetoIso_4jetsPt30_JetComb_MediumBtagging_obs.root");
22
23 //TString ControlHistoFile = "multijets_EvtSelect_NonIso_4jetsPt30_JetComb_MediumBtagging_eff.ps";
24 //TString ControlHistoFile = "multijets_EvtSelect_NonIso_VetoIso_4jetsPt30_JetComb_MediumBtagging_eff.ps";
25 TString ControlHistoFile = "allprocesses_EvtSelect_NonIso_VetoIso_4jetsPt30_JetComb_MediumBtagging_eff.ps";
26
27 //TFile *OutputFile = new TFile("multijets_EvtSelect_NonIso_VetoIso_4jetsPt30_JetComb_MediumBtagging_eff.root","RECREATE");
28 TFile *OutputFile = new TFile("allprocesses_EvtSelect_NonIso_VetoIso_4jetsPt30_JetComb_MediumBtagging_eff.root","RECREATE");
29
30 bool Print = false;
31
32 const int NbOfObs = 19;
33
34 TString Obs[NbOfObs]={"p_{T}^{#mu}[GeV/c]","H_{T}[GeV]","H_{T}[GeV]","CaloIso[GeV]","TrackIso[GeV/c]","CaloIso/p_{T}","TrackIso/p_{T}","RelIso","|d_{0}|/#sigma","M^{l#nu}_{T}[GeV/c^{2}]","p_{T}^{rel}[GeV/c]","#DeltaR(#mu,jet)","Nb of jets (p_{T}>30GeV/c)","#eta^{#mu}","Nb of b-tagged jets ","m_{qq'} [GeV/c^{2}]","|#eta_{lead. jet}|","#Chi^{2}","M^{bl#nu}_{T}[GeV/c^{2}]"};
35
36 TString ObsName[NbOfObs]={"MuonPt","HT","HTall","CaloIso","TrackIso","RelCaloIso","RelTrackIso","RelIso","d0Sig","MT","PtRel","DRMuJet","Njets","MuonEta","Nbjets","Mqq","EtaLeadJet","Chi2","LeptMT"};
37
38 TString EffName[NbOfObs]={"N_{p_{T}^{#mu}<","N_{H_{T}<","N_{H_{T}^{all}<","N_{CaloIso<","N_{TrackIso<","N_{RelCaloIso<","N_{RelTrackIso<","N_{RelIso<","N_{|d0|/#sigma<","N_{M^{l#nu}_{T}<","N_{p_{T}^{rel}<","N_{#DeltaR(#mu,jet)<","N_{Njets<","N_{#eta^{#mu}<","N_{Nbjets<","N_{m_{qq'}<","N_{|#eta_{lead.jet}|<","N_{#Chi^{2}<","N_{M^{bl#nu}_{T}<"};
39
40 TString EffNameEnd = "}/N";
41
42 const Double_t Cuts[NbOfObs] = {30,250,300,3,2,0.1,0.1,0.1,3,20,10,0.35,4,2.0,1,150,1.8,5,125}; //Cuts on MuonPt, HT, HTall, CaloIso, TrackIso, RelCaloIso, RelTrackIso, RelIso, MET, d0...
43 TString CutsValue[NbOfObs] = {"30","250","300","3","2","0.1","0.1","0.1","3","20","10","0.35","4","2.0","1","150","1.8","5","125"}; //Cuts on MuonPt, HT, HTall, CaloIso, TrackIso, RelCaloIso, RelTrackIso, RelIso, MET, d0...
44
45 const Int_t Nbins[NbOfObs] = {13,12,13,12,12,8,8,12,20,9,7,10,10,8,10,10,8,10,10};
46
47 Double_t MuonPt_bins[Nbins[0]+1] = {0,15,20,22,25,30,40,50,65,80,100,130,165,200};
48 Double_t HT_bins[Nbins[1]+1] = {0,150,200,220,250,270,300,350,400,500,600,700,1000};
49 Double_t HTall_bins[Nbins[2]+1] = {0,150,200,220,250,270,300,350,400,500,600,700,900,1200};
50 Double_t CaloIso_bins[Nbins[3]+1] = {0,1,2,3,4,8,12,16,20,24,28,32,40};
51 Double_t TrackIso_bins[Nbins[4]+1] = {0,1,2,3,4,8,12,16,20,24,28,32,40};
52 Double_t RelCaloIso_bins[Nbins[5]+1] = {0,0.1,0.2,0.3,0.5,0.8,1.3,2,3};
53 Double_t RelTrackIso_bins[Nbins[6]+1] = {0,0.1,0.2,0.3,0.5,0.8,1.3,2,3};
54 Double_t RelIso_bins[Nbins[7]+1] = {0,0.05,0.08,0.12,0.2,0.4,0.7,1.0,1.5,2.0,3.0,4.0,5.0};
55 Double_t d0Sig_bins[Nbins[8]+1] = {0,0.5,1.0,1.5,2.0,2.5,3.0,3.5,4.0,4.5,5,6,7,8,9,10,11,12,13,14,15};
56 Double_t MT_bins[Nbins[9]+1] = {0,10,20,30,40,50,70,100,150,250};
57 Double_t PtRel_bins[Nbins[10]+1] = {0,5,10,15,30,70,200,300};
58 Double_t DRMuJet_bins[Nbins[11]+1] = {0,0.1,0.2,0.3,0.4,0.5,0.8,1.5,3,4,5};
59 Double_t Njets_bins[Nbins[12]+1] = {0,1,2,3,4,5,6,7,8,9,10};
60 Double_t MuonEta_bins[Nbins[13]+1] = {-2.4,-2.1,-1.5,-0.7,0,0.7,1.5,2.1,2.4};
61 Double_t Nbjets_bins[Nbins[14]+1] = {0,1,2,3,4,5,6,7,8,9,10};
62 Double_t Nmqq_bins[Nbins[15]+1] = {0,25,50,100,125,175,200,250,300,350,400};
63 Double_t LeadJetEta_bins[Nbins[16]+1] = {-2.4,-2.1,-1.5,-0.7,0,0.7,1.5,2.1,2.4};
64 Double_t NChi2_bins[Nbins[17]+1] = {0,2.5,5,10,12,17,20,25,30,35,40};
65 Double_t Nmblv_bins[Nbins[18]+1] = {0,25,50,100,125,175,200,250,300,350,400};
66
67 Double_t **bins = new Double_t*[NbOfObs];
68 for(Int_t i=0;i<NbOfObs;i++)
69 {
70 bins[i] = new Double_t[Nbins[i]+1];
71 }
72
73 bins[0] = MuonPt_bins;
74 bins[1] = HT_bins;
75 bins[2] = HTall_bins;
76 bins[3] = CaloIso_bins;
77 bins[4] = TrackIso_bins;
78 bins[5] = RelCaloIso_bins;
79 bins[6] = RelTrackIso_bins;
80 bins[7] = RelIso_bins;
81 bins[8] = d0Sig_bins;
82 bins[9] = MT_bins;
83 bins[10] = PtRel_bins;
84 bins[11]= DRMuJet_bins;
85 bins[12]= Njets_bins;
86 bins[13]= MuonEta_bins;
87 bins[14]= Nbjets_bins;
88 bins[15]= Nmqq_bins;
89 bins[16]= LeadJetEta_bins;
90 bins[17]= NChi2_bins;
91 bins[18]= Nmblv_bins;
92 //bins[12]= Pt4thJet_bins;
93
94 std::cout<<"Initialization done"<<std::endl;
95
96 // --------------------------------------------------------------
97 // Pointers to histograms
98 // --------------------------------------------------------------
99
100 TDirectory *dir = InputFile;
101 //(TDirectory*) InputFile->Get("TSAnalyser");
102
103 TString ObsToGet = "";
104 TString ObsToGet_EffX = "";
105 TString ObsToGet_EffY = "";
106 Int_t Index = 0;
107
108 TH2F *hHisto2D[NbOfObs*(NbOfObs-1)/2];
109 TH1F *hHistoEff[NbOfObs*(NbOfObs-1)/2][2];
110
111 Double_t Eff_X = 0;
112 Double_t Eff_Y = 0;
113 Double_t S = 0;
114 Double_t Ntot = 0;
115
116 Double_t Xmin = 0;
117 Double_t Ymin = 0;
118 Double_t Xmax = 0;
119 Double_t Ymax = 0;
120
121 Int_t NbinsX = 0;
122 Int_t NbinsY = 0;
123
124 Int_t X0 = 0;
125 Int_t Y0 = 0;
126
127 Double_t error = 0;
128
129 TString YaxisName = "";
130
131 OutputFile->cd();
132
133 TCanvas c("dummy","",1);
134 c.Print(ControlHistoFile + "[");
135
136 for(Int_t i = 0; i<NbOfObs ; i++)
137 {
138 for(Int_t j = i+1; j<NbOfObs ; j++)
139 {
140 ObsToGet = "Obs_"; ObsToGet += ObsName[i]; ObsToGet +="_vs_Obs_"; ObsToGet += ObsName[j];
141
142 ObsToGet_EffX = ObsToGet; ObsToGet_EffX += "_EffX";
143 ObsToGet_EffY = ObsToGet; ObsToGet_EffY += "_EffY";
144 cout<<ObsToGet<<endl;
145
146 hHisto2D[Index] = (TH2F *) dir->Get(ObsToGet);
147
148 ////////////////// Efficiency along the X axis /////////////////////////////
149 hHistoEff[Index][0] = new TH1F(ObsToGet_EffX,ObsToGet_EffX,Nbins[i],bins[i]);
150 hHistoEff[Index][0]->SetTitle("");
151 hHistoEff[Index][0]->GetXaxis()->SetTitle(Obs[i]);
152 YaxisName = EffName[j]; YaxisName += CutsValue[j]; YaxisName += EffNameEnd;
153 hHistoEff[Index][0]->GetYaxis()->SetTitle(YaxisName);
154 hHistoEff[Index][0]->GetYaxis()->SetRangeUser(0,1.1);
155 //hHistoEff[Index][0]->GetYaxis()->SetRangeUser(0.5*hHistoEff[Index][0]->GetMinimum(),1.5*hHistoEff[Index][0]->GetMaximum());
156
157 ////////////////// Efficiency along the Y axis /////////////////////////////
158 hHistoEff[Index][1] = new TH1F(ObsToGet_EffY,ObsToGet_EffY,Nbins[j],bins[j]);
159 hHistoEff[Index][1]->SetTitle("");
160 hHistoEff[Index][1]->GetXaxis()->SetTitle(Obs[j]);
161 YaxisName = EffName[i]; YaxisName += CutsValue[i]; YaxisName += EffNameEnd;
162 hHistoEff[Index][1]->GetYaxis()->SetTitle(YaxisName);
163 hHistoEff[Index][1]->GetYaxis()->SetRangeUser(0,1.1);
164 //hHistoEff[Index][1]->GetYaxis()->SetRangeUser(0.5*hHistoEff[Index][1]->GetMinimum(),1.5*hHistoEff[Index][1]->GetMaximum());
165
166 NbinsX = hHisto2D[Index]->GetNbinsX();
167 NbinsY = hHisto2D[Index]->GetNbinsY();
168 Xmin = hHisto2D[Index]->GetXaxis()->GetXmin();Xmax = hHisto2D[Index]->GetXaxis()->GetXmax();
169 Ymin = hHisto2D[Index]->GetYaxis()->GetXmin();Ymax = hHisto2D[Index]->GetYaxis()->GetXmax();
170 X0 = NbinsX*(Cuts[i]-Xmin)/(Xmax-Xmin);
171 Y0 = NbinsY*(Cuts[j]-Ymin)/(Ymax-Ymin);
172
173 if(Print){
174 std::cout<<"----------------------------------------------------------------------------------"<<std::endl;
175 std::cout<<"------------------------------"<<ObsName[i]<<" vs "<<ObsName[j]<<"------------------------------"<<std::endl;
176 std::cout<<"----------------------------------------------------------------------------------"<<std::endl;
177 std::cout<<"--------"<<ObsName[i]<<"range = "<<Xmin<<"-"<<Xmax<<" / Cut = "<<Cuts[i]<<" / Nb of bins = "<<NbinsX<<" / X0 = "<<X0<<" --------"<<std::endl;
178 std::cout<<"--------"<<ObsName[j]<<"range = "<<Ymin<<"-"<<Ymax<<" / Cut = "<<Cuts[j]<<" / Nb of bins = "<<NbinsY<<" / Y0 = "<<Y0<<" --------"<<std::endl;
179 }
180
181 if(Print) cout<<"X Binning = ";
182 for(Int_t k=1; k <= Nbins[i]; k++)
183 {
184 if(Print) cout<<bins[i][k]<<"/";
185 if(hHisto2D[Index]->Integral((bins[i][k-1]*NbinsX/Xmax)+1,bins[i][k]*NbinsX/Xmax,1,NbinsY) > 0 )
186 {
187 S = hHisto2D[Index]->Integral((bins[i][k-1]*NbinsX/Xmax)+1,bins[i][k]*NbinsX/Xmax,1,Y0);
188 Ntot = hHisto2D[Index]->Integral((bins[i][k-1]*NbinsX/Xmax)+1,bins[i][k]*NbinsX/Xmax,1,NbinsY);
189
190 hHistoEff[Index][0]->SetBinContent(k,S/Ntot);
191 hHistoEff[Index][0]->SetBinError(k,EffError(S,Ntot));
192 }
193 else hHistoEff[Index][0]->SetBinContent(k,0);
194 }
195 if(Print) cout<<endl;
196 if(Print) cout<<"Y Binning = ";
197 for(Int_t k=1; k <= Nbins[j]; k++)
198 {
199 if(Print) cout<<bins[j][k]<<"/";
200 if( hHisto2D[Index]->Integral(1 ,NbinsX,(bins[j][k-1]*NbinsY/Ymax)+1,bins[j][k]*NbinsY/Ymax) > 0 )
201 {
202 //S = hHisto2D[Index]->Integral(X0,NbinsX,(bins[j][k-1]*NbinsY/Ymax)+1,bins[j][k]*NbinsY/Ymax);
203 S = hHisto2D[Index]->Integral(1, X0,(bins[j][k-1]*NbinsY/Ymax)+1,bins[j][k]*NbinsY/Ymax);
204 Ntot = hHisto2D[Index]->Integral(1 ,NbinsX,(bins[j][k-1]*NbinsY/Ymax)+1,bins[j][k]*NbinsY/Ymax);
205
206 hHistoEff[Index][1]->SetBinContent(k,S/Ntot);
207 hHistoEff[Index][1]->SetBinError(k,EffError(S,Ntot));
208 }
209 else hHistoEff[Index][1]->SetBinContent(k,0);
210 }
211 hHistoEff[Index][0]->Draw();
212 hHistoEff[Index][0]->Write(ObsToGet_EffX);
213 c.Print(ControlHistoFile);
214 hHistoEff[Index][1]->Draw();
215 hHistoEff[Index][1]->Write(ObsToGet_EffY);
216 c.Print(ControlHistoFile);
217 Index++;
218 }
219 }
220
221
222 //--------------------------------------------------------------------------------------------------------------------------------------------
223 //--------------------------------------------------------------------------------------------------------------------------------------------
224 //--------------------------------------------------------------------------------------------------------------------------------------------
225
226 OutputFile->Close();
227 c.Print(ControlHistoFile + "]");
228 gROOT->ProcessLine(".q");
229
230 }
231
232 Double_t EffError(Double_t S, Double_t Total)
233 {
234 Double_t Eff = S/Total;
235 return (sqrt(S*(1-Eff))/Total);
236 }
237
238