e140a1242ef4670f966185196ec19a18fb00a06e
[u/mrichter/AliRoot.git] / PWG0 / dNdEta / drawSystematics.C
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>
16
17 #endif
18
19 void SetRanges(TAxis* axis)
20 {
21   if (strcmp(axis->GetTitle(), "#eta") == 0)
22     axis->SetRangeUser(-1.7999, 1.7999);
23   if (strcmp(axis->GetTitle(), "p_{T} [GeV/c]") == 0)
24     axis->SetRangeUser(0, 9.9999);
25   if (strcmp(axis->GetTitle(), "vtx z [cm]") == 0)
26     axis->SetRangeUser(-15, 14.9999);
27   if (strcmp(axis->GetTitle(), "Ntracks") == 0)
28     axis->SetRangeUser(0, 99.9999);
29 }
30
31 void SetRanges(TH1* hist)
32 {
33   SetRanges(hist->GetXaxis());
34   SetRanges(hist->GetYaxis());
35   SetRanges(hist->GetZaxis());
36 }
37
38 void Prepare3DPlot(TH3* hist)
39 {
40   hist->GetXaxis()->SetTitleOffset(1.5);
41   hist->GetYaxis()->SetTitleOffset(1.5);
42   hist->GetZaxis()->SetTitleOffset(1.5);
43
44   hist->SetStats(kFALSE);
45 }
46
47 void Prepare2DPlot(TH2* hist)
48 {
49   hist->SetStats(kFALSE);
50   hist->GetYaxis()->SetTitleOffset(1.4);
51
52   SetRanges(hist);
53 }
54
55 void Prepare1DPlot(TH1* hist)
56 {
57   hist->SetLineWidth(2);
58   hist->SetStats(kFALSE);
59
60   SetRanges(hist);
61 }
62
63 void InitPad()
64 {
65   gPad->Range(0, 0, 1, 1);
66   gPad->SetLeftMargin(0.15);
67   //gPad->SetRightMargin(0.05);
68   //gPad->SetTopMargin(0.13);
69   //gPad->SetBottomMargin(0.1);
70
71   //gPad->SetGridx();
72   //gPad->SetGridy();
73 }
74
75 void InitPadCOLZ()
76 {
77   gPad->Range(0, 0, 1, 1);
78   gPad->SetRightMargin(0.15);
79   gPad->SetLeftMargin(0.12);
80
81   gPad->SetGridx();
82   gPad->SetGridy();
83 }
84
85 void Secondaries()
86 {
87   TFile* file = TFile::Open("systematics.root");
88
89   TH3F* secondaries = dynamic_cast<TH3F*> (file->Get("fSecondaries"));
90   if (!secondaries)
91   {
92     printf("Could not read histogram\n");
93     return;
94   }
95
96   for (Int_t ptBin=1; ptBin<=secondaries->GetNbinsZ(); ptBin++)
97   //for (Int_t ptBin = 1; ptBin<=2; ptBin++)
98   {
99     TGraph* graph = new TGraph;
100     graph->Clear();
101     graph->SetTitle(Form("%f < p_{T} < %f", secondaries->GetZaxis()->GetBinLowEdge(ptBin), secondaries->GetZaxis()->GetBinUpEdge(ptBin)));
102
103     for (Int_t cBin=1; cBin<=secondaries->GetNbinsY(); ++cBin)
104     {
105       if (secondaries->GetBinContent(0, cBin, ptBin) > 0)
106         printf("WARNING: Underflow bin not empty!");
107       if (secondaries->GetBinContent(secondaries->GetNbinsX()+1, cBin, ptBin) > 0)
108         printf("WARNING: Overflow bin not empty!");
109
110       Double_t sum = 0;
111       Double_t count = 0;
112       for (Int_t nBin=1; nBin<=secondaries->GetNbinsX(); ++nBin)
113       {
114         //printf("%f %f\n", secondaries->GetXaxis()->GetBinCenter(nBin), secondaries->GetBinContent(nBin, cBin, ptBin));
115         sum += secondaries->GetXaxis()->GetBinCenter(nBin) * secondaries->GetBinContent(nBin, cBin, ptBin);
116         count += secondaries->GetBinContent(nBin, cBin, ptBin);
117       }
118
119       printf("%f %f\n", sum, count);
120
121       if (count > 0)
122         graph->SetPoint(graph->GetN(), secondaries->GetYaxis()->GetBinCenter(cBin), sum / count);
123     }
124
125     new TCanvas;
126     graph->SetMarkerStyle(21);
127     graph->Draw("AP");
128     graph->Print();
129   }
130 }
131
132 void Composition()
133 {
134   gSystem->Load("libPWG0base");
135
136   AlidNdEtaCorrection* fdNdEtaCorrection[4];
137
138   TFile::Open("systematics.root");
139
140   for (Int_t i=0; i<4; ++i)
141   {
142     TString name;
143     name.Form("correction_%d", i);
144     fdNdEtaCorrection[i] = new AlidNdEtaCorrection(name, name);
145     fdNdEtaCorrection[i]->LoadHistograms("systematics.root", name);
146   }
147
148   //fdNdEtaCorrection[0]->GetTrack2ParticleCorrection()->GetMeasuredHistogram()->Scale(2);
149   //fdNdEtaCorrection[0]->GetTrack2ParticleCorrection()->GetGeneratedHistogram()->Scale(2);
150
151   AlidNdEtaCorrection* finalCorrection = new AlidNdEtaCorrection("dndeta_correction", "dndeta_correction");
152
153   TList* collection = new TList;
154
155   for (Int_t i=0; i<4; ++i)
156     collection->Add(fdNdEtaCorrection[i]);
157
158   finalCorrection->Merge(collection);
159
160   delete collection;
161
162   finalCorrection->Finish();
163
164   TFile* file = TFile::Open("temp.root", "RECREATE");
165   finalCorrection->SaveHistograms();
166   file->Write();
167   file->Close();
168
169   gROOT->ProcessLine(".L drawPlots.C");
170   Track2Particle1D("temp.root");
171 }
172
173 void drawSystematics()
174 {
175   Composition();
176 }