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 <TClonesArray.h>
24 #include <TParticle.h>
25 #include "AliTPythia8.h"
26 #include "AliDecayerPythia8.h"
27 #include "ParticleData.h"
29 ClassImp(AliDecayerPythia8)
31 Bool_t AliDecayerPythia8::fgInit = kFALSE;
33 AliDecayerPythia8::AliDecayerPythia8():
35 fPythia8(AliTPythia8::Instance()),
41 fPythia8->Pythia8()->readString("SoftQCD:elastic = on");
42 fPythia8->Pythia8()->init();
45 //___________________________________________________________________________
46 void AliDecayerPythia8::Decay(Int_t pdg, TLorentzVector* p)
48 // Decay a single particle
50 AppendParticle(pdg, p);
51 Int_t idPart = fPythia8->Pythia8()->event[0].id();
52 fPythia8->Pythia8()->particleData.mayDecay(idPart,kTRUE);
53 fPythia8->Pythia8()->moreDecays();
54 if (fDebug > 0) fPythia8->EventListing();
57 //___________________________________________________________________________
58 Int_t AliDecayerPythia8::ImportParticles(TClonesArray *particles)
60 //import the decay products into particles array
61 const Float_t kconvT=0.001/2.999792458e8; // mm/c to seconds conversion
62 const Float_t kconvL=1./10; // mm to cm conversion
63 int np = fPythia8->ImportParticles(particles, "All");
64 // pythia assigns decay time in mm/c, convert to seconds
65 for (int ip=1;ip<np;ip++) {
66 TParticle* prod = (TParticle*)particles->At(ip);
68 prod->SetProductionVertex(prod->Vx()*kconvL,prod->Vy()*kconvL,prod->Vz()*kconvL,kconvT*prod->T());
74 void AliDecayerPythia8::Init()
80 // fPythia->SetDecayTable();
83 // Switch on heavy flavor decays
86 Int_t heavy[14] = {411, 421, 431, 4122, 4132, 4232, 4332, 511, 521, 531, 5122, 5132, 5232, 5332};
87 // fPythia->ResetDecayTable();
88 for (j=0; j < 14; j++) {
89 if (fDecay == kNoDecayHeavy) {
90 fPythia8->ReadString(Form("%d:onMode = off", heavy[j]));
92 fPythia8->ReadString(Form("%d:onMode = on", heavy[j]));
97 //...Switch off decay of pi0, K0S, Lambda, Sigma+-, Xi0-, Omega-.
99 if (fDecay != kNeutralPion) {
100 fPythia8->ReadString("111:onMode = off");
102 fPythia8->ReadString("111:onMode = on");
105 fPythia8->ReadString("310:onMode = off");
106 fPythia8->ReadString("3112:onMode = off");
107 fPythia8->ReadString("3212:onMode = off");
108 fPythia8->ReadString("3222:onMode = off");
109 fPythia8->ReadString("3312:onMode = off");
110 fPythia8->ReadString("3322:onMode = off");
111 fPythia8->ReadString("3334:onMode = off");
112 // .. Force decay channels
116 void AliDecayerPythia8::ForceDecay()
119 // Force a particle decay mode
120 // Switch heavy flavour production off if requested
121 if (!fHeavyFlavour) SwitchOffHeavyFlavour();
123 Decay_t decay = fDecay;
124 fPythia8->ReadString("HadronLevel:Decay = on");
126 if (decay == kNoDecayHeavy) return;
134 fPythia8->ReadString("511:onMode = off");
135 fPythia8->ReadString("511:onIfAny = 13 443 100443");
137 fPythia8->ReadString("521:onMode = off");
138 fPythia8->ReadString("521:onIfAny = 13 443 100443");
140 fPythia8->ReadString("531:onMode = off");
141 fPythia8->ReadString("531:onIfAny = 13 443 100443");
143 fPythia8->ReadString("5122:onMode = off");
144 fPythia8->ReadString("5122:onIfAny = 13 443 100443");
146 fPythia8->ReadString("5132:onMode = off");
147 fPythia8->ReadString("5132:onIfAny = 13 443 100443");
149 fPythia8->ReadString("5232:onMode = off");
150 fPythia8->ReadString("5232:onIfAny = 13 443 100443");
152 fPythia8->ReadString("5332:onMode = off");
153 fPythia8->ReadString("5332:onIfAny = 13 443 100443");
155 fPythia8->ReadString("100443:onMode = off");
156 fPythia8->ReadString("100443:onIfAny = 443");
158 fPythia8->ReadString("443:onMode = off");
159 fPythia8->ReadString("443:onIfAll = 13 13");
161 fPythia8->ReadString("411:onMode = off");
162 fPythia8->ReadString("411:onIfAll = 13");
164 fPythia8->ReadString("421:onMode = off");
165 fPythia8->ReadString("421:onIfAll = 13");
167 fPythia8->ReadString("431:onMode = off");
168 fPythia8->ReadString("431:onIfAll = 13");
170 fPythia8->ReadString("4122:onMode = off");
171 fPythia8->ReadString("4122:onIfAll = 13");
173 fPythia8->ReadString("4132:onMode = off");
174 fPythia8->ReadString("4132:onIfAll = 13");
176 fPythia8->ReadString("4232:onMode = off");
177 fPythia8->ReadString("4232:onIfAll = 13");
179 fPythia8->ReadString("4332:onMode = off");
180 fPythia8->ReadString("4332:onIfAll = 13");
183 case kChiToJpsiGammaToMuonMuon:
184 // Chi_1c -> J/Psi Gamma
185 fPythia8->ReadString("20443:onMode = off");
186 fPythia8->ReadString("20443:onIfAll = 443 22");
187 // Chi_2c -> J/Psi Gamma
188 fPythia8->ReadString("445:onMode = off");
189 fPythia8->ReadString("445:onIfAll = 443 22");
191 fPythia8->ReadString("443:onMode = off");
192 fPythia8->ReadString("443:onIfAll = 13 13");
194 case kChiToJpsiGammaToElectronElectron:
195 // Chi_1c -> J/Psi Gamma
196 fPythia8->ReadString("20443:onMode = off");
197 fPythia8->ReadString("20443:onIfAll = 443 22");
198 // Chi_2c -> J/Psi Gamma
199 fPythia8->ReadString("445:onMode = off");
200 fPythia8->ReadString("445:onIfAll = 443 22");
202 fPythia8->ReadString("443:onMode = off");
203 fPythia8->ReadString("443:onIfAll = 11 11");
208 fPythia8->ReadString("511:onMode = off");
209 fPythia8->ReadString("511:onIfAny = 13");
211 fPythia8->ReadString("521:onMode = off");
212 fPythia8->ReadString("521:onIfAny = 13");
214 fPythia8->ReadString("531:onMode = off");
215 fPythia8->ReadString("531:onIfAny = 13");
217 fPythia8->ReadString("5122:onMode = off");
218 fPythia8->ReadString("5122:onIfAny = 13");
220 fPythia8->ReadString("5132:onMode = off");
221 fPythia8->ReadString("5132:onIfAny = 13");
223 fPythia8->ReadString("5232:onMode = off");
224 fPythia8->ReadString("5232:onIfAny = 13");
226 fPythia8->ReadString("5332:onMode = off");
227 fPythia8->ReadString("5332:onIfAny = 13");
231 fPythia8->ReadString("411:onMode = off");
232 fPythia8->ReadString("411:onIfAll = 13");
234 fPythia8->ReadString("421:onMode = off");
235 fPythia8->ReadString("421:onIfAll = 13");
237 fPythia8->ReadString("431:onMode = off");
238 fPythia8->ReadString("431:onIfAll = 13");
240 fPythia8->ReadString("4122:onMode = off");
241 fPythia8->ReadString("4122:onIfAll = 13");
243 fPythia8->ReadString("4132:onMode = off");
244 fPythia8->ReadString("4132:onIfAll = 13");
246 fPythia8->ReadString("4232:onMode = off");
247 fPythia8->ReadString("4232:onIfAll = 13");
249 fPythia8->ReadString("4332:onMode = off");
250 fPythia8->ReadString("4332:onIfAll = 13");
252 fPythia8->ReadString("511:onMode = off");
253 fPythia8->ReadString("511:onIfAny = 13");
255 fPythia8->ReadString("521:onMode = off");
256 fPythia8->ReadString("521:onIfAny = 13");
258 fPythia8->ReadString("531:onMode = off");
259 fPythia8->ReadString("531:onIfAny = 13");
261 fPythia8->ReadString("5122:onMode = off");
262 fPythia8->ReadString("5122:onIfAny = 13");
264 fPythia8->ReadString("5132:onMode = off");
265 fPythia8->ReadString("5132:onIfAny = 13");
267 fPythia8->ReadString("5232:onMode = off");
268 fPythia8->ReadString("5232:onIfAny = 13");
270 fPythia8->ReadString("5332:onMode = off");
271 fPythia8->ReadString("5332:onIfAny = 13");
276 fPythia8->ReadString("443:onMode = off");
277 fPythia8->ReadString("443:onIfAll = 13 13");
281 fPythia8->ReadString("113:onMode = off");
282 fPythia8->ReadString("113:onIfAll = 13 13");
284 fPythia8->ReadString("221:onMode = off");
285 fPythia8->ReadString("221:onIfAll = 13 13");
287 fPythia8->ReadString("223:onMode = off");
288 fPythia8->ReadString("223:onIfAll = 13 13");
290 fPythia8->ReadString("333:onMode = off");
291 fPythia8->ReadString("333:onIfAll = 13 13");
293 fPythia8->ReadString("443:onMode = off");
294 fPythia8->ReadString("443:onIfAll = 13 13");
296 fPythia8->ReadString("100443:onMode = off");
297 fPythia8->ReadString("100443:onIfAll = 13 13");
299 fPythia8->ReadString("553:onMode = off");
300 fPythia8->ReadString("553:onIfAll = 13 13");
302 fPythia8->ReadString("100553:onMode = off");
303 fPythia8->ReadString("100553:onIfAll = 13 13");
305 fPythia8->ReadString("200553:onMode = off");
306 fPythia8->ReadString("200553:onIfAll = 13 13");
308 case kBSemiElectronic:
310 fPythia8->ReadString("511:onMode = off");
311 fPythia8->ReadString("511:onIfAny = 11");
313 fPythia8->ReadString("521:onMode = off");
314 fPythia8->ReadString("521:onIfAny = 11");
316 fPythia8->ReadString("531:onMode = off");
317 fPythia8->ReadString("531:onIfAny = 11");
319 fPythia8->ReadString("5122:onMode = off");
320 fPythia8->ReadString("5122:onIfAny = 11");
322 fPythia8->ReadString("5132:onMode = off");
323 fPythia8->ReadString("5132:onIfAny = 11");
325 fPythia8->ReadString("5232:onMode = off");
326 fPythia8->ReadString("5232:onIfAny = 11");
328 fPythia8->ReadString("5332:onMode = off");
329 fPythia8->ReadString("5332:onIfAny = 11");
331 case kSemiElectronic:
333 fPythia8->ReadString("411:onMode = off");
334 fPythia8->ReadString("411:onIfAll = 11");
336 fPythia8->ReadString("421:onMode = off");
337 fPythia8->ReadString("421:onIfAll = 11");
339 fPythia8->ReadString("431:onMode = off");
340 fPythia8->ReadString("431:onIfAll = 11");
342 fPythia8->ReadString("4122:onMode = off");
343 fPythia8->ReadString("4122:onIfAll = 11");
345 fPythia8->ReadString("4132:onMode = off");
346 fPythia8->ReadString("4132:onIfAll = 11");
348 fPythia8->ReadString("4232:onMode = off");
349 fPythia8->ReadString("4232:onIfAll = 11");
351 fPythia8->ReadString("4332:onMode = off");
352 fPythia8->ReadString("4332:onIfAll = 11");
354 fPythia8->ReadString("511:onMode = off");
355 fPythia8->ReadString("511:onIfAny = 11");
357 fPythia8->ReadString("521:onMode = off");
358 fPythia8->ReadString("521:onIfAny = 11");
360 fPythia8->ReadString("531:onMode = off");
361 fPythia8->ReadString("531:onIfAny = 11");
363 fPythia8->ReadString("5122:onMode = off");
364 fPythia8->ReadString("5122:onIfAny = 11");
366 fPythia8->ReadString("5132:onMode = off");
367 fPythia8->ReadString("5132:onIfAny = 11");
369 fPythia8->ReadString("5232:onMode = off");
370 fPythia8->ReadString("5232:onIfAny = 11");
372 fPythia8->ReadString("5332:onMode = off");
373 fPythia8->ReadString("5332:onIfAny = 11");
377 fPythia8->ReadString("113:onMode = off");
378 fPythia8->ReadString("113:onIfAll = 11 11");
380 fPythia8->ReadString("221:onMode = off");
381 fPythia8->ReadString("221:onIfAll = 11 11");
383 fPythia8->ReadString("223:onMode = off");
384 fPythia8->ReadString("223:onIfAll = 11 11");
386 fPythia8->ReadString("333:onMode = off");
387 fPythia8->ReadString("333:onIfAll = 11 11");
389 fPythia8->ReadString("443:onMode = off");
390 fPythia8->ReadString("443:onIfAll = 11 11");
392 fPythia8->ReadString("100443:onMode = off");
393 fPythia8->ReadString("100443:onIfAll = 11 11");
395 fPythia8->ReadString("553:onMode = off");
396 fPythia8->ReadString("553:onIfAll = 11 11");
398 fPythia8->ReadString("100553:onMode = off");
399 fPythia8->ReadString("100553:onIfAll = 11 11");
401 fPythia8->ReadString("200553:onMode = off");
402 fPythia8->ReadString("200553:onIfAll = 11 11");
405 // B0 -> J/Psi (Psi') X
406 fPythia8->ReadString("511:onMode = off");
407 fPythia8->ReadString("511:onIfAny = 443 100443");
408 // B+/- -> J/Psi (Psi') X
409 fPythia8->ReadString("521:onMode = off");
410 fPythia8->ReadString("521:onIfAny = 443 100443");
411 // B_s -> J/Psi (Psi') X
412 fPythia8->ReadString("531:onMode = off");
413 fPythia8->ReadString("531:onIfAny = 443 100443");
414 // Lambda_b -> J/Psi (Psi') X
415 fPythia8->ReadString("5122:onMode = off");
416 fPythia8->ReadString("5122:onIfAny = 443 100443");
419 fPythia8->ReadString("443:onMode = off");
420 fPythia8->ReadString("443:onIfAll = 13 13");
422 fPythia8->ReadString("100443:onMode = off");
423 fPythia8->ReadString("100443:onIfAll = 13 13");
425 case kBPsiPrimeDiMuon:
427 fPythia8->ReadString("511:onMode = off");
428 fPythia8->ReadString("511:onIfAny = 100443");
430 fPythia8->ReadString("521:onMode = off");
431 fPythia8->ReadString("521:onIfAny = 100443");
433 fPythia8->ReadString("531:onMode = off");
434 fPythia8->ReadString("531:onIfAny = 100443");
435 // Lambda_b -> Psi' X
436 fPythia8->ReadString("5122:onMode = off");
437 fPythia8->ReadString("5122:onIfAny = 100443");
440 fPythia8->ReadString("100443:onMode = off");
441 fPythia8->ReadString("100443:onIfAll = 13 13");
443 case kBJpsiDiElectron:
445 fPythia8->ReadString("511:onMode = off");
446 fPythia8->ReadString("511:onIfAny = 443");
448 fPythia8->ReadString("521:onMode = off");
449 fPythia8->ReadString("521:onIfAny = 443");
451 fPythia8->ReadString("531:onMode = off");
452 fPythia8->ReadString("531:onIfAny = 443");
454 fPythia8->ReadString("5122:onMode = off");
455 fPythia8->ReadString("5122:onIfAny = 443");
458 fPythia8->ReadString("443:onMode = off");
459 fPythia8->ReadString("443:onIfAll = 11 11");
464 fPythia8->ReadString("511:onMode = off");
465 fPythia8->ReadString("511:onIfAny = 443");
467 fPythia8->ReadString("521:onMode = off");
468 fPythia8->ReadString("521:onIfAny = 443");
470 fPythia8->ReadString("531:onMode = off");
471 fPythia8->ReadString("531:onIfAny = 443");
473 fPythia8->ReadString("5122:onMode = off");
474 fPythia8->ReadString("5122:onIfAny = 443");
476 case kBPsiPrimeDiElectron:
478 fPythia8->ReadString("511:onMode = off");
479 fPythia8->ReadString("511:onIfAny = 100443");
481 fPythia8->ReadString("521:onMode = off");
482 fPythia8->ReadString("521:onIfAny = 100443");
484 fPythia8->ReadString("531:onMode = off");
485 fPythia8->ReadString("531:onIfAny = 100443");
486 // Lambda_b -> Psi' X
487 fPythia8->ReadString("5122:onMode = off");
488 fPythia8->ReadString("5122:onIfAny = 100443");
491 fPythia8->ReadString("100443:onMode = off");
492 fPythia8->ReadString("100443:onIfAll = 11 11");
496 fPythia8->ReadString("211:onMode = off");
497 fPythia8->ReadString("211:onIfAny = 13");
501 fPythia8->ReadString("321:onMode = off");
502 fPythia8->ReadString("321:onIfAny = 13");
506 fPythia8->ReadString("211:onMode = off");
507 fPythia8->ReadString("211:onIfAny = 13");
508 fPythia8->ReadString("321:onMode = off");
509 fPythia8->ReadString("321:onIfAny = 13");
513 fPythia8->ReadString("24:onMode = off");
514 fPythia8->ReadString("24:onIfAny = 13");
518 fPythia8->ReadString("24:onMode = off");
519 fPythia8->ReadString("24:onIfAny = 4");
521 case kWToCharmToMuon:
523 fPythia8->ReadString("24:onMode = off");
524 fPythia8->ReadString("24:onIfAny = 4");
526 fPythia8->ReadString("411:onMode = off");
527 fPythia8->ReadString("411:onIfAll = 13");
529 fPythia8->ReadString("421:onMode = off");
530 fPythia8->ReadString("421:onIfAll = 13");
532 fPythia8->ReadString("431:onMode = off");
533 fPythia8->ReadString("431:onIfAll = 13");
535 fPythia8->ReadString("4122:onMode = off");
536 fPythia8->ReadString("4122:onIfAll = 13");
538 fPythia8->ReadString("4132:onMode = off");
539 fPythia8->ReadString("4132:onIfAll = 13");
541 fPythia8->ReadString("4232:onMode = off");
542 fPythia8->ReadString("4232:onIfAll = 13");
544 fPythia8->ReadString("4332:onMode = off");
545 fPythia8->ReadString("4332:onIfAll = 13");
549 fPythia8->ReadString("23:onMode = off");
550 fPythia8->ReadString("23:onIfAll = 13 13");
554 fPythia8->ReadString("23:onMode = off");
555 fPythia8->ReadString("23:onIfAll = 11 11");
560 case kHadronicDWithout4Bodies:
565 fPythia8->ReadString("333:onMode = off");
566 fPythia8->ReadString("333:onIfAll = 321 321");
570 fPythia8->ReadString("3334:onMode = off");
571 fPythia8->ReadString("3334:onIfAll = 3122 321 ");
575 fPythia8->ReadString("3122:onMode = off");
576 fPythia8->ReadString("3122:onIfAll = 2212 211 ");
579 fPythia8->ReadString("5122:onMode = off");
580 fPythia8->ReadString("4122:onMode = off");
581 fPythia8->ReadString("5122:onIfAll = 4122");
582 fPythia8->ReadString("4122:onIfAll = 3122");
587 fPythia8->ReadString("HadronLevel:Decay = off");
592 case kPsiPrimeJpsiDiElectron:
600 Float_t AliDecayerPythia8::GetPartialBranchingRatio(Int_t ipart)
602 // Get the partial branching ration for the forced decay channels
604 Pythia8::Pythia* thePythia = fPythia8->Pythia8();
605 Pythia8::ParticleData & table = thePythia->particleData;
606 Pythia8::ParticleDataEntry* pd = table.particleDataEntryPtr(ipart);
608 Int_t nc = pd->sizeChannels();
611 // Loop over decay channels
612 for (Int_t ic = 0; ic < nc; ic++) {
613 Pythia8::DecayChannel& decCh = pd->channel(ic);
614 for (Int_t i = 0; i < decCh.multiplicity(); i++) {
615 br += decCh.bRatio();
622 Float_t AliDecayerPythia8::GetLifetime(Int_t kf)
624 // Return lifetime of particle
625 Pythia8::Pythia* thePythia = fPythia8->Pythia8();
626 Pythia8::ParticleData& table = thePythia->particleData;
627 Float_t tau = table.tau0(kf);
631 void AliDecayerPythia8::SwitchOffHeavyFlavour()
633 // Switch off heavy flavour production
635 // Maximum number of quark flavours used in pdf
636 fPythia8->ReadString("PDFinProcess:nQuarkIn = 3");
637 // Maximum number of flavors that can be used in showers
638 fPythia8->ReadString("SpaceShower:nQuarkIn = 3");
639 fPythia8->ReadString("TimeShower:nGammaToQuark = 3");
640 fPythia8->ReadString("TimeShower:nGluonToQuark = 3");
644 void AliDecayerPythia8::ForceHadronicD(Int_t optUse4Bodies)
647 // Force golden D decay modes
650 fPythia8->ReadString("313:onMode = off");
651 fPythia8->ReadString("313:onIfAll = 321 211");
653 fPythia8->ReadString("333:onMode = off");
654 fPythia8->ReadString("333:onIfAll = 321 321");
655 // for D0 -> rho0 pi+ k-
656 fPythia8->ReadString("113:onMode = off");
657 fPythia8->ReadString("113:onIfAll = 211 211");
658 // for Lambda_c -> Delta++ K-
659 fPythia8->ReadString("2224:onMode = off");
660 fPythia8->ReadString("2224:onIfAll = 2212 211");
661 // for Lambda_c -> Lambda(1520) K-
662 fPythia8->ReadString("3124:onMode = off");
663 fPythia8->ReadString("3124:onIfAll = 2212 321");
666 fPythia8->ReadString("411:onMode = off");
667 fPythia8->ReadString("421:onMode = off");
668 fPythia8->ReadString("431:onMode = off");
669 fPythia8->ReadString("4112:onMode = off");
670 fPythia8->ReadString("4122:onMode = off");
673 fPythia8->ReadString("411:onIfMatch = 321 211 211");
675 fPythia8->ReadString("411:onIfMatch = 313 211");
677 fPythia8->ReadString("421:onIfMatch = 321 211");
681 fPythia8->ReadString("421:onIfMatch = 321 211 211 211");
683 fPythia8->ReadString("421:onIfMatch = 321 211 113");
685 fPythia8->ReadString("421:onIfMatch = 313 211 211");
689 fPythia8->ReadString("431:onIfMatch = 321 313");
691 fPythia8->ReadString("431:onIfMatch = 333 211");
694 fPythia8->ReadString("4122:onIfMatch = 2212 313");
695 // Lambda_c -> Delta K
696 fPythia8->ReadString("4122:onIfMatch = 2224 321");
697 // Lambda_c -> Lambda(1520) pi
698 fPythia8->ReadString("4122:onIfMatch = 3124 211");
699 // Lambda_c -> p K pi
700 fPythia8->ReadString("4122:onIfMatch = 2212 321 211");
701 // Lambda_c -> Lambda pi
702 fPythia8->ReadString("4122:onIfMatch = 3122 211");
706 //___________________________________________________________________________
707 void AliDecayerPythia8::ReadDecayTable()
709 //to read a decay table (not yet implemented)
713 //___________________________________________________________________________
714 void AliDecayerPythia8::AppendParticle(Int_t pdg, TLorentzVector* p)
716 // Append a particle to the stack
717 fPythia8->Pythia8()->event.append(pdg, 11, 0, 0, p->Px(), p->Py(), p->Pz(), p->E(), p->M());
721 //___________________________________________________________________________
722 void AliDecayerPythia8::ClearEvent()
724 // Clear the event stack
725 fPythia8->Pythia8()->event.clear();