264{
265
266
267
268
269
270 if (verboseLevel > 3)
271 G4cout <<
"Calling SampleSecondaries() of G4LivermorePolarizedGammaConversionModel" <<
G4endl;
272
274
275 if(photonEnergy <= lowEnergyLimit)
276 {
279 return;
280 }
281
284
285
286
287 if(!(gammaPolarization0.
isOrthogonal(gammaDirection0, 1e-6))||(gammaPolarization0.
mag()==0))
288 {
289 gammaPolarization0 = GetRandomPolarization(gammaDirection0);
290 }
291 else
292 {
294 {
295 gammaPolarization0 = GetPerpendicularPolarization(gammaDirection0, gammaPolarization0);
296 }
297 }
298
299
300
302 G4double epsilon0Local = electron_mass_c2 / photonEnergy ;
303
304
305
306 if (photonEnergy < smallEnergy )
307 {
309 }
310 else
311 {
312
315
316 if (element == nullptr)
317 {
318 G4cout <<
"G4LivermorePolarizedGammaConversionModel::SampleSecondaries - element = 0" <<
G4endl;
319 return;
320 }
321
322
324 if (ionisation == nullptr)
325 {
326 G4cout <<
"G4LivermorePolarizedGammaConversionModel::SampleSecondaries - ionisation = 0" <<
G4endl;
327 return;
328 }
329
330
332 if (photonEnergy > 50. * MeV) fZ += 8. * (element->
GetfCoulomb());
333
334
337 G4double screenMin = std::min(4.*screenFactor,screenMax) ;
338
339
340 G4double epsilon1 = 0.5 - 0.5 * sqrt(1. - screenMin / screenMax) ;
341 G4double epsilonMin = std::max(epsilon0Local,epsilon1);
342 G4double epsilonRange = 0.5 - epsilonMin ;
343
344
347
348 G4double f10 = ScreenFunction1(screenMin) - fZ;
349 G4double f20 = ScreenFunction2(screenMin) - fZ;
350 G4double normF1 = std::max(f10 * epsilonRange * epsilonRange,0.);
351 G4double normF2 = std::max(1.5 * f20,0.);
352
353 do {
355 {
358 gReject = (ScreenFunction1(screen) - fZ) / f10 ;
359 }
360 else
361 {
364 gReject = (ScreenFunction2(screen) - fZ) / f20 ;
365 }
367 }
368
369
372
374 {
375 electronTotEnergy = (1. -
epsilon) * photonEnergy;
376 positronTotEnergy =
epsilon * photonEnergy;
377 }
378 else
379 {
380 positronTotEnergy = (1. -
epsilon) * photonEnergy;
381 electronTotEnergy =
epsilon * photonEnergy;
382 }
383
384
385
386
387 G4double Ene = electronTotEnergy/electron_mass_c2;
388
391
392 SetTheta(&cosTheta,&sinTheta,Ene);
394
395
396 phi = SetPhi(photonEnergy);
397 psi = SetPsi(photonEnergy,phi);
398 Psi = psi;
399 Phi = phi;
400
405
406 if (choice2 <= 0.5)
407 {
408
409
410 }
411 else
412 {
413 phi = -phi;
414 }
415
416 if (choice <= 0.5)
417 {
418 phie = psi;
419 phip = phie+phi;
420 }
421 else
422 {
423
424 phip = psi;
425 phie = phip + phi;
426 }
427
428
429
434
435
436
437
438
439 G4double electronKineEnergy = std::max(0.,electronTotEnergy - electron_mass_c2) ;
440
441 SystemOfRefChange(gammaDirection0,electronDirection,
442 gammaPolarization0);
443
445 electronDirection,
446 electronKineEnergy);
447
448
449 Ene = positronTotEnergy/electron_mass_c2;
450
451 cosTheta = 0.;
452 sinTheta = 0.;
453
454 SetTheta(&cosTheta,&sinTheta,Ene);
455
456
457 dirX = sinTheta*cos(phip);
458 dirY = sinTheta*sin(phip);
459 dirZ = cosTheta;
461
462 G4double positronKineEnergy = std::max(0.,positronTotEnergy - electron_mass_c2) ;
463 SystemOfRefChange(gammaDirection0,positronDirection,
464 gammaPolarization0);
465
466
468 positronDirection, positronKineEnergy);
469 fvect->push_back(particle1);
470 fvect->push_back(particle2);
471
472
475}
G4double epsilon(G4double density, G4double temperature)
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
bool isOrthogonal(const Hep3Vector &v, double epsilon=tolerance) const
double howOrthogonal(const Hep3Vector &v) const
const G4ThreeVector & GetMomentumDirection() const
G4ParticleDefinition * GetDefinition() const
G4double GetKineticEnergy() const
const G4ThreeVector & GetPolarization() const
static G4Electron * Electron()
G4double GetfCoulomb() const
G4IonisParamElm * GetIonisation() const
G4double GetlogZ3() const
void SetProposedKineticEnergy(G4double proposedKinEnergy)
static G4Positron * Positron()
const G4Element * SelectRandomAtom(const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
void ProposeTrackStatus(G4TrackStatus status)