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

#include <G4LowEnergyCompton.hh>

Inheritance diagram for G4LowEnergyCompton:
Collaboration diagram for G4LowEnergyCompton:

Public Member Functions

 G4LowEnergyCompton (const G4String &processName="LowEnCompton")
 
 ~G4LowEnergyCompton ()
 
G4bool IsApplicable (const G4ParticleDefinition &definition)
 
void BuildPhysicsTable (const G4ParticleDefinition &definition)
 
G4VParticleChangePostStepDoIt (const G4Track &track, const G4Step &step)
 
G4double DumpMeanFreePath (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
 
- Public Member Functions inherited from G4VDiscreteProcess
 G4VDiscreteProcess (const G4String &, G4ProcessType aType=fNotDefined)
 
 G4VDiscreteProcess (G4VDiscreteProcess &)
 
virtual ~G4VDiscreteProcess ()
 
virtual G4double PostStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
 
virtual G4double AlongStepGetPhysicalInteractionLength (const G4Track &, G4double, G4double, G4double &, G4GPILSelection *)
 
virtual G4double AtRestGetPhysicalInteractionLength (const G4Track &, G4ForceCondition *)
 
virtual G4VParticleChangeAtRestDoIt (const G4Track &, const G4Step &)
 
virtual G4VParticleChangeAlongStepDoIt (const G4Track &, const G4Step &)
 
- Public Member Functions inherited from G4VProcess
 G4VProcess (const G4String &aName="NoName", G4ProcessType aType=fNotDefined)
 
 G4VProcess (const G4VProcess &right)
 
virtual ~G4VProcess ()
 
G4int operator== (const G4VProcess &right) const
 
G4int operator!= (const G4VProcess &right) const
 
G4double GetCurrentInteractionLength () const
 
void SetPILfactor (G4double value)
 
G4double GetPILfactor () const
 
G4double AlongStepGPIL (const G4Track &track, G4double previousStepSize, G4double currentMinimumStep, G4double &proposedSafety, G4GPILSelection *selection)
 
G4double AtRestGPIL (const G4Track &track, G4ForceCondition *condition)
 
G4double PostStepGPIL (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
 
virtual void PreparePhysicsTable (const G4ParticleDefinition &)
 
virtual G4bool StorePhysicsTable (const G4ParticleDefinition *, const G4String &, G4bool)
 
virtual G4bool RetrievePhysicsTable (const G4ParticleDefinition *, const G4String &, G4bool)
 
const G4StringGetPhysicsTableFileName (const G4ParticleDefinition *, const G4String &directory, const G4String &tableName, G4bool ascii=false)
 
const G4StringGetProcessName () const
 
G4ProcessType GetProcessType () const
 
void SetProcessType (G4ProcessType)
 
G4int GetProcessSubType () const
 
void SetProcessSubType (G4int)
 
virtual void StartTracking (G4Track *)
 
virtual void EndTracking ()
 
virtual void SetProcessManager (const G4ProcessManager *)
 
virtual const G4ProcessManagerGetProcessManager ()
 
virtual void ResetNumberOfInteractionLengthLeft ()
 
G4double GetNumberOfInteractionLengthLeft () const
 
G4double GetTotalNumberOfInteractionLengthTraversed () const
 
G4bool isAtRestDoItIsEnabled () const
 
G4bool isAlongStepDoItIsEnabled () const
 
G4bool isPostStepDoItIsEnabled () const
 
virtual void DumpInfo () const
 
void SetVerboseLevel (G4int value)
 
G4int GetVerboseLevel () const
 
virtual void SetMasterProcess (G4VProcess *masterP)
 
const G4VProcessGetMasterProcess () const
 
virtual void BuildWorkerPhysicsTable (const G4ParticleDefinition &part)
 
virtual void PrepareWorkerPhysicsTable (const G4ParticleDefinition &)
 

Protected Member Functions

G4double GetMeanFreePath (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
 
- Protected Member Functions inherited from G4VProcess
void SubtractNumberOfInteractionLengthLeft (G4double previousStepSize)
 
void ClearNumberOfInteractionLengthLeft ()
 

Additional Inherited Members

- Static Public Member Functions inherited from G4VProcess
static const G4StringGetProcessTypeName (G4ProcessType)
 
- Protected Attributes inherited from G4VProcess
const G4ProcessManageraProcessManager
 
G4VParticleChangepParticleChange
 
G4ParticleChange aParticleChange
 
G4double theNumberOfInteractionLengthLeft
 
G4double currentInteractionLength
 
G4double theInitialNumberOfInteractionLength
 
G4String theProcessName
 
G4String thePhysicsTableFileName
 
G4ProcessType theProcessType
 
G4int theProcessSubType
 
G4double thePILfactor
 
G4bool enableAtRestDoIt
 
G4bool enableAlongStepDoIt
 
G4bool enablePostStepDoIt
 
G4int verboseLevel
 

Detailed Description

Definition at line 61 of file G4LowEnergyCompton.hh.

Constructor & Destructor Documentation

G4LowEnergyCompton::G4LowEnergyCompton ( const G4String processName = "LowEnCompton")

Definition at line 73 of file G4LowEnergyCompton.cc.

74  : G4VDiscreteProcess(processName),
75  lowEnergyLimit(250*eV),
76  highEnergyLimit(100*GeV),
77  intrinsicLowEnergyLimit(10*eV),
78  intrinsicHighEnergyLimit(100*GeV)
79 {
80  if (lowEnergyLimit < intrinsicLowEnergyLimit ||
81  highEnergyLimit > intrinsicHighEnergyLimit)
82  {
83  G4Exception("G4LowEnergyCompton::G4LowEnergyCompton()",
84  "OutOfRange", FatalException,
85  "Energy outside intrinsic process validity range!");
86  }
87 
88  crossSectionHandler = new G4RDCrossSectionHandler;
89 
90  G4RDVDataSetAlgorithm* scatterInterpolation = new G4RDLogLogInterpolation;
91  G4String scatterFile = "comp/ce-sf-";
92  scatterFunctionData = new G4RDCompositeEMDataSet(scatterInterpolation, 1., 1.);
93  scatterFunctionData->LoadData(scatterFile);
94 
95  meanFreePathTable = 0;
96 
97  rangeTest = new G4RDRangeNoTest;
98 
99  // For Doppler broadening
100  shellData.SetOccupancyData();
101 
102  if (verboseLevel > 0)
103  {
104  G4cout << GetProcessName() << " is created " << G4endl
105  << "Energy range: "
106  << lowEnergyLimit / keV << " keV - "
107  << highEnergyLimit / GeV << " GeV"
108  << G4endl;
109  }
110 }
void SetOccupancyData()
G4int verboseLevel
Definition: G4VProcess.hh:368
G4GLOB_DLL std::ostream G4cout
static constexpr double eV
Definition: G4SIunits.hh:215
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
static constexpr double GeV
Definition: G4SIunits.hh:217
#define G4endl
Definition: G4ios.hh:61
static constexpr double keV
Definition: G4SIunits.hh:216
virtual G4bool LoadData(const G4String &fileName)=0

Here is the call graph for this function:

G4LowEnergyCompton::~G4LowEnergyCompton ( )

Definition at line 112 of file G4LowEnergyCompton.cc.

113 {
114  delete meanFreePathTable;
115  delete crossSectionHandler;
116  delete scatterFunctionData;
117  delete rangeTest;
118 }

Member Function Documentation

void G4LowEnergyCompton::BuildPhysicsTable ( const G4ParticleDefinition definition)
virtual

Reimplemented from G4VProcess.

Definition at line 120 of file G4LowEnergyCompton.cc.

121 {
122 
123  crossSectionHandler->Clear();
124  G4String crossSectionFile = "comp/ce-cs-";
125  crossSectionHandler->LoadData(crossSectionFile);
126 
127  delete meanFreePathTable;
128  meanFreePathTable = crossSectionHandler->BuildMeanFreePathForMaterials();
129 
130  // For Doppler broadening
131  G4String file = "/doppler/shell-doppler";
132  shellData.LoadData(file);
133 }
G4RDVEMDataSet * BuildMeanFreePathForMaterials(const G4DataVector *energyCuts=0)
void LoadData(const G4String &fileName)
void LoadData(const G4String &dataFile)

Here is the call graph for this function:

G4double G4LowEnergyCompton::DumpMeanFreePath ( const G4Track track,
G4double  previousStepSize,
G4ForceCondition condition 
)
inline

Definition at line 76 of file G4LowEnergyCompton.hh.

79  { return GetMeanFreePath(track, previousStepSize, condition); }
G4double condition(const G4ErrorSymMatrix &m)
G4double GetMeanFreePath(const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)

Here is the call graph for this function:

G4double G4LowEnergyCompton::GetMeanFreePath ( const G4Track track,
G4double  previousStepSize,
G4ForceCondition condition 
)
protectedvirtual

Implements G4VDiscreteProcess.

Definition at line 331 of file G4LowEnergyCompton.cc.

334 {
336  G4double energy = photon->GetKineticEnergy();
337  const G4MaterialCutsCouple* couple = track.GetMaterialCutsCouple();
338  size_t materialIndex = couple->GetIndex();
339 
340  G4double meanFreePath;
341  if (energy > highEnergyLimit) meanFreePath = meanFreePathTable->FindValue(highEnergyLimit,materialIndex);
342  else if (energy < lowEnergyLimit) meanFreePath = DBL_MAX;
343  else meanFreePath = meanFreePathTable->FindValue(energy,materialIndex);
344  return meanFreePath;
345 }
G4double GetKineticEnergy() const
const G4DynamicParticle * GetDynamicParticle() const
const G4MaterialCutsCouple * GetMaterialCutsCouple() const
virtual G4double FindValue(G4double x, G4int componentId=0) const =0
G4double energy(const ThreeVector &p, const G4double m)
double G4double
Definition: G4Types.hh:76
#define DBL_MAX
Definition: templates.hh:83

Here is the call graph for this function:

Here is the caller graph for this function:

G4bool G4LowEnergyCompton::IsApplicable ( const G4ParticleDefinition definition)
virtual

Reimplemented from G4VProcess.

Definition at line 326 of file G4LowEnergyCompton.cc.

327 {
328  return ( &particle == G4Gamma::Gamma() );
329 }
static G4Gamma * Gamma()
Definition: G4Gamma.cc:86

Here is the call graph for this function:

G4VParticleChange * G4LowEnergyCompton::PostStepDoIt ( const G4Track track,
const G4Step step 
)
virtual

Reimplemented from G4VDiscreteProcess.

Definition at line 135 of file G4LowEnergyCompton.cc.

137 {
138  // The scattered gamma energy is sampled according to Klein - Nishina formula.
139  // then accepted or rejected depending on the Scattering Function multiplied
140  // by factor from Klein - Nishina formula.
141  // Expression of the angular distribution as Klein Nishina
142  // angular and energy distribution and Scattering fuctions is taken from
143  // D. E. Cullen "A simple model of photon transport" Nucl. Instr. Meth.
144  // Phys. Res. B 101 (1995). Method of sampling with form factors is different
145  // data are interpolated while in the article they are fitted.
146  // Reference to the article is from J. Stepanek New Photon, Positron
147  // and Electron Interaction Data for GEANT in Energy Range from 1 eV to 10
148  // TeV (draft).
149  // The random number techniques of Butcher & Messel are used
150  // (Nucl Phys 20(1960),15).
151 
152  aParticleChange.Initialize(aTrack);
153 
154  // Dynamic particle quantities
155  const G4DynamicParticle* incidentPhoton = aTrack.GetDynamicParticle();
156  G4double photonEnergy0 = incidentPhoton->GetKineticEnergy();
157 
158  if (photonEnergy0 <= lowEnergyLimit)
159  {
163  return G4VDiscreteProcess::PostStepDoIt(aTrack,aStep);
164  }
165 
166  G4double e0m = photonEnergy0 / electron_mass_c2 ;
167  G4ParticleMomentum photonDirection0 = incidentPhoton->GetMomentumDirection();
168  G4double epsilon0 = 1. / (1. + 2. * e0m);
169  G4double epsilon0Sq = epsilon0 * epsilon0;
170  G4double alpha1 = -std::log(epsilon0);
171  G4double alpha2 = 0.5 * (1. - epsilon0Sq);
172  G4double wlPhoton = h_Planck*c_light/photonEnergy0;
173 
174  // Select randomly one element in the current material
175  const G4MaterialCutsCouple* couple = aTrack.GetMaterialCutsCouple();
176  G4int Z = crossSectionHandler->SelectRandomAtom(couple,photonEnergy0);
177 
178  // Sample the energy of the scattered photon
180  G4double epsilonSq;
181  G4double oneCosT;
182  G4double sinT2;
183  G4double gReject;
184  do
185  {
186  if ( alpha1/(alpha1+alpha2) > G4UniformRand())
187  {
188  epsilon = std::exp(-alpha1 * G4UniformRand()); // std::pow(epsilon0,G4UniformRand())
189  epsilonSq = epsilon * epsilon;
190  }
191  else
192  {
193  epsilonSq = epsilon0Sq + (1. - epsilon0Sq) * G4UniformRand();
194  epsilon = std::sqrt(epsilonSq);
195  }
196 
197  oneCosT = (1. - epsilon) / ( epsilon * e0m);
198  sinT2 = oneCosT * (2. - oneCosT);
199  G4double x = std::sqrt(oneCosT/2.) / (wlPhoton/cm);
200  G4double scatteringFunction = scatterFunctionData->FindValue(x,Z-1);
201  gReject = (1. - epsilon * sinT2 / (1. + epsilonSq)) * scatteringFunction;
202 
203  } while(gReject < G4UniformRand()*Z);
204 
205  G4double cosTheta = 1. - oneCosT;
206  G4double sinTheta = std::sqrt (sinT2);
207  G4double phi = twopi * G4UniformRand() ;
208  G4double dirX = sinTheta * std::cos(phi);
209  G4double dirY = sinTheta * std::sin(phi);
210  G4double dirZ = cosTheta ;
211 
212  // Doppler broadening - Method based on:
213  // Y. Namito, S. Ban and H. Hirayama,
214  // "Implementation of the Doppler Broadening of a Compton-Scattered Photon Into the EGS4 Code"
215  // NIM A 349, pp. 489-494, 1994
216 
217  // Maximum number of sampling iterations
218  G4int maxDopplerIterations = 1000;
219  G4double bindingE = 0.;
220  G4double photonEoriginal = epsilon * photonEnergy0;
221  G4double photonE = -1.;
222  G4int iteration = 0;
223  G4double eMax = photonEnergy0;
224  do
225  {
226  iteration++;
227  // Select shell based on shell occupancy
228  G4int shell = shellData.SelectRandomShell(Z);
229  bindingE = shellData.BindingEnergy(Z,shell);
230 
231  eMax = photonEnergy0 - bindingE;
232 
233  // Randomly sample bound electron momentum (memento: the data set is in Atomic Units)
234  G4double pSample = profileData.RandomSelectMomentum(Z,shell);
235  // Rescale from atomic units
236  G4double pDoppler = pSample * fine_structure_const;
237  G4double pDoppler2 = pDoppler * pDoppler;
238  G4double var2 = 1. + oneCosT * e0m;
239  G4double var3 = var2*var2 - pDoppler2;
240  G4double var4 = var2 - pDoppler2 * cosTheta;
241  G4double var = var4*var4 - var3 + pDoppler2 * var3;
242  if (var > 0.)
243  {
244  G4double varSqrt = std::sqrt(var);
245  G4double scale = photonEnergy0 / var3;
246  // Random select either root
247  if (G4UniformRand() < 0.5) photonE = (var4 - varSqrt) * scale;
248  else photonE = (var4 + varSqrt) * scale;
249  }
250  else
251  {
252  photonE = -1.;
253  }
254  } while ( iteration <= maxDopplerIterations &&
255  (photonE < 0. || photonE > eMax || photonE < eMax*G4UniformRand()) );
256 
257  // End of recalculation of photon energy with Doppler broadening
258  // Revert to original if maximum number of iterations threshold has been reached
259  if (iteration >= maxDopplerIterations)
260  {
261  photonE = photonEoriginal;
262  bindingE = 0.;
263  }
264 
265  // Update G4VParticleChange for the scattered photon
266 
267  G4ThreeVector photonDirection1(dirX,dirY,dirZ);
268  photonDirection1.rotateUz(photonDirection0);
269  aParticleChange.ProposeMomentumDirection(photonDirection1);
270 
271  G4double photonEnergy1 = photonE;
272  //G4cout << "--> PHOTONENERGY1 = " << photonE/keV << G4endl;
273 
274  if (photonEnergy1 > 0.)
275  {
276  aParticleChange.ProposeEnergy(photonEnergy1) ;
277  }
278  else
279  {
282  }
283 
284  // Kinematics of the scattered electron
285  G4double eKineticEnergy = photonEnergy0 - photonEnergy1 - bindingE;
286  G4double eTotalEnergy = eKineticEnergy + electron_mass_c2;
287 
288  G4double electronE = photonEnergy0 * (1. - epsilon) + electron_mass_c2;
289  G4double electronP2 = electronE*electronE - electron_mass_c2*electron_mass_c2;
290  G4double sinThetaE = -1.;
291  G4double cosThetaE = 0.;
292  if (electronP2 > 0.)
293  {
294  cosThetaE = (eTotalEnergy + photonEnergy1 )* (1. - epsilon) / std::sqrt(electronP2);
295  sinThetaE = -1. * std::sqrt(1. - cosThetaE * cosThetaE);
296  }
297 
298  G4double eDirX = sinThetaE * std::cos(phi);
299  G4double eDirY = sinThetaE * std::sin(phi);
300  G4double eDirZ = cosThetaE;
301 
302  // Generate the electron only if with large enough range w.r.t. cuts and safety
303 
304  G4double safety = aStep.GetPostStepPoint()->GetSafety();
305 
306  if (rangeTest->Escape(G4Electron::Electron(),couple,eKineticEnergy,safety))
307  {
308  G4ThreeVector eDirection(eDirX,eDirY,eDirZ);
309  eDirection.rotateUz(photonDirection0);
310 
311  G4DynamicParticle* electron = new G4DynamicParticle (G4Electron::Electron(),eDirection,eKineticEnergy) ;
313  aParticleChange.AddSecondary(electron);
314  // Binding energy deposited locally
316  }
317  else
318  {
320  aParticleChange.ProposeLocalEnergyDeposit(eKineticEnergy + bindingE);
321  }
322 
323  return G4VDiscreteProcess::PostStepDoIt( aTrack, aStep);
324 }
G4double GetKineticEnergy() const
float h_Planck
Definition: hepunit.py:263
G4double BindingEnergy(G4int Z, G4int shellIndex) const
tuple x
Definition: test.py:50
int G4int
Definition: G4Types.hh:78
void ProposeLocalEnergyDeposit(G4double anEnergyPart)
static constexpr double twopi
Definition: G4SIunits.hh:76
#define G4UniformRand()
Definition: Randomize.hh:97
const G4ThreeVector & GetMomentumDirection() const
static constexpr double cm
Definition: G4SIunits.hh:119
float electron_mass_c2
Definition: hepunit.py:274
virtual G4double FindValue(G4double x, G4int componentId=0) const =0
G4double RandomSelectMomentum(G4int Z, G4int shellIndex) const
virtual G4bool Escape(const G4ParticleDefinition *particle, const G4MaterialCutsCouple *couple, G4double energy, G4double safety) const =0
virtual void Initialize(const G4Track &)
void SetNumberOfSecondaries(G4int totSecondaries)
G4int SelectRandomAtom(const G4MaterialCutsCouple *couple, G4double e) const
void ProposeEnergy(G4double finalEnergy)
G4ParticleChange aParticleChange
Definition: G4VProcess.hh:289
void AddSecondary(G4Track *aSecondary)
static G4Electron * Electron()
Definition: G4Electron.cc:94
void ProposeMomentumDirection(G4double Px, G4double Py, G4double Pz)
double G4double
Definition: G4Types.hh:76
void ProposeTrackStatus(G4TrackStatus status)
G4int SelectRandomShell(G4int Z) const
virtual G4VParticleChange * PostStepDoIt(const G4Track &, const G4Step &)
double epsilon(double density, double temperature)
float c_light
Definition: hepunit.py:257

Here is the call graph for this function:


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