#include <G4PhotoElectricAngularGeneratorPolarized.hh>
Inheritance diagram for G4PhotoElectricAngularGeneratorPolarized:
Public Member Functions | |
G4PhotoElectricAngularGeneratorPolarized () | |
~G4PhotoElectricAngularGeneratorPolarized () | |
virtual G4ThreeVector & | SampleDirection (const G4DynamicParticle *dp, G4double eKinEnergy, G4int shellId, const G4Material *mat=0) |
void | PrintGeneratorInformation () const |
Protected Member Functions | |
G4ThreeVector | PerpendicularVector (const G4ThreeVector &a) const |
Definition at line 52 of file G4PhotoElectricAngularGeneratorPolarized.hh.
G4PhotoElectricAngularGeneratorPolarized::G4PhotoElectricAngularGeneratorPolarized | ( | ) |
Definition at line 70 of file G4PhotoElectricAngularGeneratorPolarized.cc.
References FatalException, and G4Exception().
00071 :G4VEmAngularDistribution("AngularGenSauterGavrilaPolarized") 00072 { 00073 const G4int arrayDim = 980; 00074 00075 //minimum electron beta parameter allowed 00076 betaArray[0] = 0.02; 00077 //beta step 00078 betaArray[1] = 0.001; 00079 //maximum index array for a and c tables 00080 betaArray[2] = arrayDim - 1; 00081 00082 // read Majorant Surface Parameters. This are required in order to generate Gavrila angular photoelectron distribution 00083 for(G4int level = 0; level < 2; level++){ 00084 00085 char nameChar0[100] = "ftab0.dat"; // K-shell Majorant Surface Parameters 00086 char nameChar1[100] = "ftab1.dat"; // L-shell Majorant Surface Parameters 00087 00088 G4String filename; 00089 if(level == 0) filename = nameChar0; 00090 if(level == 1) filename = nameChar1; 00091 00092 char* path = getenv("G4LEDATA"); 00093 if (!path) 00094 { 00095 G4String excep = "G4EMDataSet - G4LEDATA environment variable not set"; 00096 G4Exception("G4PhotoElectricAngularGeneratorPolarized::G4PhotoElectricAngularGeneratorPolarized", 00097 "em0006",FatalException,"G4LEDATA environment variable not set"); 00098 return; 00099 } 00100 00101 G4String pathString(path); 00102 G4String dirFile = pathString + "/photoelectric_angular/" + filename; 00103 std::ifstream infile(dirFile); 00104 if (!infile.is_open()) 00105 { 00106 G4String excep = "data file: " + dirFile + " not found"; 00107 G4Exception("G4PhotoElectricAngularGeneratorPolarized::G4PhotoElectricAngularGeneratorPolarized", 00108 "em0003",FatalException,excep); 00109 return; 00110 } 00111 00112 // Read parameters into tables. The parameters are function of incident electron energy and shell level 00113 G4float aRead=0,cRead=0, beta=0; 00114 for(G4int i=0 ; i<arrayDim ;i++){ 00115 //fscanf(infile,"%f\t %e\t %e",&beta,&aRead,&cRead); 00116 infile >> beta >> aRead >> cRead; 00117 aMajorantSurfaceParameterTable[i][level] = aRead; 00118 cMajorantSurfaceParameterTable[i][level] = cRead; 00119 } 00120 infile.close(); 00121 } 00122 }
G4PhotoElectricAngularGeneratorPolarized::~G4PhotoElectricAngularGeneratorPolarized | ( | ) |
G4ThreeVector G4PhotoElectricAngularGeneratorPolarized::PerpendicularVector | ( | const G4ThreeVector & | a | ) | const [protected] |
Definition at line 413 of file G4PhotoElectricAngularGeneratorPolarized.cc.
00415 { 00416 G4double dx = a.x(); 00417 G4double dy = a.y(); 00418 G4double dz = a.z(); 00419 G4double x = dx < 0.0 ? -dx : dx; 00420 G4double y = dy < 0.0 ? -dy : dy; 00421 G4double z = dz < 0.0 ? -dz : dz; 00422 if (x < y) { 00423 return x < z ? G4ThreeVector(-dy,dx,0) : G4ThreeVector(0,-dz,dy); 00424 }else{ 00425 return y < z ? G4ThreeVector(dz,0,-dx) : G4ThreeVector(-dy,dx,0); 00426 } 00427 }
void G4PhotoElectricAngularGeneratorPolarized::PrintGeneratorInformation | ( | ) | const |
Definition at line 402 of file G4PhotoElectricAngularGeneratorPolarized.cc.
References G4cout, and G4endl.
00403 { 00404 G4cout << "\n" << G4endl; 00405 G4cout << "Polarized Photoelectric Angular Generator" << G4endl; 00406 G4cout << "PhotoElectric Electron Angular Generator based on the general Gavrila photoelectron angular distribution" << G4endl; 00407 G4cout << "Includes polarization effects for K and L1 atomic shells, according to Gavrilla (1959, 1961)." << G4endl; 00408 G4cout << "For higher shells the L1 cross-section is used." << G4endl; 00409 G4cout << "(see Physics Reference Manual) \n" << G4endl; 00410 }
G4ThreeVector & G4PhotoElectricAngularGeneratorPolarized::SampleDirection | ( | const G4DynamicParticle * | dp, | |
G4double | eKinEnergy, | |||
G4int | shellId, | |||
const G4Material * | mat = 0 | |||
) | [virtual] |
Implements G4VEmAngularDistribution.
Definition at line 128 of file G4PhotoElectricAngularGeneratorPolarized.cc.
References G4VEmAngularDistribution::fLocalDirection, G4DynamicParticle::GetMomentumDirection(), and G4DynamicParticle::GetPolarization().
00133 { 00134 // (shellId == 0) - K-shell - Polarized model for K-shell 00135 // (shellId > 0) - L1-shell - Polarized model for L1 and higher shells 00136 00137 // Calculate Lorentz term (gamma) and beta parameters 00138 G4double gamma = 1 + eKinEnergy/electron_mass_c2; 00139 G4double beta = std::sqrt((gamma - 1)*(gamma + 1))/gamma; 00140 00141 const G4ThreeVector& direction = dp->GetMomentumDirection(); 00142 const G4ThreeVector& polarization = dp->GetPolarization(); 00143 00144 G4double theta, phi = 0; 00145 // Majorant surface parameters 00146 // function of the outgoing electron kinetic energy 00147 G4double aBeta = 0; 00148 G4double cBeta = 0; 00149 00150 // For the outgoing kinetic energy 00151 // find the current majorant surface parameters 00152 PhotoElectronGetMajorantSurfaceAandCParameters(shellId, beta, &aBeta, &cBeta); 00153 00154 // Generate pho and theta according to the shell level 00155 // and beta parameter of the electron 00156 PhotoElectronGeneratePhiAndTheta(shellId, beta, aBeta, cBeta, &phi, &theta); 00157 00158 // Determine the rotation matrix 00159 const G4RotationMatrix rotation = 00160 PhotoElectronRotationMatrix(direction, polarization); 00161 00162 // Compute final direction of the outgoing electron 00163 fLocalDirection = PhotoElectronComputeFinalDirection(rotation, theta, phi); 00164 00165 return fLocalDirection; 00166 }