| Geant4 Cross Reference |
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 // J.L. Chuma, TRIUMF, 31-Oct-1996
27 // last modified: 19-Dec-1996
28 // modified by J.L.Chuma, 24-Jul-1997 to include total momentum
29 // inluded operator *, and some minor modifications.
30 // modified by H.P.Wellisch to add functionality needed by string models,
31 // cascade and Nucleus. (Mon Mar 16 1998)
32 // M. Kelsey 29-Aug-2011 -- Use G4Allocator model to avoid memory churn.
33
34 #ifndef G4ReactionProduct_h
35 #define G4ReactionProduct_h 1
36
37 #include "globals.hh"
38 #include "G4Allocator.hh"
39 #include "G4DynamicParticle.hh"
40 #include "G4HadProjectile.hh"
41 #include "G4HadronicException.hh"
42
43 class G4ReactionProduct;
44
45 // To support better memory management and reduced fragmentation
46 //
47 #if defined G4HADRONIC_ALLOC_EXPORT
48 extern G4DLLEXPORT G4Allocator<G4ReactionProduct>*& aRPAllocator();
49 #else
50 extern G4DLLIMPORT G4Allocator<G4ReactionProduct>*& aRPAllocator();
51 #endif
52
53 class G4ReactionProduct
54 {
55 friend G4ReactionProduct operator+(
56 const G4ReactionProduct & p1, const G4ReactionProduct &p2 );
57
58 friend G4ReactionProduct operator-(
59 const G4ReactionProduct & p1, const G4ReactionProduct &p2 );
60
61 friend G4ReactionProduct operator*(
62 const G4double aDouble, const G4ReactionProduct &p2 )
63 {
64 G4ReactionProduct result;
65 result.SetMomentum(aDouble*p2.GetMomentum());
66 result.SetMass(p2.GetMass());
67 result.SetTotalEnergy(std::sqrt(result.GetMass()*result.GetMass()+
68 result.GetMomentum()*result.GetMomentum()));
69 return result;
70 }
71
72 public:
73 G4ReactionProduct();
74
75 G4ReactionProduct(const G4ParticleDefinition *aParticleDefinition );
76
77 ~G4ReactionProduct() {}
78
79 G4ReactionProduct( const G4ReactionProduct &right );
80
81 // Override new and delete for use with G4Allocator
82 inline void* operator new(size_t) {
83 if (!aRPAllocator()) aRPAllocator() = new G4Allocator<G4ReactionProduct> ;
84 return (void *)aRPAllocator()->MallocSingle();
85 }
86 #ifdef __IBMCPP__
87 inline void* operator new(size_t, void *p) {
88 return p;
89 }
90 #endif
91 inline void operator delete(void* aReactionProduct) {
92 aRPAllocator()->FreeSingle((G4ReactionProduct*)aReactionProduct);
93 }
94
95 G4ReactionProduct &operator= ( const G4ReactionProduct &right );
96
97 G4ReactionProduct &operator= ( const G4DynamicParticle &right );
98
99 G4ReactionProduct &operator= ( const G4HadProjectile &right );
100
101 inline G4bool operator== ( const G4ReactionProduct &right ) const
102 { return ( this == (G4ReactionProduct*) &right ); }
103
104 inline G4bool operator!= ( const G4ReactionProduct &right ) const
105 { return ( this != (G4ReactionProduct*) &right ); }
106
107 inline const G4ParticleDefinition* GetDefinition() const
108 { return theParticleDefinition; }
109
110 void SetDefinition(const G4ParticleDefinition* aParticleDefinition );
111
112 void SetDefinitionAndUpdateE(const G4ParticleDefinition* aParticleDefinition );
113
114 void SetMomentum( const G4double x, const G4double y, const G4double z );
115
116 void SetMomentum( const G4double x, const G4double y );
117
118 void SetMomentum( const G4double z );
119
120 inline void SetMomentum( const G4ThreeVector &mom )
121 { momentum = mom; }
122
123 inline G4ThreeVector GetMomentum() const
124 { return momentum; }
125
126 inline G4double GetTotalMomentum() const
127 { return std::sqrt(std::abs(kineticEnergy*(totalEnergy+mass))); }
128
129 inline G4double GetTotalEnergy() const
130 { return totalEnergy; }
131
132 inline void SetKineticEnergy( const G4double en )
133 {
134 kineticEnergy = en;
135 totalEnergy = kineticEnergy + mass;
136 }
137
138 inline G4double GetKineticEnergy() const
139 { return kineticEnergy; }
140
141 inline void SetTotalEnergy( const G4double en )
142 {
143 totalEnergy = en;
144 kineticEnergy = totalEnergy - mass;
145 }
146
147 inline void SetMass( const G4double mas )
148 { mass = mas; }
149
150 inline G4double GetMass() const
151 { return mass; }
152
153 inline void SetTOF( const G4double t )
154 { timeOfFlight = t; }
155
156 inline G4double GetTOF() const
157 { return timeOfFlight; }
158
159 inline void SetSide( const G4int sid )
160 { side = sid; }
161
162 inline G4int GetSide() const
163 { return side; }
164
165 inline void SetCreatorModelID( const G4int mod )
166 { theCreatorModel = mod; }
167
168 inline G4int GetCreatorModelID() const
169 { return theCreatorModel; }
170
171 inline const G4ParticleDefinition* GetParentResonanceDef() const
172 { return theParentResonanceDef; }
173
174 inline void SetParentResonanceDef( const G4ParticleDefinition* parentDef )
175 { theParentResonanceDef = parentDef; }
176
177 inline G4int GetParentResonanceID() const { return theParentResonanceID; }
178
179 inline void SetParentResonanceID ( const G4int parentID )
180 { theParentResonanceID = parentID; }
181
182 inline void SetNewlyAdded( const G4bool f )
183 { NewlyAdded = f; }
184
185 inline G4bool GetNewlyAdded() const
186 { return NewlyAdded; }
187
188 inline void SetMayBeKilled( const G4bool f )
189 { MayBeKilled = f; }
190
191 inline G4bool GetMayBeKilled() const
192 { return MayBeKilled; }
193
194 void SetZero();
195
196 void Lorentz( const G4ReactionProduct &p1, const G4ReactionProduct &p2 );
197
198 G4double Angle( const G4ReactionProduct &p ) const;
199
200 inline void SetPositionInNucleus(G4double x, G4double y, G4double z)
201 {
202 positionInNucleus.setX(x);
203 positionInNucleus.setY(y);
204 positionInNucleus.setZ(z);
205 }
206
207 inline void SetPositionInNucleus( G4ThreeVector & aPosition )
208 {
209 positionInNucleus = aPosition;
210 }
211
212 inline G4ThreeVector GetPositionInNucleus() const { return positionInNucleus; }
213 inline G4double GetXPositionInNucleus() const { return positionInNucleus.x(); }
214 inline G4double GetYPositionInNucleus() const { return positionInNucleus.y(); }
215 inline G4double GetZPositionInNucleus() const { return positionInNucleus.z(); }
216
217 inline void SetFormationTime(G4double aTime) { formationTime = aTime; }
218
219 inline G4double GetFormationTime() const { return formationTime; }
220
221 inline void HasInitialStateParton(G4bool aFlag) { hasInitialStateParton = aFlag; }
222
223 inline G4bool HasInitialStateParton() const { return hasInitialStateParton; }
224
225 private:
226
227 const G4ParticleDefinition *theParticleDefinition;
228
229 // for use with string models and cascade.
230 G4ThreeVector positionInNucleus;
231 G4double formationTime;
232 G4bool hasInitialStateParton;
233
234 // mass is included here, since pseudo-particles are created with masses different
235 // than the standard particle masses, and we are not allowed to create particles
236 G4double mass;
237
238 G4ThreeVector momentum;
239
240 G4double totalEnergy;
241 G4double kineticEnergy;
242
243 G4double timeOfFlight;
244
245 // side refers to how the particles are distributed in the
246 // forward (+) and backward (-) hemispheres in the center of mass system
247 G4int side;
248
249 G4int theCreatorModel;
250
251 const G4ParticleDefinition* theParentResonanceDef = nullptr;
252 G4int theParentResonanceID;
253
254 // NewlyAdded refers to particles added by "nuclear excitation", or as
255 // "black track" particles, or as deuterons, tritons, and alphas
256 G4bool NewlyAdded;
257 G4bool MayBeKilled;
258 };
259
260 #endif
261
262