| Geant4 Cross Reference |
1 // 1 //
2 // ******************************************* 2 // ********************************************************************
3 // * License and Disclaimer 3 // * License and Disclaimer *
4 // * 4 // * *
5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. 9 // * include a list of copyright holders. *
10 // * 10 // * *
11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. *
17 // * 17 // * *
18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. *
24 // ******************************************* 24 // ********************************************************************
25 // 25 //
26 // 26 //
27 //-------------------------------------------- 27 //---------------------------------------------------------------------
28 // 28 //
29 // Geant4 class G4Fragment 29 // Geant4 class G4Fragment
30 // 30 //
31 // Hadronic Process: Nuclear De-excitations 31 // Hadronic Process: Nuclear De-excitations
32 // by V. Lara (May 1998) 32 // by V. Lara (May 1998)
33 // 33 //
34 // Modifications: 34 // Modifications:
35 // 03.05.2010 V.Ivanchenko General cleanup; mo 35 // 03.05.2010 V.Ivanchenko General cleanup; moved obsolete methods from
36 // inline to source 36 // inline to source
37 // 25.09.2010 M. Kelsey -- Change "setprecisio 37 // 25.09.2010 M. Kelsey -- Change "setprecision" to "setwidth" in printout,
38 // add null pointer check. 38 // add null pointer check.
39 // 27.10.2021 A.Ribon extension for hypernucle 39 // 27.10.2021 A.Ribon extension for hypernuclei.
40 40
41 #include "G4Fragment.hh" 41 #include "G4Fragment.hh"
42 #include "G4SystemOfUnits.hh" 42 #include "G4SystemOfUnits.hh"
43 #include "G4Exception.hh" 43 #include "G4Exception.hh"
44 #include "G4ios.hh" 44 #include "G4ios.hh"
45 #include <iomanip> 45 #include <iomanip>
46 46
47 G4Allocator<G4Fragment>*& pFragmentAllocator() 47 G4Allocator<G4Fragment>*& pFragmentAllocator()
48 { 48 {
49 G4ThreadLocalStatic G4Allocator<G4Fragment>* 49 G4ThreadLocalStatic G4Allocator<G4Fragment>* _instance = nullptr;
50 return _instance; 50 return _instance;
51 } 51 }
52 52
53 // Default constructor 53 // Default constructor
54 G4Fragment::G4Fragment() : 54 G4Fragment::G4Fragment() :
55 theA(0), 55 theA(0),
56 theZ(0), 56 theZ(0),
57 theL(0), 57 theL(0),
58 theExcitationEnergy(0.0), 58 theExcitationEnergy(0.0),
59 theGroundStateMass(0.0), 59 theGroundStateMass(0.0),
60 theMomentum(G4LorentzVector(0,0,0,0)), 60 theMomentum(G4LorentzVector(0,0,0,0)),
61 thePolarization(nullptr), 61 thePolarization(nullptr),
62 creatorModel(-1), 62 creatorModel(-1),
63 numberOfParticles(0), 63 numberOfParticles(0),
64 numberOfCharged(0), 64 numberOfCharged(0),
65 numberOfHoles(0), 65 numberOfHoles(0),
66 numberOfChargedHoles(0), 66 numberOfChargedHoles(0),
67 numberOfShellElectrons(0), 67 numberOfShellElectrons(0),
68 xLevel(0), 68 xLevel(0),
69 theParticleDefinition(nullptr), 69 theParticleDefinition(nullptr),
70 spin(0.0), 70 spin(0.0),
71 theCreationTime(0.0) 71 theCreationTime(0.0)
72 {} 72 {}
73 73
74 // Copy Constructor 74 // Copy Constructor
75 G4Fragment::G4Fragment(const G4Fragment &right 75 G4Fragment::G4Fragment(const G4Fragment &right) :
76 theA(right.theA), 76 theA(right.theA),
77 theZ(right.theZ), 77 theZ(right.theZ),
78 theL(right.theL), 78 theL(right.theL),
79 theExcitationEnergy(right.theExcitationEner 79 theExcitationEnergy(right.theExcitationEnergy),
80 theGroundStateMass(right.theGroundStateMass 80 theGroundStateMass(right.theGroundStateMass),
81 theMomentum(right.theMomentum), 81 theMomentum(right.theMomentum),
82 thePolarization(right.thePolarization), 82 thePolarization(right.thePolarization),
83 creatorModel(right.creatorModel), 83 creatorModel(right.creatorModel),
84 numberOfParticles(right.numberOfParticles), 84 numberOfParticles(right.numberOfParticles),
85 numberOfCharged(right.numberOfCharged), 85 numberOfCharged(right.numberOfCharged),
86 numberOfHoles(right.numberOfHoles), 86 numberOfHoles(right.numberOfHoles),
87 numberOfChargedHoles(right.numberOfChargedH 87 numberOfChargedHoles(right.numberOfChargedHoles),
88 numberOfShellElectrons(right.numberOfShellE 88 numberOfShellElectrons(right.numberOfShellElectrons),
89 xLevel(right.xLevel), 89 xLevel(right.xLevel),
90 theParticleDefinition(right.theParticleDefi 90 theParticleDefinition(right.theParticleDefinition),
91 spin(right.spin), 91 spin(right.spin),
92 theCreationTime(right.theCreationTime), 92 theCreationTime(right.theCreationTime),
93 isLongLived(right.isLongLived) 93 isLongLived(right.isLongLived)
94 {} 94 {}
95 95
96 G4Fragment::~G4Fragment() 96 G4Fragment::~G4Fragment()
97 {} 97 {}
98 98
99 G4Fragment::G4Fragment(G4int A, G4int Z, const 99 G4Fragment::G4Fragment(G4int A, G4int Z, const G4LorentzVector& aMomentum) :
100 theA(A), 100 theA(A),
101 theZ(Z), 101 theZ(Z),
102 theL(0), 102 theL(0),
103 theExcitationEnergy(0.0), 103 theExcitationEnergy(0.0),
104 theGroundStateMass(0.0), 104 theGroundStateMass(0.0),
105 theMomentum(aMomentum), 105 theMomentum(aMomentum),
106 thePolarization(nullptr), 106 thePolarization(nullptr),
107 creatorModel(-1), 107 creatorModel(-1),
108 numberOfParticles(0), 108 numberOfParticles(0),
109 numberOfCharged(0), 109 numberOfCharged(0),
110 numberOfHoles(0), 110 numberOfHoles(0),
111 numberOfChargedHoles(0), 111 numberOfChargedHoles(0),
112 numberOfShellElectrons(0), 112 numberOfShellElectrons(0),
113 xLevel(0), 113 xLevel(0),
114 theParticleDefinition(nullptr), 114 theParticleDefinition(nullptr),
115 spin(0.0), 115 spin(0.0),
116 theCreationTime(0.0) 116 theCreationTime(0.0)
117 { 117 {
118 if(theA > 0) { 118 if(theA > 0) {
119 CalculateMassAndExcitationEnergy(); 119 CalculateMassAndExcitationEnergy();
120 } 120 }
121 } 121 }
122 122
123 G4Fragment::G4Fragment(G4int A, G4int Z, G4int 123 G4Fragment::G4Fragment(G4int A, G4int Z, G4int numberOfLambdas,
124 const G4LorentzVector& 124 const G4LorentzVector& aMomentum) :
125 theA(A), 125 theA(A),
126 theZ(Z), 126 theZ(Z),
127 theL(std::max(numberOfLambdas,0)), 127 theL(std::max(numberOfLambdas,0)),
128 theExcitationEnergy(0.0), 128 theExcitationEnergy(0.0),
129 theGroundStateMass(0.0), 129 theGroundStateMass(0.0),
130 theMomentum(aMomentum), 130 theMomentum(aMomentum),
131 thePolarization(nullptr), 131 thePolarization(nullptr),
132 creatorModel(-1), 132 creatorModel(-1),
133 numberOfParticles(0), 133 numberOfParticles(0),
134 numberOfCharged(0), 134 numberOfCharged(0),
135 numberOfHoles(0), 135 numberOfHoles(0),
136 numberOfChargedHoles(0), 136 numberOfChargedHoles(0),
137 numberOfShellElectrons(0), 137 numberOfShellElectrons(0),
138 xLevel(0), 138 xLevel(0),
139 theParticleDefinition(nullptr), 139 theParticleDefinition(nullptr),
140 spin(0.0), 140 spin(0.0),
141 theCreationTime(0.0) 141 theCreationTime(0.0)
142 { 142 {
143 if(theA > 0) { 143 if(theA > 0) {
144 CalculateMassAndExcitationEnergy(); 144 CalculateMassAndExcitationEnergy();
145 } 145 }
146 } 146 }
147 147
148 // This constructor is for initialize photons 148 // This constructor is for initialize photons or electrons
149 G4Fragment::G4Fragment(const G4LorentzVector& 149 G4Fragment::G4Fragment(const G4LorentzVector& aMomentum,
150 const G4ParticleDefinition * aParti 150 const G4ParticleDefinition * aParticleDefinition) :
151 theA(0), 151 theA(0),
152 theZ(0), 152 theZ(0),
153 theL(0), 153 theL(0),
154 theExcitationEnergy(0.0), 154 theExcitationEnergy(0.0),
155 theMomentum(aMomentum), 155 theMomentum(aMomentum),
156 thePolarization(nullptr), 156 thePolarization(nullptr),
157 creatorModel(-1), 157 creatorModel(-1),
158 numberOfParticles(0), 158 numberOfParticles(0),
159 numberOfCharged(0), 159 numberOfCharged(0),
160 numberOfHoles(0), 160 numberOfHoles(0),
161 numberOfChargedHoles(0), 161 numberOfChargedHoles(0),
162 numberOfShellElectrons(0), 162 numberOfShellElectrons(0),
163 xLevel(0), 163 xLevel(0),
164 theParticleDefinition(aParticleDefinition), 164 theParticleDefinition(aParticleDefinition),
165 spin(0.0), 165 spin(0.0),
166 theCreationTime(0.0) 166 theCreationTime(0.0)
167 { 167 {
168 if(aParticleDefinition->GetPDGEncoding() != 168 if(aParticleDefinition->GetPDGEncoding() != 22 &&
169 aParticleDefinition->GetPDGEncoding() != 169 aParticleDefinition->GetPDGEncoding() != 11) {
170 G4ExceptionDescription ed; 170 G4ExceptionDescription ed;
171 ed << "Particle: " << aParticleDefinition- 171 ed << "Particle: " << aParticleDefinition->GetParticleName() << G4endl;
172 G4Exception( "G4Fragment::G4Fragment: cons 172 G4Exception( "G4Fragment::G4Fragment: constructor for gamma used for another type of particle ! ",
173 "HAD_FRAGMENT_01", FatalExcep 173 "HAD_FRAGMENT_01", FatalException, ed );
174 } 174 }
175 theGroundStateMass = aParticleDefinition->Ge 175 theGroundStateMass = aParticleDefinition->GetPDGMass();
176 } 176 }
177 177
178 void G4Fragment::CalculateMassAndExcitationEne 178 void G4Fragment::CalculateMassAndExcitationEnergy()
179 { 179 {
180 // check input 180 // check input
181 if(theZ > theA || theZ + theL > theA) { 181 if(theZ > theA || theZ + theL > theA) {
182 G4ExceptionDescription ed; 182 G4ExceptionDescription ed;
183 ed << "Fragment: Z=" << theZ << " A=" << 183 ed << "Fragment: Z=" << theZ << " A=" << theA << " nLambdas=" << theL << G4endl;
184 G4Exception( "G4Fragment::CalculateMassAnd 184 G4Exception( "G4Fragment::CalculateMassAndExcitationEnergy: inconsistent number of nucleons ! ",
185 "HAD_FRAGMENT_02", EventMustB 185 "HAD_FRAGMENT_02", EventMustBeAborted, ed );
186 } 186 }
187 // compute mass 187 // compute mass
188 theGroundStateMass = ( theL == 0 ) 188 theGroundStateMass = ( theL == 0 )
189 ? G4NucleiProperties::GetNuclearMass(theA, 189 ? G4NucleiProperties::GetNuclearMass(theA, theZ)
190 : G4HyperNucleiProperties::GetNuclearMass( 190 : G4HyperNucleiProperties::GetNuclearMass(theA, theZ, theL);
191 191
192 // excitation energy 192 // excitation energy
193 const G4double minFragExcitation = 10.*CLHEP 193 const G4double minFragExcitation = 10.*CLHEP::eV;
194 theExcitationEnergy = theMomentum.mag() - th 194 theExcitationEnergy = theMomentum.mag() - theGroundStateMass;
195 if(theExcitationEnergy < minFragExcitation) 195 if(theExcitationEnergy < minFragExcitation) {
196 if(theExcitationEnergy < -minFragExcitatio 196 if(theExcitationEnergy < -minFragExcitation) {
197 ExcitationEnergyWarning(); 197 ExcitationEnergyWarning();
198 } 198 }
199 theExcitationEnergy = 0.0; 199 theExcitationEnergy = 0.0;
200 } 200 }
201 } 201 }
202 202
203 void G4Fragment::SetExcEnergyAndMomentum(G4dou 203 void G4Fragment::SetExcEnergyAndMomentum(G4double eexc,
204 const G4LorentzVector& v) 204 const G4LorentzVector& v)
205 { 205 {
206 theExcitationEnergy = eexc; 206 theExcitationEnergy = eexc;
207 theMomentum.set(0.0, 0.0, 0.0, theGroundStat 207 theMomentum.set(0.0, 0.0, 0.0, theGroundStateMass + eexc);
208 theMomentum.boost(v.boostVector()); 208 theMomentum.boost(v.boostVector());
209 } 209 }
210 210
211 G4double G4Fragment::GetBindingEnergy() const 211 G4double G4Fragment::GetBindingEnergy() const
212 { 212 {
213 const G4double lambdaMass = 1.115683*CLHEP:: 213 const G4double lambdaMass = 1.115683*CLHEP::GeV;
214 return (theA-theZ-theL)*CLHEP::neutron_mass_ 214 return (theA-theZ-theL)*CLHEP::neutron_mass_c2
215 + theZ*CLHEP::proton_mass_c2 + theL*lambda 215 + theZ*CLHEP::proton_mass_c2 + theL*lambdaMass
216 - theGroundStateMass; 216 - theGroundStateMass;
217 } 217 }
218 218
219 G4Fragment & G4Fragment::operator=(const G4Fra 219 G4Fragment & G4Fragment::operator=(const G4Fragment &right)
220 { 220 {
221 if (this != &right) { 221 if (this != &right) {
222 theA = right.theA; 222 theA = right.theA;
223 theZ = right.theZ; 223 theZ = right.theZ;
224 theL = right.theL; 224 theL = right.theL;
225 theExcitationEnergy = right.theExcitationE 225 theExcitationEnergy = right.theExcitationEnergy;
226 theGroundStateMass = right.theGroundStateM 226 theGroundStateMass = right.theGroundStateMass;
227 theMomentum = right.theMomentum; 227 theMomentum = right.theMomentum;
228 thePolarization = right.thePolarization; 228 thePolarization = right.thePolarization;
229 creatorModel = right.creatorModel; 229 creatorModel = right.creatorModel;
230 numberOfParticles = right.numberOfParticle 230 numberOfParticles = right.numberOfParticles;
231 numberOfCharged = right.numberOfCharged; 231 numberOfCharged = right.numberOfCharged;
232 numberOfHoles = right.numberOfHoles; 232 numberOfHoles = right.numberOfHoles;
233 numberOfChargedHoles = right.numberOfCharg 233 numberOfChargedHoles = right.numberOfChargedHoles;
234 numberOfShellElectrons = right.numberOfShe 234 numberOfShellElectrons = right.numberOfShellElectrons;
235 xLevel = right.xLevel; 235 xLevel = right.xLevel;
236 theParticleDefinition = right.theParticleD 236 theParticleDefinition = right.theParticleDefinition;
237 spin = right.spin; 237 spin = right.spin;
238 theCreationTime = right.theCreationTime; 238 theCreationTime = right.theCreationTime;
239 isLongLived = right.isLongLived; 239 isLongLived = right.isLongLived;
240 } 240 }
241 return *this; 241 return *this;
242 } 242 }
243 243
244 G4bool G4Fragment::operator==(const G4Fragment 244 G4bool G4Fragment::operator==(const G4Fragment &right) const
245 { 245 {
246 return (this == (G4Fragment *) &right); 246 return (this == (G4Fragment *) &right);
247 } 247 }
248 248
249 G4bool G4Fragment::operator!=(const G4Fragment 249 G4bool G4Fragment::operator!=(const G4Fragment &right) const
250 { 250 {
251 return (this != (G4Fragment *) &right); 251 return (this != (G4Fragment *) &right);
252 } 252 }
253 253
254 std::ostream& operator << (std::ostream &out, 254 std::ostream& operator << (std::ostream &out, const G4Fragment &theFragment)
255 { 255 {
256 std::ios::fmtflags old_floatfield = out.flag 256 std::ios::fmtflags old_floatfield = out.flags();
257 out.setf(std::ios::floatfield); 257 out.setf(std::ios::floatfield);
258 258
259 out << "Fragment: A = " << std::setw(3) << t 259 out << "Fragment: A = " << std::setw(3) << theFragment.theA
260 << ", Z = " << std::setw(3) << theFragme 260 << ", Z = " << std::setw(3) << theFragment.theZ
261 << ", numberOfLambdas = " << std::setw(3 261 << ", numberOfLambdas = " << std::setw(3) << theFragment.theL ;
262 out.setf(std::ios::scientific,std::ios::floa 262 out.setf(std::ios::scientific,std::ios::floatfield);
263 263
264 // Store user's precision setting and reset 264 // Store user's precision setting and reset to (3) here: back-compatibility
265 std::streamsize floatPrec = out.precision(); 265 std::streamsize floatPrec = out.precision();
266 266
267 out << std::setprecision(3) 267 out << std::setprecision(3)
268 << ", U = " << theFragment.GetExcitation 268 << ", U = " << theFragment.GetExcitationEnergy()/CLHEP::MeV
269 << " MeV "; 269 << " MeV ";
270 if(theFragment.GetCreatorModelID() >= 0) { 270 if(theFragment.GetCreatorModelID() >= 0) {
271 out << " creatorModelID= " << theFragment. 271 out << " creatorModelID= " << theFragment.GetCreatorModelID();
272 } 272 }
273 if(theFragment.GetCreationTime() > 0.0) { 273 if(theFragment.GetCreationTime() > 0.0) {
274 out << " Time= " << theFragment.GetCreati 274 out << " Time= " << theFragment.GetCreationTime()/CLHEP::ns << " ns";
275 } 275 }
276 out << G4endl 276 out << G4endl
277 << " P = (" 277 << " P = ("
278 << theFragment.GetMomentum().x()/CLHEP:: 278 << theFragment.GetMomentum().x()/CLHEP::MeV << ","
279 << theFragment.GetMomentum().y()/CLHEP:: 279 << theFragment.GetMomentum().y()/CLHEP::MeV << ","
280 << theFragment.GetMomentum().z()/CLHEP:: 280 << theFragment.GetMomentum().z()/CLHEP::MeV
281 << ") MeV E = " 281 << ") MeV E = "
282 << theFragment.GetMomentum().t()/CLHEP:: 282 << theFragment.GetMomentum().t()/CLHEP::MeV << " MeV"
283 << G4endl; 283 << G4endl;
284 284
285 out << " #spin= " << theFragment.GetSpin( 285 out << " #spin= " << theFragment.GetSpin()
286 << " #floatLevelNo= " << theFragment. 286 << " #floatLevelNo= " << theFragment.GetFloatingLevelNumber() << " ";
287 287
288 if (theFragment.GetNumberOfExcitons() != 0) 288 if (theFragment.GetNumberOfExcitons() != 0) {
289 out << " " 289 out << " "
290 << "#Particles= " << theFragment.GetNumberOf 290 << "#Particles= " << theFragment.GetNumberOfParticles()
291 << ", #Charged= " << theFragment.GetNumberOf 291 << ", #Charged= " << theFragment.GetNumberOfCharged()
292 << ", #Holes= " << theFragment.GetNumberOf 292 << ", #Holes= " << theFragment.GetNumberOfHoles()
293 << ", #ChargedHoles= " << theFragment.GetNum 293 << ", #ChargedHoles= " << theFragment.GetNumberOfChargedHoles();
294 } 294 }
295 out << G4endl; 295 out << G4endl;
296 if(theFragment.GetNuclearPolarization()) { 296 if(theFragment.GetNuclearPolarization()) {
297 out << *(theFragment.GetNuclearPolarizatio 297 out << *(theFragment.GetNuclearPolarization());
298 } 298 }
299 //out << G4endl; 299 //out << G4endl;
300 out.setf(old_floatfield,std::ios::floatfield 300 out.setf(old_floatfield,std::ios::floatfield);
301 out.precision(floatPrec); 301 out.precision(floatPrec);
302 302
303 return out; 303 return out;
304 } 304 }
305 305
306 void G4Fragment::ExcitationEnergyWarning() 306 void G4Fragment::ExcitationEnergyWarning()
307 { 307 {
308 #ifdef G4VERBOSE 308 #ifdef G4VERBOSE
309 G4cout << "G4Fragment::CalculateExcitationEn 309 G4cout << "G4Fragment::CalculateExcitationEnergy(): WARNING "
310 << " GraundStateMass(MeV)= " << theGroundSt 310 << " GraundStateMass(MeV)= " << theGroundStateMass
311 <<G4endl; 311 <<G4endl;
312 G4cout << *this << G4endl; 312 G4cout << *this << G4endl;
313 #endif 313 #endif
314 } 314 }
315 315
316 void G4Fragment::NumberOfExitationWarning(cons 316 void G4Fragment::NumberOfExitationWarning(const G4String& value)
317 { 317 {
318 G4ExceptionDescription ed; 318 G4ExceptionDescription ed;
319 ed << "Value=" << value << G4endl; 319 ed << "Value=" << value << G4endl;
320 G4Exception( "G4Fragment::NumberOfExitationW 320 G4Exception( "G4Fragment::NumberOfExitationWarning : wrong exciton number ! ",
321 "HAD_FRAGMENT_03", FatalException, ed 321 "HAD_FRAGMENT_03", FatalException, ed );
322 } 322 }
323 323
324 void G4Fragment::SetAngularMomentum(const G4Th 324 void G4Fragment::SetAngularMomentum(const G4ThreeVector& v)
325 { 325 {
326 spin = v.mag(); 326 spin = v.mag();
327 } 327 }
328 328
329 G4ThreeVector G4Fragment::GetAngularMomentum() 329 G4ThreeVector G4Fragment::GetAngularMomentum() const
330 { 330 {
331 G4ThreeVector v(0.0,0.0,spin); 331 G4ThreeVector v(0.0,0.0,spin);
332 return v; 332 return v;
333 } 333 }
334 334