]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG0/dNdEta/drawSystematics.C
Not needed.
[u/mrichter/AliRoot.git] / PWG0 / dNdEta / drawSystematics.C
CommitLineData
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
20extern TPad* gPad;
21
22void Track2Particle1DCreatePlots(const char* fileName = "correction_map.root", const char* folderName = "dndeta_correction", Float_t upperPtLimit = 10);
10ebe68d 23
24#endif
25
8ca1a6d9 26Int_t markers[] = {20,20,21,22,23,28,29};
27Int_t colors[] = {1,2,3,4,6,8,102};
9f469bf5 28
dd367a14 29void loadlibs()
30{
31 gSystem->Load("libTree");
32 gSystem->Load("libVMC");
33
34 gSystem->Load("libSTEERBase");
35 gSystem->Load("libANALYSIS");
36 gSystem->Load("libPWG0base");
37}
38
10ebe68d 39void SetRanges(TAxis* axis)
40{
41 if (strcmp(axis->GetTitle(), "#eta") == 0)
42 axis->SetRangeUser(-1.7999, 1.7999);
43 if (strcmp(axis->GetTitle(), "p_{T} [GeV/c]") == 0)
44 axis->SetRangeUser(0, 9.9999);
45 if (strcmp(axis->GetTitle(), "vtx z [cm]") == 0)
46 axis->SetRangeUser(-15, 14.9999);
47 if (strcmp(axis->GetTitle(), "Ntracks") == 0)
48 axis->SetRangeUser(0, 99.9999);
49}
50
51void SetRanges(TH1* hist)
52{
53 SetRanges(hist->GetXaxis());
54 SetRanges(hist->GetYaxis());
55 SetRanges(hist->GetZaxis());
56}
57
58void Prepare3DPlot(TH3* hist)
59{
60 hist->GetXaxis()->SetTitleOffset(1.5);
61 hist->GetYaxis()->SetTitleOffset(1.5);
62 hist->GetZaxis()->SetTitleOffset(1.5);
63
64 hist->SetStats(kFALSE);
65}
66
67void Prepare2DPlot(TH2* hist)
68{
69 hist->SetStats(kFALSE);
70 hist->GetYaxis()->SetTitleOffset(1.4);
71
72 SetRanges(hist);
73}
74
7af955da 75void Prepare1DPlot(TH1* hist, Bool_t setRanges = kTRUE)
10ebe68d 76{
77 hist->SetLineWidth(2);
78 hist->SetStats(kFALSE);
79
7af955da 80 hist->GetXaxis()->SetTitleOffset(1.2);
81 hist->GetYaxis()->SetTitleOffset(1.2);
82
83 if (setRanges)
84 SetRanges(hist);
10ebe68d 85}
86
87void InitPad()
88{
9f469bf5 89 if (!gPad)
90 return;
91
10ebe68d 92 gPad->Range(0, 0, 1, 1);
93 gPad->SetLeftMargin(0.15);
94 //gPad->SetRightMargin(0.05);
95 //gPad->SetTopMargin(0.13);
96 //gPad->SetBottomMargin(0.1);
97
72e597d7 98 gPad->SetGridx();
99 gPad->SetGridy();
10ebe68d 100}
101
102void InitPadCOLZ()
103{
104 gPad->Range(0, 0, 1, 1);
105 gPad->SetRightMargin(0.15);
106 gPad->SetLeftMargin(0.12);
107
108 gPad->SetGridx();
109 gPad->SetGridy();
110}
111
112void Secondaries()
113{
114 TFile* file = TFile::Open("systematics.root");
115
7af955da 116 TH2F* secondaries = dynamic_cast<TH2F*> (file->Get("fSecondaries"));
10ebe68d 117 if (!secondaries)
118 {
119 printf("Could not read histogram\n");
120 return;
121 }
122
7af955da 123 TCanvas* canvas = new TCanvas("Secondaries", "Secondaries", 1000, 1000);
124 canvas->Divide(3, 3);
125 for (Int_t i=1; i<=8; i++)
10ebe68d 126 {
7af955da 127 TH1D* hist = secondaries->ProjectionY(Form("proj_%d", i), i, i);
128 hist->SetTitle(secondaries->GetXaxis()->GetBinLabel(i));
10ebe68d 129
7af955da 130 canvas->cd(i);
131 hist->Draw();
10ebe68d 132 }
133}
134
6b7fa615 135void Track2Particle1DComposition(const char** fileNames, Int_t folderCount, const char** folderNames, Float_t upperPtLimit = 9.9)
10ebe68d 136{
137 gSystem->Load("libPWG0base");
138
9f469bf5 139 TString canvasName;
140 canvasName.Form("Track2Particle1DComposition");
141 TCanvas* canvas = new TCanvas(canvasName, canvasName, 1200, 400);
142 canvas->Divide(3, 1);
10ebe68d 143
9f469bf5 144 TLegend* legend = new TLegend(0.7, 0.7, 0.95, 0.95);
10ebe68d 145
9f469bf5 146 for (Int_t i=0; i<folderCount; ++i)
10ebe68d 147 {
6b7fa615 148 Track2Particle1DCreatePlots(fileNames[i], folderNames[i], upperPtLimit);
9f469bf5 149
6b7fa615 150 TH1* corrX = dynamic_cast<TH1*> (gROOT->FindObject(Form("gene_%s_nTrackToNPart_x_div_meas_%s_nTrackToNPart_x", folderNames[i], folderNames[i])));
151 TH1* corrY = dynamic_cast<TH1*> (gROOT->FindObject(Form("gene_%s_nTrackToNPart_y_div_meas_%s_nTrackToNPart_y", folderNames[i], folderNames[i])));
152 TH1* corrZ = dynamic_cast<TH1*> (gROOT->FindObject(Form("gene_%s_nTrackToNPart_z_div_meas_%s_nTrackToNPart_z", folderNames[i], folderNames[i])));
9f469bf5 153
154 Prepare1DPlot(corrX);
155 Prepare1DPlot(corrY);
156 Prepare1DPlot(corrZ);
157
158 const char* title = "Track2Particle Correction";
159 corrX->SetTitle(title);
160 corrY->SetTitle(title);
161 corrZ->SetTitle(title);
162
163 corrZ->GetXaxis()->SetRangeUser(0, upperPtLimit);
164
165 corrX->SetLineColor(i+1);
166 corrY->SetLineColor(i+1);
167 corrZ->SetLineColor(i+1);
168
169 canvas->cd(1);
170 InitPad();
171 corrX->DrawCopy(((i>0) ? "SAME" : ""));
172
173 canvas->cd(2);
174 InitPad();
175 corrY->DrawCopy(((i>0) ? "SAME" : ""));
176
177 canvas->cd(3);
178 InitPad();
179 corrZ->DrawCopy(((i>0) ? "SAME" : ""));
180
181 legend->AddEntry(corrZ, folderNames[i]);
10ebe68d 182 }
183
9f469bf5 184 legend->Draw();
185
186 //canvas->SaveAs(Form("Track2Particle1D_%s_%d_%f.gif", fileName, gMax, upperPtLimit));
187 //canvas->SaveAs(Form("Track2Particle1D_%s_%d_%f.eps", fileName, gMax, upperPtLimit));
188}
189
6b7fa615 190TH1** DrawRatios(const char* fileName = "systematics.root")
191{
192 gSystem->Load("libPWG0base");
193
194 TFile* file = TFile::Open(fileName);
195
196 TString canvasName;
197 canvasName.Form("DrawRatios");
198 TCanvas* canvas = new TCanvas(canvasName, canvasName, 800, 400);
199 canvas->Divide(2, 1);
200 canvas->cd(1);
201
202 TH1** ptDists = new TH1*[4];
203
7af955da 204 TLegend* legend = new TLegend(0.73, 0.73, 0.98, 0.98);
6b7fa615 205
206 const char* folderNames[] = { "correction_0", "correction_1", "correction_2", "correction_3" };
207 const char* particleNames[] = { "#pi", "K", "p", "other" };
208 for (Int_t i=0; i<4; ++i)
209 {
210 AlidNdEtaCorrection* dNdEtaCorrection = new AlidNdEtaCorrection(folderNames[i], folderNames[i]);
211 dNdEtaCorrection->LoadHistograms(fileName, folderNames[i]);
212
213 TH3F* gene = dNdEtaCorrection->GetTrack2ParticleCorrection()->GetGeneratedHistogram();
214
215 gene->GetYaxis()->SetRangeUser(-0.8, 0.8);
216 gene->GetXaxis()->SetRangeUser(-10, 10);
217
218 ptDists[i] = dynamic_cast<TH1*> (gene->Project3D("z")->Clone(Form("%s_z", folderNames[i])));
219 ptDists[i]->SetTitle(";p_{T};Count");
220 if (!ptDists[i])
221 {
222 printf("Problem getting distribution %d.\n", i);
223 return 0;
224 }
225
226 ptDists[i]->GetYaxis()->SetRangeUser(1, ptDists[i]->GetMaximum()*1.1);
227 ptDists[i]->GetXaxis()->SetRangeUser(0, 9.9);
228 ptDists[i]->SetLineColor(i+1);
229 ptDists[i]->DrawCopy((i == 0) ? "" : "SAME");
230 ptDists[i]->GetYaxis()->SetRange(0, 0);
231
232 legend->AddEntry(ptDists[i], particleNames[i]);
233 }
234 gPad->SetLogy();
235
236 TH1* total = dynamic_cast<TH1*> (ptDists[0]->Clone("total"));
237
238 for (Int_t i=1; i<4; ++i)
239 total->Add(ptDists[i]);
240
241 canvas->cd(2);
242 for (Int_t i=0; i<4; ++i)
243 {
244 ptDists[i]->Divide(total);
7af955da 245 ptDists[i]->SetStats(kFALSE);
246 ptDists[i]->SetTitle(";p_{T};Fraction of total");
6b7fa615 247 ptDists[i]->GetYaxis()->SetRangeUser(0, 1);
248 ptDists[i]->Draw((i == 0) ? "" : "SAME");
249 }
7af955da 250 legend->SetFillColor(0);
6b7fa615 251 legend->Draw();
252
253 canvas->SaveAs("DrawRatios.gif");
254
7af955da 255
256 canvas = new TCanvas("PythiaRatios", "PythiaRatios", 500, 500);
257 for (Int_t i=0; i<4; ++i)
258 {
259 TH1* hist = ptDists[i]->Clone();
260 hist->GetXaxis()->SetRangeUser(0, 1.9);
261 hist->Draw((i == 0) ? "" : "SAME");
262 }
263 legend->Draw();
264
265 canvas->SaveAs("pythiaratios.eps");
266
6b7fa615 267 file->Close();
268
269 return ptDists;
270}
271
7af955da 272void DrawCompareToReal()
273{
274 gROOT->ProcessLine(".L drawPlots.C");
275
276 const char* fileNames[] = { "correction_map.root", "new_compositions.root" };
277 const char* folderNames[] = { "dndeta_correction", "PythiaRatios" };
278
279 Track2Particle1DComposition(fileNames, 2, folderNames);
280}
281
9f469bf5 282void DrawDifferentSpecies()
283{
284 gROOT->ProcessLine(".L drawPlots.C");
285
6b7fa615 286 const char* fileNames[] = { "systematics.root", "systematics.root", "systematics.root", "systematics.root" };
9f469bf5 287 const char* folderNames[] = { "correction_0", "correction_1", "correction_2", "correction_3" };
288
6b7fa615 289 Track2Particle1DComposition(fileNames, 4, folderNames);
9f469bf5 290}
10ebe68d 291
8ca1a6d9 292void DrawpiKpAndCombined()
293{
294 gROOT->ProcessLine(".L drawPlots.C");
295
296 const char* fileNames[] = { "systematics.root", "systematics.root", "systematics.root", "correction_map.root" };
297 const char* folderNames[] = { "correction_0", "correction_1", "correction_2", "dndeta_correction" };
298
299 Track2Particle1DComposition(fileNames, 4, folderNames);
300}
301
6b7fa615 302void DrawSpeciesAndCombination()
303{
304 gROOT->ProcessLine(".L drawPlots.C");
305
306 const char* fileNames[] = { "systematics.root", "systematics.root", "systematics.root", "new_compositions.root" };
307 const char* folderNames[] = { "correction_0", "correction_1", "correction_2", "PythiaRatios" };
308
309 Track2Particle1DComposition(fileNames, 4, folderNames);
310}
311
312void DrawSimulatedVsCombined()
313{
314 gROOT->ProcessLine(".L drawPlots.C");
315
316 const char* fileNames[] = { "new_compositions.root", "new_compositions.root" };
317 const char* folderNames[] = { "Pythia", "PythiaRatios" };
318
319 Track2Particle1DComposition(fileNames, 2, folderNames);
320}
321
322void DrawBoosts()
323{
324 gROOT->ProcessLine(".L drawPlots.C");
325
326 const char* fileNames[] = { "new_compositions.root", "new_compositions.root", "new_compositions.root", "new_compositions.root" };
327 const char* folderNames[] = { "PythiaRatios", "PiBoosted", "KBoosted", "pBoosted" };
328
329 Track2Particle1DComposition(fileNames, 4, folderNames);
330}
331
7af955da 332TH2F* HistogramDifferences(const char* filename, const char* folder1, const char* folder2)
6b7fa615 333{
334 gSystem->Load("libPWG0base");
335
336 AlidNdEtaCorrection* fdNdEtaCorrection[2];
337
338 TFile::Open(filename);
339
340 fdNdEtaCorrection[0] = new AlidNdEtaCorrection(folder1, folder1);
341 fdNdEtaCorrection[0]->LoadHistograms(filename, folder1);
342
343 fdNdEtaCorrection[1] = new AlidNdEtaCorrection(folder2, folder2);
344 fdNdEtaCorrection[1]->LoadHistograms(filename, folder2);
345
6b7fa615 346 TH3F* hist1 = fdNdEtaCorrection[0]->GetTrack2ParticleCorrection()->GetCorrectionHistogram();
347 TH3F* hist2 = fdNdEtaCorrection[1]->GetTrack2ParticleCorrection()->GetCorrectionHistogram();
348
7af955da 349 //TH1F* difference = new TH1F("difference", Form(";#DeltaC_{pT, z, #eta} %s / %s;Count", folder2, folder1), 1000, 0.9, 1.1);
350 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());
351
6b7fa615 352 for (Int_t x=hist1->GetXaxis()->FindBin(-10); x<=hist1->GetXaxis()->FindBin(10); ++x)
353 for (Int_t y=hist1->GetYaxis()->FindBin(-0.8); y<=hist1->GetYaxis()->FindBin(0.8); ++y)
354 for (Int_t z=hist1->GetZaxis()->FindBin(0.3); z<=hist1->GetZaxis()->FindBin(9.9); ++z)
7af955da 355 if (hist1->GetBinContent(x, y, z) != 0)
356 difference->Fill(hist2->GetBinContent(x, y, z) / hist1->GetBinContent(x, y, z), hist1->GetYaxis()->GetBinCenter(y));
357
358 difference->GetYaxis()->SetRangeUser(-0.8, 0.8);
6b7fa615 359
360 printf("Over-/Underflow bins: %d %d\n", difference->GetBinContent(0), difference->GetBinContent(difference->GetNbinsX()+1));
361
362 return difference;
363}
364
365void HistogramDifferences()
366{
7af955da 367 TH2F* KBoosted = HistogramDifferences("new_compositions.root", "PythiaRatios", "KBoosted");
368 TH2F* pBoosted = HistogramDifferences("new_compositions.root", "PythiaRatios", "pBoosted");
369 TH2F* KReduced = HistogramDifferences("new_compositions.root", "PythiaRatios", "KReduced");
370 TH2F* pReduced = HistogramDifferences("new_compositions.root", "PythiaRatios", "pReduced");
6b7fa615 371
7af955da 372 TCanvas* canvas = new TCanvas("HistogramDifferences", "HistogramDifferences", 1000, 1000);
373 canvas->Divide(2, 2);
6b7fa615 374
375 canvas->cd(1);
7af955da 376 KBoosted->GetXaxis()->SetRangeUser(-0.05, 0.05);
377 KBoosted->Draw("COLZ");
6b7fa615 378
379 canvas->cd(2);
7af955da 380 KReduced->GetXaxis()->SetRangeUser(-0.05, 0.05);
381 KReduced->Draw("COLZ");
6b7fa615 382
383 canvas->cd(3);
7af955da 384 pBoosted->GetXaxis()->SetRangeUser(-0.02, 0.02);
385 pBoosted->Draw("COLZ");
386
387 canvas->cd(4);
388 pReduced->GetXaxis()->SetRangeUser(-0.02, 0.02);
389 pReduced->Draw("COLZ");
6b7fa615 390
391 canvas->SaveAs("HistogramDifferences.gif");
7af955da 392
393 canvas = new TCanvas("HistogramDifferencesProfile", "HistogramDifferencesProfile", 1000, 500);
394 canvas->Divide(2, 1);
395
396 canvas->cd(1);
397 gPad->SetBottomMargin(0.13);
398 KBoosted->SetTitle("a) Ratio of correction maps");
399 KBoosted->SetStats(kFALSE);
400 KBoosted->GetXaxis()->SetTitleOffset(1.4);
401 KBoosted->GetXaxis()->SetLabelOffset(0.02);
402 KBoosted->Draw("COLZ");
403
404 canvas->cd(2);
405 gPad->SetGridx();
406 gPad->SetGridy();
407
408 TLegend* legend = new TLegend(0.73, 0.73, 0.98, 0.98);
409
410 for (Int_t i=0; i<4; ++i)
411 {
412 TH2F* hist = 0;
413 TString name;
414 switch (i)
415 {
416 case 0: hist = KBoosted; name = "K enhanced"; break;
417 case 1: hist = KReduced; name = "K reduced"; break;
418 case 2: hist = pBoosted; name = "p enhanced"; break;
419 case 3: hist = pReduced; name = "p reduced"; break;
420 }
421
422 TProfile* profile = hist->ProfileY();
423 profile->SetTitle("b) Mean and RMS");
424 profile->GetXaxis()->SetRange(hist->GetYaxis()->GetFirst(), hist->GetYaxis()->GetLast());
425 profile->GetXaxis()->SetTitleOffset(1.2);
426 profile->GetXaxis()->SetLabelOffset(0.02);
427 profile->GetYaxis()->SetRangeUser(0.98, 1.02);
428 profile->SetStats(kFALSE);
429 profile->SetLineColor(i+1);
430 profile->SetMarkerColor(i+1);
431 profile->DrawCopy(((i > 0) ? "SAME" : ""));
432
433
434 legend->AddEntry(profile, name);
435 }
436
437 legend->Draw();
438 canvas->SaveAs("particlecomposition_result.eps");
6b7fa615 439}
440
441
9f469bf5 442void ScalePtDependent(TH3F* hist, TF1* function)
443{
444 // assumes that pt is the third dimension of hist
445 // scales with function(pt)
10ebe68d 446
9f469bf5 447 for (Int_t z=1; z<=hist->GetNbinsZ(); ++z)
448 {
449 Double_t factor = function->Eval(hist->GetZaxis()->GetBinCenter(z));
450 printf("z = %d, pt = %f, scaling with %f\n", z, hist->GetZaxis()->GetBinCenter(z), factor);
451
452 for (Int_t x=1; x<=hist->GetNbinsX(); ++x)
453 for (Int_t y=1; y<=hist->GetNbinsY(); ++y)
454 hist->SetBinContent(x, y, z, hist->GetBinContent(x, y, z) * factor);
455 }
456}
457
6b7fa615 458void ScalePtDependent(TH3F* hist, TH1* function)
459{
460 // assumes that pt is the third dimension of hist
461 // scales with histogram(pt)
462
463 for (Int_t z=1; z<=hist->GetNbinsZ(); ++z)
464 {
465 Double_t factor = function->GetBinContent(function->GetXaxis()->FindBin(hist->GetZaxis()->GetBinCenter(z)));
466 printf("z = %d, pt = %f, scaling with %f\n", z, hist->GetZaxis()->GetBinCenter(z), factor);
467
468 for (Int_t x=1; x<=hist->GetNbinsX(); ++x)
469 for (Int_t y=1; y<=hist->GetNbinsY(); ++y)
470 hist->SetBinContent(x, y, z, hist->GetBinContent(x, y, z) * factor);
471 }
472}
473
474const char* ChangeComposition(void** correctionPointer, Int_t index)
9f469bf5 475{
476 AlidNdEtaCorrection** fdNdEtaCorrection = (AlidNdEtaCorrection**) correctionPointer;
477
478 switch (index)
479 {
6b7fa615 480 case 0: // result from pp events
481 {
482 TFile::Open("pythiaratios.root");
483
484 for (Int_t i=0; i<4; ++i)
485 {
486 TString name;
487 name.Form("correction_%d", i);
488 fdNdEtaCorrection[i] = new AlidNdEtaCorrection(name, name);
489 fdNdEtaCorrection[i]->LoadHistograms("pythiaratios.root", name);
490 }
491 }
492 return "Pythia";
493 break;
494
495 case 1: // each species rated with pythia ratios
7af955da 496 /*TH1** ptDists = DrawRatios("pythiaratios.root");
6b7fa615 497
498 for (Int_t i=0; i<3; ++i)
499 {
500 ScalePtDependent(fdNdEtaCorrection[i]->GetTrack2ParticleCorrection()->GetMeasuredHistogram(), ptDists[i]);
501 ScalePtDependent(fdNdEtaCorrection[i]->GetTrack2ParticleCorrection()->GetGeneratedHistogram(), ptDists[i]);
7af955da 502 }*/
6b7fa615 503 return "PythiaRatios";
9f469bf5 504 break;
505
7af955da 506 // one species enhanced / reduced
507 case 2: // + 50% pions
508 case 3: // - 50% pions
509 case 4: // + 50% kaons
510 case 5: // - 50% kaons
511 case 6: // + 50% protons
512 case 7: // - 50% protons
513 Int_t correctionIndex = (index - 2) / 2;
514 Double_t scaleFactor = (index % 2 == 0) ? 1.5 : 0.5;
515
516 fdNdEtaCorrection[correctionIndex]->GetTrack2ParticleCorrection()->GetMeasuredHistogram()->Scale(scaleFactor);
517 fdNdEtaCorrection[correctionIndex]->GetTrack2ParticleCorrection()->GetGeneratedHistogram()->Scale(scaleFactor);
518
519 TString* str = new TString;
520 str->Form("%s%s", (correctionIndex == 0) ? "Pi" : ((correctionIndex == 1) ? "K" : "p"), (index % 2 == 0) ? "Boosted" : "Reduced");
521 return str->Data();
522 break;
523
6b7fa615 524 // each species rated with pythia ratios
7af955da 525 case 12: // + 50% pions
526 case 13: // - 50% pions
527 case 14: // + 50% kaons
528 case 15: // - 50% kaons
529 case 16: // + 50% protons
530 case 17: // - 50% protons
6b7fa615 531 TH1** ptDists = DrawRatios("pythiaratios.root");
532 Int_t functionIndex = (index - 2) / 2;
7af955da 533 Double_t scaleFactor = (index % 2 == 0) ? 1.5 : 0.5;
6b7fa615 534 ptDists[functionIndex]->Scale(scaleFactor);
535
536 for (Int_t i=0; i<3; ++i)
537 {
538 ScalePtDependent(fdNdEtaCorrection[i]->GetTrack2ParticleCorrection()->GetMeasuredHistogram(), ptDists[i]);
539 ScalePtDependent(fdNdEtaCorrection[i]->GetTrack2ParticleCorrection()->GetGeneratedHistogram(), ptDists[i]);
540 }
541 TString* str = new TString;
542 str->Form("%s%s", (functionIndex == 0) ? "Pi" : ((functionIndex == 1) ? "K" : "p"), (index % 2 == 0) ? "Boosted" : "Reduced");
543 return str->Data();
9f469bf5 544 break;
545
6b7fa615 546 case 999:
9f469bf5 547 TF1* ptDependence = new TF1("simple", "x", 0, 100);
548 ScalePtDependent(fdNdEtaCorrection[0]->GetTrack2ParticleCorrection()->GetGeneratedHistogram(), ptDependence);
549 ScalePtDependent(fdNdEtaCorrection[0]->GetTrack2ParticleCorrection()->GetMeasuredHistogram(), ptDependence);
550 break;
551
552 }
6b7fa615 553
554 return "noname";
9f469bf5 555}
556
557void Composition()
558{
559 gSystem->Load("libPWG0base");
560
561 gSystem->Unlink("new_compositions.root");
562
6b7fa615 563 Int_t nCompositions = 8;
564 for (Int_t comp = 0; comp < nCompositions; ++comp)
9f469bf5 565 {
566 AlidNdEtaCorrection* fdNdEtaCorrection[4];
567
568 TFile::Open("systematics.root");
569
570 for (Int_t i=0; i<4; ++i)
571 {
572 TString name;
573 name.Form("correction_%d", i);
574 fdNdEtaCorrection[i] = new AlidNdEtaCorrection(name, name);
575 fdNdEtaCorrection[i]->LoadHistograms("systematics.root", name);
576 }
577
6b7fa615 578 const char* newName = ChangeComposition(fdNdEtaCorrection, comp);
9f469bf5 579
580 Double_t geneCount[5];
581 Double_t measCount[5];
582 geneCount[4] = 0;
583 measCount[4] = 0;
584
585 for (Int_t i=0; i<4; ++i)
586 {
587 geneCount[i] = fdNdEtaCorrection[i]->GetTrack2ParticleCorrection()->GetGeneratedHistogram()->Integral();
588 measCount[i] = fdNdEtaCorrection[i]->GetTrack2ParticleCorrection()->GetMeasuredHistogram()->Integral();
589
590 geneCount[4] += geneCount[i];
591 measCount[4] += measCount[i];
592
593 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]);
594 }
10ebe68d 595
9f469bf5 596 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 597
9f469bf5 598 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 599
9f469bf5 600 TList* collection = new TList;
10ebe68d 601
9f469bf5 602 for (Int_t i=0; i<4; ++i)
603 collection->Add(fdNdEtaCorrection[i]);
10ebe68d 604
6b7fa615 605 AlidNdEtaCorrection* newComposition = new AlidNdEtaCorrection(newName, newName);
9f469bf5 606 newComposition->Merge(collection);
607 newComposition->Finish();
608
609 delete collection;
610
611 TFile* file = TFile::Open("new_compositions.root", "UPDATE");
612 newComposition->SaveHistograms();
613 //file->Write();
614 file->Close();
615 }
10ebe68d 616
617 gROOT->ProcessLine(".L drawPlots.C");
9f469bf5 618
6b7fa615 619 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" };
620 const char* folderNames[] = { "Pythia", "PythiaRatios", "PiBoosted", "PiReduced", "KBoosted", "KReduced", "pBoosted", "pReduced" };
9f469bf5 621
6b7fa615 622 Track2Particle1DComposition(fileNames, nCompositions, folderNames);
9f469bf5 623}
624
9f469bf5 625
10ebe68d 626void drawSystematics()
627{
9f469bf5 628 //Secondaries();
629 //DrawDifferentSpecies();
630 //Composition();
631
632 Sigma2VertexSimulation();
27c04e01 633
10ebe68d 634}
7af955da 635
636void DrawdNdEtaDifferences()
637{
638 TH1* hists[5];
639
72e597d7 640 TLegend* legend = new TLegend(0.3, 0.73, 0.70, 0.98);
f1076daa 641 legend->SetFillColor(0);
642
7af955da 643 TCanvas* canvas = new TCanvas("DrawdNdEtaDifferences", "DrawdNdEtaDifferences", 1000, 500);
644 canvas->Divide(2, 1);
645
646 canvas->cd(1);
647
648 for (Int_t i=0; i<5; ++i)
649 {
72e597d7 650 hists[i] = 0;
7af955da 651 TFile* file = 0;
f1076daa 652 TString title;
7af955da 653
654 switch(i)
655 {
f1076daa 656 case 0 : file = TFile::Open("systematics_dndeta_reference.root"); title = "standard composition"; break;
657 case 1 : file = TFile::Open("systematics_dndeta_KBoosted.root"); title = "+ 50% kaons"; break;
658 case 2 : file = TFile::Open("systematics_dndeta_KReduced.root"); title = "- 50% kaons"; break;
659 case 3 : file = TFile::Open("systematics_dndeta_pBoosted.root"); title = "+ 50% protons"; break;
660 case 4 : file = TFile::Open("systematics_dndeta_pReduced.root"); title = "- 50% protons"; break;
7af955da 661 default: return;
662 }
663
72e597d7 664 if (file)
665 {
666 hists[i] = (TH1*) file->Get("dndeta/dndeta_dNdEta_corrected_2");
667 hists[i]->SetTitle("a)");
668
669 Prepare1DPlot(hists[i], kFALSE);
670 hists[i]->GetXaxis()->SetRangeUser(-0.7999, 0.7999);
da7acf97 671 hists[i]->GetYaxis()->SetRangeUser(6, 7);
8ca1a6d9 672 hists[i]->SetLineColor(colors[i]);
673 hists[i]->SetMarkerColor(colors[i]);
674 hists[i]->SetMarkerStyle(markers[i]);
72e597d7 675 hists[i]->GetXaxis()->SetLabelOffset(0.015);
676 hists[i]->GetYaxis()->SetTitleOffset(1.5);
677 gPad->SetLeftMargin(0.12);
678 hists[i]->DrawCopy(((i > 0) ? "SAME" : ""));
679
680 legend->AddEntry(hists[i], title);
681 hists[i]->SetTitle(title);
682 }
7af955da 683 }
f1076daa 684 legend->Draw();
7af955da 685
686 canvas->cd(2);
f1076daa 687 gPad->SetLeftMargin(0.14);
688
689 TLegend* legend2 = new TLegend(0.73, 0.73, 0.98, 0.98);
690 legend2->SetFillColor(0);
7af955da 691
692 for (Int_t i=1; i<5; ++i)
693 {
72e597d7 694 if (hists[i])
695 {
696 legend2->AddEntry(hists[i]);
697
698 hists[i]->Divide(hists[0]);
699 hists[i]->SetTitle("b)");
8ca1a6d9 700 hists[i]->SetLineColor(colors[i-1]);
701 hists[i]->SetMarkerColor(colors[i-1]);
da7acf97 702 hists[i]->GetYaxis()->SetRangeUser(0.95, 1.05);
72e597d7 703 hists[i]->GetYaxis()->SetTitle("Ratio to standard composition");
704 hists[i]->GetYaxis()->SetTitleOffset(1.8);
705 hists[i]->DrawCopy(((i > 1) ? "SAME" : ""));
706 }
7af955da 707 }
f1076daa 708
709 legend2->Draw();
710
8ca1a6d9 711 canvas->SaveAs("particlecomposition_result_detail.gif");
712
713 TCanvas* canvas2 = new TCanvas("DrawdNdEtaDifferences2", "DrawdNdEtaDifferences2", 700, 500);
714
715 for (Int_t i=1; i<5; ++i)
716 {
717 if (hists[i])
718 {
719 hists[i]->SetTitle("");
720 hists[i]->GetYaxis()->SetTitleOffset(1.1);
721 hists[i]->DrawCopy(((i > 1) ? "SAME" : ""));
722 }
723 }
724
725 legend2->Draw();
726
727 canvas2->SaveAs("particlecomposition_result.gif");
728 canvas2->SaveAs("particlecomposition_result.eps");
7af955da 729}
27c04e01 730
51f6de65 731void mergeCorrectionsWithDifferentCrosssections(Int_t correctionTarget = 3, Char_t* correctionFileName="correction_mapprocess-types.root", const char* analysisFileName = "analysis_esd_raw.root", const Char_t* outputFileName="systematics_vtxtrigger_compositions.root") {
27c04e01 732 //
733 // Function used to merge standard corrections with vertex
734 // reconstruction corrections obtained by a certain mix of ND, DD
735 // and SD events.
dd367a14 736 //
737 // the dn/deta spectrum is corrected and the ratios
738 // (standard to changed x-section) of the different dN/deta
739 // distributions are saved to a file.
740 //
51f6de65 741 // correctionTarget is of type AlidNdEtaCorrection::CorrectionType
742 // kINEL = 3
743 // kNSD = 4
27c04e01 744
dd367a14 745 loadlibs();
27c04e01 746
da7acf97 747 const Char_t* typeName[] = { "vertexreco", "trigger", "vtxtrigger" };
eaa3702a 748
51f6de65 749
750
dd367a14 751 const Char_t* changes[] = { "pythia","ddmore","ddless","sdmore","sdless", "dmore", "dless", "sdmoreddless", "sdlessddmore", "ddmore25","ddless25","sdmore25","sdless25", "dmore25", "dless25", "sdmoreddless25", "sdlessddmore25"};
51f6de65 752 Float_t scalesND[] = {1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0 };
dd367a14 753 Float_t scalesDD[] = {1.0, 1.5, 0.5, 1.0, 1.0, 1.5, 0.5, 1.5, 0.5, 1.25, 0.75, 1.0, 1.0, 1.25, 0.75, 1.25, 0.75};
754 Float_t scalesSD[] = {1.0, 1.0, 1.0, 1.5, 0.5, 1.5, 0.5, 0.5, 1.5, 1.0, 1.0, 1.25, 0.75, 1.25, 0.75, 0.75, 1.25};
eaa3702a 755 Int_t nChanges = 17;
eaa3702a 756
51f6de65 757 /*
eaa3702a 758 const Char_t* changes[] = { "pythia", "qgsm", "phojet"};
759 Float_t scalesND[] = {1.0, 1.10, 1.11};
760 Float_t scalesSD[] = {1.0, 0.69, 0.86};
761 Float_t scalesDD[] = {1.0, 0.98, 0.61};
762 Int_t nChanges = 3;
51f6de65 763 */
764
27c04e01 765 // cross section from Pythia
51f6de65 766 // 14 TeV!
27c04e01 767 Float_t sigmaND = 55.2;
768 Float_t sigmaDD = 9.78;
769 Float_t sigmaSD = 14.30;
770
dd367a14 771 // standard correction
772 TFile::Open(correctionFileName);
773 AlidNdEtaCorrection* correctionStandard = new AlidNdEtaCorrection("dndeta_correction","dndeta_correction");
774 correctionStandard->LoadHistograms();
775
776 // dont take vertexreco from this one
777 correctionStandard->GetVertexRecoCorrection()->Reset();
778 // dont take triggerbias from this one
779 correctionStandard->GetTriggerBiasCorrectionINEL()->Reset();
780 correctionStandard->GetTriggerBiasCorrectionNSD()->Reset();
781 correctionStandard->GetTriggerBiasCorrectionND()->Reset();
782
783 AlidNdEtaCorrection* corrections[100];
784 TH1F* hRatios[100];
785
786 Int_t counter = 0;
da7acf97 787 for (Int_t j=0; j<3; j++) { // j = 0 (change vtx), j = 1 (change trg), j = 2 (change both)
27c04e01 788
eaa3702a 789 for (Int_t i=0; i<nChanges; i++) {
dd367a14 790 TFile::Open(correctionFileName);
27c04e01 791
da7acf97 792 TString name;
793 name.Form("dndeta_correction_syst_%s_%s", typeName[j], changes[i]);
794 AlidNdEtaCorrection* current = new AlidNdEtaCorrection(name, name);
770a1f1d 795 current->LoadHistograms("dndeta_correction");
796 current->Reset();
da7acf97 797
dd367a14 798 name.Form("dndeta_correction_ND");
799 AlidNdEtaCorrection* dNdEtaCorrectionND = new AlidNdEtaCorrection(name,name);
800 dNdEtaCorrectionND->LoadHistograms();
801 name.Form("dndeta_correction_DD");
802 AlidNdEtaCorrection* dNdEtaCorrectionDD = new AlidNdEtaCorrection(name,name);
803 dNdEtaCorrectionDD->LoadHistograms();
804 name.Form("dndeta_correction_SD");
805 AlidNdEtaCorrection* dNdEtaCorrectionSD = new AlidNdEtaCorrection(name,name);
806 dNdEtaCorrectionSD->LoadHistograms();
da7acf97 807
808 // calculating relative
809 Float_t nd = 100 * sigmaND/(sigmaND + (scalesDD[i]*sigmaDD) + (scalesDD[i]*sigmaSD));
810 Float_t dd = 100 * (scalesDD[i]*sigmaDD)/(sigmaND + (scalesDD[i]*sigmaDD) + (scalesDD[i]*sigmaSD));
811 Float_t sd = 100 * (scalesSD[i]*sigmaSD)/(sigmaND + (scalesDD[i]*sigmaDD) + (scalesDD[i]*sigmaSD));
812
dd367a14 813 printf(Form("%s : ND=%.2f\%, DD=%.2f\%, SD=%.2f\% \n",changes[i],nd,dd,sd));
da7acf97 814 current->SetTitle(Form("ND=%.2f\%,DD=%.2f\%,SD=%.2f\%",nd,dd,sd));
815 current->SetTitle(name);
816
817 // scale
818 if (j == 0 || j == 2)
819 {
eaa3702a 820 dNdEtaCorrectionND->GetVertexRecoCorrection()->Scale(scalesND[i]);
dd367a14 821 dNdEtaCorrectionDD->GetVertexRecoCorrection()->Scale(scalesDD[i]);
822 dNdEtaCorrectionSD->GetVertexRecoCorrection()->Scale(scalesSD[i]);
da7acf97 823 }
824 if (j == 1 || j == 2)
825 {
eaa3702a 826 dNdEtaCorrectionND->GetTriggerBiasCorrectionINEL()->Scale(scalesND[i]);
dd367a14 827 dNdEtaCorrectionDD->GetTriggerBiasCorrectionINEL()->Scale(scalesDD[i]);
828 dNdEtaCorrectionSD->GetTriggerBiasCorrectionINEL()->Scale(scalesSD[i]);
829
eaa3702a 830 dNdEtaCorrectionND->GetTriggerBiasCorrectionNSD()->Scale(scalesND[i]);
dd367a14 831 dNdEtaCorrectionDD->GetTriggerBiasCorrectionNSD()->Scale(scalesDD[i]);
832 dNdEtaCorrectionSD->GetTriggerBiasCorrectionNSD()->Scale(scalesSD[i]);
833
eaa3702a 834 dNdEtaCorrectionND->GetTriggerBiasCorrectionND()->Scale(scalesND[i]);
dd367a14 835 dNdEtaCorrectionDD->GetTriggerBiasCorrectionND()->Scale(scalesDD[i]);
836 dNdEtaCorrectionSD->GetTriggerBiasCorrectionND()->Scale(scalesSD[i]);
da7acf97 837 }
838
dd367a14 839 //clear track in correction
840 dNdEtaCorrectionND->GetTrack2ParticleCorrection()->Reset();
841 dNdEtaCorrectionDD->GetTrack2ParticleCorrection()->Reset();
842 dNdEtaCorrectionSD->GetTrack2ParticleCorrection()->Reset();
da7acf97 843
844 TList collection;
845 collection.Add(correctionStandard);
dd367a14 846 collection.Add(dNdEtaCorrectionND);
847 collection.Add(dNdEtaCorrectionDD);
848 collection.Add(dNdEtaCorrectionSD);
da7acf97 849
850 current->Merge(&collection);
851 current->Finish();
852
dd367a14 853 corrections[counter] = current;
854
855 // now correct dNdeta distribution with modified correction map
856 TFile::Open(analysisFileName);
857
858 dNdEtaAnalysis* fdNdEtaAnalysis = new dNdEtaAnalysis("fdNdEtaAnalysisESD", "fdNdEtaAnalysisESD");
859 fdNdEtaAnalysis->LoadHistograms();
860
51f6de65 861 fdNdEtaAnalysis->Finish(current, 0.3, correctionTarget, Form("%d %d", j, i));
dd367a14 862
863 name = "ratio";
864 if (j==0) name.Append("_vetexReco_");
865 if (j==1) name.Append("_triggerBias_");
866 if (j==2) name.Append("_vertexReco_triggerBias_");
867 name.Append(changes[i]);
868
869 hRatios[counter] = (TH1F*) fdNdEtaAnalysis->GetdNdEtaHistogram()->Clone(name);
870
871 name.Append(Form(" (DD #times %0.2f, SD #times %0.2f)",scalesDD[i],scalesSD[i]));
872 hRatios[counter]->SetTitle(name.Data());
873 hRatios[counter]->SetYTitle("ratio (Pythia x-sections)/(changed x-sections)");
874
875 if (counter > 0)
876 hRatios[counter]->Divide(hRatios[0],hRatios[counter],1,1);
877
878 delete fdNdEtaAnalysis;
879
880 counter++;
da7acf97 881 }
27c04e01 882 }
da7acf97 883
27c04e01 884 TFile* fout = new TFile(outputFileName,"RECREATE");
885
dd367a14 886 // to make everything consistent
887 hRatios[0]->Divide(hRatios[0],hRatios[0],1,1);
888
eaa3702a 889 for (Int_t i=0; i<nChanges * 3; i++)
dd367a14 890 {
27c04e01 891 corrections[i]->SaveHistograms();
dd367a14 892 hRatios[i]->Write();
893 }
27c04e01 894
895 fout->Write();
896 fout->Close();
897}
5e08578b 898
899
7f0efb86 900DrawVertexRecoSyst() {
901 // Draws the ratio of the dN/dEta obtained with changed SD and DD
902 // cross-sections vertex reco corrections to the dN/dEta obtained
903 // from the standard pythia cross-sections
904 //
905 // The files with the vertex reco corrections for different
906 // processes (and with the other standard corrections) are expected
907 // to have the names "analysis_esd_X.root", where the Xs are defined
908 // in the array changes below.
5e08578b 909
910 Char_t* changes[] = {"pythia","ddmore","ddless","sdmore","sdless", "dmore", "dless"};
911 Char_t* descr[] = {"",
912 "#sigma_{DD}' = 1.5#times#sigma_{DD}",
913 "#sigma_{DD}' = 0.5#times#sigma_{DD}",
914 "#sigma_{SD}' = 1.5#times#sigma_{SD}",
915 "#sigma_{SD}' = 0.5#times#sigma_{SD}",
916 "#sigma_{D}' = 1.5#times#sigma_{D}",
917 "#sigma_{D}' = 0.5#times#sigma_{D}"};
918
919 Float_t scalesDD[] = {1.0, 1.5, 0.5, 1.5, 0.5, 1.5, 0.5};
920 Float_t scalesSD[] = {1.0, 1.0, 1.0, 1.0, 1.0, 1.5, 0.5};
921
fb30f6e3 922 Int_t markers[] = {20,20,21,22,23,28,29};
5e08578b 923 Int_t colors[] = {1,2,3,4,6,8,102};
924
925 // cross section from Pythia
926 Float_t sigmaND = 55.2;
927 Float_t sigmaDD = 9.78;
928 Float_t sigmaSD = 14.30;
929
5e08578b 930 TH1F* dNdEta[7];
5e08578b 931 TH1F* ratios[7];
932
933 TFile* fin;
934
5e08578b 935 for (Int_t i=0; i<7; i++) {
fb30f6e3 936 // calculating relative
7f0efb86 937 fin = TFile::Open(Form("analysis_esd_%s.root",changes[i]));
fb30f6e3 938
5e08578b 939 dNdEta[i] = (TH1F*)(fin->Get("dndeta/dndeta_dNdEta_corrected_2"))->Clone();
fb30f6e3 940
5e08578b 941 for (Int_t b=0; b<dNdEta[i]->GetNbinsX(); b++) {
942 if (TMath::Abs(dNdEta[i]->GetBinCenter(b))>0.9) {
943 dNdEta[i]->SetBinContent(b,0);
944 dNdEta[i]->SetBinError(b,0);
945 }
946 }
fb30f6e3 947
5e08578b 948 dNdEta[i]->Rebin(4);
fb30f6e3 949
5e08578b 950 dNdEta[i]->SetLineWidth(2);
951 dNdEta[i]->SetLineColor(colors[i]);
952 dNdEta[i]->SetMarkerStyle(markers[i]);
fb30f6e3 953 dNdEta[i]->SetMarkerSize(0.9);
954 dNdEta[i]->SetMarkerColor(colors[i]);
955
5e08578b 956 ratios[i] = (TH1F*)dNdEta[i]->Clone("ratio");
957 ratios[i]->Divide(ratios[i],dNdEta[0],1,1,"B");
958
959 ratios[i]->SetName(changes[i]);
960 ratios[i]->SetTitle(changes[i]);
961 }
962
963 //##########################################################
964
965 gStyle->SetOptStat(0);
966 gStyle->SetOptTitle(0);
967 gStyle->SetOptFit(0);
968
969 gStyle->SetTextSize(0.2);
970 gStyle->SetTitleSize(0.05,"xyz");
5e08578b 971 gStyle->SetLabelOffset(0.01, "xyz");
972
973
974 gStyle->SetTitleOffset(1.2, "y");
975 gStyle->SetTitleOffset(1.2, "x");
976 gStyle->SetEndErrorSize(0.0);
977
978 //##############################################
979
980 //making canvas and pads
da7acf97 981 TCanvas *c = new TCanvas(Form("vertex_reco_syst_%s", plot), Form("vertex_reco_syst_%s", plot),600,500);
5e08578b 982
983 TPad* p1 = new TPad("pad1","", 0, 0.0, 1.0, 1.0, 0, 0, 0);
984
985 p1->SetBottomMargin(0.15);
986 p1->SetTopMargin(0.03);
987 p1->SetLeftMargin(0.15);
988 p1->SetRightMargin(0.03);
da7acf97 989
5e08578b 990 p1->SetGridx();
991 p1->SetGridy();
992
993 p1->Draw();
994 p1->cd();
5e08578b 995
996
997 TH2F* null = new TH2F("","",100,-1.5,1.5,100,0.9601,1.099);
998 null->SetXTitle("#eta");
999 null->GetXaxis()->CenterTitle(kTRUE);
1000 null->SetYTitle("dN/d#eta / dN/d#eta_{pythia}");
1001 null->GetYaxis()->CenterTitle(kTRUE);
1002 null->Draw();
8ca1a6d9 1003
5e08578b 1004 for (Int_t i=1; i<7; i++)
1005 ratios[i]->Draw("same");
fb30f6e3 1006
1007 TLegend* legend = new TLegend(0.6, 0.6, 0.95, 0.95);
1008 legend->SetFillColor(0);
1009
5e08578b 1010 TLatex* text[7];
1011 for (Int_t i=1; i<7; i++) {
fb30f6e3 1012 legend->AddEntry(dNdEta[i], descr[i]);
5e08578b 1013 }
fb30f6e3 1014
1015 legend->Draw();
da7acf97 1016 //text(0.2,0.88,"Effect of changing",0.045,1,kTRUE);
1017 //text(0.2,0.83,"relative cross-sections",0.045,1,kTRUE);
1018 //text(0.2,0.78,"(vertex reconstruction corr.)",0.043,13,kTRUE);
1019
1020 c->SaveAs(Form("%s.gif", c->GetName()));
fb30f6e3 1021 c->SaveAs(Form("%s.eps", c->GetName()));
5e08578b 1022}
1023
1024
1025
1026DrawTriggerEfficiency(Char_t* fileName) {
1027
1028 gStyle->SetOptStat(0);
1029 gStyle->SetOptTitle(0);
1030 gStyle->SetOptFit(0);
1031
1032 gStyle->SetTextSize(0.04);
1033 gStyle->SetTitleSize(0.05,"xyz");
1034 //gStyle->SetTitleFont(133, "xyz");
1035 //gStyle->SetLabelFont(133, "xyz");
1036 //gStyle->SetLabelSize(17, "xyz");
1037 gStyle->SetLabelOffset(0.01, "xyz");
1038
1039 gStyle->SetTitleOffset(1.1, "y");
1040 gStyle->SetTitleOffset(1.1, "x");
1041 gStyle->SetEndErrorSize(0.0);
1042
1043 //##############################################
1044
1045 //making canvas and pads
1046 TCanvas *c = new TCanvas("trigger_eff", "",600,500);
1047
1048 TPad* p1 = new TPad("pad1","", 0, 0.0, 1.0, 1.0, 0, 0, 0);
1049
1050 p1->SetBottomMargin(0.15);
1051 p1->SetTopMargin(0.03);
1052 p1->SetLeftMargin(0.15);
1053 p1->SetRightMargin(0.03);
1054
1055 p1->SetGridx();
1056 p1->SetGridy();
1057
1058 p1->Draw();
1059 p1->cd();
1060
1061 gSystem->Load("libPWG0base");
1062
1063 AlidNdEtaCorrection* corrections[4];
1064 AliCorrectionMatrix2D* triggerBiasCorrection[4];
1065
1066 TH1F* hTriggerEffInv[4];
1067 TH1F* hTriggerEff[4];
1068
1069 Char_t* names[] = {"triggerBiasND", "triggerBiasDD", "triggerBiasSD", "triggerBiasINEL"};
1070
1071 for (Int_t i=0; i<4; i++) {
1072 corrections[i] = new AlidNdEtaCorrection(names[i], names[i]);
1073 corrections[i]->LoadHistograms(fileName, names[i]);
1074
1075 triggerBiasCorrection[i] = corrections[i]->GetTriggerBiasCorrectionINEL();
1076
1077
1078 hTriggerEffInv[i] = triggerBiasCorrection[i]->Get1DCorrection();
1079 hTriggerEff[i] = (TH1F*)hTriggerEffInv[i]->Clone();
1080
1081 for (Int_t b=0; b<=hTriggerEff[i]->GetNbinsX(); b++) {
1082 hTriggerEff[i]->SetBinContent(b,1);
1083 hTriggerEff[i]->SetBinError(b,0);
1084 }
1085 hTriggerEff[i]->Divide(hTriggerEff[i],hTriggerEffInv[i]);
1086 hTriggerEff[i]->Scale(100);
1087 }
1088
1089 Int_t colors[] = {2,3,4,1};
1090 Float_t effs[4];
1091 for (Int_t i=0; i<4; i++) {
1092 hTriggerEff[i]->Fit("pol0","","",-20,20);
1093 effs[i] = ((TF1*)hTriggerEff[i]->GetListOfFunctions()->At(0))->GetParameter(0);
1094 ((TF1*)hTriggerEff[i]->GetListOfFunctions()->At(0))->SetLineWidth(1);
1095 ((TF1*)hTriggerEff[i]->GetListOfFunctions()->At(0))->SetLineStyle(2);
1096 ((TF1*)hTriggerEff[i]->GetListOfFunctions()->At(0))->SetLineColor(colors[i]);
1097 cout << effs[i] << endl;
1098 }
1099
1100
1101 Char_t* text[] = {"ND", "DD", "SD", "INEL"};
1102 TLatex* latex[4];
1103
1104 TH2F* null = new TH2F("","",100,-25,35,100,0,110);
1105 null->SetXTitle("Vertex z [cm]");
1106 null->GetXaxis()->CenterTitle(kTRUE);
1107 null->SetYTitle("Trigger efficiency [%]");
1108 null->GetYaxis()->CenterTitle(kTRUE);
1109 null->Draw();
1110
1111
1112 for (Int_t i=0; i<4; i++) {
1113 hTriggerEff[i]->SetLineWidth(2);
1114 hTriggerEff[i]->SetLineColor(colors[i]);
1115
1116 hTriggerEff[i]->Draw("same");
1117
1118 latex[i] = new TLatex(22,effs[i]-1.5, Form("%s (%0.1f)",text[i],effs[i]));
1119 latex[i]->SetTextColor(colors[i]);
1120 latex[i]->Draw();
1121 }
1122
1123}
1124
1125
1126DrawSpectraPID(Char_t* fileName) {
1127
1128 gSystem->Load("libPWG0base");
1129
1130 Char_t* names[] = {"correction_0", "correction_1", "correction_2", "correction_3"};
1131 AlidNdEtaCorrection* corrections[4];
1132 AliCorrectionMatrix3D* trackToPartCorrection[4];
1133
1134 TH1F* measuredPt[4];
1135 TH1F* generatedPt[4];
1136 TH1F* ratioPt[4];
1137
1138 for (Int_t i=0; i<4; i++) {
1139 corrections[i] = new AlidNdEtaCorrection(names[i], names[i]);
1140 corrections[i]->LoadHistograms(fileName, names[i]);
1141
1142 trackToPartCorrection[i] = corrections[i]->GetTrack2ParticleCorrection();
1143
1144 Int_t binX1 = (TH1F*)trackToPartCorrection[i]->GetMeasuredHistogram()->GetXaxis()->FindBin(-10);
1145 Int_t binX2 = (TH1F*)trackToPartCorrection[i]->GetMeasuredHistogram()->GetXaxis()->FindBin(10);
1146 Int_t binY1 = (TH1F*)trackToPartCorrection[i]->GetMeasuredHistogram()->GetYaxis()->FindBin(-1);
1147 Int_t binY2 = (TH1F*)trackToPartCorrection[i]->GetMeasuredHistogram()->GetYaxis()->FindBin(1);
1148
1149 measuredPt[i] = (TH1F*)trackToPartCorrection[i]->GetMeasuredHistogram()->ProjectionZ(Form("m_%d",i),binX1,binX2,binY1,binY2);
1150 generatedPt[i] = (TH1F*)trackToPartCorrection[i]->GetGeneratedHistogram()->ProjectionZ(Form("g_%d",i),binX1,binX2,binY1,binY2);
1151 ratioPt[i] = (TH1F*)generatedPt[i]->Clone(Form("r_%d",i));
1152 ratioPt[i]->Divide(measuredPt[i], generatedPt[i], 1,1,"B");
1153 }
1154
1155 ratioPt[0]->Draw();
1156
1157 for (Int_t i=0; i<3; i++) {
1158 ratioPt[i]->SetLineColor(i+1);
1159 ratioPt[i]->SetLineWidth(2);
1160
1161 ratioPt[i]->Draw("same");
1162
1163 }
1164
1165 return;
1166 measuredPt[0]->SetLineColor(2);
1167 measuredPt[0]->SetLineWidth(5);
1168
1169 measuredPt[0]->Draw();
1170 generatedPt[0]->Draw("same");
1171}
da7acf97 1172
8277513e 1173void changePtSpectrum(const char* fileName = "analysis_mc.root")
da7acf97 1174{
8277513e 1175 TFile* file = TFile::Open(fileName);
da7acf97 1176 TH1F* hist = dynamic_cast<TH1F*> (file->Get("dndeta_check_pt"));
1177
8ca1a6d9 1178 //hist->Rebin(3);
1179 //hist->Scale(1.0/3);
da7acf97 1180
1181 TH1F* clone1 = dynamic_cast<TH1F*> (hist->Clone("clone1"));
1182 TH1F* clone2 = dynamic_cast<TH1F*> (hist->Clone("clone2"));
1183
1184 TH1F* scale1 = dynamic_cast<TH1F*> (hist->Clone("scale1"));
1185 TH1F* scale2 = dynamic_cast<TH1F*> (hist->Clone("scale2"));
1186
1187 Float_t ptCutOff = 0.3;
1188
1189 for (Int_t i=1; i <= hist->GetNbinsX(); ++i)
1190 {
1191 if (hist->GetBinCenter(i) > ptCutOff)
1192 {
1193 scale1->SetBinContent(i, 1);
1194 scale2->SetBinContent(i, 1);
1195 }
1196 else
1197 {
1198 // 90 % at pt = 0, 0% at pt = ptcutoff
1199 scale1->SetBinContent(i, 1 - (ptCutOff - hist->GetBinCenter(i)) / ptCutOff * 0.3);
1200
1201 // 110% at pt = 0, ...
1202 scale2->SetBinContent(i, 1 + (ptCutOff - hist->GetBinCenter(i)) / ptCutOff * 0.3);
1203 }
8ca1a6d9 1204 scale1->SetBinError(i, 0);
1205 scale2->SetBinError(i, 0);
da7acf97 1206 }
1207
1208 new TCanvas;
8ca1a6d9 1209
da7acf97 1210 scale1->Draw();
1211 scale2->SetLineColor(kRed);
1212 scale2->Draw("SAME");
1213
1214 clone1->Multiply(scale1);
1215 clone2->Multiply(scale2);
1216
1217 Prepare1DPlot(hist);
1218 Prepare1DPlot(clone1);
1219 Prepare1DPlot(clone2);
1220
8ca1a6d9 1221 /*hist->SetMarkerStyle(markers[0]);
1222 clone1->SetMarkerStyle(markers[0]);
1223 clone2->SetMarkerStyle(markers[0]);*/
1224
8277513e 1225 hist->SetTitle(";p_{T} in GeV/c;dN_{ch}/dp_{T} in c/GeV");
da7acf97 1226 hist->GetXaxis()->SetRangeUser(0, 0.7);
8ca1a6d9 1227 hist->GetYaxis()->SetRangeUser(0.01, clone2->GetMaximum() * 1.1);
da7acf97 1228 hist->GetYaxis()->SetTitleOffset(1);
1229
1230 TCanvas* canvas = new TCanvas;
8ca1a6d9 1231 gPad->SetBottomMargin(0.12);
1232 hist->Draw("H");
da7acf97 1233 clone1->SetLineColor(kRed);
8ca1a6d9 1234 clone1->Draw("HSAME");
da7acf97 1235 clone2->SetLineColor(kBlue);
8ca1a6d9 1236 clone2->Draw("HSAME");
1237 hist->Draw("HSAME");
da7acf97 1238
1239 Float_t fraction = hist->Integral(hist->GetXaxis()->FindBin(ptCutOff), hist->GetNbinsX()) / hist->Integral(1, hist->GetNbinsX());
1240 Float_t fraction1 = clone1->Integral(clone1->GetXaxis()->FindBin(ptCutOff), clone1->GetNbinsX()) / clone1->Integral(1, clone1->GetNbinsX());
1241 Float_t fraction2 = clone2->Integral(clone2->GetXaxis()->FindBin(ptCutOff), clone2->GetNbinsX()) / clone2->Integral(1, clone2->GetNbinsX());
1242
1243 printf("%f %f %f\n", fraction, fraction1, fraction2);
1244 printf("Rel. %f %f\n", fraction1 / fraction, fraction2 / fraction);
1245
1246 canvas->SaveAs("changePtSpectrum.gif");
8ca1a6d9 1247 canvas->SaveAs("changePtSpectrum.eps");
da7acf97 1248}
fb30f6e3 1249
1250void FractionBelowPt()
1251{
1252 gSystem->Load("libPWG0base");
1253
1254 AlidNdEtaCorrection* fdNdEtaCorrection[4];
1255
1256 TFile::Open("systematics.root");
1257
1258 for (Int_t i=0; i<4; ++i)
1259 {
1260 TString name;
1261 name.Form("correction_%d", i);
1262 fdNdEtaCorrection[i] = new AlidNdEtaCorrection(name, name);
1263 fdNdEtaCorrection[i]->LoadHistograms("systematics.root", name);
1264 }
1265
1266 Double_t geneCount[5];
1267 Double_t measCount[5];
1268 geneCount[4] = 0;
1269 measCount[4] = 0;
1270
1271 for (Int_t i=0; i<4; ++i)
1272 {
1273 TH3F* hist = (TH3F*) fdNdEtaCorrection[i]->GetTrack2ParticleCorrection()->GetGeneratedHistogram();
1274 geneCount[i] = hist->Integral(hist->GetXaxis()->FindBin(-10), hist->GetXaxis()->FindBin(10),
1275 hist->GetYaxis()->FindBin(-0.8), hist->GetYaxis()->FindBin(0.8),
1276 1, hist->GetZaxis()->FindBin(0.3));
1277
1278 hist = (TH3F*) fdNdEtaCorrection[i]->GetTrack2ParticleCorrection()->GetMeasuredHistogram();
1279 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));
1280
1281 geneCount[4] += geneCount[i];
1282 measCount[4] += measCount[i];
1283
1284 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]);
1285 }
1286
1287 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]);
1288
1289 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 1290}
8ca1a6d9 1291
1292
1293mergeCorrectionsMisalignment(Char_t* alignedFile = "correction_map_aligned.root",
1294 Char_t* misalignedFile = "correction_map_misaligned.root",
1295 Char_t* outputFileName="correction_map_misaligned_single.root")
1296{
1297 //
1298 // from the aligned and misaligned corrections, 3 new corrections are created
1299 // in these new corrections only one of the corrections (track2particle, vertex, trigger)
1300 // is taken from the misaligned input to allow study of the effect on the different
1301 // corrections
1302
1303 gSystem->Load("libPWG0base");
1304
1305 const Char_t* typeName[] = { "track2particle", "vertex", "trigger" };
1306
1307 AlidNdEtaCorrection* corrections[3];
1308 for (Int_t j=0; j<3; j++) { // j = 0 (track2particle), j = 1 (vertex), j = 2 (trigger)
1309 AlidNdEtaCorrection* alignedCorrection = new AlidNdEtaCorrection("dndeta_correction","dndeta_correction");
1310 alignedCorrection->LoadHistograms(alignedFile);
1311
1312 AlidNdEtaCorrection* misalignedCorrection = new AlidNdEtaCorrection("dndeta_correction","dndeta_correction");
1313 misalignedCorrection->LoadHistograms(misalignedFile);
1314
1315 TString name;
1316 name.Form("dndeta_correction_alignment_%s", typeName[j]);
1317 AlidNdEtaCorrection* current = new AlidNdEtaCorrection(name, name);
1318
1319 switch (j)
1320 {
1321 case 0:
1322 alignedCorrection->GetTrack2ParticleCorrection()->Reset();
1323 misalignedCorrection->GetTriggerBiasCorrectionNSD()->Reset();
1324 misalignedCorrection->GetTriggerBiasCorrectionND()->Reset();
1325 misalignedCorrection->GetTriggerBiasCorrectionINEL()->Reset();
1326 misalignedCorrection->GetVertexRecoCorrection()->Reset();
1327 break;
1328
1329 case 1:
1330 misalignedCorrection->GetTrack2ParticleCorrection()->Reset();
1331 misalignedCorrection->GetTriggerBiasCorrectionNSD()->Reset();
1332 misalignedCorrection->GetTriggerBiasCorrectionND()->Reset();
1333 misalignedCorrection->GetTriggerBiasCorrectionINEL()->Reset();
1334 alignedCorrection->GetVertexRecoCorrection()->Reset();
1335 break;
1336
1337 case 2:
1338 misalignedCorrection->GetTrack2ParticleCorrection()->Reset();
1339 alignedCorrection->GetTriggerBiasCorrectionNSD()->Reset();
1340 alignedCorrection->GetTriggerBiasCorrectionND()->Reset();
1341 alignedCorrection->GetTriggerBiasCorrectionINEL()->Reset();
1342 misalignedCorrection->GetVertexRecoCorrection()->Reset();
1343 break;
1344
1345 default:
1346 return;
1347 }
1348
1349 TList collection;
1350 collection.Add(misalignedCorrection);
1351 collection.Add(alignedCorrection);
1352
1353 current->Merge(&collection);
1354 current->Finish();
1355
1356 corrections[j] = current;
1357 }
1358
1359 TFile* fout = new TFile(outputFileName, "RECREATE");
1360
1361 for (Int_t i=0; i<3; i++)
1362 corrections[i]->SaveHistograms();
1363
1364 fout->Write();
1365 fout->Close();
1366}
1367
1368
1369void DrawdNdEtaDifferencesAlignment()
1370{
1371 TH1* hists[5];
1372
1373 TLegend* legend = new TLegend(0.3, 0.73, 0.70, 0.98);
1374 legend->SetFillColor(0);
1375
1376 TCanvas* canvas = new TCanvas("DrawdNdEtaDifferences", "DrawdNdEtaDifferences", 1000, 500);
1377 canvas->Divide(2, 1);
1378
1379 canvas->cd(1);
1380
1381 for (Int_t i=0; i<5; ++i)
1382 {
1383 hists[i] = 0;
1384 TFile* file = 0;
1385 TString title;
1386
1387 switch(i)
1388 {
1389 case 0 : file = TFile::Open("systematics_misalignment_aligned.root"); title = "aligned"; break;
1390 case 1 : file = TFile::Open("systematics_misalignment_misaligned.root"); title = "fully misaligned"; break;
1391 case 2 : file = TFile::Open("systematics_misalignment_track2particle.root"); title = "only track2particle"; break;
1392 case 3 : file = TFile::Open("systematics_misalignment_vertex.root"); title = "only vertex rec."; break;
1393 case 4 : file = TFile::Open("systematics_misalignment_trigger.root"); title = "only trigger bias"; break;
1394 default: return;
1395 }
1396
1397 if (file)
1398 {
1399 hists[i] = (TH1*) file->Get("dndeta/dndeta_dNdEta_corrected_2");
1400 hists[i]->SetTitle("");
1401
1402 Prepare1DPlot(hists[i], kFALSE);
1403 hists[i]->GetXaxis()->SetRangeUser(-0.7999, 0.7999);
1404 hists[i]->GetYaxis()->SetRangeUser(6, 7);
1405 hists[i]->SetLineWidth(1);
1406 hists[i]->SetLineColor(colors[i]);
1407 hists[i]->SetMarkerColor(colors[i]);
1408 hists[i]->SetMarkerStyle(markers[i]);
1409 hists[i]->GetXaxis()->SetLabelOffset(0.015);
1410 hists[i]->GetYaxis()->SetTitleOffset(1.5);
1411 gPad->SetLeftMargin(0.12);
1412 hists[i]->DrawCopy(((i > 0) ? "SAME" : ""));
1413
1414 legend->AddEntry(hists[i], title);
1415 hists[i]->SetTitle(title);
1416 }
1417 }
1418 legend->Draw();
1419
1420 canvas->cd(2);
1421 gPad->SetLeftMargin(0.14);
1422
1423 TLegend* legend2 = new TLegend(0.63, 0.73, 0.98, 0.98);
1424 legend2->SetFillColor(0);
1425
1426 for (Int_t i=1; i<5; ++i)
1427 {
1428 if (hists[i])
1429 {
1430 legend2->AddEntry(hists[i]);
1431
1432 hists[i]->Divide(hists[0]);
1433 hists[i]->SetTitle("b)");
1434 hists[i]->GetYaxis()->SetRangeUser(0.9, 1.1);
1435 hists[i]->GetYaxis()->SetTitle("Ratio to standard composition");
1436 hists[i]->GetYaxis()->SetTitleOffset(1.8);
1437 hists[i]->DrawCopy(((i > 1) ? "SAME" : ""));
1438 }
1439 }
1440
1441 legend2->Draw();
1442
1443 canvas->SaveAs("misalignment_result.eps");
1444 canvas->SaveAs("misalignment_result.gif");
1445}
1446
1447void EvaluateMultiplicityEffect()
1448{
1449 gSystem->Load("libPWG0base");
1450
1451 AlidNdEtaCorrection* fdNdEtaCorrectionLow[4];
1452 AlidNdEtaCorrection* fdNdEtaCorrectionHigh[4];
1453
1454 TFile::Open("systematics-low-multiplicity.root");
1455
1456 for (Int_t i=0; i<4; ++i)
1457 {
1458 TString name;
1459 name.Form("correction_%d", i);
1460 fdNdEtaCorrectionLow[i] = new AlidNdEtaCorrection(name, name);
1461 fdNdEtaCorrectionLow[i]->LoadHistograms("systematics-low-multiplicity.root", name);
1462 }
1463
1464 TList list;
1465 for (Int_t i=1; i<4; ++i)
1466 list.Add(fdNdEtaCorrectionLow[i]);
1467
1468 fdNdEtaCorrectionLow[0]->Merge(&list);
1469 fdNdEtaCorrectionLow[0]->Finish();
1470
1471 TFile::Open("systematics-high-multiplicity.root");
1472
1473 for (Int_t i=0; i<4; ++i)
1474 {
1475 TString name;
1476 name.Form("correction_%d", i);
1477 fdNdEtaCorrectionHigh[i] = new AlidNdEtaCorrection(name, name);
1478 fdNdEtaCorrectionHigh[i]->LoadHistograms("systematics-high-multiplicity.root", name);
1479 }
1480
1481 TList list2;
1482 for (Int_t i=1; i<4; ++i)
1483 list2.Add(fdNdEtaCorrectionHigh[i]);
1484
1485 fdNdEtaCorrectionHigh[0]->Merge(&list2);
1486 fdNdEtaCorrectionHigh[0]->Finish();
1487
1488 TH1F* outputLow = new TH1F("Track2ParticleLow", "Track2Particle at low multiplicity", 200, 0, 2);
1489 TH1F* outputHigh = new TH1F("Track2ParticleHigh", "Track2Particle at high multiplicity", 200, 0, 2);
1490
1491 TH3F* hist = fdNdEtaCorrectionLow[0]->GetTrack2ParticleCorrection()->GetCorrectionHistogram();
1492 TH3F* hist2 = fdNdEtaCorrectionHigh[0]->GetTrack2ParticleCorrection()->GetCorrectionHistogram();
1493 for (Int_t x=hist->GetXaxis()->FindBin(-10); x<=hist->GetXaxis()->FindBin(10); ++x)
1494 for (Int_t y=hist->GetYaxis()->FindBin(-0.8); y<=hist->GetYaxis()->FindBin(0.8); ++y)
1495 for (Int_t z=hist->GetZaxis()->FindBin(0.3); z<=hist->GetZaxis()->FindBin(9.9); ++z)
1496 //for (Int_t z=1; z<=hist->GetNbinsZ(); ++z)
1497 {
1498 if (hist->GetBinContent(x, y, z) > 0)
1499 outputLow->Fill(hist->GetBinContent(x, y, z));
1500 //if (hist->GetBinContent(x, y, z) == 1)
1501 // printf("z = %f, eta = %f, pt = %f: %f %f %f\n", hist->GetXaxis()->GetBinCenter(x), hist->GetYaxis()->GetBinCenter(y), hist->GetZaxis()->GetBinCenter(z), hist->GetBinContent(x, y, z), fdNdEtaCorrectionLow[0]->GetTrack2ParticleCorrection()->GetGeneratedHistogram()->GetBinContent(x, y, z), fdNdEtaCorrectionLow[0]->GetTrack2ParticleCorrection()->GetMeasuredHistogram()->GetBinContent(x, y, z));
1502
1503 if (hist2->GetBinContent(x, y, z) > 0)
1504 outputHigh->Fill(hist2->GetBinContent(x, y, z));
1505 }
1506
1507 TCanvas* canvas = new TCanvas("EvaluateMultiplicityEffect", "EvaluateMultiplicityEffect", 1000, 500);
1508 canvas->Divide(2, 1);
1509
1510 canvas->cd(1);
1511 outputLow->Draw();
1512 outputLow->Fit("gaus", "0");
1513 outputLow->GetFunction("gaus")->SetLineColor(2);
1514 outputLow->GetFunction("gaus")->DrawCopy("SAME");
1515
1516 canvas->cd(2);
1517 outputHigh->Draw();
1518 outputHigh->Fit("gaus", "0");
1519 outputHigh->GetFunction("gaus")->DrawCopy("SAME");
1520
1521 canvas->SaveAs(Form("%s.gif", canvas->GetName()));
1522}
1523
1524void PlotErrors(Float_t xmin, Float_t xmax, Float_t ymin, Float_t ymax, Float_t zmin, Float_t zmax, const char* correctionMapFile = "correction_map.root", const char* correctionMapFolder = "dndeta_correction")
1525{
1526 gSystem->Load("libPWG0base");
1527
1528 AlidNdEtaCorrection* dNdEtaCorrection = new AlidNdEtaCorrection(correctionMapFolder, correctionMapFolder);
1529 dNdEtaCorrection->LoadHistograms(correctionMapFile, correctionMapFolder);
1530
1531 dNdEtaCorrection->GetTrack2ParticleCorrection()->PlotBinErrors(xmin, xmax, ymin, ymax, zmin, zmax)->Draw();
1532}
9e952c39 1533
1534
1535void runStudy(const char* baseCorrectionMapFile = "correction_map.root", const char* baseCorrectionMapFolder = "dndeta_correction", const char* changedCorrectionMapFile = "correction_map.root", const char* changedCorrectionMapFolder = "dndeta_correction", const char* dataFile = "analysis_esd_raw.root", const char* output = "analysis_esd_syst.root")
1536{
1537 gSystem->Load("libPWG0base");
1538
1539 TFile* file = TFile::Open(output, "RECREATE");
1540
1541 const Int_t max = 5;
1542 dNdEtaAnalysis* fdNdEtaAnalysis[5];
1543
1544 new TCanvas;
1545 TLegend* legend = new TLegend(0.63, 0.73, 0.98, 0.98);
1546 legend->SetFillColor(0);
1547
1548 for (Int_t i = 0; i < max; ++i)
1549 {
1550 TFile::Open(baseCorrectionMapFile);
1551 AlidNdEtaCorrection* baseCorrection = new AlidNdEtaCorrection(baseCorrectionMapFolder, baseCorrectionMapFolder);
1552 baseCorrection->LoadHistograms();
1553
1554 AlidNdEtaCorrection::CorrectionType correctionType = AlidNdEtaCorrection::kNone;
1555 const char* name = 0;
1556
1557 TFile::Open(changedCorrectionMapFile);
1558 switch (i)
1559 {
1560 case 0 :
1561 name = "default";
1562 break;
1563
1564 case 1 :
1565 baseCorrection->GetTrack2ParticleCorrection()->LoadHistograms(Form("%s/Track2Particle", changedCorrectionMapFolder));
1566 name = "Track2Particle";
1567 break;
1568
1569 case 2 :
1570 baseCorrection->GetVertexRecoCorrection()->LoadHistograms(Form("%s/VertexReconstruction", changedCorrectionMapFolder));
1571 name = "VertexReco";
1572 break;
1573
1574 case 3 :
1575 baseCorrection->GetTriggerBiasCorrectionINEL()->LoadHistograms(Form("%s/TriggerBias_MBToINEL", changedCorrectionMapFolder));
1576 name = "TriggerBias_MBToINEL";
1577 break;
1578
1579 case 4 :
1580 baseCorrection->LoadHistograms(changedCorrectionMapFolder);
1581 name = "all";
1582 break;
1583
1584 default: return;
1585 }
1586
1587 TFile::Open(dataFile);
1588 fdNdEtaAnalysis[i] = new dNdEtaAnalysis(name, name);
1589 fdNdEtaAnalysis[i]->LoadHistograms("dndeta");
1590
1591 fdNdEtaAnalysis[i]->Finish(baseCorrection, 0.3, AlidNdEtaCorrection::kINEL);
1592 file->cd();
1593 fdNdEtaAnalysis[i]->SaveHistograms();
1594
1595 TH1* hist = fdNdEtaAnalysis[i]->GetdNdEtaHistogram(0);
1596 hist->SetLineColor(colors[i]);
1597 hist->SetMarkerColor(colors[i]);
1598 hist->SetMarkerStyle(markers[i]+1);
1599 hist->DrawCopy((i == 0) ? "" : "SAME");
1600 legend->AddEntry(hist, name);
1601 }
1602
1603 legend->Draw();
1604}