HBTAN/AliHBTFunction.cxx

   1 #include "AliHBTFunction.h"
   2
   3 /* $Id: */
   4
   5 //--------------------------------------------------------------------
   6 //AliHBTFunction
   7 //Author: Piotr Krzysztof Skowronski
   8 //Piotr.Skowronski@cern.ch
   9 //Base classes for HBT functions
  10 /*
  11  OnePairFctn       Function                  TwoPairFctn
  12    | \    \        /  |   \                      /|\
  13    |  \    \      /   |    \                    / | \
  14    |   \    \   1D   2D    3D                  /  |  \
  15    |    \    \  / \   |\   | \________________/__ |   \
  16    |     \    \/   \  | \__|_____________    /   \|    \
  17    |      \    \    \_|____|__           \  /     |     \
  18    |       \ /  \___  |    |  \           \/      |\     \
  19    |        /       \ |    |   \          /\      | \     \
  20    |       / \       \|    |    \        /  \     |  \     \
  21    |      /   \      /\    |     \      /    \    |   \     |
  22    |     /     \    /  \   |      \    /      \   |    \    |
  23    |    /       \  /    \  |       \  /        \  |     \   |
  24    |   /         \/      \ |        \/          \ |      \  |
  25  OnePair1D   OnePair2D  OnePair3D  TwoPair1D  TwoPair2D TwoPair3D
  26
  27
  28  four particle functions are intendent to be resolution functions:
  29  it is mecessary to have simulated particle pair corresponding to given
  30  recontructed track pair in order to calculate function simualted value
  31  and recontructed value, to be further histogrammed
  32 */
  33 //---------------------------------------------------------------------
  34
  35 #include <TH2.h>
  36 #include <TH3.h>
  37
  38 /******************************************************************/
  39 /******************************************************************/
  40
  41 ClassImp( AliHBTFunction )
  42
  43 AliHBTFunction::AliHBTFunction():
  44   fPairCut(new AliHBTEmptyPairCut())   //dummy cut
  45 {
  46 //Default constructor
  47 }
  48 /******************************************************************/
  49 AliHBTFunction::AliHBTFunction(const char* name,const char* title):
  50   TNamed(name,title),
  51   fPairCut(new AliHBTEmptyPairCut()) //dummy cut
  52 {
  53 //Constructor
  54 }
  55 /******************************************************************/
  56
  57 AliHBTFunction::~AliHBTFunction()
  58 {
  59 //destructor
  60   delete fPairCut;
  61 }
  62 /******************************************************************/
  63
  64 void AliHBTFunction::WriteFunction()
  65 {
  66 //writes result of the function to file
  67    if (GetNumerator()) GetNumerator()->Write();
  68    if (GetDenominator()) GetDenominator()->Write();
  69    TH1* res = GetResult();
  70    if (res) res->Write();
  71 }
  72 /******************************************************************/
  73
  74 TH1* AliHBTFunction::GetRatio(Double_t normfactor)
  75  {
  76  //returns ratio of numerator and denominator
  77  //
  78    if (gDebug>0) Info("GetRatio","Norm. Factor is %f for %s",normfactor,GetName());
  79
  80    if (normfactor == 0.0)
  81     {
  82       Error("GetRatio","Scaling Factor is 0. Null poiner returned");
  83       return 0x0;
  84     }
  85    TString str = fName + " ratio";
  86    TH1 *result = (TH1*)GetNumerator()->Clone(str.Data());
  87
  88    result->SetTitle(str.Data());
  89
  90    result->Divide(GetNumerator(),GetDenominator(),normfactor);
  91
  92    return result;
  93
  94  }
  95 /******************************************************************/
  96 void AliHBTFunction::SetPairCut(AliHBTPairCut* cut)
  97 {
  98 //Sets new Pair Cut. Old one is deleted
  99 //Note that it is created new object instead of simple pointer set
 100 //I do not want to have pointer
 101 //to object created somewhere else
 102 //because in that case I could not believe that
 103 //it would always exist (sb could delete it)
 104 //so we have always own copy
 105
 106  if(!cut)
 107    {
 108      Error("AliHBTFunction::SetPairCut","argument is NULL");
 109      return;
 110    }
 111  delete fPairCut;
 112  fPairCut = (AliHBTPairCut*)cut->Clone();
 113
 114 }
 115
 116 /******************************************************************/
 117
 118 void AliHBTFunction::Rename(const Char_t * name)
 119  {
 120  //renames the function and histograms
 121   SetName(name);
 122   SetTitle(name);
 123
 124   TString numstr = fName + " Numerator";  //title and name of the
 125                                            //numerator histogram
 126   TString denstr = fName + " Denominator";//title and name of the
 127                                            //denominator histogram
 128
 129   GetNumerator()->SetName(numstr.Data());
 130   GetNumerator()->SetTitle(numstr.Data());
 131
 132   GetDenominator()->SetName(denstr.Data());
 133   GetDenominator()->SetTitle(denstr.Data());
 134
 135  }
 136
 137 void AliHBTFunction::Rename(const Char_t * name, const Char_t * title)
 138  {
 139  //renames and retitle the function and histograms
 140
 141   SetName(name);
 142   SetTitle(title);
 143
 144   TString numstrn = fName + " Numerator";  //name of the
 145                                            //numerator histogram
 146
 147   TString numstrt = fTitle + " Numerator";  //title of the
 148                                            //numerator histogram
 149
 150   TString denstrn = fName + " Denominator";//name of the
 151                                            //denominator histogram
 152
 153   TString denstrt = fTitle + " Denominator";//title of the
 154                                            //denominator histogram
 155
 156
 157   GetNumerator()->SetName(numstrn.Data());
 158   GetNumerator()->SetTitle(numstrt.Data());
 159
 160   GetDenominator()->SetName(denstrn.Data());
 161   GetDenominator()->SetTitle(denstrt.Data());
 162
 163
 164  }
 165 /******************************************************************/
 166
 167 void AliHBTFunction::InitFunction()
 168 {
 169 //Iniotializes fctn.: Resets histograms
 170 //In case histograms are not created in ctor, builds with default parameters
 171   if ( !(GetNumerator()&&GetDenominator()) ) BuildHistos();
 172   GetNumerator()->Reset();
 173   GetDenominator()->Reset();
 174 }
 175
 176 /******************************************************************/
 177 /******************************************************************/
 178 /******************************************************************/
 179
 180 ClassImp( AliHBTOnePairFctn )
 181
 182 /******************************************************************/
 183 /******************************************************************/
 184 /******************************************************************/
 185
 186 ClassImp( AliHBTTwoPairFctn)
 187
 188 /******************************************************************/
 189 /******************************************************************/
 190 /******************************************************************/
 191
 192 ClassImp( AliHBTFunction1D )
 193
 194 const Int_t AliHBTFunction1D::fgkDefaultNBins = 100;
 195 const Float_t AliHBTFunction1D::fgkDefaultMin = 0.0;
 196 const Float_t AliHBTFunction1D::fgkDefaultMax = 0.15;
 197 const UInt_t AliHBTFunction1D::fgkDefaultNBinsToScale = 30;
 198
 199
 200 AliHBTFunction1D::AliHBTFunction1D():
 201  fNumerator(0x0),
 202  fDenominator(0x0),
 203  fNBinsToScale(fgkDefaultNBinsToScale)
 204 {//default constructor
 205 }
 206 /******************************************************************/
 207
 208 AliHBTFunction1D::AliHBTFunction1D(Int_t nbins, Float_t maxXval, Float_t minXval):
 209  fNumerator(0x0),
 210  fDenominator(0x0),
 211  fNBinsToScale(fgkDefaultNBinsToScale)
 212 {
 213  //Constructor of Two Part One Dimentional Function
 214  // nbins: number of bins in histograms - default 100
 215  // maxXval and minXval: range of histgram(s) default 0 - 0.15 (GeV)
 216  BuildHistos(nbins,maxXval,minXval);
 217 }
 218 /******************************************************************/
 219 AliHBTFunction1D::AliHBTFunction1D(const Char_t *name, const Char_t *title):
 220  AliHBTFunction(name,title),
 221  fNumerator(0x0),
 222  fDenominator(0x0),
 223  fNBinsToScale(fgkDefaultNBinsToScale)
 224 {//constructor
 225 }
 226 /******************************************************************/
 227 AliHBTFunction1D::AliHBTFunction1D(const Char_t *name, const Char_t *title,
 228                                    Int_t nbins, Float_t maxXval, Float_t minXval):
 229  AliHBTFunction(name,title),
 230  fNumerator(0x0),
 231  fDenominator(0x0),
 232  fNBinsToScale(fgkDefaultNBinsToScale)
 233 {
 234 //constructor
 235   BuildHistos(nbins,maxXval,minXval);
 236 }
 237 /******************************************************************/
 238
 239 AliHBTFunction1D::~AliHBTFunction1D()
 240 {
 241 //destructor
 242   delete fNumerator;
 243   delete fDenominator;
 244 }
 245 /******************************************************************/
 246 void AliHBTFunction1D::BuildHistos()
 247 {
 248 //builds histograms with default settings
 249  BuildHistos(fgkDefaultNBins,fgkDefaultMax,fgkDefaultMin);
 250 }
 251
 252 /******************************************************************/
 253
 254 void AliHBTFunction1D::BuildHistos(Int_t nbins, Float_t max, Float_t min)
 255 {
 256 //builds numarator and denominator hitograms
 257   TString numstr = fName + " Numerator";  //title and name of the
 258                                           //numerator histogram
 259   TString denstr = fName + " Denominator";//title and name of the
 260                                           //denominator histogram
 261   fNumerator   = new TH1D(numstr.Data(),numstr.Data(),nbins,min,max);
 262   fDenominator = new TH1D(denstr.Data(),denstr.Data(),nbins,min,max);
 263
 264   fNumerator->Sumw2();
 265   fDenominator->Sumw2();
 266 }
 267 /******************************************************************/
 268
 269 Double_t AliHBTFunction1D::Scale()
 270 {
 271  //Calculates the factor that should be used to scale
 272  //quatience of fNumerator and fDenominator to 1 at tail
 273  return Scale(fNumerator,fDenominator);
 274 }
 275 /******************************************************************/
 276
 277 Double_t AliHBTFunction1D::Scale(TH1D* num,TH1D* den)
 278 {
 279  //Calculates the factor that should be used to scale
 280  //quatience of num and den to 1 at tail
 281
 282   if (gDebug>0) Info("Scale","Enetered Scale()");
 283   if(!num)
 284    {
 285      Error("Scale","No numerator");
 286      return 0.0;
 287    }
 288   if(!den)
 289    {
 290      Error("Scale","No denominator");
 291      return 0.0;
 292    }
 293
 294   if(fNBinsToScale < 1)
 295    {
 296     return 0.0;
 297     Error("Scale","Number of bins for scaling is smaller thnan 1");
 298    }
 299   UInt_t nbins = num->GetNbinsX();
 300   if (fNBinsToScale > nbins)
 301    {
 302     Error("Scale","Number of bins for scaling is bigger thnan number of bins in histograms");
 303     return 0.0;
 304    }
 305   if (gDebug>0) Info("Scale","No errors detected");
 306
 307   Double_t ratio;
 308   Double_t sum = 0;
 309   Int_t n = 0;
 310
 311   Int_t offset = nbins - fNBinsToScale - 1;
 312
 313   for (UInt_t i = offset; i< nbins; i++)
 314    {
 315     if ( num->GetBinContent(i) > 0.0 )
 316      {
 317        ratio = den->GetBinContent(i)/num->GetBinContent(i);
 318        sum += ratio;
 319        n++;
 320      }
 321    }
 322
 323   if(gDebug > 0) Info("Scale","sum=%f fNBinsToScale=%d n=%d",sum,fNBinsToScale,n);
 324
 325   if (n == 0) return 0.0;
 326   Double_t ret = sum/((Double_t)n);
 327
 328   if(gDebug > 0) Info("Scale","returning %f",ret);
 329   return ret;
 330 }
 331
 332 /******************************************************************/
 333 /******************************************************************/
 334 /******************************************************************/
 335
 336 //____________________
 337 ///////////////////////////////////////////////////////
 338 //                                                   //
 339 // AliHBTFunction2D                                  //
 340 //                                                   //
 341 // Base Calss for 2-dimensinal Functions             //
 342 //                                                   //
 343 // Piotr.Skowronski@cern.ch                          //
 344 // http://alisoft.cern.ch/people/skowron/analyzer    //
 345 //                                                   //
 346 ///////////////////////////////////////////////////////
 347
 348 ClassImp( AliHBTFunction2D )
 349
 350 const Int_t AliHBTFunction2D::fgkDefaultNBinsX = 200;//default number of Bins in X axis in histograms
 351 const Float_t AliHBTFunction2D::fgkDefaultMinX = 0.0;//Default min value of X axis in histograms
 352 const Float_t AliHBTFunction2D::fgkDefaultMaxX = 1.5;//Default max value of X axis inhistograms
 353
 354 const Int_t AliHBTFunction2D::fgkDefaultNBinsY = 200;//default number of Bins in histograms
 355 const Float_t AliHBTFunction2D::fgkDefaultMinY = -0.15;//Default min value of histograms
 356 const Float_t AliHBTFunction2D::fgkDefaultMaxY =  0.15;//Default max value of histograms
 357
 358 const UInt_t AliHBTFunction2D::fgkDefaultNBinsToScaleX = 30;//Default number of X bins used for scaling to tale
 359 const UInt_t AliHBTFunction2D::fgkDefaultNBinsToScaleY = 30;//Default number of bins used for scaling to tale
 360
 361 /******************************************************************/
 362 AliHBTFunction2D::AliHBTFunction2D():
 363  fNumerator(0x0),
 364  fDenominator(0x0),
 365  fNBinsToScaleX(fgkDefaultNBinsToScaleX),
 366  fNBinsToScaleY(fgkDefaultNBinsToScaleY)
 367 {//default constructor
 368 }
 369 /******************************************************************/
 370 AliHBTFunction2D::AliHBTFunction2D(const Char_t *name, const Char_t *title):
 371  AliHBTFunction(name,title),
 372  fNumerator(0x0),
 373  fDenominator(0x0),
 374  fNBinsToScaleX(fgkDefaultNBinsToScaleX),
 375  fNBinsToScaleY(fgkDefaultNBinsToScaleY)
 376 {//constructor
 377 }
 378 /******************************************************************/
 379
 380 AliHBTFunction2D::AliHBTFunction2D(Int_t nXbins, Double_t maxXval, Double_t minXval,
 381                                    Int_t nYbins, Double_t maxYval, Double_t minYval):
 382  fNumerator(0x0),
 383  fDenominator(0x0),
 384  fNBinsToScaleX(fgkDefaultNBinsToScaleX),
 385  fNBinsToScaleY(fgkDefaultNBinsToScaleY)
 386 {
 387   BuildHistos(nXbins,maxXval,minXval,nYbins,maxYval,minYval);
 388 }
 389 /******************************************************************/
 390
 391 AliHBTFunction2D::AliHBTFunction2D(const Char_t *name, const Char_t *title,
 392                                    Int_t nXbins, Double_t maxXval, Double_t minXval,
 393                                    Int_t nYbins, Double_t maxYval, Double_t minYval):
 394  AliHBTFunction(name,title),
 395  fNumerator(0x0),
 396  fDenominator(0x0),
 397  fNBinsToScaleX(fgkDefaultNBinsToScaleX),
 398  fNBinsToScaleY(fgkDefaultNBinsToScaleY)
 399 {
 400   BuildHistos(nXbins,maxXval,minXval,nYbins,maxYval,minYval);
 401 }
 402 /******************************************************************/
 403
 404 AliHBTFunction2D::~AliHBTFunction2D()
 405 {
 406   delete fNumerator;
 407   delete fDenominator;
 408 }
 409 /******************************************************************/
 410
 411 void AliHBTFunction2D::BuildHistos()
 412 {
 413 //Creates default histograms
 414   BuildHistos(fgkDefaultNBinsX,fgkDefaultMaxX,fgkDefaultMinX,
 415               fgkDefaultNBinsY,fgkDefaultMaxY,fgkDefaultMinY);
 416 }
 417 /******************************************************************/
 418
 419 void AliHBTFunction2D::BuildHistos(Int_t nxbins, Float_t xmax, Float_t xmin,
 420                                    Int_t nybins, Float_t ymax, Float_t ymin)
 421 {
 422   //Builds numerator and denominator histograms (2d-case)
 423  TString numstr = fName + " Numerator";  //title and name of the
 424                                            //numerator histogram
 425  TString denstr = fName + " Denominator";//title and name of the
 426                                            //denominator histogram
 427
 428  fNumerator   = new TH2D(numstr.Data(),numstr.Data(),
 429                          nxbins,xmin,xmax,nybins,ymin,ymax);
 430
 431  fDenominator = new TH2D(denstr.Data(),denstr.Data(),
 432                          nxbins,xmin,xmax,nybins,ymin,ymax);
 433
 434  fNumerator->Sumw2();
 435  fDenominator->Sumw2();
 436 }
 437 /******************************************************************/
 438
 439 void AliHBTFunction2D::SetNumberOfBinsToScale(UInt_t xn, UInt_t yn)
 440 {
 441 //defines area used for scaling factor calculation
 442   fNBinsToScaleX = xn;
 443   fNBinsToScaleY = yn;
 444 }
 445 /******************************************************************/
 446
 447 Double_t AliHBTFunction2D::Scale()
 448 {
 449   // Calculates the factor that should be used to scale
 450   // quatience of fNumerator and fDenominator to 1 at
 451   // given region
 452   if (gDebug>0) Info("Scale","Enetered Scale()");
 453   if(!fNumerator)
 454    {
 455      Error("Scale","No numerator");
 456      return 0.0;
 457    }
 458   if(!fDenominator)
 459    {
 460      Error("Scale","No denominator");
 461      return 0.0;
 462    }
 463
 464   if( (fNBinsToScaleX < 1) || (fNBinsToScaleY < 1) )
 465    {
 466     return 0.0;
 467     Error("Scale","Number of bins for scaling is smaller thnan 1");
 468    }
 469   UInt_t nbinsX = fNumerator->GetNbinsX();
 470   if (fNBinsToScaleX > nbinsX)
 471    {
 472     Error("Scale","Number of X bins for scaling is bigger thnan number of bins in histograms");
 473     return 0.0;
 474    }
 475
 476   UInt_t nbinsY = fNumerator->GetNbinsX();
 477   if (fNBinsToScaleY > nbinsY)
 478    {
 479     Error("Scale","Number of Y bins for scaling is bigger thnan number of bins in histograms");
 480     return 0.0;
 481    }
 482
 483   if (gDebug>0) Info("Scale","No errors detected");
 484
 485   Int_t offsetX = nbinsX - fNBinsToScaleX - 1; //bin that we start loop over bins in axis X
 486   Int_t offsetY = nbinsY - fNBinsToScaleY - 1; //bin that we start loop over bins in axis X
 487
 488   Double_t ratio;
 489   Double_t sum = 0;
 490   Int_t n = 0;
 491
 492   for (UInt_t j = offsetY; j< nbinsY; j++)
 493     for (UInt_t i = offsetX; i< nbinsX; i++)
 494      {
 495       if ( fNumerator->GetBinContent(i,j) > 0.0 )
 496        {
 497          ratio = fDenominator->GetBinContent(i,j)/fNumerator->GetBinContent(i,j);
 498          sum += ratio;
 499          n++;
 500        }
 501      }
 502
 503   if(gDebug > 0) Info("Scale","sum=%f fNBinsToScaleX=%d fNBinsToScaleY=%d n=%d",sum,fNBinsToScaleX,fNBinsToScaleY,n);
 504
 505   if (n == 0) return 0.0;
 506   Double_t ret = sum/((Double_t)n);
 507
 508   if(gDebug > 0) Info("Scale","returning %f",ret);
 509   return ret;
 510 }
 511
 512 /******************************************************************/
 513 /******************************************************************/
 514 /******************************************************************/
 515
 516 //____________________
 517 ///////////////////////////////////////////////////////
 518 //                                                   //
 519 // AliHBTFunction3D                                  //
 520 //                                                   //
 521 // Base Calss for 3-dimensinal Functions             //
 522 //                                                   //
 523 // Piotr.Skowronski@cern.ch                          //
 524 // http://alisoft.cern.ch/people/skowron/analyzer    //
 525 //                                                   //
 526 ///////////////////////////////////////////////////////
 527
 528 ClassImp( AliHBTFunction3D)
 529
 530 const Int_t AliHBTFunction3D::fgkDefaultNBinsX = 200;//default number of Bins in X axis in histograms
 531 const Float_t AliHBTFunction3D::fgkDefaultMinX = 0.0;//Default min value of X axis in histograms
 532 const Float_t AliHBTFunction3D::fgkDefaultMaxX = 1.5;//Default max value of X axis inhistograms
 533
 534 const Int_t AliHBTFunction3D::fgkDefaultNBinsY = 200;//default number of Bins in histograms
 535 const Float_t AliHBTFunction3D::fgkDefaultMinY = -0.15;//Default min value of histograms
 536 const Float_t AliHBTFunction3D::fgkDefaultMaxY =  0.15;//Default max value of histograms
 537
 538 const Int_t AliHBTFunction3D::fgkDefaultNBinsZ = 200;//default number of Bins in histograms
 539 const Float_t AliHBTFunction3D::fgkDefaultMinZ = -0.15;//Default min value of histograms
 540 const Float_t AliHBTFunction3D::fgkDefaultMaxZ =  0.15;//Default max value of histograms
 541
 542 const UInt_t AliHBTFunction3D::fgkDefaultNBinsToScaleX = 30;//Default number of X bins used for scaling to tale
 543 const UInt_t AliHBTFunction3D::fgkDefaultNBinsToScaleY = 30;//Default number of bins used for scaling to tale
 544 const UInt_t AliHBTFunction3D::fgkDefaultNBinsToScaleZ = 30;//Default number of bins used for scaling to tale
 545
 546 AliHBTFunction3D::AliHBTFunction3D():
 547  fNumerator(0x0),
 548  fDenominator(0x0),
 549  fNBinsToScaleX(fgkDefaultNBinsToScaleX),
 550  fNBinsToScaleY(fgkDefaultNBinsToScaleY),
 551  fNBinsToScaleZ(fgkDefaultNBinsToScaleZ)
 552 {
 553  //constructor
 554 }
 555 /******************************************************************/
 556
 557 AliHBTFunction3D::AliHBTFunction3D(const Char_t *name, const Char_t *title):
 558  AliHBTFunction(name,title),
 559  fNumerator(0x0),
 560  fDenominator(0x0),
 561  fNBinsToScaleX(fgkDefaultNBinsToScaleX),
 562  fNBinsToScaleY(fgkDefaultNBinsToScaleY),
 563  fNBinsToScaleZ(fgkDefaultNBinsToScaleZ)
 564 {
 565  //constructor
 566 }
 567 /******************************************************************/
 568
 569 AliHBTFunction3D::AliHBTFunction3D(Int_t nXbins, Double_t maxXval, Double_t minXval,
 570                   Int_t nYbins, Double_t maxYval, Double_t minYval,
 571                   Int_t nZbins, Double_t maxZval, Double_t minZval):
 572  fNumerator(0x0),
 573  fDenominator(0x0),
 574  fNBinsToScaleX(fgkDefaultNBinsToScaleX),
 575  fNBinsToScaleY(fgkDefaultNBinsToScaleY),
 576  fNBinsToScaleZ(fgkDefaultNBinsToScaleZ)
 577 {
 578  //constructor
 579  BuildHistos( nXbins,maxXval,minXval,nYbins,maxYval,minYval,nZbins,maxZval,minZval);
 580 }
 581 /******************************************************************/
 582
 583 AliHBTFunction3D::AliHBTFunction3D(const Char_t *name, const Char_t *title,
 584                  Int_t nXbins, Double_t maxXval, Double_t minXval,
 585                  Int_t nYbins, Double_t maxYval, Double_t minYval,
 586                  Int_t nZbins, Double_t maxZval, Double_t minZval):
 587  AliHBTFunction(name,title),
 588  fNumerator(0x0),
 589  fDenominator(0x0),
 590  fNBinsToScaleX(fgkDefaultNBinsToScaleX),
 591  fNBinsToScaleY(fgkDefaultNBinsToScaleY),
 592  fNBinsToScaleZ(fgkDefaultNBinsToScaleZ)
 593 {
 594  //constructor
 595  BuildHistos( nXbins,maxXval,minXval,nYbins,maxYval,minYval,nZbins,maxZval,minZval);
 596 }
 597 /******************************************************************/
 598
 599
 600 AliHBTFunction3D::~AliHBTFunction3D()
 601 {
 602   delete fNumerator;
 603   delete fDenominator;
 604 }
 605 /******************************************************************/
 606
 607 void AliHBTFunction3D::BuildHistos()
 608 {
 609 //Creates default histograms
 610   BuildHistos(fgkDefaultNBinsX,fgkDefaultMaxX,fgkDefaultMinX,
 611               fgkDefaultNBinsY,fgkDefaultMaxY,fgkDefaultMinY,
 612               fgkDefaultNBinsZ,fgkDefaultMaxZ,fgkDefaultMinZ);
 613 }
 614 /******************************************************************/
 615
 616 void AliHBTFunction3D::BuildHistos(Int_t nxbins, Float_t xmax, Float_t xmin,
 617                                    Int_t nybins, Float_t ymax, Float_t ymin,
 618                        Int_t nzbins, Float_t zmax, Float_t zmin)
 619 {
 620   //Builds numerator and denominator histograms (3d-case)
 621    TString numstr = fName + " Numerator";  //title and name of the
 622                                            //numerator histogram
 623    TString denstr = fName + " Denominator";//title and name of the
 624                                            //denominator histogram
 625
 626    fNumerator   = new TH3D(numstr.Data(),numstr.Data(),
 627                            nxbins,xmin,xmin,nybins,ymin,ymin,nzbins,zmin,zmin);
 628
 629    fDenominator = new TH3D(denstr.Data(),denstr.Data(),
 630                            nxbins,xmin,xmin,nybins,ymin,ymin,nzbins,zmin,zmin);
 631
 632    fNumerator->Sumw2();
 633    fDenominator->Sumw2();
 634 }
 635
 636 /******************************************************************/
 637
 638 Double_t AliHBTFunction3D::Scale()
 639 {
 640   // Calculates the factor that should be used to scale
 641   // quatience of fNumerator and fDenominator to 1 at
 642   // given volume
 643   if (gDebug>0) Info("Scale","Enetered Scale()");
 644   if(!fNumerator)
 645    {
 646      Error("Scale","No numerator");
 647      return 0.0;
 648    }
 649   if(!fDenominator)
 650    {
 651      Error("Scale","No denominator");
 652      return 0.0;
 653    }
 654
 655   if( (fNBinsToScaleX < 1) || (fNBinsToScaleY < 1) || (fNBinsToScaleZ < 1))
 656    {
 657     return 0.0;
 658     Error("Scale","Number of bins for scaling is smaller thnan 1");
 659    }
 660   UInt_t nbinsX = fNumerator->GetNbinsX();
 661   if (fNBinsToScaleX > nbinsX)
 662    {
 663     Error("Scale","Number of X bins for scaling is bigger thnan number of bins in histograms");
 664     return 0.0;
 665    }
 666
 667   UInt_t nbinsY = fNumerator->GetNbinsX();
 668   if (fNBinsToScaleY > nbinsY)
 669    {
 670     Error("Scale","Number of Y bins for scaling is bigger thnan number of bins in histograms");
 671     return 0.0;
 672    }
 673
 674   UInt_t nbinsZ = fNumerator->GetNbinsZ();
 675   if (fNBinsToScaleZ > nbinsZ)
 676    {
 677     Error("Scale","Number of Z bins for scaling is bigger thnan number of bins in histograms");
 678     return 0.0;
 679    }
 680
 681   if (gDebug>0) Info("Scale","No errors detected");
 682
 683   Int_t offsetX = nbinsX - fNBinsToScaleX - 1; //bin that we start loop over bins in axis X
 684   Int_t offsetY = nbinsY - fNBinsToScaleY - 1; //bin that we start loop over bins in axis Y
 685   Int_t offsetZ = nbinsZ - fNBinsToScaleZ - 1; //bin that we start loop over bins in axis Z
 686
 687   Double_t ratio;
 688   Double_t sum = 0;
 689   Int_t n = 0;
 690
 691   for (UInt_t k = offsetZ; k<nbinsZ; k++)
 692     for (UInt_t j = offsetY; j<nbinsY; j++)
 693       for (UInt_t i = offsetX; i<nbinsX; i++)
 694        {
 695         if ( fNumerator->GetBinContent(i,j,k) > 0.0 )
 696          {
 697            ratio = fDenominator->GetBinContent(i,j,k)/fNumerator->GetBinContent(i,j,k);
 698            sum += ratio;
 699            n++;
 700          }
 701        }
 702
 703   if(gDebug > 0)
 704     Info("Scale","sum=%f fNBinsToScaleX=%d fNBinsToScaleY=%d fNBinsToScaleZ=%d n=%d",
 705           sum,fNBinsToScaleX,fNBinsToScaleY,fNBinsToScaleZ,n);
 706
 707   if (n == 0) return 0.0;
 708   Double_t ret = sum/((Double_t)n);
 709
 710   if(gDebug > 0) Info("Scale","returning %f",ret);
 711   return ret;
 712 }
 713 /******************************************************************/
 714
 715 void AliHBTFunction3D::SetNumberOfBinsToScale(UInt_t xn, UInt_t yn,UInt_t zn)
 716 {
 717 //sets up the volume to be used for scaling to tail
 718  fNBinsToScaleX = xn;
 719  fNBinsToScaleY = yn;
 720  fNBinsToScaleZ = zn;
 721 }
 722
 723 /******************************************************************/
 724 /******************************************************************/
 725 /******************************************************************/
 726 /******************************************************************/
 727 /******************************************************************/
 728 /******************************************************************/
 729 /******************************************************************/
 730 /******************************************************************/
 731 /******************************************************************/
 732
 733 //____________________
 734 ///////////////////////////////////////////////////////
 735 //                                                   //
 736 // AliHBTOnePairFctn1D                               //
 737 //                                                   //
 738 // Base Calss for 1-dimensinal Functions that need   //
 739 // one pair to fill function                         //
 740 //                                                   //
 741 // Piotr.Skowronski@cern.ch                          //
 742 // http://alisoft.cern.ch/people/skowron/analyzer    //
 743 //                                                   //
 744 ///////////////////////////////////////////////////////
 745
 746 ClassImp( AliHBTOnePairFctn1D )
 747 /******************************************************************/
 748
 749 AliHBTOnePairFctn1D::AliHBTOnePairFctn1D(Int_t nbins, Float_t maxXval, Float_t minXval):
 750  AliHBTFunction1D(nbins,maxXval,minXval)
 751 {
 752 }
 753 /******************************************************************/
 754
 755 AliHBTOnePairFctn1D::AliHBTOnePairFctn1D(const Char_t *name, const Char_t *title):
 756  AliHBTFunction1D(name,title)
 757 {
 758 }
 759 /******************************************************************/
 760
 761 AliHBTOnePairFctn1D::AliHBTOnePairFctn1D(const Char_t *name, const Char_t *title,
 762                                          Int_t nbins, Float_t maxXval, Float_t minXval):
 763  AliHBTFunction1D(name,title,nbins,maxXval,minXval)
 764 {
 765 }
 766 /******************************************************************/
 767
 768 void AliHBTOnePairFctn1D::ProcessSameEventParticles(AliHBTPair* pair)
 769 {
 770  //Fills the numerator using pair from the same event
 771    pair = CheckPair(pair);
 772    if(pair) fNumerator->Fill(GetValue(pair));
 773 }
 774 /******************************************************************/
 775 void AliHBTOnePairFctn1D::ProcessDiffEventParticles(AliHBTPair* pair)
 776  {
 777   //Fills the denumerator using mixed pairs
 778    pair = CheckPair(pair);
 779    if(pair) fDenominator->Fill(GetValue(pair));
 780   }
 781
 782 /******************************************************************/
 783 /******************************************************************/
 784 /******************************************************************/
 785
 786 //____________________
 787 ///////////////////////////////////////////////////////
 788 //                                                   //
 789 // AliHBTOnePairFctn2D                               //
 790 //                                                   //
 791 // Base Calss for 2-dimensinal Functions that need   //
 792 // one pair to fill function                         //
 793 //                                                   //
 794 // Piotr.Skowronski@cern.ch                          //
 795 // http://alisoft.cern.ch/people/skowron/analyzer    //
 796 //                                                   //
 797 ///////////////////////////////////////////////////////
 798
 799 ClassImp( AliHBTOnePairFctn2D )
 800 /******************************************************************/
 801
 802 AliHBTOnePairFctn2D::AliHBTOnePairFctn2D(const Char_t *name, const Char_t *title):
 803  AliHBTFunction2D(name,title)
 804 {
 805 }
 806 /******************************************************************/
 807
 808 AliHBTOnePairFctn2D::AliHBTOnePairFctn2D(Int_t nXbins, Double_t maxXval, Double_t minXval,
 809                       Int_t nYbins, Double_t maxYval, Double_t minYval):
 810  AliHBTFunction2D(nXbins, maxXval, minXval, nYbins, maxYval, minYval)
 811 {
 812 }
 813 /******************************************************************/
 814
 815 AliHBTOnePairFctn2D::AliHBTOnePairFctn2D(const Char_t *name, const Char_t *title,
 816                       Int_t nXbins, Double_t maxXval, Double_t minXval,
 817                       Int_t nYbins, Double_t maxYval, Double_t minYval):
 818  AliHBTFunction2D(name,title, nXbins, maxXval, minXval, nYbins, maxYval, minYval)
 819 {
 820 }
 821 /******************************************************************/
 822
 823 void AliHBTOnePairFctn2D::ProcessSameEventParticles(AliHBTPair* pair)
 824 {
 825   // Fills the numerator using pairs from the same event
 826   pair = CheckPair(pair);
 827   if(pair)
 828    {
 829      Double_t x,y;
 830      GetValues(pair,x,y);
 831      fNumerator->Fill(x,y);
 832    }
 833 }
 834 /******************************************************************/
 835
 836 void AliHBTOnePairFctn2D::ProcessDiffEventParticles(AliHBTPair* pair)
 837 {
 838   // Fills the denumerator using mixed pairs
 839   pair = CheckPair(pair);
 840   if(pair)
 841    {
 842      Double_t x,y;
 843      GetValues(pair,x,y);
 844      fDenominator->Fill(x,y);
 845    }
 846 }
 847 /******************************************************************/
 848 /******************************************************************/
 849 /******************************************************************/
 850 /******************************************************************/
 851 //____________________
 852 ///////////////////////////////////////////////////////
 853 //                                                   //
 854 // AliHBTOnePairFctn3D                               //
 855 //                                                   //
 856 // Base Calss for 3-dimensinal Functions that need   //
 857 // one pair to fill function                         //
 858 //                                                   //
 859 // Piotr.Skowronski@cern.ch                          //
 860 // http://alisoft.cern.ch/people/skowron/analyzer    //
 861 //                                                   //
 862 ///////////////////////////////////////////////////////
 863 ClassImp( AliHBTOnePairFctn3D)
 864
 865 /******************************************************************/
 866
 867 AliHBTOnePairFctn3D::AliHBTOnePairFctn3D(const Char_t *name, const Char_t *title):
 868  AliHBTFunction3D(name,title)
 869 {
 870 }
 871 /******************************************************************/
 872
 873 AliHBTOnePairFctn3D::AliHBTOnePairFctn3D(Int_t nXbins, Double_t maxXval, Double_t minXval,
 874                       Int_t nYbins, Double_t maxYval, Double_t minYval,
 875                       Int_t nZbins, Double_t maxZval, Double_t minZval):
 876  AliHBTFunction3D(nXbins, maxXval, minXval, nYbins, maxYval, minYval, nZbins, maxZval, minZval)
 877 {
 878 }
 879 /******************************************************************/
 880
 881 AliHBTOnePairFctn3D::AliHBTOnePairFctn3D(const Char_t *name, const Char_t *title,
 882                       Int_t nXbins, Double_t maxXval, Double_t minXval,
 883                       Int_t nYbins, Double_t maxYval, Double_t minYval,
 884                       Int_t nZbins, Double_t maxZval, Double_t minZval):
 885  AliHBTFunction3D(name,title,nXbins, maxXval, minXval, nYbins, maxYval, minYval, nZbins, maxZval, minZval)
 886 {
 887 }
 888 /******************************************************************/
 889
 890 void AliHBTOnePairFctn3D::ProcessSameEventParticles(AliHBTPair* pair)
 891 {
 892 //Reacts on pair coming from same event (real pairs)
 893 //and fills numerator histogram
 894   pair  = CheckPair(pair);
 895   if( pair )
 896    {
 897      Double_t x,y,z;
 898      GetValues(pair,x,y,z);
 899      fNumerator->Fill(x,y,z);
 900    }
 901 }
 902 /******************************************************************/
 903
 904 void AliHBTOnePairFctn3D::ProcessDiffEventParticles(AliHBTPair* pair)
 905 {
 906 //Reacts on pair coming from different events (mixed pairs)
 907 //and fills denominator histogram
 908   pair  = CheckPair(pair);
 909   if( pair )
 910    {
 911      Double_t x,y,z;
 912      GetValues(pair,x,y,z);
 913      fDenominator->Fill(x,y,z);
 914    }
 915 }
 916
 917 /******************************************************************/
 918 /******************************************************************/
 919 /******************************************************************/
 920
 921 //____________________
 922 ///////////////////////////////////////////////////////
 923 //                                                   //
 924 // AliHBTTwoPairFctn1D                               //
 925 //                                                   //
 926 // Base Calss for 1-dimensinal Functions that need   //
 927 // two pair (simulated and reconstructed)            //
 928 // to fill function                                  //
 929 //                                                   //
 930 // Piotr.Skowronski@cern.ch                          //
 931 // http://alisoft.cern.ch/people/skowron/analyzer    //
 932 //                                                   //
 933 ///////////////////////////////////////////////////////
 934 ClassImp(AliHBTTwoPairFctn1D)
 935 /******************************************************************/
 936
 937 AliHBTTwoPairFctn1D::AliHBTTwoPairFctn1D(Int_t nbins, Float_t maxXval, Float_t minXval):
 938  AliHBTFunction1D(nbins,maxXval,minXval)
 939 {
 940 }
 941 /******************************************************************/
 942
 943 AliHBTTwoPairFctn1D::AliHBTTwoPairFctn1D(const Char_t *name, const Char_t *title):
 944  AliHBTFunction1D(name,title)
 945 {
 946 }
 947 /******************************************************************/
 948
 949 AliHBTTwoPairFctn1D::AliHBTTwoPairFctn1D(const Char_t *name, const Char_t *title,
 950                                          Int_t nbins, Float_t maxXval, Float_t minXval):
 951  AliHBTFunction1D(name,title,nbins,maxXval,minXval)
 952 {
 953 }
 954 /******************************************************************/
 955
 956 void AliHBTTwoPairFctn1D::ProcessSameEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair)
 957 {
 958   // Fills the numerator using pairs from the same event
 959   partpair  = CheckPair(partpair);
 960   if( partpair )
 961    {
 962      Double_t x = GetValue(trackpair,partpair);
 963      fNumerator->Fill(x);
 964    }
 965 }
 966 /******************************************************************/
 967
 968 void AliHBTTwoPairFctn1D::ProcessDiffEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair)
 969 {
 970   // Fills the denumerator usin mixed pairs
 971   partpair  = CheckPair(partpair);
 972   if( partpair )
 973    {
 974      Double_t x = GetValue(trackpair,partpair);
 975      fDenominator->Fill(x);
 976    }
 977 }
 978 /******************************************************************/
 979 /******************************************************************/
 980 /******************************************************************/
 981
 982 //____________________
 983 ///////////////////////////////////////////////////////
 984 //                                                   //
 985 // AliHBTTwoPairFctn2D                               //
 986 //                                                   //
 987 // Base Calss for 2-dimensinal Functions that need   //
 988 // two pair (simulated and reconstructed)            //
 989 // to fill function                                  //
 990 //                                                   //
 991 // Piotr.Skowronski@cern.ch                          //
 992 // http://alisoft.cern.ch/people/skowron/analyzer    //
 993 //                                                   //
 994 ///////////////////////////////////////////////////////
 995
 996 ClassImp(AliHBTTwoPairFctn2D)
 997
 998 /******************************************************************/
 999
