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

#include <G4eplusTo2GammaOKVIModel.hh>

+ Inheritance diagram for G4eplusTo2GammaOKVIModel:

Public Member Functions

 G4eplusTo2GammaOKVIModel (const G4ParticleDefinition *p=nullptr, const G4String &nam="eplus2ggOKVI")
 
 ~G4eplusTo2GammaOKVIModel () override
 
void Initialise (const G4ParticleDefinition *, const G4DataVector &) override
 
G4double ComputeCrossSectionPerElectron (G4double kinEnergy)
 
G4double ComputeCrossSectionPerAtom (const G4ParticleDefinition *, G4double kinEnergy, G4double Z, G4double A=0., G4double cutEnergy=0., G4double maxEnergy=DBL_MAX) final
 
G4double CrossSectionPerVolume (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX) final
 
void SampleSecondaries (std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin=0.0, G4double maxEnergy=DBL_MAX) final
 
void SetDelta (G4double val)
 
G4eplusTo2GammaOKVIModeloperator= (const G4eplusTo2GammaOKVIModel &right)=delete
 
 G4eplusTo2GammaOKVIModel (const G4eplusTo2GammaOKVIModel &)=delete
 
- Public Member Functions inherited from G4VEmModel
 G4VEmModel (const G4String &nam)
 
virtual ~G4VEmModel ()
 
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 GetPartialCrossSection (const G4Material *, G4int level, const G4ParticleDefinition *, G4double kineticEnergy)
 
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 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 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)
 
void SetLPMFlag (G4bool)
 
G4VEmModeloperator= (const G4VEmModel &right)=delete
 
 G4VEmModel (const G4VEmModel &)=delete
 

Additional Inherited Members

- Protected Member Functions inherited from G4VEmModel
G4ParticleChangeForLossGetParticleChangeForLoss ()
 
G4ParticleChangeForGammaGetParticleChangeForGamma ()
 
virtual G4double MaxSecondaryEnergy (const G4ParticleDefinition *, G4double kineticEnergy)
 
const G4MaterialCutsCoupleCurrentCouple () const
 
void SetCurrentElement (const G4Element *)
 
- 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
 
std::size_t currentCoupleIndex = 0
 
std::size_t basedCoupleIndex = 0
 
G4bool lossFlucFlag = true
 

Detailed Description

Definition at line 56 of file G4eplusTo2GammaOKVIModel.hh.

Constructor & Destructor Documentation

◆ G4eplusTo2GammaOKVIModel() [1/2]

G4eplusTo2GammaOKVIModel::G4eplusTo2GammaOKVIModel ( const G4ParticleDefinition * p = nullptr,
const G4String & nam = "eplus2ggOKVI" )
explicit

Definition at line 69 of file G4eplusTo2GammaOKVIModel.cc.

71 : G4VEmModel(nam),
72 fDelta(0.001),
73 fGammaTh(MeV)
74{
75 theGamma = G4Gamma::Gamma();
76 fParticleChange = nullptr;
77 fCuts = nullptr;
78 f3GModel = new G4eplusTo3GammaOKVIModel();
79 SetTripletModel(f3GModel);
80}
static G4Gamma * Gamma()
Definition G4Gamma.cc:81
void SetTripletModel(G4VEmModel *)
G4VEmModel(const G4String &nam)
Definition G4VEmModel.cc:67

◆ ~G4eplusTo2GammaOKVIModel()

G4eplusTo2GammaOKVIModel::~G4eplusTo2GammaOKVIModel ( )
overridedefault

◆ G4eplusTo2GammaOKVIModel() [2/2]

G4eplusTo2GammaOKVIModel::G4eplusTo2GammaOKVIModel ( const G4eplusTo2GammaOKVIModel & )
delete

Member Function Documentation

◆ ComputeCrossSectionPerAtom()

G4double G4eplusTo2GammaOKVIModel::ComputeCrossSectionPerAtom ( const G4ParticleDefinition * ,
G4double kinEnergy,
G4double Z,
G4double A = 0.,
G4double cutEnergy = 0.,
G4double maxEnergy = DBL_MAX )
finalvirtual

Reimplemented from G4VEmModel.

