Geant4 11.3.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4BOptrForceCollision Class Reference

#include <G4BOptrForceCollision.hh>

+ Inheritance diagram for G4BOptrForceCollision:

Public Member Functions

 G4BOptrForceCollision (const G4String &particleToForce, const G4String &name="ForceCollision")
 
 G4BOptrForceCollision (const G4ParticleDefinition *particleToForce, const G4String &name="ForceCollision")
 
 ~G4BOptrForceCollision ()
 
virtual void Configure () final
 
virtual void ConfigureForWorker () final
 
virtual void StartRun () final
 
virtual void StartTracking (const G4Track *track) final
 
virtual void ExitBiasing (const G4Track *, const G4BiasingProcessInterface *) final
 
virtual void EndTracking () final
 
void OperationApplied (const G4BiasingProcessInterface *callingProcess, G4BiasingAppliedCase biasingCase, G4VBiasingOperation *operationApplied, const G4VParticleChange *particleChangeProduced) final
 
void OperationApplied (const G4BiasingProcessInterface *callingProcess, G4BiasingAppliedCase biasingCase, G4VBiasingOperation *occurenceOperationApplied, G4double weightForOccurenceInteraction, G4VBiasingOperation *finalStateOperationApplied, const G4VParticleChange *particleChangeProduced) final
 
- Public Member Functions inherited from G4VBiasingOperator
 G4VBiasingOperator (const G4String &name)
 
virtual ~G4VBiasingOperator ()=default
 
const G4StringGetName () const
 
void AttachTo (const G4LogicalVolume *)
 
G4BiasingAppliedCase GetPreviousBiasingAppliedCase () const
 
G4VBiasingOperationGetProposedOccurenceBiasingOperation (const G4Track *track, const G4BiasingProcessInterface *callingProcess)
 
G4VBiasingOperationGetProposedFinalStateBiasingOperation (const G4Track *track, const G4BiasingProcessInterface *callingProcess)
 
G4VBiasingOperationGetProposedNonPhysicsBiasingOperation (const G4Track *track, const G4BiasingProcessInterface *callingProcess)
 
void ExitingBiasing (const G4Track *track, const G4BiasingProcessInterface *callingProcess)
 
void ReportOperationApplied (const G4BiasingProcessInterface *callingProcess, G4BiasingAppliedCase biasingCase, G4VBiasingOperation *operationApplied, const G4VParticleChange *particleChangeProduced)
 
void ReportOperationApplied (const G4BiasingProcessInterface *callingProcess, G4BiasingAppliedCase biasingCase, G4VBiasingOperation *occurenceOperationApplied, G4double weightForOccurenceInteraction, G4VBiasingOperation *finalStateOperationApplied, const G4VParticleChange *particleChangeProduced)
 
const G4VBiasingOperationGetPreviousNonPhysicsAppliedOperation ()
 

Additional Inherited Members

- Static Public Member Functions inherited from G4VBiasingOperator
static const std::vector< G4VBiasingOperator * > & GetBiasingOperators ()
 
static G4VBiasingOperatorGetBiasingOperator (const G4LogicalVolume *)
 

Detailed Description

Definition at line 55 of file G4BOptrForceCollision.hh.

Constructor & Destructor Documentation

◆ G4BOptrForceCollision() [1/2]

G4BOptrForceCollision::G4BOptrForceCollision ( const G4String & particleToForce,
const G4String & name = "ForceCollision" )

Definition at line 49 of file G4BOptrForceCollision.cc.

52 : G4VBiasingOperator(name),
53 fForceCollisionModelID(G4PhysicsModelCatalog::GetModelID("model_GenBiasForceCollision"))
54{
55 fSharedForceInteractionOperation = new G4BOptnForceCommonTruncatedExp("SharedForceInteraction");
56 fCloningOperation = new G4BOptnCloning("Cloning");
57 fParticleToBias = G4ParticleTable::GetParticleTable()->FindParticle(particleName);
58
59 if ( fParticleToBias == nullptr )
60 {
62 ed << " Particle `" << particleName << "' not found !" << G4endl;
63 G4Exception(" G4BOptrForceCollision::G4BOptrForceCollision(...)",
64 "BIAS.GEN.07", JustWarning, ed);
65 }
66}
@ JustWarning
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
std::ostringstream G4ExceptionDescription
#define G4endl
Definition G4ios.hh:67
G4ParticleDefinition * FindParticle(G4int PDGEncoding)
static G4ParticleTable * GetParticleTable()
static G4int GetModelID(const G4int modelIndex)
G4VBiasingOperator(const G4String &name)

