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     Error("Scale","Number of bins for scaling is smaller thnan 1");
 297     return 0.0;
 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     Error("Scale","Number of bins for scaling is smaller thnan 1");
 467     return 0.0;
 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     Error("Scale","Number of bins for scaling is smaller thnan 1");
 658     return 0.0;
 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   //constructor
 753 }
 754 /******************************************************************/
 755
 756 AliHBTOnePairFctn1D::AliHBTOnePairFctn1D(const Char_t *name, const Char_t *title):
 757  AliHBTFunction1D(name,title)
 758 {
 759
 760 }
 761 /******************************************************************/
 762
 763 AliHBTOnePairFctn1D::AliHBTOnePairFctn1D(const Char_t *name, const Char_t *title,
 764                                          Int_t nbins, Float_t maxXval, Float_t minXval):
 765  AliHBTFunction1D(name,title,nbins,maxXval,minXval)
 766 {
 767   //constructor
 768 }
 769 /******************************************************************/
 770
 771 void AliHBTOnePairFctn1D::ProcessSameEventParticles(AliHBTPair* pair)
 772 {
 773  //Fills the numerator using pair from the same event
 774    pair = CheckPair(pair);
 775    if(pair) fNumerator->Fill(GetValue(pair));
 776 }
 777 /******************************************************************/
 778 void AliHBTOnePairFctn1D::ProcessDiffEventParticles(AliHBTPair* pair)
 779  {
 780   //Fills the denumerator using mixed pairs
 781    pair = CheckPair(pair);
 782    if(pair) fDenominator->Fill(GetValue(pair));
 783   }
 784
 785 /******************************************************************/
 786 /******************************************************************/
 787 /******************************************************************/
 788
 789 //____________________
 790 ///////////////////////////////////////////////////////
 791 //                                                   //
 792 // AliHBTOnePairFctn2D                               //
 793 //                                                   //
 794 // Base Calss for 2-dimensinal Functions that need   //
 795 // one pair to fill function                         //
 796 //                                                   //
 797 // Piotr.Skowronski@cern.ch                          //
 798 // http://alisoft.cern.ch/people/skowron/analyzer    //
 799 //                                                   //
 800 ///////////////////////////////////////////////////////
 801
 802 ClassImp( AliHBTOnePairFctn2D )
 803 /******************************************************************/
 804
 805 AliHBTOnePairFctn2D::AliHBTOnePairFctn2D(const Char_t *name, const Char_t *title):
 806  AliHBTFunction2D(name,title)
 807 {
 808   //constructor
 809 }
 810 /******************************************************************/
 811
 812 AliHBTOnePairFctn2D::AliHBTOnePairFctn2D(Int_t nXbins, Double_t maxXval, Double_t minXval,
 813                       Int_t nYbins, Double_t maxYval, Double_t minYval):
 814  AliHBTFunction2D(nXbins, maxXval, minXval, nYbins, maxYval, minYval)
 815 {
 816   //constructor
 817 }
 818 /******************************************************************/
 819
 820 AliHBTOnePairFctn2D::AliHBTOnePairFctn2D(const Char_t *name, const Char_t *title,
 821                       Int_t nXbins, Double_t maxXval, Double_t minXval,
 822                       Int_t nYbins, Double_t maxYval, Double_t minYval):
 823  AliHBTFunction2D(name,title, nXbins, maxXval, minXval, nYbins, maxYval, minYval)
 824 {
 825   //constructor
 826 }
 827 /******************************************************************/
 828
 829 void AliHBTOnePairFctn2D::ProcessSameEventParticles(AliHBTPair* pair)
 830 {
 831   // Fills the numerator using pairs from the same event
 832   pair = CheckPair(pair);
 833   if(pair)
 834    {
 835      Double_t x,y;
 836      GetValues(pair,x,y);
 837      fNumerator->Fill(x,y);
 838    }
 839 }
 840 /******************************************************************/
 841
 842 void AliHBTOnePairFctn2D::ProcessDiffEventParticles(AliHBTPair* pair)
 843 {
 844   // Fills the denumerator using mixed pairs
 845   pair = CheckPair(pair);
 846   if(pair)
 847    {
 848      Double_t x,y;
 849      GetValues(pair,x,y);
 850      fDenominator->Fill(x,y);
 851    }
 852 }
 853 /******************************************************************/
 854 /******************************************************************/
 855 /******************************************************************/
 856 /******************************************************************/
 857 //____________________
 858 ///////////////////////////////////////////////////////
 859 //                                                   //
 860 // AliHBTOnePairFctn3D                               //
 861 //                                                   //
 862 // Base Calss for 3-dimensinal Functions that need   //
 863 // one pair to fill function                         //
 864 //                                                   //
 865 // Piotr.Skowronski@cern.ch                          //
 866 // http://alisoft.cern.ch/people/skowron/analyzer    //
 867 //                                                   //
 868 ///////////////////////////////////////////////////////
 869 ClassImp( AliHBTOnePairFctn3D)
 870
 871 /******************************************************************/
 872
 873 AliHBTOnePairFctn3D::AliHBTOnePairFctn3D(const Char_t *name, const Char_t *title):
 874  AliHBTFunction3D(name,title)
 875 {
 876   //constructor
 877 }
 878 /******************************************************************/
 879
 880 AliHBTOnePairFctn3D::AliHBTOnePairFctn3D(Int_t nXbins, Double_t maxXval, Double_t minXval,
 881                       Int_t nYbins, Double_t maxYval, Double_t minYval,
 882                       Int_t nZbins, Double_t maxZval, Double_t minZval):
 883  AliHBTFunction3D(nXbins, maxXval, minXval, nYbins, maxYval, minYval, nZbins, maxZval, minZval)
 884 {
 885   //constructor
 886 }
 887 /******************************************************************/
 888
 889 AliHBTOnePairFctn3D::AliHBTOnePairFctn3D(const Char_t *name, const Char_t *title,
 890                       Int_t nXbins, Double_t maxXval, Double_t minXval,
 891                       Int_t nYbins, Double_t maxYval, Double_t minYval,
 892                       Int_t nZbins, Double_t maxZval, Double_t minZval):
 893  AliHBTFunction3D(name,title,nXbins, maxXval, minXval, nYbins, maxYval, minYval, nZbins, maxZval, minZval)
 894 {
 895   //constructor
 896 }
 897 /******************************************************************/
 898
 899 void AliHBTOnePairFctn3D::ProcessSameEventParticles(AliHBTPair* pair)
 900 {
 901 //Reacts on pair coming from same event (real pairs)
 902 //and fills numerator histogram
 903   pair  = CheckPair(pair);
 904   if( pair )
 905    {
 906      Double_t x,y,z;
 907      GetValues(pair,x,y,z);
 908      fNumerator->Fill(x,y,z);
 909    }
 910 }
 911 /******************************************************************/
 912
 913 void AliHBTOnePairFctn3D::ProcessDiffEventParticles(AliHBTPair* pair)
 914 {
 915 //Reacts on pair coming from different events (mixed pairs)
 916 //and fills denominator histogram
 917   pair  = CheckPair(pair);
 918   if( pair )
 919    {
 920      Double_t x,y,z;
 921      GetValues(pair,x,y,z);
 922      fDenominator->Fill(x,y,z);
 923    }
 924 }
 925
 926 /******************************************************************/
 927 /******************************************************************/
 928 /******************************************************************/
 929
 930 //____________________
 931 ///////////////////////////////////////////////////////
 932 //                                                   //
 933 // AliHBTTwoPairFctn1D                               //
 934 //                                                   //
 935 // Base Calss for 1-dimensinal Functions that need   //
 936 // two pair (simulated and reconstructed)            //
 937 // to fill function                                  //
 938 //                                                   //
 939 // Piotr.Skowronski@cern.ch                          //
 940 // http://alisoft.cern.ch/people/skowron/analyzer    //
 941 //                                                   //
 942 ///////////////////////////////////////////////////////
 943 ClassImp(AliHBTTwoPairFctn1D)
 944 /******************************************************************/
 945
 946 AliHBTTwoPairFctn1D::AliHBTTwoPairFctn1D(Int_t nbins, Float_t maxXval, Float_t minXval):
 947  AliHBTFunction1D(nbins,maxXval,minXval)
 948 {
 949   //constructor
 950 }
 951 /******************************************************************/
 952
 953 AliHBTTwoPairFctn1D::AliHBTTwoPairFctn1D(const Char_t *name, const Char_t *title):
 954  AliHBTFunction1D(name,title)
 955 {
 956   //constructor
 957 }
 958 /******************************************************************/
 959
 960 AliHBTTwoPairFctn1D::AliHBTTwoPairFctn1D(const Char_t *name, const Char_t *title,
 961                                          Int_t nbins, Float_t maxXval, Float_t minXval):
 962  AliHBTFunction1D(name,title,nbins,maxXval,minXval)
 963 {
 964   //constructor
 965 }
 966 /******************************************************************/
 967
 968 void AliHBTTwoPairFctn1D::ProcessSameEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair)
 969 {
 970   // Fills the numerator using pairs from the same event
 971   partpair  = CheckPair(partpair);
 972   if( partpair )
 973    {
 974      Double_t x = GetValue(trackpair,partpair);
 975      fNumerator->Fill(x);
 976    }
 977 }
 978 /******************************************************************/
 979
 980 void AliHBTTwoPairFctn1D::ProcessDiffEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair)
 981 {
 982   // Fills the denumerator usin mixed pairs
 983   partpair  = CheckPair(partpair);
 984   if( partpair )
 985    {
 986      Double_t x = GetValue(trackpair,partpair);
 987      fDenominator->Fill(x);
 988    }
 989 }
 990 /******************************************************************/
 991 /******************************************************************/
 992 /******************************************************************/
 993
 994 //____________________
 995 ///////////////////////////////////////////////////////
 996 //                                                   //
 997 // AliHBTTwoPairFctn2D                               //
 998 //                                                   //
 999 // Base Calss for 2-dimensinal Functions that need   //
