Public Member Functions | Protected Member Functions | Protected Attributes | Private Member Functions | Private Attributes | Friends | List of all members
genie::HAIntranuke Class Reference

#include "/cvmfs/nova.opensciencegrid.org/externals/genie/v3_00_06_p01/Linux64bit+2.6-2.12-e17-debug/GENIE-Generator/src/Physics/HadronTransport/HAIntranuke.h"

Inheritance diagram for genie::HAIntranuke:
genie::Intranuke genie::EventRecordVisitorI genie::Algorithm

Public Member Functions

 HAIntranuke ()
 
 HAIntranuke (string config)
 
 ~HAIntranuke ()
 
void ProcessEventRecord (GHepRecord *event_rec) const
 
virtual void Configure (string param_set)
 
void Configure (const Registry &config)
 
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...
 

Protected Member Functions

void TransportHadrons (GHepRecord *ev) const
 
void GenerateVertex (GHepRecord *ev) const
 
bool NeedsRescattering (const GHepParticle *p) const
 
bool CanRescatter (const GHepParticle *p) const
 
bool IsInNucleus (const GHepParticle *p) const
 
void SetTrackingRadius (const GHepParticle *p) const
 
double GenerateStep (GHepRecord *ev, GHepParticle *p) const
 
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 >
bool GetParamVect (const std::string &comm_name, std::vector< T > &v, unsigned int max, 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...
 

Protected Attributes

double fTrackingRadius
 tracking radius for the nucleus in the current event More...
 
TGenPhaseSpace fGenPhaseSpace
 a phase space generator More...
 
INukeHadroDatafHadroData
 a collection of h+N,h+A data & calculations More...
 
AlgFactoryfAlgf
 algorithm factory instance More...
 
const NuclearModelIfNuclmodel
 nuclear model used to generate fermi momentum More...
 
int fRemnA
 remnant nucleus A More...
 
int fRemnZ
 remnant nucleus Z More...
 
TLorentzVector fRemnP4
 P4 of remnant system. More...
 
GEvGenMode_t fGMode
 event generation mode (lepton+A, hadron+A, ...) More...
 
double fR0
 effective nuclear size param More...
 
double fNR
 param multiplying the nuclear radius, determining how far to track hadrons beyond the "nuclear boundary" More...
 
double fNucRmvE
 binding energy to subtract from cascade nucleons More...
 
double fDelRPion
 factor by which Pion Compton wavelength gets multiplied to become nuclear size enhancement More...
 
double fDelRNucleon
 factor by which Nucleon Compton wavelength gets multiplied to become nuclear size enhancement More...
 
double fHadStep
 step size for intranuclear hadron transport More...
 
double fNucAbsFac
 absorption xsec correction factor (hN Mode) More...
 
double fNucCEXFac
 charge exchange xsec correction factor (hN Mode) More...
 
double fEPreEq
 threshold for pre-equilibrium reaction More...
 
double fFermiFac
 testing parameter to modify fermi momentum More...
 
double fFermiMomentum
 whether or not particle collision is pauli blocked More...
 
bool fDoFermi
 whether or not to do fermi mom. More...
 
bool fDoMassDiff
 whether or not to do mass diff. mode More...
 
bool fDoCompoundNucleus
 whether or not to do compound nucleus considerations More...
 
double fPionMFPScale
 tweaking factors for tuning More...
 
double fPionFracCExScale
 
double fPionFracElasScale
 
double fPionFracInelScale
 
double fPionFracAbsScale
 
double fPionFracPiProdScale
 
double fNucleonMFPScale
 
double fNucleonFracCExScale
 
double fNucleonFracElasScale
 
double fNucleonFracInelScale
 
double fNucleonFracAbsScale
 
double fNucleonFracPiProdScale
 
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...
 

Private Member Functions

void LoadConfig (void)
 
void SimulateHadronicFinalState (GHepRecord *ev, GHepParticle *p) const
 
void SimulateHadronicFinalStateKinematics (GHepRecord *ev, GHepParticle *p) const
 
INukeFateHA_t HadronFateHA (const GHepParticle *p) const
 
void Inelastic (GHepRecord *ev, GHepParticle *p, INukeFateHA_t fate) const
 
void ElasHA (GHepRecord *ev, GHepParticle *p, INukeFateHA_t fate) const
 
void InelasticHA (GHepRecord *ev, GHepParticle *p, INukeFateHA_t fate) const
 
double PiBounce (void) const
 
double PnBounce (void) const
 
bool HandleCompoundNucleus (GHepRecord *ev, GHepParticle *p, int mom) const
 

Private Attributes

unsigned int fNumIterations
 

Friends

class IntranukeTester
 

Detailed Description

Definition at line 51 of file HAIntranuke.h.

Constructor & Destructor Documentation

HAIntranuke::HAIntranuke ( )

Definition at line 91 of file HAIntranuke.cxx.

91  :
92 Intranuke("genie::HAIntranuke")
93 {
94 
95 }
HAIntranuke::HAIntranuke ( string  config)

Definition at line 97 of file HAIntranuke.cxx.

97  :
98 Intranuke("genie::HAIntranuke",config)
99 {
100 
101 }
Definition: config.py:1
HAIntranuke::~HAIntranuke ( )

Definition at line 103 of file HAIntranuke.cxx.

104 {
105 
106 }

Member Function Documentation

int Algorithm::AddLowRegistry ( Registry rp,
bool  owns = true 
)
protectedinherited

add registry with lowest priority, also update ownership

Definition at line 601 of file Algorithm.cxx.

Referenced by genie::EventGenerator::Configure().

601  {
602 
603  fConfVect.push_back( rp ) ;
604  fOwnerships.push_back( own ) ;
605 
606  if ( fConfig ) {
607  delete fConfig ;
608  fConfig = 0 ;
609  }
610 
611  return fConfVect.size() ;
612 
613 }
vector< Registry * > fConfVect
Definition: Algorithm.h:161
vector< bool > fOwnerships
ownership for every registry in fConfVect
Definition: Algorithm.h:164
Registry * fConfig
Summary configuration derived from fConvVect, not necessarily allocated.
Definition: Algorithm.h:194
int Algorithm::AddTopRegisties ( const vector< Registry * > &  rs,
bool  owns = false 
)
protectedinherited

Add registries with top priority, also udated Ownerships.

Definition at line 653 of file Algorithm.cxx.

653  {
654 
655  fConfVect.insert( fConfVect.begin(), rs.begin(), rs.end() ) ;
656 
657  fOwnerships.insert( fOwnerships.begin(), rs.size(), own ) ;
658 
659  if ( fConfig ) {
660  delete fConfig ;
661  fConfig = 0 ;
662  }
663 
664  return fConfVect.size() ;
665 
666 }
vector< Registry * > fConfVect
Definition: Algorithm.h:161
vector< bool > fOwnerships
ownership for every registry in fConfVect
Definition: Algorithm.h:164
Registry * fConfig
Summary configuration derived from fConvVect, not necessarily allocated.
Definition: Algorithm.h:194
int Algorithm::AddTopRegistry ( Registry rp,
bool  owns = true 
)
protectedinherited

add registry with top priority, also update ownership

Definition at line 585 of file Algorithm.cxx.

Referenced by genie::EventGeneratorListAssembler::AssembleGeneratorList().

585  {
586 
587  fConfVect.insert( fConfVect.begin(), rp ) ;
588  fOwnerships.insert( fOwnerships.begin(), own ) ;
589 
590  if ( fConfig ) {
591  delete fConfig ;
592  fConfig = 0 ;
593  }
594 
595  return fConfVect.size() ;
596 
597 }
vector< Registry * > fConfVect
Definition: Algorithm.h:161
vector< bool > fOwnerships
ownership for every registry in fConfVect
Definition: Algorithm.h:164
Registry * fConfig
Summary configuration derived from fConvVect, not necessarily allocated.
Definition: Algorithm.h:194
void Algorithm::AdoptConfig ( void  )
inherited

Clone the configuration registry looked up from the configuration pool and take its ownership

Definition at line 394 of file Algorithm.cxx.

References Configure(), GetConfig(), LOG, and pNOTICE.

Referenced by genie::Algorithm::AllowReconfig().

394  {
395 
396  LOG("Algorithm", pNOTICE)
397  << this->Id().Key() << " is taking ownership of its configuration";
398 
399  // if(fOwnsConfig) {
400  // LOG("Algorithm", pWARN)
401  // << this->Id().Key() << " already owns its configuration!";
402  // return;
403  // }
404 
405  this->Configure( GetConfig() );
406 }
virtual const Registry & GetConfig(void) const
Definition: Algorithm.cxx:254
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
virtual void Configure(const Registry &config)
Definition: Algorithm.cxx:70
virtual const AlgId & Id(void) const
Get algorithm ID.
Definition: Algorithm.h:98
#define pNOTICE
Definition: Messenger.h:62
string Key(void) const
Definition: AlgId.h:47
void Algorithm::AdoptSubstructure ( void  )
inherited

Take ownership of the algorithms subtructure (sub-algorithms,...) by copying them from the AlgFactory pool to the local pool Also bring all the configuration variables to the top level config Registry. This can be used to group together a series of algorithms & their configurations and extract (a clone of) this group from the shared pools. Having a series of algorithms/configurations behaving as a monolithic block, with a single point of configuration (the top level) is to be used when bits & pieces of GENIE are used in isolation for data fitting or reweighting

Definition at line 408 of file Algorithm.cxx.

References genie::AlgFactory::AdoptAlgorithm(), genie::Algorithm::AdoptSubstructure(), GetConfig(), genie::AlgFactory::Instance(), genie::kRgAlg, LOG, pDEBUG, pNOTICE, and genie::RegistryItemI::TypeInfo().

Referenced by genie::Algorithm::AdoptSubstructure(), genie::Algorithm::AllowReconfig(), main(), and testReconfigInOwnedModules().

409 {
410 // Take ownership of the algorithms subtructure (sub-algorithms,..) by copying
411 // them from the AlgFactory pool to the local pool. Also bring all the
412 // configuration variables to the top level. See the header for more details.
413 // A secial naming convention is required for configuration parameter keys
414 // for parameters belonging to sub-algorithms (or sub-algorithms of these
415 // sub-algorithms and so on...).
416 // The convention is: "sub-alg-key/sub-sub-alg-key/.../original name"
417 // This is a recursive method: The AdoptSubtructure()of all sub-algorithms is
418 // invoked.
419 //
420  LOG("Algorithm", pNOTICE)
421  << "Algorithm: " << this->Id().Key() << " is adopting its substructure";
422 
423 // Registry deep_config;
424 // deep_config.UnLock();
425 // deep_config.SetName(this->Id().Key());
426 
427  // deep_config.SetName(this->Id().Config() + ";D");
428  // fID.SetConfig(this->Id().Config() + ";D");
429 
430  if(fOwnsSubstruc) this->DeleteSubstructure();
431 
432  fOwnedSubAlgMp = new AlgMap;
433  fOwnsSubstruc = true;
434 
435  AlgFactory * algf = AlgFactory::Instance();
436 
437  const RgIMap & rgmap = GetConfig().GetItemMap();
438 
439  RgIMapConstIter iter = rgmap.begin();
440  for( ; iter != rgmap.end(); ++iter) {
441 
442  RgKey reg_key = iter->first;
443  RegistryItemI * ri = iter->second;
444 
445  if(ri->TypeInfo() == kRgAlg) {
446  LOG("Algorithm", pDEBUG)
447  << "Found sub-algorithm pointed to by " << reg_key;
448  RgAlg reg_alg = fConfig->GetAlg(reg_key);
449  AlgId id(reg_alg);
450 
451  LOG("Algorithm", pDEBUG) << "Adopting sub-algorithm = " << id.Key();
452  Algorithm * subalg = algf->AdoptAlgorithm(id.Name(),id.Config());
453  subalg->AdoptSubstructure();
454 
455  LOG("Algorithm", pDEBUG) << "Adding sub-algorithm to local pool";
456  AlgMapPair key_alg_pair(reg_key, subalg);
457  fOwnedSubAlgMp->insert(key_alg_pair);
458 
459  }
460 
461  }
462 
463 
464  if ( fConfig ) {
465  delete fConfig ;
466  fConfig = 0 ;
467  }
468 
469 }
::xsd::cxx::tree::id< char, ncname > id
Definition: Database.h:165
void DeleteSubstructure(void)
Definition: Algorithm.cxx:496
AlgMap * fOwnedSubAlgMp
local pool for owned sub-algs (taken out of the factory pool)
Definition: Algorithm.h:167
bool fOwnsSubstruc
true if it owns its substructure (sub-algs,...)
Definition: Algorithm.h:155
Algorithm abstract base class.
Definition: Algorithm.h:54
map< string, Algorithm * > AlgMap
Definition: Algorithm.h:49
Registry item pABC.
Definition: RegistryItemI.h:30
virtual const Registry & GetConfig(void) const
Definition: Algorithm.cxx:254
virtual RgType_t TypeInfo(void) const =0
map< RgKey, RegistryItemI * >::const_iterator RgIMapConstIter
Definition: Registry.h:50
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
const RgIMap & GetItemMap(void) const
Definition: Registry.h:162
void AdoptSubstructure(void)
Definition: Algorithm.cxx:408
pair< string, Algorithm * > AlgMapPair
Definition: Algorithm.h:52
Algorithm * AdoptAlgorithm(const AlgId &algid) const
Definition: AlgFactory.cxx:127
Algorithm ID (algorithm name + configuration set name)
Definition: AlgId.h:35
virtual const AlgId & Id(void) const
Get algorithm ID.
Definition: Algorithm.h:98
static AlgFactory * Instance()
Definition: AlgFactory.cxx:75
string RgKey
Registry * fConfig
Summary configuration derived from fConvVect, not necessarily allocated.
Definition: Algorithm.h:194
#define pNOTICE
Definition: Messenger.h:62
The GENIE Algorithm Factory.
Definition: AlgFactory.h:40
string Key(void) const
Definition: AlgId.h:47
RgAlg GetAlg(RgKey key) const
Definition: Registry.cxx:503
#define pDEBUG
Definition: Messenger.h:64
map< RgKey, RegistryItemI * > RgIMap
Definition: Registry.h:46
virtual bool genie::Algorithm::AllowReconfig ( void  ) const
inlinevirtualinherited
bool Intranuke::CanRescatter ( const GHepParticle p) const
protectedinherited

Definition at line 231 of file Intranuke.cxx.

References ana::assert(), genie::kPdgKP, genie::kPdgNeutron, genie::kPdgPi0, genie::kPdgPiM, genie::kPdgPiP, genie::kPdgProton, and genie::GHepParticle::Pdg().

Referenced by genie::Intranuke::TransportHadrons().

232 {
233 // checks whether a particle that needs to be rescattered, can in fact be
234 // rescattered by this cascade MC
235 
236  assert(p);
237  return ( p->Pdg() == kPdgPiP ||
238  p->Pdg() == kPdgPiM ||
239  p->Pdg() == kPdgPi0 ||
240  p->Pdg() == kPdgProton ||
241  p->Pdg() == kPdgNeutron ||
242  // p->Pdg() == kPdgGamma ||
243  p->Pdg() == kPdgKP //||
244  // p->Pdg() == kPdgKM
245  );
246 }
int Pdg(void) const
Definition: GHepParticle.h:64
const int kPdgKP
Definition: PDGCodes.h:149
const int kPdgPiP
Definition: PDGCodes.h:135
const int kPdgPi0
Definition: PDGCodes.h:137
assert(nhit_max >=nhit_nbins)
const int kPdgPiM
Definition: PDGCodes.h:136
const int kPdgProton
Definition: PDGCodes.h:65
const int kPdgNeutron
Definition: PDGCodes.h:67
AlgCmp_t Algorithm::Compare ( const Algorithm alg) const
virtualinherited

Compare with input algorithm.

Definition at line 294 of file Algorithm.cxx.

References genie::AlgId::Config(), genie::Algorithm::Id(), genie::kAlgCmpDiffAlg, genie::kAlgCmpDiffConfig, genie::kAlgCmpIdentical, genie::kAlgCmpUnknown, and genie::AlgId::Name().

Referenced by genie::Algorithm::AllowReconfig().

295 {
296 // Compares itself with the input algorithm
297 
298  string alg1 = this->Id().Name();
299  string config1 = this->Id().Config();
300  string alg2 = algo->Id().Name();
301  string config2 = algo->Id().Config();
302 
303  if(alg1 == alg2)
304  {
305  if(config1 == config2) return kAlgCmpIdentical;
306  else return kAlgCmpDiffConfig;
307  }
308  else return kAlgCmpDiffAlg;
309 
310  return kAlgCmpUnknown;
311 }
string Name(void) const
Definition: AlgId.h:45
virtual const AlgId & Id(void) const
Get algorithm ID.
Definition: Algorithm.h:98
string Config(void) const
Definition: AlgId.h:46
void HAIntranuke::Configure ( string  config)
virtual

Configure the algorithm from the AlgoConfigPool based on param_set string given in input An algorithm contains a vector of registries coming from different xml configuration files, which are loaded according a very precise prioriy This methods will load a number registries in order of priority: 1) "Tunable" parameter set from CommonParametes. This is loaded with the highest prioriry and it is designed to be used for tuning procedure Usage not expected from the user. 2) For every string defined in "CommonParame" the corresponding parameter set will be loaded from CommonParameter.xml 3) parameter set specified by the config string and defined in the xml file of the algorithm 4) if config is not "Default" also the Default parameter set from the same xml file will be loaded Effectively this avoids the repetion of a parameter when it is not changed in the requested configuration

Reimplemented from genie::Algorithm.

Definition at line 1379 of file HAIntranuke.cxx.

References genie::Algorithm::Configure(), and LoadConfig().

1379  {
1380 
1381  Algorithm::Configure(param_set) ;
1382  this -> LoadConfig() ;
1383 }
virtual void Configure(const Registry &config)
Definition: Algorithm.cxx:70
void Intranuke::Configure ( const Registry config)
virtualinherited

Configure the algorithm with an external registry The registry is merged with the top level registry if it is owned, Otherwise a copy of it is added with the highest priority

Reimplemented from genie::Algorithm.

Definition at line 422 of file Intranuke.cxx.

References genie::Algorithm::Configure(), and genie::Intranuke::LoadConfig().

423 {
424  Algorithm::Configure(config);
425  this->LoadConfig();
426 }
virtual void Configure(const Registry &config)
Definition: Algorithm.cxx:70
virtual void LoadConfig(void)=0
void Algorithm::DeleteConfig ( void  )
protectedinherited

Definition at line 471 of file Algorithm.cxx.

References MECModelEnuComparisons::i.

Referenced by genie::Algorithm::AllowReconfig().

472 {
473  // there is nothing to delete if the configuration is not owned but is
474  // rather looked up from the configuration pool
475  //
476 
477  for ( unsigned int i = 0 ; i < fConfVect.size() ; ++i ) {
478  if ( fOwnerships[i] ) {
479  delete fConfVect[i] ;
480  }
481  }
482 
483  fConfVect.clear() ;
484  fOwnerships.clear() ;
485 
486  // delete owned configuration registry
487 
488  if(fConfig) {
489  delete fConfig;
490  fConfig=0;
491  }
492 
493 }
vector< Registry * > fConfVect
Definition: Algorithm.h:161
vector< bool > fOwnerships
ownership for every registry in fConfVect
Definition: Algorithm.h:164
Registry * fConfig
Summary configuration derived from fConvVect, not necessarily allocated.
Definition: Algorithm.h:194
void Algorithm::DeleteSubstructure ( void  )
protectedinherited

Definition at line 496 of file Algorithm.cxx.

Referenced by genie::Algorithm::AllowReconfig().

497 {
498  // there is nothing to delete if the sub-algorithms are not owned but rather
499  // taken from the AlgFactory's pool
500  //
501  if(!fOwnsSubstruc) return;
502 
503  // delete local algorithm pool
504  //
505  AlgMapIter iter = fOwnedSubAlgMp->begin();
506  for( ; iter != fOwnedSubAlgMp->end(); ++iter) {
507  Algorithm * alg = iter->second;
508  if(alg) {
509  delete alg;
510  alg=0;
511  }
512  }
513  delete fOwnedSubAlgMp;
514  fOwnedSubAlgMp = 0;
515 }
AlgMap * fOwnedSubAlgMp
local pool for owned sub-algs (taken out of the factory pool)
Definition: Algorithm.h:167
bool fOwnsSubstruc
true if it owns its substructure (sub-algs,...)
Definition: Algorithm.h:155
Algorithm abstract base class.
Definition: Algorithm.h:54
map< string, Algorithm * >::iterator AlgMapIter
Definition: Algorithm.h:50
void HAIntranuke::ElasHA ( GHepRecord ev,
GHepParticle p,
INukeFateHA_t  fate 
) const
private

Definition at line 468 of file HAIntranuke.cxx.

References genie::GHepParticle::A(), genie::GHepRecord::AddParticle(), genie::INukeHadroFates::AsString(), genie::PDGLibrary::Find(), genie::Intranuke::fRemnA, genie::Intranuke::fRemnP4, genie::Intranuke::fRemnZ, genie::PDGLibrary::Instance(), genie::kIHAFtElas, genie::kIStStableFinalState, genie::kPdgNeutron, genie::kPdgProton, LOG, genie::utils::res::Mass(), genie::GHepParticle::Mass(), genie::GHepParticle::Name(), genie::GHepParticle::P4(), pDEBUG, genie::GHepParticle::Pdg(), PiBounce(), pINFO, PnBounce(), pNOTICE, pWARN, genie::GHepParticle::SetMomentum(), genie::exceptions::INukeException::SetReason(), genie::GHepParticle::SetStatus(), genie::GHepRecord::TargetNucleus(), and genie::utils::intranuke::TwoBodyKinematics().

Referenced by SimulateHadronicFinalStateKinematics().