◆ G4BOptrForceCollision() [2/2]

G4BOptrForceCollision::G4BOptrForceCollision ( const G4ParticleDefinition * particleToForce,
const G4String & name = "ForceCollision" )

Definition at line 68 of file G4BOptrForceCollision.cc.

71 : G4VBiasingOperator(name),
72 fForceCollisionModelID(G4PhysicsModelCatalog::GetModelID("model_GenBiasForceCollision"))
73{
74 fSharedForceInteractionOperation = new G4BOptnForceCommonTruncatedExp("SharedForceInteraction");
75 fCloningOperation = new G4BOptnCloning("Cloning");
76 fParticleToBias = particle;
77}

◆ ~G4BOptrForceCollision()

G4BOptrForceCollision::~G4BOptrForceCollision ( )

Definition at line 79 of file G4BOptrForceCollision.cc.

80{
81 for ( auto it = fFreeFlightOperations.cbegin();
82 it != fFreeFlightOperations.cend(); ++it )
83 {
84 delete (*it).second;
85 }
86 delete fSharedForceInteractionOperation;
87 delete fCloningOperation;
88}

Member Function Documentation

◆ Configure()

void G4BOptrForceCollision::Configure ( )
finalvirtual

Reimplemented from G4VBiasingOperator.

Definition at line 90 of file G4BOptrForceCollision.cc.

91{
92 // -- build free flight operations:
94}
virtual void ConfigureForWorker() final

◆ ConfigureForWorker()

void G4BOptrForceCollision::ConfigureForWorker ( )
finalvirtual

Reimplemented from G4VBiasingOperator.

Definition at line 96 of file G4BOptrForceCollision.cc.

97{
98 // -- start by remembering processes under biasing,
99 // and create needed biasing operations:
100 if ( fSetup )
101 {
102 const G4ProcessManager* processManager = fParticleToBias->GetProcessManager();
103 const G4BiasingProcessSharedData* interfaceProcessSharedData = G4BiasingProcessInterface::GetSharedData( processManager );
104 if ( interfaceProcessSharedData ) // -- sharedData tested, as is can happen a user attaches an operator
105 // -- to a volume but without defining BiasingProcessInterface processes.
106 {
107 for ( std::size_t i = 0 ; i < (interfaceProcessSharedData->GetPhysicsBiasingProcessInterfaces()).size(); ++i )
108 {
109 const G4BiasingProcessInterface* wrapperProcess =
110 (interfaceProcessSharedData->GetPhysicsBiasingProcessInterfaces())[i];
111 const G4String& operationName = "FreeFlight-"+wrapperProcess->GetWrappedProcess()->GetProcessName();
112 fFreeFlightOperations[wrapperProcess] = new G4BOptnForceFreeFlight(operationName);
113 }
114 }
115 fSetup = false;
116 }
117}
const G4BiasingProcessSharedData * GetSharedData() const
const std::vector< const G4BiasingProcessInterface * > & GetPhysicsBiasingProcessInterfaces() const
const G4String & GetProcessName() const

Referenced by Configure().

◆ EndTracking()

void G4BOptrForceCollision::EndTracking ( )
finalvirtual

Reimplemented from G4VBiasingOperator.

Definition at line 307 of file G4BOptrForceCollision.cc.

308{
309 // -- check for consistency, operator should have cleaned the track:
310 if ( fCurrentTrackData != nullptr )
311 {
312 if ( !fCurrentTrackData->IsFreeFromBiasing() )
313 {
314 if ( (fCurrentTrack->GetTrackStatus() == fStopAndKill)
315 || (fCurrentTrack->GetTrackStatus() == fKillTrackAndSecondaries) )
316 {
318 ed << "Current track deleted while under biasing by "
319 << GetName() << ". Will result in inconsistencies.";
320 G4Exception(" G4BOptrForceCollision::EndTracking()",
321 "BIAS.GEN.18", JustWarning, ed);
322 }
323 }
324 }
325}
@ fKillTrackAndSecondaries
@ fStopAndKill
const G4String & GetName() const

◆ ExitBiasing()

virtual void G4BOptrForceCollision::ExitBiasing ( const G4Track * ,
const G4BiasingProcessInterface *  )
inlinefinalvirtual

Reimplemented from G4VBiasingOperator.

Definition at line 69 of file G4BOptrForceCollision.hh.

69{};

◆ OperationApplied() [1/2]

