Loading...
Searching...
No Matches
G4eplusTo3GammaOKVIModel Class Reference
#include <G4eplusTo3GammaOKVIModel.hh>
Inheritance diagram for G4eplusTo3GammaOKVIModel:
Additional Inherited Members | |
| Protected Member Functions inherited from G4VEmModel | |
| G4ParticleChangeForLoss * | GetParticleChangeForLoss () |
| G4ParticleChangeForGamma * | GetParticleChangeForGamma () |
| virtual G4double | MaxSecondaryEnergy (const G4ParticleDefinition *, G4double kineticEnergy) |
| const G4MaterialCutsCouple * | CurrentCouple () const |
| void | SetCurrentElement (const G4Element *) |
| Protected Attributes inherited from G4VEmModel | |
| G4ElementData * | fElementData |
| G4VParticleChange * | pParticleChange |
| G4PhysicsTable * | xSectionTable |
| const G4Material * | pBaseMaterial |
| const std::vector< G4double > * | theDensityFactor |
| const std::vector< G4int > * | theDensityIdx |
| size_t | idxTable |
| G4bool | lossFlucFlag |
| G4double | inveplus |
| G4double | pFactor |
Detailed Description
Definition at line 53 of file G4eplusTo3GammaOKVIModel.hh.
Constructor & Destructor Documentation
◆ G4eplusTo3GammaOKVIModel()
|
explicit |
Definition at line 60 of file G4eplusTo3GammaOKVIModel.cc.
62 : G4VEmModel(nam), fDelta(0.001)
63{
64 theGamma = G4Gamma::Gamma();
65 fParticleChange = nullptr;
66}
Definition: G4VEmModel.hh:103
◆ ~G4eplusTo3GammaOKVIModel()
|
virtual |
Definition at line 70 of file G4eplusTo3GammaOKVIModel.cc.
71{}
Member Function Documentation
◆ ComputeCrossSectionPerAtom()
|
finalvirtual |
Reimplemented from G4VEmModel.
Definition at line 203 of file G4eplusTo3GammaOKVIModel.cc.
207{
208 // Calculates the cross section per atom of annihilation into two photons
209
210
211
213 return cross;
214}
G4double ComputeCrossSectionPerElectron(G4double kinEnergy)
Definition: G4eplusTo3GammaOKVIModel.cc:182
◆ ComputeCrossSectionPerElectron()
Definition at line 182 of file G4eplusTo3GammaOKVIModel.cc.
183{
184 // Calculates the cross section per electron of annihilation into 3 photons
185 // from the Heilter formula.
186
187 G4double ekin = std::max(eV,kinEnergy);
188 G4double tau = ekin/electron_mass_c2;
191 G4double bg2 = tau * (tau+2.0);
192 G4double bg = sqrt(bg2);
193
195 - (gam+3.)/(sqrt(gam*gam - 1.));
196
198 return cross;
199}
Referenced by ComputeCrossSectionPerAtom(), CrossSectionPerVolume(), G4eplusTo2GammaOKVIModel::Initialise(), and G4eplusTo2GammaOKVIModel::SampleSecondaries().
◆ ComputeF()
| G4double G4eplusTo3GammaOKVIModel::ComputeF | ( | G4double | fr1, |
| G4double | fr2, | ||
| G4double | fr3, | ||
| G4double | kinEnergy | ||
| ) |
Definition at line 86 of file G4eplusTo3GammaOKVIModel.cc.
88{
89 G4double ekin = std::max(eV,kinEnergy);
90 G4double tau = ekin/electron_mass_c2;
93 G4double bg2 = tau * (tau+2.0);
94 G4double bg = sqrt(bg2);
95
97 - (gam+3.)/(sqrt(gam*gam - 1.)) + 1.;
98 G4double border;
99
100 if(ekin < 500*MeV) {
101 border = 1. - (electron_mass_c2)/(2*(ekin + electron_mass_c2));
102 } else {
103 border = 1. - (100*electron_mass_c2)/(2*(ekin + electron_mass_c2));
104 }
105
106 border = std::min(border, 0.9999);
107
108 if (fr1>border) { fr1 = border; }
109 if (fr2>border) { fr2 = border; }
110 if (fr3>border) { fr3 = border; }
111
113 G4double fr2s = fr2*fr2;
114 G4double fr3s = fr3*fr3;
115
116 G4double aa = (1.-fr1)*(1.-fr2);
117 G4double ab = fr3s + (fr1-fr2)*(fr1-fr2);
118 G4double add= ((1.-fr1)*(1.-fr1) + (1.-fr2)*(1.-fr2))/(fr3s*aa);
119
120 G4double fres = -rho*(1./fr1s + 1./fr2s)
121 + (ab/(2.*(fr1*fr2*aa)))*(G4Log(2.*gam*aa/(fr1*fr2)))
122 + (ab/(2.*fr1*fr2*(1-fr3)))*G4Log(2.*gam*(1.-fr3)/(fr1*fr2)) - add;
123
124 return fres;
125}
Referenced by ComputeFS().
◆ ComputeF0()
Definition at line 130 of file G4eplusTo3GammaOKVIModel.cc.
132{
133 G4double tau = 0.0;
136 G4double bg2 = tau * (tau+2.0);
137 G4double bg = sqrt(bg2);
138
140 - (gam+3.)/(sqrt(gam*gam - 1.)) + 1.;
141 G4double border = 0.5;
142
143 if (fr1>border) { fr1 = border; }
144 if (fr2>border) { fr2 = border; }
145 if (fr3>border) { fr3 = border; }
146
148 G4double fr2s = fr2*fr2;
149 G4double fr3s = fr3*fr3;
150
151 G4double aa = (1.-fr1)*(1.-fr2);
152 G4double ab = fr3s + (fr1-fr2)*(fr1-fr2);
153 G4double add= ((1.-fr1)*(1.-fr1) + (1.-fr2)*(1.-fr2))/(fr3s*aa);
154
155 G4double fres = -rho*(1./fr1s + 1./fr2s)
156 + (ab/(2.*(fr1*fr2*aa)))*(G4Log(2.*gam*aa/(fr1*fr2)))
157 + (ab/(2.*fr1*fr2*(1-fr3)))*G4Log(2.*gam*(1.-fr3)/(fr1*fr2)) - add;
158
159 return fres;
160}
◆ ComputeFS()
| G4double G4eplusTo3GammaOKVIModel::ComputeFS | ( | G4double | fr1, |
| G4double | fr2, | ||
| G4double | fr3, | ||
| G4double | kinEnergy | ||
| ) |
Definition at line 163 of file G4eplusTo3GammaOKVIModel.cc.
165{
166 G4double ekin = std::max(eV,kinEnergy);
167 G4double tau = ekin/electron_mass_c2;
169
171 ComputeF(fr3,fr1,fr2,ekin) +
172 ComputeF(fr2,fr3,fr1,ekin));
173
175
176 return dcross;
177}
G4double ComputeF(G4double fr1, G4double fr2, G4double fr3, G4double kinEnergy)
Definition: G4eplusTo3GammaOKVIModel.cc:86
Referenced by SampleSecondaries().
◆ CrossSectionPerVolume()
|
finalvirtual |
Reimplemented from G4VEmModel.
Definition at line 218 of file G4eplusTo3GammaOKVIModel.cc.
223{
224 // Calculates the cross section per volume of annihilation into two photons
225
228 return cross;
229}
◆ Initialise()
|
finalvirtual |
Implements G4VEmModel.
Definition at line 75 of file G4eplusTo3GammaOKVIModel.cc.
77{
78 // here particle change is set for the triplet model
79 if(fParticleChange) { return; }
80 fParticleChange = GetParticleChangeForGamma();
81}
G4ParticleChangeForGamma * GetParticleChangeForGamma()
Definition: G4VEmModel.cc:133
Referenced by G4eplusTo2GammaOKVIModel::Initialise().
◆ SampleSecondaries()
|
finalvirtual |
Implements G4VEmModel.
Definition at line 237 of file G4eplusTo3GammaOKVIModel.cc.
241{
242
244 G4DynamicParticle *aGamma1, *aGamma2;
246 G4double border;
247
248 if(posiKinEnergy < 500*MeV) {
249 border = 1. - (electron_mass_c2)/(2*(posiKinEnergy + electron_mass_c2));
250 } else {
251 border = 1. - (100*electron_mass_c2)/(2*(posiKinEnergy + electron_mass_c2));
252 }
253 border = std::min(border, 0.9999);
254
255 CLHEP::HepRandomEngine* rndmEngine = G4Random::getTheEngine();
256
257 // Case at rest
258 if(posiKinEnergy == 0.0) {
260 G4double sint = sqrt((1. - cost)*(1. + cost));
262 G4ThreeVector dir(sint*cos(phi), sint*sin(phi), cost);
263 phi = twopi * rndmEngine->flat();
264 G4double cosphi = cos(phi);
265 G4double sinphi = sin(phi);
266 G4ThreeVector pol(cosphi, sinphi, 0.0);
267 pol.rotateUz(dir);
269 aGamma1->SetPolarization(pol.x(),pol.y(),pol.z());
271 pol.set(-sinphi, cosphi, 0.0);
272 pol.rotateUz(dir);
273 aGamma2->SetPolarization(pol.x(),pol.y(),pol.z());
274
275 } else {
276
278
279 // (A.A.) LIMITS FOR 1st GAMMA
280 G4double xmin = 0.01;
282
283 G4double d1, d0, x1, x2, dmax, x2min;
284
285 // (A.A.) sampling of x1 x2 x3 (whole cycle of rejection)
286 do {
287 x1 = 1/((1/xmin) - ((1/xmin)-(1/xmax))*rndmEngine->flat());
288 dmax = ComputeFS(posiKinEnergy, x1,1.-x1,border);
289 x2min = 1.-x1;
290 x2 = 1 - rndmEngine->flat()*(1-x2min);
291 d1 = dmax*rndmEngine->flat();
292 d0 = ComputeFS(posiKinEnergy,x1,x2,2-x1-x2);
293 }
294 while(d0 < d1);
295
296 G4double x3 = 2 - x1 - x2;
297
298 //
299 // angles between Gammas
300 //
301
302 G4double psi13 = 2*asin(sqrt(std::abs((x1+x3-1)/(x1*x3))));
303 G4double psi12 = 2*asin(sqrt(std::abs((x1+x2-1)/(x1*x2))));
304
305 // sin^t
306
307 //G4double phi = twopi * rndmEngine->flat();
308 //G4double psi = acos(x3); // Angle of the plane
309
310 //
311 // kinematic of the created pair
312 //
313
314 G4double TotalAvailableEnergy = posiKinEnergy + 2.0*electron_mass_c2;
315
316 G4double phot1Energy = 0.5*x1*TotalAvailableEnergy;
317 G4double phot2Energy = 0.5*x2*TotalAvailableEnergy;
318 G4double phot3Energy = 0.5*x3*TotalAvailableEnergy;
319
320
321 // DIRECTIONS
322
323 // The azimuthal angles of ql and q3 with respect to some plane
324 // through the beam axis are generated at random.
325
326 G4ThreeVector phot1Direction(0, 0, 1);
327 G4ThreeVector phot2Direction(0, sin(psi12), cos(psi12));
328 G4ThreeVector phot3Direction(0, sin(psi13), cos(psi13));
329
330 phot1Direction.rotateUz(posiDirection);
331 phot2Direction.rotateUz(posiDirection);
332 phot3Direction.rotateUz(posiDirection);
333
337
338
339 //POLARIZATION - ???
340 /*
341
342
343 phi = twopi * rndmEngine->flat();
344 G4double cosphi = cos(phi);
345 G4double sinphi = sin(phi);
346 G4ThreeVector pol(cosphi, sinphi, 0.0);
347 pol.rotateUz(phot1Direction);
348 aGamma1->SetPolarization(pol.x(),pol.y(),pol.z());
349
350 G4double phot2Energy =(1.-epsil)*TotalAvailableEnergy;
351 G4double posiP= sqrt(posiKinEnergy*(posiKinEnergy+2.*electron_mass_c2));
352 G4ThreeVector dir = posiDirection*posiP - phot1Direction*phot1Energy;
353 G4ThreeVector phot2Direction = dir.unit();
354
355 // create G4DynamicParticle object for the particle2
356 aGamma2 = new G4DynamicParticle (theGamma,phot2Direction, phot2Energy);
357
358 //!!! likely problematic direction to be checked
359 pol.set(-sinphi, cosphi, 0.0);
360 pol.rotateUz(phot1Direction);
361 cost = pol*phot2Direction;
362 pol -= cost*phot2Direction;
363 pol = pol.unit();
364 aGamma2->SetPolarization(pol.x(),pol.y(),pol.z());
365
366 */
367
368 }
369 /*
370 G4cout << "Annihilation in fly: e0= " << posiKinEnergy
371 << " m= " << electron_mass_c2
372 << " e1= " << phot1Energy
373 << " e2= " << phot2Energy << " dir= " << dir
374 << " -> " << phot1Direction << " "
375 << phot2Direction << G4endl;
376 */
377
378 vdp->push_back(aGamma1);
379 vdp->push_back(aGamma2);
380 if(aGamma3 != nullptr) { vdp->push_back(aGamma3); }
381
382 // kill primary positron
383 fParticleChange->SetProposedKineticEnergy(0.0);
385}
Definition: ThreeVector.h:36
Definition: RandomEngine.h:53
virtual double flat()=0
Definition: G4DynamicParticle.hh:65
void SetPolarization(const G4ThreeVector &)
const G4ThreeVector & GetMomentumDirection() const
G4double GetKineticEnergy() const
void SetProposedKineticEnergy(G4double proposedKinEnergy)
Definition: G4ParticleChangeForGamma.hh:119
void ProposeTrackStatus(G4TrackStatus status)
G4double ComputeFS(G4double fr1, G4double fr2, G4double fr3, G4double kinEnergy)
Definition: G4eplusTo3GammaOKVIModel.cc:163
Referenced by G4eplusTo2GammaOKVIModel::SampleSecondaries().
◆ SetDelta()
|
inline |
Definition at line 96 of file G4eplusTo3GammaOKVIModel.hh.
96{ if(val > 0.0) { fDelta = val; } };
Referenced by G4eplusTo2GammaOKVIModel::Initialise(), and G4eplusTo2GammaOKVIModel::SampleSecondaries().
The documentation for this class was generated from the following files:
