Geant4  10.03.p01
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
G4AblaInterface Class Reference

#include <G4AblaInterface.hh>

Inheritance diagram for G4AblaInterface:
Collaboration diagram for G4AblaInterface:

Public Member Functions

 G4AblaInterface ()
 
virtual ~G4AblaInterface ()
 
virtual G4ReactionProductVectorDeExcite (G4Fragment &aFragment)
 
virtual G4HadFinalStateApplyYourself (G4HadProjectile const &, G4Nucleus &)
 
virtual void ModelDescription (std::ostream &outFile) const
 
virtual void DeExciteModelDescription (std::ostream &outFile) const
 
- Public Member Functions inherited from G4VPreCompoundModel
 G4VPreCompoundModel (G4ExcitationHandler *ptr=nullptr, const G4String &modelName="PrecompoundModel")
 
virtual ~G4VPreCompoundModel ()
 
void SetExcitationHandler (G4ExcitationHandler *ptr)
 
G4ExcitationHandlerGetExcitationHandler () const
 
- Public Member Functions inherited from G4HadronicInteraction
 G4HadronicInteraction (const G4String &modelName="HadronicModel")
 
virtual ~G4HadronicInteraction ()
 
virtual G4double SampleInvariantT (const G4ParticleDefinition *p, G4double plab, G4int Z, G4int A)
 
virtual G4bool IsApplicable (const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)
 
G4double GetMinEnergy () const
 
G4double GetMinEnergy (const G4Material *aMaterial, const G4Element *anElement) const
 
void SetMinEnergy (G4double anEnergy)
 
void SetMinEnergy (G4double anEnergy, const G4Element *anElement)
 
void SetMinEnergy (G4double anEnergy, const G4Material *aMaterial)
 
G4double GetMaxEnergy () const
 
G4double GetMaxEnergy (const G4Material *aMaterial, const G4Element *anElement) const
 
void SetMaxEnergy (const G4double anEnergy)
 
void SetMaxEnergy (G4double anEnergy, const G4Element *anElement)
 
void SetMaxEnergy (G4double anEnergy, const G4Material *aMaterial)
 
G4int GetVerboseLevel () const
 
void SetVerboseLevel (G4int value)
 
const G4StringGetModelName () const
 
void DeActivateFor (const G4Material *aMaterial)
 
void ActivateFor (const G4Material *aMaterial)
 
void DeActivateFor (const G4Element *anElement)
 
void ActivateFor (const G4Element *anElement)
 
G4bool IsBlocked (const G4Material *aMaterial) const
 
G4bool IsBlocked (const G4Element *anElement) const
 
void SetRecoilEnergyThreshold (G4double val)
 
G4double GetRecoilEnergyThreshold () const
 
virtual const std::pair
< G4double, G4double
GetFatalEnergyCheckLevels () const
 
virtual std::pair< G4double,
G4double
GetEnergyMomentumCheckLevels () const
 
void SetEnergyMomentumCheckLevels (G4double relativeLevel, G4double absoluteLevel)
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &)
 
virtual void InitialiseModel ()
 

Additional Inherited Members

- Protected Member Functions inherited from G4HadronicInteraction
void SetModelName (const G4String &nam)
 
G4bool IsBlocked () const
 
void Block ()
 
- Protected Attributes inherited from G4HadronicInteraction
G4HadFinalState theParticleChange
 
G4int verboseLevel
 
G4double theMinEnergy
 
G4double theMaxEnergy
 
G4bool isBlocked
 

Detailed Description

Definition at line 49 of file G4AblaInterface.hh.

Constructor & Destructor Documentation

G4AblaInterface::G4AblaInterface ( )

Definition at line 49 of file G4AblaInterface.cc.

49  :
50  G4VPreCompoundModel(NULL, "ABLA"),
51  ablaResult(new G4VarNtp),
52  volant(new G4Volant),
53  theABLAModel(new G4Abla(volant, ablaResult)),
54  eventNumber(0)
55 {
56  theABLAModel->initEvapora();
57 }
Definition: G4Abla.hh:54
G4VPreCompoundModel(G4ExcitationHandler *ptr=nullptr, const G4String &modelName="PrecompoundModel")
void initEvapora()
Definition: G4Abla.cc:682

Here is the call graph for this function:

G4AblaInterface::~G4AblaInterface ( )
virtual

Definition at line 59 of file G4AblaInterface.cc.

59  {
60  delete volant;
61  delete ablaResult;
62  delete theABLAModel;
63 }

Member Function Documentation

virtual G4HadFinalState* G4AblaInterface::ApplyYourself ( G4HadProjectile const &  ,
G4Nucleus  
)
inlinevirtual

Implements G4HadronicInteraction.

Definition at line 56 of file G4AblaInterface.hh.

