Using the event range in the selectors mode as well. To be used in CAF reco in order...
[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() :
fb494025 42 AliCFFrame(),
43 fSumW2(kFALSE),
db6722a5 44 fData(0x0)
45{
46 // default constructor
47}
48//____________________________________________________________________
49AliCFGridSparse::AliCFGridSparse(const Char_t* name, const Char_t* title) :
fb494025 50 AliCFFrame(name,title),
51 fSumW2(kFALSE),
db6722a5 52 fData(0x0)
53{
54 // default constructor
55}
56//____________________________________________________________________
fb494025 57AliCFGridSparse::AliCFGridSparse(const Char_t* name, const Char_t* title, Int_t nVarIn, const Int_t * nBinIn) :
58 AliCFFrame(name,title),
59 fSumW2(kFALSE),
db6722a5 60 fData(0x0)
61{
62 //
63 // main constructor
64 //
65
fb494025 66 fData = new THnSparseF(name,title,nVarIn,nBinIn);
db6722a5 67}
fb494025 68
db6722a5 69//____________________________________________________________________
fb494025 70AliCFGridSparse::~AliCFGridSparse()
db6722a5 71{
72 //
fb494025 73 // destructor
db6722a5 74 //
fb494025 75 if (fData) delete fData;
db6722a5 76}
77
78//____________________________________________________________________
fb494025 79AliCFGridSparse::AliCFGridSparse(const AliCFGridSparse& c) :
80 AliCFFrame(c),
81 fSumW2(kFALSE),
82 fData(0x0)
db6722a5 83{
84 //
fb494025 85 // copy constructor
db6722a5 86 //
fb494025 87 ((AliCFGridSparse &)c).Copy(*this);
db6722a5 88}
89
90//____________________________________________________________________
fb494025 91AliCFGridSparse& AliCFGridSparse::operator=(const AliCFGridSparse &c)
db6722a5 92{
93 //
94 // assigment operator
95 //
fb494025 96 if (this != &c) c.Copy(*this);
db6722a5 97 return *this;
98}
99
100//____________________________________________________________________
fb494025 101void AliCFGridSparse::SetBinLimits(Int_t ivar, const Double_t *array)
db6722a5 102{
103 //
104 // setting the arrays containing the bin limits
105 //
106 fData->SetBinEdges(ivar, array);
db6722a5 107}
108
109//____________________________________________________________________
fb494025 110void AliCFGridSparse::Fill(const Double_t *var, Double_t weight)
db6722a5 111{
112 //
113 // Fill the grid,
114 // given a set of values of the input variable,
115 // with weight (by default w=1)
116 //
117 fData->Fill(var,weight);
118}
119
120//___________________________________________________________________
121TH1D *AliCFGridSparse::Project(Int_t ivar) const
122{
123 //
124 // Make a 1D projection along variable ivar
125 //
126
db6722a5 127 TH1D *hist=fData->Projection(ivar);
fb494025 128 hist->SetXTitle(GetVarTitle(ivar));
91e32eeb 129 hist->SetName(Form("%s_proj-%s",GetName(),GetVarTitle(ivar)));
fb494025 130 hist->SetTitle(Form("%s: projection on %s",GetTitle(),GetVarTitle(ivar)));
db6722a5 131 return hist;
db6722a5 132}
133//___________________________________________________________________
134TH2D *AliCFGridSparse::Project(Int_t ivar1, Int_t ivar2) const
135{
136 //
137 // Make a 2D projection along variables ivar1 & ivar2
138 //
139
318f64b1 140 TH2D *hist=fData->Projection(ivar2,ivar1); //notice inverted axis (THnSparse uses TH3 2d-projection convention...)
fb494025 141 hist->SetXTitle(GetVarTitle(ivar1));
142 hist->SetYTitle(GetVarTitle(ivar2));
91e32eeb 143 hist->SetName(Form("%s_proj-%s,%s",GetName(),GetVarTitle(ivar1),GetVarTitle(ivar2)));
fb494025 144 hist->SetTitle(Form("%s: projection on %s-%s",GetTitle(),GetVarTitle(ivar1),GetVarTitle(ivar2)));
db6722a5 145 return hist;
db6722a5 146}
147//___________________________________________________________________
148TH3D *AliCFGridSparse::Project(Int_t ivar1, Int_t ivar2, Int_t ivar3) const
149{
150 //
151 // Make a 3D projection along variables ivar1 & ivar2 & ivar3
152 //
db6722a5 153
154 TH3D *hist=fData->Projection(ivar1,ivar2,ivar3);
fb494025 155 hist->SetXTitle(GetVarTitle(ivar1));
156 hist->SetYTitle(GetVarTitle(ivar2));
157 hist->SetZTitle(GetVarTitle(ivar3));
91e32eeb 158 hist->SetName(Form("%s_proj-%s,%s,%s",GetName(),GetVarTitle(ivar1),GetVarTitle(ivar2),GetVarTitle(ivar3)));
fb494025 159 hist->SetTitle(Form("%s: projection on %s-%s-%s",GetTitle(),GetVarTitle(ivar1),GetVarTitle(ivar2),GetVarTitle(ivar3)));
db6722a5 160 return hist;
db6722a5 161}
162
9105291d 163//___________________________________________________________________
fb494025 164AliCFGridSparse* AliCFGridSparse::Project(Int_t nVars, const Int_t* vars, const Double_t* varMin, const Double_t* varMax) const
9105291d 165{
38b1447f 166 //
167 // projects the grid on the nVars dimensions defined in vars.
168 // axis ranges can be defined in arrays varMin, varMax
169 //
9105291d 170
171 // binning for new grid
172 Int_t* bins = new Int_t[nVars];
173 for (Int_t iVar=0; iVar<nVars; iVar++) {
fb494025 174 bins[iVar] = GetNBins(vars[iVar]);
9105291d 175 }
176
177 // create new grid sparse
178 AliCFGridSparse* out = new AliCFGridSparse(fName,fTitle,nVars,bins);
179
180 //set the range in the THnSparse to project
181 THnSparse* clone = ((THnSparse*)fData->Clone());
38b1447f 182 if (varMin && varMax) {
fb494025 183 for (Int_t iAxis=0; iAxis<GetNVar(); iAxis++) {
38b1447f 184 clone->GetAxis(iAxis)->SetRangeUser(varMin[iAxis],varMax[iAxis]);
185 }
9105291d 186 }
38b1447f 187 else AliInfo("Keeping same axis ranges");
188
9105291d 189 out->SetGrid(clone->Projection(nVars,vars));
190 return out;
191}
192
fb494025 193
db6722a5 194//____________________________________________________________________
fb494025 195Float_t AliCFGridSparse::GetBinCenter(Int_t ivar, Int_t ibin) const
db6722a5 196{
197 //
fb494025 198 // Returns the center of specified bin for variable axis ivar
199 //
200
201 return (Float_t) fData->GetAxis(ivar)->GetBinCenter(ibin);
db6722a5 202}
203
204//____________________________________________________________________
fb494025 205Float_t AliCFGridSparse::GetBinSize(Int_t ivar, Int_t ibin) const
db6722a5 206{
207 //
fb494025 208 // Returns the size of specified bin for variable axis ivar
209 //
210
211 return (Float_t) fData->GetAxis(ivar)->GetBinUpEdge(ibin) - fData->GetAxis(ivar)->GetBinLowEdge(ibin);
db6722a5 212}
213
db6722a5 214//____________________________________________________________________
215Float_t AliCFGridSparse::GetEntries() const
216{
217 //
218 // total entries (including overflows and underflows)
219 //
220
221 return fData->GetEntries();
222}
223
224//____________________________________________________________________
fb494025 225Float_t AliCFGridSparse::GetElement(Long_t index) const
db6722a5 226{
227 //
fb494025 228 // Returns content of grid element index
db6722a5 229 //
db6722a5 230
fb494025 231 return fData->GetBinContent(index);
db6722a5 232}
233//____________________________________________________________________
fb494025 234Float_t AliCFGridSparse::GetElement(const Int_t *bin) const
db6722a5 235{
236 //
237 // Get the content in a bin corresponding to a set of bin indexes
238 //
239 return fData->GetBinContent(bin);
240
241}
242//____________________________________________________________________
fb494025 243Float_t AliCFGridSparse::GetElement(const Double_t *var) const
db6722a5 244{
245 //
246 // Get the content in a bin corresponding to a set of input variables
247 //
248
fb494025 249 Long_t index = fData->GetBin(var,kFALSE);
250 if (index<0) return 0.;
251 return fData->GetBinContent(index);
db6722a5 252}
253
254//____________________________________________________________________
fb494025 255Float_t AliCFGridSparse::GetElementError(Long_t index) const
db6722a5 256{
257 //
fb494025 258 // Returns the error on the content
db6722a5 259 //
260
fb494025 261 return fData->GetBinContent(index);
db6722a5 262}
263//____________________________________________________________________
fb494025 264Float_t AliCFGridSparse::GetElementError(const Int_t *bin) const
db6722a5 265{
266 //
267 // Get the error in a bin corresponding to a set of bin indexes
268 //
269 return fData->GetBinError(bin);
270
271}
272//____________________________________________________________________
fb494025 273Float_t AliCFGridSparse::GetElementError(const Double_t *var) const
db6722a5 274{
275 //
276 // Get the error in a bin corresponding to a set of input variables
277 //
278
279 Long_t index=fData->GetBin(var,kFALSE); //this is the THnSparse index (do not allocate new cells if content is empy)
fb494025 280 if (index<0) return 0.;
281 return fData->GetBinError(index);
db6722a5 282}
283
db6722a5 284//____________________________________________________________________
fb494025 285void AliCFGridSparse::SetElement(Long_t index, Float_t val)
db6722a5 286{
287 //
fb494025 288 // Sets grid element value
db6722a5 289 //
fb494025 290 Int_t* bin = new Int_t[GetNVar()];
291 fData->GetBinContent(index,bin); //affects the bin coordinates
292 SetElement(bin,val);
293 delete [] bin ;
db6722a5 294}
fb494025 295
db6722a5 296//____________________________________________________________________
fb494025 297void AliCFGridSparse::SetElement(const Int_t *bin, Float_t val)
db6722a5 298{
299 //
300 // Sets grid element of bin indeces bin to val
301 //
302 fData->SetBinContent(bin,val);
303}
304//____________________________________________________________________
fb494025 305void AliCFGridSparse::SetElement(const Double_t *var, Float_t val)
db6722a5 306{
307 //
308 // Set the content in a bin to value val corresponding to a set of input variables
309 //
fb494025 310 Long_t index=fData->GetBin(var,kTRUE); //THnSparse index: allocate the cell
311 Int_t *bin = new Int_t[GetNVar()];
db6722a5 312 fData->GetBinContent(index,bin); //trick to access the array of bins
fb494025 313 SetElement(bin,val);
db6722a5 314 delete [] bin;
db6722a5 315}
316
317//____________________________________________________________________
fb494025 318void AliCFGridSparse::SetElementError(Long_t index, Float_t val)
db6722a5 319{
320 //
321 // Sets grid element iel error to val (linear indexing) in AliCFFrame
322 //
fb494025 323 Int_t *bin = new Int_t[GetNVar()];
324 fData->GetBinContent(index,bin);
325 SetElementError(bin,val);
db6722a5 326 delete [] bin;
327}
fb494025 328
db6722a5 329//____________________________________________________________________
fb494025 330void AliCFGridSparse::SetElementError(const Int_t *bin, Float_t val)
db6722a5 331{
332 //
333 // Sets grid element error of bin indeces bin to val
334 //
335 fData->SetBinError(bin,val);
336}
337//____________________________________________________________________
fb494025 338void AliCFGridSparse::SetElementError(const Double_t *var, Float_t val)
db6722a5 339{
340 //
341 // Set the error in a bin to value val corresponding to a set of input variables
342 //
343 Long_t index=fData->GetBin(var); //THnSparse index
fb494025 344 Int_t *bin = new Int_t[GetNVar()];
db6722a5 345 fData->GetBinContent(index,bin); //trick to access the array of bins
fb494025 346 SetElementError(bin,val);
db6722a5 347 delete [] bin;
348}
349
350//____________________________________________________________________
351void AliCFGridSparse::SumW2()
352{
353 //
354 //set calculation of the squared sum of the weighted entries
355 //
356 if(!fSumW2){
357 fData->CalculateErrors(kTRUE);
358 }
db6722a5 359 fSumW2=kTRUE;
360}
361
362//____________________________________________________________________
fb494025 363void AliCFGridSparse::Add(const AliCFGridSparse* aGrid, Double_t c)
db6722a5 364{
365 //
366 //add aGrid to the current one
367 //
368
fb494025 369 if (aGrid->GetNVar() != GetNVar()){
370 AliError("Different number of variables, cannot add the grids");
db6722a5 371 return;
372 }
db6722a5 373
fb494025 374 if (!fSumW2 && aGrid->GetSumW2()) SumW2();
375 fData->Add(aGrid->GetGrid(),c);
db6722a5 376}
377
378//____________________________________________________________________
fb494025 379void AliCFGridSparse::Add(const AliCFGridSparse* aGrid1, const AliCFGridSparse* aGrid2, Double_t c1,Double_t c2)
db6722a5 380{
381 //
382 //Add aGrid1 and aGrid2 and deposit the result into the current one
383 //
384
fb494025 385 if (GetNVar() != aGrid1->GetNVar() || GetNVar() != aGrid2->GetNVar()) {
db6722a5 386 AliInfo("Different number of variables, cannot add the grids");
387 return;
388 }
db6722a5 389
fb494025 390 if (!fSumW2 && (aGrid1->GetSumW2() || aGrid2->GetSumW2())) SumW2();
db6722a5 391
392 fData->Reset();
fb494025 393 fData->Add(aGrid1->GetGrid(),c1);
394 fData->Add(aGrid2->GetGrid(),c2);
db6722a5 395}
396
397//____________________________________________________________________
fb494025 398void AliCFGridSparse::Multiply(const AliCFGridSparse* aGrid, Double_t c)
db6722a5 399{
400 //
401 // Multiply aGrid to the current one
402 //
403
fb494025 404 if (aGrid->GetNVar() != GetNVar()) {
405 AliError("Different number of variables, cannot multiply the grids");
db6722a5 406 return;
407 }
db6722a5 408
fb494025 409 if(!fSumW2 && aGrid->GetSumW2()) SumW2();
410 THnSparse *h = aGrid->GetGrid();
db6722a5 411 fData->Multiply(h);
412 fData->Scale(c);
db6722a5 413}
414
415//____________________________________________________________________
fb494025 416void AliCFGridSparse::Multiply(const AliCFGridSparse* aGrid1, const AliCFGridSparse* aGrid2, Double_t c1,Double_t c2)
db6722a5 417{
418 //
419 //Multiply aGrid1 and aGrid2 and deposit the result into the current one
420 //
421
fb494025 422 if (GetNVar() != aGrid1->GetNVar() || GetNVar() != aGrid2->GetNVar()) {
423 AliError("Different number of variables, cannot multiply the grids");
db6722a5 424 return;
fb494025 425 }
db6722a5 426
fb494025 427 if(!fSumW2 && (aGrid1->GetSumW2() || aGrid2->GetSumW2())) SumW2();
db6722a5 428
429 fData->Reset();
fb494025 430 THnSparse *h1 = aGrid1->GetGrid();
431 THnSparse *h2 = aGrid2->GetGrid();
db6722a5 432 h2->Multiply(h1);
433 h2->Scale(c1*c2);
434 fData->Add(h2);
435}
436
db6722a5 437//____________________________________________________________________
fb494025 438void AliCFGridSparse::Divide(const AliCFGridSparse* aGrid, Double_t c)
db6722a5 439{
440 //
441 // Divide aGrid to the current one
442 //
443
fb494025 444 if (aGrid->GetNVar() != GetNVar()) {
445 AliError("Different number of variables, cannot divide the grids");
db6722a5 446 return;
447 }
db6722a5 448
fb494025 449 if (!fSumW2 && aGrid->GetSumW2()) SumW2();
db6722a5 450
d7c1cb54 451 THnSparse *h1 = aGrid->GetGrid();
452 THnSparse *h2 = (THnSparse*)fData->Clone();
453 fData->Divide(h2,h1);
db6722a5 454 fData->Scale(c);
db6722a5 455}
456
457//____________________________________________________________________
fb494025 458void AliCFGridSparse::Divide(const AliCFGridSparse* aGrid1, const AliCFGridSparse* aGrid2, Double_t c1,Double_t c2, Option_t *option)
db6722a5 459{
460 //
461 //Divide aGrid1 and aGrid2 and deposit the result into the current one
fb494025 462 //binomial errors are supported
db6722a5 463 //
464
fb494025 465 if (GetNVar() != aGrid1->GetNVar() || GetNVar() != aGrid2->GetNVar()) {
466 AliError("Different number of variables, cannot divide the grids");
db6722a5 467 return;
468 }
db6722a5 469
fb494025 470 if (!fSumW2 && (aGrid1->GetSumW2() || aGrid2->GetSumW2())) SumW2();
db6722a5 471
fb494025 472 THnSparse *h1= aGrid1->GetGrid();
473 THnSparse *h2= aGrid2->GetGrid();
db6722a5 474 fData->Divide(h1,h2,c1,c2,option);
475}
476
477
478//____________________________________________________________________
7411edfd 479void AliCFGridSparse::Rebin(const Int_t* group)
480{
481 //
482 // rebin the grid according to Rebin() as in THnSparse
483 // Please notice that the original number of bins on
484 // a given axis has to be divisible by the rebin group.
485 //
486
fb494025 487 for(Int_t i=0;i<GetNVar();i++){
488 if (group[i]!=1) AliInfo(Form(" merging bins along dimension %i in groups of %i bins", i,group[i]));
7411edfd 489 }
490
491 THnSparse *rebinned =fData->Rebin(group);
492 fData->Reset();
493 fData = rebinned;
fb494025 494}
495//____________________________________________________________________
496void AliCFGridSparse::Scale(Long_t index, const Double_t *fact)
497{
498 //
499 //scale content of a certain cell by (positive) fact (with error)
500 //
7411edfd 501
fb494025 502 if (GetElement(index)==0 || fact[0]==0) return;
7411edfd 503
fb494025 504 Double_t in[2], out[2];
505 in[0]=GetElement(index);
506 in[1]=GetElementError(index);
507 GetScaledValues(fact,in,out);
508 SetElement(index,out[0]);
509 if (fSumW2) SetElementError(index,out[1]);
510}
511//____________________________________________________________________
512void AliCFGridSparse::Scale(const Int_t *bin, const Double_t *fact)
513{
514 //
515 //scale content of a certain cell by (positive) fact (with error)
516 //
517 if(GetElement(bin)==0 || fact[0]==0)return;
7411edfd 518
fb494025 519 Double_t in[2], out[2];
520 in[0]=GetElement(bin);
521 in[1]=GetElementError(bin);
522 GetScaledValues(fact,in,out);
523 SetElement(bin,out[0]);
524 if(fSumW2)SetElementError(bin,out[1]);
525
526}
527//____________________________________________________________________
528void AliCFGridSparse::Scale(const Double_t *var, const Double_t *fact)
529{
530 //
531 //scale content of a certain cell by (positive) fact (with error)
532 //
533 if(GetElement(var)==0 || fact[0]==0)return;
7411edfd 534
fb494025 535 Double_t in[2], out[2];
536 in[0]=GetElement(var);
537 in[1]=GetElementError(var);
538 GetScaledValues(fact,in,out);
539 SetElement(var,out[0]);
540 if(fSumW2)SetElementError(var,out[1]);
541
542}
543//____________________________________________________________________
544void AliCFGridSparse::Scale(const Double_t* fact)
545{
546 //
547 //scale contents of the whole grid by fact
548 //
549
550 for (Long_t iel=0; iel<GetNFilledBins(); iel++) {
551 Scale(iel,fact);
7411edfd 552 }
fb494025 553}
554//____________________________________________________________________
555Long_t AliCFGridSparse::GetEmptyBins() const {
556 //
557 // Get empty bins
558 //
7411edfd 559
fb494025 560 return (GetNBinsTotal() - GetNFilledBins()) ;
561}
7411edfd 562
fb494025 563//_____________________________________________________________________
564// Int_t AliCFGridSparse::GetEmptyBins( Double_t *varMin, Double_t* varMax ) const
565// {
566// //
567// // Get empty bins in a range specified by varMin and varMax
568// //
569
570// Int_t *indexMin = new Int_t[GetNVar()];
571// Int_t *indexMax = new Int_t[GetNVar()];
572
573// //Find out the min and max bins
574
575// for (Int_t i=0;i<GetNVar();i++) {
576// Double_t xmin=varMin[i]; // the min values
577// Double_t xmax=varMax[i]; // the min values
578// Int_t nbins = fData->GetAxis(i)->GetNbins()+1;
579// Double_t *bins=new Double_t[nbins];
580// for(Int_t ibin =0;ibin<nbins;ibin++){
581// bins[ibin] = fVarBinLimits[ibin+fOffset[i]];
582// }
583// indexMin[i] = TMath::BinarySearch(nbins,bins,xmin);
584// indexMax[i] = TMath::BinarySearch(nbins,bins,xmax);
585// delete [] bins;
586// }
587
588// Int_t val=0;
589// for(Int_t i=0;i<fNDim;i++){
590// for (Int_t j=0;j<GetNVar();j++)fIndex[j]=GetBinIndex(j,i);
591// Bool_t isIn=kTRUE;
592// for (Int_t j=0;j<GetNVar();j++){
593// if(!(fIndex[j]>=indexMin[j] && fIndex[j]<=indexMax[j]))isIn=kFALSE;
594// }
595// if(isIn && GetElement(i)<=0)val++;
596// }
597// AliInfo(Form(" the empty bins = %i ",val));
598
599// delete [] indexMin;
600// delete [] indexMax;
601// return val;
602// }
603
604//____________________________________________________________________
605Int_t AliCFGridSparse::CheckStats(Double_t thr) const
606{
607 //
608 // Count the cells below a certain threshold
609 //
610 Int_t ncellsLow=0;
611 for (Int_t i=0; i<GetNBinsTotal(); i++) {
612 if (GetElement(i)<thr) ncellsLow++;
613 }
614 return ncellsLow;
615}
7411edfd 616
fb494025 617//_____________________________________________________________________
618Double_t AliCFGridSparse::GetIntegral() const
619{
620 //
621 // Get full Integral
622 //
623 return fData->ComputeIntegral();
624}
7411edfd 625
fb494025 626//_____________________________________________________________________
627// Double_t AliCFGridSparse::GetIntegral(Int_t *binMin, Int_t* binMax ) const
628// {
629// //
630// // Get Integral in a range of bin indeces (extremes included)
631// //
632
633// Double_t val=0;
634
635// for(Int_t i=0;i<GetNVar();i++){
636// if(binMin[i]<1)binMin[i]=1;
637// if(binMax[i]>fNVarBins[i])binMax[i]=fNVarBins[i];
638// if((binMin[i]>binMax[i])){
639// AliInfo(Form(" Bin indices in variable %i in reverse order, please check!", i));
640// return val;
641// }
642// }
643// val=GetSum(0,binMin,binMax);
644
645// return val;
646// }
647
648//_____________________________________________________________________
649// Double_t AliCFGridSparse::GetIntegral(Double_t *varMin, Double_t* varMax ) const
650// {
651// //
652// // Get Integral in a range (extremes included)
653// //
654
655// Int_t *indexMin=new Int_t[GetNVar()];
656// Int_t *indexMax=new Int_t[GetNVar()];
657
658// //Find out the min and max bins
659
660// for(Int_t i=0;i<GetNVar();i++){
661// Double_t xmin=varMin[i]; // the min values
662// Double_t xmax=varMax[i]; // the min values
663// Int_t nbins=fNVarBins[i]+1;
664// Double_t *bins=new Double_t[nbins];
665// for(Int_t ibin =0;ibin<nbins;ibin++){
666// bins[ibin] = fVarBinLimits[ibin+fOffset[i]];
667// }
668// indexMin[i] = TMath::BinarySearch(nbins,bins,xmin);
669// indexMax[i] = TMath::BinarySearch(nbins,bins,xmax);
670// delete [] bins;
671// }
672
673// //move to the TH/THnSparse convention in N-dim bin numbering
674// for(Int_t i=0;i<GetNVar();i++){
675// indexMin[i]+=1;
676// indexMax[i]+=1;
677// }
678
679// Double_t val=GetIntegral(indexMin,indexMax);
680
681// delete [] indexMin;
682// delete [] indexMax;
683
684// return val;
685// }
686
687// //_____________________________________________________________________
688// Double_t AliCFGridSparse::GetSum(Int_t ivar, Int_t *binMin, Int_t* binMax) const
689// {
690// //
691// // recursively add over nested loops....
692// //
693
694// static Double_t val;
695// if (ivar==0) val=0.;
696// for(Int_t ibin=binMin[ivar]-1 ; ibin<=binMax[ivar]-1 ; ibin++){
697// //-1 is to move from TH/ThnSparse N-dim bin convention to one in AliCFFrame
698// fIndex[ivar]=ibin;
699// if(ivar<GetNVar()-1) {
700// val=GetSum(ivar+1,binMin,binMax);
701// }
702// else {
703// Int_t iel=GetBinIndex(fIndex);
704// val+=GetElement(iel);
705// }
706// }
707
708// return val;
709// }
710
711//____________________________________________________________________
712Long64_t AliCFGridSparse::Merge(TCollection* list)
713{
714 //
715 // Merge a list of AliCFGridSparse with this (needed for PROOF).
716 // Returns the number of merged objects (including this).
717 //
718
719 if (!list)
720 return 0;
7411edfd 721
fb494025 722 if (list->IsEmpty())
723 return 1;
724
725 TIterator* iter = list->MakeIterator();
726 TObject* obj;
727
728 Int_t count = 0;
729 while ((obj = iter->Next())) {
730 AliCFGridSparse* entry = dynamic_cast<AliCFGridSparse*> (obj);
731 if (entry == 0)
732 continue;
733 this->Add(entry);
734 count++;
735 }
736
737 return count+1;
7411edfd 738}
fb494025 739
740//____________________________________________________________________
741void AliCFGridSparse::GetScaledValues(const Double_t *fact, const Double_t *in, Double_t *out) const{
742 //
743 // scale input *in and its error by (positive) fact (with error)
744 // and erite it to *out
745 //
746 out[0]=in[0]*fact[0];
747 out[1]=TMath::Sqrt(in[1]*in[1]/in[0]/in[0]
748 +fact[1]*fact[1]/fact[0]/fact[0])*out[0];
749
750}
751
7411edfd 752//____________________________________________________________________
db6722a5 753void AliCFGridSparse::Copy(TObject& c) const
754{
755 //
756 // copy function
757 //
fb494025 758 AliCFFrame::Copy(c);
db6722a5 759 AliCFGridSparse& target = (AliCFGridSparse &) c;
fb494025 760 target.fSumW2 = fSumW2 ;
761 if (fData) {
762 target.fData = (THnSparse*)fData->Clone();
763 }
db6722a5 764}
765
c8df672e 766//____________________________________________________________________
fb494025 767TH1D* AliCFGridSparse::Slice(Int_t iVar, const Double_t *varMin, const Double_t *varMax) const
c8df672e 768{
769 //
770 // return a slice (1D-projection) on variable iVar while axis ranges are defined with varMin,varMax
771 // arrays varMin and varMax contain the min and max values of each variable.
fb494025 772 // therefore varMin and varMax must have their dimensions equal to GetNVar()
c8df672e 773 //
774
775 THnSparse* clone = (THnSparse*)fData->Clone();
fb494025 776 for (Int_t iAxis=0; iAxis<GetNVar(); iAxis++) {
c8df672e 777 clone->GetAxis(iAxis)->SetRangeUser(varMin[iAxis],varMax[iAxis]);
778 }
779 return clone->Projection(iVar);
780}
781
782//____________________________________________________________________
fb494025 783TH2D* AliCFGridSparse::Slice(Int_t iVar1, Int_t iVar2, const Double_t *varMin, const Double_t *varMax) const
c8df672e 784{
785 //
786 // return a slice (2D-projection) on variables iVar1 and iVar2 while axis ranges are defined with varMin,varMax
787 // arrays varMin and varMax contain the min and max values of each variable.
fb494025 788 // therefore varMin and varMax must have their dimensions equal to GetNVar()
c8df672e 789 //
790
791 THnSparse* clone = (THnSparse*)fData->Clone();
fb494025 792 for (Int_t iAxis=0; iAxis<GetNVar(); iAxis++) {
c8df672e 793 clone->GetAxis(iAxis)->SetRangeUser(varMin[iAxis],varMax[iAxis]);
794 }
795 return clone->Projection(iVar1,iVar2);
796}
797
798//____________________________________________________________________
fb494025 799TH3D* AliCFGridSparse::Slice(Int_t iVar1, Int_t iVar2, Int_t iVar3, const Double_t *varMin, const Double_t *varMax) const
c8df672e 800{
801 //
802 // return a slice (3D-projection) on variables iVar1, iVar2 and iVar3 while axis ranges are defined with varMin,varMax
803 // arrays varMin and varMax contain the min and max values of each variable.
fb494025 804 // therefore varMin and varMax must have their dimensions equal to GetNVar()
c8df672e 805 //
806
807 THnSparse* clone = (THnSparse*)fData->Clone();
fb494025 808 for (Int_t iAxis=0; iAxis<GetNVar(); iAxis++) {
c8df672e 809 clone->GetAxis(iAxis)->SetRangeUser(varMin[iAxis],varMax[iAxis]);
810 }
811 return clone->Projection(iVar1,iVar2,iVar3);
812}
813
814//____________________________________________________________________
9105291d 815void AliCFGridSparse::SetRangeUser(Int_t iVar, Double_t varMin, Double_t varMax) {
816 //
817 // set range of axis iVar.
818 //
819 fData->GetAxis(iVar)->SetRangeUser(varMin,varMax);
820 AliWarning(Form("THnSparse axis %d range has been modified",iVar));
821}
822
823//____________________________________________________________________
fb494025 824void AliCFGridSparse::SetRangeUser(const Double_t *varMin, const Double_t *varMax) {
c8df672e 825 //
fb494025 826 // set range of every axis. varMin and varMax must be of dimension GetNVar()
c8df672e 827 //
fb494025 828 for (Int_t iAxis=0; iAxis<GetNVar() ; iAxis++) { // set new range for every axis
9105291d 829 SetRangeUser(iAxis,varMin[iAxis],varMax[iAxis]);
c8df672e 830 }
9105291d 831 AliWarning("THnSparse axes ranges have been modified");
832}
833
834//____________________________________________________________________
835void AliCFGridSparse::UseAxisRange(Bool_t b) const {
fb494025 836 for (Int_t iAxis=0; iAxis<GetNVar(); iAxis++) fData->GetAxis(iAxis)->SetBit(TAxis::kAxisRange,b);
837}
838
839//____________________________________________________________________
840Float_t AliCFGridSparse::GetOverFlows(Int_t ivar, Bool_t exclusive) const
841{
842 //
843 // Returns overflows in variable ivar
844 // Set 'exclusive' to true for an exclusive check on variable ivar
845 //
846 Int_t* bin = new Int_t[GetNVar()];
847 memset(bin, 0, sizeof(Int_t) * GetNVar());
848 Float_t ovfl=0.;
849 for (Long64_t i = 0; i < fData->GetNbins(); i++) {
850 Double_t v = fData->GetBinContent(i, bin);
851 Bool_t add=kTRUE;
852 if (exclusive) {
853 for(Int_t j=0;j<GetNVar();j++){
854 if(ivar==j)continue;
855 if((bin[j]==0) || (bin[j]==GetNBins(j)+1))add=kFALSE;
856 }
857 }
858 if(bin[ivar]==GetNBins(ivar)+1 && add) ovfl+=v;
859 }
860
861 delete[] bin;
862 return ovfl;
863}
864
865//____________________________________________________________________
866Float_t AliCFGridSparse::GetUnderFlows(Int_t ivar, Bool_t exclusive) const
867{
868 //
869 // Returns exclusive overflows in variable ivar
870 // Set 'exclusive' to true for an exclusive check on variable ivar
871 //
872 Int_t* bin = new Int_t[GetNVar()];
873 memset(bin, 0, sizeof(Int_t) * GetNVar());
874 Float_t unfl=0.;
875 for (Long64_t i = 0; i < fData->GetNbins(); i++) {
876 Double_t v = fData->GetBinContent(i, bin);
877 Bool_t add=kTRUE;
878 if (exclusive) {
879 for(Int_t j=0;j<GetNVar();j++){
880 if(ivar==j)continue;
881 if((bin[j]==0) || (bin[j]==GetNBins(j)+1))add=kFALSE;
882 }
883 }
884 if(bin[ivar]==0 && add) unfl+=v;
885 }
886
887 delete[] bin;
888 return unfl;
c8df672e 889}
fb494025 890