#include <TClonesArray.h>
#include <TMath.h>
#include <TRandom.h>
#include "AliH2F.h"
ClassImp(AliH2F)
AliH2F::AliH2F():TH2F()
{
}
AliH2F::AliH2F(const Text_t *name,const Text_t *title,
Int_t nbinsx,Axis_t xlow,Axis_t xup
,Int_t nbinsy,Axis_t ylow,Axis_t yup):
TH2F(name,title,nbinsx,xlow,xup
,nbinsy,ylow,yup)
{
}
AliH2F::~AliH2F()
{
}
AliH2F::AliH2F(const AliH2F &his) :
TH2F(his)
{
}
AliH2F & AliH2F::operator = (const AliH2F & )
{
return *this;
}
void AliH2F::ClearSpectrum()
{
Int_t dimx = fXaxis.GetNbins();
Int_t dimy = fYaxis.GetNbins();
for (Int_t i = 0 ;i<dimx;i++)
for (Int_t j = 0 ;j<dimy;j++)
{
SetBinContent(GetBin(i,j),0);
SetBinError(GetBin(i,j),0);
}
}
void AliH2F::AddNoise(Float_t sn)
{
Int_t dimx = fXaxis.GetNbins();
Int_t dimy = fYaxis.GetNbins();
for (Int_t i = 0 ;i<dimx;i++)
for (Int_t j = 0 ;j<dimy;j++)
{
Float_t noise = gRandom->Gaus(0,sn);
Float_t oldv =GetBinContent(GetBin(i,j));
Float_t olds =GetBinError(GetBin(i,j));
if (noise >0)
{
SetBinContent(GetBin(i,j),noise+oldv);
SetBinError(GetBin(i,j),TMath::Sqrt((noise*noise+olds*olds)));
}
}
}
void AliH2F::AddGauss(Float_t x, Float_t y,
Float_t sx, Float_t sy, Float_t max)
{
Int_t dimx = fXaxis.GetNbins();
Int_t dimy = fYaxis.GetNbins();
Float_t dx =(GetXaxis()->GetXmax()-GetXaxis()->GetXmin())/Float_t(dimx);
Float_t dy =(GetYaxis()->GetXmax()-GetYaxis()->GetXmin())/Float_t(dimy);
sx/=dx;
sy/=dy;
for (Int_t i = 0 ;i<dimx;i++)
for (Int_t j = 0 ;j<dimy;j++)
{
Float_t x2 =GetXaxis()->GetBinCenter(i+1);
Float_t y2 =GetYaxis()->GetBinCenter(j+1);
Float_t dx2 = (x2-x)*(x2-x);
Float_t dy2 = (y2-y)*(y2-y);
Float_t amp =max*exp(-(dx2/(2*sx*sx)+dy2/(2*sy*sy)));
Fill(x2,y2,amp);
}
}
void AliH2F::ClearUnderTh(Int_t threshold)
{
Int_t dimx = fXaxis.GetNbins();
Int_t dimy = fYaxis.GetNbins();
for (Int_t i = 0 ;i<=dimx;i++)
for (Int_t j = 0 ;j<=dimy;j++)
{
Float_t oldv =GetBinContent(GetBin(i,j));
if (oldv <threshold)
SetBinContent(GetBin(i,j),0);
}
}
void AliH2F::Round()
{
Int_t dimx = fXaxis.GetNbins();
Int_t dimy = fYaxis.GetNbins();
for (Int_t i = 0 ;i<=dimx;i++)
for (Int_t j = 0 ;j<=dimy;j++)
{
Float_t oldv =GetBinContent(GetBin(i,j));
oldv=(Int_t)oldv;
SetBinContent(GetBin(i,j),oldv);
}
}
AliH2F *AliH2F::GetSubrange2d(Float_t xmin, Float_t xmax,
Float_t ymin, Float_t ymax)
{
if (xmax<=xmin) {
xmin=fXaxis.GetXmin();
xmax=fXaxis.GetXmax();
}
if (ymax<=ymin) {
ymin=fYaxis.GetXmin();
ymax=fYaxis.GetXmax();
}
Int_t nx = Int_t((xmax-xmin)/(fXaxis.GetXmax()-fXaxis.GetXmin()) *
Float_t(fXaxis.GetNbins()));
Int_t ny = Int_t((ymax-ymin)/(fYaxis.GetXmax()-fYaxis.GetXmin()) *
Float_t(fYaxis.GetNbins()));
TString t1 = fName ;
TString t2 = fTitle ;
t1+="_subrange";
t2+="_subrange";
const Text_t *ktt1 = t1;
const Text_t *ktt2 = t2;
AliH2F * sub = new AliH2F(ktt1,ktt2,nx,xmin,xmax,ny,ymin,ymax);
Int_t i1 = Int_t( Float_t(fXaxis.GetNbins())*(xmin-fXaxis.GetXmin())/
(fXaxis.GetXmax()-fXaxis.GetXmin()) ) ;
Int_t i2 = Int_t( Float_t(fYaxis.GetNbins())*(ymin-fYaxis.GetXmin())/
(fYaxis.GetXmax()-fYaxis.GetXmin()) ) ;
for (Int_t i=0;i<nx;i++)
for (Int_t j=0;j<ny;j++)
{
Int_t index1 = GetBin(i1+i,i2+j);
Float_t val = GetBinContent(index1);
sub->SetBinContent(GetBin(i,j),val);
}
return sub;
}
TH1F *AliH2F::GetAmplitudes(Float_t zmin, Float_t zmax, Float_t th, Float_t xmin, Float_t xmax,
Float_t ymin, Float_t ymax)
{
if (xmax<=xmin) {
xmin=fXaxis.GetXmin();
xmax=fXaxis.GetXmax();
}
if (ymax<=ymin) {
ymin=fYaxis.GetXmin();
ymax=fYaxis.GetXmax();
}
Int_t nx = Int_t((xmax-xmin)/(fXaxis.GetXmax()-fXaxis.GetXmin()) *
Float_t(fXaxis.GetNbins()));
Int_t ny = Int_t((ymax-ymin)/(fYaxis.GetXmax()-fYaxis.GetXmin()) *
Float_t(fYaxis.GetNbins()));
TString t1 = fName ;
TString t2 = fTitle ;
t1+="_amplitudes";
t2+="_amplitudes";
const Text_t *ktt1 = t1;
const Text_t *ktt2 = t2;
TH1F * h = new TH1F(ktt1,ktt2,100,zmin,zmax);
Int_t i1 = Int_t( Float_t(fXaxis.GetNbins())*(xmin-fXaxis.GetXmin())/
(fXaxis.GetXmax()-fXaxis.GetXmin()) ) ;
Int_t i2 = Int_t( Float_t(fYaxis.GetNbins())*(ymin-fYaxis.GetXmin())/
(fYaxis.GetXmax()-fYaxis.GetXmin()) ) ;
for (Int_t i=0;i<nx;i++)
for (Int_t j=0;j<ny;j++)
{
Int_t index1 = GetBin(i1+i,i2+j);
Float_t val = GetBinContent(index1);
if (val>th) h->Fill(val);
}
return h;
}
Float_t AliH2F::GetOccupancy(Float_t th , Float_t xmin, Float_t xmax,
Float_t ymin, Float_t ymax)
{
if (xmax<=xmin) {
xmin=fXaxis.GetXmin();
xmax=fXaxis.GetXmax();
}
if (ymax<=ymin) {
ymin=fYaxis.GetXmin();
ymax=fYaxis.GetXmax();
}
Int_t nx = Int_t((xmax-xmin)/(fXaxis.GetXmax()-fXaxis.GetXmin()) *
Float_t(fXaxis.GetNbins()));
Int_t ny = Int_t((ymax-ymin)/(fYaxis.GetXmax()-fYaxis.GetXmin()) *
Float_t(fYaxis.GetNbins()));
Int_t over =0;
Int_t i1 = Int_t( Float_t(fXaxis.GetNbins())*(xmin-fXaxis.GetXmin())/
(fXaxis.GetXmax()-fXaxis.GetXmin()) ) ;
Int_t i2 = Int_t( Float_t(fYaxis.GetNbins())*(ymin-fYaxis.GetXmin())/
(fYaxis.GetXmax()-fYaxis.GetXmin()) ) ;
for (Int_t i=0;i<nx;i++)
for (Int_t j=0;j<ny;j++)
{
Int_t index1 = GetBin(i1+i,i2+j);
Float_t val = GetBinContent(index1);
if (val>th) over++;
}
Int_t all = nx*ny;
if (all>0) return Float_t(over)/Float_t(all);
else
return 0;
}