[u/mrichter/AliRoot.git] / CORRFW / AliCFGridSparse.cxx

/**************************************************************************
 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
 *                                                                        *
 * Author: The ALICE Off-line Project.                                    *
 * Contributors are mentioned in the code where appropriate.              *
 *                                                                        *
 * Permission to use, copy, modify and distribute this software and its   *
 * documentation strictly for non-commercial purposes is hereby granted   *
 * without fee, provided that the above copyright notice appears in all   *
 * copies and that both the copyright notice and this permission notice   *
 * appear in the supporting documentation. The authors make no claims     *
 * about the suitability of this software for any purpose. It is          *
 * provided "as is" without express or implied warranty.                  *
 **************************************************************************/
//--------------------------------------------------------------------//
//                                                                    //
// AliCFGridSparse Class                                              //
// Class to accumulate data on an N-dimensional grid, to be used      //
// as input to get corrections for Reconstruction & Trigger efficiency// 
// Based on root THnSparse                                            //
// -- Author : S.Arcelli                                              //
// Still to be done:                                                  //
// --Interpolate among bins in a range                                // 
//--------------------------------------------------------------------//
//
//
#include "AliCFGridSparse.h"
#include "THnSparse.h"
#include "AliLog.h"
#include "TMath.h"
#include "TROOT.h"
#include "TH1D.h"
#include "TH2D.h"
#include "TH3D.h"
#include "TAxis.h"

//____________________________________________________________________
ClassImp(AliCFGridSparse)

//____________________________________________________________________
AliCFGridSparse::AliCFGridSparse() : 
  AliCFVGrid(),
  fData(0x0)
{
  // default constructor
}
//____________________________________________________________________
AliCFGridSparse::AliCFGridSparse(const Char_t* name, const Char_t* title) : 
  AliCFVGrid(name,title),
  fData(0x0)
{
  // default constructor
}
//____________________________________________________________________
AliCFGridSparse::AliCFGridSparse(const Char_t* name, const Char_t* title, const Int_t nVarIn, const Int_t * nBinIn, const Double_t *binLimitsIn) :  
  AliCFVGrid(name,title,nVarIn,nBinIn,binLimitsIn),
  fData(0x0)
{
  //
  // main constructor
  //

  fData=new THnSparseF(name,title,fNVar,fNVarBins);
  
  if(binLimitsIn){
    for(Int_t ivar=0;ivar<fNVar;ivar++){    
      Int_t nbins=fNVarBins[ivar]+1;
      Double_t *array= new Double_t[nbins];
      for(Int_t i=0;i<nbins;i++){
	array[i]=fVarBinLimits[fOffset[ivar]+i];
      } 
      fData->SetBinEdges(ivar, array);
      delete array;
    }
  }
}
//____________________________________________________________________
AliCFGridSparse::AliCFGridSparse(const AliCFGridSparse& c) : 
  AliCFVGrid(c),
  fData(c.fData)
{
  //
  // copy constructor
  //
  ((AliCFGridSparse &)c).Copy(*this);
}

//____________________________________________________________________
AliCFGridSparse::~AliCFGridSparse()
{
  //
  // destructor
  //
  if(fData) delete fData;
}

//____________________________________________________________________
AliCFGridSparse &AliCFGridSparse::operator=(const AliCFGridSparse &c)
{
  //
  // assigment operator
  //
  if (this != &c)
    ((AliCFGridSparse &) c).Copy(*this);
  return *this;
} 

//____________________________________________________________________
void AliCFGridSparse::SetBinLimits(Int_t ivar, Double_t *array)
{
  //
  // setting the arrays containing the bin limits 
  //
  fData->SetBinEdges(ivar, array);
  //then fill the appropriate array in ALICFFrame, to be able to use 
  //the getter, in case....
  Int_t nbins=fNVarBins[ivar]+1;
  for(Int_t i=0;i<nbins;i++){
    fVarBinLimits[fOffset[ivar]+i] =array[i];
  } 
} 

//____________________________________________________________________
void AliCFGridSparse::Fill(Double_t *var, Double_t weight)
{
  //
  // Fill the grid,
  // given a set of values of the input variable, 
  // with weight (by default w=1)
  //
  fData->Fill(var,weight);
}

//___________________________________________________________________
TH1D *AliCFGridSparse::Project(Int_t ivar) const
{
  //
  // Make a 1D projection along variable ivar 
  //

  //exclude overflows in hidden dimesniosn by default (used in projections)
  if(fExclOffEntriesInProj){
    for(Int_t i=0;i<fNVar;i++){
      if(i!=ivar)
	fData->GetAxis(i)->SetBit(TAxis::kAxisRange);
    }
  }
  TH1D *hist=fData->Projection(ivar);
  return hist;
}
//___________________________________________________________________
TH2D *AliCFGridSparse::Project(Int_t ivar1, Int_t ivar2) const
{
  //
  // Make a 2D projection along variables ivar1 & ivar2 
  //

  //exclude overflows by default (used in projections)
  if(fExclOffEntriesInProj){
    for(Int_t i=0;i<fNVar;i++){
      if(i!=ivar1 && i!=ivar2)
	fData->GetAxis(i)->SetBit(TAxis::kAxisRange);
    }
  }
  TH2D *hist=fData->Projection(ivar2,ivar1); //notice inverted axis (THnSparse uses TH3 2d-projection convention...)
  return hist;

}
//___________________________________________________________________
TH3D *AliCFGridSparse::Project(Int_t ivar1, Int_t ivar2, Int_t ivar3) const
{
  //
  // Make a 3D projection along variables ivar1 & ivar2 & ivar3 
  //
  //exclude overflows by default (used in projections)
  if(fExclOffEntriesInProj){
    for(Int_t i=0;i<fNVar;i++){
      if(i!=ivar1 && i!=ivar2 && i!=ivar3)
	fData->GetAxis(i)->SetBit(TAxis::kAxisRange);
    }
  }

  TH3D *hist=fData->Projection(ivar1,ivar2,ivar3); 
  return hist;

}

