Geant4  10.02.p03
G4Pythia6Decayer Class Reference

#include <G4Pythia6Decayer.hh>

Inheritance diagram for G4Pythia6Decayer:
Collaboration diagram for G4Pythia6Decayer:

Public Member Functions

 G4Pythia6Decayer ()
 
virtual ~G4Pythia6Decayer ()
 
virtual G4DecayProductsImportDecayProducts (const G4Track &track)
 
void ForceDecayType (EDecayType decayType)
 
void SetVerboseLevel (G4int verboseLevel)
 
- Public Member Functions inherited from G4VExtDecayer
 G4VExtDecayer (const G4String &name="")
 
virtual ~G4VExtDecayer ()
 
const G4StringGetName () const
 

Private Member Functions

 G4Pythia6Decayer (const G4Pythia6Decayer &right)
 Not implemented. More...
 
G4Pythia6Decayeroperator= (const G4Pythia6Decayer &right)
 Not implemented. More...
 
G4ParticleDefinitionGetParticleDefinition (const Pythia6Particle *p, G4bool warn=true) const
 
G4DynamicParticleCreateDynamicParticle (const Pythia6Particle *p) const
 
G4ThreeVector GetParticlePosition (const Pythia6Particle *particle) const
 
G4ThreeVector GetParticleMomentum (const Pythia6Particle *particle) const
 
G4int CountProducts (G4int channel, G4int particle)
 
void ForceParticleDecay (G4int particle, G4int product, G4int mult)
 
void ForceParticleDecay (G4int particle, G4int *products, G4int *mult, G4int npart)
 
void ForceHadronicD ()
 
void ForceOmega ()
 
void ForceDecay (EDecayType decayType)
 
void Decay (G4int pdg, const CLHEP::HepLorentzVector &p)
 
G4int ImportParticles (ParticleVector *particles)
 

Private Attributes

G4Pythia6DecayerMessenger fMessenger
 command messenger More...
 
G4int fVerboseLevel
 verbose level More...
 
EDecayType fDecayType
 selected decay type More...
 
ParticleVectorfDecayProductsArray
 array of decay products More...
 

Static Private Attributes

static const EDecayType fgkDefaultDecayType = kAll
 default decay type More...
 

Additional Inherited Members

- Protected Attributes inherited from G4VExtDecayer
G4String decayerName
 

Detailed Description

Pythia6 decayer

Implements the G4VExtDecayer abstract class using the Pythia6 interface. According to TPythia6Decayer class in Root: http://root.cern.ch/ see http://root.cern.ch/root/License.html

Definition at line 53 of file G4Pythia6Decayer.hh.

Constructor & Destructor Documentation

◆ G4Pythia6Decayer() [1/2]

G4Pythia6Decayer::G4Pythia6Decayer ( )

Standard constructor

Definition at line 55 of file G4Pythia6Decayer.cc.

56  : G4VExtDecayer("G4Pythia6Decayer"),
57  fMessenger(this),
58  fVerboseLevel(0),
61 {
63 
65 
67 }
static const EDecayType fgkDefaultDecayType
default decay type
ParticleVector * fDecayProductsArray
array of decay products
std::vector< Pythia6Particle * > ParticleVector
Definition: Pythia6.hh:149
EDecayType fDecayType
selected decay type
G4VExtDecayer(const G4String &name="")
void ForceDecay(EDecayType decayType)
G4int fVerboseLevel
verbose level
G4Pythia6DecayerMessenger fMessenger
command messenger
Here is the call graph for this function:
Here is the caller graph for this function:

◆ ~G4Pythia6Decayer()

G4Pythia6Decayer::~G4Pythia6Decayer ( )
virtual

Destructor

Definition at line 71 of file G4Pythia6Decayer.cc.

72 {
74 
75  delete fDecayProductsArray;
76 }
ParticleVector * fDecayProductsArray
array of decay products
Here is the call graph for this function:

◆ G4Pythia6Decayer() [2/2]

G4Pythia6Decayer::G4Pythia6Decayer ( const G4Pythia6Decayer right)
private

Not implemented.

Member Function Documentation

◆ CountProducts()

