GENIEGenerator
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Public Member Functions | Private Member Functions | Private Attributes | List of all members
genie::PythiaDecayer Class Reference

Interface to PYTHIA particle decayer.
The PythiaDecayer is a concrete implementation of the Decayer interface. More...

#include <PythiaDecayer.h>

Inheritance diagram for genie::PythiaDecayer:
Inheritance graph
[legend]
Collaboration diagram for genie::PythiaDecayer:
Collaboration graph
[legend]

Public Member Functions

 PythiaDecayer ()
 
 PythiaDecayer (string config)
 
virtual ~PythiaDecayer ()
 
void ProcessEventRecord (GHepRecord *event) const
 

Private Member Functions

void Initialize (void) const
 
bool IsHandled (int pdgc) const
 
void InhibitDecay (int pdgc, TDecayChannel *ch=0) const
 
void UnInhibitDecay (int pdgc, TDecayChannel *ch=0) const
 
bool Decay (int ip, GHepRecord *event) const
 
double SumOfBranchingRatios (int kc) const
 
int FindPythiaDecayChannel (int kc, TDecayChannel *ch) const
 
bool MatchDecayChannels (int ic, TDecayChannel *ch) const
 

Private Attributes

TPythia6 * fPythia
 PYTHIA6 wrapper class. More...
 
double fWeight
 

Additional Inherited Members

- Protected Member Functions inherited from genie::Decayer
 Decayer ()
 
 Decayer (string name)
 
 Decayer (string name, string config)
 
virtual void LoadConfig (void)
 
virtual bool ToBeDecayed (int pdgc, GHepStatus_t ist) const
 
virtual bool IsUnstable (int pdgc) const
 
virtual ~Decayer ()
 
void Configure (const Registry &config)
 
void Configure (string config)
 
- Protected Member Functions inherited from genie::EventRecordVisitorI
 EventRecordVisitorI ()
 
 EventRecordVisitorI (string name)
 
 EventRecordVisitorI (string name, string config)
 
virtual ~EventRecordVisitorI ()
 
- Protected Member Functions inherited from genie::Algorithm
 Algorithm ()
 
 Algorithm (string name)
 
 Algorithm (string name, string config)
 
void Initialize (void)
 
void DeleteConfig (void)
 
void DeleteSubstructure (void)
 
RegistryExtractLocalConfig (const Registry &in) const
 
RegistryExtractLowerConfig (const Registry &in, const string &alg_key) const
 Split an incoming configuration Registry into a block valid for the sub-algo identified by alg_key. More...
 
template<class T >
bool GetParam (const RgKey &name, T &p, bool is_top_call=true) const
 
template<class T >
bool GetParamDef (const RgKey &name, T &p, const T &def) const
 
template<class T >
int GetParamVect (const std::string &comm_name, std::vector< T > &v, bool is_top_call=true) const
 Handle to load vectors of parameters. More...
 
int GetParamVectKeys (const std::string &comm_name, std::vector< RgKey > &k, bool is_top_call=true) const
 
template<class T >
int GetParamMat (const std::string &comm_name, TMatrixT< T > &mat, bool is_top_call=true) const
 Handle to load matrix of parameters. More...
 
template<class T >
int GetParamMatSym (const std::string &comm_name, TMatrixTSym< T > &mat, bool is_top_call=true) const
 
int GetParamMatKeys (const std::string &comm_name, std::vector< RgKey > &k, bool is_top_call=true) const
 
int AddTopRegistry (Registry *rp, bool owns=true)
 add registry with top priority, also update ownership More...
 
int AddLowRegistry (Registry *rp, bool owns=true)
 add registry with lowest priority, also update ownership More...
 
int MergeTopRegistry (const Registry &r)
 
int AddTopRegisties (const vector< Registry * > &rs, bool owns=false)
 Add registries with top priority, also udated Ownerships. More...
 
virtual ~Algorithm ()
 
virtual void FindConfig (void)
 
