Geant4
9.6.p02
Main Page
Related Pages
Modules
Namespaces
Classes
Files
File List
File Members
All
Classes
Namespaces
Files
Functions
Variables
Typedefs
Enumerations
Enumerator
Friends
Macros
Groups
Pages
geant4_9_6_p02
source
processes
electromagnetic
adjoint
src
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
// $Id: G4AdjointAlongStepWeightCorrection.cc 69844 2013-05-16 09:19:33Z gcosmo $
27
//
28
#include "
G4AdjointAlongStepWeightCorrection.hh
"
29
#include "
G4Step.hh
"
30
#include "
G4ParticleDefinition.hh
"
31
#include "
G4VParticleChange.hh
"
32
#include "
G4AdjointCSManager.hh
"
33
35
//
36
37
G4AdjointAlongStepWeightCorrection::G4AdjointAlongStepWeightCorrection
(
const
G4String
&
name
,
38
G4ProcessType
type):
G4VContinuousProcess
(name, type)
39
{
fParticleChange
=
new
G4ParticleChange
();
40
currentMaterialIndex=0;
41
preStepKinEnergy=1.;
42
currentCouple=0;
43
}
44
46
//
47
48
G4AdjointAlongStepWeightCorrection::~G4AdjointAlongStepWeightCorrection
()
49
{
delete
fParticleChange
;
50
}
52
//
53
void
G4AdjointAlongStepWeightCorrection::PreparePhysicsTable
(
54
const
G4ParticleDefinition
& )
55
{
56
;
57
}
59
//
60
61
void
G4AdjointAlongStepWeightCorrection::BuildPhysicsTable
(
const
G4ParticleDefinition
& )
62
{;
63
}
65
//
66
G4VParticleChange
*
G4AdjointAlongStepWeightCorrection::AlongStepDoIt
(
const
G4Track
& track,
67
const
G4Step
& step)
68
{
69
70
fParticleChange
->
Initialize
(track);
71
72
// Get the actual (true) Step length
73
//----------------------------------
74
G4double
length = step.
GetStepLength
();
75
76
77
G4double
Tkin = step.
GetPostStepPoint
()->
GetKineticEnergy
();
78
G4ParticleDefinition
* thePartDef=
const_cast<
G4ParticleDefinition
*
>
(track.
GetDynamicParticle
()->
GetDefinition
());
79
G4double
weight_correction=
G4AdjointCSManager::GetAdjointCSManager
()->
GetContinuousWeightCorrection
(thePartDef,
80
preStepKinEnergy,Tkin, currentCouple,length);
81
82
83
84
85
//Caution!!!
86
// It is important to select the weight of the post_step_point
87
// as the current weight and not the weight of the track, as t
88
// the weight of the track is changed after having applied all
89
// the along_step_do_it.
90
91
// G4double new_weight=weight_correction*track.GetWeight(); //old
92
G4double
new_weight=weight_correction*step.
GetPostStepPoint
()->
GetWeight
();
93
94
95
96
//if (weight_correction >2.) new_weight=1.e-300;
97
98
99
//The following test check for zero weight.
100
//This happens after weight correction of gamma for photo electric effect.
101
//When the new weight is 0 it will be later on consider as nan by G4.
102
//Therefore we do put a lower limit of 1.e-300. for new_weight
103
//Correction by L.Desorgher on 15 July 2009
104
if
(new_weight==0 || (new_weight<=0 && new_weight>0)){
105
//G4cout<<new_weight<<'\t'<<weight_correction<<'\t'<<track.GetWeight()<<G4endl;
106
new_weight=1.e-300;
107
}
108
109
//G4cout<<new_weight<<'\t'<<weight_correction<<'\t'<<track.GetWeight()<<G4endl;
110
fParticleChange
->
SetParentWeightByProcess
(
false
);
111
fParticleChange
->
SetSecondaryWeightByProcess
(
false
);
112
fParticleChange
->
ProposeParentWeight
(new_weight);
113
114
115
return
fParticleChange
;
116
117
}
119
//
120
G4double
G4AdjointAlongStepWeightCorrection::GetContinuousStepLimit
(
const
G4Track
& track,
121
G4double
,
G4double
,
G4double
& )
122
{
123
G4double
x
=
DBL_MAX
;
124
DefineMaterial(track.
GetMaterialCutsCouple
());
125
preStepKinEnergy = track.
GetKineticEnergy
();
126
return
x
;
127
}
Generated on Sat May 25 2013 14:33:29 for Geant4 by
1.8.4