G4int G4Pythia6Decayer::CountProducts ( G4int  channel,
G4int  particle 
)
private

Count number of decay products

Definition at line 160 of file G4Pythia6Decayer.cc.

161 {
163 
164  G4int np = 0;
165  for ( G4int i=1; i<=5; i++ )
166  if ( std::abs(Pythia6::Instance()->GetKFDP(channel,i) ) == particle )
167  np++;
168  return np;
169 }
static Pythia6 * Instance()
Definition: Pythia6.cc:117
int G4int
Definition: G4Types.hh:78
Here is the call graph for this function:
Here is the caller graph for this function:

◆ CreateDynamicParticle()

G4DynamicParticle * G4Pythia6Decayer::CreateDynamicParticle ( const Pythia6Particle p) const
private

Create G4DynamicParticle.

Definition at line 112 of file G4Pythia6Decayer.cc.

113 {
115 
116  // get particle properties
117  const G4ParticleDefinition* particleDefinition
118  = GetParticleDefinition(particle);
119  if ( ! particleDefinition ) return 0;
120 
121  G4ThreeVector momentum = GetParticleMomentum(particle);
122 
123  // create G4DynamicParticle
124  G4DynamicParticle* dynamicParticle
125  = new G4DynamicParticle(particleDefinition, momentum);
126 
127  return dynamicParticle;
128 }
G4ThreeVector GetParticleMomentum(const Pythia6Particle *particle) const
G4ParticleDefinition * GetParticleDefinition(const Pythia6Particle *p, G4bool warn=true) const
Here is the call graph for this function:
Here is the caller graph for this function:

◆ Decay()

void G4Pythia6Decayer::Decay ( G4int  pdg,
const CLHEP::HepLorentzVector p 
)
private

Decay a particle of type IDPART (PDG code) and momentum P.

Definition at line 517 of file G4Pythia6Decayer.cc.

518 {
520 
521  Pythia6::Instance()->Py1ent(0, pdg, p.e(), p.theta(), p.phi());
522 }
double theta() const
static Pythia6 * Instance()
Definition: Pythia6.cc:117
void Py1ent(int line, int kf, double pe, double theta, double phi)
Definition: Pythia6.cc:179
Here is the call graph for this function:
Here is the caller graph for this function:

◆ ForceDecay()

void G4Pythia6Decayer::ForceDecay ( EDecayType  decayType)
private

Force a particle decay mode

Definition at line 315 of file G4Pythia6Decayer.cc.