//____________________________________________________________________
Float_t AliCFGridSparse::GetOverFlows(Int_t ivar) const
{
  //
  // Returns exclusive overflows in variable ivar
  //
  Int_t* bin = new Int_t[fNVar];
  memset(bin, 0, sizeof(Int_t) * fNVar);
  Float_t ovfl=0.;
  for (Long64_t i = 0; i < fData->GetNbins(); ++i) {
    Double_t v = fData->GetBinContent(i, bin);
    Bool_t add=kTRUE;
    for(Int_t j=0;j<fNVar;j++){
      if(ivar==j)continue;
      if((bin[j]==0) || (bin[j]==fNVarBins[j]+1))add=kFALSE;
    }
    if(bin[ivar]==fNVarBins[ivar]+1 && add) ovfl+=v;
  }

  delete[] bin;
  return ovfl;
}

//____________________________________________________________________
Float_t AliCFGridSparse::GetUnderFlows(Int_t ivar) const
{
  //
  // Returns exclusive overflows in variable ivar
  //
  Int_t* bin = new Int_t[fNVar];
  memset(bin, 0, sizeof(Int_t) * fNVar);
  Float_t unfl=0.;
  for (Long64_t i = 0; i < fData->GetNbins(); ++i) {
    Double_t v = fData->GetBinContent(i, bin);
    Bool_t add=kTRUE;
    for(Int_t j=0;j<fNVar;j++){
      if(ivar==j)continue;
      if((bin[j]==0) || (bin[j]==fNVarBins[j]+1))add=kFALSE;
    }
    if(bin[ivar]==0 && add) unfl+=v;
  }

  delete[] bin;
  return unfl;
}


//____________________________________________________________________
Float_t AliCFGridSparse::GetEntries() const
{
  //
  // total entries (including overflows and underflows)
  //

  return fData->GetEntries();
}

//____________________________________________________________________
Float_t AliCFGridSparse::GetElement(Int_t index) const
{
  //
  // Returns content of grid element index according to the
  // linear indexing in AliCFFrame
  //
  Int_t *bin = new Int_t[fNVar];
  GetBinIndex(index, bin);
  for(Int_t i=0;i<fNVar;i++)fIndex[i]=bin[i]+1; //consistency with AliCFGrid
  Float_t val=GetElement(fIndex);
  delete [] bin;
  return val; 
  
}
//____________________________________________________________________
Float_t AliCFGridSparse::GetElement(Int_t *bin) const
{
  //
  // Get the content in a bin corresponding to a set of bin indexes
  //
  return fData->GetBinContent(bin);

}  
//____________________________________________________________________
Float_t AliCFGridSparse::GetElement(Double_t *var) const
{
  //
  // Get the content in a bin corresponding to a set of input variables
  //

  Long_t index=fData->GetBin(var,kFALSE); //this is the THnSparse index (do not allocate new cells if content is empty)
  if(index<0){
    return 0.;
  }else{
    return fData->GetBinContent(index);
  }
} 

//____________________________________________________________________
Float_t AliCFGridSparse::GetElementError(Int_t index) const
{
  //
  // Returns the error on the content of a bin according to a linear
  // indexing in AliCFFrame
  //

  Int_t *bin = new Int_t[fNVar];
  GetBinIndex(index, bin);
  for(Int_t i=0;i<fNVar;i++)fIndex[i]=bin[i]+1; //consistency with AliCFGrid
  Float_t val=GetElementError(fIndex);
  delete [] bin;
  return val;

}
//____________________________________________________________________
Float_t AliCFGridSparse::GetElementError(Int_t *bin) const
{
 //
  // Get the error in a bin corresponding to a set of bin indexes
  //
  return fData->GetBinError(bin);

}  
//____________________________________________________________________
Float_t AliCFGridSparse::GetElementError(Double_t *var) const
{
  //
  // Get the error in a bin corresponding to a set of input variables
  //

  Long_t index=fData->GetBin(var,kFALSE); //this is the THnSparse index (do not allocate new cells if content is empy)
  if(index<0){
    return 0.;
  }else{
    return fData->GetBinError(index);
  }
} 