56  {
57  return NULL;
58  }
G4ReactionProductVector * G4AblaInterface::DeExcite ( G4Fragment aFragment)
virtual

Implements G4VPreCompoundModel.

Definition at line 65 of file G4AblaInterface.cc.

65  {
66  volant->clear();
67  ablaResult->clear();
68 
69  const G4int ARem = aFragment.GetA_asInt();
70  const G4int ZRem = aFragment.GetZ_asInt();
71  const G4double nuclearMass = aFragment.GetGroundStateMass() / MeV;
72  const G4double eStarRem = aFragment.GetExcitationEnergy() / MeV;
73  const G4double jRem = aFragment.GetAngularMomentum().mag() / hbar_Planck;
74  const G4LorentzVector &pRem = aFragment.GetMomentum();
75  const G4double eTotRem = pRem.e();
76  const G4double eKinRem = (eTotRem - pRem.invariantMass()) / MeV;
77  const G4double pxRem = pRem.x() / MeV;
78  const G4double pyRem = pRem.y() / MeV;
79  const G4double pzRem = pRem.z() / MeV;
80 
81  eventNumber++;
82 
83  theABLAModel->breakItUp(ARem, ZRem,
84  nuclearMass,
85  eStarRem,
86  jRem,
87  eKinRem,
88  pxRem,
89  pyRem,
90  pzRem,
91  eventNumber);
92 
94 
95  for(int j = 0; j < ablaResult->ntrack; ++j) { // Copy ABLA result to the EventInfo
96  G4ReactionProduct *product = toG4Particle(ablaResult->avv[j],
97  ablaResult->zvv[j],
98  ablaResult->enerj[j],
99  ablaResult->plab[j]*std::sin(ablaResult->tetlab[j]*pi/180.0)*std::cos(ablaResult->philab[j]*pi/180.0),
100  ablaResult->plab[j]*std::sin(ablaResult->tetlab[j]*pi/180.0)*std::sin(ablaResult->philab[j]*pi/180.0),
101  ablaResult->plab[j]*std::cos(ablaResult->tetlab[j]*pi/180.0));
102  if(product)
103  result->push_back(product);
104  }
105  return result;
106 }
G4double G4ParticleHPJENDLHEData::G4double result
void clear()
G4double plab[VARNTPSIZE]
G4int avv[VARNTPSIZE]
G4ThreeVector GetAngularMomentum() const
Definition: G4Fragment.cc:272
int G4int
Definition: G4Types.hh:78
G4double enerj[VARNTPSIZE]
std::vector< G4ReactionProduct * > G4ReactionProductVector
void breakItUp(G4int nucleusA, G4int nucleusZ, G4double nucleusMass, G4double excitationEnergy, G4double angularMomentum, G4double recoilEnergy, G4double momX, G4double momY, G4double momZ, G4int eventnumber)
Definition: G4Abla.cc:102
G4int GetA_asInt() const
Definition: G4Fragment.hh:266
const G4LorentzVector & GetMomentum() const
Definition: G4Fragment.hh:307
G4double tetlab[VARNTPSIZE]
G4double GetGroundStateMass() const
Definition: G4Fragment.hh:288
G4int GetZ_asInt() const
Definition: G4Fragment.hh:271
void clear()
G4double philab[VARNTPSIZE]
double invariantMass(const HepLorentzVector &w) const
G4int zvv[VARNTPSIZE]
static constexpr double MeV
Definition: G4SIunits.hh:214
static constexpr double pi
Definition: G4SIunits.hh:75
double G4double
Definition: G4Types.hh:76
static constexpr double hbar_Planck
double mag() const
G4double GetExcitationEnergy() const
Definition: G4Fragment.hh:283

Here is the call graph for this function:

void G4AblaInterface::DeExciteModelDescription ( std::ostream &  outFile) const
virtual

Implements G4VPreCompoundModel.

Definition at line 147 of file G4AblaInterface.cc.

147  {
148  outFile << "ABLA V3 is a statistical model for nuclear de-excitation. It simulates\n"
149  << "evaporation of neutrons, protons and alpha particles, as well as fission\n"
150  << "where applicable. The code included in Geant4 is a C++ translation of the\n"
151  << "original Fortran code. More details about the physics are available in the\n"
152  << "the Geant4 Physics Reference Manual and in the reference articles.\n\n"
153  << "References: A.R. Junghans et al., Nucl. Phys. A629 (1998) 635;\n"
154  << " J. Benlliure et al., Nucl. Phys. A628 (1998) 458.\n\n"; }
void G4AblaInterface::ModelDescription ( std::ostream &  outFile) const
virtual

Reimplemented from G4HadronicInteraction.

Definition at line 143 of file G4AblaInterface.cc.

143  {
144  outFile << "ABLA V3 does not provide an implementation of the ApplyYourself method!\n\n";
145 }

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