virtual const RegistryGetConfig (void) const
 
RegistryGetOwnedConfig (void)
 
virtual const AlgIdId (void) const
 Get algorithm ID. More...
 
virtual AlgStatus_t GetStatus (void) const
 Get algorithm status. More...
 
virtual bool AllowReconfig (void) const
 
virtual AlgCmp_t Compare (const Algorithm *alg) const
 Compare with input algorithm. More...
 
virtual void SetId (const AlgId &id)
 Set algorithm ID. More...
 
virtual void SetId (string name, string config)
 
const AlgorithmSubAlg (const RgKey &registry_key) const
 
void AdoptConfig (void)
 
void AdoptSubstructure (void)
 
virtual void Print (ostream &stream) const
 Print algorithm info. More...
 
- Static Protected Member Functions inherited from genie::Algorithm
static string BuildParamVectKey (const std::string &comm_name, unsigned int i)
 
static string BuildParamVectSizeKey (const std::string &comm_name)
 
static string BuildParamMatKey (const std::string &comm_name, unsigned int i, unsigned int j)
 
static string BuildParamMatRowSizeKey (const std::string &comm_name)
 
static string BuildParamMatColSizeKey (const std::string &comm_name)
 
- Protected Attributes inherited from genie::Decayer
bool fGenerateWeighted
 generate weighted or unweighted decays? More...
 
bool fRunBefHadroTransp
 is invoked before or after FSI? More...
 
PDGCodeList fParticlesToDecay
 list of particles to be decayed More...
 
PDGCodeList fParticlesNotToDecay
 list of particles for which decay is inhibited More...
 
- Protected Attributes inherited from genie::Algorithm
bool fAllowReconfig
 
bool fOwnsSubstruc
 true if it owns its substructure (sub-algs,...) More...
 
AlgId fID
 algorithm name and configuration set More...
 
vector< Registry * > fConfVect
 
vector< bool > fOwnerships
 ownership for every registry in fConfVect More...
 
AlgStatus_t fStatus
 algorithm execution status More...
 
AlgMapfOwnedSubAlgMp
 local pool for owned sub-algs (taken out of the factory pool) More...
 

Detailed Description

Interface to PYTHIA particle decayer.
The PythiaDecayer is a concrete implementation of the Decayer interface.

Author
Costas Andreopoulos <c.andreopoulos cern.ch> University of Liverpool
Created:
June 20, 2004
License:
Copyright (c) 2003-2024, The GENIE Collaboration For the full text of the license visit http://copyright.genie-mc.org

Definition at line 30 of file PythiaDecayer.h.

Constructor & Destructor Documentation

PythiaDecayer::PythiaDecayer ( )

Definition at line 45 of file PythiaDecayer.cxx.

References Initialize().

45  :
46 Decayer("genie::PythiaDecayer")
47 {
48  this->Initialize();
49 }
void Initialize(void) const
PythiaDecayer::PythiaDecayer ( string  config)

Definition at line 51 of file PythiaDecayer.cxx.

References Initialize().

51  :
52 Decayer("genie::PythiaDecayer", config)
53 {
54  this->Initialize();
55 }
void Initialize(void) const
PythiaDecayer::~PythiaDecayer ( )
virtual

Definition at line 57 of file PythiaDecayer.cxx.

58 {
59 
60 }

Member Function Documentation

bool PythiaDecayer::Decay ( int  ip,
GHepRecord event 
) const
private

Definition at line 95 of file PythiaDecayer.cxx.

References genie::GHepParticle::Energy(), genie::units::fm, fPythia, genie::Decayer::fRunBefHadroTransp, fWeight, genie::PDGLibrary::Instance(), genie::pdg::IsHadron(), genie::kIStDecayedState, genie::kIStHadronInTheNucleus, genie::kIStNucleonTarget, genie::kIStStableFinalState, LOG, genie::units::m, genie::GHepParticle::Mass(), genie::units::mm, genie::GHepParticle::P4(), genie::GHepParticle::Pdg(), pINFO, pNOTICE, pWARN, genie::GHepParticle::Px(), genie::GHepParticle::Py(), py1ent_(), genie::GHepParticle::Pz(), genie::units::s, genie::GHepParticle::SetStatus(), SLOG, genie::GHepParticle::Status(), SumOfBranchingRatios(), and genie::GHepParticle::X4().