316 {
318 
319  Pythia6::Instance()->SetMSTJ(21,2);
320 
321  if ( fDecayType == kNoDecayHeavy ) return;
322 
323  //
324  // select mode
325  G4int products[3];
326  G4int mult[3];
327 
328  switch ( decayType ) {
329 
330  case kHardMuons:
331  products[0] = 13;
332  products[1] = 443;
333  products[2] = 100443;
334  mult[0] = 1;
335  mult[1] = 1;
336  mult[2] = 1;
337  ForceParticleDecay( 511, products, mult, 3);
338  ForceParticleDecay( 521, products, mult, 3);
339  ForceParticleDecay( 531, products, mult, 3);
340  ForceParticleDecay( 5122, products, mult, 3);
341  ForceParticleDecay( 5132, products, mult, 3);
342  ForceParticleDecay( 5232, products, mult, 3);
343  ForceParticleDecay( 5332, products, mult, 3);
344  ForceParticleDecay( 100443, 443, 1); // Psi' -> J/Psi X
345  ForceParticleDecay( 443, 13, 2); // J/Psi -> mu+ mu-
346 
347  ForceParticleDecay( 411,13,1); // D+/-
348  ForceParticleDecay( 421,13,1); // D0
349  ForceParticleDecay( 431,13,1); // D_s
350  ForceParticleDecay( 4122,13,1); // Lambda_c
351  ForceParticleDecay( 4132,13,1); // Xsi_c
352  ForceParticleDecay( 4232,13,1); // Sigma_c
353  ForceParticleDecay( 4332,13,1); // Omega_c
354  break;
355 
356  case kSemiMuonic:
357  ForceParticleDecay( 411,13,1); // D+/-
358  ForceParticleDecay( 421,13,1); // D0
359  ForceParticleDecay( 431,13,1); // D_s
360  ForceParticleDecay( 4122,13,1); // Lambda_c
361  ForceParticleDecay( 4132,13,1); // Xsi_c
362  ForceParticleDecay( 4232,13,1); // Sigma_c
363  ForceParticleDecay( 4332,13,1); // Omega_c
364  ForceParticleDecay( 511,13,1); // B0
365  ForceParticleDecay( 521,13,1); // B+/-
366  ForceParticleDecay( 531,13,1); // B_s
367  ForceParticleDecay( 5122,13,1); // Lambda_b
368  ForceParticleDecay( 5132,13,1); // Xsi_b
369  ForceParticleDecay( 5232,13,1); // Sigma_b
370  ForceParticleDecay( 5332,13,1); // Omega_b
371  break;
372 
373  case kDiMuon:
374  ForceParticleDecay( 113,13,2); // rho
375  ForceParticleDecay( 221,13,2); // eta
376  ForceParticleDecay( 223,13,2); // omega
377  ForceParticleDecay( 333,13,2); // phi
378  ForceParticleDecay( 443,13,2); // J/Psi
379  ForceParticleDecay(100443,13,2);// Psi'
380  ForceParticleDecay( 553,13,2); // Upsilon
381  ForceParticleDecay(100553,13,2);// Upsilon'
382  ForceParticleDecay(200553,13,2);// Upsilon''
383  break;
384 
385  case kSemiElectronic:
386  ForceParticleDecay( 411,11,1); // D+/-
387  ForceParticleDecay( 421,11,1); // D0
388  ForceParticleDecay( 431,11,1); // D_s
389  ForceParticleDecay( 4122,11,1); // Lambda_c
390  ForceParticleDecay( 4132,11,1); // Xsi_c
391  ForceParticleDecay( 4232,11,1); // Sigma_c
392  ForceParticleDecay( 4332,11,1); // Omega_c
393  ForceParticleDecay( 511,11,1); // B0
394  ForceParticleDecay( 521,11,1); // B+/-
395  ForceParticleDecay( 531,11,1); // B_s
396  ForceParticleDecay( 5122,11,1); // Lambda_b
397  ForceParticleDecay( 5132,11,1); // Xsi_b
398  ForceParticleDecay( 5232,11,1); // Sigma_b
399  ForceParticleDecay( 5332,11,1); // Omega_b
400  break;
401 
402  case kDiElectron:
403  ForceParticleDecay( 113,11,2); // rho
404  ForceParticleDecay( 333,11,2); // phi
405  ForceParticleDecay( 221,11,2); // eta
406  ForceParticleDecay( 223,11,2); // omega
407  ForceParticleDecay( 443,11,2); // J/Psi
408  ForceParticleDecay(100443,11,2);// Psi'
409  ForceParticleDecay( 553,11,2); // Upsilon
410  ForceParticleDecay(100553,11,2);// Upsilon'
411  ForceParticleDecay(200553,11,2);// Upsilon''
412  break;
413 
414  case kBJpsiDiMuon:
415 
416  products[0] = 443;
417  products[1] = 100443;
418  mult[0] = 1;
419  mult[1] = 1;
420 
421  ForceParticleDecay( 511, products, mult, 2); // B0 -> J/Psi (Psi') X
422  ForceParticleDecay( 521, products, mult, 2); // B+/- -> J/Psi (Psi') X
423  ForceParticleDecay( 531, products, mult, 2); // B_s -> J/Psi (Psi') X
424  ForceParticleDecay( 5122, products, mult, 2); // Lambda_b -> J/Psi (Psi')X
425  ForceParticleDecay( 100443, 443, 1); // Psi' -> J/Psi X
426  ForceParticleDecay( 443,13,2); // J/Psi -> mu+ mu-
427  break;
428 
429  case kBPsiPrimeDiMuon:
430  ForceParticleDecay( 511,100443,1); // B0
431  ForceParticleDecay( 521,100443,1); // B+/-
432  ForceParticleDecay( 531,100443,1); // B_s
433  ForceParticleDecay( 5122,100443,1); // Lambda_b
434  ForceParticleDecay(100443,13,2); // Psi'
435  break;
436 
437  case kBJpsiDiElectron:
438  ForceParticleDecay( 511,443,1); // B0
439  ForceParticleDecay( 521,443,1); // B+/-
440  ForceParticleDecay( 531,443,1); // B_s
441  ForceParticleDecay( 5122,443,1); // Lambda_b
442  ForceParticleDecay( 443,11,2); // J/Psi
443  break;
444 
445  case kBJpsi:
446  ForceParticleDecay( 511,443,1); // B0
447  ForceParticleDecay( 521,443,1); // B+/-
448  ForceParticleDecay( 531,443,1); // B_s
449  ForceParticleDecay( 5122,443,1); // Lambda_b
450  break;
451 
453  ForceParticleDecay( 511,100443,1); // B0
454  ForceParticleDecay( 521,100443,1); // B+/-
455  ForceParticleDecay( 531,100443,1); // B_s
456  ForceParticleDecay( 5122,100443,1); // Lambda_b
457  ForceParticleDecay(100443,11,2); // Psi'
458  break;
459 
460  case kPiToMu:
461  ForceParticleDecay(211,13,1); // pi->mu
462  break;
463 
464  case kKaToMu:
465  ForceParticleDecay(321,13,1); // K->mu
466  break;
467 
468  case kWToMuon:
469  ForceParticleDecay( 24, 13,1); // W -> mu
470  break;
471 
472  case kWToCharm:
473  ForceParticleDecay( 24, 4,1); // W -> c
474  break;
475 
476  case kWToCharmToMuon:
477  ForceParticleDecay( 24, 4,1); // W -> c
478  ForceParticleDecay( 411,13,1); // D+/- -> mu
479  ForceParticleDecay( 421,13,1); // D0 -> mu
480  ForceParticleDecay( 431,13,1); // D_s -> mu
481  ForceParticleDecay( 4122,13,1); // Lambda_c
482  ForceParticleDecay( 4132,13,1); // Xsi_c
483  ForceParticleDecay( 4232,13,1); // Sigma_c
484  ForceParticleDecay( 4332,13,1); // Omega_c
485  break;
486 
487  case kZDiMuon:
488  ForceParticleDecay( 23, 13,2); // Z -> mu+ mu-
489  break;
490 
491  case kHadronicD:
492  ForceHadronicD();
493  break;
494 
495  case kPhiKK:
496  ForceParticleDecay(333,321,2); // Phi->K+K-
497  break;
498 
499  case kOmega:
500  ForceOmega();
501 
502  case kAll:
503  break;
504 
505  case kNoDecay:
506  Pythia6::Instance()->SetMSTJ(21,0);
507  break;
508 
509  case kNoDecayHeavy: break;
510 
511  case kMaxDecay: break;
512  }
513 }
void SetMSTJ(int i, int m)
Definition: Pythia6.hh:182
static Pythia6 * Instance()
Definition: Pythia6.cc:117
int G4int
Definition: G4Types.hh:78
EDecayType fDecayType
selected decay type
void ForceParticleDecay(G4int particle, G4int product, G4int mult)
Here is the call graph for this function:
Here is the caller graph for this function:

