Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4RDPhotoElectricAngularGeneratorPolarized Class Reference

#include <G4RDPhotoElectricAngularGeneratorPolarized.hh>

Inheritance diagram for G4RDPhotoElectricAngularGeneratorPolarized:
Collaboration diagram for G4RDPhotoElectricAngularGeneratorPolarized:

Public Member Functions

 G4RDPhotoElectricAngularGeneratorPolarized (const G4String &name)
 
 ~G4RDPhotoElectricAngularGeneratorPolarized ()
 
G4ThreeVector GetPhotoElectronDirection (const G4ThreeVector &direction, const G4double kineticEnergy, const G4ThreeVector &polarization, const G4int shellId) const
 
void PrintGeneratorInformation () const
 
- Public Member Functions inherited from G4RDVPhotoElectricAngularDistribution
 G4RDVPhotoElectricAngularDistribution (const G4String &name)
 
virtual ~G4RDVPhotoElectricAngularDistribution ()
 

Protected Member Functions

G4ThreeVector SetPerpendicularVector (const G4ThreeVector &a) const
 

Detailed Description

Constructor & Destructor Documentation

G4RDPhotoElectricAngularGeneratorPolarized::G4RDPhotoElectricAngularGeneratorPolarized ( const G4String name)

Definition at line 71 of file G4RDPhotoElectricAngularGeneratorPolarized.cc.

72 {
73  const G4int arrayDim = 980;
74 
75  //minimum electron beta parameter allowed
76  betaArray[0] = 0.02;
77  //beta step
78  betaArray[1] = 0.001;
79  //maximum index array for a and c tables
80  betaArray[2] = arrayDim - 1;
81 
82  // read Majorant Surface Parameters. This are required in order to generate Gavrila angular photoelectron distribution
83  for(G4int level = 0; level < 2; level++){
84 
85  char nameChar0[100] = "ftab0.dat"; // K-shell Majorant Surface Parameters
86  char nameChar1[100] = "ftab1.dat"; // L-shell Majorant Surface Parameters
87 
88  G4String filename;
89  if(level == 0) filename = nameChar0;
90  if(level == 1) filename = nameChar1;
91 
92  char* path = getenv("G4LEDATA");
93  if (!path)
94  {
95  G4String excep = "G4LEDATA environment variable not set!";
96  G4Exception("G4RDPhotoElectricAngularGeneratorPolarized()",
97  "InvalidSetup", FatalException, excep);
98  }
99 
100  G4String pathString(path);
101  G4String dirFile = pathString + "/photoelectric_angular/" + filename;
102  FILE *infile;
103  infile = fopen(dirFile,"r");
104  if (infile == 0)
105  {
106  G4String excep = "Data file: " + dirFile + " not found";
107  G4Exception("G4RDPhotoElectricAngularGeneratorPolarized()",
108  "DataNotFound", FatalException, excep);
109  }
110 
111  // Read parameters into tables. The parameters are function of incident electron energy and shell level
112  G4float aRead,cRead, beta;
113  for(G4int i=0 ; i<arrayDim ;i++){
114  fscanf(infile,"%f\t %e\t %e",&beta,&aRead,&cRead);
115  aMajorantSurfaceParameterTable[i][level] = aRead;
116  cMajorantSurfaceParameterTable[i][level] = cRead;
117  }
118  fclose(infile);
119 
120  }
121 }
float G4float
Definition: G4Types.hh:77
int G4int
Definition: G4Types.hh:78
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41

Here is the call graph for this function:

G4RDPhotoElectricAngularGeneratorPolarized::~G4RDPhotoElectricAngularGeneratorPolarized ( )

Definition at line 125 of file G4RDPhotoElectricAngularGeneratorPolarized.cc.

126 {;}

Member Function Documentation