Referenced by ProcessEventRecord().

96 {
97  fWeight = 1.; // reset previous decay weight
98 
99  // Get particle to be decayed
100  GHepParticle * decay_particle = event->Particle(decay_particle_id);
101  if(!decay_particle) return 0;
102 
103  // Get the particle 4-momentum, 4-position and PDG code
104  TLorentzVector decay_particle_p4 = *(decay_particle->P4());
105  int decay_particle_pdg_code = decay_particle->Pdg();
106 
107  // Convert to PYTHIA6 particle code and check whether decay is inhibited
108  int kc = fPythia->Pycomp(decay_particle_pdg_code);
109  int mdcy = fPythia->GetMDCY(kc, 1);
110  if(mdcy == 0) {
111  LOG("Pythia6Decay", pNOTICE)
112  << (PDGLibrary::Instance())->Find(decay_particle_pdg_code)->GetName()
113  << " decays are inhibited!";
114  return false;
115  }
116 
117  // Get sub of BRs and compute weight if decay channels were inhibited
118  double sumbr = this->SumOfBranchingRatios(kc);
119  if(sumbr <= 0) {
120  LOG("Pythia6Decay", pWARN)
121  << "The sum of enabled "
122  << (PDGLibrary::Instance())->Find(decay_particle_pdg_code)->GetName()
123  << " decay channel branching ratios is non-positive!";
124  return false;
125  }
126  fWeight = 1./sumbr; // update weight to account for inhibited channels
127 
128  // Run PYTHIA6 decay
129  int ip = 0;
130  double E = decay_particle_p4.Energy();
131  double theta = decay_particle_p4.Theta();
132  double phi = decay_particle_p4.Phi();
133  fPythia->SetMSTJ(22,1);
134  py1ent_(&ip, &decay_particle_pdg_code, &E, &theta, &phi);
135 
136  // Get decay products
137  fPythia->GetPrimaries();
138  TClonesArray * impl = (TClonesArray *) fPythia->ImportParticles("All");
139  if(!impl) {
140  LOG("Pythia6Decay", pWARN)
141  << "TPythia6::ImportParticles() returns a null list!";
142  return false;
143  }
144 
145  // Copy the PYTHIA6 container to the GENIE event record
146 
147  // Check whether the interaction is off a nuclear target or free nucleon
148  // Depending on whether this module is run before or after the hadron
149  // transport module it would affect the daughters status code
150  GHepParticle * target_nucleus = event->TargetNucleus();
151  bool in_nucleus = (target_nucleus!=0);
152 
153  // the values of space coordinates from pythia are in mm.
154  // our conventions want it in fm
155  constexpr double space_scale = units::mm / units::fm ;
156 
157  // the values of time coordinate from pythia is in mm/c.
158  // our conventions want it in ys
159  constexpr double time_scale = 1e21 * units::m / units::s ;
160 
161  TMCParticle * p = 0;
162  TIter particle_iter(impl);
163  while( (p = (TMCParticle *) particle_iter.Next()) ) {
164  // Convert from TMCParticle to GHepParticle
166  p->GetKF(), // pdg
167  GHepStatus_t(p->GetKS()), // status
168  p->GetParent(), // first parent
169  0, // second parent
170  p->GetFirstChild(), // first daughter
171  p->GetLastChild(), // second daughter
172  p->GetPx(), // px
173  p->GetPy(), // py
174  p->GetPz(), // pz
175  p->GetEnergy(), // e
176  p->GetVx() * space_scale , // x
177  p->GetVy() * space_scale , // y
178  p->GetVz() * space_scale , // z
179  p->GetTime() * time_scale // t
180  );
181 
182  if(mcp.Status()==kIStNucleonTarget) continue; // mother particle, already in GHEP
183 
184  int daughter_pdg_code = mcp.Pdg();
185  SLOG("Pythia6Decay", pINFO)
186  << "Adding daughter particle wit PDG code = "
187  << daughter_pdg_code << ", m = " << mcp.Mass()
188  << " GeV, E = " << mcp.Energy() << " GeV)";
189 
190  bool is_hadron = pdg::IsHadron(daughter_pdg_code);
191  bool hadron_in_nuc = (in_nucleus && is_hadron && fRunBefHadroTransp);
192 
193  GHepStatus_t daughter_status_code = (hadron_in_nuc) ?
195 
196  TLorentzVector daughter_p4(
197  mcp.Px(),mcp.Py(),mcp.Pz(),mcp.Energy());
198 
199  event->AddParticle(
200  daughter_pdg_code, daughter_status_code,
201  decay_particle_id,-1,-1,-1,
202  daughter_p4, * mcp.X4() );
203  }
204 
205  // Update the event weight for each weighted particle decay
206  double weight = event->Weight() * fWeight;
207  event->SetWeight(weight);
208 
209  // Mark input particle as a 'decayed state' & add its daughter links
210  decay_particle->SetStatus(kIStDecayedState);
211 
212  return true;
213 }
const TLorentzVector * P4(void) const
Definition: GHepParticle.h:78
static constexpr double s
Definition: Units.h:95
double Mass(void) const
Mass that corresponds to the PDG code.
double Pz(void) const
Get Pz.
Definition: GHepParticle.h:90
GHepStatus_t Status(void) const
Definition: GHepParticle.h:64
double Energy(void) const
Get energy.
Definition: GHepParticle.h:92
double Px(void) const
Get Px.
Definition: GHepParticle.h:88
int Pdg(void) const
Definition: GHepParticle.h:63
bool IsHadron(int pdgc)
Definition: PDGUtils.cxx:392
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
bool fRunBefHadroTransp
is invoked before or after FSI?
Definition: Decayer.h:57
TPythia6 * fPythia
PYTHIA6 wrapper class.
Definition: PythiaDecayer.h:51
#define pINFO
Definition: Messenger.h:62
#define pWARN
Definition: Messenger.h:60
double SumOfBranchingRatios(int kc) const
void py1ent_(int *, int *, double *, double *, double *)
static PDGLibrary * Instance(void)
Definition: PDGLibrary.cxx:68
const TLorentzVector * X4(void) const
Definition: GHepParticle.h:79
void SetStatus(GHepStatus_t s)
Definition: GHepParticle.h:126
static constexpr double mm
Definition: Units.h:65
static constexpr double fm
Definition: Units.h:75
#define pNOTICE
Definition: Messenger.h:61
#define SLOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a short string (using the FUNCTION and...
Definition: Messenger.h:84
static constexpr double m
Definition: Units.h:71
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:39
enum genie::EGHepStatus GHepStatus_t
double Py(void) const
Get Py.
Definition: GHepParticle.h:89
int PythiaDecayer::FindPythiaDecayChannel ( int  kc,
TDecayChannel *  ch 
) const
private

