ECCE @ EIC Software
Reference for
ECCE @ EIC
simulation and reconstruction software on GitHub
Home page
Related Pages
Modules
Namespaces
Classes
Files
External Links
File List
File Members
All
Classes
Namespaces
Files
Functions
Variables
Typedefs
Enumerations
Enumerator
Friends
Macros
Groups
Pages
G4InuclNuclei.cc
Go to the documentation of this file.
Or view
the newest version in sPHENIX GitHub for file G4InuclNuclei.cc
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
//
27
// 20100301 M. Kelsey -- Add function to create unphysical nuclei for use
28
// as temporary final-state fragments.
29
// 20100319 M. Kelsey -- Add information message to makeNuclearFragment().
30
// Use new GetBindingEnergy() function instead of bindingEnergy().
31
// 20100622 M. Kelsey -- Use local "bindingEnergy()" function to call through.
32
// 20100627 M. Kelsey -- Test for non-physical fragments and abort job.
33
// 20100630 M. Kelsey -- Use excitation energy in G4Ions
34
// 20100714 M. Kelsey -- Use G4DynamicParticle::theDynamicalMass to deal with
35
// excitation energy without instantianting "infinite" G4PartDefns.
36
// 20100719 M. Kelsey -- Change excitation energy without altering momentum
37
// 20100906 M. Kelsey -- Add fill() functions to rewrite contents
38
// 20100910 M. Kelsey -- Add clearExitonConfiguration() to fill() functions
39
// 20100914 M. Kelsey -- Make printout symmetric with G4InuclElemPart,
40
// migrate to integer A and Z
41
// 20100924 M. Kelsey -- Add constructor to copy G4Fragment input, and output
42
// functions to create G4Fragment
43
// 20110214 M. Kelsey -- Replace integer "model" with enum
44
// 20110308 M. Kelsey -- Follow new G4Fragment interface for hole types
45
// 20110427 M. Kelsey -- Remove PDG-code warning
46
// 20110721 M. Kelsey -- Follow base-class ctor change to pass model directly
47
// 20110829 M. Kelsey -- Add constructor to copy G4V3DNucleus input
48
// 20110919 M. Kelsey -- Special case: Allow fill(A=0,Z=0) to make dummy
49
// 20110922 M. Kelsey -- Add stream argument to printParticle() => print()
50
// 20121009 M. Kelsey -- Add report of excitons if non-empty
51
// 20130314 M. Kelsey -- Use G4IonList typedef for fragment map, encapsulate
52
// it in a static function with mutexes.
53
// 20130620 M. Kelsey -- Address Coverity #37503, check self in op=()
54
// 20140523 M. Kelsey -- Avoid FPE in setExcitationEnergy() for zero Ekin
55
// 20150608 M. Kelsey -- Label all while loops as terminating.
56
57
#include "
G4InuclNuclei.hh
"
58
#include "
G4AutoLock.hh
"
59
#include "
G4Fragment.hh
"
60
#include "
G4HadronicException.hh
"
61
#include "
G4InuclSpecialFunctions.hh
"
62
#include "
G4IonTable.hh
"
63
#include "
G4Ions.hh
"
64
#include "
G4NucleiProperties.hh
"
65
#include "
G4Nucleon.hh
"
66
#include "
G4ParticleDefinition.hh
"
67
#include "
G4ParticleTable.hh
"
68
#include "
G4SystemOfUnits.hh
"
69
#include "
G4Threading.hh
"
70
#include "
G4V3DNucleus.hh
"
71
72
#include <assert.h>
73
#include <sstream>
74
#include <map>
75
76
using namespace
G4InuclSpecialFunctions;
77
78
79
// Convert contents from (via constructor) and to G4Fragment
80
81
G4InuclNuclei::G4InuclNuclei
(
const
G4Fragment
& aFragment,
82
G4InuclParticle::Model
model
)
83
:
G4InuclParticle
() {
84
copy
(aFragment, model);
85
}
86
87
void
G4InuclNuclei::copy
(
const
G4Fragment
& aFragment,
Model
model
) {
88
fill
(aFragment.
GetMomentum
()/
GeV
, aFragment.
GetA_asInt
(),
89
aFragment.
GetZ_asInt
(), aFragment.
GetExcitationEnergy
(),
model
);
90
91
// Exciton configuration must be set by hand
92
theExitonConfiguration
.
protonQuasiParticles
= aFragment.
GetNumberOfCharged
();
93
94
theExitonConfiguration
.
neutronQuasiParticles
=
95
aFragment.
GetNumberOfParticles
() - aFragment.
GetNumberOfCharged
();
96
97
theExitonConfiguration
.
protonHoles
= aFragment.
GetNumberOfChargedHoles
();
98
99
theExitonConfiguration
.
neutronHoles
=
100
aFragment.
GetNumberOfHoles
() -
theExitonConfiguration
.
protonHoles
;
101
}
102
103
104
// FIXME: Should we have a local buffer and return by const-reference instead?
105
G4Fragment
G4InuclNuclei::makeG4Fragment
()
const
{
106
G4Fragment
frag(
getA
(),
getZ
(),
getMomentum
()*
GeV
);
// From Bertini units
107
108
// Note: exciton configuration has to be set piece by piece
109
frag.
SetNumberOfHoles
(
theExitonConfiguration
.
protonHoles
110
+
theExitonConfiguration
.
neutronHoles
,
111
theExitonConfiguration
.
protonHoles
);
112
113
frag.
SetNumberOfExcitedParticle
(
theExitonConfiguration
.
protonQuasiParticles
114
+
theExitonConfiguration
.
neutronQuasiParticles
,
115
theExitonConfiguration
.
protonQuasiParticles
);
116
117
return
frag;
118
}
119
120
G4InuclNuclei::operator
G4Fragment
()
const
{
121
return
makeG4Fragment();
122
}
123
124
125
// Convert contents from (via constructor) G4V3DNucleus
126
127
G4InuclNuclei::G4InuclNuclei
(
G4V3DNucleus
* a3DNucleus,
128
G4InuclParticle::Model
model
)
129
:
G4InuclParticle
() {
130
copy
(a3DNucleus, model);
131
}
132
133
void
G4InuclNuclei::copy
(
G4V3DNucleus
* a3DNucleus,
Model
model
) {
134
if
(!a3DNucleus)
return
;
// Null pointer means no action
135
136
fill
(0., a3DNucleus->
GetMassNumber
(), a3DNucleus->
GetCharge
(), 0.,
model
);
137
138
// Convert every hit nucleon into an exciton hole
139
if
(a3DNucleus->
StartLoop
()) {
140
G4Nucleon
* nucl = 0;
141
142
/* Loop checking 08.06.2015 MHK */
143
while
((nucl = a3DNucleus->
GetNextNucleon
())) {
144
if
(nucl->
AreYouHit
()) {
// Found previously interacted nucleon
145
if
(nucl->
GetParticleType
() ==
G4Proton::Definition
())
146
theExitonConfiguration
.
protonHoles
++;
147
148
if
(nucl->
GetParticleType
() ==
G4Neutron::Definition
())
149
theExitonConfiguration
.
neutronHoles
++;
150
}
151
}
152
}
153
}
154
155
156
// Overwrite data structure (avoids creating/copying temporaries)
157
158
void
G4InuclNuclei::fill
(
const
G4LorentzVector
&
mom
,
G4int
a
,
G4int
z
,
159
G4double
exc,
G4InuclParticle::Model
model
) {
160
setDefinition
(
makeDefinition
(a,z));
161
setMomentum
(mom);
162
setExitationEnergy
(exc);
163
clearExitonConfiguration
();
164
setModel
(model);
165
}
166
167
void
G4InuclNuclei::fill
(
G4double
ekin,
G4int
a
,
G4int
z
,
G4double
exc,
168
G4InuclParticle::Model
model
) {
169
setDefinition
(
makeDefinition
(a,z));
170
setKineticEnergy
(ekin);
171
setExitationEnergy
(exc);
172
clearExitonConfiguration
();
173
setModel
(model);
174
}
175
176
void
G4InuclNuclei::clear
() {
177
setDefinition
(0);
178
clearExitonConfiguration
();
179
setModel
(
G4InuclParticle::DefaultModel
);
180
}
181
182
183
// Change excitation energy while keeping momentum vector constant
184
185
void
G4InuclNuclei::setExitationEnergy
(
G4double
e
) {
186
G4double
ekin =
getKineticEnergy
();
// Current kinetic energy
187
188
G4double
emass =
getNucleiMass
() + e*
MeV
/
GeV
;
// From Bertini to G4 units
189
190
// Safety check -- if zero energy, don't do computation
191
G4double
ekin_new = (ekin == 0.) ? 0.
192
: std::sqrt(emass*emass + ekin*(2.*
getMass
()+ekin)) - emass;
193
194
setMass
(emass);
// Momentum is computed from mass and Ekin
195
setKineticEnergy
(ekin_new);
196
}
197
198
199
// Convert nuclear configuration to standard GEANT4 pointer
200
201
// WARNING: Opposite conventions! G4InuclNuclei uses (A,Z) everywhere, while
202
// G4ParticleTable::GetIon() uses (Z,A)!
203
204
G4ParticleDefinition
*
G4InuclNuclei::makeDefinition
(
G4int
a
,
G4int
z
) {
205
// SPECIAL CASE: (0,0) means create dummy without definition
206
if
(0 == a && 0 == z)
return
0;
207
208
G4ParticleTable
* pTable =
G4ParticleTable::GetParticleTable
();
209
G4ParticleDefinition
*pd = pTable->
GetIonTable
()->
GetIon
(z, a, 0);
210
211
// SPECIAL CASE: Non-physical nuclear fragment, for final-state return
212
if
(!pd) pd =
makeNuclearFragment
(a,z);
213
214
return
pd;
// This could return a null pointer if above fails
215
}
216
217
218
// Shared buffer of nuclear fragments created below, to avoid memory leaks
219
220
namespace
{
221
static
std::map<G4int,G4ParticleDefinition*> fragmentList;
222
G4Mutex
fragListMutex =
G4MUTEX_INITIALIZER
;
223
}
224
225
// Creates a non-physical pseudo-nucleus, for return as final-state fragment
226
// from G4IntraNuclearCascader
227
228
G4ParticleDefinition
*
229
G4InuclNuclei::makeNuclearFragment
(
G4int
a
,
G4int
z
) {
230
if
(a<=0 || z<0 || a<z) {
231
G4cerr
<<
" >>> G4InuclNuclei::makeNuclearFragment() called with"
232
<<
" impossible arguments A="
<< a <<
" Z="
<< z <<
G4endl
;
233
throw
G4HadronicException
(__FILE__, __LINE__,
234
"G4InuclNuclei impossible A/Z arguments"
);
235
}
236
237
G4int
code
=
G4IonTable::GetNucleusEncoding
(z, a);
238
239
// Use local lookup table (see above) to maintain singletons
240
// NOTE: G4ParticleDefinitions don't need to be explicitly deleted
241
// (see comments in G4IonTable.cc::~G4IonTable)
242
243
G4AutoLock
fragListLock(&fragListMutex);
244
if
(fragmentList.find(code) != fragmentList.end())
return
fragmentList[code];
245
fragListLock.
unlock
();
246
247
// Name string follows format in G4IonTable.cc::GetIonName(Z,A,E)
248
std::stringstream zstr, astr;
249
zstr <<
z
;
250
astr <<
a
;
251
252
G4String
name
=
"Z"
+ zstr.str() +
"A"
+ astr.str();
253
254
G4double
mass
=
getNucleiMass
(a,z) *
GeV
/
MeV
;
// From Bertini to GEANT4 units
255
256
// Arguments for constructor are as follows
257
// name mass width charge
258
// 2*spin parity C-conjugation
259
// 2*Isospin 2*Isospin3 G-parity
260
// type lepton number baryon number PDG encoding
261
// stable lifetime decay table
262
// shortlived subType anti_encoding Excitation-energy
263
264
G4Ions
* fragPD =
new
G4Ions
(name, mass, 0., z*
eplus
,
265
0, +1, 0,
266
0, 0, 0,
267
"nucleus"
, 0, a, code,
268
true
, 0., 0,
269
true
,
"generic"
, 0, 0.);
270
fragPD->
SetAntiPDGEncoding
(0);
271
272
fragListLock.
lock
();
// Protect before saving new fragment
273
return
(fragmentList[code] = fragPD);
// Store in table for next lookup
274
}
275
276
G4double
G4InuclNuclei::getNucleiMass
(
G4int
a
,
G4int
z
,
G4double
exc) {
277
// Simple minded mass calculation use constants in CLHEP (all in MeV)
278
G4double
mass
=
G4NucleiProperties::GetNuclearMass
(a,z) + exc;
279
280
return
mass*
MeV
/
GeV
;
// Convert from GEANT4 to Bertini units
281
}
282
283
// Assignment operator for use with std::sort()
284
G4InuclNuclei
&
G4InuclNuclei::operator=
(
const
G4InuclNuclei
&
right
) {
285
if
(
this
!= &right) {
286
theExitonConfiguration
= right.
theExitonConfiguration
;
287
G4InuclParticle::operator=
(right);
288
}
289
return
*
this
;
290
}
291
292
// Dump particle properties for diagnostics
293
294
void
G4InuclNuclei::print
(std::ostream& os)
const
{
295
G4InuclParticle::print
(os);
296
os <<
G4endl
<<
" Nucleus: "
<<
getDefinition
()->
GetParticleName
()
297
<<
" A "
<<
getA
() <<
" Z "
<<
getZ
() <<
" mass "
<<
getMass
()
298
<<
" Eex (MeV) "
<<
getExitationEnergy
();
299
300
if
(!
theExitonConfiguration
.
empty
())
301
os <<
G4endl
<<
" "
<<
theExitonConfiguration
;
302
}
geant4
tree
geant4-10.6-release
source
processes
hadronic
models
cascade
cascade
src
G4InuclNuclei.cc
Built by
Jin Huang
. updated:
Wed Jun 29 2022 17:25:40
using
1.8.2 with
ECCE GitHub integration