Geant4  10.03.p02
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
GB04BOptnBremSplitting.cc
Go to the documentation of this file.
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 // $Id$
27 //
30 
33 
35 
36 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
37 
39 : G4VBiasingOperation(name),
40  fSplittingFactor(1),
41  fParticleChange()
42 {
43 }
44 
45 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
46 
48 {
49 }
50 
54  const G4Track* track,
55  const G4Step* step,
56  G4bool& )
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 }
144 
145 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
const XML_Char * name
Definition: expat.h:151
G4int GetNumberOfSecondaries() const
G4Track * GetSecondary(G4int anIndex) const
const G4ThreeVector & GetProposedMomentumDirection() const
virtual G4VParticleChange * ApplyFinalStateBiasing(const G4BiasingProcessInterface *, const G4Track *, const G4Step *, G4bool &)
G4double GetProposedKineticEnergy() const
int G4int
Definition: G4Types.hh:78
void SetWeight(G4double aValue)
G4VProcess * GetWrappedProcess() const
void SetSecondaryWeightByProcess(G4bool)
bool G4bool
Definition: G4Types.hh:79
Definition: G4Step.hh:76
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
Definition of the GB04BOptnBremSplitting class.