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

#include <G4PenelopeBremsstrahlungAngular.hh>

+ Inheritance diagram for G4PenelopeBremsstrahlungAngular:

Public Member Functions

 G4PenelopeBremsstrahlungAngular ()
 
 ~G4PenelopeBremsstrahlungAngular ()
 
void Initialize ()
 The Initialize() method forces the cleaning of tables.
 
G4double PolarAngle (const G4double initial_energy, const G4double final_energy, const G4int Z)
 
G4ThreeVectorSampleDirection (const G4DynamicParticle *dp, G4double out_energy, G4int Z, const G4Material *mat=0)
 
void SetVerbosityLevel (G4int vl)
 Set/Get Verbosity level.
 
G4int GetVerbosityLevel ()
 
- Public Member Functions inherited from G4VEmAngularDistribution
 G4VEmAngularDistribution (const G4String &name)
 
virtual ~G4VEmAngularDistribution ()
 
virtual G4ThreeVectorSampleDirection (const G4DynamicParticle *dp, G4double finalTotalEnergy, G4int Z, const G4Material *)=0
 
const G4StringGetName () const
 

Additional Inherited Members

- Protected Attributes inherited from G4VEmAngularDistribution
G4ThreeVector fLocalDirection
 

Detailed Description

Definition at line 55 of file G4PenelopeBremsstrahlungAngular.hh.

Constructor & Destructor Documentation

◆ G4PenelopeBremsstrahlungAngular()

G4PenelopeBremsstrahlungAngular::G4PenelopeBremsstrahlungAngular ( )

Definition at line 57 of file G4PenelopeBremsstrahlungAngular.cc.

57 :
58 G4VEmAngularDistribution("Penelope"), theEffectiveZSq(0),
59 theLorentzTables1(0),theLorentzTables2(0)
60
61{
62 dataRead = false;
63 verbosityLevel = 0;
64}

◆ ~G4PenelopeBremsstrahlungAngular()

G4PenelopeBremsstrahlungAngular::~G4PenelopeBremsstrahlungAngular ( )

Definition at line 68 of file G4PenelopeBremsstrahlungAngular.cc.

69{
70 ClearTables();
71}

Member Function Documentation

◆ GetVerbosityLevel()

G4int G4PenelopeBremsstrahlungAngular::GetVerbosityLevel ( )
inline

Definition at line 80 of file G4PenelopeBremsstrahlungAngular.hh.

80{return verbosityLevel;};

◆ Initialize()

void G4PenelopeBremsstrahlungAngular::Initialize ( )

The Initialize() method forces the cleaning of tables.

Definition at line 75 of file G4PenelopeBremsstrahlungAngular.cc.

76{
77 ClearTables();
78}

Referenced by G4PenelopeBremsstrahlungModel::Initialise().

◆ PolarAngle()

G4double G4PenelopeBremsstrahlungAngular::PolarAngle ( const G4double  initial_energy,
const G4double  final_energy,
const G4int  Z 
)

Old interface, backwards compatibility. Will not work in this case it will produce a G4Exception().

Definition at line 439 of file G4PenelopeBremsstrahlungAngular.cc.

442{
443 G4cout << "WARNING: G4PenelopeBremsstrahlungAngular() does NOT support PolarAngle()" << G4endl;
444 G4cout << "Please use the alternative interface SampleDirection()" << G4endl;
445 G4Exception("G4PenelopeBremsstrahlungAngular::PolarAngle()",
446 "em0005",FatalException,"Unsupported interface");
447 return 0;
448}
@ FatalException
#define G4endl
Definition: G4ios.hh:52
G4DLLIMPORT std::ostream G4cout
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41

◆ SampleDirection()

G4ThreeVector & G4PenelopeBremsstrahlungAngular::SampleDirection ( const G4DynamicParticle dp,
G4double  out_energy,
G4int  Z,
const G4Material mat = 0 
)
virtual

Samples the direction of the outgoing photon (in global coordinates). Forces the calculation of tables, if they are not available

Implements G4VEmAngularDistribution.

Definition at line 309 of file G4PenelopeBremsstrahlungAngular.cc.

