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