Definition at line 321 of file PythiaDecayer.cxx.

References fPythia, LOG, MatchDecayChannels(), pINFO, pNOTICE, and pWARN.

Referenced by InhibitDecay(), and UnInhibitDecay().

322 {
323  if(!dc) return -1;
324 
325  int first_channel = fPythia->GetMDCY(kc,2);
326  int last_channel = fPythia->GetMDCY(kc,2) + fPythia->GetMDCY(kc,3) - 1;
327 
328  bool found_match = false;
329 
330  // loop over pythia decay channels
331  for(int ichannel = first_channel;
332  ichannel < last_channel; ichannel++) {
333 
334  // does the current pythia channel matches the input TDecayChannel?
335  LOG("Pythia6Decay", pINFO)
336  << "\nComparing PYTHIA's channel = " << ichannel
337  << " with TDecayChannel = " << dc->Number();
338 
339  found_match = this->MatchDecayChannels(ichannel,dc);
340  if(found_match) {
341  LOG("Pythia6Decay", pNOTICE)
342  << " ** TDecayChannel id = " << dc->Number()
343  << " corresponds to PYTHIA6 channel id = " << ichannel;
344  return ichannel;
345  }//match?
346  }//loop pythia decay ch.
347 
348  LOG("Pythia6Decay", pWARN)
349  << " ** No PYTHIA6 decay channel match found for "
350  << "TDecayChannel id = " << dc->Number();
351 
352  return -1;
353 }
bool MatchDecayChannels(int ic, TDecayChannel *ch) const
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
TPythia6 * fPythia
PYTHIA6 wrapper class.
Definition: PythiaDecayer.h:51
#define pINFO
Definition: Messenger.h:62
#define pWARN
Definition: Messenger.h:60
#define pNOTICE
Definition: Messenger.h:61
void PythiaDecayer::InhibitDecay ( int  pdgc,
TDecayChannel *  ch = 0 
) const
privatevirtual

