]>
Commit | Line | Data |
---|---|---|
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 | ||
16 | /* $Id$ */ | |
17 | ||
18 | ///////////////////////////////////////////////////////////////////////////// | |
19 | // // | |
20 | // Implementation of AliGenlib to collect parametrisations used for // | |
21 | // GSI simulations. // | |
22 | // It is an extension of AliMUONLib providing in addition the option // | |
23 | // for different parametrisations of pt, y and ip for every particle type // | |
24 | // // | |
25 | // Responsible: Andres.Sandoval@cern.ch // | |
26 | // // | |
27 | ///////////////////////////////////////////////////////////////////////////// | |
28 | ||
29 | #include "TMath.h" | |
30 | #include "TRandom.h" | |
31 | #include "TString.h" | |
32 | #include "AliGenGSIlib.h" | |
33 | ||
34 | ||
35 | ClassImp(AliGenGSIlib) | |
36 | ||
37 | //========================================================================== | |
38 | // | |
39 | // Definition of Particle Distributions | |
40 | // | |
41 | //========================================================================== | |
42 | // | |
43 | // Upsilon | |
44 | // | |
45 | //-------------------------------------------------------------------------- | |
46 | // | |
47 | // upsilon particle composition | |
48 | // | |
49 | //-------------------------------------------------------------------------- | |
50 | Int_t AliGenGSIlib::IpUpsilon(TRandom *) | |
51 | { | |
52 | // Return upsilon pdg code | |
53 | ||
54 | return 553; | |
55 | ||
56 | } | |
57 | Double_t AliGenGSIlib::PtUpsilonFlat( const Double_t *px, const Double_t */*dummy*/ ) | |
58 | { | |
59 | //-------------------------------------------------------------------------- | |
60 | // | |
61 | // upsilon pt-distribution FLAT | |
62 | // | |
63 | //____________________________________________________________-------------- | |
64 | ||
65 | const Double_t kptmin = 0.0; | |
66 | const Double_t kptmax = 15.0; | |
67 | Double_t x=*px; | |
68 | Double_t weight = 0.; | |
69 | ||
70 | if ((x > kptmin) && (x < kptmax)) weight = 1.; | |
71 | ||
72 | return weight; | |
73 | ||
74 | } | |
75 | Double_t AliGenGSIlib::YUpsilonFlat(const Double_t *py, const Double_t */*dummy*/) | |
76 | { | |
77 | //-------------------------------------------------------------------------- | |
78 | // | |
79 | // upsilon y-distribution FLAT | |
80 | // | |
81 | //-------------------------------------------------------------------------- | |
82 | ||
83 | const Double_t ky0 = 1.5; | |
84 | const Double_t kb=1.; | |
85 | Double_t yu; | |
86 | Double_t y=TMath::Abs(*py); | |
87 | ||
88 | if (y < ky0) | |
89 | yu=kb; | |
90 | else | |
91 | yu = 0.; | |
92 | ||
93 | return yu; | |
94 | ||
95 | } | |
96 | Double_t AliGenGSIlib::PtUpsilonRitman( const Double_t *px, const Double_t */*dummy*/ ) | |
97 | { | |
98 | //-------------------------------------------------------------------------- | |
99 | // | |
100 | // upsilon pt-distribution RITMAN | |
101 | // | |
102 | //-------------------------------------------------------------------------- | |
103 | ||
104 | const Double_t kpt0 = 4.7; | |
105 | const Double_t kxn = 3.5; | |
106 | Double_t x=*px; | |
107 | ||
108 | Double_t pass1 = 1.+((x*x)/(kpt0*kpt0)); | |
109 | ||
110 | return x/TMath::Power(pass1,kxn); | |
111 | ||
112 | } | |
113 | Double_t AliGenGSIlib::YUpsilonRitman(const Double_t *py, const Double_t */*dummy*/) | |
114 | { | |
115 | //-------------------------------------------------------------------------- | |
116 | // | |
117 | // upsilon y-distribution RITMAN | |
118 | // | |
119 | //-------------------------------------------------------------------------- | |
120 | ||
121 | const Double_t ky0 = 3.; | |
122 | const Double_t kb=1.; | |
123 | Double_t yu; | |
124 | Double_t y=TMath::Abs(*py); | |
125 | ||
126 | if (y < ky0) | |
127 | yu=kb; | |
128 | else | |
129 | yu=kb*TMath::Exp(-(y-ky0)*(y-ky0)/2); | |
130 | ||
131 | return yu; | |
132 | ||
133 | } | |
134 | Double_t AliGenGSIlib::PtUpsilonKarel( const Double_t */*px*/, const Double_t */*dummy*/ ) | |
135 | { | |
136 | //-------------------------------------------------------------------------- | |
137 | // | |
138 | // upsilon pt-distribution kAREL | |
139 | // | |
140 | //-------------------------------------------------------------------------- | |
141 | // to implement | |
142 | ||
143 | return 0.1; | |
144 | ||
145 | } | |
146 | Double_t AliGenGSIlib::YUpsilonKarel(const Double_t */*py*/, const Double_t */*dummy*/) | |
147 | { | |
148 | //-------------------------------------------------------------------------- | |
149 | // | |
150 | // upsilon y-distribution KAREL | |
151 | // | |
152 | //-------------------------------------------------------------------------- | |
153 | ||
154 | //to implement | |
155 | ||
156 | return 0.2; | |
157 | ||
158 | } | |
159 | Double_t AliGenGSIlib::PtUpsilonMUON( const Double_t *px, const Double_t */*dummy*/ ) | |
160 | { | |
161 | //-------------------------------------------------------------------------- | |
162 | // | |
163 | // upsilon pt-distribution MUONlib | |
164 | // | |
165 | //-------------------------------------------------------------------------- | |
166 | ||
167 | const Double_t kpt0 = 5.3; | |
168 | const Double_t kxn = 2.5; | |
169 | Double_t x=*px; | |
170 | ||
171 | Double_t pass1 = 1.+(x/kpt0)*(x/kpt0); | |
172 | ||
173 | return x/TMath::Power(pass1,kxn); | |
174 | ||
175 | } | |
176 | Double_t AliGenGSIlib::YUpsilonMUON(const Double_t *py, const Double_t */*dummy*/) | |
177 | { | |
178 | //-------------------------------------------------------------------------- | |
179 | // | |
180 | // upsilon y-distribution MUONlib | |
181 | // | |
182 | //-------------------------------------------------------------------------- | |
183 | ||
184 | const Double_t ky0 = 3.; | |
185 | const Double_t kb=1.; | |
186 | Double_t yu; | |
187 | Double_t y=TMath::Abs(*py); | |
188 | ||
189 | if (y < ky0) | |
190 | yu=kb; | |
191 | else | |
192 | yu=kb*TMath::Exp(-(y-ky0)*(y-ky0)/2); | |
193 | ||
194 | return yu; | |
195 | ||
196 | } | |
197 | //-------------------------------------------------------------------------- | |
198 | // | |
199 | // J/Psi | |
200 | // | |
201 | Int_t AliGenGSIlib::IpJpsi(TRandom *) | |
202 | { | |
203 | //-------------------------------------------------------------------------- | |
204 | // | |
205 | // J/Psi particle composition | |
206 | // | |
207 | //-------------------------------------------------------------------------- | |
208 | ||
209 | return 443; | |
210 | ||
211 | } | |
212 | Double_t AliGenGSIlib::PtJpsiFlat( const Double_t *px, const Double_t */*dummy*/ ) | |
213 | { | |
214 | //-------------------------------------------------------------------------- | |
215 | // | |
216 | // J/Psi pt-distribution FLAT | |
217 | // | |
218 | //-------------------------------------------------------------------------- | |
219 | ||
220 | const Double_t kptmin = 0.0; | |
221 | const Double_t kptmax = 15.0; | |
222 | Double_t x=*px; | |
223 | Double_t weight = 0.; | |
224 | ||
225 | if ((x > kptmin) && (x < kptmax)) weight = 1.; | |
226 | ||
227 | return weight; | |
228 | ||
229 | } | |
230 | Double_t AliGenGSIlib::YJpsiFlat(const Double_t *py, const Double_t */*dummy*/) | |
231 | { | |
232 | //-------------------------------------------------------------------------- | |
233 | // | |
234 | // J/Psi y-distribution FLAT | |
235 | // | |
236 | //-------------------------------------------------------------------------- | |
237 | ||
238 | const Double_t ky0 = 1.5; | |
239 | const Double_t kb=1.; | |
240 | Double_t yu; | |
241 | Double_t y=TMath::Abs(*py); | |
242 | ||
243 | if (y < ky0) | |
244 | yu=kb; | |
245 | else | |
246 | yu = 0.; | |
247 | ||
248 | return yu; | |
249 | ||
250 | } | |
251 | Double_t AliGenGSIlib::PtJpsiMUON( const Double_t *px, const Double_t */*dummy*/ ) | |
252 | { | |
253 | //-------------------------------------------------------------------------- | |
254 | // | |
255 | // J/Psi pt-distribution MUONlib | |
256 | // | |
257 | //-------------------------------------------------------------------------- | |
258 | ||
259 | const Double_t kpt0 = 4.; | |
260 | const Double_t kxn = 3.6; | |
261 | Double_t x=*px; | |
262 | ||
263 | Double_t pass1 = 1.+(x/kpt0)*(x/kpt0); | |
264 | return x/TMath::Power(pass1,kxn); | |
265 | ||
266 | } | |
267 | Double_t AliGenGSIlib::PtJpsiRitman( const Double_t *px, const Double_t */*dummy*/ ) | |
268 | { | |
269 | //-------------------------------------------------------------------------- | |
270 | // | |
271 | // J/Psi pt-distribution Ritman | |
272 | // | |
273 | //-------------------------------------------------------------------------- | |
274 | ||
275 | const Double_t kpt0 = 2.3; | |
276 | const Double_t kxn = 3.5; | |
277 | Double_t x=*px; | |
278 | ||
279 | Double_t pass1 = 1.+(x/kpt0)*(x/kpt0); | |
280 | ||
281 | return x/TMath::Power(pass1,kxn); | |
282 | ||
283 | } | |
284 | Double_t AliGenGSIlib::YJpsiMUON(const Double_t *py, const Double_t */*dummy*/) | |
285 | { | |
286 | //-------------------------------------------------------------------------- | |
287 | // | |
288 | // J/Psi y-distribution MUONlib | |
289 | // | |
290 | //-------------------------------------------------------------------------- | |
291 | ||
292 | const Double_t ky0 = 4.; | |
293 | const Double_t kb=1.; | |
294 | Double_t yj; | |
295 | Double_t y=TMath::Abs(*py); | |
296 | ||
297 | if (y < ky0) | |
298 | yj=kb; | |
299 | else | |
300 | yj=kb*TMath::Exp(-(y-ky0)*(y-ky0)/2); | |
301 | return yj; | |
302 | ||
303 | } | |
304 | //-------------------------------------------------------------------------- | |
305 | // | |
306 | // J/Psi pt-distribution by Sergei | |
307 | // | |
308 | //-------------------------------------------------------------------------- | |
309 | //Double_t AliGenGSIlib::PtJpsi( Double_t *px, Double_t */*dummy*/ ) | |
310 | //{ | |
311 | ||
312 | // return x = gRandom->Rndm()*10.; | |
313 | ||
314 | //} | |
315 | //-------------------------------------------------------------------------- | |
316 | // | |
317 | // J/Psi y-distribution by Sergei | |
318 | // | |
319 | //-------------------------------------------------------------------------- | |
320 | /*Double_t AliGenGSIlib::YJpsi(Double_t *py, Double_t *dummy) | |
321 | { | |
322 | ||
323 | const Double_t ky0 = 4.; | |
324 | const Double_t kb=1.; | |
325 | Double_t yj; | |
326 | Double_t y=TMath::Abs(*py); | |
327 | // | |
328 | if (y < ky0) | |
329 | yj=kb; | |
330 | else | |
331 | yj=kb*TMath::Exp(-(y-ky0)*(y-ky0)/2); | |
332 | return yj; | |
333 | ||
334 | } | |
335 | */ | |
336 | //-------------------------------------------------------------------------- | |
337 | // | |
338 | // Charm | |
339 | // | |
340 | //-------------------------------------------------------------------------- | |
341 | Int_t AliGenGSIlib::IpCharm(TRandom *ran) | |
342 | { | |
343 | // | |
344 | // charm particle composition | |
345 | // | |
346 | //-------------------------------------------------------------------------- | |
347 | ||
348 | Float_t random; | |
349 | Int_t ip; | |
350 | // 411,421,431,4122 | |
351 | random = ran->Rndm(); | |
352 | if (random < 0.5) { | |
353 | ip=411; | |
354 | } else if (random < 0.75) { | |
355 | ip=421; | |
356 | } else if (random < 0.90) { | |
357 | ip=431; | |
358 | } else { | |
359 | ip=4122; | |
360 | } | |
361 | if (ran->Rndm() < 0.5) {ip=-ip;} | |
362 | ||
363 | return ip; | |
364 | ||
365 | } | |
366 | Double_t AliGenGSIlib::PtCharmFlat( const Double_t *px, const Double_t */*dummy*/) | |
367 | { | |
368 | //-------------------------------------------------------------------------- | |
369 | // | |
370 | // charm pt-distribution, FLAT | |
371 | // | |
372 | //-------------------------------------------------------------------------- | |
373 | ||
374 | Double_t x=*px; | |
375 | ||
376 | if (x>10.) x = 0.; | |
377 | else x=1.; | |
378 | return x ; | |
379 | ||
380 | } | |
381 | Double_t AliGenGSIlib::PtCharmGSI( const Double_t *px, const Double_t */*dummy*/) | |
382 | { | |
383 | //-------------------------------------------------------------------------- | |
384 | // | |
385 | // charm pt-distribution, from Dariuzs Miskowiec | |
386 | // | |
387 | //-------------------------------------------------------------------------- | |
388 | ||
389 | //Taken from PYTHIA with MRS D-' (3031 from PDFLIB), K=3.0 | |
390 | const Double_t kp1 = 1.3; | |
391 | const Double_t kp2 = 0.39; | |
392 | const Double_t kp3 = 0.018; | |
393 | const Double_t kp4 = 0.91; | |
394 | Double_t x=*px; | |
395 | ||
396 | Double_t pass1 =TMath::Exp(-x/kp2) ; | |
397 | Double_t pass2 =TMath::Exp(-x/kp4) ; | |
398 | return TMath::Power(x,kp1) * (pass1 + kp3 * pass2); | |
399 | ||
400 | } | |
401 | Double_t AliGenGSIlib::PtCharmMUON( const Double_t *px, const Double_t */*dummy*/) | |
402 | { | |
403 | //-------------------------------------------------------------------------- | |
404 | // | |
405 | // charm pt-distribution, from MUONlib | |
406 | // | |
407 | //-------------------------------------------------------------------------- | |
408 | ||
409 | const Double_t kpt0 = 4.08; | |
410 | const Double_t kxn = 9.40; | |
411 | Double_t x=*px; | |
412 | ||
413 | Double_t pass1 = 1.+(x/kpt0)*(x/kpt0); | |
414 | ||
415 | return x/TMath::Power(pass1,kxn); | |
416 | ||
417 | } | |
418 | Double_t AliGenGSIlib::YCharm( const Double_t *px, const Double_t */*dummy*/) | |
419 | { | |
420 | //-------------------------------------------------------------------------- | |
421 | // | |
422 | // charm y-distribution | |
423 | // | |
424 | //-------------------------------------------------------------------------- | |
425 | ||
426 | Double_t *dum=0; | |
427 | ||
428 | return YJpsiMUON(px,dum); | |
429 | ||
430 | } | |
431 | //-------------------------------------------------------------------------- | |
432 | // | |
433 | // Beauty | |
434 | // | |
435 | //-------------------------------------------------------------------------- | |
436 | Int_t AliGenGSIlib::IpBeauty(TRandom *ran) | |
437 | { | |
438 | // | |
439 | // beauty particle composition | |
440 | // | |
441 | //-------------------------------------------------------------------------- | |
442 | ||
443 | Float_t random; | |
444 | Int_t ip; | |
445 | random = ran->Rndm(); | |
446 | if (random < 0.5) { | |
447 | ip=511; | |
448 | } else if (random < 0.75) { | |
449 | ip=521; | |
450 | } else if (random < 0.90) { | |
451 | ip=531; | |
452 | } else { | |
453 | ip=5122; | |
454 | } | |
455 | if (ran->Rndm() < 0.5) {ip=-ip;} | |
456 | ||
457 | return ip; | |
458 | } | |
459 | Double_t AliGenGSIlib::PtBeautyFlat( const Double_t *px, const Double_t */*dummy*/) | |
460 | { | |
461 | //-------------------------------------------------------------------------- | |
462 | // | |
463 | // beauty pt-distribution, FLAT | |
464 | // | |
465 | //-------------------------------------------------------------------------- | |
466 | ||
467 | Double_t x=*px; | |
468 | ||
469 | if (x>10.) x=0.; | |
470 | else x = 1.; | |
471 | return x ; | |
472 | ||
473 | } | |
474 | Double_t AliGenGSIlib::PtBeautyGSI( const Double_t *px, const Double_t */*dummy*/) | |
475 | { | |
476 | //-------------------------------------------------------------------------- | |
477 | // | |
478 | // | |
479 | // beauty pt-distribution, from D. Miskowiec | |
480 | // | |
481 | //-------------------------------------------------------------------------- | |
482 | ||
483 | //Taken from PYTHIA with MRS D-' (3031 from PDFLIB), K=3.0 | |
484 | const Double_t kp1 = 1.3; | |
485 | const Double_t kp2 = 1.78; | |
486 | const Double_t kp3 = 0.0096; | |
487 | const Double_t kp4 = 4.16; | |
488 | Double_t x=*px; | |
489 | ||
490 | Double_t pass1 =TMath::Exp(-x/kp2) ; | |
491 | Double_t pass2 =TMath::Exp(-x/kp4) ; | |
492 | ||
493 | return TMath::Power(x,kp1) * (pass1 + kp3 * pass2); | |
494 | ||
495 | } | |
496 | Double_t AliGenGSIlib::PtBeautyMUON( const Double_t *px, const Double_t */*dummy*/) | |
497 | { | |
498 | //-------------------------------------------------------------------------- | |
499 | // | |
500 | // beauty pt-distribution, from MUONlib | |
501 | // | |
502 | //-------------------------------------------------------------------------- | |
503 | ||
504 | const Double_t kpt0 = 4.; | |
505 | const Double_t kxn = 3.6; | |
506 | Double_t x=*px; | |
507 | ||
508 | Double_t pass1 = 1.+(x/kpt0)*(x/kpt0); | |
509 | ||
510 | return x/TMath::Power(pass1,kxn); | |
511 | ||
512 | } | |
513 | Double_t AliGenGSIlib::YBeauty( const Double_t *px, const Double_t */*dummy*/) | |
514 | { | |
515 | //-------------------------------------------------------------------------- | |
516 | // | |
517 | // beauty y-distribution | |
518 | // | |
519 | //-------------------------------------------------------------------------- | |
520 | ||
521 | Double_t *dum=0; | |
522 | ||
523 | return YJpsiMUON(px,dum); | |
524 | ||
525 | } | |
526 | //-------------------------------------------------------------------------- | |
527 | // | |
528 | // Eta | |
529 | // | |
530 | //-------------------------------------------------------------------------- | |
531 | Int_t AliGenGSIlib::IpEta(TRandom *) | |
532 | { | |
533 | // | |
534 | // eta particle composition | |
535 | // | |
536 | //-------------------------------------------------------------------------- | |
537 | ||
538 | return 221; | |
539 | ||
540 | } | |
541 | Double_t AliGenGSIlib::PtEtaPHOS( const Double_t *px, const Double_t */*dummy*/ ) | |
542 | { | |
543 | //-------------------------------------------------------------------------- | |
544 | // | |
545 | // eta pt-distribution | |
546 | // | |
547 | //____________________________________________________________-------------- | |
548 | ||
549 | return PtScal(*px,2); // 2==> Eta in the PtScal function | |
550 | ||
551 | } | |
552 | Double_t AliGenGSIlib::YEtaPHOS(const Double_t *py, const Double_t */*dummy*/) | |
553 | { | |
554 | //-------------------------------------------------------------------------- | |
555 | // | |
556 | // eta y-distribution | |
557 | // | |
558 | //-------------------------------------------------------------------------- | |
559 | ||
560 | const Double_t ka = 1000.; | |
561 | const Double_t kdy = 4.; | |
562 | ||
563 | Double_t y=TMath::Abs(*py); | |
564 | ||
565 | Double_t ex = y*y/(2*kdy*kdy); | |
566 | ||
567 | return ka*TMath::Exp(-ex); | |
568 | ||
569 | } | |
570 | //-------------------------------------------------------------------------- | |
571 | // | |
572 | // Etaprime | |
573 | // | |
574 | //-------------------------------------------------------------------------- | |
575 | Int_t AliGenGSIlib::IpEtaprime(TRandom *) | |
576 | { | |
577 | // | |
578 | // etaprime particle composition | |
579 | // | |
580 | //-------------------------------------------------------------------------- | |
581 | ||
582 | return 331; | |
583 | ||
584 | } | |
585 | Double_t AliGenGSIlib::PtEtaprimePHOS( const Double_t *px, const Double_t */*dummy*/ ) | |
586 | { | |
587 | //-------------------------------------------------------------------------- | |
588 | // | |
589 | // etaprime pt-distribution | |
590 | // | |
591 | //____________________________________________________________-------------- | |
592 | ||
593 | return PtScal(*px,4); // 4==> Etaprime in the PtScal function | |
594 | ||
595 | } | |
596 | Double_t AliGenGSIlib::YEtaprimePHOS(const Double_t *py, const Double_t */*dummy*/) | |
597 | { | |
598 | //-------------------------------------------------------------------------- | |
599 | // | |
600 | // etaprime y-distribution | |
601 | // | |
602 | //-------------------------------------------------------------------------- | |
603 | ||
604 | const Double_t ka = 1000.; | |
605 | const Double_t kdy = 4.; | |
606 | ||
607 | Double_t y=TMath::Abs(*py); | |
608 | ||
609 | Double_t ex = y*y/(2*kdy*kdy); | |
610 | ||
611 | return ka*TMath::Exp(-ex); | |
612 | ||
613 | } | |
614 | //-------------------------------------------------------------------------- | |
615 | // | |
616 | // omega | |
617 | // | |
618 | //-------------------------------------------------------------------------- | |
619 | Int_t AliGenGSIlib::IpOmega(TRandom *) | |
620 | { | |
621 | // | |
622 | // omega particle composition | |
623 | // | |
624 | //-------------------------------------------------------------------------- | |
625 | ||
626 | return 223; | |
627 | ||
628 | } | |
629 | Double_t AliGenGSIlib::PtOmega( const Double_t *px, const Double_t */*dummy*/ ) | |
630 | { | |
631 | //-------------------------------------------------------------------------- | |
632 | // | |
633 | // omega pt-distribution | |
634 | // | |
635 | //____________________________________________________________-------------- | |
636 | ||
637 | return PtScal(*px,3); // 3==> Omega in the PtScal function | |
638 | ||
639 | } | |
640 | Double_t AliGenGSIlib::YOmega(const Double_t *py, const Double_t */*dummy*/) | |
641 | { | |
642 | //-------------------------------------------------------------------------- | |
643 | // | |
644 | // omega y-distribution | |
645 | // | |
646 | //-------------------------------------------------------------------------- | |
647 | ||
648 | const Double_t ka = 1000.; | |
649 | const Double_t kdy = 4.; | |
650 | ||
651 | ||
652 | Double_t y=TMath::Abs(*py); | |
653 | ||
654 | Double_t ex = y*y/(2*kdy*kdy); | |
655 | ||
656 | return ka*TMath::Exp(-ex); | |
657 | ||
658 | } | |
659 | //-------------------------------------------------------------------------- | |
660 | // | |
661 | // Rho | |
662 | // | |
663 | //-------------------------------------------------------------------------- | |
664 | ||
665 | Int_t AliGenGSIlib::IpRho(TRandom *) | |
666 | { | |
667 | // | |
668 | // rho particle composition | |
669 | // | |
670 | //-------------------------------------------------------------------------- | |
671 | ||
672 | return 113; | |
673 | ||
674 | } | |
675 | Double_t AliGenGSIlib::PtRho( const Double_t *px, const Double_t */*dummy*/ ) | |
676 | { | |
677 | //-------------------------------------------------------------------------- | |
678 | // | |
679 | // rho pt-distribution | |
680 | // | |
681 | //____________________________________________________________-------------- | |
682 | ||
683 | return PtScal(*px,10); // 10==> Rho in the PtScal function | |
684 | ||
685 | } | |
686 | Double_t AliGenGSIlib::YRho(const Double_t *py, const Double_t */*dummy*/) | |
687 | { | |
688 | //-------------------------------------------------------------------------- | |
689 | // | |
690 | // rho y-distribution | |
691 | // | |
692 | //-------------------------------------------------------------------------- | |
693 | ||
694 | const Double_t ka = 1000.; | |
695 | const Double_t kdy = 4.; | |
696 | ||
697 | ||
698 | Double_t y=TMath::Abs(*py); | |
699 | ||
700 | Double_t ex = y*y/(2*kdy*kdy); | |
701 | ||
702 | return ka*TMath::Exp(-ex); | |
703 | ||
704 | } | |
705 | //-------------------------------------------------------------------------- | |
706 | // | |
707 | // Pion | |
708 | // | |
709 | //-------------------------------------------------------------------------- | |
710 | Int_t AliGenGSIlib::IpPionPHOS(TRandom *ran) | |
711 | { | |
712 | // | |
713 | // particle composition pi+, pi0, pi- | |
714 | // | |
715 | //-------------------------------------------------------------------------- | |
716 | ||
717 | Float_t random = ran->Rndm(); | |
718 | ||
719 | if ( (3.*random) < 1. ) | |
720 | { | |
721 | return 211 ; | |
722 | } | |
723 | else | |
724 | { | |
725 | if ( (3.*random) >= 2.) | |
726 | { | |
727 | return -211 ; | |
728 | } | |
729 | else | |
730 | { | |
731 | return 111 ; | |
732 | } | |
733 | } | |
734 | } | |
735 | Double_t AliGenGSIlib::PtPion( const Double_t *px, const Double_t */*dummy*/ ) | |
736 | { | |
737 | //-------------------------------------------------------------------------- | |
738 | // | |
739 | // pion pt-distribution | |
740 | // | |
741 | // Pion transverse momentum distribtuion as in AliGenMUONlib class, | |
742 | // version 3.01 of aliroot: | |
743 | // PT-PARAMETERIZATION CDF, PRL 61(88) 1819 | |
744 | // POWER LAW FOR PT > 500 MEV | |
745 | // MT SCALING BELOW (T=160 MEV) | |
746 | // | |
747 | //____________________________________________________________-------------- | |
748 | ||
749 | const Double_t kp0 = 1.3; | |
750 | const Double_t kxn = 8.28; | |
751 | const Double_t kxlim=0.5; | |
752 | const Double_t kt=0.160; | |
753 | const Double_t kxmpi=0.139; | |
754 | const Double_t kb=1.; | |
755 | Double_t y, y1, kxmpi2, ynorm, a; | |
756 | Double_t x=*px; | |
757 | ||
758 | y1=TMath::Power(kp0/(kp0+kxlim),kxn); | |
759 | kxmpi2=kxmpi*kxmpi; | |
760 | ynorm=kb*(TMath::Exp(-sqrt(kxlim*kxlim+kxmpi2)/kt)); | |
761 | a=ynorm/y1; | |
762 | if (x > kxlim) | |
763 | y=a*TMath::Power(kp0/(kp0+x),kxn); | |
764 | else | |
765 | y=kb*TMath::Exp(-sqrt(x*x+kxmpi2)/kt); | |
766 | return y*x; | |
767 | ||
768 | } | |
769 | Double_t AliGenGSIlib::YPion(const Double_t *py, const Double_t */*dummy*/) | |
770 | { | |
771 | //-------------------------------------------------------------------------- | |
772 | // | |
773 | // pion y-distribution | |
774 | // | |
775 | //-------------------------------------------------------------------------- | |
776 | ||
777 | const Double_t ka = 7000.; | |
778 | const Double_t kdy = 4.; | |
779 | ||
780 | Double_t y=TMath::Abs(*py); | |
781 | ||
782 | Double_t ex = y*y/(2*kdy*kdy); | |
783 | ||
784 | return ka*TMath::Exp(-ex); | |
785 | ||
786 | } | |
787 | Int_t AliGenGSIlib::IpKaonPHOS(TRandom *ran) | |
788 | { | |
789 | //-------------------------------------------------------------------------- | |
790 | // | |
791 | // | |
792 | // Kaon | |
793 | //-------------------------------------------------------------------------- | |
794 | // | |
795 | // kaon particle composition K+, K-, Ko_short, Ko_long | |
796 | // | |
797 | //-------------------------------------------------------------------------- | |
798 | ||
799 | Float_t random = ran->Rndm(); | |
800 | Float_t random2 = ran->Rndm(); | |
801 | if (random2 < 0.5) | |
802 | { | |
803 | if (random < 0.5) { | |
804 | return 321; // K+ | |
805 | } else { | |
806 | return -321; // K- | |
807 | } | |
808 | } | |
809 | else | |
810 | { | |
811 | if (random < 0.5) { | |
812 | return 130; // K^0 short | |
813 | } else { | |
814 | return 310; // K^0 long | |
815 | } | |
816 | } | |
817 | } | |
818 | Double_t AliGenGSIlib::PtKaonPHOS( const Double_t *px, const Double_t */*dummy*/ ) | |
819 | { | |
820 | //-------------------------------------------------------------------------- | |
821 | // | |
822 | // kaon pt-distribution | |
823 | // | |
824 | //____________________________________________________________-------------- | |
825 | ||
826 | return PtScal(*px,1); // 1==> Kaon in the PtScal function | |
827 | ||
828 | } | |
829 | Double_t AliGenGSIlib::YKaonPHOS(const Double_t *py, const Double_t */*dummy*/) | |
830 | { | |
831 | //-------------------------------------------------------------------------- | |
832 | // | |
833 | // kaon y-distribution | |
834 | // | |
835 | //-------------------------------------------------------------------------- | |
836 | ||
837 | const Double_t ka = 1000.; | |
838 | const Double_t kdy = 4.; | |
839 | ||
840 | Double_t y=TMath::Abs(*py); | |
841 | ||
842 | Double_t ex = y*y/(2*kdy*kdy); | |
843 | ||
844 | return ka*TMath::Exp(-ex); | |
845 | ||
846 | } | |
847 | //-------------------------------------------------------------------------- | |
848 | // | |
849 | // Phi | |
850 | // | |
851 | Int_t AliGenGSIlib::IpPhi(TRandom *) | |
852 | { | |
853 | //-------------------------------------------------------------------------- | |
854 | // | |
855 | // particle composition | |
856 | // | |
857 | //-------------------------------------------------------------------------- | |
858 | ||
859 | return 333; | |
860 | ||
861 | } | |
862 | Double_t AliGenGSIlib::PtPhiPHOS( const Double_t *px, const Double_t */*dummy*/ ) | |
863 | { | |
864 | //-------------------------------------------------------------------------- | |
865 | // | |
866 | // phi pt-distribution | |
867 | // | |
868 | //____________________________________________________________-------------- | |
869 | ||
870 | return PtScal(*px,5); // 5==> Phi in the PtScal function | |
871 | ||
872 | } | |
873 | Double_t AliGenGSIlib::YPhiPHOS(const Double_t *py, const Double_t */*dummy*/) | |
874 | { | |
875 | //-------------------------------------------------------------------------- | |
876 | // | |
877 | // phi y-distribution | |
878 | // | |
879 | //-------------------------------------------------------------------------- | |
880 | ||
881 | const Double_t ka = 1000.; | |
882 | const Double_t kdy = 4.; | |
883 | ||
884 | ||
885 | Double_t y=TMath::Abs(*py); | |
886 | ||
887 | Double_t ex = y*y/(2*kdy*kdy); | |
888 | ||
889 | return ka*TMath::Exp(-ex); | |
890 | ||
891 | } | |
892 | Int_t AliGenGSIlib::IpBaryons(TRandom *ran) | |
893 | { | |
894 | //-------------------------------------------------------------------------- | |
895 | // | |
896 | // Baryons | |
897 | // | |
898 | //-------------------------------------------------------------------------- | |
899 | // | |
900 | // baryons particle composition p, pbar, n, nbar | |
901 | // | |
902 | //-------------------------------------------------------------------------- | |
903 | ||
904 | Float_t random = ran->Rndm(); | |
905 | Float_t random2 = ran->Rndm(); | |
906 | if (random2 < 0.5) | |
907 | { | |
908 | if (random < 0.5) { | |
909 | return 2212; // p | |
910 | } else { | |
911 | return -2212; // pbar | |
912 | } | |
913 | } | |
914 | else | |
915 | { | |
916 | if (random < 0.5) { | |
917 | return 2112; // n | |
918 | } else { | |
919 | return -2112; // n bar | |
920 | } | |
921 | } | |
922 | } | |
923 | Double_t AliGenGSIlib::PtBaryons( const Double_t *px, const Double_t */*dummy*/ ) | |
924 | { | |
925 | //-------------------------------------------------------------------------- | |
926 | // | |
927 | // baryons pt-distribution | |
928 | // | |
929 | //____________________________________________________________-------------- | |
930 | ||
931 | return PtScal(*px,6); // 6==> Baryon in the PtScal function | |
932 | ||
933 | } | |
934 | Double_t AliGenGSIlib::YBaryons(const Double_t *py, const Double_t */*dummy*/) | |
935 | { | |
936 | //-------------------------------------------------------------------------- | |
937 | // | |
938 | // baryons y-distribution | |
939 | // | |
940 | //-------------------------------------------------------------------------- | |
941 | ||
942 | const Double_t ka = 1000.; | |
943 | const Double_t kdy = 4.; | |
944 | ||
945 | Double_t y=TMath::Abs(*py); | |
946 | ||
947 | Double_t ex = y*y/(2*kdy*kdy); | |
948 | ||
949 | return ka*TMath::Exp(-ex); | |
950 | ||
951 | } | |
952 | //============================================================= | |
953 | // | |
954 | // Mt-scaling as in AliGenPHOSlib | |
955 | // | |
956 | //============================================================= | |
957 | // | |
958 | Double_t AliGenGSIlib::PtScal(Double_t pt, Int_t np) | |
959 | { | |
960 | // Function for the calculation of the Pt distribution for a | |
961 | // given particle np, from the pion Pt distribution using the | |
962 | // mt scaling. | |
963 | ||
964 | // It was taken from AliGenPHOSlib aliroot version 3.04, which | |
965 | // is an update of the one in AliGenMUONlib aliroot version 3.01 | |
966 | // with an extension for Baryons but supressing Rhos | |
967 | // np = 1=>Pions 2=>Kaons 3=>Etas 4=>Omegas 5=>ETA' 6=>PHI | |
968 | // 7=>BARYONS-BARYONBARS | |
969 | ||
970 | // The present adds the Rhos | |
971 | ||
972 | // MASS 0=>PI, 1=>K, 2=>ETA, 3=>OMEGA, 4=>ETA', 5=>PHI | |
973 | // 6=>BARYON-BARYONBAR, 10==>RHO | |
974 | ||
975 | const Double_t khm[11] = {0.1396, 0.494, 0.547, 0.782, 0.957, 1.02, | |
976 | 0.938, 0. , 0., 0., 0.769}; | |
977 | ||
978 | // VALUE MESON/PI AT 5 GEV | |
979 | ||
980 | const Double_t kfmax[11]={1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1.}; | |
981 | ||
982 | Double_t f5=TMath::Power(((sqrt(100.018215)+2.)/(sqrt(100.+khm[np]*khm[np])+2.0)),12.3); | |
983 | Double_t kfmax2=f5/kfmax[np]; | |
984 | // PIONS | |
985 | Double_t ptpion=100.*PtPion(&pt, (Double_t*) 0); | |
986 | Double_t fmtscal=TMath::Power(((sqrt(pt*pt+0.018215)+2.)/ | |
987 | (sqrt(pt*pt+khm[np]*khm[np])+2.0)),12.3)/ kfmax2; | |
988 | return fmtscal*ptpion; | |
989 | } | |
990 | ||
991 | //========================================================================== | |
992 | // | |
993 | // Set Getters | |
994 | // | |
995 | //========================================================================== | |
996 | ||
997 | typedef Double_t (*GenFunc) (const Double_t*, const Double_t*); | |
998 | ||
999 | typedef Int_t (*GenFuncIp) (TRandom *); | |
1000 | ||
1001 | GenFunc AliGenGSIlib::GetPt(Int_t param, const char * tname) const | |
1002 | { | |
1003 | // Return pointer to pT parameterisation | |
1004 | GenFunc func=0; | |
1005 | TString sname(tname); | |
1006 | ||
1007 | switch (param) | |
1008 | { | |
1009 | case kUpsilon: | |
1010 | if (sname=="FLAT"){ | |
1011 | func= PtUpsilonFlat; | |
1012 | break; | |
1013 | } | |
1014 | if (sname=="MUON"){ | |
1015 | func= PtUpsilonMUON; | |
1016 | break; | |
1017 | } | |
1018 | if (sname=="RITMAN"){ | |
1019 | func=PtUpsilonRitman; | |
1020 | break; | |
1021 | } | |
1022 | if (sname=="KAREL"){ | |
1023 | func=PtUpsilonKarel; | |
1024 | break; | |
1025 | } | |
1026 | func=0; | |
1027 | printf("<AliGenGSIlib::GetPt> unknown parametrisation\n"); | |
1028 | break; | |
1029 | ||
1030 | case kJPsi: | |
1031 | if (sname=="FLAT"){ | |
1032 | func= PtJpsiFlat; | |
1033 | break; | |
1034 | } | |
1035 | if (sname=="MUON"){ | |
1036 | func= PtJpsiMUON; | |
1037 | break; | |
1038 | } | |
1039 | // if (sname=="SERGEI"){ | |
1040 | // func= PtJpsi; | |
1041 | // break; | |
1042 | // } | |
1043 | func=0; | |
1044 | printf("<AliGenGSIlib::GetPt> unknown parametrisation\n"); | |
1045 | break; | |
1046 | ||
1047 | case kCharm: | |
1048 | if (sname=="FLAT"){ | |
1049 | func= PtCharmFlat; | |
1050 | break; | |
1051 | } | |
1052 | ||
1053 | if (sname=="MUON"){ | |
1054 | func= PtCharmMUON; | |
1055 | break; | |
1056 | } | |
1057 | ||
1058 | if (sname=="GSI"){ | |
1059 | func= PtCharmGSI; | |
1060 | break; | |
1061 | } | |
1062 | func=0; | |
1063 | printf("<AliGenGSIlib::GetPt> unknown parametrisation\n"); | |
1064 | break; | |
1065 | ||
1066 | case kBeauty: | |
1067 | if (sname=="FLAT"){ | |
1068 | func= PtBeautyFlat; | |
1069 | break; | |
1070 | } | |
1071 | if (sname=="MUON"){ | |
1072 | func= PtBeautyMUON; | |
1073 | break; | |
1074 | } | |
1075 | if (sname=="GSI"){ | |
1076 | func= PtBeautyGSI; | |
1077 | break; | |
1078 | } | |
1079 | func=0; | |
1080 | printf("<AliGenGSIlib::GetPt> unknown parametrisation\n"); | |
1081 | break; | |
1082 | ||
1083 | ||
1084 | case kEta: | |
1085 | func=PtEtaPHOS; | |
1086 | break; | |
1087 | ||
1088 | case kEtaprime: | |
1089 | func=PtEtaprimePHOS; | |
1090 | break; | |
1091 | ||
1092 | case kOmega: | |
1093 | func=PtOmega; | |
1094 | break; | |
1095 | ||
1096 | case kRho: | |
1097 | func=PtRho; | |
1098 | break; | |
1099 | ||
1100 | case kKaon: | |
1101 | func=PtKaonPHOS; | |
1102 | break; | |
1103 | ||
1104 | case kPion: | |
1105 | func=PtPion; | |
1106 | break; | |
1107 | ||
1108 | case kPhi: | |
1109 | func=PtPhiPHOS; | |
1110 | break; | |
1111 | ||
1112 | // case kLambda: | |
1113 | // func=PtLambda; | |
1114 | // break; | |
1115 | ||
1116 | case kBaryons: | |
1117 | func=PtBaryons; | |
1118 | break; | |
1119 | ||
1120 | default: | |
1121 | func=0; | |
1122 | printf("<AliGenGSIlib::GetPt> unknown parametrisation\n"); | |
1123 | } | |
1124 | return func; | |
1125 | } | |
1126 | ||
1127 | GenFunc AliGenGSIlib::GetY(Int_t param, const char * tname) const | |
1128 | { | |
1129 | // Return pointer to y- parameterisation | |
1130 | GenFunc func=0; | |
1131 | TString sname(tname); | |
1132 | ||
1133 | switch (param) | |
1134 | { | |
1135 | case kUpsilon: | |
1136 | if (sname=="FLAT"){ | |
1137 | func= YUpsilonFlat; | |
1138 | break; | |
1139 | } | |
1140 | if (sname=="MUON"){ | |
1141 | func= YUpsilonMUON; | |
1142 | break; | |
1143 | } | |
1144 | if (sname=="RITMAN"){ | |
1145 | func=YUpsilonRitman; | |
1146 | break; | |
1147 | } | |
1148 | if (sname=="KAREL"){ | |
1149 | func=YUpsilonKarel; | |
1150 | break; | |
1151 | } | |
1152 | func=0; | |
1153 | printf("<AliGenGSIlib::GetY> unknown parametrisation\n"); | |
1154 | break; | |
1155 | ||
1156 | case kJPsi: | |
1157 | if (sname=="FLAT"){ | |
1158 | func= YJpsiFlat; | |
1159 | break; | |
1160 | } | |
1161 | if (sname=="MUON"){ | |
1162 | func= YJpsiMUON; | |
1163 | break; | |
1164 | } | |
1165 | // if (sname=="SERGEI"){ | |
1166 | // func= YJpsi; | |
1167 | // break; | |
1168 | // } | |
1169 | ||
1170 | func=0; | |
1171 | printf("<AliGenGSIlib::GetY> unknown parametrisation\n"); | |
1172 | break; | |
1173 | ||
1174 | case kCharm: | |
1175 | func= YCharm; | |
1176 | break; | |
1177 | ||
1178 | case kBeauty: | |
1179 | func= YBeauty; | |
1180 | break; | |
1181 | ||
1182 | case kEta: | |
1183 | func=YEtaPHOS; | |
1184 | break; | |
1185 | ||
1186 | case kEtaprime: | |
1187 | func=YEtaprimePHOS; | |
1188 | break; | |
1189 | ||
1190 | case kOmega: | |
1191 | func=YOmega; | |
1192 | break; | |
1193 | ||
1194 | case kRho: | |
1195 | func=YRho; | |
1196 | break; | |
1197 | ||
1198 | case kKaon: | |
1199 | func=YKaonPHOS; | |
1200 | break; | |
1201 | ||
1202 | case kPion: | |
1203 | func=YPion; | |
1204 | break; | |
1205 | ||
1206 | case kPhi: | |
1207 | func=YPhiPHOS; | |
1208 | break; | |
1209 | ||
1210 | // case kLambda: | |
1211 | // func=YLambda; | |
1212 | // break; | |
1213 | ||
1214 | case kBaryons: | |
1215 | func=YBaryons; | |
1216 | break; | |
1217 | ||
1218 | default: | |
1219 | func=0; | |
1220 | printf("<AliGenGSIlib::GetY> unknown parametrisation\n"); | |
1221 | } | |
1222 | return func; | |
1223 | } | |
1224 | ||
1225 | GenFuncIp AliGenGSIlib::GetIp(Int_t param, const char * tname) const | |
1226 | { | |
1227 | // Return pointer to particle type parameterisation | |
1228 | GenFuncIp func=0; | |
1229 | TString sname(tname); | |
1230 | ||
1231 | switch (param) | |
1232 | { | |
1233 | case kUpsilon: | |
1234 | func= IpUpsilon; | |
1235 | break; | |
1236 | ||
1237 | case kJPsi: | |
1238 | func= IpJpsi; | |
1239 | break; | |
1240 | ||
1241 | case kCharm: | |
1242 | func= IpCharm; | |
1243 | break; | |
1244 | ||
1245 | case kBeauty: | |
1246 | func= IpBeauty; | |
1247 | break; | |
1248 | ||
1249 | case kEta: | |
1250 | func=IpEta; | |
1251 | break; | |
1252 | ||
1253 | case kEtaprime: | |
1254 | func=IpEtaprime; | |
1255 | break; | |
1256 | ||
1257 | case kOmega: | |
1258 | func=IpOmega; | |
1259 | break; | |
1260 | ||
1261 | case kRho: | |
1262 | func=IpRho; | |
1263 | break; | |
1264 | ||
1265 | case kKaon: | |
1266 | func=IpKaonPHOS; | |
1267 | break; | |
1268 | ||
1269 | case kPion: | |
1270 | func=IpPionPHOS; | |
1271 | break; | |
1272 | ||
1273 | case kPhi: | |
1274 | func=IpPhi; | |
1275 | break; | |
1276 | ||
1277 | // case kLambda: | |
1278 | // func=IpLambda; | |
1279 | // break; | |
1280 | ||
1281 | case kBaryons: | |
1282 | func=IpBaryons; | |
1283 | break; | |
1284 | ||
1285 | default: | |
1286 | func=0; | |
1287 | printf("<AliGenGSIlib::GetIp> unknown parametrisation\n"); | |
1288 | } | |
1289 | return func; | |
1290 | } | |
1291 | ||
1292 | ||
1293 | ||
1294 | ||
1295 | ||
1296 | ||
1297 | ||
1298 | ||
1299 | ||
1300 | ||
1301 | ||
1302 |