Geant4_10
G4CascadeDeexciteBase.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 // $Id$
27 //
28 // Semi-concrete base class for de-excitation modules, analogous to
29 // G4CascadeColliderBase.
30 //
31 // 20130806 M. Kelsey -- Per A. Dotti, move zero vector to file scope to
32 // address thread-collision problem.
33 
34 #include "G4CascadeDeexciteBase.hh"
35 #include "G4CascadeCheckBalance.hh"
36 #include "G4CascadeParameters.hh"
37 #include "G4CollisionOutput.hh"
38 #include "G4Fragment.hh"
39 #include "G4InuclNuclei.hh"
41 #include "G4SystemOfUnits.hh"
42 #include <vector>
43 
44 using namespace G4InuclSpecialFunctions;
45 
46 
47 // Constructor and destructor
48 
50  : G4VCascadeDeexcitation(name), balance(0), A(0), Z(0), EEXS(0) {
52  balance = new G4CascadeCheckBalance(name);
53 }
54 
56  delete balance;
57 }
58 
61  if (balance) balance->setVerboseLevel(verbose);
62 }
63 
64 
65 // Copy pertinent information from G4Fragment for modules
66 
68  A = target.GetA_asInt();
69  Z = target.GetZ_asInt();
70  PEX = target.GetMomentum()/GeV; // Convert from G4 to Bertini units
71  EEXS = target.GetExcitationEnergy();
72 }
73 
74 
75 // Create (fill) new G4Fragment with proper momentum/energy handling
76 
77 namespace {
78  static const G4LorentzVector zero(0.,0.,0.,0.); // File scope avoids churn
79 }
80 
81 const G4Fragment&
83  return makeFragment(zero, fragA, fragZ, EX);
84 }
85 
86 const G4Fragment&
88  G4int fragZ, G4double EX) {
89  if (verboseLevel>2) {
90  G4cout << " >>> " << theName << "::makeFragment " << mom << " " << fragA
91  << " " << fragZ << " " << EX << G4endl;
92  }
93 
94  // Adjust four-momentum so that mass is nucleus + excitation
95  G4double mass =
96  G4InuclNuclei::getNucleiMass(fragA,fragZ) + EX/GeV;
97  mom.setVectM(mom.vect(), mass);
98 
99  // Overwrite previous fragment contents, zeroing out excitons
100  aFragment.SetZandA_asInt(fragZ, fragA);
101  aFragment.SetMomentum(mom*GeV); // Bertini uses GeV!
104 
105  return aFragment;
106 }
107 
108 // Decide wether nuclear fragment is candidate for G4BigBanger
109 
111  return explosion(fragment.GetA_asInt(), fragment.GetZ_asInt(),
112  fragment.GetExcitationEnergy()); // in MeV
113 }
114 
116  G4double excitation) const {
117  if (verboseLevel) G4cout << " >>> " << theName << "::explosion ?" << G4endl;
118 
119  const G4int a_cut = 20;
120  const G4double be_cut = 3.0;
121 
122  // Neutron balls, or small fragments with high excitations can explode
123  return ((fragA <= a_cut || fragZ==0) &&
124  (excitation >= be_cut * bindingEnergy(fragA,fragZ))
125  );
126 }
127 
128 
129 // Validate output for energy, momentum conservation, etc.
130 
132  G4CollisionOutput& output) {
133  if (!balance) return true; // Skip checks unless requested
134 
135  if (verboseLevel > 1)
136  G4cout << " >>> " << theName << "::validateOutput" << G4endl;
137 
139  balance->collide(target, output);
140  return balance->okay(); // Returns false if violations
141 }
142 
144  const std::vector<G4InuclElementaryParticle>& particles) {
145  if (!balance) return true; // Skip checks unless requested
146 
147  if (verboseLevel > 1)
148  G4cout << " >>> " << theName << "::validateOutput" << G4endl;
149 
151  balance->collide(target, particles);
152  return balance->okay(); // Returns false if violations
153 }
154 
156  const std::vector<G4InuclNuclei>& fragments) {
157  if (!balance) return true; // Skip checks unless requested
158 
159  if (verboseLevel > 1)
160  G4cout << " >>> " << theName << "::validateOutput" << G4endl;
161 
163  balance->collide(target, fragments);
164  return balance->okay(); // Returns false if violations
165 }
virtual void setVerboseLevel(G4int verbose=0)
void collide(G4InuclParticle *bullet, G4InuclParticle *target, G4CollisionOutput &output)
void SetNumberOfHoles(G4int valueTot, G4int valueP=0)
Definition: G4Fragment.hh:335
const XML_Char * name
Definition: expat.h:151
static G4bool checkConservation()
const XML_Char * target
Definition: expat.h:268
int G4int
Definition: G4Types.hh:78
virtual void setVerboseLevel(G4int verbose=0)
void setVectM(const Hep3Vector &spatial, double mass)
void getTargetData(const G4Fragment &target)
void SetNumberOfExcitedParticle(G4int valueTot, G4int valueP)
Definition: G4Fragment.hh:316
Hep3Vector vect() const
G4GLOB_DLL std::ostream G4cout
G4int GetA_asInt() const
Definition: G4Fragment.hh:218
Float_t Z
Definition: plot.C:39
bool G4bool
Definition: G4Types.hh:79
const G4LorentzVector & GetMomentum() const
Definition: G4Fragment.hh:251
void SetMomentum(const G4LorentzVector &value)
Definition: G4Fragment.hh:256
G4double getNucleiMass() const
virtual G4bool explosion(const G4Fragment &target) const
G4CascadeCheckBalance * balance
const G4Fragment & makeFragment(G4LorentzVector mom, G4int A, G4int Z, G4double EX=0.)
void SetZandA_asInt(G4int Znew, G4int Anew)
Definition: G4Fragment.hh:228
G4int GetZ_asInt() const
Definition: G4Fragment.hh:223
#define G4endl
Definition: G4ios.hh:61
virtual G4bool validateOutput(const G4Fragment &target, G4CollisionOutput &output)
G4CascadeDeexciteBase(const char *name)
double G4double
Definition: G4Types.hh:76
G4double bindingEnergy(G4int A, G4int Z)
G4double GetExcitationEnergy() const
Definition: G4Fragment.hh:235