Memory leaks corrected (Ch. Klein-Boesing)
[u/mrichter/AliRoot.git] / PWG0 / AliCorrectionMatrix3D.cxx
CommitLineData
bf21645b 1/* $Id$ */
2
3// ------------------------------------------------------
4//
5// Class to handle 3d-corrections.
6//
7// ------------------------------------------------------
8//
9
10#include <TH3F.h>
06e4b91b 11#include <TH2F.h>
f10a1859 12#include <TH1F.h>
06e4b91b 13#include <TString.h>
bf21645b 14
15#include <AliLog.h>
16
17#include "AliCorrectionMatrix3D.h"
06e4b91b 18#include "AliCorrectionMatrix2D.h"
847489f7 19#include "AliPWG0Helper.h"
bf21645b 20
21//____________________________________________________________________
22ClassImp(AliCorrectionMatrix3D)
23
24//____________________________________________________________________
25AliCorrectionMatrix3D::AliCorrectionMatrix3D() :
26 AliCorrectionMatrix()
27{
28 // default constructor
29}
30
31//____________________________________________________________________
32AliCorrectionMatrix3D::AliCorrectionMatrix3D(const AliCorrectionMatrix3D& c)
33 : AliCorrectionMatrix(c)
34{
35 // copy constructor
36 ((AliCorrectionMatrix3D &)c).Copy(*this);
37}
38
61385583 39//____________________________________________________________________
40AliCorrectionMatrix3D &AliCorrectionMatrix3D::operator=(const AliCorrectionMatrix3D &c)
41{
42 // assigment operator
43
44 if (this != &c)
45 ((AliCorrectionMatrix3D &) c).Copy(*this);
46
47 return *this;
48}
49
bf21645b 50//____________________________________________________________________
51AliCorrectionMatrix3D::AliCorrectionMatrix3D(const Char_t* name, const Char_t* title,
52 Int_t nBinX, Float_t Xmin, Float_t Xmax,
53 Int_t nBinY, Float_t Ymin, Float_t Ymax,
54 Int_t nBinZ, Float_t Zmin, Float_t Zmax)
55 : AliCorrectionMatrix(name, title)
56{
57 //
58 // constructor
59 //
60
083a636e 61 Float_t* binLimitsX = new Float_t[nBinX+1];
62 for (Int_t i=0; i<=nBinX; ++i)
63 binLimitsX[i] = Xmin + (Xmax - Xmin) / nBinX * i;
29771dc8 64
083a636e 65 Float_t* binLimitsY = new Float_t[nBinY+1];
66 for (Int_t i=0; i<=nBinY; ++i)
67 binLimitsY[i] = Ymin + (Ymax - Ymin) / nBinY * i;
29771dc8 68
083a636e 69 Float_t* binLimitsZ = new Float_t[nBinZ+1];
70 for (Int_t i=0; i<=nBinZ; ++i)
71 binLimitsZ[i] = Zmin + (Zmax - Zmin) / nBinZ * i;
bf21645b 72
083a636e 73 CreateHists(nBinX, binLimitsX, nBinY, binLimitsY, nBinZ, binLimitsZ);
74
75 delete[] binLimitsX;
76 delete[] binLimitsY;
77 delete[] binLimitsZ;
bf21645b 78}
79
1afae8ff 80AliCorrectionMatrix3D::AliCorrectionMatrix3D(const Char_t* name, const Char_t* title,
81 Int_t nBinX, Float_t Xmin, Float_t Xmax,
82 Int_t nBinY, Float_t Ymin, Float_t Ymax,
83 Int_t nBinZ, const Float_t* zbins)
84 : AliCorrectionMatrix(name, title)
85{
86 // constructor with variable bin sizes
87
88 Float_t* binLimitsX = new Float_t[nBinX+1];
89 for (Int_t i=0; i<=nBinX; ++i)
90 binLimitsX[i] = Xmin + (Xmax - Xmin) / nBinX * i;
91
92 Float_t* binLimitsY = new Float_t[nBinY+1];
93 for (Int_t i=0; i<=nBinY; ++i)
94 binLimitsY[i] = Ymin + (Ymax - Ymin) / nBinY * i;
95
083a636e 96 CreateHists(nBinX, binLimitsX, nBinY, binLimitsY, nBinZ, zbins);
29771dc8 97
083a636e 98 delete[] binLimitsX;
99 delete[] binLimitsY;
100}
29771dc8 101
06e4b91b 102AliCorrectionMatrix3D::AliCorrectionMatrix3D(const Char_t* name, const Char_t* title, TH3F* hBinning)
083a636e 103 : AliCorrectionMatrix(name, title)
104{
06e4b91b 105 // constructor with variable bin sizes (uses binning of hBinning)
083a636e 106
06e4b91b 107 // do not add this hists to the directory
108 Bool_t oldStatus = TH1::AddDirectoryStatus();
109 TH1::AddDirectory(kFALSE);
083a636e 110
06e4b91b 111 fhMeas = (TH3F*)hBinning->Clone("measured");
112 fhGene = (TH3F*)hBinning->Clone("generated");
113 fhCorr = (TH3F*)hBinning->Clone("correction");
1afae8ff 114
06e4b91b 115 fhMeas->SetTitle(Form("%s measured",title) );
116 fhGene->SetTitle(Form("%s generated",title ));
117 fhCorr->SetTitle(Form("%s correction",title));
118
119 fhMeas->Reset();
120 fhGene->Reset();
121 fhCorr->Reset();
122
123 TH1::AddDirectory(oldStatus);
124
125 fhMeas->Sumw2();
126 fhGene->Sumw2();
127 fhCorr->Sumw2();
083a636e 128}
129
130//____________________________________________________________________
131void AliCorrectionMatrix3D::CreateHists(Int_t nBinX, const Float_t* binLimitsX,
132 Int_t nBinY, const Float_t* binLimitsY,
133 Int_t nBinZ, const Float_t* binLimitsZ)
134{
135 // create the histograms
136
137 // do not add this hists to the directory
138 Bool_t oldStatus = TH1::AddDirectoryStatus();
139 TH1::AddDirectory(kFALSE);
140
141 fhMeas = new TH3F("measured", Form("%s measured",GetTitle()), nBinX, binLimitsX, nBinY, binLimitsY, nBinZ, binLimitsZ);
142 fhGene = new TH3F("generated", Form("%s generated",GetTitle()), nBinX, binLimitsX, nBinY, binLimitsY, nBinZ, binLimitsZ);
143 fhCorr = new TH3F("correction", Form("%s correction",GetTitle()), nBinX, binLimitsX, nBinY, binLimitsY, nBinZ, binLimitsZ);
1afae8ff 144
145 fhMeas->Sumw2();
146 fhGene->Sumw2();
147 fhCorr->Sumw2();
083a636e 148
149 TH1::AddDirectory(oldStatus);
1afae8ff 150}
151
083a636e 152
bf21645b 153//____________________________________________________________________
154AliCorrectionMatrix3D::~AliCorrectionMatrix3D()
155{
156 //
157 // destructor
158 //
159
160 // histograms already deleted in base class
161}
162
163//____________________________________________________________________
164TH3F* AliCorrectionMatrix3D::GetGeneratedHistogram()
165{
166 // return generated histogram casted to correct type
167 return dynamic_cast<TH3F*> (fhGene);
168}
169
170//____________________________________________________________________
171TH3F* AliCorrectionMatrix3D::GetMeasuredHistogram()
172{
173 // return measured histogram casted to correct type
174 return dynamic_cast<TH3F*> (fhMeas);
175}
176
177//____________________________________________________________________
178TH3F* AliCorrectionMatrix3D::GetCorrectionHistogram()
179{
180 // return correction histogram casted to correct type
181 return dynamic_cast<TH3F*> (fhCorr);
182}
183
06e4b91b 184//____________________________________________________________________
185AliCorrectionMatrix2D* AliCorrectionMatrix3D::Get2DCorrection(Char_t* opt, Float_t aMin, Float_t aMax)
186{
187 // returns a 2D projection of this correction
188
189 TString option = opt;
190
191 if (aMin<aMax) {
192 if (option.Contains("xy") || option.Contains("yx")) {
193 Int_t bMin = fhMeas->GetZaxis()->FindBin(aMin);
194 Int_t bMax = fhMeas->GetZaxis()->FindBin(aMax);
195 fhMeas->GetZaxis()->SetRange(bMin, bMax);
196 }
197 else if (option.Contains("xz") || option.Contains("zx")) {
198 Int_t bMin = fhMeas->GetYaxis()->FindBin(aMin);
199 Int_t bMax = fhMeas->GetYaxis()->FindBin(aMax);
200 fhMeas->GetYaxis()->SetRange(bMin, bMax);
201 }
202 else if (option.Contains("yz") || option.Contains("zy")) {
203 Int_t bMin = fhMeas->GetXaxis()->FindBin(aMin);
204 Int_t bMax = fhMeas->GetXaxis()->FindBin(aMax);
205 fhMeas->GetXaxis()->SetRange(bMin, bMax);
206 }
207 else {
208 AliDebug(AliLog::kWarning, Form("WARNING: unknown projection option %s", opt));
209 return 0;
210 }
211 }
212 AliCorrectionMatrix2D* corr2D = new AliCorrectionMatrix2D(Form("%s_%s",GetName(),opt),Form("%s projection %s",GetName(),opt),100,0,100,100,0,100);
213
214 TH2F* meas = (TH2F*) ((TH3F*)fhMeas)->Project3D(opt);
215 TH2F* gene = (TH2F*) ((TH3F*)fhGene)->Project3D(opt);
216
217 TH2F* corr = (TH2F*)gene->Clone("corr");
218 corr->Reset();
219
220 corr2D->SetGeneratedHistogram(gene);
221 corr2D->SetMeasuredHistogram(meas);
222 corr2D->SetCorrectionHistogram(corr);
223
224 corr2D->Divide();
225
226 // unzoom
227 fhMeas->GetXaxis()->UnZoom();
228 fhMeas->GetYaxis()->UnZoom();
229 fhMeas->GetZaxis()->UnZoom();
230
231 fhGene->GetXaxis()->UnZoom();
232 fhGene->GetYaxis()->UnZoom();
233 fhGene->GetZaxis()->UnZoom();
234
235 return corr2D;
236}
237
238//____________________________________________________________________
239TH1F* AliCorrectionMatrix3D::Get1DCorrectionHistogram(Char_t* opt, Float_t aMin1, Float_t aMax1, Float_t aMin2, Float_t aMax2)
240{
241 // returns a 1D projection of this correction
242 AliDebug(AliLog::kWarning, Form("WARNING: test"));
243
244 AliCorrectionMatrix2D* corr2D;
245 if (strcmp(opt,"x")==0) {
246 corr2D = Get2DCorrection("xy",aMin1,aMax1);
247 return corr2D->Get1DCorrectionHistogram("x",aMin2,aMax2);
248 }
249 if (strcmp(opt,"y")==0) {
250 corr2D = Get2DCorrection("xy",aMin1,aMax1);
251 return corr2D->Get1DCorrectionHistogram("y",aMin2,aMax2);
252 }
253 if (strcmp(opt,"z")==0) {
254 corr2D = Get2DCorrection("yz",aMin1,aMax1);
255 return corr2D->Get1DCorrectionHistogram("x",aMin2,aMax2);
256 }
257 AliDebug(AliLog::kWarning, Form("WARNING: unknown projection option %s (should be x,y or z)", opt));
258
259 return 0;
260}
261
262
263
bf21645b 264//____________________________________________________________________
265void AliCorrectionMatrix3D::FillMeas(Float_t ax, Float_t ay, Float_t az)
266{
267 // add value to measured histogram
268 GetMeasuredHistogram()->Fill(ax, ay, az);
269}
270
271//____________________________________________________________________
272void AliCorrectionMatrix3D::FillGene(Float_t ax, Float_t ay, Float_t az)
273{
274 // add value to generated histogram
275 GetGeneratedHistogram()->Fill(ax, ay, az);
276}
277
278//____________________________________________________________________
279Float_t AliCorrectionMatrix3D::GetCorrection(Float_t ax, Float_t ay, Float_t az) const
280{
281 // returns a value of the correction map
282 return fhCorr->GetBinContent(fhCorr->FindBin(ax, ay, az));
283}
284
285//____________________________________________________________________
286//void AliCorrectionMatrix3D::RemoveEdges(Float_t cut, Int_t nBinsXedge, Int_t nBinsYedge, Int_t nBinsZedge)
287void AliCorrectionMatrix3D::RemoveEdges(Float_t, Int_t, Int_t, Int_t)
288{
289 // so what do we do here...
290}
291
292//____________________________________________________________________
293void AliCorrectionMatrix3D::SaveHistograms()
294{
295 //
296 // saves the histograms
297 //
298
299 AliCorrectionMatrix::SaveHistograms();
300
92d2d8ad 301 if (GetGeneratedHistogram() && GetMeasuredHistogram())
29771dc8 302 {
303 gDirectory->cd(GetName());
304
305 AliPWG0Helper::CreateDividedProjections(GetGeneratedHistogram(), GetMeasuredHistogram(), 0, kFALSE, kTRUE);
306
307 gDirectory->cd("..");
308 }
bf21645b 309}
f10a1859 310
311//____________________________________________________________________
312Int_t AliCorrectionMatrix3D::CheckEmptyBins(Float_t xmin, Float_t xmax, Float_t ymin, Float_t ymax, Float_t zmin, Float_t zmax, Bool_t quiet)
313{
314 //
315 // counts the number of empty Bins in a given region
316 //
317
318 TH3F* hist = GetGeneratedHistogram();
319 if (!hist)
320 return -1;
321
322 Int_t emptyBins = 0;
323 for (Int_t x=hist->GetXaxis()->FindBin(xmin); x<=hist->GetXaxis()->FindBin(xmax); ++x)
324 for (Int_t y=hist->GetYaxis()->FindBin(ymin); y<=hist->GetYaxis()->FindBin(ymax); ++y)
325 for (Int_t z=hist->GetZaxis()->FindBin(zmin); z<=hist->GetZaxis()->FindBin(zmax); ++z)
326 if (hist->GetBinContent(x, y, z) == 0)
327 {
328 if (!quiet)
329 printf("Empty bin in %s at vtx = %f, eta = %f, pt = %f\n", GetName(), hist->GetXaxis()->GetBinCenter(x), hist->GetYaxis()->GetBinCenter(y), hist->GetZaxis()->GetBinCenter(z));
330 ++emptyBins;
331 }
332
333 return emptyBins;
334}
335
336//____________________________________________________________________
337TH1F* AliCorrectionMatrix3D::PlotBinErrors(Float_t xmin, Float_t xmax, Float_t ymin, Float_t ymax, Float_t zmin, Float_t zmax)
338{
339 //
340 // makes a 1d plots of the relative bin errors
341 //
342
343 TH3F* hist = GetCorrectionHistogram();
344 if (!hist)
345 return 0;
346
347 TH1F* target = new TH1F("relerrors", "relerrors", 100, 0, 10);
348
349 for (Int_t x=hist->GetXaxis()->FindBin(xmin); x<=hist->GetXaxis()->FindBin(xmax); ++x)
350 for (Int_t y=hist->GetYaxis()->FindBin(ymin); y<=hist->GetYaxis()->FindBin(ymax); ++y)
351 for (Int_t z=hist->GetZaxis()->FindBin(zmin); z<=hist->GetZaxis()->FindBin(zmax); ++z)
352 if (hist->GetBinContent(x, y, z) != 0)
353 target->Fill(100 * hist->GetBinError(x, y, z) / hist->GetBinContent(x, y, z));
354
355 return target;
356}
357