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

#include <G4ICRU73QOModel.hh>

+ Inheritance diagram for G4ICRU73QOModel:

Public Member Functions

 G4ICRU73QOModel (const G4ParticleDefinition *p=nullptr, const G4String &nam="ICRU73QO")
 
 ~G4ICRU73QOModel ()=default
 
void Initialise (const G4ParticleDefinition *, const G4DataVector &) override
 
G4double ComputeCrossSectionPerElectron (const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy, G4double maxEnergy)
 
G4double ComputeCrossSectionPerAtom (const G4ParticleDefinition *, G4double kineticEnergy, G4double Z, G4double A, G4double cutEnergy, G4double maxEnergy) override
 
G4double CrossSectionPerVolume (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy, G4double maxEnergy) override
 
G4double ComputeDEDXPerVolume (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double) override
 
void SampleSecondaries (std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin, G4double maxEnergy) override
 
void CorrectionsAlongStep (const G4MaterialCutsCouple *, const G4DynamicParticle *, const G4double &length, G4double &eloss) override
 
G4ICRU73QOModeloperator= (const G4ICRU73QOModel &right)=delete
 
 G4ICRU73QOModel (const G4ICRU73QOModel &)=delete
 
- Public Member Functions inherited from G4VEmModel
 G4VEmModel (const G4String &nam)
 
virtual ~G4VEmModel ()
 
virtual void Initialise (const G4ParticleDefinition *, const G4DataVector &)=0
 
virtual void SampleSecondaries (std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin=0.0, G4double tmax=DBL_MAX)=0
 
virtual void InitialiseLocal (const G4ParticleDefinition *, G4VEmModel *masterModel)
 
virtual void InitialiseForMaterial (const G4ParticleDefinition *, const G4Material *)
 
virtual void InitialiseForElement (const G4ParticleDefinition *, G4int Z)
 
virtual G4double ComputeDEDXPerVolume (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=DBL_MAX)
 
