[u/mrichter/AliRoot.git] / CONTAINERS / AliH2F.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.                  *
 **************************************************************************/

/* $Id$ */

//----------------------------------------------------------------------------
//  Author:   Marian Ivanov
//
//  Implementation of class AliH2F
//
//-----------------------------------------------------------------------------

#include "AliH2F.h"
#include "TClonesArray.h"
#include "TRandom.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 & /*his*/) 
{
  //
  return *this;
}

/*
TClonesArray * AliH2F::FindPeaks(Float_t threshold, Float_t noise)
{
  //find peaks and write it in form of AliTPCcluster to array
    
  //firstly we need to create object for cluster finding
  //and fill it with contents of histogram
  AliTPCClusterFinder cfinder;
  cfinder.SetThreshold(threshold);
  cfinder.SetNoise(noise);
  cfinder.GetHisto(this);
  return cfinder.FindPeaks3();
}
*/

void AliH2F::ClearSpectrum()
{
  //clera histogram
  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++) 
      {
	SetCellContent(i,j,0);
	SetCellError(i,j,0);
      }
}


void AliH2F::AddNoise(Float_t sn)
{
  // add gauss noise with sigma 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  =GetCellContent(i,j);
	Float_t olds  =GetCellError(i,j);
	if (noise >0)
	  {
	    SetCellContent(i,j,noise+oldv);
	    SetCellError(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)
{  
  //transform to histogram coordinata  
  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);  
  //  x=(x-GetXaxis()->GetXmin())/dx;
  //y=(y-GetYaxis()->GetXmin())/dy;
  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)));
	//Float_t oldv  =GetCellContent(i+1,j+1);
	//	SetCellContent(i+1,j+1,amp+oldv);
	Fill(x2,y2,amp);
      }
}

void AliH2F::ClearUnderTh(Int_t threshold)
{
  //clear histogram for bin under 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  =GetCellContent(i,j);
        if (oldv <threshold)
	  SetCellContent(i,j,0);
      }
}

void AliH2F::Round()
{
  //round float to integer 
  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  =GetCellContent(i,j);
        oldv=(Int_t)oldv;
	SetCellContent(i,j,oldv);
      }
}


AliH2F *AliH2F::GetSubrange2d(Float_t xmin, Float_t xmax, 
				      Float_t ymin, Float_t ymax)
{
  //this function return pointer to the new created 
  //histogram which is subhistogram of the 
  //calculate number
  //subhistogram range must be inside histogram

  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);
	//        Int_t index2 = sub->GetBin(i,j);
        Float_t val = GetBinContent(index1);
	//        sub->SetBinContent(index2,val);
	//        Float_t err = GetBinError(index1);
        //sub->SetBinError(index2,GetBinError(index1));
        sub->SetCellContent(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)
{
  //this function return pointer to the new created 
  //histogram which is subhistogram of the 
  //calculate number
  //subhistogram range must be inside histogram
 
  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)
{
  //this function return pointer to the new created 
  //histogram which is subhistogram of the 
  //calculate number
  //subhistogram range must be inside histogram
 
  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;
}
Commit	Line	Data
cc80f89e	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
5270561d	16	/* $Id$ */
cc80f89e	17
	18	//----------------------------------------------------------------------------
	19	// Author: Marian Ivanov
	20	//
	21	// Implementation of class AliH2F
	22	//
	23	//-----------------------------------------------------------------------------
	24
	25	#include "AliH2F.h"
	26	#include "TClonesArray.h"
