Geant4.10
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Public Member Functions
G4BOptnForceCommonTruncatedExp Class Reference

#include <G4BOptnForceCommonTruncatedExp.hh>

Inheritance diagram for G4BOptnForceCommonTruncatedExp:
G4VBiasingOperation

Public Member Functions

 G4BOptnForceCommonTruncatedExp (G4String name)
 
virtual ~G4BOptnForceCommonTruncatedExp ()
 
virtual const
G4VBiasingInteractionLaw
ProvideOccurenceBiasingInteractionLaw (const G4BiasingProcessInterface *)
 
virtual G4ForceCondition ProposeForceCondition (const G4ForceCondition processCondition)
 
virtual G4bool DenyProcessPostStepDoIt (const G4BiasingProcessInterface *, const G4Track *, const G4Step *, G4double &)
 
virtual G4double ProposeAlongStepLimit (const G4BiasingProcessInterface *)
 
virtual G4GPILSelection ProposeGPILSelection (const G4GPILSelection processSelection)
 
virtual G4VParticleChangeApplyFinalStateBiasing (const G4BiasingProcessInterface *, const G4Track *, const G4Step *)
 
virtual G4double DistanceToApplyOperation (const G4Track *, G4double, G4ForceCondition *)
 
virtual G4VParticleChangeGenerateBiasingFinalState (const G4Track *, const G4Step *)
 
G4ILawCommonTruncatedExpGetCommonTruncatedExpLaw ()
 
void Initialize (const G4Track *)
 
void UpdateForStep (const G4Step *)
 
void Sample ()
 
const G4ThreeVectorGetInitialMomentum () const
 
G4double GetMaximumDistance () const
 
void ChooseProcessToApply ()
 
const G4VProcessGetProcessToApply () const
 
void AddCrossSection (const G4VProcess *, G4double)
 
G4double GetTriggeredProcessXSfraction ()
 
void PostStepInteractionOccured (const G4VProcess *)
 
void SetInteractionOccured (G4bool b)
 
G4bool GetInteractionOccured () const
 
- Public Member Functions inherited from G4VBiasingOperation
 G4VBiasingOperation (G4String name)
 
virtual ~G4VBiasingOperation ()
 
virtual void AlongMoveBy (const G4BiasingProcessInterface *, const G4Step *, G4double)
 
const G4StringGetName () const
 
std::size_t GetUniqueID () const
 

Detailed Description

Definition at line 58 of file G4BOptnForceCommonTruncatedExp.hh.

Constructor & Destructor Documentation

G4BOptnForceCommonTruncatedExp::G4BOptnForceCommonTruncatedExp ( G4String  name)

Definition at line 33 of file G4BOptnForceCommonTruncatedExp.cc.

References G4ILawCommonTruncatedExp::SetOperation().

34  : G4VBiasingOperation(name),
35  fInteractionOccured(false)
36 {
37  fCommonTruncatedExpLaw = new G4ILawCommonTruncatedExp("ExpLawForOperation"+name);
38  fCommonTruncatedExpLaw->SetOperation( this );
39  fTotalCrossSection = 0.0;
40 }
void SetOperation(G4BOptnForceCommonTruncatedExp *operation)
G4VBiasingOperation(G4String name)
G4BOptnForceCommonTruncatedExp::~G4BOptnForceCommonTruncatedExp ( )
virtual

Definition at line 42 of file G4BOptnForceCommonTruncatedExp.cc.

43 {}

Member Function Documentation

void G4BOptnForceCommonTruncatedExp::AddCrossSection ( const G4VProcess process,
G4double  crossSection 
)

Definition at line 90 of file G4BOptnForceCommonTruncatedExp.cc.

References G4ILawCommonTruncatedExp::SetForceCrossSection(), and G4ILawCommonTruncatedExp::SetNumberOfSharing().

91 {
92  fTotalCrossSection += crossSection;
93  fCrossSections[process] = crossSection;
94  fNumberOfSharing = fCrossSections.size();
95  fCommonTruncatedExpLaw->SetNumberOfSharing ( fNumberOfSharing );
96  fCommonTruncatedExpLaw->SetForceCrossSection( fTotalCrossSection );
97 }
virtual G4VParticleChange* G4BOptnForceCommonTruncatedExp::ApplyFinalStateBiasing ( const G4BiasingProcessInterface ,
const G4Track ,
const G4Step  
)
inlinevirtual

Implements G4VBiasingOperation.

Definition at line 76 of file G4BOptnForceCommonTruncatedExp.hh.

78  {return 0;}
void G4BOptnForceCommonTruncatedExp::ChooseProcessToApply ( )

Definition at line 141 of file G4BOptnForceCommonTruncatedExp.cc.

References G4UniformRand.

