]>
Commit | Line | Data |
---|---|---|
bc92a3aa | 1 | //Christine Nattrass, University of Tennessee at Knoxville |
2 | //This macro is to calculate the contributions to Et from various particles based on Levy fits to ALICE data at 900 GeV | |
3 | //It uses d^2N/dpTdy from the papers, weighs it by Et, and integrates over all pT. | |
4 | //A=0 for mesons | |
5 | //A=1 for antinucleons since they would anihilate in the calorimeter if they were observed by a calorimeter | |
6 | //A=-1 for nucleons since they would not anihilate and their mass would not be measured | |
7 | //At the end this is used to calculate the neutral energy correction | |
8 | class AliAnalysisLevyPtModified{ | |
9 | public: | |
10 | virtual ~AliAnalysisLevyPtModified(){;}; | |
11 | ||
12 | ||
13 | Double_t Evaluate(Double_t *pt, Double_t *par) | |
14 | { | |
15 | Double_t ldNdy = par[0]; | |
16 | Double_t l2pi = 2*TMath::Pi(); | |
17 | Double_t lTemp = par[1]; | |
18 | Double_t lPower = par[2]; | |
19 | Double_t lMass = par[3]; | |
20 | Double_t A = par[4]; | |
21 | Double_t Et = TMath::Sqrt(pt[0]*pt[0]+lMass*lMass)+A*lMass; | |
22 | ||
23 | Double_t lBigCoef = ((lPower-1)*(lPower-2)) / (l2pi*lPower*lTemp*(lPower*lTemp+lMass*(lPower-2))); | |
24 | Double_t lInPower = 1 + (TMath::Sqrt(pt[0]*pt[0]+lMass*lMass)-lMass) / (lPower*lTemp); | |
25 | ||
26 | return ldNdy *Et* pt[0] * lBigCoef * TMath::Power(lInPower,(-1)*lPower); | |
27 | } | |
28 | ClassDef(AliAnalysisLevyPtModified, 1); | |
29 | }; | |
df8b9b1f | 30 | class AliAnalysisLevyPtModifiedBaryonEnhanced{ |
31 | public: | |
32 | virtual ~AliAnalysisLevyPtModifiedBaryonEnhanced(){;}; | |
33 | ||
34 | ||
35 | Double_t Evaluate(Double_t *pt, Double_t *par) | |
36 | { | |
37 | Double_t ldNdy = par[0]; | |
38 | Double_t l2pi = 2*TMath::Pi(); | |
39 | Double_t lTemp = par[1]; | |
40 | Double_t lPower = par[2]; | |
41 | Double_t lMass = par[3]; | |
42 | Double_t A = par[4]; | |
43 | Double_t lMassEt = par[3];//only used for calculating Et | |
44 | Double_t lBary0 = par[6]; | |
45 | Double_t lBary1 = par[7]; | |
46 | Double_t lBary2 = par[8]; | |
47 | Double_t lBary3 = par[9]; | |
48 | Double_t lBary4 = par[10]; | |
49 | Double_t lBary5 = par[11]; | |
50 | //this is the Et we would calculate if we had identified the particle as having a mass lMassEt | |
51 | Double_t Et = TMath::Sqrt(pt[0]*pt[0]+lMassEt*lMassEt)+A*lMassEt; | |
52 | ||
53 | Double_t lBigCoef = ((lPower-1)*(lPower-2)) / (l2pi*lPower*lTemp*(lPower*lTemp+lMass*(lPower-2))); | |
54 | Double_t lInPower = 1 + (TMath::Sqrt(pt[0]*pt[0]+lMass*lMass)-lMass) / (lPower*lTemp); | |
55 | //this is the baryon enhancement factor | |
56 | Double_t baryonEnhancement = lBary0*pow(pt[0],lBary1)*exp(-pow(pt[0]/lBary2,lBary3))/(lBary4+lBary5*pt[0]); | |
57 | if(baryonEnhancement<1.0) baryonEnhancement = 1.0; | |
58 | //This is the density of particles times the Et weight | |
59 | return ldNdy *Et* baryonEnhancement*pt[0] * lBigCoef * TMath::Power(lInPower,(-1)*lPower); } | |
60 | ClassDef(AliAnalysisLevyPtModifiedBaryonEnhanced, 1); | |
61 | }; | |
62 | ||
63 | class AliAnalysisLevyPtModifiedStrangeness{ | |
64 | public: | |
65 | virtual ~AliAnalysisLevyPtModifiedStrangeness(){;}; | |
66 | ||
67 | ||
68 | Double_t Evaluate(Double_t *pt, Double_t *par) | |
69 | { | |
70 | Double_t ldNdy = par[0]; | |
71 | Double_t l2pi = 2*TMath::Pi(); | |
72 | Double_t lTemp = par[1]; | |
73 | Double_t lPower = par[2]; | |
74 | Double_t lMass = par[3]; | |
75 | Double_t A = par[4]; | |
76 | Double_t lMassEt = par[3];//only used for calculating Et | |
77 | //this is the Et we would calculate if we had identified the particle as having a mass lMassEt | |
78 | Double_t Et = TMath::Sqrt(pt[0]*pt[0]+lMassEt*lMassEt)+A*lMassEt; | |
79 | ||
80 | Double_t lBigCoef = ((lPower-1)*(lPower-2)) / (l2pi*lPower*lTemp*(lPower*lTemp+lMass*(lPower-2))); | |
81 | Double_t lInPower = 1 + (TMath::Sqrt(pt[0]*pt[0]+lMass*lMass)-lMass) / (lPower*lTemp); | |
82 | Double_t strangeness = (0.4333333*pt[0]-0.0666666)/(0.2333333333*pt[0]-0.01666666); | |
83 | if(strangeness<1.0 || pt[0]<0.1) strangeness = 1.0; | |
84 | //This is the density of particles times the Et weight | |
85 | return ldNdy *strangeness* Et* pt[0] * lBigCoef * TMath::Power(lInPower,(-1)*lPower); | |
86 | } | |
87 | ClassDef(AliAnalysisLevyPtModifiedStrangeness, 1); | |
88 | }; | |
89 | class AliAnalysisLevyPtModifiedStrangenessBaryonEnhanced{ | |
90 | public: | |
91 | virtual ~AliAnalysisLevyPtModifiedStrangenessBaryonEnhanced(){;}; | |
92 | ||
93 | ||
94 | Double_t Evaluate(Double_t *pt, Double_t *par) | |
95 | { | |
96 | Double_t ldNdy = par[0]; | |
97 | Double_t l2pi = 2*TMath::Pi(); | |
98 | Double_t lTemp = par[1]; | |
99 | Double_t lPower = par[2]; | |
100 | Double_t lMass = par[3]; | |
101 | Double_t A = par[4]; | |
102 | Double_t lMassEt = par[3];//only used for calculating Et | |
103 | Double_t lBary0 = par[6]; | |
104 | Double_t lBary1 = par[7]; | |
105 | Double_t lBary2 = par[8]; | |
106 | Double_t lBary3 = par[9]; | |
107 | Double_t lBary4 = par[10]; | |
108 | Double_t lBary5 = par[11]; | |
109 | //this is the Et we would calculate if we had identified the particle as having a mass lMassEt | |
110 | Double_t Et = TMath::Sqrt(pt[0]*pt[0]+lMassEt*lMassEt)+A*lMassEt; | |
111 | ||
112 | Double_t lBigCoef = ((lPower-1)*(lPower-2)) / (l2pi*lPower*lTemp*(lPower*lTemp+lMass*(lPower-2))); | |
113 | Double_t lInPower = 1 + (TMath::Sqrt(pt[0]*pt[0]+lMass*lMass)-lMass) / (lPower*lTemp); | |
114 | //this is the baryon enhancement factor | |
115 | Double_t baryonEnhancement = lBary0*pow(pt[0],lBary1)*exp(-pow(pt[0]/lBary2,lBary3))/(lBary4+lBary5*pt[0]); | |
116 | if(baryonEnhancement<1.0) baryonEnhancement = 1.0; | |
117 | Double_t strangeness = (0.4333333*pt[0]-0.0666666)/(0.2333333333*pt[0]-0.01666666); | |
118 | if(strangeness<1.0 || pt[0]<0.1) strangeness = 1.0; | |
119 | //This is the density of particles times the Et weight | |
120 | return ldNdy *Et*strangeness* baryonEnhancement*pt[0] * lBigCoef * TMath::Power(lInPower,(-1)*lPower); } | |
121 | ClassDef(AliAnalysisLevyPtModifiedStrangenessBaryonEnhanced, 1); | |
122 | }; | |
123 | ||
bc92a3aa | 124 | |
836ce14a | 125 | void CorrNeutralLevyFit(bool hadronic = false){ |
126 | ||
127 | float factor = 0.0; | |
128 | float factorerr = 0.0; | |
129 | if(hadronic){ | |
130 | factor = 0.548; | |
131 | factorerr = 0.003; | |
132 | } | |
bc92a3aa | 133 | |
134 | // particle & $\frac{dN}{dy}$ & T (GeV) & n & $\frac{dE_T}{dy}$& a &\ET \\ \hline | |
135 | // $\pi^{+}+\pi^{-}$ & 2.977 $\pm$ 0.15 & 0.126 $\pm$ 0.001 & 7.82 $\pm$ 0.1 & & & \\ | |
136 | AliAnalysisLevyPtModified *function = new AliAnalysisLevyPtModified(); | |
137 | fPion = new TF1("fPion",function, &AliAnalysisLevyPtModified::Evaluate,0,50,5,"AliAnalysisLevyPtModified","Evaluate"); | |
138 | fPion->SetParameter(0,2.977);//dN/dy | |
139 | fPion->SetParameter(1,0.126);//T | |
140 | fPion->SetParameter(2,7.82);//n | |
141 | fPion->SetParameter(3,0.13957);//mass | |
142 | fPion->SetParameter(4,0);//A=0 for pions | |
143 | float integralPion = fPion->Integral(0,50); | |
144 | float integralErrPion = integralPion*0.15/2.977; | |
145 | float myerrorPionT = 0.0; | |
146 | float myerrorPionn = 0.0; | |
147 | float tmpPion; | |
148 | float T = 0.126; | |
149 | float Terr = 0.001; | |
150 | float n = 7.82; | |
151 | float nerr = 0.1; | |
152 | fPion->SetParameter(1,T+Terr);//T | |
153 | tmpPion = fPion->Integral(0,50); | |
154 | myerrorPionT = TMath::Abs(integralPion-tmpPion); | |
155 | fPion->SetParameter(1,T-Terr);//T | |
156 | tmpPion = fPion->Integral(0,50); | |
157 | if(TMath::Abs(integralPion-tmpPion)>myerrorPionT) myerrorPionT = TMath::Abs(integralPion-tmpPion); | |
158 | fPion->SetParameter(1,T);//T | |
159 | fPion->SetParameter(2,n+nerr);//n | |
160 | tmpPion = fPion->Integral(0,50); | |
bc92a3aa | 161 | myerrorPionn = TMath::Abs(integralPion-tmpPion); |
162 | fPion->SetParameter(2,n-nerr);//n | |
163 | tmpPion = fPion->Integral(0,50); | |
bc92a3aa | 164 | if(TMath::Abs(integralPion-tmpPion)>myerrorPionn) myerrorPionn = TMath::Abs(integralPion-tmpPion); |
bc92a3aa | 165 | //This isn't strictly correct because the errors on the parameters should be correlated but it's close |
836ce14a | 166 | //To get the correct error one would have to fit the spectra data to get the covariance matrix... |
bc92a3aa | 167 | integralErrPion = TMath::Sqrt(TMath::Power(integralErrPion,2)+TMath::Power(myerrorPionT,2)+TMath::Power(myerrorPionn,2)); |
168 | cout<<"Pion Et = "<<integralPion<<"$\\pm$"<<integralErrPion<<endl; | |
169 | ||
170 | // particle & $\frac{dN}{dy}$ & T (GeV) & n & $\frac{dE_T}{dy}$& a &\ET \\ \hline | |
171 | // $K^{+}+K^{-}$ & 0.366 $\pm$ 0.03 & 0.160 $\pm$ 0.006 & 6.087 $\pm$ 0.4 & & & \\ | |
172 | fKaon = new TF1("fKaon",function, &AliAnalysisLevyPtModified::Evaluate,0,50,5,"AliAnalysisLevyPtModified","Evaluate"); | |
173 | fKaon->SetParameter(0,0.366);//dN/dy | |
174 | fKaon->SetParameter(1,0.160);//T | |
175 | fKaon->SetParameter(2,6.087);//n | |
176 | fKaon->SetParameter(3,0.493677);//mass | |
177 | fKaon->SetParameter(4,0);//A=0 for kaons | |
178 | float integralKaon = fKaon->Integral(0,50); | |
179 | float integralErrKaon = integralKaon*0.03/0.366; | |
df8b9b1f | 180 | fKaonStrange = new TF1("fKaonStrange",function, &AliAnalysisLevyPtModifiedStrangeness::Evaluate,0,50,5,"AliAnalysisLevyPtModifiedStrangeness","Evaluate"); |
181 | for(int i=0; i<5;i++){fKaonStrange->SetParameter(i,fKaon->GetParameter(i));} | |
182 | float integralKaonStrange = fKaonStrange->Integral(0,50); | |
bc92a3aa | 183 | float myerrorKaonT = 0.0; |
184 | float myerrorKaonn = 0.0; | |
185 | float tmpKaon; | |
186 | float T = 0.160; | |
187 | float Terr = 0.006; | |
188 | float n = 6.087; | |
189 | float nerr = 0.4; | |
190 | fKaon->SetParameter(1,T+Terr);//T | |
191 | tmpKaon = fKaon->Integral(0,50); | |
192 | myerrorKaonT = TMath::Abs(integralKaon-tmpKaon); | |
193 | fKaon->SetParameter(1,T-Terr);//T | |
194 | tmpKaon = fKaon->Integral(0,50); | |
195 | if(TMath::Abs(integralKaon-tmpKaon)>myerrorKaonT) myerrorKaonT = TMath::Abs(integralKaon-tmpKaon); | |
196 | fKaon->SetParameter(1,T);//T | |
197 | fKaon->SetParameter(2,n+nerr);//n | |
198 | tmpKaon = fKaon->Integral(0,50); | |
bc92a3aa | 199 | myerrorKaonn = TMath::Abs(integralKaon-tmpKaon); |
200 | fKaon->SetParameter(2,n-nerr);//n | |
201 | tmpKaon = fKaon->Integral(0,50); | |
bc92a3aa | 202 | if(TMath::Abs(integralKaon-tmpKaon)>myerrorKaonn) myerrorKaonn = TMath::Abs(integralKaon-tmpKaon); |
bc92a3aa | 203 | //This isn't strictly correct because the errors on the parameters should be correlated but it's close |
204 | integralErrKaon = TMath::Sqrt(TMath::Power(integralErrKaon,2)+TMath::Power(myerrorKaonT,2)+TMath::Power(myerrorKaonn,2)); | |
205 | cout<<"Kaon Et = "<<integralKaon<<"$\\pm$"<<integralErrKaon<<endl; | |
df8b9b1f | 206 | cout<<"Kaon Et Strange = "<<integralKaonStrange<<endl; |
bc92a3aa | 207 | |
208 | // particle & $\frac{dN}{dy}$ & T (GeV) & n & $\frac{dE_T}{dy}$& a &\ET \\ \hline | |
209 | // $p + \bar{p}$ & 0.157 $\pm$ 0.012 & 0.196 $\pm$ 0.009 & 8.6 $\pm$ 1.1 & & & \\ | |
210 | fProton = new TF1("fProton",function, &AliAnalysisLevyPtModified::Evaluate,0,50,5,"AliAnalysisLevyPtModified","Evaluate"); | |
211 | fProton->SetParameter(0,0.157/2.0);//dN/dy | |
212 | fProton->SetParameter(1,0.196);//T | |
213 | fProton->SetParameter(2,8.6);//n | |
214 | fProton->SetParameter(3,0.938272);//mass | |
215 | fProton->SetParameter(4,-1);//A=0 for protons | |
df8b9b1f | 216 | fProtonEnhanced = new TF1("fProtonEnhanced",function, &AliAnalysisLevyPtModifiedBaryonEnhanced::Evaluate,0,50,12,"AliAnalysisLevyPtModifiedBaryonEnhanced","Evaluate"); |
217 | for(int i=0;i<6;i++){fProtonEnhanced->SetParameter(i,fProton->GetParameter(i));}//set all of the spectra parameters to their normal values | |
218 | //and now set the baryon enhancement parameters | |
219 | fProtonEnhanced->SetParameter(6,0.900878*1.2); | |
220 | fProtonEnhanced->SetParameter(7,1.38882); | |
221 | fProtonEnhanced->SetParameter(8,2.6361); | |
222 | fProtonEnhanced->SetParameter(9,1.37751); | |
223 | fProtonEnhanced->SetParameter(10,0.5); | |
224 | fProtonEnhanced->SetParameter(11,-0.03); | |
225 | float integralProtonEnhanced = fProtonEnhanced->Integral(0,50); | |
bc92a3aa | 226 | float integralProton = fProton->Integral(0,50); |
227 | float integralErrProton = integralProton*0.012/0.157; | |
228 | float myerrorProtonT = 0.0; | |
229 | float myerrorProtonn = 0.0; | |
230 | float tmpProton; | |
231 | float T = 0.196; | |
232 | float Terr = 0.009; | |
233 | float n = 8.6; | |
234 | float nerr = 1.1; | |
235 | fProton->SetParameter(1,T+Terr);//T | |
236 | tmpProton = fProton->Integral(0,50); | |
237 | myerrorProtonT = TMath::Abs(integralProton-tmpProton); | |
238 | fProton->SetParameter(1,T-Terr);//T | |
239 | tmpProton = fProton->Integral(0,50); | |
240 | if(TMath::Abs(integralProton-tmpProton)>myerrorProtonT) myerrorProtonT = TMath::Abs(integralProton-tmpProton); | |
241 | fProton->SetParameter(1,T);//T | |
242 | fProton->SetParameter(2,n+nerr);//n | |
243 | tmpProton = fProton->Integral(0,50); | |
bc92a3aa | 244 | myerrorProtonn = TMath::Abs(integralProton-tmpProton); |
245 | fProton->SetParameter(2,n-nerr);//n | |
246 | tmpProton = fProton->Integral(0,50); | |
bc92a3aa | 247 | if(TMath::Abs(integralProton-tmpProton)>myerrorProtonn) myerrorProtonn = TMath::Abs(integralProton-tmpProton); |
bc92a3aa | 248 | //This isn't strictly correct because the errors on the parameters should be correlated but it's close |
249 | integralErrProton = TMath::Sqrt(TMath::Power(integralErrProton,2)+TMath::Power(myerrorProtonT,2)+TMath::Power(myerrorProtonn,2)); | |
250 | cout<<"Proton Et = "<<integralProton<<"$\\pm$"<<integralErrProton<<endl; | |
251 | ||
252 | ||
253 | ||
df8b9b1f | 254 | |
255 | ||
bc92a3aa | 256 | //Antiprotons... |
df8b9b1f | 257 | fProton->SetParameter(0,0.157/2.0);//dN/dy |
258 | fProton->SetParameter(1,0.196);//T | |
259 | fProton->SetParameter(2,8.6);//n | |
260 | fProton->SetParameter(3,0.938272);//mass | |
bc92a3aa | 261 | fProton->SetParameter(2,n);//n |
262 | fProton->SetParameter(4,1);//A=0 for protons | |
263 | float integralAntiProton = fProton->Integral(0,50); | |
264 | float integralErrAntiProton = integralAntiProton*0.012/0.157; | |
df8b9b1f | 265 | fAntiProtonEnhanced = new TF1("fAntiProtonEnhanced",function, &AliAnalysisLevyPtModifiedBaryonEnhanced::Evaluate,0,50,12,"AliAnalysisLevyPtModifiedBaryonEnhanced","Evaluate"); |
266 | for(int i=0;i<6;i++){fAntiProtonEnhanced->SetParameter(i,fProton->GetParameter(i));}//set all of the spectra parameters to their normal values | |
267 | fAntiProtonEnhanced->SetParameter(2,n);//n | |
268 | fAntiProtonEnhanced->SetParameter(4,1);//A=0 for protons | |
269 | //and now set the baryon enhancement parameters | |
270 | fAntiProtonEnhanced->SetParameter(6,0.900878*1.2); | |
271 | fAntiProtonEnhanced->SetParameter(7,1.38882); | |
272 | fAntiProtonEnhanced->SetParameter(8,2.6361); | |
273 | fAntiProtonEnhanced->SetParameter(9,1.37751); | |
274 | fAntiProtonEnhanced->SetParameter(10,0.5); | |
275 | fAntiProtonEnhanced->SetParameter(11,-0.03); | |
276 | float integralAntiProtonEnhanced = fAntiProtonEnhanced->Integral(0,50); | |
bc92a3aa | 277 | float myerrorAntiProtonT = 0.0; |
278 | float myerrorAntiProtonn = 0.0; | |
279 | float tmpAntiProton; | |
280 | float T = 0.196; | |
281 | float Terr = 0.009; | |
282 | float n = 8.6; | |
283 | float nerr = 1.1; | |
284 | fProton->SetParameter(1,T+Terr);//T | |
285 | tmpAntiProton = fProton->Integral(0,50); | |
286 | myerrorAntiProtonT = TMath::Abs(integralAntiProton-tmpAntiProton); | |
287 | fProton->SetParameter(1,T-Terr);//T | |
288 | tmpAntiProton = fProton->Integral(0,50); | |
289 | if(TMath::Abs(integralAntiProton-tmpAntiProton)>myerrorAntiProtonT) myerrorAntiProtonT = TMath::Abs(integralAntiProton-tmpAntiProton); | |
290 | fProton->SetParameter(1,T);//T | |
291 | fProton->SetParameter(2,n+nerr);//n | |
292 | tmpAntiProton = fProton->Integral(0,50); | |
bc92a3aa | 293 | myerrorAntiProtonn = TMath::Abs(integralAntiProton-tmpAntiProton); |
294 | fProton->SetParameter(2,n-nerr);//n | |
295 | tmpAntiProton = fProton->Integral(0,50); | |
bc92a3aa | 296 | if(TMath::Abs(integralAntiProton-tmpAntiProton)>myerrorAntiProtonn) myerrorAntiProtonn = TMath::Abs(integralAntiProton-tmpAntiProton); |
bc92a3aa | 297 | //This isn't strictly correct because the errors on the parameters should be correlated but it's close |
298 | integralErrAntiProton = TMath::Sqrt(TMath::Power(integralErrAntiProton,2)+TMath::Power(myerrorAntiProtonT,2)+TMath::Power(myerrorAntiProtonn,2)); | |
299 | cout<<"AntiProton Et = "<<integralAntiProton<<"$\\pm$"<<integralErrAntiProton<<endl; | |
df8b9b1f | 300 | cout<<"AntiProton enhanced "<<integralAntiProtonEnhanced<<", "<<integralAntiProton/integralAntiProtonEnhanced*100.0<<endl; |
301 | cout<<"Proton enhanced "<<integralProtonEnhanced<<", "<<integralProton/integralProtonEnhanced*100.0<<endl; | |
302 | cout<<"Proton and antiproton enhanced "<< (integralProtonEnhanced+integralAntiProtonEnhanced) <<", " | |
303 | << (integralProtonEnhanced+integralAntiProtonEnhanced) / (integralProton+integralAntiProton)*100.0 | |
304 | <<endl; | |
bc92a3aa | 305 | |
306 | ||
307 | // particle & $\frac{dN}{dy}$ & T (GeV) & n & $\frac{dE_T}{dy}$& a &\ET \\ \hline | |
308 | // \kzeroshort &0.184 $\pm$ 0.006 & 0.168 $\pm$ 0.005 & 6.6 $\pm$ 0.3 & & & \\ | |
309 | fK0S = new TF1("fK0S",function, &AliAnalysisLevyPtModified::Evaluate,0,50,5,"AliAnalysisLevyPtModified","Evaluate"); | |
310 | fK0S->SetParameter(0,0.184);//dN/dy | |
311 | fK0S->SetParameter(1,0.168);//T | |
312 | fK0S->SetParameter(2,6.6);//n | |
313 | fK0S->SetParameter(3,.497614);//mass | |
314 | fK0S->SetParameter(4,0);//A=0 for kaons | |
315 | float integralK0S = fK0S->Integral(0,50); | |
316 | float integralErrK0S = integralK0S*0.006/0.184; | |
df8b9b1f | 317 | fK0Strange = new TF1("fK0Strange",function, &AliAnalysisLevyPtModifiedStrangeness::Evaluate,0,50,5,"AliAnalysisLevyPtModifiedStrangeness","Evaluate"); |
318 | for(int i=0; i<5;i++){fK0Strange->SetParameter(i,fK0S->GetParameter(i));} | |
319 | float integralK0SStrange = fK0Strange->Integral(0,50); | |
bc92a3aa | 320 | float myerrorK0ST = 0.0; |
321 | float myerrorK0Sn = 0.0; | |
322 | float tmpK0S; | |
323 | float T = 0.184; | |
324 | float Terr = 0.006; | |
325 | float n = 6.6; | |
326 | float nerr = 0.3; | |
327 | fK0S->SetParameter(1,T+Terr);//T | |
328 | tmpK0S = fK0S->Integral(0,50); | |
329 | myerrorK0ST = TMath::Abs(integralK0S-tmpK0S); | |
330 | fK0S->SetParameter(1,T-Terr);//T | |
331 | tmpK0S = fK0S->Integral(0,50); | |
332 | if(TMath::Abs(integralK0S-tmpK0S)>myerrorK0ST) myerrorK0ST = TMath::Abs(integralK0S-tmpK0S); | |
333 | fK0S->SetParameter(1,T);//T | |
334 | fK0S->SetParameter(2,n+nerr);//n | |
335 | tmpK0S = fK0S->Integral(0,50); | |
bc92a3aa | 336 | myerrorK0Sn = TMath::Abs(integralK0S-tmpK0S); |
337 | fK0S->SetParameter(2,n-nerr);//n | |
338 | tmpK0S = fK0S->Integral(0,50); | |
bc92a3aa | 339 | if(TMath::Abs(integralK0S-tmpK0S)>myerrorK0Sn) myerrorK0Sn = TMath::Abs(integralK0S-tmpK0S); |
bc92a3aa | 340 | //This isn't strictly correct because the errors on the parameters should be correlated but it's close |
341 | integralErrK0S = TMath::Sqrt(TMath::Power(integralErrK0S,2)+TMath::Power(myerrorK0ST,2)+TMath::Power(myerrorK0Sn,2)); | |
342 | cout<<"K0S Et = "<<integralK0S<<"$\\pm$"<<integralErrK0S<<endl; | |
df8b9b1f | 343 | cout<<"K0S Et Strange = "<<integralK0SStrange<<endl; |
bc92a3aa | 344 | |
345 | // particle & $\frac{dN}{dy}$ & T (GeV) & n & $\frac{dE_T}{dy}$& a &\ET \\ \hline | |
346 | // \lam &0.048 $\pm$ 0.004 & 0.229 $\pm$ 0.015 & 10.8 $\pm$ 2.0 & & & \\ | |
347 | fLambda = new TF1("fLambda",function, &AliAnalysisLevyPtModified::Evaluate,0,50,5,"AliAnalysisLevyPtModified","Evaluate"); | |
348 | fLambda->SetParameter(0,0.048);//dN/dy | |
349 | fLambda->SetParameter(1,0.229);//T | |
350 | fLambda->SetParameter(2,10.8);//n | |
351 | fLambda->SetParameter(3,1.115683);//mass | |
352 | fLambda->SetParameter(4,0);//A=0 for kaons | |
353 | float integralLambda = fLambda->Integral(0,50); | |
354 | float integralErrLambda = integralLambda*0.004/0.048; | |
df8b9b1f | 355 | fLambdaEnhanced = new TF1("fLambdaEnhanced",function, &AliAnalysisLevyPtModifiedBaryonEnhanced::Evaluate,0,50,12,"AliAnalysisLevyPtModifiedBaryonEnhanced","Evaluate"); |
356 | for(int i=0;i<6;i++){fLambdaEnhanced->SetParameter(i,fLambda->GetParameter(i));}//set all of the spectra parameters to their normal values | |
357 | //and now set the baryon enhancement parameters | |
358 | fLambdaEnhanced->SetParameter(6,0.900878); | |
359 | fLambdaEnhanced->SetParameter(7,1.38882); | |
360 | fLambdaEnhanced->SetParameter(8,2.6361); | |
361 | fLambdaEnhanced->SetParameter(9,1.37751); | |
362 | fLambdaEnhanced->SetParameter(10,0.5); | |
363 | fLambdaEnhanced->SetParameter(11,-0.03); | |
364 | float integralLambdaEnhanced = fLambdaEnhanced->Integral(0,50); | |
365 | // fLambdaEnhanced->Draw(); | |
366 | // return; | |
367 | fLambdaStrange = new TF1("fLambdaStrange",function, &AliAnalysisLevyPtModifiedStrangeness::Evaluate,0,50,5,"AliAnalysisLevyPtModifiedStrangeness","Evaluate"); | |
368 | for(int i=0; i<5;i++){fLambdaStrange->SetParameter(i,fLambda->GetParameter(i));} | |
369 | //for(int i=0;i<12;i++){cout<<"Lambda a"<<i<<" "<<fLambdaStrange->GetParameter(i)<<" b"<<i<<" "<<fLambda->GetParameter(i)<<endl;} | |
370 | float integralLambdaStrange = fLambdaStrange->Integral(0,50); | |
371 | //cout<<"test "<<integralLambdaStrange<<endl; | |
372 | fLambdaStrangeEnhanced = new TF1("fLambdaStrangeEnhanced",function, &AliAnalysisLevyPtModifiedStrangenessBaryonEnhanced::Evaluate,0,50,12,"AliAnalysisLevyPtModifiedStrangenessBaryonEnhanced","Evaluate"); | |
373 | for(int i=0; i<12;i++){fLambdaStrangeEnhanced->SetParameter(i,fLambdaEnhanced->GetParameter(i));} | |
374 | // //and now set the baryon enhancement parameters | |
375 | // fLambdaStrangeEnhanced->SetParameter(6,0.900878); | |
376 | // fLambdaStrangeEnhanced->SetParameter(7,1.38882); | |
377 | // fLambdaStrangeEnhanced->SetParameter(8,2.6361); | |
378 | // fLambdaStrangeEnhanced->SetParameter(9,1.37751); | |
379 | // fLambdaStrangeEnhanced->SetParameter(10,0.5); | |
380 | // fLambdaStrangeEnhanced->SetParameter(11,-0.03); | |
381 | //for(int i=0;i<12;i++){cout<<"a"<<i<<" "<<fLambdaStrangeEnhanced->GetParameter(i)<<" b"<<i<<" "<<fLambdaEnhanced->GetParameter(i)<<" c"<<i<<" "<<fLambda->GetParameter(i)<<endl;} | |
382 | float integralLambdaStrangeEnhanced = fLambdaStrangeEnhanced->Integral(0,50); | |
383 | cout<<"Lambda enhanced "<<integralLambdaEnhanced<<", "; | |
384 | if(integralLambdaEnhanced>0.0) cout<<integralLambda/integralLambdaEnhanced*100.0; | |
385 | cout<<endl; | |
bc92a3aa | 386 | float myerrorLambdaT = 0.0; |
387 | float myerrorLambdan = 0.0; | |
388 | float tmpLambda; | |
389 | float T = 0.229; | |
390 | float Terr = 0.015; | |
391 | float n = 10.8; | |
392 | float nerr = 2.0; | |
393 | fLambda->SetParameter(1,T+Terr);//T | |
394 | tmpLambda = fLambda->Integral(0,50); | |
395 | myerrorLambdaT = TMath::Abs(integralLambda-tmpLambda); | |
396 | fLambda->SetParameter(1,T-Terr);//T | |
397 | tmpLambda = fLambda->Integral(0,50); | |
398 | if(TMath::Abs(integralLambda-tmpLambda)>myerrorLambdaT) myerrorLambdaT = TMath::Abs(integralLambda-tmpLambda); | |
399 | fLambda->SetParameter(1,T);//T | |
400 | fLambda->SetParameter(2,n+nerr);//n | |
401 | tmpLambda = fLambda->Integral(0,50); | |
bc92a3aa | 402 | myerrorLambdan = TMath::Abs(integralLambda-tmpLambda); |
403 | fLambda->SetParameter(2,n-nerr);//n | |
404 | tmpLambda = fLambda->Integral(0,50); | |
bc92a3aa | 405 | if(TMath::Abs(integralLambda-tmpLambda)>myerrorLambdan) myerrorLambdan = TMath::Abs(integralLambda-tmpLambda); |
bc92a3aa | 406 | //This isn't strictly correct because the errors on the parameters should be correlated but it's close |
407 | integralErrLambda = TMath::Sqrt(TMath::Power(integralErrLambda,2)+TMath::Power(myerrorLambdaT,2)+TMath::Power(myerrorLambdan,2)); | |
408 | cout<<"Lambda Et = "<<integralLambda<<"$\\pm$"<<integralErrLambda<<endl; | |
df8b9b1f | 409 | cout<<"Lambda Et Strange = "<<integralLambdaStrange<<endl; |
410 | cout<<"Lambda Et Strange Enhanced = "<<integralLambdaStrangeEnhanced<<endl; | |
bc92a3aa | 411 | |
412 | ||
413 | // particle & $\frac{dN}{dy}$ & T (GeV) & n & $\frac{dE_T}{dy}$& a &\ET \\ \hline | |
414 | // \alam & 0.047 $\pm$ 0.005 & 0.210 $\pm$ 0.015 & 9.2 $\pm$ 1.4 & & & \\ & & \\ | |
415 | fAntiLambda = new TF1("fAntiLambda",function, &AliAnalysisLevyPtModified::Evaluate,0,50,5,"AliAnalysisLevyPtModified","Evaluate"); | |
416 | fAntiLambda->SetParameter(0,0.047);//dN/dy | |
417 | fAntiLambda->SetParameter(1,0.210);//T | |
418 | fAntiLambda->SetParameter(2,9.2);//n | |
419 | fAntiLambda->SetParameter(3,1.115683);//mass | |
420 | fAntiLambda->SetParameter(4,0);//A=0 for kaons | |
421 | float integralAntiLambda = fAntiLambda->Integral(0,50); | |
422 | float integralErrAntiLambda = integralAntiLambda*0.005/0.047; | |
df8b9b1f | 423 | fAntiLambdaEnhanced = new TF1("fAntiLambdaEnhanced",function, &AliAnalysisLevyPtModifiedBaryonEnhanced::Evaluate,0,50,12,"AliAnalysisLevyPtModifiedBaryonEnhanced","Evaluate"); |
424 | for(int i=0;i<6;i++){fAntiLambdaEnhanced->SetParameter(i,fAntiLambda->GetParameter(i));}//set all of the spectra parameters to their normal values | |
425 | //and now set the baryon enhancement parameters | |
426 | fAntiLambdaEnhanced->SetParameter(6,0.900878*1.2); | |
427 | fAntiLambdaEnhanced->SetParameter(7,1.38882); | |
428 | fAntiLambdaEnhanced->SetParameter(8,2.6361); | |
429 | fAntiLambdaEnhanced->SetParameter(9,1.37751); | |
430 | fAntiLambdaEnhanced->SetParameter(10,0.5); | |
431 | fAntiLambdaEnhanced->SetParameter(11,-0.03); | |
432 | float integralAntiLambdaEnhanced = fAntiLambdaEnhanced->Integral(0,50); | |
433 | cout<<"AntiLambda enhanced "<<integralAntiLambdaEnhanced<<", "<<integralAntiLambda/integralAntiLambdaEnhanced*100.0<<endl; | |
434 | fAntiLambdaStrange = new TF1("fAntiLambdaStrange",function, &AliAnalysisLevyPtModifiedStrangeness::Evaluate,0,50,5,"AliAnalysisLevyPtModifiedStrangeness","Evaluate"); | |
435 | for(int i=0; i<5;i++){fAntiLambdaStrange->SetParameter(i,fAntiLambda->GetParameter(i));} | |
436 | float integralAntiLambdaStrange = fAntiLambdaStrange->Integral(0,50); | |
437 | ||
438 | fAntiLambdaStrangeEnhanced = new TF1("fAntiLambdaStrangeEnhanced",function, &AliAnalysisLevyPtModifiedStrangenessBaryonEnhanced::Evaluate,0,50,12,"AliAnalysisLevyPtModifiedStrangenessBaryonEnhanced","Evaluate"); | |
439 | for(int i=0; i<12;i++){fAntiLambdaStrangeEnhanced->SetParameter(i,fAntiLambdaEnhanced->GetParameter(i));} | |
440 | float integralAntiLambdaStrangeEnhanced = fAntiLambdaStrangeEnhanced->Integral(0,50); | |
441 | ||
bc92a3aa | 442 | float myerrorAntiLambdaT = 0.0; |
443 | float myerrorAntiLambdan = 0.0; | |
444 | float tmpAntiLambda; | |
445 | float T = 0.210; | |
446 | float Terr = 0.015; | |
447 | float n = 9.2; | |
448 | float nerr = 1.4; | |
449 | fAntiLambda->SetParameter(1,T+Terr);//T | |
450 | tmpAntiLambda = fAntiLambda->Integral(0,50); | |
451 | myerrorAntiLambdaT = TMath::Abs(integralAntiLambda-tmpAntiLambda); | |
452 | fAntiLambda->SetParameter(1,T-Terr);//T | |
453 | tmpAntiLambda = fAntiLambda->Integral(0,50); | |
454 | if(TMath::Abs(integralAntiLambda-tmpAntiLambda)>myerrorAntiLambdaT) myerrorAntiLambdaT = TMath::Abs(integralAntiLambda-tmpAntiLambda); | |
455 | fAntiLambda->SetParameter(1,T);//T | |
456 | fAntiLambda->SetParameter(2,n+nerr);//n | |
457 | tmpAntiLambda = fAntiLambda->Integral(0,50); | |
bc92a3aa | 458 | myerrorAntiLambdan = TMath::Abs(integralAntiLambda-tmpAntiLambda); |
459 | fAntiLambda->SetParameter(2,n-nerr);//n | |
460 | tmpAntiLambda = fAntiLambda->Integral(0,50); | |
bc92a3aa | 461 | if(TMath::Abs(integralAntiLambda-tmpAntiLambda)>myerrorAntiLambdan) myerrorAntiLambdan = TMath::Abs(integralAntiLambda-tmpAntiLambda); |
bc92a3aa | 462 | //This isn't strictly correct because the errors on the parameters should be correlated but it's close |
463 | integralErrAntiLambda = TMath::Sqrt(TMath::Power(integralErrAntiLambda,2)+TMath::Power(myerrorAntiLambdaT,2)+TMath::Power(myerrorAntiLambdan,2)); | |
464 | cout<<"AntiLambda Et = "<<integralAntiLambda<<"$\\pm$"<<integralErrAntiLambda<<endl; | |
df8b9b1f | 465 | cout<<"AntiLambda Et Strange = "<<integralAntiLambdaStrange<<endl; |
466 | cout<<"AntiLambda Et Strange Enhanced = "<<integralAntiLambdaStrangeEnhanced<<endl; | |
bc92a3aa | 467 | |
468 | //now we apply various assumptions | |
469 | float integralK0L = integralK0S; | |
470 | float integralErrK0L = integralErrK0S; | |
471 | float integralNeutron = integralProton; | |
472 | float integralErrNeutron = integralErrProton; | |
473 | ||
836ce14a | 474 | float totalEt = (1.0+factor)*integralPion+integralKaon+2.0*integralProton+2.0*integralAntiProton+2.0*integralK0S+integralLambda+integralAntiLambda; |
bc92a3aa | 475 | float measuredEt = integralPion+integralKaon+1.0*integralProton+1.0*integralAntiProton; |
476 | float fneutral = measuredEt/totalEt; | |
477 | cout<<"fneutral = "<<fneutral<<endl; | |
df8b9b1f | 478 | |
479 | float totalEtEnhanced = (1.0+factor)*integralPion+integralKaon+2.0*integralProtonEnhanced+2.0*integralAntiProtonEnhanced+2.0*integralK0S+integralLambdaEnhanced+integralAntiLambdaEnhanced; | |
480 | float measuredEtEnhanced = integralPion+integralKaon+1.0*integralProtonEnhanced+1.0*integralAntiProtonEnhanced; | |
481 | float fneutralEnhanced = measuredEtEnhanced/totalEtEnhanced; | |
482 | cout<<"fneutralEnhanced = "<<fneutralEnhanced<<endl; | |
483 | ||
484 | ||
485 | float totalEtStrange = (1.0+factor)*integralPion+integralKaonStrange+2.0*integralProton+2.0*integralAntiProton+2.0*integralK0SStrange+integralLambdaStrange+integralAntiLambdaStrange; | |
486 | float measuredEtStrange = integralPion+integralKaonStrange+1.0*integralProton+1.0*integralAntiProton; | |
487 | float fneutralStrange = measuredEtStrange/totalEtStrange; | |
488 | cout<<"fneutralStrange = "<<fneutralStrange<<endl; | |
489 | ||
490 | ||
491 | ||
492 | float totalEtStrangeEnhanced = (1.0+factor)*integralPion+integralKaonStrange+2.0*integralProtonEnhanced+2.0*integralAntiProtonEnhanced+2.0*integralK0SStrange+integralLambdaStrangeEnhanced+integralAntiLambdaStrangeEnhanced; | |
493 | float measuredEtStrangeEnhanced = integralPion+integralKaonStrange+1.0*integralProtonEnhanced+1.0*integralAntiProtonEnhanced; | |
494 | float fneutralStrangeEnhanced = measuredEtStrangeEnhanced/totalEtStrangeEnhanced; | |
495 | cout<<"fneutralStrangeEnhanced = "<<fneutralStrangeEnhanced<<endl; | |
496 | ||
bc92a3aa | 497 | //this is from ugly derivative taking for error propagation |
498 | //form 1: pions, kaons, protons | |
499 | //for f=(xa+A)/(ya+B) | |
500 | //df/da=(xf-yf^2)/(xa+A) | |
501 | //x=y=1; xa+A=measuredEt | |
502 | float errPion = (fneutral-fneutral*fneutral)/measuredEt*integralErrPion; | |
836ce14a | 503 | if(hadronic){ |
504 | //then we have x=1, y = 1+factor | |
505 | //the error on the extra bit of et from pi0s | |
506 | float extraerror = integralPion*factorerr; | |
507 | integralErrPion = TMath::Sqrt(TMath::Power(integralErrPion,2)+TMath::Power(extraerror,2)); | |
508 | errPion = (fneutral-(1.0+factor)*fneutral*fneutral)/measuredEt*integralErrPion; | |
509 | } | |
bc92a3aa | 510 | float errKaon = (fneutral-fneutral*fneutral)/measuredEt*integralErrKaon; |
511 | //x=1,y=2 | |
512 | float errProton = (fneutral-2.0*fneutral*fneutral)/measuredEt*integralErrProton; | |
513 | float errAntiProton = (fneutral-2.0*fneutral*fneutral)/measuredEt*integralErrAntiProton; | |
514 | //form 2: K0S, lambda, antilambda | |
515 | //for f=A/(B+xa) | |
516 | //df/da=-xf^2/A | |
517 | //x=1,A=measuredEt | |
518 | float errLambda = fneutral*fneutral/measuredEt*integralErrLambda; | |
519 | float errAntiLambda = fneutral*fneutral/measuredEt*integralErrAntiLambda; | |
520 | //x=2,A=measuredEt | |
521 | float errK0S = 2.0*fneutral*fneutral/measuredEt*integralErrK0S; | |
522 | cout<<"Errors Pion "<<errPion<<" Kaon "<<errKaon<<" Proton "<<errProton<<" Lambda "<<errLambda<<" AntiLambda "<<errAntiLambda<<" K0S "<<errK0S<<endl; | |
523 | float totalErr = TMath::Sqrt( TMath::Power(errPion,2) + TMath::Power(errKaon,2) + TMath::Power(errProton,2)+ TMath::Power(errAntiProton,2) + TMath::Power(errLambda,2) + TMath::Power(errAntiLambda,2) + TMath::Power(errK0S,2) ); | |
524 | ||
525 | cout<<"fneutral = "<<fneutral<<"$\\pm$"<<totalErr<<endl; | |
526 | cout<<"1/fneutral = "<<1.0/fneutral<<"$\\pm$"<<1.0/fneutral*totalErr/fneutral<<endl; | |
527 | cout<<"Percentage error "<<totalErr/fneutral*100.0<<endl; | |
df8b9b1f | 528 | |
529 | float fneutralPb = (fneutralStrangeEnhanced+fneutral)/2.0; | |
530 | float fneutralPbErr = TMath::Sqrt( TMath::Power((fneutralStrangeEnhanced-fneutral)/2.0,2.0) + TMath::Power(totalErr,2.0)); | |
531 | cout<<"fneutral strange enhanced "<<fneutralPb<<"$\\pm$"<<fneutralPbErr<<endl; | |
532 | cout<<"1/fneutral strange enhanced = "<<1.0/fneutralPb<<"$\\pm$"<<1.0/fneutralPb*fneutralPbErr/fneutralPb<<endl; | |
533 | cout<<"Percentage error "<<fneutralPbErr/fneutralPb*100.0<<endl; | |
534 | ||
bc92a3aa | 535 | } |
536 |