CRAB added
[u/mrichter/AliRoot.git] / HBTAN / AliHBTWeightsPID.cxx
1 /* $Id $ */
2
3 //-----------------------------------------------------------
4 //This class introduces the weights calculated according 
5 //with functions of efficiency of identification (TPC+TOF) 
6 //(calculated by B.V. Batyunia).
7 //Author: Ludmila Malinina, JINR (malinina@sunhe.jinr.ru)
8 //-----------------------------------------------------------
9
10 #include "AliHBTWeightsPID.h"
11 #include "AliHBTPair.h"
12 #include "AliHBTParticle.h"
13 #include <TRandom.h>
14 #include <TMath.h>
15 #include <TH1.h>
16 #include <TF1.h>
17
18 ClassImp(AliHBTWeightsPID)  
19  
20 AliHBTWeightsPID* AliHBTWeightsPID::fgWeightsPID=NULL; 
21
22 AliHBTWeightsPID::AliHBTWeightsPID()
23 {
24   //ctor
25   //initial parameters of model
26
27   fPtK   = new TH1F("fPtK"," pt of K+- ",40,0,4);
28   fPtKefftpc   = new TH1F("fPtEfficKpos"," pt of K+ after efficiency cut ",40,0,4);
29   fPtKefftpcboth   = new TH1F("fPtEfficKneg"," pt of K- after efficiency cut ",40,0,4);
30
31   // efficiency functions (from pt) for K+- in the TPC
32   fEffic1pol = new TF1("fEffic1pol4","pol4",100.,480.); 
33   // for -0.9 < eta < -0.7
34   // for  0.7 < eta <  0.9
35   fEffic1pol->SetParameter(0,-0.266362);
36   fEffic1pol->SetParameter(1,0.00565461);
37   fEffic1pol->SetParameter(2,-4.06686e-05);
38   fEffic1pol->SetParameter(3,1.39387e-07);
39   fEffic1pol->SetParameter(4,-1.59674e-10);
40
41   fEffic2pol = new TF1("fEffic2pol4","pol4",100.,540.); 
42
43   fEffic2pol->SetParameter(0,-0.324881);
44   fEffic2pol->SetParameter(1,0.00565381);
45   fEffic2pol->SetParameter(2,-3.23633e-05);
46   fEffic2pol->SetParameter(3,9.72523e-08);
47   fEffic2pol->SetParameter(4,-1.01013e-10);
48                                                                                
49   fEffic3pol = new TF1("fEffic3pol4","pol4",100.,585.); 
50   // for -0.5 < eta < -0.3
51   // for  0.3 < eta <  0.5
52   fEffic3pol->SetParameter(0,-0.306572);
53   fEffic3pol->SetParameter(1,0.00557472);
54   fEffic3pol->SetParameter(2,-3.33752e-05);
55   fEffic3pol->SetParameter(3,9.83241e-08);
56   fEffic3pol->SetParameter(4,-9.5827e-11);
57   fEffic4pol = new TF1("fEffic4pol4","pol4",100.,600.);
58   // for -0.3 < eta < 0.3
59   fEffic4pol->SetParameter(0,-0.168648);
60   fEffic4pol->SetParameter(1,0.00252021);
61   fEffic4pol->SetParameter(2,-1.09113e-05);
62   fEffic4pol->SetParameter(3,3.34871e-08);
63   fEffic4pol->SetParameter(4,-3.31691e-11);
64
65   // efficiency functions (from pt) for K+- in the TOF
66   fEffic1polTOF = new TF1("fEffic1pol4TOF","pol4",0.2,2.0); 
67   // for -0.9 < eta < -0.7
68   // for  0.7 < eta <  0.9
69   fEffic1polTOF->SetParameter(0,-0.165885);
70   fEffic1polTOF->SetParameter(1,0.717459);
71   fEffic1polTOF->SetParameter(2,-0.457131);
72   fEffic1polTOF->SetParameter(3,0.284753);
73   fEffic1polTOF->SetParameter(4,-0.105215);
74                                                                                
75   fEffic2polTOF = new TF1("fEffic2pol4TOF","pol4",0.2,2.4); 
76   // for -0.7 < eta < -0.5
77   // for  0.5 < eta <  0.7
78
79   fEffic2polTOF->SetParameter(0,-0.165947);
80   fEffic2polTOF->SetParameter(1,0.702475);
81   fEffic2polTOF->SetParameter(2,-0.300313);
82   fEffic2polTOF->SetParameter(3,0.127047);
83   fEffic2polTOF->SetParameter(4,-0.0489395);
84                                                                                
85
86   fEffic3polTOF = new TF1("fEffic3pol4TOF","pol4",0.2,2.4); 
87   // for -0.5 < eta < -0.3
88   // for  0.3 < eta <  0.5
89   fEffic3polTOF->SetParameter(0,-0.339516);
90   fEffic3polTOF->SetParameter(1,1.56942);
91   fEffic3polTOF->SetParameter(2,-1.43132);
92   fEffic3polTOF->SetParameter(3,0.727148);
93   fEffic3polTOF->SetParameter(4,-0.158444);
94   fEffic4polTOF = new TF1("fEffic4pol4TOF","pol4",0.2,2.6);
95   // for -0.3 < eta < 0.3
96   fEffic4polTOF->SetParameter(0,-0.243435);
97   fEffic4polTOF->SetParameter(1,1.00928);
98   fEffic4polTOF->SetParameter(2,-0.594597);
99   fEffic4polTOF->SetParameter(3,0.212601);
100   fEffic4polTOF->SetParameter(4,-0.0453419);
101
102 }                                            
103
104 AliHBTWeightsPID* AliHBTWeightsPID::Instance()
105 {
106   //Creates an instance of the class 
107   //or returns pointer to already existing one
108   if (fgWeightsPID) {
109     return fgWeightsPID;
110   } else {
111     fgWeightsPID = new AliHBTWeightsPID();
112     return fgWeightsPID;
113   }
114 }
115
116
117 Double_t AliHBTWeightsPID::GetWeightPID(const AliHBTPair* trackpair)
118 {
119   //Calculates the weight of "trackpair"
120   AliHBTParticle *track1 = trackpair->Particle1();
121   AliHBTParticle *track2 = trackpair->Particle2();
122     
123   Double_t pt1=track1->Pt();
124   Double_t eta1=track1->Eta();
125
126   fPtK->Fill(pt1); 
127
128   if(TMath::Abs(eta1) > 0.7 && TMath::Abs(eta1) < 0.9 && pt1 < 0.48) {
129     fEfficTPC1 = fEffic1pol->Eval(pt1*1000.);
130   }else if(TMath::Abs(eta1) > 0.5 && TMath::Abs(eta1) < 0.7 && pt1 < 0.54) {
131     fEfficTPC1 = fEffic2pol->Eval(pt1*1000.);
132   }else if(TMath::Abs(eta1) > 0.3 && TMath::Abs(eta1) < 0.5 && pt1 < 0.585) {
133     fEfficTPC1 = fEffic3pol->Eval(pt1*1000.);
134   }else if(eta1 > -0.3 && eta1 < 0.3 && pt1 < 0.6) {
135     fEfficTPC1 = fEffic4pol->Eval(pt1*1000.);
136   }
137   
138   // TOF efficiency
139   
140   if(pt1 > 0.2) {
141     if(TMath::Abs(eta1) > 0.7 && TMath::Abs(eta1) < 0.9 && pt1 < 2.0) {
142       fEfficTOF1 = fEffic1polTOF->Eval(pt1);
143     }else if(TMath::Abs(eta1) > 0.5 && TMath::Abs(eta1) < 0.7 && pt1 < 2.4) {
144       fEfficTOF1 = fEffic2polTOF->Eval(pt1);
145     }else if(TMath::Abs(eta1) > 0.3 && TMath::Abs(eta1) < 0.5 && pt1 < 2.4) {
146       fEfficTOF1 = fEffic3polTOF->Eval(pt1);
147     }else if(eta1 > -0.3 && eta1 < 0.3 && pt1 < 2.6) {
148       fEfficTOF1 = fEffic4polTOF->Eval(pt1);
149     }
150   }
151   
152   Double_t rndmtpc=gRandom->Rndm();
153   Double_t rndmtof=gRandom->Rndm();
154   Double_t weightPID1=1.;
155   if(fEfficTPC1 < rndmtpc && fEfficTOF1 < rndmtof) { weightPID1=0.;}
156   fPtKefftpc->Fill(pt1,weightPID1); 
157   
158   
159   Double_t pt2=track2->Pt();
160   Double_t eta2=track2->Eta();
161   
162   // TPC efficiency
163
164   if(TMath::Abs(eta2) > 0.7 && TMath::Abs(eta2) < 0.9 && pt2 < 0.48) {
165     fEfficTPC1 = fEffic1pol->Eval(pt2*1000.);
166   }else if(TMath::Abs(eta2) > 0.5 && TMath::Abs(eta2) < 0.7 && pt2 < 0.54) {
167     fEfficTPC1 = fEffic2pol->Eval(pt2*1000.);
168   }else if(TMath::Abs(eta2) > 0.3 && TMath::Abs(eta2) < 0.5 && pt2 < 0.585) {
169     fEfficTPC1 = fEffic3pol->Eval(pt2*1000.);
170   }else if(eta2 > -0.3 && eta2 < 0.3 && pt2 < 0.6) {
171     fEfficTPC1 = fEffic4pol->Eval(pt2*1000.);
172   }
173   
174   // TOF efficiency
175   
176   if(pt2 > 0.2) {
177     if(TMath::Abs(eta2) > 0.7 && TMath::Abs(eta2) < 0.9 && pt2 < 2.0) {
178       fEfficTOF1 = fEffic1polTOF->Eval(pt2);
179     }else if(TMath::Abs(eta2) > 0.5 && TMath::Abs(eta2) < 0.7 && pt2 < 2.4) {
180       fEfficTOF1 = fEffic2polTOF->Eval(pt2);
181     }else if(TMath::Abs(eta2) > 0.3 && TMath::Abs(eta2) < 0.5 && pt2 < 2.4) {
182       fEfficTOF1 = fEffic3polTOF->Eval(pt2);
183     }else if(eta2 > -0.3 && eta2 < 0.3 && pt2 < 2.6) {
184       fEfficTOF1 = fEffic4polTOF->Eval(pt2);
185     }
186   }
187   
188   rndmtpc=gRandom->Rndm();
189   rndmtof=gRandom->Rndm();
190   
191   Double_t weightPID2=1.;
192   if(fEfficTPC1 < rndmtpc && fEfficTOF1 < rndmtof) { weightPID2=0.;}
193   
194   Double_t weightPID=weightPID1*weightPID2;
195   fPtKefftpcboth->Fill(pt1,weightPID); 
196   
197   return weightPID;
198 }