◆ ForceDecayType()

void G4Pythia6Decayer::ForceDecayType ( EDecayType  decayType)

Force a given decay type

Definition at line 618 of file G4Pythia6Decayer.cc.

619 {
621 
622  // Do nothing if the decay type is not different from current one
623  if ( decayType == fDecayType ) return;
624 
625  fDecayType = decayType;
627 }
EDecayType fDecayType
selected decay type
void ForceDecay(EDecayType decayType)
Here is the call graph for this function:
Here is the caller graph for this function:

◆ ForceHadronicD()

void G4Pythia6Decayer::ForceHadronicD ( )
private

Force golden D decay modes

Definition at line 226 of file G4Pythia6Decayer.cc.

227 {
229 
230  const G4int kNHadrons = 4;
231  G4int channel;
232  G4int hadron[kNHadrons] = {411, 421, 431, 4112};
233 
234  // for D+ -> K0* (-> K- pi+) pi+
235  G4int iKstar0 = 313;
236  G4int iKstarbar0 = -313;
237  G4int iKPlus = 321;
238  G4int iKMinus = -321;
239  G4int iPiPlus = 211;
240  G4int iPiMinus = -211;
241 
242  G4int products[2] = {iKPlus, iPiMinus}, mult[2] = {1, 1};
243  ForceParticleDecay(iKstar0, products, mult, 2);
244 
245  // for Ds -> Phi pi+
246  G4int iPhi = 333;
247  ForceParticleDecay(iPhi,iKPlus,2); // Phi->K+K-
248 
249  G4int decayP1[kNHadrons][3] = {
250  {iKMinus, iPiPlus, iPiPlus},
251  {iKMinus, iPiPlus, 0 },
252  {iKPlus , iKstarbar0, 0 },
253  {-1 , -1 , -1 }
254  };
255  G4int decayP2[kNHadrons][3] = {
256  {iKstarbar0, iPiPlus, 0 },
257  {-1 , -1 , -1 },
258  {iPhi , iPiPlus, 0 },
259  {-1 , -1 , -1 }
260  };
261 
262  Pythia6* pythia6 = Pythia6::Instance();
263  for ( G4int ihadron = 0; ihadron < kNHadrons; ihadron++ ) {
264  G4int kc = pythia6->Pycomp(hadron[ihadron]);
265  pythia6->SetMDCY(kc,1,1);
266  G4int ifirst = pythia6->GetMDCY(kc,2);
267  G4int ilast = ifirst + pythia6->GetMDCY(kc,3)-1;
268 
269  for (channel = ifirst; channel <= ilast; channel++) {
270  if ((pythia6->GetKFDP(channel,1) == decayP1[ihadron][0] &&
271  pythia6->GetKFDP(channel,2) == decayP1[ihadron][1] &&
272  pythia6->GetKFDP(channel,3) == decayP1[ihadron][2] &&
273  pythia6->GetKFDP(channel,4) == 0) ||
274  (pythia6->GetKFDP(channel,1) == decayP2[ihadron][0] &&
275  pythia6->GetKFDP(channel,2) == decayP2[ihadron][1] &&
276  pythia6->GetKFDP(channel,3) == decayP2[ihadron][2] &&
277  pythia6->GetKFDP(channel,4) == 0)) {
278  pythia6->SetMDME(channel,1,1);
279  } else {
280  pythia6->SetMDME(channel,1,0);
281  } // selected channel ?
282  } // decay channels
283  } // hadrons
284 }
static Pythia6 * Instance()
Definition: Pythia6.cc:117
int Pycomp(int kf)
Definition: Pythia6.cc:170
int GetMDCY(int i, int j)
Definition: Pythia6.hh:186
int G4int
Definition: G4Types.hh:78
int GetKFDP(int i, int j)
Definition: Pythia6.hh:187
void ForceParticleDecay(G4int particle, G4int product, G4int mult)
void SetMDCY(int i, int j, int m)
Definition: Pythia6.hh:188
void SetMDME(int i, int j, int m)
Definition: Pythia6.hh:189
Here is the call graph for this function:
Here is the caller graph for this function:

