137 for ( std::vector<G4PhysicsTable*>::iterator it =
fAngleBank.begin();
139 if ( (*it) ) (*it)->clearAndDestroy();
163 for(jEl = 0 ; jEl < numOfEl; ++jEl)
173 G4cout<<
"G4NuclNuclDiffuseElastic::Initialise() the element: "
174 <<(*theElementTable)[jEl]->GetName()<<
G4endl;
223 if (iZ == 1 && iA == 1) theDef =
theProton;
226 else if (iZ == 2 && iA == 3) theDef =
G4He3::He3();
227 else if (iZ == 2 && iA == 4) theDef =
theAlpha;
241 G4double cost = 1 - 0.5*std::fabs(tMand)/ptot2;
243 if( cost >= 1.0 ) cost = 1.0;
244 else if( cost <= -1.0) cost = -1.0;
246 G4double thetaCMS = std::acos(cost);
278 if( z && (kRt > kRtC) )
310 if (iZ == 1 && iA == 1) theDef =
theProton;
313 else if (iZ == 2 && iA == 3) theDef =
G4He3::He3();
314 else if (iZ == 2 && iA == 4) theDef =
theAlpha;
328 G4double cost = 1 - 0.5*std::fabs(tMand)/ptot2;
330 if( cost >= 1.0 ) cost = 1.0;
331 else if( cost <= -1.0) cost = -1.0;
333 G4double thetaCMS = std::acos(cost);
360 if (iZ == 1 && iA == 1) theDef =
theProton;
363 else if (iZ == 2 && iA == 3) theDef =
G4He3::He3();
364 else if (iZ == 2 && iA == 4) theDef =
theAlpha;
378 G4double cost = 1 - 0.5*std::fabs(tMand)/ptot2;
380 if( cost >= 1.0 ) cost = 1.0;
381 else if( cost <= -1.0) cost = -1.0;
383 G4double thetaCMS = std::acos(cost);
403 G4double sigma, bzero, bzero2, bonebyarg, bonebyarg2, damp, damp2;
416 bzero2 = bzero*bzero;
420 bonebyarg2 = bonebyarg*bonebyarg;
435 diffuse = 0.63*
fermi;
464 sigma += mode2k2*bone2 + e2dk3t*bzero*bone;
465 sigma += kr2*bonebyarg2;
483 G4double sigma, bzero, bzero2, bonebyarg, bonebyarg2, damp, damp2;
496 bzero2 = bzero*bzero;
500 bonebyarg2 = bonebyarg*bonebyarg;
504 diffuse = 0.63*
fermi;
513 diffuse = 0.63*
fermi;
527 G4double sinHalfTheta = std::sin(0.5*theta);
528 G4double sinHalfTheta2 = sinHalfTheta*sinHalfTheta;
554 sigma += mode2k2*bone2;
555 sigma += e2dk3t*bzero*bone;
558 sigma += kr2*bonebyarg2;
576 theta = std::sqrt(alpha);
580 G4double sigma, bzero, bzero2, bonebyarg, bonebyarg2, damp, damp2;
593 bzero2 = bzero*bzero;
597 bonebyarg2 = bonebyarg*bonebyarg;
601 diffuse = 0.63*
fermi;
610 diffuse = 0.63*
fermi;
625 G4double sinHalfTheta2 = sinHalfTheta*sinHalfTheta;
651 sigma += mode2k2*bone2;
652 sigma += e2dk3t*bzero*bone;
655 sigma += kr2*bonebyarg2;
713 G4double t = 2*p*p*( 1 - std::cos(theta) );
737 if (thetaMax >
pi) thetaMax =
pi;
746 for(i = 1; i <= iMax; i++)
748 theta1 = (i-1)*thetaMax/iMax;
749 theta2 = i*thetaMax/iMax;
754 result = 0.5*(theta1 + theta2);
758 if (i > iMax ) result = 0.5*(theta1 + theta2);
764 if(result < 0.) result = 0.;
765 if(result > thetaMax) result = thetaMax;
785 G4double totElab = std::sqrt(m1*m1+p*p);
859 G4int iMomentum, iAngle;
881 G4double kinE = std::sqrt(momentum*momentum + m1*m1) - m1;
883 for( iMomentum = 0; iMomentum <
fEnergyBin; iMomentum++)
887 if( kinE < fEnergyVector->GetLowEdgeEnergy(iMomentum) )
break;
892 if ( iMomentum >= fEnergyBin ) iMomentum = fEnergyBin-1;
893 if ( iMomentum < 0 ) iMomentum = 0;
896 if (iMomentum == fEnergyBin -1 || iMomentum == 0 )
902 for(iAngle = 0; iAngle <
fAngleBin-1; iAngle++)
904 if( position < (*(*
fAngleTable)(iMomentum))(iAngle) )
break;
906 if (iAngle >= fAngleBin-1) iAngle = fAngleBin-2;
921 for(iAngle = 0; iAngle <
fAngleBin-1; iAngle++)
924 if( position > (*(*
fAngleTable)(iMomentum))(iAngle) )
break;
926 if (iAngle >= fAngleBin-1) iAngle = fAngleBin-2;
944 for(iAngle = 0; iAngle < fAngleBin-1; iAngle++)
947 if( position > (*(*
fAngleTable)(iMomentum))(iAngle) )
break;
949 if (iAngle >= fAngleBin-1) iAngle = fAngleBin-2;
963 randAngle = W1*theta1 + W2*theta2;
990 G4cout<<
"G4NuclNuclDiffuseElastic::Initialise() the element with Z = "
991 <<Z<<
"; and A = "<<A<<
G4endl;
1011 G4double alpha1, alpha2, alphaMax, alphaCoulomb, delta = 0.,
sum = 0.;
1026 partMom = std::sqrt( kinE*(kinE + 2*m1) );
1032 if(alphaMax >
pi) alphaMax =
pi;
1061 alpha1 = alphaCoulomb + delth*(j-1);
1063 alpha2 = alpha1 + delth;
1092 randAngle = (*fAngleTable)(iMomentum)->GetLowEdgeEnergy(iAngle);
1099 iAngle = (*fAngleTable)(iMomentum)->GetVectorLength() - 1;
1101 y1 = (*(*fAngleTable)(iMomentum))(iAngle-1);
1102 y2 = (*(*fAngleTable)(iMomentum))(iAngle);
1104 x1 = (*fAngleTable)(iMomentum)->GetLowEdgeEnergy(iAngle-1);
1105 x2 = (*fAngleTable)(iMomentum)->GetLowEdgeEnergy(iAngle);
1107 if ( x1 == x2 ) randAngle =
x2;
1113 randAngle = x1 + ( position -
y1 )*( x2 - x1 )/( y2 -
y1 );
1152 t =
SampleT( theParticle, ptot, A);
1155 if(!(t < 0.0 || t >= 0.0))
1159 G4cout <<
"G4NuclNuclDiffuseElastic:WARNING: A = " << A
1160 <<
" mom(GeV)= " << plab/
GeV
1161 <<
" S-wave will be sampled"
1168 G4cout <<
" t= " << t <<
" tmax= " << tmax
1169 <<
" ptot= " << ptot <<
G4endl;
1182 else if( cost <= -1.0)
1189 sint = std::sqrt((1.0-cost)*(1.0+cost));
1193 G4cout <<
"cos(t)=" << cost <<
" std::sin(t)=" << sint <<
G4endl;
1195 G4ThreeVector v1(sint*std::cos(phi),sint*std::sin(phi),cost);
1236 G4double cost = std::cos(thetaCMS);
1244 else if( cost <= -1.0)
1251 sint = std::sqrt((1.0-cost)*(1.0+cost));
1255 G4cout <<
"cos(tcms)=" << cost <<
" std::sin(tcms)=" << sint <<
G4endl;
1257 G4ThreeVector v1(sint*std::cos(phi),sint*std::sin(phi),cost);
1297 G4double cost = std::cos(thetaLab);
1305 else if( cost <= -1.0)
1312 sint = std::sqrt((1.0-cost)*(1.0+cost));
1316 G4cout <<
"cos(tlab)=" << cost <<
" std::sin(tlab)=" << sint <<
G4endl;
1318 G4ThreeVector v1(sint*std::cos(phi),sint*std::sin(phi),cost);
1346 G4cout<<
"G4NuclNuclDiffuseElastic::TestAngleTable() init the element with Z = "
1347 <<Z<<
"; and A = "<<A<<
G4endl;
1356 G4double alpha1=0., alpha2=0., alphaMax=0., alphaCoulomb=0.;
1357 G4double deltaL10 = 0., deltaL96 = 0., deltaAG = 0.;
1358 G4double sumL10 = 0.,sumL96 = 0.,sumAG = 0.;
1372 alphaMax = kRmax*kRmax/kR2;
1374 if (alphaMax > 4.) alphaMax = 4.;
1376 alphaCoulomb = kRcoul*kRcoul/kR2;
1381 fBeta = a/std::sqrt(1+a*a);
1397 alpha1 = alphaMax*(j-1)/fAngleBin;
1398 alpha2 = alphaMax*( j )/fAngleBin;
1400 if( ( alpha2 > alphaCoulomb ) &&
z )
fAddCoulomb =
true;
1405 alpha1, alpha2,epsilon);
1415 <<sumL10<<
"\t"<<sumL96<<
"\t"<<sumAG<<
G4endl;
1417 angleVector->
PutValue( j-1 , alpha1, sumL10 );
1445 if ( n < 0 ) legPol = 0.;
1446 else if( n == 0 ) legPol = 1.;
1447 else if( n == 1 ) legPol =
x;
1448 else if( n == 2 ) legPol = (3.*x*x-1.)/2.;
1449 else if( n == 3 ) legPol = (5.*x*x*x-3.*
x)/2.;
1450 else if( n == 4 ) legPol = (35.*x*x*x*x-30.*x*x+3.)/8.;
1451 else if( n == 5 ) legPol = (63.*x*x*x*x*x-70.*x*x*x+15.*
x)/8.;
1452 else if( n == 6 ) legPol = (231.*x*x*x*x*x*x-315.*x*x*x*x+105.*x*x-5.)/16.;
1457 legPol = std::sqrt( 2./(n*
CLHEP::pi*std::sin(theta+epsilon)) )*std::sin( (n+0.5)*theta+0.25*
CLHEP::pi );
1469 G4double n2, cofn, shny, chny, fn, gn;
1488 for( n = 1; n <= nMax; n++)
1492 cofn =
G4Exp(-0.5*n2)/(n2+twox2);
1494 chny = std::cosh(n*y);
1495 shny = std::sinh(n*y);
1497 fn = twox - twoxcos2xy*chny + n*sin2xy*shny;
1498 gn = twoxsin2xy*chny + n*cos2xy*shny;
1509 if(std::abs(x) < 0.0001)
1516 outRe +=
GetErf(x) + cof1*(1-cos2xy)/twox;
1517 outIm += cof1*sin2xy/twox;
1567 order /= std::sqrt(2.);
1570 G4complex a0 = 0.5*(1. + 4.*(1.+im*u2)*cosHalfThetaR2/3.)/sinThetaR;
1571 G4complex a1 = 0.5*(1. + 2.*(1.+im*u2m2p3)*cosHalfThetaR2)/sinThetaR;
1572 G4complex out = gamma*(1. - a1*dTheta) - a0;
1595 order /= std::sqrt(2.);
1597 G4complex a0 = 0.5*(1. + 4.*(1.+im*u2)*cosHalfThetaR2/3.)/sinThetaR;
1598 G4complex a1 = 0.5*(1. + 2.*(1.+im*u2m2p3)*cosHalfThetaR2)/sinThetaR;
1599 G4complex out = -gamma*(1. - a1*dTheta) - a0;
1637 G4double sindTheta = std::sin(dTheta);
1638 G4double persqrt2 = std::sqrt(0.5);
1671 for( n = 0; n <
fMaxL; n++)
1678 shiftN = std::exp( -0.5*(1.-im*
fEtaRatio)*T12b ) - 1.;
1694 G4double T12b,
a, aTemp, b2, sinThetaH = std::sin(0.5*theta);
1695 G4double sinThetaH2 = sinThetaH*sinThetaH;
1704 for( n = 1; n <
fMaxL; n++)
1706 T12b = aTemp*
G4Exp(-b2/n)/
n;
1746 fBeta = a/std::sqrt(1+a*a);
1781 fBeta = a/std::sqrt(1+a*a);
1821 if( pN < 0. ) pN = 0.;
1824 if( tN < 0. ) tN = 0.;
1843 fBeta = a/std::sqrt(1+a*a);
1874 G4double proj_energy = proj_mass + pTkin;
1875 G4double proj_momentum = std::sqrt(pTkin*(pTkin+2*proj_mass));
1890 if( proj_momentum >= 1.2 )
1894 else if( proj_momentum >= 0.6 )
1901 fEtaRatio = 15.5*proj_momentum/(27*proj_momentum*proj_momentum*proj_momentum+2);
1907 if( proj_momentum >= 10. )
1917 A0 = 100. - B0*
G4Log(3.0e7);
1919 xsection = A0 + B0*
G4Log(proj_energy) - 11
1921 0.93827*0.93827,-0.165);
1926 if(pParticle == tParticle)
1928 if( proj_momentum < 0.73 )
1932 else if( proj_momentum < 1.05 )
1934 hnXsc = 23 + 40*(
G4Log(proj_momentum/0.73))*
1935 (
G4Log(proj_momentum/0.73));
1946 if( proj_momentum < 0.8 )
1950 else if( proj_momentum < 1.4 )
1976 G4double sindTheta = std::sin(dTheta);
1983 order = std::abs(order);
1993 out = 1. + 0.5*( (0.5-cosFresnel)*(0.5-cosFresnel)+(0.5-sinFresnel)*(0.5-sinFresnel) )*prof2;
1994 out += ( cosFresnel + sinFresnel - 1. )*prof;
1998 out = 0.5*( (0.5-cosFresnel)*(0.5-cosFresnel)+(0.5-sinFresnel)*(0.5-sinFresnel) )*prof2;
2011 const G4double cof[6] = { 76.18009172947146, -86.50532032941677,
2012 24.01409824083091, -1.231739572450155,
2013 0.1208650973866179e-2, -0.5395239384953e-5 } ;
2017 tmp -= (z + 0.5) * std::log(tmp);
2020 for ( j = 0; j <= 5; j++ )
2025 return -tmp + std::log(2.5066282746310005*ser);
2035 G4double modvalue, value2, fact1, fact2, arg, shift, bessel;
2037 modvalue = std::fabs(value);
2041 value2 = value*
value;
2043 fact1 = 57568490574.0 + value2*(-13362590354.0
2044 + value2*( 651619640.7
2045 + value2*(-11214424.18
2046 + value2*( 77392.33017
2047 + value2*(-184.9052456 ) ) ) ) );
2049 fact2 = 57568490411.0 + value2*( 1029532985.0
2050 + value2*( 9494680.718
2051 + value2*(59272.64853
2052 + value2*(267.8532712
2053 + value2*1.0 ) ) ) );
2055 bessel = fact1/fact2;
2063 shift = modvalue-0.785398164;
2065 fact1 = 1.0 + value2*(-0.1098628627e-2
2066 + value2*(0.2734510407e-4
2067 + value2*(-0.2073370639e-5
2068 + value2*0.2093887211e-6 ) ) );
2070 fact2 = -0.1562499995e-1 + value2*(0.1430488765e-3
2071 + value2*(-0.6911147651e-5
2072 + value2*(0.7621095161e-6
2073 - value2*0.934945152e-7 ) ) );
2075 bessel = std::sqrt(0.636619772/modvalue)*(std::cos(shift)*fact1 - arg*std::sin(shift)*fact2 );
2087 G4double modvalue, value2, fact1, fact2, arg, shift, bessel;
2089 modvalue = std::fabs(value);
2091 if ( modvalue < 8.0 )
2093 value2 = value*
value;
2095 fact1 = value*(72362614232.0 + value2*(-7895059235.0
2096 + value2*( 242396853.1
2097 + value2*(-2972611.439
2098 + value2*( 15704.48260
2099 + value2*(-30.16036606 ) ) ) ) ) );
2101 fact2 = 144725228442.0 + value2*(2300535178.0
2102 + value2*(18583304.74
2103 + value2*(99447.43394
2104 + value2*(376.9991397
2105 + value2*1.0 ) ) ) );
2106 bessel = fact1/fact2;
2114 shift = modvalue - 2.356194491;
2116 fact1 = 1.0 + value2*( 0.183105e-2
2117 + value2*(-0.3516396496e-4
2118 + value2*(0.2457520174e-5
2119 + value2*(-0.240337019e-6 ) ) ) );
2121 fact2 = 0.04687499995 + value2*(-0.2002690873e-3
2122 + value2*( 0.8449199096e-5
2123 + value2*(-0.88228987e-6
2124 + value2*0.105787412e-6 ) ) );
2126 bessel = std::sqrt( 0.636619772/modvalue)*(std::cos(shift)*fact1 - arg*std::sin(shift)*fact2);
2128 if (value < 0.0) bessel = -bessel;
const G4ParticleDefinition * thePionPlus
G4double CalculateNuclearRad(G4double A)
static G4PionMinus * PionMinus()
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
G4complex GetErfcInt(G4complex z)
G4double GetDiffElasticProb(G4double theta)
G4complex AmplitudeGG(G4double theta)
G4complex GetErfInt(G4complex z)
void SetMinEnergy(G4double anEnergy)
std::vector< ExP01TrackerHit * > a
void SetMaxEnergy(const G4double anEnergy)
static G4ParticleTable * GetParticleTable()
static constexpr double MeV
G4double GetLowEdgeEnergy(size_t binNumber) const
G4double CalculateParticleBeta(const G4ParticleDefinition *particle, G4double momentum)
G4double GetDiffElasticSumProbA(G4double alpha)
Float_t y1[n_points_granero]
G4ParticleDefinition * theAlpha
G4complex GetErfComp(G4complex z, G4int nMax)
G4complex GammaMore(G4double theta)
void InitialiseOnFly(G4double Z, G4double A)
const G4ParticleDefinition * fParticle
static constexpr double keV
Float_t x1[n_points_granero]
G4complex CoulombAmplitude(G4double theta)
static constexpr double hbarc
G4double Legendre96(T &typeT, F f, G4double a, G4double b)
G4double fRutherfordTheta
static G4PionPlus * PionPlus()
G4double GetDiffElasticSumProb(G4double theta)
G4IonTable * GetIonTable() const
G4double lowEnergyLimitHE
G4double GetPDGCharge() const
G4double GetTotalMomentum() const
G4double BesselJone(G4double z)
void InitDynParameters(const G4ParticleDefinition *theParticle, G4double partMom)
G4complex GammaLogarithm(G4complex xx)
G4double GetInvCoulombElasticXsc(const G4ParticleDefinition *particle, G4double tMand, G4double momentum, G4double A, G4double Z)
static G4Proton * Proton()
void insertAt(size_t, G4PhysicsVector *)
G4ParticleDefinition * theDeuteron
G4double GetExpCos(G4double x)
G4double GetDiffuseElasticSumXsc(const G4ParticleDefinition *particle, G4double theta, G4double momentum, G4double A, G4double Z)
std::vector< G4PhysicsTable * > fAngleBank
G4double SampleThetaCMS(const G4ParticleDefinition *aParticle, G4double p, G4double A)
G4double CalculateAm(G4double momentum, G4double n, G4double Z)
G4double GetAtomicMassAmu(const G4String &symb) const
G4double GetPDGMass() const
G4double G4Log(G4double x)
std::vector< G4String > fElementNameVector
G4double GetSint(G4double x)
static G4Deuteron * Deuteron()
G4double lowEnergyRecoilLimit
G4complex GammaLess(G4double theta)
Float_t y2[n_points_geant4]
G4double GetRatioGen(G4double theta)
G4ParticleDefinition * theNeutron
G4double AdaptiveGauss(T &typeT, F f, G4double a, G4double b, G4double e)
G4PhysicsTable * fAngleTable
virtual G4double SampleInvariantT(const G4ParticleDefinition *p, G4double plab, G4int Z, G4int A)
static constexpr double TeV
G4double SampleCoulombMuCMS(const G4ParticleDefinition *aParticle, G4double p)
G4double GetHadronNucleonXscNS(G4ParticleDefinition *pParticle, G4double pTkin, G4ParticleDefinition *tParticle)
std::complex< G4double > G4complex
G4double ThetaCMStoThetaLab(const G4DynamicParticle *aParticle, G4double tmass, G4double thetaCMS)
static constexpr double fermi
static G4Pow * GetInstance()
G4NuclNuclDiffuseElastic()
void TestAngleTable(const G4ParticleDefinition *theParticle, G4double partMom, G4double Z, G4double A)
G4ParticleDefinition * GetDefinition() const
const XML_Char int const XML_Char * value
void InitParametersGla(const G4DynamicParticle *aParticle, G4double partMom, G4double Z, G4double A)
G4double powA(G4double A, G4double y) const
G4double CalculateZommerfeld(G4double beta, G4double Z1, G4double Z2)
G4double ThetaLabToThetaCMS(const G4DynamicParticle *aParticle, G4double tmass, G4double thetaLab)
G4int GetBaryonNumber() const
G4double GetFresnelIntegrandXsc(G4double alpha)
G4double SampleTableT(const G4ParticleDefinition *aParticle, G4double p, G4double Z, G4double A)
G4double SampleT(const G4ParticleDefinition *aParticle, G4double p, G4double A)
virtual ~G4NuclNuclDiffuseElastic()
std::vector< G4Element * > G4ElementTable
static const G4double alpha
static constexpr double twopi
G4ParticleDefinition * theProton
G4LorentzVector Get4Momentum() const
double A(double temperature)
G4double DampFactor(G4double z)
static G4double GetNuclearMass(const G4double A, const G4double Z)
G4double CalcMandelstamS(const G4double mp, const G4double mt, const G4double Plab)
void CalculateRutherfordAnglePar()
G4double CalculateCoulombPhase(G4int n)
ThreeVector shoot(const G4int Ap, const G4int Af)
G4double GetErf(G4double x)
const G4LorentzVector & Get4Momentum() const
G4double G4ParticleHPJENDLHEData::G4double result
G4PhysicsLogVector * fEnergyVector
G4double GetDiffuseElasticXsc(const G4ParticleDefinition *particle, G4double theta, G4double momentum, G4double A)
G4double GetInvElasticSumXsc(const G4ParticleDefinition *particle, G4double tMand, G4double momentum, G4double A, G4double Z)
void CalculateCoulombPhaseZero()
double epsilon(double density, double temperature)
static constexpr double twopi
G4complex AmplitudeNear(G4double theta)
G4complex PhaseNear(G4double theta)
const G4ParticleDefinition * GetDefinition() const
void InitParameters(const G4ParticleDefinition *theParticle, G4double partMom, G4double Z, G4double A)
static G4Triton * Triton()
G4double GetExpSin(G4double x)
G4double fNuclearRadiusCof
static G4Neutron * Neutron()
G4double IntegralElasticProb(const G4ParticleDefinition *particle, G4double theta, G4double momentum, G4double A)
G4double GetInvElasticXsc(const G4ParticleDefinition *particle, G4double theta, G4double momentum, G4double A, G4double Z)
G4double ProfileNear(G4double theta)
G4double GetCint(G4double x)
G4double lowestEnergyLimit
G4double fCofAlphaCoulomb
G4double BesselJzero(G4double z)
static constexpr double degree
G4double GetKineticEnergy() const
G4GLOB_DLL std::ostream G4cout
G4double fRutherfordRatio
Hep3Vector boostVector() const
std::vector< G4double > fElementNumberVector
G4double fHalfRutThetaTg2
G4double SampleTableThetaCMS(const G4ParticleDefinition *aParticle, G4double p, G4double Z, G4double A)
G4double GetCoulombElasticXsc(const G4ParticleDefinition *particle, G4double theta, G4double momentum, G4double Z)
static size_t GetNumberOfElements()
void PutValue(size_t index, G4double energy, G4double dataValue)
G4double GetTotalMomentum() const
static constexpr double pi
G4ParticleDefinition * GetIon(G4int Z, G4int A, G4int lvl=0)
G4double Profile(G4double theta)
G4double SampleThetaLab(const G4HadProjectile *aParticle, G4double tmass, G4double A)
Float_t x2[n_points_geant4]
const G4ParticleDefinition * thePionMinus
G4double fNuclearRadiusSquare
G4double GetScatteringAngle(G4int iMomentum, G4int iAngle, G4double position)
static constexpr double GeV
G4double Legendre10(T &typeT, F f, G4double a, G4double b)
static constexpr double GeV
G4double GetIntegrandFunction(G4double theta)
static constexpr double pi
static constexpr double millibarn
G4double BesselOneByArg(G4double z)
G4double GetLegendrePol(G4int n, G4double x)
HepLorentzVector & boost(double, double, double)
G4complex AmplitudeGla(G4double theta)
static G4NistManager * Instance()
static G4ElementTable * GetElementTable()
G4complex AmplitudeSim(G4double theta)