Geant4  10.00.p02
G4Fragment.hh
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: G4Fragment.hh 79213 2014-02-20 14:43:12Z gcosmo $
27 //
28 //---------------------------------------------------------------------
29 //
30 // Geant4 header G4Fragment
31 //
32 // by V. Lara (May 1998)
33 //
34 // Modifications:
35 // 03.05.2010 V.Ivanchenko General cleanup of inline functions: objects
36 // are accessed by reference; remove double return
37 // tolerance of excitation energy at modent it is computed;
38 // safe computation of excitation for exotic fragments
39 // 18.05.2010 V.Ivanchenko added member theGroundStateMass and inline
40 // method which allowing to compute this value once and use
41 // many times
42 // 26.09.2010 V.Ivanchenko added number of protons, neutrons, proton holes
43 // and neutron holes as members of the class and Get/Set methods;
44 // removed not needed 'const'; removed old debug staff and unused
45 // private methods; add comments and reorder methods for
46 // better reading
47 
48 #ifndef G4Fragment_h
49 #define G4Fragment_h 1
50 
51 #include <vector>
52 
53 #include "globals.hh"
54 #include "G4Allocator.hh"
55 #include "G4LorentzVector.hh"
56 #include "G4ThreeVector.hh"
57 #include "G4NucleiProperties.hh"
58 #include "Randomize.hh"
59 #include "G4Proton.hh"
60 #include "G4Neutron.hh"
61 #include "G4HadTmpUtil.hh"
62 
64 
65 class G4Fragment;
66 typedef std::vector<G4Fragment*> G4FragmentVector;
67 
68 class G4Fragment
69 {
70 public:
71 
72  // ============= CONSTRUCTORS ==================
73 
74  // Default constructor - obsolete
75  G4Fragment();
76 
77  // Destructor
78  ~G4Fragment();
79 
80  // Copy constructor
81  G4Fragment(const G4Fragment &right);
82 
83  // A,Z and 4-momentum - main constructor for fragment
84  G4Fragment(G4int A, G4int Z, const G4LorentzVector& aMomentum);
85 
86  // 4-momentum and pointer to G4particleDefinition (for gammas, e-)
87  G4Fragment(const G4LorentzVector& aMomentum,
88  G4ParticleDefinition* aParticleDefinition);
89 
90  // ============= OPERATORS ==================
91 
92  G4Fragment & operator=(const G4Fragment &right);
93  G4bool operator==(const G4Fragment &right) const;
94  G4bool operator!=(const G4Fragment &right) const;
95 
96  friend std::ostream& operator<<(std::ostream&, const G4Fragment*);
97  friend std::ostream& operator<<(std::ostream&, const G4Fragment&);
98 
99  // new/delete operators are overloded to use G4Allocator
100  inline void *operator new(size_t);
101  inline void operator delete(void *aFragment);
102 
103  // ============= GENERAL METHODS ==================
104 
105  inline G4int GetZ_asInt() const;
106  inline G4int GetA_asInt() const;
107  inline void SetZandA_asInt(G4int Znew, G4int Anew);
108 
109  inline G4double GetExcitationEnergy() const;
110 
111  inline G4double GetGroundStateMass() const;
112 
113  inline G4double GetBindingEnergy() const;
114 
115  inline const G4LorentzVector& GetMomentum() const;
116  inline void SetMomentum(const G4LorentzVector& value);
117 
118  inline const G4ThreeVector& GetAngularMomentum() const;
119  inline void SetAngularMomentum(const G4ThreeVector& value);
120 
121  // computation of mass for any Z and A
122  inline G4double ComputeGroundStateMass(G4int Z, G4int A) const;
123 
124  // obsolete methods
125  inline G4double GetZ() const;
126  inline G4double GetA() const;
127  inline void SetZ(G4double value);
128  inline void SetA(G4double value);
129 
130  // ============= METHODS FOR PRE-COMPOUND MODEL ===============
131 
132  inline G4int GetNumberOfExcitons() const;
133 
134  inline G4int GetNumberOfParticles() const;
135  inline G4int GetNumberOfCharged() const;
136  inline void SetNumberOfExcitedParticle(G4int valueTot, G4int valueP);
137 
138  inline G4int GetNumberOfHoles() const;
139  inline G4int GetNumberOfChargedHoles() const;
140  inline void SetNumberOfHoles(G4int valueTot, G4int valueP=0);
141 
142  // these methods will be removed in future
143  inline void SetNumberOfParticles(G4int value);
144  inline void SetNumberOfCharged(G4int value);
145 
146  // ============= METHODS FOR PHOTON EVAPORATION ===============
147 
148  inline G4int GetNumberOfElectrons() const;
149  inline void SetNumberOfElectrons(G4int value);
150 
153 
154  inline G4double GetCreationTime() const;
155  inline void SetCreationTime(G4double time);
156 
157  inline G4bool IsStable() const;
158  inline void SetStable(G4bool val);
159 
160  // ============= PRIVATE METHODS ==============================
161 
162 private:
163 
165 
166  void NumberOfExitationWarning(const G4String&);
167 
168  inline void CalculateExcitationEnergy();
169 
170  inline void CalculateGroundStateMass();
171 
172  // ============= DATA MEMBERS ==================
173 
175 
177 
179 
181 
183 
185 
186  // Exciton model data members
187 
189 
191 
193 
195 
196  // Gamma evaporation data members
197 
199 
201 
203 
205 
206 };
207 
208 // ============= INLINE METHOD IMPLEMENTATIONS ===================
209 
210 #if defined G4HADRONIC_ALLOC_EXPORT
212 #else
214 #endif
215 
216 inline void * G4Fragment::operator new(size_t)
217 {
219  return (void*) pFragmentAllocator->MallocSingle();
220 }
221 
222 inline void G4Fragment::operator delete(void * aFragment)
223 {
224  pFragmentAllocator->FreeSingle((G4Fragment *) aFragment);
225 }
226 
228 {
231 }
232 
234 {
236 }
237 
239 {
240  return theA;
241 }
242 
244 {
245  return theZ;
246 }
247 
248 inline void G4Fragment::SetZandA_asInt(G4int Znew, G4int Anew)
249 {
250  theZ = Znew;
251  theA = Anew;
253 }
254 
256 {
257  return theExcitationEnergy;
258 }
259 
261 {
262  return theGroundStateMass;
263 }
264 
266 {
267  return (theA-theZ)*CLHEP::neutron_mass_c2 + theZ*CLHEP::proton_mass_c2
269 }
270 
272 {
273  return theMomentum;
274 }
275 
276 inline void G4Fragment::SetMomentum(const G4LorentzVector& value)
277 {
278  theMomentum = value;
280 }
281 
283 {
284  return theAngularMomentum;
285 }
286 
288 {
289  theAngularMomentum = value;
290 }
291 
292 inline G4double
294 {
296 }
297 
299 {
300  return G4double(theZ);
301 }
302 
304 {
305  return G4double(theA);
306 }
307 
308 inline void G4Fragment::SetZ(const G4double value)
309 {
310  theZ = G4lrint(value);
312 }
313 
314 inline void G4Fragment::SetA(const G4double value)
315 {
316  theA = G4lrint(value);
318 }
319 
321 {
323 }
324 
326 {
327  return numberOfParticles;
328 }
329 
331 {
332  return numberOfCharged;
333 }
334 
335 inline
337 {
338  numberOfParticles = valueTot;
339  numberOfCharged = valueP;
340  if(valueTot < valueP) {
341  NumberOfExitationWarning("SetNumberOfExcitedParticle");
342  }
343 }
344 
346 {
347  return numberOfHoles;
348 }
349 
351 {
352  return numberOfChargedHoles;
353 }
354 
355 inline void G4Fragment::SetNumberOfHoles(G4int valueTot, G4int valueP)
356 {
357  numberOfHoles = valueTot;
358  numberOfChargedHoles = valueP;
359  if(valueTot < valueP) {
360  NumberOfExitationWarning("SetNumberOfHoles");
361  }
362 }
363 
365 {
366  numberOfParticles = value;
367 }
368 
370 {
371  numberOfCharged = value;
372  if(value > numberOfParticles) {
373  NumberOfExitationWarning("SetNumberOfCharged");
374  }
375 }
376 
378 {
379  return numberOfShellElectrons;
380 }
381 
383 {
384  numberOfShellElectrons = value;
385 }
386 
387 inline
389 {
390  return theParticleDefinition;
391 }
392 
394 {
396 }
397 
399 {
400  return theCreationTime;
401 }
402 
404 {
405  theCreationTime = time;
406 }
407 
409 {
410  return isStable;
411 }
412 
414 {
415  isStable = val;
416 }
417 
418 #endif
419 
420 
G4int theZ
Definition: G4Fragment.hh:176
G4double theCreationTime
Definition: G4Fragment.hh:202
G4double theGroundStateMass
Definition: G4Fragment.hh:180
static G4double GetNuclearMass(const G4double A, const G4double Z)
G4int numberOfCharged
Definition: G4Fragment.hh:190
#define G4DLLEXPORT
Definition: G4Types.hh:62
friend std::ostream & operator<<(std::ostream &, const G4Fragment *)
Definition: G4Fragment.cc:175
CLHEP::Hep3Vector G4ThreeVector
void ExcitationEnergyWarning()
Definition: G4Fragment.cc:225
G4double GetA() const
Definition: G4Fragment.hh:303
#define G4DLLIMPORT
Definition: G4Types.hh:63
const G4ThreeVector & GetAngularMomentum() const
Definition: G4Fragment.hh:282
G4int GetNumberOfElectrons() const
Definition: G4Fragment.hh:377
void SetAngularMomentum(const G4ThreeVector &value)
Definition: G4Fragment.hh:287
void SetNumberOfHoles(G4int valueTot, G4int valueP=0)
Definition: G4Fragment.hh:355
void SetNumberOfElectrons(G4int value)
Definition: G4Fragment.hh:382
G4int GetNumberOfParticles() const
Definition: G4Fragment.hh:325
#define G4ThreadLocal
Definition: tls.hh:52
int G4int
Definition: G4Types.hh:78
void SetA(G4double value)
Definition: G4Fragment.hh:314
void SetParticleDefinition(G4ParticleDefinition *p)
Definition: G4Fragment.hh:393
G4int numberOfShellElectrons
Definition: G4Fragment.hh:198
G4int GetNumberOfHoles() const
Definition: G4Fragment.hh:345
void SetNumberOfExcitedParticle(G4int valueTot, G4int valueP)
Definition: G4Fragment.hh:336
G4int GetA_asInt() const
Definition: G4Fragment.hh:238
G4bool operator==(const G4Fragment &right) const
Definition: G4Fragment.cc:165
G4double GetCreationTime() const
Definition: G4Fragment.hh:398
G4DLLIMPORT G4ThreadLocal G4Allocator< G4Fragment > * pFragmentAllocator
Definition: G4Fragment.cc:50
bool G4bool
Definition: G4Types.hh:79
const G4LorentzVector & GetMomentum() const
Definition: G4Fragment.hh:271
void SetMomentum(const G4LorentzVector &value)
Definition: G4Fragment.hh:276
std::vector< G4Fragment * > G4FragmentVector
Definition: G4Fragment.hh:65
G4double GetGroundStateMass() const
Definition: G4Fragment.hh:260
void CalculateExcitationEnergy()
Definition: G4Fragment.hh:227
void SetNumberOfParticles(G4int value)
Definition: G4Fragment.hh:364
G4int theA
Definition: G4Fragment.hh:174
G4ThreeVector theAngularMomentum
Definition: G4Fragment.hh:184
static const G4double A[nN]
G4LorentzVector theMomentum
Definition: G4Fragment.hh:182
G4double GetBindingEnergy() const
Definition: G4Fragment.hh:265
G4int numberOfHoles
Definition: G4Fragment.hh:192
G4int numberOfParticles
Definition: G4Fragment.hh:188
void SetCreationTime(G4double time)
Definition: G4Fragment.hh:403
G4bool IsStable() const
Definition: G4Fragment.hh:408
G4int GetNumberOfExcitons() const
Definition: G4Fragment.hh:320
int G4lrint(double ad)
Definition: templates.hh:163
G4ParticleDefinition * theParticleDefinition
Definition: G4Fragment.hh:200
void SetZandA_asInt(G4int Znew, G4int Anew)
Definition: G4Fragment.hh:248
G4bool operator!=(const G4Fragment &right) const
Definition: G4Fragment.cc:170
G4int GetZ_asInt() const
Definition: G4Fragment.hh:243
G4ParticleDefinition * GetParticleDefinition() const
Definition: G4Fragment.hh:388
void SetStable(G4bool val)
Definition: G4Fragment.hh:413
void SetZ(G4double value)
Definition: G4Fragment.hh:308
void NumberOfExitationWarning(const G4String &)
Definition: G4Fragment.cc:244
G4int GetNumberOfChargedHoles() const
Definition: G4Fragment.hh:350
void SetNumberOfCharged(G4int value)
Definition: G4Fragment.hh:369
double G4double
Definition: G4Types.hh:76
void CalculateGroundStateMass()
Definition: G4Fragment.hh:233
G4Fragment & operator=(const G4Fragment &right)
Definition: G4Fragment.cc:144
G4int GetNumberOfCharged() const
Definition: G4Fragment.hh:330
G4double ComputeGroundStateMass(G4int Z, G4int A) const
Definition: G4Fragment.hh:293
G4bool isStable
Definition: G4Fragment.hh:204
G4int numberOfChargedHoles
Definition: G4Fragment.hh:194
G4double GetExcitationEnergy() const
Definition: G4Fragment.hh:255
G4double GetZ() const
Definition: G4Fragment.hh:298
G4double theExcitationEnergy
Definition: G4Fragment.hh:178
CLHEP::HepLorentzVector G4LorentzVector