◆ ForceOmega()

void G4Pythia6Decayer::ForceOmega ( )
private

Force Omega -> Lambda K- Decay

Definition at line 288 of file G4Pythia6Decayer.cc.

289 {
291 
292  Pythia6* pythia6 = Pythia6::Instance();
293 
294  G4int iLambda0 = 3122;
295  G4int iKMinus = -321;
296 
297  G4int kc = pythia6->Pycomp(3334);
298  pythia6->SetMDCY(kc,1,1);
299  G4int ifirst = pythia6->GetMDCY(kc,2);
300  G4int ilast = ifirst + pythia6->GetMDCY(kc,3)-1;
301 
302  for (G4int channel = ifirst; channel <= ilast; channel++) {
303  if (pythia6->GetKFDP(channel,1) == iLambda0 &&
304  pythia6->GetKFDP(channel,2) == iKMinus &&
305  pythia6->GetKFDP(channel,3) == 0)
306  pythia6->SetMDME(channel,1,1);
307  else
308  pythia6->SetMDME(channel,1,0);
309  // selected channel ?
310  } // decay channels
311 }
static Pythia6 * Instance()
Definition: Pythia6.cc:117
int Pycomp(int kf)
Definition: Pythia6.cc:170
int GetMDCY(int i, int j)
Definition: Pythia6.hh:186
int G4int
Definition: G4Types.hh:78
int GetKFDP(int i, int j)
Definition: Pythia6.hh:187
void SetMDCY(int i, int j, int m)
Definition: Pythia6.hh:188
void SetMDME(int i, int j, int m)
Definition: Pythia6.hh:189
Here is the call graph for this function:
Here is the caller graph for this function:

◆ ForceParticleDecay() [1/2]

void G4Pythia6Decayer::ForceParticleDecay ( G4int  particle,
G4int  product,
G4int  mult 
)
private

Force decay of particle into products with multiplicity mult

Definition at line 174 of file G4Pythia6Decayer.cc.

175 {
177 
178  Pythia6* pythia6 = Pythia6::Instance();
179 
180  G4int kc = pythia6->Pycomp(particle);
181  pythia6->SetMDCY(kc,1,1);
182 
183  G4int ifirst = pythia6->GetMDCY(kc,2);
184  G4int ilast = ifirst + pythia6->GetMDCY(kc,3)-1;
185 
186  //
187  // Loop over decay channels
188  for (G4int channel= ifirst; channel <= ilast; channel++) {
189  if (CountProducts(channel,product) >= mult) {
190  pythia6->SetMDME(channel,1,1);
191  } else {
192  pythia6->SetMDME(channel,1,0);
193  }
194  }
195 }
static Pythia6 * Instance()
Definition: Pythia6.cc:117
int Pycomp(int kf)
Definition: Pythia6.cc:170
int GetMDCY(int i, int j)
Definition: Pythia6.hh:186
int G4int
Definition: G4Types.hh:78
G4int CountProducts(G4int channel, G4int particle)
void SetMDCY(int i, int j, int m)
Definition: Pythia6.hh:188
void SetMDME(int i, int j, int m)
Definition: Pythia6.hh:189
Here is the call graph for this function:
Here is the caller graph for this function:

◆ ForceParticleDecay() [2/2]

void G4Pythia6Decayer::ForceParticleDecay ( G4int  particle,
G4int products,
G4int mult,
G4int  npart 
)
private

Force decay of particle into products with multiplicity mult

Definition at line 199 of file G4Pythia6Decayer.cc.

201 {
203 
204  Pythia6* pythia6 = Pythia6::Instance();
205 
206  G4int kc = pythia6->Pycomp(particle);
207  pythia6->SetMDCY(kc,1,1);
208  G4int ifirst = pythia6->GetMDCY(kc,2);
209  G4int ilast = ifirst+pythia6->GetMDCY(kc,3)-1;
210  //
211  // Loop over decay channels
212  for (G4int channel = ifirst; channel <= ilast; channel++) {
213  G4int nprod = 0;
214  for (G4int i = 0; i < npart; i++)
215  nprod += (CountProducts(channel, products[i]) >= mult[i]);
216  if (nprod)
217  pythia6->SetMDME(channel,1,1);
218  else {
219  pythia6->SetMDME(channel,1,0);
220  }
221  }
222 }
static Pythia6 * Instance()
Definition: Pythia6.cc:117
int Pycomp(int kf)
Definition: Pythia6.cc:170
int GetMDCY(int i, int j)
Definition: Pythia6.hh:186
int G4int
Definition: G4Types.hh:78
G4int CountProducts(G4int channel, G4int particle)
Int_t npart
void SetMDCY(int i, int j, int m)
Definition: Pythia6.hh:188
void SetMDME(int i, int j, int m)
Definition: Pythia6.hh:189
Here is the call graph for this function:

◆ GetParticleDefinition()