1000 AliHBTTwoPairFctn2D::AliHBTTwoPairFctn2D(const Char_t *name, const Char_t *title):
1001  AliHBTFunction2D(name,title)
1002 {
1003 }
1004 /******************************************************************/
1005
1006 AliHBTTwoPairFctn2D::AliHBTTwoPairFctn2D(Int_t nXbins, Double_t maxXval, Double_t minXval,
1007                       Int_t nYbins, Double_t maxYval, Double_t minYval):
1008  AliHBTFunction2D(nXbins, maxXval, minXval, nYbins, maxYval, minYval)
1009 {
1010 }
1011 /******************************************************************/
1012
1013 AliHBTTwoPairFctn2D::AliHBTTwoPairFctn2D(const Char_t *name, const Char_t *title,
1014                       Int_t nXbins, Double_t maxXval, Double_t minXval,
1015                       Int_t nYbins, Double_t maxYval, Double_t minYval):
1016  AliHBTFunction2D(name,title,nXbins, maxXval, minXval, nYbins, maxYval, minYval)
1017 {
1018 }
1019 /******************************************************************/
1020
1021 void AliHBTTwoPairFctn2D::ProcessSameEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair)
1022 {
1023 //processes pair of particles coming from a same events (real pair)
1024   partpair  = CheckPair(partpair);  //check cuts
1025   if( partpair )
1026    {
1027      Double_t x,y;
1028      GetValues(trackpair,partpair,x,y);
1029      fNumerator->Fill(x,y);
1030    }
1031 }
1032 /******************************************************************/
1033
1034 void AliHBTTwoPairFctn2D::ProcessDiffEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair)
1035 {
1036 //processes pair of particles coming from a different events (mixed pair)
1037   partpair  = CheckPair(partpair);
1038   if( partpair )
1039    {
1040      Double_t x,y;
1041      GetValues(trackpair,partpair,x,y);
1042      fDenominator->Fill(x,y);
1043    }
1044 }
1045
1046 /******************************************************************/
1047 /******************************************************************/
1048 /******************************************************************/
1049
1050 //____________________
1051 ///////////////////////////////////////////////////////
1052 //                                                   //
1053 // AliHBTTwoPairFctn3D                               //
1054 //                                                   //
1055 // Base Calss for 3-dimensinal Functions that need   //
1056 // two pair (simulated and reconstructed)            //
1057 // to fill function                                  //
1058 //                                                   //
1059 // Piotr.Skowronski@cern.ch                          //
1060 // http://alisoft.cern.ch/people/skowron/analyzer    //
1061 //                                                   //
1062 ///////////////////////////////////////////////////////
1063
1064 ClassImp(AliHBTTwoPairFctn3D)
1065
1066 /******************************************************************/
1067
1068 AliHBTTwoPairFctn3D::AliHBTTwoPairFctn3D(const Char_t *name, const Char_t *title):
1069  AliHBTFunction3D(name,title)
1070 {
1071 }
1072 /******************************************************************/
1073
1074 AliHBTTwoPairFctn3D::AliHBTTwoPairFctn3D(Int_t nXbins, Double_t maxXval, Double_t minXval,
1075                       Int_t nYbins, Double_t maxYval, Double_t minYval,
1076                       Int_t nZbins, Double_t maxZval, Double_t minZval):
1077  AliHBTFunction3D(nXbins, maxXval, minXval, nYbins, maxYval, minYval, nZbins, maxZval, minZval)
1078 {
1079 }
1080 /******************************************************************/
1081
1082 AliHBTTwoPairFctn3D::AliHBTTwoPairFctn3D(const Char_t *name, const Char_t *title,
1083                       Int_t nXbins, Double_t maxXval, Double_t minXval,
1084                       Int_t nYbins, Double_t maxYval, Double_t minYval,
1085                       Int_t nZbins, Double_t maxZval, Double_t minZval):
1086  AliHBTFunction3D(name,title,nXbins, maxXval, minXval, nYbins, maxYval, minYval, nZbins, maxZval, minZval)
1087 {
1088 }
1089 /******************************************************************/
1090
1091 void AliHBTTwoPairFctn3D::ProcessSameEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair)
1092 {
1093   // Fills th numerator using pairs from the same event
1094   partpair  = CheckPair(partpair);
1095   if( partpair )
1096    {
1097      Double_t x,y,z;
1098      GetValues(trackpair,partpair,x,y,z);
1099      fNumerator->Fill(x,y,z);
1100    }
1101 }
1102 /******************************************************************/
1103
1104 void AliHBTTwoPairFctn3D::ProcessDiffEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair)
1105 {
1106   // Fills the denumerator using mixed pairs
1107   partpair  = CheckPair(partpair);
1108   if( partpair )
1109    {
1110      Double_t x,y,z;
1111      GetValues(trackpair,partpair,x,y,z);
1112      fDenominator->Fill(x,y,z);
1113    }
1114
1115 }
1116
1117
1118 /******************************************************************/
1119 /******************************************************************/
1120 /******************************************************************/
1121