Definition at line 149 of file G4eplusTo2GammaOKVIModel.cc.

153{
154 // Calculates the cross section per atom of annihilation into two photons
155 G4double cross = Z*fCrossSection->Value(kineticEnergy);
156 return cross;
157}
double G4double
Definition G4Types.hh:83
G4double Value(const G4double energy, std::size_t &lastidx) const

◆ ComputeCrossSectionPerElectron()

G4double G4eplusTo2GammaOKVIModel::ComputeCrossSectionPerElectron ( G4double kinEnergy)

Definition at line 126 of file G4eplusTo2GammaOKVIModel.cc.

127{
128 // Calculates the cross section per electron of annihilation into two
129 // photons from the Heilter formula with the radiation correction to 3 gamma
130 // annihilation channel. (A.A.) rho is changed
131
132 G4double ekin = std::max(eV,kinEnergy);
133 G4double tau = ekin/electron_mass_c2;
134 G4double gam = tau + 1.0;
135 G4double gamma2 = gam*gam;
136 G4double bg2 = tau * (tau+2.0);
137 G4double bg = sqrt(bg2);
138 G4double rho = (gamma2+4.*gam+1.)*G4Log(gam+bg)/(gamma2-1.)
139 - (gam+3.)/(sqrt(gam*gam - 1.));
140
141 static const G4double pir2 = pi*classic_electr_radius*classic_electr_radius;
142 G4double cross = (pir2*rho + alpha_rcl2*2.*G4Log(fDelta)*rho*rho)/(gam+1.);
143
144 return cross;
145}
G4double G4Log(G4double x)
Definition G4Log.hh:227

Referenced by Initialise().

◆ CrossSectionPerVolume()

G4double G4eplusTo2GammaOKVIModel::CrossSectionPerVolume ( const G4Material * material,
const G4ParticleDefinition * ,
G4double kineticEnergy,
G4double cutEnergy = 0.0,
G4double maxEnergy = DBL_MAX )
finalvirtual

Reimplemented from G4VEmModel.

Definition at line 161 of file G4eplusTo2GammaOKVIModel.cc.

166{
167 // Calculates the cross section per volume of annihilation into two photons
168 G4double eDensity = material->GetElectronDensity();
169 G4double cross = eDensity*fCrossSection->Value(kineticEnergy);
170 return cross;
171}
G4double GetElectronDensity() const

◆ Initialise()

void G4eplusTo2GammaOKVIModel::Initialise ( const G4ParticleDefinition * p,
const G4DataVector & cuts )
overridevirtual

Implements G4VEmModel.

Definition at line 88 of file G4eplusTo2GammaOKVIModel.cc.

90{
91 f3GModel->Initialise(p, cuts);
92 fCuts = &cuts;
94 f3GModel->SetDelta(fDelta);
95
96 if(IsMaster()) {
97 if(!fCrossSection) {
98 G4double emin = 10*eV;
99 G4double emax = 100*TeV;
100 G4int nbins = 20*G4lrint(std::log10(emax/emin));
101 fCrossSection = new G4PhysicsLogVector(emin, emax, nbins, true);
102 fCrossSection3G = new G4PhysicsLogVector(emin, emax, nbins, true);
103 f3GProbability = new G4PhysicsLogVector(emin, emax, nbins, true);
104 for(G4int i=0; i<= nbins; ++i) {
105 G4double e = fCrossSection->Energy(i);
107 G4double cs3 = f3GModel->ComputeCrossSectionPerElectron(e);
108 cs2 += cs3;
109 fCrossSection->PutValue(i, cs2);
110 fCrossSection3G->PutValue(i, cs3);
111 f3GProbability->PutValue(i, cs3/cs2);
112 }
113 fCrossSection->FillSecondDerivatives();
114 fCrossSection3G->FillSecondDerivatives();
115 f3GProbability->FillSecondDerivatives();
116 }
117 }
118 // here particle change is set for the triplet model
119 if(fParticleChange) { return; }
120 fParticleChange = GetParticleChangeForGamma();
121}
int G4int
Definition G4Types.hh:85
static G4EmParameters * Instance()
G4double LowestTripletEnergy() const
void PutValue(const std::size_t index, const G4double value)
G4double Energy(const std::size_t index) const
void FillSecondDerivatives(const G4SplineType=G4SplineType::Base, const G4double dir1=0.0, const G4double dir2=0.0)
G4ParticleChangeForGamma * GetParticleChangeForGamma()
G4bool IsMaster() const
G4double ComputeCrossSectionPerElectron(G4double kinEnergy)
G4double ComputeCrossSectionPerElectron(G4double kinEnergy)
void Initialise(const G4ParticleDefinition *, const G4DataVector &) final
int G4lrint(double ad)
Definition templates.hh:134

