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