cc80f89e	27	#include "TRandom.h"
cc80f89e	28
	29
	30	ClassImp(AliH2F)
	31	//***********************************************************************
	32	//***********************************************************************
	33	//***********************************************************************
	34	//***********************************************************************
	35	AliH2F::AliH2F():TH2F()
	36	{
73042f01	37	//
73042f01	38
cc80f89e	39	}
	40	AliH2F::AliH2F(const Text_t name,const Text_t title,
	41	Int_t nbinsx,Axis_t xlow,Axis_t xup
	42	,Int_t nbinsy,Axis_t ylow,Axis_t yup):
	43	TH2F(name,title,nbinsx,xlow,xup
	44	,nbinsy,ylow,yup)
	45	{
73042f01	46	//
73042f01	47
cc80f89e	48	}
	49
	50	AliH2F::~AliH2F()
	51	{
73042f01	52	//
cc80f89e	53	}
cc80f89e	54
88cb7938	55	AliH2F::AliH2F(const AliH2F &his) :
88cb7938	56	TH2F(his)
cc80f89e	57	{
73042f01	58	//
73042f01	59
cc80f89e	60	}
cc80f89e	61
88cb7938	62	AliH2F & AliH2F::operator = (const AliH2F & /his/)
cc80f89e	63	{
73042f01	64	//
73042f01	65	return *this;
cc80f89e	66	}
cc80f89e	67
88cb7938	68	/*
cc80f89e	69	TClonesArray * AliH2F::FindPeaks(Float_t threshold, Float_t noise)
cc80f89e	70	{
88cb7938	71	//find peaks and write it in form of AliTPCcluster to array
	72
	73	//firstly we need to create object for cluster finding
	74	//and fill it with contents of histogram
	75	AliTPCClusterFinder cfinder;
	76	cfinder.SetThreshold(threshold);
	77	cfinder.SetNoise(noise);
	78	cfinder.GetHisto(this);
	79	return cfinder.FindPeaks3();
cc80f89e	80	}
88cb7938	81	*/
cc80f89e	82
	83	void AliH2F::ClearSpectrum()
	84	{
73042f01	85	//clera histogram
cc80f89e	86	Int_t dimx = fXaxis.GetNbins();
	87	Int_t dimy = fYaxis.GetNbins();
	88	for (Int_t i = 0 ;i<dimx;i++)
	89	for (Int_t j = 0 ;j<dimy;j++)
	90	{
	91	SetCellContent(i,j,0);
	92	SetCellError(i,j,0);
	93	}
	94	}
	95
	96
	97	void AliH2F::AddNoise(Float_t sn)
	98	{
73042f01	99	// add gauss noise with sigma sn
cc80f89e	100	Int_t dimx = fXaxis.GetNbins();
	101	Int_t dimy = fYaxis.GetNbins();
	102	for (Int_t i = 0 ;i<dimx;i++)
	103	for (Int_t j = 0 ;j<dimy;j++)
	104	{
	105	Float_t noise = gRandom->Gaus(0,sn);
	106	Float_t oldv =GetCellContent(i,j);
	107	Float_t olds =GetCellError(i,j);
	108	if (noise >0)
	109	{
	110	SetCellContent(i,j,noise+oldv);
	111	SetCellError(i,j,TMath::Sqrt((noisenoise+oldsolds)));
	112	}
	113	}
	114	}
	115	void AliH2F::AddGauss(Float_t x, Float_t y,
	116	Float_t sx, Float_t sy, Float_t max)
