Unused variable i removed
[u/mrichter/AliRoot.git] / PWG1 / AliTreeDraw.cxx
CommitLineData
c92725b7 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16
17///////////////////////////////////////////////////////////////////////////
18/*
19
20Origin: marian.ivanov@cern.ch
21Frequenlty used function for visualization
22marian.ivanov@cern.ch
23*/
24
25#if !defined(__CINT__) || defined(__MAKECINT__)
26#include <stdio.h>
27#include <string.h>
28//ROOT includes
29#include "TROOT.h"
30#include "Rtypes.h"
31#include "TFile.h"
32#include "TTree.h"
33#include "TChain.h"
34#include "TCut.h"
35#include "TString.h"
36#include "TBenchmark.h"
37#include "TStopwatch.h"
38#include "TParticle.h"
39#include "TSystem.h"
40#include "TTimer.h"
41#include "TVector3.h"
42#include "TH1F.h"
43#include "TH2F.h"
44#include "TCanvas.h"
45#include "TPad.h"
46#include "TF1.h"
47#include "TView.h"
48#include "TView3D.h"
49#include "TGeometry.h"
50#include "TPolyLine3D.h"
51#include "TPolyMarker3D.h"
52
53//ALIROOT includes
54#include "AliTrackPointArray.h"
55#include "AliTreeDraw.h"
56
57#endif
58
59//
60// Class for visualization and some statistacal analysis using tree
61// To be used in comparisons
62// and calib viewers based on tree
63
64
65ClassImp(AliTreeDraw)
66
67
68AliTreeDraw::AliTreeDraw():
69 fTree(0),
70 fRes(0),
71 fMean(0),
72 fPoints(0){
73 //
74 // default constructor
75 //
76}
77
78void AliTreeDraw::ClearHisto(){
79 //
80 //
81 delete fRes;
82 delete fMean;
83 fRes=0;
84 fMean=0;
85}
86
87
88
89TH1F * AliTreeDraw::DrawXY(const char * chx, const char *chy, const char* selection,
90 const char * quality, Int_t nbins, Float_t minx, Float_t maxx, Float_t miny, Float_t maxy, Int_t nBinsRes)
91{
92 //
93 Double_t* bins = CreateLogBins(nbins, minx, maxx);
94 TH2F* hRes2 = new TH2F("hRes2", "residuals", nbins, minx, maxx, nBinsRes, miny, maxy);
95 char cut[1000];
96 sprintf(cut,"%s&&%s",selection,quality);
97 char expression[1000];
98 sprintf(expression,"%s:%s>>hRes2",chy,chx);
99 fTree->Draw(expression, cut, "groff");
100 TH1F* hMean=0;
101 TH1F* hRes = CreateResHisto(hRes2, &hMean);
102 AliLabelAxes(hRes, chx, chy);
103 //
104 delete hRes2;
105 delete[] bins;
106 ClearHisto();
107 fRes = hRes;
108 fMean = hMean;
109 return hRes;
110}
111
112
113
114TH1F * AliTreeDraw::DrawLogXY(const char * chx, const char *chy, const char* selection,
115 const char * quality, Int_t nbins, Float_t minx, Float_t maxx, Float_t miny, Float_t maxy, Int_t nBinsRes)
116{
117 //
118 //
119 //
120 Double_t* bins = CreateLogBins(nbins, minx, maxx);
121 TH2F* hRes2 = new TH2F("hRes2", "residuals", nbins, bins, nBinsRes, miny, maxy);
122 char cut[1000];
123 sprintf(cut,"%s&&%s",selection,quality);
124 char expression[1000];
125 sprintf(expression,"%s:%s>>hRes2",chy,chx);
126 fTree->Draw(expression, cut, "groff");
127 TH1F* hMean=0;
128 TH1F* hRes = CreateResHisto(hRes2, &hMean);
129 AliLabelAxes(hRes, chx, chy);
130 //
131 delete hRes2;
132 delete[] bins;
133 ClearHisto();
134 fRes = hRes;
135 fMean = hMean;
136 return hRes;
137}
138
139///////////////////////////////////////////////////////////////////////////////////
140///////////////////////////////////////////////////////////////////////////////////
141TH1F * AliTreeDraw::Eff(const char *variable, const char* selection, const char * quality,
142 Int_t nbins, Float_t min, Float_t max)
143{
144 //
145 //
146 TH1F* hGen = new TH1F("hGen", "gen. tracks", nbins, min, max);
147 TH1F* hRec = new TH1F("hRec", "rec. tracks", nbins, min, max);
148 char inputGen[1000];
149 sprintf(inputGen,"%s>>hGen", variable);
150 fTree->Draw(inputGen, selection, "groff");
151 char selectionRec[256];
152 sprintf(selectionRec, "%s && %s", selection, quality);
153 char inputRec[1000];
154 sprintf(inputRec,"%s>>hRec", variable);
155 fTree->Draw(inputRec, selectionRec, "groff");
156 //
157 TH1F* hEff = CreateEffHisto(hGen, hRec);
158 AliLabelAxes(hEff, variable, "#epsilon [%]");
159 fRes = hEff;
160 delete hRec;
161 delete hGen;
162 return hEff;
163}
164
165
166
167///////////////////////////////////////////////////////////////////////////////////
168///////////////////////////////////////////////////////////////////////////////////
169TH1F * AliTreeDraw::EffLog(const char *variable, const char* selection, const char * quality,
170 Int_t nbins, Float_t min, Float_t max)
171{
172 //
173 //
174 Double_t* bins = CreateLogBins(nbins, min, max);
175 TH1F* hGen = new TH1F("hGen", "gen. tracks", nbins, bins);
176 TH1F* hRec = new TH1F("hRec", "rec. tracks", nbins, bins);
177 char inputGen[1000];
178 sprintf(inputGen,"%s>>hGen", variable);
179 fTree->Draw(inputGen, selection, "groff");
180 char selectionRec[256];
181 sprintf(selectionRec, "%s && %s", selection, quality);
182 char inputRec[1000];
183 sprintf(inputRec,"%s>>hRec", variable);
184 fTree->Draw(inputRec, selectionRec, "groff");
185 //
186 TH1F* hEff = CreateEffHisto(hGen, hRec);
187 AliLabelAxes(hEff, variable, "#epsilon [%]");
188 fRes = hEff;
189 delete hRec;
190 delete hGen;
191 delete[] bins;
192 return hEff;
193}
194
195
196///////////////////////////////////////////////////////////////////////////////////
197///////////////////////////////////////////////////////////////////////////////////
198
199Double_t* AliTreeDraw::CreateLogBins(Int_t nBins, Double_t xMin, Double_t xMax)
200{
201 Double_t* bins = new Double_t[nBins+1];
202 bins[0] = xMin;
203 Double_t factor = pow(xMax/xMin, 1./nBins);
204 for (Int_t i = 1; i <= nBins; i++)
205 bins[i] = factor * bins[i-1];
206 return bins;
207}
208
209
210
211
212void AliTreeDraw::AliLabelAxes(TH1* histo, const char* xAxisTitle, const char* yAxisTitle)
213{
214 //
215 histo->GetXaxis()->SetTitle(xAxisTitle);
216 histo->GetYaxis()->SetTitle(yAxisTitle);
217}
218
219
220TH1F* AliTreeDraw::CreateEffHisto(TH1F* hGen, TH1F* hRec)
221{
222 //
223 Int_t nBins = hGen->GetNbinsX();
224 TH1F* hEff = (TH1F*) hGen->Clone("hEff");
225 hEff->SetTitle("");
226 hEff->SetStats(kFALSE);
227 hEff->SetMinimum(0.);
228 hEff->SetMaximum(110.);
229 //
230 for (Int_t iBin = 0; iBin <= nBins; iBin++) {
231 Double_t nGen = hGen->GetBinContent(iBin);
232 Double_t nRec = hRec->GetBinContent(iBin);
233 if (nGen > 0) {
234 Double_t eff = nRec/nGen;
235 hEff->SetBinContent(iBin, 100. * eff);
236 Double_t error = sqrt((eff*(1.-eff)+0.01) / nGen);
237 // if (error == 0) error = sqrt(0.1/nGen);
238 //
239 if (error == 0) error = 0.0001;
240 hEff->SetBinError(iBin, 100. * error);
241 } else {
242 hEff->SetBinContent(iBin, 100. * 0.5);
243 hEff->SetBinError(iBin, 100. * 0.5);
244 }
245 }
246 return hEff;
247}
248
249
250
251TH1F* AliTreeDraw::CreateResHisto(TH2F* hRes2, TH1F **phMean, Bool_t drawBinFits,
252 Bool_t overflowBinFits)
253{
254 TVirtualPad* currentPad = gPad;
255 TAxis* axis = hRes2->GetXaxis();
256 Int_t nBins = axis->GetNbins();
257 TH1F* hRes, *hMean;
258 if (axis->GetXbins()->GetSize()){
259 hRes = new TH1F("hRes", "", nBins, axis->GetXbins()->GetArray());
260 hMean = new TH1F("hMean", "", nBins, axis->GetXbins()->GetArray());
261 }
262 else{
263 hRes = new TH1F("hRes", "", nBins, axis->GetXmin(), axis->GetXmax());
264 hMean = new TH1F("hMean", "", nBins, axis->GetXmin(), axis->GetXmax());
265
266 }
267 hRes->SetStats(false);
268 hRes->SetOption("E");
269 hRes->SetMinimum(0.);
270 //
271 hMean->SetStats(false);
272 hMean->SetOption("E");
273
274 // create the fit function
275 TF1 * fitFunc = new TF1("G","[0]*exp(-(x-[1])*(x-[1])/(2.*[2]*[2]))",-3,3);
276
277 fitFunc->SetLineWidth(2);
278 fitFunc->SetFillStyle(0);
279 // create canvas for fits
280 TCanvas* canBinFits = NULL;
281 Int_t nPads = (overflowBinFits) ? nBins+2 : nBins;
282 Int_t nx = Int_t(sqrt(nPads-1.));// + 1;
283 Int_t ny = (nPads-1) / nx + 1;
284 if (drawBinFits) {
285 canBinFits = (TCanvas*)gROOT->FindObject("canBinFits");
286 if (canBinFits) delete canBinFits;
287 canBinFits = new TCanvas("canBinFits", "fits of bins", 200, 100, 500, 700);
288 canBinFits->Divide(nx, ny);
289 }
290
291 // loop over x bins and fit projection
292 Int_t dBin = ((overflowBinFits) ? 1 : 0);
293 for (Int_t bin = 1-dBin; bin <= nBins+dBin; bin++) {
294 if (drawBinFits) canBinFits->cd(bin + dBin);
295 TH1D* hBin = hRes2->ProjectionY("hBin", bin, bin);
296 //
297 if (hBin->GetEntries() > 5) {
298 fitFunc->SetParameters(hBin->GetMaximum(),hBin->GetMean(),hBin->GetRMS());
299 hBin->Fit(fitFunc,"s");
300 Double_t sigma = TMath::Abs(fitFunc->GetParameter(2));
301
302 if (sigma > 0.){
303 hRes->SetBinContent(bin, TMath::Abs(fitFunc->GetParameter(2)));
304 hMean->SetBinContent(bin, fitFunc->GetParameter(1));
305 }
306 else{
307 hRes->SetBinContent(bin, 0.);
308 hMean->SetBinContent(bin,0);
309 }
310 hRes->SetBinError(bin, fitFunc->GetParError(2));
311 hMean->SetBinError(bin, fitFunc->GetParError(1));
312
313 //
314 //
315
316 } else {
317 hRes->SetBinContent(bin, 0.);
318 hRes->SetBinError(bin, 0.);
319 hMean->SetBinContent(bin, 0.);
320 hMean->SetBinError(bin, 0.);
321 }
322
323
324 if (drawBinFits) {
325 char name[256];
326 if (bin == 0) {
327 sprintf(name, "%s < %.4g", axis->GetTitle(), axis->GetBinUpEdge(bin));
328 } else if (bin == nBins+1) {
329 sprintf(name, "%.4g < %s", axis->GetBinLowEdge(bin), axis->GetTitle());
330 } else {
331 sprintf(name, "%.4g < %s < %.4g", axis->GetBinLowEdge(bin),
332 axis->GetTitle(), axis->GetBinUpEdge(bin));
333 }
334 canBinFits->cd(bin + dBin);
335 hBin->SetTitle(name);
336 hBin->SetStats(kTRUE);
337 hBin->DrawCopy("E");
338 canBinFits->Update();
339 canBinFits->Modified();
340 canBinFits->Update();
341 }
342
343 delete hBin;
344 }
345
346 delete fitFunc;
347 currentPad->cd();
348 *phMean = hMean;
349 return hRes;
350}
351
352
353
354
355void AliTreeDraw::GetPoints3D(const char * label, const char * chpoints,
356 const char* selection, TTree * tpoints, Int_t color,Float_t rmin){
357 //
358 // load selected points from tree
359 //
360 if (!fPoints) fPoints = new TObjArray;
361 if (tpoints->GetIndex()==0) tpoints->BuildIndex("fLabel","Label");
362 TBranch * br = tpoints->GetBranch(chpoints);
363 if (!br) return;
364 AliTrackPointArray * points = new AliTrackPointArray;
365 br->SetAddress(&points);
366
367 Int_t npoints = fTree->Draw(label,selection);
368 Float_t xyz[30000];
369 rmin*=rmin;
370 for (Int_t i=0;i<npoints;i++){
371 Int_t index = (Int_t)fTree->GetV1()[i];
372 tpoints->GetEntryWithIndex(index,index);
373 if (points->GetNPoints()<2) continue;
374 Int_t goodpoints=0;
375 for (Int_t i=0;i<points->GetNPoints(); i++){
376 xyz[goodpoints*3] = points->GetX()[i];
377 xyz[goodpoints*3+1] = points->GetY()[i];
378 xyz[goodpoints*3+2] = points->GetZ()[i];
379 if ( points->GetX()[i]*points->GetX()[i]+points->GetY()[i]*points->GetY()[i]>rmin) goodpoints++;
380 }
381 TPolyMarker3D * marker = new TPolyMarker3D(goodpoints,xyz);
382 marker->SetMarkerColor(color);
383 marker->SetMarkerStyle(1);
384 fPoints->AddLast(marker);
385 }
386}
387