8c6f68376d32e8d3073f1f09b3bed361cb415441
[u/mrichter/AliRoot.git] / CORRFW / AliCFEffGrid.cxx
1 /* $Id$ */
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 // AliCFEffGrid Class                                              //
19 // Class to handle efficiency grids                                   // 
20 //                                                                    //
21 // -- Author : S.Arcelli                                              //
22 //                                                                    //
23 //                                                                    //
24 //                                                                    //
25 //--------------------------------------------------------------------//
26 //
27 //
28 #include "TMath.h"
29 #include "AliLog.h"
30 #include "AliCFEffGrid.h"
31 #include "TH1D.h"
32 #include "TH2D.h"
33 #include "TH3D.h"
34
35 //____________________________________________________________________
36 ClassImp(AliCFEffGrid)
37
38 //____________________________________________________________________
39 AliCFEffGrid::AliCFEffGrid() : 
40   AliCFGridSparse(),
41   fContainer(0x0),
42   fSelNum(-1),
43   fSelDen(-1)
44 {
45   //
46   // default constructor
47   //
48 }
49
50 //____________________________________________________________________
51 AliCFEffGrid::AliCFEffGrid(const Char_t* name, const Char_t* title, const Int_t nVarIn, const Int_t * nBinIn, const Double_t *binLimitsIn) :  
52   AliCFGridSparse(name,title,nVarIn,nBinIn,binLimitsIn),
53   fContainer(0x0),
54   fSelNum(-1),
55   fSelDen(-1)
56 {
57   //
58   // ctor
59   //
60   SumW2();
61 }
62 //____________________________________________________________________
63 AliCFEffGrid::AliCFEffGrid(const Char_t* name, const Char_t* title, const AliCFContainer &c) :  
64   AliCFGridSparse(name,title,c.GetNVar(),c.GetNBins(),c.GetBinLimits()),
65   fContainer(NULL),
66   fSelNum(-1),
67   fSelDen(-1)
68 {
69   //
70   // main constructor
71   //
72   SumW2();
73   //assign the container;
74   fContainer=&c;
75 }
76 //____________________________________________________________________
77 AliCFEffGrid::AliCFEffGrid(const AliCFEffGrid& eff) : AliCFGridSparse(),
78   fContainer(0x0),
79   fSelNum(-1),
80   fSelDen(-1)
81 {
82   //
83   // copy constructor
84   //
85   ((AliCFEffGrid &)eff).Copy(*this);
86 }
87
88 //____________________________________________________________________
89 AliCFEffGrid::~AliCFEffGrid()
90 {
91   //
92   // destructor
93   //
94 }
95
96 //____________________________________________________________________
97 AliCFEffGrid &AliCFEffGrid::operator=(const AliCFEffGrid &eff)
98 {
99   //
100   // assigment operator
101   //
102   if (this != &eff)
103     ((AliCFEffGrid &) eff).Copy(*this);
104   return *this;
105
106 //____________________________________________________________________
107
108 void AliCFEffGrid::CalculateEfficiency(Int_t istep1,Int_t istep2, Option_t *option)
109 {
110   //
111   // Calculate the efficiency matrix and its error between selection
112   // Steps istep1 and istep2
113   //
114
115   fSelNum=istep1;
116   fSelDen=istep2;
117   AliCFVGrid *num=fContainer->GetGrid(fSelNum);
118   AliCFVGrid *den=fContainer->GetGrid(fSelDen);
119   num->SumW2();
120   den->SumW2();
121   this->SumW2();
122   this->Divide(num,den,1.,1.,option);
123
124   Int_t nEmptyBinsNum=0;
125   Int_t nEmptyBinsNumAndDen=0;
126   for(Int_t iel=0;iel<fNDim;iel++){
127     if(den->GetElement(iel)>0){
128       if(num->GetElement(iel)==0)nEmptyBinsNum++; //num==0,den!=0
129     }
130     else{
131       nEmptyBinsNumAndDen++;
132     }
133   }    
134   // Some monitoring printout:
135   AliInfo(Form("Efficiency calculated for steps %i and %i: %i empty bins in the numerator && !denominator and %i empty bins in numerator && denominator were found.",fSelNum,fSelDen,nEmptyBinsNumAndDen,nEmptyBinsNum));
136   AliInfo(Form("The correction map contains %i empty bins ",nEmptyBinsNum+nEmptyBinsNumAndDen));
137
138 //_____________________________________________________________________
139 Double_t AliCFEffGrid::GetAverage() const 
140 {
141   //
142   // Get the average efficiency 
143   //
144
145   Double_t val=0;
146   Double_t valnum=0;
147   Double_t valden=0;
148   for(Int_t i=0;i<fNDim;i++){
149     valnum+=fContainer->GetGrid(fSelNum)->GetElement(i);
150     valden+=fContainer->GetGrid(fSelDen)->GetElement(i);
151   }
152   if(valden>0)val=valnum/valden;
153   AliInfo(Form(" The Average Efficiency = %f ",val)); 
154
155   return val;
156
157 //_____________________________________________________________________
158 Double_t AliCFEffGrid::GetAverage(Double_t *varMin, Double_t* varMax ) const 
159 {
160   //
161   // Get ave efficiency in a range
162   //
163
164
165   Double_t val=0;
166   Int_t *indexMin = new Int_t[fNVar];
167   Int_t *indexMax = new Int_t[fNVar];
168   Int_t *index    = new Int_t[fNVar];
169   
170   //Find out the min and max bins
171   
172   for(Int_t i=0;i<fNVar;i++){
173     Double_t xmin=varMin[i]; // the min values  
174     Double_t xmax=varMax[i]; // the max values  
175     Int_t nbins=fNVarBins[i]+1;
176     Double_t *bins=new Double_t[nbins];
177     for(Int_t ibin =0;ibin<nbins;ibin++){
178       bins[ibin] = fVarBinLimits[ibin+fOffset[i]];
179     }
180     indexMin[i] = TMath::BinarySearch(nbins,bins,xmin);
181     indexMax[i] = TMath::BinarySearch(nbins,bins,xmax);
182     if(xmax>=bins[nbins-1]){
183       indexMax[i]=indexMax[i]-1;
184     }  
185     delete [] bins;
186   }
187   
188   Double_t valnum=0;
189   Double_t valden=0;
190   for(Int_t i=0;i<fNDim;i++){
191     for (Int_t j=0;j<fNVar;j++)index[j]=GetBinIndex(j,i);
192     Bool_t isIn=kTRUE;
193     for (Int_t j=0;j<fNVar;j++){
194       if(!(index[j]>=indexMin[j] && index[j]<=indexMax[j]))isIn=kFALSE;   
195     }
196     if(isIn){
197       valnum+=fContainer->GetGrid(fSelNum)->GetElement(i);
198       valden+=fContainer->GetGrid(fSelDen)->GetElement(i);
199     }
200   } 
201   delete [] index;
202   delete [] indexMin;
203   delete [] indexMax;
204   if(valden>0)val=valnum/valden;
205   AliInfo(Form(" the Average Efficiency = %f ",val)); 
206   return val;
207
208 //____________________________________________________________________
209 void AliCFEffGrid::Copy(TObject& eff) const
210 {
211   //
212   // copy function
213   //
214   Copy(eff);
215   AliCFEffGrid& target = (AliCFEffGrid &) eff;
216   
217   target.fSelNum=fSelNum; 
218   target.fSelDen=fSelDen; 
219   if(fContainer)
220     target.fContainer=fContainer;
221 }
222 //___________________________________________________________________
223 TH1D *AliCFEffGrid::Project(Int_t ivar) const
224 {
225   //
226   // Make a 1D projection along variable ivar 
227   //
228  
229   TH1D *proj1D=0;
230   Int_t nbins =fNVarBins[ivar];
231   Double_t *bins = new Double_t[nbins+1];    
232   for(Int_t ibin =0;ibin<nbins+1;ibin++){
233     bins[ibin] = fVarBinLimits[ibin+fOffset[ivar]];
234   }
235   
236   char pname[30];
237   sprintf(pname,"%s%s%i%i%s%i",GetName(),"_SelStep",fSelNum,fSelDen,"_proj1D_var", ivar);
238   char htitle[30];
239   sprintf(htitle,"%s%s%i%i%s%i",GetName(),"_SelStep",fSelNum,fSelDen,"_proj1D_var", ivar);
240   
241   if(!proj1D){
242     proj1D =new TH1D(pname,htitle, nbins, bins);
243   }  
244   
245   //unset the use of axis range to be able to divide the histos
246   fContainer->GetGrid(fSelNum)->UseAxisRange(kFALSE);
247   fContainer->GetGrid(fSelDen)->UseAxisRange(kFALSE);
248
249   proj1D->Sumw2();
250   proj1D->Divide(fContainer->GetGrid(fSelNum)->Project(ivar),fContainer->GetGrid(fSelDen)->Project(ivar),1.,1.,"B");
251
252   fContainer->GetGrid(fSelNum)->UseAxisRange(kTRUE);
253   fContainer->GetGrid(fSelDen)->UseAxisRange(kTRUE);
254   proj1D->GetXaxis()->SetRange(
255                                ((AliCFGridSparse*)GetNum())->GetGrid()->GetAxis(ivar)->GetFirst(),
256                                ((AliCFGridSparse*)GetNum())->GetGrid()->GetAxis(ivar)->GetLast()
257                                );
258
259   delete [] bins; 
260   return proj1D;
261
262 //___________________________________________________________________
263 TH2D *AliCFEffGrid::Project(Int_t ivar1,Int_t ivar2) const
264 {
265   //
266   // Make a 2D projection along variable ivar1,ivar2 
267   //
268  
269   TH2D *proj2D=0;
270
271   Int_t nbins1 =fNVarBins[ivar1];
272   Double_t *bins1 = new Double_t[nbins1+1];    
273   Int_t nbins2 =fNVarBins[ivar2];
274   Double_t *bins2 = new Double_t[nbins2+1];    
275   for(Int_t ibin1 =0;ibin1<nbins1+1;ibin1++){
276     bins1[ibin1] = fVarBinLimits[ibin1+fOffset[ivar1]];
277   }
278   for(Int_t ibin2 =0;ibin2<nbins2+1;ibin2++){
279     bins2[ibin2] = fVarBinLimits[ibin2+fOffset[ivar2]];
280   }
281   
282   char pname[30];
283   sprintf(pname,"%s%s%i%i%s%i%i",GetName(),"_SelStep",fSelNum,fSelDen,"_proj2D_var", ivar1,ivar2);
284   char htitle[30];
285   sprintf(htitle,"%s%s%i%i%s%i%i",GetName(),"_SelStep",fSelNum,fSelDen,"_proj2D_var",ivar1,ivar2);
286   
287   if(!proj2D){
288     proj2D =new TH2D(pname,htitle, nbins1,bins1,nbins2,bins2);
289   }  
290   
291   //unset the use of axis range to be able to divide the histos
292   fContainer->GetGrid(fSelNum)->UseAxisRange(kFALSE);
293   fContainer->GetGrid(fSelDen)->UseAxisRange(kFALSE);
294   
295   proj2D->Sumw2();
296   proj2D->Divide(fContainer->GetGrid(fSelNum)->Project(ivar1,ivar2),fContainer->GetGrid(fSelDen)->Project(ivar1,ivar2),1.,1.,"B");
297   
298   fContainer->GetGrid(fSelNum)->UseAxisRange(kTRUE);
299   fContainer->GetGrid(fSelDen)->UseAxisRange(kTRUE);
300   proj2D->GetXaxis()->SetRange(
301                                ((AliCFGridSparse*)GetNum())->GetGrid()->GetAxis(ivar1)->GetFirst(),
302                                ((AliCFGridSparse*)GetNum())->GetGrid()->GetAxis(ivar1)->GetLast()
303                                );
304   proj2D->GetYaxis()->SetRange(
305                                ((AliCFGridSparse*)GetNum())->GetGrid()->GetAxis(ivar2)->GetFirst(),
306                                ((AliCFGridSparse*)GetNum())->GetGrid()->GetAxis(ivar2)->GetLast()
307                                );
308
309   delete [] bins1;
310   delete [] bins2; 
311   return proj2D;
312
313 //___________________________________________________________________
314 TH3D *AliCFEffGrid::Project(Int_t ivar1, Int_t ivar2, Int_t ivar3) const
315 {
316   //
317   // Make a 3D projection along variable ivar1,ivar2,ivar3 
318   //
319
320   TH3D *proj3D=0;
321
322   Int_t nbins1 =fNVarBins[ivar1];
323   Int_t nbins2 =fNVarBins[ivar2];
324   Int_t nbins3 =fNVarBins[ivar3];
325   
326   Double_t *bins1 = new Double_t[nbins1+1];         
327   Double_t *bins2 = new Double_t[nbins2+1];     
328   Double_t *bins3 = new Double_t[nbins3+1];     
329   
330   for(Int_t ibin =0;ibin<nbins1+1;ibin++){
331     bins1[ibin] = fVarBinLimits[ibin+fOffset[ivar1]];
332   }
333   for(Int_t ibin =0;ibin<nbins2+1;ibin++){
334     bins2[ibin] = fVarBinLimits[ibin+fOffset[ivar2]];
335   }
336   for(Int_t ibin =0;ibin<nbins3+1;ibin++){
337     bins3[ibin] = fVarBinLimits[ibin+fOffset[ivar3]];
338   }
339   
340   char pname[30];
341   sprintf(pname,"%s%s%i%i%s%i%i%i",GetName(),"_SelStep",fSelNum,fSelDen,"_proj3D_var",ivar1,ivar2,ivar3);
342   char htitle[30];
343   sprintf(htitle,"%s%s%i%i%s%i%i%i",GetName(),"_SelStep",fSelNum,fSelDen,"_proj3D_var",ivar1,ivar2,ivar3);
344    
345   
346   if(!proj3D){
347     proj3D =new TH3D(pname,htitle, nbins1, bins1,nbins2,bins2,nbins3,bins3);
348   }  
349   
350   //unset the use of axis range to be able to divide the histos
351   fContainer->GetGrid(fSelNum)->UseAxisRange(kFALSE);
352   fContainer->GetGrid(fSelDen)->UseAxisRange(kFALSE);
353
354   proj3D->Sumw2();
355   proj3D->Divide(fContainer->GetGrid(fSelNum)->Project(ivar1,ivar2,ivar3),fContainer->GetGrid(fSelDen)->Project(ivar1,ivar2,ivar3),1.,1.,"B");
356   
357   fContainer->GetGrid(fSelNum)->UseAxisRange(kTRUE);
358   fContainer->GetGrid(fSelDen)->UseAxisRange(kTRUE);
359
360   proj3D->GetXaxis()->SetRange(
361                                ((AliCFGridSparse*)GetNum())->GetGrid()->GetAxis(ivar1)->GetFirst(),
362                                ((AliCFGridSparse*)GetNum())->GetGrid()->GetAxis(ivar1)->GetLast()
363                                );
364   proj3D->GetYaxis()->SetRange(
365                                ((AliCFGridSparse*)GetNum())->GetGrid()->GetAxis(ivar2)->GetFirst(),
366                                ((AliCFGridSparse*)GetNum())->GetGrid()->GetAxis(ivar2)->GetLast()
367                                );
368   proj3D->GetZaxis()->SetRange(
369                                ((AliCFGridSparse*)GetNum())->GetGrid()->GetAxis(ivar3)->GetFirst(),
370                                ((AliCFGridSparse*)GetNum())->GetGrid()->GetAxis(ivar3)->GetLast()
371                                );
372
373   delete [] bins1;
374   delete [] bins2;
375   delete [] bins3; 
376   
377   return proj3D;
378
379