73 lowEnergyLimit(250*
eV),
74 highEnergyLimit(100*
GeV),
75 intrinsicLowEnergyLimit(10*
eV),
76 intrinsicHighEnergyLimit(100*
GeV)
81 G4Exception(
"G4LowEnergyCompton::G4LowEnergyCompton()",
83 "Energy outside intrinsic process validity range!");
89 G4String scatterFile =
"comp/ce-sf-";
122 G4String crossSectionFile =
"comp/ce-cs-";
168 G4double alpha1 = -std::log(epsilon0);
169 G4double alpha2 = 0.5 * (1. - epsilon0Sq);
191 epsilonSq = epsilon0Sq + (1. - epsilon0Sq) *
G4UniformRand();
192 epsilon = std::sqrt(epsilonSq);
195 oneCosT = (1. -
epsilon) / ( epsilon * e0m);
196 sinT2 = oneCosT * (2. - oneCosT);
197 G4double x = std::sqrt(oneCosT/2.) / (wlPhoton/
cm);
199 gReject = (1. - epsilon * sinT2 / (1. + epsilonSq)) * scatteringFunction;
204 G4double sinTheta = std::sqrt (sinT2);
206 G4double dirX = sinTheta * std::cos(phi);
207 G4double dirY = sinTheta * std::sin(phi);
216 G4int maxDopplerIterations = 1000;
218 G4double photonEoriginal = epsilon * photonEnergy0;
229 eMax = photonEnergy0 - bindingE;
235 G4double pDoppler2 = pDoppler * pDoppler;
237 G4double var3 = var2*var2 - pDoppler2;
238 G4double var4 = var2 - pDoppler2 * cosTheta;
239 G4double var = var4*var4 - var3 + pDoppler2 * var3;
245 if (
G4UniformRand() < 0.5) photonE = (var4 - varSqrt) * scale;
246 else photonE = (var4 + varSqrt) * scale;
252 }
while ( iteration <= maxDopplerIterations &&
253 (photonE < 0. || photonE > eMax || photonE < eMax*
G4UniformRand()) );
257 if (iteration >= maxDopplerIterations)
259 photonE = photonEoriginal;
266 photonDirection1.
rotateUz(photonDirection0);
272 if (photonEnergy1 > 0.)
283 G4double eKineticEnergy = photonEnergy0 - photonEnergy1 - bindingE;
286 G4double electronE = photonEnergy0 * (1. -
epsilon) + electron_mass_c2;
292 cosThetaE = (eTotalEnergy + photonEnergy1 )* (1. - epsilon) / std::sqrt(electronP2);
293 sinThetaE = -1. * std::sqrt(1. - cosThetaE * cosThetaE);
296 G4double eDirX = sinThetaE * std::cos(phi);
297 G4double eDirY = sinThetaE * std::sin(phi);
307 eDirection.
rotateUz(photonDirection0);
336 size_t materialIndex = couple->
GetIndex();