In Two Pair Functions only particle pair is checked for cut pass. All cout's substitu...
[u/mrichter/AliRoot.git] / HBTAN / AliHBTFunction.cxx
CommitLineData
1b446896 1#include "AliHBTFunction.h"
2/******************************************************************/
3/*
4Piotr Krzysztof Skowronski
5Piotr.Skowronski@cern.ch
6Base classes for HBT functions
7
8 function
9 / \
10 / \
11 / \
12 / \
13 / \
14 / \
15 / \
27b3fe5d 16 one pair two pair
1b446896 17 / | \ / | \
18 / | \ / | \
19 1D 2D 3D 1D 2D 3D
20
21 four particle functions are intendent to be resolution functions:
22 it is mecessary to have simulated particle pair corresponding to given
23 recontructed track pair in order to calculate function simualted value
24 and recontructed value to be further histogrammed
25
26*/
27/******************************************************************/
28/******************************************************************/
29
30ClassImp( AliHBTFunction )
31
32AliHBTFunction::AliHBTFunction()
33{
ba95ae3f 34 fPairCut = new AliHBTEmptyPairCut(); //dummy cut
35}
36/******************************************************************/
37AliHBTFunction::AliHBTFunction(const char* name,const char* title):TNamed(name,title)
38{
39 fPairCut = new AliHBTEmptyPairCut(); //dummy cut
1b446896 40}
2674b9ff 41/******************************************************************/
42
43AliHBTFunction::~AliHBTFunction()
44 {
45 if (fPairCut) delete fPairCut;
46 }
47/******************************************************************/
1b446896 48
49void AliHBTFunction::
50Write()
51 {
52 if (GetNumerator()) GetNumerator()->Write();
53 if (GetDenominator()) GetDenominator()->Write();
54 TH1* res = GetResult();
976183fd 55 if (res) res->Write();
1b446896 56 }
57/******************************************************************/
58
266c51ca 59TH1* AliHBTFunction::GetRatio(Double_t normfactor)
1b446896 60 {
266c51ca 61 if (gDebug>0) Info("GetRatio","Norm. Factor is %f for %s",normfactor,GetName());
976183fd 62
63 if (normfactor == 0.0)
64 {
65 Error("GetRatio","Scaling Factor is 0. Null poiner returned");
66 return 0x0;
67 }
1b446896 68 TString str = fName + " ratio";
69 TH1 *result = (TH1*)GetNumerator()->Clone(str.Data());
70
71 result->SetTitle(str.Data());
1b446896 72
73 result->Divide(GetNumerator(),GetDenominator(),normfactor);
74
75 return result;
76
77 }
78/******************************************************************/
79void AliHBTFunction::SetPairCut(AliHBTPairCut* cut)
80{
81//Sets new Pair Cut. Old one is deleted
82//Note that it is created new object instead of simple pointer set
83//I do not want to have pointer
84//to object created somewhere else
85//because in that case I could not believe that
86//it would always exist (sb could delete it)
87//so we have always own copy
88
89 if(!cut)
90 {
91 Error("AliHBTFunction::SetPairCut","argument is NULL");
92 return;
93 }
94 delete fPairCut;
95 fPairCut = (AliHBTPairCut*)cut->Clone();
96
97}
98
99/******************************************************************/
100
101void AliHBTFunction::
102Rename(const Char_t * name)
103 {
104 //renames the function and histograms
105 SetName(name);
106 SetTitle(name);
107
108 TString numstr = fName + " Numerator"; //title and name of the
109 //numerator histogram
110 TString denstr = fName + " Denominator";//title and name of the
111 //denominator histogram
112
113 GetNumerator()->SetName(numstr.Data());
114 GetNumerator()->SetTitle(numstr.Data());
115
116 GetDenominator()->SetName(denstr.Data());
117 GetDenominator()->SetTitle(denstr.Data());
118
119 }
120
121void AliHBTFunction::
122Rename(const Char_t * name, const Char_t * title)
123 {
124 //renames and retitle the function and histograms
125
126 SetName(name);
127 SetTitle(title);
128
129 TString numstrn = fName + " Numerator"; //name of the
130 //numerator histogram
131
132 TString numstrt = fTitle + " Numerator"; //title of the
133 //numerator histogram
134
135 TString denstrn = fName + " Denominator";//name of the
136 //denominator histogram
137
138 TString denstrt = fTitle + " Denominator";//title of the
139 //denominator histogram
140
141
142 GetNumerator()->SetName(numstrn.Data());
143 GetNumerator()->SetTitle(numstrt.Data());
144
145 GetDenominator()->SetName(denstrn.Data());
146 GetDenominator()->SetTitle(denstrt.Data());
147
148
149 }
150
151/******************************************************************/
152/******************************************************************/
153/******************************************************************/
154
27b3fe5d 155ClassImp( AliHBTOnePairFctn )
1b446896 156
157/******************************************************************/
158/******************************************************************/
159/******************************************************************/
160
27b3fe5d 161ClassImp( AliHBTTwoPairFctn)
1b446896 162
163/******************************************************************/
164/******************************************************************/
165/******************************************************************/
166
27b3fe5d 167ClassImp( AliHBTOnePairFctn1D )
ba95ae3f 168AliHBTOnePairFctn1D::AliHBTOnePairFctn1D():fNBinsToScale(30)
169 {
170 fNumerator = 0x0;
171 fDenominator = 0x0;
172 }
1b446896 173
27b3fe5d 174AliHBTOnePairFctn1D::
175AliHBTOnePairFctn1D(Int_t nbins, Double_t maxXval, Double_t minXval)
1b446896 176 {
177 //Constructor of Two Part One Dimentional Function
178 // nbins: number of bins in histograms - default 100
179 // maxXval and minXval: range of histgram(s) default 0 - 0.15 (GeV)
180
181
182 TString numstr = fName + " Numerator"; //title and name of the
183 //numerator histogram
184 TString denstr = fName + " Denominator";//title and name of the
185 //denominator histogram
186
187 fNumerator = new TH1D(numstr.Data(),numstr.Data(),nbins,minXval,maxXval);
188 fDenominator = new TH1D(denstr.Data(),denstr.Data(),nbins,minXval,maxXval);
189
190 fNumerator->Sumw2();
191 fDenominator->Sumw2();
192
976183fd 193 fNBinsToScale = 30;
194
1b446896 195 }
ba95ae3f 196AliHBTOnePairFctn1D::
197AliHBTOnePairFctn1D(const Char_t *name, const Char_t *title,
198 Int_t nbins, Double_t maxXval, Double_t minXval)
199 :AliHBTOnePairFctn(name,title)
200{
201 TString numstr = fName + " Numerator"; //title and name of the
202 //numerator histogram
203 TString denstr = fName + " Denominator";//title and name of the
204 //denominator histogram
205
206 fNumerator = new TH1D(numstr.Data(),numstr.Data(),nbins,minXval,maxXval);
207 fDenominator = new TH1D(denstr.Data(),denstr.Data(),nbins,minXval,maxXval);
208
209 fNumerator->Sumw2();
210 fDenominator->Sumw2();
211
212 fNBinsToScale = 30;
213}
1b446896 214/******************************************************************/
27b3fe5d 215AliHBTOnePairFctn1D::~AliHBTOnePairFctn1D()
1b446896 216{
217 delete fNumerator;
218 delete fDenominator;
219}
976183fd 220/******************************************************************/
221
27b3fe5d 222void AliHBTOnePairFctn1D::ProcessSameEventParticles(AliHBTPair* pair)
976183fd 223{
224 //Fills the numerator
225 pair = CheckPair(pair);
226 if(pair) fNumerator->Fill(GetValue(pair));
227}
228/******************************************************************/
27b3fe5d 229void AliHBTOnePairFctn1D::ProcessDiffEventParticles(AliHBTPair* pair)
976183fd 230 {
231 //fills denumerator
232 pair = CheckPair(pair);
233 if(pair) fDenominator->Fill(GetValue(pair));
234
235 }
236/******************************************************************/
27b3fe5d 237Double_t AliHBTOnePairFctn1D::Scale()
976183fd 238{
266c51ca 239 if (gDebug>0) Info("Scale","Enetered Scale()");
976183fd 240 if(!fNumerator)
241 {
242 Error("Scale","No numerator");
243 return 0.0;
244 }
245 if(!fDenominator)
246 {
247 Error("Scale","No denominator");
248 return 0.0;
249 }
250
251 if(fNBinsToScale < 1)
252 {
253 return 0.0;
254 Error("Scale","Number of bins for scaling is smaller thnan 1");
255 }
256 Int_t nbins = fNumerator->GetNbinsX();
257 if (fNBinsToScale > nbins)
258 {
259 Error("Scale","Number of bins for scaling is bigger thnan number of bins in histograms");
260 return 0.0;
261 }
266c51ca 262 if (gDebug>0) Info("Scale","No errors detected");
976183fd 263
2674b9ff 264 Double_t ratio;
265 Double_t sum = 0;
266 Int_t N = 0;
267
976183fd 268 Int_t offset = nbins - fNBinsToScale - 1;
269 Int_t i;
270 for ( i = offset; i< nbins; i++)
271 {
2674b9ff 272 if ( fNumerator->GetBinContent(i) > 0.0 )
976183fd 273 {
2674b9ff 274 ratio = fDenominator->GetBinContent(i)/fNumerator->GetBinContent(i);
275 sum += ratio;
276 N++;
976183fd 277 }
278 }
2674b9ff 279
266c51ca 280 if(gDebug > 0) Info("Scale","sum=%f fNBinsToScale=%d N=%d",sum,fNBinsToScale,N);
2674b9ff 281
282 if (N == 0) return 0.0;
283 Double_t ret = sum/((Double_t)N);
976183fd 284
266c51ca 285 if(gDebug > 0) Info("Scale","returning %f",ret);
2674b9ff 286 return ret;
976183fd 287}
288
1b446896 289/******************************************************************/
290/******************************************************************/
291/******************************************************************/
292
27b3fe5d 293ClassImp( AliHBTOnePairFctn2D )
1b446896 294
27b3fe5d 295AliHBTOnePairFctn2D::
296AliHBTOnePairFctn2D(Int_t nXbins, Double_t maxXval, Double_t minXval ,
1b446896 297 Int_t nYbins, Double_t maxYval, Double_t minYval)
298
299{
300 TString numstr = fName + " Numerator"; //title and name of the
301 //numerator histogram
302 TString denstr = fName + " Denominator";//title and name of the
303 //denominator histogram
304
305 fNumerator = new TH2D(numstr.Data(),numstr.Data(),
306 nXbins,minXval,maxXval,
307 nYbins,minYval,maxYval);
308
309 fDenominator = new TH2D(denstr.Data(),denstr.Data(),
310 nXbins,minXval,maxXval,
311 nYbins,minYval,maxYval);
312
313 fNumerator->Sumw2();
314 fDenominator->Sumw2();
315
316}
27b3fe5d 317AliHBTOnePairFctn2D::~AliHBTOnePairFctn2D()
1b446896 318{
319 delete fNumerator;
320 delete fDenominator;
321}
27b3fe5d 322void AliHBTOnePairFctn2D::ProcessSameEventParticles(AliHBTPair* pair)
1b446896 323{
324 pair = CheckPair(pair);
325 if(pair)
326 {
327 Double_t x,y;
328 GetValues(pair,x,y);
19b27ae5 329 fNumerator->Fill(x,y);
1b446896 330 }
331}
332
27b3fe5d 333void AliHBTOnePairFctn2D::ProcessDiffEventParticles(AliHBTPair* pair)
1b446896 334{
335 pair = CheckPair(pair);
336 if(pair)
337 {
338 Double_t x,y;
339 GetValues(pair,x,y);
19b27ae5 340 fDenominator->Fill(x,y);
1b446896 341 }
342
343}
344
345
346/******************************************************************/
347/******************************************************************/
348/******************************************************************/
349
27b3fe5d 350ClassImp( AliHBTOnePairFctn3D)
1b446896 351
27b3fe5d 352AliHBTOnePairFctn3D::
353AliHBTOnePairFctn3D(Int_t nXbins, Double_t maxXval, Double_t minXval,
1b446896 354 Int_t nYbins, Double_t maxYval, Double_t minYval,
355 Int_t nZbins, Double_t maxZval, Double_t minZval)
356
357{
358 TString numstr = fName + " Numerator"; //title and name of the
359 //numerator histogram
360 TString denstr = fName + " Denominator";//title and name of the
361 //denominator histogram
362
363 fNumerator = new TH3D(numstr.Data(),numstr.Data(),
364 nXbins,minXval,maxXval,
365 nYbins,minYval,maxYval,
366 nZbins,minZval,maxZval);
367
368 fDenominator = new TH3D(denstr.Data(),denstr.Data(),
369 nXbins,minXval,maxXval,
370 nYbins,minYval,maxYval,
371 nZbins,minZval,maxZval);
372
373 fNumerator->Sumw2();
374 fDenominator->Sumw2();
375
376}
ba95ae3f 377/******************************************************************/
7836ee94 378
5a9432a9 379AliHBTOnePairFctn3D::~AliHBTOnePairFctn3D()
7836ee94 380{
381 delete fNumerator;
382 delete fDenominator;
383}
ba95ae3f 384/******************************************************************/
7836ee94 385
7836ee94 386
ba95ae3f 387/******************************************************************/
388/******************************************************************/
389/******************************************************************/
390ClassImp( AliHBTTwoPairFctn1D)
391
392AliHBTTwoPairFctn1D::
393AliHBTTwoPairFctn1D(Int_t nbins, Double_t maxval, Double_t minval)
394 {
395 TString numstr = fName + " Numerator"; //title and name of the
396 //numerator histogram
397 TString denstr = fName + " Denominator";//title and name of the
398 //denominator histogram
399
400 fNumerator = new TH1D(numstr.Data(),numstr.Data(),
401 nbins,minval,maxval);
402
403 fDenominator = new TH1D(denstr.Data(),denstr.Data(),
404 nbins,minval,maxval);
405
406 fNumerator->Sumw2();
407 fDenominator->Sumw2();
8e0e9b48 408 fNBinsToScale = 30;
ba95ae3f 409 }
410
411AliHBTTwoPairFctn1D::
412AliHBTTwoPairFctn1D(const Char_t* name, const Char_t* title,
413 Int_t nbins, Double_t maxval, Double_t minval)
414 :AliHBTTwoPairFctn(name,title)
415 {
416 TString numstr = fName + " Numerator"; //title and name of the
417 //numerator histogram
418 TString denstr = fName + " Denominator";//title and name of the
419 //denominator histogram
420
421 fNumerator = new TH1D(numstr.Data(),numstr.Data(),
422 nbins,minval,maxval);
423
424 fDenominator = new TH1D(denstr.Data(),denstr.Data(),
425 nbins,minval,maxval);
426
427 fNumerator->Sumw2();
428 fDenominator->Sumw2();
8e0e9b48 429 fNBinsToScale = 30;
ba95ae3f 430 }
431
432
433/******************************************************************/
434AliHBTTwoPairFctn1D::~AliHBTTwoPairFctn1D()
435{
436 delete fNumerator;
437 delete fDenominator;
438}
266c51ca 439void AliHBTTwoPairFctn1D::ProcessSameEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair)
ba95ae3f 440{
441 partpair = CheckPair(partpair);
266c51ca 442 if( partpair )
ba95ae3f 443 {
444 Double_t x = GetValue(trackpair,partpair);
445 fNumerator->Fill(x);
446 }
447}
448/******************************************************************/
449
4e871a8a 450void AliHBTTwoPairFctn1D::ProcessDiffEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair)
ba95ae3f 451{
452 partpair = CheckPair(partpair);
266c51ca 453 if( partpair )
ba95ae3f 454 {
455 Double_t x = GetValue(trackpair,partpair);
456 fDenominator->Fill(x);
457 }
ba95ae3f 458}
4e871a8a 459/******************************************************************/
460Double_t AliHBTTwoPairFctn1D::Scale()
461{
266c51ca 462 if (gDebug>0) Info("Scale","Enetered Scale()");
4e871a8a 463 if(!fNumerator)
464 {
465 Error("Scale","No numerator");
466 return 0.0;
467 }
468 if(!fDenominator)
469 {
470 Error("Scale","No denominator");
471 return 0.0;
472 }
473
474 if(fNBinsToScale < 1)
475 {
476 return 0.0;
477 Error("Scale","Number of bins for scaling is smaller thnan 1");
478 }
479 Int_t nbins = fNumerator->GetNbinsX();
480 if (fNBinsToScale > nbins)
481 {
482 Error("Scale","Number of bins for scaling is bigger thnan number of bins in histograms");
483 return 0.0;
484 }
266c51ca 485 if (gDebug>0) Info("Scale","No errors detected");
4e871a8a 486
487 Double_t ratio;
488 Double_t sum = 0;
489 Int_t N = 0;
490
491 Int_t offset = nbins - fNBinsToScale - 1;
492 Int_t i;
493 for ( i = offset; i< nbins; i++)
494 {
495 if ( fNumerator->GetBinContent(i) > 0.0 )
496 {
497 ratio = fDenominator->GetBinContent(i)/fNumerator->GetBinContent(i);
498 sum += ratio;
499 N++;
500 }
501 }
502
266c51ca 503 if(gDebug > 0) Info("Scale","sum=%f fNBinsToScale=%d N=%d",sum,fNBinsToScale,N);
4e871a8a 504
505 if (N == 0) return 0.0;
506 Double_t ret = sum/((Double_t)N);
507
266c51ca 508 if(gDebug > 0) Info("Scale","returning %f",ret);
4e871a8a 509 return ret;
510}
ba95ae3f 511
1b446896 512/******************************************************************/
513/******************************************************************/
514/******************************************************************/
27b3fe5d 515ClassImp( AliHBTTwoPairFctn2D)
1b446896 516
517
27b3fe5d 518AliHBTTwoPairFctn2D::
519AliHBTTwoPairFctn2D(Int_t nXbins, Double_t maxXval, Double_t minXval ,
1b446896 520 Int_t nYbins, Double_t maxYval, Double_t minYval)
521
522{
523 TString numstr = fName + " Numerator"; //title and name of the
524 //numerator histogram
525 TString denstr = fName + " Denominator";//title and name of the
526 //denominator histogram
527
528 fNumerator = new TH2D(numstr.Data(),numstr.Data(),
529 nXbins,minXval,maxXval,
530 nYbins,minYval,maxYval);
531
532 fDenominator = new TH2D(denstr.Data(),denstr.Data(),
533 nXbins,minXval,maxXval,
534 nYbins,minYval,maxYval);
535
536 fNumerator->Sumw2();
537 fDenominator->Sumw2();
538
539}
27b3fe5d 540AliHBTTwoPairFctn2D::~AliHBTTwoPairFctn2D()
1b446896 541{
542 delete fNumerator;
543 delete fDenominator;
544}
27b3fe5d 545void AliHBTTwoPairFctn2D::
1b446896 546ProcessSameEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair)
547{
548 partpair = CheckPair(partpair);
266c51ca 549 if( partpair )
1b446896 550 {
551 Double_t x,y;
552 GetValues(trackpair,partpair,x,y);
19b27ae5 553 fNumerator->Fill(x,y);
1b446896 554 }
555}
556
27b3fe5d 557void AliHBTTwoPairFctn2D::
1b446896 558ProcessDiffEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair)
559{
560 partpair = CheckPair(partpair);
266c51ca 561 if( partpair )
1b446896 562 {
563 Double_t x,y;
564 GetValues(trackpair,partpair,x,y);
19b27ae5 565 fDenominator->Fill(x,y);
1b446896 566 }
567
568}
569
ba95ae3f 570/******************************************************************/
571/******************************************************************/
572/******************************************************************/
573ClassImp(AliHBTTwoPairFctn3D)
574
575void AliHBTTwoPairFctn3D::
576ProcessSameEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair)
577{
578 partpair = CheckPair(partpair);
266c51ca 579 if( partpair )
ba95ae3f 580 {
581 Double_t x,y,z;
582 GetValues(trackpair,partpair,x,y,z);
19b27ae5 583 fNumerator->Fill(x,y,z);
ba95ae3f 584 }
585}
586
587void AliHBTTwoPairFctn3D::
588ProcessDiffEventParticles(AliHBTPair* trackpair, AliHBTPair* partpair)
589{
590 partpair = CheckPair(partpair);
266c51ca 591 if( partpair )
ba95ae3f 592 {
593 Double_t x,y,z;
594 GetValues(trackpair,partpair,x,y,z);
19b27ae5 595 fDenominator->Fill(x,y,z);
ba95ae3f 596 }
597
598}
599