EVGEN/AliGenGSIlib.cxx

   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( Double_t *px, 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 (kptmin < x < kptmax) weight = 1.;
  71
  72   return weight;
  73
  74 }
  75 Double_t AliGenGSIlib::YUpsilonFlat(Double_t *py, 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( Double_t *px, 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(Double_t *py, 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( Double_t */*px*/, 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(Double_t */*py*/, 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( Double_t *px, 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(Double_t *py, 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( Double_t *px, 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 (kptmin < x < kptmax) weight = 1.;
 226
 227   return weight;
 228
 229 }
 230 Double_t AliGenGSIlib::YJpsiFlat(Double_t *py, 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( Double_t *px, 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( Double_t *px, 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(Double_t *py, 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( Double_t *px, 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( Double_t *px, 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( Double_t *px, 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( Double_t *px, 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( Double_t *px, 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( Double_t *px, 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( Double_t *px, 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( Double_t *px, 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( Double_t *px, Double_t */*dummy*/ )
 542 {
 543 //--------------------------------------------------------------------------
 544 //
 545 //                  eta pt-distribution
 546 //
 547 //____________________________________________________________--------------
 548
 549   return PtScal(*px,3);  //  3==> Eta in the PtScal function
 550
 551 }
 552 Double_t AliGenGSIlib::YEtaPHOS(Double_t *py, 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( Double_t *px, Double_t */*dummy*/ )
 586 {
 587 //--------------------------------------------------------------------------
 588 //
 589 //                 etaprime pt-distribution
 590 //
 591 //____________________________________________________________--------------
 592
 593   return PtScal(*px,5);  //  5==> Etaprime in the PtScal function
 594
 595 }
 596 Double_t AliGenGSIlib::YEtaprimePHOS(Double_t *py, 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( Double_t *px, Double_t */*dummy*/ )
 630 {
 631 //--------------------------------------------------------------------------
 632 //
 633 //                  omega pt-distribution
 634 //
 635 //____________________________________________________________--------------
 636
 637   return PtScal(*px,4);  //  4==> Omega in the PtScal function
 638
 639 }
 640 Double_t AliGenGSIlib::YOmega(Double_t *py, 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( Double_t *px, Double_t */*dummy*/ )
 676 {
 677 //--------------------------------------------------------------------------
 678 //
 679 //                  rho pt-distribution
 680 //
 681 //____________________________________________________________--------------
 682
 683   return PtScal(*px,11);  //  11==> Rho in the PtScal function
 684
 685 }
 686 Double_t AliGenGSIlib::YRho(Double_t *py, 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( Double_t *px, 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(Double_t *py, 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( Double_t *px, Double_t */*dummy*/ )
 819 {
 820 //--------------------------------------------------------------------------
 821 //
 822 //                   kaon pt-distribution
 823 //
 824 //____________________________________________________________--------------
 825
 826   return PtScal(*px,2);  //  2==> Kaon in the PtScal function
 827
 828 }
 829 Double_t AliGenGSIlib::YKaonPHOS(Double_t *py, 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( Double_t *px, Double_t */*dummy*/ )
 863 {
 864 //--------------------------------------------------------------------------
 865 //
 866 //                   phi pt-distribution
 867 //
 868 //____________________________________________________________--------------
 869
 870   return PtScal(*px,6);  //  6==> Phi in the PtScal function
 871
 872 }
 873 Double_t AliGenGSIlib::YPhiPHOS(Double_t *py, 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( Double_t *px, Double_t */*dummy*/ )
 924 {
 925 //--------------------------------------------------------------------------
 926 //
 927 //                  baryons pt-distribution
 928 //
 929 //____________________________________________________________--------------
 930
 931   return PtScal(*px,7);  //  7==> Baryon in the PtScal function
 932
 933 }
 934 Double_t AliGenGSIlib::YBaryons(Double_t *py, 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   1=>PI, 2=>K, 3=>ETA, 4=>OMEGA, 5=>ETA', 6=>PHI
 973 //        7=>BARYON-BARYONBAR, 11==>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   np--;
 983   Double_t f5=TMath::Power(((sqrt(100.018215)+2.)/(sqrt(100.+khm[np]*khm[np])+2.0)),12.3);
 984   Double_t kfmax2=f5/kfmax[np];
 985   // PIONS
 986   Double_t ptpion=100.*PtPion(&pt, (Double_t*) 0);
 987   Double_t fmtscal=TMath::Power(((sqrt(pt*pt+0.018215)+2.)/
 988                                  (sqrt(pt*pt+khm[np]*khm[np])+2.0)),12.3)/ kfmax2;
 989   return fmtscal*ptpion;
 990 }
 991
 992 //==========================================================================
 993 //
 994 //                     Set Getters
 995 //
 996 //==========================================================================
 997
 998 typedef Double_t (*GenFunc) (Double_t*,  Double_t*);
 999
1000 typedef Int_t (*GenFuncIp) (TRandom *);
1001
1002 GenFunc AliGenGSIlib::GetPt(Int_t param, const char * tname) const
1003 {
1004 // Return pointer to pT parameterisation
1005    GenFunc func=0;
1006    TString sname(tname);
1007
1008    switch (param)
1009     {
1010     case kUpsilon:
1011       if (sname=="FLAT"){
1012         func= PtUpsilonFlat;
1013         break;
1014       }
1015       if (sname=="MUON"){
1016         func= PtUpsilonMUON;
1017         break;
1018       }
1019       if (sname=="RITMAN"){
1020         func=PtUpsilonRitman;
1021         break;
1022       }
1023       if (sname=="KAREL"){
1024         func=PtUpsilonKarel;
1025         break;
1026       }
1027       func=0;
1028       printf("<AliGenGSIlib::GetPt> unknown parametrisation\n");
1029       break;
1030
1031     case kJPsi:
1032       if (sname=="FLAT"){
1033         func= PtJpsiFlat;
1034         break;
1035       }
1036       if (sname=="MUON"){
1037         func= PtJpsiMUON;
1038         break;
1039       }
1040       //      if (sname=="SERGEI"){
1041       //        func= PtJpsi;
1042       //        break;
1043       //     }
1044       func=0;
1045       printf("<AliGenGSIlib::GetPt> unknown parametrisation\n");
1046       break;
1047
1048     case kCharm:
1049       if (sname=="FLAT"){
1050         func= PtCharmFlat;
1051         break;
1052       }
1053
1054       if (sname=="MUON"){
1055         func= PtCharmMUON;
1056         break;
1057       }
1058
1059       if (sname=="GSI"){
1060         func= PtCharmGSI;
1061         break;
1062       }
1063       func=0;
1064       printf("<AliGenGSIlib::GetPt> unknown parametrisation\n");
1065       break;
1066
1067     case kBeauty:
1068       if (sname=="FLAT"){
1069         func= PtBeautyFlat;
1070         break;
1071       }
1072       if (sname=="MUON"){
1073         func= PtBeautyMUON;
1074         break;
1075       }
1076       if (sname=="GSI"){
1077         func= PtBeautyGSI;
1078         break;
1079       }
1080       func=0;
1081       printf("<AliGenGSIlib::GetPt> unknown parametrisation\n");
1082       break;
1083
1084
1085     case kEta:
1086       func=PtEtaPHOS;
1087       break;
1088
1089     case kEtaprime:
1090       func=PtEtaprimePHOS;
1091       break;
1092
1093     case kOmega:
1094       func=PtOmega;
1095       break;
1096
1097     case kRho:
1098       func=PtRho;
1099       break;
1100
1101     case kKaon:
1102       func=PtKaonPHOS;
1103       break;
1104
1105     case kPion:
1106       func=PtPion;
1107       break;
1108
1109     case kPhi:
1110       func=PtPhiPHOS;
1111       break;
1112
1113          //    case kLambda:
1114          //         func=PtLambda;
1115          //         break;
1116
1117     case kBaryons:
1118       func=PtBaryons;
1119       break;
1120
1121     default:
1122       func=0;
1123       printf("<AliGenGSIlib::GetPt> unknown parametrisation\n");
1124     }
1125    return func;
1126 }
1127
1128 GenFunc AliGenGSIlib::GetY(Int_t param, const char * tname) const
1129 {
1130 // Return pointer to y- parameterisation
1131    GenFunc func=0;
1132    TString sname(tname);
1133
1134    switch (param)
1135     {
1136     case kUpsilon:
1137       if (sname=="FLAT"){
1138         func= YUpsilonFlat;
1139         break;
1140       }
1141       if (sname=="MUON"){
1142         func= YUpsilonMUON;
1143         break;
1144       }
1145       if (sname=="RITMAN"){
1146         func=YUpsilonRitman;
1147         break;
1148       }
1149       if (sname=="KAREL"){
1150         func=YUpsilonKarel;
1151         break;
1152       }
1153       func=0;
1154       printf("<AliGenGSIlib::GetY> unknown parametrisation\n");
1155       break;
1156
1157     case kJPsi:
1158       if (sname=="FLAT"){
1159         func= YJpsiFlat;
1160         break;
1161       }
1162       if (sname=="MUON"){
1163         func= YJpsiMUON;
1164         break;
1165       }
1166       //      if (sname=="SERGEI"){
1167       //        func= YJpsi;
1168       //        break;
1169       //      }
1170
1171       func=0;
1172       printf("<AliGenGSIlib::GetY> unknown parametrisation\n");
1173       break;
1174
1175     case kCharm:
1176         func= YCharm;
1177         break;
1178
1179     case kBeauty:
1180         func= YBeauty;
1181         break;
1182
1183     case kEta:
1184          func=YEtaPHOS;
1185          break;
1186
1187     case kEtaprime:
1188          func=YEtaprimePHOS;
1189          break;
1190
1191     case kOmega:
1192          func=YOmega;
1193          break;
1194
1195     case kRho:
1196          func=YRho;
1197          break;
1198
1199     case kKaon:
1200          func=YKaonPHOS;
1201          break;
1202
1203     case kPion:
1204          func=YPion;
1205          break;
1206
1207     case kPhi:
1208          func=YPhiPHOS;
1209          break;
1210
1211          //    case kLambda:
1212          //         func=YLambda;
1213          //         break;
1214
1215     case kBaryons:
1216          func=YBaryons;
1217          break;
1218
1219     default:
1220         func=0;
1221         printf("<AliGenGSIlib::GetY> unknown parametrisation\n");
1222     }
1223     return func;
1224 }
1225
1226 GenFuncIp AliGenGSIlib::GetIp(Int_t param, const char * tname) const
1227 {
1228 // Return pointer to particle type parameterisation
1229    GenFuncIp func=0;
1230    TString sname(tname);
1231
1232    switch (param)
1233     {
1234     case kUpsilon:
1235         func= IpUpsilon;
1236         break;
1237
1238     case kJPsi:
1239         func= IpJpsi;
1240         break;
1241
1242     case kCharm:
1243         func= IpCharm;
1244         break;
1245
1246     case kBeauty:
1247         func= IpBeauty;
1248         break;
1249
1250     case kEta:
1251          func=IpEta;
1252          break;
1253
1254     case kEtaprime:
1255          func=IpEtaprime;
1256          break;
1257
1258     case kOmega:
1259          func=IpOmega;
1260          break;
1261
1262     case kRho:
1263          func=IpRho;
1264          break;
1265
1266     case kKaon:
1267          func=IpKaonPHOS;
1268          break;
1269
1270     case kPion:
1271          func=IpPionPHOS;
1272          break;
1273
1274     case kPhi:
1275          func=IpPhi;
1276          break;
1277
1278          //    case kLambda:
1279          //         func=IpLambda;
1280          //         break;
1281
1282     case kBaryons:
1283          func=IpBaryons;
1284          break;
1285
1286     default:
1287         func=0;
1288         printf("<AliGenGSIlib::GetIp> unknown parametrisation\n");
1289     }
1290     return func;
1291 }
1292
1293
1294
1295
1296
1297
1298
1299
1300
1301
1302
1303