Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
GB04BOptnBremSplitting Class Reference

#include <GB04BOptnBremSplitting.hh>

Inheritance diagram for GB04BOptnBremSplitting:
Collaboration diagram for GB04BOptnBremSplitting:

Public Member Functions

 GB04BOptnBremSplitting (G4String name)
 
virtual ~GB04BOptnBremSplitting ()
 
virtual const
G4VBiasingInteractionLaw
ProvideOccurenceBiasingInteractionLaw (const G4BiasingProcessInterface *, G4ForceCondition &)
 
virtual G4VParticleChangeApplyFinalStateBiasing (const G4BiasingProcessInterface *, const G4Track *, const G4Step *, G4bool &)
 
virtual G4double DistanceToApplyOperation (const G4Track *, G4double, G4ForceCondition *)
 
virtual G4VParticleChangeGenerateBiasingFinalState (const G4Track *, const G4Step *)
 
void SetSplittingFactor (G4int splittingFactor)
 
G4int GetSplittingFactor () const
 
- Public Member Functions inherited from G4VBiasingOperation
 G4VBiasingOperation (G4String name)
 
virtual ~G4VBiasingOperation ()
 
virtual G4double ProposeAlongStepLimit (const G4BiasingProcessInterface *)
 
virtual G4GPILSelection ProposeGPILSelection (const G4GPILSelection wrappedProcessSelection)
 
virtual void AlongMoveBy (const G4BiasingProcessInterface *, const G4Step *, G4double)
 
const G4StringGetName () const
 
std::size_t GetUniqueID () const
 

Detailed Description

Definition at line 48 of file GB04BOptnBremSplitting.hh.

Constructor & Destructor Documentation

GB04BOptnBremSplitting::GB04BOptnBremSplitting ( G4String  name)

Definition at line 38 of file GB04BOptnBremSplitting.cc.

39 : G4VBiasingOperation(name),
40  fSplittingFactor(1),
41  fParticleChange()
42 {
43 }
G4VBiasingOperation(G4String name)
GB04BOptnBremSplitting::~GB04BOptnBremSplitting ( )
virtual

Definition at line 47 of file GB04BOptnBremSplitting.cc.

48 {
49 }

Member Function Documentation

G4VParticleChange * GB04BOptnBremSplitting::ApplyFinalStateBiasing ( const G4BiasingProcessInterface callingProcess,
const G4Track track,
const G4Step step,
G4bool  
)
virtual

Implements G4VBiasingOperation.

Definition at line 53 of file GB04BOptnBremSplitting.cc.

