use conventional TH numbering (from 1 to N) in bin indeces for GetBinSize(...) and...
[u/mrichter/AliRoot.git] / CORRFW / AliCFGridSparse.cxx
CommitLineData
db6722a5 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// AliCFGridSparse Class //
18// Class to accumulate data on an N-dimensional grid, to be used //
19// as input to get corrections for Reconstruction & Trigger efficiency//
20// Based on root THnSparse //
21// -- Author : S.Arcelli //
22// Still to be done: //
23// --Interpolate among bins in a range //
24//--------------------------------------------------------------------//
25//
26//
27#include "AliCFGridSparse.h"
28#include "THnSparse.h"
29#include "AliLog.h"
30#include "TMath.h"
31#include "TROOT.h"
32#include "TH1D.h"
33#include "TH2D.h"
34#include "TH3D.h"
35#include "TAxis.h"
36
37//____________________________________________________________________
38ClassImp(AliCFGridSparse)
39
40//____________________________________________________________________
41AliCFGridSparse::AliCFGridSparse() :
42 AliCFVGrid(),
db6722a5 43 fData(0x0)
44{
45 // default constructor
46}
47//____________________________________________________________________
48AliCFGridSparse::AliCFGridSparse(const Char_t* name, const Char_t* title) :
49 AliCFVGrid(name,title),
db6722a5 50 fData(0x0)
51{
52 // default constructor
53}
54//____________________________________________________________________
55AliCFGridSparse::AliCFGridSparse(const Char_t* name, const Char_t* title, const Int_t nVarIn, const Int_t * nBinIn, const Double_t *binLimitsIn) :
56 AliCFVGrid(name,title,nVarIn,nBinIn,binLimitsIn),
db6722a5 57 fData(0x0)
58{
59 //
60 // main constructor
61 //
62
63 fData=new THnSparseF(name,title,fNVar,fNVarBins);
64
65 if(binLimitsIn){
66 for(Int_t ivar=0;ivar<fNVar;ivar++){
67 Int_t nbins=fNVarBins[ivar]+1;
68 Double_t *array= new Double_t[nbins];
69 for(Int_t i=0;i<nbins;i++){
70 array[i]=fVarBinLimits[fOffset[ivar]+i];
71 }
72 fData->SetBinEdges(ivar, array);
73 delete array;
74 }
75 }
76}
77//____________________________________________________________________
78AliCFGridSparse::AliCFGridSparse(const AliCFGridSparse& c) :
79 AliCFVGrid(c),
db6722a5 80 fData(c.fData)
81{
82 //
83 // copy constructor
84 //
85 ((AliCFGridSparse &)c).Copy(*this);
86}
87
88//____________________________________________________________________
89AliCFGridSparse::~AliCFGridSparse()
90{
91 //
92 // destructor
93 //
94 if(fData) delete fData;
95}
96
97//____________________________________________________________________
98AliCFGridSparse &AliCFGridSparse::operator=(const AliCFGridSparse &c)
99{
100 //
101 // assigment operator
102 //
103 if (this != &c)
104 ((AliCFGridSparse &) c).Copy(*this);
105 return *this;
106}
107
108//____________________________________________________________________
109void AliCFGridSparse::SetBinLimits(Int_t ivar, Double_t *array)
110{
111 //
112 // setting the arrays containing the bin limits
113 //
114 fData->SetBinEdges(ivar, array);
115 //then fill the appropriate array in ALICFFrame, to be able to use
116 //the getter, in case....
117 Int_t nbins=fNVarBins[ivar]+1;
118 for(Int_t i=0;i<nbins;i++){
119 fVarBinLimits[fOffset[ivar]+i] =array[i];
120 }
121}
122
123//____________________________________________________________________
124void AliCFGridSparse::Fill(Double_t *var, Double_t weight)
125{
126 //
127 // Fill the grid,
128 // given a set of values of the input variable,
129 // with weight (by default w=1)
130 //
131 fData->Fill(var,weight);
132}
133
134//___________________________________________________________________
135TH1D *AliCFGridSparse::Project(Int_t ivar) const
136{
137 //
138 // Make a 1D projection along variable ivar
139 //
140
141 //exclude overflows by default (used in projections)
318f64b1 142 if(fExclOffEntriesInProj){
db6722a5 143 for(Int_t i=0;i<fNVar;i++){
144 fData->GetAxis(i)->SetBit(TAxis::kAxisRange);
145 }
146 }
318f64b1 147
148
db6722a5 149 TH1D *hist=fData->Projection(ivar);
db6722a5 150
db6722a5 151 return hist;
152
153}
154//___________________________________________________________________
155TH2D *AliCFGridSparse::Project(Int_t ivar1, Int_t ivar2) const
156{
157 //
158 // Make a 2D projection along variables ivar1 & ivar2
159 //
160
161 //exclude overflows by default (used in projections)
318f64b1 162 if(fExclOffEntriesInProj){
db6722a5 163 for(Int_t i=0;i<fNVar;i++){
164 fData->GetAxis(i)->SetBit(TAxis::kAxisRange);
165 }
166 }
318f64b1 167 TH2D *hist=fData->Projection(ivar2,ivar1); //notice inverted axis (THnSparse uses TH3 2d-projection convention...)
db6722a5 168 return hist;
169
170}
171//___________________________________________________________________
172TH3D *AliCFGridSparse::Project(Int_t ivar1, Int_t ivar2, Int_t ivar3) const
173{
174 //
175 // Make a 3D projection along variables ivar1 & ivar2 & ivar3
176 //
177 //exclude overflows by default (used in projections)
318f64b1 178 if(fExclOffEntriesInProj){
db6722a5 179 for(Int_t i=0;i<fNVar;i++){
180 fData->GetAxis(i)->SetBit(TAxis::kAxisRange);
181 }
182 }
183
184 TH3D *hist=fData->Projection(ivar1,ivar2,ivar3);
db6722a5 185 return hist;
186
187}
188
189//____________________________________________________________________
190Float_t AliCFGridSparse::GetOverFlows(Int_t ivar) const
191{
192 //
318f64b1 193 // Returns exclusive overflows in variable ivar
db6722a5 194 //
195 Int_t* bin = new Int_t[fNDim];
196 memset(bin, 0, sizeof(Int_t) * fNDim);
197 Float_t ovfl=0.;
198 for (Long64_t i = 0; i < fData->GetNbins(); ++i) {
199 Double_t v = fData->GetBinContent(i, bin);
200 Bool_t add=kTRUE;
201 for(Int_t j=0;j<fNVar;j++){
202 if(ivar==j)continue;
203 if((bin[j]==0) || (bin[j]==fNVarBins[j]+1))add=kFALSE;
204 }
205 if(bin[ivar]==fNVarBins[ivar]+1 && add) ovfl+=v;
206 }
207
208 delete[] bin;
209 return ovfl;
210}
211
212//____________________________________________________________________
213Float_t AliCFGridSparse::GetUnderFlows(Int_t ivar) const
214{
215 //
318f64b1 216 // Returns exclusive overflows in variable ivar
db6722a5 217 //
218 Int_t* bin = new Int_t[fNDim];
219 memset(bin, 0, sizeof(Int_t) * fNDim);
220 Float_t unfl=0.;
221 for (Long64_t i = 0; i < fData->GetNbins(); ++i) {
222 Double_t v = fData->GetBinContent(i, bin);
223 Bool_t add=kTRUE;
224 for(Int_t j=0;j<fNVar;j++){
225 if(ivar==j)continue;
226 if((bin[j]==0) || (bin[j]==fNVarBins[j]+1))add=kFALSE;
227 }
228 if(bin[ivar]==0 && add) unfl+=v;
229 }
230
231 delete[] bin;
232 return unfl;
233}
234
235
236//____________________________________________________________________
237Float_t AliCFGridSparse::GetEntries() const
238{
239 //
240 // total entries (including overflows and underflows)
241 //
242
243 return fData->GetEntries();
244}
245
246//____________________________________________________________________
247Float_t AliCFGridSparse::GetElement(Int_t index) const
248{
249 //
250 // Returns content of grid element index according to the
251 // linear indexing in AliCFFrame
252 //
253 Int_t *bin = new Int_t[fNVar];
254 GetBinIndex(index, bin);
255 for(Int_t i=0;i<fNVar;i++)fIndex[i]=bin[i]+1; //consistency with AliCFGrid
256 Float_t val=GetElement(fIndex);
257 delete [] bin;
258 return val;
259
260}
261//____________________________________________________________________
262Float_t AliCFGridSparse::GetElement(Int_t *bin) const
263{
264 //
265 // Get the content in a bin corresponding to a set of bin indexes
266 //
267 return fData->GetBinContent(bin);
268
269}
270//____________________________________________________________________
271Float_t AliCFGridSparse::GetElement(Double_t *var) const
272{
273 //
274 // Get the content in a bin corresponding to a set of input variables
275 //
276
277 Long_t index=fData->GetBin(var,kFALSE); //this is the THnSparse index (do not allocate new cells if content is empty)
278 if(index<0){
279 return 0.;
280 }else{
281 return fData->GetBinContent(index);
282 }
283}
284
285//____________________________________________________________________
286Float_t AliCFGridSparse::GetElementError(Int_t index) const
287{
288 //
289 // Returns the error on the content of a bin according to a linear
290 // indexing in AliCFFrame
291 //
292
293 Int_t *bin = new Int_t[fNVar];
294 GetBinIndex(index, bin);
295 for(Int_t i=0;i<fNVar;i++)fIndex[i]=bin[i]+1; //consistency with AliCFGrid
296 Float_t val=GetElementError(fIndex);
297 delete [] bin;
298 return val;
299
300}
301//____________________________________________________________________
302Float_t AliCFGridSparse::GetElementError(Int_t *bin) const
303{
304 //
305 // Get the error in a bin corresponding to a set of bin indexes
306 //
307 return fData->GetBinError(bin);
308
309}
310//____________________________________________________________________
311Float_t AliCFGridSparse::GetElementError(Double_t *var) const
312{
313 //
314 // Get the error in a bin corresponding to a set of input variables
315 //
316
317 Long_t index=fData->GetBin(var,kFALSE); //this is the THnSparse index (do not allocate new cells if content is empy)
318 if(index<0){
319 return 0.;
320 }else{
321 return fData->GetBinError(index);
322 }
323}
324
325
326//____________________________________________________________________
327void AliCFGridSparse::SetElement(Int_t index, Float_t val)
328{
329 //
330 // Sets grid element iel to val (linear indexing) in AliCFFrame
331 //
332 Int_t *bin = new Int_t[fNVar];
333 GetBinIndex(index, bin);
334 for(Int_t i=0;i<fNVar;i++)fIndex[i]=bin[i]+1;
335 SetElement(fIndex,val);
336 delete [] bin;
337}
338//____________________________________________________________________
339void AliCFGridSparse::SetElement(Int_t *bin, Float_t val)
340{
341 //
342 // Sets grid element of bin indeces bin to val
343 //
344 fData->SetBinContent(bin,val);
345}
346//____________________________________________________________________
347void AliCFGridSparse::SetElement(Double_t *var, Float_t val)
348{
349 //
350 // Set the content in a bin to value val corresponding to a set of input variables
351 //
352 Long_t index=fData->GetBin(var); //THnSparse index: allocate the cell
353 Int_t *bin = new Int_t[fNVar];
354 fData->GetBinContent(index,bin); //trick to access the array of bins
355 fData->SetBinContent(bin,val);
356 delete [] bin;
357
358}
359
360//____________________________________________________________________
361void AliCFGridSparse::SetElementError(Int_t index, Float_t val)
362{
363 //
364 // Sets grid element iel error to val (linear indexing) in AliCFFrame
365 //
366 Int_t *bin = new Int_t[fNVar];
367 GetBinIndex(index, bin);
368 for(Int_t i=0;i<fNVar;i++)fIndex[i]=bin[i]+1;
369 SetElementError(fIndex,val);
370 delete [] bin;
371}
372//____________________________________________________________________
373void AliCFGridSparse::SetElementError(Int_t *bin, Float_t val)
374{
375 //
376 // Sets grid element error of bin indeces bin to val
377 //
378 fData->SetBinError(bin,val);
379}
380//____________________________________________________________________
381void AliCFGridSparse::SetElementError(Double_t *var, Float_t val)
382{
383 //
384 // Set the error in a bin to value val corresponding to a set of input variables
385 //
386 Long_t index=fData->GetBin(var); //THnSparse index
387 Int_t *bin = new Int_t[fNVar];
388 fData->GetBinContent(index,bin); //trick to access the array of bins
389 fData->SetBinError(bin,val);
390 delete [] bin;
391}
392
393//____________________________________________________________________
394void AliCFGridSparse::SumW2()
395{
396 //
397 //set calculation of the squared sum of the weighted entries
398 //
399 if(!fSumW2){
400 fData->CalculateErrors(kTRUE);
401 }
402
403 fSumW2=kTRUE;
404}
405
406//____________________________________________________________________
407void AliCFGridSparse::Add(AliCFVGrid* aGrid, Double_t c)
408{
409 //
410 //add aGrid to the current one
411 //
412
413
414 if(aGrid->GetNVar()!=fNVar){
415 AliInfo("Different number of variables, cannot add the grids");
416 return;
417 }
418 if(aGrid->GetNDim()!=fNDim){
419 AliInfo("Different number of dimensions, cannot add the grids!");
420 return;
421 }
422
423 if(!fSumW2 && aGrid->GetSumW2())SumW2();
424
425
426 fData->Add(((AliCFGridSparse*)aGrid)->GetGrid(),c);
427
428}
429
430//____________________________________________________________________
431void AliCFGridSparse::Add(AliCFVGrid* aGrid1, AliCFVGrid* aGrid2, Double_t c1,Double_t c2)
432{
433 //
434 //Add aGrid1 and aGrid2 and deposit the result into the current one
435 //
436
437 if(fNVar!=aGrid1->GetNVar()|| fNVar!=aGrid2->GetNVar()){
438 AliInfo("Different number of variables, cannot add the grids");
439 return;
440 }
441 if(fNDim!=aGrid1->GetNDim()|| fNDim!=aGrid2->GetNDim()){
442 AliInfo("Different number of dimensions, cannot add the grids!");
443 return;
444 }
445
446 if(!fSumW2 && (aGrid1->GetSumW2() || aGrid2->GetSumW2()))SumW2();
447
448
449 fData->Reset();
450 fData->Add(((AliCFGridSparse*)aGrid1)->GetGrid(),c1);
451 fData->Add(((AliCFGridSparse*)aGrid2)->GetGrid(),c2);
452
453}
454
455//____________________________________________________________________
456void AliCFGridSparse::Multiply(AliCFVGrid* aGrid, Double_t c)
457{
458 //
459 // Multiply aGrid to the current one
460 //
461
462
463 if(aGrid->GetNVar()!=fNVar){
464 AliInfo("Different number of variables, cannot multiply the grids");
465 return;
466 }
467 if(aGrid->GetNDim()!=fNDim){
468 AliInfo("Different number of dimensions, cannot multiply the grids!");
469 return;
470 }
471
472 if(!fSumW2 && aGrid->GetSumW2())SumW2();
473
474 THnSparse *h= ((AliCFGridSparse*)aGrid)->GetGrid();
475
476 fData->Multiply(h);
477 fData->Scale(c);
478
479}
480
481//____________________________________________________________________
482void AliCFGridSparse::Multiply(AliCFVGrid* aGrid1, AliCFVGrid* aGrid2, Double_t c1,Double_t c2)
483{
484 //
485 //Multiply aGrid1 and aGrid2 and deposit the result into the current one
486 //
487
488 if(fNVar!=aGrid1->GetNVar()|| fNVar!=aGrid2->GetNVar()){
489 AliInfo("Different number of variables, cannot multiply the grids");
490 return;
491 }
492 if(fNDim!=aGrid1->GetNDim()|| fNDim!=aGrid2->GetNDim()){
493 AliInfo("Different number of dimensions, cannot multiply the grids!");
494 return;
495 }
496
497 if(!fSumW2 && (aGrid1->GetSumW2() || aGrid2->GetSumW2()))SumW2();
498
499
500 fData->Reset();
501 THnSparse *h1= ((AliCFGridSparse*)aGrid1)->GetGrid();
502 THnSparse *h2= ((AliCFGridSparse*)aGrid2)->GetGrid();
503 h2->Multiply(h1);
504 h2->Scale(c1*c2);
505 fData->Add(h2);
506}
507
508
509
510//____________________________________________________________________
511void AliCFGridSparse::Divide(AliCFVGrid* aGrid, Double_t c)
512{
513 //
514 // Divide aGrid to the current one
515 //
516
517
518 if(aGrid->GetNVar()!=fNVar){
519 AliInfo("Different number of variables, cannot divide the grids");
520 return;
521 }
522 if(aGrid->GetNDim()!=fNDim){
523 AliInfo("Different number of dimensions, cannot divide the grids!");
524 return;
525 }
526
527 if(!fSumW2 && aGrid->GetSumW2())SumW2();
528
529 THnSparse *h= ((AliCFGridSparse*)aGrid)->GetGrid();
530
531 fData->Divide(h);
532 fData->Scale(c);
533
534}
535
536//____________________________________________________________________
537void AliCFGridSparse::Divide(AliCFVGrid* aGrid1, AliCFVGrid* aGrid2, Double_t c1,Double_t c2, Option_t *option)
538{
539 //
540 //Divide aGrid1 and aGrid2 and deposit the result into the current one
541 //bynomial errors are supported
542 //
543
544 if(fNVar!=aGrid1->GetNVar()|| fNVar!=aGrid2->GetNVar()){
545 AliInfo("Different number of variables, cannot divide the grids");
546 return;
547 }
548 if(fNDim!=aGrid1->GetNDim()|| fNDim!=aGrid2->GetNDim()){
549 AliInfo("Different number of dimensions, cannot divide the grids!");
550 return;
551 }
552
553 if(!fSumW2 && (aGrid1->GetSumW2() || aGrid2->GetSumW2()))SumW2();
554
555
556 THnSparse *h1= ((AliCFGridSparse*)aGrid1)->GetGrid();
557 THnSparse *h2= ((AliCFGridSparse*)aGrid2)->GetGrid();
558 fData->Divide(h1,h2,c1,c2,option);
559}
560
561
562//____________________________________________________________________
563void AliCFGridSparse::Copy(TObject& c) const
564{
565 //
566 // copy function
567 //
568 AliCFGridSparse& target = (AliCFGridSparse &) c;
569
570 if(fData)target.fData = fData;
571}
572