G4ParticleDefinition * G4Pythia6Decayer::GetParticleDefinition ( const Pythia6Particle p,
G4bool  warn = true 
) const
private

Return G4 particle definition for given TParticle

Definition at line 86 of file G4Pythia6Decayer.cc.

87 {
89 
90  // get particle definition from G4ParticleTable
91  G4int pdgEncoding = particle->fKF;
92  G4ParticleTable* particleTable
94  G4ParticleDefinition* particleDefinition = 0;
95  if (pdgEncoding != 0)
96  particleDefinition = particleTable->FindParticle(pdgEncoding);
97 
98  if ( particleDefinition == 0 && warn) {
99  G4cerr
100  << "G4Pythia6Decayer: GetParticleDefinition: " << std::endl
101  << "G4ParticleTable::FindParticle() for particle with PDG = "
102  << pdgEncoding
103  << " failed." << std::endl;
104  }
105 
106  return particleDefinition;
107 }
G4ParticleDefinition * FindParticle(G4int PDGEncoding)
int G4int
Definition: G4Types.hh:78
static G4ParticleTable * GetParticleTable()
G4GLOB_DLL std::ostream G4cerr
Here is the call graph for this function:
Here is the caller graph for this function:

◆ GetParticleMomentum()

G4ThreeVector G4Pythia6Decayer::GetParticleMomentum ( const Pythia6Particle particle) const
private

Return particle momentum.

Definition at line 146 of file G4Pythia6Decayer.cc.

148 {
150 
151  G4ThreeVector momentum
152  = G4ThreeVector(particle->fPx * GeV,
153  particle->fPy * GeV,
154  particle->fPz * GeV);
155  return momentum;
156 }
CLHEP::Hep3Vector G4ThreeVector
static const double GeV
Definition: G4SIunits.hh:214
Here is the caller graph for this function:

◆ GetParticlePosition()

G4ThreeVector G4Pythia6Decayer::GetParticlePosition ( const Pythia6Particle particle) const
private

Return particle vertex position.

Definition at line 132 of file G4Pythia6Decayer.cc.

134 {
136 
138  = G4ThreeVector(particle->fVx * cm,
139  particle->fVy * cm,
140  particle->fVz * cm);
141  return position;
142 }
static const double cm
Definition: G4SIunits.hh:118
CLHEP::Hep3Vector G4ThreeVector
#define position
Definition: xmlparse.cc:622
Here is the caller graph for this function:

◆ ImportDecayProducts()

G4DecayProducts * G4Pythia6Decayer::ImportDecayProducts ( const G4Track &  track)
virtual

Import decay products

Implements G4VExtDecayer.

Definition at line 539 of file G4Pythia6Decayer.cc.

