HBTAN/AliHBTFunction.cxx

   1 #include "AliHBTFunction.h"
   2 /******************************************************************/
   3 /*
   4 Piotr Krzysztof Skowronski
   5 Piotr.Skowronski@cern.ch
   6 Base classes for HBT functions
   7
   8              function
   9               /    \
  10              /      \
  11             /        \
  12            /          \
  13           /            \
  14          /              \
  15         /                \
  16     two part          four part
  17     /   |   \         /   |   \
  18    /    |    \       /    |    \
  19   1D   2D    3D     1D   2D    3D
  20
  21  four particle functions are intendent to be resolution functions:
  22  it is mecessary to have simulated particle pair corresponding to given
  23  recontructed track pair in order to calculate function simualted value
  24  and recontructed value to be further histogrammed
  25
  26 */
  27 /******************************************************************/
  28 /******************************************************************/
  29
  30 #include <iostream.h>
  31 ClassImp( AliHBTFunction )
  32
  33 AliHBTFunction::AliHBTFunction()
  34 {
  35
  36  fPairCut = new AliHBTEmptyPairCut(); //dummy cut
  37 }
  38 /******************************************************************/
  39
  40 AliHBTFunction::~AliHBTFunction()
  41  {
  42   if (fPairCut) delete fPairCut;
  43  }
  44 /******************************************************************/
  45
  46 void AliHBTFunction::
  47 Write()
  48  {
  49    if (GetNumerator()) GetNumerator()->Write();
  50    if (GetDenominator()) GetDenominator()->Write();
  51    TH1* res = GetResult();
  52    if (res) res->Write();
  53  }
  54 /******************************************************************/
  55
  56 TH1* AliHBTFunction::
  57 GetRatio(Double_t normfactor)
  58  {
  59    if (gDebug>0) cout<<"Mormfactor is "<<normfactor<<" for "<<fName<<endl;
  60
  61    if (normfactor == 0.0)
  62     {
  63       Error("GetRatio","Scaling Factor is 0. Null poiner returned");
  64       return 0x0;
  65     }
  66    TString str = fName + " ratio";
  67    TH1 *result = (TH1*)GetNumerator()->Clone(str.Data());
  68
  69    result->SetTitle(str.Data());
  70    //result->Sumw2();
  71
  72    result->Divide(GetNumerator(),GetDenominator(),normfactor);
  73
  74    return result;
  75
  76  }
  77 /******************************************************************/
  78 void AliHBTFunction::SetPairCut(AliHBTPairCut* cut)
  79 {
  80 //Sets new Pair Cut. Old one is deleted
  81 //Note that it is created new object instead of simple pointer set
  82 //I do not want to have pointer
  83 //to object created somewhere else
  84 //because in that case I could not believe that
  85 //it would always exist (sb could delete it)
  86 //so we have always own copy
  87
  88  if(!cut)
  89    {
  90      Error("AliHBTFunction::SetPairCut","argument is NULL");
  91      return;
  92    }
  93  delete fPairCut;
  94  fPairCut = (AliHBTPairCut*)cut->Clone();
  95
  96 }
  97
  98 /******************************************************************/
  99
 100 void AliHBTFunction::
 101 Rename(const Char_t * name)
 102  {
 103  //renames the function and histograms
 104   SetName(name);
 105   SetTitle(name);
 106
 107   TString numstr = fName + " Numerator";  //title and name of the
 108                                            //numerator histogram
 109   TString denstr = fName + " Denominator";//title and name of the
 110                                            //denominator histogram
 111
 112   GetNumerator()->SetName(numstr.Data());
 113   GetNumerator()->SetTitle(numstr.Data());
 114
 115   GetDenominator()->SetName(denstr.Data());
 116   GetDenominator()->SetTitle(denstr.Data());
 117
 118  }
 119
 120 void AliHBTFunction::
 121 Rename(const Char_t * name, const Char_t * title)
 122  {
 123  //renames and retitle the function and histograms
 124
 125   SetName(name);
 126   SetTitle(title);
 127
 128   TString numstrn = fName + " Numerator";  //name of the
 129                                            //numerator histogram
 130
 131   TString numstrt = fTitle + " Numerator";  //title of the
 132                                            //numerator histogram
 133
 134   TString denstrn = fName + " Denominator";//name of the
 135                                            //denominator histogram
 136
 137   TString denstrt = fTitle + " Denominator";//title of the
 138                                            //denominator histogram
 139
 140
 141   GetNumerator()->SetName(numstrn.Data());
 142   GetNumerator()->SetTitle(numstrt.Data());
 143
 144   GetDenominator()->SetName(denstrn.Data());
 145   GetDenominator()->SetTitle(denstrt.Data());
 146
 147
 148  }
 149
 150 /******************************************************************/
 151 /******************************************************************/
 152 /******************************************************************/
 153
 154 ClassImp( AliHBTTwoPartFctn )
 155
 156 /******************************************************************/
 157 /******************************************************************/
 158 /******************************************************************/
 159
 160 ClassImp( AliHBTFourPartFctn)
 161
 162 /******************************************************************/
 163 /******************************************************************/
 164 /******************************************************************/
 165
 166 ClassImp( AliHBTTwoPartFctn1D )
 167
 168 AliHBTTwoPartFctn1D::
 169 AliHBTTwoPartFctn1D(Int_t nbins, Double_t maxXval, Double_t minXval)
 170  {
 171  //Constructor of Two Part One Dimentional Function
 172  // nbins: number of bins in histograms - default 100
 173  // maxXval and minXval: range of histgram(s) default 0 - 0.15 (GeV)
 174
 175
 176    TString numstr = fName + " Numerator";  //title and name of the
 177                                            //numerator histogram
 178    TString denstr = fName + " Denominator";//title and name of the
 179                                            //denominator histogram
 180
 181    fNumerator   = new TH1D(numstr.Data(),numstr.Data(),nbins,minXval,maxXval);
 182    fDenominator = new TH1D(denstr.Data(),denstr.Data(),nbins,minXval,maxXval);
 183
 184    fNumerator->Sumw2();
 185    fDenominator->Sumw2();
 186
 187    fNBinsToScale = 30;
 188
 189  }
 190 /******************************************************************/
 191 AliHBTTwoPartFctn1D::~AliHBTTwoPartFctn1D()
 192 {
 193   delete fNumerator;
 194   delete fDenominator;
 195 }
 196 /******************************************************************/
 197
 198 void AliHBTTwoPartFctn1D::ProcessSameEventParticles(AliHBTPair* pair)
 199 {
 200  //Fills the numerator
 201    pair = CheckPair(pair);
 202    if(pair) fNumerator->Fill(GetValue(pair));
 203 }
 204 /******************************************************************/
 205 void AliHBTTwoPartFctn1D::ProcessDiffEventParticles(AliHBTPair* pair)
 206  {
 207   //fills denumerator
 208    pair = CheckPair(pair);
 209    if(pair) fDenominator->Fill(GetValue(pair));
 210
 211   }
 212 /******************************************************************/
 213 Double_t AliHBTTwoPartFctn1D::Scale()
 214 {
 215   if (gDebug>0) cout<<"Enetered Scale()"<<endl;
 216   if(!fNumerator)
 217    {
 218      Error("Scale","No numerator");
 219      return 0.0;
 220    }
 221   if(!fDenominator)
 222    {
 223      Error("Scale","No denominator");
 224      return 0.0;
 225    }
 226
 227   if(fNBinsToScale < 1)
 228    {
 229     return 0.0;
 230     Error("Scale","Number of bins for scaling is smaller thnan 1");
 231    }
 232   Int_t nbins = fNumerator->GetNbinsX();
 233   if (fNBinsToScale > nbins)
 234    {
 235     Error("Scale","Number of bins for scaling is bigger thnan number of bins in histograms");
 236     return 0.0;
 237    }
 238   if (gDebug>0) cout<<"No errors detected"<<endl;
 239
 240   Double_t ratio;
 241   Double_t sum = 0;
 242   Int_t N = 0;
 243
 244   Int_t offset = nbins - fNBinsToScale - 1;
 245   Int_t i;
 246   for ( i = offset; i< nbins; i++)
 247    {
 248     if ( fNumerator->GetBinContent(i) > 0.0 )
 249      {
 250        ratio = fDenominator->GetBinContent(i)/fNumerator->GetBinContent(i);
 251        sum += ratio;
 252        N++;
 253      }
 254    }
 255
 256   if(gDebug > 0) cout<<"sum="<<sum<<" fNBinsToScale="<<fNBinsToScale<<" N="<<N<<endl;
 257
 258   if (N == 0) return 0.0;
 259   Double_t ret = sum/((Double_t)N);
 260
 261   if(gDebug > 0) cout<<"Scale() returning "<<ret<<endl;
 262   return ret;
 263 }
 264
 265 /******************************************************************/
 266 /******************************************************************/
 267 /******************************************************************/
 268
 269 ClassImp( AliHBTTwoPartFctn2D )
 270
 271 AliHBTTwoPartFctn2D::
 272 AliHBTTwoPartFctn2D(Int_t nXbins, Double_t maxXval, Double_t minXval ,
 273                     Int_t nYbins, Double_t maxYval, Double_t minYval)
 274
 275 {
 276    TString numstr = fName + " Numerator";  //title and name of the
 277                                            //numerator histogram
 278    TString denstr = fName + " Denominator";//title and name of the
 279                                            //denominator histogram
 280
 281    fNumerator   = new TH2D(numstr.Data(),numstr.Data(),
 282                            nXbins,minXval,maxXval,
 283                nYbins,minYval,maxYval);
 284
 285    fDenominator = new TH2D(denstr.Data(),denstr.Data(),
 286                            nXbins,minXval,maxXval,
 287                nYbins,minYval,maxYval);
 288
 289    fNumerator->Sumw2();
 290    fDenominator->Sumw2();
 291
 292 }
 293 AliHBTTwoPartFctn2D::~AliHBTTwoPartFctn2D()
 294 {
 295   delete fNumerator;
 296   delete fDenominator;
 297 }
 298 void AliHBTTwoPartFctn2D::ProcessSameEventParticles(AliHBTPair* pair)
 299 {
 300   pair = CheckPair(pair);
 301   if(pair)
 302    {
 303      Double_t x,y;
 304      GetValues(pair,x,y);
 305      fNumerator->Fill(y,x);
 306    }
 307 }
 308
 309 void AliHBTTwoPartFctn2D::ProcessDiffEventParticles(AliHBTPair* pair)
 310 {
 311   pair = CheckPair(pair);
 312   if(pair)
 313    {
 314      Double_t x,y;
 315      GetValues(pair,x,y);
 316      fDenominator->Fill(y,x);
 317    }
 318
 319 }
 320
 321
 322 /******************************************************************/
 323 /******************************************************************/
 324 /******************************************************************/
 325
 326 ClassImp( AliHBTTwoPartFctn3D)
 327
 328 AliHBTTwoPartFctn3D::
 329 AliHBTTwoPartFctn3D(Int_t nXbins, Double_t maxXval, Double_t minXval,
 330                     Int_t nYbins, Double_t maxYval, Double_t minYval,
 331                     Int_t nZbins, Double_t maxZval, Double_t minZval)
 332
 333 {
 334    TString numstr = fName + " Numerator";  //title and name of the
 335                                            //numerator histogram
 336    TString denstr = fName + " Denominator";//title and name of the
 337                                            //denominator histogram
 338
 339    fNumerator   = new TH3D(numstr.Data(),numstr.Data(),
 340                            nXbins,minXval,maxXval,
 341                nYbins,minYval,maxYval,
 342                nZbins,minZval,maxZval);
 343
 344    fDenominator = new TH3D(denstr.Data(),denstr.Data(),
 345                            nXbins,minXval,maxXval,
 346                nYbins,minYval,maxYval,
 347                nZbins,minZval,maxZval);
 348
 349    fNumerator->Sumw2();
 350    fDenominator->Sumw2();
 351
 352 }
 353
 354
 355 AliHBTTwoPartFctn3D::~AliHBTTwoPartFctn3D()
 356 {
 357   delete fNumerator;
 358   delete fDenominator;
 359 }
 360
 361
 362 /******************************************************************/
 363 /******************************************************************/
 364 /******************************************************************/
 365 ClassImp( AliHBTFourPartFctn2D)
 366
 367
 368 AliHBTFourPartFctn2D::
 369 AliHBTFourPartFctn2D(Int_t nXbins, Double_t maxXval, Double_t minXval ,
 370                     Int_t nYbins, Double_t maxYval, Double_t minYval)
 371
 372 {
 373    TString numstr = fName + " Numerator";  //title and name of the
 374                                            //numerator histogram
 375    TString denstr = fName + " Denominator";//title and name of the
 376                                            //denominator histogram
 377
 378    fNumerator   = new TH2D(numstr.Data(),numstr.Data(),
 379                            nXbins,minXval,maxXval,
 380                nYbins,minYval,maxYval);
 381
 382    fDenominator = new TH2D(denstr.Data(),denstr.Data(),
 383                            nXbins,minXval,maxXval,
 384                nYbins,minYval,maxYval);
 385
 386    fNumerator->Sumw2();
 387    fDenominator->Sumw2();
 388
 389 }
 390 AliHBTFourPartFctn2D::~AliHBTFourPartFctn2D()
 391 {
 392   delete fNumerator;
 393   delete fDenominator;
 394 }
 395 void AliHBTFourPartFctn2D::
 396 ProcessSameEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair)
 397 {
 398   partpair  = CheckPair(partpair);
 399   trackpair = CheckPair(trackpair);
 400   if( partpair && trackpair)
 401    {
 402      Double_t x,y;
 403      GetValues(trackpair,partpair,x,y);
 404      fNumerator->Fill(y,x);
 405    }
 406 }
 407
 408 void AliHBTFourPartFctn2D::
 409 ProcessDiffEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair)
 410 {
 411   partpair  = CheckPair(partpair);
 412   trackpair = CheckPair(trackpair);
 413   if( partpair && trackpair)
 414    {
 415      Double_t x,y;
 416      GetValues(trackpair,partpair,x,y);
 417      fDenominator->Fill(y,x);
 418    }
 419
 420 }
 421