134{
139
142
143 fParticleIsLooping = false ;
144
145
146
147
148
149
150
151
152
153
155
156
157
163
164#ifdef G4DEBUG_TRANSPORT
165 if( fVerboseLevel > 1 ) {
166 G4cout <<
"G4CoupledTransportation::AlongStepGPIL> called in volume "
168 }
169#endif
170
171
172
173
174
175
176 G4ThreeVector OriginShift = startPosition - fPreviousSftOrigin ;
178 startMassSafety = 0.0;
179 startFullSafety= 0.0;
180
181
182
183 if( MagSqShift <
sqr(fPreviousFullSafety) ) {
184 G4double mag_shift= std::sqrt(MagSqShift);
185 startMassSafety = std::max( (fPreviousMassSafety - mag_shift), 0.0);
186 startFullSafety = std::max( (fPreviousFullSafety - mag_shift), 0.0);
187
188
189
190
191
192 }
193
194
195
199
200 fMassGeometryLimitedStep = false ;
201 fAnyGeometryLimitedStep = false;
202
203
204
205
206
207
208
209
210
212 G4bool fieldExertsForce = false ;
213
216
218 if( fieldMgr != 0 )
219 {
220
222
223
224
225
227
228 if( ptrField != 0)
229 {
231
232 if( (particleCharge != 0.0)
233 || (fUseMagneticMoment && (magneticMoment != 0.0) )
234 || (gravityOn && (restMass != 0.0))
235 )
236 {
237 fieldExertsForce = true;
238 }
239 }
240 }
242
244 momentumMagnitude,
245 restMass ) ;
246
247
251 0.0,
253 restMass,
254 0.0,
257 &spin ) ;
259
261
264
265 fMassGeometryLimitedStep = false ;
266 fAnyGeometryLimitedStep = false ;
267 if( currentMinimumStep > 0 ) {
269
270
271
272 lengthAlongCurve = fPathFinder->
ComputeStep( theFieldTrack,
273 currentMinimumStep,
274 fNavigatorId,
275 stepNo,
276 newMassSafety,
277 limitedStep,
278 endTrackState,
279 currentVolume ) ;
280
281
283
284
285
287 fMassGeometryLimitedStep = true ;
288 }
289
291
292
293 if( fMassGeometryLimitedStep && !fAnyGeometryLimitedStep ){
294 G4cerr <<
" Error in determining geometries limiting the step" <<
G4endl;
295 G4cerr <<
" Limiting: mass=" << fMassGeometryLimitedStep
296 <<
" any= " << fAnyGeometryLimitedStep <<
G4endl;
297 G4Exception(
"G4CoupledTransportation::AlongStepGetPhysicalInteractionLength()",
298 "PathFinderConfused",
300 "Incompatible conditions - was limited by a geometry?");
301 }
302
303
304
305
306
307
308
309
310 geometryStepLength = std::min( lengthAlongCurve, currentMinimumStep);
311
312
313
314 fMomentumChanged = true ;
315 fTransportEndMomentumDir = endTrackState.GetMomentumDir() ;
316
317
318 fPreviousSftOrigin = startPosition ;
319 fPreviousMassSafety = newMassSafety ;
320 fPreviousFullSafety = newFullSafety ;
321
322
323#ifdef G4DEBUG_TRANSPORT
324 if( fVerboseLevel > 1 ){
325 G4cout <<
"G4Transport:CompStep> "
326 << " called the pathfinder for a new step at " << startPosition
327 <<
" and obtained step = " << lengthAlongCurve <<
G4endl;
328 G4cout <<
" New safety (preStep) = " << newMassSafety
329 <<
" versus precalculated = " << startMassSafety <<
G4endl;
330 }
331#endif
332
333
334 startMassSafety = newMassSafety ;
335 startFullSafety = newFullSafety ;
336
337
338 fTransportEndPosition = endTrackState.GetPosition() ;
339 fTransportEndKineticEnergy = endTrackState.GetKineticEnergy() ;
340 } else {
341 geometryStepLength = lengthAlongCurve= 0.0 ;
342 fMomentumChanged = false ;
343
344
347
348 fTransportEndPosition = startPosition;
349
350
351 if( startMassSafety == 0.0 ) {
352 fMassGeometryLimitedStep = true ;
353 fAnyGeometryLimitedStep = true;
354 }
355
356 }
357
358
359 if( !fieldExertsForce )
360 {
361 fParticleIsLooping = false ;
362 fMomentumChanged = false ;
363 fEndGlobalTimeComputed = false ;
364
365 }
366 else
367 {
368
369#ifdef G4DEBUG_TRANSPORT
370 if( fVerboseLevel > 1 ){
371 G4cout <<
" G4CT::CS End Position = " << fTransportEndPosition <<
G4endl;
372 G4cout <<
" G4CT::CS End Direction = " << fTransportEndMomentumDir <<
G4endl;
373 }
374#endif
375
377 {
378
379
380
381 fCandidateEndGlobalTime = endTrackState.GetLabTimeOfFlight();
382 fEndGlobalTimeComputed = true;
383
384
385
386
387 }
388 else
389 {
390
391
392
393 fEndGlobalTimeComputed = false;
394
395
396
398 G4double endEnergy= fTransportEndKineticEnergy;
399
400 static G4int no_inexact_steps=0;
401 G4double absEdiff = std::fabs(startEnergy- endEnergy);
402 if( absEdiff > perMillion * endEnergy ) {
403 no_inexact_steps++;
404
405 }
406#ifdef G4VERBOSE
407 if( (fVerboseLevel > 1) && ( absEdiff > perThousand * endEnergy) ){
409 }
410#endif
411
412
413
414
416 }
417 }
418
419 endpointDistance = (fTransportEndPosition - startPosition).mag() ;
421
422 fTransportEndSpin = endTrackState.GetSpin();
423
424
425 safetyProposal= startFullSafety;
426
427
428
429
430 if( (startFullSafety < endpointDistance )
431 && ( particleCharge != 0.0 ) )
432
433 {
436
437
438
439
440
441
443
444
447
448
449 fPreviousMassSafety = endMassSafety ;
450 fPreviousFullSafety = endFullSafety;
451 fPreviousSftOrigin = fTransportEndPosition ;
452
453
454
455 safetyProposal = endFullSafety + endpointDistance;
456
457
458
459
460
461#ifdef G4DEBUG_TRANSPORT
462 int prec=
G4cout.precision(12) ;
463 G4cout <<
"***Transportation::AlongStepGPIL ** " <<
G4endl ;
464 G4cout <<
" Revised Safety at endpoint " << fTransportEndPosition
465 << " give safety values: Mass= " << endMassSafety
466 <<
" All= " << endFullSafety <<
G4endl ;
467 G4cout <<
" Adding endpoint distance " << endpointDistance
468 <<
" to obtain pseudo-safety= " << safetyProposal <<
G4endl ;
470 }
471 else{
472 int prec=
G4cout.precision(12) ;
473 G4cout <<
"***Transportation::AlongStepGPIL ** " <<
G4endl ;
474 G4cout <<
" Quick Safety estimate at endpoint " << fTransportEndPosition
475 << " gives safety endpoint value = " << startFullSafety - endpointDistance
476 << " using start-point value " << startFullSafety
477 <<
" and endpointDistance " << endpointDistance <<
G4endl;
479#endif
480 }
481
482 proposedSafetyForStart= safetyProposal;
484
485 return geometryStepLength ;
486}
CLHEP::Hep3Vector G4ThreeVector
G4DLLIMPORT 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
void SetChargeAndMoments(G4double charge, G4double magnetic_dipole_moment=DBL_MAX, G4double electric_dipole_moment=DBL_MAX, G4double magnetic_charge=DBL_MAX)
G4bool IsGravityActive() const
G4double GetPDGMass() 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
void SetChargeMomentumMass(G4double charge, G4double momentum, G4double pMass)
void SetCurrentSafety(G4double val, const G4ThreeVector &pos)
G4VPhysicalVolume * GetVolume() const
const G4ThreeVector & GetPosition() const
const G4ThreeVector & GetMomentumDirection() const
const G4ThreeVector & GetPolarization() const
void ProposeTrueStepLength(G4double truePathLength)
const G4String & GetName() const
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)