]>
Commit | Line | Data |
---|---|---|
85bfac17 | 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 | /* $Id: AliTHn.cxx 20164 2007-08-14 15:31:50Z morsch $ */ | |
17 | ||
18 | // | |
19 | // this storage container is optimized for small memory usage | |
20 | // under/over flow bins do not exist | |
21 | // sumw2 structure is float only and only create when the weight != 1 | |
22 | // all histogram functionality (projections, axis, ranges, etc) are taken from THnSparse by propagating | |
23 | // the information up into the THnSparse structure (in the ananalysis *after* data processing and merging) | |
24 | // | |
25 | // the derivation from THnSparse is obviously against many OO rules. correct would be a common baseclass of THnSparse and THn. | |
26 | // | |
27 | // Author: Jan Fiete Grosse-Oetringhaus | |
28 | ||
29 | #include "AliTHn.h" | |
30 | #include "TList.h" | |
31 | #include "TCollection.h" | |
32 | #include "AliLog.h" | |
33 | #include "TArrayF.h" | |
34 | #include "THnSparse.h" | |
35 | #include "TMath.h" | |
36 | ||
37 | ClassImp(AliTHn) | |
38 | ||
39 | AliTHn::AliTHn() : | |
40 | AliCFContainer(), | |
41 | fNBins(0), | |
42 | fNVars(0), | |
43 | fNSteps(0), | |
44 | fValues(0), | |
eed401dc | 45 | fSumw2(0), |
46 | axisCache(0) | |
85bfac17 | 47 | { |
48 | // Constructor | |
49 | } | |
50 | ||
51 | AliTHn::AliTHn(const Char_t* name, const Char_t* title,const Int_t nSelStep, const Int_t nVarIn, const Int_t* nBinIn) : | |
52 | AliCFContainer(name, title, nSelStep, nVarIn, nBinIn), | |
53 | fNBins(0), | |
54 | fNVars(nVarIn), | |
55 | fNSteps(nSelStep), | |
56 | fValues(0), | |
eed401dc | 57 | fSumw2(0), |
58 | axisCache(0) | |
85bfac17 | 59 | { |
60 | // Constructor | |
61 | ||
62 | fNBins = 1; | |
63 | for (Int_t i=0; i<fNVars; i++) | |
64 | fNBins *= nBinIn[i]; | |
65 | ||
66 | Init(); | |
67 | } | |
68 | ||
69 | void AliTHn::Init() | |
70 | { | |
71 | // initialize | |
72 | ||
73 | fValues = new TArrayF*[fNSteps]; | |
74 | fSumw2 = new TArrayF*[fNSteps]; | |
75 | ||
76 | for (Int_t i=0; i<fNSteps; i++) | |
77 | { | |
78 | fValues[i] = 0; | |
79 | fSumw2[i] = 0; | |
80 | } | |
81 | } | |
82 | ||
83 | AliTHn::AliTHn(const AliTHn &c) : | |
84 | AliCFContainer(), | |
85 | fNBins(0), | |
86 | fNVars(0), | |
87 | fNSteps(0), | |
88 | fValues(0), | |
eed401dc | 89 | fSumw2(0), |
90 | axisCache(0) | |
85bfac17 | 91 | { |
92 | // | |
93 | // AliTHn copy constructor | |
94 | // | |
95 | ||
96 | ((AliTHn &) c).Copy(*this); | |
97 | } | |
98 | ||
99 | AliTHn::~AliTHn() | |
100 | { | |
101 | // Destructor | |
102 | ||
103 | DeleteContainers(); | |
104 | ||
105 | if (fValues) | |
106 | { | |
eed401dc | 107 | delete[] fValues; |
85bfac17 | 108 | fValues = 0; |
109 | } | |
110 | ||
111 | if (fSumw2) | |
112 | { | |
eed401dc | 113 | delete[] fSumw2; |
85bfac17 | 114 | fSumw2 = 0; |
115 | } | |
eed401dc | 116 | |
117 | if (axisCache) | |
118 | { | |
119 | delete[] axisCache; | |
120 | axisCache = 0; | |
121 | } | |
85bfac17 | 122 | } |
123 | ||
124 | void AliTHn::DeleteContainers() | |
125 | { | |
126 | // delete data containers | |
127 | ||
128 | for (Int_t i=0; i<fNSteps; i++) | |
129 | { | |
130 | if (fValues && fValues[i]) | |
131 | { | |
132 | delete fValues[i]; | |
133 | fValues[i] = 0; | |
134 | } | |
135 | ||
136 | if (fSumw2 && fSumw2[i]) | |
137 | { | |
138 | delete fSumw2[i]; | |
139 | fSumw2[i] = 0; | |
140 | } | |
141 | } | |
142 | } | |
143 | ||
144 | //____________________________________________________________________ | |
145 | AliTHn &AliTHn::operator=(const AliTHn &c) | |
146 | { | |
147 | // assigment operator | |
148 | ||
149 | if (this != &c) | |
150 | ((AliTHn &) c).Copy(*this); | |
151 | ||
152 | return *this; | |
153 | } | |
154 | ||
155 | //____________________________________________________________________ | |
156 | void AliTHn::Copy(TObject& c) const | |
157 | { | |
158 | // copy function | |
159 | ||
160 | AliTHn& target = (AliTHn &) c; | |
161 | ||
162 | AliCFContainer::Copy(target); | |
163 | ||
164 | target.fNSteps = fNSteps; | |
165 | target.fNBins = fNBins; | |
166 | target.fNVars = fNVars; | |
167 | ||
168 | target.Init(); | |
169 | ||
170 | for (Int_t i=0; i<fNSteps; i++) | |
171 | { | |
172 | if (fValues[i]) | |
173 | target.fValues[i] = new TArrayF(*(fValues[i])); | |
174 | else | |
175 | target.fValues[i] = 0; | |
176 | ||
177 | if (fSumw2[i]) | |
178 | target.fSumw2[i] = new TArrayF(*(fSumw2[i])); | |
179 | else | |
180 | target.fSumw2[i] = 0; | |
181 | } | |
182 | } | |
183 | ||
184 | //____________________________________________________________________ | |
185 | Long64_t AliTHn::Merge(TCollection* list) | |
186 | { | |
187 | // Merge a list of AliTHn objects with this (needed for | |
188 | // PROOF). | |
189 | // Returns the number of merged objects (including this). | |
190 | ||
191 | if (!list) | |
192 | return 0; | |
193 | ||
194 | if (list->IsEmpty()) | |
195 | return 1; | |
196 | ||
197 | AliCFContainer::Merge(list); | |
198 | ||
199 | TIterator* iter = list->MakeIterator(); | |
200 | TObject* obj; | |
201 | ||
202 | Int_t count = 0; | |
203 | while ((obj = iter->Next())) { | |
204 | ||
205 | AliTHn* entry = dynamic_cast<AliTHn*> (obj); | |
206 | if (entry == 0) | |
207 | continue; | |
208 | ||
209 | for (Int_t i=0; i<fNSteps; i++) | |
210 | { | |
211 | if (entry->fValues[i]) | |
212 | { | |
213 | if (!fValues[i]) | |
214 | fValues[i] = new TArrayF(fNBins); | |
215 | ||
216 | for (Long64_t l = 0; l<fNBins; l++) | |
217 | fValues[i]->GetArray()[l] += entry->fValues[i]->GetArray()[l]; | |
218 | } | |
219 | ||
220 | if (entry->fSumw2[i]) | |
221 | { | |
222 | if (!fSumw2[i]) | |
223 | fSumw2[i] = new TArrayF(fNBins); | |
224 | ||
225 | for (Long64_t l = 0; l<fNBins; l++) | |
226 | fSumw2[i]->GetArray()[l] += entry->fSumw2[i]->GetArray()[l]; | |
227 | } | |
228 | } | |
229 | ||
230 | count++; | |
231 | } | |
232 | ||
233 | return count+1; | |
234 | } | |
235 | ||
236 | void AliTHn::Fill(const Double_t *var, Int_t istep, Double_t weight) | |
237 | { | |
238 | // fills an entry | |
239 | ||
eed401dc | 240 | // fill axis cache |
241 | if (!axisCache) | |
242 | { | |
243 | axisCache = new TAxis*[fNVars]; | |
244 | for (Int_t i=0; i<fNVars; i++) | |
245 | axisCache[i] = GetAxis(i, 0); | |
246 | } | |
247 | ||
85bfac17 | 248 | // calculate global bin index |
249 | Long64_t bin = 0; | |
250 | for (Int_t i=0; i<fNVars; i++) | |
251 | { | |
eed401dc | 252 | bin *= axisCache[i]->GetNbins(); |
85bfac17 | 253 | |
eed401dc | 254 | Int_t tmpBin = axisCache[i]->FindBin(var[i]); |
85bfac17 | 255 | // Printf("%d", tmpBin); |
256 | // under/overflow not supported | |
eed401dc | 257 | if (tmpBin < 1 || tmpBin > axisCache[i]->GetNbins()) |
85bfac17 | 258 | return; |
259 | ||
260 | // bins start from 0 here | |
261 | bin += tmpBin - 1; | |
262 | // Printf("%lld", bin); | |
263 | } | |
264 | ||
265 | if (!fValues[istep]) | |
266 | { | |
267 | fValues[istep] = new TArrayF(fNBins); | |
268 | AliInfo(Form("Created values container for step %d", istep)); | |
269 | } | |
270 | ||
271 | if (weight != 1) | |
272 | { | |
273 | // initialize with already filled entries (which have been filled with weight == 1), in this case fSumw2 := fValues | |
274 | if (!fSumw2[istep]) | |
275 | { | |
276 | fSumw2[istep] = new TArrayF(*fValues[istep]); | |
277 | AliInfo(Form("Created sumw2 container for step %d", istep)); | |
278 | } | |
279 | } | |
280 | ||
281 | fValues[istep]->GetArray()[bin] += weight; | |
282 | if (fSumw2[istep]) | |
283 | fSumw2[istep]->GetArray()[bin] += weight * weight; | |
284 | ||
285 | // Printf("%f", fValues[istep][bin]); | |
286 | ||
287 | // debug | |
288 | // AliCFContainer::Fill(var, istep, weight); | |
289 | } | |
290 | ||
291 | Long64_t AliTHn::GetGlobalBinIndex(const Int_t* binIdx) | |
292 | { | |
293 | // calculates global bin index | |
294 | // binIdx contains TAxis bin indexes | |
295 | // here bin count starts at 0 because we do not have over/underflow bins | |
296 | ||
297 | Long64_t bin = 0; | |
298 | for (Int_t i=0; i<fNVars; i++) | |
299 | { | |
300 | bin *= GetAxis(i, 0)->GetNbins(); | |
301 | bin += binIdx[i] - 1; | |
302 | } | |
303 | ||
304 | return bin; | |
305 | } | |
306 | ||
307 | void AliTHn::FillParent() | |
308 | { | |
309 | // fills the information stored in the buffer in this class into the baseclass containers | |
310 | ||
311 | for (Int_t i=0; i<fNSteps; i++) | |
312 | { | |
313 | if (!fValues[i]) | |
314 | continue; | |
315 | ||
316 | Float_t* source = fValues[i]->GetArray(); | |
317 | // if fSumw2 is not stored, the sqrt of the number of bin entries in source is filled below; otherwise we use fSumw2 | |
318 | Float_t* sourceSumw2 = source; | |
319 | if (fSumw2[i]) | |
320 | sourceSumw2 = fSumw2[i]->GetArray(); | |
321 | ||
322 | THnSparse* target = GetGrid(i)->GetGrid(); | |
323 | ||
324 | Int_t* binIdx = new Int_t[fNVars]; | |
c32a0ca9 | 325 | Int_t* nBins = new Int_t[fNVars]; |
85bfac17 | 326 | for (Int_t j=0; j<fNVars; j++) |
c32a0ca9 | 327 | { |
85bfac17 | 328 | binIdx[j] = 1; |
c32a0ca9 | 329 | nBins[j] = target->GetAxis(j)->GetNbins(); |
330 | } | |
85bfac17 | 331 | |
332 | Long64_t count = 0; | |
333 | ||
334 | while (1) | |
335 | { | |
336 | // for (Int_t j=0; j<fNVars; j++) | |
337 | // printf("%d ", binIdx[j]); | |
338 | ||
339 | Long64_t globalBin = GetGlobalBinIndex(binIdx); | |
340 | // Printf(" --> %lld", globalBin); | |
341 | ||
342 | if (source[globalBin] != 0) | |
343 | { | |
344 | target->SetBinContent(binIdx, source[globalBin]); | |
345 | target->SetBinError(binIdx, TMath::Sqrt(sourceSumw2[globalBin])); | |
346 | ||
347 | count++; | |
348 | } | |
349 | ||
350 | binIdx[fNVars-1]++; | |
351 | ||
352 | for (Int_t j=fNVars-1; j>0; j--) | |
353 | { | |
c32a0ca9 | 354 | if (binIdx[j] > nBins[j]) |
85bfac17 | 355 | { |
356 | binIdx[j] = 1; | |
357 | binIdx[j-1]++; | |
358 | } | |
359 | } | |
360 | ||
c32a0ca9 | 361 | if (binIdx[0] > nBins[0]) |
85bfac17 | 362 | break; |
363 | } | |
364 | ||
365 | AliInfo(Form("Step %d: copied %lld entries out of %lld bins", i, count, GetGlobalBinIndex(binIdx))); | |
c32a0ca9 | 366 | |
367 | delete[] binIdx; | |
368 | delete[] nBins; | |
369 | } | |
370 | } | |
371 | ||
372 | void AliTHn::ReduceAxis() | |
373 | { | |
374 | // "removes" one axis by summing over the axis and putting the entry to bin 1 | |
375 | // TODO presently only implemented for the last axis | |
376 | ||
377 | Int_t axis = fNVars-1; | |
378 | ||
379 | for (Int_t i=0; i<fNSteps; i++) | |
380 | { | |
381 | if (!fValues[i]) | |
382 | continue; | |
383 | ||
384 | Float_t* source = fValues[i]->GetArray(); | |
385 | Float_t* sourceSumw2 = 0; | |
386 | if (fSumw2[i]) | |
387 | sourceSumw2 = fSumw2[i]->GetArray(); | |
388 | ||
389 | THnSparse* target = GetGrid(i)->GetGrid(); | |
390 | ||
391 | Int_t* binIdx = new Int_t[fNVars]; | |
392 | Int_t* nBins = new Int_t[fNVars]; | |
393 | for (Int_t j=0; j<fNVars; j++) | |
394 | { | |
395 | binIdx[j] = 1; | |
396 | nBins[j] = target->GetAxis(j)->GetNbins(); | |
397 | } | |
398 | ||
399 | Long64_t count = 0; | |
400 | ||
401 | while (1) | |
402 | { | |
403 | // sum over axis <axis> | |
404 | Float_t sumValues = 0; | |
405 | Float_t sumSumw2 = 0; | |
406 | for (Int_t j=1; j<=nBins[axis]; j++) | |
407 | { | |
408 | binIdx[axis] = j; | |
409 | Long64_t globalBin = GetGlobalBinIndex(binIdx); | |
410 | sumValues += source[globalBin]; | |
411 | source[globalBin] = 0; | |
412 | ||
413 | if (sourceSumw2) | |
414 | { | |
415 | sumSumw2 += sourceSumw2[globalBin]; | |
416 | sourceSumw2[globalBin] = 0; | |
417 | } | |
418 | } | |
419 | binIdx[axis] = 1; | |
420 | ||
421 | Long64_t globalBin = GetGlobalBinIndex(binIdx); | |
422 | source[globalBin] = sumValues; | |
423 | if (sourceSumw2) | |
424 | sourceSumw2[globalBin] = sumSumw2; | |
425 | ||
426 | count++; | |
427 | ||
428 | // next bin | |
429 | binIdx[fNVars-2]++; | |
430 | ||
431 | for (Int_t j=fNVars-2; j>0; j--) | |
432 | { | |
433 | if (binIdx[j] > nBins[j]) | |
434 | { | |
435 | binIdx[j] = 1; | |
436 | binIdx[j-1]++; | |
437 | } | |
438 | } | |
439 | ||
440 | if (binIdx[0] > nBins[0]) | |
441 | break; | |
442 | } | |
443 | ||
444 | AliInfo(Form("Step %d: reduced %lld bins to %lld entries", i, GetGlobalBinIndex(binIdx), count)); | |
445 | ||
446 | delete[] binIdx; | |
447 | delete[] nBins; | |
85bfac17 | 448 | } |
449 | } |