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