void G4BOptrForceCollision::OperationApplied ( const G4BiasingProcessInterface * callingProcess,
G4BiasingAppliedCase biasingCase,
G4VBiasingOperation * occurenceOperationApplied,
G4double weightForOccurenceInteraction,
G4VBiasingOperation * finalStateOperationApplied,
const G4VParticleChange * particleChangeProduced )
finalvirtual

Reimplemented from G4VBiasingOperator.

Definition at line 389 of file G4BOptrForceCollision.cc.

396{
397 if ( fCurrentTrackData->fForceCollisionState == ForceCollisionState::toBeForced )
398 {
399 if ( finalStateOperationApplied != fSharedForceInteractionOperation )
400 {
402 ed << " Internal inconsistency : please submit bug report. " << G4endl;
403 G4Exception(" G4BOptrForceCollision::OperationApplied(...)",
404 "BIAS.GEN.20.5", JustWarning, ed);
405 }
406 if ( fSharedForceInteractionOperation->GetInteractionOccured() )
407 fCurrentTrackData->Reset(); // -- off biasing for this track
408 }
409 else
410 {
412 ed << " Internal inconsistency : please submit bug report. " << G4endl;
413 G4Exception(" G4BOptrForceCollision::OperationApplied(...)",
414 "BIAS.GEN.20.6", JustWarning, ed);
415 }
416}

◆ OperationApplied() [2/2]

void G4BOptrForceCollision::OperationApplied ( const G4BiasingProcessInterface * callingProcess,
G4BiasingAppliedCase biasingCase,
G4VBiasingOperation * operationApplied,
const G4VParticleChange * particleChangeProduced )
finalvirtual

Reimplemented from G4VBiasingOperator.

Definition at line 327 of file G4BOptrForceCollision.cc.

332{
333 if ( fCurrentTrackData == nullptr )
334 {
335 if ( BAC != BAC_None )
336 {
338 ed << " Internal inconsistency : please submit bug report. " << G4endl;
339 G4Exception(" G4BOptrForceCollision::OperationApplied(...)",
340 "BIAS.GEN.20.1", JustWarning, ed);
341 }
342 return;
343 }
344
345 if ( fCurrentTrackData->fForceCollisionState == ForceCollisionState::toBeCloned )
346 {
347 fCurrentTrackData->fForceCollisionState = ForceCollisionState::toBeFreeFlight;
348 auto cloneData = new G4BOptrForceCollisionTrackData( this );
349 cloneData->fForceCollisionState = ForceCollisionState::toBeForced;
350 fCloningOperation->GetCloneTrack()->SetAuxiliaryTrackInformation(fForceCollisionModelID, cloneData);
351 }
352 else if ( fCurrentTrackData->fForceCollisionState == ForceCollisionState::toBeFreeFlight )
353 {
354 if ( fFreeFlightOperations[callingProcess]->OperationComplete() )
355 fCurrentTrackData->Reset(); // -- off biasing for this track
356 }
357 else if ( fCurrentTrackData->fForceCollisionState == ForceCollisionState::toBeForced )
358 {
359 if ( operationApplied != fSharedForceInteractionOperation )
360 {
362 ed << " Internal inconsistency : please submit bug report. " << G4endl;
363 G4Exception(" G4BOptrForceCollision::OperationApplied(...)",
364 "BIAS.GEN.20.2", JustWarning, ed);
365 }
366 if ( fSharedForceInteractionOperation->GetInteractionOccured() )
367 {
368 if ( operationApplied != fSharedForceInteractionOperation )
369 {
371 ed << " Internal inconsistency : please submit bug report. " << G4endl;
372 G4Exception(" G4BOptrForceCollision::OperationApplied(...)",
373 "BIAS.GEN.20.3", JustWarning, ed);
374 }
375 }
376 }
377 else
378 {
379 if ( fCurrentTrackData->fForceCollisionState != ForceCollisionState::free )
380 {
382 ed << " Internal inconsistency : please submit bug report. " << G4endl;
383 G4Exception(" G4BOptrForceCollision::OperationApplied(...)",
384 "BIAS.GEN.20.4", JustWarning, ed);
385 }
386 }
387}

◆ StartRun()

void G4BOptrForceCollision::StartRun ( )
finalvirtual

Reimplemented from G4VBiasingOperator.

Definition at line 119 of file G4BOptrForceCollision.cc.

120{
121}

◆ StartTracking()

void G4BOptrForceCollision::StartTracking ( const G4Track * track)
finalvirtual

Reimplemented from G4VBiasingOperator.

Definition at line 301 of file G4BOptrForceCollision.cc.

302{
303 fCurrentTrack = track;
304 fCurrentTrackData = nullptr;
305}

The documentation for this class was generated from the following files: