Geant4
9.6.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4NeutronKiller.hh
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
//
28
//---------------------------------------------------------------------------
29
//
30
// ClassName: G4NeutronKiller
31
//
32
// Description: The process to kill neutrons to save CPU
33
//
34
// Author: V.Ivanchenko 26/09/00 for HARP software
35
//
36
//----------------------------------------------------------------------------
37
//
38
// Class description:
39
//
40
// G4NeutronKiller allows to remove unwanted neutrons from simulation in
41
// order to improve CPU performance. There are two parameters:
42
//
43
// 1) low energy threshold for neutron kinetic energy (default 0)
44
// 2) time limit for neutron track (default DBL_MAX)
45
//
46
// These parameters can be changed via Set methods or by UI commands:
47
// /physics_engine/neutron/energyCut
48
// /physics_engine/neutron/timeLimit
49
//
50
// If a neutron track is killed no energy deposition is added to the step
51
//
52
53
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
54
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
55
56
#ifndef NeutronKiller_h
57
#define NeutronKiller_h 1
58
59
#include "
globals.hh
"
60
#include "
G4VDiscreteProcess.hh
"
61
#include "
G4ParticleDefinition.hh
"
62
#include "
G4Step.hh
"
63
#include "
G4Track.hh
"
64
65
class
G4NeutronKillerMessenger
;
66
67
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
68
69
class
G4NeutronKiller
:
public
G4VDiscreteProcess
70
{
71
public
:
72
73
G4NeutronKiller
(
const
G4String
& processName =
"nKiller"
,
74
G4ProcessType
aType =
fGeneral
);
75
76
virtual
~G4NeutronKiller
();
77
78
G4bool
IsApplicable
(
const
G4ParticleDefinition
&);
79
80
void
SetTimeLimit
(
G4double
);
81
82
void
SetKinEnergyLimit
(
G4double
);
83
84
G4double
PostStepGetPhysicalInteractionLength
(
const
G4Track
& track,
85
G4double
previousStepSize,
86
G4ForceCondition
*
condition
);
87
88
G4VParticleChange
*
PostStepDoIt
(
const
G4Track
&,
const
G4Step
&);
89
90
G4double
GetMeanFreePath
(
const
G4Track
&,
G4double
,
G4ForceCondition
*);
91
92
private
:
93
94
// hide assignment operator as private
95
G4NeutronKiller
(
const
G4NeutronKiller
&);
96
G4NeutronKiller
& operator = (
const
G4NeutronKiller
&right);
97
98
G4double
kinEnergyThreshold;
99
G4double
timeThreshold;
100
101
G4NeutronKillerMessenger
* pMess;
102
};
103
104
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
105
106
inline
G4double
G4NeutronKiller::PostStepGetPhysicalInteractionLength
(
107
const
G4Track
& aTrack,
108
G4double
,
G4ForceCondition
*
condition
)
109
{
110
// condition is set to "Not Forced"
111
*
condition
=
NotForced
;
112
113
G4double
limit =
DBL_MAX
;
114
if
(aTrack.
GetGlobalTime
() > timeThreshold ||
115
aTrack.
GetKineticEnergy
() < kinEnergyThreshold) limit = 0.0;
116
return
limit;
117
}
118
119
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
120
121
inline
G4double
G4NeutronKiller::GetMeanFreePath
(
const
G4Track
&,
G4double
,
122
G4ForceCondition
*)
123
{
124
return
DBL_MAX
;
125
}
126
127
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
128
129
inline
G4VParticleChange
*
G4NeutronKiller::PostStepDoIt
(
const
G4Track
& aTrack,
130
const
G4Step
&)
131
{
132
pParticleChange
->
Initialize
(aTrack);
133
pParticleChange
->
ProposeTrackStatus
(
fStopAndKill
);
134
return
pParticleChange
;
135
}
136
137
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
138
139
#endif
140
condition
G4double condition(const G4ErrorSymMatrix &m)
G4ForceCondition
G4ForceCondition
Definition:
G4ForceCondition.hh:51
NotForced
@ NotForced
Definition:
G4ForceCondition.hh:56
G4ParticleDefinition.hh
G4ProcessType
G4ProcessType
Definition:
G4ProcessType.hh:44
fGeneral
@ fGeneral
Definition:
G4ProcessType.hh:52
G4Step.hh
fStopAndKill
@ fStopAndKill
Definition:
G4TrackStatus.hh:56
G4Track.hh
G4double
double G4double
Definition:
G4Types.hh:64
G4bool
bool G4bool
Definition:
G4Types.hh:67
G4VDiscreteProcess.hh
G4NeutronKillerMessenger
Definition:
G4NeutronKillerMessenger.hh:54
G4NeutronKiller
Definition:
G4NeutronKiller.hh:70
G4NeutronKiller::SetTimeLimit
void SetTimeLimit(G4double)
Definition:
G4NeutronKiller.cc:76
G4NeutronKiller::~G4NeutronKiller
virtual ~G4NeutronKiller()
Definition:
G4NeutronKiller.cc:62
G4NeutronKiller::PostStepDoIt
G4VParticleChange * PostStepDoIt(const G4Track &, const G4Step &)
Definition:
G4NeutronKiller.hh:129
G4NeutronKiller::PostStepGetPhysicalInteractionLength
G4double PostStepGetPhysicalInteractionLength(const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
Definition:
G4NeutronKiller.hh:106
G4NeutronKiller::IsApplicable
G4bool IsApplicable(const G4ParticleDefinition &)
Definition:
G4NeutronKiller.cc:69
G4NeutronKiller::GetMeanFreePath
G4double GetMeanFreePath(const G4Track &, G4double, G4ForceCondition *)
Definition:
G4NeutronKiller.hh:121
G4NeutronKiller::SetKinEnergyLimit
void SetKinEnergyLimit(G4double)
Definition:
G4NeutronKiller.cc:86
G4ParticleDefinition
Definition:
G4ParticleDefinition.hh:68
G4Step
Definition:
G4Step.hh:78
G4String
Definition:
G4String.hh:105
G4Track
Definition:
G4Track.hh:72
G4Track::GetGlobalTime
G4double GetGlobalTime() const
G4Track::GetKineticEnergy
G4double GetKineticEnergy() const
G4VDiscreteProcess
Definition:
G4VDiscreteProcess.hh:59
G4VParticleChange
Definition:
G4VParticleChange.hh:95
G4VParticleChange::ProposeTrackStatus
void ProposeTrackStatus(G4TrackStatus status)
G4VParticleChange::Initialize
virtual void Initialize(const G4Track &)
G4VProcess::pParticleChange
G4VParticleChange * pParticleChange
Definition:
G4VProcess.hh:283
globals.hh
DBL_MAX
#define DBL_MAX
Definition:
templates.hh:83
geant4-v9.6.0
source
processes
transportation
include
G4NeutronKiller.hh
Generated by
1.9.6