Added projection with axis range defined from bin numbers
authorrvernet <rvernet@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 2 Oct 2009 13:27:15 +0000 (13:27 +0000)
committerrvernet <rvernet@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 2 Oct 2009 13:27:15 +0000 (13:27 +0000)
CORRFW/AliCFContainer.cxx
CORRFW/AliCFContainer.h
CORRFW/AliCFEffGrid.cxx
CORRFW/AliCFEffGrid.h
CORRFW/AliCFGridSparse.cxx
CORRFW/AliCFGridSparse.h

index 7c18e54..ff47cbb 100644 (file)
@@ -182,10 +182,11 @@ TH3D *AliCFContainer::Project(Int_t ivar1, Int_t ivar2, Int_t ivar3, Int_t istep
   return fGrid[istep]->Project(ivar1,ivar2,ivar3);
 }
 //___________________________________________________________________
-TH1D *AliCFContainer::ShowSlice(Int_t ivar, Double_t *varMin, Double_t* varMax, Int_t istep) const
+TH1D *AliCFContainer::ShowSlice(Int_t ivar, const Double_t *varMin, const Double_t* varMax, Int_t istep, Bool_t useBins) const
 {
   //
   // Make a slice along variable ivar at selection level istep in range [varMin,varMax]
+  // If useBins=true, varMin and varMax are taken as bin numbers
   //
   if(istep >= fNStep || istep < 0){
     AliError("Non-existent selection step, return NULL");
@@ -195,13 +196,14 @@ TH1D *AliCFContainer::ShowSlice(Int_t ivar, Double_t *varMin, Double_t* varMax,
     AliError("Non-existent variable, return NULL");
     return 0x0;
   }
-  return (TH1D*)fGrid[istep]->Slice(ivar,varMin,varMax);
+  return (TH1D*)fGrid[istep]->Slice(ivar,varMin,varMax,useBins);
 }
 //___________________________________________________________________
-TH2D *AliCFContainer::ShowSlice(Int_t ivar1, Int_t ivar2, Double_t *varMin, Double_t* varMax, Int_t istep) const
+TH2D *AliCFContainer::ShowSlice(Int_t ivar1, Int_t ivar2, const Double_t *varMin, const Double_t* varMax, Int_t istep, Bool_t useBins) const
 {
   //
   // Make a slice along variables ivar1 and ivar2 at selection level istep in range [varMin,varMax]
+  // If useBins=true, varMin and varMax are taken as bin numbers
   //
   if(istep >= fNStep || istep < 0){
     AliError("Non-existent selection step, return NULL");
@@ -211,13 +213,14 @@ TH2D *AliCFContainer::ShowSlice(Int_t ivar1, Int_t ivar2, Double_t *varMin, Doub
     AliError("Non-existent variable, return NULL");
     return 0x0;
   }
-  return (TH2D*)fGrid[istep]->Slice(ivar1,ivar2,varMin,varMax);
+  return (TH2D*)fGrid[istep]->Slice(ivar1,ivar2,varMin,varMax,useBins);
 }
 //___________________________________________________________________
-TH3D *AliCFContainer::ShowSlice(Int_t ivar1, Int_t ivar2, Int_t ivar3, Double_t *varMin, Double_t* varMax, Int_t istep) const
+TH3D *AliCFContainer::ShowSlice(Int_t ivar1, Int_t ivar2, Int_t ivar3, const Double_t *varMin, const Double_t* varMax, Int_t istep, Bool_t useBins) const
 {
   //
   // Make a slice along variables ivar1, ivar2and ivar3 at selection level istep in range [varMin,varMax]
+  // If useBins=true, varMin and varMax are taken as bin numbers
   //
   if(istep >= fNStep || istep < 0){
     AliError("Non-existent selection step, return NULL");
@@ -227,30 +230,32 @@ TH3D *AliCFContainer::ShowSlice(Int_t ivar1, Int_t ivar2, Int_t ivar3, Double_t
     AliError("Non-existent variable, return NULL");
     return 0x0;
   }
-  return (TH3D*)fGrid[istep]->Slice(ivar1,ivar2,ivar3,varMin,varMax);
+  return (TH3D*)fGrid[istep]->Slice(ivar1,ivar2,ivar3,varMin,varMax,useBins);
 }
 //____________________________________________________________________
-AliCFContainer* AliCFContainer::MakeSlice(Int_t nVars, const Int_t* vars, const Double_t* varMin, const Double_t* varMax) const
+AliCFContainer* AliCFContainer::MakeSlice(Int_t nVars, const Int_t* vars, const Double_t* varMin, const Double_t* varMax, Bool_t useBins) const
 {
   //
   // Makes a slice along the "nVars" variables defined in the array "vars[nVars]" for all the container steps.
   // The ranges of ALL the container variables must be defined in the array varMin[GetNVar()] and varMax[GetNVar()]
   // The function returns a new container of nVars variables.
+  // If useBins=true, varMin and varMax are taken as bin numbers
   //
   Int_t* steps = new Int_t[fNStep];
   for (Int_t iStep=0;iStep<fNStep;iStep++) steps[iStep]=iStep;
-  AliCFContainer* out = MakeSlice(nVars,vars,varMin,varMax,fNStep,steps);
+  AliCFContainer* out = MakeSlice(nVars,vars,varMin,varMax,fNStep,steps,useBins);
   delete [] steps ;
   return out;
 }
 
 //____________________________________________________________________
-AliCFContainer* AliCFContainer::MakeSlice(Int_t nVars, const Int_t* vars, const Double_t* varMin, const Double_t* varMax, Int_t nSteps, const Int_t* steps) const
+AliCFContainer* AliCFContainer::MakeSlice(Int_t nVars, const Int_t* vars, const Double_t* varMin, const Double_t* varMax, Int_t nSteps, const Int_t* steps, Bool_t useBins) const
 {
   //
   // Makes a slice along the "nVars" variables defined in the array "vars[nVars]" for the given "nSteps" defined in "steps[nSteps]".
   // The ranges of ALL the container variables must be defined in the array varMin[GetNVar()] and varMax[GetNVar()]
   // The function returns a new container of nVars variables.
+  // If useBins=true, varMin and varMax are taken as bin numbers
   //
 
   if (nVars < 1 || nVars > GetNVar())   AliError("Bad number of dimensions required for the slice");
@@ -260,7 +265,7 @@ AliCFContainer* AliCFContainer::MakeSlice(Int_t nVars, const Int_t* vars, const
 
   // create the output grids
   AliCFGridSparse** grids = new AliCFGridSparse*[nSteps] ;
-  for (Int_t iStep=0; iStep<nSteps; iStep++) grids[iStep] = fGrid[steps[iStep]]->Project(nVars,vars,varMin,varMax);
+  for (Int_t iStep=0; iStep<nSteps; iStep++) grids[iStep] = fGrid[steps[iStep]]->Project(nVars,vars,varMin,varMax,useBins);
   
   TAxis ** axis = new TAxis*[nVars];
   for (Int_t iVar=0; iVar<nVars; iVar++) axis[iVar] = ((AliCFGridSparse*)grids[0])->GetGrid()->GetAxis(iVar); //same axis for every grid
@@ -286,7 +291,7 @@ AliCFContainer* AliCFContainer::MakeSlice(Int_t nVars, const Int_t* vars, const
   //set grid for the given steps
   for (Int_t iStep=0; iStep<nSteps; iStep++) out->SetGrid(iStep,grids[iStep]);
 
-  delete bins;
+  delete [] bins;
   delete [] axis ;
   return out;
 }
index 83fe538..aa79143 100644 (file)
@@ -82,11 +82,13 @@ class AliCFContainer : public AliCFFrame
   virtual TH1D* Project( Int_t ivar, Int_t istep) const;
   virtual TH2D* Project( Int_t ivar1, Int_t ivar2, Int_t istep) const;
   virtual TH3D* Project( Int_t ivar1, Int_t ivar2,Int_t ivar3, Int_t istep) const;
-  virtual TH1D* ShowSlice( Int_t ivar, Double_t *varMin, Double_t *varMax, Int_t istep) const;
-  virtual TH2D* ShowSlice( Int_t ivar1, Int_t ivar2, Double_t *varMin, Double_t *varMax, Int_t istep) const;
-  virtual TH3D* ShowSlice( Int_t ivar1, Int_t ivar2, Int_t ivar3, Double_t *varMin, Double_t *varMax, Int_t istep) const;
-  virtual AliCFContainer*  MakeSlice (Int_t nVars, const Int_t* vars, const Double_t* varMin, const Double_t* varMax) const ;
-  virtual AliCFContainer*  MakeSlice (Int_t nVars, const Int_t* vars, const Double_t* varMin, const Double_t* varMax, Int_t nStep, const Int_t* steps) const ;
+
+  virtual TH1D* ShowSlice(Int_t ivar, const Double_t *varMin, const Double_t *varMax, Int_t istep, Bool_t useBins=0) const ;
+  virtual TH2D* ShowSlice(Int_t ivar1, Int_t ivar2, const Double_t *varMin, const Double_t *varMax, Int_t istep, Bool_t useBins=0) const ;
+  virtual TH3D* ShowSlice(Int_t ivar1, Int_t ivar2, Int_t ivar3, const Double_t *varMin, const Double_t *varMax, Int_t istep, Bool_t useBins=0) const ;
+  virtual AliCFContainer* MakeSlice(Int_t nVars, const Int_t* vars, const Double_t* varMin, const Double_t* varMax, Bool_t useBins=0) const ;
+  virtual AliCFContainer* MakeSlice(Int_t nVars, const Int_t* vars, const Double_t* varMin, const Double_t* varMax, Int_t nStep, const Int_t* steps, Bool_t useBins=0) const ;
+
   virtual void  SetRangeUser(Int_t ivar, Double_t varMin, Double_t varMax, Int_t istep) ;
   virtual void  SetRangeUser(Double_t* varMin, Double_t* varMax, Int_t istep) ;
   virtual void  SetGrid(Int_t step, AliCFGridSparse* grid) {fGrid[step]=grid;}
index dfd7ae6..c187757 100644 (file)
@@ -296,14 +296,15 @@ TH3D *AliCFEffGrid::Project(Int_t ivar1, Int_t ivar2, Int_t ivar3) const
   return h;
 } 
 //___________________________________________________________________
-AliCFEffGrid* AliCFEffGrid::MakeSlice(Int_t nVars, const Int_t* vars, const Double_t* varMin, const Double_t* varMax, Int_t numStep, Int_t denStep) const {
+AliCFEffGrid* AliCFEffGrid::MakeSlice(Int_t nVars, const Int_t* vars, const Double_t* varMin, const Double_t* varMax, Int_t numStep, Int_t denStep, Bool_t useBins) const {
   //
   // Makes a slice along the "nVars" variables defined in the array "vars[nVars]" for all the container steps.
   // The ranges of ALL the container variables must be defined in the array varMin[fNVar] and varMax[fNVar].
-  // This function returns the effiency relative to this new 'sliced' container, between steps defined in numStep and denStep
+  // This function returns the efficiency relative to this new 'sliced' container, between steps defined in numStep and denStep
+  // If useBins=true, varMin and varMax are taken as bin numbers
   //
   
-  AliCFContainer* cont = fContainer->MakeSlice(nVars,vars,varMin,varMax);
+  AliCFContainer* cont = fContainer->MakeSlice(nVars,vars,varMin,varMax,useBins);
   AliCFEffGrid  * eff  = new AliCFEffGrid(Form("%s_sliced",GetName()), Form("%s_sliced",GetTitle()), *cont);
   eff->CalculateEfficiency(numStep,denStep);
   return eff;
index 7bc2d65..7030b73 100644 (file)
@@ -33,8 +33,8 @@ class AliCFEffGrid : public AliCFGridSparse
   virtual TH1D* Project( Int_t ivar) const;
   virtual TH2D* Project( Int_t ivar1, Int_t ivar2) const;
   virtual TH3D* Project( Int_t ivar1, Int_t ivar2,Int_t ivar3) const;
-  virtual AliCFEffGrid* Project(Int_t,const Int_t*, const Double_t*, const Double_t*) const {AliWarning("Function not to be used"); return 0x0;}
-  virtual AliCFEffGrid* MakeSlice(Int_t nVars, const Int_t* vars, const Double_t* varMin, const Double_t* varMax, Int_t numStep, Int_t denStep) const;
+  virtual AliCFEffGrid* Project(Int_t,const Int_t*, const Double_t*, const Double_t*,Bool_t) const {AliWarning("Function not to be used"); return 0x0;}
+  virtual AliCFEffGrid* MakeSlice(Int_t nVars, const Int_t* vars, const Double_t* varMin, const Double_t* varMax, Int_t numStep, Int_t denStep, Bool_t useBins=0) const;
 
   //Efficiency calculation
   virtual void  CalculateEfficiency(Int_t istep1, Int_t istep2, Option_t *option ="B" /*binomial*/);
@@ -42,10 +42,9 @@ class AliCFEffGrid : public AliCFGridSparse
   virtual AliCFGridSparse*  GetDen() const {return fContainer->GetGrid(fSelDen);};
   virtual void  SetContainer(const AliCFContainer &c) {fContainer=&c;};
 
-/*   //basic operations */
-/*   virtual void Copy(TObject& eff) const; */
-  
+  //basic operations
+  /*   virtual void Copy(TObject& eff) const; */
+
  private:
   const AliCFContainer *fContainer; //pointer to the input AliContainer
   Int_t fSelNum;                    //numerator selection step
index 3641f66..d286149 100755 (executable)
@@ -161,11 +161,12 @@ TH3D *AliCFGridSparse::Project(Int_t ivar1, Int_t ivar2, Int_t ivar3) const
 }
 
 //___________________________________________________________________
-AliCFGridSparse* AliCFGridSparse::Project(Int_t nVars, const Int_t* vars, const Double_t* varMin, const Double_t* varMax) const
+AliCFGridSparse* AliCFGridSparse::Project(Int_t nVars, const Int_t* vars, const Double_t* varMin, const Double_t* varMax, Bool_t useBins) const
 {
   //
   // projects the grid on the nVars dimensions defined in vars.
   // axis ranges can be defined in arrays varMin, varMax
+  // If useBins=true, varMin and varMax are taken as bin numbers
   //
 
   // binning for new grid
@@ -181,7 +182,8 @@ AliCFGridSparse* AliCFGridSparse::Project(Int_t nVars, const Int_t* vars, const
   THnSparse* clone = ((THnSparse*)fData->Clone());
   if (varMin && varMax) {
     for (Int_t iAxis=0; iAxis<GetNVar(); iAxis++) {
-      clone->GetAxis(iAxis)->SetRangeUser(varMin[iAxis],varMax[iAxis]);
+      if (useBins)  clone->GetAxis(iAxis)->SetRange((Int_t)varMin[iAxis],(Int_t)varMax[iAxis]);
+      else          clone->GetAxis(iAxis)->SetRangeUser(varMin[iAxis],varMax[iAxis]);
     }
   }
   else AliInfo("Keeping same axis ranges");
@@ -764,49 +766,52 @@ void AliCFGridSparse::Copy(TObject& c) const
 }
 
 //____________________________________________________________________
-TH1D* AliCFGridSparse::Slice(Int_t iVar, const Double_t *varMin, const Double_t *varMax) const
+TH1D* AliCFGridSparse::Slice(Int_t iVar, const Double_t *varMin, const Double_t *varMax, Bool_t useBins) const
 {
   //
   // return a slice (1D-projection) on variable iVar while axis ranges are defined with varMin,varMax
   // arrays varMin and varMax contain the min and max values of each variable.
   // therefore varMin and varMax must have their dimensions equal to GetNVar()
-  //
+  // If useBins=true, varMin and varMax are taken as bin numbers
   
   THnSparse* clone = (THnSparse*)fData->Clone();
   for (Int_t iAxis=0; iAxis<GetNVar(); iAxis++) {
-    clone->GetAxis(iAxis)->SetRangeUser(varMin[iAxis],varMax[iAxis]);
+    if (useBins)  clone->GetAxis(iAxis)->SetRange((Int_t)varMin[iAxis],(Int_t)varMax[iAxis]);
+    else          clone->GetAxis(iAxis)->SetRangeUser(varMin[iAxis],varMax[iAxis]);
   }
   return clone->Projection(iVar);
 }
 
 //____________________________________________________________________
-TH2D* AliCFGridSparse::Slice(Int_t iVar1, Int_t iVar2, const Double_t *varMin, const Double_t *varMax) const
+TH2D* AliCFGridSparse::Slice(Int_t iVar1, Int_t iVar2, const Double_t *varMin, const Double_t *varMax, Bool_t useBins) const
 {
   //
   // return a slice (2D-projection) on variables iVar1 and iVar2 while axis ranges are defined with varMin,varMax
   // arrays varMin and varMax contain the min and max values of each variable.
   // therefore varMin and varMax must have their dimensions equal to GetNVar()
-  //
+  // If useBins=true, varMin and varMax are taken as bin numbers
   
   THnSparse* clone = (THnSparse*)fData->Clone();
   for (Int_t iAxis=0; iAxis<GetNVar(); iAxis++) {
-    clone->GetAxis(iAxis)->SetRangeUser(varMin[iAxis],varMax[iAxis]);
+    if (useBins)  clone->GetAxis(iAxis)->SetRange((Int_t)varMin[iAxis],(Int_t)varMax[iAxis]);
+    else          clone->GetAxis(iAxis)->SetRangeUser(varMin[iAxis],varMax[iAxis]);
   }
   return clone->Projection(iVar1,iVar2);
 }
 
 //____________________________________________________________________
-TH3D* AliCFGridSparse::Slice(Int_t iVar1, Int_t iVar2, Int_t iVar3, const Double_t *varMin, const Double_t *varMax) const
+TH3D* AliCFGridSparse::Slice(Int_t iVar1, Int_t iVar2, Int_t iVar3, const Double_t *varMin, const Double_t *varMax, Bool_t useBins) const
 {
   //
   // return a slice (3D-projection) on variables iVar1, iVar2 and iVar3 while axis ranges are defined with varMin,varMax
   // arrays varMin and varMax contain the min and max values of each variable.
   // therefore varMin and varMax must have their dimensions equal to GetNVar()
-  //
+  // If useBins=true, varMin and varMax are taken as bin numbers
 
   THnSparse* clone = (THnSparse*)fData->Clone();
   for (Int_t iAxis=0; iAxis<GetNVar(); iAxis++) {
-    clone->GetAxis(iAxis)->SetRangeUser(varMin[iAxis],varMax[iAxis]);
+    if (useBins)  clone->GetAxis(iAxis)->SetRange((Int_t)varMin[iAxis],(Int_t)varMax[iAxis]);
+    else          clone->GetAxis(iAxis)->SetRangeUser(varMin[iAxis],varMax[iAxis]);
   }
   return clone->Projection(iVar1,iVar2,iVar3);
 }
index ec19e7e..1ebd098 100755 (executable)
@@ -72,10 +72,10 @@ class AliCFGridSparse : public AliCFFrame
   virtual TH1D*            Project( Int_t ivar) const;
   virtual TH2D*            Project( Int_t ivar1, Int_t ivar2) const;
   virtual TH3D*            Project( Int_t ivar1, Int_t ivar2,Int_t ivar3) const;
-  virtual AliCFGridSparse* Project(Int_t nVars, const Int_t* vars, const Double_t* varMin, const Double_t* varMax) const ;
-  virtual TH1D*            Slice(Int_t ivar, const Double_t* varMin, const Double_t* varMax) const ; 
-  virtual TH2D*            Slice(Int_t ivar1, Int_t ivar2, const Double_t *varMin, const Double_t *varMax) const ;
-  virtual TH3D*            Slice(Int_t ivar1, Int_t ivar2, Int_t ivar3, const Double_t *varMin, const Double_t *varMax) const ;
+  virtual AliCFGridSparse* Project(Int_t nVars, const Int_t* vars, const Double_t* varMin, const Double_t* varMax, Bool_t useBins=0) const ;
+  virtual TH1D*            Slice(Int_t ivar, const Double_t* varMin, const Double_t* varMax, Bool_t useBins=0) const ; 
+  virtual TH2D*            Slice(Int_t ivar1, Int_t ivar2, const Double_t *varMin, const Double_t *varMax, Bool_t useBins=0) const ;
+  virtual TH3D*            Slice(Int_t ivar1, Int_t ivar2, Int_t ivar3, const Double_t *varMin, const Double_t *varMax, Bool_t useBins=0) const ;
   virtual void             SetRangeUser(Int_t iVar, Double_t varMin, Double_t varMax) ;
   virtual void             SetRangeUser(const Double_t* varMin, const Double_t* varMax) ;
   virtual void             UseAxisRange(Bool_t b) const ;