1000 // two pair (simulated and reconstructed)            //
1001 // to fill function                                  //
1002 //                                                   //
1003 // Piotr.Skowronski@cern.ch                          //
1004 // http://alisoft.cern.ch/people/skowron/analyzer    //
1005 //                                                   //
1006 ///////////////////////////////////////////////////////
1007
1008 ClassImp(AliHBTTwoPairFctn2D)
1009
1010 /******************************************************************/
1011
1012 AliHBTTwoPairFctn2D::AliHBTTwoPairFctn2D(const Char_t *name, const Char_t *title):
1013  AliHBTFunction2D(name,title)
1014 {
1015   //constructor
1016 }
1017 /******************************************************************/
1018
1019 AliHBTTwoPairFctn2D::AliHBTTwoPairFctn2D(Int_t nXbins, Double_t maxXval, Double_t minXval,
1020                       Int_t nYbins, Double_t maxYval, Double_t minYval):
1021  AliHBTFunction2D(nXbins, maxXval, minXval, nYbins, maxYval, minYval)
1022 {
1023   //constructor
1024 }
1025 /******************************************************************/
1026
1027 AliHBTTwoPairFctn2D::AliHBTTwoPairFctn2D(const Char_t *name, const Char_t *title,
1028                       Int_t nXbins, Double_t maxXval, Double_t minXval,
1029                       Int_t nYbins, Double_t maxYval, Double_t minYval):
1030  AliHBTFunction2D(name,title,nXbins, maxXval, minXval, nYbins, maxYval, minYval)
1031 {
1032   //constructor
1033 }
1034 /******************************************************************/
1035
1036 void AliHBTTwoPairFctn2D::ProcessSameEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair)
1037 {
1038 //processes pair of particles coming from a same events (real pair)
1039   partpair  = CheckPair(partpair);  //check cuts
1040   if( partpair )
1041    {
1042      Double_t x,y;
1043      GetValues(trackpair,partpair,x,y);
1044      fNumerator->Fill(x,y);
1045    }
1046 }
1047 /******************************************************************/
1048
1049 void AliHBTTwoPairFctn2D::ProcessDiffEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair)
1050 {
1051 //processes pair of particles coming from a different events (mixed pair)
1052   partpair  = CheckPair(partpair);
1053   if( partpair )
1054    {
1055      Double_t x,y;
1056      GetValues(trackpair,partpair,x,y);
1057      fDenominator->Fill(x,y);
1058    }
1059 }
1060
1061 /******************************************************************/
1062 /******************************************************************/
1063 /******************************************************************/
1064
1065 //____________________
1066 ///////////////////////////////////////////////////////
1067 //                                                   //
1068 // AliHBTTwoPairFctn3D                               //
1069 //                                                   //
1070 // Base Calss for 3-dimensinal Functions that need   //
1071 // two pair (simulated and reconstructed)            //
1072 // to fill function                                  //
1073 //                                                   //
1074 // Piotr.Skowronski@cern.ch                          //
1075 // http://alisoft.cern.ch/people/skowron/analyzer    //
1076 //                                                   //
1077 ///////////////////////////////////////////////////////
1078
1079 ClassImp(AliHBTTwoPairFctn3D)
1080
1081 /******************************************************************/
1082
1083 AliHBTTwoPairFctn3D::AliHBTTwoPairFctn3D(const Char_t *name, const Char_t *title):
1084  AliHBTFunction3D(name,title)
1085 {
1086   //constructor
1087 }
1088 /******************************************************************/
1089
1090 AliHBTTwoPairFctn3D::AliHBTTwoPairFctn3D(Int_t nXbins, Double_t maxXval, Double_t minXval,
1091                       Int_t nYbins, Double_t maxYval, Double_t minYval,
1092                       Int_t nZbins, Double_t maxZval, Double_t minZval):
1093  AliHBTFunction3D(nXbins, maxXval, minXval, nYbins, maxYval, minYval, nZbins, maxZval, minZval)
1094 {
1095   //constructor
1096 }
1097 /******************************************************************/
1098
1099 AliHBTTwoPairFctn3D::AliHBTTwoPairFctn3D(const Char_t *name, const Char_t *title,
1100                       Int_t nXbins, Double_t maxXval, Double_t minXval,
1101                       Int_t nYbins, Double_t maxYval, Double_t minYval,
1102                       Int_t nZbins, Double_t maxZval, Double_t minZval):
1103  AliHBTFunction3D(name,title,nXbins, maxXval, minXval, nYbins, maxYval, minYval, nZbins, maxZval, minZval)
1104 {
1105   //constructor
1106 }
1107 /******************************************************************/
1108
1109 void AliHBTTwoPairFctn3D::ProcessSameEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair)
1110 {
1111   // Fills th numerator using pairs from the same event
1112   partpair  = CheckPair(partpair);
1113   if( partpair )
1114    {
1115      Double_t x,y,z;
1116      GetValues(trackpair,partpair,x,y,z);
1117      fNumerator->Fill(x,y,z);
1118    }
1119 }
1120 /******************************************************************/
1121
1122 void AliHBTTwoPairFctn3D::ProcessDiffEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair)
1123 {
1124   // Fills the denumerator using mixed pairs
1125   partpair  = CheckPair(partpair);
1126   if( partpair )
1127    {
1128      Double_t x,y,z;
1129      GetValues(trackpair,partpair,x,y,z);
1130      fDenominator->Fill(x,y,z);
1131    }
1132
1133 }
1134
1135 /******************************************************************/
1136 /******************************************************************/
1137 /******************************************************************/
1138