170{
175
178
179 fParticleIsLooping = false ;
180
181
182
183
184
185
186
187
188
189
191
192 fFirstStepInMassVolume = fNewTrack | fMassGeometryLimitedStep ;
193 fFirstStepInAnyVolume = fNewTrack | fAnyGeometryLimitedStep ;
194
195#ifdef G4DEBUG_TRANSPORT
196 G4cout <<
" CoupledTransport::AlongStep GPIL: "
197 << " 1st-step: any= " <<fFirstStepInAnyVolume << " ( geom= "
198 << fAnyGeometryLimitedStep << " ) "
199 << " mass= " << fFirstStepInMassVolume << " ( geom= "
200 << fMassGeometryLimitedStep << " ) "
201 <<
" newTrack= " << fNewTrack <<
G4endl;
202#endif
203
204
205 fNewTrack = false;
206
207
208
214
215#ifdef G4DEBUG_TRANSPORT
217 {
218 G4cout <<
"G4CoupledTransportation::AlongStepGPIL> called in volume "
220 }
221#endif
222
223
224
225
226
227
228 G4ThreeVector OriginShift = startPosition - fPreviousSftOrigin ;
230 startMassSafety = 0.0;
231 startFullSafety= 0.0;
232
233
234
235 if( MagSqShift <
sqr(fPreviousFullSafety) )
236 {
237 G4double mag_shift= std::sqrt(MagSqShift);
238 startMassSafety = std::max( (fPreviousMassSafety - mag_shift), 0.0);
239 startFullSafety = std::max( (fPreviousFullSafety - mag_shift), 0.0);
240
241
242
243
244
245
246 }
247
248
249
253
254 fMassGeometryLimitedStep = false ;
255 fAnyGeometryLimitedStep = false;
256
257
258
259
260
261
262
264 G4bool fieldExertsForce = false ;
265
266 const G4Field* ptrField=
nullptr;
267
269 G4bool eligibleEM = (particleCharge != 0.0)
270 || ( fUseMagneticMoment && (magneticMoment != 0.0) );
271 G4bool eligibleGrav = fUseGravity && (restMass != 0.0) ;
272
273 if( (fieldMgr!=nullptr) && (eligibleEM||eligibleGrav) )
274 {
275
276
278
279
280
281
282
284
285 if( ptrField != nullptr)
286 {
287 fieldExertsForce = eligibleEM
289 }
290 }
292
293 if( fieldExertsForce )
294 {
296
298 magneticMoment,
300 if( equationOfMotion )
301 {
302 equationOfMotion->SetChargeMomentumMass( chargeState,
303 momentumMagnitude,
304 restMass );
305 }
306#ifdef G4DEBUG_TRANSPORT
307 else
308 {
309 G4cerr <<
" ERROR in G4CoupledTransportation> "
310 <<
"Cannot find valid Equation of motion: " <<
G4endl;
311 <<
" Unable to pass Charge, Momentum and Mass " <<
G4endl;
312 }
313#endif
314 }
315
321 restMass,
322 particleCharge,
323 polarizationVec,
325 0.0,
328
331
332 fMassGeometryLimitedStep = false ;
333 fAnyGeometryLimitedStep = false ;
334 if( currentMinimumStep > 0 )
335 {
337
338
339
340 lengthAlongCurve = fPathFinder->
ComputeStep( aFieldTrack,
341 currentMinimumStep,
342 fNavigatorId,
343 stepNo,
344 newMassSafety,
345 limitedStep,
346 endTrackState,
347 currentVolume ) ;
348
350
351
353 {
354 fMassGeometryLimitedStep = true ;
355 }
356
358
359#ifdef G4DEBUG_TRANSPORT
360 if( fMassGeometryLimitedStep && !fAnyGeometryLimitedStep )
361 {
362 std::ostringstream message;
363 message <<
" ERROR in determining geometries limiting the step" <<
G4endl;
364 message << " Limiting: mass=" << fMassGeometryLimitedStep
365 <<
" any= " << fAnyGeometryLimitedStep <<
G4endl;
366 message << "Incompatible conditions - by which geometries was it limited ?";
367 G4Exception(
"G4CoupledTransportation::AlongStepGetPhysicalInteractionLength()",
369 }
370#endif
371
372 geometryStepLength = std::min( lengthAlongCurve, currentMinimumStep);
373
374
375
376 fMomentumChanged = true ;
377 fTransportEndMomentumDir = endTrackState.GetMomentumDir() ;
378
379
381 fPreviousMassSafety = newMassSafety ;
382 fPreviousFullSafety = newFullSafety ;
383
384
385#ifdef G4DEBUG_TRANSPORT
387 {
388 G4cout <<
"G4Transport:CompStep> "
389 << " called the pathfinder for a new step at " << startPosition
390 <<
" and obtained step = " << lengthAlongCurve <<
G4endl;
391 G4cout <<
" New safety (preStep) = " << newMassSafety
392 <<
" versus precalculated = " << startMassSafety <<
G4endl;
393 }
394#endif
395
396
397 startMassSafety = newMassSafety ;
398 startFullSafety = newFullSafety ;
399
400
401 fTransportEndPosition = endTrackState.GetPosition() ;
402 fTransportEndKineticEnergy = endTrackState.GetKineticEnergy() ;
403 }
404 else
405 {
406 geometryStepLength = lengthAlongCurve= 0.0 ;
407 fMomentumChanged = false ;
408
409
412
413 fTransportEndPosition = startPosition;
414
415 endTrackState= aFieldTrack;
416
417
418
419 if( startMassSafety == 0.0 )
420 {
421 fMassGeometryLimitedStep = true ;
422 fAnyGeometryLimitedStep = true;
423 }
424
425 }
426
427
428 if( !fieldExertsForce )
429 {
430 fParticleIsLooping = false ;
431 fMomentumChanged = false ;
432 fEndGlobalTimeComputed = false ;
433 }
434 else
435 {
437
438#ifdef G4DEBUG_TRANSPORT
440 {
441 G4cout <<
" G4CT::CS End Position = "
442 << fTransportEndPosition <<
G4endl;
443 G4cout <<
" G4CT::CS End Direction = "
444 << fTransportEndMomentumDir <<
G4endl;
445 }
446#endif
448 {
449
450
451
452 fCandidateEndGlobalTime = endTrackState.GetLabTimeOfFlight();
453 fEndGlobalTimeComputed = true;
454
455
456
457
458 }
459 else
460 {
461
462
463
464 fEndGlobalTimeComputed = false;
465
466
467
469 G4double endEnergy= fTransportEndKineticEnergy;
470
472 G4double absEdiff = std::fabs(startEnergy- endEnergy);
473 if( absEdiff > perMillion * endEnergy )
474 {
475 no_inexact_steps++;
476
477 }
478#ifdef G4VERBOSE
479 if( (
verboseLevel > 1) && ( absEdiff > perThousand * endEnergy) )
480 {
482 }
483#endif
484
485
486
488 }
489 }
490
491 fEndpointDistance = (fTransportEndPosition - startPosition).mag() ;
492 fTransportEndSpin = endTrackState.GetSpin();
493
494
495
496 safetyProposal= startFullSafety;
497
498
499
500
501 if( (startFullSafety < fEndpointDistance )
502 && ( particleCharge != 0.0 ) )
503
504 {
507
508
509
510
511
512
514
515
518
519
520 fPreviousMassSafety = endMassSafety ;
521 fPreviousFullSafety = endFullSafety;
523
524
525
526 safetyProposal = endFullSafety + fEndpointDistance;
527
528
529
530#ifdef G4DEBUG_TRANSPORT
532 G4cout <<
"***CoupledTransportation::AlongStepGPIL ** " <<
G4endl ;
533 G4cout <<
" Revised Safety at endpoint " << fTransportEndPosition
534 << " give safety values: Mass= " << endMassSafety
535 <<
" All= " << endFullSafety <<
G4endl ;
536 G4cout <<
" Adding endpoint distance " << fEndpointDistance
537 <<
" to obtain pseudo-safety= " << safetyProposal <<
G4endl ;
539 }
540 else
541 {
543 G4cout <<
"***CoupledTransportation::AlongStepGPIL ** " <<
G4endl ;
544 G4cout <<
" Quick Safety estimate at endpoint "
545 << fTransportEndPosition
546 << " gives safety endpoint value = "
547 << startFullSafety - fEndpointDistance
548 << " using start-point value " << startFullSafety
549 <<
" and endpointDistance " << fEndpointDistance <<
G4endl;
551#endif
552 }
553
554 proposedSafetyForStart= safetyProposal;
556
557 return geometryStepLength ;
558}
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
#define fPreviousSftOrigin
CLHEP::Hep3Vector G4ThreeVector
G4GLOB_DLL std::ostream G4cerr
void ReportInexactEnergy(G4double startEnergy, G4double endEnergy)
const G4ThreeVector & GetMomentumDirection() const
G4double GetCharge() const
G4ParticleDefinition * GetDefinition() const
G4double GetMagneticMoment() const
G4double GetTotalMomentum() const
G4bool DoesFieldChangeEnergy() const
virtual void ConfigureForTrack(const G4Track *)
const G4Field * GetDetectorField() const
G4bool IsGravityActive() const
G4double GetPDGMagneticMoment() const
G4double GetPDGSpin() const
G4double ComputeSafety(const G4ThreeVector &globalPoint)
unsigned int GetNumberGeometriesLimitingStep() const
G4double ObtainSafety(G4int navId, G4ThreeVector &globalCenterPoint)
G4double ComputeStep(const G4FieldTrack &pFieldTrack, G4double pCurrentProposedStepLength, G4int navigatorId, G4int stepNo, G4double &pNewSafety, ELimited &limitedStep, G4FieldTrack &EndState, G4VPhysicalVolume *currentVolume)
G4double GetCurrentSafety() const
G4FieldManager * FindAndSetFieldManager(G4VPhysicalVolume *pCurrentPhysVol)
G4FieldManager * GetCurrentFieldManager()
G4bool IsParticleLooping() const
G4EquationOfMotion * GetCurrentEquationOfMotion()
void SetCurrentSafety(G4double val, const G4ThreeVector &pos)
G4VPhysicalVolume * GetVolume() const
const G4ThreeVector & GetPosition() const
const G4ThreeVector & GetMomentumDirection() const
G4double GetKineticEnergy() const
const G4ThreeVector & GetPolarization() const
void ProposeTrueStepLength(G4double truePathLength)
const G4String & GetName() const