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