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 **************************************************************************/
18 // Implementation of AliDecayer using Pythia8
19 // Author: andreas.morsch@cern.ch
22 #include <TLorentzVector.h>
23 #include "AliTPythia8.h"
24 #include "AliDecayerPythia8.h"
25 #include "ParticleData.h"
27 ClassImp(AliDecayerPythia8)
29 Bool_t AliDecayerPythia8::fgInit = kFALSE;
31 AliDecayerPythia8::AliDecayerPythia8():
33 fPythia8(new AliTPythia8()),
39 fPythia8->Pythia8()->readString("SoftQCD:elastic = on");
40 fPythia8->Pythia8()->init();
43 //___________________________________________________________________________
44 void AliDecayerPythia8::Decay(Int_t pdg, TLorentzVector* p)
46 // Decay a single particle
48 AppendParticle(pdg, p);
49 Int_t idPart = fPythia8->Pythia8()->event[0].id();
50 fPythia8->Pythia8()->particleData.mayDecay(idPart,kTRUE);
51 fPythia8->Pythia8()->moreDecays();
52 if (fDebug > 0) fPythia8->EventListing();
55 //___________________________________________________________________________
56 Int_t AliDecayerPythia8::ImportParticles(TClonesArray *particles)
58 //import the decay products into particles array
59 return (fPythia8->ImportParticles(particles, "All"));
63 void AliDecayerPythia8::Init()
69 // fPythia->SetDecayTable();
72 // Switch on heavy flavor decays
75 Int_t heavy[14] = {411, 421, 431, 4122, 4132, 4232, 4332, 511, 521, 531, 5122, 5132, 5232, 5332};
76 // fPythia->ResetDecayTable();
77 for (j=0; j < 14; j++) {
78 if (fDecay == kNoDecayHeavy) {
79 AliTPythia8::Instance()->ReadString(Form("%d:onMode = off", heavy[j]));
81 AliTPythia8::Instance()->ReadString(Form("%d:onMode = on", heavy[j]));
86 //...Switch off decay of pi0, K0S, Lambda, Sigma+-, Xi0-, Omega-.
88 if (fDecay != kNeutralPion) {
89 AliTPythia8::Instance()->ReadString("111:onMode = off");
91 AliTPythia8::Instance()->ReadString("111:onMode = on");
94 AliTPythia8::Instance()->ReadString("310:onMode = off");
95 AliTPythia8::Instance()->ReadString("3112:onMode = off");
96 AliTPythia8::Instance()->ReadString("3212:onMode = off");
97 AliTPythia8::Instance()->ReadString("3222:onMode = off");
98 AliTPythia8::Instance()->ReadString("3312:onMode = off");
99 AliTPythia8::Instance()->ReadString("3322:onMode = off");
100 AliTPythia8::Instance()->ReadString("3334:onMode = off");
101 // .. Force decay channels
105 void AliDecayerPythia8::ForceDecay()
108 // Force a particle decay mode
109 // Switch heavy flavour production off if requested
110 if (!fHeavyFlavour) SwitchOffHeavyFlavour();
112 Decay_t decay = fDecay;
113 AliTPythia8::Instance()->ReadString("HadronLevel:Decay = on");
115 if (decay == kNoDecayHeavy) return;
123 AliTPythia8::Instance()->ReadString("511:onMode = off");
124 AliTPythia8::Instance()->ReadString("511:onIfAny = 13 443 100443");
126 AliTPythia8::Instance()->ReadString("521:onMode = off");
127 AliTPythia8::Instance()->ReadString("521:onIfAny = 13 443 100443");
129 AliTPythia8::Instance()->ReadString("531:onMode = off");
130 AliTPythia8::Instance()->ReadString("531:onIfAny = 13 443 100443");
132 AliTPythia8::Instance()->ReadString("5122:onMode = off");
133 AliTPythia8::Instance()->ReadString("5122:onIfAny = 13 443 100443");
135 AliTPythia8::Instance()->ReadString("5132:onMode = off");
136 AliTPythia8::Instance()->ReadString("5132:onIfAny = 13 443 100443");
138 AliTPythia8::Instance()->ReadString("5232:onMode = off");
139 AliTPythia8::Instance()->ReadString("5232:onIfAny = 13 443 100443");
141 AliTPythia8::Instance()->ReadString("5332:onMode = off");
142 AliTPythia8::Instance()->ReadString("5332:onIfAny = 13 443 100443");
144 AliTPythia8::Instance()->ReadString("100443:onMode = off");
145 AliTPythia8::Instance()->ReadString("100443:onIfAny = 443");
147 AliTPythia8::Instance()->ReadString("443:onMode = off");
148 AliTPythia8::Instance()->ReadString("443:onIfAll = 13 13");
150 AliTPythia8::Instance()->ReadString("411:onMode = off");
151 AliTPythia8::Instance()->ReadString("411:onIfAll = 13");
153 AliTPythia8::Instance()->ReadString("421:onMode = off");
154 AliTPythia8::Instance()->ReadString("421:onIfAll = 13");
156 AliTPythia8::Instance()->ReadString("431:onMode = off");
157 AliTPythia8::Instance()->ReadString("431:onIfAll = 13");
159 AliTPythia8::Instance()->ReadString("4122:onMode = off");
160 AliTPythia8::Instance()->ReadString("4122:onIfAll = 13");
162 AliTPythia8::Instance()->ReadString("4132:onMode = off");
163 AliTPythia8::Instance()->ReadString("4132:onIfAll = 13");
165 AliTPythia8::Instance()->ReadString("4232:onMode = off");
166 AliTPythia8::Instance()->ReadString("4232:onIfAll = 13");
168 AliTPythia8::Instance()->ReadString("4332:onMode = off");
169 AliTPythia8::Instance()->ReadString("4332:onIfAll = 13");
172 case kChiToJpsiGammaToMuonMuon:
173 // Chi_1c -> J/Psi Gamma
174 AliTPythia8::Instance()->ReadString("20443:onMode = off");
175 AliTPythia8::Instance()->ReadString("20443:onIfAll = 443 22");
176 // Chi_2c -> J/Psi Gamma
177 AliTPythia8::Instance()->ReadString("445:onMode = off");
178 AliTPythia8::Instance()->ReadString("445:onIfAll = 443 22");
180 AliTPythia8::Instance()->ReadString("443:onMode = off");
181 AliTPythia8::Instance()->ReadString("443:onIfAll = 13 13");
183 case kChiToJpsiGammaToElectronElectron:
184 // Chi_1c -> J/Psi Gamma
185 AliTPythia8::Instance()->ReadString("20443:onMode = off");
186 AliTPythia8::Instance()->ReadString("20443:onIfAll = 443 22");
187 // Chi_2c -> J/Psi Gamma
188 AliTPythia8::Instance()->ReadString("445:onMode = off");
189 AliTPythia8::Instance()->ReadString("445:onIfAll = 443 22");
191 AliTPythia8::Instance()->ReadString("443:onMode = off");
192 AliTPythia8::Instance()->ReadString("443:onIfAll = 11 11");
197 AliTPythia8::Instance()->ReadString("511:onMode = off");
198 AliTPythia8::Instance()->ReadString("511:onIfAny = 13");
200 AliTPythia8::Instance()->ReadString("521:onMode = off");
201 AliTPythia8::Instance()->ReadString("521:onIfAny = 13");
203 AliTPythia8::Instance()->ReadString("531:onMode = off");
204 AliTPythia8::Instance()->ReadString("531:onIfAny = 13");
206 AliTPythia8::Instance()->ReadString("5122:onMode = off");
207 AliTPythia8::Instance()->ReadString("5122:onIfAny = 13");
209 AliTPythia8::Instance()->ReadString("5132:onMode = off");
210 AliTPythia8::Instance()->ReadString("5132:onIfAny = 13");
212 AliTPythia8::Instance()->ReadString("5232:onMode = off");
213 AliTPythia8::Instance()->ReadString("5232:onIfAny = 13");
215 AliTPythia8::Instance()->ReadString("5332:onMode = off");
216 AliTPythia8::Instance()->ReadString("5332:onIfAny = 13");
220 AliTPythia8::Instance()->ReadString("411:onMode = off");
221 AliTPythia8::Instance()->ReadString("411:onIfAll = 13");
223 AliTPythia8::Instance()->ReadString("421:onMode = off");
224 AliTPythia8::Instance()->ReadString("421:onIfAll = 13");
226 AliTPythia8::Instance()->ReadString("431:onMode = off");
227 AliTPythia8::Instance()->ReadString("431:onIfAll = 13");
229 AliTPythia8::Instance()->ReadString("4122:onMode = off");
230 AliTPythia8::Instance()->ReadString("4122:onIfAll = 13");
232 AliTPythia8::Instance()->ReadString("4132:onMode = off");
233 AliTPythia8::Instance()->ReadString("4132:onIfAll = 13");
235 AliTPythia8::Instance()->ReadString("4232:onMode = off");
236 AliTPythia8::Instance()->ReadString("4232:onIfAll = 13");
238 AliTPythia8::Instance()->ReadString("4332:onMode = off");
239 AliTPythia8::Instance()->ReadString("4332:onIfAll = 13");
241 AliTPythia8::Instance()->ReadString("511:onMode = off");
242 AliTPythia8::Instance()->ReadString("511:onIfAny = 13");
244 AliTPythia8::Instance()->ReadString("521:onMode = off");
245 AliTPythia8::Instance()->ReadString("521:onIfAny = 13");
247 AliTPythia8::Instance()->ReadString("531:onMode = off");
248 AliTPythia8::Instance()->ReadString("531:onIfAny = 13");
250 AliTPythia8::Instance()->ReadString("5122:onMode = off");
251 AliTPythia8::Instance()->ReadString("5122:onIfAny = 13");
253 AliTPythia8::Instance()->ReadString("5132:onMode = off");
254 AliTPythia8::Instance()->ReadString("5132:onIfAny = 13");
256 AliTPythia8::Instance()->ReadString("5232:onMode = off");
257 AliTPythia8::Instance()->ReadString("5232:onIfAny = 13");
259 AliTPythia8::Instance()->ReadString("5332:onMode = off");
260 AliTPythia8::Instance()->ReadString("5332:onIfAny = 13");
265 AliTPythia8::Instance()->ReadString("113:onMode = off");
266 AliTPythia8::Instance()->ReadString("113:onIfAll = 13 13");
268 AliTPythia8::Instance()->ReadString("221:onMode = off");
269 AliTPythia8::Instance()->ReadString("221:onIfAll = 13 13");
271 AliTPythia8::Instance()->ReadString("223:onMode = off");
272 AliTPythia8::Instance()->ReadString("223:onIfAll = 13 13");
274 AliTPythia8::Instance()->ReadString("333:onMode = off");
275 AliTPythia8::Instance()->ReadString("333:onIfAll = 13 13");
277 AliTPythia8::Instance()->ReadString("443:onMode = off");
278 AliTPythia8::Instance()->ReadString("443:onIfAll = 13 13");
280 AliTPythia8::Instance()->ReadString("100443:onMode = off");
281 AliTPythia8::Instance()->ReadString("100443:onIfAll = 13 13");
283 AliTPythia8::Instance()->ReadString("553:onMode = off");
284 AliTPythia8::Instance()->ReadString("553:onIfAll = 13 13");
286 AliTPythia8::Instance()->ReadString("100553:onMode = off");
287 AliTPythia8::Instance()->ReadString("100553:onIfAll = 13 13");
289 AliTPythia8::Instance()->ReadString("200553:onMode = off");
290 AliTPythia8::Instance()->ReadString("200553:onIfAll = 13 13");
292 case kBSemiElectronic:
294 AliTPythia8::Instance()->ReadString("511:onMode = off");
295 AliTPythia8::Instance()->ReadString("511:onIfAny = 11");
297 AliTPythia8::Instance()->ReadString("521:onMode = off");
298 AliTPythia8::Instance()->ReadString("521:onIfAny = 11");
300 AliTPythia8::Instance()->ReadString("531:onMode = off");
301 AliTPythia8::Instance()->ReadString("531:onIfAny = 11");
303 AliTPythia8::Instance()->ReadString("5122:onMode = off");
304 AliTPythia8::Instance()->ReadString("5122:onIfAny = 11");
306 AliTPythia8::Instance()->ReadString("5132:onMode = off");
307 AliTPythia8::Instance()->ReadString("5132:onIfAny = 11");
309 AliTPythia8::Instance()->ReadString("5232:onMode = off");
310 AliTPythia8::Instance()->ReadString("5232:onIfAny = 11");
312 AliTPythia8::Instance()->ReadString("5332:onMode = off");
313 AliTPythia8::Instance()->ReadString("5332:onIfAny = 11");
315 case kSemiElectronic:
317 AliTPythia8::Instance()->ReadString("411:onMode = off");
318 AliTPythia8::Instance()->ReadString("411:onIfAll = 11");
320 AliTPythia8::Instance()->ReadString("421:onMode = off");
321 AliTPythia8::Instance()->ReadString("421:onIfAll = 11");
323 AliTPythia8::Instance()->ReadString("431:onMode = off");
324 AliTPythia8::Instance()->ReadString("431:onIfAll = 11");
326 AliTPythia8::Instance()->ReadString("4122:onMode = off");
327 AliTPythia8::Instance()->ReadString("4122:onIfAll = 11");
329 AliTPythia8::Instance()->ReadString("4132:onMode = off");
330 AliTPythia8::Instance()->ReadString("4132:onIfAll = 11");
332 AliTPythia8::Instance()->ReadString("4232:onMode = off");
333 AliTPythia8::Instance()->ReadString("4232:onIfAll = 11");
335 AliTPythia8::Instance()->ReadString("4332:onMode = off");
336 AliTPythia8::Instance()->ReadString("4332:onIfAll = 11");
338 AliTPythia8::Instance()->ReadString("511:onMode = off");
339 AliTPythia8::Instance()->ReadString("511:onIfAny = 11");
341 AliTPythia8::Instance()->ReadString("521:onMode = off");
342 AliTPythia8::Instance()->ReadString("521:onIfAny = 11");
344 AliTPythia8::Instance()->ReadString("531:onMode = off");
345 AliTPythia8::Instance()->ReadString("531:onIfAny = 11");
347 AliTPythia8::Instance()->ReadString("5122:onMode = off");
348 AliTPythia8::Instance()->ReadString("5122:onIfAny = 11");
350 AliTPythia8::Instance()->ReadString("5132:onMode = off");
351 AliTPythia8::Instance()->ReadString("5132:onIfAny = 11");
353 AliTPythia8::Instance()->ReadString("5232:onMode = off");
354 AliTPythia8::Instance()->ReadString("5232:onIfAny = 11");
356 AliTPythia8::Instance()->ReadString("5332:onMode = off");
357 AliTPythia8::Instance()->ReadString("5332:onIfAny = 11");
361 AliTPythia8::Instance()->ReadString("113:onMode = off");
362 AliTPythia8::Instance()->ReadString("113:onIfAll = 11 11");
364 AliTPythia8::Instance()->ReadString("221:onMode = off");
365 AliTPythia8::Instance()->ReadString("221:onIfAll = 11 11");
367 AliTPythia8::Instance()->ReadString("223:onMode = off");
368 AliTPythia8::Instance()->ReadString("223:onIfAll = 11 11");
370 AliTPythia8::Instance()->ReadString("333:onMode = off");
371 AliTPythia8::Instance()->ReadString("333:onIfAll = 11 11");
373 AliTPythia8::Instance()->ReadString("443:onMode = off");
374 AliTPythia8::Instance()->ReadString("443:onIfAll = 11 11");
376 AliTPythia8::Instance()->ReadString("100443:onMode = off");
377 AliTPythia8::Instance()->ReadString("100443:onIfAll = 11 11");
379 AliTPythia8::Instance()->ReadString("553:onMode = off");
380 AliTPythia8::Instance()->ReadString("553:onIfAll = 11 11");
382 AliTPythia8::Instance()->ReadString("100553:onMode = off");
383 AliTPythia8::Instance()->ReadString("100553:onIfAll = 11 11");
385 AliTPythia8::Instance()->ReadString("200553:onMode = off");
386 AliTPythia8::Instance()->ReadString("200553:onIfAll = 11 11");
389 // B0 -> J/Psi (Psi') X
390 AliTPythia8::Instance()->ReadString("511:onMode = off");
391 AliTPythia8::Instance()->ReadString("511:onIfAny = 443 100443");
392 // B+/- -> J/Psi (Psi') X
393 AliTPythia8::Instance()->ReadString("521:onMode = off");
394 AliTPythia8::Instance()->ReadString("521:onIfAny = 443 100443");
395 // B_s -> J/Psi (Psi') X
396 AliTPythia8::Instance()->ReadString("531:onMode = off");
397 AliTPythia8::Instance()->ReadString("531:onIfAny = 443 100443");
398 // Lambda_b -> J/Psi (Psi') X
399 AliTPythia8::Instance()->ReadString("5122:onMode = off");
400 AliTPythia8::Instance()->ReadString("5122:onIfAny = 443 100443");
403 AliTPythia8::Instance()->ReadString("443:onMode = off");
404 AliTPythia8::Instance()->ReadString("443:onIfAll = 13 13");
406 AliTPythia8::Instance()->ReadString("100443:onMode = off");
407 AliTPythia8::Instance()->ReadString("100443:onIfAll = 13 13");
409 case kBPsiPrimeDiMuon:
411 AliTPythia8::Instance()->ReadString("511:onMode = off");
412 AliTPythia8::Instance()->ReadString("511:onIfAny = 100443");
414 AliTPythia8::Instance()->ReadString("521:onMode = off");
415 AliTPythia8::Instance()->ReadString("521:onIfAny = 100443");
417 AliTPythia8::Instance()->ReadString("531:onMode = off");
418 AliTPythia8::Instance()->ReadString("531:onIfAny = 100443");
419 // Lambda_b -> Psi' X
420 AliTPythia8::Instance()->ReadString("5122:onMode = off");
421 AliTPythia8::Instance()->ReadString("5122:onIfAny = 100443");
424 AliTPythia8::Instance()->ReadString("100443:onMode = off");
425 AliTPythia8::Instance()->ReadString("100443:onIfAll = 13 13");
427 case kBJpsiDiElectron:
429 AliTPythia8::Instance()->ReadString("511:onMode = off");
430 AliTPythia8::Instance()->ReadString("511:onIfAny = 443");
432 AliTPythia8::Instance()->ReadString("521:onMode = off");
433 AliTPythia8::Instance()->ReadString("521:onIfAny = 443");
435 AliTPythia8::Instance()->ReadString("531:onMode = off");
436 AliTPythia8::Instance()->ReadString("531:onIfAny = 443");
438 AliTPythia8::Instance()->ReadString("5122:onMode = off");
439 AliTPythia8::Instance()->ReadString("5122:onIfAny = 443");
442 AliTPythia8::Instance()->ReadString("443:onMode = off");
443 AliTPythia8::Instance()->ReadString("443:onIfAll = 11 11");
448 AliTPythia8::Instance()->ReadString("511:onMode = off");
449 AliTPythia8::Instance()->ReadString("511:onIfAny = 443");
451 AliTPythia8::Instance()->ReadString("521:onMode = off");
452 AliTPythia8::Instance()->ReadString("521:onIfAny = 443");
454 AliTPythia8::Instance()->ReadString("531:onMode = off");
455 AliTPythia8::Instance()->ReadString("531:onIfAny = 443");
457 AliTPythia8::Instance()->ReadString("5122:onMode = off");
458 AliTPythia8::Instance()->ReadString("5122:onIfAny = 443");
460 case kBPsiPrimeDiElectron:
462 AliTPythia8::Instance()->ReadString("511:onMode = off");
463 AliTPythia8::Instance()->ReadString("511:onIfAny = 100443");
465 AliTPythia8::Instance()->ReadString("521:onMode = off");
466 AliTPythia8::Instance()->ReadString("521:onIfAny = 100443");
468 AliTPythia8::Instance()->ReadString("531:onMode = off");
469 AliTPythia8::Instance()->ReadString("531:onIfAny = 100443");
470 // Lambda_b -> Psi' X
471 AliTPythia8::Instance()->ReadString("5122:onMode = off");
472 AliTPythia8::Instance()->ReadString("5122:onIfAny = 100443");
475 AliTPythia8::Instance()->ReadString("100443:onMode = off");
476 AliTPythia8::Instance()->ReadString("100443:onIfAll = 11 11");
480 AliTPythia8::Instance()->ReadString("211:onMode = off");
481 AliTPythia8::Instance()->ReadString("211:onIfAny = 13");
485 AliTPythia8::Instance()->ReadString("321:onMode = off");
486 AliTPythia8::Instance()->ReadString("321:onIfAny = 13");
490 AliTPythia8::Instance()->ReadString("211:onMode = off");
491 AliTPythia8::Instance()->ReadString("211:onIfAny = 13");
492 AliTPythia8::Instance()->ReadString("321:onMode = off");
493 AliTPythia8::Instance()->ReadString("321:onIfAny = 13");
497 AliTPythia8::Instance()->ReadString("24:onMode = off");
498 AliTPythia8::Instance()->ReadString("24:onIfAny = 13");
502 AliTPythia8::Instance()->ReadString("24:onMode = off");
503 AliTPythia8::Instance()->ReadString("24:onIfAny = 4");
505 case kWToCharmToMuon:
507 AliTPythia8::Instance()->ReadString("24:onMode = off");
508 AliTPythia8::Instance()->ReadString("24:onIfAny = 4");
510 AliTPythia8::Instance()->ReadString("411:onMode = off");
511 AliTPythia8::Instance()->ReadString("411:onIfAll = 13");
513 AliTPythia8::Instance()->ReadString("421:onMode = off");
514 AliTPythia8::Instance()->ReadString("421:onIfAll = 13");
516 AliTPythia8::Instance()->ReadString("431:onMode = off");
517 AliTPythia8::Instance()->ReadString("431:onIfAll = 13");
519 AliTPythia8::Instance()->ReadString("4122:onMode = off");
520 AliTPythia8::Instance()->ReadString("4122:onIfAll = 13");
522 AliTPythia8::Instance()->ReadString("4132:onMode = off");
523 AliTPythia8::Instance()->ReadString("4132:onIfAll = 13");
525 AliTPythia8::Instance()->ReadString("4232:onMode = off");
526 AliTPythia8::Instance()->ReadString("4232:onIfAll = 13");
528 AliTPythia8::Instance()->ReadString("4332:onMode = off");
529 AliTPythia8::Instance()->ReadString("4332:onIfAll = 13");
533 AliTPythia8::Instance()->ReadString("23:onMode = off");
534 AliTPythia8::Instance()->ReadString("23:onIfAll = 13 13");
538 AliTPythia8::Instance()->ReadString("23:onMode = off");
539 AliTPythia8::Instance()->ReadString("23:onIfAll = 11 11");
544 case kHadronicDWithout4Bodies:
549 AliTPythia8::Instance()->ReadString("333:onMode = off");
550 AliTPythia8::Instance()->ReadString("333:onIfAll = 321 321");
554 AliTPythia8::Instance()->ReadString("3334:onMode = off");
555 AliTPythia8::Instance()->ReadString("3334:onIfAll = 3122 321 ");
558 AliTPythia8::Instance()->ReadString("3122:onMode = off");
559 AliTPythia8::Instance()->ReadString("3122:onIfAll = 2212 211 ");
563 AliTPythia8::Instance()->ReadString("HadronLevel:Decay = off");
568 case kPsiPrimeJpsiDiElectron:
576 Float_t AliDecayerPythia8::GetPartialBranchingRatio(Int_t ipart)
578 // Get the partial branching ration for the forced decay channels
580 Pythia8::Pythia* thePythia = AliTPythia8::Instance()->Pythia8();
581 Pythia8::ParticleData & table = thePythia->particleData;
582 Pythia8::ParticleDataEntry* pd = table.particleDataEntryPtr(ipart);
584 Int_t nc = pd->sizeChannels();
587 // Loop over decay channels
588 for (Int_t ic = 0; ic < nc; ic++) {
589 Pythia8::DecayChannel& decCh = pd->channel(ic);
590 for (Int_t i = 0; i < decCh.multiplicity(); i++) {
591 br += decCh.bRatio();
598 Float_t AliDecayerPythia8::GetLifetime(Int_t kf)
600 // Return lifetime of particle
601 Pythia8::Pythia* thePythia = AliTPythia8::Instance()->Pythia8();
602 Pythia8::ParticleData& table = thePythia->particleData;
603 Float_t tau = table.tau0(kf);
607 void AliDecayerPythia8::SwitchOffHeavyFlavour()
609 // Switch off heavy flavour production
611 // Maximum number of quark flavours used in pdf
612 AliTPythia8::Instance()->ReadString("PDFinProcess:nQuarkIn = 3");
613 // Maximum number of flavors that can be used in showers
614 AliTPythia8::Instance()->ReadString("SpaceShower:nQuarkIn = 3");
615 AliTPythia8::Instance()->ReadString("TimeShower:nGammaToQuark = 3");
616 AliTPythia8::Instance()->ReadString("TimeShower:nGluonToQuark = 3");
620 void AliDecayerPythia8::ForceHadronicD(Int_t optUse4Bodies)
623 // Force golden D decay modes
626 AliTPythia8::Instance()->ReadString("313:onMode = off");
627 AliTPythia8::Instance()->ReadString("313:onIfAll = 321 211");
629 AliTPythia8::Instance()->ReadString("333:onMode = off");
630 AliTPythia8::Instance()->ReadString("333:onIfAll = 321 321");
631 // for D0 -> rho0 pi+ k-
632 AliTPythia8::Instance()->ReadString("113:onMode = off");
633 AliTPythia8::Instance()->ReadString("113:onIfAll = 211 211");
634 // for Lambda_c -> Delta++ K-
635 AliTPythia8::Instance()->ReadString("2224:onMode = off");
636 AliTPythia8::Instance()->ReadString("2224:onIfAll = 2212 211");
637 // for Lambda_c -> Lambda(1520) K-
638 AliTPythia8::Instance()->ReadString("3124:onMode = off");
639 AliTPythia8::Instance()->ReadString("3124:onIfAll = 2212 321");
642 AliTPythia8::Instance()->ReadString("411:onMode = off");
643 AliTPythia8::Instance()->ReadString("421:onMode = off");
644 AliTPythia8::Instance()->ReadString("431:onMode = off");
645 AliTPythia8::Instance()->ReadString("4112:onMode = off");
646 AliTPythia8::Instance()->ReadString("4122:onMode = off");
649 AliTPythia8::Instance()->ReadString("411:onIfMatch = 321 211 211");
651 AliTPythia8::Instance()->ReadString("411:onIfMatch = 313 211");
653 AliTPythia8::Instance()->ReadString("421:onIfMatch = 321 211");
657 AliTPythia8::Instance()->ReadString("421:onIfMatch = 321 211 211 211");
659 AliTPythia8::Instance()->ReadString("421:onIfMatch = 321 211 113");
661 AliTPythia8::Instance()->ReadString("421:onIfMatch = 313 211 211");
665 AliTPythia8::Instance()->ReadString("431:onIfMatch = 321 313");
667 AliTPythia8::Instance()->ReadString("431:onIfMatch = 333 211");
670 AliTPythia8::Instance()->ReadString("4122:onIfMatch = 2212 313");
671 // Lambda_c -> Delta K
672 AliTPythia8::Instance()->ReadString("4122:onIfMatch = 2224 321");
673 // Lambda_c -> Lambda(1520) pi
674 AliTPythia8::Instance()->ReadString("4122:onIfMatch = 3124 211");
675 // Lambda_c -> p K pi
676 AliTPythia8::Instance()->ReadString("4122:onIfMatch = 2212 321 211");
677 // Lambda_c -> Lambda pi
678 AliTPythia8::Instance()->ReadString("4122:onIfMatch = 3122 211");
682 //___________________________________________________________________________
683 void AliDecayerPythia8::ReadDecayTable()
685 //to read a decay table (not yet implemented)
689 //___________________________________________________________________________
690 void AliDecayerPythia8::AppendParticle(Int_t pdg, TLorentzVector* p)
692 // Append a particle to the stack
693 fPythia8->Pythia8()->event.append(pdg, 11, 0, 0, p->Px(), p->Py(), p->Pz(), p->E(), p->M());
697 //___________________________________________________________________________
698 void AliDecayerPythia8::ClearEvent()
700 // Clear the event stack
701 fPythia8->Pythia8()->event.clear();