Improved eloss fitting - see NIM B1, 16
[u/mrichter/AliRoot.git] / PWG2 / FORWARD / analysis2 / scripts / DrawFits.C
CommitLineData
c389303e 1#include <TFile.h>
2#include <THStack.h>
3#include <TList.h>
4#include <TError.h>
5#include <TCanvas.h>
6#include <TPad.h>
7#include <TStyle.h>
8#include <TF1.h>
9#include <TLegend.h>
10#include <TMath.h>
11
12//____________________________________________________________________
13// Global
14TList* fitter = 0;
15TCanvas* canvas = 0;
16const char* pdfName = "FitResults.pdf";
17
18//____________________________________________________________________
19TList* OpenFile(const char* fname)
20{
21 TFile* file = TFile::Open(fname, "READ");
22 if (!file) {
23 Error("DrawFits", "Couldn't open %s", fname);
24 return 0;
25 }
26
27 TList* forward =
28 static_cast<TList*>(file->Get("PWG2forwardDnDeta/Forward"));
29 if (!forward) {
30 Error("DrawFits", "Couldn't get forward list from %s", fname);
31 return 0;
32 }
33
34 fitter = static_cast<TList*>(forward->FindObject("fmdEnergyFitter"));
35 if (!fitter) {
36 Error("DrawFits", "Couldn't get fitter folder");
37 return 0;
38 }
39 return fitter;
40}
41//____________________________________________________________________
42TList* CheckFitter(const char* fname="AnalysisResults.root")
43{
44 if (!fitter) return OpenFile(fname);
45 return fitter;
46}
47
48//____________________________________________________________________
49TCanvas* CheckCanvas()
50{
51 if (canvas) return canvas;
52
53 gStyle->SetOptFit(111111);
54 gStyle->SetTitleFillColor(0);
55 gStyle->SetTitleBorderSize(0);
56 gStyle->SetTitleX(0);
57 gStyle->SetTitleY(.9);
58 gStyle->SetTitleW(.4);
59 gStyle->SetTitleH(.1);
60 gStyle->SetStatW(0.2);
61 gStyle->SetStatH(0.2);
62 gStyle->SetStatColor(0);
63 gStyle->SetStatBorderSize(1);
64 gStyle->SetOptTitle(0);
65 gStyle->SetOptFit(1111);
66 gStyle->SetStatW(0.3);
67 gStyle->SetStatH(0.15);
68 gStyle->SetStatColor(0);
69 gStyle->SetStatBorderSize(1);
70
71 canvas = new TCanvas("c", "C", Int_t(800 / TMath::Sqrt(2)), 800);
72 canvas->SetFillColor(0);
73 canvas->SetBorderMode(0);
74 canvas->SetBorderSize(0);
75 canvas->SetBottomMargin(0.15);
76 return canvas;
77}
78
79//____________________________________________________________________
80void DrawSummary(const char* fname="AnalysisResults.root")
81{
82 if (!CheckFitter(fname)) {
83 Error("DrawFits", "File not opened");
84 return;
85 }
86 if (!CheckCanvas()) {
87 Error("DrawFits", "No canvas");
88 return;
89 }
90 canvas->Clear();
91
92 THStack* chi2nu;
93 THStack* c;
94 THStack* delta;
95 THStack* xi;
96 THStack* sigma;
97 THStack* sigman;
98 THStack* n;
99 TList stacks;
100 stacks.Add(chi2nu = static_cast<THStack*>(fitter->FindObject("chi2")));
101 stacks.Add(c = static_cast<THStack*>(fitter->FindObject("c")));
102 stacks.Add(delta = static_cast<THStack*>(fitter->FindObject("delta")));
103 stacks.Add(xi = static_cast<THStack*>(fitter->FindObject("xi")));
104 stacks.Add(sigma = static_cast<THStack*>(fitter->FindObject("sigma")));
105 stacks.Add(sigman = static_cast<THStack*>(fitter->FindObject("sigman")));
106 stacks.Add(n = static_cast<THStack*>(fitter->FindObject("n")));
107 Int_t baseA = stacks.GetEntries()+1;
108 Int_t i=2;
109 while (true) {
110 TObject* o = fitter->FindObject(Form("a%d",i++));
111 if (!o) break;
112 Info("DrawFits", "Adding %s", o->GetName());
113 stacks.Add(o);
114 }
115 // stacks.ls();
116 Int_t nMax = stacks.GetEntries();
117 for (Int_t i = nMax-1; i >= baseA; i--) {
118 THStack* stack = static_cast<THStack*>(stacks.At(i));
119 TIter nextH(stack->GetHists());
120 TH1* hist = 0;
121 Bool_t hasData = kFALSE;
122 while ((hist = static_cast<TH1*>(nextH())))
123 if (hist->Integral() > 0) hasData = kTRUE;
124 if (!hasData) nMax--;
125 }
126
127 canvas->SetRightMargin(0.05);
128 canvas->SetTopMargin(0.05);
129 canvas->Divide(2, (nMax+1)/2, 0.1, 0, 0);
130
131 TIter next(&stacks);
132 THStack* stack = 0;
133 i = 0;
134 Int_t b = 1;
135 while ((stack = static_cast<THStack*>(next()))) {
136 if (i > nMax) break;
137 TVirtualPad* p = canvas->cd(1+i/5 + 2*(i%5));
138 p->SetLeftMargin(.15);
139 p->SetFillColor(0);
140 p->SetFillStyle(0);
141 p->SetGridx();
142 stack->Draw("nostack");
143 stack->GetHistogram()->SetYTitle(stack->GetTitle());
144 stack->GetHistogram()->SetXTitle("#eta");
145
146 TAxis* yaxis = stack->GetHistogram()->GetYaxis();
147 if (i == 0) yaxis->SetRangeUser(0,20); // Chi2
148 if (i == 1) stack->SetMaximum(1); // c
149 if (i == 2) stack->SetMaximum(1); // delta
150 if (i == 3) stack->SetMaximum(0.1); // xi
151 if (i == 4 || i == 5) stack->SetMaximum(0.5); // sigma{,n}
152 if (i == 7) stack->SetMaximum(0.5); // a
153 yaxis->SetTitleSize(0.15);
154 yaxis->SetLabelSize(0.08);
155 yaxis->SetTitleOffset(0.35);
156 yaxis->SetNdivisions(5);
157
158 TAxis* xaxis = stack->GetHistogram()->GetXaxis();
159 xaxis->SetTitleSize(0.15);
160 xaxis->SetLabelSize(0.08);
161 xaxis->SetTitleOffset(0.35);
162 xaxis->SetNdivisions(10);
163
164 // Redraw
165 stack->Draw("nostack");
166 i++;
167 if (i >= 5) b = 2;
168 p->cd();
169 }
170 canvas->Print(pdfName, "Title:Fit summary");
171}
172
173//____________________________________________________________________
174void DrawRings(const char* fname="AnalysisResults.root")
175{
176 if (!CheckFitter(fname)) {
177 Error("DrawFits", "File not opened");
178 return;
179 }
180 if (!CheckCanvas()) {
181 Error("DrawFits", "No canvas");
182 return;
183 }
184 canvas->Clear();
185
186 canvas->Clear();
187 canvas->Divide(1, 5,0,0);
188
189 const char* dets[] = { "FMD1I", "FMD2I", "FMD2O", "FMD3I", "FMD3O", 0 };
190 for (Int_t i = 0; i < 5; i++) {
191 TVirtualPad* p = canvas->cd(i+1);
192 p->SetGridx();
193 p->SetFillColor(0);
194 p->SetFillStyle(0);
195 TList* d = static_cast<TList*>(fitter->FindObject(dets[i]));
196 if (!d) {
197 Warning("DrawFits", "List %s not found", dets[i]);
198 continue;
199 }
200 TH1* edist = static_cast<TH1*>(d->FindObject(Form("%s_edist", dets[i])));
201 if (!edist) {
202 Warning("DrawFits", "Histogram %s_edist not found", dets[i]);
203 continue;
204 }
205 edist->Draw();
206 TF1* f = 0;
207 TIter nextF(edist->GetListOfFunctions());
208 while ((f = static_cast<TF1*>(nextF()))) {
209 Double_t chi2 = f->GetChisquare();
210 Int_t ndf = f->GetNDF();
211 Printf("%s %s:\n Range: %f-%f\n"
212 "chi^2/ndf= %f / %d = %f",
213 edist->GetName(), f->GetName(),
214 f->GetXmin(), f->GetXmax(), chi2, ndf,
215 (ndf > 0) ? chi2/ndf : 0);
216 for (Int_t j = 0; j < f->GetNpar(); j++) {
217 Printf(" %-20s : %9.4f +/- %9.4f",
218 f->GetParName(j), f->GetParameter(j), f->GetParError(j));
219 }
220 }
221 p->cd();
222 }
223 canvas->cd();
224 canvas->Print(pdfName, "Title:Fit to rings");
225}
226
227//____________________________________________________________________
228void DrawEtaBins(const char* fname="AnalysisResults.root")
229{
230 if (!CheckFitter(fname)) {
231 Error("DrawFits", "File not opened");
232 return;
233 }
234 if (!CheckCanvas()) {
235 Error("DrawFits", "No canvas");
236 return;
237 }
238 canvas->Clear();
239 canvas->Divide(2,2,0,0);
240
241 for (UShort_t d=1; d<=3; d++) {
242 UShort_t nQ = (d == 1 ? 1 : 2);
243 for (UShort_t q = 0; q < nQ; q++) {
244 Char_t r = (q == 0 ? 'I' : 'O');
245
246 TList* ring =
247 static_cast<TList*>(fitter->FindObject(Form("FMD%d%c",d,r)));
248 if (!ring) {
249 Error("PrintFits", "Couldn't get ring FMD%d%c", d,r);
250 continue;
251 }
252 TList* edists = static_cast<TList*>(ring->FindObject("EDists"));
253 if (!edists) {
254 Error("PrintFits", "Couldn't get EDists list for FMD%d%c", d,r);
255 continue;
256 }
257 TIter next(edists);
258 TH1* dist = 0;
259 Int_t i = 0;
260 Int_t j = 1;
261 while ((dist = static_cast<TH1*>(next()))) {
262 if (i == 4) {
263 i = 0;
264 j++;
265 canvas->Print(pdfName, Form("Title:FMD%d%c page %2d", d,r,j));
266 }
267 TVirtualPad* p = canvas->cd(++i);
268 p->SetFillColor(kWhite);
269 p->SetFillStyle(0);
270 p->SetBorderSize(0);
271 p->SetLogy();
272 dist->SetMaximum(15);
273 dist->Draw();
274
275 }
276 if (i != 0) {
277 i++;
278 for (; i <= 4; i++) {
279 TVirtualPad* p = canvas->cd(i);
280 p->Clear();
281 p->SetFillColor(kMagenta-3);
282 p->SetFillStyle(0);
283 p->SetBorderSize(0);
284 }
285 canvas->Print(pdfName, Form("FMD%d%c page %2d", d,r,j++));
286 }
287 }
288 }
289}
290
291//____________________________________________________________________
292void
293DrawFits(const char* fname="AnalysisResults.root")
294{
295 if (!CheckCanvas()) {
296 Error("DrawFits", "No canvas");
297 return;
298 }
299 canvas->Print(Form("%s[", pdfName));
300 DrawSummary(fname);
301 DrawRings(fname);
302 DrawEtaBins(fname);
303 canvas->Print(Form("%s]", pdfName));
304}