//____________________________________________________________________
void AliCFGridSparse::SetElement(Int_t index, Float_t val)
{
  //
  // Sets grid element iel to val (linear indexing) in AliCFFrame
  //
  Int_t *bin = new Int_t[fNVar];
  GetBinIndex(index, bin);
  for(Int_t i=0;i<fNVar;i++)fIndex[i]=bin[i]+1;
  SetElement(fIndex,val);
  delete [] bin;
}
//____________________________________________________________________
void AliCFGridSparse::SetElement(Int_t *bin, Float_t val)
{
  //
  // Sets grid element of bin indeces bin to val
  //
  fData->SetBinContent(bin,val);
}
//____________________________________________________________________
void AliCFGridSparse::SetElement(Double_t *var, Float_t val) 
{
  //
  // Set the content in a bin to value val corresponding to a set of input variables
  //
  Long_t index=fData->GetBin(var); //THnSparse index: allocate the cell
  Int_t *bin = new Int_t[fNVar];
  fData->GetBinContent(index,bin); //trick to access the array of bins
  fData->SetBinContent(bin,val);
  delete [] bin;

}

//____________________________________________________________________
void AliCFGridSparse::SetElementError(Int_t index, Float_t val)
{
  //
  // Sets grid element iel error to val (linear indexing) in AliCFFrame
  //
  Int_t *bin = new Int_t[fNVar];
  GetBinIndex(index, bin);
  for(Int_t i=0;i<fNVar;i++)fIndex[i]=bin[i]+1;
  SetElementError(fIndex,val);
  delete [] bin;
}
//____________________________________________________________________
void AliCFGridSparse::SetElementError(Int_t *bin, Float_t val)
{
  //
  // Sets grid element error of bin indeces bin to val
  //
  fData->SetBinError(bin,val);
}
//____________________________________________________________________
void AliCFGridSparse::SetElementError(Double_t *var, Float_t val) 
{
  //
  // Set the error in a bin to value val corresponding to a set of input variables
  //
  Long_t index=fData->GetBin(var); //THnSparse index
  Int_t *bin = new Int_t[fNVar];
  fData->GetBinContent(index,bin); //trick to access the array of bins
  fData->SetBinError(bin,val);
  delete [] bin;
}

//____________________________________________________________________
void AliCFGridSparse::SumW2()
{
  //
  //set calculation of the squared sum of the weighted entries
  //
  if(!fSumW2){
    fData->CalculateErrors(kTRUE); 
  }

  fSumW2=kTRUE;
}

//____________________________________________________________________
void AliCFGridSparse::Add(AliCFVGrid* aGrid, Double_t c)
{
  //
  //add aGrid to the current one
  //


  if(aGrid->GetNVar()!=fNVar){
    AliInfo("Different number of variables, cannot add the grids");
    return;
  } 
  if(aGrid->GetNDim()!=fNDim){
    AliInfo("Different number of dimensions, cannot add the grids!");
    return;
  } 
  
  if(!fSumW2  && aGrid->GetSumW2())SumW2();


  fData->Add(((AliCFGridSparse*)aGrid)->GetGrid(),c);

}

//____________________________________________________________________
void AliCFGridSparse::Add(AliCFVGrid* aGrid1, AliCFVGrid* aGrid2, Double_t c1,Double_t c2)
{
  //
  //Add aGrid1 and aGrid2 and deposit the result into the current one
  //

  if(fNVar!=aGrid1->GetNVar()|| fNVar!=aGrid2->GetNVar()){
    AliInfo("Different number of variables, cannot add the grids");
    return;
  } 
  if(fNDim!=aGrid1->GetNDim()|| fNDim!=aGrid2->GetNDim()){
    AliInfo("Different number of dimensions, cannot add the grids!");
    return;
  } 
  
  if(!fSumW2  && (aGrid1->GetSumW2() || aGrid2->GetSumW2()))SumW2();


  fData->Reset();
  fData->Add(((AliCFGridSparse*)aGrid1)->GetGrid(),c1);
  fData->Add(((AliCFGridSparse*)aGrid2)->GetGrid(),c2);

}

//____________________________________________________________________
void AliCFGridSparse::Multiply(AliCFVGrid* aGrid, Double_t c)
{
  //
  // Multiply aGrid to the current one
  //


  if(aGrid->GetNVar()!=fNVar){
    AliInfo("Different number of variables, cannot multiply the grids");
    return;
  } 
  if(aGrid->GetNDim()!=fNDim){
    AliInfo("Different number of dimensions, cannot multiply the grids!");
    return;
  } 
  
  if(!fSumW2  && aGrid->GetSumW2())SumW2();

  THnSparse *h= ((AliCFGridSparse*)aGrid)->GetGrid();

  fData->Multiply(h);
  fData->Scale(c);

}

//____________________________________________________________________
void AliCFGridSparse::Multiply(AliCFVGrid* aGrid1, AliCFVGrid* aGrid2, Double_t c1,Double_t c2)
{
  //
  //Multiply aGrid1 and aGrid2 and deposit the result into the current one
  //

  if(fNVar!=aGrid1->GetNVar()|| fNVar!=aGrid2->GetNVar()){
    AliInfo("Different number of variables, cannot multiply the grids");
    return;
  } 
  if(fNDim!=aGrid1->GetNDim()|| fNDim!=aGrid2->GetNDim()){
    AliInfo("Different number of dimensions, cannot multiply the grids!");
    return;
  } 
  
  if(!fSumW2  && (aGrid1->GetSumW2() || aGrid2->GetSumW2()))SumW2();


  fData->Reset();
  THnSparse *h1= ((AliCFGridSparse*)aGrid1)->GetGrid();
  THnSparse *h2= ((AliCFGridSparse*)aGrid2)->GetGrid();
  h2->Multiply(h1);
  h2->Scale(c1*c2);
  fData->Add(h2);
}