313{
314 if (!material)
315 {
316 G4Exception("G4PenelopeBremsstrahlungAngular::SampleDirection()",
317 "em2040",FatalException,"The pointer to G4Material* is NULL");
318 return fLocalDirection;
319 }
320
321 G4double Zmat = GetEffectiveZ(material);
322 if (verbosityLevel > 0)
323 {
324 G4cout << "Effective <Z> for material : " << material->GetName() <<
325 " = " << Zmat << G4endl;
326 }
327
328 G4double ePrimary = dp->GetKineticEnergy();
329
330 G4double beta = std::sqrt(ePrimary*(ePrimary+2*electron_mass_c2))/
331 (ePrimary+electron_mass_c2);
332 G4double cdt = 0;
333 G4double sinTheta = 0;
334 G4double phi = 0;
335
336 //Use a pure dipole distribution for energy above 500 keV
337 if (ePrimary > 500*keV)
338 {
339 cdt = 2.0*G4UniformRand() - 1.0;
340 if (G4UniformRand() > 0.75)
341 {
342 if (cdt<0)
343 cdt = -1.0*std::pow(-cdt,1./3.);
344 else
345 cdt = std::pow(cdt,1./3.);
346 }
347 cdt = (cdt+beta)/(1.0+beta*cdt);
348 //Get primary kinematics
349 sinTheta = std::sqrt(1. - cdt*cdt);
350 phi = twopi * G4UniformRand();
351 fLocalDirection.set(sinTheta* std::cos(phi),
352 sinTheta* std::sin(phi),
353 cdt);
354 //rotate
356 //return
357 return fLocalDirection;
358 }
359
360 //Else, retrieve tables and go through the full thing
361 if (!theLorentzTables1)
362 theLorentzTables1 = new std::map<G4double,G4PhysicsTable*>;
363 if (!theLorentzTables2)
364 theLorentzTables2 = new std::map<G4double,G4PhysicsTable*>;
365
366 //Check if tables exist for the given Zmat
367 if (!(theLorentzTables1->count(Zmat)))
368 PrepareInterpolationTables(Zmat);
369
370 if (!(theLorentzTables1->count(Zmat)) || !(theLorentzTables2->count(Zmat)))
371 {
373 ed << "Unable to retrieve Lorentz tables for Z= " << Zmat << G4endl;
374 G4Exception("G4PenelopeBremsstrahlungAngular::SampleDirection()",
375 "em2006",FatalException,ed);
376 }
377
378 //retrieve actual tables
379 G4PhysicsTable* theTable1 = theLorentzTables1->find(Zmat)->second;
380 G4PhysicsTable* theTable2 = theLorentzTables2->find(Zmat)->second;
381
382 G4double RK=20.0*eGamma/ePrimary;
383 G4int ik=std::min((G4int) RK,19);
384
385 G4double P10=0,P11=0,P1=0;
386 G4double P20=0,P21=0,P2=0;
387
388 //First coefficient
389 G4PhysicsFreeVector* v1 = (G4PhysicsFreeVector*) (*theTable1)[ik];
390 G4PhysicsFreeVector* v2 = (G4PhysicsFreeVector*) (*theTable1)[ik+1];
391 P10 = v1->Value(beta);
392 P11 = v2->Value(beta);
393 P1=P10+(RK-(G4double) ik)*(P11-P10);
394
395 //Second coefficient
396 G4PhysicsFreeVector* v3 = (G4PhysicsFreeVector*) (*theTable2)[ik];
397 G4PhysicsFreeVector* v4 = (G4PhysicsFreeVector*) (*theTable2)[ik+1];
398 P20=v3->Value(beta);
399 P21=v4->Value(beta);
400 P2=P20+(RK-(G4double) ik)*(P21-P20);
401
402 //Sampling from the Lorenz-trasformed dipole distributions
403 P1=std::min(std::exp(P1)/beta,1.0);
404 G4double betap = std::min(std::max(beta*(1.0+P2/beta),0.0),0.9999);
405
406 G4double testf=0;
407
408 if (G4UniformRand() < P1)
409 {
410 do{
411 cdt = 2.0*G4UniformRand()-1.0;
412 testf=2.0*G4UniformRand()-(1.0+cdt*cdt);
413 }while(testf>0);
414 }
415 else
416 {
417 do{
418 cdt = 2.0*G4UniformRand()-1.0;
419 testf=G4UniformRand()-(1.0-cdt*cdt);
420 }while(testf>0);
421 }
422 cdt = (cdt+betap)/(1.0+betap*cdt);
423
424 //Get primary kinematics
425 sinTheta = std::sqrt(1. - cdt*cdt);
426 phi = twopi * G4UniformRand();
427 fLocalDirection.set(sinTheta* std::cos(phi),
428 sinTheta* std::sin(phi),
429 cdt);
430 //rotate
432 //return
433 return fLocalDirection;
434
435}
double G4double
Definition: G4Types.hh:64
int G4int
Definition: G4Types.hh:66
#define G4UniformRand()
Definition: Randomize.hh:53
void set(double x, double y, double z)
Hep3Vector & rotateUz(const Hep3Vector &)
Definition: ThreeVector.cc:72
const G4ThreeVector & GetMomentumDirection() const
G4double GetKineticEnergy() const
G4double Value(G4double theEnergy)
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76

Referenced by G4PenelopeBremsstrahlungModel::SampleSecondaries().

◆ SetVerbosityLevel()

void G4PenelopeBremsstrahlungAngular::SetVerbosityLevel ( G4int  vl)
inline

Set/Get Verbosity level.

Definition at line 79 of file G4PenelopeBremsstrahlungAngular.hh.

79{verbosityLevel = vl;};

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