virtual G4double CrossSectionPerVolume (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
virtual G4double GetPartialCrossSection (const G4Material *, G4int level, const G4ParticleDefinition *, G4double kineticEnergy)
 
virtual G4double ComputeCrossSectionPerAtom (const G4ParticleDefinition *, G4double kinEnergy, G4double Z, G4double A=0., G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
virtual G4double ComputeCrossSectionPerShell (const G4ParticleDefinition *, G4int Z, G4int shellIdx, G4double kinEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
virtual G4double ChargeSquareRatio (const G4Track &)
 
virtual G4double GetChargeSquareRatio (const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
 
virtual G4double GetParticleCharge (const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
 
virtual void StartTracking (G4Track *)
 
virtual void CorrectionsAlongStep (const G4MaterialCutsCouple *, const G4DynamicParticle *, const G4double &length, G4double &eloss)
 
virtual G4double Value (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy)
 
virtual G4double MinPrimaryEnergy (const G4Material *, const G4ParticleDefinition *, G4double cut=0.0)
 
virtual G4double MinEnergyCut (const G4ParticleDefinition *, const G4MaterialCutsCouple *)
 
virtual void SetupForMaterial (const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
 
virtual void DefineForRegion (const G4Region *)
 
virtual void FillNumberOfSecondaries (G4int &numberOfTriplets, G4int &numberOfRecoil)
 
virtual void ModelDescription (std::ostream &outFile) const
 
void InitialiseElementSelectors (const G4ParticleDefinition *, const G4DataVector &)
 
std::vector< G4EmElementSelector * > * GetElementSelectors ()
 
void SetElementSelectors (std::vector< G4EmElementSelector * > *)
 
G4double ComputeDEDX (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=DBL_MAX)
 
G4double CrossSection (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
G4double ComputeMeanFreePath (const G4ParticleDefinition *, G4double kineticEnergy, const G4Material *, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
G4double ComputeCrossSectionPerAtom (const G4ParticleDefinition *, const G4Element *, G4double kinEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
const G4ElementSelectRandomAtom (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
const G4ElementSelectTargetAtom (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double logKineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
const G4ElementSelectRandomAtom (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
const G4ElementGetCurrentElement (const G4Material *mat=nullptr) const
 
G4int SelectRandomAtomNumber (const G4Material *) const
 
const G4IsotopeGetCurrentIsotope (const G4Element *elm=nullptr) const
 
G4int SelectIsotopeNumber (const G4Element *) const
 
void SetParticleChange (G4VParticleChange *, G4VEmFluctuationModel *f=nullptr)
 
void SetCrossSectionTable (G4PhysicsTable *, G4bool isLocal)
 
G4ElementDataGetElementData ()
 
G4PhysicsTableGetCrossSectionTable ()
 
G4VEmFluctuationModelGetModelOfFluctuations ()
 
G4VEmAngularDistributionGetAngularDistribution ()
 
G4VEmModelGetTripletModel ()
 
void SetTripletModel (G4VEmModel *)
 
void SetAngularDistribution (G4VEmAngularDistribution *)
 
G4double HighEnergyLimit () const
 
G4double LowEnergyLimit () const
 
G4double HighEnergyActivationLimit () const
 
G4double LowEnergyActivationLimit () const
 
G4double PolarAngleLimit () const
 
G4double SecondaryThreshold () const
 
G4bool LPMFlag () const
 
G4bool DeexcitationFlag () const
 
G4bool ForceBuildTableFlag () const
 
G4bool UseAngularGeneratorFlag () const
 
void SetAngularGeneratorFlag (G4bool)
 
void SetHighEnergyLimit (G4double)
 
void SetLowEnergyLimit (G4double)
 
void SetActivationHighEnergyLimit (G4double)
 
void SetActivationLowEnergyLimit (G4double)
 
G4bool IsActive (G4double kinEnergy) const
 
void SetPolarAngleLimit (G4double)
 
void SetSecondaryThreshold (G4double)
 
void SetLPMFlag (G4bool val)
 
void SetDeexcitationFlag (G4bool val)
 
void SetForceBuildTable (G4bool val)
 
void SetFluctuationFlag (G4bool val)
 
void SetMasterThread (G4bool val)
 
G4bool IsMaster () const
 
void SetUseBaseMaterials (G4bool val)
 
G4bool UseBaseMaterials () const
 
G4double MaxSecondaryKinEnergy (const G4DynamicParticle *dynParticle)
 
const G4StringGetName () const
 
void SetCurrentCouple (const G4MaterialCutsCouple *)
 
G4bool IsLocked () const
 
void SetLocked (G4bool)
 
G4VEmModeloperator= (const G4VEmModel &right)=delete
 
 G4VEmModel (const G4VEmModel &)=delete
 

Protected Member Functions

G4double MaxSecondaryEnergy (const G4ParticleDefinition *, G4double kinEnergy) final
 
- Protected Member Functions inherited from G4VEmModel
G4ParticleChangeForLossGetParticleChangeForLoss ()
 
G4ParticleChangeForGammaGetParticleChangeForGamma ()
 
virtual G4double MaxSecondaryEnergy (const G4ParticleDefinition *, G4double kineticEnergy)
 
const G4MaterialCutsCoupleCurrentCouple () const
 
void SetCurrentElement (const G4Element *)
 

Additional Inherited Members

- Protected Attributes inherited from G4VEmModel
G4ElementDatafElementData = nullptr
 
G4VParticleChangepParticleChange = nullptr
 
G4PhysicsTablexSectionTable = nullptr
 
const G4MaterialpBaseMaterial = nullptr
 
const std::vector< G4double > * theDensityFactor = nullptr
 
const std::vector< G4int > * theDensityIdx = nullptr
 
G4double inveplus
 
G4double pFactor = 1.0
 
size_t currentCoupleIndex = 0
 
size_t basedCoupleIndex = 0
 
G4bool lossFlucFlag = true
 

Detailed Description

Definition at line 70 of file G4ICRU73QOModel.hh.

Constructor & Destructor Documentation

◆ G4ICRU73QOModel() [1/2]

G4ICRU73QOModel::G4ICRU73QOModel ( const G4ParticleDefinition p = nullptr,
const G4String nam = "ICRU73QO" 
)
explicit

Definition at line 64 of file G4ICRU73QOModel.cc.

65 : G4VEmModel(nam),
66 particle(nullptr),
67 isInitialised(false)
68{
69 mass = charge = chargeSquare = massRate = ratio = 0.0;
70 if(p) { SetParticle(p); }
71 SetHighEnergyLimit(10.0*MeV);
72
73 lowestKinEnergy = 5.0*keV;
74
75 sizeL0 = 67;
76 sizeL1 = 22;
77 sizeL2 = 14;
78
79 theElectron = G4Electron::Electron();
80
81 for (G4int i = 0; i < 100; ++i)
82 {
83 indexZ[i] = -1;
84 }
85 for(G4int i = 0; i < NQOELEM; ++i)
86 {
87 if(ZElementAvailable[i] > 0) {
88 indexZ[ZElementAvailable[i]] = i;
89 }
90 }
91 fParticleChange = nullptr;
92 denEffData = nullptr;
93}
int G4int
Definition: G4Types.hh:85
static G4Electron * Electron()
Definition: G4Electron.cc:93
void SetHighEnergyLimit(G4double)
Definition: G4VEmModel.hh:746

◆ ~G4ICRU73QOModel()

G4ICRU73QOModel::~G4ICRU73QOModel ( )
default

◆ G4ICRU73QOModel() [2/2]

G4ICRU73QOModel::G4ICRU73QOModel ( const G4ICRU73QOModel )
delete

Member Function Documentation

◆ ComputeCrossSectionPerAtom()

G4double G4ICRU73QOModel::ComputeCrossSectionPerAtom ( const G4ParticleDefinition p,
G4double  kineticEnergy,
G4double  Z,
G4double  A,
G4double  cutEnergy,
G4double  maxEnergy 
)
overridevirtual

Reimplemented from G4VEmModel.

Definition at line 147 of file G4ICRU73QOModel.cc.

153{
155 (p,kineticEnergy,cutEnergy,maxEnergy);
156 return cross;
157}
double G4double
Definition: G4Types.hh:83
const G4int Z[17]
G4double ComputeCrossSectionPerElectron(const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy, G4double maxEnergy)

◆ ComputeCrossSectionPerElectron()

G4double G4ICRU73QOModel::ComputeCrossSectionPerElectron ( const G4ParticleDefinition p,
G4double  kineticEnergy,
G4double  cutEnergy,
G4double  maxEnergy 
)

Definition at line 121 of file G4ICRU73QOModel.cc.

126{
127 G4double cross = 0.0;
128 const G4double tmax = MaxSecondaryEnergy(p, kineticEnergy);
129 const G4double maxEnergy = std::min(tmax, maxKinEnergy);
130 const G4double cutEnergy = std::max(cut, lowestKinEnergy*massRate);
131 if(cutEnergy < maxEnergy) {
132
133 const G4double energy = kineticEnergy + mass;
134 const G4double energy2 = energy*energy;
135 const G4double beta2 = kineticEnergy*(kineticEnergy + 2.0*mass)/energy2;
136 cross = (maxEnergy - cutEnergy)/(cutEnergy*maxEnergy)
137 - beta2*G4Log(maxEnergy/cutEnergy)/tmax;
138
139 cross *= CLHEP::twopi_mc2_rcl2*chargeSquare/beta2;
140 }
141
142 return cross;
143}
G4double G4Log(G4double x)
Definition: G4Log.hh:227
G4double MaxSecondaryEnergy(const G4ParticleDefinition *, G4double kinEnergy) final
G4double energy(const ThreeVector &p, const G4double m)

Referenced by ComputeCrossSectionPerAtom(), and CrossSectionPerVolume().

◆ ComputeDEDXPerVolume()

G4double G4ICRU73QOModel::ComputeDEDXPerVolume ( const G4Material material,
const G4ParticleDefinition p,
G4double  kineticEnergy,
G4double  cut 
)
overridevirtual

Reimplemented from G4VEmModel.

Definition at line 176 of file G4ICRU73QOModel.cc.

180{
181 SetParticle(p);
182 const G4double tmax = MaxSecondaryEnergy(p, kineticEnergy);
183 const G4double tkin = kineticEnergy/massRate;
184 const G4double cutEnergy = std::max(cut, lowestKinEnergy*massRate);
185 G4double dedx = 0.0;
186 if(tkin > lowestKinEnergy) { dedx = DEDX(material, tkin); }
187 else { dedx = DEDX(material, lowestKinEnergy)*sqrt(tkin/lowestKinEnergy); }
188
189 if (cutEnergy < tmax) {
190
191 const G4double tau = kineticEnergy/mass;
192 const G4double x = cutEnergy/tmax;
193 dedx += (G4Log(x)*(tau + 1.)*(tau + 1.)/(tau * (tau + 2.0)) + 1.0 - x) *
194 CLHEP::twopi_mc2_rcl2 *chargeSquare * material->GetElectronDensity();
195 }
196 dedx = std::max(dedx, 0.0);
197 return dedx;
198}
G4double GetElectronDensity() const
Definition: G4Material.hh:212

Referenced by G4ICRU73NoDeltaModel::ComputeDEDXPerVolume().

◆ CorrectionsAlongStep()

void G4ICRU73QOModel::CorrectionsAlongStep ( const G4MaterialCutsCouple ,
const G4DynamicParticle ,
const G4double length,
G4double eloss 
)
overridevirtual

Reimplemented from G4VEmModel.

Definition at line 401 of file G4ICRU73QOModel.cc.

405{}

◆ CrossSectionPerVolume()

G4double G4ICRU73QOModel::CrossSectionPerVolume ( const G4Material material,
const G4ParticleDefinition p,
G4double  kineticEnergy,
G4double  cutEnergy,
G4double  maxEnergy 
)
overridevirtual

Reimplemented from G4VEmModel.

Definition at line 161 of file G4ICRU73QOModel.cc.

167{
168 G4double eDensity = material->GetElectronDensity();
170 (p,kineticEnergy,cutEnergy,maxEnergy);
171 return cross;
172}

◆ Initialise()

void G4ICRU73QOModel::Initialise ( const G4ParticleDefinition p,
const G4DataVector  
)
overridevirtual

Implements G4VEmModel.

Definition at line 97 of file G4ICRU73QOModel.cc.

99{
100 if(p != particle) SetParticle(p);
101
102 // always false before the run
103 SetDeexcitationFlag(false);
104
105 if(!isInitialised) {
106 isInitialised = true;
107
110 }
111
112 G4String pname = particle->GetParticleName();
113 fParticleChange = GetParticleChangeForLoss();
115 denEffData = (*mtab)[0]->GetIonisation()->GetDensityEffectData();
116 }
117}
std::vector< G4Material * > G4MaterialTable
static G4MaterialTable * GetMaterialTable()
Definition: G4Material.cc:677
const G4String & GetParticleName() const
G4VEmAngularDistribution * GetAngularDistribution()
Definition: G4VEmModel.hh:600
void SetDeexcitationFlag(G4bool val)
Definition: G4VEmModel.hh:802
void SetAngularDistribution(G4VEmAngularDistribution *)
Definition: G4VEmModel.hh:607
G4bool UseAngularGeneratorFlag() const
Definition: G4VEmModel.hh:697
G4ParticleChangeForLoss * GetParticleChangeForLoss()
Definition: G4VEmModel.cc:109

◆ MaxSecondaryEnergy()

G4double G4ICRU73QOModel::MaxSecondaryEnergy ( const G4ParticleDefinition pd,
G4double  kinEnergy 
)
finalprotectedvirtual

Reimplemented from G4VEmModel.

Definition at line 486 of file G4ICRU73QOModel.cc.

488{
489 if(pd != particle) { SetParticle(pd); }
490 G4double tau = kinEnergy/mass;
491 G4double tmax = 2.0*electron_mass_c2*tau*(tau + 2.) /
492 (1. + 2.0*(tau + 1.)*ratio + ratio*ratio);
493 return tmax;
494}

Referenced by ComputeCrossSectionPerElectron(), and ComputeDEDXPerVolume().

◆ operator=()

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

◆ SampleSecondaries()

void G4ICRU73QOModel::SampleSecondaries ( std::vector< G4DynamicParticle * > *  vdp,
const G4MaterialCutsCouple couple,
const G4DynamicParticle dp,
G4double  tmin,
G4double  maxEnergy 
)
overridevirtual

Implements G4VEmModel.

Definition at line 409 of file G4ICRU73QOModel.cc.

414{
415 const G4double tmax = MaxSecondaryKinEnergy(dp);
416 const G4double xmax = std::min(tmax, maxEnergy);
417 const G4double xmin = std::max(minEnergy, lowestKinEnergy*massRate);
418 if(xmin >= xmax) { return; }
419
420 G4double kineticEnergy = dp->GetKineticEnergy();
421 const G4double energy = kineticEnergy + mass;
422 const G4double energy2 = energy*energy;
423 const G4double beta2 = kineticEnergy*(kineticEnergy + 2.0*mass)/energy2;
424 G4double grej = 1.0;
425 G4double deltaKinEnergy, f;
426
427 G4ThreeVector direction = dp->GetMomentumDirection();
428
429 // sampling follows ...
430 do {
432 deltaKinEnergy = xmin*xmax/(xmin*(1.0 - x) + xmax*x);
433
434 f = 1.0 - beta2*deltaKinEnergy/tmax;
435
436 if(f > grej) {
437 G4cout << "G4ICRU73QOModel::SampleSecondary Warning! "
438 << "Majorant " << grej << " < "
439 << f << " for e= " << deltaKinEnergy
440 << G4endl;
441 }
442
443 // Loop checking, 03-Aug-2015, Vladimir Ivanchenko
444 } while( grej*G4UniformRand() >= f );
445
446 G4ThreeVector deltaDirection;
447
449 const G4Material* mat = couple->GetMaterial();
451
452 deltaDirection =
453 GetAngularDistribution()->SampleDirection(dp, deltaKinEnergy, Z, mat);
454
455 } else {
456
457 G4double deltaMomentum =
458 sqrt(deltaKinEnergy * (deltaKinEnergy + 2.0*electron_mass_c2));
459 G4double totMomentum = energy*sqrt(beta2);
460 G4double cost = deltaKinEnergy * (energy + electron_mass_c2) /
461 (deltaMomentum * totMomentum);
462 if(cost > 1.0) { cost = 1.0; }
463 G4double sint = sqrt((1.0 - cost)*(1.0 + cost));
464
465 G4double phi = twopi * G4UniformRand() ;
466
467 deltaDirection.set(sint*cos(phi),sint*sin(phi), cost) ;
468 deltaDirection.rotateUz(direction);
469 }
470 // create G4DynamicParticle object for delta ray
471 auto delta = new G4DynamicParticle(theElectron,deltaDirection,deltaKinEnergy);
472
473 // Change kinematics of primary particle
474 kineticEnergy -= deltaKinEnergy;
475 G4ThreeVector finalP = dp->GetMomentum() - delta->GetMomentum();
476 finalP = finalP.unit();
477
478 fParticleChange->SetProposedKineticEnergy(kineticEnergy);
479 fParticleChange->SetProposedMomentumDirection(finalP);
480
481 vdp->push_back(delta);
482}
#define G4endl
Definition: G4ios.hh:57
G4GLOB_DLL std::ostream G4cout
#define G4UniformRand()
Definition: Randomize.hh:52
Hep3Vector unit() const
void set(double x, double y, double z)
Hep3Vector & rotateUz(const Hep3Vector &)
Definition: ThreeVector.cc:33
const G4ThreeVector & GetMomentumDirection() const
G4double GetKineticEnergy() const
G4ThreeVector GetMomentum() const
const G4Material * GetMaterial() const
void SetProposedKineticEnergy(G4double proposedKinEnergy)
void SetProposedMomentumDirection(const G4ThreeVector &dir)
virtual G4ThreeVector & SampleDirection(const G4DynamicParticle *dp, G4double finalTotalEnergy, G4int Z, const G4Material *)=0
G4int SelectRandomAtomNumber(const G4Material *) const
Definition: G4VEmModel.cc:253
G4double MaxSecondaryKinEnergy(const G4DynamicParticle *dynParticle)
Definition: G4VEmModel.hh:501

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