1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
21 #include "AliPythia.h"
26 # define lu1ent lu1ent_
29 # define lu1ent LU1ENT
30 # define type_of_call _stdcall
33 extern "C" void type_of_call
34 lu1ent(Int_t&, Int_t&, Float_t&, Float_t&, Float_t&);
38 //_____________________________________________________________________________
40 Int_t AliPythia::fgInit=0;
42 AliPythia::AliPythia()
44 for (Int_t i=0; i<501; i++) {
50 void AliPythia::Lu1Ent(Int_t flag, Int_t idpart,
51 Float_t mom, Float_t theta,Float_t phi)
53 printf("%d %d %f %f %f\n",flag, idpart, mom, theta, phi);
54 lu1ent(flag, idpart, mom, theta, phi);
57 void AliPythia::DecayParticle(Int_t idpart,
58 Float_t mom, Float_t theta,Float_t phi)
60 Lu1Ent(0, idpart, mom, theta, phi);
64 void AliPythia::ProcInit(Process_t process, Float_t energy, StrucFunc_t strucfunc)
68 fStrucFunc = strucfunc;
70 SetMDCY(LuComp(111),1,0);
71 // select structure function
73 SetMSTP(51,strucfunc);
75 // Pythia initialisation for selected processes//
79 for (Int_t i=1; i<= 200; i++) {
82 // select charm production
126 case beauty_unforced:
136 // Minimum Bias pp-Collisions
138 // Tuning of parameters descibed in G. Ciapetti and A. Di Ciaccio
139 // Proc. of the LHC Workshop, Aachen 1990, Vol. II p. 155
141 // select Pythia min. bias model
147 // Multiple interactions switched on
150 // Low-pT cut-off for hard scattering
152 // model for subsequent non-hardest interaction
153 // 90% gg->gg 10% gg->qq
155 // 90% of gluon interactions have minimum string length
160 Initialize("CMS","p","p",fEcms);
163 Int_t AliPythia::CountProducts(Int_t channel, Int_t particle)
166 for (Int_t i=1; i<=5; i++) {
167 if (TMath::Abs(GetKFDP(channel,i)) == particle) np++;
172 void AliPythia::AllowAllDecays()
175 for (i=1; i<= 2000; i++) {
179 for (i=0; i<501; i++){
184 void AliPythia::ForceParticleDecay(Int_t particle, Int_t product, Int_t mult)
187 // force decay of particle into products with multiplicity mult
189 Int_t kc=LuComp(particle);
191 Int_t ifirst=GetMDCY(kc,2);
192 Int_t ilast=ifirst+GetMDCY(kc,3)-1;
195 // Loop over decay channels
196 for (Int_t channel=ifirst; channel<=ilast;channel++) {
197 if (CountProducts(channel,product) >= mult) {
198 SetMDME(channel,1,1);
200 SetMDME(channel,1,0);
201 fBraPart[kc]-=GetBRAT(channel);
206 void AliPythia::ForceDecay(Decay_t decay)
218 if (fProcess==charm || fProcess == charm_unforced) {
219 ForceParticleDecay( 411,13,1); // D+/-
220 ForceParticleDecay( 421,13,1); // D0
221 ForceParticleDecay( 431,13,1); // D_s
222 ForceParticleDecay( 4122,13,1); // Lambda_c
224 if (fProcess==beauty || fProcess == beauty_unforced) {
225 ForceParticleDecay( 511,13,1); // B0
226 ForceParticleDecay( 521,13,1); // B+/-
227 ForceParticleDecay( 531,13,1); // B_s
228 ForceParticleDecay( 5122,13,1); // Lambda_b
232 ForceParticleDecay( 41,13,2); // phi
233 ForceParticleDecay( 443,13,2); // J/Psi
234 ForceParticleDecay(30443,13,2); // Psi'
235 ForceParticleDecay( 553,13,2); // Upsilon
236 ForceParticleDecay(30553,13,2); // Upsilon'
240 ForceParticleDecay( 411,11,1); // D+/-
241 ForceParticleDecay( 421,11,1); // D0
242 ForceParticleDecay( 431,11,1); // D_s
243 ForceParticleDecay( 4122,11,1); // Lambda_c
245 ForceParticleDecay( 511,11,1); // B0
246 ForceParticleDecay( 521,11,1); // B+/-
247 ForceParticleDecay( 531,11,1); // B_s
248 ForceParticleDecay( 5122,11,1); // Lambda_b
252 ForceParticleDecay( 41,11,2); // phi
253 ForceParticleDecay( 443,11,2); // J/Psi
254 ForceParticleDecay(30443,11,2); // Psi'
255 ForceParticleDecay( 553,11,2); // Upsilon
256 ForceParticleDecay(30553,11,2); // Upsilon'
259 ForceParticleDecay( 511,443,1); // B0
260 ForceParticleDecay( 521,443,1); // B+/-
261 ForceParticleDecay( 531,443,1); // B_s
262 ForceParticleDecay( 5122,443,1); // Lambda_b
263 ForceParticleDecay( 443,13,2); // J/Psi
266 ForceParticleDecay( 511,30443,1); // B0
267 ForceParticleDecay( 521,30443,1); // B+/-
268 ForceParticleDecay( 531,30443,1); // B_s
269 ForceParticleDecay( 5122,30443,1); // Lambda_b
270 ForceParticleDecay(30443,13,2); // Psi'
272 case b_jpsi_dielectron:
273 ForceParticleDecay( 511,443,1); // B0
274 ForceParticleDecay( 521,443,1); // B+/-
275 ForceParticleDecay( 531,443,1); // B_s
276 ForceParticleDecay( 5122,443,1); // Lambda_b
277 ForceParticleDecay( 443,11,2); // J/Psi
279 case b_psip_dielectron:
280 ForceParticleDecay( 511,30443,1); // B0
281 ForceParticleDecay( 521,30443,1); // B+/-
282 ForceParticleDecay( 531,30443,1); // B_s
283 ForceParticleDecay( 5122,30443,1); // Lambda_b
284 ForceParticleDecay(30443,11,2); // Psi'
287 ForceParticleDecay(211,13,1); // pi->mu
290 ForceParticleDecay(321,13,1); // K->mu
296 void AliPythia::DefineParticles()
306 // Some particles cloned for rare decays
308 // phi-> mu+mu- and phi -> e+e-
309 // clone the original phi
314 SetKCHG(nkc,i,GetKCHG(kc,i));
318 SetPMAS(nkc,i,GetPMAS(kc,i));
320 SetCHAF(nkc,GetCHAF(kc));
331 SetBRAT(993,2.99e-4);
350 // phi clone for dilepton decay-channel
354 gMC->Gspart(113,"Phi",3,mass,0,tlife);
361 gMC->Gspart(114,"J/Psi",3,mass,0,tlife);
368 gMC->Gspart(115,"Psi'",3,mass,0,tlife);
375 gMC->Gspart(116,"Upsilon",3,mass,0,tlife);
382 gMC->Gspart(117,"Upsilon'",3,mass,0,tlife);
389 gMC->Gspart(118,"Upsilon''",3,mass,0,tlife);
399 gMC->Gspart(119,"D^+",3,mass, 1,tlife);
400 gMC->Gspart(120,"D^-",3,mass,-1,tlife);
407 gMC->Gspart(121,"D^0",3,mass,0,tlife);
408 gMC->Gspart(122,"D^0bar",3,mass,0,tlife);
415 gMC->Gspart(123,"D_s^+",3,mass, 1,tlife);
416 gMC->Gspart(124,"D_s^-",3,mass,-1,tlife);
423 gMC->Gspart(125,"Lambda_c+",3,mass, 1,tlife);
424 gMC->Gspart(126,"Lambda_c-",3,mass,-1,tlife);
433 gMC->Gspart(127,"B^0",3,mass, 0,tlife);
434 gMC->Gspart(128,"B^0bar",3,mass, 0,tlife);
441 gMC->Gspart(129,"B^+",3,mass, 1,tlife);
442 gMC->Gspart(130,"B^-",3,mass,-1,tlife);
449 gMC->Gspart(131,"B_s",3,mass, 0,tlife);
450 gMC->Gspart(132,"B_s^bar",3,mass,0,tlife);
457 gMC->Gspart(133,"Lambda_b",3,mass, 0,tlife);
458 gMC->Gspart(134,"Lambda_b^bar",3,mass,0,tlife);
462 // set up correspondance between standard GEANT particle codes
465 kc=LuComp(22); // gamma
469 kc=LuComp(11); // positron
473 kc=LuComp(12); // neutrino
477 kc=LuComp(13); // muon
481 kc=LuComp(111); // pi0
485 kc=LuComp(211); // pi+
489 kc=LuComp(130); // K0 short
493 kc=LuComp(321); // K+/-
497 kc=LuComp(2112); // neutron/anti-neutron
501 kc=LuComp(2212); // proton/anti-proton
505 kc=LuComp(310); // K0 short
509 kc=LuComp(221); // eta
513 kc=LuComp(3122); // lambda
517 kc=LuComp(3222); // sigma+/antisigma+
521 kc=LuComp(3212); // sigma0/antisigma0
525 kc=LuComp(3112); // sigma-/antisigma-
529 kc=LuComp(3322); // xsi0-/antixsi0
533 kc=LuComp(3312); // xsi-/antixsi+
537 kc=LuComp(3334); // omega/antiomega
544 Int_t AliPythia::GetGeantCode(Int_t kf)
546 Int_t kc=LuComp(TMath::Abs(kf));
547 return (kf > 0) ? fGPCode[kc][0] : fGPCode[kc][1];
550 Float_t AliPythia::GetBraPart(Int_t kf)
552 Int_t kc=LuComp(TMath::Abs(kf));