]>
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) : | |
03631860 | 84 | AliCFContainer(c), |
85 | fNBins(c.fNBins), | |
86 | fNVars(c.fNVars), | |
87 | fNSteps(c.fNSteps), | |
88 | fValues(new TArrayF*[c.fNSteps]), | |
89 | fSumw2(new TArrayF*[c.fNSteps]), | |
eed401dc | 90 | axisCache(0) |
85bfac17 | 91 | { |
92 | // | |
93 | // AliTHn copy constructor | |
94 | // | |
95 | ||
03631860 | 96 | memset(fValues,0,fNSteps*sizeof(TArrayF*)); |
97 | memset(fSumw2,0,fNSteps*sizeof(TArrayF*)); | |
98 | ||
99 | for (Int_t i=0; i<fNSteps; i++) { | |
100 | if (c.fValues[i]) fValues[i] = new TArrayF(*(c.fValues[i])); | |
101 | if (c.fSumw2[i]) fSumw2[i] = new TArrayF(*(c.fSumw2[i])); | |
102 | } | |
103 | ||
85bfac17 | 104 | } |
105 | ||
106 | AliTHn::~AliTHn() | |
107 | { | |
108 | // Destructor | |
109 | ||
110 | DeleteContainers(); | |
111 | ||
03631860 | 112 | delete[] fValues; |
113 | delete[] fSumw2; | |
114 | delete[] axisCache; | |
85bfac17 | 115 | |
85bfac17 | 116 | } |
117 | ||
118 | void AliTHn::DeleteContainers() | |
119 | { | |
120 | // delete data containers | |
121 | ||
122 | for (Int_t i=0; i<fNSteps; i++) | |
123 | { | |
124 | if (fValues && fValues[i]) | |
125 | { | |
126 | delete fValues[i]; | |
127 | fValues[i] = 0; | |
128 | } | |
129 | ||
130 | if (fSumw2 && fSumw2[i]) | |
131 | { | |
132 | delete fSumw2[i]; | |
133 | fSumw2[i] = 0; | |
134 | } | |
135 | } | |
136 | } | |
137 | ||
138 | //____________________________________________________________________ | |
139 | AliTHn &AliTHn::operator=(const AliTHn &c) | |
140 | { | |
141 | // assigment operator | |
142 | ||
03631860 | 143 | if (this != &c) { |
144 | AliCFContainer::operator=(c); | |
145 | fNBins=c.fNBins; | |
146 | fNVars=c.fNVars; | |
03631860 | 147 | if(fNSteps) { |
148 | for(Int_t i=0; i< fNSteps; ++i) { | |
149 | delete fValues[i]; | |
150 | delete fSumw2[i]; | |
151 | } | |
152 | delete [] fValues; | |
153 | delete [] fSumw2; | |
3e1096ac | 154 | } |
155 | fNSteps=c.fNSteps; | |
156 | if(fNSteps) { | |
03631860 | 157 | fValues=new TArrayF*[fNSteps]; |
158 | fSumw2=new TArrayF*[fNSteps]; | |
159 | memset(fValues,0,fNSteps*sizeof(TArrayF*)); | |
160 | memset(fSumw2,0,fNSteps*sizeof(TArrayF*)); | |
161 | ||
162 | for (Int_t i=0; i<fNSteps; i++) { | |
163 | if (c.fValues[i]) fValues[i] = new TArrayF(*(c.fValues[i])); | |
164 | if (c.fSumw2[i]) fSumw2[i] = new TArrayF(*(c.fSumw2[i])); | |
165 | } | |
3e1096ac | 166 | } else { |
167 | fValues = 0; | |
168 | fSumw2 = 0; | |
03631860 | 169 | } |
170 | delete [] axisCache; | |
171 | axisCache = 0; | |
172 | } | |
85bfac17 | 173 | return *this; |
174 | } | |
175 | ||
176 | //____________________________________________________________________ | |
177 | void AliTHn::Copy(TObject& c) const | |
178 | { | |
179 | // copy function | |
180 | ||
181 | AliTHn& target = (AliTHn &) c; | |
182 | ||
183 | AliCFContainer::Copy(target); | |
184 | ||
185 | target.fNSteps = fNSteps; | |
186 | target.fNBins = fNBins; | |
187 | target.fNVars = fNVars; | |
188 | ||
189 | target.Init(); | |
190 | ||
191 | for (Int_t i=0; i<fNSteps; i++) | |
192 | { | |
193 | if (fValues[i]) | |
194 | target.fValues[i] = new TArrayF(*(fValues[i])); | |
195 | else | |
196 | target.fValues[i] = 0; | |
197 | ||
198 | if (fSumw2[i]) | |
199 | target.fSumw2[i] = new TArrayF(*(fSumw2[i])); | |
200 | else | |
201 | target.fSumw2[i] = 0; | |
202 | } | |
203 | } | |
204 | ||
205 | //____________________________________________________________________ | |
206 | Long64_t AliTHn::Merge(TCollection* list) | |
207 | { | |
208 | // Merge a list of AliTHn objects with this (needed for | |
209 | // PROOF). | |
210 | // Returns the number of merged objects (including this). | |
211 | ||
212 | if (!list) | |
213 | return 0; | |
214 | ||
215 | if (list->IsEmpty()) | |
216 | return 1; | |
217 | ||
218 | AliCFContainer::Merge(list); | |
219 | ||
220 | TIterator* iter = list->MakeIterator(); | |
221 | TObject* obj; | |
222 | ||
223 | Int_t count = 0; | |
224 | while ((obj = iter->Next())) { | |
225 | ||
226 | AliTHn* entry = dynamic_cast<AliTHn*> (obj); | |
227 | if (entry == 0) | |
228 | continue; | |
229 | ||
230 | for (Int_t i=0; i<fNSteps; i++) | |
231 | { | |
232 | if (entry->fValues[i]) | |
233 | { | |
234 | if (!fValues[i]) | |
235 | fValues[i] = new TArrayF(fNBins); | |
236 | ||
237 | for (Long64_t l = 0; l<fNBins; l++) | |
238 | fValues[i]->GetArray()[l] += entry->fValues[i]->GetArray()[l]; | |
239 | } | |
240 | ||
241 | if (entry->fSumw2[i]) | |
242 | { | |
243 | if (!fSumw2[i]) | |
244 | fSumw2[i] = new TArrayF(fNBins); | |
245 | ||
246 | for (Long64_t l = 0; l<fNBins; l++) | |
247 | fSumw2[i]->GetArray()[l] += entry->fSumw2[i]->GetArray()[l]; | |
248 | } | |
249 | } | |
250 | ||
251 | count++; | |
252 | } | |
253 | ||
254 | return count+1; | |
255 | } | |
256 | ||
257 | void AliTHn::Fill(const Double_t *var, Int_t istep, Double_t weight) | |
258 | { | |
259 | // fills an entry | |
260 | ||
eed401dc | 261 | // fill axis cache |
262 | if (!axisCache) | |
263 | { | |
264 | axisCache = new TAxis*[fNVars]; | |
265 | for (Int_t i=0; i<fNVars; i++) | |
266 | axisCache[i] = GetAxis(i, 0); | |
267 | } | |
268 | ||
85bfac17 | 269 | // calculate global bin index |
270 | Long64_t bin = 0; | |
271 | for (Int_t i=0; i<fNVars; i++) | |
272 | { | |
eed401dc | 273 | bin *= axisCache[i]->GetNbins(); |
85bfac17 | 274 | |
eed401dc | 275 | Int_t tmpBin = axisCache[i]->FindBin(var[i]); |
85bfac17 | 276 | // Printf("%d", tmpBin); |
277 | // under/overflow not supported | |
eed401dc | 278 | if (tmpBin < 1 || tmpBin > axisCache[i]->GetNbins()) |
85bfac17 | 279 | return; |
280 | ||
281 | // bins start from 0 here | |
282 | bin += tmpBin - 1; | |
283 | // Printf("%lld", bin); | |
284 | } | |
285 | ||
286 | if (!fValues[istep]) | |
287 | { | |
288 | fValues[istep] = new TArrayF(fNBins); | |
289 | AliInfo(Form("Created values container for step %d", istep)); | |
290 | } | |
291 | ||
292 | if (weight != 1) | |
293 | { | |
294 | // initialize with already filled entries (which have been filled with weight == 1), in this case fSumw2 := fValues | |
295 | if (!fSumw2[istep]) | |
296 | { | |
297 | fSumw2[istep] = new TArrayF(*fValues[istep]); | |
298 | AliInfo(Form("Created sumw2 container for step %d", istep)); | |
299 | } | |
300 | } | |
301 | ||
302 | fValues[istep]->GetArray()[bin] += weight; | |
303 | if (fSumw2[istep]) | |
304 | fSumw2[istep]->GetArray()[bin] += weight * weight; | |
305 | ||
306 | // Printf("%f", fValues[istep][bin]); | |
307 | ||
308 | // debug | |
309 | // AliCFContainer::Fill(var, istep, weight); | |
310 | } | |
311 | ||
312 | Long64_t AliTHn::GetGlobalBinIndex(const Int_t* binIdx) | |
313 | { | |
314 | // calculates global bin index | |
315 | // binIdx contains TAxis bin indexes | |
316 | // here bin count starts at 0 because we do not have over/underflow bins | |
317 | ||
318 | Long64_t bin = 0; | |
319 | for (Int_t i=0; i<fNVars; i++) | |
320 | { | |
321 | bin *= GetAxis(i, 0)->GetNbins(); | |
322 | bin += binIdx[i] - 1; | |
323 | } | |
324 | ||
325 | return bin; | |
326 | } | |
327 | ||
1bba939a | 328 | void AliTHn::FillContainer(AliCFContainer* cont) |
85bfac17 | 329 | { |
1bba939a | 330 | // fills the information stored in the buffer in this class into the container <cont> |
85bfac17 | 331 | |
332 | for (Int_t i=0; i<fNSteps; i++) | |
333 | { | |
334 | if (!fValues[i]) | |
335 | continue; | |
336 | ||
337 | Float_t* source = fValues[i]->GetArray(); | |
338 | // if fSumw2 is not stored, the sqrt of the number of bin entries in source is filled below; otherwise we use fSumw2 | |
339 | Float_t* sourceSumw2 = source; | |
340 | if (fSumw2[i]) | |
341 | sourceSumw2 = fSumw2[i]->GetArray(); | |
342 | ||
1bba939a | 343 | THnSparse* target = cont->GetGrid(i)->GetGrid(); |
85bfac17 | 344 | |
345 | Int_t* binIdx = new Int_t[fNVars]; | |
c32a0ca9 | 346 | Int_t* nBins = new Int_t[fNVars]; |
85bfac17 | 347 | for (Int_t j=0; j<fNVars; j++) |
c32a0ca9 | 348 | { |
85bfac17 | 349 | binIdx[j] = 1; |
c32a0ca9 | 350 | nBins[j] = target->GetAxis(j)->GetNbins(); |
351 | } | |
85bfac17 | 352 | |
353 | Long64_t count = 0; | |
354 | ||
355 | while (1) | |
356 | { | |
357 | // for (Int_t j=0; j<fNVars; j++) | |
358 | // printf("%d ", binIdx[j]); | |
359 | ||
360 | Long64_t globalBin = GetGlobalBinIndex(binIdx); | |
361 | // Printf(" --> %lld", globalBin); | |
362 | ||
363 | if (source[globalBin] != 0) | |
364 | { | |
365 | target->SetBinContent(binIdx, source[globalBin]); | |
366 | target->SetBinError(binIdx, TMath::Sqrt(sourceSumw2[globalBin])); | |
367 | ||
368 | count++; | |
369 | } | |
370 | ||
371 | binIdx[fNVars-1]++; | |
372 | ||
373 | for (Int_t j=fNVars-1; j>0; j--) | |
374 | { | |
c32a0ca9 | 375 | if (binIdx[j] > nBins[j]) |
85bfac17 | 376 | { |
377 | binIdx[j] = 1; | |
378 | binIdx[j-1]++; | |
379 | } | |
380 | } | |
381 | ||
c32a0ca9 | 382 | if (binIdx[0] > nBins[0]) |
85bfac17 | 383 | break; |
384 | } | |
385 | ||
386 | AliInfo(Form("Step %d: copied %lld entries out of %lld bins", i, count, GetGlobalBinIndex(binIdx))); | |
c32a0ca9 | 387 | |
388 | delete[] binIdx; | |
389 | delete[] nBins; | |
1bba939a | 390 | } |
391 | } | |
392 | ||
393 | void AliTHn::FillParent() | |
394 | { | |
395 | // fills the information stored in the buffer in this class into the baseclass containers | |
396 | ||
397 | FillContainer(this); | |
c32a0ca9 | 398 | } |
399 | ||
400 | void AliTHn::ReduceAxis() | |
401 | { | |
402 | // "removes" one axis by summing over the axis and putting the entry to bin 1 | |
403 | // TODO presently only implemented for the last axis | |
404 | ||
405 | Int_t axis = fNVars-1; | |
406 | ||
407 | for (Int_t i=0; i<fNSteps; i++) | |
408 | { | |
409 | if (!fValues[i]) | |
410 | continue; | |
411 | ||
412 | Float_t* source = fValues[i]->GetArray(); | |
413 | Float_t* sourceSumw2 = 0; | |
414 | if (fSumw2[i]) | |
415 | sourceSumw2 = fSumw2[i]->GetArray(); | |
416 | ||
417 | THnSparse* target = GetGrid(i)->GetGrid(); | |
418 | ||
419 | Int_t* binIdx = new Int_t[fNVars]; | |
420 | Int_t* nBins = new Int_t[fNVars]; | |
421 | for (Int_t j=0; j<fNVars; j++) | |
422 | { | |
423 | binIdx[j] = 1; | |
424 | nBins[j] = target->GetAxis(j)->GetNbins(); | |
425 | } | |
426 | ||
427 | Long64_t count = 0; | |
428 | ||
429 | while (1) | |
430 | { | |
431 | // sum over axis <axis> | |
432 | Float_t sumValues = 0; | |
433 | Float_t sumSumw2 = 0; | |
434 | for (Int_t j=1; j<=nBins[axis]; j++) | |
435 | { | |
436 | binIdx[axis] = j; | |
437 | Long64_t globalBin = GetGlobalBinIndex(binIdx); | |
438 | sumValues += source[globalBin]; | |
439 | source[globalBin] = 0; | |
440 | ||
441 | if (sourceSumw2) | |
442 | { | |
443 | sumSumw2 += sourceSumw2[globalBin]; | |
444 | sourceSumw2[globalBin] = 0; | |
445 | } | |
446 | } | |
447 | binIdx[axis] = 1; | |
448 | ||
449 | Long64_t globalBin = GetGlobalBinIndex(binIdx); | |
450 | source[globalBin] = sumValues; | |
451 | if (sourceSumw2) | |
452 | sourceSumw2[globalBin] = sumSumw2; | |
453 | ||
454 | count++; | |
455 | ||
456 | // next bin | |
457 | binIdx[fNVars-2]++; | |
458 | ||
459 | for (Int_t j=fNVars-2; j>0; j--) | |
460 | { | |
461 | if (binIdx[j] > nBins[j]) | |
462 | { | |
463 | binIdx[j] = 1; | |
464 | binIdx[j-1]++; | |
465 | } | |
466 | } | |
467 | ||
468 | if (binIdx[0] > nBins[0]) | |
469 | break; | |
470 | } | |
471 | ||
472 | AliInfo(Form("Step %d: reduced %lld bins to %lld entries", i, GetGlobalBinIndex(binIdx), count)); | |
473 | ||
474 | delete[] binIdx; | |
475 | delete[] nBins; | |
85bfac17 | 476 | } |
477 | } |