Formatting changes.
[u/mrichter/AliRoot.git] / TPC / AliH2F.cxx
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 /* $Id$ */
17
18 //----------------------------------------------------------------------------
19 //  Author:   Marian Ivanov
20 //
21 //  Implementation of class AliH2F
22 //
23 //-----------------------------------------------------------------------------
24
25 #include "AliH2F.h"
26 #include "TClonesArray.h"
27 #include "TRandom.h"
28
29
30 ClassImp(AliH2F)
31 //***********************************************************************
32 //***********************************************************************
33 //***********************************************************************
34 //***********************************************************************
35 AliH2F::AliH2F():TH2F() 
36 {
37   //
38  
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 {
46   //
47   
48 }
49      
50 AliH2F::~AliH2F() 
51 {
52   //
53 }
54
55 AliH2F::AliH2F(const AliH2F &his) :
56   TH2F(his)
57 {
58   //
59   
60 }
61
62 AliH2F & AliH2F::operator = (const AliH2F & /*his*/) 
63 {
64   //
65   return *this;
66 }
67
68 /*
69 TClonesArray * AliH2F::FindPeaks(Float_t threshold, Float_t noise)
70 {
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();
80 }
81 */
82
83 void AliH2F::ClearSpectrum()
84 {
85   //clera histogram
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 {
99   // add gauss noise with sigma sn
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((noise*noise+olds*olds)));
112           }
113       }
114 }
115 void AliH2F::AddGauss(Float_t x, Float_t y, 
116                           Float_t sx, Float_t sy, Float_t max)
117 {  
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 =max*exp(-(dx2/(2*sx*sx)+dy2/(2*sy*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 {
145   //clear histogram for bin under threshold
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 {
159   //round float to integer 
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";
198   const Text_t *ktt1 = t1;
199   const Text_t *ktt2 = t2;
200   
201   AliH2F * sub = new AliH2F(ktt1,ktt2,nx,xmin,xmax,ny,ymin,ymax); 
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";
245   const  Text_t *ktt1 = t1;
246   const Text_t *ktt2 = t2;
247   
248   TH1F * h = new TH1F(ktt1,ktt2,100,zmin,zmax); 
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 }