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