//____________________________________________________________________
void AliCFGridSparse::Divide(AliCFVGrid* aGrid, Double_t c)
{
  //
  // Divide aGrid to the current one
  //


  if(aGrid->GetNVar()!=fNVar){
    AliInfo("Different number of variables, cannot divide the grids");
    return;
  } 
  if(aGrid->GetNDim()!=fNDim){
    AliInfo("Different number of dimensions, cannot divide the grids!");
    return;
  } 
  
  if(!fSumW2  && aGrid->GetSumW2())SumW2();

  THnSparse *h= ((AliCFGridSparse*)aGrid)->GetGrid();

  fData->Divide(h);
  fData->Scale(c);

}

//____________________________________________________________________
void AliCFGridSparse::Divide(AliCFVGrid* aGrid1, AliCFVGrid* aGrid2, Double_t c1,Double_t c2, Option_t *option)
{
  //
  //Divide aGrid1 and aGrid2 and deposit the result into the current one
  //bynomial errors are supported
  //

  if(fNVar!=aGrid1->GetNVar()|| fNVar!=aGrid2->GetNVar()){
    AliInfo("Different number of variables, cannot divide the grids");
    return;
  } 
  if(fNDim!=aGrid1->GetNDim()|| fNDim!=aGrid2->GetNDim()){
    AliInfo("Different number of dimensions, cannot divide the grids!");
    return;
  } 
  
  if(!fSumW2  && (aGrid1->GetSumW2() || aGrid2->GetSumW2()))SumW2();


  THnSparse *h1= ((AliCFGridSparse*)aGrid1)->GetGrid();
  THnSparse *h2= ((AliCFGridSparse*)aGrid2)->GetGrid();
  fData->Divide(h1,h2,c1,c2,option);
}


//____________________________________________________________________
void AliCFGridSparse::Rebin(const Int_t* group)
{
  //
  // rebin the grid according to Rebin() as in THnSparse
  // Please notice that the original number of bins on
  // a given axis has to be divisible by the rebin group.
  //

  for(Int_t i=0;i<fNVar;i++){
    if(group[i]!=1)AliInfo(Form(" merging bins along dimension %i in groups of %i bins", i,group[i]));
  }

  THnSparse *rebinned =fData->Rebin(group);
  fData->Reset();
  fData = rebinned;

  //redefine the needed stuff

  Int_t ndimTot=1;
  Int_t nbinTot=0;

  //number of bins in each dimension, auxiliary variables

  for(Int_t ivar=0;ivar<fNVar;ivar++){
    Int_t nbins = fData->GetAxis(ivar)->GetNbins();
    fNVarBins[ivar]=nbins;
    ndimTot*=fNVarBins[ivar];
    nbinTot+=(fNVarBins[ivar]+1);
    Int_t offset=0;
    for(Int_t i =0;i<ivar;i++)offset+=(fNVarBins[i]+1);      
    fOffset[ivar]=offset;
    Int_t prod=1;
    for(Int_t i=0;i<ivar;i++)prod*=fNVarBins[i];
    fProduct[ivar]=prod;
  }

  fNDim=ndimTot;

  //now the array of bin limits

  delete fVarBinLimits;
  fNVarBinLimits=nbinTot;
  fVarBinLimits=new Double_t[fNVarBinLimits];

  for(Int_t ivar=0;ivar<fNVar;ivar++){
    Double_t low = fData->GetAxis(ivar)->GetXmin();
    Double_t high = fData->GetAxis(ivar)->GetXmax();    
    const TArrayD *xbins = fData->GetAxis(ivar)->GetXbins();
    if (xbins->fN == 0){
      for(Int_t ibin=0;ibin<=fNVarBins[ivar];ibin++){
	fVarBinLimits[ibin+fOffset[ivar]] = low + ibin*(high-low)/((Double_t) fNVarBins[ivar]);
      }
    }
    else{
      
      for(Int_t ibin=0;ibin<=fNVarBins[ivar];ibin++) {
	fVarBinLimits[ibin+fOffset[ivar]] = xbins->At(ibin);
      }
    }
  }   
  
}
//____________________________________________________________________
void AliCFGridSparse::Copy(TObject& c) const
{
  //
  // copy function
  //
  AliCFGridSparse& target = (AliCFGridSparse &) c;

  if(fData)target.fData = fData;
}
Commit	Line	Data
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	//____________________________________________________________________
	38	ClassImp(AliCFGridSparse)
	39
	40	//____________________________________________________________________
	41	AliCFGridSparse::AliCFGridSparse() :
	42	AliCFVGrid(),
db6722a5	43	fData(0x0)
	44	{
	45	// default constructor
	46	}
	47	//____________________________________________________________________
	48	AliCFGridSparse::AliCFGridSparse(const Char_t* name, const Char_t* title) :
	49	AliCFVGrid(name,title),
db6722a5	50	fData(0x0)
	51	{
	52	// default constructor
	53	}
	54	//____________________________________________________________________
	55	AliCFGridSparse::AliCFGridSparse(const Char_t* name, const Char_t* title, const Int_t nVarIn, const Int_t * nBinIn, const Double_t *binLimitsIn) :
	56	AliCFVGrid(name,title,nVarIn,nBinIn,binLimitsIn),
