]> git.uio.no Git - u/mrichter/AliRoot.git/blame - CORRFW/AliCFGrid.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
added low and high flux parameters to the array
[u/mrichter/AliRoot.git] / CORRFW / AliCFGrid.cxx
CommitLineData
563113d0 1/* $Id$ */
1e9dad92 2/**************************************************************************
3 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * *
5 * Author: The ALICE Off-line Project. *
6 * Contributors are mentioned in the code where appropriate. *
7 * *
8 * Permission to use, copy, modify and distribute this software and its *
9 * documentation strictly for non-commercial purposes is hereby granted *
10 * without fee, provided that the above copyright notice appears in all *
11 * copies and that both the copyright notice and this permission notice *
12 * appear in the supporting documentation. The authors make no claims *
13 * about the suitability of this software for any purpose. It is *
14 * provided "as is" without express or implied warranty. *
15 **************************************************************************/
16//---------------------------------------------------------------------//
17// //
18// AliCFGrid Class //
19// Class to accumulate data on an N-dimensional grid, to be used //
20// as input to get corrections for Reconstruction & Trigger efficiency //
21// The class uses a one-dimensional array of floats to store the grid //
22// --Author : S.Arcelli //
23// Still to be done: //
24// --Implement methods to merge cells //
25// --Interpolate among bins in a range //
365c0ea4 26// This implementation will be eventually replaced by AliCFGridSparse //
1e9dad92 27//---------------------------------------------------------------------//
563113d0 28//
29//
9f6be3a2 30#include "AliLog.h"
563113d0 31#include "AliCFGrid.h"
32#include "TMath.h"
33#include "TROOT.h"
1e9dad92 34#include "TH1D.h"
35#include "TH2D.h"
36#include "TH3D.h"
563113d0 37
38//____________________________________________________________________
39ClassImp(AliCFGrid)
40
41//____________________________________________________________________
42AliCFGrid::AliCFGrid() :
1e9dad92 43 AliCFVGrid(),
563113d0 44 fNentriesTot(0),
45 fNunfl(0x0),
46 fNovfl(0x0),
47 fData(0x0),
48 fErr2(0x0)
49{
50 // default constructor
51}
52//____________________________________________________________________
53AliCFGrid::AliCFGrid(const Char_t* name, const Char_t* title) :
1e9dad92 54 AliCFVGrid(name,title),
563113d0 55 fNentriesTot(0),
56 fNunfl(0x0),
57 fNovfl(0x0),
58 fData(0x0),
59 fErr2(0x0)
60{
61 // default constructor
62}
63
64//____________________________________________________________________
1e9dad92 65AliCFGrid::AliCFGrid(const Char_t* name, const Char_t* title, const Int_t nVarIn, const Int_t * nBinIn, const Double_t *binLimitsIn) :
66 AliCFVGrid(name,title,nVarIn,nBinIn,binLimitsIn),
563113d0 67 fNentriesTot(0),
68 fNunfl(0x0),
69 fNovfl(0x0),
70 fData(0x0),
71 fErr2(0x0)
72{
73 //
74 // main constructor
75 //
1e9dad92 76
77 //The over/underflows
563113d0 78 fNunfl=new Float_t[fNVar];
79 fNovfl= new Float_t[fNVar];
80
81
563113d0 82 //Initialization
563113d0 83 fNentriesTot =0;
84 for(Int_t j=0;j<fNVar;j++){
85 fNunfl[j] =0;
86 fNovfl[j] =0;
87 }
1e9dad92 88
89
90 // the grid
91
318f64b1 92 fData = new Float_t[fNDim];
1e9dad92 93
94 //Initialization
318f64b1 95 for(Int_t j=0;j<fNDim;j++){
96 fData[j] =0;
97 }
1e9dad92 98
563113d0 99}
100
101//____________________________________________________________________
102AliCFGrid::AliCFGrid(const AliCFGrid& c) :
1e9dad92 103 AliCFVGrid(c),
563113d0 104 fNentriesTot(0),
105 fNunfl(0x0),
106 fNovfl(0x0),
107 fData(0x0),
108 fErr2(0x0)
109{
110 //
111 // copy constructor
112 //
113 ((AliCFGrid &)c).Copy(*this);
114}
115
116//____________________________________________________________________
117AliCFGrid::~AliCFGrid()
118{
119 //
120 // destructor
121 //
1e9dad92 122 if(fNunfl)delete fNunfl;
123 if(fNovfl)delete fNovfl;
124 if(fData)delete fData;
563113d0 125 if(fSumW2)delete fErr2;
563113d0 126
127}
563113d0 128//____________________________________________________________________
129AliCFGrid &AliCFGrid::operator=(const AliCFGrid &c)
130{
131 //
132 // assigment operator
133 //
134 if (this != &c)
135 ((AliCFGrid &) c).Copy(*this);
136 return *this;
137}
138//____________________________________________________________________
139Float_t AliCFGrid::GetElement(Int_t bin) const
140{
141 //
142 // Returns grid element bin
143 //
144 if(bin>=fNDim){
145 AliInfo(Form(" element index outside the grid, return -1"));
146 return -1.;
147 }
148 return fData[bin];
149}
150//____________________________________________________________________
151Float_t AliCFGrid::GetElement(Int_t *bin) const
152{
153 //
154 // Get the content in a bin corresponding to a set of bin indexes
155 //
1e9dad92 156 //-1 is to move from TH/ThnSparse N-dim bin convention to one in AliCFFrame
157 for(Int_t i=0;i<fNVar; i++)fIndex[i]=bin[i]-1;
158 Int_t ind =GetBinIndex(fIndex);
159 return GetElement(ind);
563113d0 160}
161//____________________________________________________________________
1e9dad92 162Float_t AliCFGrid::GetElement(Double_t *var) const
563113d0 163{
164 //
165 // Get the content in a bin corresponding to a set of input variables
166 //
167 Int_t unfl=0;
168 Int_t ovfl=0;
169 for(Int_t i=0;i<fNVar;i++){
170 Int_t nbins=fNVarBins[i]+1;
365c0ea4 171 Double_t *bins=new Double_t[nbins];
563113d0 172 for(Int_t ibin =0;ibin<nbins;ibin++){
173 bins[ibin] = fVarBinLimits[ibin+fOffset[i]];
174 }
175
176 fIndex[i] = TMath::BinarySearch(nbins,bins,var[i]);
177
178 //underflows
179
180 if(var[i] < bins[0]){
181 unfl=1;
182 }
183
184 //overflows
185
186 if(var[i] > bins[nbins-1]){
187 ovfl=1;
188 }
189 delete [] bins;
190 }
191
1e9dad92 192
193 //move to the TH/THnSparse convention in N-dim bin numbering
194 for(Int_t i=0;i<fNVar; i++)fIndex[i]+=1;
195
563113d0 196 if(!(ovfl==1 || unfl==1)){
197 return GetElement(fIndex);
198 }
199 else{
200 AliInfo(Form(" input variables outside the grid, return -1"));
201 return -1.;
202 }
203}
204//____________________________________________________________________
205Float_t AliCFGrid::GetElementError(Int_t iel) const
206{
207 //
1e9dad92 208 // Return the error on grid element iel
563113d0 209 //
210 if(iel>=fNDim){
211 AliInfo(Form(" element index outside the grid, return -1"));
212 return -1.;
213 }
1e9dad92 214 if(fSumW2)return TMath::Sqrt(fErr2[iel]);
215 return TMath::Sqrt(fData[iel]);
563113d0 216}
217//____________________________________________________________________
218Float_t AliCFGrid::GetElementError(Int_t *bin) const
219{
220 //
1e9dad92 221 // Get the error in a bin corresponding to a set of bin indeces
563113d0 222 //
1e9dad92 223 //-1 is to move from TH/ThnSparse N-dim bin convention to one in AliCFFrame
224 for(Int_t i=0;i<fNVar; i++)fIndex[i]=bin[i]-1;
225 Int_t ind =GetBinIndex(fIndex);
226 return GetElementError(ind);
563113d0 227
228}
229//____________________________________________________________________
1e9dad92 230Float_t AliCFGrid::GetElementError(Double_t *var) const
563113d0 231{
232 //
1e9dad92 233 // Get the error in a bin corresponding to a set of input variables
563113d0 234 //
235 Int_t unfl=0;
236 Int_t ovfl=0;
237 for(Int_t i=0;i<fNVar;i++){
238 Int_t nbins=fNVarBins[i]+1;
365c0ea4 239 Double_t *bins=new Double_t[nbins];
563113d0 240 for(Int_t ibin =0;ibin<nbins;ibin++){
241 bins[ibin] = fVarBinLimits[ibin+fOffset[i]];
242 }
243
244 fIndex[i] = TMath::BinarySearch(nbins,bins,var[i]);
245 //underflows
246
247 if(var[i] < bins[0]){
248 unfl=1;
249 }
250
251 //overflows
252
253 if(var[i] > bins[nbins-1]){
254 ovfl=1;
255 }
256 delete [] bins;
257 }
258
1e9dad92 259 //move to the TH/THnSparse convention in N-dim bin numbering
260 for(Int_t i=0;i<fNVar; i++)fIndex[i]+=1;
563113d0 261 if(!(ovfl==1 || unfl==1)){
262 return GetElementError(fIndex);
263 }
264 else{
265 AliInfo(Form(" input variables outside the grid, return -1"));
266 return -1.;
267 }
268}
269//____________________________________________________________________
270void AliCFGrid::SetElement(Int_t iel, Float_t val)
271{
272 //
273 // Sets grid element iel to val
274 //
275 if(iel>=fNDim){
276 AliInfo(Form(" element index outside the grid, no value set"));
277 }else {
278 fData[iel]=val;
279 }
280}
281//____________________________________________________________________
282void AliCFGrid::SetElement(Int_t *bin, Float_t val)
283{
284 //
285 // Sets grid element of bin indeces bin to val
286 //
1e9dad92 287 //-1 is to move from TH/ThnSparse N-dim bin convention to one in AliCFFrame
288 for(Int_t i=0;i<fNVar; i++)fIndex[i]=bin[i]-1;
289 Int_t ind =GetBinIndex(fIndex);
290 SetElement(ind,val);
563113d0 291}
292//____________________________________________________________________
1e9dad92 293void AliCFGrid::SetElement(Double_t *var, Float_t val)
563113d0 294{
295 //
296 // Set the content in a bin to value val corresponding to a set of input variables
297 //
298 Int_t unfl=0;
299 Int_t ovfl=0;
300 for(Int_t i=0;i<fNVar;i++){
301 Int_t nbins=fNVarBins[i]+1;
365c0ea4 302 Double_t *bins=new Double_t[nbins];
563113d0 303 for(Int_t ibin =0;ibin<nbins;ibin++){
304 bins[ibin] = fVarBinLimits[ibin+fOffset[i]];
305 }
306
307 fIndex[i] = TMath::BinarySearch(nbins,bins,var[i]);
308 //underflows
309
310 if(var[i] < bins[0]){
311 unfl=1;
312 }
313
314 //overflows
315
316 if(var[i] > bins[nbins-1]){
317 ovfl=1;
318 }
319 delete [] bins;
320 }
321
1e9dad92 322 //move to the TH/THnSparse convention in N-dim bin numbering
323 for(Int_t i=0;i<fNVar; i++)fIndex[i]+=1;
563113d0 324 if(!(ovfl==1 || unfl==1)){
325 SetElement(fIndex,val);
326 }
327 else{
328 AliInfo(Form(" input variables outside the grid, no value set"));
329 }
330}
331//____________________________________________________________________
332void AliCFGrid::SetElementError(Int_t iel, Float_t val)
333{
334 //
1e9dad92 335 // Set squared error on grid element iel to val*val
563113d0 336 //
337 if(iel>=fNDim){
338 AliInfo(Form(" element index outside the grid, no value set"));
339 return;
340 }
341 if(!fErr2)SumW2();
1e9dad92 342 fErr2[iel]=val*val;
563113d0 343}
344//____________________________________________________________________
345void AliCFGrid::SetElementError(Int_t *bin, Float_t val)
346{
347 //
348 // Set squared error to val on grid element of bin indeces bin
349 //
1e9dad92 350 //-1 is to move from TH/ThnSparse N-dim bin convention to one in AliCFFrame
351 for(Int_t i=0;i<fNVar; i++)fIndex[i]=bin[i]-1;
352 Int_t ind =GetBinIndex(fIndex);
353 SetElementError(ind,val);
563113d0 354}
355//____________________________________________________________________
1e9dad92 356void AliCFGrid::SetElementError(Double_t *var, Float_t val)
563113d0 357{
358 //
359 // Set squared error to val in a bin corresponding to a set of input variables
360 //
361 Int_t unfl=0;
362 Int_t ovfl=0;
363 for(Int_t i=0;i<fNVar;i++){
364 Int_t nbins=fNVarBins[i]+1;
365c0ea4 365 Double_t *bins=new Double_t[nbins];
563113d0 366 for(Int_t ibin =0;ibin<nbins;ibin++){
367 bins[ibin] = fVarBinLimits[ibin+fOffset[i]];
368 }
369
370 fIndex[i] = TMath::BinarySearch(nbins,bins,var[i]);
371 //underflows
372
373 if(var[i] < bins[0]){
374 unfl=1;
375 }
376
377 //overflows
378
379 if(var[i] > bins[nbins-1]){
380 ovfl=1;
381 }
382 delete [] bins;
383 }
384
1e9dad92 385 //move to the TH/THnSparse convention in N-dim bin numbering
386 for(Int_t i=0;i<fNVar; i++)fIndex[i]+=1;
387
563113d0 388 if(!(ovfl==1 || unfl==1)){
389 SetElementError(fIndex,val);
390 }
391 else{
392 AliInfo(Form(" input variables outside the grid, no value set"));
393 }
394}
395//____________________________________________________________________
1e9dad92 396void AliCFGrid::Fill(Double_t *var, Double_t weight)
563113d0 397{
1e9dad92 398
563113d0 399 //
1e9dad92 400 // Fill the grid,
401 // given a set of values of the input variable,
402 // with weight (by default w=1)
563113d0 403 //
1e9dad92 404
405
406 Int_t isunfl=0;
407 Int_t isovfl=0;
408 Int_t *unfl=new Int_t[fNVar];
409 Int_t *ovfl=new Int_t[fNVar];
410
411 for(Int_t i=0;i<fNVar;i++){
412 unfl[i]=0;
413 ovfl[i]=0;
563113d0 414 }
1e9dad92 415
563113d0 416 for(Int_t i=0;i<fNVar;i++){
417 Int_t nbins=fNVarBins[i]+1;
365c0ea4 418 Double_t *bins=new Double_t[nbins];
563113d0 419 for(Int_t ibin =0;ibin<nbins;ibin++){
420 bins[ibin] = fVarBinLimits[ibin+fOffset[i]];
421 }
422
423 fIndex[i] = TMath::BinarySearch(nbins,bins,var[i]);
424 //underflows
425
426 if(var[i] < bins[0]){
1e9dad92 427 unfl[i]=1;
428 isunfl=1;
563113d0 429 }
430
431 //overflows
432
433 if(var[i] > bins[nbins-1]){
1e9dad92 434 ovfl[i]=1;
435 isovfl=1;
563113d0 436 }
437 delete [] bins;
438 }
439
1e9dad92 440 //exclusive under/overflows
563113d0 441
563113d0 442 for(Int_t i=0;i<fNVar;i++){
1e9dad92 443 Bool_t add=kTRUE;
444 for(Int_t j=0;j<fNVar;j++){
445 if(i==j)continue;
446 if(!(unfl[j]==0 && ovfl[j]==0))add=kFALSE;
563113d0 447 }
1e9dad92 448 if(add && unfl[i]==1)fNunfl[i]++;
449 if(add && ovfl[i]==1)fNovfl[i]++;
563113d0 450 }
451
1e9dad92 452 delete [] unfl;
453 delete [] ovfl;
454
563113d0 455 // Total number of entries, overflows and underflows
456
457 fNentriesTot++;
563113d0 458
459 //if not ovfl/unfl, fill the element
1e9dad92 460
461 if(!(isovfl==1 || isunfl==1)){
563113d0 462 Int_t ind =GetBinIndex(fIndex);
463 fData[ind]+=weight;
464 if(fSumW2)fErr2[ind]+=(weight*weight);
465 }
466}
467//____________________________________________________________________
468Float_t AliCFGrid::GetOverFlows( Int_t ivar) const {
469 //
470 // Get overflows in variable var
471 //
472 return fNovfl[ivar];
473}
474//____________________________________________________________________
563113d0 475Float_t AliCFGrid::GetUnderFlows( Int_t ivar) const {
476 //
477 // Get overflows in variable var
478 //
479 return fNunfl[ivar];
480}
481//____________________________________________________________________
563113d0 482Float_t AliCFGrid::GetEntries( ) const {
483 //
484 // Get total entries (in grid + over/underflows)
485 //
486 return fNentriesTot;
487}
563113d0 488//___________________________________________________________________
1e9dad92 489TH1D *AliCFGrid::Project(Int_t ivar) const
563113d0 490{
491 //
492 // Make a 1D projection along variable ivar
493
494
495 Int_t nbins =fNVarBins[ivar];
496 Float_t *bins = new Float_t[nbins+1];
497 for (Int_t i=0;i<=fNVar;i++){
498 }
499 for(Int_t ibin =0;ibin<nbins+1;ibin++){
500 bins[ibin] = fVarBinLimits[ibin+fOffset[ivar]];
501 }
502
503 char pname[40];
504 sprintf(pname,"%s%s_%i",GetName(),"_proj1D_var", ivar);
505 char htitle[40];
506 sprintf(htitle,"%s%s_%i",GetName(),"_proj1D_var", ivar);
507
1e9dad92 508 TH1D *proj1D=0;
563113d0 509
510 //check if a projection with identical name exist
511 TObject *obj = gROOT->FindObject(pname);
1e9dad92 512 if (obj && obj->InheritsFrom("TH1D")) {
513 proj1D = (TH1D*)obj;
563113d0 514 proj1D->Reset();
515 }
516
517 if(!proj1D){
1e9dad92 518 proj1D =new TH1D(pname,htitle, nbins, bins);
563113d0 519 }
520
521 delete [] bins;
522 Float_t sum=0;
523 Float_t *data= new Float_t[nbins];
524 Float_t *err= new Float_t[nbins];
525
526 for(Int_t ibin=0;ibin<nbins;ibin++)data[ibin]=0;
527 for(Int_t ibin=0;ibin<nbins;ibin++)err[ibin]=0;
528 for(Int_t iel=0;iel<fNDim;iel++){
529 data[GetBinIndex(ivar,iel)]+=fData[iel];
530 if(fSumW2)err[GetBinIndex(ivar,iel)]+=fErr2[iel];
531 }
532
533 for(Int_t ibin =0;ibin<nbins;ibin++){
534 proj1D->SetBinContent(ibin+1,data[ibin]);
535 proj1D->SetBinError(ibin+1,TMath::Sqrt(data[ibin]));
536 if(fSumW2)proj1D->SetBinError(ibin+1,TMath::Sqrt(err[ibin]));
537 sum+=data[ibin];
538 }
539
540 delete [] data;
541 delete [] err;
542 proj1D->SetBinContent(nbins+1,GetOverFlows(ivar));
543 proj1D->SetBinContent(0,GetUnderFlows(ivar));
318f64b1 544 proj1D->SetEntries(fNentriesTot);
563113d0 545 return proj1D;
546}
547
548//___________________________________________________________________
1e9dad92 549TH2D *AliCFGrid::Project(Int_t ivar1, Int_t ivar2) const
563113d0 550{
551 //
552 // Make a 2D projection along variable ivar
553
554 Int_t nbins1 =fNVarBins[ivar1];
555 Int_t nbins2 =fNVarBins[ivar2];
556
557 Float_t *bins1 = new Float_t[nbins1+1];
558 Float_t *bins2 = new Float_t[nbins2+1];
559
560 for(Int_t ibin =0;ibin<nbins1+1;ibin++){
561 bins1[ibin] = fVarBinLimits[ibin+fOffset[ivar1]];
562 }
563 for(Int_t ibin =0;ibin<nbins2+1;ibin++){
564 bins2[ibin] = fVarBinLimits[ibin+fOffset[ivar2]];
565 }
566
567 char pname[40];
568 sprintf(pname,"%s%s_%i_%i",GetName(),"_proj2D_var",ivar1,ivar2);
569 char htitle[40];
570 sprintf(htitle,"%s%s_%i_%i",GetName(),"_proj2D_var",ivar1,ivar2);
571
1e9dad92 572 TH2D *proj2D=0;
563113d0 573
574 //check if a projection with identical name exist
575 TObject *obj = gROOT->FindObject(pname);
1e9dad92 576 if (obj && obj->InheritsFrom("TH2D")) {
577 proj2D = (TH2D*)obj;
563113d0 578 proj2D->Reset();
579 }
580
581 if(!proj2D){
1e9dad92 582 proj2D =new TH2D(pname,htitle, nbins1, bins1,nbins2,bins2);
563113d0 583 }
584
585 delete [] bins1;
586 delete [] bins2;
587
588
589 Float_t sum=0;
590 Float_t **data=new Float_t*[nbins1];
591 Float_t *data2=new Float_t[nbins1*nbins2];
592 Float_t **err=new Float_t*[nbins1];
593 Float_t *err2=new Float_t[nbins1*nbins2];
594 for(Int_t i=0;i<nbins1;i++)data[i] = data2+i*nbins2;
595 for(Int_t i=0;i<nbins1;i++)err[i] = err2+i*nbins2;
596
597 for(Int_t ibin1 =0;ibin1<nbins1;ibin1++){
598 for(Int_t ibin2 =0;ibin2<nbins2;ibin2++){
599 data[ibin1][ibin2]=0;
600 err[ibin1][ibin2]=0;
601 }
602 }
603
604 for(Int_t iel=0;iel<fNDim;iel++){
605 data[GetBinIndex(ivar1,iel)][GetBinIndex(ivar2,iel)]+=fData[iel];
606 if(fSumW2)err[GetBinIndex(ivar1,iel)][GetBinIndex(ivar2,iel)]+=fErr2[iel];
607 }
608
609 for(Int_t ibin1 =0;ibin1<nbins1;ibin1++){
610 for(Int_t ibin2 =0;ibin2<nbins2;ibin2++){
611 proj2D->SetBinContent(ibin1+1,ibin2+1,data[ibin1][ibin2]);
612 proj2D->SetBinError(ibin1+1,ibin2+1,TMath::Sqrt(data[ibin1][ibin2]));
613 if(fSumW2)proj2D->SetBinError(ibin1+1,ibin2+1,TMath::Sqrt(err[ibin1][ibin2]));
614 sum+=data[ibin1][ibin2];
615 }
616
617 }
618 delete data;
619 delete data2;
620 delete err;
621 delete err2;
622 proj2D->SetBinContent(0,nbins2/2,GetUnderFlows(ivar1));
623 proj2D->SetBinContent(nbins1+1,nbins2/2,GetOverFlows(ivar1));
624 proj2D->SetBinContent(nbins1/2,0,GetUnderFlows(ivar2));
625 proj2D->SetBinContent(nbins1/2,nbins2+1,GetOverFlows(ivar2));
318f64b1 626 proj2D->SetEntries(fNentriesTot);
563113d0 627 return proj2D;
628}
629//___________________________________________________________________
1e9dad92 630TH3D *AliCFGrid::Project(Int_t ivar1, Int_t ivar2, Int_t ivar3) const
563113d0 631{
632 //
633 // Make a 3D projection along variable ivar
634
635 Int_t nbins1 =fNVarBins[ivar1];
636 Int_t nbins2 =fNVarBins[ivar2];
637 Int_t nbins3 =fNVarBins[ivar3];
638
639 Float_t *bins1 = new Float_t[nbins1+1];
640 Float_t *bins2 = new Float_t[nbins2+1];
641 Float_t *bins3 = new Float_t[nbins3+1];
642
643 for(Int_t ibin =0;ibin<nbins1+1;ibin++){
644 bins1[ibin] = fVarBinLimits[ibin+fOffset[ivar1]];
645 }
646 for(Int_t ibin =0;ibin<nbins2+1;ibin++){
647 bins2[ibin] = fVarBinLimits[ibin+fOffset[ivar2]];
648 }
649 for(Int_t ibin =0;ibin<nbins3+1;ibin++){
650 bins3[ibin] = fVarBinLimits[ibin+fOffset[ivar3]];
651 }
652
653 char pname[40];
654 sprintf(pname,"%s%s_%i_%i_%i",GetName(),"_proj3D_var",ivar1,ivar2,ivar3);
655 char htitle[40];
656 sprintf(htitle,"%s%s_%i_%i_%i",GetName(),"_proj3D_var",ivar1,ivar2,ivar3);
657
1e9dad92 658 TH3D *proj3D=0;
563113d0 659
660 //check if a projection with identical name exist
661 TObject *obj = gROOT->FindObject(pname);
1e9dad92 662 if (obj && obj->InheritsFrom("TH3D")) {
663 proj3D = (TH3D*)obj;
563113d0 664 proj3D->Reset();
665 }
666
667 if(!proj3D){
1e9dad92 668 proj3D =new TH3D(pname,htitle, nbins1,bins1,nbins2,bins2,nbins3,bins3);
563113d0 669 }
670
671 delete [] bins1;
672 delete [] bins2;
673 delete [] bins3;
674
675
676 Float_t sum=0;
677 Float_t ***data=new Float_t**[nbins1];
678 Float_t **data2=new Float_t*[nbins1*nbins2];
679 Float_t *data3=new Float_t[nbins1*nbins2*nbins3];
680 Float_t ***err=new Float_t**[nbins1];
681 Float_t **err2=new Float_t*[nbins1*nbins2];
682 Float_t *err3=new Float_t[nbins1*nbins2*nbins3];
683 for(Int_t i=0;i<nbins1;i++)data[i] = data2+i*nbins2;
684 for(Int_t i=0;i<nbins1;i++)err[i] = err2+i*nbins2;
685 for(Int_t i=0;i<nbins1;i++){
686 for(Int_t j=0;j<nbins2;j++){
687 data[i][j] = data3+i*nbins2*nbins3+j*nbins3;
688 err[i][j] = err3+i*nbins2*nbins3+j*nbins3;
689 }
690 }
691 for(Int_t ibin1 =0;ibin1<nbins1;ibin1++){
692 for(Int_t ibin2 =0;ibin2<nbins2;ibin2++){
693 for(Int_t ibin3 =0;ibin3<nbins3;ibin3++){
694 data[ibin1][ibin2][ibin3]=0;
695 err[ibin1][ibin2][ibin3]=0;
696 }
697 }
698 }
699
700 for(Int_t iel=0;iel<fNDim;iel++){
701 data[GetBinIndex(ivar1,iel)][GetBinIndex(ivar2,iel)][GetBinIndex(ivar3,iel)]+=fData[iel];
702 if(fSumW2)err[GetBinIndex(ivar1,iel)][GetBinIndex(ivar2,iel)][GetBinIndex(ivar3,iel)]+=fErr2[iel];
703 }
704
705 for(Int_t ibin1 =0;ibin1<nbins1;ibin1++){
706 for(Int_t ibin2 =0;ibin2<nbins2;ibin2++){
707 for(Int_t ibin3 =0;ibin3<nbins3;ibin3++){
708 proj3D->SetBinContent(ibin1+1,ibin2+1,ibin3+1,data[ibin1][ibin2][ibin3]);
709 proj3D->SetBinError(ibin1+1,ibin2+1,ibin3+1,TMath::Sqrt(data[ibin1][ibin2][ibin3]));
710 if(fSumW2)proj3D->SetBinError(ibin1+1,ibin2+1,ibin3+1,TMath::Sqrt(err[ibin1][ibin2][ibin3]));
711 sum+=data[ibin1][ibin2][ibin3];
712 }
713 }
714 }
715
716 delete data;
717 delete data2;
718 delete data3;
719 delete err;
720 delete err2;
721 delete err3;
318f64b1 722
723 proj3D->SetEntries(fNentriesTot);
563113d0 724 return proj3D;
725}
726
727//___________________________________________________________________
1e9dad92 728TH1D *AliCFGrid::Slice(Int_t ivar, Double_t *varMin, Double_t* varMax) const
563113d0 729{
730 //
731 // Make a slice along variable ivar in range [varMin,varMax]
732
733
734 Int_t nbins =fNVarBins[ivar];
735 Float_t *bins = new Float_t[nbins+1];
736 for (Int_t i=0;i<=fNVar;i++){
737 }
738 for(Int_t ibin =0;ibin<nbins+1;ibin++){
739 bins[ibin] = fVarBinLimits[ibin+fOffset[ivar]];
740 }
741
742 char pname[40];
743 sprintf(pname,"%s%s_%i",GetName(),"_proj1D_var", ivar);
744 char htitle[40];
745 sprintf(htitle,"%s%s_%i",GetName(),"_proj1D_var", ivar);
746
1e9dad92 747 TH1D *proj1D=0;
563113d0 748
749 //check if a projection with identical name exist
750 TObject *obj = gROOT->FindObject(pname);
1e9dad92 751 if (obj && obj->InheritsFrom("TH1D")) {
752 proj1D = (TH1D*)obj;
563113d0 753 proj1D->Reset();
754 }
755
756 if(!proj1D){
1e9dad92 757 proj1D =new TH1D(pname,htitle, nbins, bins);
563113d0 758 }
759
760 delete [] bins;
761
762
763 Int_t *indexMin=new Int_t[fNVar];
764 Int_t *indexMax=new Int_t[fNVar];
765
766
767 //Find out the min and max bins
768
769 for(Int_t i=0;i<fNVar;i++){
365c0ea4 770 Double_t xmin=varMin[i]; // the min values
771 Double_t xmax=varMax[i]; // the max values
cdd7aeae 772 Int_t nBins=fNVarBins[i]+1;
773 Double_t *Bins=new Double_t[nBins];
774 for(Int_t ibin =0;ibin<nBins;ibin++){
775 Bins[ibin] = fVarBinLimits[ibin+fOffset[i]];
563113d0 776 }
cdd7aeae 777 indexMin[i] = TMath::BinarySearch(nBins,Bins,xmin);
778 indexMax[i] = TMath::BinarySearch(nBins,Bins,xmax);
779 delete [] Bins;
563113d0 780 }
781
782 Float_t sum=0;
783 Float_t *data= new Float_t[nbins];
784 for(Int_t ibin=0;ibin<nbins;ibin++)data[ibin]=0;
785
786 Int_t *index= new Int_t[fNVar];
787 Int_t ielmin=GetBinIndex(indexMin);
788 Int_t ielmax=GetBinIndex(indexMax);
365c0ea4 789
790
563113d0 791 for(Int_t iel=ielmin;iel<=ielmax;iel++){
792 GetBinIndex(iel,index);
793 Bool_t isIn=kTRUE;
794 for (Int_t j=0;j<fNVar;j++){
795 if(!(index[j]>=indexMin[j] && index[j]<=indexMax[j]))isIn=kFALSE;
796 break;
797 }
798 if(isIn)data[GetBinIndex(ivar,iel)]+=fData[iel];
799 }
800
801 delete [] index;
365c0ea4 802 delete [] indexMin;
803 delete [] indexMax;
563113d0 804
805
806 for(Int_t ibin =0;ibin<nbins;ibin++){
807 proj1D->SetBinContent(ibin+1,data[ibin]);
808 proj1D->SetBinError(ibin+1,TMath::Sqrt(data[ibin]));
809 sum+=data[ibin];
810 }
811
812 delete [] data;
813
814 proj1D->SetEntries(sum);
815 return proj1D;
816}
817
563113d0 818
819//____________________________________________________________________
1e9dad92 820void AliCFGrid::Add(AliCFVGrid* aGrid, Double_t c)
563113d0 821{
822 //
823 //add aGrid to the current one
824 //
825
826 if(aGrid->GetNVar()!=fNVar){
827 AliInfo("Different number of variables, cannot add the grids");
828 return;
829 }
830 if(aGrid->GetNDim()!=fNDim){
831 AliInfo("Different number of dimensions, cannot add the grids!");
832 return;
833 }
834
835 if(!fSumW2 && aGrid->GetSumW2())SumW2();
836
837 for(Int_t iel=0;iel<fNDim;iel++){
838 fData[iel]+=(c*aGrid->GetElement(iel));
839 if(fSumW2){
1e9dad92 840 Float_t err=aGrid->GetElementError(iel);
563113d0 841 fErr2[iel]+=c*c*err*err;
842 }
843 }
844
845 //Add entries, overflows and underflows
846
847 fNentriesTot+= c*aGrid->GetEntries();
563113d0 848 for(Int_t j=0;j<fNVar;j++){
849 fNunfl[j]+= c*aGrid->GetUnderFlows(j);
db6722a5 850 fNovfl[j]+= c*aGrid->GetOverFlows(j);
563113d0 851 }
852}
853//____________________________________________________________________
1e9dad92 854void AliCFGrid::Add(AliCFVGrid* aGrid1, AliCFVGrid* aGrid2, Double_t c1,Double_t c2)
563113d0 855{
856 //
857 //add aGrid1 and aGrid2
858 //
859
860 if(fNVar!=aGrid1->GetNVar()|| fNVar!=aGrid2->GetNVar()){
861 AliInfo("Different number of variables, cannot add the grids");
862 return;
863 }
864 if(fNDim!=aGrid1->GetNDim()|| fNDim!=aGrid2->GetNDim()){
865 AliInfo("Different number of dimensions, cannot add the grids!");
866 return;
867 }
868
869 if(!fSumW2 && (aGrid1->GetSumW2() || aGrid2->GetSumW2()))SumW2();
870
871 Float_t cont1,cont2,err1,err2;
872
873 for(Int_t iel=0;iel<fNDim;iel++){
874 cont1=aGrid1->GetElement(iel);
875 cont2=aGrid2->GetElement(iel);
876 SetElement(iel,c1*cont1+c2*cont2);
877 if(fSumW2){
1e9dad92 878 err1=aGrid1->GetElementError(iel);
879 err2=aGrid2->GetElementError(iel);
880 SetElementError(iel,TMath::Sqrt(c1*c1*err1*err1+c2*c2*err2*err2));
563113d0 881 }
882 }
883
884 //Add entries, overflows and underflows
885
886 fNentriesTot= c1*aGrid1->GetEntries()+c2*aGrid2->GetEntries();
563113d0 887 for(Int_t j=0;j<fNVar;j++){
888 fNunfl[j]= c1*aGrid1->GetUnderFlows(j)+c2*aGrid2->GetUnderFlows(j);
db6722a5 889 fNovfl[j]= c1*aGrid1->GetOverFlows(j)+c2*aGrid2->GetOverFlows(j);
563113d0 890 }
891}
892//____________________________________________________________________
1e9dad92 893void AliCFGrid::Multiply(AliCFVGrid* aGrid, Double_t c)
563113d0 894{
895 //
896 //multiply grid aGrid by the current one
897 //
898
899 if(aGrid->GetNVar()!=fNVar){
900 AliInfo("Different number of variables, cannot multiply the grids");
901 return;
902 }
903 if(aGrid->GetNDim()!=fNDim){
904 AliInfo("Different number of dimensions, cannot multiply the grids!");
905 return;
906 }
907
908 if(!fSumW2 && aGrid->GetSumW2())SumW2();
909
910 Float_t cont1,cont2,err1,err2;
911
912 for(Int_t iel=0;iel<fNDim;iel++){
913 cont1=GetElement(iel);
914 cont2=c*aGrid->GetElement(iel);
915 SetElement(iel,cont1*cont2);
916 if(fSumW2){
1e9dad92 917 err1=GetElementError(iel);
918 err2=aGrid->GetElementError(iel);
919 SetElementError(iel,TMath::Sqrt(c*c*(cont2*cont2*err1*err1+cont1*cont1*err2*err2)));
563113d0 920 }
921 }
922
923 //Set entries to the number of bins, preserve original overflows and underflows
924
925 fNentriesTot=fNDim;
563113d0 926 for(Int_t j=0;j<fNVar;j++){
927 fNunfl[j]= GetUnderFlows(j);
db6722a5 928 fNovfl[j]= GetOverFlows(j);
563113d0 929 }
930}
931//____________________________________________________________________
1e9dad92 932void AliCFGrid::Multiply(AliCFVGrid* aGrid1, AliCFVGrid* aGrid2, Double_t c1, Double_t c2)
563113d0 933{
934 //
935 //multiply grids aGrid1 and aGrid2
936 //
937
938 if(fNVar!=aGrid1->GetNVar()|| fNVar!=aGrid2->GetNVar()){
939 AliInfo("Different number of variables, cannot multiply the grids");
940 return;
941 }
942 if(fNDim!=aGrid1->GetNDim()|| fNDim!=aGrid2->GetNDim()){
943 AliInfo("Different number of dimensions, cannot multiply the grids!");
944 return;
945 }
946
947 if(!fSumW2 && (aGrid1->GetSumW2() || aGrid2->GetSumW2()))SumW2();
948
949 Float_t cont1,cont2,err1,err2;
950 for(Int_t iel=0;iel<fNDim;iel++){
951 cont1=c1*aGrid1->GetElement(iel);
952 cont2=c2*aGrid2->GetElement(iel);
953 SetElement(iel,cont1*cont2);
954 if(fSumW2){
1e9dad92 955 err1=aGrid1->GetElementError(iel);
956 err2=aGrid2->GetElementError(iel);
957 SetElementError(iel,TMath::Sqrt(c1*c1*c2*c2*(cont2*cont2*err1*err1+cont1*cont1*err2*err2)));
563113d0 958 }
959 }
960
961 //Set entries to the number of bins, preserve original overflows and underflows
962
963 fNentriesTot=fNDim;
563113d0 964 for(Int_t j=0;j<fNVar;j++){
965 fNunfl[j]= GetUnderFlows(j);
db6722a5 966 fNovfl[j]= GetOverFlows(j);
563113d0 967 }
968}
969//____________________________________________________________________
1e9dad92 970void AliCFGrid::Divide(AliCFVGrid* aGrid, Double_t c)
563113d0 971{
972 //
973 //divide current grid by grid aGrid
974 //
975
563113d0 976 if(aGrid->GetNVar()!=fNVar){
977 AliInfo("Different number of variables, cannot divide the grids");
978 return;
979 }
980 if(aGrid->GetNDim()!=fNDim){
981 AliInfo("Different number of dimensions, cannot divide the grids!");
982 return;
983 }
984 if(!c){AliInfo(Form("c is %f, cannot divide!",c)); return;}
985
986 if(!fSumW2 && aGrid->GetSumW2())SumW2();
987
1e9dad92 988 Float_t cont1,cont2,err1,err2,den;
563113d0 989 for(Int_t iel=0;iel<fNDim;iel++){
990 cont1=GetElement(iel);
991 cont2=aGrid->GetElement(iel);
992 if(cont2)SetElement(iel,cont1/(c*cont2));
993 else SetElement(iel,0);
994 if(fSumW2){
1e9dad92 995 err1=GetElementError(iel);
996 err2=aGrid->GetElementError(iel);
563113d0 997 if(!cont2){SetElementError(iel,0.); continue;}
1e9dad92 998 den=cont2*cont2*cont2*c*c;
999 SetElementError(iel,TMath::Sqrt((cont2*cont2*err1*err1+cont1*cont1*err2*err2)/den));
563113d0 1000 }
1001 }
1e9dad92 1002
563113d0 1003 //Set entries to the number of bins, preserve original overflows and underflows
1004
1005 fNentriesTot=fNDim;
563113d0 1006 for(Int_t j=0;j<fNVar;j++){
1007 fNunfl[j]= GetUnderFlows(j);
db6722a5 1008 fNovfl[j]= GetOverFlows(j);
563113d0 1009 }
1010}
1011//____________________________________________________________________
1e9dad92 1012void AliCFGrid::Divide(AliCFVGrid* aGrid1, AliCFVGrid* aGrid2, Double_t c1,Double_t c2, Option_t *option)
563113d0 1013{
1014 //
1015 //divide grids aGrid1,aGrid2
1016 //
1017
1018 TString opt = option;
1019 opt.ToUpper();
1020
1021 if(fNVar!=aGrid1->GetNVar()|| fNVar!=aGrid2->GetNVar()){
1022 AliInfo("Different number of variables, cannot divide the grids");
1023 return;
1024 }
1025 if(fNDim!=aGrid1->GetNDim()|| fNDim!=aGrid2->GetNDim()){
1026 AliInfo("Different number of dimensions, cannot divide the grids!");
1027 return;
1028 }
1029 if(!c2){AliInfo(Form("c2 is %f, cannot divide!",c2)); return;}
1030
1031 if(!fSumW2 && (aGrid1->GetSumW2() || aGrid2->GetSumW2()))SumW2();
1032
1033 Float_t cont1,cont2,err1,err2,r,den;
1034
1035 for(Int_t iel=0;iel<fNDim;iel++){
1036 cont1=aGrid1->GetElement(iel);
1037 cont2=aGrid2->GetElement(iel);
1038 if(cont2)SetElement(iel,c1*cont1/(c2*cont2));
1039 else SetElement(iel,0);
1040 if(fSumW2){
1e9dad92 1041 err1=aGrid1->GetElementError(iel);
1042 err2=aGrid2->GetElementError(iel);
563113d0 1043 if(!cont2){SetElementError(iel,0.); continue;}
1044 if (opt.Contains("B")){
1045 if(cont1!=cont2){
1046 r=cont1/cont2;
1e9dad92 1047 SetElementError(iel,TMath::Sqrt(TMath::Abs(((1.-2.*r)*err1*err1+r*r*err2*err2)/(cont2*cont2))));
563113d0 1048 }else{
1049 SetElementError(iel,0.);
1050 }
1051 }else{
1052 den=cont2*cont2*cont2*cont2*c2*c2;
1e9dad92 1053 SetElementError(iel,TMath::Sqrt(c1*c1*(cont2*cont2*err1*err1+cont1*cont1*err2*err2)/den));
563113d0 1054 }
1055 }
1056 }
1057
1058 //Set entries to the number of bins, preserve original overflows and underflows
1059
1060 fNentriesTot=fNDim;
563113d0 1061 for(Int_t j=0;j<fNVar;j++){
1062 fNunfl[j]= GetUnderFlows(j);
db6722a5 1063 fNovfl[j]= GetOverFlows(j);
563113d0 1064 }
1065}
1066//____________________________________________________________________
7411edfd 1067void AliCFGrid::Rebin(const Int_t* group)
1068{
1069 //
1070 // Not yet implemented
1071 //
1072 for(Int_t i=0;i<fNVar;i++){
1073 if(group[i]!=1)AliInfo(Form(" merging bins along dimension %i in groups of %i bins", i,group[i]));
1074 }
365c0ea4 1075 AliInfo(Form("This method was so far not implemented for AliCFGrid, but it is available for AliCFGridSparse"));
7411edfd 1076
1077}
1078//____________________________________________________________________
563113d0 1079void AliCFGrid::SumW2()
1080{
1081 //
1082 //set calculation of the squared sum of the weighted entries
1083 //
1084 if(!fSumW2){
1085 fErr2=new Float_t [fNDim];
1086 //init....
1087 for(Int_t iel=0;iel<fNDim;iel++){
1088 fErr2[iel]=fData[iel];
1089 }
1090 }
1091
1092 fSumW2=kTRUE;
1093}
563113d0 1094//____________________________________________________________________
1095void AliCFGrid::Copy(TObject& c) const
1096{
1097 //
1098 // copy function
1099 //
1100 AliCFGrid& target = (AliCFGrid &) c;
1101
563113d0 1102 target.fNentriesTot = fNentriesTot;
563113d0 1103 if (fNunfl)
1104 target.fNunfl = fNunfl;
1105 if (fNovfl)
1106 target.fNovfl = fNovfl;
563113d0 1107 if (fData)
1108 target.fData = fData;
1109 if (fErr2)
1110 target.fErr2 = fErr2;
1111
1112}
318f64b1 1113//____________________________________________________________________
1114void AliCFGrid::SetExcludeOffEntriesInProj(Bool_t in)
1115{
1116 //
1117 // require under/overflows in 'hidden dimensions' to be excluded
1118 // or included, when performing projections.
1119 // For AliCFGrid implementation, only option = kTRUE is available
1120
46308b66 1121 if(!in){
1122 AliInfo(Form("This option is not available for AliCFGrid, Under/Overflows in hidden dimensions are always excluded"));
318f64b1 1123 return;
1124 }
1125
1126 fExclOffEntriesInProj=in;
1127}
1128//____________________________________________________________________
1129Bool_t AliCFGrid::GetExcludeOffEntriesInProj( ) const
1130{
1131 //
1132 // return flag saying whether under/overflows are excluded in projections
1133 //
1134
1135 return fExclOffEntriesInProj;
1136}