Implements genie::Decayer.

Definition at line 237 of file PythiaDecayer.cxx.

References FindPythiaDecayChannel(), fPythia, IsHandled(), LOG, and pINFO.

238 {
239  if(! this->IsHandled(pdg_code)) return;
240 
241  int kc = fPythia->Pycomp(pdg_code);
242 
243  if(!dc) {
244  LOG("Pythia6Decay", pINFO)
245  << "Switching OFF ALL decay channels for particle = " << pdg_code;
246  fPythia->SetMDCY(kc, 1,0);
247  return;
248  }
249 
250  LOG("Pythia6Decay", pINFO)
251  << "Switching OFF decay channel = " << dc->Number()
252  << " for particle = " << pdg_code;
253 
254  int ichannel = this->FindPythiaDecayChannel(kc, dc);
255  if(ichannel != -1) {
256  fPythia->SetMDME(ichannel,1,0); // switch-off
257  }
258 }
bool IsHandled(int pdgc) const
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
TPythia6 * fPythia
PYTHIA6 wrapper class.
Definition: PythiaDecayer.h:51
#define pINFO
Definition: Messenger.h:62
int FindPythiaDecayChannel(int kc, TDecayChannel *ch) const
void PythiaDecayer::Initialize ( void  ) const
private

Definition at line 215 of file PythiaDecayer.cxx.

References fPythia, fWeight, and genie::RandomGen::Instance().

Referenced by PythiaDecayer().

216 {
217  fPythia = TPythia6::Instance();
218  fWeight = 1.;
219 
220  // sync GENIE/PYTHIA6 seeds
222 }
static RandomGen * Instance()
Access instance.
Definition: RandomGen.cxx:71
TPythia6 * fPythia
PYTHIA6 wrapper class.
Definition: PythiaDecayer.h:51
bool PythiaDecayer::IsHandled ( int  pdgc) const
privatevirtual

Implements genie::Decayer.

Definition at line 224 of file PythiaDecayer.cxx.

References genie::utils::res::IsBaryonResonance(), LOG, and pDEBUG.

Referenced by InhibitDecay(), ProcessEventRecord(), and UnInhibitDecay().

225 {
226 // does not handle requests to decay baryon resonances
227 
228  bool is_handled = (!utils::res::IsBaryonResonance(pdg_code));
229 
230  LOG("Pythia6Decay", pDEBUG)
231  << "Can decay particle with PDG code = " << pdg_code
232  << "? " << ((is_handled)? "Yes" : "No");
233 
234  return is_handled;
235 }
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
bool IsBaryonResonance(int pdgc)
is input a baryon resonance?
#define pDEBUG
Definition: Messenger.h:63
bool PythiaDecayer::MatchDecayChannels ( int  ic,
TDecayChannel *  ch 
) const
private