db6722a5	57	fData(0x0)
	58	{
	59	//
	60	// main constructor
	61	//
	62
	63	fData=new THnSparseF(name,title,fNVar,fNVarBins);
	64
	65	if(binLimitsIn){
	66	for(Int_t ivar=0;ivar<fNVar;ivar++){
	67	Int_t nbins=fNVarBins[ivar]+1;
	68	Double_t *array= new Double_t[nbins];
	69	for(Int_t i=0;i<nbins;i++){
	70	array[i]=fVarBinLimits[fOffset[ivar]+i];
	71	}
	72	fData->SetBinEdges(ivar, array);
	73	delete array;
	74	}
	75	}
	76	}
	77	//____________________________________________________________________
	78	AliCFGridSparse::AliCFGridSparse(const AliCFGridSparse& c) :
	79	AliCFVGrid(c),
db6722a5	80	fData(c.fData)
	81	{
	82	//
	83	// copy constructor
	84	//
	85	((AliCFGridSparse &)c).Copy(*this);
	86	}
	87
	88	//____________________________________________________________________
	89	AliCFGridSparse::~AliCFGridSparse()
	90	{
	91	//
	92	// destructor
	93	//
	94	if(fData) delete fData;
	95	}
	96
	97	//____________________________________________________________________
	98	AliCFGridSparse &AliCFGridSparse::operator=(const AliCFGridSparse &c)
	99	{
	100	//
	101	// assigment operator
	102	//
	103	if (this != &c)
	104	((AliCFGridSparse &) c).Copy(*this);
	105	return *this;
	106	}
	107
	108	//____________________________________________________________________
	109	void AliCFGridSparse::SetBinLimits(Int_t ivar, Double_t *array)
	110	{
	111	//
	112	// setting the arrays containing the bin limits
	113	//
	114	fData->SetBinEdges(ivar, array);
	115	//then fill the appropriate array in ALICFFrame, to be able to use
	116	//the getter, in case....
	117	Int_t nbins=fNVarBins[ivar]+1;
	118	for(Int_t i=0;i<nbins;i++){
	119	fVarBinLimits[fOffset[ivar]+i] =array[i];
	120	}
	121	}
	122
	123	//____________________________________________________________________
	124	void AliCFGridSparse::Fill(Double_t *var, Double_t weight)
	125	{
	126	//
	127	// Fill the grid,
	128	// given a set of values of the input variable,
	129	// with weight (by default w=1)
	130	//
	131	fData->Fill(var,weight);
	132	}
	133
	134	//___________________________________________________________________
	135	TH1D *AliCFGridSparse::Project(Int_t ivar) const
	136	{
	137	//
	138	// Make a 1D projection along variable ivar
	139	//
	140
7411edfd	141	//exclude overflows in hidden dimesniosn by default (used in projections)
318f64b1	142	if(fExclOffEntriesInProj){
db6722a5	143	for(Int_t i=0;i<fNVar;i++){
7411edfd	144	if(i!=ivar)
7411edfd	145	fData->GetAxis(i)->SetBit(TAxis::kAxisRange);
db6722a5	146	}
	147	}
	148	TH1D *hist=fData->Projection(ivar);
db6722a5	149	return hist;
db6722a5	150	}
	151	//___________________________________________________________________
	152	TH2D *AliCFGridSparse::Project(Int_t ivar1, Int_t ivar2) const
	153	{
	154	//
	155	// Make a 2D projection along variables ivar1 & ivar2
	156	//
	157
	158	//exclude overflows by default (used in projections)
318f64b1	159	if(fExclOffEntriesInProj){
db6722a5	160	for(Int_t i=0;i<fNVar;i++){
7411edfd	161	if(i!=ivar1 && i!=ivar2)
7411edfd	162	fData->GetAxis(i)->SetBit(TAxis::kAxisRange);
db6722a5	163	}
db6722a5	164	}
318f64b1	165	TH2D *hist=fData->Projection(ivar2,ivar1); //notice inverted axis (THnSparse uses TH3 2d-projection convention...)
db6722a5	166	return hist;
	167
	168	}
	169	//___________________________________________________________________
	170	TH3D *AliCFGridSparse::Project(Int_t ivar1, Int_t ivar2, Int_t ivar3) const
	171	{
	172	//
	173	// Make a 3D projection along variables ivar1 & ivar2 & ivar3
	174	//
	175	//exclude overflows by default (used in projections)
318f64b1	176	if(fExclOffEntriesInProj){
db6722a5	177	for(Int_t i=0;i<fNVar;i++){
7411edfd	178	if(i!=ivar1 && i!=ivar2 && i!=ivar3)
7411edfd	179	fData->GetAxis(i)->SetBit(TAxis::kAxisRange);
db6722a5	180	}
	181	}
	182
	183	TH3D *hist=fData->Projection(ivar1,ivar2,ivar3);
db6722a5	184	return hist;
	185
	186	}
	187
	188	//____________________________________________________________________
	189	Float_t AliCFGridSparse::GetOverFlows(Int_t ivar) const
	190	{
	191	//
318f64b1	192	// Returns exclusive overflows in variable ivar
db6722a5	193	//
7411edfd	194	Int_t* bin = new Int_t[fNVar];
7411edfd	195	memset(bin, 0, sizeof(Int_t) * fNVar);
db6722a5	196	Float_t ovfl=0.;
	197	for (Long64_t i = 0; i < fData->GetNbins(); ++i) {
	198	Double_t v = fData->GetBinContent(i, bin);
	199	Bool_t add=kTRUE;
	200	for(Int_t j=0;j<fNVar;j++){
	201	if(ivar==j)continue;
	202	if((bin[j]==0) \|\| (bin[j]==fNVarBins[j]+1))add=kFALSE;
	203	}
	204	if(bin[ivar]==fNVarBins[ivar]+1 && add) ovfl+=v;
	205	}
	206
	207	delete[] bin;
	208	return ovfl;
	209	}
	210
	211	//____________________________________________________________________
	212	Float_t AliCFGridSparse::GetUnderFlows(Int_t ivar) const
	213	{
	214	//
318f64b1	215	// Returns exclusive overflows in variable ivar
db6722a5	216	//
7411edfd	217	Int_t* bin = new Int_t[fNVar];
7411edfd	218	memset(bin, 0, sizeof(Int_t) * fNVar);
db6722a5	219	Float_t unfl=0.;
	220	for (Long64_t i = 0; i < fData->GetNbins(); ++i) {
	221	Double_t v = fData->GetBinContent(i, bin);
	222	Bool_t add=kTRUE;
	223	for(Int_t j=0;j<fNVar;j++){
	224	if(ivar==j)continue;
	225	if((bin[j]==0) \|\| (bin[j]==fNVarBins[j]+1))add=kFALSE;
	226	}
	227	if(bin[ivar]==0 && add) unfl+=v;
	228	}
	229
	230	delete[] bin;
	231	return unfl;
	232	}
	233
	234
	235	//____________________________________________________________________
	236	Float_t AliCFGridSparse::GetEntries() const
	237	{
	238	//
	239	// total entries (including overflows and underflows)
	240	//
	241
	242	return fData->GetEntries();
	243	}
	244
	245	//____________________________________________________________________
	246	Float_t AliCFGridSparse::GetElement(Int_t index) const
	247	{
	248	//
	249	// Returns content of grid element index according to the
	250	// linear indexing in AliCFFrame
	251	//
	252	Int_t *bin = new Int_t[fNVar];
	253	GetBinIndex(index, bin);
	254	for(Int_t i=0;i<fNVar;i++)fIndex[i]=bin[i]+1; //consistency with AliCFGrid
	255	Float_t val=GetElement(fIndex);
	256	delete [] bin;
	257	return val;
	258
	259	}
	260	//____________________________________________________________________
	261	Float_t AliCFGridSparse::GetElement(Int_t *bin) const
	262	{
	263	//
	264	// Get the content in a bin corresponding to a set of bin indexes
	265	//
	266	return fData->GetBinContent(bin);
	267
	268	}
	269	//____________________________________________________________________
	270	Float_t AliCFGridSparse::GetElement(Double_t *var) const
	271	{
	272	//
	273	// Get the content in a bin corresponding to a set of input variables
	274	//
	275
	276	Long_t index=fData->GetBin(var,kFALSE); //this is the THnSparse index (do not allocate new cells if content is empty)
	277	if(index<0){
	278	return 0.;
	279	}else{
	280	return fData->GetBinContent(index);
	281	}
	282	}
