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