Referenced by Sample().

142 {
143  G4double sigmaRand = G4UniformRand() * fTotalCrossSection;
144  G4double sigmaSelect = 0.0;
145  for ( std::map< const G4VProcess*, G4double>::const_iterator it = fCrossSections.begin();
146  it != fCrossSections.end();
147  it++)
148  {
149  sigmaSelect += (*it).second;
150  if ( sigmaRand <= sigmaSelect )
151  {
152  fProcessToApply = (*it).first;
153  break;
154  }
155  }
156 }
#define G4UniformRand()
Definition: Randomize.hh:87
double G4double
Definition: G4Types.hh:76
G4bool G4BOptnForceCommonTruncatedExp::DenyProcessPostStepDoIt ( const G4BiasingProcessInterface callingProcess,
const G4Track ,
const G4Step step,
G4double  
)
virtual

Reimplemented from G4VBiasingOperation.

Definition at line 67 of file G4BOptnForceCommonTruncatedExp.cc.

References G4BiasingProcessInterface::GetAlongStepGPIL(), G4BiasingProcessInterface::GetPostStepGPIL(), and G4BiasingProcessInterface::GetWrappedProcess().

69 {
70  G4double processGPIL = callingProcess->GetPostStepGPIL() < callingProcess->GetAlongStepGPIL() ?
71  callingProcess->GetPostStepGPIL() : callingProcess->GetAlongStepGPIL() ;
72  if ( processGPIL <= step->GetStepLength() )
73  {
74  G4bool denyProcess = (callingProcess->GetWrappedProcess() != fProcessToApply );
75  fInteractionOccured = (fInteractionOccured || !denyProcess);
76  return denyProcess;
77  }
78  else return true;
79 }
G4VProcess * GetWrappedProcess() const
bool G4bool
Definition: G4Types.hh:79
double G4double
Definition: G4Types.hh:76
virtual G4double G4BOptnForceCommonTruncatedExp::DistanceToApplyOperation ( const G4Track ,
G4double  ,
G4ForceCondition  
)
inlinevirtual

Implements G4VBiasingOperation.

Definition at line 79 of file G4BOptnForceCommonTruncatedExp.hh.

References DBL_MAX.

81  {return DBL_MAX;}
#define DBL_MAX
Definition: templates.hh:83
virtual G4VParticleChange* G4BOptnForceCommonTruncatedExp::GenerateBiasingFinalState ( const G4Track ,
const G4Step  
)
inlinevirtual

Implements G4VBiasingOperation.

Definition at line 82 of file G4BOptnForceCommonTruncatedExp.hh.

83  {return 0;}
G4ILawCommonTruncatedExp* G4BOptnForceCommonTruncatedExp::GetCommonTruncatedExpLaw ( )
inline

Definition at line 89 of file G4BOptnForceCommonTruncatedExp.hh.

90  {
91  return fCommonTruncatedExpLaw;
92  }
const G4ThreeVector& G4BOptnForceCommonTruncatedExp::GetInitialMomentum ( ) const
inline

Definition at line 96 of file G4BOptnForceCommonTruncatedExp.hh.

96 { return fInitialMomentum; }
G4bool G4BOptnForceCommonTruncatedExp::GetInteractionOccured ( ) const
inline

Definition at line 104 of file G4BOptnForceCommonTruncatedExp.hh.

104 { return fInteractionOccured; }
G4double G4BOptnForceCommonTruncatedExp::GetMaximumDistance ( ) const
inline

Definition at line 97 of file G4BOptnForceCommonTruncatedExp.hh.

97 { return fMaximumDistance;}
const G4VProcess* G4BOptnForceCommonTruncatedExp::GetProcessToApply ( ) const
inline

Definition at line 99 of file G4BOptnForceCommonTruncatedExp.hh.

99 { return fProcessToApply; }
G4double G4BOptnForceCommonTruncatedExp::GetTriggeredProcessXSfraction ( )

Definition at line 81 of file G4BOptnForceCommonTruncatedExp.cc.

Referenced by G4ILawCommonTruncatedExp::ComputeEffectiveCrossSectionAt().

82 {
83  return fCrossSections[fProcessToApply] / fTotalCrossSection;
84 }
void G4BOptnForceCommonTruncatedExp::Initialize ( const G4Track track)

Definition at line 100 of file G4BOptnForceCommonTruncatedExp.cc.

References DBL_MIN, G4VSolid::DistanceToOut(), G4VPhysicalVolume::GetLogicalVolume(), G4Track::GetMomentum(), G4Track::GetMomentumDirection(), G4Track::GetPosition(), G4LogicalVolume::GetSolid(), G4TransportationManager::GetTransportationManager(), G4Track::GetVolume(), G4ILawCommonTruncatedExp::reset(), and G4ILawCommonTruncatedExp::SetMaximumDistance().

