| Geant4 Cross Reference |
1 // 1 //
2 // ******************************************* 2 // ********************************************************************
3 // * License and Disclaimer << 3 // * DISCLAIMER *
4 // * 4 // * *
5 // * The Geant4 software is copyright of th << 5 // * The following disclaimer summarizes all the specific disclaimers *
6 // * the Geant4 Collaboration. It is provided << 6 // * of contributors to this software. The specific disclaimers,which *
7 // * conditions of the Geant4 Software License << 7 // * govern, are listed with their locations in: *
8 // * LICENSE and available at http://cern.ch/ << 8 // * http://cern.ch/geant4/license *
9 // * include a list of copyright holders. <<
10 // * 9 // * *
11 // * Neither the authors of this software syst 10 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing fin 11 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warran 12 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assum 13 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file << 14 // * use. *
16 // * for the full disclaimer and the limitatio <<
17 // * 15 // * *
18 // * This code implementation is the result << 16 // * This code implementation is the intellectual property of the *
19 // * technical work of the GEANT4 collaboratio << 17 // * GEANT4 collaboration. *
20 // * By using, copying, modifying or distri << 18 // * By copying, distributing or modifying the Program (or any work *
21 // * any work based on the software) you ag << 19 // * based on the Program) you indicate your acceptance of this *
22 // * use in resulting scientific publicati << 20 // * statement, and all its terms. *
23 // * acceptance of all terms of the Geant4 Sof <<
24 // ******************************************* 21 // ********************************************************************
25 // 22 //
26 // 23 //
>> 24 // $Id: G4FastStep.cc,v 1.10 2002/11/20 17:05:14 gcosmo Exp $
>> 25 // GEANT4 tag $Name: geant4-05-01 $
27 // 26 //
28 //-------------------------------------------- 27 //---------------------------------------------------------------
29 // 28 //
30 // G4FastStep.cc 29 // G4FastStep.cc
31 // 30 //
32 // Description: 31 // Description:
33 // Encapsulates a G4ParticleChange and insu 32 // Encapsulates a G4ParticleChange and insure friendly interface
34 // methods to manage the primary/secondarie << 33 // methods to manage the primary/secondaries final state for
35 // Fast Simulation Models. 34 // Fast Simulation Models.
36 // 35 //
37 // History: 36 // History:
38 // Oct 97: Verderi && MoraDeFreitas - First 37 // Oct 97: Verderi && MoraDeFreitas - First Implementation.
39 // Apr 98: MoraDeFreitas - G4FastStep becom 38 // Apr 98: MoraDeFreitas - G4FastStep becomes the G4ParticleChange
40 // for the Fast Simulatio 39 // for the Fast Simulation Process.
41 // 40 //
42 //-------------------------------------------- 41 //---------------------------------------------------------------
43 42
44 #include "G4FastStep.hh" 43 #include "G4FastStep.hh"
45 <<
46 #include "G4DynamicParticle.hh" <<
47 #include "G4Step.hh" <<
48 #include "G4SystemOfUnits.hh" <<
49 #include "G4Track.hh" 44 #include "G4Track.hh"
>> 45 #include "G4Step.hh"
50 #include "G4TrackFastVector.hh" 46 #include "G4TrackFastVector.hh"
51 #include "G4UnitsTable.hh" << 47 #include "G4DynamicParticle.hh"
52 48
53 void G4FastStep::Initialize(const G4FastTrack& 49 void G4FastStep::Initialize(const G4FastTrack& fastTrack)
54 { 50 {
55 // keeps the fastTrack reference 51 // keeps the fastTrack reference
56 fFastTrack = &fastTrack; << 52 fFastTrack=&fastTrack;
57 53
58 // currentTrack will be used to Initialize t 54 // currentTrack will be used to Initialize the other data members
59 const G4Track& currentTrack = *(fFastTrack-> 55 const G4Track& currentTrack = *(fFastTrack->GetPrimaryTrack());
60 56
61 // use base class's method at first 57 // use base class's method at first
62 G4VParticleChange::Initialize(currentTrack); 58 G4VParticleChange::Initialize(currentTrack);
63 59
64 // set Energy/Momentum etc. equal to those o 60 // set Energy/Momentum etc. equal to those of the parent particle
65 const G4DynamicParticle* pParticle = current << 61 const G4DynamicParticle* pParticle = currentTrack.GetDynamicParticle();
66 theEnergyChange = pParticle->GetKineticEnerg << 62 theEnergyChange = pParticle->GetKineticEnergy();
67 theMomentumChange = pParticle->GetMomentumDi << 63 theMomentumChange = pParticle->GetMomentumDirection();
68 thePolarizationChange = pParticle->GetPolari << 64 thePolarizationChange = pParticle->GetPolarization();
69 theProperTimeChange = pParticle->GetProperTi << 65 theProperTimeChange = pParticle->GetProperTime();
70 66
71 // set Position/Time etc. equal to those of 67 // set Position/Time etc. equal to those of the parent track
72 thePositionChange = currentTrack.GetPosition << 68 thePositionChange = currentTrack.GetPosition();
73 theTimeChange = currentTrack.GetGlobalTime() << 69 theTimeChange = currentTrack.GetGlobalTime();
74 70
75 // switch off stepping hit invokation by def 71 // switch off stepping hit invokation by default:
76 theSteppingControlFlag = AvoidHitInvocation; 72 theSteppingControlFlag = AvoidHitInvocation;
77 73
78 // event biasing weigth: 74 // event biasing weigth:
79 theWeightChange = currentTrack.GetWeight(); << 75 theWeightChange = currentTrack.GetWeight();
80 } << 76 }
81 77
82 //---------------------------------------- 78 //----------------------------------------
83 // -- Set the StopAndKilled signal 79 // -- Set the StopAndKilled signal
84 // -- and put kinetic energy to 0.0. in the 80 // -- and put kinetic energy to 0.0. in the
85 // -- G4ParticleChange. 81 // -- G4ParticleChange.
86 //---------------------------------------- 82 //----------------------------------------
87 void G4FastStep::KillPrimaryTrack() 83 void G4FastStep::KillPrimaryTrack()
88 { 84 {
89 ProposePrimaryTrackFinalKineticEnergy(0.); << 85 SetPrimaryTrackFinalKineticEnergy(0.) ;
90 ProposeTrackStatus(fStopAndKill); << 86 SetStatusChange(fStopAndKill) ;
91 } 87 }
92 88
93 //-------------------- 89 //--------------------
94 // 90 //
95 //-------------------- 91 //--------------------
96 void G4FastStep::ProposePrimaryTrackFinalPosit << 92 void
97 << 93 G4FastStep::
>> 94 SetPrimaryTrackFinalPosition(const G4ThreeVector &position,
>> 95 G4bool localCoordinates)
98 { 96 {
99 // Compute the position coordinate in global 97 // Compute the position coordinate in global
100 // reference system if needed ... 98 // reference system if needed ...
101 G4ThreeVector globalPosition = position; 99 G4ThreeVector globalPosition = position;
102 if (localCoordinates) << 100 if (localCoordinates)
103 globalPosition = fFastTrack->GetInverseAff << 101 globalPosition = fFastTrack->GetInverseAffineTransformation()->
>> 102 TransformPoint(position);
104 // ...and feed the globalPosition: 103 // ...and feed the globalPosition:
105 thePositionChange = globalPosition; 104 thePositionChange = globalPosition;
106 } 105 }
107 106
108 void G4FastStep::SetPrimaryTrackFinalPosition( <<
109 <<
110 { <<
111 ProposePrimaryTrackFinalPosition(position, l <<
112 } <<
113 <<
114 //-------------------- 107 //--------------------
115 // 108 //
116 //-------------------- 109 //--------------------
117 void G4FastStep::ProposePrimaryTrackFinalMomen << 110 void
118 << 111 G4FastStep::
>> 112 SetPrimaryTrackFinalMomentum(const G4ThreeVector &momentum,
>> 113 G4bool localCoordinates)
119 { 114 {
120 // Compute the momentum in global reference 115 // Compute the momentum in global reference
121 // system if needed ... 116 // system if needed ...
122 G4ThreeVector globalMomentum = momentum; 117 G4ThreeVector globalMomentum = momentum;
123 if (localCoordinates) 118 if (localCoordinates)
124 globalMomentum = fFastTrack->GetInverseAff << 119 globalMomentum = fFastTrack->GetInverseAffineTransformation()->
>> 120 TransformAxis(momentum);
125 // ...and feed the globalMomentum (ensuring 121 // ...and feed the globalMomentum (ensuring unitarity)
126 SetMomentumChange(globalMomentum.unit()); 122 SetMomentumChange(globalMomentum.unit());
127 } 123 }
128 124
129 void G4FastStep::SetPrimaryTrackFinalMomentum( <<
130 <<
131 { <<
132 ProposePrimaryTrackFinalMomentumDirection(mo <<
133 } <<
134 125
135 //-------------------- 126 //--------------------
136 // 127 //
137 //-------------------- 128 //--------------------
138 void G4FastStep::ProposePrimaryTrackFinalKinet << 129 void
139 << 130 G4FastStep::
140 << 131 SetPrimaryTrackFinalKineticEnergyAndDirection(G4double kineticEnergy,
>> 132 const G4ThreeVector &direction,
>> 133 G4bool localCoordinates)
141 { 134 {
142 // Compute global direction if needed... 135 // Compute global direction if needed...
143 G4ThreeVector globalDirection = direction; 136 G4ThreeVector globalDirection = direction;
144 if (localCoordinates) 137 if (localCoordinates)
145 globalDirection = fFastTrack->GetInverseAf << 138 globalDirection =fFastTrack->GetInverseAffineTransformation()->
>> 139 TransformAxis(direction);
146 // ...and feed the globalMomentum (ensuring 140 // ...and feed the globalMomentum (ensuring unitarity)
147 SetMomentumChange(globalDirection.unit()); 141 SetMomentumChange(globalDirection.unit());
148 ProposePrimaryTrackFinalKineticEnergy(kineti << 142 SetPrimaryTrackFinalKineticEnergy(kineticEnergy);
149 } <<
150 <<
151 void G4FastStep::SetPrimaryTrackFinalKineticEn <<
152 <<
153 <<
154 { <<
155 ProposePrimaryTrackFinalKineticEnergyAndDire <<
156 } 143 }
157 144
158 //-------------------- 145 //--------------------
159 // 146 //
160 //-------------------- 147 //--------------------
161 void G4FastStep::ProposePrimaryTrackFinalPolar << 148 void
162 << 149 G4FastStep::
>> 150 SetPrimaryTrackFinalPolarization(const G4ThreeVector &polarization,
>> 151 G4bool localCoordinates)
163 { 152 {
164 // Compute polarization in global system if 153 // Compute polarization in global system if needed:
165 G4ThreeVector globalPolarization(polarizatio 154 G4ThreeVector globalPolarization(polarization);
166 if (localCoordinates) 155 if (localCoordinates)
167 globalPolarization = << 156 globalPolarization = fFastTrack->GetInverseAffineTransformation()->
168 fFastTrack->GetInverseAffineTransformati << 157 TransformAxis(globalPolarization);
169 // Feed the particle globalPolarization: 158 // Feed the particle globalPolarization:
170 thePolarizationChange = globalPolarization; 159 thePolarizationChange = globalPolarization;
171 } 160 }
172 161
173 void G4FastStep::SetPrimaryTrackFinalPolarizat <<
174 <<
175 { <<
176 ProposePrimaryTrackFinalPolarization(polariz <<
177 } <<
178 <<
179 //-------------------- 162 //--------------------
180 // 163 //
181 //-------------------- 164 //--------------------
182 G4Track* G4FastStep::CreateSecondaryTrack(cons << 165 G4Track* G4FastStep::
183 G4Th << 166 CreateSecondaryTrack(const G4DynamicParticle& dynamics,
184 G4do << 167 G4ThreeVector polarization,
>> 168 G4ThreeVector position,
>> 169 G4double time,
>> 170 G4bool localCoordinates )
185 { 171 {
186 G4DynamicParticle dummyDynamics(dynamics); 172 G4DynamicParticle dummyDynamics(dynamics);
187 << 173
188 // ----------------------------------------- 174 // ------------------------------------------
189 // Add the polarization to the dummyDynamics 175 // Add the polarization to the dummyDynamics:
190 // ----------------------------------------- 176 // ------------------------------------------
191 dummyDynamics.SetPolarization(polarization.x << 177 dummyDynamics.SetPolarization(polarization.x(),
192 << 178 polarization.y(),
>> 179 polarization.z());
>> 180
193 return CreateSecondaryTrack(dummyDynamics, p 181 return CreateSecondaryTrack(dummyDynamics, position, time, localCoordinates);
194 } 182 }
195 183
196 //-------------------- 184 //--------------------
197 // 185 //
198 //-------------------- 186 //--------------------
199 G4Track* G4FastStep::CreateSecondaryTrack(cons << 187 G4Track* G4FastStep::
200 G4do << 188 CreateSecondaryTrack(const G4DynamicParticle& dynamics,
>> 189 G4ThreeVector position,
>> 190 G4double time,
>> 191 G4bool localCoordinates )
201 { 192 {
202 // ---------------------------------------- 193 // ----------------------------------------
203 // Quantities in global coordinates system. 194 // Quantities in global coordinates system.
204 // << 195 //
205 // The allocated globalDynamics is deleted 196 // The allocated globalDynamics is deleted
206 // by the destructor of the G4Track. 197 // by the destructor of the G4Track.
207 // ---------------------------------------- 198 // ----------------------------------------
208 auto globalDynamics = new G4DynamicParticle( << 199 G4DynamicParticle* globalDynamics =
>> 200 new G4DynamicParticle(dynamics);
209 G4ThreeVector globalPosition(position); 201 G4ThreeVector globalPosition(position);
210 << 202
211 // ----------------------------------- 203 // -----------------------------------
212 // Convert to global system if needed: 204 // Convert to global system if needed:
213 // ----------------------------------- 205 // -----------------------------------
214 if (localCoordinates) { << 206 if (localCoordinates)
215 // -- Momentum Direction: << 207 {
216 globalDynamics->SetMomentumDirection( << 208 // -- Momentum Direction:
217 fFastTrack->GetInverseAffineTransformati << 209 globalDynamics->SetMomentumDirection(fFastTrack->
218 globalDynamics->GetMomentumDirection() << 210 GetInverseAffineTransformation()->
219 // -- Polarization: << 211 TransformAxis(globalDynamics->
220 G4ThreeVector globalPolarization; << 212 GetMomentumDirection()));
221 globalPolarization = fFastTrack->GetInvers << 213 // -- Polarization:
222 globalDynamics->GetPolarization()); << 214 G4ThreeVector globalPolarization;
223 globalDynamics->SetPolarization(globalPola << 215 globalPolarization = fFastTrack->GetInverseAffineTransformation()->
224 globalPola << 216 TransformAxis(globalDynamics->GetPolarization());
225 << 217 globalDynamics->SetPolarization(
226 // -- Position: << 218 globalPolarization.x(),
227 globalPosition = fFastTrack->GetInverseAff << 219 globalPolarization.y(),
228 } << 220 globalPolarization.z()
229 << 221 );
>> 222
>> 223 // -- Position:
>> 224 globalPosition = fFastTrack->GetInverseAffineTransformation()->
>> 225 TransformPoint(globalPosition);
>> 226 }
>> 227
230 //------------------------------------- 228 //-------------------------------------
231 // Create the G4Track of the secondary: 229 // Create the G4Track of the secondary:
232 //------------------------------------- 230 //-------------------------------------
233 auto secondary = new G4Track(globalDynamics, << 231 G4Track* secondary = new G4Track(
234 << 232 globalDynamics,
>> 233 time,
>> 234 globalPosition
>> 235 );
>> 236
235 //------------------------------- 237 //-------------------------------
236 // and feed the changes: 238 // and feed the changes:
237 //------------------------------- 239 //-------------------------------
238 AddSecondary(secondary); 240 AddSecondary(secondary);
239 << 241
240 //-------------------------------------- 242 //--------------------------------------
241 // returns the pointer on the secondary: 243 // returns the pointer on the secondary:
242 //-------------------------------------- 244 //--------------------------------------
243 return secondary; 245 return secondary;
244 } 246 }
245 247
246 // G4FastStep should never be Initialized in t 248 // G4FastStep should never be Initialized in this way
247 // but we must define it to avoid warnings. 249 // but we must define it to avoid warnings.
248 void G4FastStep::Initialize(const G4Track&) << 250 void G4FastStep::Initialize(const G4Track&) {
>> 251 G4Exception("G4FastStep::Initialize(const G4Track&) should never be called,\nyou must use instead the G4FastStep::Initialize(const G4FastTrack&)\nmethod!");
>> 252 }
>> 253
>> 254 G4FastStep::G4FastStep():G4VParticleChange()
>> 255 {
>> 256 if (verboseLevel>2) {
>> 257 G4cerr << "G4FastStep::G4FastStep() " << G4endl;
>> 258 }
>> 259 }
>> 260
>> 261 G4FastStep::~G4FastStep()
>> 262 {
>> 263 if (verboseLevel>2) {
>> 264 G4cerr << "G4FastStep::~G4FastStep() " << G4endl;
>> 265 }
>> 266 }
>> 267
>> 268 // copy and assignment operators are implemented as "shallow copy"
>> 269 G4FastStep::G4FastStep(const G4FastStep &right)
>> 270 {
>> 271 *this = right;
>> 272 }
>> 273
>> 274
>> 275 G4FastStep & G4FastStep::operator=(const G4FastStep &right)
>> 276 {
>> 277 if (this != &right)
>> 278 {
>> 279 G4VParticleChange::operator=(right);
>> 280 theListOfSecondaries = right.theListOfSecondaries;
>> 281 theSizeOftheListOfSecondaries = right.theSizeOftheListOfSecondaries;
>> 282 theNumberOfSecondaries = right.theNumberOfSecondaries;
>> 283 theStatusChange = right.theStatusChange;
>> 284 theMomentumChange = right.theMomentumChange;
>> 285 thePolarizationChange = right.thePolarizationChange;
>> 286 thePositionChange = right.thePositionChange;
>> 287 theTimeChange = right.theTimeChange;
>> 288 theEnergyChange = right.theEnergyChange;
>> 289 theTrueStepLength = right.theTrueStepLength;
>> 290 theLocalEnergyDeposit = right.theLocalEnergyDeposit;
>> 291 theSteppingControlFlag = right.theSteppingControlFlag;
>> 292 theWeightChange = right.theWeightChange;
>> 293 }
>> 294 return *this;
>> 295 }
>> 296
>> 297
>> 298
>> 299
>> 300
>> 301 G4bool G4FastStep::operator==(const G4FastStep &right) const
249 { 302 {
250 G4ExceptionDescription tellWhatIsWrong; << 303 return ((G4VParticleChange *)this == (G4VParticleChange *) &right);
251 tellWhatIsWrong << "G4FastStep can be initia << 304 }
252 G4Exception("G4FastStep::Initialize(const G4 << 305
253 tellWhatIsWrong); << 306 G4bool G4FastStep::operator!=(const G4FastStep &right) const
>> 307 {
>> 308 return ((G4VParticleChange *)this != (G4VParticleChange *) &right);
254 } 309 }
255 310
256 //-------------------------------------------- 311 //----------------------------------------------------------------
257 // methods for updating G4Step << 312 // methods for updating G4Step
258 // 313 //
259 314
260 G4Step* G4FastStep::UpdateStepForPostStep(G4St 315 G4Step* G4FastStep::UpdateStepForPostStep(G4Step* pStep)
261 { << 316 {
262 // A physics process always calculates the f 317 // A physics process always calculates the final state of the particle
263 318
264 // Take note that the return type of GetMome 319 // Take note that the return type of GetMomentumChange is a
265 // pointer to G4ParticleMometum. Also it is << 320 // pointer to G4ParticleMometum. Also it is a normalized
266 // momentum vector. 321 // momentum vector.
267 322
268 // G4StepPoint* pPreStepPoint = pStep->Get << 323 // G4StepPoint* pPreStepPoint = pStep->GetPreStepPoint();
269 G4StepPoint* pPostStepPoint = pStep->GetPost << 324 G4StepPoint* pPostStepPoint = pStep->GetPostStepPoint();
270 G4Track* aTrack = pStep->GetTrack(); << 325 G4Track* aTrack = pStep->GetTrack();
271 // G4double mass = aTrack->GetDynamicPa 326 // G4double mass = aTrack->GetDynamicParticle()->GetMass();
272 << 327
273 // update kinetic energy and momentum direct 328 // update kinetic energy and momentum direction
274 pPostStepPoint->SetMomentumDirection(theMome 329 pPostStepPoint->SetMomentumDirection(theMomentumChange);
275 pPostStepPoint->SetKineticEnergy(theEnergyCh << 330 pPostStepPoint->SetKineticEnergy( theEnergyChange );
276 <<
277 // update polarization <<
278 pPostStepPoint->SetPolarization(thePolarizat <<
279 331
>> 332 // update polarization
>> 333 pPostStepPoint->SetPolarization( thePolarizationChange );
>> 334
280 // update position and time 335 // update position and time
281 pPostStepPoint->SetPosition(thePositionChang << 336 pPostStepPoint->SetPosition( thePositionChange );
282 pPostStepPoint->SetGlobalTime(theTimeChange) << 337 pPostStepPoint->SetGlobalTime( theTimeChange );
283 pPostStepPoint->AddLocalTime(theTimeChange - << 338 pPostStepPoint->AddLocalTime( theTimeChange
284 pPostStepPoint->SetProperTime(theProperTimeC << 339 - aTrack->GetGlobalTime());
>> 340 pPostStepPoint->SetProperTime( theProperTimeChange );
285 341
286 // update weight 342 // update weight
287 pPostStepPoint->SetWeight(theWeightChange); << 343 pPostStepPoint->SetWeight( theWeightChange );
288 344
289 if (debugFlag) CheckIt(*aTrack); 345 if (debugFlag) CheckIt(*aTrack);
290 346
291 // Update the G4Step specific attributes << 347
>> 348 // Update the G4Step specific attributes
292 return UpdateStepInfo(pStep); 349 return UpdateStepInfo(pStep);
293 } 350 }
294 351
295 G4Step* G4FastStep::UpdateStepForAtRest(G4Step 352 G4Step* G4FastStep::UpdateStepForAtRest(G4Step* pStep)
296 { << 353 {
297 // A physics process always calculates the f 354 // A physics process always calculates the final state of the particle
298 355
299 // G4StepPoint* pPreStepPoint = pStep->GetP << 356 // G4StepPoint* pPreStepPoint = pStep->GetPreStepPoint();
300 G4StepPoint* pPostStepPoint = pStep->GetPost << 357 G4StepPoint* pPostStepPoint = pStep->GetPostStepPoint();
301 G4Track* aTrack = pStep->GetTrack(); << 358 G4Track* aTrack = pStep->GetTrack();
302 // G4double mass = aTrack->GetDynamicPar 359 // G4double mass = aTrack->GetDynamicParticle()->GetMass();
303 << 360
304 // update kinetic energy and momentum direct 361 // update kinetic energy and momentum direction
305 pPostStepPoint->SetMomentumDirection(theMome 362 pPostStepPoint->SetMomentumDirection(theMomentumChange);
306 pPostStepPoint->SetKineticEnergy(theEnergyCh << 363 pPostStepPoint->SetKineticEnergy( theEnergyChange );
307 364
308 // update polarization 365 // update polarization
309 pPostStepPoint->SetPolarization(thePolarizat << 366 pPostStepPoint->SetPolarization( thePolarizationChange );
310 << 367
311 // update position and time 368 // update position and time
312 pPostStepPoint->SetPosition(thePositionChang << 369 pPostStepPoint->SetPosition( thePositionChange );
313 pPostStepPoint->SetGlobalTime(theTimeChange) << 370 pPostStepPoint->SetGlobalTime( theTimeChange );
314 pPostStepPoint->AddLocalTime(theTimeChange - << 371 pPostStepPoint->AddLocalTime( theTimeChange
315 pPostStepPoint->SetProperTime(theProperTimeC << 372 - aTrack->GetGlobalTime());
>> 373 pPostStepPoint->SetProperTime( theProperTimeChange );
316 374
317 // update weight 375 // update weight
318 pPostStepPoint->SetWeight(theWeightChange); << 376 pPostStepPoint->SetWeight( theWeightChange );
319 377
320 if (debugFlag) CheckIt(*aTrack); 378 if (debugFlag) CheckIt(*aTrack);
321 379
322 // Update the G4Step specific attributes << 380 // Update the G4Step specific attributes
323 return UpdateStepInfo(pStep); 381 return UpdateStepInfo(pStep);
324 } 382 }
325 383
326 //-------------------------------------------- 384 //----------------------------------------------------------------
327 // methods for printing messages << 385 // methods for printing messages
328 // 386 //
329 387
330 void G4FastStep::DumpInfo() const 388 void G4FastStep::DumpInfo() const
331 { 389 {
332 // use base-class DumpInfo << 390 // use base-class DumpInfo
333 G4VParticleChange::DumpInfo(); 391 G4VParticleChange::DumpInfo();
334 392
335 G4cout << " Position - x (mm) : " < <<
336 << G4endl; <<
337 G4cout << " Position - y (mm) : " < <<
338 << G4endl; <<
339 G4cout << " Position - z (mm) : " < <<
340 << G4endl; <<
341 G4cout << " Time (ns) : " < <<
342 G4cout << " Proper Time (ns) : " < <<
343 G4long olprc = G4cout.precision(3); <<
344 G4cout << " Momentum Direct - x : " < <<
345 G4cout << " Momentum Direct - y : " < <<
346 G4cout << " Momentum Direct - z : " < <<
347 G4cout.precision(olprc); <<
348 G4cout << " Kinetic Energy (MeV): " < <<
349 G4cout.precision(3); 393 G4cout.precision(3);
350 G4cout << " Polarization - x : " < << 394 G4cout << " Position - x (mm) : "
351 << G4endl; << 395 << G4std::setw(20) << thePositionChange.x()/mm
352 G4cout << " Polarization - y : " < << 396 << G4endl;
353 << G4endl; << 397 G4cout << " Position - y (mm) : "
354 G4cout << " Polarization - z : " < << 398 << G4std::setw(20) << thePositionChange.y()/mm
355 << G4endl; << 399 << G4endl;
356 G4cout.precision(olprc); << 400 G4cout << " Position - z (mm) : "
>> 401 << G4std::setw(20) << thePositionChange.z()/mm
>> 402 << G4endl;
>> 403 G4cout << " Time (ns) : "
>> 404 << G4std::setw(20) << theTimeChange/ns
>> 405 << G4endl;
>> 406 G4cout << " Proper Time (ns) : "
>> 407 << G4std::setw(20) << theProperTimeChange/ns
>> 408 << G4endl;
>> 409 G4cout << " Momentum Direct - x : "
>> 410 << G4std::setw(20) << theMomentumChange.x()
>> 411 << G4endl;
>> 412 G4cout << " Momentum Direct - y : "
>> 413 << G4std::setw(20) << theMomentumChange.y()
>> 414 << G4endl;
>> 415 G4cout << " Momentum Direct - z : "
>> 416 << G4std::setw(20) << theMomentumChange.z()
>> 417 << G4endl;
>> 418 G4cout << " Kinetic Energy (MeV): "
>> 419 << G4std::setw(20) << theEnergyChange/MeV
>> 420 << G4endl;
>> 421 G4cout << " Polarization - x : "
>> 422 << G4std::setw(20) << thePolarizationChange.x()
>> 423 << G4endl;
>> 424 G4cout << " Polarization - y : "
>> 425 << G4std::setw(20) << thePolarizationChange.y()
>> 426 << G4endl;
>> 427 G4cout << " Polarization - z : "
>> 428 << G4std::setw(20) << thePolarizationChange.z()
>> 429 << G4endl;
357 } 430 }
358 431
359 G4bool G4FastStep::CheckIt(const G4Track& aTra 432 G4bool G4FastStep::CheckIt(const G4Track& aTrack)
360 { 433 {
361 // 434 //
362 // In the G4FastStep::CheckIt 435 // In the G4FastStep::CheckIt
363 // We only check a bit 436 // We only check a bit
364 // << 437 //
365 // If the user violates the energy, 438 // If the user violates the energy,
366 // We don't care, we agree. 439 // We don't care, we agree.
367 // 440 //
368 // But for theMomentumDirectionChange, 441 // But for theMomentumDirectionChange,
369 // We do pay attention. 442 // We do pay attention.
370 // And if too large is its range, 443 // And if too large is its range,
371 // We issue an Exception. 444 // We issue an Exception.
372 // 445 //
373 // 446 //
374 // It means, the G4FastStep::CheckIt issues 447 // It means, the G4FastStep::CheckIt issues an exception only for the
375 // theMomentumDirectionChange which should b 448 // theMomentumDirectionChange which should be an unit vector
376 // and it corrects it because it could cause 449 // and it corrects it because it could cause problems for the ulterior
377 // tracking.For the rest, only warning are i 450 // tracking.For the rest, only warning are issued.
378 451
379 G4bool itsOK = true; << 452 G4bool itsOK = true;
380 G4bool exitWithError = false; << 453 G4bool exitWithError = false;
381 G4double accuracy; << 454 G4double accuracy;
382 << 455
383 // Energy should not be larger than the init 456 // Energy should not be larger than the initial value
384 accuracy = (theEnergyChange - aTrack.GetKine << 457 accuracy = ( theEnergyChange - aTrack.GetKineticEnergy())/MeV;
385 if (accuracy > GetAccuracyForWarning()) { 458 if (accuracy > GetAccuracyForWarning()) {
386 G4ExceptionDescription ed; << 459 G4cout << " G4FastStep::CheckIt : ";
387 ed << "The energy becomes larger than the << 460 G4cout << "the energy becomes larger than the initial value !!" << G4endl;
388 << G4endl; << 461 G4cout << " Difference: " << accuracy << "[MeV] " <<G4endl;
389 G4Exception("G4FastStep::CheckIt(const G4T <<
390 itsOK = false; 462 itsOK = false;
391 if (accuracy > GetAccuracyForException()) << 463 if (accuracy > GetAccuracyForException()) exitWithError = true;
392 exitWithError = true; <<
393 } <<
394 } 464 }
395 465
396 G4bool itsOKforMomentum = true; 466 G4bool itsOKforMomentum = true;
397 if (theEnergyChange > 0.) { << 467 if ( theEnergyChange >0.) {
398 accuracy = std::abs(theMomentumChange.mag2 << 468 accuracy = abs(theMomentumChange.mag2()-1.0);
399 if (accuracy > GetAccuracyForWarning()) { 469 if (accuracy > GetAccuracyForWarning()) {
400 G4ExceptionDescription ed; << 470 G4cout << " G4FastStep::CheckIt : ";
401 ed << "The Momentum Change is not a unit << 471 G4cout << "the Momentum Change is not unit vector !!"
402 G4Exception("G4FastStep::CheckIt(const G << 472 << " Difference: " << accuracy << G4endl;
403 itsOK = itsOKforMomentum = false; 473 itsOK = itsOKforMomentum = false;
404 if (accuracy > GetAccuracyForException() << 474 if (accuracy > GetAccuracyForException()) exitWithError = true;
405 exitWithError = true; <<
406 } <<
407 } 475 }
408 } 476 }
409 << 477
410 accuracy = (aTrack.GetGlobalTime() - theTime << 478 accuracy = (aTrack.GetGlobalTime()- theTimeChange)/ns;
411 if (accuracy > GetAccuracyForWarning()) { 479 if (accuracy > GetAccuracyForWarning()) {
412 G4ExceptionDescription ed; << 480 G4cout << " G4FastStep::CheckIt : ";
413 ed << "The global time is getting backward << 481 G4cout << "the global time goes back !!"
414 G4Exception("G4FastStep::CheckIt(const G4T << 482 << " Difference: " << accuracy << "[ns] " <<G4endl;
415 itsOK = false; 483 itsOK = false;
416 } 484 }
417 << 485
418 accuracy = (aTrack.GetProperTime() - theProp << 486 accuracy = (aTrack.GetProperTime() - theProperTimeChange )/ns;
419 if (accuracy > GetAccuracyForWarning()) { << 487 if (accuracy) {
420 G4ExceptionDescription ed; << 488 G4cout << " G4FastStep::CheckIt : ";
421 ed << "The proper time is getting backward << 489 G4cout << "the proper time goes back !!"
422 G4Exception("G4FastStep::CheckIt(const G4T << 490 << " Difference: " << accuracy << "[ns] " <<G4endl;
423 itsOK = false; 491 itsOK = false;
424 } 492 }
425 << 493
426 if (!itsOK) { << 494 if (!itsOK) {
427 G4cout << "ERROR - G4FastStep::CheckIt() " << 495 G4cout << " G4FastStep::CheckIt " <<G4endl;
428 G4cout << " Pointer : " << this << << 496 G4cout << " pointer : " << this <<G4endl ;
429 DumpInfo(); 497 DumpInfo();
430 } 498 }
431 << 499
432 // Exit with error 500 // Exit with error
433 if (exitWithError) { << 501 if (exitWithError) G4Exception("G4ParticleChange::CheckIt");
434 G4ExceptionDescription ed; <<
435 ed << "An inaccuracy in G4FastStep is beyo <<
436 G4Exception("G4FastStep::CheckIt(const G4T <<
437 } <<
438 502
439 // correction for Momentum only. << 503 //correction for Momentum only.
440 if (!itsOKforMomentum) { 504 if (!itsOKforMomentum) {
441 G4double vmag = theMomentumChange.mag(); 505 G4double vmag = theMomentumChange.mag();
442 theMomentumChange = (1. / vmag) * theMomen << 506 theMomentumChange = (1./vmag)*theMomentumChange;
443 } 507 }
444 508
445 itsOK = (itsOK) && G4VParticleChange::CheckI << 509 itsOK = (itsOK) && G4VParticleChange::CheckIt(aTrack);
446 return itsOK; 510 return itsOK;
447 } 511 }
448 512