◆ operator=()

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

◆ SampleSecondaries()

void G4eplusTo2GammaOKVIModel::SampleSecondaries ( std::vector< G4DynamicParticle * > * vdp,
const G4MaterialCutsCouple * mcc,
const G4DynamicParticle * dp,
G4double tmin = 0.0,
G4double maxEnergy = DBL_MAX )
finalvirtual

!! likely problematic direction to be checked

Implements G4VEmModel.

Definition at line 179 of file G4eplusTo2GammaOKVIModel.cc.

183{
184 G4double posiKinEnergy = dp->GetKineticEnergy();
185 CLHEP::HepRandomEngine* rndmEngine = G4Random::getTheEngine();
186
187 if(rndmEngine->flat() < f3GProbability->Value(posiKinEnergy)) {
188 G4double cutd = std::max(fGammaTh,(*fCuts)[mcc->GetIndex()])
189 /(posiKinEnergy + electron_mass_c2);
190 // check cut to avoid production of 3d gamma below
191 if(cutd > fDelta) {
192 G4double cs30 = fCrossSection3G->Value(posiKinEnergy);
193 f3GModel->SetDelta(cutd);
194 G4double cs3 = f3GModel->ComputeCrossSectionPerElectron(posiKinEnergy);
195 if(rndmEngine->flat()*cs30 < cs3) {
196 f3GModel->SampleSecondaries(vdp, mcc, dp);
197 return;
198 }
199 } else {
200 f3GModel->SampleSecondaries(vdp, mcc, dp);
201 return;
202 }
203 }
204
205 G4DynamicParticle *aGamma1, *aGamma2;
206
207 // Case at rest
208 if(posiKinEnergy == 0.0) {
209 G4double cost = 2.*rndmEngine->flat()-1.;
210 G4double sint = sqrt((1. - cost)*(1. + cost));
211 G4double phi = twopi * rndmEngine->flat();
212 G4ThreeVector dir(sint*cos(phi), sint*sin(phi), cost);
213 phi = twopi * rndmEngine->flat();
214 G4double cosphi = cos(phi);
215 G4double sinphi = sin(phi);
216 G4ThreeVector pol(cosphi, sinphi, 0.0);
217 pol.rotateUz(dir);
218 aGamma1 = new G4DynamicParticle(theGamma, dir, electron_mass_c2);
219 aGamma1->SetPolarization(pol.x(),pol.y(),pol.z());
220 aGamma2 = new G4DynamicParticle(theGamma,-dir, electron_mass_c2);
221 pol.set(-sinphi, cosphi, 0.0);
222 pol.rotateUz(dir);
223 aGamma2->SetPolarization(pol.x(),pol.y(),pol.z());
224
225 } else {
226
227 G4ThreeVector posiDirection = dp->GetMomentumDirection();
228
229 G4double tau = posiKinEnergy/electron_mass_c2;
230 G4double gam = tau + 1.0;
231 G4double tau2 = tau + 2.0;
232 G4double sqgrate = sqrt(tau/tau2)*0.5;
233 G4double sqg2m1 = sqrt(tau*tau2);
234
235 // limits of the energy sampling
236 G4double epsilmin = 0.5 - sqgrate;
237 G4double epsilmax = 0.5 + sqgrate;
238 G4double epsilqot = epsilmax/epsilmin;
239
240 //
241 // sample the energy rate of the created gammas
242 //
243 G4double epsil, greject;
244
245 do {
246 epsil = epsilmin*G4Exp(G4Log(epsilqot)*rndmEngine->flat());
247 greject = 1. - epsil + (2.*gam*epsil-1.)/(epsil*tau2*tau2);
248 // Loop checking, 03-Aug-2015, Vladimir Ivanchenko
249 } while( greject < rndmEngine->flat());
250
251 //
252 // scattered Gamma angles. ( Z - axis along the parent positron)
253 //
254
255 G4double cost = (epsil*tau2-1.)/(epsil*sqg2m1);
256 if(std::abs(cost) > 1.0) {
257 G4cout << "### G4eplusTo2GammaOKVIModel WARNING cost= " << cost
258 << " positron Ekin(MeV)= " << posiKinEnergy
259 << " gamma epsil= " << epsil
260 << G4endl;
261 if(cost > 1.0) cost = 1.0;
262 else cost = -1.0;
263 }
264 G4double sint = sqrt((1.+cost)*(1.-cost));
265 G4double phi = twopi * rndmEngine->flat();
266
267 //
268 // kinematic of the created pair
269 //
270
271 G4double TotalAvailableEnergy = posiKinEnergy + 2.0*electron_mass_c2;
272 G4double phot1Energy = epsil*TotalAvailableEnergy;
273
274 G4ThreeVector phot1Direction(sint*cos(phi), sint*sin(phi), cost);
275 phot1Direction.rotateUz(posiDirection);
276 aGamma1 = new G4DynamicParticle (theGamma,phot1Direction, phot1Energy);
277 phi = twopi * rndmEngine->flat();
278 G4double cosphi = cos(phi);
279 G4double sinphi = sin(phi);
280 G4ThreeVector pol(cosphi, sinphi, 0.0);
281 pol.rotateUz(phot1Direction);
282 aGamma1->SetPolarization(pol.x(),pol.y(),pol.z());
283
284 G4double phot2Energy =(1.-epsil)*TotalAvailableEnergy;
285 G4double posiP= sqrt(posiKinEnergy*(posiKinEnergy+2.*electron_mass_c2));
286 G4ThreeVector dir = posiDirection*posiP - phot1Direction*phot1Energy;
287 G4ThreeVector phot2Direction = dir.unit();
288
289 // create G4DynamicParticle object for the particle2
290 aGamma2 = new G4DynamicParticle (theGamma,phot2Direction, phot2Energy);
291
292 //!!! likely problematic direction to be checked
293 pol.set(-sinphi, cosphi, 0.0);
294 pol.rotateUz(phot1Direction);
295 cost = pol*phot2Direction;
296 pol -= cost*phot2Direction;
297 pol = pol.unit();
298 aGamma2->SetPolarization(pol.x(),pol.y(),pol.z());
299 }
300 /*
301 G4cout << "Annihilation in fly: e0= " << posiKinEnergy
302 << " m= " << electron_mass_c2
303 << " e1= " << phot1Energy
304 << " e2= " << phot2Energy << " dir= " << dir
305 << " -> " << phot1Direction << " "
306 << phot2Direction << G4endl;
307 */
308
309 vdp->push_back(aGamma1);
310 vdp->push_back(aGamma2);
311
312 // kill primary positron
313 fParticleChange->SetProposedKineticEnergy(0.0);
314 fParticleChange->ProposeTrackStatus(fStopAndKill);
315}
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
Definition G4Exp.hh:180
@ fStopAndKill
#define G4endl
Definition G4ios.hh:67
G4GLOB_DLL std::ostream G4cout
Hep3Vector unit() const
virtual double flat()=0
void SetPolarization(const G4ThreeVector &)
const G4ThreeVector & GetMomentumDirection() const
G4double GetKineticEnergy() const
void SetProposedKineticEnergy(G4double proposedKinEnergy)
void ProposeTrackStatus(G4TrackStatus status)
void SampleSecondaries(std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin=0.0, G4double maxEnergy=DBL_MAX) final
double flat()

◆ SetDelta()

void G4eplusTo2GammaOKVIModel::SetDelta ( G4double val)
inline

Definition at line 90 of file G4eplusTo2GammaOKVIModel.hh.

90{ if(val > 0.0) { fDelta = val; } };

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