Definition at line 355 of file PythiaDecayer.cxx.

References fPythia, LOG, and pDEBUG.

Referenced by FindPythiaDecayChannel().

356 {
357  // num. of daughters in the input TDecayChannel & the input PYTHIA ichannel
358  int nd = dc->NDaughters();
359 
360  int py_nd = 0;
361  for (int i = 1; i <= 5; i++) {
362  if(fPythia->GetKFDP(ichannel,i) != 0) py_nd++;
363  }
364 
365  LOG("Pythia6Decay", pDEBUG)
366  << "NDaughters: PYTHIA = " << py_nd << ", ROOT's TDecayChannel = " << nd;
367 
368  if(nd != py_nd) return false;
369 
370  //
371  // if the two channels have the same num. of daughters, then compare them
372  //
373 
374  // store decay daughters for the input TDecayChannel
375  vector<int> dc_daughter(nd);
376  int k=0;
377  for( ; k < nd; k++) {
378  dc_daughter[k] = dc->DaughterPdgCode(k);
379  }
380  // store decay daughters for the input PYTHIA's ichannel
381  vector<int> py_daughter(nd);
382  k=0;
383  for(int i = 1; i <= 5; i++) {
384  if(fPythia->GetKFDP(ichannel,i) == 0) continue;
385  py_daughter[k] = fPythia->GetKFDP(ichannel,i);
386  k++;
387  }
388 
389  // sort both daughter lists
390  sort( dc_daughter.begin(), dc_daughter.end() );
391  sort( py_daughter.begin(), py_daughter.end() );
392 
393  // compare
394  for(int i = 0; i < nd; i++) {
395  if(dc_daughter[i] != py_daughter[i]) return false;
396  }
397  return true;
398 }
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
TPythia6 * fPythia
PYTHIA6 wrapper class.
Definition: PythiaDecayer.h:51
#define pDEBUG
Definition: Messenger.h:63
void PythiaDecayer::ProcessEventRecord ( GHepRecord event) const
virtual

Implements genie::EventRecordVisitorI.

Definition at line 62 of file PythiaDecayer.cxx.

References Decay(), genie::Decayer::fRunBefHadroTransp, IsHandled(), LOG, genie::GHepParticle::Name(), pDEBUG, genie::GHepParticle::Pdg(), pINFO, pNOTICE, genie::GHepParticle::Status(), and genie::Decayer::ToBeDecayed().

63 {
64  LOG("ResonanceDecay", pINFO)
65  << "Running PYTHIA6 particle decayer "
66  << ((fRunBefHadroTransp) ? "*before*" : "*after*") << " FSI";
67 
68  // Loop over particles, find unstable ones and decay them
69  TObjArrayIter piter(event);
70  GHepParticle * p = 0;
71  int ipos = -1;
72 
73  while( (p = (GHepParticle *) piter.Next()) ) {
74  ipos++;
75 
76  LOG("Pythia6Decay", pDEBUG) << "Checking: " << p->Name();
77 
78  int pdg_code = p->Pdg();
79  GHepStatus_t status_code = p->Status();
80 
81  if(!this->IsHandled (pdg_code)) continue;
82  if(!this->ToBeDecayed(pdg_code, status_code)) continue;
83 
84  LOG("Pythia6Decay", pNOTICE)
85  << "Decaying unstable particle: " << p->Name();
86 
87  bool decayed = this->Decay(ipos, event);
88  assert(decayed); // handle this more graciously and throw an exception
89  }
90 
91  LOG("Pythia6Decay", pNOTICE)
92  << "Done finding & decaying unstable particles";
93 }
virtual bool ToBeDecayed(int pdgc, GHepStatus_t ist) const
Definition: Decayer.cxx:51
bool IsHandled(int pdgc) const
GHepStatus_t Status(void) const
Definition: GHepParticle.h:64
int Pdg(void) const
Definition: GHepParticle.h:63
string Name(void) const
Name that corresponds to the PDG code.
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
bool fRunBefHadroTransp
is invoked before or after FSI?
Definition: Decayer.h:57
bool Decay(int ip, GHepRecord *event) const
#define pINFO
Definition: Messenger.h:62
#define pNOTICE
Definition: Messenger.h:61
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:39
enum genie::EGHepStatus GHepStatus_t
#define pDEBUG
Definition: Messenger.h:63
double PythiaDecayer::SumOfBranchingRatios ( int  kc) const
private

