Geant4  10.02
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 92198 2015-08-21 10:00:50Z 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 "globals.hh"
52 #include "G4Allocator.hh"
53 #include "G4LorentzVector.hh"
54 #include "G4ThreeVector.hh"
55 #include "G4NuclearPolarization.hh"
56 #include "G4NucleiProperties.hh"
57 #include "G4Proton.hh"
58 #include "G4Neutron.hh"
59 #include <vector>
60 
62 
63 class G4Fragment;
64 typedef std::vector<G4Fragment*> G4FragmentVector;
65 
66 class G4Fragment
67 {
68 public:
69 
70  // ============= CONSTRUCTORS ==================
71 
72  // Default constructor - obsolete
73  G4Fragment();
74 
75  // Destructor
76  ~G4Fragment();
77 
78  // Copy constructor
79  G4Fragment(const G4Fragment &right);
80 
81  // A,Z and 4-momentum - main constructor for fragment
82  G4Fragment(G4int A, G4int Z, const G4LorentzVector& aMomentum);
83 
84  // 4-momentum and pointer to G4particleDefinition (for gammas, e-)
85  G4Fragment(const G4LorentzVector& aMomentum,
86  const G4ParticleDefinition* aParticleDefinition);
87 
88  // ============= OPERATORS ==================
89 
90  G4Fragment & operator=(const G4Fragment &right);
91  G4bool operator==(const G4Fragment &right) const;
92  G4bool operator!=(const G4Fragment &right) const;
93 
94  friend std::ostream& operator<<(std::ostream&, const G4Fragment*);
95  friend std::ostream& operator<<(std::ostream&, const G4Fragment&);
96 
97  // new/delete operators are overloded to use G4Allocator
98  inline void *operator new(size_t);
99  inline void operator delete(void *aFragment);
100 
101  // ============= GENERAL METHODS ==================
102 
103  inline G4int GetZ_asInt() const;
104  inline G4int GetA_asInt() const;
105  inline void SetZandA_asInt(G4int Znew, G4int Anew);
106 
107  inline G4double GetExcitationEnergy() const;
108 
109  inline G4double GetGroundStateMass() const;
110 
111  inline G4double GetBindingEnergy() const;
112 
113  inline const G4LorentzVector& GetMomentum() const;
114  inline void SetMomentum(const G4LorentzVector& value);
115 
116  // computation of mass for any Z and A
117  inline G4double ComputeGroundStateMass(G4int Z, G4int A) const;
118 
119  inline G4int GetCreatorModelType() const;
120  inline void SetCreatorModelType(G4int value);
121 
122  // obsolete methods
123 
124  inline G4double GetZ() const;
125  inline G4double GetA() const;
126  inline void SetZ(G4double value);
127  inline void SetA(G4double value);
128 
129  // ============= METHODS FOR PRE-COMPOUND MODEL ===============
130 
131  inline G4int GetNumberOfExcitons() const;
132 
133  inline G4int GetNumberOfParticles() const;
134  inline G4int GetNumberOfCharged() const;
135  inline void SetNumberOfExcitedParticle(G4int valueTot, G4int valueP);
136 
137  inline G4int GetNumberOfHoles() const;
138  inline G4int GetNumberOfChargedHoles() const;
139  inline void SetNumberOfHoles(G4int valueTot, G4int valueP=0);
140 
141  // these methods will be removed in future
142  inline void SetNumberOfParticles(G4int value);
143  inline void SetNumberOfCharged(G4int value);
144 
145  // ============= METHODS FOR PHOTON EVAPORATION ===============
146 
147  inline G4int GetNumberOfElectrons() const;
148  inline void SetNumberOfElectrons(G4int value);
149 
150  inline const G4ParticleDefinition * GetParticleDefinition() const;
151  inline void SetParticleDefinition(const G4ParticleDefinition * p);
152 
153  inline G4double GetCreationTime() const;
154  inline void SetCreationTime(G4double time);
155 
158 
161 
162  // ============= PRIVATE METHODS ==============================
163 
164 private:
165 
167 
168  void NumberOfExitationWarning(const G4String&);
169 
170  inline void CalculateExcitationEnergy();
171 
172  inline void CalculateGroundStateMass();
173 
174  // ============= DATA MEMBERS ==================
175 
177 
179 
181 
183 
185 
186  // Nuclear polarisation by default is nullptr
188 
189  // creator model type
191 
192  // Exciton model data members
197 
198  // Gamma evaporation data members
200 
202 
204 };
205 
206 // ============= INLINE METHOD IMPLEMENTATIONS ===================
207 
209 {
210  if(p != thePolarization) {
211  delete thePolarization;
212  thePolarization = p;
213  }
214 }
215 
216 #if defined G4HADRONIC_ALLOC_EXPORT
218 #else
220 #endif
221 
222 inline void * G4Fragment::operator new(size_t)
223 {
225  return (void*) pFragmentAllocator->MallocSingle();
226 }
227 
228 inline void G4Fragment::operator delete(void * aFragment)
229 {
230  ((G4Fragment *)aFragment)->SetNuclearPolarization(0);
231  pFragmentAllocator->FreeSingle((G4Fragment *) aFragment);
232 }
233 
235 {
238 }
239 
240 inline G4double
242 {
244 }
245 
247 {
249 }
250 
252 {
253  return theA;
254 }
255 
257 {
258  return theZ;
259 }
260 
261 inline void G4Fragment::SetZandA_asInt(G4int Znew, G4int Anew)
262 {
263  theZ = Znew;
264  theA = Anew;
266 }
267 
269 {
270  return theExcitationEnergy;
271 }
272 
274 {
275  return theGroundStateMass;
276 }
277 
279 {
280  return (theA-theZ)*CLHEP::neutron_mass_c2 + theZ*CLHEP::proton_mass_c2
282 }
283 
285 {
286  return theMomentum;
287 }
288 
289 inline void G4Fragment::SetMomentum(const G4LorentzVector& value)
290 {
291  theMomentum = value;
293 }
294 
296 {
297  return G4double(theZ);
298 }
299 
301 {
302  return G4double(theA);
303 }
304 
305 inline void G4Fragment::SetZ(const G4double value)
306 {
307  theZ = G4lrint(value);
309 }
310 
311 inline void G4Fragment::SetA(const G4double value)
312 {
313  theA = G4lrint(value);
315 }
316 
318 {
320 }
321 
323 {
324  return numberOfParticles;
325 }
326 
328 {
329  return numberOfCharged;
330 }
331 
332 inline
334 {
335  numberOfParticles = valueTot;
336  numberOfCharged = valueP;
337  if(valueTot < valueP) {
338  NumberOfExitationWarning("SetNumberOfExcitedParticle");
339  }
340 }
341 
343 {
344  return numberOfHoles;
345 }
346 
348 {
349  return numberOfChargedHoles;
350 }
351 
352 inline void G4Fragment::SetNumberOfHoles(G4int valueTot, G4int valueP)
353 {
354  numberOfHoles = valueTot;
355  numberOfChargedHoles = valueP;
356  if(valueTot < valueP) {
357  NumberOfExitationWarning("SetNumberOfHoles");
358  }
359 }
360 
362 {
363  numberOfParticles = value;
364 }
365 
367 {
368  numberOfCharged = value;
369  if(value > numberOfParticles) {
370  NumberOfExitationWarning("SetNumberOfCharged");
371  }
372 }
373 
375 {
376  return numberOfShellElectrons;
377 }
378 
380 {
381  numberOfShellElectrons = value;
382 }
383 
385 {
386  return creatorModel;
387 }
388 
390 {
391  creatorModel = value;
392 }
393 
394 inline
396 {
397  return theParticleDefinition;
398 }
399 
401 {
403 }
404 
406 {
407  return theCreationTime;
408 }
409 
411 {
412  theCreationTime = time;
413 }
414 
416 {
417  return thePolarization;
418 }
419 
420 #endif
421 
422 
G4int theZ
Definition: G4Fragment.hh:178
G4double theCreationTime
Definition: G4Fragment.hh:203
G4double theGroundStateMass
Definition: G4Fragment.hh:182
static G4double GetNuclearMass(const G4double A, const G4double Z)
G4int numberOfCharged
Definition: G4Fragment.hh:194
G4NuclearPolarization * thePolarization
Definition: G4Fragment.hh:187
#define G4DLLEXPORT
Definition: G4Types.hh:62
friend std::ostream & operator<<(std::ostream &, const G4Fragment *)
Definition: G4Fragment.cc:176
CLHEP::Hep3Vector G4ThreeVector
void ExcitationEnergyWarning()
Definition: G4Fragment.cc:232
G4double GetA() const
Definition: G4Fragment.hh:300
#define G4DLLIMPORT
Definition: G4Types.hh:63
G4int GetNumberOfElectrons() const
Definition: G4Fragment.hh:374
const G4ParticleDefinition * GetParticleDefinition() const
Definition: G4Fragment.hh:395
void SetNumberOfHoles(G4int valueTot, G4int valueP=0)
Definition: G4Fragment.hh:352
void SetNumberOfElectrons(G4int value)
Definition: G4Fragment.hh:379
G4int GetNumberOfParticles() const
Definition: G4Fragment.hh:322
G4ThreeVector GetAngularMomentum() const
Definition: G4Fragment.cc:269
#define G4ThreadLocal
Definition: tls.hh:89
int G4int
Definition: G4Types.hh:78
void SetA(G4double value)
Definition: G4Fragment.hh:311
void SetCreatorModelType(G4int value)
Definition: G4Fragment.hh:389
G4int numberOfShellElectrons
Definition: G4Fragment.hh:199
G4int GetNumberOfHoles() const
Definition: G4Fragment.hh:342
void SetParticleDefinition(const G4ParticleDefinition *p)
Definition: G4Fragment.hh:400
void SetNumberOfExcitedParticle(G4int valueTot, G4int valueP)
Definition: G4Fragment.hh:333
double A(double temperature)
G4int GetA_asInt() const
Definition: G4Fragment.hh:251
G4bool operator==(const G4Fragment &right) const
Definition: G4Fragment.cc:166
void SetNuclearPolarization(G4NuclearPolarization *)
Definition: G4Fragment.hh:208
G4double GetCreationTime() const
Definition: G4Fragment.hh:405
G4DLLIMPORT G4ThreadLocal G4Allocator< G4Fragment > * pFragmentAllocator
Definition: G4Fragment.cc:48
bool G4bool
Definition: G4Types.hh:79
const G4LorentzVector & GetMomentum() const
Definition: G4Fragment.hh:284
void SetMomentum(const G4LorentzVector &value)
Definition: G4Fragment.hh:289
std::vector< G4Fragment * > G4FragmentVector
Definition: G4Fragment.hh:63
G4double GetGroundStateMass() const
Definition: G4Fragment.hh:273
void CalculateExcitationEnergy()
Definition: G4Fragment.hh:234
void SetNumberOfParticles(G4int value)
Definition: G4Fragment.hh:361
G4int theA
Definition: G4Fragment.hh:176
G4LorentzVector theMomentum
Definition: G4Fragment.hh:184
G4double GetBindingEnergy() const
Definition: G4Fragment.hh:278
G4int numberOfHoles
Definition: G4Fragment.hh:195
G4int creatorModel
Definition: G4Fragment.hh:190
G4int numberOfParticles
Definition: G4Fragment.hh:193
void SetCreationTime(G4double time)
Definition: G4Fragment.hh:410
G4int GetNumberOfExcitons() const
Definition: G4Fragment.hh:317
int G4lrint(double ad)
Definition: templates.hh:163
void SetZandA_asInt(G4int Znew, G4int Anew)
Definition: G4Fragment.hh:261
G4bool operator!=(const G4Fragment &right) const
Definition: G4Fragment.cc:171
G4int GetZ_asInt() const
Definition: G4Fragment.hh:256
void SetAngularMomentum(G4ThreeVector &)
Definition: G4Fragment.cc:260
void SetZ(G4double value)
Definition: G4Fragment.hh:305
void NumberOfExitationWarning(const G4String &)
Definition: G4Fragment.cc:251
G4int GetNumberOfChargedHoles() const
Definition: G4Fragment.hh:347
void SetNumberOfCharged(G4int value)
Definition: G4Fragment.hh:366
double G4double
Definition: G4Types.hh:76
void CalculateGroundStateMass()
Definition: G4Fragment.hh:246
G4Fragment & operator=(const G4Fragment &right)
Definition: G4Fragment.cc:142
G4int GetCreatorModelType() const
Definition: G4Fragment.hh:384
G4int GetNumberOfCharged() const
Definition: G4Fragment.hh:327
G4double ComputeGroundStateMass(G4int Z, G4int A) const
Definition: G4Fragment.hh:241
G4int numberOfChargedHoles
Definition: G4Fragment.hh:196
G4double GetExcitationEnergy() const
Definition: G4Fragment.hh:268
G4NuclearPolarization * GetNuclearPolarization() const
Definition: G4Fragment.hh:415
G4double GetZ() const
Definition: G4Fragment.hh:295
const G4ParticleDefinition * theParticleDefinition
Definition: G4Fragment.hh:201
G4double theExcitationEnergy
Definition: G4Fragment.hh:180
CLHEP::HepLorentzVector G4LorentzVector