Geant4 11.3.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4DNABornAngle Class Reference

#include <G4DNABornAngle.hh>

+ Inheritance diagram for G4DNABornAngle:

Public Member Functions

 G4DNABornAngle (const G4String &name="")
 
 ~G4DNABornAngle () override
 
G4ThreeVectorSampleDirection (const G4DynamicParticle *dp, G4double kinEnergyFinal, G4int Z, const G4Material *mat=nullptr) override
 
G4ThreeVectorSampleDirectionForShell (const G4DynamicParticle *dp, G4double kinEnergyFinal, G4int Z, G4int shellIdx, const G4Material *mat=nullptr) override
 
void PrintGeneratorInformation () const override
 
G4DNABornAngleoperator= (const G4DNABornAngle &right)=delete
 
 G4DNABornAngle (const G4DNABornAngle &)=delete
 
- Public Member Functions inherited from G4VEmAngularDistribution
 G4VEmAngularDistribution (const G4String &name)
 
virtual ~G4VEmAngularDistribution ()
 
virtual void SamplePairDirections (const G4DynamicParticle *dp, G4double elecKinEnergy, G4double posiKinEnergy, G4ThreeVector &dirElectron, G4ThreeVector &dirPositron, G4int Z=0, const G4Material *mat=nullptr)
 
const G4StringGetName () const
 
G4VEmAngularDistributionoperator= (const G4VEmAngularDistribution &right)=delete
 
 G4VEmAngularDistribution (const G4VEmAngularDistribution &)=delete
 

Additional Inherited Members

- Protected Attributes inherited from G4VEmAngularDistribution
G4ThreeVector fLocalDirection
 
G4bool fPolarisation
 

Detailed Description

Definition at line 56 of file G4DNABornAngle.hh.

Constructor & Destructor Documentation

◆ G4DNABornAngle() [1/2]

G4DNABornAngle::G4DNABornAngle ( const G4String & name = "")

Definition at line 58 of file G4DNABornAngle.cc.

59 : G4VEmAngularDistribution("deltaBorn")
60{
61 fElectron = G4Electron::Electron();
62}
static G4Electron * Electron()
Definition G4Electron.cc:91
G4VEmAngularDistribution(const G4String &name)

Referenced by G4DNABornAngle(), and operator=().

◆ ~G4DNABornAngle()

G4DNABornAngle::~G4DNABornAngle ( )
overridedefault

◆ G4DNABornAngle() [2/2]

G4DNABornAngle::G4DNABornAngle ( const G4DNABornAngle & )
delete

Member Function Documentation

◆ operator=()

G4DNABornAngle & G4DNABornAngle::operator= ( const G4DNABornAngle & right)
delete

◆ PrintGeneratorInformation()

void G4DNABornAngle::PrintGeneratorInformation ( ) const
overridevirtual

Reimplemented from G4VEmAngularDistribution.

Definition at line 119 of file G4DNABornAngle.cc.

120{}

◆ SampleDirection()

G4ThreeVector & G4DNABornAngle::SampleDirection ( const G4DynamicParticle * dp,
G4double kinEnergyFinal,
G4int Z,
const G4Material * mat = nullptr )
overridevirtual

Implements G4VEmAngularDistribution.

Definition at line 112 of file G4DNABornAngle.cc.

115{
116 return SampleDirectionForShell(dp, secEkin, Z, 0, mat);
117}
G4ThreeVector & SampleDirectionForShell(const G4DynamicParticle *dp, G4double kinEnergyFinal, G4int Z, G4int shellIdx, const G4Material *mat=nullptr) override

◆ SampleDirectionForShell()

G4ThreeVector & G4DNABornAngle::SampleDirectionForShell ( const G4DynamicParticle * dp,
G4double kinEnergyFinal,
G4int Z,
G4int shellIdx,
const G4Material * mat = nullptr )
overridevirtual

Reimplemented from G4VEmAngularDistribution.

Definition at line 68 of file G4DNABornAngle.cc.

72{
73 G4double k = dp->GetKineticEnergy();
74 G4double cosTheta = 1.0;
75
76 if (dp->GetDefinition() == fElectron)
77 {
78 if (secKinetic < 50.*eV) cosTheta = (2.*G4UniformRand())-1.;
79 else if (secKinetic <= 200.*eV)
80 {
81 if (G4UniformRand() <= 0.1) cosTheta = (2.*G4UniformRand())-1.;
82 else cosTheta = G4UniformRand()*(std::sqrt(2.)/2);
83 }
84 else
85 {
86 G4double sin2O =
87 (1.-secKinetic/k)/(1.+secKinetic/(2.*electron_mass_c2));
88 cosTheta = std::sqrt(1.-sin2O);
89 }
90 }
91 else
92 {
93 G4double tau = k/dp->GetDefinition()->GetPDGMass();
94 G4double maxSecKinetic = 2.* electron_mass_c2*tau*(tau + 2.0);
95
96 // Restriction below 100 eV from Emfietzoglou (2000)
97
98 if (secKinetic>100*eV) cosTheta = std::min(std::sqrt(secKinetic / maxSecKinetic), 1.0);
99 else cosTheta = (2.*G4UniformRand())-1.;
100 }
101
102 G4double sint = sqrt((1.0 - cosTheta)*(1.0 + cosTheta));
103 G4double phi = twopi*G4UniformRand();
104
105 fLocalDirection.set(sint*cos(phi), sint*sin(phi), cosTheta);
107
108 return fLocalDirection;
109}
double G4double
Definition G4Types.hh:83
#define G4UniformRand()
Definition Randomize.hh:52
const G4ThreeVector & GetMomentumDirection() const
G4ParticleDefinition * GetDefinition() const
G4double GetKineticEnergy() const

Referenced by SampleDirection().


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