283
284	//____________________________________________________________________
285	Float_t AliCFGridSparse::GetElementError(Int_t index) const
286	{
287	//
288	// Returns the error on the content of a bin according to a linear
289	// indexing in AliCFFrame
290	//
291
292	Int_t *bin = new Int_t[fNVar];
293	GetBinIndex(index, bin);
294	for(Int_t i=0;i<fNVar;i++)fIndex[i]=bin[i]+1; //consistency with AliCFGrid
295	Float_t val=GetElementError(fIndex);
296	delete [] bin;
297	return val;
298
299	}
300	//____________________________________________________________________
301	Float_t AliCFGridSparse::GetElementError(Int_t *bin) const
302	{
303	//
304	// Get the error in a bin corresponding to a set of bin indexes
305	//
306	return fData->GetBinError(bin);
307
308	}
309	//____________________________________________________________________
310	Float_t AliCFGridSparse::GetElementError(Double_t *var) const
311	{
312	//
313	// Get the error in a bin corresponding to a set of input variables
314	//
315
316	Long_t index=fData->GetBin(var,kFALSE); //this is the THnSparse index (do not allocate new cells if content is empy)
317	if(index<0){
318	return 0.;
319	}else{
320	return fData->GetBinError(index);
321	}
322	}
323
324
325	//____________________________________________________________________
326	void AliCFGridSparse::SetElement(Int_t index, Float_t val)
327	{
328	//
329	// Sets grid element iel to val (linear indexing) in AliCFFrame
330	//
331	Int_t *bin = new Int_t[fNVar];
332	GetBinIndex(index, bin);
333	for(Int_t i=0;i<fNVar;i++)fIndex[i]=bin[i]+1;
334	SetElement(fIndex,val);
335	delete [] bin;
336	}
337	//____________________________________________________________________
338	void AliCFGridSparse::SetElement(Int_t *bin, Float_t val)
339	{
340	//
341	// Sets grid element of bin indeces bin to val
342	//
343	fData->SetBinContent(bin,val);
344	}
345	//____________________________________________________________________
346	void AliCFGridSparse::SetElement(Double_t *var, Float_t val)
347	{
348	//
349	// Set the content in a bin to value val corresponding to a set of input variables
350	//
351	Long_t index=fData->GetBin(var); //THnSparse index: allocate the cell
352	Int_t *bin = new Int_t[fNVar];
353	fData->GetBinContent(index,bin); //trick to access the array of bins
354	fData->SetBinContent(bin,val);
355	delete [] bin;
356
357	}
358
359	//____________________________________________________________________
360	void AliCFGridSparse::SetElementError(Int_t index, Float_t val)
361	{
362	//
363	// Sets grid element iel error to val (linear indexing) in AliCFFrame
364	//
365	Int_t *bin = new Int_t[fNVar];
366	GetBinIndex(index, bin);
367	for(Int_t i=0;i<fNVar;i++)fIndex[i]=bin[i]+1;
368	SetElementError(fIndex,val);
369	delete [] bin;
370	}
371	//____________________________________________________________________
372	void AliCFGridSparse::SetElementError(Int_t *bin, Float_t val)
373	{
374	//
375	// Sets grid element error of bin indeces bin to val
376	//
377	fData->SetBinError(bin,val);
378	}
379	//____________________________________________________________________
380	void AliCFGridSparse::SetElementError(Double_t *var, Float_t val)
381	{
382	//
383	// Set the error in a bin to value val corresponding to a set of input variables
384	//
385	Long_t index=fData->GetBin(var); //THnSparse index
386	Int_t *bin = new Int_t[fNVar];
387	fData->GetBinContent(index,bin); //trick to access the array of bins
388	fData->SetBinError(bin,val);
389	delete [] bin;
390	}
391
392	//____________________________________________________________________
393	void AliCFGridSparse::SumW2()
394	{
395	//
396	//set calculation of the squared sum of the weighted entries
397	//
398	if(!fSumW2){
399	fData->CalculateErrors(kTRUE);
400	}
401
402	fSumW2=kTRUE;
403	}
404
405	//____________________________________________________________________
406	void AliCFGridSparse::Add(AliCFVGrid* aGrid, Double_t c)
407	{
408	//
409	//add aGrid to the current one
410	//
411
412
413	if(aGrid->GetNVar()!=fNVar){
414	AliInfo("Different number of variables, cannot add the grids");
415	return;
416	}
417	if(aGrid->GetNDim()!=fNDim){
418	AliInfo("Different number of dimensions, cannot add the grids!");
419	return;
420	}
421
422	if(!fSumW2 && aGrid->GetSumW2())SumW2();
423
424
425	fData->Add(((AliCFGridSparse*)aGrid)->GetGrid(),c);
426
427	}
428
429	//____________________________________________________________________
430	void AliCFGridSparse::Add(AliCFVGrid* aGrid1, AliCFVGrid* aGrid2, Double_t c1,Double_t c2)
431	{
432	//
433	//Add aGrid1 and aGrid2 and deposit the result into the current one
434	//
435
436	if(fNVar!=aGrid1->GetNVar()\|\| fNVar!=aGrid2->GetNVar()){
437	AliInfo("Different number of variables, cannot add the grids");
438	return;
439	}
440	if(fNDim!=aGrid1->GetNDim()\|\| fNDim!=aGrid2->GetNDim()){
441	AliInfo("Different number of dimensions, cannot add the grids!");
442	return;
443	}
444
445	if(!fSumW2 && (aGrid1->GetSumW2() \|\| aGrid2->GetSumW2()))SumW2();
446
447
448	fData->Reset();
449	fData->Add(((AliCFGridSparse*)aGrid1)->GetGrid(),c1);
450	fData->Add(((AliCFGridSparse*)aGrid2)->GetGrid(),c2);
451
452	}
453
454	//____________________________________________________________________
455	void AliCFGridSparse::Multiply(AliCFVGrid* aGrid, Double_t c)
456	{
457	//
458	// Multiply aGrid to the current one
459	//
460
461
462	if(aGrid->GetNVar()!=fNVar){
463	AliInfo("Different number of variables, cannot multiply the grids");
464	return;
465	}
466	if(aGrid->GetNDim()!=fNDim){
467	AliInfo("Different number of dimensions, cannot multiply the grids!");
468	return;
469	}
470
471	if(!fSumW2 && aGrid->GetSumW2())SumW2();
472
473	THnSparse h= ((AliCFGridSparse)aGrid)->GetGrid();
474
475	fData->Multiply(h);
476	fData->Scale(c);
477
478	}
479
480	//____________________________________________________________________
481	void AliCFGridSparse::Multiply(AliCFVGrid* aGrid1, AliCFVGrid* aGrid2, Double_t c1,Double_t c2)
482	{
483	//
484	//Multiply aGrid1 and aGrid2 and deposit the result into the current one
485	//
486
487	if(fNVar!=aGrid1->GetNVar()\|\| fNVar!=aGrid2->GetNVar()){
488	AliInfo("Different number of variables, cannot multiply the grids");
489	return;
490	}
491	if(fNDim!=aGrid1->GetNDim()\|\| fNDim!=aGrid2->GetNDim()){
492	AliInfo("Different number of dimensions, cannot multiply the grids!");
493	return;
494	}
495
496	if(!fSumW2 && (aGrid1->GetSumW2() \|\| aGrid2->GetSumW2()))SumW2();
497
498
499	fData->Reset();
500	THnSparse h1= ((AliCFGridSparse)aGrid1)->GetGrid();
501	THnSparse h2= ((AliCFGridSparse)aGrid2)->GetGrid();
502	h2->Multiply(h1);
503	h2->Scale(c1*c2);
504	fData->Add(h2);
505	}
506
507
508
509	//____________________________________________________________________
510	void AliCFGridSparse::Divide(AliCFVGrid* aGrid, Double_t c)
511	{
512	//
513	// Divide aGrid to the current one
514	//
515
516
517	if(aGrid->GetNVar()!=fNVar){
518	AliInfo("Different number of variables, cannot divide the grids");
519	return;
520	}
521	if(aGrid->GetNDim()!=fNDim){
522	AliInfo("Different number of dimensions, cannot divide the grids!");
523	return;
524	}
525
526	if(!fSumW2 && aGrid->GetSumW2())SumW2();
527
528	THnSparse h= ((AliCFGridSparse)aGrid)->GetGrid();
529
530	fData->Divide(h);
531	fData->Scale(c);
532
533	}
534
535	//____________________________________________________________________
536	void AliCFGridSparse::Divide(AliCFVGrid* aGrid1, AliCFVGrid* aGrid2, Double_t c1,Double_t c2, Option_t *option)
537	{
538	//
539	//Divide aGrid1 and aGrid2 and deposit the result into the current one
540	//bynomial errors are supported
541	//
542
543	if(fNVar!=aGrid1->GetNVar()\|\| fNVar!=aGrid2->GetNVar()){
544	AliInfo("Different number of variables, cannot divide the grids");
545	return;
546	}
547	if(fNDim!=aGrid1->GetNDim()\|\| fNDim!=aGrid2->GetNDim()){
548	AliInfo("Different number of dimensions, cannot divide the grids!");
549	return;
550	}
551
552	if(!fSumW2 && (aGrid1->GetSumW2() \|\| aGrid2->GetSumW2()))SumW2();
553
554
555	THnSparse h1= ((AliCFGridSparse)aGrid1)->GetGrid();
556	THnSparse h2= ((AliCFGridSparse)aGrid2)->GetGrid();
557	fData->Divide(h1,h2,c1,c2,option);
558	}
559
560
7411edfd	561	//____________________________________________________________________
	562	void AliCFGridSparse::Rebin(const Int_t* group)
	563	{
	564	//
	565	// rebin the grid according to Rebin() as in THnSparse
	566	// Please notice that the original number of bins on
	567	// a given axis has to be divisible by the rebin group.
	568	//
	569
	570	for(Int_t i=0;i<fNVar;i++){
	571	if(group[i]!=1)AliInfo(Form(" merging bins along dimension %i in groups of %i bins", i,group[i]));
	572	}
	573
	574	THnSparse *rebinned =fData->Rebin(group);
	575	fData->Reset();
	576	fData = rebinned;
	577
	578	//redefine the needed stuff
	579
	580	Int_t ndimTot=1;
	581	Int_t nbinTot=0;
	582
	583	//number of bins in each dimension, auxiliary variables
	584
	585	for(Int_t ivar=0;ivar<fNVar;ivar++){
	586	Int_t nbins = fData->GetAxis(ivar)->GetNbins();
	587	fNVarBins[ivar]=nbins;
	588	ndimTot*=fNVarBins[ivar];
	589	nbinTot+=(fNVarBins[ivar]+1);
	590	Int_t offset=0;
	591	for(Int_t i =0;i<ivar;i++)offset+=(fNVarBins[i]+1);
	592	fOffset[ivar]=offset;
	593	Int_t prod=1;
	594	for(Int_t i=0;i<ivar;i++)prod*=fNVarBins[i];
	595	fProduct[ivar]=prod;
	596	}
	597
	598	fNDim=ndimTot;
	599
	600	//now the array of bin limits
	601
	602	delete fVarBinLimits;
	603	fNVarBinLimits=nbinTot;
	604	fVarBinLimits=new Double_t[fNVarBinLimits];
	605
	606	for(Int_t ivar=0;ivar<fNVar;ivar++){
	607	Double_t low = fData->GetAxis(ivar)->GetXmin();
	608	Double_t high = fData->GetAxis(ivar)->GetXmax();
	609	const TArrayD *xbins = fData->GetAxis(ivar)->GetXbins();
	610	if (xbins->fN == 0){
	611	for(Int_t ibin=0;ibin<=fNVarBins[ivar];ibin++){
	612	fVarBinLimits[ibin+fOffset[ivar]] = low + ibin*(high-low)/((Double_t) fNVarBins[ivar]);
	613	}
	614	}
	615	else{
	616
	617	for(Int_t ibin=0;ibin<=fNVarBins[ivar];ibin++) {
	618	fVarBinLimits[ibin+fOffset[ivar]] = xbins->At(ibin);
	619	}
	620	}
	621	}
	622
	623	}
db6722a5	624	//____________________________________________________________________
	625	void AliCFGridSparse::Copy(TObject& c) const
	626	{
	627	//
	628	// copy function
	629	//
	630	AliCFGridSparse& target = (AliCFGridSparse &) c;
	631
	632	if(fData)target.fData = fData;
	633	}
	634