Fork me on GitHub

source: svn/trunk/routines/resolutions.C@ 94

Last change on this file since 94 was 91, checked in by severine ovyn, 16 years ago

for resol plots
File size: 9.3 KB
Line 
1#include "TROOT.h"
2#include "TFile.h"
3#include "TTree.h"
4#include "TCanvas.h"
5#include "TProfile.h"
6#include "TF1.h"
7#include "TGraph.h"
8
9#include "interface/FuncDef.h"
10
11void JetResol()
12{
13 setTDRStyle();
14 gROOT->Reset();
15
16 TFile *f1 = new TFile("JET.root","read");
17 TTree *Analyze = (TTree*)f1->Get("Analysis");
18
19 const Int_t numBin=7;
20 double bins[numBin]={0,20,40,60,80,100,220};
21 TProfile *ETSoverET = new TProfile("ETSoverET","Jet resolution: E_{T}^{rec}/E_{T}^{mc}",(numBin-1),bins,-10,10);
22
23 double rms[numBin];
24 double mean[numBin];
25
26 TCanvas *c1 = new TCanvas("c1","JET resol",0,0,1000,650);
27 c1->cd(); int frame=0;
28 c1->Divide(6,2);
29
30 float x[numBin-1];
31 float y[numBin-1];
32 float ex[numBin-1];
33 float ey[numBin-1];
34
35 float finval=0;//valeur a remplir puis a fitter
36
37 for ( int i=0; i<(numBin-1); i++) // premiÚre bin : i ==1 et pas i == 0
38 {
39 TAxis *xaxis = ETSoverET->GetXaxis();
40 float binCenter = xaxis->GetBinCenter(i+1);
41 int binMin = (int)xaxis->GetBinLowEdge(i+1);
42 int binMax = (int)xaxis->GetBinUpEdge(i+1);
43 char tempMin[500];
44 if(i==0)binMin=5;
45 sprintf(tempMin,"JetPTResol.PT > %d",binMin);
46 string mystringMin(tempMin);
47 char tempMax[500];
48 sprintf(tempMax,"JetPTResol.PT < %d",binMax);
49 string mystringMax(tempMax);
50 char tempName[500];
51 sprintf(tempName,"(JetPTResol.SmearedPT)>>hETSoverET%d",i);
52 string mystringName(tempName);
53
54 c1->cd(++frame);
55 GaussValues(Analyze,tempName,rms[i],mean[i],mystringMin,mystringMax);
56 x[i]=binCenter;
57 finval=rms[i]/mean[i];
58 y[i]=(finval*100);
59 ex[i]=0;
60 ey[i]=0;
61 }
62
63 TCanvas *c2 = new TCanvas("c2","JET resol",100,100,600,450);
64 c2->cd();
65
66 TF1 *fitfun = new TF1("user","sqrt(pow([0],2)+pow([1]/sqrt(x),2)+pow([2]/x,2))",1,170);
67
68 TGraphErrors *gr11 = new TGraphErrors((numBin-1),x,y,ex,ey);
69 gr11->Draw("AP");
70 gr11->SetTitle("");
71 gr11->GetXaxis()->SetTitle("E_{T}^{MC} [GeV]");
72 gr11->GetYaxis()->SetRangeUser(0,50);
73 gr11->GetYaxis()->SetTitle("#sigma(E_{T}^{rec}/E_{T}^{MC})_{fit}/<E_{T}^{rec}/E_{T}^{MC}>_{fit}");
74
75 Double_t* params = fitfun->GetParameters();
76
77 gr11->Fit("user","QR");
78 gr11->Fit("user","QRI");
79 gr11->Fit("user","QRI");
80 gr11->Fit("user","QRI");
81 char tempResol[500];
82 if(params[2]==0.000000)sprintf(tempResol,"#frac{#sigma(E_{T}^{rec}/E_{T}^{MC})}{<E_{T}^{rec}/E_{T}^{MC}>} = #frac{%f}{#sqrt{E_{T}^{MC}}} #oplus %f",params[1]/100,params[0]/100);
83 else sprintf(tempResol,"#frac{#sigma(E_{T}^{rec}/E_{T}^{MC})}{<E_{T}^{rec}/E_{T}^{MC}>} = #frac{%f}{#sqrt{E_{T}^{MC}}} #oplus #frac{%f}{E_{T}^{MC}} #oplus %f",params[1]/100,params[2]/100,params[0]/100);
84
85 TPaveText *leg1 = MakeTPave(0.4,0.6,0.8,0.65,tempResol);
86 leg1->Draw();
87
88 TPaveText *events = MakeTPave(0.2,0.75,0.35,0.8,"Events: pp #rightarrow gg ");
89 events->Draw();
90
91
92 TPaveText *Delphes = MakeTPave(0.2,0.15,0.35,0.2,"MG/ME + Delphes");
93 Delphes->Draw();
94
95
96 delete fitfun;
97}
98
99void ElecResol()
100{
101
102 setTDRStyle();
103 gROOT->Reset();
104
105 TFile *f1 = new TFile("lept.root","read");
106 TTree *Analyze = (TTree*)f1->Get("Analysis");
107
108 const Int_t numBin=6;
109 double bins[numBin]={0,20,40,60,80,140};
110 TProfile *ETSminusET = new TProfile("ETSminusET","Electron resolution: E_{T}^{rec}/E_{T}^{mc}",(numBin-1),bins,-10,10);
111
112 double rms[numBin];
113 double mean[numBin];
114
115 TCanvas *c3 = new TCanvas("c3","ELEC resol",0,0,1000,650);
116 c3->cd(); int frame=0;
117 c3->Divide(6,2);
118
119 float x[numBin-1];
120 float y[numBin-1];
121 float ex[numBin-1];
122 float ey[numBin-1];
123
124 float finval=0;//valeur a remplir puis a fitter
125
126 for ( int i=0; i<(numBin-1); i++) // premiÚre bin : i ==1 et pas i == 0
127 {
128 TAxis *xaxis = ETSminusET->GetXaxis();
129 float binCenter = xaxis->GetBinCenter(i+1);
130 int binMin = (int)xaxis->GetBinLowEdge(i+1);
131 int binMax = (int)xaxis->GetBinUpEdge(i+1);
132 char tempMin[500];
133 if(i==0)binMin=5;
134 sprintf(tempMin,"ElecEResol.E > %d",binMin);
135 string mystringMin(tempMin);
136 char tempMax[500];
137 sprintf(tempMax,"ElecEResol.E < %d",binMax);
138 string mystringMax(tempMax);
139 char tempName[500];
140 sprintf(tempName,"(ElecEResol.E-ElecEResol.SmearedE)>>hETSoverET%d",i);
141 string mystringName(tempName);
142
143 c3->cd(++frame);
144 GaussValues(Analyze,tempName,rms[i],mean[i],mystringMin,mystringMax);
145 x[i]=binCenter;
146 finval=rms[i]/binCenter;
147 y[i]=(finval*100);
148 ex[i]=0;
149 ey[i]=0;
150 }
151
152 TCanvas *c4 = new TCanvas("c4","ELEC resol",100,100,600,450);
153 c4->cd();
154
155 TF1 *fitfun = new TF1("user","sqrt(pow([0],2)+pow([1]/sqrt(x),2)+pow([2]/x,2))",1,400);
156
157 TGraphErrors *gr11 = new TGraphErrors((numBin-1),x,y,ex,ey);
158 gr11->Draw("AP");
159 gr11->SetTitle("");
160 gr11->GetXaxis()->SetTitle("E [GeV]");
161 gr11->GetYaxis()->SetRangeUser(0,5);
162 gr11->GetYaxis()->SetTitle("#sigma/E");
163
164 Double_t* params = fitfun->GetParameters();
165
166 gr11->Fit("user","QR");
167 gr11->Fit("user","QRI");
168 gr11->Fit("user","QRI");
169 gr11->Fit("user","QRI");
170 char tempResol[500];
171 sprintf(tempResol,"#frac{#sigma}{E} = #frac{%f}{#sqrt{E}} #oplus #frac{%f}{E} #oplus %f",params[1]/100,params[2]/100,params[0]/100);
172
173 TPaveText *leg1 = MakeTPave(0.4,0.6,0.8,0.65,tempResol);
174 leg1->Draw();
175
176 TPaveText *Delphes = MakeTPave(0.2,0.15,0.35,0.2,"MG/ME + Delphes");
177 Delphes->Draw();
178
179 TPaveText *events = MakeTPave(0.2,0.75,0.35,0.8,"Events: WHq'#rightarrow W#tau#tauq'#rightarrowjjl#tauq', m_{H}=150 GeV ");
180 events->Draw();
181
182 delete fitfun;
183}
184
185void TauJetInfo()
186{
187
188 setTDRStyle();
189 gROOT->Reset();
190
191 TFile *f1 = new TFile("Etmis.root","read");
192 TTree *Analyze = (TTree*)f1->Get("Analysis");
193
194 TCanvas *ct1 = new TCanvas("ct1","Tau information",100,100,600,450);
195 ct1->cd();
196
197 TH1F *tauEnergy =MakeNormTH1F(20,0.8,1,Analyze,"TauJetPTResol.EnergieCen>>tauEnergy",1, 0, 1,2,false);
198 tauEnergy->Draw();
199 tauEnergy->SetTitle("");
200 tauEnergy->GetYaxis()->SetTitle("Fraction of events");
201 tauEnergy->GetXaxis()->SetTitle("C_{#tau}");
202
203
204 TPaveText *Delphes1 = MakeTPave(0.3,0.85,0.45,0.9,"MG/ME + Delphes");
205 Delphes1->Draw();
206 TPaveText *ctau = MakeTPave(0.3,0.75,0.45,0.8,"C_{#tau} = #frac{#sum E^{towers} (#DeltaR = 0.15)}{E^{jet}}");
207 ctau->Draw();
208 TPaveText *events1 = MakeTPave(0.3,0.65,0.45,0.7,"Events: WHq'#rightarrow W#tau#tauq'#rightarrowjjl#tauq', m_{H}=150 GeV ");
209 events1->Draw();
210
211
212 TCanvas *ct2 = new TCanvas("ct2","Tau information",100,100,600,450);
213 ct2->cd();
214 TH1F *NumTrack =MakeNormTH1F(6,0,6,Analyze,"TauJetPTResol.NumTrack>>NumTrack",1, 0, 1,2,false);
215 NumTrack->Draw();
216 NumTrack->SetTitle("");
217 NumTrack->GetYaxis()->SetTitle("Fraction of events");
218 NumTrack->GetXaxis()->SetTitle("N^{tracks}");
219
220 TPaveText *Delphes = MakeTPave(0.6,0.85,0.85,0.9,"MG/ME + Delphes");
221 Delphes->Draw();
222 TPaveText *numtracks = MakeTPave(0.6,0.75,0.85,0.8,"#DeltaR < 0.4, p_{T}^{track} > 2 GeV");
223 numtracks->Draw();
224 TPaveText *events = MakeTPave(0.6,0.65,0.85,0.7,"Events: WHq'#rightarrow W#tau#tauq'#rightarrowjjl#tauq', m_{H}=150 GeV ");
225 events->Draw();
226
227}
228
229void ETmisResol()
230{
231
232 setTDRStyle();
233 gROOT->Reset();
234
235 TFile *f1 = new TFile("Etmis.root","read");
236 TTree *Analyze = (TTree*)f1->Get("Analysis");
237
238 TF1 *fitfun = new TF1("user","[0]*sqrt(x)",0,600);
239 const Int_t numBin=7;
240 double bins[numBin]={200,250,300,350,400,450,500};
241 TProfile *ETSoverET = new TProfile("ETSoverET","ETmis resol",(numBin-1),bins,-1000,1000);
242
243
244 double rms[numBin-1];
245
246 TCanvas *c5 = new TCanvas("c5","PTmis resol",0,0,1000,650);
247 c5->cd(); int frame=0;
248 c5->Divide(6,2);
249
250 double x[numBin];
251 double y[numBin];
252
253 for ( int i=0; i<(numBin-1); i++) // premiÚre bin : i ==1 et pas i == 0
254 {
255 TAxis *xaxis = ETSoverET->GetXaxis();
256 float binCenter = xaxis->GetBinCenter(i+1);
257 int binMin = (int)xaxis->GetBinLowEdge(i+1);
258 int binMax = (int)xaxis->GetBinUpEdge(i+1);
259 char tempMin[500];
260 sprintf(tempMin,"ETmisResol.SEt>%d",binMin);
261 string mystringMin(tempMin);
262 char tempMax[500];
263 sprintf(tempMax,"ETmisResol.SEt<%d",binMax);
264 string mystringMax(tempMax);
265
266 char tempName[500];
267 sprintf(tempName,"ETmisResol.ExSmeare>>hETSoverET%d",i);
268 string mystringName(tempName);
269
270 c5->cd(++frame);
271 GaussValuesETmis(Analyze,tempName,rms[i],mystringMin,mystringMax);
272 x[i]=binCenter;
273 y[i]=rms[i];
274
275 }
276
277 TCanvas *c6 = new TCanvas("c6","ETmis resolution",100,100,600,450);
278 c6->cd();
279
280 x[numBin]=0;
281 y[numBin]=0;
282 TGraph *gr11 = new TGraph((numBin),x,rms);
283 gr11->Draw("AP");
284 gr11->GetXaxis()->SetTitle("Offline sum of E_{T} [GeV]");
285 gr11->GetYaxis()->SetTitle("Resolution of x-component of MET [GeV]");
286 gr11->Fit("user","RQ");
287 gr11->Fit("user","RQ");
288 gr11->Fit("user","RQ");
289 gr11->Fit("user","RQ");
290 gr11->GetYaxis()->SetRangeUser(0,30);
291 gr11->GetXaxis()->SetRangeUser(1,600);
292
293 Double_t* params = fitfun->GetParameters();
294
295 char tempResol[500];
296 sprintf(tempResol,"%f * #sqrt{E_{T}}",params[0]);
297
298 TPaveText *leg1 = MakeTPave(0.4,0.6,0.8,0.65,tempResol);
299 leg1->Draw();
300
301 TPaveText *leg2 = MakeTPave(0.2,0.8,0.8,0.85,"WHq'#rightarrow W#tau#tauq'#rightarrowjjl#tauq', m_{H}=150 GeV ");
302 leg2->Draw();
303
304 TPaveText *Delphes = MakeTPave(0.2,0.15,0.35,0.2,"MG/ME + Delphes");
305 Delphes->Draw();
306
307
308 delete fitfun;
309}
310
311void General()
312{
313 JetResol();
314 ElecResol();
315 ETmisResol();
316 TauJetInfo();
317
318}
319
320
Note: See TracBrowser for help on using the repository browser.