101 {
102  fCrossSections.clear();
103  fTotalCrossSection = 0.0;
104  fNumberOfSharing = 0;
105  fProcessToApply = 0;
106  fInitialMomentum = track->GetMomentum();
107 
108  fCommonTruncatedExpLaw->reset();
109  G4VSolid* currentSolid = track->GetVolume()->GetLogicalVolume()->GetSolid();
111  GetNavigatorForTracking()->
112  GetGlobalToLocalTransform()).TransformPoint(track->GetPosition());
114  GetNavigatorForTracking()->
115  GetGlobalToLocalTransform()).TransformAxis(track->GetMomentumDirection());
116  fMaximumDistance = currentSolid->DistanceToOut(localPosition, localDirection);
117  if ( fMaximumDistance <= DBL_MIN ) fMaximumDistance = 0.0;
118  fCommonTruncatedExpLaw->SetMaximumDistance( fMaximumDistance );
119 }
const G4ThreeVector & GetPosition() const
static G4TransportationManager * GetTransportationManager()
G4ThreeVector GetMomentum() const
const G4ThreeVector & GetMomentumDirection() const
G4LogicalVolume * GetLogicalVolume() const
#define DBL_MIN
Definition: templates.hh:75
G4VPhysicalVolume * GetVolume() const
virtual G4double DistanceToOut(const G4ThreeVector &p, const G4ThreeVector &v, const G4bool calcNorm=false, G4bool *validNorm=0, G4ThreeVector *n=0) const =0
G4VSolid * GetSolid() const
void G4BOptnForceCommonTruncatedExp::PostStepInteractionOccured ( const G4VProcess )

Definition at line 86 of file G4BOptnForceCommonTruncatedExp.cc.

87 {
88 }
G4double G4BOptnForceCommonTruncatedExp::ProposeAlongStepLimit ( const G4BiasingProcessInterface callingProcess)
virtual

Reimplemented from G4VBiasingOperation.

Definition at line 55 of file G4BOptnForceCommonTruncatedExp.cc.

References DBL_MAX, G4ILawCommonTruncatedExp::GetInteractionDistance(), and G4BiasingProcessInterface::GetWrappedProcess().

56 {
57  if ( callingProcess->GetWrappedProcess() == fProcessToApply )
58  return fCommonTruncatedExpLaw->GetInteractionDistance();
59  else return DBL_MAX;
60 }
G4VProcess * GetWrappedProcess() const
G4double GetInteractionDistance() const
#define DBL_MAX
Definition: templates.hh:83
G4ForceCondition G4BOptnForceCommonTruncatedExp::ProposeForceCondition ( const G4ForceCondition  processCondition)
virtual

Reimplemented from G4VBiasingOperation.

Definition at line 50 of file G4BOptnForceCommonTruncatedExp.cc.

References Forced.

51 {
52  return Forced;
53 }
G4GPILSelection G4BOptnForceCommonTruncatedExp::ProposeGPILSelection ( const G4GPILSelection  processSelection)
virtual

Reimplemented from G4VBiasingOperation.

Definition at line 62 of file G4BOptnForceCommonTruncatedExp.cc.

References CandidateForSelection.

const G4VBiasingInteractionLaw * G4BOptnForceCommonTruncatedExp::ProvideOccurenceBiasingInteractionLaw ( const G4BiasingProcessInterface )
virtual

Implements G4VBiasingOperation.

Definition at line 45 of file G4BOptnForceCommonTruncatedExp.cc.

46 {
47  return fCommonTruncatedExpLaw;
48 }
void G4BOptnForceCommonTruncatedExp::Sample ( )
void G4BOptnForceCommonTruncatedExp::SetInteractionOccured ( G4bool  b)
inline

Definition at line 103 of file G4BOptnForceCommonTruncatedExp.hh.

References test::b.

103 { fInteractionOccured = b; }
void G4BOptnForceCommonTruncatedExp::UpdateForStep ( const G4Step step)

Definition at line 122 of file G4BOptnForceCommonTruncatedExp.cc.

References G4ILawCommonTruncatedExp::GetMaximumDistance(), G4Step::GetStepLength(), and G4VBiasingInteractionLaw::UpdateForStep().

123 {
124  fCrossSections.clear();
125  fTotalCrossSection = 0.0;
126  fNumberOfSharing = 0;
127  fProcessToApply = 0;
128 
129  fCommonTruncatedExpLaw->UpdateForStep( step->GetStepLength() );
130  fMaximumDistance = fCommonTruncatedExpLaw->GetMaximumDistance();
131 }
G4double GetStepLength() const
G4double UpdateForStep(G4double truePathLength)

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