Geant4  10.03.p02
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4AblaInterface.cc
Go to the documentation of this file.
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 // ABLAXX statistical de-excitation model
27 // Pekka Kaitaniemi, HIP (translation)
28 // Christelle Schmidt, IPNL (fission code)
29 // Davide Mancusi, CEA (contact person INCL/ABLA)
30 // Aatos Heikkinen, HIP (project coordination)
31 //
32 #define ABLAXX_IN_GEANT4_MODE 1
33 
34 #include "globals.hh"
35 
36 #ifdef ABLAXX_IN_GEANT4_MODE
37 
38 #include "G4AblaInterface.hh"
39 #include "G4ParticleDefinition.hh"
41 #include "G4ReactionProduct.hh"
42 #include "G4DynamicParticle.hh"
43 #include "G4IonTable.hh"
44 #include "G4SystemOfUnits.hh"
45 #include "G4PhysicalConstants.hh"
46 #include <iostream>
47 #include <cmath>
48 
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 }
58 
60  delete volant;
61  delete ablaResult;
62  delete theABLAModel;
63 }
64 
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 }
107 
108 G4ParticleDefinition *G4AblaInterface::toG4ParticleDefinition(G4int A, G4int Z) const {
109  if (A == 1 && Z == 1) return G4Proton::Proton();
110  else if(A == 1 && Z == 0) return G4Neutron::Neutron();
111  else if(A == 0 && Z == 1) return G4PionPlus::PionPlus();
112  else if(A == 0 && Z == -1) return G4PionMinus::PionMinus();
113  else if(A == 0 && Z == 0) return G4PionZero::PionZero();
114  else if(A == 2 && Z == 1) return G4Deuteron::Deuteron();
115  else if(A == 3 && Z == 1) return G4Triton::Triton();
116  else if(A == 3 && Z == 2) return G4He3::He3();
117  else if(A == 4 && Z == 2) return G4Alpha::Alpha();
118  else if(A > 0 && Z > 0 && A >= Z) { // Returns ground state ion definition
119  return G4IonTable::GetIonTable()->GetIon(Z, A, 0);
120  } else { // Error, unrecognized particle
121  G4cout << "Can't convert particle with A=" << A << ", Z=" << Z << " to G4ParticleDefinition, trouble ahead" << G4endl;
122  return 0;
123  }
124 }
125 
126 G4ReactionProduct *G4AblaInterface::toG4Particle(G4int A, G4int Z,
127  G4double kinE,
128  G4double px,
129  G4double py, G4double pz) const {
130  G4ParticleDefinition *def = toG4ParticleDefinition(A, Z);
131  if(def == 0) { // Check if we have a valid particle definition
132  return 0;
133  }
134  const double energy = kinE * MeV;
135  const G4ThreeVector momentum(px, py, pz);
136  const G4ThreeVector momentumDirection = momentum.unit();
137  G4DynamicParticle p(def, momentumDirection, energy);
138  G4ReactionProduct *r = new G4ReactionProduct(def);
139  (*r) = p;
140  return r;
141 }
142 
143 void G4AblaInterface::ModelDescription(std::ostream& outFile) const {
144  outFile << "ABLA V3 does not provide an implementation of the ApplyYourself method!\n\n";
145 }
146 
147 void G4AblaInterface::DeExciteModelDescription(std::ostream& outFile) const {
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"; }
155 
156 #endif // ABLAXX_IN_GEANT4_MODE
G4double G4ParticleHPJENDLHEData::G4double result
Definition: G4Abla.hh:54
void clear()
const char * p
Definition: xmltok.h:285
G4ParticleDefinition * GetIon(G4int Z, G4int A, G4int lvl=0)
Definition: G4IonTable.cc:503
G4double plab[VARNTPSIZE]
virtual G4ReactionProductVector * DeExcite(G4Fragment &aFragment)
G4int avv[VARNTPSIZE]
G4ThreeVector GetAngularMomentum() const
Definition: G4Fragment.cc:271
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
G4GLOB_DLL std::ostream G4cout
double A(double temperature)
G4int GetA_asInt() const
Definition: G4Fragment.hh:266
virtual void ModelDescription(std::ostream &outFile) const
virtual ~G4AblaInterface()
const G4LorentzVector & GetMomentum() const
Definition: G4Fragment.hh:307
static G4Triton * Triton()
Definition: G4Triton.cc:95
G4double tetlab[VARNTPSIZE]
static G4Proton * Proton()
Definition: G4Proton.cc:93
static G4PionPlus * PionPlus()
Definition: G4PionPlus.cc:98
G4double GetGroundStateMass() const
Definition: G4Fragment.hh:288
static G4Neutron * Neutron()
Definition: G4Neutron.cc:104
static G4PionZero * PionZero()
Definition: G4PionZero.cc:108
static G4Deuteron * Deuteron()
Definition: G4Deuteron.cc:94
static G4IonTable * GetIonTable()
Definition: G4IonTable.hh:78
virtual void DeExciteModelDescription(std::ostream &outFile) const
void initEvapora()
Definition: G4Abla.cc:682
G4double energy(const ThreeVector &p, const G4double m)
Hep3Vector unit() const
static G4PionMinus * PionMinus()
Definition: G4PionMinus.cc:98
G4int GetZ_asInt() const
Definition: G4Fragment.hh:271
void clear()
G4double philab[VARNTPSIZE]
#define G4endl
Definition: G4ios.hh:61
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
static G4Alpha * Alpha()
Definition: G4Alpha.cc:89
double G4double
Definition: G4Types.hh:76
static G4He3 * He3()
Definition: G4He3.cc:94
double mag() const
G4double GetExcitationEnergy() const
Definition: G4Fragment.hh:283