Correcting some trivial warnings on Alpha
[u/mrichter/AliRoot.git] / 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