540 {
542 
543  // get particle momentum
544  G4ThreeVector momentum = track.GetMomentum();
545  G4double etot = track.GetDynamicParticle()->GetTotalEnergy();;
547  p[0] = momentum.x() / GeV;
548  p[1] = momentum.y() / GeV;
549  p[2] = momentum.z() / GeV;
550  p[3] = etot / GeV;
551 
552  // get particle PDG
553  // ask G4Pythia6Decayer to get PDG encoding
554  // (in order to get PDG from extended TDatabasePDG
555  // in case the standard PDG code is not defined)
556  G4ParticleDefinition* particleDef = track.GetDefinition();
557  G4int pdgEncoding = particleDef->GetPDGEncoding();
558 
559  // let Pythia6Decayer decay the particle
560  // and import the decay products
561  Decay(pdgEncoding, p);
562  G4int nofParticles = ImportParticles(fDecayProductsArray);
563 
564  if ( fVerboseLevel > 0 ) {
565  G4cout << "nofParticles: " << nofParticles << G4endl;
566  }
567 
568  // convert decay products Pythia6Particle type
569  // to G4DecayProducts
570  G4DecayProducts* decayProducts
571  = new G4DecayProducts(*(track.GetDynamicParticle()));
572 
573  G4int counter = 0;
574  for (G4int i=0; i<nofParticles; i++) {
575 
576  // get particle from ParticleVector
577  Pythia6Particle* particle = (*fDecayProductsArray)[i];
578 
579  G4int status = particle->fKS;
580  G4int pdg = particle->fKF;
581  if ( status>0 && status<11 &&
582  std::abs(pdg)!=12 && std::abs(pdg)!=14 && std::abs(pdg)!=16 ) {
583  // pass to tracking final particles only;
584  // skip neutrinos
585 
586  if ( fVerboseLevel > 0 ) {
587  G4cout << " " << i << "th particle PDG: " << pdg << " ";
588  }
589 
590  // create G4DynamicParticle
591  G4DynamicParticle* dynamicParticle
592  = CreateDynamicParticle(particle);
593 
594  if (dynamicParticle) {
595 
596  if ( fVerboseLevel > 0 ) {
597  G4cout << " G4 particle name: "
598  << dynamicParticle->GetDefinition()->GetParticleName()
599  << G4endl;
600  }
601 
602  // add dynamicParticle to decayProducts
603  decayProducts->PushProducts(dynamicParticle);
604 
605  counter++;
606  }
607  }
608  }
609  if ( fVerboseLevel > 0 ) {
610  G4cout << "nofParticles for tracking: " << counter << G4endl;
611  }
612 
613  return decayProducts;
614 }
void Decay(G4int pdg, const CLHEP::HepLorentzVector &p)
G4int PushProducts(G4DynamicParticle *aParticle)
G4int ImportParticles(ParticleVector *particles)
Structure for Pythia6 particle properties.
Definition: Pythia6.hh:119
ParticleVector * fDecayProductsArray
array of decay products
int G4int
Definition: G4Types.hh:78
G4DynamicParticle * CreateDynamicParticle(const Pythia6Particle *p) const
const G4String & GetParticleName() const
G4GLOB_DLL std::ostream G4cout
static const double GeV
Definition: G4SIunits.hh:214
double x() const
double y() const
double z() const
G4ParticleDefinition * GetDefinition() const
#define G4endl
Definition: G4ios.hh:61
G4int fVerboseLevel
verbose level
double G4double
Definition: G4Types.hh:76
Here is the call graph for this function:

◆ ImportParticles()

G4int G4Pythia6Decayer::ImportParticles ( ParticleVector particles)
private

Get the decay products into the passed PARTICLES vector

Definition at line 526 of file G4Pythia6Decayer.cc.

527 {
529 
530  return Pythia6::Instance()->ImportParticles(particles,"All");
531 }
static Pythia6 * Instance()
Definition: Pythia6.cc:117
ParticleVector * ImportParticles()
Here is the call graph for this function:
Here is the caller graph for this function:

◆ operator=()

G4Pythia6Decayer& G4Pythia6Decayer::operator= ( const G4Pythia6Decayer right)
private

Not implemented.

Here is the caller graph for this function:

◆ SetVerboseLevel()

void G4Pythia6Decayer::SetVerboseLevel ( G4int  verboseLevel)
inline

Definition at line 63 of file G4Pythia6Decayer.hh.

63 { fVerboseLevel = verboseLevel; }
G4int fVerboseLevel
verbose level
Here is the call graph for this function:
Here is the caller graph for this function:

Member Data Documentation

◆ fDecayProductsArray

ParticleVector* G4Pythia6Decayer::fDecayProductsArray
private

array of decay products

Definition at line 95 of file G4Pythia6Decayer.hh.

◆ fDecayType

EDecayType G4Pythia6Decayer::fDecayType
private

selected decay type

Definition at line 94 of file G4Pythia6Decayer.hh.

◆ fgkDefaultDecayType

const EDecayType G4Pythia6Decayer::fgkDefaultDecayType = kAll
staticprivate

default decay type

Definition at line 90 of file G4Pythia6Decayer.hh.

◆ fMessenger

G4Pythia6DecayerMessenger G4Pythia6Decayer::fMessenger
private

command messenger

Definition at line 92 of file G4Pythia6Decayer.hh.

◆ fVerboseLevel

G4int G4Pythia6Decayer::fVerboseLevel
private

verbose level

Definition at line 93 of file G4Pythia6Decayer.hh.


The documentation for this class was generated from the following files: