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: geant4-08-03-patch-01 $ 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