470 {
471  // scatters particle within nucleus, copy of hN code meant to run only once
472  // in hA mode
473 
474 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
475  LOG("HAIntranuke", pDEBUG)
476  << "ElasHA() is invoked for a : " << p->Name()
477  << " whose fate is : " << INukeHadroFates::AsString(fate);
478 #endif
479 
480  if(fate!=kIHAFtElas)
481  {
482  LOG("HAIntranuke", pWARN)
483  << "ElasHA() cannot handle fate: " << INukeHadroFates::AsString(fate);
484  return;
485  }
486 
487  // check remnants
488  if(fRemnA<0 || fRemnZ<0) // best to stop it here and not try again.
489  {
490  LOG("HAIntranuke", pWARN) << "Invalid Nucleus! : (A,Z) = ("<<fRemnA<<','<<fRemnZ<<')';
492  ev->AddParticle(*p);
493  return;
494  }
495 
496  // vars for incoming particle, target, and scattered pdg codes
497  int pcode = p->Pdg();
498  double Mp = p->Mass();
499  double Mt = 0.;
500  if (ev->TargetNucleus()->A()==fRemnA)
501  { Mt = PDGLibrary::Instance()->Find(ev->TargetNucleus()->Pdg())->Mass(); }
502  else
503  {
504  Mt = fRemnP4.M();
505  }
506  TLorentzVector t4PpL = *p->P4();
507  TLorentzVector t4PtL = fRemnP4;
508  double C3CM = 0.0;
509 
510  // calculate scattering angle
511  if(pcode==kPdgNeutron||pcode==kPdgProton) C3CM = TMath::Cos(this->PnBounce());
512  else C3CM = TMath::Cos(this->PiBounce());
513 
514  // calculate final 4 momentum of probe
515  TLorentzVector t4P3L, t4P4L;
516 
517  if (!utils::intranuke::TwoBodyKinematics(Mp,Mt,t4PpL,t4PtL,t4P3L,t4P4L,C3CM,fRemnP4))
518  {
519  LOG("HAIntranuke", pNOTICE) << "ElasHA() failed";
521  exception.SetReason("TwoBodyKinematics failed in ElasHA");
522  throw exception;
523  }
524 
525  // Update probe particle
526  p->SetMomentum(t4P3L);
528 
529  // Update Remnant nucleus
530  fRemnP4 = t4P4L;
531  LOG("HAIntranuke",pINFO)
532  << "C3cm = " << C3CM;
533  LOG("HAIntranuke",pINFO)
534  << "|p3| = " << t4P3L.Vect().Mag() << ", E3 = " << t4P3L.E() << ",Mp = " << Mp;
535  LOG("HAIntranuke",pINFO)
536  << "|p4| = " << fRemnP4.Vect().Mag() << ", E4 = " << fRemnP4.E() << ",Mt = " << Mt;
537 
538  ev->AddParticle(*p);
539 
540 }
An exception thrown by SimulateHadronState for kinematics problems. TwoBodyCollision/Kinematics used ...
TLorentzVector fRemnP4
P4 of remnant system.
Definition: Intranuke.h:98
const TLorentzVector * P4(void) const
Definition: GHepParticle.h:79
double PnBounce(void) const
double PiBounce(void) const
::xsd::cxx::tree::exception< char > exception
Definition: Database.h:225
double Mass(Resonance_t res)
resonance mass (GeV)
void SetMomentum(const TLorentzVector &p4)
double Mass(void) const
Mass that corresponds to the PDG code.
int Pdg(void) const
Definition: GHepParticle.h:64
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:97
bool TwoBodyKinematics(double M3, double M4, TLorentzVector tP1L, TLorentzVector tP2L, TLorentzVector &tP3L, TLorentzVector &tP4L, double C3CM, TLorentzVector &RemnP4, double bindE=0)
Definition: INukeUtils.cxx:754
virtual GHepParticle * TargetNucleus(void) const
Definition: GHepRecord.cxx:330
#define pINFO
Definition: Messenger.h:63
int fRemnA
remnant nucleus A
Definition: Intranuke.h:96
#define pWARN
Definition: Messenger.h:61
static PDGLibrary * Instance(void)
Definition: PDGLibrary.cxx:43
void SetStatus(GHepStatus_t s)
Definition: GHepParticle.h:127
static string AsString(INukeFateHN_t fate)
TParticlePDG * Find(int pdgc)
Definition: PDGLibrary.cxx:61
int fRemnZ
remnant nucleus Z
Definition: Intranuke.h:97
virtual void AddParticle(const GHepParticle &p)
Definition: GHepRecord.cxx:535
const int kPdgProton
Definition: PDGCodes.h:65
int A(void) const
#define pNOTICE
Definition: Messenger.h:62
const int kPdgNeutron
Definition: PDGCodes.h:67
#define pDEBUG
Definition: Messenger.h:64
Registry * Algorithm::ExtractLocalConfig ( const Registry in) const
protectedinherited

Split an incoming configuration Registry into a block valid for this algorithm Ownership of the returned registry belongs to the algo

Definition at line 518 of file Algorithm.cxx.

References genie::RegistryItemI::Clone(), genie::Registry::GetItemMap(), genie::Registry::Name(), and confusionMatrixTree::out.

Referenced by genie::Algorithm::AllowReconfig().

518  {
519 
520  const RgIMap & rgmap = in.GetItemMap();
521  Registry * out = new Registry( in.Name(), false );
522 
523  for( RgIMapConstIter reg_iter = rgmap.begin();
524  reg_iter != rgmap.end(); ++reg_iter ) {
525 
526  RgKey reg_key = reg_iter->first;
527  if( reg_key.find( '/' ) != string::npos) continue;
528 
529  // at this point
530  // this key is referred to the local algorithm
531  // it has to be copied in out;
532 
533  RegistryItemI * ri = reg_iter->second;
534  RgIMapPair key_item_pair( reg_key, ri->Clone() );
535  out -> Set(key_item_pair);
536 
537  }
538 
539  if ( out -> NEntries() <= 0 ) {
540  delete out ;
541  out = 0 ;
542  }
543 
544  return out ;
545 }
Registry item pABC.
Definition: RegistryItemI.h:30
string Name(void) const
get the registry name
Definition: Registry.cxx:612
map< RgKey, RegistryItemI * >::const_iterator RgIMapConstIter
Definition: Registry.h:50
const RgIMap & GetItemMap(void) const
Definition: Registry.h:162
pair< RgKey, RegistryItemI * > RgIMapPair
Definition: Registry.h:47
string RgKey
A registry. Provides the container for algorithm configuration parameters.
Definition: Registry.h:66
virtual RegistryItemI * Clone(void) const =0
map< RgKey, RegistryItemI * > RgIMap
Definition: Registry.h:46
Registry * Algorithm::ExtractLowerConfig ( const Registry in,
const string alg_key 
) const
protectedinherited

Split an incoming configuration Registry into a block valid for the sub-algo identified by alg_key.

Definition at line 549 of file Algorithm.cxx.

References genie::RegistryItemI::Clone(), genie::Registry::GetItemMap(), genie::Registry::Name(), and confusionMatrixTree::out.

Referenced by genie::Algorithm::AllowReconfig().

549  {
550 
551  const RgIMap & rgmap = in.GetItemMap();
552  Registry * out = new Registry( in.Name(), false );
553 
554  for( RgIMapConstIter reg_iter = rgmap.begin();
555  reg_iter != rgmap.end(); ++reg_iter ) {
556 
557  RgKey reg_key = reg_iter->first;
558  if( reg_key.find(alg_key+"/") == string::npos) continue;
559 
560  // at this point
561  // this key is referred to the sub-algorithm
562  // indicated by alg_key: it has to be copied in out;
563 
564  int new_key_start = reg_key.find_first_of('/')+1;
565  RgKey new_reg_key = reg_key.substr( new_key_start, reg_key.length() );
566 
567  RegistryItemI * ri = reg_iter->second;
568  RgIMapPair key_item_pair(new_reg_key, ri->Clone());
569  out -> Set(key_item_pair);
570 
571  }
572 
573  if ( out -> NEntries() <= 0 ) {
574  delete out ;
575  out = 0 ;
576  }
577 
578  return out ;
579 
580 }
Registry item pABC.
Definition: RegistryItemI.h:30
string Name(void) const
get the registry name
Definition: Registry.cxx:612
map< RgKey, RegistryItemI * >::const_iterator RgIMapConstIter
Definition: Registry.h:50
const RgIMap & GetItemMap(void) const
Definition: Registry.h:162
pair< RgKey, RegistryItemI * > RgIMapPair
Definition: Registry.h:47
string RgKey
A registry. Provides the container for algorithm configuration parameters.
Definition: Registry.h:66
virtual RegistryItemI * Clone(void) const =0
map< RgKey, RegistryItemI * > RgIMap
Definition: Registry.h:46
void Algorithm::FindConfig ( void  )
virtualinherited

Lookup configuration from the config pool Similar logic from void Configure(string)

Definition at line 135 of file Algorithm.cxx.

References gen_flatrecord::config, exit(), genie::AlgConfigPool::FindRegistry(), genie::Registry::GetItemMap(), genie::Registry::GetString(), MECModelEnuComparisons::i, genie::AlgConfigPool::Instance(), it, genie::Registry::ItemIsLocal(), parse_dependency_file_t::list, LOG, pDEBUG, pFATAL, time_estimates::pool, pWARN, moon_position_table_new3::second, genie::utils::str::Split(), string, and APDHVSetting::temp.

136 {
137 // Finds its configration Registry from the ConfigPool and gets a pointer to
138 // it. If the Registry comes from the ConfigPool then the Algorithm does not
139 // own its configuration (the ConfigPool singleton keeps the ownership and the
140 // responsibility to -eventually- delete all the Registries it instantiates
141 // by parsing the XML config files).
142 
143  DeleteConfig() ;
144 
146 
147  Registry * config = 0 ;
148 
149  // load the Default config if config is not the default
150  if ( fID.Config() != "Default" ) {
151  config = pool -> FindRegistry( fID.Name(), "Default" );
152  if ( config ) {
153  if ( config -> NEntries() > 0 ) {
154  AddTopRegistry( config, false ) ;
155  LOG("Algorithm", pDEBUG) << "\n" << *config;
156  }
157  }
158  }
159 
160  // Load the right config
161  config = pool->FindRegistry(this);
162 
163  if(!config)
164  // notify & keep whatever config Registry was used before.
165  LOG("Algorithm", pWARN)
166  << "No Configuration available for "
167  << this->Id().Key() << " at the ConfigPool";
168  else {
169  if ( config -> NEntries() > 0 ) {
170  AddTopRegistry( config, false ) ;
171  LOG("Algorithm", pDEBUG) << "\n" << config;
172  }
173  }
174 
175  const string common_key_root = "Common" ;
176  std::map<string, string> common_lists;
177 
178  // Load Common Parameters if key that start with "Common" is found
179  for ( unsigned int i = 0 ; i < fConfVect.size() ; ++i ) {
180  const Registry & temp = * fConfVect[i] ;
181  for ( RgIMapConstIter it = temp.GetItemMap().begin() ; it != temp.GetItemMap().end() ; ++it ) {
182 
183  // check if it is a "Common" entry
184  if ( it -> first.find( common_key_root ) == 0 ) {
185  // retrieve the type of the common entry
186  std::string type = it -> first.substr(common_key_root.size() ) ;
187 
188  if ( temp.ItemIsLocal( it -> first ) ) {
189 
190  string temp_list = temp.GetString( it -> first ) ;
191  if ( temp_list.length() > 0 ) {
192  common_lists[type] = temp_list ;
193  }
194  }
195  }
196 
197  }
198 
199  } // loop over the local registries
200 
201 
202  for ( std::map<string, string>::const_iterator it = common_lists.begin() ;
203  it != common_lists.end() ; ++it ) {
204 
205  vector<string> list = str::Split( it -> second , "," ) ;
206 
207  for ( unsigned int i = 0; i < list.size(); ++i ) {
208 
209  config = pool -> CommonList( it -> first, list[i] ) ;
210 
211  if ( ! config ) {
212  LOG("Algorithm", pFATAL)
213  << "No Commom parameters available for " << it -> first << " list "
214  << list[i] << " at the ConfigPool";
215 
216  exit( 78 ) ;
217  }
218  else {
219  AddLowRegistry( config, false ) ;
220  LOG("Algorithm", pDEBUG) << "Loading "
221  << it -> first << " registry "
222  << list[i] << " \n" << config;
223  }
224 
225  }
226 
227  }
228 
229 
230  // Load Tunable from CommonParameters
231  // only if the option is specified in RunOpt
232  config = pool -> CommonList( "Param", "Tunable" ) ;
233  if ( config ) {
234  if ( config -> NEntries() > 0 ) {
235  AddTopRegistry( config, false ) ;
236  LOG("Algorithm", pDEBUG) << "Loading Tunable registry \n" << config;
237  }
238  }
239  else {
240  // notify & keep whatever config Registry was used before.
241  LOG("Algorithm", pWARN)
242  << "No Tunable parameter set available at the ConfigPool";
243  }
244 
245  if ( fConfig ) {
246  delete fConfig ;
247  fConfig = 0 ;
248  }
249 
250 }
set< int >::iterator it
A singleton class holding all configuration registries built while parsing all loaded XML configurati...
Definition: AlgConfigPool.h:41
#define pFATAL
Definition: Messenger.h:57
Definition: config.py:1
AlgId fID
algorithm name and configuration set
Definition: Algorithm.h:156
map< RgKey, RegistryItemI * >::const_iterator RgIMapConstIter
Definition: Registry.h:50
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
const RgIMap & GetItemMap(void) const
Definition: Registry.h:162
string Name(void) const
Definition: AlgId.h:45
bool ItemIsLocal(RgKey key) const
local or global?
Definition: Registry.cxx:193
int AddTopRegistry(Registry *rp, bool owns=true)
add registry with top priority, also update ownership
Definition: Algorithm.cxx:585
#define pWARN
Definition: Messenger.h:61
void DeleteConfig(void)
Definition: Algorithm.cxx:471
virtual const AlgId & Id(void) const
Get algorithm ID.
Definition: Algorithm.h:98
RgStr GetString(RgKey key) const
Definition: Registry.cxx:496
vector< Registry * > fConfVect
Definition: Algorithm.h:161
vector< string > Split(string input, string delim)
Definition: StringUtils.cxx:42
A registry. Provides the container for algorithm configuration parameters.
Definition: Registry.h:66
exit(0)
Registry * fConfig
Summary configuration derived from fConvVect, not necessarily allocated.
Definition: Algorithm.h:194
Registry * FindRegistry(string key) const
int AddLowRegistry(Registry *rp, bool owns=true)
add registry with lowest priority, also update ownership
Definition: Algorithm.cxx:601
string Key(void) const
Definition: AlgId.h:47
static AlgConfigPool * Instance()
#define pDEBUG
Definition: Messenger.h:64
string Config(void) const
Definition: AlgId.h:46
enum BeamMode string
double Intranuke::GenerateStep ( GHepRecord ev,
GHepParticle p 
) const
protectedinherited

Definition at line 393 of file Intranuke.cxx.

References d, genie::Intranuke::fDelRNucleon, genie::Intranuke::fDelRPion, genie::Intranuke::fNucleonMFPScale, genie::Intranuke::fPionMFPScale, genie::Intranuke::fRemnA, genie::Intranuke::fRemnZ, genie::RandomGen::Instance(), genie::kPdgNeutron, genie::kPdgPi0, genie::kPdgPiM, genie::kPdgPiP, genie::kPdgProton, CLHEP::L, LOG, genie::utils::intranuke::MeanFreePath(), genie::GHepParticle::P4(), pDEBUG, genie::GHepParticle::Pdg(), generate_hists::rnd, genie::RandomGen::RndFsi(), scale, and genie::GHepParticle::X4().

Referenced by genie::Intranuke::TransportHadrons().

