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.14 2004/12/07 08:22:33 gcosmo Exp $ >> 25 // GEANT4 tag $Name: geant4-07-00-cand-03 $ 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 ProposeTrackStatus(fStopAndKill) ; 91 } 87 } 92 88 93 //-------------------- 89 //-------------------- 94 // 90 // 95 //-------------------- 91 //-------------------- 96 void G4FastStep::ProposePrimaryTrackFinalPosit << 92 void 97 << 93 G4FastStep:: >> 94 ProposePrimaryTrackFinalPosition(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( << 107 void 109 << 108 G4FastStep:: >> 109 SetPrimaryTrackFinalPosition(const G4ThreeVector &position, >> 110 G4bool localCoordinates) 110 { 111 { 111 ProposePrimaryTrackFinalPosition(position, l 112 ProposePrimaryTrackFinalPosition(position, localCoordinates); 112 } 113 } 113 114 114 //-------------------- 115 //-------------------- 115 // 116 // 116 //-------------------- 117 //-------------------- 117 void G4FastStep::ProposePrimaryTrackFinalMomen << 118 void 118 << 119 G4FastStep:: >> 120 ProposePrimaryTrackFinalMomentumDirection(const G4ThreeVector &momentum, >> 121 G4bool localCoordinates) 119 { 122 { 120 // Compute the momentum in global reference 123 // Compute the momentum in global reference 121 // system if needed ... 124 // system if needed ... 122 G4ThreeVector globalMomentum = momentum; 125 G4ThreeVector globalMomentum = momentum; 123 if (localCoordinates) 126 if (localCoordinates) 124 globalMomentum = fFastTrack->GetInverseAff << 127 globalMomentum = fFastTrack->GetInverseAffineTransformation()-> >> 128 TransformAxis(momentum); 125 // ...and feed the globalMomentum (ensuring 129 // ...and feed the globalMomentum (ensuring unitarity) 126 SetMomentumChange(globalMomentum.unit()); 130 SetMomentumChange(globalMomentum.unit()); 127 } 131 } 128 132 129 void G4FastStep::SetPrimaryTrackFinalMomentum( << 133 void 130 << 134 G4FastStep:: >> 135 SetPrimaryTrackFinalMomentum(const G4ThreeVector &momentum, >> 136 G4bool localCoordinates) 131 { 137 { 132 ProposePrimaryTrackFinalMomentumDirection(mo 138 ProposePrimaryTrackFinalMomentumDirection(momentum, localCoordinates); 133 } 139 } 134 140 135 //-------------------- 141 //-------------------- 136 // 142 // 137 //-------------------- 143 //-------------------- 138 void G4FastStep::ProposePrimaryTrackFinalKinet << 144 void 139 << 145 G4FastStep:: 140 << 146 ProposePrimaryTrackFinalKineticEnergyAndDirection(G4double kineticEnergy, >> 147 const G4ThreeVector &direction, >> 148 G4bool localCoordinates) 141 { 149 { 142 // Compute global direction if needed... 150 // Compute global direction if needed... 143 G4ThreeVector globalDirection = direction; 151 G4ThreeVector globalDirection = direction; 144 if (localCoordinates) 152 if (localCoordinates) 145 globalDirection = fFastTrack->GetInverseAf << 153 globalDirection =fFastTrack->GetInverseAffineTransformation()-> >> 154 TransformAxis(direction); 146 // ...and feed the globalMomentum (ensuring 155 // ...and feed the globalMomentum (ensuring unitarity) 147 SetMomentumChange(globalDirection.unit()); 156 SetMomentumChange(globalDirection.unit()); 148 ProposePrimaryTrackFinalKineticEnergy(kineti << 157 SetPrimaryTrackFinalKineticEnergy(kineticEnergy); 149 } 158 } 150 159 151 void G4FastStep::SetPrimaryTrackFinalKineticEn << 160 void 152 << 161 G4FastStep:: 153 << 162 SetPrimaryTrackFinalKineticEnergyAndDirection(G4double kineticEnergy, >> 163 const G4ThreeVector &direction, >> 164 G4bool localCoordinates) 154 { 165 { 155 ProposePrimaryTrackFinalKineticEnergyAndDire 166 ProposePrimaryTrackFinalKineticEnergyAndDirection(kineticEnergy, direction, localCoordinates); 156 } 167 } 157 168 158 //-------------------- 169 //-------------------- 159 // 170 // 160 //-------------------- 171 //-------------------- 161 void G4FastStep::ProposePrimaryTrackFinalPolar << 172 void 162 << 173 G4FastStep:: >> 174 ProposePrimaryTrackFinalPolarization(const G4ThreeVector &polarization, >> 175 G4bool localCoordinates) 163 { 176 { 164 // Compute polarization in global system if 177 // Compute polarization in global system if needed: 165 G4ThreeVector globalPolarization(polarizatio 178 G4ThreeVector globalPolarization(polarization); 166 if (localCoordinates) 179 if (localCoordinates) 167 globalPolarization = << 180 globalPolarization = fFastTrack->GetInverseAffineTransformation()-> 168 fFastTrack->GetInverseAffineTransformati << 181 TransformAxis(globalPolarization); 169 // Feed the particle globalPolarization: 182 // Feed the particle globalPolarization: 170 thePolarizationChange = globalPolarization; 183 thePolarizationChange = globalPolarization; 171 } 184 } 172 185 173 void G4FastStep::SetPrimaryTrackFinalPolarizat << 186 void 174 << 187 G4FastStep:: >> 188 SetPrimaryTrackFinalPolarization(const G4ThreeVector &polarization, >> 189 G4bool localCoordinates) 175 { 190 { 176 ProposePrimaryTrackFinalPolarization(polariz 191 ProposePrimaryTrackFinalPolarization(polarization, localCoordinates); 177 } 192 } 178 193 179 //-------------------- 194 //-------------------- 180 // 195 // 181 //-------------------- 196 //-------------------- 182 G4Track* G4FastStep::CreateSecondaryTrack(cons << 197 G4Track* G4FastStep:: 183 G4Th << 198 CreateSecondaryTrack(const G4DynamicParticle& dynamics, 184 G4do << 199 G4ThreeVector polarization, >> 200 G4ThreeVector position, >> 201 G4double time, >> 202 G4bool localCoordinates ) 185 { 203 { 186 G4DynamicParticle dummyDynamics(dynamics); 204 G4DynamicParticle dummyDynamics(dynamics); 187 << 205 188 // ----------------------------------------- 206 // ------------------------------------------ 189 // Add the polarization to the dummyDynamics 207 // Add the polarization to the dummyDynamics: 190 // ----------------------------------------- 208 // ------------------------------------------ 191 dummyDynamics.SetPolarization(polarization.x << 209 dummyDynamics.SetPolarization(polarization.x(), 192 << 210 polarization.y(), >> 211 polarization.z()); >> 212 193 return CreateSecondaryTrack(dummyDynamics, p 213 return CreateSecondaryTrack(dummyDynamics, position, time, localCoordinates); 194 } 214 } 195 215 196 //-------------------- 216 //-------------------- 197 // 217 // 198 //-------------------- 218 //-------------------- 199 G4Track* G4FastStep::CreateSecondaryTrack(cons << 219 G4Track* G4FastStep:: 200 G4do << 220 CreateSecondaryTrack(const G4DynamicParticle& dynamics, >> 221 G4ThreeVector position, >> 222 G4double time, >> 223 G4bool localCoordinates ) 201 { 224 { 202 // ---------------------------------------- 225 // ---------------------------------------- 203 // Quantities in global coordinates system. 226 // Quantities in global coordinates system. 204 // << 227 // 205 // The allocated globalDynamics is deleted 228 // The allocated globalDynamics is deleted 206 // by the destructor of the G4Track. 229 // by the destructor of the G4Track. 207 // ---------------------------------------- 230 // ---------------------------------------- 208 auto globalDynamics = new G4DynamicParticle( << 231 G4DynamicParticle* globalDynamics = >> 232 new G4DynamicParticle(dynamics); 209 G4ThreeVector globalPosition(position); 233 G4ThreeVector globalPosition(position); 210 << 234 211 // ----------------------------------- 235 // ----------------------------------- 212 // Convert to global system if needed: 236 // Convert to global system if needed: 213 // ----------------------------------- 237 // ----------------------------------- 214 if (localCoordinates) { << 238 if (localCoordinates) 215 // -- Momentum Direction: << 239 { 216 globalDynamics->SetMomentumDirection( << 240 // -- Momentum Direction: 217 fFastTrack->GetInverseAffineTransformati << 241 globalDynamics->SetMomentumDirection(fFastTrack-> 218 globalDynamics->GetMomentumDirection() << 242 GetInverseAffineTransformation()-> 219 // -- Polarization: << 243 TransformAxis(globalDynamics-> 220 G4ThreeVector globalPolarization; << 244 GetMomentumDirection())); 221 globalPolarization = fFastTrack->GetInvers << 245 // -- Polarization: 222 globalDynamics->GetPolarization()); << 246 G4ThreeVector globalPolarization; 223 globalDynamics->SetPolarization(globalPola << 247 globalPolarization = fFastTrack->GetInverseAffineTransformation()-> 224 globalPola << 248 TransformAxis(globalDynamics->GetPolarization()); 225 << 249 globalDynamics->SetPolarization( 226 // -- Position: << 250 globalPolarization.x(), 227 globalPosition = fFastTrack->GetInverseAff << 251 globalPolarization.y(), 228 } << 252 globalPolarization.z() 229 << 253 ); >> 254 >> 255 // -- Position: >> 256 globalPosition = fFastTrack->GetInverseAffineTransformation()-> >> 257 TransformPoint(globalPosition); >> 258 } >> 259 230 //------------------------------------- 260 //------------------------------------- 231 // Create the G4Track of the secondary: 261 // Create the G4Track of the secondary: 232 //------------------------------------- 262 //------------------------------------- 233 auto secondary = new G4Track(globalDynamics, << 263 G4Track* secondary = new G4Track( 234 << 264 globalDynamics, >> 265 time, >> 266 globalPosition >> 267 ); >> 268 235 //------------------------------- 269 //------------------------------- 236 // and feed the changes: 270 // and feed the changes: 237 //------------------------------- 271 //------------------------------- 238 AddSecondary(secondary); 272 AddSecondary(secondary); 239 << 273 240 //-------------------------------------- 274 //-------------------------------------- 241 // returns the pointer on the secondary: 275 // returns the pointer on the secondary: 242 //-------------------------------------- 276 //-------------------------------------- 243 return secondary; 277 return secondary; 244 } 278 } 245 279 246 // G4FastStep should never be Initialized in t 280 // G4FastStep should never be Initialized in this way 247 // but we must define it to avoid warnings. 281 // but we must define it to avoid warnings. 248 void G4FastStep::Initialize(const G4Track&) 282 void G4FastStep::Initialize(const G4Track&) 249 { 283 { 250 G4ExceptionDescription tellWhatIsWrong; << 284 G4Exception("G4FastStep::Initialize(const G4Track&) should never be called,\nyou must use instead the G4FastStep::Initialize(const G4FastTrack&)\nmethod!"); 251 tellWhatIsWrong << "G4FastStep can be initia << 285 } 252 G4Exception("G4FastStep::Initialize(const G4 << 286 253 tellWhatIsWrong); << 287 G4FastStep::G4FastStep() >> 288 : G4VParticleChange() >> 289 { >> 290 if (verboseLevel>2) >> 291 { >> 292 G4cerr << "G4FastStep::G4FastStep() " << G4endl; >> 293 } >> 294 } >> 295 >> 296 G4FastStep::~G4FastStep() >> 297 { >> 298 if (verboseLevel>2) >> 299 { >> 300 G4cerr << "G4FastStep::~G4FastStep() " << G4endl; >> 301 } >> 302 } >> 303 >> 304 // copy and assignment operators are implemented as "shallow copy" >> 305 G4FastStep::G4FastStep(const G4FastStep &right) >> 306 : G4VParticleChange() >> 307 { >> 308 *this = right; >> 309 } >> 310 >> 311 >> 312 G4FastStep & G4FastStep::operator=(const G4FastStep &right) >> 313 { >> 314 if (this != &right) >> 315 { >> 316 G4VParticleChange::operator=(right); >> 317 theListOfSecondaries = right.theListOfSecondaries; >> 318 theSizeOftheListOfSecondaries = right.theSizeOftheListOfSecondaries; >> 319 theNumberOfSecondaries = right.theNumberOfSecondaries; >> 320 theStatusChange = right.theStatusChange; >> 321 theMomentumChange = right.theMomentumChange; >> 322 thePolarizationChange = right.thePolarizationChange; >> 323 thePositionChange = right.thePositionChange; >> 324 theTimeChange = right.theTimeChange; >> 325 theEnergyChange = right.theEnergyChange; >> 326 theTrueStepLength = right.theTrueStepLength; >> 327 theLocalEnergyDeposit = right.theLocalEnergyDeposit; >> 328 theSteppingControlFlag = right.theSteppingControlFlag; >> 329 theWeightChange = right.theWeightChange; >> 330 } >> 331 return *this; >> 332 } >> 333 >> 334 >> 335 >> 336 >> 337 >> 338 G4bool G4FastStep::operator==(const G4FastStep &right) const >> 339 { >> 340 return ((G4VParticleChange *)this == (G4VParticleChange *) &right); >> 341 } >> 342 >> 343 G4bool G4FastStep::operator!=(const G4FastStep &right) const >> 344 { >> 345 return ((G4VParticleChange *)this != (G4VParticleChange *) &right); 254 } 346 } 255 347 256 //-------------------------------------------- 348 //---------------------------------------------------------------- 257 // methods for updating G4Step << 349 // methods for updating G4Step 258 // 350 // 259 351 260 G4Step* G4FastStep::UpdateStepForPostStep(G4St 352 G4Step* G4FastStep::UpdateStepForPostStep(G4Step* pStep) 261 { << 353 { 262 // A physics process always calculates the f 354 // A physics process always calculates the final state of the particle 263 355 264 // Take note that the return type of GetMome 356 // Take note that the return type of GetMomentumChange is a 265 // pointer to G4ParticleMometum. Also it is << 357 // pointer to G4ParticleMometum. Also it is a normalized 266 // momentum vector. 358 // momentum vector. 267 359 268 // G4StepPoint* pPreStepPoint = pStep->Get << 360 // G4StepPoint* pPreStepPoint = pStep->GetPreStepPoint(); 269 G4StepPoint* pPostStepPoint = pStep->GetPost << 361 G4StepPoint* pPostStepPoint = pStep->GetPostStepPoint(); 270 G4Track* aTrack = pStep->GetTrack(); << 362 G4Track* aTrack = pStep->GetTrack(); 271 // G4double mass = aTrack->GetDynamicPa 363 // G4double mass = aTrack->GetDynamicParticle()->GetMass(); 272 << 364 273 // update kinetic energy and momentum direct 365 // update kinetic energy and momentum direction 274 pPostStepPoint->SetMomentumDirection(theMome 366 pPostStepPoint->SetMomentumDirection(theMomentumChange); 275 pPostStepPoint->SetKineticEnergy(theEnergyCh << 367 pPostStepPoint->SetKineticEnergy( theEnergyChange ); 276 << 277 // update polarization << 278 pPostStepPoint->SetPolarization(thePolarizat << 279 368 >> 369 // update polarization >> 370 pPostStepPoint->SetPolarization( thePolarizationChange ); >> 371 280 // update position and time 372 // update position and time 281 pPostStepPoint->SetPosition(thePositionChang << 373 pPostStepPoint->SetPosition( thePositionChange ); 282 pPostStepPoint->SetGlobalTime(theTimeChange) << 374 pPostStepPoint->SetGlobalTime( theTimeChange ); 283 pPostStepPoint->AddLocalTime(theTimeChange - << 375 pPostStepPoint->AddLocalTime( theTimeChange 284 pPostStepPoint->SetProperTime(theProperTimeC << 376 - aTrack->GetGlobalTime()); >> 377 pPostStepPoint->SetProperTime( theProperTimeChange ); 285 378 286 // update weight 379 // update weight 287 pPostStepPoint->SetWeight(theWeightChange); << 380 pPostStepPoint->SetWeight( theWeightChange ); 288 381 289 if (debugFlag) CheckIt(*aTrack); 382 if (debugFlag) CheckIt(*aTrack); 290 383 291 // Update the G4Step specific attributes << 384 >> 385 // Update the G4Step specific attributes 292 return UpdateStepInfo(pStep); 386 return UpdateStepInfo(pStep); 293 } 387 } 294 388 295 G4Step* G4FastStep::UpdateStepForAtRest(G4Step 389 G4Step* G4FastStep::UpdateStepForAtRest(G4Step* pStep) 296 { << 390 { 297 // A physics process always calculates the f 391 // A physics process always calculates the final state of the particle 298 392 299 // G4StepPoint* pPreStepPoint = pStep->GetP << 393 // G4StepPoint* pPreStepPoint = pStep->GetPreStepPoint(); 300 G4StepPoint* pPostStepPoint = pStep->GetPost << 394 G4StepPoint* pPostStepPoint = pStep->GetPostStepPoint(); 301 G4Track* aTrack = pStep->GetTrack(); << 395 G4Track* aTrack = pStep->GetTrack(); 302 // G4double mass = aTrack->GetDynamicPar 396 // G4double mass = aTrack->GetDynamicParticle()->GetMass(); 303 << 397 304 // update kinetic energy and momentum direct 398 // update kinetic energy and momentum direction 305 pPostStepPoint->SetMomentumDirection(theMome 399 pPostStepPoint->SetMomentumDirection(theMomentumChange); 306 pPostStepPoint->SetKineticEnergy(theEnergyCh << 400 pPostStepPoint->SetKineticEnergy( theEnergyChange ); 307 401 308 // update polarization 402 // update polarization 309 pPostStepPoint->SetPolarization(thePolarizat << 403 pPostStepPoint->SetPolarization( thePolarizationChange ); 310 << 404 311 // update position and time 405 // update position and time 312 pPostStepPoint->SetPosition(thePositionChang << 406 pPostStepPoint->SetPosition( thePositionChange ); 313 pPostStepPoint->SetGlobalTime(theTimeChange) << 407 pPostStepPoint->SetGlobalTime( theTimeChange ); 314 pPostStepPoint->AddLocalTime(theTimeChange - << 408 pPostStepPoint->AddLocalTime( theTimeChange 315 pPostStepPoint->SetProperTime(theProperTimeC << 409 - aTrack->GetGlobalTime()); >> 410 pPostStepPoint->SetProperTime( theProperTimeChange ); 316 411 317 // update weight 412 // update weight 318 pPostStepPoint->SetWeight(theWeightChange); << 413 pPostStepPoint->SetWeight( theWeightChange ); 319 414 320 if (debugFlag) CheckIt(*aTrack); 415 if (debugFlag) CheckIt(*aTrack); 321 416 322 // Update the G4Step specific attributes << 417 // Update the G4Step specific attributes 323 return UpdateStepInfo(pStep); 418 return UpdateStepInfo(pStep); 324 } 419 } 325 420 326 //-------------------------------------------- 421 //---------------------------------------------------------------- 327 // methods for printing messages << 422 // methods for printing messages 328 // 423 // 329 424 330 void G4FastStep::DumpInfo() const 425 void G4FastStep::DumpInfo() const 331 { 426 { 332 // use base-class DumpInfo << 427 // use base-class DumpInfo 333 G4VParticleChange::DumpInfo(); 428 G4VParticleChange::DumpInfo(); 334 429 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); 430 G4cout.precision(3); 350 G4cout << " Polarization - x : " < << 431 G4cout << " Position - x (mm) : " 351 << G4endl; << 432 << std::setw(20) << thePositionChange.x()/mm 352 G4cout << " Polarization - y : " < << 433 << G4endl; 353 << G4endl; << 434 G4cout << " Position - y (mm) : " 354 G4cout << " Polarization - z : " < << 435 << std::setw(20) << thePositionChange.y()/mm 355 << G4endl; << 436 << G4endl; 356 G4cout.precision(olprc); << 437 G4cout << " Position - z (mm) : " >> 438 << std::setw(20) << thePositionChange.z()/mm >> 439 << G4endl; >> 440 G4cout << " Time (ns) : " >> 441 << std::setw(20) << theTimeChange/ns >> 442 << G4endl; >> 443 G4cout << " Proper Time (ns) : " >> 444 << std::setw(20) << theProperTimeChange/ns >> 445 << G4endl; >> 446 G4cout << " Momentum Direct - x : " >> 447 << std::setw(20) << theMomentumChange.x() >> 448 << G4endl; >> 449 G4cout << " Momentum Direct - y : " >> 450 << std::setw(20) << theMomentumChange.y() >> 451 << G4endl; >> 452 G4cout << " Momentum Direct - z : " >> 453 << std::setw(20) << theMomentumChange.z() >> 454 << G4endl; >> 455 G4cout << " Kinetic Energy (MeV): " >> 456 << std::setw(20) << theEnergyChange/MeV >> 457 << G4endl; >> 458 G4cout << " Polarization - x : " >> 459 << std::setw(20) << thePolarizationChange.x() >> 460 << G4endl; >> 461 G4cout << " Polarization - y : " >> 462 << std::setw(20) << thePolarizationChange.y() >> 463 << G4endl; >> 464 G4cout << " Polarization - z : " >> 465 << std::setw(20) << thePolarizationChange.z() >> 466 << G4endl; 357 } 467 } 358 468 359 G4bool G4FastStep::CheckIt(const G4Track& aTra 469 G4bool G4FastStep::CheckIt(const G4Track& aTrack) 360 { 470 { 361 // 471 // 362 // In the G4FastStep::CheckIt 472 // In the G4FastStep::CheckIt 363 // We only check a bit 473 // We only check a bit 364 // << 474 // 365 // If the user violates the energy, 475 // If the user violates the energy, 366 // We don't care, we agree. 476 // We don't care, we agree. 367 // 477 // 368 // But for theMomentumDirectionChange, 478 // But for theMomentumDirectionChange, 369 // We do pay attention. 479 // We do pay attention. 370 // And if too large is its range, 480 // And if too large is its range, 371 // We issue an Exception. 481 // We issue an Exception. 372 // 482 // 373 // 483 // 374 // It means, the G4FastStep::CheckIt issues 484 // It means, the G4FastStep::CheckIt issues an exception only for the 375 // theMomentumDirectionChange which should b 485 // theMomentumDirectionChange which should be an unit vector 376 // and it corrects it because it could cause 486 // and it corrects it because it could cause problems for the ulterior 377 // tracking.For the rest, only warning are i 487 // tracking.For the rest, only warning are issued. 378 488 379 G4bool itsOK = true; << 489 G4bool itsOK = true; 380 G4bool exitWithError = false; << 490 G4bool exitWithError = false; 381 G4double accuracy; << 491 G4double accuracy; 382 << 492 383 // Energy should not be larger than the init 493 // Energy should not be larger than the initial value 384 accuracy = (theEnergyChange - aTrack.GetKine << 494 accuracy = ( theEnergyChange - aTrack.GetKineticEnergy())/MeV; 385 if (accuracy > GetAccuracyForWarning()) { 495 if (accuracy > GetAccuracyForWarning()) { 386 G4ExceptionDescription ed; << 496 G4cout << " G4FastStep::CheckIt : "; 387 ed << "The energy becomes larger than the << 497 G4cout << "the energy becomes larger than the initial value !!" << G4endl; 388 << G4endl; << 498 G4cout << " Difference: " << accuracy << "[MeV] " <<G4endl; 389 G4Exception("G4FastStep::CheckIt(const G4T << 390 itsOK = false; 499 itsOK = false; 391 if (accuracy > GetAccuracyForException()) << 500 if (accuracy > GetAccuracyForException()) exitWithError = true; 392 exitWithError = true; << 393 } << 394 } 501 } 395 502 396 G4bool itsOKforMomentum = true; 503 G4bool itsOKforMomentum = true; 397 if (theEnergyChange > 0.) { << 504 if ( theEnergyChange >0.) { 398 accuracy = std::abs(theMomentumChange.mag2 << 505 accuracy = std::abs(theMomentumChange.mag2()-1.0); 399 if (accuracy > GetAccuracyForWarning()) { 506 if (accuracy > GetAccuracyForWarning()) { 400 G4ExceptionDescription ed; << 507 G4cout << " G4FastStep::CheckIt : "; 401 ed << "The Momentum Change is not a unit << 508 G4cout << "the Momentum Change is not unit vector !!" 402 G4Exception("G4FastStep::CheckIt(const G << 509 << " Difference: " << accuracy << G4endl; 403 itsOK = itsOKforMomentum = false; 510 itsOK = itsOKforMomentum = false; 404 if (accuracy > GetAccuracyForException() << 511 if (accuracy > GetAccuracyForException()) exitWithError = true; 405 exitWithError = true; << 406 } << 407 } 512 } 408 } 513 } 409 << 514 410 accuracy = (aTrack.GetGlobalTime() - theTime << 515 accuracy = (aTrack.GetGlobalTime()- theTimeChange)/ns; 411 if (accuracy > GetAccuracyForWarning()) { 516 if (accuracy > GetAccuracyForWarning()) { 412 G4ExceptionDescription ed; << 517 G4cout << " G4FastStep::CheckIt : "; 413 ed << "The global time is getting backward << 518 G4cout << "the global time goes back !!" 414 G4Exception("G4FastStep::CheckIt(const G4T << 519 << " Difference: " << accuracy << "[ns] " <<G4endl; 415 itsOK = false; 520 itsOK = false; 416 } 521 } 417 << 522 418 accuracy = (aTrack.GetProperTime() - theProp << 523 accuracy = (aTrack.GetProperTime() - theProperTimeChange )/ns; 419 if (accuracy > GetAccuracyForWarning()) { << 524 if (accuracy) { 420 G4ExceptionDescription ed; << 525 G4cout << " G4FastStep::CheckIt : "; 421 ed << "The proper time is getting backward << 526 G4cout << "the proper time goes back !!" 422 G4Exception("G4FastStep::CheckIt(const G4T << 527 << " Difference: " << accuracy << "[ns] " <<G4endl; 423 itsOK = false; 528 itsOK = false; 424 } 529 } 425 << 530 426 if (!itsOK) { << 531 if (!itsOK) { 427 G4cout << "ERROR - G4FastStep::CheckIt() " << 532 G4cout << " G4FastStep::CheckIt " <<G4endl; 428 G4cout << " Pointer : " << this << << 533 G4cout << " pointer : " << this <<G4endl ; 429 DumpInfo(); 534 DumpInfo(); 430 } 535 } 431 << 536 432 // Exit with error 537 // Exit with error 433 if (exitWithError) { << 538 if (exitWithError) G4Exception("G4ParticleChange::CheckIt"); 434 G4ExceptionDescription ed; << 435 ed << "An inaccuracy in G4FastStep is beyo << 436 G4Exception("G4FastStep::CheckIt(const G4T << 437 } << 438 539 439 // correction for Momentum only. << 540 //correction for Momentum only. 440 if (!itsOKforMomentum) { 541 if (!itsOKforMomentum) { 441 G4double vmag = theMomentumChange.mag(); 542 G4double vmag = theMomentumChange.mag(); 442 theMomentumChange = (1. / vmag) * theMomen << 543 theMomentumChange = (1./vmag)*theMomentumChange; 443 } 544 } 444 545 445 itsOK = (itsOK) && G4VParticleChange::CheckI << 546 itsOK = (itsOK) && G4VParticleChange::CheckIt(aTrack); 446 return itsOK; 547 return itsOK; 447 } 548 } 448 549