Definition at line 292 of file PythiaDecayer.cxx.

References fPythia, LOG, and pINFO.

Referenced by Decay().

293 {
294 // Sum of branching ratios for enabled channels
295 //
296  double sumbr=0.;
297 
298  int first_channel = fPythia->GetMDCY(kc,2);
299  int last_channel = fPythia->GetMDCY(kc,2) + fPythia->GetMDCY(kc,3) - 1;
300 
301  bool has_inhibited_channels=false;
302 
303  // loop over pythia decay channels
304  for(int ichannel = first_channel;
305  ichannel < last_channel; ichannel++) {
306 
307  bool enabled = (fPythia->GetMDME(ichannel,1) == 1);
308  if (!enabled) {
309  has_inhibited_channels = true;
310  } else {
311  sumbr += fPythia->GetBRAT(ichannel);
312  }
313  }
314 
315  if(!has_inhibited_channels) return 1.;
316  LOG("Pythia6Decay", pINFO) << "Sum{BR} = " << sumbr;
317 
318  return sumbr;
319 }
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
TPythia6 * fPythia
PYTHIA6 wrapper class.
Definition: PythiaDecayer.h:51
#define pINFO
Definition: Messenger.h:62
void PythiaDecayer::UnInhibitDecay ( int  pdgc,
TDecayChannel *  ch = 0 
) const
privatevirtual

dc

Implements genie::Decayer.

Definition at line 260 of file PythiaDecayer.cxx.

References FindPythiaDecayChannel(), fPythia, IsHandled(), LOG, and pINFO.

261 {
262  if(! this->IsHandled(pdg_code)) return;
263 
264  int kc = fPythia->Pycomp(pdg_code);
265 
266  if(!dc) {
267  LOG("Pythia6Decay", pINFO)
268  << "Switching ON all PYTHIA decay channels for particle = " << pdg_code;
269 
270  fPythia->SetMDCY(kc, 1,1);
271 
272  int first_channel = fPythia->GetMDCY(kc,2);
273  int last_channel = fPythia->GetMDCY(kc,2) + fPythia->GetMDCY(kc,3) - 1;
274 
275  for(int ichannel = first_channel;
276  ichannel < last_channel; ichannel++) {
277  fPythia->SetMDME(ichannel,1,1); // switch-on
278  }
279  return;
280  }//!dc
281 
282  LOG("Pythia6Decay", pINFO)
283  << "Switching OFF decay channel = " << dc->Number()
284  << " for particle = " << pdg_code;
285 
286  int ichannel = this->FindPythiaDecayChannel(kc, dc);
287  if(ichannel != -1) {
288  fPythia->SetMDME(ichannel,1,1); // switch-on
289  }
290 }
bool IsHandled(int pdgc) const
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
TPythia6 * fPythia
PYTHIA6 wrapper class.
Definition: PythiaDecayer.h:51
#define pINFO
Definition: Messenger.h:62
int FindPythiaDecayChannel(int kc, TDecayChannel *ch) const

Member Data Documentation

TPythia6* genie::PythiaDecayer::fPythia
mutableprivate
double genie::PythiaDecayer::fWeight
mutableprivate

Definition at line 52 of file PythiaDecayer.h.

Referenced by Decay(), and Initialize().


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