/************************************************************************** * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * * * * Author: The ALICE Off-line Project. * * Contributors are mentioned in the code where appropriate. * * * * Permission to use, copy, modify and distribute this software and its * * documentation strictly for non-commercial purposes is hereby granted * * without fee, provided that the above copyright notice appears in all * * copies and that both the copyright notice and this permission notice * * appear in the supporting documentation. The authors make no claims * * about the suitability of this software for any purpose. It is * * provided "as is" without express or implied warranty. * **************************************************************************/ //====================================================================== // AliGenSTRANGElib class contains parameterizations of the // kaon, phi and hyperon (Lambda, Anti-Lambda, Xi, anti-Xi, Omega, // anti-Omega) for the PPR study of the strange particle production. // These parameterizations are used by the // AliGenParam class: // AliGenParam(npar, param, AliGenSTRANGElib::GetPt(param), // AliGenSTRANGElib::GetY(param), // AliGenSTRANGElib::GetIp(param) ) // param represents the particle to be simulated. // ????????? // Pt distributions are calculated from the transverse mass scaling // with Pions, using the PtScal function taken from AliGenMUONlib // version aliroot 3.01 // // Rocco CALIANDRO. Rosa Anna FINI, Tiziano VIRGILI // Rocco.Caliandro@cern.ch Rosanna.Fini@ba.infn.it, // Tiziano.Virgili@roma1.infn.it //====================================================================== /* $Id$ */ #include #include #include "AliGenSTRANGElib.h" ClassImp(AliGenSTRANGElib) //============================================================= // Double_t AliGenSTRANGElib::PtScal(Double_t pt, Int_t np) { // Mt-scaling // Function for the calculation of the Pt distribution for a // given particle np, from the pion Pt distribution using the // mt scaling. This function was taken from AliGenMUONlib // aliroot version 3.01, and was extended for hyperons. // np = 1=>Pions 2=>Kaons 3=>Etas 4=>Omegas 5=>ETA' 6=>PHI // 7=>BARYONS-BARYONBARS // 8=>Lambda-antiLambda // 9=>Xi-antiXi // 10=>Omega-antiOmega // MASS SCALING RESPECT TO PIONS // MASS 0=>PI, 1=>K, 2=>ETA,3=>OMEGA,4=>ETA',5=>PHI const Double_t khm[11] = {0.1396, 0.494,0.547, 0.782, 0.957, 1.02, // MASS 6=>BARYON-BARYONBAR 0.938, // MASS 7=>Lambda-antiLambda 1.1157, // MASS 8=>Xi-antiXi 1.3213, // MASS 9=>Omega-antiOmega 1.6725, // MASS 10=>Lambda(1520) 1.5195}; // VALUE MESON/PI AT 5 GEV const Double_t kfmax[11]={1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1.}; Double_t f5=TMath::Power(((sqrt(100.018215)+2.)/(sqrt(100.+khm[np]*khm[np])+2.0)),12.3); Double_t kfmax2=f5/kfmax[np]; // PIONS Double_t ptpion=100.*PtPion(&pt, (Double_t*) 0); Double_t fmtscal=TMath::Power(((sqrt(pt*pt+0.018215)+2.)/ (sqrt(pt*pt+khm[np]*khm[np])+2.0)),12.3)/ kfmax2; return fmtscal*ptpion; } //============================================================= // Double_t AliGenSTRANGElib::PtPion(const Double_t *px, const Double_t *) { // Pion transverse momentum distribtuion taken // from AliGenMUONlib class, version 3.01 of aliroot // PT-PARAMETERIZATION CDF, PRL 61(88) 1819 // POWER LAW FOR PT > 500 MEV // MT SCALING BELOW (T=160 MEV) // const Double_t kp0 = 1.3; const Double_t kxn = 8.28; const Double_t kxlim=0.5; const Double_t kt=0.160; const Double_t kxmpi=0.139; const Double_t kb=1.; Double_t y, y1, kxmpi2, ynorm, a; Double_t x=*px; // y1=TMath::Power(kp0/(kp0+kxlim),kxn); kxmpi2=kxmpi*kxmpi; ynorm=kb*(TMath::Exp(-sqrt(kxlim*kxlim+kxmpi2)/kt)); a=ynorm/y1; if (x > kxlim) y=a*TMath::Power(kp0/(kp0+x),kxn); else y=kb*TMath::Exp(-sqrt(x*x+kxmpi2)/kt); return y*x; } // End Scaling //============================================================================ // K A O N Double_t AliGenSTRANGElib::PtKaon( const Double_t *px, const Double_t *) { // kaon // pt-distribution //____________________________________________________________ return PtScal(*px,1); // 1==> Kaon in the PtScal function } Double_t AliGenSTRANGElib::YKaon( const Double_t *py, const Double_t *) { // y-distribution //____________________________________________________________ const Double_t ka = 1000.; const Double_t kdy = 4.*4; Double_t y=TMath::Abs(*py); // Double_t ex = y*y/(2*kdy*kdy); return ka*TMath::Exp(-ex); } Int_t AliGenSTRANGElib::IpKaon(TRandom *ran) { // particle composition // Float_t random = ran->Rndm(); Float_t random2 = ran->Rndm(); if (random2 < 0.5) { if (random < 0.5) { return 321; // K+ } else { return -321; // K- } } else { if (random < 0.5) { return 130; // K^0 short } else { return 310; // K^0 long } } } // End Kaons //============================================================================ //============================================================================ // P H I Double_t AliGenSTRANGElib::PtPhi( const Double_t *px, const Double_t *) { // phi // pt-distribution //____________________________________________________________ return PtScal(*px,5); // 5==> Phi in the PtScal function } Double_t AliGenSTRANGElib::YPhi( const Double_t *py, const Double_t *) { // y-distribution //____________________________________________________________ const Double_t ka = 1000.; const Double_t kdy = 4.*4; Double_t y=TMath::Abs(*py); // Double_t ex = y*y/(2*kdy*kdy); return ka*TMath::Exp(-ex); } Int_t AliGenSTRANGElib::IpPhi(TRandom *) { // particle composition // return 333; // Phi } // End Phis //=================================================================== //============================================================================ // Lambda Double_t AliGenSTRANGElib::PtLambda( const Double_t *px, const Double_t *) { // Lambda // pt-distribution //____________________________________________________________ return PtScal(*px,7); // 7==> Lambda-antiLambda in the PtScal function } Double_t AliGenSTRANGElib::YLambda( const Double_t *py, const Double_t *) { // y-distribution //____________________________________________________________ const Double_t ka = 1000.; const Double_t kdy = 4.*4; Double_t y=TMath::Abs(*py); // Double_t ex = y*y/(2*kdy*kdy); return ka*TMath::Exp(-ex); } Int_t AliGenSTRANGElib::IpLambda(TRandom *ran) { // particle composition // generation of fixed type of particle // Float_t random = ran->Rndm(); if (random < 0.5) { return 3122; // Lambda } else { return -3122; // Anti-Lambda } } // End Lambda //============================================================================ // XIminus Double_t AliGenSTRANGElib::PtXiMinus( const Double_t *px, const Double_t *) { // Xi // pt-distribution //____________________________________________________________ return PtScal(*px,8); // 8==> Xi-antiXi in the PtScal function } Double_t AliGenSTRANGElib::YXiMinus( const Double_t *py, const Double_t *) { // y-distribution //____________________________________________________________ const Double_t ka = 1000.; const Double_t kdy = 4.*4; Double_t y=TMath::Abs(*py); // Double_t ex = y*y/(2*kdy*kdy); return ka*TMath::Exp(-ex); } Int_t AliGenSTRANGElib::IpXiMinus(TRandom *ran) { // particle composition // generation of fixed type of particle // Float_t random = ran->Rndm(); if (random < 0.5) { return 3312; // Xi- } else { return -3312; // Xi+ } } // End Ximinus //============================================================================ // Omegaminus Double_t AliGenSTRANGElib::PtOmegaMinus( const Double_t *px, const Double_t *) { // Omega // pt-distribution //____________________________________________________________ return PtScal(*px,9); // 9==> Omega-antiOmega in the PtScal function } Double_t AliGenSTRANGElib::YOmegaMinus( const Double_t *py, const Double_t *) { // y-distribution //____________________________________________________________ const Double_t ka = 1000.; const Double_t kdy = 4.*4; Double_t y=TMath::Abs(*py); // Double_t ex = y*y/(2*kdy*kdy); return ka*TMath::Exp(-ex); } Int_t AliGenSTRANGElib::IpOmegaMinus(TRandom * ran) { // particle composition // generation of fixed type of particle // Float_t random = ran->Rndm(); if (random < 0.5) { return 3334; // Omega- } else { return -3334; // Omega+ } } // End Omegaminus //============================================================================ // Lambda(1520) Double_t AliGenSTRANGElib::PtLambda1520( const Double_t *px, const Double_t *) { // Lambda(1520) // pt-distribution //____________________________________________________________ return PtScal(*px,10); // 10=> Lambda(1520) in the PtScal function } Double_t AliGenSTRANGElib::YLambda1520( const Double_t *py, const Double_t *) { // y-distribution //____________________________________________________________ const Double_t ka = 1000.; const Double_t kdy = 4.; Double_t y=TMath::Abs(*py); // Double_t ex = y*y/(2*kdy*kdy); return ka*TMath::Exp(-ex); } Int_t AliGenSTRANGElib::IpLambda1520(TRandom * ran) { // particle composition // generation of fixed type of particle // Float_t random = ran->Rndm(); if (random < 0.5) { return 3124; // Lambda(1520) } else { return -3124; // antiLambda(1520) } } // End Lambda(1520) //============================================================================ typedef Double_t (*GenFunc) (const Double_t*, const Double_t*); GenFunc AliGenSTRANGElib::GetPt(Int_t param, const char* /*tname*/) const { // Return pinter to pT parameterisation GenFunc func; switch (param) { case kKaon: func=PtKaon; break; case kPhi: func=PtPhi; break; case kLambda: func=PtLambda; break; case kXiMinus: func=PtXiMinus; break; case kOmegaMinus: func=PtOmegaMinus; break; case kLambda1520: func=PtLambda1520; break; default: func=0; printf(" unknown parametrisationn"); } return func; } GenFunc AliGenSTRANGElib::GetY(Int_t param, const char* /*tname*/) const { // Return pointer to Y parameterisation GenFunc func; switch (param) { case kKaon: func=YKaon; break; case kPhi: func=YPhi; break; case kLambda: func=YLambda; break; case kXiMinus: func=YXiMinus; break; case kOmegaMinus: func=YOmegaMinus; break; case kLambda1520: func=YLambda1520; break; default: func=0; printf(" unknown parametrisationn"); } return func; } typedef Int_t (*GenFuncIp) (TRandom *); GenFuncIp AliGenSTRANGElib::GetIp(Int_t param, const char* /*tname*/) const { // Return pointer to particle composition GenFuncIp func; switch (param) { case kKaon: func=IpKaon; break; case kPhi: func=IpPhi; break; case kLambda: func=IpLambda; break; case kXiMinus: func=IpXiMinus; break; case kOmegaMinus: func=IpOmegaMinus; break; case kLambda1520: func=IpLambda1520; break; default: func=0; printf(" unknown parametrisationn"); } return func; }