G4ThreeVector G4RDPhotoElectricAngularGeneratorPolarized::GetPhotoElectronDirection ( const G4ThreeVector direction,
const G4double  kineticEnergy,
const G4ThreeVector polarization,
const G4int  shellId 
) const
virtual

Implements G4RDVPhotoElectricAngularDistribution.

Definition at line 130 of file G4RDPhotoElectricAngularGeneratorPolarized.cc.

132 {
133  // Calculate Lorentz term (gamma) and beta parameters
134  G4double gamma = 1. + eKineticEnergy/electron_mass_c2;
135  G4double beta = std::sqrt(gamma*gamma-1.)/gamma;
136 
137  G4double theta, phi = 0;
138  G4double aBeta = 0; // Majorant surface parameter (function of the outgoing electron kinetic energy)
139  G4double cBeta = 0; // Majorant surface parameter (function of the outgoing electron kinetic energy)
140 
141  G4int shellLevel = 0;
142  if(shellId < 2) shellLevel = 0; // K-shell // Polarized model for K-shell
143  if(shellId >= 2) shellLevel = 1; // L1-shell // Polarized model for L1 and higher shells
144 
145  // For the outgoing kinetic energy find the current majorant surface parameters
146  PhotoElectronGetMajorantSurfaceAandCParameters( shellLevel, beta, &aBeta, &cBeta);
147 
148  // Generate pho and theta according to the shell level and beta parameter of the electron
149  PhotoElectronGeneratePhiAndTheta(shellLevel, beta, aBeta, cBeta, &phi, &theta);
150 
151  // Determine the rotation matrix
152  G4RotationMatrix rotation = PhotoElectronRotationMatrix(direction, polarization);
153 
154  // Compute final direction of the outgoing electron
155  G4ThreeVector final_direction = PhotoElectronComputeFinalDirection(rotation, theta, phi);
156 
157  return final_direction;
158 }
int G4int
Definition: G4Types.hh:78
float electron_mass_c2
Definition: hepunit.py:274
double G4double
Definition: G4Types.hh:76
void G4RDPhotoElectricAngularGeneratorPolarized::PrintGeneratorInformation ( ) const
virtual

Implements G4RDVPhotoElectricAngularDistribution.

Definition at line 397 of file G4RDPhotoElectricAngularGeneratorPolarized.cc.

398 {
399  G4cout << "\n" << G4endl;
400  G4cout << "Polarized Photoelectric Angular Generator" << G4endl;
401  G4cout << "PhotoElectric Electron Angular Generator based on the general Gavrila photoelectron angular distribution" << G4endl;
402  G4cout << "Includes polarization effects for K and L1 atomic shells, according to Gavrilla (1959, 1961)." << G4endl;
403  G4cout << "For higher shells the L1 cross-section is used." << G4endl;
404  G4cout << "(see Physics Reference Manual) \n" << G4endl;
405 }
G4GLOB_DLL std::ostream G4cout
#define G4endl
Definition: G4ios.hh:61
G4ThreeVector G4RDPhotoElectricAngularGeneratorPolarized::SetPerpendicularVector ( const G4ThreeVector a) const
protected

Definition at line 407 of file G4RDPhotoElectricAngularGeneratorPolarized.cc.

408 {
409  G4double dx = a.x();
410  G4double dy = a.y();
411  G4double dz = a.z();
412  G4double x = dx < 0.0 ? -dx : dx;
413  G4double y = dy < 0.0 ? -dy : dy;
414  G4double z = dz < 0.0 ? -dz : dz;
415  if (x < y) {
416  return x < z ? G4ThreeVector(-dy,dx,0) : G4ThreeVector(0,-dz,dy);
417  }else{
418  return y < z ? G4ThreeVector(dz,0,-dx) : G4ThreeVector(-dy,dx,0);
419  }
420 }
CLHEP::Hep3Vector G4ThreeVector
double x() const
tuple x
Definition: test.py:50
double z() const
double y() const
tuple z
Definition: test.py:28
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:


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