HLT/hough/AliL3HistogramAdaptive.cxx

   1 // @(#) $Id$
   2
   3 // Author: Anders Vestbo <mailto:vestbo@fi.uib.no>
   4 //*-- Copyright &copy ALICE HLT Group
   5
   6 #include "AliL3StandardIncludes.h"
   7 #include "AliL3Logging.h"
   8 #include "AliL3HistogramAdaptive.h"
   9 #include "AliL3Transform.h"
  10 #include "AliL3Track.h"
  11
  12 #if GCCVERSION == 3
  13 using namespace std;
  14 #endif
  15
  16 //_____________________________________________________________
  17 // AliL3HistogramAdaptive
  18 //
  19 // 2D histogram class adapted for kappa and psi as used in the Circle Hough Transform.
  20 // The bins in kappa is not linear, but has a width which is specified by argument
  21 // ptres in the constructor. This gives the relative pt resolution which should
  22 // be kept throughout the kappa range.
  23
  24 ClassImp(AliL3HistogramAdaptive)
  25
  26 AliL3HistogramAdaptive::AliL3HistogramAdaptive() : AliL3Histogram()
  27 {
  28   fKappaBins=0;
  29 }
  30
  31
  32 AliL3HistogramAdaptive::AliL3HistogramAdaptive(Char_t *name,Double_t minpt,Double_t maxpt,Double_t ptres,
  33                                                Int_t nybins,Double_t ymin,Double_t ymax)
  34 {
  35   strcpy(fName,name);
  36
  37   fPtres = ptres;
  38   fXmin = -1*AliL3Transform::GetBFact()*AliL3Transform::GetBField()/minpt;
  39   fXmax = AliL3Transform::GetBFact()*AliL3Transform::GetBField()/minpt;
  40
  41   fMinPt = minpt;
  42   fMaxPt = maxpt;
  43   fNxbins = InitKappaBins();
  44   fNybins = nybins;
  45
  46   fYmin = ymin;
  47   fYmax = ymax;
  48   fFirstXbin=1;
  49   fFirstYbin=1;
  50   fLastXbin = fNxbins;
  51   fLastYbin = fNybins;
  52   fNcells = (fNxbins+2)*(fNybins+2);
  53
  54   fThreshold=0;
  55   fContent = new Int_t[fNcells];
  56   Reset();
  57 }
  58
  59 AliL3HistogramAdaptive::~AliL3HistogramAdaptive()
  60 {
  61   if(fKappaBins)
  62     delete [] fKappaBins;
  63 }
  64
  65 Int_t AliL3HistogramAdaptive::InitKappaBins()
  66 {
  67   //Here a LUT for the kappa values created. This has to be done since
  68   //the binwidth in kappa is not constant, but change according to the
  69   //set relative resolution in pt.
  70   //Since the kappa values are symmetric about origo, the size of the
  71   //LUT is half of the total number of bins in kappa direction.
  72
  73   Double_t pt = fMinPt,delta_pt,local_pt;
  74   Int_t bin=0;
  75
  76   while(pt < fMaxPt)
  77     {
  78       local_pt = pt;
  79       delta_pt = fPtres*local_pt;
  80       pt += delta_pt;
  81       bin++;
  82     }
  83   fKappaBins = new Double_t[bin+1];
  84   pt=fMinPt;
  85   bin=0;
  86   fKappaBins[bin] = AliL3Transform::GetBFact()*AliL3Transform::GetBField()/fMinPt;
  87   while(pt < fMaxPt)
  88     {
  89       local_pt = pt;
  90       delta_pt = fPtres*local_pt;
  91       pt += delta_pt;
  92       bin++;
  93       fKappaBins[bin] = AliL3Transform::GetBFact()*AliL3Transform::GetBField()/pt;
  94     }
  95   return (bin+1)*2; //Both negative and positive kappa.
  96 }
  97
  98
  99 void AliL3HistogramAdaptive::Fill(Double_t x,Double_t y,Int_t weight)
 100 {
 101   Int_t bin = FindBin(x,y);
 102   if(bin < 0)
 103     return;
 104   AddBinContent(bin,weight);
 105
 106 }
 107
 108 Int_t AliL3HistogramAdaptive::FindBin(Double_t x,Double_t y)
 109 {
 110
 111   Int_t xbin = FindXbin(x);
 112   Int_t ybin = FindYbin(y);
 113
 114   if(xbin < 0)
 115     return -1;
 116   return GetBin(xbin,ybin);
 117 }
 118
 119 Int_t AliL3HistogramAdaptive::FindXbin(Double_t x)
 120 {
 121
 122   if(x < fXmin || x > fXmax || fabs(x) < fKappaBins[(fNxbins/2-1)])
 123     return -1;
 124
 125   //Remember that kappa value is decreasing with bin number!
 126   //Also, the bin numbering starts at 1 and ends at fNxbins,
 127   //so the corresponding elements in the LUT is bin - 1.
 128
 129   Int_t bin=0;
 130   while(bin < fNxbins/2)
 131     {
 132       if(fabs(x) <= fKappaBins[bin] && fabs(x) > fKappaBins[bin+1])
 133         break;
 134       bin++;
 135     }
 136   if(x < 0)
 137     return bin + 1;
 138   else
 139     return fNxbins - bin;
 140
 141 }
 142
 143 Int_t AliL3HistogramAdaptive::FindYbin(Double_t y)
 144 {
 145   if(y < fYmin || y > fYmax)
 146     return 0;
 147
 148   return 1 + (Int_t)(fNybins*(y-fYmin)/(fYmax-fYmin));
 149 }
 150
 151 Double_t AliL3HistogramAdaptive::GetBinCenterX(Int_t xbin)
 152 {
 153   if(xbin < fFirstXbin || xbin > fLastXbin)
 154     {
 155       LOG(AliL3Log::kWarning,"AliL3HistogramAdaptive::GetBinCenterX","Bin-value")
 156         <<"XBinvalue out of range "<<xbin<<ENDLOG;
 157       return 0;
 158     }
 159
 160   //The bin numbers go from 1 to fNxbins, so the corresponding
 161   //element in the LUT is xbin - 1. This is the reason why we
 162   //substract a 1 here:
 163
 164   Int_t bin = xbin;
 165   bin -= 1;
 166   if(bin >= fNxbins/2)
 167     bin = fNxbins - 1 - bin;
 168
 169   //Remember again that the kappa-values are _decreasing_ with bin number.
 170
 171   Double_t binwidth = fKappaBins[bin] - fKappaBins[bin+1];
 172   Double_t kappa = fKappaBins[bin] - 0.5*binwidth;
 173   if(xbin < fNxbins/2)
 174     return -1.*kappa;
 175   else
 176     return kappa;
 177
 178 }
 179
 180 Double_t AliL3HistogramAdaptive::GetBinCenterY(Int_t ybin)
 181 {
 182   if(ybin < fFirstYbin || ybin > fLastYbin)
 183     {
 184       LOG(AliL3Log::kError,"AliL3HistogramAdaptive::GetBinCenterY","ybin")
 185         <<"Bin-value out of range "<<ybin<<ENDLOG;
 186       return -1;
 187     }
 188   Double_t binwidth = (fYmax - fYmin) / fNybins;
 189   return fYmin + (ybin-0.5) * binwidth;
 190
 191 }
 192
 193
 194 void AliL3HistogramAdaptive::Draw(Char_t *option)
 195 {
 196 #ifdef use_root
 197   if(!fRootHisto)
 198     CreateRootHisto();
 199
 200   Double_t kappa,psi;
 201   Int_t content,bin;
 202   for(Int_t i=fFirstXbin; i<=fLastXbin; i++)
 203     {
 204       kappa = GetBinCenterX(i);
 205       for(Int_t j=fFirstYbin; j<=fLastYbin; j++)
 206         {
 207           psi = GetBinCenterY(j);
 208           bin = GetBin(i,j);
 209           content = GetBinContent(bin);
 210           fRootHisto->Fill(kappa,psi,content);
 211         }
 212     }
 213   fRootHisto->Draw(option);
 214   return;
 215 #endif
 216   cerr<<"AliL3HistogramAdaptive::Draw : You need to compile with ROOT in order to draw histogram"<<endl;
 217 }
 218
 219 void AliL3HistogramAdaptive::Print()
 220 {
 221   cout<<"Printing content of histogram "<<fName<<endl;
 222   for(Int_t i=0; i<fNcells; i++)
 223     {
 224       if(GetBinContent(i)==0) continue;
 225       cout<<"Bin "<<i<<": "<<GetBinContent(i)<<endl;
 226     }
 227
 228 }