Geant4  10.03.p02
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
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 102724 2017-02-20 13:00:39Z 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  inline void SetExcEnergyAndMomentum(G4double eexc, const G4LorentzVector&);
109 
110  inline G4double GetGroundStateMass() const;
111 
112  inline G4double GetBindingEnergy() const;
113 
114  inline const G4LorentzVector& GetMomentum() const;
115  inline void SetMomentum(const G4LorentzVector& value);
116 
117  // computation of mass for any Z and A
118  inline G4double ComputeGroundStateMass(G4int Z, G4int A) const;
119 
120  // extra methods
121  inline G4double GetSpin() const;
122  inline void SetSpin(G4double value);
123 
124  inline G4int GetCreatorModelType() const;
125  inline void SetCreatorModelType(G4int value);
126 
127  // obsolete methods
128 
129  inline G4double GetZ() const;
130  inline G4double GetA() const;
131  inline void SetZ(G4double value);
132  inline void SetA(G4double value);
133 
134  // ============= METHODS FOR PRE-COMPOUND MODEL ===============
135 
136  inline G4int GetNumberOfExcitons() const;
137 
138  inline G4int GetNumberOfParticles() const;
139  inline G4int GetNumberOfCharged() const;
140  inline void SetNumberOfExcitedParticle(G4int valueTot, G4int valueP);
141 
142  inline G4int GetNumberOfHoles() const;
143  inline G4int GetNumberOfChargedHoles() const;
144  inline void SetNumberOfHoles(G4int valueTot, G4int valueP=0);
145 
146  // these methods will be removed in future
147  inline void SetNumberOfParticles(G4int value);
148  inline void SetNumberOfCharged(G4int value);
149 
150  // ============= METHODS FOR PHOTON EVAPORATION ===============
151 
152  inline G4int GetNumberOfElectrons() const;
153  inline void SetNumberOfElectrons(G4int value);
154 
155  inline G4int GetFloatingLevelNumber() const;
156  inline void SetFloatingLevelNumber(G4int value);
157 
158  inline const G4ParticleDefinition * GetParticleDefinition() const;
159  inline void SetParticleDefinition(const G4ParticleDefinition * p);
160 
161  inline G4double GetCreationTime() const;
162  inline void SetCreationTime(G4double time);
163 
166 
167  void SetAngularMomentum(const G4ThreeVector&);
169 
170  // ============= PRIVATE METHODS ==============================
171 
172 private:
173 
174  void ExcitationEnergyWarning();
175 
176  void NumberOfExitationWarning(const G4String&);
177 
178  inline void CalculateExcitationEnergy();
179 
180  inline void CalculateGroundStateMass();
181 
182  // ============= DATA MEMBERS ==================
183 
184  G4int theA;
185 
186  G4int theZ;
187 
188  G4double theExcitationEnergy;
189 
190  G4double theGroundStateMass;
191 
192  G4LorentzVector theMomentum;
193 
194  // Nuclear polarisation by default is nullptr
195  G4NuclearPolarization* thePolarization;
196 
197  // creator model type
198  G4int creatorModel;
199 
200  // Exciton model data members
201  G4int numberOfParticles;
202  G4int numberOfCharged;
203  G4int numberOfHoles;
204  G4int numberOfChargedHoles;
205 
206  // Gamma evaporation data members
207  G4int numberOfShellElectrons;
208  G4int xLevel;
209 
210  const G4ParticleDefinition* theParticleDefinition;
211 
212  G4double spin;
213  G4double theCreationTime;
214 
215  static const G4double minFragExcitation;
216 };
217 
218 // ============= INLINE METHOD IMPLEMENTATIONS ===================
219 
221 {
222  if(p != thePolarization) {
223  delete thePolarization;
224  thePolarization = p;
225  }
226 }
227 
228 #if defined G4HADRONIC_ALLOC_EXPORT
230 #else
232 #endif
233 
234 inline void * G4Fragment::operator new(size_t)
235 {
237  return (void*) pFragmentAllocator->MallocSingle();
238 }
239 
240 inline void G4Fragment::operator delete(void * aFragment)
241 {
242  ((G4Fragment *)aFragment)->SetNuclearPolarization(nullptr);
243  pFragmentAllocator->FreeSingle((G4Fragment *) aFragment);
244 }
245 
246 inline void G4Fragment::CalculateExcitationEnergy()
247 {
248  theExcitationEnergy = theMomentum.mag() - theGroundStateMass;
249  if(theExcitationEnergy < minFragExcitation) {
250  if(theExcitationEnergy < -minFragExcitation) { ExcitationEnergyWarning(); }
251  theExcitationEnergy = 0.0;
252  }
253 }
254 
255 inline G4double
257 {
259 }
260 
261 inline void G4Fragment::CalculateGroundStateMass()
262 {
263  theGroundStateMass = G4NucleiProperties::GetNuclearMass(theA, theZ);
264 }
265 
267 {
268  return theA;
269 }
270 
272 {
273  return theZ;
274 }
275 
276 inline void G4Fragment::SetZandA_asInt(G4int Znew, G4int Anew)
277 {
278  theZ = Znew;
279  theA = Anew;
280  CalculateGroundStateMass();
281 }
282 
284 {
285  return theExcitationEnergy;
286 }
287 
289 {
290  return theGroundStateMass;
291 }
292 
294  const G4LorentzVector& v)
295 {
296  theExcitationEnergy = eexc;
297  theMomentum.set(0.0, 0.0, 0.0, theGroundStateMass + eexc);
298  theMomentum.boost(v.boostVector());
299 }
300 
302 {
303  return (theA-theZ)*CLHEP::neutron_mass_c2 + theZ*CLHEP::proton_mass_c2
304  - theGroundStateMass;
305 }
306 
308 {
309  return theMomentum;
310 }
311 
313 {
314  theMomentum = value;
315  CalculateExcitationEnergy();
316 }
317 
319 {
320  return G4double(theZ);
321 }
322 
324 {
325  return G4double(theA);
326 }
327 
328 inline void G4Fragment::SetZ(const G4double value)
329 {
330  theZ = G4lrint(value);
331  CalculateGroundStateMass();
332 }
333 
334 inline void G4Fragment::SetA(const G4double value)
335 {
336  theA = G4lrint(value);
337  CalculateGroundStateMass();
338 }
339 
341 {
342  return numberOfParticles + numberOfHoles;
343 }
344 
346 {
347  return numberOfParticles;
348 }
349 
351 {
352  return numberOfCharged;
353 }
354 
355 inline
357 {
358  numberOfParticles = valueTot;
359  numberOfCharged = valueP;
360  if(valueTot < valueP) {
361  NumberOfExitationWarning("SetNumberOfExcitedParticle");
362  }
363 }
364 
366 {
367  return numberOfHoles;
368 }
369 
371 {
372  return numberOfChargedHoles;
373 }
374 
375 inline void G4Fragment::SetNumberOfHoles(G4int valueTot, G4int valueP)
376 {
377  numberOfHoles = valueTot;
378  numberOfChargedHoles = valueP;
379  if(valueTot < valueP) {
380  NumberOfExitationWarning("SetNumberOfHoles");
381  }
382 }
383 
385 {
386  numberOfParticles = value;
387 }
388 
390 {
391  numberOfCharged = value;
392  if(value > numberOfParticles) {
393  NumberOfExitationWarning("SetNumberOfCharged");
394  }
395 }
396 
398 {
399  return numberOfShellElectrons;
400 }
401 
403 {
404  numberOfShellElectrons = value;
405 }
406 
408 {
409  return creatorModel;
410 }
411 
413 {
414  creatorModel = value;
415 }
416 
418 {
419  return spin;
420 }
421 
423 {
424  spin = value;
425 }
426 
428 {
429  return xLevel;
430 }
431 
433 {
434  xLevel = value;
435 }
436 
437 inline
439 {
440  return theParticleDefinition;
441 }
442 
444 {
445  theParticleDefinition = p;
446 }
447 
449 {
450  return theCreationTime;
451 }
452 
454 {
455  theCreationTime = time;
456 }
457 
459 {
460  return thePolarization;
461 }
462 
463 #endif
464 
465 
Hep3Vector boostVector() const
G4int GetFloatingLevelNumber() const
Definition: G4Fragment.hh:427
static G4double GetNuclearMass(const G4double A, const G4double Z)
#define G4DLLEXPORT
Definition: G4Types.hh:62
static constexpr double proton_mass_c2
friend std::ostream & operator<<(std::ostream &, const G4Fragment *)
Definition: G4Fragment.cc:187
G4double GetA() const
Definition: G4Fragment.hh:323
#define G4DLLIMPORT
Definition: G4Types.hh:63
G4int GetNumberOfElectrons() const
Definition: G4Fragment.hh:397
const char * p
Definition: xmltok.h:285
const G4ParticleDefinition * GetParticleDefinition() const
Definition: G4Fragment.hh:438
void SetNumberOfHoles(G4int valueTot, G4int valueP=0)
Definition: G4Fragment.hh:375
void SetNumberOfElectrons(G4int value)
Definition: G4Fragment.hh:402
G4int GetNumberOfParticles() const
Definition: G4Fragment.hh:345
G4ThreeVector GetAngularMomentum() const
Definition: G4Fragment.cc:271
#define G4ThreadLocal
Definition: tls.hh:89
int G4int
Definition: G4Types.hh:78
void SetA(G4double value)
Definition: G4Fragment.hh:334
void SetCreatorModelType(G4int value)
Definition: G4Fragment.hh:412
G4int GetNumberOfHoles() const
Definition: G4Fragment.hh:365
void SetParticleDefinition(const G4ParticleDefinition *p)
Definition: G4Fragment.hh:443
void SetNumberOfExcitedParticle(G4int valueTot, G4int valueP)
Definition: G4Fragment.hh:356
double A(double temperature)
G4int GetA_asInt() const
Definition: G4Fragment.hh:266
G4bool operator==(const G4Fragment &right) const
Definition: G4Fragment.cc:177
void SetNuclearPolarization(G4NuclearPolarization *)
Definition: G4Fragment.hh:220
const XML_Char int const XML_Char * value
Definition: expat.h:331
G4double GetCreationTime() const
Definition: G4Fragment.hh:448
G4DLLIMPORT G4ThreadLocal G4Allocator< G4Fragment > * pFragmentAllocator
Definition: G4Fragment.cc:46
double mag() const
bool G4bool
Definition: G4Types.hh:79
const G4LorentzVector & GetMomentum() const
Definition: G4Fragment.hh:307
HepLorentzVector & boost(double, double, double)
void SetMomentum(const G4LorentzVector &value)
Definition: G4Fragment.hh:312
static constexpr double neutron_mass_c2
std::vector< G4Fragment * > G4FragmentVector
Definition: G4Fragment.hh:63
G4double GetGroundStateMass() const
Definition: G4Fragment.hh:288
void SetNumberOfParticles(G4int value)
Definition: G4Fragment.hh:384
tuple v
Definition: test.py:18
G4double GetBindingEnergy() const
Definition: G4Fragment.hh:301
G4double GetSpin() const
Definition: G4Fragment.hh:417
void SetCreationTime(G4double time)
Definition: G4Fragment.hh:453
void set(double x, double y, double z, double t)
G4int GetNumberOfExcitons() const
Definition: G4Fragment.hh:340
int G4lrint(double ad)
Definition: templates.hh:163
void SetZandA_asInt(G4int Znew, G4int Anew)
Definition: G4Fragment.hh:276
G4bool operator!=(const G4Fragment &right) const
Definition: G4Fragment.cc:182
G4int GetZ_asInt() const
Definition: G4Fragment.hh:271
void SetFloatingLevelNumber(G4int value)
Definition: G4Fragment.hh:432
void SetExcEnergyAndMomentum(G4double eexc, const G4LorentzVector &)
Definition: G4Fragment.hh:293
void SetZ(G4double value)
Definition: G4Fragment.hh:328
void SetSpin(G4double value)
Definition: G4Fragment.hh:422
G4int GetNumberOfChargedHoles() const
Definition: G4Fragment.hh:370
void SetNumberOfCharged(G4int value)
Definition: G4Fragment.hh:389
void SetAngularMomentum(const G4ThreeVector &)
Definition: G4Fragment.cc:266
double G4double
Definition: G4Types.hh:76
G4Fragment & operator=(const G4Fragment &right)
Definition: G4Fragment.cc:150
G4int GetCreatorModelType() const
Definition: G4Fragment.hh:407
G4int GetNumberOfCharged() const
Definition: G4Fragment.hh:350
G4double ComputeGroundStateMass(G4int Z, G4int A) const
Definition: G4Fragment.hh:256
G4double GetExcitationEnergy() const
Definition: G4Fragment.hh:283
G4NuclearPolarization * GetNuclearPolarization() const
Definition: G4Fragment.hh:458
G4double GetZ() const
Definition: G4Fragment.hh:318