57 {
58 
59  // -- Collect brem. process (wrapped process) final state:
60  G4VParticleChange* processFinalState =
61  callingProcess->GetWrappedProcess()->PostStepDoIt(*track, *step);
62 
63  // -- if no splitting requested, let the brem. process to return directly its
64  // -- generated final state:
65  if ( fSplittingFactor == 1 ) return processFinalState;
66 
67  // -- a special case here: the brem. process corrects for cross-section change
68  // -- over the step due to energy loss by sometimes "abandoning" the interaction,
69  // -- returning an unchanged incoming electron/positron.
70  // -- We respect this correction, and if no secondary is produced, its means this
71  // -- case is happening:
72  if ( processFinalState->GetNumberOfSecondaries() == 0 ) return processFinalState;
73 
74  // -- Now start the biasing:
75  // -- - the electron state will be taken as the first one produced by the brem.
76  // -- process, hence the one stored in above processFinalState particle change.
77  // -- This state will be stored in our fParticleChange object.
78  // -- - the photon accompagnying the electron will be stored also this way.
79  // -- - we will then do fSplittingFactor - 1 call to the brem. process to collect
80  // -- fSplittingFactor - 1 additionnal gammas. All these will be stored in our
81  // -- fParticleChange object.
82 
83  // -- We called the brem. process above. Its concrete particle change is indeed
84  // -- a "G4ParticleChangeForLoss" object. We cast this particle change to access
85  // -- methods of the concrete G4ParticleChangeForLoss type:
86  G4ParticleChangeForLoss* actualParticleChange =
87  ( G4ParticleChangeForLoss* ) processFinalState ;
88 
89  fParticleChange.Initialize(*track);
90 
91  // -- Store electron final state:
92  fParticleChange.
93  ProposeTrackStatus ( actualParticleChange->GetTrackStatus() );
94  fParticleChange.
95  ProposeEnergy ( actualParticleChange->GetProposedKineticEnergy() );
96  fParticleChange.
97  ProposeMomentumDirection( actualParticleChange->GetProposedMomentumDirection() );
98 
99  // -- Now deal with the gamma's:
100  // -- their common weight:
101  G4double gammaWeight = track->GetWeight() / fSplittingFactor;
102 
103  // -- inform we will have fSplittingFactor gamma's:
104  fParticleChange.SetNumberOfSecondaries( fSplittingFactor );
105 
106  // -- inform we take care of secondaries weight (otherwise these
107  // -- secondaries are by default given the primary weight).
108  fParticleChange.SetSecondaryWeightByProcess(true);
109 
110  // -- Store first gamma:
111  G4Track* gammaTrack = actualParticleChange->GetSecondary(0);
112  gammaTrack->SetWeight( gammaWeight );
113  fParticleChange.AddSecondary( gammaTrack );
114  // -- and clean-up the brem. process particle change:
115  actualParticleChange->Clear();
116 
117  // -- now start the fSplittingFactor-1 calls to the brem. process to store each
118  // -- related gamma:
119  G4int nCalls = 1;
120  while ( nCalls < fSplittingFactor )
121  {
122  // ( note: we don't need to cast to actual type here, as methods for accessing
123  // secondary particles are from base class G4VParticleChange )
124  processFinalState = callingProcess->GetWrappedProcess()->PostStepDoIt(*track, *step);
125  if ( processFinalState->GetNumberOfSecondaries() == 1 )
126  {
127  gammaTrack = processFinalState->GetSecondary(0);
128  gammaTrack->SetWeight( gammaWeight );
129  fParticleChange.AddSecondary( gammaTrack );
130  nCalls++;
131  }
132  // -- very rare special case: we ignore for now.
133  else if ( processFinalState->GetNumberOfSecondaries() > 1 )
134  {
135  for ( G4int i = 0 ; i < processFinalState->GetNumberOfSecondaries() ; i++)
136  delete processFinalState->GetSecondary(i);
137  }
138  processFinalState->Clear();
139  }
140 
141  // -- we are done:
142  return &fParticleChange;
143 }
G4int GetNumberOfSecondaries() const
G4Track * GetSecondary(G4int anIndex) const
const G4ThreeVector & GetProposedMomentumDirection() const
G4double GetProposedKineticEnergy() const
int G4int
Definition: G4Types.hh:78
void SetWeight(G4double aValue)
G4VProcess * GetWrappedProcess() const
void SetSecondaryWeightByProcess(G4bool)
virtual void Initialize(const G4Track &)
void SetNumberOfSecondaries(G4int totSecondaries)
void AddSecondary(G4Track *aSecondary)
G4double GetWeight() const
G4TrackStatus GetTrackStatus() const
double G4double
Definition: G4Types.hh:76
virtual G4VParticleChange * PostStepDoIt(const G4Track &track, const G4Step &stepData)=0

Here is the call graph for this function:

virtual G4double GB04BOptnBremSplitting::DistanceToApplyOperation ( const G4Track ,
G4double  ,
G4ForceCondition  
)
inlinevirtual

Implements G4VBiasingOperation.

Definition at line 72 of file GB04BOptnBremSplitting.hh.

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

Implements G4VBiasingOperation.

Definition at line 76 of file GB04BOptnBremSplitting.hh.

78  {return 0;}
G4int GB04BOptnBremSplitting::GetSplittingFactor ( ) const
inline

Definition at line 88 of file GB04BOptnBremSplitting.hh.

88 { return fSplittingFactor; }
virtual const G4VBiasingInteractionLaw* GB04BOptnBremSplitting::ProvideOccurenceBiasingInteractionLaw ( const G4BiasingProcessInterface ,
G4ForceCondition  
)
inlinevirtual

Implements G4VBiasingOperation.

Definition at line 61 of file GB04BOptnBremSplitting.hh.

63  { return 0; }
void GB04BOptnBremSplitting::SetSplittingFactor ( G4int  splittingFactor)
inline

Definition at line 86 of file GB04BOptnBremSplitting.hh.

87  { fSplittingFactor = splittingFactor; }

Here is the caller graph for this function:


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