73042f01	117	{
cc80f89e	118	//transform to histogram coordinata
	119	Int_t dimx = fXaxis.GetNbins();
	120	Int_t dimy = fYaxis.GetNbins();
	121	Float_t dx =(GetXaxis()->GetXmax()-GetXaxis()->GetXmin())/Float_t(dimx);
	122	Float_t dy =(GetYaxis()->GetXmax()-GetYaxis()->GetXmin())/Float_t(dimy);
	123	// x=(x-GetXaxis()->GetXmin())/dx;
	124	//y=(y-GetYaxis()->GetXmin())/dy;
	125	sx/=dx;
	126	sy/=dy;
	127
	128
	129	for (Int_t i = 0 ;i<dimx;i++)
	130	for (Int_t j = 0 ;j<dimy;j++)
	131	{
	132	Float_t x2 =GetXaxis()->GetBinCenter(i+1);
	133	Float_t y2 =GetYaxis()->GetBinCenter(j+1);
	134	Float_t dx2 = (x2-x)*(x2-x);
	135	Float_t dy2 = (y2-y)*(y2-y);
	136	Float_t amp =maxexp(-(dx2/(2sxsx)+dy2/(2sy*sy)));
	137	//Float_t oldv =GetCellContent(i+1,j+1);
	138	// SetCellContent(i+1,j+1,amp+oldv);
	139	Fill(x2,y2,amp);
	140	}
	141	}
	142
	143	void AliH2F::ClearUnderTh(Int_t threshold)
	144	{
73042f01	145	//clear histogram for bin under threshold
cc80f89e	146	Int_t dimx = fXaxis.GetNbins();
	147	Int_t dimy = fYaxis.GetNbins();
	148	for (Int_t i = 0 ;i<=dimx;i++)
	149	for (Int_t j = 0 ;j<=dimy;j++)
	150	{
	151	Float_t oldv =GetCellContent(i,j);
	152	if (oldv <threshold)
	153	SetCellContent(i,j,0);
	154	}
	155	}
	156
	157	void AliH2F::Round()
	158	{
73042f01	159	//round float to integer
cc80f89e	160	Int_t dimx = fXaxis.GetNbins();
	161	Int_t dimy = fYaxis.GetNbins();
	162	for (Int_t i = 0 ;i<=dimx;i++)
	163	for (Int_t j = 0 ;j<=dimy;j++)
	164	{
	165	Float_t oldv =GetCellContent(i,j);
	166	oldv=(Int_t)oldv;
	167	SetCellContent(i,j,oldv);
	168	}
	169	}
	170
	171
	172
	173	AliH2F *AliH2F::GetSubrange2d(Float_t xmin, Float_t xmax,
	174	Float_t ymin, Float_t ymax)
	175	{
	176	//this function return pointer to the new created
	177	//histogram which is subhistogram of the
	178	//calculate number
	179	//subhistogram range must be inside histogram
	180
	181	if (xmax<=xmin) {
	182	xmin=fXaxis.GetXmin();
	183	xmax=fXaxis.GetXmax();
	184	}
	185	if (ymax<=ymin) {
	186	ymin=fYaxis.GetXmin();
	187	ymax=fYaxis.GetXmax();
	188	}
	189
	190	Int_t nx = Int_t((xmax-xmin)/(fXaxis.GetXmax()-fXaxis.GetXmin()) *
	191	Float_t(fXaxis.GetNbins()));
	192	Int_t ny = Int_t((ymax-ymin)/(fYaxis.GetXmax()-fYaxis.GetXmin()) *
	193	Float_t(fYaxis.GetNbins()));
	194	TString t1 = fName ;
	195	TString t2 = fTitle ;
	196	t1+="_subrange";
	197	t2+="_subrange";
73042f01	198	const Text_t *ktt1 = t1;
73042f01	199	const Text_t *ktt2 = t2;
cc80f89e	200
73042f01	201	AliH2F * sub = new AliH2F(ktt1,ktt2,nx,xmin,xmax,ny,ymin,ymax);
cc80f89e	202
	203	Int_t i1 = Int_t( Float_t(fXaxis.GetNbins())*(xmin-fXaxis.GetXmin())/
	204	(fXaxis.GetXmax()-fXaxis.GetXmin()) ) ;
	205	Int_t i2 = Int_t( Float_t(fYaxis.GetNbins())*(ymin-fYaxis.GetXmin())/
	206	(fYaxis.GetXmax()-fYaxis.GetXmin()) ) ;
	207	for (Int_t i=0;i<nx;i++)
	208	for (Int_t j=0;j<ny;j++)
	209	{
	210	Int_t index1 = GetBin(i1+i,i2+j);
	211	// Int_t index2 = sub->GetBin(i,j);
	212	Float_t val = GetBinContent(index1);
	213	// sub->SetBinContent(index2,val);
	214	// Float_t err = GetBinError(index1);
	215	//sub->SetBinError(index2,GetBinError(index1));
	216	sub->SetCellContent(i,j,val);
	217	}
	218	return sub;
	219	}
	220
	221	TH1F *AliH2F::GetAmplitudes(Float_t zmin, Float_t zmax, Float_t th, Float_t xmin, Float_t xmax,
	222	Float_t ymin, Float_t ymax)
	223	{
	224	//this function return pointer to the new created
	225	//histogram which is subhistogram of the
	226	//calculate number
	227	//subhistogram range must be inside histogram
	228
	229	if (xmax<=xmin) {
	230	xmin=fXaxis.GetXmin();
	231	xmax=fXaxis.GetXmax();
	232	}
	233	if (ymax<=ymin) {
	234	ymin=fYaxis.GetXmin();
	235	ymax=fYaxis.GetXmax();
	236	}
	237	Int_t nx = Int_t((xmax-xmin)/(fXaxis.GetXmax()-fXaxis.GetXmin()) *
	238	Float_t(fXaxis.GetNbins()));
	239	Int_t ny = Int_t((ymax-ymin)/(fYaxis.GetXmax()-fYaxis.GetXmin()) *
	240	Float_t(fYaxis.GetNbins()));
	241	TString t1 = fName ;
	242	TString t2 = fTitle ;
	243	t1+="_amplitudes";
	244	t2+="_amplitudes";
73042f01	245	const Text_t *ktt1 = t1;
73042f01	246	const Text_t *ktt2 = t2;
cc80f89e	247
73042f01	248	TH1F * h = new TH1F(ktt1,ktt2,100,zmin,zmax);
cc80f89e	249
	250	Int_t i1 = Int_t( Float_t(fXaxis.GetNbins())*(xmin-fXaxis.GetXmin())/
	251	(fXaxis.GetXmax()-fXaxis.GetXmin()) ) ;
	252	Int_t i2 = Int_t( Float_t(fYaxis.GetNbins())*(ymin-fYaxis.GetXmin())/
	253	(fYaxis.GetXmax()-fYaxis.GetXmin()) ) ;
	254	for (Int_t i=0;i<nx;i++)
	255	for (Int_t j=0;j<ny;j++)
	256	{
	257	Int_t index1 = GetBin(i1+i,i2+j);
	258	Float_t val = GetBinContent(index1);
	259	if (val>th) h->Fill(val);
	260	}
	261	return h;
	262	}
	263
	264	Float_t AliH2F::GetOccupancy(Float_t th , Float_t xmin, Float_t xmax,
	265	Float_t ymin, Float_t ymax)
	266	{
	267	//this function return pointer to the new created
	268	//histogram which is subhistogram of the
	269	//calculate number
	270	//subhistogram range must be inside histogram
	271
	272	if (xmax<=xmin) {
	273	xmin=fXaxis.GetXmin();
	274	xmax=fXaxis.GetXmax();
	275	}
	276	if (ymax<=ymin) {
	277	ymin=fYaxis.GetXmin();
	278	ymax=fYaxis.GetXmax();
	279	}
	280	Int_t nx = Int_t((xmax-xmin)/(fXaxis.GetXmax()-fXaxis.GetXmin()) *
	281	Float_t(fXaxis.GetNbins()));
	282	Int_t ny = Int_t((ymax-ymin)/(fYaxis.GetXmax()-fYaxis.GetXmin()) *
	283	Float_t(fYaxis.GetNbins()));
	284
	285	Int_t over =0;
	286	Int_t i1 = Int_t( Float_t(fXaxis.GetNbins())*(xmin-fXaxis.GetXmin())/
	287	(fXaxis.GetXmax()-fXaxis.GetXmin()) ) ;
	288	Int_t i2 = Int_t( Float_t(fYaxis.GetNbins())*(ymin-fYaxis.GetXmin())/
	289	(fYaxis.GetXmax()-fYaxis.GetXmin()) ) ;
	290	for (Int_t i=0;i<nx;i++)
	291	for (Int_t j=0;j<ny;j++)
	292	{
	293	Int_t index1 = GetBin(i1+i,i2+j);
	294	Float_t val = GetBinContent(index1);
	295	if (val>th) over++;
	296	}
	297	Int_t all = nx*ny;
	298	if (all>0) return Float_t(over)/Float_t(all);
	299	else
	300	return 0;
	301	}