]>
Commit | Line | Data |
---|---|---|
5b3a5a5d | 1 | /************************************************************************** |
2 | * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * | |
3 | * * | |
4 | * Author: The ALICE Off-line Project. * | |
5 | * Contributors are mentioned in the code where appropriate. * | |
6 | * * | |
7 | * Permission to use, copy, modify and distribute this software and its * | |
8 | * documentation strictly for non-commercial purposes is hereby granted * | |
9 | * without fee, provided that the above copyright notice appears in all * | |
10 | * copies and that both the copyright notice and this permission notice * | |
11 | * appear in the supporting documentation. The authors make no claims * | |
12 | * about the suitability of this software for any purpose. It is * | |
13 | * provided "as is" without express or implied warranty. * | |
14 | **************************************************************************/ | |
15 | ||
88cb7938 | 16 | /* $Id$ */ |
5b3a5a5d | 17 | |
18 | // from AliRoot | |
19 | #include "AliFastGlauber.h" | |
20 | // from root | |
21 | #include <TH1F.h> | |
22 | #include <TF1.h> | |
23 | #include <TF2.h> | |
24 | #include <TCanvas.h> | |
25 | #include <TRandom.h> | |
26 | #include <TFile.h> | |
27 | ||
28 | ClassImp(AliFastGlauber) | |
29 | ||
30 | TF1* AliFastGlauber::fWSb = NULL; | |
31 | TF2* AliFastGlauber::fWSbz = NULL; | |
32 | TF1* AliFastGlauber::fWSz = NULL; | |
33 | TF1* AliFastGlauber::fWSta = NULL; | |
34 | TF2* AliFastGlauber::fWStarfi = NULL; | |
35 | TF1* AliFastGlauber::fWStaa = NULL; | |
36 | TF1* AliFastGlauber::fWSgeo = NULL; | |
37 | TF1* AliFastGlauber::fWSbinary = NULL; | |
38 | TF1* AliFastGlauber::fWSN = NULL; | |
39 | Float_t AliFastGlauber::fbMax = 0.; | |
40 | ||
41 | AliFastGlauber::AliFastGlauber() | |
42 | { | |
43 | // Default Constructor | |
44 | // | |
45 | // Defaults for Pb | |
46 | SetWoodSaxonParameters(6.624, 0.549, 0.00, 7.69e-4); | |
47 | SetHardCrossSection(); | |
48 | SetMaxImpact(); | |
49 | } | |
50 | ||
51 | void AliFastGlauber::Init(Int_t mode) | |
52 | { | |
53 | // Initialisation | |
54 | // | |
55 | // Wood-Saxon | |
56 | // | |
57 | fWSb = new TF1("WSb", WSb, 0, fbMax, 4); | |
58 | fWSb->SetParameter(0, fWSr0); | |
59 | fWSb->SetParameter(1, fWSd); | |
60 | fWSb->SetParameter(2, fWSw); | |
61 | fWSb->SetParameter(3, fWSn); | |
62 | ||
63 | fWSbz = new TF2("WSbz", WSbz, 0, fbMax, 4); | |
64 | fWSbz->SetParameter(0, fWSr0); | |
65 | fWSbz->SetParameter(1, fWSd); | |
66 | fWSbz->SetParameter(2, fWSw); | |
67 | fWSbz->SetParameter(3, fWSn); | |
68 | ||
69 | fWSz = new TF1("WSz", WSz, 0, fbMax, 5); | |
70 | fWSz->SetParameter(0, fWSr0); | |
71 | fWSz->SetParameter(1, fWSd); | |
72 | fWSz->SetParameter(2, fWSw); | |
73 | fWSz->SetParameter(3, fWSn); | |
74 | ||
75 | // | |
76 | // Thickness | |
77 | // | |
78 | fWSta = new TF1("WSta", WSta, 0., fbMax, 0); | |
79 | ||
80 | // | |
81 | // Overlap Kernel | |
82 | // | |
83 | fWStarfi = new TF2("WStarfi", WStarfi, 0., fbMax, 0., TMath::Pi(), 1); | |
84 | fWStarfi->SetParameter(0, 0.); | |
85 | fWStarfi->SetNpx(200); | |
86 | fWStarfi->SetNpy(20); | |
87 | // | |
88 | // Overlap | |
89 | // | |
90 | if (! mode) { | |
91 | fWStaa = new TF1("WStaa", WStaa, 0., fbMax, 0); | |
92 | fWStaa->SetNpx(100); | |
93 | } else { | |
94 | TFile* f = new TFile("$(ALICE_ROOT)/FASTSIM/data/glauberPbPb.root"); | |
95 | fWStaa = (TF1*) f->Get("WStaa"); | |
96 | } | |
97 | ||
98 | // | |
99 | // Geometrical Cross-Section | |
100 | // | |
101 | fWSgeo = new TF1("WSgeo", WSgeo, 0., fbMax, 0); | |
102 | fWSgeo->SetNpx(100); | |
103 | // | |
104 | // Hard cross section (~ binary collisions) | |
105 | // | |
106 | fWSbinary = new TF1("WSbinary", WSbinary, 0., fbMax, 1); | |
107 | fWSbinary->SetParameter(0, fSigmaHard); // mb | |
108 | fWSbinary->SetNpx(100); | |
109 | // | |
110 | // Hard collisions per event | |
111 | // | |
112 | fWSN = new TF1("WSN", WSN, 0., fbMax, 1); | |
113 | fWSN->SetNpx(100); | |
114 | } | |
115 | ||
116 | void AliFastGlauber::DrawWSb() | |
117 | { | |
118 | // | |
119 | // Draw Wood-Saxon Nuclear Density Function | |
120 | // | |
121 | TCanvas *c1 = new TCanvas("c1","Wood Saxon",400,10,600,700); | |
122 | c1->cd(); | |
123 | fWSb->Draw(); | |
124 | } | |
125 | ||
126 | void AliFastGlauber::DrawOverlap() | |
127 | { | |
128 | // | |
129 | // Draw Overlap Function | |
130 | // | |
131 | TCanvas *c2 = new TCanvas("c2","Overlap",400,10,600,700); | |
132 | c2->cd(); | |
133 | fWStaa->Draw(); | |
134 | } | |
135 | ||
136 | void AliFastGlauber::DrawThickness() | |
137 | { | |
138 | // | |
139 | // Draw Thickness Function | |
140 | // | |
141 | TCanvas *c3 = new TCanvas("c3","Thickness",400,10,600,700); | |
142 | c3->cd(); | |
143 | fWSta->Draw(); | |
144 | } | |
145 | ||
146 | void AliFastGlauber::DrawGeo() | |
147 | { | |
148 | // | |
149 | // Draw Geometrical Cross-Section | |
150 | // | |
151 | TCanvas *c3 = new TCanvas("c3","Geometrical Cross-Section",400,10,600,700); | |
152 | c3->cd(); | |
153 | fWSgeo->Draw(); | |
154 | } | |
155 | ||
156 | void AliFastGlauber::DrawBinary() | |
157 | { | |
158 | // | |
159 | // Draw Wood-Saxon Nuclear Density Function | |
160 | // | |
161 | TCanvas *c4 = new TCanvas("c4","Binary Cross-Section",400,10,600,700); | |
162 | c4->cd(); | |
163 | fWSbinary->Draw(); | |
164 | } | |
165 | ||
166 | void AliFastGlauber::DrawN() | |
167 | { | |
168 | // | |
169 | // Draw Binaries per event | |
170 | // | |
171 | TCanvas *c5 = new TCanvas("c5","Binaries per event",400,10,600,700); | |
172 | c5->cd(); | |
173 | fWSN->Draw(); | |
174 | } | |
175 | ||
176 | void AliFastGlauber::DrawKernel() | |
177 | { | |
178 | // | |
179 | // Draw Kernel | |
180 | // | |
181 | TCanvas *c6 = new TCanvas("c6","Kernel",400,10,600,700); | |
182 | c6->cd(); | |
183 | fWStarfi->Draw(); | |
184 | } | |
185 | ||
186 | Double_t AliFastGlauber::WSb(Double_t* x, Double_t* par) | |
187 | { | |
188 | // | |
189 | // Woods-Saxon Parameterisation | |
190 | // as a function of radius | |
191 | // | |
192 | Double_t xx = x[0]; | |
193 | Double_t r0 = par[0]; | |
194 | Double_t d = par[1]; | |
195 | Double_t w = par[2]; | |
196 | Double_t n = par[3]; | |
197 | ||
198 | Double_t y = n * (1.+w*(xx/r0)*(xx/r0))/(1.+TMath::Exp((xx-r0)/d)); | |
199 | return y; | |
200 | } | |
201 | ||
202 | Double_t AliFastGlauber::WSbz(Double_t* x, Double_t* par) | |
203 | { | |
204 | // | |
205 | // Wood Saxon Parameterisation | |
206 | // as a function of z and b | |
207 | // | |
208 | Double_t bb = x[0]; | |
209 | Double_t zz = x[1]; | |
210 | Double_t r0 = par[0]; | |
211 | Double_t d = par[1]; | |
212 | Double_t w = par[2]; | |
213 | Double_t n = par[3]; | |
214 | Double_t xx = TMath::Sqrt(bb*bb+zz*zz); | |
215 | Double_t y = n * (1.+w*(xx/r0)*(xx/r0))/(1.+TMath::Exp((xx-r0)/d)); | |
216 | return y; | |
217 | } | |
218 | ||
219 | Double_t AliFastGlauber::WSz(Double_t* x, Double_t* par) | |
220 | { | |
221 | // | |
222 | // Wood Saxon Parameterisation | |
223 | // as a function of z for fixed b | |
224 | // | |
225 | Double_t bb = par[4]; | |
226 | Double_t zz = x[0]; | |
227 | Double_t r0 = par[0]; | |
228 | Double_t d = par[1]; | |
229 | Double_t w = par[2]; | |
230 | Double_t n = par[3]; | |
231 | Double_t xx = TMath::Sqrt(bb*bb+zz*zz); | |
232 | Double_t y = n * (1.+w*(xx/r0)*(xx/r0))/(1.+TMath::Exp((xx-r0)/d)); | |
233 | // if (y < 1.e-6) y = 0; | |
234 | return y; | |
235 | } | |
236 | ||
237 | Double_t AliFastGlauber::WSta(Double_t* x, Double_t* par) | |
238 | { | |
239 | // | |
240 | // Thickness function | |
241 | // | |
242 | Double_t b = x[0]; | |
243 | fWSz->SetParameter(4, b); | |
244 | Double_t y = 2. * fWSz->Integral(0., fbMax); | |
245 | return y; | |
246 | } | |
247 | ||
248 | ||
249 | ||
250 | Double_t AliFastGlauber::WStarfi(Double_t* x, Double_t* par) | |
251 | { | |
252 | // | |
253 | // Kernel for overlap function | |
254 | // | |
255 | Double_t b = par[0]; | |
256 | Double_t r1 = x[0]; | |
257 | Double_t phi = x[1]; | |
258 | Double_t r2 = TMath::Sqrt(r1 * r1 + b * b - 2. * r1 * b * TMath::Cos(phi)); | |
259 | Double_t y = r1 * fWSta->Eval(r1) * fWSta->Eval(r2); | |
260 | return y; | |
261 | } | |
262 | ||
263 | ||
264 | Double_t AliFastGlauber::WStaa(Double_t* x, Double_t* par) | |
265 | { | |
266 | // | |
267 | // Overlap function | |
268 | // | |
269 | Double_t b = x[0]; | |
270 | fWStarfi->SetParameter(0, b); | |
271 | /* | |
272 | Double_t al[2]; | |
273 | Double_t bl[2]; | |
274 | al[0] = 0.; | |
275 | al[1] = 0.; | |
276 | bl[0] = 6.6; | |
277 | bl[1] = TMath::Pi(); | |
278 | Double_t err; | |
279 | ||
280 | Double_t y = 2. * fWStarfi->IntegralMultiple(2, al, bl, 0.001, err); | |
281 | printf("WStaa: %f %f %f\n", b, y, err); | |
282 | */ | |
283 | // | |
284 | // MC Integration | |
285 | // | |
286 | Double_t y = 0; | |
287 | for (Int_t i = 0; i < 100000; i++) | |
288 | { | |
289 | Double_t phi = TMath::Pi() * gRandom->Rndm(); | |
290 | Double_t b1 = fbMax * gRandom->Rndm(); | |
291 | y += fWStarfi->Eval(b1, phi); | |
292 | } | |
293 | y *= 2. * 0.1 * 208. * 208. * TMath::Pi() * fbMax / 100000.; | |
294 | return y; | |
295 | } | |
296 | ||
297 | Double_t AliFastGlauber::WSgeo(Double_t* x, Double_t* par) | |
298 | { | |
299 | // | |
300 | // Geometrical Cross-Section | |
301 | // | |
302 | Double_t b = x[0]; | |
303 | Double_t taa = fWStaa->Eval(b); | |
304 | const Double_t sigma = 55.6; // mbarn | |
305 | ||
306 | Double_t y = 2. * TMath::Pi() * b * (1. - TMath::Exp(- sigma * taa)); // fm | |
307 | return y; | |
308 | } | |
309 | ||
310 | ||
311 | Double_t AliFastGlauber::WSbinary(Double_t* x, Double_t* par) | |
312 | { | |
313 | // | |
c2140715 | 314 | // Number of binary collisions |
5b3a5a5d | 315 | // |
316 | Double_t b = x[0]; | |
317 | Double_t sigma = par[0]; | |
318 | Double_t taa = fWStaa->Eval(b); | |
319 | ||
320 | Double_t y = 2. * TMath::Pi() * b * sigma * taa; // fm | |
321 | return y; | |
322 | } | |
323 | ||
324 | Double_t AliFastGlauber::WSN(Double_t* x, Double_t* par) | |
325 | { | |
326 | // | |
c2140715 | 327 | // Number of hard processes per event |
5b3a5a5d | 328 | // |
329 | Double_t b = x[0]; | |
88cb7938 | 330 | Double_t y = fWSbinary->Eval(b)/fWSgeo->Eval(b); |
5b3a5a5d | 331 | return y; |
332 | } | |
333 | ||
334 | void AliFastGlauber::SimulateTrigger(Int_t n) | |
335 | { | |
336 | // | |
337 | // Simulates Trigger | |
338 | // | |
339 | TH1F* mbtH = new TH1F("mbtH", "MB Trigger b-Distribution", 100, 0., 20.); | |
340 | TH1F* hdtH = new TH1F("hdtH", "Hard Trigger b-Distribution", 100, 0., 20.); | |
341 | TH1F* mbmH = new TH1F("mbmH", "MB Trigger Multiplicity Distribution", 100, 0., 8000.); | |
342 | TH1F* hdmH = new TH1F("hdmH", "Hard Trigger Multiplicity Distribution", 100, 0., 8000.); | |
343 | ||
344 | mbtH->SetXTitle("b [fm]"); | |
345 | hdtH->SetXTitle("b [fm]"); | |
346 | mbmH->SetXTitle("Multiplicity"); | |
347 | hdmH->SetXTitle("Multiplicity"); | |
348 | ||
349 | TCanvas *c0 = new TCanvas("c0","Trigger Simulation",400,10,600,700); | |
350 | c0->Divide(2,1); | |
351 | TCanvas *c1 = new TCanvas("c1","Trigger Simulation",400,10,600,700); | |
352 | c1->Divide(1,2); | |
353 | ||
354 | // | |
355 | // | |
356 | Init(1); | |
357 | for (Int_t iev = 0; iev < n; iev++) | |
358 | { | |
359 | Float_t b, p, mult; | |
360 | GetRandom(b, p, mult); | |
361 | mbtH->Fill(b,1.); | |
362 | hdtH->Fill(b, p); | |
363 | mbmH->Fill(mult, 1.); | |
364 | hdmH->Fill(mult, p); | |
365 | ||
366 | c0->cd(1); | |
367 | mbtH->Draw(); | |
368 | c0->cd(2); | |
369 | hdtH->Draw(); | |
370 | c0->Update(); | |
371 | ||
372 | c1->cd(1); | |
373 | mbmH->Draw(); | |
374 | c1->cd(2); | |
375 | hdmH->Draw(); | |
376 | c1->Update(); | |
377 | } | |
378 | } | |
379 | ||
380 | void AliFastGlauber::GetRandom(Float_t& b, Float_t& p, Float_t& mult) | |
381 | { | |
382 | // | |
383 | // Gives back a random impact parameter, hard trigger probability and multiplicity | |
384 | // | |
385 | b = fWSgeo->GetRandom(); | |
386 | Float_t mu = fWSN->Eval(b); | |
387 | p = 1.-TMath::Exp(-mu); | |
88cb7938 | 388 | mult = 6000./fWSN->Eval(1.) * mu; |
5b3a5a5d | 389 | } |
390 | ||
c2140715 | 391 | void AliFastGlauber::GetRandom(Int_t& bin, Bool_t& hard) |
392 | { | |
393 | // | |
394 | // Gives back a random impact parameter bin, and hard trigger decission | |
395 | // | |
396 | Float_t b = fWSgeo->GetRandom(); | |
397 | Float_t mu = fWSN->Eval(b) * fSigmaHard; | |
398 | Float_t p = 1.-TMath::Exp(-mu); | |
399 | if (b < 5.) { | |
400 | bin = 1; | |
401 | } else if (b < 8.6) { | |
402 | bin = 2; | |
403 | } else if (b < 11.2) { | |
404 | bin = 3; | |
405 | } else if (b < 13.2) { | |
406 | bin = 4; | |
204e011f | 407 | } else if (b < 15.0) { |
c2140715 | 408 | bin = 5; |
204e011f | 409 | } else { |
410 | bin = 6; | |
c2140715 | 411 | } |
412 | ||
413 | hard = kFALSE; | |
414 | ||
415 | Float_t r = gRandom->Rndm(); | |
416 | ||
417 | if (r < p) hard = kTRUE; | |
418 | } | |
419 | ||
420 | ||
5b3a5a5d | 421 | Float_t AliFastGlauber::GetRandomImpactParameter(Float_t bmin, Float_t bmax) |
422 | { | |
423 | // | |
424 | // Gives back a random impact parameter in the range bmin .. bmax | |
425 | // | |
426 | ||
427 | Float_t b = -1.; | |
428 | while(b < bmin || b > bmax) | |
429 | b = fWSgeo->GetRandom(); | |
430 | return b; | |
431 | } | |
432 | ||
433 | Double_t AliFastGlauber::CrossSection(Double_t b1, Double_t b2) | |
434 | { | |
435 | // | |
436 | // Return cross-section integrated from b1 to b2 | |
437 | // | |
438 | ||
439 | return fWSgeo->Integral(b1, b2)/100.; | |
440 | } | |
441 | ||
442 | Double_t AliFastGlauber::FractionOfHardCrossSection(Double_t b1, Double_t b2) | |
443 | { | |
444 | // | |
445 | // Return raction of hard cross-section integrated from b1 to b2 | |
446 | // | |
447 | ||
448 | return fWSbinary->Integral(b1, b2)/fWSbinary->Integral(0., 100.); | |
449 | } | |
450 | ||
451 | ||
452 | Double_t AliFastGlauber::Binaries(Double_t b) | |
453 | { | |
454 | // | |
455 | // Return number of binary collisions normalized to 1 at b=0 | |
456 | // | |
457 | ||
270181a5 | 458 | return fWSN->Eval(b)/fWSN->Eval(0.); |
5b3a5a5d | 459 | } |