]>
Commit | Line | Data |
---|---|---|
10ebe68d | 1 | /* $Id$ */ |
2 | ||
3 | #if !defined(__CINT__) || defined(__MAKECINT__) | |
4 | ||
5 | #include "AliPWG0Helper.h" | |
6 | #include "dNdEtaAnalysis.h" | |
7 | #include "AlidNdEtaCorrection.h" | |
8 | ||
9 | #include <TCanvas.h> | |
10 | #include <TFile.h> | |
11 | #include <TH1.h> | |
12 | #include <TH2F.h> | |
13 | #include <TH3F.h> | |
14 | #include <TLine.h> | |
15 | #include <TSystem.h> | |
6b7fa615 | 16 | #include <TLegend.h> |
17 | #include <TPad.h> | |
18 | #include <TF1.h> | |
19 | ||
20 | extern TPad* gPad; | |
21 | ||
22 | void Track2Particle1DCreatePlots(const char* fileName = "correction_map.root", const char* folderName = "dndeta_correction", Float_t upperPtLimit = 10); | |
10ebe68d | 23 | |
24 | #endif | |
25 | ||
9f469bf5 | 26 | |
10ebe68d | 27 | void SetRanges(TAxis* axis) |
28 | { | |
29 | if (strcmp(axis->GetTitle(), "#eta") == 0) | |
30 | axis->SetRangeUser(-1.7999, 1.7999); | |
31 | if (strcmp(axis->GetTitle(), "p_{T} [GeV/c]") == 0) | |
32 | axis->SetRangeUser(0, 9.9999); | |
33 | if (strcmp(axis->GetTitle(), "vtx z [cm]") == 0) | |
34 | axis->SetRangeUser(-15, 14.9999); | |
35 | if (strcmp(axis->GetTitle(), "Ntracks") == 0) | |
36 | axis->SetRangeUser(0, 99.9999); | |
37 | } | |
38 | ||
39 | void SetRanges(TH1* hist) | |
40 | { | |
41 | SetRanges(hist->GetXaxis()); | |
42 | SetRanges(hist->GetYaxis()); | |
43 | SetRanges(hist->GetZaxis()); | |
44 | } | |
45 | ||
46 | void Prepare3DPlot(TH3* hist) | |
47 | { | |
48 | hist->GetXaxis()->SetTitleOffset(1.5); | |
49 | hist->GetYaxis()->SetTitleOffset(1.5); | |
50 | hist->GetZaxis()->SetTitleOffset(1.5); | |
51 | ||
52 | hist->SetStats(kFALSE); | |
53 | } | |
54 | ||
55 | void Prepare2DPlot(TH2* hist) | |
56 | { | |
57 | hist->SetStats(kFALSE); | |
58 | hist->GetYaxis()->SetTitleOffset(1.4); | |
59 | ||
60 | SetRanges(hist); | |
61 | } | |
62 | ||
7af955da | 63 | void Prepare1DPlot(TH1* hist, Bool_t setRanges = kTRUE) |
10ebe68d | 64 | { |
65 | hist->SetLineWidth(2); | |
66 | hist->SetStats(kFALSE); | |
67 | ||
7af955da | 68 | hist->GetXaxis()->SetTitleOffset(1.2); |
69 | hist->GetYaxis()->SetTitleOffset(1.2); | |
70 | ||
71 | if (setRanges) | |
72 | SetRanges(hist); | |
10ebe68d | 73 | } |
74 | ||
75 | void InitPad() | |
76 | { | |
9f469bf5 | 77 | if (!gPad) |
78 | return; | |
79 | ||
10ebe68d | 80 | gPad->Range(0, 0, 1, 1); |
81 | gPad->SetLeftMargin(0.15); | |
82 | //gPad->SetRightMargin(0.05); | |
83 | //gPad->SetTopMargin(0.13); | |
84 | //gPad->SetBottomMargin(0.1); | |
85 | ||
72e597d7 | 86 | gPad->SetGridx(); |
87 | gPad->SetGridy(); | |
10ebe68d | 88 | } |
89 | ||
90 | void InitPadCOLZ() | |
91 | { | |
92 | gPad->Range(0, 0, 1, 1); | |
93 | gPad->SetRightMargin(0.15); | |
94 | gPad->SetLeftMargin(0.12); | |
95 | ||
96 | gPad->SetGridx(); | |
97 | gPad->SetGridy(); | |
98 | } | |
99 | ||
100 | void Secondaries() | |
101 | { | |
102 | TFile* file = TFile::Open("systematics.root"); | |
103 | ||
7af955da | 104 | TH2F* secondaries = dynamic_cast<TH2F*> (file->Get("fSecondaries")); |
10ebe68d | 105 | if (!secondaries) |
106 | { | |
107 | printf("Could not read histogram\n"); | |
108 | return; | |
109 | } | |
110 | ||
7af955da | 111 | TCanvas* canvas = new TCanvas("Secondaries", "Secondaries", 1000, 1000); |
112 | canvas->Divide(3, 3); | |
113 | for (Int_t i=1; i<=8; i++) | |
10ebe68d | 114 | { |
7af955da | 115 | TH1D* hist = secondaries->ProjectionY(Form("proj_%d", i), i, i); |
116 | hist->SetTitle(secondaries->GetXaxis()->GetBinLabel(i)); | |
10ebe68d | 117 | |
7af955da | 118 | canvas->cd(i); |
119 | hist->Draw(); | |
10ebe68d | 120 | } |
121 | } | |
122 | ||
6b7fa615 | 123 | void Track2Particle1DComposition(const char** fileNames, Int_t folderCount, const char** folderNames, Float_t upperPtLimit = 9.9) |
10ebe68d | 124 | { |
125 | gSystem->Load("libPWG0base"); | |
126 | ||
9f469bf5 | 127 | TString canvasName; |
128 | canvasName.Form("Track2Particle1DComposition"); | |
129 | TCanvas* canvas = new TCanvas(canvasName, canvasName, 1200, 400); | |
130 | canvas->Divide(3, 1); | |
10ebe68d | 131 | |
9f469bf5 | 132 | TLegend* legend = new TLegend(0.7, 0.7, 0.95, 0.95); |
10ebe68d | 133 | |
9f469bf5 | 134 | for (Int_t i=0; i<folderCount; ++i) |
10ebe68d | 135 | { |
6b7fa615 | 136 | Track2Particle1DCreatePlots(fileNames[i], folderNames[i], upperPtLimit); |
9f469bf5 | 137 | |
6b7fa615 | 138 | TH1* corrX = dynamic_cast<TH1*> (gROOT->FindObject(Form("gene_%s_nTrackToNPart_x_div_meas_%s_nTrackToNPart_x", folderNames[i], folderNames[i]))); |
139 | TH1* corrY = dynamic_cast<TH1*> (gROOT->FindObject(Form("gene_%s_nTrackToNPart_y_div_meas_%s_nTrackToNPart_y", folderNames[i], folderNames[i]))); | |
140 | TH1* corrZ = dynamic_cast<TH1*> (gROOT->FindObject(Form("gene_%s_nTrackToNPart_z_div_meas_%s_nTrackToNPart_z", folderNames[i], folderNames[i]))); | |
9f469bf5 | 141 | |
142 | Prepare1DPlot(corrX); | |
143 | Prepare1DPlot(corrY); | |
144 | Prepare1DPlot(corrZ); | |
145 | ||
146 | const char* title = "Track2Particle Correction"; | |
147 | corrX->SetTitle(title); | |
148 | corrY->SetTitle(title); | |
149 | corrZ->SetTitle(title); | |
150 | ||
151 | corrZ->GetXaxis()->SetRangeUser(0, upperPtLimit); | |
152 | ||
153 | corrX->SetLineColor(i+1); | |
154 | corrY->SetLineColor(i+1); | |
155 | corrZ->SetLineColor(i+1); | |
156 | ||
157 | canvas->cd(1); | |
158 | InitPad(); | |
159 | corrX->DrawCopy(((i>0) ? "SAME" : "")); | |
160 | ||
161 | canvas->cd(2); | |
162 | InitPad(); | |
163 | corrY->DrawCopy(((i>0) ? "SAME" : "")); | |
164 | ||
165 | canvas->cd(3); | |
166 | InitPad(); | |
167 | corrZ->DrawCopy(((i>0) ? "SAME" : "")); | |
168 | ||
169 | legend->AddEntry(corrZ, folderNames[i]); | |
10ebe68d | 170 | } |
171 | ||
9f469bf5 | 172 | legend->Draw(); |
173 | ||
174 | //canvas->SaveAs(Form("Track2Particle1D_%s_%d_%f.gif", fileName, gMax, upperPtLimit)); | |
175 | //canvas->SaveAs(Form("Track2Particle1D_%s_%d_%f.eps", fileName, gMax, upperPtLimit)); | |
176 | } | |
177 | ||
6b7fa615 | 178 | TH1** DrawRatios(const char* fileName = "systematics.root") |
179 | { | |
180 | gSystem->Load("libPWG0base"); | |
181 | ||
182 | TFile* file = TFile::Open(fileName); | |
183 | ||
184 | TString canvasName; | |
185 | canvasName.Form("DrawRatios"); | |
186 | TCanvas* canvas = new TCanvas(canvasName, canvasName, 800, 400); | |
187 | canvas->Divide(2, 1); | |
188 | canvas->cd(1); | |
189 | ||
190 | TH1** ptDists = new TH1*[4]; | |
191 | ||
7af955da | 192 | TLegend* legend = new TLegend(0.73, 0.73, 0.98, 0.98); |
6b7fa615 | 193 | |
194 | const char* folderNames[] = { "correction_0", "correction_1", "correction_2", "correction_3" }; | |
195 | const char* particleNames[] = { "#pi", "K", "p", "other" }; | |
196 | for (Int_t i=0; i<4; ++i) | |
197 | { | |
198 | AlidNdEtaCorrection* dNdEtaCorrection = new AlidNdEtaCorrection(folderNames[i], folderNames[i]); | |
199 | dNdEtaCorrection->LoadHistograms(fileName, folderNames[i]); | |
200 | ||
201 | TH3F* gene = dNdEtaCorrection->GetTrack2ParticleCorrection()->GetGeneratedHistogram(); | |
202 | ||
203 | gene->GetYaxis()->SetRangeUser(-0.8, 0.8); | |
204 | gene->GetXaxis()->SetRangeUser(-10, 10); | |
205 | ||
206 | ptDists[i] = dynamic_cast<TH1*> (gene->Project3D("z")->Clone(Form("%s_z", folderNames[i]))); | |
207 | ptDists[i]->SetTitle(";p_{T};Count"); | |
208 | if (!ptDists[i]) | |
209 | { | |
210 | printf("Problem getting distribution %d.\n", i); | |
211 | return 0; | |
212 | } | |
213 | ||
214 | ptDists[i]->GetYaxis()->SetRangeUser(1, ptDists[i]->GetMaximum()*1.1); | |
215 | ptDists[i]->GetXaxis()->SetRangeUser(0, 9.9); | |
216 | ptDists[i]->SetLineColor(i+1); | |
217 | ptDists[i]->DrawCopy((i == 0) ? "" : "SAME"); | |
218 | ptDists[i]->GetYaxis()->SetRange(0, 0); | |
219 | ||
220 | legend->AddEntry(ptDists[i], particleNames[i]); | |
221 | } | |
222 | gPad->SetLogy(); | |
223 | ||
224 | TH1* total = dynamic_cast<TH1*> (ptDists[0]->Clone("total")); | |
225 | ||
226 | for (Int_t i=1; i<4; ++i) | |
227 | total->Add(ptDists[i]); | |
228 | ||
229 | canvas->cd(2); | |
230 | for (Int_t i=0; i<4; ++i) | |
231 | { | |
232 | ptDists[i]->Divide(total); | |
7af955da | 233 | ptDists[i]->SetStats(kFALSE); |
234 | ptDists[i]->SetTitle(";p_{T};Fraction of total"); | |
6b7fa615 | 235 | ptDists[i]->GetYaxis()->SetRangeUser(0, 1); |
236 | ptDists[i]->Draw((i == 0) ? "" : "SAME"); | |
237 | } | |
7af955da | 238 | legend->SetFillColor(0); |
6b7fa615 | 239 | legend->Draw(); |
240 | ||
241 | canvas->SaveAs("DrawRatios.gif"); | |
242 | ||
7af955da | 243 | |
244 | canvas = new TCanvas("PythiaRatios", "PythiaRatios", 500, 500); | |
245 | for (Int_t i=0; i<4; ++i) | |
246 | { | |
247 | TH1* hist = ptDists[i]->Clone(); | |
248 | hist->GetXaxis()->SetRangeUser(0, 1.9); | |
249 | hist->Draw((i == 0) ? "" : "SAME"); | |
250 | } | |
251 | legend->Draw(); | |
252 | ||
253 | canvas->SaveAs("pythiaratios.eps"); | |
254 | ||
6b7fa615 | 255 | file->Close(); |
256 | ||
257 | return ptDists; | |
258 | } | |
259 | ||
7af955da | 260 | void DrawCompareToReal() |
261 | { | |
262 | gROOT->ProcessLine(".L drawPlots.C"); | |
263 | ||
264 | const char* fileNames[] = { "correction_map.root", "new_compositions.root" }; | |
265 | const char* folderNames[] = { "dndeta_correction", "PythiaRatios" }; | |
266 | ||
267 | Track2Particle1DComposition(fileNames, 2, folderNames); | |
268 | } | |
269 | ||
9f469bf5 | 270 | void DrawDifferentSpecies() |
271 | { | |
272 | gROOT->ProcessLine(".L drawPlots.C"); | |
273 | ||
6b7fa615 | 274 | const char* fileNames[] = { "systematics.root", "systematics.root", "systematics.root", "systematics.root" }; |
9f469bf5 | 275 | const char* folderNames[] = { "correction_0", "correction_1", "correction_2", "correction_3" }; |
276 | ||
6b7fa615 | 277 | Track2Particle1DComposition(fileNames, 4, folderNames); |
9f469bf5 | 278 | } |
10ebe68d | 279 | |
6b7fa615 | 280 | void DrawSpeciesAndCombination() |
281 | { | |
282 | gROOT->ProcessLine(".L drawPlots.C"); | |
283 | ||
284 | const char* fileNames[] = { "systematics.root", "systematics.root", "systematics.root", "new_compositions.root" }; | |
285 | const char* folderNames[] = { "correction_0", "correction_1", "correction_2", "PythiaRatios" }; | |
286 | ||
287 | Track2Particle1DComposition(fileNames, 4, folderNames); | |
288 | } | |
289 | ||
290 | void DrawSimulatedVsCombined() | |
291 | { | |
292 | gROOT->ProcessLine(".L drawPlots.C"); | |
293 | ||
294 | const char* fileNames[] = { "new_compositions.root", "new_compositions.root" }; | |
295 | const char* folderNames[] = { "Pythia", "PythiaRatios" }; | |
296 | ||
297 | Track2Particle1DComposition(fileNames, 2, folderNames); | |
298 | } | |
299 | ||
300 | void DrawBoosts() | |
301 | { | |
302 | gROOT->ProcessLine(".L drawPlots.C"); | |
303 | ||
304 | const char* fileNames[] = { "new_compositions.root", "new_compositions.root", "new_compositions.root", "new_compositions.root" }; | |
305 | const char* folderNames[] = { "PythiaRatios", "PiBoosted", "KBoosted", "pBoosted" }; | |
306 | ||
307 | Track2Particle1DComposition(fileNames, 4, folderNames); | |
308 | } | |
309 | ||
7af955da | 310 | TH2F* HistogramDifferences(const char* filename, const char* folder1, const char* folder2) |
6b7fa615 | 311 | { |
312 | gSystem->Load("libPWG0base"); | |
313 | ||
314 | AlidNdEtaCorrection* fdNdEtaCorrection[2]; | |
315 | ||
316 | TFile::Open(filename); | |
317 | ||
318 | fdNdEtaCorrection[0] = new AlidNdEtaCorrection(folder1, folder1); | |
319 | fdNdEtaCorrection[0]->LoadHistograms(filename, folder1); | |
320 | ||
321 | fdNdEtaCorrection[1] = new AlidNdEtaCorrection(folder2, folder2); | |
322 | fdNdEtaCorrection[1]->LoadHistograms(filename, folder2); | |
323 | ||
6b7fa615 | 324 | TH3F* hist1 = fdNdEtaCorrection[0]->GetTrack2ParticleCorrection()->GetCorrectionHistogram(); |
325 | TH3F* hist2 = fdNdEtaCorrection[1]->GetTrack2ParticleCorrection()->GetCorrectionHistogram(); | |
326 | ||
7af955da | 327 | //TH1F* difference = new TH1F("difference", Form(";#DeltaC_{pT, z, #eta} %s / %s;Count", folder2, folder1), 1000, 0.9, 1.1); |
328 | TH2F* difference = new TH2F(Form("difference_%s_%s", folder1, folder2), Form(";#Sigma (C_{pT, z} %s / C_{pT, z} %s);#eta;Count", folder2, folder1), 100, 0.9, 1.1, hist1->GetYaxis()->GetNbins(), hist1->GetYaxis()->GetXmin(), hist1->GetYaxis()->GetXmax()); | |
329 | ||
6b7fa615 | 330 | for (Int_t x=hist1->GetXaxis()->FindBin(-10); x<=hist1->GetXaxis()->FindBin(10); ++x) |
331 | for (Int_t y=hist1->GetYaxis()->FindBin(-0.8); y<=hist1->GetYaxis()->FindBin(0.8); ++y) | |
332 | for (Int_t z=hist1->GetZaxis()->FindBin(0.3); z<=hist1->GetZaxis()->FindBin(9.9); ++z) | |
7af955da | 333 | if (hist1->GetBinContent(x, y, z) != 0) |
334 | difference->Fill(hist2->GetBinContent(x, y, z) / hist1->GetBinContent(x, y, z), hist1->GetYaxis()->GetBinCenter(y)); | |
335 | ||
336 | difference->GetYaxis()->SetRangeUser(-0.8, 0.8); | |
6b7fa615 | 337 | |
338 | printf("Over-/Underflow bins: %d %d\n", difference->GetBinContent(0), difference->GetBinContent(difference->GetNbinsX()+1)); | |
339 | ||
340 | return difference; | |
341 | } | |
342 | ||
343 | void HistogramDifferences() | |
344 | { | |
7af955da | 345 | TH2F* KBoosted = HistogramDifferences("new_compositions.root", "PythiaRatios", "KBoosted"); |
346 | TH2F* pBoosted = HistogramDifferences("new_compositions.root", "PythiaRatios", "pBoosted"); | |
347 | TH2F* KReduced = HistogramDifferences("new_compositions.root", "PythiaRatios", "KReduced"); | |
348 | TH2F* pReduced = HistogramDifferences("new_compositions.root", "PythiaRatios", "pReduced"); | |
6b7fa615 | 349 | |
7af955da | 350 | TCanvas* canvas = new TCanvas("HistogramDifferences", "HistogramDifferences", 1000, 1000); |
351 | canvas->Divide(2, 2); | |
6b7fa615 | 352 | |
353 | canvas->cd(1); | |
7af955da | 354 | KBoosted->GetXaxis()->SetRangeUser(-0.05, 0.05); |
355 | KBoosted->Draw("COLZ"); | |
6b7fa615 | 356 | |
357 | canvas->cd(2); | |
7af955da | 358 | KReduced->GetXaxis()->SetRangeUser(-0.05, 0.05); |
359 | KReduced->Draw("COLZ"); | |
6b7fa615 | 360 | |
361 | canvas->cd(3); | |
7af955da | 362 | pBoosted->GetXaxis()->SetRangeUser(-0.02, 0.02); |
363 | pBoosted->Draw("COLZ"); | |
364 | ||
365 | canvas->cd(4); | |
366 | pReduced->GetXaxis()->SetRangeUser(-0.02, 0.02); | |
367 | pReduced->Draw("COLZ"); | |
6b7fa615 | 368 | |
369 | canvas->SaveAs("HistogramDifferences.gif"); | |
7af955da | 370 | |
371 | canvas = new TCanvas("HistogramDifferencesProfile", "HistogramDifferencesProfile", 1000, 500); | |
372 | canvas->Divide(2, 1); | |
373 | ||
374 | canvas->cd(1); | |
375 | gPad->SetBottomMargin(0.13); | |
376 | KBoosted->SetTitle("a) Ratio of correction maps"); | |
377 | KBoosted->SetStats(kFALSE); | |
378 | KBoosted->GetXaxis()->SetTitleOffset(1.4); | |
379 | KBoosted->GetXaxis()->SetLabelOffset(0.02); | |
380 | KBoosted->Draw("COLZ"); | |
381 | ||
382 | canvas->cd(2); | |
383 | gPad->SetGridx(); | |
384 | gPad->SetGridy(); | |
385 | ||
386 | TLegend* legend = new TLegend(0.73, 0.73, 0.98, 0.98); | |
387 | ||
388 | for (Int_t i=0; i<4; ++i) | |
389 | { | |
390 | TH2F* hist = 0; | |
391 | TString name; | |
392 | switch (i) | |
393 | { | |
394 | case 0: hist = KBoosted; name = "K enhanced"; break; | |
395 | case 1: hist = KReduced; name = "K reduced"; break; | |
396 | case 2: hist = pBoosted; name = "p enhanced"; break; | |
397 | case 3: hist = pReduced; name = "p reduced"; break; | |
398 | } | |
399 | ||
400 | TProfile* profile = hist->ProfileY(); | |
401 | profile->SetTitle("b) Mean and RMS"); | |
402 | profile->GetXaxis()->SetRange(hist->GetYaxis()->GetFirst(), hist->GetYaxis()->GetLast()); | |
403 | profile->GetXaxis()->SetTitleOffset(1.2); | |
404 | profile->GetXaxis()->SetLabelOffset(0.02); | |
405 | profile->GetYaxis()->SetRangeUser(0.98, 1.02); | |
406 | profile->SetStats(kFALSE); | |
407 | profile->SetLineColor(i+1); | |
408 | profile->SetMarkerColor(i+1); | |
409 | profile->DrawCopy(((i > 0) ? "SAME" : "")); | |
410 | ||
411 | ||
412 | legend->AddEntry(profile, name); | |
413 | } | |
414 | ||
415 | legend->Draw(); | |
416 | canvas->SaveAs("particlecomposition_result.eps"); | |
6b7fa615 | 417 | } |
418 | ||
419 | ||
9f469bf5 | 420 | void ScalePtDependent(TH3F* hist, TF1* function) |
421 | { | |
422 | // assumes that pt is the third dimension of hist | |
423 | // scales with function(pt) | |
10ebe68d | 424 | |
9f469bf5 | 425 | for (Int_t z=1; z<=hist->GetNbinsZ(); ++z) |
426 | { | |
427 | Double_t factor = function->Eval(hist->GetZaxis()->GetBinCenter(z)); | |
428 | printf("z = %d, pt = %f, scaling with %f\n", z, hist->GetZaxis()->GetBinCenter(z), factor); | |
429 | ||
430 | for (Int_t x=1; x<=hist->GetNbinsX(); ++x) | |
431 | for (Int_t y=1; y<=hist->GetNbinsY(); ++y) | |
432 | hist->SetBinContent(x, y, z, hist->GetBinContent(x, y, z) * factor); | |
433 | } | |
434 | } | |
435 | ||
6b7fa615 | 436 | void ScalePtDependent(TH3F* hist, TH1* function) |
437 | { | |
438 | // assumes that pt is the third dimension of hist | |
439 | // scales with histogram(pt) | |
440 | ||
441 | for (Int_t z=1; z<=hist->GetNbinsZ(); ++z) | |
442 | { | |
443 | Double_t factor = function->GetBinContent(function->GetXaxis()->FindBin(hist->GetZaxis()->GetBinCenter(z))); | |
444 | printf("z = %d, pt = %f, scaling with %f\n", z, hist->GetZaxis()->GetBinCenter(z), factor); | |
445 | ||
446 | for (Int_t x=1; x<=hist->GetNbinsX(); ++x) | |
447 | for (Int_t y=1; y<=hist->GetNbinsY(); ++y) | |
448 | hist->SetBinContent(x, y, z, hist->GetBinContent(x, y, z) * factor); | |
449 | } | |
450 | } | |
451 | ||
452 | const char* ChangeComposition(void** correctionPointer, Int_t index) | |
9f469bf5 | 453 | { |
454 | AlidNdEtaCorrection** fdNdEtaCorrection = (AlidNdEtaCorrection**) correctionPointer; | |
455 | ||
456 | switch (index) | |
457 | { | |
6b7fa615 | 458 | case 0: // result from pp events |
459 | { | |
460 | TFile::Open("pythiaratios.root"); | |
461 | ||
462 | for (Int_t i=0; i<4; ++i) | |
463 | { | |
464 | TString name; | |
465 | name.Form("correction_%d", i); | |
466 | fdNdEtaCorrection[i] = new AlidNdEtaCorrection(name, name); | |
467 | fdNdEtaCorrection[i]->LoadHistograms("pythiaratios.root", name); | |
468 | } | |
469 | } | |
470 | return "Pythia"; | |
471 | break; | |
472 | ||
473 | case 1: // each species rated with pythia ratios | |
7af955da | 474 | /*TH1** ptDists = DrawRatios("pythiaratios.root"); |
6b7fa615 | 475 | |
476 | for (Int_t i=0; i<3; ++i) | |
477 | { | |
478 | ScalePtDependent(fdNdEtaCorrection[i]->GetTrack2ParticleCorrection()->GetMeasuredHistogram(), ptDists[i]); | |
479 | ScalePtDependent(fdNdEtaCorrection[i]->GetTrack2ParticleCorrection()->GetGeneratedHistogram(), ptDists[i]); | |
7af955da | 480 | }*/ |
6b7fa615 | 481 | return "PythiaRatios"; |
9f469bf5 | 482 | break; |
483 | ||
7af955da | 484 | // one species enhanced / reduced |
485 | case 2: // + 50% pions | |
486 | case 3: // - 50% pions | |
487 | case 4: // + 50% kaons | |
488 | case 5: // - 50% kaons | |
489 | case 6: // + 50% protons | |
490 | case 7: // - 50% protons | |
491 | Int_t correctionIndex = (index - 2) / 2; | |
492 | Double_t scaleFactor = (index % 2 == 0) ? 1.5 : 0.5; | |
493 | ||
494 | fdNdEtaCorrection[correctionIndex]->GetTrack2ParticleCorrection()->GetMeasuredHistogram()->Scale(scaleFactor); | |
495 | fdNdEtaCorrection[correctionIndex]->GetTrack2ParticleCorrection()->GetGeneratedHistogram()->Scale(scaleFactor); | |
496 | ||
497 | TString* str = new TString; | |
498 | str->Form("%s%s", (correctionIndex == 0) ? "Pi" : ((correctionIndex == 1) ? "K" : "p"), (index % 2 == 0) ? "Boosted" : "Reduced"); | |
499 | return str->Data(); | |
500 | break; | |
501 | ||
6b7fa615 | 502 | // each species rated with pythia ratios |
7af955da | 503 | case 12: // + 50% pions |
504 | case 13: // - 50% pions | |
505 | case 14: // + 50% kaons | |
506 | case 15: // - 50% kaons | |
507 | case 16: // + 50% protons | |
508 | case 17: // - 50% protons | |
6b7fa615 | 509 | TH1** ptDists = DrawRatios("pythiaratios.root"); |
510 | Int_t functionIndex = (index - 2) / 2; | |
7af955da | 511 | Double_t scaleFactor = (index % 2 == 0) ? 1.5 : 0.5; |
6b7fa615 | 512 | ptDists[functionIndex]->Scale(scaleFactor); |
513 | ||
514 | for (Int_t i=0; i<3; ++i) | |
515 | { | |
516 | ScalePtDependent(fdNdEtaCorrection[i]->GetTrack2ParticleCorrection()->GetMeasuredHistogram(), ptDists[i]); | |
517 | ScalePtDependent(fdNdEtaCorrection[i]->GetTrack2ParticleCorrection()->GetGeneratedHistogram(), ptDists[i]); | |
518 | } | |
519 | TString* str = new TString; | |
520 | str->Form("%s%s", (functionIndex == 0) ? "Pi" : ((functionIndex == 1) ? "K" : "p"), (index % 2 == 0) ? "Boosted" : "Reduced"); | |
521 | return str->Data(); | |
9f469bf5 | 522 | break; |
523 | ||
6b7fa615 | 524 | case 999: |
9f469bf5 | 525 | TF1* ptDependence = new TF1("simple", "x", 0, 100); |
526 | ScalePtDependent(fdNdEtaCorrection[0]->GetTrack2ParticleCorrection()->GetGeneratedHistogram(), ptDependence); | |
527 | ScalePtDependent(fdNdEtaCorrection[0]->GetTrack2ParticleCorrection()->GetMeasuredHistogram(), ptDependence); | |
528 | break; | |
529 | ||
530 | } | |
6b7fa615 | 531 | |
532 | return "noname"; | |
9f469bf5 | 533 | } |
534 | ||
535 | void Composition() | |
536 | { | |
537 | gSystem->Load("libPWG0base"); | |
538 | ||
539 | gSystem->Unlink("new_compositions.root"); | |
540 | ||
6b7fa615 | 541 | Int_t nCompositions = 8; |
542 | for (Int_t comp = 0; comp < nCompositions; ++comp) | |
9f469bf5 | 543 | { |
544 | AlidNdEtaCorrection* fdNdEtaCorrection[4]; | |
545 | ||
546 | TFile::Open("systematics.root"); | |
547 | ||
548 | for (Int_t i=0; i<4; ++i) | |
549 | { | |
550 | TString name; | |
551 | name.Form("correction_%d", i); | |
552 | fdNdEtaCorrection[i] = new AlidNdEtaCorrection(name, name); | |
553 | fdNdEtaCorrection[i]->LoadHistograms("systematics.root", name); | |
554 | } | |
555 | ||
6b7fa615 | 556 | const char* newName = ChangeComposition(fdNdEtaCorrection, comp); |
9f469bf5 | 557 | |
558 | Double_t geneCount[5]; | |
559 | Double_t measCount[5]; | |
560 | geneCount[4] = 0; | |
561 | measCount[4] = 0; | |
562 | ||
563 | for (Int_t i=0; i<4; ++i) | |
564 | { | |
565 | geneCount[i] = fdNdEtaCorrection[i]->GetTrack2ParticleCorrection()->GetGeneratedHistogram()->Integral(); | |
566 | measCount[i] = fdNdEtaCorrection[i]->GetTrack2ParticleCorrection()->GetMeasuredHistogram()->Integral(); | |
567 | ||
568 | geneCount[4] += geneCount[i]; | |
569 | measCount[4] += measCount[i]; | |
570 | ||
571 | printf("Particle %s: %d gene, %d meas\n", ((i == 0) ? "pi" : (i == 1) ? "K" : (i == 2) ? "p" : "others"), (Int_t) geneCount[i], (Int_t) measCount[i]); | |
572 | } | |
10ebe68d | 573 | |
9f469bf5 | 574 | printf("Generated ratios are: %f pi, %f K, %f p, %f others\n", geneCount[0] / geneCount[4], geneCount[1] / geneCount[4], geneCount[2] / geneCount[4], geneCount[3] / geneCount[4]); |
10ebe68d | 575 | |
9f469bf5 | 576 | printf("Reconstructed ratios are: %f pi, %f K, %f p, %f others\n", measCount[0] / measCount[4], measCount[1] / measCount[4], measCount[2] / measCount[4], measCount[3] / measCount[4]); |
10ebe68d | 577 | |
9f469bf5 | 578 | TList* collection = new TList; |
10ebe68d | 579 | |
9f469bf5 | 580 | for (Int_t i=0; i<4; ++i) |
581 | collection->Add(fdNdEtaCorrection[i]); | |
10ebe68d | 582 | |
6b7fa615 | 583 | AlidNdEtaCorrection* newComposition = new AlidNdEtaCorrection(newName, newName); |
9f469bf5 | 584 | newComposition->Merge(collection); |
585 | newComposition->Finish(); | |
586 | ||
587 | delete collection; | |
588 | ||
589 | TFile* file = TFile::Open("new_compositions.root", "UPDATE"); | |
590 | newComposition->SaveHistograms(); | |
591 | //file->Write(); | |
592 | file->Close(); | |
593 | } | |
10ebe68d | 594 | |
595 | gROOT->ProcessLine(".L drawPlots.C"); | |
9f469bf5 | 596 | |
6b7fa615 | 597 | const char* fileNames[] = {"new_compositions.root", "new_compositions.root", "new_compositions.root", "new_compositions.root", "new_compositions.root", "new_compositions.root", "new_compositions.root", "new_compositions.root" }; |
598 | const char* folderNames[] = { "Pythia", "PythiaRatios", "PiBoosted", "PiReduced", "KBoosted", "KReduced", "pBoosted", "pReduced" }; | |
9f469bf5 | 599 | |
6b7fa615 | 600 | Track2Particle1DComposition(fileNames, nCompositions, folderNames); |
9f469bf5 | 601 | } |
602 | ||
603 | Double_t ConvSigma1To2D(Double_t sigma) | |
604 | { | |
605 | return TMath::Sqrt( - TMath::Log( 1 - TMath::Erf(sigma / TMath::Sqrt(2)) ) * 2); | |
606 | } | |
607 | ||
608 | Double_t ConvDistance1DTo2D(Double_t distance) | |
609 | { | |
610 | return TMath::ErfInverse(1 - TMath::Exp(-distance * distance / 2)) * TMath::Sqrt(2); | |
611 | } | |
612 | ||
613 | Double_t Sigma2VertexCount(TH2F* tracks, Double_t nSigma) | |
614 | { | |
615 | Double_t count = 0; | |
616 | ||
617 | //nSigma = ConvSigma1To2D(nSigma); | |
618 | ||
619 | for (Int_t x=1; x<=tracks->GetNbinsX(); ++x) | |
620 | for (Int_t y=1; y<=tracks->GetNbinsY(); ++y) | |
621 | { | |
622 | Double_t impactX = tracks->GetXaxis()->GetBinCenter(x); | |
623 | Double_t impactY = tracks->GetYaxis()->GetBinCenter(y); | |
624 | ||
625 | Float_t d = TMath::Sqrt(impactX*impactX + impactY*impactY); | |
626 | ||
627 | d = ConvDistance1DTo2D(d); | |
628 | ||
629 | if (d < nSigma) | |
630 | count += tracks->GetBinContent(x, y); | |
631 | } | |
632 | ||
633 | return count; | |
634 | } | |
635 | ||
636 | TH2F* Sigma2VertexGaussianTracksHist() | |
637 | { | |
638 | TH2F* tracks = new TH2F("Sigma2Vertex_tracks", "Sigma2Vertex_tracks", 200, -5, 5, 200, -5, 5); | |
639 | ||
640 | TF2* gaussian2D = new TF2("gaussian2D", "xgausn(0) * ygausn(3)", -5, 5, -5, 5); | |
641 | gaussian2D->SetParameters(1, 0, 1, 1, 0, 1); | |
642 | ||
643 | for (Int_t x=1; x<=tracks->GetNbinsX(); ++x) | |
644 | for (Int_t y=1; y<=tracks->GetNbinsY(); ++y) | |
645 | tracks->SetBinContent(x, y, gaussian2D->Eval(tracks->GetXaxis()->GetBinCenter(x), tracks->GetYaxis()->GetBinCenter(y))); | |
646 | ||
647 | //normalize | |
648 | tracks->Scale(1.0 / tracks->Integral()); | |
649 | ||
650 | return tracks; | |
651 | } | |
652 | ||
4c6b34a8 | 653 | TH1F* Sigma2VertexGaussian() |
9f469bf5 | 654 | { |
655 | TH2F* tracks = Sigma2VertexGaussianTracksHist(); | |
656 | ||
657 | TCanvas* canvas = new TCanvas("Sigma2VertexGaussian", "Sigma2VertexGaussian", 1000, 1000); | |
658 | canvas->Divide(2, 2); | |
659 | ||
660 | canvas->cd(1); | |
661 | tracks->Draw("COLZ"); | |
662 | ||
72e597d7 | 663 | TH1F* ratio = new TH1F("Sigma2Vertex_ratio", "Sigma2Vertex_ratio;n sigma;included", 50, 0.05, 5.05); |
664 | for (Double_t nSigma = 0.1; nSigma < 5.05; nSigma += 0.1) | |
9f469bf5 | 665 | ratio->Fill(nSigma, Sigma2VertexCount(tracks, nSigma)); |
666 | ratio->SetMarkerStyle(21); | |
667 | ||
668 | canvas->cd(2); | |
4c6b34a8 | 669 | ratio->DrawCopy("P"); |
9f469bf5 | 670 | |
72e597d7 | 671 | TH1F* ratio2 = new TH1F("Sigma2Vertex_ratio2", "Sigma2Vertex_ratio2;nSigma;% included 3 sigma / % included n sigma", 50, 0.05, 5.05); |
9f469bf5 | 672 | Double_t sigma3 = Sigma2VertexCount(tracks, 3); |
72e597d7 | 673 | for (Double_t nSigma = 0.1; nSigma < 5.05; nSigma += 0.1) |
674 | ratio2->Fill(nSigma, sigma3 / ratio->GetBinContent(ratio->FindBin(nSigma))); | |
9f469bf5 | 675 | ratio2->SetMarkerStyle(21); |
676 | ||
677 | canvas->cd(3); | |
4c6b34a8 | 678 | ratio2->DrawCopy("P"); |
9f469bf5 | 679 | |
680 | canvas->SaveAs("Sigma2Vertex.eps"); | |
4c6b34a8 | 681 | |
682 | return ratio2; | |
9f469bf5 | 683 | } |
684 | ||
4c6b34a8 | 685 | TH1F* Sigma2VertexSimulation() |
9f469bf5 | 686 | { |
687 | TFile* file = TFile::Open("systematics.root"); | |
688 | ||
689 | TH1F* sigmavertex = dynamic_cast<TH1F*> (file->Get("fSigmaVertex")); | |
690 | if (!sigmavertex) | |
691 | { | |
692 | printf("Could not read histogram\n"); | |
693 | return; | |
694 | } | |
695 | ||
5c495d37 | 696 | TH1F* ratio = new TH1F("sigmavertexsimulation_ratio", "sigmavertexsimulation_ratio;N#sigma;% included in 3 #sigma / % included in N#sigma", sigmavertex->GetNbinsX(), sigmavertex->GetXaxis()->GetXmin(), sigmavertex->GetXaxis()->GetXmax()); |
9f469bf5 | 697 | |
698 | for (Int_t i=1; i<=sigmavertex->GetNbinsX(); ++i) | |
699 | ratio->SetBinContent(i, sigmavertex->GetBinContent(sigmavertex->GetXaxis()->FindBin(3)) / sigmavertex->GetBinContent(i)); | |
700 | ||
701 | TCanvas* canvas = new TCanvas("Sigma2VertexSimulation", "Sigma2VertexSimulation", 1000, 500); | |
702 | canvas->Divide(2, 1); | |
703 | ||
704 | canvas->cd(1); | |
705 | sigmavertex->SetMarkerStyle(21); | |
706 | sigmavertex->Draw("P"); | |
707 | ||
708 | canvas->cd(2); | |
709 | ratio->SetMarkerStyle(21); | |
4c6b34a8 | 710 | ratio->DrawCopy("P"); |
711 | ||
712 | return ratio; | |
10ebe68d | 713 | } |
714 | ||
4c6b34a8 | 715 | void Sigma2VertexCompare() |
716 | { | |
717 | TH1F* ratio1 = Sigma2VertexGaussian(); | |
718 | TH1F* ratio2 = Sigma2VertexSimulation(); | |
719 | ||
720 | ratio1->SetStats(kFALSE); | |
721 | ratio2->SetStats(kFALSE); | |
722 | ||
723 | ratio1->SetMarkerStyle(0); | |
724 | ratio2->SetMarkerStyle(0); | |
725 | ||
0bd1f8a0 | 726 | ratio1->SetLineWidth(2); |
727 | ratio2->SetLineWidth(2); | |
728 | ||
729 | TLegend* legend = new TLegend(0.7, 0.8, 0.95, 0.95); | |
730 | legend->SetFillColor(0); | |
4c6b34a8 | 731 | legend->AddEntry(ratio1, "Gaussian"); |
732 | legend->AddEntry(ratio2, "Simulation"); | |
733 | ||
0bd1f8a0 | 734 | ratio2->SetTitle(""); |
735 | ratio2->GetYaxis()->SetTitleOffset(1.5); | |
736 | ratio2->GetXaxis()->SetRangeUser(2, 4); | |
4c6b34a8 | 737 | |
738 | TCanvas* canvas = new TCanvas("Sigma2VertexCompare", "Sigma2VertexCompare", 500, 500); | |
739 | InitPad(); | |
740 | ||
5c495d37 | 741 | ratio2->SetMarkerStyle(21); |
742 | ratio1->SetMarkerStyle(22); | |
743 | ||
744 | ratio2->GetYaxis()->SetRangeUser(0.8, 1.2); | |
4c6b34a8 | 745 | ratio2->SetLineColor(kRed); |
5c495d37 | 746 | ratio2->SetMarkerColor(kRed); |
747 | ratio2->Draw("PL"); | |
748 | ratio1->Draw("SAMEPL"); | |
4c6b34a8 | 749 | |
750 | legend->Draw(); | |
0bd1f8a0 | 751 | |
752 | canvas->SaveAs("Sigma2VertexCompare.eps"); | |
4c6b34a8 | 753 | } |
9f469bf5 | 754 | |
10ebe68d | 755 | void drawSystematics() |
756 | { | |
9f469bf5 | 757 | //Secondaries(); |
758 | //DrawDifferentSpecies(); | |
759 | //Composition(); | |
760 | ||
761 | Sigma2VertexSimulation(); | |
27c04e01 | 762 | |
10ebe68d | 763 | } |
7af955da | 764 | |
765 | void DrawdNdEtaDifferences() | |
766 | { | |
767 | TH1* hists[5]; | |
768 | ||
72e597d7 | 769 | TLegend* legend = new TLegend(0.3, 0.73, 0.70, 0.98); |
f1076daa | 770 | legend->SetFillColor(0); |
771 | ||
7af955da | 772 | TCanvas* canvas = new TCanvas("DrawdNdEtaDifferences", "DrawdNdEtaDifferences", 1000, 500); |
773 | canvas->Divide(2, 1); | |
774 | ||
775 | canvas->cd(1); | |
776 | ||
777 | for (Int_t i=0; i<5; ++i) | |
778 | { | |
72e597d7 | 779 | hists[i] = 0; |
7af955da | 780 | TFile* file = 0; |
f1076daa | 781 | TString title; |
7af955da | 782 | |
783 | switch(i) | |
784 | { | |
f1076daa | 785 | case 0 : file = TFile::Open("systematics_dndeta_reference.root"); title = "standard composition"; break; |
786 | case 1 : file = TFile::Open("systematics_dndeta_KBoosted.root"); title = "+ 50% kaons"; break; | |
787 | case 2 : file = TFile::Open("systematics_dndeta_KReduced.root"); title = "- 50% kaons"; break; | |
788 | case 3 : file = TFile::Open("systematics_dndeta_pBoosted.root"); title = "+ 50% protons"; break; | |
789 | case 4 : file = TFile::Open("systematics_dndeta_pReduced.root"); title = "- 50% protons"; break; | |
7af955da | 790 | default: return; |
791 | } | |
792 | ||
72e597d7 | 793 | if (file) |
794 | { | |
795 | hists[i] = (TH1*) file->Get("dndeta/dndeta_dNdEta_corrected_2"); | |
796 | hists[i]->SetTitle("a)"); | |
797 | ||
798 | Prepare1DPlot(hists[i], kFALSE); | |
799 | hists[i]->GetXaxis()->SetRangeUser(-0.7999, 0.7999); | |
da7acf97 | 800 | hists[i]->GetYaxis()->SetRangeUser(6, 7); |
72e597d7 | 801 | hists[i]->SetLineColor(i+1); |
802 | hists[i]->SetMarkerColor(i+1); | |
803 | hists[i]->GetXaxis()->SetLabelOffset(0.015); | |
804 | hists[i]->GetYaxis()->SetTitleOffset(1.5); | |
805 | gPad->SetLeftMargin(0.12); | |
806 | hists[i]->DrawCopy(((i > 0) ? "SAME" : "")); | |
807 | ||
808 | legend->AddEntry(hists[i], title); | |
809 | hists[i]->SetTitle(title); | |
810 | } | |
7af955da | 811 | } |
f1076daa | 812 | legend->Draw(); |
7af955da | 813 | |
814 | canvas->cd(2); | |
f1076daa | 815 | gPad->SetLeftMargin(0.14); |
816 | ||
817 | TLegend* legend2 = new TLegend(0.73, 0.73, 0.98, 0.98); | |
818 | legend2->SetFillColor(0); | |
7af955da | 819 | |
820 | for (Int_t i=1; i<5; ++i) | |
821 | { | |
72e597d7 | 822 | if (hists[i]) |
823 | { | |
824 | legend2->AddEntry(hists[i]); | |
825 | ||
826 | hists[i]->Divide(hists[0]); | |
827 | hists[i]->SetTitle("b)"); | |
da7acf97 | 828 | hists[i]->GetYaxis()->SetRangeUser(0.95, 1.05); |
72e597d7 | 829 | hists[i]->GetYaxis()->SetTitle("Ratio to standard composition"); |
830 | hists[i]->GetYaxis()->SetTitleOffset(1.8); | |
831 | hists[i]->DrawCopy(((i > 1) ? "SAME" : "")); | |
832 | } | |
7af955da | 833 | } |
f1076daa | 834 | |
835 | legend2->Draw(); | |
836 | ||
837 | canvas->SaveAs("particlecomposition_result.eps"); | |
838 | canvas->SaveAs("particlecomposition_result.gif"); | |
7af955da | 839 | } |
27c04e01 | 840 | |
da7acf97 | 841 | mergeCorrections4SystematicStudies(Char_t* standardCorrectionFileName="correction_map.root", |
27c04e01 | 842 | Char_t* systematicCorrectionFileName="systematics.root", |
da7acf97 | 843 | Char_t* outputFileName="systematics_vtxtrigger_compositions.root") { |
27c04e01 | 844 | // |
845 | // Function used to merge standard corrections with vertex | |
846 | // reconstruction corrections obtained by a certain mix of ND, DD | |
847 | // and SD events. | |
848 | // | |
849 | ||
850 | gSystem->Load("libPWG0base"); | |
851 | ||
da7acf97 | 852 | const Char_t* typeName[] = { "vertexreco", "trigger", "vtxtrigger" }; |
853 | const Char_t* changes[] = {"pythia","ddmore","ddless","sdmore","sdless", "dmore", "dless"}; | |
854 | Float_t scalesDD[] = {1.0, 1.5, 0.5, 1.0, 1.0, 1.5, 0.5}; | |
855 | Float_t scalesSD[] = {1.0, 1.0, 1.0, 1.5, 0.5, 1.5, 0.5}; | |
27c04e01 | 856 | |
857 | // cross section from Pythia | |
858 | Float_t sigmaND = 55.2; | |
859 | Float_t sigmaDD = 9.78; | |
860 | Float_t sigmaSD = 14.30; | |
861 | ||
da7acf97 | 862 | AlidNdEtaCorrection* corrections[21]; |
863 | for (Int_t j=0; j<3; j++) { // j = 0 (change vtx), j = 1 (change trg), j = 2 (change both) | |
864 | AlidNdEtaCorrection* correctionStandard = new AlidNdEtaCorrection("dndeta_correction","dndeta_correction"); | |
865 | correctionStandard->LoadHistograms(standardCorrectionFileName); | |
27c04e01 | 866 | |
da7acf97 | 867 | // dont take vertexreco from this one |
868 | correctionStandard->GetVertexRecoCorrection()->Reset(); | |
869 | // dont take triggerbias from this one | |
870 | correctionStandard->GetTriggerBiasCorrectionINEL()->Reset(); | |
27c04e01 | 871 | |
da7acf97 | 872 | for (Int_t i=0; i<7; i++) { |
873 | TString name; | |
874 | name.Form("dndeta_correction_syst_%s_%s", typeName[j], changes[i]); | |
875 | AlidNdEtaCorrection* current = new AlidNdEtaCorrection(name, name); | |
876 | ||
877 | name.Form("vertexRecoND"); | |
878 | AlidNdEtaCorrection* dNdEtaCorrectionVtxND = new AlidNdEtaCorrection(name,name); | |
879 | dNdEtaCorrectionVtxND->LoadHistograms(systematicCorrectionFileName, name); | |
880 | name.Form("vertexRecoDD"); | |
881 | AlidNdEtaCorrection* dNdEtaCorrectionVtxDD = new AlidNdEtaCorrection(name,name); | |
882 | dNdEtaCorrectionVtxDD->LoadHistograms(systematicCorrectionFileName, name); | |
883 | name.Form("vertexRecoSD"); | |
884 | AlidNdEtaCorrection* dNdEtaCorrectionVtxSD = new AlidNdEtaCorrection(name,name); | |
885 | dNdEtaCorrectionVtxSD->LoadHistograms(systematicCorrectionFileName, name); | |
886 | ||
887 | name.Form("triggerBiasND"); | |
888 | AlidNdEtaCorrection* dNdEtaCorrectionTriggerND = new AlidNdEtaCorrection(name,name); | |
889 | dNdEtaCorrectionTriggerND->LoadHistograms(systematicCorrectionFileName, name); | |
890 | name.Form("triggerBiasDD"); | |
891 | AlidNdEtaCorrection* dNdEtaCorrectionTriggerDD = new AlidNdEtaCorrection(name,name); | |
892 | dNdEtaCorrectionTriggerDD->LoadHistograms(systematicCorrectionFileName, name); | |
893 | name.Form("triggerBiasSD"); | |
894 | AlidNdEtaCorrection* dNdEtaCorrectionTriggerSD = new AlidNdEtaCorrection(name,name); | |
895 | dNdEtaCorrectionTriggerSD->LoadHistograms(systematicCorrectionFileName, name); | |
896 | ||
897 | // calculating relative | |
898 | Float_t nd = 100 * sigmaND/(sigmaND + (scalesDD[i]*sigmaDD) + (scalesDD[i]*sigmaSD)); | |
899 | Float_t dd = 100 * (scalesDD[i]*sigmaDD)/(sigmaND + (scalesDD[i]*sigmaDD) + (scalesDD[i]*sigmaSD)); | |
900 | Float_t sd = 100 * (scalesSD[i]*sigmaSD)/(sigmaND + (scalesDD[i]*sigmaDD) + (scalesDD[i]*sigmaSD)); | |
901 | ||
902 | printf(Form("%s : ND=%.1f\%, DD=%.1f\%, SD=%.1f\% \n",changes[i],nd,dd,sd)); | |
903 | current->SetTitle(Form("ND=%.2f\%,DD=%.2f\%,SD=%.2f\%",nd,dd,sd)); | |
904 | current->SetTitle(name); | |
905 | ||
906 | // scale | |
907 | if (j == 0 || j == 2) | |
908 | { | |
909 | dNdEtaCorrectionVtxDD->GetVertexRecoCorrection()->GetMeasuredHistogram()->Scale(scalesDD[i]); | |
910 | dNdEtaCorrectionVtxSD->GetVertexRecoCorrection()->GetMeasuredHistogram()->Scale(scalesSD[i]); | |
911 | dNdEtaCorrectionVtxDD->GetVertexRecoCorrection()->GetGeneratedHistogram()->Scale(scalesDD[i]); | |
912 | dNdEtaCorrectionVtxSD->GetVertexRecoCorrection()->GetGeneratedHistogram()->Scale(scalesSD[i]); | |
913 | } | |
914 | if (j == 1 || j == 2) | |
915 | { | |
916 | dNdEtaCorrectionTriggerDD->GetTriggerBiasCorrectionINEL()->GetMeasuredHistogram()->Scale(scalesDD[i]); | |
917 | dNdEtaCorrectionTriggerSD->GetTriggerBiasCorrectionINEL()->GetMeasuredHistogram()->Scale(scalesSD[i]); | |
918 | dNdEtaCorrectionTriggerDD->GetTriggerBiasCorrectionINEL()->GetGeneratedHistogram()->Scale(scalesDD[i]); | |
919 | dNdEtaCorrectionTriggerSD->GetTriggerBiasCorrectionINEL()->GetGeneratedHistogram()->Scale(scalesSD[i]); | |
920 | } | |
921 | ||
922 | //clear track, trigger in Vtx correction | |
923 | dNdEtaCorrectionVtxND->GetTrack2ParticleCorrection()->Reset(); | |
924 | dNdEtaCorrectionVtxND->GetTriggerBiasCorrectionNSD()->Reset(); | |
925 | dNdEtaCorrectionVtxND->GetTriggerBiasCorrectionND()->Reset(); | |
926 | dNdEtaCorrectionVtxND->GetTriggerBiasCorrectionINEL()->Reset(); | |
927 | dNdEtaCorrectionVtxSD->GetTrack2ParticleCorrection()->Reset(); | |
928 | dNdEtaCorrectionVtxSD->GetTriggerBiasCorrectionNSD()->Reset(); | |
929 | dNdEtaCorrectionVtxSD->GetTriggerBiasCorrectionND()->Reset(); | |
930 | dNdEtaCorrectionVtxSD->GetTriggerBiasCorrectionINEL()->Reset(); | |
931 | dNdEtaCorrectionVtxDD->GetTrack2ParticleCorrection()->Reset(); | |
932 | dNdEtaCorrectionVtxDD->GetTriggerBiasCorrectionNSD()->Reset(); | |
933 | dNdEtaCorrectionVtxDD->GetTriggerBiasCorrectionND()->Reset(); | |
934 | dNdEtaCorrectionVtxDD->GetTriggerBiasCorrectionINEL()->Reset(); | |
935 | ||
936 | //clear track, vertexreco in trigger correction | |
937 | dNdEtaCorrectionTriggerND->GetTrack2ParticleCorrection()->Reset(); | |
938 | dNdEtaCorrectionTriggerND->GetTriggerBiasCorrectionNSD()->Reset(); | |
939 | dNdEtaCorrectionTriggerND->GetTriggerBiasCorrectionND()->Reset(); | |
940 | dNdEtaCorrectionTriggerND->GetVertexRecoCorrection()->Reset(); | |
941 | dNdEtaCorrectionTriggerSD->GetTrack2ParticleCorrection()->Reset(); | |
942 | dNdEtaCorrectionTriggerSD->GetTriggerBiasCorrectionNSD()->Reset(); | |
943 | dNdEtaCorrectionTriggerSD->GetTriggerBiasCorrectionND()->Reset(); | |
944 | dNdEtaCorrectionTriggerSD->GetVertexRecoCorrection()->Reset(); | |
945 | dNdEtaCorrectionTriggerDD->GetTrack2ParticleCorrection()->Reset(); | |
946 | dNdEtaCorrectionTriggerDD->GetTriggerBiasCorrectionNSD()->Reset(); | |
947 | dNdEtaCorrectionTriggerDD->GetTriggerBiasCorrectionND()->Reset(); | |
948 | dNdEtaCorrectionTriggerDD->GetVertexRecoCorrection()->Reset(); | |
949 | ||
950 | TList collection; | |
951 | collection.Add(correctionStandard); | |
952 | collection.Add(dNdEtaCorrectionVtxND); | |
953 | collection.Add(dNdEtaCorrectionVtxDD); | |
954 | collection.Add(dNdEtaCorrectionVtxSD); | |
955 | collection.Add(dNdEtaCorrectionTriggerND); | |
956 | collection.Add(dNdEtaCorrectionTriggerDD); | |
957 | collection.Add(dNdEtaCorrectionTriggerSD); | |
958 | ||
959 | current->Merge(&collection); | |
960 | current->Finish(); | |
961 | ||
962 | corrections[i+j*7] = current; | |
963 | } | |
27c04e01 | 964 | } |
da7acf97 | 965 | |
27c04e01 | 966 | TFile* fout = new TFile(outputFileName,"RECREATE"); |
967 | ||
da7acf97 | 968 | for (Int_t i=0; i<21; i++) |
27c04e01 | 969 | corrections[i]->SaveHistograms(); |
970 | ||
971 | fout->Write(); | |
972 | fout->Close(); | |
973 | } | |
5e08578b | 974 | |
975 | ||
7f0efb86 | 976 | DrawVertexRecoSyst() { |
977 | // Draws the ratio of the dN/dEta obtained with changed SD and DD | |
978 | // cross-sections vertex reco corrections to the dN/dEta obtained | |
979 | // from the standard pythia cross-sections | |
980 | // | |
981 | // The files with the vertex reco corrections for different | |
982 | // processes (and with the other standard corrections) are expected | |
983 | // to have the names "analysis_esd_X.root", where the Xs are defined | |
984 | // in the array changes below. | |
5e08578b | 985 | |
986 | Char_t* changes[] = {"pythia","ddmore","ddless","sdmore","sdless", "dmore", "dless"}; | |
987 | Char_t* descr[] = {"", | |
988 | "#sigma_{DD}' = 1.5#times#sigma_{DD}", | |
989 | "#sigma_{DD}' = 0.5#times#sigma_{DD}", | |
990 | "#sigma_{SD}' = 1.5#times#sigma_{SD}", | |
991 | "#sigma_{SD}' = 0.5#times#sigma_{SD}", | |
992 | "#sigma_{D}' = 1.5#times#sigma_{D}", | |
993 | "#sigma_{D}' = 0.5#times#sigma_{D}"}; | |
994 | ||
995 | Float_t scalesDD[] = {1.0, 1.5, 0.5, 1.5, 0.5, 1.5, 0.5}; | |
996 | Float_t scalesSD[] = {1.0, 1.0, 1.0, 1.0, 1.0, 1.5, 0.5}; | |
997 | ||
fb30f6e3 | 998 | Int_t markers[] = {20,20,21,22,23,28,29}; |
5e08578b | 999 | Int_t colors[] = {1,2,3,4,6,8,102}; |
1000 | ||
1001 | // cross section from Pythia | |
1002 | Float_t sigmaND = 55.2; | |
1003 | Float_t sigmaDD = 9.78; | |
1004 | Float_t sigmaSD = 14.30; | |
1005 | ||
5e08578b | 1006 | TH1F* dNdEta[7]; |
5e08578b | 1007 | TH1F* ratios[7]; |
1008 | ||
1009 | TFile* fin; | |
1010 | ||
5e08578b | 1011 | for (Int_t i=0; i<7; i++) { |
fb30f6e3 | 1012 | // calculating relative |
7f0efb86 | 1013 | fin = TFile::Open(Form("analysis_esd_%s.root",changes[i])); |
fb30f6e3 | 1014 | |
5e08578b | 1015 | dNdEta[i] = (TH1F*)(fin->Get("dndeta/dndeta_dNdEta_corrected_2"))->Clone(); |
fb30f6e3 | 1016 | |
5e08578b | 1017 | for (Int_t b=0; b<dNdEta[i]->GetNbinsX(); b++) { |
1018 | if (TMath::Abs(dNdEta[i]->GetBinCenter(b))>0.9) { | |
1019 | dNdEta[i]->SetBinContent(b,0); | |
1020 | dNdEta[i]->SetBinError(b,0); | |
1021 | } | |
1022 | } | |
fb30f6e3 | 1023 | |
5e08578b | 1024 | dNdEta[i]->Rebin(4); |
fb30f6e3 | 1025 | |
5e08578b | 1026 | dNdEta[i]->SetLineWidth(2); |
1027 | dNdEta[i]->SetLineColor(colors[i]); | |
1028 | dNdEta[i]->SetMarkerStyle(markers[i]); | |
fb30f6e3 | 1029 | dNdEta[i]->SetMarkerSize(0.9); |
1030 | dNdEta[i]->SetMarkerColor(colors[i]); | |
1031 | ||
5e08578b | 1032 | ratios[i] = (TH1F*)dNdEta[i]->Clone("ratio"); |
1033 | ratios[i]->Divide(ratios[i],dNdEta[0],1,1,"B"); | |
1034 | ||
1035 | ratios[i]->SetName(changes[i]); | |
1036 | ratios[i]->SetTitle(changes[i]); | |
1037 | } | |
1038 | ||
1039 | //########################################################## | |
1040 | ||
1041 | gStyle->SetOptStat(0); | |
1042 | gStyle->SetOptTitle(0); | |
1043 | gStyle->SetOptFit(0); | |
1044 | ||
1045 | gStyle->SetTextSize(0.2); | |
1046 | gStyle->SetTitleSize(0.05,"xyz"); | |
5e08578b | 1047 | gStyle->SetLabelOffset(0.01, "xyz"); |
1048 | ||
1049 | ||
1050 | gStyle->SetTitleOffset(1.2, "y"); | |
1051 | gStyle->SetTitleOffset(1.2, "x"); | |
1052 | gStyle->SetEndErrorSize(0.0); | |
1053 | ||
1054 | //############################################## | |
1055 | ||
1056 | //making canvas and pads | |
da7acf97 | 1057 | TCanvas *c = new TCanvas(Form("vertex_reco_syst_%s", plot), Form("vertex_reco_syst_%s", plot),600,500); |
5e08578b | 1058 | |
1059 | TPad* p1 = new TPad("pad1","", 0, 0.0, 1.0, 1.0, 0, 0, 0); | |
1060 | ||
1061 | p1->SetBottomMargin(0.15); | |
1062 | p1->SetTopMargin(0.03); | |
1063 | p1->SetLeftMargin(0.15); | |
1064 | p1->SetRightMargin(0.03); | |
da7acf97 | 1065 | |
5e08578b | 1066 | p1->SetGridx(); |
1067 | p1->SetGridy(); | |
1068 | ||
1069 | p1->Draw(); | |
1070 | p1->cd(); | |
5e08578b | 1071 | |
1072 | ||
1073 | TH2F* null = new TH2F("","",100,-1.5,1.5,100,0.9601,1.099); | |
1074 | null->SetXTitle("#eta"); | |
1075 | null->GetXaxis()->CenterTitle(kTRUE); | |
1076 | null->SetYTitle("dN/d#eta / dN/d#eta_{pythia}"); | |
1077 | null->GetYaxis()->CenterTitle(kTRUE); | |
1078 | null->Draw(); | |
1079 | ||
1080 | for (Int_t i=1; i<7; i++) | |
1081 | ratios[i]->Draw("same"); | |
fb30f6e3 | 1082 | |
1083 | TLegend* legend = new TLegend(0.6, 0.6, 0.95, 0.95); | |
1084 | legend->SetFillColor(0); | |
1085 | ||
5e08578b | 1086 | TLatex* text[7]; |
1087 | for (Int_t i=1; i<7; i++) { | |
fb30f6e3 | 1088 | legend->AddEntry(dNdEta[i], descr[i]); |
5e08578b | 1089 | } |
fb30f6e3 | 1090 | |
1091 | legend->Draw(); | |
da7acf97 | 1092 | //text(0.2,0.88,"Effect of changing",0.045,1,kTRUE); |
1093 | //text(0.2,0.83,"relative cross-sections",0.045,1,kTRUE); | |
1094 | //text(0.2,0.78,"(vertex reconstruction corr.)",0.043,13,kTRUE); | |
1095 | ||
1096 | c->SaveAs(Form("%s.gif", c->GetName())); | |
fb30f6e3 | 1097 | c->SaveAs(Form("%s.eps", c->GetName())); |
5e08578b | 1098 | } |
1099 | ||
1100 | ||
1101 | ||
1102 | DrawTriggerEfficiency(Char_t* fileName) { | |
1103 | ||
1104 | gStyle->SetOptStat(0); | |
1105 | gStyle->SetOptTitle(0); | |
1106 | gStyle->SetOptFit(0); | |
1107 | ||
1108 | gStyle->SetTextSize(0.04); | |
1109 | gStyle->SetTitleSize(0.05,"xyz"); | |
1110 | //gStyle->SetTitleFont(133, "xyz"); | |
1111 | //gStyle->SetLabelFont(133, "xyz"); | |
1112 | //gStyle->SetLabelSize(17, "xyz"); | |
1113 | gStyle->SetLabelOffset(0.01, "xyz"); | |
1114 | ||
1115 | gStyle->SetTitleOffset(1.1, "y"); | |
1116 | gStyle->SetTitleOffset(1.1, "x"); | |
1117 | gStyle->SetEndErrorSize(0.0); | |
1118 | ||
1119 | //############################################## | |
1120 | ||
1121 | //making canvas and pads | |
1122 | TCanvas *c = new TCanvas("trigger_eff", "",600,500); | |
1123 | ||
1124 | TPad* p1 = new TPad("pad1","", 0, 0.0, 1.0, 1.0, 0, 0, 0); | |
1125 | ||
1126 | p1->SetBottomMargin(0.15); | |
1127 | p1->SetTopMargin(0.03); | |
1128 | p1->SetLeftMargin(0.15); | |
1129 | p1->SetRightMargin(0.03); | |
1130 | ||
1131 | p1->SetGridx(); | |
1132 | p1->SetGridy(); | |
1133 | ||
1134 | p1->Draw(); | |
1135 | p1->cd(); | |
1136 | ||
1137 | gSystem->Load("libPWG0base"); | |
1138 | ||
1139 | AlidNdEtaCorrection* corrections[4]; | |
1140 | AliCorrectionMatrix2D* triggerBiasCorrection[4]; | |
1141 | ||
1142 | TH1F* hTriggerEffInv[4]; | |
1143 | TH1F* hTriggerEff[4]; | |
1144 | ||
1145 | Char_t* names[] = {"triggerBiasND", "triggerBiasDD", "triggerBiasSD", "triggerBiasINEL"}; | |
1146 | ||
1147 | for (Int_t i=0; i<4; i++) { | |
1148 | corrections[i] = new AlidNdEtaCorrection(names[i], names[i]); | |
1149 | corrections[i]->LoadHistograms(fileName, names[i]); | |
1150 | ||
1151 | triggerBiasCorrection[i] = corrections[i]->GetTriggerBiasCorrectionINEL(); | |
1152 | ||
1153 | ||
1154 | hTriggerEffInv[i] = triggerBiasCorrection[i]->Get1DCorrection(); | |
1155 | hTriggerEff[i] = (TH1F*)hTriggerEffInv[i]->Clone(); | |
1156 | ||
1157 | for (Int_t b=0; b<=hTriggerEff[i]->GetNbinsX(); b++) { | |
1158 | hTriggerEff[i]->SetBinContent(b,1); | |
1159 | hTriggerEff[i]->SetBinError(b,0); | |
1160 | } | |
1161 | hTriggerEff[i]->Divide(hTriggerEff[i],hTriggerEffInv[i]); | |
1162 | hTriggerEff[i]->Scale(100); | |
1163 | } | |
1164 | ||
1165 | Int_t colors[] = {2,3,4,1}; | |
1166 | Float_t effs[4]; | |
1167 | for (Int_t i=0; i<4; i++) { | |
1168 | hTriggerEff[i]->Fit("pol0","","",-20,20); | |
1169 | effs[i] = ((TF1*)hTriggerEff[i]->GetListOfFunctions()->At(0))->GetParameter(0); | |
1170 | ((TF1*)hTriggerEff[i]->GetListOfFunctions()->At(0))->SetLineWidth(1); | |
1171 | ((TF1*)hTriggerEff[i]->GetListOfFunctions()->At(0))->SetLineStyle(2); | |
1172 | ((TF1*)hTriggerEff[i]->GetListOfFunctions()->At(0))->SetLineColor(colors[i]); | |
1173 | cout << effs[i] << endl; | |
1174 | } | |
1175 | ||
1176 | ||
1177 | Char_t* text[] = {"ND", "DD", "SD", "INEL"}; | |
1178 | TLatex* latex[4]; | |
1179 | ||
1180 | TH2F* null = new TH2F("","",100,-25,35,100,0,110); | |
1181 | null->SetXTitle("Vertex z [cm]"); | |
1182 | null->GetXaxis()->CenterTitle(kTRUE); | |
1183 | null->SetYTitle("Trigger efficiency [%]"); | |
1184 | null->GetYaxis()->CenterTitle(kTRUE); | |
1185 | null->Draw(); | |
1186 | ||
1187 | ||
1188 | for (Int_t i=0; i<4; i++) { | |
1189 | hTriggerEff[i]->SetLineWidth(2); | |
1190 | hTriggerEff[i]->SetLineColor(colors[i]); | |
1191 | ||
1192 | hTriggerEff[i]->Draw("same"); | |
1193 | ||
1194 | latex[i] = new TLatex(22,effs[i]-1.5, Form("%s (%0.1f)",text[i],effs[i])); | |
1195 | latex[i]->SetTextColor(colors[i]); | |
1196 | latex[i]->Draw(); | |
1197 | } | |
1198 | ||
1199 | } | |
1200 | ||
1201 | ||
1202 | DrawSpectraPID(Char_t* fileName) { | |
1203 | ||
1204 | gSystem->Load("libPWG0base"); | |
1205 | ||
1206 | Char_t* names[] = {"correction_0", "correction_1", "correction_2", "correction_3"}; | |
1207 | AlidNdEtaCorrection* corrections[4]; | |
1208 | AliCorrectionMatrix3D* trackToPartCorrection[4]; | |
1209 | ||
1210 | TH1F* measuredPt[4]; | |
1211 | TH1F* generatedPt[4]; | |
1212 | TH1F* ratioPt[4]; | |
1213 | ||
1214 | for (Int_t i=0; i<4; i++) { | |
1215 | corrections[i] = new AlidNdEtaCorrection(names[i], names[i]); | |
1216 | corrections[i]->LoadHistograms(fileName, names[i]); | |
1217 | ||
1218 | trackToPartCorrection[i] = corrections[i]->GetTrack2ParticleCorrection(); | |
1219 | ||
1220 | Int_t binX1 = (TH1F*)trackToPartCorrection[i]->GetMeasuredHistogram()->GetXaxis()->FindBin(-10); | |
1221 | Int_t binX2 = (TH1F*)trackToPartCorrection[i]->GetMeasuredHistogram()->GetXaxis()->FindBin(10); | |
1222 | Int_t binY1 = (TH1F*)trackToPartCorrection[i]->GetMeasuredHistogram()->GetYaxis()->FindBin(-1); | |
1223 | Int_t binY2 = (TH1F*)trackToPartCorrection[i]->GetMeasuredHistogram()->GetYaxis()->FindBin(1); | |
1224 | ||
1225 | measuredPt[i] = (TH1F*)trackToPartCorrection[i]->GetMeasuredHistogram()->ProjectionZ(Form("m_%d",i),binX1,binX2,binY1,binY2); | |
1226 | generatedPt[i] = (TH1F*)trackToPartCorrection[i]->GetGeneratedHistogram()->ProjectionZ(Form("g_%d",i),binX1,binX2,binY1,binY2); | |
1227 | ratioPt[i] = (TH1F*)generatedPt[i]->Clone(Form("r_%d",i)); | |
1228 | ratioPt[i]->Divide(measuredPt[i], generatedPt[i], 1,1,"B"); | |
1229 | } | |
1230 | ||
1231 | ratioPt[0]->Draw(); | |
1232 | ||
1233 | for (Int_t i=0; i<3; i++) { | |
1234 | ratioPt[i]->SetLineColor(i+1); | |
1235 | ratioPt[i]->SetLineWidth(2); | |
1236 | ||
1237 | ratioPt[i]->Draw("same"); | |
1238 | ||
1239 | } | |
1240 | ||
1241 | return; | |
1242 | measuredPt[0]->SetLineColor(2); | |
1243 | measuredPt[0]->SetLineWidth(5); | |
1244 | ||
1245 | measuredPt[0]->Draw(); | |
1246 | generatedPt[0]->Draw("same"); | |
1247 | } | |
da7acf97 | 1248 | |
1249 | void changePtSpectrum() | |
1250 | { | |
1251 | TFile* file = TFile::Open("pt.root"); | |
1252 | TH1F* hist = dynamic_cast<TH1F*> (file->Get("dndeta_check_pt")); | |
1253 | ||
1254 | hist->Rebin(3); | |
1255 | hist->Scale(1.0/3); | |
1256 | ||
1257 | TH1F* clone1 = dynamic_cast<TH1F*> (hist->Clone("clone1")); | |
1258 | TH1F* clone2 = dynamic_cast<TH1F*> (hist->Clone("clone2")); | |
1259 | ||
1260 | TH1F* scale1 = dynamic_cast<TH1F*> (hist->Clone("scale1")); | |
1261 | TH1F* scale2 = dynamic_cast<TH1F*> (hist->Clone("scale2")); | |
1262 | ||
1263 | Float_t ptCutOff = 0.3; | |
1264 | ||
1265 | for (Int_t i=1; i <= hist->GetNbinsX(); ++i) | |
1266 | { | |
1267 | if (hist->GetBinCenter(i) > ptCutOff) | |
1268 | { | |
1269 | scale1->SetBinContent(i, 1); | |
1270 | scale2->SetBinContent(i, 1); | |
1271 | } | |
1272 | else | |
1273 | { | |
1274 | // 90 % at pt = 0, 0% at pt = ptcutoff | |
1275 | scale1->SetBinContent(i, 1 - (ptCutOff - hist->GetBinCenter(i)) / ptCutOff * 0.3); | |
1276 | ||
1277 | // 110% at pt = 0, ... | |
1278 | scale2->SetBinContent(i, 1 + (ptCutOff - hist->GetBinCenter(i)) / ptCutOff * 0.3); | |
1279 | } | |
1280 | } | |
1281 | ||
1282 | new TCanvas; | |
1283 | scale1->Draw(); | |
1284 | scale2->SetLineColor(kRed); | |
1285 | scale2->Draw("SAME"); | |
1286 | ||
1287 | clone1->Multiply(scale1); | |
1288 | clone2->Multiply(scale2); | |
1289 | ||
1290 | Prepare1DPlot(hist); | |
1291 | Prepare1DPlot(clone1); | |
1292 | Prepare1DPlot(clone2); | |
1293 | ||
1294 | hist->SetTitle(";p_{T} in GeV/c;dN/dp_{T} in c/GeV"); | |
1295 | hist->GetXaxis()->SetRangeUser(0, 0.7); | |
1296 | hist->GetYaxis()->SetRangeUser(0, clone2->GetMaximum() * 1.1); | |
1297 | hist->GetYaxis()->SetTitleOffset(1); | |
1298 | ||
1299 | TCanvas* canvas = new TCanvas; | |
1300 | hist->Draw(); | |
1301 | clone1->SetLineColor(kRed); | |
1302 | clone1->Draw("SAME"); | |
1303 | clone2->SetLineColor(kBlue); | |
1304 | clone2->Draw("SAME"); | |
1305 | ||
1306 | Float_t fraction = hist->Integral(hist->GetXaxis()->FindBin(ptCutOff), hist->GetNbinsX()) / hist->Integral(1, hist->GetNbinsX()); | |
1307 | Float_t fraction1 = clone1->Integral(clone1->GetXaxis()->FindBin(ptCutOff), clone1->GetNbinsX()) / clone1->Integral(1, clone1->GetNbinsX()); | |
1308 | Float_t fraction2 = clone2->Integral(clone2->GetXaxis()->FindBin(ptCutOff), clone2->GetNbinsX()) / clone2->Integral(1, clone2->GetNbinsX()); | |
1309 | ||
1310 | printf("%f %f %f\n", fraction, fraction1, fraction2); | |
1311 | printf("Rel. %f %f\n", fraction1 / fraction, fraction2 / fraction); | |
1312 | ||
1313 | canvas->SaveAs("changePtSpectrum.gif"); | |
1314 | } | |
fb30f6e3 | 1315 | |
1316 | void FractionBelowPt() | |
1317 | { | |
1318 | gSystem->Load("libPWG0base"); | |
1319 | ||
1320 | AlidNdEtaCorrection* fdNdEtaCorrection[4]; | |
1321 | ||
1322 | TFile::Open("systematics.root"); | |
1323 | ||
1324 | for (Int_t i=0; i<4; ++i) | |
1325 | { | |
1326 | TString name; | |
1327 | name.Form("correction_%d", i); | |
1328 | fdNdEtaCorrection[i] = new AlidNdEtaCorrection(name, name); | |
1329 | fdNdEtaCorrection[i]->LoadHistograms("systematics.root", name); | |
1330 | } | |
1331 | ||
1332 | Double_t geneCount[5]; | |
1333 | Double_t measCount[5]; | |
1334 | geneCount[4] = 0; | |
1335 | measCount[4] = 0; | |
1336 | ||
1337 | for (Int_t i=0; i<4; ++i) | |
1338 | { | |
1339 | TH3F* hist = (TH3F*) fdNdEtaCorrection[i]->GetTrack2ParticleCorrection()->GetGeneratedHistogram(); | |
1340 | geneCount[i] = hist->Integral(hist->GetXaxis()->FindBin(-10), hist->GetXaxis()->FindBin(10), | |
1341 | hist->GetYaxis()->FindBin(-0.8), hist->GetYaxis()->FindBin(0.8), | |
1342 | 1, hist->GetZaxis()->FindBin(0.3)); | |
1343 | ||
1344 | hist = (TH3F*) fdNdEtaCorrection[i]->GetTrack2ParticleCorrection()->GetMeasuredHistogram(); | |
1345 | measCount[i] = hist->Integral(hist->GetXaxis()->FindBin(-10), hist->GetXaxis()->FindBin(10), hist->GetYaxis()->FindBin(-0.8), hist->GetYaxis()->FindBin(0.8), 1, hist->GetZaxis()->FindBin(0.3)); | |
1346 | ||
1347 | geneCount[4] += geneCount[i]; | |
1348 | measCount[4] += measCount[i]; | |
1349 | ||
1350 | printf("Particle %s: %d gene, %d meas\n", ((i == 0) ? "pi" : (i == 1) ? "K" : (i == 2) ? "p" : "others"), (Int_t) geneCount[i], (Int_t) measCount[i]); | |
1351 | } | |
1352 | ||
1353 | printf("Generated ratios are: %f pi, %f K, %f p, %f others\n", geneCount[0] / geneCount[4], geneCount[1] / geneCount[4], geneCount[2] / geneCount[4], geneCount[3] / geneCount[4]); | |
1354 | ||
1355 | printf("Reconstructed ratios are: %f pi, %f K, %f p, %f others\n", measCount[0] / measCount[4], measCount[1] / measCount[4], measCount[2] / measCount[4], measCount[3] / measCount[4]); | |
7f0efb86 | 1356 | } |