Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
TPythia8Decayer.cxx
Go to the documentation of this file.
1// @(#)root/pythia8:$Name$:$Id$
2// Author: Andreas Morsch 04/07/2008
3
4/**************************************************************************
5 * Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
6 * *
7 * Author: The ALICE Off-line Project. *
8 * Contributors are mentioned in the code where appropriate. *
9 * *
10 * Permission to use, copy, modify and distribute this software and its *
11 * documentation strictly for non-commercial purposes is hereby granted *
12 * without fee, provided that the above copyright notice appears in all *
13 * copies and that both the copyright notice and this permission notice *
14 * appear in the supporting documentation. The authors make no claims *
15 * about the suitability of this software for any purpose. It is *
16 * provided "as is" without express or implied warranty. *
17 **************************************************************************/
18
19/** \class TPythia8Decayer
20 \ingroup pythia8
21
22This class implements the TVirtualMCDecayer interface using TPythia8.
23
24Author: Andreas Morsch 04/07/2008
25*/
26
27#include "TLorentzVector.h"
28#include "TPythia8.h"
29#include "TPythia8Decayer.h"
30
32
33////////////////////////////////////////////////////////////////////////////////
34///constructor
35
37 fPythia8(new TPythia8()),
38 fDebug(0)
39{
40 fPythia8->Pythia8()->readString("SoftQCD:elastic = on");
41 fPythia8->Pythia8()->init();
42}
43
44////////////////////////////////////////////////////////////////////////////////
45/// Initialize the decayer
46
48{
49}
50
51////////////////////////////////////////////////////////////////////////////////
52/// Decay a single particle
53
55{
56 ClearEvent();
57 AppendParticle(pdg, p);
58 Int_t idPart = fPythia8->Pythia8()->event[0].id();
59 fPythia8->Pythia8()->particleData.mayDecay(idPart,kTRUE);
60 fPythia8->Pythia8()->moreDecays();
61 if (fDebug > 0) fPythia8->EventListing();
62}
63
64////////////////////////////////////////////////////////////////////////////////
65///import the decay products into particles array
66
68{
69 return (fPythia8->ImportParticles(particles, "All"));
70}
71
72////////////////////////////////////////////////////////////////////////////////
73/// Set forced decay mode
74
76{
77 printf("SetForceDecay not yet implemented !\n");
78}
79////////////////////////////////////////////////////////////////////////////////
80/// ForceDecay not yet implemented
81
83{
84 printf("ForceDecay not yet implemented !\n");
85}
86////////////////////////////////////////////////////////////////////////////////
87
89{
90 return 0.0;
91}
92////////////////////////////////////////////////////////////////////////////////
93///return lifetime in seconds of teh particle with PDG number pdg
94
96{
97 return (fPythia8->Pythia8()->particleData.tau0(pdg) * 3.3333e-12) ;
98}
99
100////////////////////////////////////////////////////////////////////////////////
101///to read a decay table (not yet implemented)
102
104{
105}
106
107
108////////////////////////////////////////////////////////////////////////////////
109/// Append a particle to the stack
110
112{
113 fPythia8->Pythia8()->event.append(pdg, 11, 0, 0, p->Px(), p->Py(), p->Pz(), p->E(), p->M());
114}
115
116
117////////////////////////////////////////////////////////////////////////////////
118/// Clear the event stack
119
121{
122 fPythia8->Pythia8()->event.clear();
123}
124
float Float_t
Definition RtypesCore.h:57
const Bool_t kTRUE
Definition RtypesCore.h:100
#define ClassImp(name)
Definition Rtypes.h:364
An array of clone (identical) objects.
Double_t Px() const
Double_t M() const
Double_t Pz() const
Double_t Py() const
Double_t E() const
This class implements the TVirtualMCDecayer interface using TPythia8.
TPythia8 * fPythia8
virtual Float_t GetLifetime(Int_t kf)
return lifetime in seconds of teh particle with PDG number pdg
virtual void Decay(Int_t pdg, TLorentzVector *p)
Decay a single particle.
virtual void Init()
Initialize the decayer.
TPythia8Decayer()
constructor
virtual Int_t ImportParticles(TClonesArray *particles)
import the decay products into particles array
virtual void SetForceDecay(Int_t type)
Set forced decay mode.
void AppendParticle(Int_t pdg, TLorentzVector *p)
Append a particle to the stack.
virtual void ReadDecayTable()
to read a decay table (not yet implemented)
virtual void ForceDecay()
ForceDecay not yet implemented.
void ClearEvent()
Clear the event stack.
virtual Float_t GetPartialBranchingRatio(Int_t ipart)
Get the partial branching ratio for a particle of type IPART (a PDG code).
TPythia8 is an interface class to C++ version of Pythia 8.1 event generators, written by T....
Definition TPythia8.h:77
Pythia8::Pythia * Pythia8()
Definition TPythia8.h:89
void EventListing() const
Event listing.
Definition TPythia8.cxx:364
virtual Int_t ImportParticles(TClonesArray *particles, Option_t *option="")
Import particles from Pythia stack.
Definition TPythia8.cxx:194