394 {
395 // Generate a step (in fermis) for particle p in the input event.
396 // Computes the mean free path L and generate an 'interaction' distance d
397 // from an exp(-d/L) distribution
398 
400 
401  int pdgc = p->Pdg();
402 
403  double scale = 1.;
404  if (pdgc==kPdgPiP || pdgc==kPdgPiM || pdgc==kPdgPi0) {
405  scale = fPionMFPScale;
406  }
407  else if (pdgc==kPdgProton || pdgc==kPdgNeutron) {
408  scale = fNucleonMFPScale;
409  }
410 
411  double L = utils::intranuke::MeanFreePath(pdgc, *p->X4(), *p->P4(), fRemnA,
413  L *= scale;
414  double d = -1.*L * TMath::Log(rnd->RndFsi().Rndm());
415 
416  LOG("Intranuke", pDEBUG)
417  << "Mean free path = " << L << " fm / "
418  << "Generated path length = " << d << " fm";
419  return d;
420 }
TRandom3 & RndFsi(void) const
rnd number generator used by intranuclear cascade monte carlos
Definition: RandomGen.h:60
static RandomGen * Instance()
Access instance.
Definition: RandomGen.cxx:79
const TLorentzVector * P4(void) const
Definition: GHepParticle.h:79
A singleton holding random number generator classes. All random number generation in GENIE should tak...
Definition: RandomGen.h:30
int Pdg(void) const
Definition: GHepParticle.h:64
double fDelRNucleon
factor by which Nucleon Compton wavelength gets multiplied to become nuclear size enhancement ...
Definition: Intranuke.h:106
double fNucleonMFPScale
Definition: Intranuke.h:123
Double_t scale
Definition: plot.C:25
double fDelRPion
factor by which Pion Compton wavelength gets multiplied to become nuclear size enhancement ...
Definition: Intranuke.h:105
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
static constexpr double L
Float_t d
Definition: plot.C:236
const int kPdgPiP
Definition: PDGCodes.h:135
const int kPdgPi0
Definition: PDGCodes.h:137
int fRemnA
remnant nucleus A
Definition: Intranuke.h:96
double MeanFreePath(int pdgc, const TLorentzVector &x4, const TLorentzVector &p4, double A, double Z, double nRpi=0.5, double nRnuc=1.0)
Mean free path (pions, nucleons)
Definition: INukeUtils.cxx:70
const TLorentzVector * X4(void) const
Definition: GHepParticle.h:80
int fRemnZ
remnant nucleus Z
Definition: Intranuke.h:97
const int kPdgPiM
Definition: PDGCodes.h:136
const int kPdgProton
Definition: PDGCodes.h:65
double fPionMFPScale
tweaking factors for tuning
Definition: Intranuke.h:117
const int kPdgNeutron
Definition: PDGCodes.h:67
#define pDEBUG
Definition: Messenger.h:64
void Intranuke::GenerateVertex ( GHepRecord ev) const
protectedinherited

Definition at line 151 of file Intranuke.cxx.

References ana::assert(), epsilon, genie::Intranuke::fTrackingRadius, genie::RandomGen::Instance(), genie::pdg::IsIon(), genie::pdg::IsPseudoParticle(), LOG, genie::GHepParticle::P4(), genie::GHepParticle::Pdg(), pNOTICE, genie::GHepRecord::Probe(), generate_hists::rnd, genie::RandomGen::RndFsi(), genie::GHepParticle::SetPosition(), ana::Sqrt(), genie::GHepRecord::TargetNucleus(), genie::utils::print::Vec3AsString(), submit_syst::x, and submit_syst::y.

Referenced by genie::Intranuke::ProcessEventRecord().

152 {
153 // Sets a vertex in the nucleus periphery
154 // Called onlt in hadron/photon-nucleus interactions.
155 
156  GHepParticle * nucltgt = evrec->TargetNucleus();
157  assert(nucltgt);
158 
160  TVector3 vtx(999999.,999999.,999999.);
161 
162  // *** For h+A events (test mode):
163  // Assume a hadron beam with uniform intensity across an area,
164  // so we need to choose events uniformly within that area.
165  double x=999999., y=999999., epsilon = 0.001;
166  double R2 = TMath::Power(fTrackingRadius,2.);
167  double rp2 = TMath::Power(x,2.) + TMath::Power(y,2.);
168  while(rp2 > R2-epsilon) {
169  x = (fTrackingRadius-epsilon) * rnd->RndFsi().Rndm();
170  y = -fTrackingRadius + 2*fTrackingRadius * rnd->RndFsi().Rndm();
171  y -= ((y>0) ? epsilon : -epsilon);
172  rp2 = TMath::Power(x,2.) + TMath::Power(y,2.);
173  }
174  vtx.SetXYZ(x,y, -1.*TMath::Sqrt(TMath::Max(0.,R2-rp2)) + epsilon);
175 
176  // get the actual unit vector along the incoming hadron direction
177  TVector3 direction = evrec->Probe()->P4()->Vect().Unit();
178 
179  // rotate the vtx position
180  vtx.RotateUz(direction);
181 
182  LOG("Intranuke", pNOTICE)
183  << "Generated vtx @ R = " << vtx.Mag() << " fm / "
184  << print::Vec3AsString(&vtx);
185 
186  TObjArrayIter piter(evrec);
187  GHepParticle * p = 0;
188  while( (p = (GHepParticle *) piter.Next()) )
189  {
190  if(pdg::IsPseudoParticle(p->Pdg())) continue;
191  if(pdg::IsIon (p->Pdg())) continue;
192 
193  p->SetPosition(vtx.x(), vtx.y(), vtx.z(), 0.);
194  }
195 }
TRandom3 & RndFsi(void) const
rnd number generator used by intranuclear cascade monte carlos
Definition: RandomGen.h:60
static RandomGen * Instance()
Access instance.
Definition: RandomGen.cxx:79
const char * p
Definition: xmltok.h:285
double fTrackingRadius
tracking radius for the nucleus in the current event
Definition: Intranuke.h:91
A singleton holding random number generator classes. All random number generation in GENIE should tak...
Definition: RandomGen.h:30
int Pdg(void) const
Definition: GHepParticle.h:64
void SetPosition(const TLorentzVector &v4)
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
bool IsIon(int pdgc)
Definition: PDGUtils.cxx:40
Var Sqrt(const Var &v)
Use to take sqrt of a var.
Definition: Var.cxx:326
bool IsPseudoParticle(int pdgc)
Definition: PDGUtils.cxx:25
double epsilon
assert(nhit_max >=nhit_nbins)
#define pNOTICE
Definition: Messenger.h:62
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:40
string Vec3AsString(const TVector3 *vec)
Definition: PrintUtils.cxx:87
const Registry & Algorithm::GetConfig ( void  ) const
virtualinherited

Get configuration registry Evaluate the summary of the configuration and returns it The summary of a configuration is a merge of all the registries known to the algorithm (see Configure() methods) but every parameter is appearing only once and in case of repetitions, only the parameter from the registry with the highest prioriry is considered.

Definition at line 254 of file Algorithm.cxx.

References febshutoff_auto::end, genie::Algorithm::GetConfig(), MECModelEnuComparisons::i, LOG, pDEBUG, r(), and moon_position_table_new3::second.

Referenced by genie::EventGeneratorListAssembler::AssembleGeneratorList(), GetAlgorithms(), genie::Algorithm::GetConfig(), genie::GRV98LO::GRV98LO(), genie::NewQELXSec::Integrate(), genie::LHAPDF5::LHAPDF5(), genie::IBDXSecMap::LoadConfig(), genie::Decayer::LoadConfig(), genie::PythiaHadronization::LoadConfig(), genie::FGMBodekRitchie::LoadConfig(), genie::NuclearModelMap::LoadConfig(), genie::SmithMonizUtils::LoadConfig(), main(), genie::AlgFactory::Print(), TestPythiaTauDecays(), testReconfigInOwnedModules(), and genie::P33PaschosLalakulichPXSec::XSec().

254  {
255 
256  if ( fConfig ) return * fConfig ;
257 
258  const_cast<Algorithm*>( this ) -> fConfig = new Registry( fID.Key() + "_summary", false ) ;
259 
260  // loop and append
261  // understand the append mechanism
262  for ( unsigned int i = 0 ; i < fConfVect.size(); ++i ) {
263  fConfig -> Append( * fConfVect[i] ) ;
264  }
265 
266  if ( fOwnsSubstruc ) {
267 
268  for ( AlgMapConstIter iter = fOwnedSubAlgMp -> begin() ;
269  iter != fOwnedSubAlgMp -> end() ; ++iter ) {
270 
271  Algorithm * subalg = iter -> second ;
272 
273  LOG("Algorithm", pDEBUG) << "Appending config from " << iter -> first << " -> " << subalg -> Id() ;
274  const Registry & r = subalg->GetConfig();
275  RgKey prefix = iter -> first + "/";
276  fConfig -> Append(r,prefix);
277 
278  }
279 
280  } //if owned substructure
281 
282  return * fConfig ;
283 }
AlgMap * fOwnedSubAlgMp
local pool for owned sub-algs (taken out of the factory pool)
Definition: Algorithm.h:167
bool fOwnsSubstruc
true if it owns its substructure (sub-algs,...)
Definition: Algorithm.h:155
Algorithm abstract base class.
Definition: Algorithm.h:54
AlgId fID
algorithm name and configuration set
Definition: Algorithm.h:156
virtual const Registry & GetConfig(void) const
Definition: Algorithm.cxx:254
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
const XML_Char * prefix
Definition: expat.h:380
map< string, Algorithm * >::const_iterator AlgMapConstIter
Definition: Algorithm.h:51
virtual const AlgId & Id(void) const
Get algorithm ID.
Definition: Algorithm.h:98
vector< Registry * > fConfVect
Definition: Algorithm.h:161
string RgKey
A registry. Provides the container for algorithm configuration parameters.
Definition: Registry.h:66
TRandom3 r(0)
Registry * fConfig
Summary configuration derived from fConvVect, not necessarily allocated.
Definition: Algorithm.h:194
string Key(void) const
Definition: AlgId.h:47
#define pDEBUG
Definition: Messenger.h:64
Registry * Algorithm::GetOwnedConfig ( void  )
inherited

Returns the pointer of the summary registry, see previous method Gives access to the summary so it could be changed. The usage of this method is deprecated as it is mantained only for back compatibility. If you need to add or chage a parter (or more), use the AddTopRegistry() instead

Definition at line 287 of file Algorithm.cxx.

References GetConfig().

Referenced by genie::TransverseEnhancementFFModel::LoadConfig(), and genie::EffectiveSF::LoadConfig().

288 {
289 
290  GetConfig() ;
291  return fConfig;
292 }
virtual const Registry & GetConfig(void) const
Definition: Algorithm.cxx:254
Registry * fConfig
Summary configuration derived from fConvVect, not necessarily allocated.
Definition: Algorithm.h:194
template<class T >
bool genie::Algorithm::GetParam ( const RgKey name,
T p,
bool  is_top_call = true 
) const
protectedinherited

Ideal access to a parameter value from the vector of registries Returns true if the value is found and the parameters is set

Referenced by genie::CollinsSpillerFragm::BuildFunction(), genie::PetersonFragm::BuildFunction(), genie::COHXSec::LoadConfig(), genie::DISXSec::LoadConfig(), genie::INukeDeltaPropg::LoadConfig(), genie::HadronTransporter::LoadConfig(), genie::DFRKinematicsGenerator::LoadConfig(), genie::RSHelicityAmplModelNCn::LoadConfig(), genie::RSHelicityAmplModelNCp::LoadConfig(), genie::BaryonResonanceDecayer::LoadConfig(), genie::RESKinematicsGenerator::LoadConfig(), genie::DMDISXSec::LoadConfig(), genie::DipoleAxialFormFactorModel::LoadConfig(), genie::DipoleELFormFactorsModel::LoadConfig(), genie::COHKinematicsGenerator::LoadConfig(), genie::VertexGenerator::LoadConfig(), genie::H3AMNuGammaPXSec::LoadConfig(), genie::RSPPResonanceSelector::LoadConfig(), genie::IBDXSecMap::LoadConfig(), genie::Decayer::LoadConfig(), genie::DISHadronicSystemGenerator::LoadConfig(), genie::EmpiricalMECPXSec2015::LoadConfig(), genie::COHElasticPXSec::LoadConfig(), genie::KuzminNaumov2016AxialFormFactorModel::LoadConfig(), genie::SlowRsclCharmDISPXSecLO::LoadConfig(), genie::UnstableParticleDecayer::LoadConfig(), genie::AhrensNCELPXSec::LoadConfig(), genie::AlamSimoAtharVacasSKPXSec2014::LoadConfig(), genie::PythiaHadronization::LoadConfig(), genie::ReinDFRPXSec::LoadConfig(), genie::QPMDISPXSec::LoadConfig(), genie::DFRXSec::LoadConfig(), genie::RosenbluthPXSec::LoadConfig(), genie::StrumiaVissaniIBDPXSec::LoadConfig(), genie::MECGenerator::LoadConfig(), genie::BYPDF::LoadConfig(), genie::KNOPythiaHadronization::LoadConfig(), genie::P33PaschosLalakulichPXSec::LoadConfig(), genie::AhrensDMELPXSec::LoadConfig(), genie::MECXSec::LoadConfig(), genie::AivazisCharmPXSecLO::LoadConfig(), genie::BergerSehgalFMCOHPiPXSec2015::LoadConfig(), genie::ZExpAxialFormFactorModel::LoadConfig(), genie::QPMDMDISPXSec::LoadConfig(), genie::BBA05ELFormFactorsModel::LoadConfig(), genie::BergerSehgalCOHPiPXSec2015::LoadConfig(), genie::BBA03ELFormFactorsModel::LoadConfig(), genie::LwlynSmithQELCCPXSec::LoadConfig(), genie::ReinSehgalRESXSec::LoadConfig(), genie::PrimaryLeptonGenerator::LoadConfig(), genie::NuElectronPXSec::LoadConfig(), genie::FGMBodekRitchie::LoadConfig(), genie::SpectralFunc1d::LoadConfig(), genie::ReinSehgalCOHPiPXSec::LoadConfig(), genie::PaisQELLambdaPXSec::LoadConfig(), genie::OutgoingDarkGenerator::LoadConfig(), genie::ReinSehgalRESXSecFast::LoadConfig(), genie::CharmHadronization::LoadConfig(), genie::NievesSimoVacasMECPXSec2016::LoadConfig(), genie::LHAPDF6::LoadConfig(), genie::EventGenerator::LoadConfig(), genie::NuclearModelMap::LoadConfig(), genie::ReinSehgalSPPXSec::LoadConfig(), genie::ReinSehgalRESPXSec::LoadConfig(), genie::LwlynSmithFF::LoadConfig(), genie::QPMDISStrucFuncBase::LoadConfig(), genie::SmithMonizQELCCPXSec::LoadConfig(), genie::BBA07ELFormFactorsModel::LoadConfig(), LoadConfig(), genie::NievesQELCCPXSec::LoadConfig(), genie::HAIntranuke2018::LoadConfig(), genie::LocalFGM::LoadConfig(), genie::HNIntranuke2018::LoadConfig(), genie::BSKLNBaseRESPXSec2014::LoadConfig(), genie::EffectiveSF::LoadConfig(), genie::ReinSehgalSPPPXSec::LoadConfig(), genie::KNOHadronization::LoadConfig(), genie::SmithMonizUtils::LoadConfig(), genie::MECInteractionListGenerator::LoadConfigData(), genie::PhysInteractionSelector::LoadConfigData(), genie::RESInteractionListGenerator::LoadConfigData(), genie::PauliBlocker::LoadModelType(), genie::BYStrucFunc::ReadBYParams(), and genie::LHAPDF5::SetPDFSetFromConfig().

template<class T >
bool genie::Algorithm::GetParamDef ( const RgKey name,
T p,
const T def 
) const
protectedinherited

Ideal access to a parameter value from the vector of registries, With default value. Returns true if the value is set from the registries, false if the value is the default

Referenced by genie::IMDXSec::LoadConfig(), genie::COHXSec::LoadConfig(), genie::DISXSec::LoadConfig(), genie::RESXSec::LoadConfig(), genie::DFRKinematicsGenerator::LoadConfig(), genie::COHXSecAR::LoadConfig(), genie::SKKinematicsGenerator::LoadConfig(), genie::COHElKinematicsGenerator::LoadConfig(), genie::NuEKinematicsGenerator::LoadConfig(), genie::QELXSec::LoadConfig(), genie::RESKinematicsGenerator::LoadConfig(), genie::DMDISXSec::LoadConfig(), genie::BaryonResonanceDecayer::LoadConfig(), genie::COHKinematicsGenerator::LoadConfig(), genie::IBDKinematicsGenerator::LoadConfig(), genie::NuEInteractionListGenerator::LoadConfig(), genie::QELKinematicsGenerator::LoadConfig(), genie::DMELXSec::LoadConfig(), genie::DISHadronicSystemGenerator::LoadConfig(), genie::DISKinematicsGenerator::LoadConfig(), genie::NucBindEnergyAggregator::LoadConfig(), genie::DMELKinematicsGenerator::LoadConfig(), genie::DMDISKinematicsGenerator::LoadConfig(), genie::QPMDISPXSec::LoadConfig(), genie::DFRXSec::LoadConfig(), genie::FermiMover::LoadConfig(), genie::AlamSimoAtharVacasSKXSec::LoadConfig(), genie::AhrensDMELPXSec::LoadConfig(), genie::NuElectronXSec::LoadConfig(), genie::QELEventGenerator::LoadConfig(), genie::P33PaschosLalakulichPXSec::LoadConfig(), genie::MECXSec::LoadConfig(), genie::QPMDMDISPXSec::LoadConfig(), genie::LwlynSmithQELCCPXSec::LoadConfig(), genie::ReinSehgalRESXSec::LoadConfig(), genie::FGMBodekRitchie::LoadConfig(), genie::ReinSehgalRESXSecFast::LoadConfig(), genie::KovalenkoQELCharmPXSec::LoadConfig(), genie::SmithMonizQELCCXSec::LoadConfig(), genie::ReinSehgalSPPXSec::LoadConfig(), genie::ReinSehgalRESPXSec::LoadConfig(), genie::QELEventGeneratorSM::LoadConfig(), genie::QPMDISStrucFuncBase::LoadConfig(), genie::SmithMonizQELCCPXSec::LoadConfig(), genie::NievesQELCCPXSec::LoadConfig(), LoadConfig(), genie::LocalFGM::LoadConfig(), genie::HNIntranuke2018::LoadConfig(), genie::HAIntranuke2018::LoadConfig(), genie::BSKLNBaseRESPXSec2014::LoadConfig(), genie::EffectiveSF::LoadConfig(), genie::KNOHadronization::LoadConfig(), genie::NewQELXSec::LoadConfig(), genie::DFRInteractionListGenerator::LoadConfigData(), genie::MECInteractionListGenerator::LoadConfigData(), genie::QELInteractionListGenerator::LoadConfigData(), genie::COHInteractionListGenerator::LoadConfigData(), genie::RESInteractionListGenerator::LoadConfigData(), genie::SKInteractionListGenerator::LoadConfigData(), genie::DMELInteractionListGenerator::LoadConfigData(), genie::RSPPInteractionListGenerator::LoadConfigData(), genie::DISInteractionListGenerator::LoadConfigData(), and genie::DMDISInteractionListGenerator::LoadConfigData().

template<class T >
bool genie::Algorithm::GetParamVect ( const std::string comm_name,
std::vector< T > &  v,
unsigned int  max,
bool  is_top_call = true 
) const
protectedinherited

Handle to load vectors of parameters It looks for different registry item with name comm_name0, comm_name1, etc...

virtual AlgStatus_t genie::Algorithm::GetStatus ( void  ) const
inlinevirtualinherited

Get algorithm status.

Definition at line 101 of file Algorithm.h.

References genie::Algorithm::fStatus.

101 { return fStatus; }
AlgStatus_t fStatus
algorithm execution status
Definition: Algorithm.h:166
INukeFateHA_t HAIntranuke::HadronFateHA ( const GHepParticle p) const
private

Definition at line 212 of file HAIntranuke.cxx.

References genie::INukeHadroFates::AsString(), genie::Intranuke::fHadroData, genie::Intranuke::fNucleonFracAbsScale, genie::Intranuke::fNucleonFracCExScale, genie::Intranuke::fNucleonFracElasScale, genie::Intranuke::fNucleonFracInelScale, genie::Intranuke::fNucleonFracPiProdScale, genie::Intranuke::fPionFracAbsScale, genie::Intranuke::fPionFracCExScale, genie::Intranuke::fPionFracElasScale, genie::Intranuke::fPionFracInelScale, genie::Intranuke::fPionFracPiProdScale, genie::INukeHadroData::Frac(), genie::RandomGen::Instance(), genie::kIHAFtAbs, genie::kIHAFtCEx, genie::kIHAFtElas, genie::kIHAFtInelas, genie::kIHAFtPiProd, genie::kIHAFtUndefined, genie::GHepParticle::KinE(), genie::kPdgKM, genie::kPdgKP, genie::kPdgNeutron, genie::kPdgPi0, genie::kPdgPiM, genie::kPdgPiP, genie::kPdgProton, genie::controls::kRjMaxIterations, LOG, genie::units::MeV, genie::GHepParticle::Name(), pDEBUG, genie::GHepParticle::Pdg(), pINFO, pWARN, r(), generate_hists::rnd, and genie::RandomGen::RndFsi().

Referenced by SimulateHadronicFinalState().

213 {
214 // Select a hadron fate in HA mode
215 //
217 
218  // get pdgc code & kinetic energy in MeV
219  int pdgc = p->Pdg();
220  double ke = p->KinE() / units::MeV;
221 
222  LOG("HAIntranuke", pINFO)
223  << "Selecting hA fate for " << p->Name() << " with KE = " << ke << " MeV";
224 
225  // try to generate a hadron fate
226  unsigned int iter = 0;
227  while(iter++ < kRjMaxIterations) {
228 
229  // handle pions
230  //
231  if (pdgc==kPdgPiP || pdgc==kPdgPiM || pdgc==kPdgPi0) {
232 
233  double frac_cex = fHadroData->Frac(pdgc, kIHAFtCEx, ke);
234  double frac_elas = fHadroData->Frac(pdgc, kIHAFtElas, ke);
235  double frac_inel = fHadroData->Frac(pdgc, kIHAFtInelas, ke);
236  double frac_abs = fHadroData->Frac(pdgc, kIHAFtAbs, ke);
237  double frac_piprod = fHadroData->Frac(pdgc, kIHAFtPiProd, ke);
238 
239  LOG("HAIntranuke", pINFO)
240  << "\n frac{" << INukeHadroFates::AsString(kIHAFtCEx) << "} = " << frac_cex
241  << "\n frac{" << INukeHadroFates::AsString(kIHAFtElas) << "} = " << frac_elas
242  << "\n frac{" << INukeHadroFates::AsString(kIHAFtInelas) << "} = " << frac_inel
243  << "\n frac{" << INukeHadroFates::AsString(kIHAFtAbs) << "} = " << frac_abs
244  << "\n frac{" << INukeHadroFates::AsString(kIHAFtPiProd) << "} = " << frac_piprod;
245 
246  // apply external tweaks to fractions
247  frac_cex *= fPionFracCExScale;
248  frac_elas *= fPionFracElasScale;
249  frac_inel *= fPionFracInelScale;
250  frac_abs *= fPionFracAbsScale;
251  frac_piprod *= fPionFracPiProdScale;
252  double frac_rescale = 1./(frac_cex + frac_elas + frac_inel + frac_abs + frac_piprod);
253  frac_cex *= frac_rescale;
254  frac_elas *= frac_rescale;
255  frac_inel *= frac_rescale;
256  frac_abs *= frac_rescale;
257  frac_piprod *= frac_rescale;
258 
259  LOG("HAIntranuke", pINFO)
260  << "\n frac{" << INukeHadroFates::AsString(kIHAFtCEx) << "} = " << frac_cex
261  << "\n frac{" << INukeHadroFates::AsString(kIHAFtElas) << "} = " << frac_elas
262  << "\n frac{" << INukeHadroFates::AsString(kIHAFtInelas) << "} = " << frac_inel
263  << "\n frac{" << INukeHadroFates::AsString(kIHAFtAbs) << "} = " << frac_abs
264  << "\n frac{" << INukeHadroFates::AsString(kIHAFtPiProd) << "} = " << frac_piprod;
265 
266  // compute total fraction (can be <1 if fates have been switched off)
267  double tf = frac_cex +
268  frac_elas +
269  frac_inel +
270  frac_abs +
271  frac_piprod;
272 
273  double r = tf * rnd->RndFsi().Rndm();
274 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
275  LOG("HAIntranuke", pDEBUG) << "r = " << r << " (max = " << tf << ")";
276 #endif
277  double cf=0; // current fraction
278  if(r < (cf += frac_cex )) return kIHAFtCEx; // cex
279  if(r < (cf += frac_elas )) return kIHAFtElas; // elas
280  if(r < (cf += frac_inel )) return kIHAFtInelas; // inelas
281  if(r < (cf += frac_abs )) return kIHAFtAbs; // abs
282  if(r < (cf += frac_piprod )) return kIHAFtPiProd; // pi prod
283 
284  LOG("HAIntranuke", pWARN)
285  << "No selection after going through all fates! "
286  << "Total fraction = " << tf << " (r = " << r << ")";
287  }
288 
289  // handle nucleons
290  else if (pdgc==kPdgProton || pdgc==kPdgNeutron) {
291 
292  double frac_cex = fHadroData->Frac(pdgc, kIHAFtCEx, ke);
293  double frac_elas = fHadroData->Frac(pdgc, kIHAFtElas, ke);
294  double frac_inel = fHadroData->Frac(pdgc, kIHAFtInelas, ke);
295  double frac_abs = fHadroData->Frac(pdgc, kIHAFtAbs, ke);
296  double frac_pipro = fHadroData->Frac(pdgc, kIHAFtPiProd, ke);
297 
298  // apply external tweaks to fractions
299  frac_cex *= fNucleonFracCExScale;
300  frac_elas *= fNucleonFracElasScale;
301  frac_inel *= fNucleonFracInelScale;
302  frac_abs *= fNucleonFracAbsScale;
303  frac_pipro *= fNucleonFracPiProdScale;
304  double frac_rescale = 1./(frac_cex + frac_elas + frac_inel + frac_abs + frac_pipro);
305  frac_cex *= frac_rescale;
306  frac_elas *= frac_rescale;
307  frac_inel *= frac_rescale;
308  frac_abs *= frac_rescale;
309  frac_pipro *= frac_rescale;
310 
311  LOG("HAIntranuke", pDEBUG)
312  << "\n frac{" << INukeHadroFates::AsString(kIHAFtCEx) << "} = " << frac_cex
313  << "\n frac{" << INukeHadroFates::AsString(kIHAFtElas) << "} = " << frac_elas
314  << "\n frac{" << INukeHadroFates::AsString(kIHAFtInelas) << "} = " << frac_inel
315  << "\n frac{" << INukeHadroFates::AsString(kIHAFtAbs) << "} = " << frac_abs
316  << "\n frac{" << INukeHadroFates::AsString(kIHAFtPiProd) << "} = " << frac_pipro;
317 
318  // compute total fraction (can be <1 if fates have been switched off)
319  double tf = frac_cex +
320  frac_elas +
321  frac_inel +
322  frac_abs +
323  frac_pipro;
324 
325  double r = tf * rnd->RndFsi().Rndm();
326 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
327  LOG("HAIntranuke", pDEBUG) << "r = " << r << " (max = " << tf << ")";
328 #endif
329  double cf=0; // current fraction
330  if(r < (cf += frac_cex )) return kIHAFtCEx; // cex
331  if(r < (cf += frac_elas )) return kIHAFtElas; // elas
332  if(r < (cf += frac_inel )) return kIHAFtInelas; // inelas
333  if(r < (cf += frac_abs )) return kIHAFtAbs; // abs
334  if(r < (cf += frac_pipro )) return kIHAFtPiProd; // pi prod
335 
336  LOG("HAIntranuke", pWARN)
337  << "No selection after going through all fates! "
338  << "Total fraction = " << tf << " (r = " << r << ")";
339  }
340  // handle kaons
341  else if (pdgc==kPdgKP || pdgc==kPdgKM) {
342  double frac_inel = fHadroData->Frac(pdgc, kIHAFtInelas, ke);
343  double frac_abs = fHadroData->Frac(pdgc, kIHAFtAbs, ke);
344  LOG("HAIntranuke", pDEBUG)
345  << "\n frac{" << INukeHadroFates::AsString(kIHAFtInelas) << "} = " << frac_inel
346  << "\n frac{" << INukeHadroFates::AsString(kIHAFtAbs) << "} = " << frac_abs;
347  // compute total fraction (can be <1 if fates have been switched off)
348  double tf = frac_inel +
349  frac_abs;
350  double r = tf * rnd->RndFsi().Rndm();
351 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
352  LOG("HAIntranuke", pDEBUG) << "r = " << r << " (max = " << tf << ")";
353 #endif
354  double cf=0; // current fraction
355  if(r < (cf += frac_inel )) return kIHAFtInelas; // inelas
356  if(r < (cf += frac_abs )) return kIHAFtAbs; // abs
357  }
358  }//iterations
359 
360  return kIHAFtUndefined;
361 }
TRandom3 & RndFsi(void) const
rnd number generator used by intranuclear cascade monte carlos
Definition: RandomGen.h:60
double fNucleonFracInelScale
Definition: Intranuke.h:126
INukeHadroData * fHadroData
a collection of h+N,h+A data & calculations
Definition: Intranuke.h:93
static RandomGen * Instance()
Access instance.
Definition: RandomGen.cxx:79
static const double MeV
Definition: Units.h:130
double fNucleonFracCExScale
Definition: Intranuke.h:124
double fNucleonFracAbsScale
Definition: Intranuke.h:127
Timing fit.
A singleton holding random number generator classes. All random number generation in GENIE should tak...
Definition: RandomGen.h:30
double Frac(int hpdgc, INukeFateHA_t fate, double ke) const
int Pdg(void) const
Definition: GHepParticle.h:64
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:97
double fPionFracPiProdScale
Definition: Intranuke.h:122
const int kPdgKM
Definition: PDGCodes.h:150
double fPionFracInelScale
Definition: Intranuke.h:120
const int kPdgKP
Definition: PDGCodes.h:149
const int kPdgPiP
Definition: PDGCodes.h:135
const int kPdgPi0
Definition: PDGCodes.h:137
#define pINFO
Definition: Messenger.h:63
#define pWARN
Definition: Messenger.h:61
double KinE(bool mass_from_pdg=false) const
Get kinetic energy.
double fPionFracAbsScale
Definition: Intranuke.h:121
static string AsString(INukeFateHN_t fate)
double fPionFracCExScale
Definition: Intranuke.h:118
double fPionFracElasScale
Definition: Intranuke.h:119
TRandom3 r(0)
const int kPdgPiM
Definition: PDGCodes.h:136
const int kPdgProton
Definition: PDGCodes.h:65
double fNucleonFracPiProdScale
Definition: Intranuke.h:128
static const unsigned int kRjMaxIterations
Definition: Controls.h:27
const int kPdgNeutron
Definition: PDGCodes.h:67
double fNucleonFracElasScale
Definition: Intranuke.h:125
#define pDEBUG
Definition: Messenger.h:64
bool HAIntranuke::HandleCompoundNucleus ( GHepRecord ev,
GHepParticle p,
int  mom 
) const
privatevirtual

Implements genie::Intranuke.

Definition at line 1372 of file HAIntranuke.cxx.

1374 {
1375  // only relevant for hN mode
1376  return false;
1377 }
virtual const AlgId& genie::Algorithm::Id ( void  ) const
inlinevirtualinherited

Get algorithm ID.

Definition at line 98 of file Algorithm.h.

References genie::Algorithm::fID.

Referenced by genie::KineGeneratorWithCache::AccessCacheBranch(), genie::QELEventGeneratorSM::AccessCacheBranch2(), genie::QELEventGeneratorSM::AccessCacheBranchDiffv(), genie::InteractionListAssembler::AssembleInteractionList(), genie::XSecAlgorithmMap::BuildMap(), genie::InteractionGeneratorMap::BuildMap(), genie::XSecSplineList::BuildSplineKey(), genie::DISXSec::CacheBranchName(), genie::ReinSehgalRESXSecWithCache::CacheBranchName(), genie::DMDISXSec::CacheBranchName(), genie::ReinSehgalRESXSecWithCacheFast::CacheBranchName(), genie::Algorithm::Compare(), genie::RESKinematicsGenerator::ComputeMaxXSec(), genie::COHElKinematicsGenerator::ComputeMaxXSec(), genie::SKKinematicsGenerator::ComputeMaxXSec(), genie::COHKinematicsGenerator::ComputeMaxXSec(), genie::Algorithm::Configure(), genie::GEVGDriver::CreateSplines(), genie::QPMDISPXSec::DISRESJoinSuppressionFactor(), genie::QPMDMDISPXSec::DMDISRESJoinSuppressionFactor(), genie::AlgConfigPool::FindRegistry(), genie::AlgFactory::ForceReconfiguration(), genie::GEVGDriver::GenerateEvent(), GetAlgorithms(), genie::LwlynSmithQELCCPXSec::Integral(), genie::NievesQELCCPXSec::Integral(), genie::COHXSec::Integrate(), genie::QPMDISPXSec::LoadConfig(), genie::QPMDMDISPXSec::LoadConfig(), genie::EventGenerator::LoadConfig(), genie::EventGeneratorListAssembler::LoadGenerator(), main(), genie::COHKinematicsGenerator::MaxXSec_AlvarezRuso(), genie::XSecAlgorithmMap::Print(), genie::InteractionGeneratorMap::Print(), genie::AlgFactory::Print(), genie::COHHadronicSystemGenerator::ProcessEventRecord(), genie::COHPrimaryLeptonGenerator::ProcessEventRecord(), genie::COHKinematicsGenerator::ProcessEventRecord(), genie::MECGenerator::ProcessEventRecord(), genie::EventGenerator::ProcessEventRecord(), genie::KNOPythiaHadronization::SelectHadronizer(), TestPythiaTauDecays(), and genie::GEVGDriver::UseSplines().

98 { return fID; }
AlgId fID
algorithm name and configuration set
Definition: Algorithm.h:156
void HAIntranuke::Inelastic ( GHepRecord ev,
GHepParticle p,
INukeFateHA_t  fate 
) const
private

Definition at line 704 of file HAIntranuke.cxx.

References genie::GHepRecord::AddParticle(), genie::INukeHadroFates::AsString(), demo5::c1, demo5::c2, chisquared::c3, genie::Intranuke::fDoFermi, genie::Intranuke::fFermiFac, genie::Intranuke::fFermiMomentum, genie::Intranuke::fHadroData, genie::PDGLibrary::Find(), genie::GHepParticle::FirstMother(), genie::Intranuke::fNuclmodel, genie::Intranuke::fNucRmvE, util::frac(), genie::Intranuke::fRemnA, genie::Intranuke::fRemnP4, genie::Intranuke::fRemnZ, genie::NuclearModelI::GenerateNucleon(), MECModelEnuComparisons::i, genie::PDGLibrary::Instance(), genie::RandomGen::Instance(), makeTrainCVSamples::int, genie::INukeHadroData::IntBounce(), genie::pdg::IsKaon(), genie::pdg::IsNeutron(), genie::pdg::IsNeutronOrProton(), genie::pdg::IsPion(), genie::pdg::IsProton(), genie::kIHAFtAbs, genie::kIHAFtPiProd, genie::kIHNFtAbs, genie::kIMdHA, genie::GHepParticle::KinE(), genie::kIStNucleonClusterTarget, genie::kIStStableFinalState, genie::kPdgCompNuclCluster, genie::kPdgKM, genie::kPdgKP, genie::kPdgNeutron, genie::kPdgPi0, genie::kPdgPiM, genie::kPdgPiP, genie::kPdgProton, kPi, genie::GHepParticle::LastMother(), parse_dependency_file_t::list, LOG, cet::sqlite::max(), genie::units::MeV, genie::NuclearModelI::Momentum3(), genie::GHepParticle::Name(), nd, ns, PandAna.Demos.pi0_spectra::p0, plot_validation_datamc::p1, plot_validation_datamc::p2, make_associated_cosmic_defs::p3, genie::GHepParticle::P4(), make_associated_cosmic_defs::p4, pDEBUG, genie::GHepParticle::Pdg(), genie::utils::intranuke::PhaseSpaceDecay(), pINFO, genie::utils::intranuke::PionProduction(), pNOTICE, genie::PDGCodeList::push_back(), pWARN, generate_hists::rnd, genie::RandomGen::RndFsi(), genie::GHepParticle::SetFirstMother(), genie::GHepParticle::SetLastMother(), genie::GHepParticle::SetMomentum(), genie::GHepParticle::SetPdgCode(), genie::exceptions::INukeException::SetReason(), genie::GHepParticle::SetStatus(), gen_flatrecord::size, ana::Sqrt(), log::success(), getGoodRuns4SAM::t1, t2, genie::GHepRecord::TargetNucleus(), genie::utils::intranuke::TwoBodyKinematics(), x1, submit_syst::x2, and genie::GHepParticle::X4().

Referenced by SimulateHadronicFinalStateKinematics().

706 {
707 
708  // Aaron Meyer (05/25/10)
709  //
710  // Called to handle all absorption and pi production reactions
711  //
712  // Nucleons -> Reaction approximated by exponential decay in p+n (sum) space,
713  // gaussian in p-n (difference) space
714  // -fit to hN simulations p C, Fe, Pb at 200, 800 MeV
715  // -get n from isospin, np-nn smaller by 2
716  // Pions -> Reaction approximated with a modified gaussian in p+n space,
717  // normal gaussian in p-n space
718  // -based on fits to multiplicity distributions of hN model
719  // for pi+ C, Fe, Pb at 250, 500 MeV
720  // -fit sum and diff of nn, np to Gaussian
721  // -get pi0 from isospin, np-nn smaller by 2
722  // -get pi- from isospin, np-nn smaller by 4
723  // -add 2-body absorption to better match McKeown data
724  // Kaons -> no guidance, use same code as pions.
725  //
726  // Normally distributed random number generated using Box-Muller transformation
727  //
728  // Pion production reactions rescatter pions in nucleus, otherwise unchanged from
729  // older versions of GENIE
730  //
731 
732 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
733  LOG("HAIntranuke", pDEBUG)
734  << "Inelastic() is invoked for a : " << p->Name()
735  << " whose fate is : " << INukeHadroFates::AsString(fate);
736 #endif
737 
738  bool allow_dup = true;
739  PDGCodeList list(allow_dup); // list of final state particles
740 
741  // only absorption/pipro fates allowed
742  if (fate == kIHAFtPiProd) {
743 
744  GHepParticle s1(*p);
745  GHepParticle s2(*p);
746  GHepParticle s3(*p);
747 
750 
751  if (success){
752  LOG ("HAIntranuke",pINFO) << " successful pion production fate";
753  // set status of particles and conserve charge/baryon number
754  s1.SetStatus(kIStStableFinalState); //should be redundant
755  // if (pdg::IsPion(s2->Pdg())) s2->SetStatus(kIStHadronInTheNucleus);
756  s2.SetStatus(kIStStableFinalState);
757  // if (pdg::IsPion(s3->Pdg())) s3->SetStatus(kIStHadronInTheNucleus);
758  s3.SetStatus(kIStStableFinalState);
759 
760  ev->AddParticle(s1);
761  ev->AddParticle(s2);
762  ev->AddParticle(s3);
763 
764  return;
765  }
766  else {
767  LOG("HAIntranuke", pNOTICE) << "Error: could not create pion production final state";
769  exception.SetReason("PionProduction kinematics failed - retry kinematics");
770  throw exception;
771  }
772  }
773 
774  else if (fate==kIHAFtAbs)
775 // tuned for pions - mixture of 2-body and many-body
776 // use same for kaons as there is no guidance
777  {
778  // Instances for reference
779  PDGLibrary * pLib = PDGLibrary::Instance();
781 
782  double ke = p->KinE() / units::MeV;
783  int pdgc = p->Pdg();
784 
785  if (fRemnA<2)
786  {
787  LOG("HAIntranuke", pNOTICE) << "stop propagation - could not create absorption final state: too few particles";
789  ev->AddParticle(*p);
790  return;
791  }
792  if (fRemnZ<1 && (pdgc==kPdgPiM || pdgc==kPdgKM))
793  {
794  LOG("HAIntranuke", pNOTICE) << "stop propagation - could not create absorption final state: Pi- or K- cannot be absorbed by only neutrons";
796  ev->AddParticle(*p);
797  return;
798  }
799  if (fRemnA-fRemnZ<1 && (pdgc==kPdgPiP || pdgc==kPdgKP))
800  {
801  LOG("HAIntranuke", pINFO) << "stop propagation - could not create absorption final state: Pi+ or K+ cannot be absorbed by only protons";
803  ev->AddParticle(*p);
804  return;
805  }
806 
807  // for now, empirical split between multi-nucleon absorption and pi d -> N N
808  //
809  // added 03/21/11 - Aaron Meyer
810  //
811  if (pdg::IsPion(pdgc) && rnd->RndFsi().Rndm()<1.14*(.903-0.00189*fRemnA)*(1.35-0.00467*ke))
812  { // pi d -> N N, probability determined empirically with McKeown data
813 
814  INukeFateHN_t fate_hN=kIHNFtAbs;
815  int t1code,t2code,scode,s2code;
816  double ppcnt = (double) fRemnZ / (double) fRemnA; // % of protons
817 
818  // choose target nucleon
819  // -- fates weighted by values from Engel, Mosel...
820  if (pdgc==kPdgPiP) {
821  double Prob_pipd_pp=2.*ppcnt*(1.-ppcnt);
822  double Prob_pipnn_pn=.083*(1.-ppcnt)*(1.-ppcnt);
823  if (rnd->RndFsi().Rndm()*(Prob_pipd_pp+Prob_pipnn_pn)<Prob_pipd_pp){
824  t1code=kPdgNeutron; t2code=kPdgProton;
825  scode=kPdgProton; s2code=kPdgProton;}
826  else{
827  t1code=kPdgNeutron; t2code=kPdgNeutron;
828  scode=kPdgProton; s2code=kPdgNeutron;}
829  }
830  if (pdgc==kPdgPiM) {
831  double Prob_pimd_nn=2.*ppcnt*(1.-ppcnt);
832  double Prob_pimpp_pn=.083*ppcnt*ppcnt;
833  if (rnd->RndFsi().Rndm()*(Prob_pimd_nn+Prob_pimpp_pn)<Prob_pimd_nn){
834  t1code=kPdgProton; t2code=kPdgNeutron;
835  scode=kPdgNeutron; s2code=kPdgNeutron; }
836  else{
837  t1code=kPdgProton; t2code=kPdgProton;
838  scode=kPdgProton; s2code=kPdgNeutron;}
839  }
840  else { // pi0
841  double Prob_pi0d_pn=0.88*ppcnt*(1.-ppcnt); // 2 * .44
842  double Prob_pi0pp_pp=.14*ppcnt*ppcnt;
843  double Prob_pi0nn_nn=.14*(1.-ppcnt)*(1.-ppcnt);
844  if (rnd->RndFsi().Rndm()*(Prob_pi0d_pn+Prob_pi0pp_pp+Prob_pi0nn_nn)<Prob_pi0d_pn){
845  t1code=kPdgNeutron; t2code=kPdgProton;
846  scode=kPdgNeutron; s2code=kPdgProton; }
847  else if (rnd->RndFsi().Rndm()*(Prob_pi0d_pn+Prob_pi0pp_pp+Prob_pi0nn_nn)<(Prob_pi0d_pn+Prob_pi0pp_pp)){
848  t1code=kPdgProton; t2code=kPdgProton;
849  scode=kPdgProton; s2code=kPdgProton; }
850  else {
851  t1code=kPdgNeutron; t2code=kPdgNeutron;
852  scode=kPdgNeutron; s2code=kPdgNeutron; }
853  }
854  LOG("HAIntranuke",pINFO) << "choose 2 body absorption, probe, fs = " << pdgc <<" "<< scode <<" "<<s2code;
855  // assign proper masses
856  //double M1 = pLib->Find(pdgc) ->Mass();
857  double M2_1 = pLib->Find(t1code)->Mass();
858  double M2_2 = pLib->Find(t2code)->Mass();
859  //double M2 = M2_1 + M2_2;
860  double M3 = pLib->Find(scode) ->Mass();
861  double M4 = pLib->Find(s2code)->Mass();
862 
863  // handle fermi momentum
864  double E2_1L, E2_2L;
865  TVector3 tP2_1L, tP2_2L;
866  //TLorentzVector dNucl_P4;
867  Target target(ev->TargetNucleus()->Pdg());
868  if(fDoFermi)
869  {
870  target.SetHitNucPdg(t1code);
872  tP2_1L=fFermiFac * fNuclmodel->Momentum3();
873  E2_1L = TMath::Sqrt(tP2_1L.Mag2() + M2_1*M2_1);
874 
875  target.SetHitNucPdg(t2code);
877  tP2_2L=fFermiFac * fNuclmodel->Momentum3();
878  E2_2L = TMath::Sqrt(tP2_2L.Mag2() + M2_2*M2_2);
879 
880  //dNucl_P4=TLorentzVector(tP2_1L+tP2_2L,E2_1L+E2_2L);
881  }
882  else
883  {
884  tP2_1L.SetXYZ(0.0, 0.0, 0.0);
885  E2_1L = M2_1;
886 
887  tP2_2L.SetXYZ(0.0, 0.0, 0.0);
888  E2_2L = M2_2;
889  }
890  TLorentzVector dNucl_P4=TLorentzVector(tP2_1L+tP2_2L,E2_1L+E2_2L);
891 
892  double E2L = E2_1L + E2_2L;
893 
894  // adjust p to reflect scattering
895  // get random scattering angle
896  double C3CM = fHadroData->IntBounce(p,t1code,scode,fate_hN);
897  if (C3CM<-1.)
898  {
899  LOG("HAIntranuke", pWARN) << "Inelastic() failed: IntBounce returned bad angle - try again";
901  exception.SetReason("unphysical angle for hN scattering");
902  throw exception;
903  return;
904  }
905 
906  TLorentzVector t4P1L,t4P2L,t4P3L,t4P4L;
907  t4P1L=*p->P4();
908  t4P2L=TLorentzVector(TVector3(tP2_1L+tP2_2L),E2L);
909  double bindE=0.075; // set to fit McKeown data
910  //double bindE=0.0;
911  if (utils::intranuke::TwoBodyKinematics(M3,M4,t4P1L,t4P2L,t4P3L,t4P4L,C3CM,fRemnP4,bindE))
912  {
913  if (pdgc==kPdgPiP || pdgc==kPdgKP) fRemnZ++;
914  if (pdgc==kPdgPiM || pdgc==kPdgKM) fRemnZ--;
915  if (t1code==kPdgProton) fRemnZ--;
916  if (t2code==kPdgProton) fRemnZ--;
917  fRemnA-=2;
918 
919  fRemnP4-=dNucl_P4;
920 
921  // create t particles w/ appropriate momenta, code, and status
922  // Set target's mom to be the mom of the hadron that was cloned
923  GHepParticle* t1 = new GHepParticle(*p);
924  GHepParticle* t2 = new GHepParticle(*p);
925  t1->SetFirstMother(p->FirstMother());
926  t1->SetLastMother(p->LastMother());
927  t2->SetFirstMother(p->FirstMother());
928  t2->SetLastMother(p->LastMother());
929 
930  // adjust p to reflect scattering
931  t1->SetPdgCode(scode);
932  t1->SetMomentum(t4P3L);
933 
934  t2->SetPdgCode(s2code);
935  t2->SetMomentum(t4P4L);
936 
939 
940  ev->AddParticle(*t1);
941  ev->AddParticle(*t2);
942  delete t1;
943  delete t2;
944 
945  return;
946  }
947  else
948  {
949  LOG("HAIntranuke", pNOTICE) << "Inelastic in hA failed calling TwoBodyKineamtics";
951  exception.SetReason("Pion absorption kinematics through TwoBodyKinematics failed");
952  throw exception;
953 
954  }
955 
956  } // end pi d -> N N
957  else // multi-nucleon
958  {
959 
960  // declare some parameters for double gaussian and determine values chosen
961  // parameters for proton and pi+, others come from isospin transformations
962 
963  double ns0=0; // mean - sum of nucleons
964  double nd0=0; // mean - difference of nucleons
965  double Sig_ns=0; // std dev - sum
966  double Sig_nd=0; // std dev - diff
967  double gam_ns=0; // exponential decay rate (for nucleons)
968 
969  if ( pdg::IsNeutronOrProton (pdgc) ) // nucleon probe
970  {
971  // antisymmetric about Z=N
972  if (fRemnA-fRemnZ > fRemnZ)
973  nd0 = 135.227 * TMath::Exp(-7.124*(fRemnA-fRemnZ)/double(fRemnA)) - 2.762;
974  else
975  nd0 = -135.227 * TMath::Exp(-7.124* fRemnZ /double(fRemnA)) + 4.914;
976 
977  Sig_nd = 2.034 + fRemnA * 0.007846;
978 
979  double c1 = 0.041 + ke * 0.0001525;
980  double c2 = -0.003444 - ke * 0.00002324;
981 //change last factor from 30 to 15 so that gam_ns always larger than 0
982 //add check to be certain
983  double c3 = 0.064 - ke * 0.000015;
984  gam_ns = c1 * TMath::Exp(c2*fRemnA) + c3;
985  if(gam_ns<0.002) gam_ns = 0.002;
986  //gam_ns = 10.;
987  LOG("HAIntranuke", pINFO) << "nucleon absorption";
988  LOG("HAIntranuke", pINFO) << "--> mean diff distr = " << nd0 << ", stand dev = " << Sig_nd;
989  // LOG("HAIntranuke", pINFO) << "--> mean sum distr = " << ns0 << ", Stand dev = " << Sig_ns;
990  LOG("HAIntranuke", pINFO) << "--> gam_ns = " << gam_ns;
991  }
992  else if ( pdgc==kPdgPiP || pdgc==kPdgPi0 || pdgc==kPdgPiM) //pion probe
993  {
994  ns0 = .0001*(1.+ke/250.) * (fRemnA-10)*(fRemnA-10) + 3.5;
995  nd0 = (1.+ke/250.) - ((fRemnA/200.)*(1. + 2.*ke/250.));
996  Sig_ns = (10. + 4. * ke/250.)*TMath::Power(fRemnA/250.,0.9); //(1. - TMath::Exp(-0.02*fRemnA));
997  Sig_nd = 4*(1 - TMath::Exp(-0.03*ke));
998  LOG("HAIntranuke", pINFO) << "pion absorption";
999  LOG("HAIntranuke", pINFO) << "--> mean diff distr = " << nd0 << ", stand dev = " << Sig_nd;
1000  LOG("HAIntranuke", pINFO) << "--> mean sum distr = " << ns0 << ", Stand dev = " << Sig_ns;
1001  }
1002  else if (pdgc==kPdgKP || pdgc==kPdgKM) // kaon probe
1003  {
1004  ns0 = (rnd->RndFsi().Rndm()>0.5?3:2);
1005  nd0 = 1.;
1006  Sig_ns = 0.1;
1007  Sig_nd = 0.1;
1008  LOG("HAIntranuke", pINFO) << "kaon absorption - set ns, nd later";
1009  // LOG("HAIntranuke", pINFO) << "--> mean diff distr = " << nd0 << ", stand dev = " << Sig_nd;
1010  // LOG("HAIntranuke", pINFO) << "--> mean sum distr = " << ns0 << ", Stand dev = " << Sig_ns;
1011  }
1012  else
1013  {
1014  LOG("HAIntranuke", pWARN) << "Inelastic() cannot handle absorption reaction for " << p->Name();
1015  }
1016 
1017  // account for different isospin
1018  if (pdgc==kPdgPi0 || pdgc==kPdgNeutron) nd0-=2.;
1019  if (pdgc==kPdgPiM) nd0-=4.;
1020 
1021  int iter=0; // counter
1022  int np=0,nn=0; // # of p, # of n
1023  bool not_done=true;
1024  double u1 = 0, u2 = 0;
1025 
1026  while (not_done)
1027  {
1028  // infinite loop check
1029  if (iter>=100) {
1030  LOG("HAIntranuke", pNOTICE) << "Error: could not choose absorption final state";
1031  LOG("HAIntranuke", pNOTICE) << "--> mean diff distr = " << nd0 << ", stand dev = " << Sig_nd;
1032  LOG("HAIntranuke", pNOTICE) << "--> mean sum distr = " << ns0 << ", Stand dev = " << Sig_ns;
1033  LOG("HAIntranuke", pNOTICE) << "--> gam_ns = " << gam_ns;
1034  LOG("HAIntranuke", pNOTICE) << "--> A = " << fRemnA << ", Z = " << fRemnZ << ", Energy = " << ke;
1036  exception.SetReason("Absorption choice of # of p,n failed");
1037  throw exception;
1038  }
1039  //here??
1040 
1041  // Box-Muller transform
1042  // Takes two uniform distribution random variables on (0,1]
1043  // Creates two normally distributed random variables on (0,inf)
1044 
1045  u1 = rnd->RndFsi().Rndm(); // uniform random variable 1
1046  u2 = rnd->RndFsi().Rndm(); // " " 2
1047  if (u1==0) u1 = rnd->RndFsi().Rndm();
1048  if (u2==0) u2 = rnd->RndFsi().Rndm(); // Just in case
1049 
1050  // normally distributed random variable
1051  double x2 = TMath::Sqrt(-2*TMath::Log(u1))*TMath::Sin(2*kPi*u2);
1052 
1053  double ns = 0;
1054 
1055  if ( pdg::IsNeutronOrProton (pdgc) ) //nucleon probe
1056  {
1057  ns = -TMath::Log(rnd->RndFsi().Rndm())/gam_ns; // exponential random variable
1058  }
1059  if ( pdg::IsKaon (pdgc) ) //charged kaon probe - either 2 or 3 nucleons to stay simple
1060  {
1061  ns = (rnd->RndFsi().Rndm()<0.5?2:3);
1062  }
1063  else if ( pdgc==kPdgPiP || pdgc==kPdgPi0 || pdgc==kPdgPiM) //pion probe
1064  {
1065  // Pion fit for sum takes for xs*exp((xs-x0)^2 / 2*sig_xs0)
1066  // Find random variable by first finding gaussian random variable
1067  // then throwing the value against a linear P.D.F.
1068  //
1069  // max is the maximum value allowed for the random variable (10 std + mean)
1070  // minimum allowed value is 0
1071 
1072  double max = ns0 + Sig_ns * 10;
1073  if(max>fRemnA) max=fRemnA;
1074  double x1 = 0;
1075  bool not_found = true;
1076  int iter2 = 0;
1077 
1078  while (not_found)
1079  {
1080  // infinite loop check
1081  if (iter2>=100)
1082  {
1083  LOG("HAIntranuke", pNOTICE) << "Error: stuck in random variable loop for ns";
1084  LOG("HAIntranuke", pNOTICE) << "--> mean of sum parent distr = " << ns0 << ", Stand dev = " << Sig_ns;
1085  LOG("HAIntranuke", pNOTICE) << "--> A = " << fRemnA << ", Z = " << fRemnZ << ", Energy = " << ke;
1086 
1088  exception.SetReason("Random number generator for choice of #p,n final state failed - unusual - redo kinematics");
1089  throw exception;
1090  }
1091 
1092  // calculate exponential random variable
1093  u1 = rnd->RndFsi().Rndm();
1094  u2 = rnd->RndFsi().Rndm();
1095  if (u1==0) u1 = rnd->RndFsi().Rndm();
1096  if (u2==0) u2 = rnd->RndFsi().Rndm();
1097  x1 = TMath::Sqrt(-2*TMath::Log(u1))*TMath::Cos(2*kPi*u2);
1098 
1099  ns = ns0 + Sig_ns * x1;
1100  if ( ns>max || ns<0 ) {iter2++; continue;}
1101  else if ( rnd->RndFsi().Rndm() > (ns/max) ) {iter2++; continue;}
1102  else {
1103  // accept this sum value
1104  not_found=false;
1105  }
1106  } //while(not_found)
1107  }//else pion
1108 
1109  double nd = nd0 + Sig_nd * x2; // difference (p-n) for both pion, nucleon probe
1110  if (pdgc==kPdgKP) // special for KP
1111  { if (ns==2) nd=0;
1112  if (ns>2) nd=1; }
1113 
1114  np = int((ns+nd)/2.+.5); // Addition of .5 for rounding correction
1115  nn = int((ns-nd)/2.+.5);
1116 
1117  LOG("HAIntranuke", pINFO) << "ns = "<<ns<<", nd = "<<nd<<", np = "<<np<<", nn = "<<nn;
1118  //LOG("HAIntranuke", pNOTICE) << "RemA = "<<fRemnA<<", RemZ = "<<fRemnZ<<", probe = "<<pdgc;
1119 
1120  /*if ((ns+nd)/2. < 0 || (ns-nd)/2. < 0) {iter++; continue;}
1121  else */
1122  //check for problems befor executing phase space 'decay'
1123  if (np < 0 || nn < 0 ) {iter++; continue;}
1124  else if (np + nn < 2. ) {iter++; continue;}
1125  else if ((np + nn == 2.) && pdg::IsNeutronOrProton (pdgc)) {iter++; continue;}
1126  else if (np > fRemnZ + ((pdg::IsProton(pdgc) || pdgc==kPdgPiP || pdgc==kPdgKP)?1:0)
1127  - ((pdgc==kPdgPiM || pdgc==kPdgKM)?1:0)) {iter++; continue;}
1128  else if (nn > fRemnA-fRemnZ + ((pdg::IsNeutron(pdgc)||pdgc==kPdgPiM||pdgc==kPdgKM)?1:0)
1129  - ((pdgc==kPdgPiP||pdgc==kPdgKP)?1:0)) {iter++; continue;}
1130  else {
1131  not_done=false; //success
1132  LOG("HAIntranuke",pINFO) << "success, iter = " << iter << " np, nn = " << np << " " << nn;
1133  if (np+nn>86) // too many particles, scale down
1134  {
1135  double frac = 85./double(np+nn);
1136  np = int(np*frac);
1137  nn = int(nn*frac);
1138  }
1139 
1140  if ( (np==fRemnZ +((pdg::IsProton (pdgc)||pdgc==kPdgPiP||pdgc==kPdgKP)?1:0)-(pdgc==kPdgPiM||pdgc==kPdgKM?1:0))
1141  &&(nn==fRemnA-fRemnZ+((pdg::IsNeutron(pdgc)||pdgc==kPdgPiM||pdgc==kPdgKM)?1:0)-(pdgc==kPdgPiP||pdgc==kPdgKP?1:0)) )
1142  { // leave at least one nucleon in the nucleus to prevent excess momentum
1143  if (rnd->RndFsi().Rndm()<np/(double)(np+nn)) np--;
1144  else nn--;
1145  }
1146 
1147  LOG("HAIntranuke", pNOTICE) << "Final state chosen; # protons : "
1148  << np << ", # neutrons : " << nn;
1149  }
1150  } //while(not_done)
1151 
1152  // change remnants to reflect probe
1153  if ( pdgc==kPdgProton || pdgc==kPdgPiP || pdgc==kPdgKP) fRemnZ++;
1154  if ( pdgc==kPdgPiM || pdgc==kPdgKM) fRemnZ--;
1155  if ( pdg::IsNeutronOrProton (pdgc) ) fRemnA++;
1156 
1157  // PhaseSpaceDecay forbids anything over 18 particles
1158  //
1159  // If over 18 particles, split into 5 groups and perform decay on each group
1160  // Anything over 85 particles reduced to 85 in previous step
1161  //
1162  // 4 of the nucleons are used as "probes" as well as original probe,
1163  // with each one sharing 1/5 of the total incident momentum
1164  //
1165  // Note: choosing 5 groups and distributing momentum evenly was arbitrary
1166  // Needs to be revised later
1167 
1168  if ((np+nn)>18)
1169  {
1170  // code lists
1171  PDGCodeList list0(allow_dup);
1172  PDGCodeList list1(allow_dup);
1173  PDGCodeList list2(allow_dup);
1174  PDGCodeList list3(allow_dup);
1175  PDGCodeList list4(allow_dup);
1176  PDGCodeList* listar[5] = {&list0, &list1, &list2, &list3, &list4};
1177 
1178  //set up HadronClusters
1179  // simple for now, each (of 5) hadron cluster has 1/5 of mom and KE
1180 
1181  double probM = pLib->Find(pdgc) ->Mass();
1182  TVector3 pP3 = p->P4()->Vect() * (1./5.);
1183  double probKE = p->P4()->E() -probM;
1184  double clusKE = probKE * (1./5.);
1185  TLorentzVector clusP4(pP3,clusKE); //no mass
1186 
1187  TLorentzVector X4(*p->X4());
1189 
1190  int mom = p->FirstMother();
1191 
1192  GHepParticle * p0 = new GHepParticle(kPdgCompNuclCluster,ist, mom,-1,-1,-1,clusP4,X4);
1193  GHepParticle * p1 = new GHepParticle(kPdgCompNuclCluster,ist, mom,-1,-1,-1,clusP4,X4);
1194  GHepParticle * p2 = new GHepParticle(kPdgCompNuclCluster,ist, mom,-1,-1,-1,clusP4,X4);
1195  GHepParticle * p3 = new GHepParticle(kPdgCompNuclCluster,ist, mom,-1,-1,-1,clusP4,X4);
1196  GHepParticle * p4 = new GHepParticle(kPdgCompNuclCluster,ist, mom,-1,-1,-1,clusP4,X4);
1197 
1198  // To conserve 4-momenta
1199  // fRemnP4 -= probP4 + protP4*np_p + neutP4*(4-np_p) - *p->P4();
1200  fRemnP4 -= 5.*clusP4 - *p->P4();
1201 
1202  for (int i=0;i<(np+nn);i++)
1203  {
1204  if (i<np)
1205  {
1206  listar[i%5]->push_back(kPdgProton);
1207  fRemnZ--;
1208  }
1209  else listar[i%5]->push_back(kPdgNeutron);
1210  fRemnA--;
1211  }
1212  for (int i=0;i<5;i++)
1213  {
1214  LOG("HAIntranuke", pINFO) << "List" << i << " size: " << listar[i]->size();
1215  if (listar[i]->size() <2)
1216  {
1218  exception.SetReason("too few particles for Phase Space decay - try again");
1219  throw exception;
1220  }
1221  }
1222 
1223  // commented out to better fit with absorption reactions
1224  // Add the fermi energy of the nucleons to the phase space
1225  /*if(fDoFermi)
1226  {
1227  GHepParticle* p_ar[5] = {cl, p1, p2, p3, p4};
1228  for (int i=0;i<5;i++)
1229  {
1230  Target target(ev->TargetNucleus()->Pdg());
1231  TVector3 pBuf = p_ar[i]->P4()->Vect();
1232  double mBuf = p_ar[i]->Mass();
1233  double eBuf = TMath::Sqrt(pBuf.Mag2() + mBuf*mBuf);
1234  TLorentzVector tSum(pBuf,eBuf);
1235  double mSum = 0.0;
1236  vector<int>::const_iterator pdg_iter;
1237  for(pdg_iter=++(listar[i]->begin());pdg_iter!=listar[i]->end();++pdg_iter)
1238  {
1239  target.SetHitNucPdg(*pdg_iter);
1240  fNuclmodel->GenerateNucleon(target);
1241  mBuf = pLib->Find(*pdg_iter)->Mass();
1242  mSum += mBuf;
1243  pBuf = fFermiFac * fNuclmodel->Momentum3();
1244  eBuf = TMath::Sqrt(pBuf.Mag2() + mBuf*mBuf);
1245  tSum += TLorentzVector(pBuf,eBuf);
1246  fRemnP4 -= TLorentzVector(pBuf,eBuf-mBuf);
1247  }
1248  TLorentzVector dP4 = tSum + TLorentzVector(TVector3(0,0,0),-mSum);
1249  p_ar[i]->SetMomentum(dP4);
1250  }
1251  }*/
1252 
1253  bool success1 = utils::intranuke::PhaseSpaceDecay(ev,p0,*listar[0],fRemnP4,fNucRmvE,kIMdHA);
1254  bool success2 = utils::intranuke::PhaseSpaceDecay(ev,p1,*listar[1],fRemnP4,fNucRmvE,kIMdHA);
1255  bool success3 = utils::intranuke::PhaseSpaceDecay(ev,p2,*listar[2],fRemnP4,fNucRmvE,kIMdHA);
1256  bool success4 = utils::intranuke::PhaseSpaceDecay(ev,p3,*listar[3],fRemnP4,fNucRmvE,kIMdHA);
1257  bool success5 = utils::intranuke::PhaseSpaceDecay(ev,p4,*listar[4],fRemnP4,fNucRmvE,kIMdHA);
1258  if(success1 && success2 && success3 && success4 && success5)
1259  {
1260  LOG("HAIntranuke", pINFO)<<"Successful many-body absorption - n>=18";
1261  }
1262  else
1263  {
1264  // try to recover
1265  LOG("HAIntranuke", pWARN) << "PhaseSpace decay fails for HadrCluster- recovery likely incorrect - rethrow event";
1267  ev->AddParticle(*p);
1268  fRemnA+=np+nn;
1269  fRemnZ+=np;
1270  if ( pdgc==kPdgProton || pdgc==kPdgPiP ) fRemnZ--;
1271  if ( pdgc==kPdgPiM ) fRemnZ++;
1272  if ( pdg::IsNeutronOrProton (pdgc) ) fRemnA--;
1273  /* exceptions::INukeException exception;
1274  exception.SetReason("Phase space generation of absorption final state failed");
1275  throw exception;
1276  */
1277  }
1278 
1279  // delete cl;
1280  delete p0;
1281  delete p1;
1282  delete p2;
1283  delete p3;
1284  delete p4;
1285 
1286  }
1287  else // less than 18 particles pion
1288  {
1289  if (pdgc==kPdgKP) list.push_back(kPdgKP); //normally conserve strangeness
1290  if (pdgc==kPdgKM) list.push_back(kPdgKM);
1291  for (int i=0;i<np;i++)
1292  {
1293  list.push_back(kPdgProton);
1294  fRemnA--;
1295  fRemnZ--;
1296  }
1297  for (int i=0;i<nn;i++)
1298  {
1299  list.push_back(kPdgNeutron);
1300  fRemnA--;
1301  }
1302 
1303  // Library instance for reference
1304  //PDGLibrary * pLib = PDGLibrary::Instance();
1305 
1306  // commented out to better fit with absorption reactions
1307  // Add the fermi energy of the nucleons to the phase space
1308  /*if(fDoFermi)
1309  {
1310  Target target(ev->TargetNucleus()->Pdg());
1311  TVector3 pBuf = p->P4()->Vect();
1312  double mBuf = p->Mass();
1313  double eBuf = TMath::Sqrt(pBuf.Mag2() + mBuf*mBuf);
1314  TLorentzVector tSum(pBuf,eBuf);
1315  double mSum = 0.0;
1316  vector<int>::const_iterator pdg_iter;
1317  for(pdg_iter=++(list.begin());pdg_iter!=list.end();++pdg_iter)
1318  {
1319  target.SetHitNucPdg(*pdg_iter);
1320  fNuclmodel->GenerateNucleon(target);
1321  mBuf = pLib->Find(*pdg_iter)->Mass();
1322  mSum += mBuf;
1323  pBuf = fFermiFac * fNuclmodel->Momentum3();
1324  eBuf = TMath::Sqrt(pBuf.Mag2() + mBuf*mBuf);
1325  tSum += TLorentzVector(pBuf,eBuf);
1326  fRemnP4 -= TLorentzVector(pBuf,eBuf-mBuf);
1327  }
1328  TLorentzVector dP4 = tSum + TLorentzVector(TVector3(0,0,0),-mSum);
1329  p->SetMomentum(dP4);
1330  }*/
1331 
1332  LOG("HAIntranuke", pDEBUG)
1333  << "Remnant nucleus (A,Z) = (" << fRemnA << ", " << fRemnZ << ")";
1334  LOG("HAIntranuke", pINFO) << " list size: " << np+nn;
1335  if (np+nn <2)
1336  {
1338  exception.SetReason("too few particles for Phase Space decay - try again");
1339  throw exception;
1340  }
1341  // GHepParticle * cl = new GHepParticle(*p);
1342  // cl->SetPdgCode(kPdgDecayNuclCluster);
1344  if (success)
1345  {
1346  LOG ("HAIntranuke",pINFO) << "Successful many-body absorption, n<=18";
1347  }
1348  else {
1349  // recover
1351  ev->AddParticle(*p);
1352  fRemnA+=np+nn;
1353  fRemnZ+=np;
1354  if ( pdgc==kPdgProton || pdgc==kPdgPiP ) fRemnZ--;
1355  if ( pdgc==kPdgPiM ) fRemnZ++;
1356  if ( pdg::IsNeutronOrProton (pdgc) ) fRemnA--;
1358  exception.SetReason("Phase space generation of absorption final state failed");
1359  throw exception;
1360  }
1361  }
1362  } // end multi-nucleon FS
1363  }
1364  else // not absorption/pipro
1365  {
1366  LOG("HAIntranuke", pWARN)
1367  << "Inelastic() can not handle fate: " << INukeHadroFates::AsString(fate);
1368  return;
1369  }
1370 }
double fFermiMomentum
whether or not particle collision is pauli blocked
Definition: Intranuke.h:112
const double kPi
An exception thrown by SimulateHadronState for kinematics problems. TwoBodyCollision/Kinematics used ...
void SetFirstMother(int m)
Definition: GHepParticle.h:133
bool IsPion(int pdgc)
Definition: PDGUtils.cxx:289
TRandom3 & RndFsi(void) const
rnd number generator used by intranuclear cascade monte carlos
Definition: RandomGen.h:60
const XML_Char * target
Definition: expat.h:268
INukeHadroData * fHadroData
a collection of h+N,h+A data & calculations
Definition: Intranuke.h:93
TLorentzVector fRemnP4
P4 of remnant system.
Definition: Intranuke.h:98
static RandomGen * Instance()
Access instance.
Definition: RandomGen.cxx:79
const TLorentzVector * P4(void) const
Definition: GHepParticle.h:79
Float_t x1[n_points_granero]
Definition: compare.C:5
enum genie::EGHepStatus GHepStatus_t
static const double MeV
Definition: Units.h:130
::xsd::cxx::tree::exception< char > exception
Definition: Database.h:225
void SetMomentum(const TLorentzVector &p4)
bool IsKaon(int pdgc)
Definition: PDGUtils.cxx:294
A singleton holding random number generator classes. All random number generation in GENIE should tak...
Definition: RandomGen.h:30
const TVector3 & Momentum3(void) const
Definition: NuclearModelI.h:67
A list of PDG codes.
Definition: PDGCodeList.h:33
int LastMother(void) const
Definition: GHepParticle.h:68
int Pdg(void) const
Definition: GHepParticle.h:64
int FirstMother(void) const
Definition: GHepParticle.h:67
c2
Definition: demo5.py:33
string Name(void) const
Name that corresponds to the PDG code.
bool IsNeutron(int pdgc)
Definition: PDGUtils.cxx:304
const int kPdgCompNuclCluster
Definition: PDGCodes.h:194
bool IsProton(int pdgc)
Definition: PDGUtils.cxx:299
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
enum genie::EINukeFateHN_t INukeFateHN_t
bool PhaseSpaceDecay(GHepRecord *ev, GHepParticle *p, const PDGCodeList &pdgv, TLorentzVector &RemnP4, double NucRmvE, EINukeMode mode=kIMdHA)
general phase space decay method
const int kPdgKM
Definition: PDGCodes.h:150
bool PionProduction(GHepRecord *ev, GHepParticle *p, GHepParticle *s1, GHepParticle *s2, GHepParticle *s3, int &RemnA, int &RemnZ, TLorentzVector &RemnP4, bool DoFermi, double FermiFac, double FermiMomentum, const NuclearModelI *Nuclmodel)
A Neutrino Interaction Target. Is a transparent encapsulation of quite different physical systems suc...
Definition: Target.h:41
const int kPdgKP
Definition: PDGCodes.h:149
bool TwoBodyKinematics(double M3, double M4, TLorentzVector tP1L, TLorentzVector tP2L, TLorentzVector &tP3L, TLorentzVector &tP4L, double C3CM, TLorentzVector &RemnP4, double bindE=0)
Definition: INukeUtils.cxx:754
const int kPdgPiP
Definition: PDGCodes.h:135
const int kPdgPi0
Definition: PDGCodes.h:137
virtual GHepParticle * TargetNucleus(void) const
Definition: GHepRecord.cxx:330
def success(message)
Definition: log.py:5
double frac(double x)
Fractional part.
#define pINFO
Definition: Messenger.h:63
int fRemnA
remnant nucleus A
Definition: Intranuke.h:96
double t2
static const double ns
Module that plots metrics from reconstructed cosmic ray data.
#define pWARN
Definition: Messenger.h:61
double fNucRmvE
binding energy to subtract from cascade nucleons
Definition: Intranuke.h:104
double IntBounce(const GHepParticle *p, int target, int s1, INukeFateHN_t fate)
double KinE(bool mass_from_pdg=false) const
Get kinetic energy.
static PDGLibrary * Instance(void)
Definition: PDGLibrary.cxx:43
void SetLastMother(int m)
Definition: GHepParticle.h:134
Singleton class to load & serve a TDatabasePDG.
Definition: PDGLibrary.h:30
Var Sqrt(const Var &v)
Use to take sqrt of a var.
Definition: Var.cxx:326
const TLorentzVector * X4(void) const
Definition: GHepParticle.h:80
void SetStatus(GHepStatus_t s)
Definition: GHepParticle.h:127
static string AsString(INukeFateHN_t fate)
bool IsNeutronOrProton(int pdgc)
Definition: PDGUtils.cxx:314
TParticlePDG * Find(int pdgc)
Definition: PDGLibrary.cxx:61
int fRemnZ
remnant nucleus Z
Definition: Intranuke.h:97
const int kPdgPiM
Definition: PDGCodes.h:136
bool fDoFermi
whether or not to do fermi mom.
Definition: Intranuke.h:113
virtual bool GenerateNucleon(const Target &) const =0
c1
Definition: demo5.py:24
virtual void AddParticle(const GHepParticle &p)
Definition: GHepRecord.cxx:535
double fFermiFac
testing parameter to modify fermi momentum
Definition: Intranuke.h:111
const int kPdgProton
Definition: PDGCodes.h:65
#define pNOTICE
Definition: Messenger.h:62
T max(sqlite3 *const db, std::string const &table_name, std::string const &column_name)
Definition: statistics.h:68
const int kPdgNeutron
Definition: PDGCodes.h:67
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:40
const NuclearModelI * fNuclmodel
nuclear model used to generate fermi momentum
Definition: Intranuke.h:95
void push_back(int pdg_code)
Definition: PDGCodeList.cxx:67
#define pDEBUG
Definition: Messenger.h:64
void SetPdgCode(int c)
void HAIntranuke::InelasticHA ( GHepRecord ev,
GHepParticle p,
INukeFateHA_t  fate 
) const
private

Definition at line 542 of file HAIntranuke.cxx.

References genie::GHepRecord::AddParticle(), genie::INukeHadroFates::AsString(), genie::GHepParticle::E(), genie::Intranuke::fDoFermi, genie::Intranuke::fFermiFac, genie::Intranuke::fHadroData, genie::Intranuke::fNuclmodel, genie::Intranuke::fRemnA, genie::Intranuke::fRemnP4, genie::Intranuke::fRemnZ, genie::NuclearModelI::GenerateNucleon(), genie::RandomGen::Instance(), genie::INukeHadroData::IntBounce(), genie::kIHAFtCEx, genie::kIHAFtInelas, genie::kIHNFtCEx, genie::kIHNFtElas, genie::kIMdHA, genie::GHepParticle::KinE(), genie::kIStStableFinalState, genie::kPdgNeutron, genie::kPdgPi0, genie::kPdgPiM, genie::kPdgPiP, genie::kPdgProton, LOG, Mag2(), genie::GHepParticle::Mass(), genie::NuclearModelI::Momentum3(), genie::GHepParticle::Name(), genie::GHepParticle::P4(), pDEBUG, genie::GHepParticle::Pdg(), pINFO, pNOTICE, genie::GHepRecord::Probe(), pWARN, genie::GHepParticle::Px(), genie::GHepParticle::Py(), genie::GHepParticle::Pz(), generate_hists::rnd, genie::RandomGen::RndFsi(), genie::GHepParticle::SetMomentum(), genie::GHepParticle::SetPdgCode(), genie::exceptions::INukeException::SetReason(), genie::GHepParticle::SetStatus(), ana::Sqrt(), log::success(), confusionMatrixTree::t, genie::GHepRecord::TargetNucleus(), and genie::utils::intranuke::TwoBodyCollision().

Referenced by SimulateHadronicFinalStateKinematics().

544 {
545  // charge exch and inelastic - scatters particle within nucleus, hA version
546  // each are treated as quasielastic, particle scatters off single nucleon
547 
548 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
549  LOG("HAIntranuke", pDEBUG)
550  << "InelasticHA() is invoked for a : " << p->Name()
551  << " whose fate is : " << INukeHadroFates::AsString(fate);
552 #endif
553  if(ev->Probe() ) {
554  LOG("HAIntranuke", pINFO) << " probe KE = " << ev->Probe()->KinE();
555  }
556  if(fate!=kIHAFtCEx && fate!=kIHAFtInelas)
557  {
558  LOG("HAIntranuke", pWARN)
559  << "InelasticHA() cannot handle fate: " << INukeHadroFates::AsString(fate);
560  return;
561  }
562 
563  // Random number generator
565 
566  // vars for incoming particle, target, and scattered pdg codes
567  int pcode = p->Pdg();
568  int tcode, scode, s2code;
569  double ppcnt = (double) fRemnZ / (double) fRemnA; // % of protons
570 
571  // Select a hadron fate in HN mode
572  INukeFateHN_t h_fate;
573  if (fate == kIHAFtCEx) h_fate = kIHNFtCEx;
574  else h_fate = kIHNFtElas;
575 
576  // Select a target randomly, weighted to #
577  // -- Unless, of course, the fate is CEx,
578  // -- in which case the target may be deterministic
579  // Also assign scattered particle code
580  if(fate==kIHAFtCEx)
581  {
582  if(pcode==kPdgPiP) {tcode = kPdgNeutron; scode = kPdgPi0; s2code = kPdgProton;}
583  else if(pcode==kPdgPiM) {tcode = kPdgProton; scode = kPdgPi0; s2code = kPdgNeutron;}
584  else if(pcode==kPdgPi0)
585  {
586  // for pi0
587  tcode = (rnd->RndFsi().Rndm()<=ppcnt)?(kPdgProton) :(kPdgNeutron);
588  scode = (tcode == kPdgProton) ?(kPdgPiP) :(kPdgPiM);
589  s2code = (tcode == kPdgProton) ?(kPdgNeutron):(kPdgProton);
590  }
591  else if(pcode==kPdgProton) {tcode = kPdgNeutron; scode = kPdgNeutron; s2code = kPdgProton;}
592  else if(pcode==kPdgNeutron){tcode = kPdgProton; scode = kPdgProton; s2code = kPdgNeutron;}
593  else
594  { LOG("HAIntranuke", pWARN) << "InelasticHA() cannot handle fate: "
596  << " for particle " << p->Name();
597  return;
598  }
599  }
600  else
601  {
602  tcode = (rnd->RndFsi().Rndm()<=ppcnt)?(kPdgProton):(kPdgNeutron);
603  // if(pcode == kPdgKP || pcode == kPdgKM) tcode = kPdgProton;
604  scode = pcode;
605  s2code = tcode;
606  }
607 
608  // check remnants
609  if ( fRemnA < 1 ) //we've blown nucleus apart, no need to retry anything - exit
610  {
611  LOG("HAIntranuke",pNOTICE) << "InelasticHA() stops : not enough nucleons";
613  ev->AddParticle(*p);
614  return;
615  }
616  else if ( fRemnZ + (((pcode==kPdgProton)||(pcode==kPdgPiP))?1:0) - (pcode==kPdgPiM?1:0)
617  < ((( scode==kPdgProton)||( scode==kPdgPiP)) ?1:0) - (scode ==kPdgPiM ?1:0)
618  + (((s2code==kPdgProton)||(s2code==kPdgPiP)) ?1:0) - (s2code==kPdgPiM ?1:0) )
619  {
620  LOG("HAIntranuke",pWARN) << "InelasticHA() failed : too few protons in nucleus";
622  ev->AddParticle(*p);
623  return; // another extreme case, best strategy is to exit and go to next event
624  }
625 
626  GHepParticle t(*p);
627  t.SetPdgCode(tcode);
628 
629  // set up fermi target
630  Target target(ev->TargetNucleus()->Pdg());
631  double tM = t.Mass();
632 
633  // handle fermi momentum
634  if(fDoFermi)
635  {
636  target.SetHitNucPdg(tcode);
638  TVector3 tP3 = fFermiFac * fNuclmodel->Momentum3();
639  double tE = TMath::Sqrt(tP3.Mag2()+ tM*tM);
640  t.SetMomentum(TLorentzVector(tP3,tE));
641  }
642  else
643  {
644  t.SetMomentum(0,0,0,tM);
645  }
646 
647  GHepParticle * cl = new GHepParticle(*p); // clone particle, to run IntBounce at proper energy
648  // calculate energy and momentum using invariant mass
649  double pM = p->Mass();
650  double E_p = ((*p->P4() + *t.P4()).Mag2() - tM*tM - pM*pM)/(2.0*tM);
651  double P_p = TMath::Sqrt(E_p*E_p - pM*pM);
652  cl->SetMomentum(TLorentzVector(P_p,0,0,E_p));
653  // momentum doesn't have to be in right direction, only magnitude
654  double C3CM = fHadroData->IntBounce(cl,tcode,scode,h_fate);
655  delete cl;
656  if (C3CM<-1.) // hope this doesn't occur too often - unphysical but we just pass it on
657  {
658  LOG("HAIntranuke", pWARN) << "unphysical angle chosen in InelasicHA - put particle outside nucleus";
660  ev->AddParticle(*p);
661  return;
662  }
663  double KE1L = p->KinE();
664  double KE2L = t.KinE();
665  LOG("HAIntranuke",pINFO)
666  << " KE1L = " << KE1L << " " << KE1L << " KE2L = " << KE2L;
667  GHepParticle cl1(*p);
668  GHepParticle cl2(t);
669  bool success = utils::intranuke::TwoBodyCollision(ev,pcode,tcode,scode,s2code,C3CM,
670  &cl1,&cl2,fRemnA,fRemnZ,fRemnP4,kIMdHA);
671  if(success)
672  {
673  double P3L = TMath::Sqrt(cl1.Px()*cl1.Px() + cl1.Py()*cl1.Py() + cl1.Pz()*cl1.Pz());
674  double P4L = TMath::Sqrt(cl2.Px()*cl2.Px() + cl2.Py()*cl2.Py() + cl2.Pz()*cl2.Pz());
675  double E3L = cl1.KinE();
676  double E4L = cl2.KinE();
677  LOG ("HAIntranuke",pINFO) << "Successful quasielastic scattering or charge exchange";
678  LOG("HAIntranuke",pINFO)
679  << "C3CM = " << C3CM << "\n P3L, E3L = "
680  << P3L << " " << E3L << " P4L, E4L = "<< P4L << " " << E4L ;
681  if(ev->Probe() ) { LOG("HAIntranuke",pINFO)
682  << "P4L = " << P4L << " ;E4L= " << E4L << "\n probe KE = " << ev->Probe()->KinE() << "\n";
683  if (ev->Probe() && (E3L>ev->Probe()->E()||E4L>ev->Probe()->E())) //is this redundant?
684  {
686  exception.SetReason("TwoBodyCollison gives KE> probe KE in hA simulation");
687  throw exception;
688  }
689  }
690  // always get here, even nucleon decay
691  ev->AddParticle(cl1);
692  ev->AddParticle(cl2);
693  LOG("HAIntranuke", pDEBUG) << "Nucleus : (A,Z) = ("<<fRemnA<<','<<fRemnZ<<')';
694 
695  } else
696  {
698  exception.SetReason("TwoBodyCollison failed in hA simulation");
699  throw exception;
700  }
701 
702 }
An exception thrown by SimulateHadronState for kinematics problems. TwoBodyCollision/Kinematics used ...
TRandom3 & RndFsi(void) const
rnd number generator used by intranuclear cascade monte carlos
Definition: RandomGen.h:60
const XML_Char * target
Definition: expat.h:268
INukeHadroData * fHadroData
a collection of h+N,h+A data & calculations
Definition: Intranuke.h:93
double E(void) const
Get energy.
Definition: GHepParticle.h:92
TLorentzVector fRemnP4
P4 of remnant system.
Definition: Intranuke.h:98
static RandomGen * Instance()
Access instance.
Definition: RandomGen.cxx:79
const TLorentzVector * P4(void) const
Definition: GHepParticle.h:79
::xsd::cxx::tree::exception< char > exception
Definition: Database.h:225
void SetMomentum(const TLorentzVector &p4)
double Mass(void) const
Mass that corresponds to the PDG code.
A singleton holding random number generator classes. All random number generation in GENIE should tak...
Definition: RandomGen.h:30
const TVector3 & Momentum3(void) const
Definition: NuclearModelI.h:67
virtual GHepParticle * Probe(void) const
Definition: GHepRecord.cxx:321
int Pdg(void) const
Definition: GHepParticle.h:64
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:97
enum genie::EINukeFateHN_t INukeFateHN_t
A Neutrino Interaction Target. Is a transparent encapsulation of quite different physical systems suc...
Definition: Target.h:41
const int kPdgPiP
Definition: PDGCodes.h:135
const int kPdgPi0
Definition: PDGCodes.h:137
virtual GHepParticle * TargetNucleus(void) const
Definition: GHepRecord.cxx:330
def success(message)
Definition: log.py:5
#define pINFO
Definition: Messenger.h:63
int fRemnA
remnant nucleus A
Definition: Intranuke.h:96
#define pWARN
Definition: Messenger.h:61
float Mag2() const
double IntBounce(const GHepParticle *p, int target, int s1, INukeFateHN_t fate)
double KinE(bool mass_from_pdg=false) const
Get kinetic energy.
Var Sqrt(const Var &v)
Use to take sqrt of a var.
Definition: Var.cxx:326
void SetStatus(GHepStatus_t s)
Definition: GHepParticle.h:127
static string AsString(INukeFateHN_t fate)
int fRemnZ
remnant nucleus Z
Definition: Intranuke.h:97
const int kPdgPiM
Definition: PDGCodes.h:136
bool fDoFermi
whether or not to do fermi mom.
Definition: Intranuke.h:113
virtual bool GenerateNucleon(const Target &) const =0
virtual void AddParticle(const GHepParticle &p)
Definition: GHepRecord.cxx:535
double fFermiFac
testing parameter to modify fermi momentum
Definition: Intranuke.h:111
const int kPdgProton
Definition: PDGCodes.h:65
#define pNOTICE
Definition: Messenger.h:62
const int kPdgNeutron
Definition: PDGCodes.h:67
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:40
bool TwoBodyCollision(GHepRecord *ev, int pcode, int tcode, int scode, int s2code, double C3CM, GHepParticle *p, GHepParticle *t, int &RemnA, int &RemnZ, TLorentzVector &RemnP4, EINukeMode mode=kIMdHA)
Intranuke utility functions.
Definition: INukeUtils.cxx:625
const NuclearModelI * fNuclmodel
nuclear model used to generate fermi momentum
Definition: Intranuke.h:95
#define pDEBUG
Definition: Messenger.h:64
void Algorithm::Initialize ( void  )
protectedinherited

Definition at line 343 of file Algorithm.cxx.

Referenced by genie::Algorithm::AllowReconfig().

344 {
345 // Algorithm initialization
346 //
347  fAllowReconfig = true;
348  fOwnsSubstruc = false;
349  fConfig = 0;
350  fOwnedSubAlgMp = 0;
351 }
AlgMap * fOwnedSubAlgMp
local pool for owned sub-algs (taken out of the factory pool)
Definition: Algorithm.h:167
bool fOwnsSubstruc
true if it owns its substructure (sub-algs,...)
Definition: Algorithm.h:155
Registry * fConfig
Summary configuration derived from fConvVect, not necessarily allocated.
Definition: Algorithm.h:194
bool Intranuke::IsInNucleus ( const GHepParticle p) const
protectedinherited

Definition at line 248 of file Intranuke.cxx.

References genie::Intranuke::fHadStep, genie::Intranuke::fTrackingRadius, and genie::GHepParticle::X4().

Referenced by genie::Intranuke::TransportHadrons().

249 {
250 // check whether the input particle is still within the nucleus
251 //
252  return (p->X4()->Vect().Mag() < fTrackingRadius + fHadStep);
253 }
double fHadStep
step size for intranuclear hadron transport
Definition: Intranuke.h:107
double fTrackingRadius
tracking radius for the nucleus in the current event
Definition: Intranuke.h:91
const TLorentzVector * X4(void) const
Definition: GHepParticle.h:80
void HAIntranuke::LoadConfig ( void  )
privatevirtual

Implements genie::Intranuke.

Definition at line 1385 of file HAIntranuke.cxx.

References genie::INukeMode::AsString(), genie::Intranuke::fDelRNucleon, genie::Intranuke::fDelRPion, genie::Intranuke::fDoCompoundNucleus, genie::Intranuke::fDoFermi, genie::Intranuke::fEPreEq, genie::Intranuke::fFermiFac, genie::Intranuke::fFermiMomentum, genie::Intranuke::fHadroData, genie::Intranuke::fHadStep, genie::Intranuke::fNR, genie::Intranuke::fNucAbsFac, genie::Intranuke::fNucCEXFac, genie::Intranuke::fNucleonFracAbsScale, genie::Intranuke::fNucleonFracCExScale, genie::Intranuke::fNucleonFracElasScale, genie::Intranuke::fNucleonFracInelScale, genie::Intranuke::fNucleonFracPiProdScale, genie::Intranuke::fNucleonMFPScale, genie::Intranuke::fNuclmodel, genie::Intranuke::fNucRmvE, genie::Intranuke::fPionFracAbsScale, genie::Intranuke::fPionFracCExScale, genie::Intranuke::fPionFracElasScale, genie::Intranuke::fPionFracInelScale, genie::Intranuke::fPionFracPiProdScale, genie::Intranuke::fPionMFPScale, genie::Intranuke::fR0, genie::Algorithm::GetParam(), genie::Algorithm::GetParamDef(), genie::INukeHadroData::Instance(), genie::kIMdHA, LOG, pINFO, and genie::Algorithm::SubAlg().

Referenced by Configure().

1386 {
1387  // load hadronic cross sections
1389 
1390  // fermi momentum setup
1391  fNuclmodel = dynamic_cast<const NuclearModelI *>( this -> SubAlg("NuclearModel") ) ;
1392 
1393  // other intranuke config params
1394  GetParam( "NUCL-R0", fR0 ); // fm
1395  GetParam( "NUCL-NR", fNR );
1396 
1397  GetParam( "INUKE-NucRemovalE", fNucRmvE ); // GeV
1398  GetParam( "INUKE-HadStep", fHadStep ) ;
1399  GetParam( "INUKE-NucAbsFac", fNucAbsFac ) ;
1400  GetParam( "INUKE-NucCEXFac", fNucCEXFac ) ;
1401  GetParam( "INUKE-Energy_Pre_Eq", fEPreEq ) ;
1402  GetParam( "INUKE-FermiFac", fFermiFac ) ;
1403  GetParam( "INUKE-FermiMomentum", fFermiMomentum ) ;
1404 
1405  GetParam( "INUKE-DoCompoundNucleus", fDoCompoundNucleus ) ;
1406  GetParam( "INUKE-DoFermi", fDoFermi ) ;
1407 
1408  GetParam( "HAINUKE-DelRPion", fDelRPion ) ;
1409  GetParam( "HAINUKE-DelRNucleon", fDelRNucleon ) ;
1410 
1411  GetParamDef( "FSI-Pion-MFPScale", fPionMFPScale, 1.0 ) ;
1412  GetParamDef( "FSI-Pion-FracCExScale", fPionFracCExScale, 1.0 ) ;
1413  GetParamDef( "FSI-Pion-FracAbsScale", fPionFracAbsScale, 1.0 ) ;
1414  GetParamDef( "FSI-Pion-FracElasScale", fPionFracElasScale, 1.0 ) ;
1415  GetParamDef( "FSI-Pion-FracInelScale", fPionFracInelScale, 1.0 ) ;
1416  GetParamDef( "FSI-Pion-FracPiProdScale", fPionFracPiProdScale, 1.0 ) ;
1417  GetParamDef( "FSI-Nucleon-MFPScale", fNucleonMFPScale, 1.0 ) ;
1418  GetParamDef( "FSI-Nucleon-FracCExScale", fNucleonFracCExScale , 1.0 ) ;
1419  GetParamDef( "FSI-Nucleon-FracElasScale", fNucleonFracElasScale, 1.0 ) ;
1420  GetParamDef( "FSI-Nucleon-FracInelScale", fNucleonFracInelScale, 1.0 ) ;
1421  GetParamDef( "FSI-Nucleon-FracAbsScale", fNucleonFracAbsScale, 1.0 ) ;
1422  GetParamDef( "FSI-Nucleon-FracPiProdScale", fNucleonFracPiProdScale, 1.0 ) ;
1423 
1424  // report
1425  LOG("HAIntranuke", pINFO) << "Settings for INTRANUKE mode: " << INukeMode::AsString(kIMdHA);
1426  LOG("HAIntranuke", pINFO) << "R0 = " << fR0 << " fermi";
1427  LOG("HAIntranuke", pINFO) << "NR = " << fNR;
1428  LOG("HAIntranuke", pINFO) << "DelRPion = " << fDelRPion;
1429  LOG("HAIntranuke", pINFO) << "DelRNucleon = " << fDelRNucleon;
1430  LOG("HAIntranuke", pINFO) << "HadStep = " << fHadStep << " fermi";
1431  LOG("HAIntranuke", pINFO) << "EPreEq = " << fHadStep << " fermi";
1432  LOG("HAIntranuke", pINFO) << "NucAbsFac = " << fNucAbsFac;
1433  LOG("HAIntranuke", pINFO) << "NucCEXFac = " << fNucCEXFac;
1434  LOG("HAIntranuke", pINFO) << "FermiFac = " << fFermiFac;
1435  LOG("HAIntranuke", pINFO) << "FermiMomtm = " << fFermiMomentum;
1436  LOG("HAIntranuke", pINFO) << "DoFermi? = " << ((fDoFermi)?(true):(false));
1437  LOG("HAIntranuke", pINFO) << "DoCmpndNuc? = " << ((fDoCompoundNucleus)?(true):(false));
1438 }
double fFermiMomentum
whether or not particle collision is pauli blocked
Definition: Intranuke.h:112
static string AsString(INukeMode_t mode)
Definition: INukeMode.h:41
double fNucleonFracInelScale
Definition: Intranuke.h:126
INukeHadroData * fHadroData
a collection of h+N,h+A data & calculations
Definition: Intranuke.h:93
double fHadStep
step size for intranuclear hadron transport
Definition: Intranuke.h:107
double fNucleonFracCExScale
Definition: Intranuke.h:124
Pure abstract base class. Defines the NuclearModelI interface to be implemented by any physics model ...
Definition: NuclearModelI.h:42
double fNucleonFracAbsScale
Definition: Intranuke.h:127
double fDelRNucleon
factor by which Nucleon Compton wavelength gets multiplied to become nuclear size enhancement ...
Definition: Intranuke.h:106
double fNucleonMFPScale
Definition: Intranuke.h:123
double fDelRPion
factor by which Pion Compton wavelength gets multiplied to become nuclear size enhancement ...
Definition: Intranuke.h:105
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
double fPionFracPiProdScale
Definition: Intranuke.h:122
double fPionFracInelScale
Definition: Intranuke.h:120
static INukeHadroData * Instance(void)
bool fDoCompoundNucleus
whether or not to do compound nucleus considerations
Definition: Intranuke.h:115
double fR0
effective nuclear size param
Definition: Intranuke.h:102
#define pINFO
Definition: Messenger.h:63
double fNucRmvE
binding energy to subtract from cascade nucleons
Definition: Intranuke.h:104
double fNucAbsFac
absorption xsec correction factor (hN Mode)
Definition: Intranuke.h:108
double fEPreEq
threshold for pre-equilibrium reaction
Definition: Intranuke.h:110
double fPionFracAbsScale
Definition: Intranuke.h:121
double fPionFracCExScale
Definition: Intranuke.h:118
double fPionFracElasScale
Definition: Intranuke.h:119
double fNucCEXFac
charge exchange xsec correction factor (hN Mode)
Definition: Intranuke.h:109
bool fDoFermi
whether or not to do fermi mom.
Definition: Intranuke.h:113
double fFermiFac
testing parameter to modify fermi momentum
Definition: Intranuke.h:111
bool GetParamDef(const RgKey &name, T &p, const T &def) const
bool GetParam(const RgKey &name, T &p, bool is_top_call=true) const
double fNucleonFracPiProdScale
Definition: Intranuke.h:128
double fPionMFPScale
tweaking factors for tuning
Definition: Intranuke.h:117
double fNucleonFracElasScale
Definition: Intranuke.h:125
const NuclearModelI * fNuclmodel
nuclear model used to generate fermi momentum
Definition: Intranuke.h:95
double fNR
param multiplying the nuclear radius, determining how far to track hadrons beyond the "nuclear bounda...
Definition: Intranuke.h:103
const Algorithm * SubAlg(const RgKey &registry_key) const
Definition: Algorithm.cxx:353
int Algorithm::MergeTopRegistry ( const Registry r)
protectedinherited

Merge with top level registry if first reg of the vector is owned Otherwise an owned copy is added as a top registry

Definition at line 618 of file Algorithm.cxx.

618  {
619 
620  if ( fOwnerships.empty() ) {
621 
622  // this algorithm is not configured right now, the incoming registry is the only configuration
623  Registry * p = new Registry( r ) ;
624  AddTopRegistry( p ) ;
625 
626  return 1 ;
627  }
628 
629  if ( fOwnerships[0] ) {
630  //the top registry is owned: it can be changed with no consequences for other algorithms
631  fConfVect[0] -> Merge( r ) ;
632  }
633  else {
634  // The top registry is not owned so it cannot be changed
635  // The registry will be added with top priority
636 
637  Registry * p = new Registry( r ) ;
638  AddTopRegistry( p ) ;
639  }
640 
641  // The configuration has changed so the summary is not updated anymore and must be deleted
642  if ( fConfig ) {
643  delete fConfig ;
644  fConfig = 0 ;
645  }
646 
647  return fConfVect.size() ;
648 }
const char * p
Definition: xmltok.h:285
int AddTopRegistry(Registry *rp, bool owns=true)
add registry with top priority, also update ownership
Definition: Algorithm.cxx:585
vector< Registry * > fConfVect
Definition: Algorithm.h:161
A registry. Provides the container for algorithm configuration parameters.
Definition: Registry.h:66
vector< bool > fOwnerships
ownership for every registry in fConfVect
Definition: Algorithm.h:164
Registry * fConfig
Summary configuration derived from fConvVect, not necessarily allocated.
Definition: Algorithm.h:194
bool Intranuke::NeedsRescattering ( const GHepParticle p) const
protectedinherited

Definition at line 212 of file Intranuke.cxx.

References ana::assert(), genie::Intranuke::fGMode, genie::pdg::IsIon(), genie::kGMdHadronNucleus, genie::kGMdPhotonNucleus, genie::kIStHadronInTheNucleus, genie::kIStInitialState, genie::GHepParticle::Pdg(), and genie::GHepParticle::Status().

Referenced by genie::Intranuke::TransportHadrons().

213 {
214 // checks whether the particle should be rescattered
215 
216  assert(p);
217 
218  if(fGMode == kGMdHadronNucleus ||
220  // hadron/photon-nucleus scattering propagate the incoming particle
221  return (
223  && !pdg::IsIon(p->Pdg()));
224  }
225  else {
226  // attempt to rescatter anything marked as 'hadron in the nucleus'
227  return (p->Status() == kIStHadronInTheNucleus);
228  }
229 }
GHepStatus_t Status(void) const
Definition: GHepParticle.h:65
int Pdg(void) const
Definition: GHepParticle.h:64
GEvGenMode_t fGMode
event generation mode (lepton+A, hadron+A, ...)
Definition: Intranuke.h:99
bool IsIon(int pdgc)
Definition: PDGUtils.cxx:40
assert(nhit_max >=nhit_nbins)
double HAIntranuke::PiBounce ( void  ) const
private

Definition at line 363 of file HAIntranuke.cxx.

References MECModelEnuComparisons::i, genie::RandomGen::Instance(), calib::j, LOG, pNOTICE, r(), generate_hists::rnd, genie::RandomGen::RndFsi(), and chisquared::theta.

Referenced by ElasHA().

364 {
365 // [adapted from neugen3 intranuke_bounce.F]
366 // [is a fortran stub / difficult to understand - needs to be improved]
367 //
368 // Generates theta in radians for elastic pion-nucleus scattering/
369 // Lookup table is based on Fig 17 of Freedman, Miller and Henley, Nucl.Phys.
370 // A389, 457 (1982)
371 //
372  const int nprob = 25;
373  double dintor = 0.0174533;
374  double denom = 47979.453;
375  double rprob[nprob] = {
376  5000., 4200., 3000., 2600., 2100., 1800., 1200., 750., 500., 230., 120.,
377  35., 9., 3., 11., 18., 29., 27., 20., 14., 10., 6., 2., 0.14, 0.19 };
378 
379  double angles[nprob];
380  for(int i=0; i<nprob; i++) angles[i] = 2.5*i;
381 
383  double r = rnd->RndFsi().Rndm();
384 
385  double xsum = 0.;
386  double theta = 0.;
387  double binl = 0.;
388  double binh = 0.;
389  int tj = 0;
390  for(int i=0; i<60; i++) {
391  theta = i+0.5;
392  for(int j=0; j < nprob-1; j++) {
393  binl = angles[j];
394  binh = angles[j+1];
395  tj=j;
396  if(binl<=theta && binh>=theta) break;
397  tj=0;
398  }//j
399  int itj = tj;
400  double tfract = (theta-binl)/2.5;
401  double delp = rprob[itj+1] - rprob[itj];
402  xsum += (rprob[itj] + tfract*delp)/denom;
403  if(xsum>r) break;
404  theta = 0.;
405  }//i
406 
407  theta *= dintor;
408 
409  LOG("HAIntranuke", pNOTICE)
410  << "Generated pi+A elastic scattering angle = " << theta << " radians";
411 
412  return theta;
413 }
TRandom3 & RndFsi(void) const
rnd number generator used by intranuclear cascade monte carlos
Definition: RandomGen.h:60
static RandomGen * Instance()
Access instance.
Definition: RandomGen.cxx:79
A singleton holding random number generator classes. All random number generation in GENIE should tak...
Definition: RandomGen.h:30
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
const double j
Definition: BetheBloch.cxx:29
TRandom3 r(0)
#define pNOTICE
Definition: Messenger.h:62
double HAIntranuke::PnBounce ( void  ) const
private

Definition at line 415 of file HAIntranuke.cxx.

References MECModelEnuComparisons::i, genie::RandomGen::Instance(), calib::j, LOG, pNOTICE, r(), generate_hists::rnd, genie::RandomGen::RndFsi(), and chisquared::theta.

Referenced by ElasHA().

416 {
417 // [adapted from neugen3 intranuke_pnbounce.F]
418 // [is a fortran stub / difficult to understand - needs to be improved]
419 //
420 // Generates theta in radians for elastic nucleon-nucleus scattering.
421 // Use 800 MeV p+O16 as template in same (highly simplified) spirit as pi+A
422 // from table in Adams et al., PRL 1979. Guess value at 0-2 deg based on Ni
423 // data.
424 //
425  const int nprob = 20;
426  double dintor = 0.0174533;
427  double denom = 11967.0;
428  double rprob[nprob] = {
429  2400., 2350., 2200., 2000., 1728., 1261., 713., 312., 106., 35.,
430  6., 5., 10., 12., 11., 9., 6., 1., 1., 1. };
431 
432  double angles[nprob];
433  for(int i=0; i<nprob; i++) angles[i] = 1.0*i;
434 
436  double r = rnd->RndFsi().Rndm();
437 
438  double xsum = 0.;
439  double theta = 0.;
440  double binl = 0.;
441  double binh = 0.;
442  int tj = 0;
443  for(int i=0; i< nprob; i++) {
444  theta = i+0.5;
445  for(int j=0; j < nprob-1; j++) {
446  binl = angles[j];
447  binh = angles[j+1];
448  tj=j;
449  if(binl<=theta && binh>=theta) break;
450  tj=0;
451  }//j
452  int itj = tj;
453  double tfract = (theta-binl)/2.5;
454  double delp = rprob[itj+1] - rprob[itj];
455  xsum += (rprob[itj] + tfract*delp)/denom;
456  if(xsum>r) break;
457  theta = 0.;
458  }//i
459 
460  theta *= dintor;
461 
462  LOG("HAIntranuke", pNOTICE)
463  << "Generated N+A elastic scattering angle = " << theta << " radians";
464 
465  return theta;
466 }
TRandom3 & RndFsi(void) const
rnd number generator used by intranuclear cascade monte carlos
Definition: RandomGen.h:60
static RandomGen * Instance()
Access instance.
Definition: RandomGen.cxx:79
A singleton holding random number generator classes. All random number generation in GENIE should tak...
Definition: RandomGen.h:30
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
const double j
Definition: BetheBloch.cxx:29
TRandom3 r(0)
#define pNOTICE
Definition: Messenger.h:62
void Algorithm::Print ( ostream &  stream) const
virtualinherited

Print algorithm info.

Definition at line 323 of file Algorithm.cxx.

References GetConfig(), and r().

Referenced by genie::Algorithm::AllowReconfig(), and genie::operator<<().

324 {
325  // print algorithm name & parameter-set
326  stream << "\nAlgorithm Key: " << this->fID.Key();
327  stream << " - Owns Substruc: " << ((fOwnsSubstruc) ? "[true]" : "[false]");
328 
329  // print algorithm configuration
330  const Registry & r = this->GetConfig();
331  stream << r;
332 
333  if(fOwnsSubstruc) {
334  AlgMapConstIter iter = fOwnedSubAlgMp->begin();
335  for(; iter!=fOwnedSubAlgMp->end(); ++iter) {
336  Algorithm * alg = iter->second;
337  stream << "<Next algorithm is owned by : " << this->fID.Key() << ">";
338  stream << *alg;
339  }
340  }
341 }
AlgMap * fOwnedSubAlgMp
local pool for owned sub-algs (taken out of the factory pool)
Definition: Algorithm.h:167
bool fOwnsSubstruc
true if it owns its substructure (sub-algs,...)
Definition: Algorithm.h:155
Algorithm abstract base class.
Definition: Algorithm.h:54
AlgId fID
algorithm name and configuration set
Definition: Algorithm.h:156
virtual const Registry & GetConfig(void) const
Definition: Algorithm.cxx:254
map< string, Algorithm * >::const_iterator AlgMapConstIter
Definition: Algorithm.h:51
A registry. Provides the container for algorithm configuration parameters.
Definition: Registry.h:66
TRandom3 r(0)
string Key(void) const
Definition: AlgId.h:47
void HAIntranuke::ProcessEventRecord ( GHepRecord event_rec) const
virtual

Reimplemented from genie::Intranuke.

Definition at line 108 of file HAIntranuke.cxx.

References LOG, pINFO, pNOTICE, and genie::Intranuke::ProcessEventRecord().

109 {
110  LOG("HAIntranuke", pNOTICE)
111  << "************ Running HA MODE INTRANUKE ************";
112 
114 
115  LOG("HAIntranuke", pINFO) << "Done with this event";
116 }
virtual void ProcessEventRecord(GHepRecord *event_rec) const
Definition: Intranuke.cxx:114
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
#define pINFO
Definition: Messenger.h:63
#define pNOTICE
Definition: Messenger.h:62
void Algorithm::SetId ( const AlgId id)
virtualinherited

Set algorithm ID.

Definition at line 313 of file Algorithm.cxx.

Referenced by genie::Algorithm::AllowReconfig().

314 {
315  fID.Copy(id);
316 }
AlgId fID
algorithm name and configuration set
Definition: Algorithm.h:156
void Copy(const AlgId &id)
Definition: AlgId.cxx:78
void Algorithm::SetId ( string  name,
string  config 
)
virtualinherited

Definition at line 318 of file Algorithm.cxx.

319 {
320  fID.SetId(name, config);
321 }
const XML_Char * name
Definition: expat.h:151
Definition: config.py:1
AlgId fID
algorithm name and configuration set
Definition: Algorithm.h:156
void SetId(string name, string config="")
Definition: AlgId.cxx:70
void Intranuke::SetTrackingRadius ( const GHepParticle p) const
protectedinherited

Definition at line 197 of file Intranuke.cxx.

References genie::units::A, genie::GHepParticle::A(), ana::assert(), genie::Intranuke::fNR, genie::Intranuke::fR0, genie::Intranuke::fTrackingRadius, genie::pdg::IsIon(), LOG, genie::GHepParticle::Pdg(), and pNOTICE.

Referenced by genie::Intranuke::ProcessEventRecord().

198 {
199  assert(p && pdg::IsIon(p->Pdg()));
200  double A = p->A();
201  fTrackingRadius = fR0 * TMath::Power(A, 1./3.);
202 
203  // multiply that by some input factor so that hadrons are tracked
204  // beyond the nuclear 'boundary' since the nuclear density distribution
205  // is not zero there
206  fTrackingRadius *= fNR;
207 
208  LOG("Intranuke", pNOTICE)
209  << "Setting tracking radius to R = " << fTrackingRadius;
210 }
double fTrackingRadius
tracking radius for the nucleus in the current event
Definition: Intranuke.h:91
int Pdg(void) const
Definition: GHepParticle.h:64
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
double fR0
effective nuclear size param
Definition: Intranuke.h:102
static const double A
Definition: Units.h:82
bool IsIon(int pdgc)
Definition: PDGUtils.cxx:40
assert(nhit_max >=nhit_nbins)
int A(void) const
#define pNOTICE
Definition: Messenger.h:62
double fNR
param multiplying the nuclear radius, determining how far to track hadrons beyond the "nuclear bounda...
Definition: Intranuke.h:103
void HAIntranuke::SimulateHadronicFinalState ( GHepRecord ev,
GHepParticle p 
) const
privatevirtual

Implements genie::Intranuke.

Definition at line 118 of file HAIntranuke.cxx.

References genie::GHepRecord::AddParticle(), genie::INukeHadroFates::AsString(), genie::GHepParticle::FirstMother(), fNumIterations, HadronFateHA(), genie::kIHAFtUndefined, genie::kIStStableFinalState, genie::kPdgGamma, genie::kPdgKM, genie::kPdgKP, genie::kPdgNeutron, genie::kPdgPi0, genie::kPdgPiM, genie::kPdgPiP, genie::kPdgProton, LOG, genie::GHepParticle::Name(), genie::GHepRecord::Particle(), genie::GHepParticle::Pdg(), pERROR, pNOTICE, genie::GHepParticle::SetStatus(), and SimulateHadronicFinalStateKinematics().

Referenced by SimulateHadronicFinalStateKinematics().

120 {
121 // Simulate a hadron interaction for the input particle p in HA mode
122 //
123  // check inputs
124  if(!p || !ev) {
125  LOG("HAIntranuke", pERROR) << "** Null input!";
126  return;
127  }
128 
129  // get particle code and check whether this particle can be handled
130  int pdgc = p->Pdg();
131  bool is_gamma = (pdgc==kPdgGamma);
132  bool is_pion = (pdgc==kPdgPiP || pdgc==kPdgPiM || pdgc==kPdgPi0);
133  bool is_kaon = (pdgc==kPdgKP || pdgc==kPdgKM);
134  bool is_baryon = (pdgc==kPdgProton || pdgc==kPdgNeutron);
135  bool is_handled = (is_baryon || is_pion || is_kaon || is_gamma);
136  if(!is_handled) {
137  LOG("HAIntranuke", pERROR) << "** Can not handle particle: " << p->Name();
138  return;
139  }
140 
141  // select a fate for the input particle
142  INukeFateHA_t fate = this->HadronFateHA(p);
143 
144  // store the fate
145  ev->Particle(p->FirstMother())->SetRescatterCode((int)fate);
146 
147  if(fate == kIHAFtUndefined) {
148  LOG("HAIntranuke", pERROR) << "** Couldn't select a fate";
150  ev->AddParticle(*p);
151  return;
152  }
153  LOG("HAIntranuke", pNOTICE)
154  << "Selected "<< p->Name() << " fate: "<< INukeHadroFates::AsString(fate);
155 
156  // try to generate kinematics - repeat till is done (should seldom need >2)
157  fNumIterations = 0;
159 }
virtual GHepParticle * Particle(int position) const
Definition: GHepRecord.cxx:148
#define pERROR
Definition: Messenger.h:60
unsigned int fNumIterations
Definition: HAIntranuke.h:79
int Pdg(void) const
Definition: GHepParticle.h:64
int FirstMother(void) const
Definition: GHepParticle.h:67
string Name(void) const
Name that corresponds to the PDG code.
void SimulateHadronicFinalStateKinematics(GHepRecord *ev, GHepParticle *p) const
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
const int kPdgKM
Definition: PDGCodes.h:150
const int kPdgGamma
Definition: PDGCodes.h:166
const int kPdgKP
Definition: PDGCodes.h:149
const int kPdgPiP
Definition: PDGCodes.h:135
const int kPdgPi0
Definition: PDGCodes.h:137
INukeFateHA_t HadronFateHA(const GHepParticle *p) const
void SetStatus(GHepStatus_t s)
Definition: GHepParticle.h:127
static string AsString(INukeFateHN_t fate)
const int kPdgPiM
Definition: PDGCodes.h:136
virtual void AddParticle(const GHepParticle &p)
Definition: GHepRecord.cxx:535
const int kPdgProton
Definition: PDGCodes.h:65
#define pNOTICE
Definition: Messenger.h:62
const int kPdgNeutron
Definition: PDGCodes.h:67
enum genie::EINukeFateHA_t INukeFateHA_t
void HAIntranuke::SimulateHadronicFinalStateKinematics ( GHepRecord ev,
GHepParticle p 
) const
private

Definition at line 161 of file HAIntranuke.cxx.

References genie::INukeHadroFates::AsString(), ElasHA(), genie::GHepParticle::FirstMother(), fNumIterations, Inelastic(), InelasticHA(), genie::kIHAFtAbs, genie::kIHAFtCEx, genie::kIHAFtElas, genie::kIHAFtInelas, genie::kIHAFtPiProd, LOG, genie::GHepParticle::Name(), genie::GHepRecord::Particle(), pINFO, pNOTICE, and SimulateHadronicFinalState().

Referenced by SimulateHadronicFinalState().

163 {
164  // get stored fate
165  INukeFateHA_t fate = (INukeFateHA_t)
166  ev->Particle(p->FirstMother())->RescatterCode();
167 
168  LOG("HAIntranuke", pINFO)
169  << "Generating kinematics for " << p->Name()
170  << " fate: "<< INukeHadroFates::AsString(fate);
171 
172  // try to generate kinematics for the selected fate
173 
174  try
175  {
176  fNumIterations++;
177  if (fate == kIHAFtElas)
178  {
179  this->ElasHA(ev,p,fate);
180  }
181  else
182  if (fate == kIHAFtInelas || fate == kIHAFtCEx)
183  {
184  this->InelasticHA(ev,p,fate);
185  }
186  else if (fate == kIHAFtAbs || fate == kIHAFtPiProd)
187  {
188  this->Inelastic(ev,p,fate);
189  }
190  }
192  {
193  LOG("HAIntranuke", pNOTICE)
194  << exception;
195  if(fNumIterations <= 100) {
196  LOG("HAIntranuke", pNOTICE)
197  << "Failed attempt to generate kinematics for "
198  << p->Name() << " fate: " << INukeHadroFates::AsString(fate)
199  << " - After " << fNumIterations << " tries, still retrying...";
201  } else {
202  LOG("HAIntranuke", pNOTICE)
203  << "Failed attempt to generate kinematics for "
204  << p->Name() << " fate: " << INukeHadroFates::AsString(fate)
205  << " after " << fNumIterations-1
206  << " attempts. Trying a new fate...";
207  this->SimulateHadronicFinalState(ev,p);
208  }
209  }
210 }
An exception thrown by SimulateHadronState for kinematics problems. TwoBodyCollision/Kinematics used ...
virtual GHepParticle * Particle(int position) const
Definition: GHepRecord.cxx:148
::xsd::cxx::tree::exception< char > exception
Definition: Database.h:225
unsigned int fNumIterations
Definition: HAIntranuke.h:79
void InelasticHA(GHepRecord *ev, GHepParticle *p, INukeFateHA_t fate) const
int FirstMother(void) const
Definition: GHepParticle.h:67
string Name(void) const
Name that corresponds to the PDG code.
void SimulateHadronicFinalStateKinematics(GHepRecord *ev, GHepParticle *p) const
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
void Inelastic(GHepRecord *ev, GHepParticle *p, INukeFateHA_t fate) const
#define pINFO
Definition: Messenger.h:63
static string AsString(INukeFateHN_t fate)
void ElasHA(GHepRecord *ev, GHepParticle *p, INukeFateHA_t fate) const
#define pNOTICE
Definition: Messenger.h:62
enum genie::EINukeFateHA_t INukeFateHA_t
void SimulateHadronicFinalState(GHepRecord *ev, GHepParticle *p) const
const Algorithm * Algorithm::SubAlg ( const RgKey registry_key) const
inherited

Access the sub-algorithm pointed to by the input key, either from the local pool or from AlgFactory's pool

Definition at line 353 of file Algorithm.cxx.

References ana::assert(), genie::AlgFactory::GetAlgorithm(), genie::AlgFactory::Instance(), LOG, pERROR, and pINFO.

Referenced by genie::Algorithm::AllowReconfig(), genie::utils::gsl::FullQELdXSec::FullQELdXSec(), genie::NewQELXSec::Integrate(), genie::HadronTransporter::LoadConfig(), genie::NucleonDecayPrimaryVtxGenerator::LoadConfig(), genie::IBDXSecMap::LoadConfig(), genie::DISHadronicSystemGenerator::LoadConfig(), genie::EmpiricalMECPXSec2015::LoadConfig(), genie::COHElasticPXSec::LoadConfig(), genie::SlowRsclCharmDISPXSecLO::LoadConfig(), genie::AhrensNCELPXSec::LoadConfig(), genie::AlamSimoAtharVacasSKPXSec2014::LoadConfig(), genie::UnstableParticleDecayer::LoadConfig(), genie::ReinDFRPXSec::LoadConfig(), genie::PythiaHadronization::LoadConfig(), genie::BYPDF::LoadConfig(), genie::QPMDISPXSec::LoadConfig(), genie::AlvarezRusoCOHPiPXSec::LoadConfig(), genie::RosenbluthPXSec::LoadConfig(), genie::StrumiaVissaniIBDPXSec::LoadConfig(), genie::MECGenerator::LoadConfig(), genie::NNBarOscPrimaryVtxGenerator::LoadConfig(), genie::FermiMover::LoadConfig(), genie::AhrensDMELPXSec::LoadConfig(), genie::IMDAnnihilationPXSec::LoadConfig(), genie::QELEventGenerator::LoadConfig(), genie::KNOPythiaHadronization::LoadConfig(), genie::AivazisCharmPXSecLO::LoadConfig(), genie::RESHadronicSystemGenerator::LoadConfig(), genie::P33PaschosLalakulichPXSec::LoadConfig(), genie::BergerSehgalFMCOHPiPXSec2015::LoadConfig(), genie::QPMDMDISPXSec::LoadConfig(), genie::BergerSehgalCOHPiPXSec2015::LoadConfig(), genie::LwlynSmithQELCCPXSec::LoadConfig(), genie::NuElectronPXSec::LoadConfig(), genie::PaisQELLambdaPXSec::LoadConfig(), genie::ReinSehgalCOHPiPXSec::LoadConfig(), genie::NievesSimoVacasMECPXSec2016::LoadConfig(), genie::CharmHadronization::LoadConfig(), genie::KovalenkoQELCharmPXSec::LoadConfig(), genie::EventGenerator::LoadConfig(), genie::NuclearModelMap::LoadConfig(), genie::SmithMonizQELCCXSec::LoadConfig(), genie::BardinIMDRadCorPXSec::LoadConfig(), genie::LwlynSmithFF::LoadConfig(), genie::ReinSehgalRESPXSec::LoadConfig(), genie::QELEventGeneratorSM::LoadConfig(), genie::MartiniEricsonChanfrayMarteauMECPXSec2016::LoadConfig(), genie::SmithMonizQELCCPXSec::LoadConfig(), genie::QPMDISStrucFuncBase::LoadConfig(), LoadConfig(), genie::NievesQELCCPXSec::LoadConfig(), genie::HAIntranuke2018::LoadConfig(), genie::HNIntranuke2018::LoadConfig(), genie::BSKLNBaseRESPXSec2014::LoadConfig(), genie::ReinSehgalSPPPXSec::LoadConfig(), genie::KNOHadronization::LoadConfig(), and genie::EventGeneratorListAssembler::LoadGenerator().

354 {
355 // Returns the sub-algorithm pointed to this algorithm's XML config file using
356 // the the values of the key.
357 // This method asserts the existence of these keys in the XML config.
358 // Note: Since only 1 parameter is used, the key value should contain both the
359 // algorithm name and its configuration set according to the usual scheme:
360 // namespace::algorithm_name/configuration_set
361 //
362  LOG("Algorithm", pINFO)
363  << "Fetching sub-alg within alg: " << this->Id().Key()
364  << " pointed to by key: " << registry_key;
365 
366  //-- if the algorithm owns its substructure:
367  // return the sub-algorithm from the local pool
368  //
369  if(fOwnsSubstruc) {
370  AlgMapConstIter iter = fOwnedSubAlgMp->find(registry_key);
371  if(iter!=fOwnedSubAlgMp->end()) return iter->second;
372  LOG("Algorithm", pERROR)
373  << "Owned sub-alg pointed to by key: " << registry_key
374  << " was not found within alg: " << this->Id().Key();
375  return 0;
376  }
377 
378  //-- if the algorithm does not own its substructure:
379  // return the sub-algorithm from the AlgFactory's pool
380  RgAlg alg ;
381  GetParam( registry_key, alg ) ;
382 
383  LOG("Algorithm", pINFO)
384  << "Registry key: " << registry_key << " points to algorithm: " << alg;
385 
386  // retrieve the Algorithm object from the the Algorithm factory
387  AlgFactory * algf = AlgFactory::Instance();
388  const Algorithm * algbase = algf->GetAlgorithm(alg.name, alg.config);
389  assert(algbase);
390 
391  return algbase;
392 }
#define pERROR
Definition: Messenger.h:60
AlgMap * fOwnedSubAlgMp
local pool for owned sub-algs (taken out of the factory pool)
Definition: Algorithm.h:167
bool fOwnsSubstruc
true if it owns its substructure (sub-algs,...)
Definition: Algorithm.h:155
Algorithm abstract base class.
Definition: Algorithm.h:54
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
const Algorithm * GetAlgorithm(const AlgId &algid)
Definition: AlgFactory.cxx:86
#define pINFO
Definition: Messenger.h:63
map< string, Algorithm * >::const_iterator AlgMapConstIter
Definition: Algorithm.h:51
virtual const AlgId & Id(void) const
Get algorithm ID.
Definition: Algorithm.h:98
static AlgFactory * Instance()
Definition: AlgFactory.cxx:75
assert(nhit_max >=nhit_nbins)
bool GetParam(const RgKey &name, T &p, bool is_top_call=true) const
The GENIE Algorithm Factory.
Definition: AlgFactory.h:40
string Key(void) const
Definition: AlgId.h:47
void Intranuke::TransportHadrons ( GHepRecord ev) const
protectedinherited

Definition at line 255 of file Intranuke.cxx.

References genie::GHepParticle::A(), genie::GHepRecord::AddParticle(), genie::Intranuke::CanRescatter(), d, genie::Intranuke::fGMode, genie::Intranuke::fHadStep, genie::GHepParticle::FirstMother(), genie::Intranuke::fRemnA, genie::Intranuke::fRemnP4, genie::Intranuke::fRemnZ, genie::Intranuke::fTrackingRadius, genie::Intranuke::GenerateStep(), genie::Intranuke::HandleCompoundNucleus(), genie::Intranuke::IsInNucleus(), genie::kGMdDarkMatterNucleus, genie::kGMdHadronNucleus, genie::kGMdLeptonNucleus, genie::kGMdNeutronOsc, genie::kGMdNucleonDecay, genie::kGMdPhotonNucleus, genie::GHepParticle::KinE(), genie::kIStFinalStateNuclearRemnant, genie::kIStIntermediateState, genie::kIStStableFinalState, genie::kPdgHadronicBlob, LOG, genie::GHepParticle::Name(), genie::Intranuke::NeedsRescattering(), genie::GHepParticle::P4(), genie::GHepRecord::Particle(), pERROR, pNOTICE, genie::GHepRecord::RemnantNucleusPosition(), genie::GHepParticle::SetFirstMother(), genie::GHepParticle::SetStatus(), genie::Intranuke::SimulateHadronicFinalState(), genie::utils::intranuke::StepParticle(), genie::GHepRecord::TargetNucleusPosition(), genie::GHepParticle::X4(), and genie::GHepParticle::Z().

Referenced by genie::Intranuke::ProcessEventRecord().

256 {
257 // transport all hadrons outside the nucleus
258 
259  int inucl = -1;
260  fRemnA = -1;
261  fRemnZ = -1;
262 
263  // Get 'nuclear environment' at the beginning of hadron transport
264  // and keep track of the remnant nucleus A,Z
265 
266  if(fGMode == kGMdHadronNucleus ||
268  {
269  inucl = evrec->TargetNucleusPosition();
270  }
271  else
272  if(fGMode == kGMdLeptonNucleus ||
276  {
277  inucl = evrec->RemnantNucleusPosition();
278  }
279 
280  LOG("Intranuke", pNOTICE)
281  << "Propagating hadrons within nucleus found in position = " << inucl;
282  GHepParticle * nucl = evrec->Particle(inucl);
283  if(!nucl) {
284  LOG("Intranuke", pERROR)
285  << "No nucleus found in position = " << inucl;
286  LOG("Intranuke", pERROR)
287  << *evrec;
288  return;
289  }
290 
291  fRemnA = nucl->A();
292  fRemnZ = nucl->Z();
293 
294  LOG("Intranuke", pNOTICE)
295  << "Nucleus (A,Z) = (" << fRemnA << ", " << fRemnZ << ")";
296 
297  const TLorentzVector & p4nucl = *(nucl->P4());
298  fRemnP4 = p4nucl;
299 
300  // Loop over GHEP and run intranuclear rescattering on handled particles
301  TObjArrayIter piter(evrec);
302  GHepParticle * p = 0;
303  int icurr = -1;
304 
305  while( (p = (GHepParticle *) piter.Next()) )
306  {
307  icurr++;
308 
309  // Check whether the particle needs rescattering, otherwise skip it
310  if( ! this->NeedsRescattering(p) ) continue;
311 
312  if(this->HandleCompoundNucleus(evrec,p,icurr)) continue;
313 
314  LOG("Intranuke", pNOTICE)
315  << " >> Stepping a " << p->Name()
316  << " with kinetic E = " << p->KinE() << " GeV";
317 
318  // Rescatter a clone, not the original particle
319  GHepParticle * sp = new GHepParticle(*p);
320 
321  // Set clone's mom to be the hadron that was cloned
322  sp->SetFirstMother(icurr);
323 
324  // Check whether the particle can be rescattered
325  if(!this->CanRescatter(sp)) {
326 
327  // if I can't rescatter it, I will just take it out of the nucleus
328  LOG("Intranuke", pNOTICE)
329  << "... Current version can't rescatter a " << sp->Name();
330  sp->SetFirstMother(icurr);
332  evrec->AddParticle(*sp);
333  delete sp;
334  continue; // <-- skip to next GHEP entry
335  }
336 
337  // Start stepping particle out of the nucleus
338  bool has_interacted = false;
339  while ( this-> IsInNucleus(sp) )
340  {
341  // advance the hadron by a step
343 
344  // check whether it interacts
345  double d = this->GenerateStep(evrec,sp);
346  has_interacted = (d<fHadStep);
347  if(has_interacted) break;
348  }//stepping
349 
350  if(has_interacted && fRemnA>0) {
351  // the particle interacts - simulate the hadronic interaction
352  LOG("Intranuke", pNOTICE)
353  << "Particle has interacted at location: "
354  << sp->X4()->Vect().Mag() << " / nucl rad= " << fTrackingRadius;
355  this->SimulateHadronicFinalState(evrec,sp);
356  } else if(has_interacted && fRemnA<=0) {
357  // nothing left to interact with!
358  LOG("Intranuke", pNOTICE)
359  << "*** Nothing left to interact with, escaping.";
361  evrec->AddParticle(*sp);
362  evrec->Particle(sp->FirstMother())->SetRescatterCode(1);
363  } else {
364  // the exits the nucleus without interacting - Done with it!
365  LOG("Intranuke", pNOTICE)
366  << "*** Hadron escaped the nucleus! Done with it.";
368  evrec->AddParticle(*sp);
369  evrec->Particle(sp->FirstMother())->SetRescatterCode(1);
370  }
371  delete sp;
372 
373  // Current snapshot
374  //LOG("Intranuke", pINFO) << "Current event record snapshot: " << *evrec;
375 
376  }// GHEP entries
377 
378  // Add remnant nucleus - that 'hadronic blob' has all the remaining hadronic
379  // 4p not put explicitly into the simulated particles
380  TLorentzVector v4(0.,0.,0.,0.);
381  GHepParticle remnant_nucleus(
383  evrec->AddParticle(remnant_nucleus);
384  // Mark the initial remnant nucleus as an intermediate state
385  // Don't do that in the hadron/photon-nucleus scatterig mode since the initial
386  // remnant nucleus and the target nucleus coincide.
387  if(fGMode != kGMdHadronNucleus &&
389  evrec->Particle(inucl)->SetStatus(kIStIntermediateState);
390  }
391 }
int Z(void) const
void SetFirstMother(int m)
Definition: GHepParticle.h:133
#define pERROR
Definition: Messenger.h:60
TLorentzVector fRemnP4
P4 of remnant system.
Definition: Intranuke.h:98
const TLorentzVector * P4(void) const
Definition: GHepParticle.h:79
const int kPdgHadronicBlob
Definition: PDGCodes.h:188
const char * p
Definition: xmltok.h:285
double fHadStep
step size for intranuclear hadron transport
Definition: Intranuke.h:107
double fTrackingRadius
tracking radius for the nucleus in the current event
Definition: Intranuke.h:91
int FirstMother(void) const
Definition: GHepParticle.h:67
string Name(void) const
Name that corresponds to the PDG code.
GEvGenMode_t fGMode
event generation mode (lepton+A, hadron+A, ...)
Definition: Intranuke.h:99
virtual void SimulateHadronicFinalState(GHepRecord *ev, GHepParticle *p) const =0
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
bool NeedsRescattering(const GHepParticle *p) const
Definition: Intranuke.cxx:212
Float_t d
Definition: plot.C:236
int fRemnA
remnant nucleus A
Definition: Intranuke.h:96
bool CanRescatter(const GHepParticle *p) const
Definition: Intranuke.cxx:231
virtual bool HandleCompoundNucleus(GHepRecord *ev, GHepParticle *p, int mom) const =0
double KinE(bool mass_from_pdg=false) const
Get kinetic energy.
void StepParticle(GHepParticle *p, double step, double nuclear_radius=-1.)
Step particle.
Definition: INukeUtils.cxx:341
const TLorentzVector * X4(void) const
Definition: GHepParticle.h:80
void SetStatus(GHepStatus_t s)
Definition: GHepParticle.h:127
int fRemnZ
remnant nucleus Z
Definition: Intranuke.h:97
double GenerateStep(GHepRecord *ev, GHepParticle *p) const
Definition: Intranuke.cxx:393
int A(void) const
#define pNOTICE
Definition: Messenger.h:62
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:40
bool IsInNucleus(const GHepParticle *p) const
Definition: Intranuke.cxx:248

Friends And Related Function Documentation

friend class IntranukeTester
friend

Definition at line 53 of file HAIntranuke.h.

Member Data Documentation

AlgFactory* genie::Intranuke::fAlgf
protectedinherited

algorithm factory instance

Definition at line 94 of file Intranuke.h.

bool genie::Algorithm::fAllowReconfig
protectedinherited
vector<Registry*> genie::Algorithm::fConfVect
protectedinherited

ideally these members should go private Registry will be access only through the GetParam method configurations registries from various sources the order of the vector is the precedence in case of repeated parameters position 0 -> Highest precedence

Definition at line 161 of file Algorithm.h.

double genie::Intranuke::fDelRNucleon
protectedinherited

factor by which Nucleon Compton wavelength gets multiplied to become nuclear size enhancement

Definition at line 106 of file Intranuke.h.

Referenced by genie::Intranuke::GenerateStep(), and LoadConfig().

double genie::Intranuke::fDelRPion
protectedinherited

factor by which Pion Compton wavelength gets multiplied to become nuclear size enhancement

Definition at line 105 of file Intranuke.h.

Referenced by genie::Intranuke::GenerateStep(), and LoadConfig().

bool genie::Intranuke::fDoCompoundNucleus
protectedinherited

whether or not to do compound nucleus considerations

Definition at line 115 of file Intranuke.h.

Referenced by LoadConfig().

bool genie::Intranuke::fDoFermi
protectedinherited

whether or not to do fermi mom.

Definition at line 113 of file Intranuke.h.

Referenced by Inelastic(), InelasticHA(), and LoadConfig().

bool genie::Intranuke::fDoMassDiff
protectedinherited

whether or not to do mass diff. mode

Definition at line 114 of file Intranuke.h.

double genie::Intranuke::fEPreEq
protectedinherited

threshold for pre-equilibrium reaction

Definition at line 110 of file Intranuke.h.

Referenced by LoadConfig().

double genie::Intranuke::fFermiFac
protectedinherited

testing parameter to modify fermi momentum

Definition at line 111 of file Intranuke.h.

Referenced by Inelastic(), InelasticHA(), and LoadConfig().

double genie::Intranuke::fFermiMomentum
protectedinherited

whether or not particle collision is pauli blocked

Definition at line 112 of file Intranuke.h.

Referenced by Inelastic(), and LoadConfig().

TGenPhaseSpace genie::Intranuke::fGenPhaseSpace
mutableprotectedinherited

a phase space generator

Definition at line 92 of file Intranuke.h.

GEvGenMode_t genie::Intranuke::fGMode
mutableprotectedinherited

event generation mode (lepton+A, hadron+A, ...)

Definition at line 99 of file Intranuke.h.

Referenced by genie::Intranuke::NeedsRescattering(), genie::Intranuke::ProcessEventRecord(), and genie::Intranuke::TransportHadrons().

INukeHadroData* genie::Intranuke::fHadroData
protectedinherited

a collection of h+N,h+A data & calculations

Definition at line 93 of file Intranuke.h.

Referenced by HadronFateHA(), Inelastic(), InelasticHA(), and LoadConfig().

double genie::Intranuke::fHadStep
protectedinherited

step size for intranuclear hadron transport

Definition at line 107 of file Intranuke.h.

Referenced by genie::Intranuke::IsInNucleus(), LoadConfig(), and genie::Intranuke::TransportHadrons().

AlgId genie::Algorithm::fID
protectedinherited

algorithm name and configuration set

Definition at line 156 of file Algorithm.h.

Referenced by genie::Algorithm::Id().

double genie::Intranuke::fNR
protectedinherited

param multiplying the nuclear radius, determining how far to track hadrons beyond the "nuclear boundary"

Definition at line 103 of file Intranuke.h.

Referenced by LoadConfig(), and genie::Intranuke::SetTrackingRadius().

double genie::Intranuke::fNucAbsFac
protectedinherited

absorption xsec correction factor (hN Mode)

Definition at line 108 of file Intranuke.h.

Referenced by LoadConfig().

double genie::Intranuke::fNucCEXFac
protectedinherited

charge exchange xsec correction factor (hN Mode)

Definition at line 109 of file Intranuke.h.

Referenced by LoadConfig().

double genie::Intranuke::fNucleonFracAbsScale
protectedinherited

Definition at line 127 of file Intranuke.h.

Referenced by HadronFateHA(), and LoadConfig().

double genie::Intranuke::fNucleonFracCExScale
protectedinherited

Definition at line 124 of file Intranuke.h.

Referenced by HadronFateHA(), and LoadConfig().

double genie::Intranuke::fNucleonFracElasScale
protectedinherited

Definition at line 125 of file Intranuke.h.

Referenced by HadronFateHA(), and LoadConfig().

double genie::Intranuke::fNucleonFracInelScale
protectedinherited

Definition at line 126 of file Intranuke.h.

Referenced by HadronFateHA(), and LoadConfig().

double genie::Intranuke::fNucleonFracPiProdScale
protectedinherited

Definition at line 128 of file Intranuke.h.

Referenced by HadronFateHA(), and LoadConfig().

double genie::Intranuke::fNucleonMFPScale
protectedinherited

Definition at line 123 of file Intranuke.h.

Referenced by genie::Intranuke::GenerateStep(), and LoadConfig().

const NuclearModelI* genie::Intranuke::fNuclmodel
protectedinherited

nuclear model used to generate fermi momentum

Definition at line 95 of file Intranuke.h.

Referenced by Inelastic(), InelasticHA(), and LoadConfig().

double genie::Intranuke::fNucRmvE
protectedinherited

binding energy to subtract from cascade nucleons

Definition at line 104 of file Intranuke.h.

Referenced by Inelastic(), and LoadConfig().

unsigned int genie::HAIntranuke::fNumIterations
mutableprivate
AlgMap* genie::Algorithm::fOwnedSubAlgMp
protectedinherited

local pool for owned sub-algs (taken out of the factory pool)

Definition at line 167 of file Algorithm.h.

vector<bool> genie::Algorithm::fOwnerships
protectedinherited

ownership for every registry in fConfVect

Definition at line 164 of file Algorithm.h.

bool genie::Algorithm::fOwnsSubstruc
protectedinherited

true if it owns its substructure (sub-algs,...)

Definition at line 155 of file Algorithm.h.

double genie::Intranuke::fPionFracAbsScale
protectedinherited

Definition at line 121 of file Intranuke.h.

Referenced by HadronFateHA(), and LoadConfig().

double genie::Intranuke::fPionFracCExScale
protectedinherited

Definition at line 118 of file Intranuke.h.

Referenced by HadronFateHA(), and LoadConfig().

double genie::Intranuke::fPionFracElasScale
protectedinherited

Definition at line 119 of file Intranuke.h.

Referenced by HadronFateHA(), and LoadConfig().

double genie::Intranuke::fPionFracInelScale
protectedinherited

Definition at line 120 of file Intranuke.h.

Referenced by HadronFateHA(), and LoadConfig().

double genie::Intranuke::fPionFracPiProdScale
protectedinherited

Definition at line 122 of file Intranuke.h.

Referenced by HadronFateHA(), and LoadConfig().

double genie::Intranuke::fPionMFPScale
protectedinherited

tweaking factors for tuning

Definition at line 117 of file Intranuke.h.

Referenced by genie::Intranuke::GenerateStep(), and LoadConfig().

double genie::Intranuke::fR0
protectedinherited

effective nuclear size param

Definition at line 102 of file Intranuke.h.

Referenced by LoadConfig(), and genie::Intranuke::SetTrackingRadius().

int genie::Intranuke::fRemnA
mutableprotectedinherited

remnant nucleus A

Definition at line 96 of file Intranuke.h.

Referenced by ElasHA(), genie::Intranuke::GenerateStep(), Inelastic(), InelasticHA(), and genie::Intranuke::TransportHadrons().

TLorentzVector genie::Intranuke::fRemnP4
mutableprotectedinherited

P4 of remnant system.

Definition at line 98 of file Intranuke.h.

Referenced by ElasHA(), Inelastic(), InelasticHA(), and genie::Intranuke::TransportHadrons().

int genie::Intranuke::fRemnZ
mutableprotectedinherited

remnant nucleus Z

Definition at line 97 of file Intranuke.h.

Referenced by ElasHA(), genie::Intranuke::GenerateStep(), Inelastic(), InelasticHA(), and genie::Intranuke::TransportHadrons().

AlgStatus_t genie::Algorithm::fStatus
protectedinherited

algorithm execution status

Definition at line 166 of file Algorithm.h.

Referenced by genie::Algorithm::GetStatus().

double genie::Intranuke::fTrackingRadius
mutableprotectedinherited

tracking radius for the nucleus in the current event

Definition at line 91 of file Intranuke.h.

Referenced by genie::Intranuke::GenerateVertex(), genie::Intranuke::IsInNucleus(), genie::Intranuke::SetTrackingRadius(), and genie::Intranuke::TransportHadrons().


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