Geant4 11.1.1
Toolkit for the simulation of the passage of particles through matter
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
G4AdjointAlongStepWeightCorrection.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
28
29#include "G4AdjointCSManager.hh"
30#include "G4ParticleChange.hh"
32#include "G4Step.hh"
33#include "G4VParticleChange.hh"
34
35///////////////////////////////////////////////////////
37 const G4String& name, G4ProcessType type)
38 : G4VContinuousProcess(name, type)
39{
40 fParticleChange = new G4ParticleChange();
42}
43
44///////////////////////////////////////////////////////
46{
47 delete fParticleChange;
48}
49
50///////////////////////////////////////////////////////
52 std::ostream& out) const
53{
54 out <<
55 "Continuous processes act on adjoint particles to continuously correct their "
56 "weight during the adjoint reverse tracking. This process is needed when "
57 "the adjoint cross sections are not scaled such that the total adjoint cross "
58 "section matches the total forward cross section. By default the mode where "
59 "the total adjoint cross section is equal to the total forward cross section "
60 "is used and therefore this along step weightcorrection factor is 1. However "
61 "in some cases (some energy ranges) the total forward cross section or the "
62 "total adjoint cross section can be zero. In this case the along step weight "
63 "correction is needed and is given by exp(-(Sigma_tot_adj-Sigma_tot_fwd).dx)"
64 "\n";
65}
66
67///////////////////////////////////////////////////////
69 const G4Track& track, const G4Step& step)
70{
71 fParticleChange->Initialize(track);
72
73 // Get the actual (true) Step length
74 G4double length = step.GetStepLength();
75
77 G4ParticleDefinition* thePartDef = const_cast<G4ParticleDefinition*>(
79 G4double weight_correction = fCSManager->GetContinuousWeightCorrection(
80 thePartDef, fPreStepKinEnergy, Tkin, fCurrentCouple, length);
81
82 // Caution!!!
83 // It is important to select the weight of the post_step_point as the current
84 // weight and not the weight of the track, as the weight of the track is
85 // changed after having applied all the along_step_do_it.
86 G4double new_weight =
87 weight_correction * step.GetPostStepPoint()->GetWeight();
88
89 // The following test check for zero weight.
90 // This happens after weight correction of gamma for photo electric effect.
91 // When the new weight is 0 it will be later on considered as NaN by G4.
92 // Therefore we put a lower limit of 1.e-300. for new_weight
93 if(new_weight == 0. || (new_weight <= 0. && new_weight > 0.))
94 {
95 new_weight = 1.e-300;
96 }
97
98 fParticleChange->SetParentWeightByProcess(false);
99 fParticleChange->SetSecondaryWeightByProcess(false);
100 fParticleChange->ProposeParentWeight(new_weight);
101
102 return fParticleChange;
103}
104
105///////////////////////////////////////////////////////
107 const G4Track& track, G4double, G4double, G4double&)
108{
109 DefineMaterial(track.GetMaterialCutsCouple());
110 fPreStepKinEnergy = track.GetKineticEnergy();
111 return DBL_MAX;
112}
G4ProcessType
double G4double
Definition: G4Types.hh:83
G4double GetContinuousStepLimit(const G4Track &track, G4double previousStepSize, G4double currentMinimumStep, G4double &currentSafety) override
void ProcessDescription(std::ostream &) const override
G4VParticleChange * AlongStepDoIt(const G4Track &, const G4Step &) override
G4AdjointAlongStepWeightCorrection(const G4String &name="ContinuousWeightCorrection", G4ProcessType type=fElectromagnetic)
G4double GetContinuousWeightCorrection(G4ParticleDefinition *aPartDef, G4double PreStepEkin, G4double AfterStepEkin, const G4MaterialCutsCouple *aCouple, G4double step_length)
static G4AdjointCSManager * GetAdjointCSManager()
G4ParticleDefinition * GetDefinition() const
void Initialize(const G4Track &) override
G4double GetKineticEnergy() const
G4double GetWeight() const
Definition: G4Step.hh:62
G4double GetStepLength() const
G4StepPoint * GetPostStepPoint() const
const G4DynamicParticle * GetDynamicParticle() const
G4double GetKineticEnergy() const
const G4MaterialCutsCouple * GetMaterialCutsCouple() const
void SetSecondaryWeightByProcess(G4bool)
void SetParentWeightByProcess(G4bool)
void ProposeParentWeight(G4double finalWeight)
#define DBL_MAX
Definition: templates.hh:62