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