79 return 4*costheta2/(2*costheta2-klk2+2*std::sqrt(costheta2*(costheta2-klk2)));
81 return 4*costheta2/klk2;
91 (2*costhetas2+klks2+2*std::sqrt(costhetas2*(costhetas2+klks2)));
105 if ((std::fabs(thetai-thetao)<AngCut) && (std::fabs(phio)<AngCut))
111 G4double sinthetai=std::sin(thetai);
112 G4double sinthetao=std::sin(thetao);
114 (sinthetai*sinthetai+sinthetao*sinthetao-
115 2.*sinthetai*sinthetao*std::cos(phio));
120 return b2*w2/twopi*std::exp(-mu_squared*w2/2);
136 if ((std::fabs(thetarefract-thetaSo)<AngCut) && (std::fabs(phiSo)<AngCut))
140 G4double sinthetai=std::sin(thetai);
141 G4double sinthetaSo=std::sin(thetaSo);
144 mu_squared=k*k*sinthetai*sinthetai+kS*kS*sinthetaSo*sinthetaSo-
145 2.*k*kS*sinthetai*sinthetaSo*std::cos(phiSo);
150 return b2*w2/twopi*std::exp(-mu_squared*w2/2);
166 G4double a_max_theta_o, max_theta_o=theta_i, a_max_phi_o, max_phi_o=0.;
178 G4double costheta_i=std::cos(theta_i);
179 G4double costheta_i_squared=costheta_i*costheta_i;
182 G4double kl4d4=neutron_mass_c2/hbarc_squared*neutron_mass_c2/
183 hbarc_squared*fermipot*fermipot;
191 G4double k2=2*neutron_mass_c2*E/hbarc_squared;
197 for (theta_o=0.*
degree; theta_o<=90.*
degree+1e-6; theta_o+=ang_steptheta)
199 costheta_o_squared=std::cos(theta_o)*std::cos(theta_o);
203 for (phi_o=-180.*
degree; phi_o<=180.*
degree+1e-6; phi_o+=ang_stepphi)
207 Intens=kl4d4/costheta_i*
S2(costheta_i_squared,klk2)*
208 S2(costheta_o_squared,klk2)*
209 Fmu(k2,theta_i,theta_o,phi_o,b2,w2,AngCut)*std::sin(theta_o);
226 wkeit+=Intens*ang_steptheta*ang_stepphi;
238 while ((ang_stepphi>=AngCut*AngCut) || (ang_steptheta>=AngCut*AngCut))
240 a_max_theta_o=max_theta_o;
241 a_max_phi_o=max_phi_o;
249 for (theta_o=a_max_theta_o-ang_steptheta;
250 theta_o<=a_max_theta_o-ang_steptheta+1e-6;
251 theta_o+=ang_steptheta)
254 costheta_o_squared=std::cos(theta_o)*std::cos(theta_o);
255 for (phi_o=a_max_phi_o-ang_stepphi;
256 phi_o<=a_max_phi_o+ang_stepphi+1e-6;
260 Intens=kl4d4/costheta_i*
S2(costheta_i_squared, klk2)*
261 S2(costheta_o_squared,klk2)*
262 Fmu(k2,theta_i,theta_o,phi_o,b2,w2,AngCut)*std::sin(theta_o);
285 G4double a_max_thetas_o, max_thetas_o = theta_i;
286 G4double a_max_phis_o, max_phis_o = 0.;
292 G4double costheta_i=std::cos(theta_i);
293 G4double costheta_i_squared=costheta_i*costheta_i;
298 if (E*costheta_i_squared < fermipot)
return wkeit;
301 G4double kl4d4=neutron_mass_c2/hbarc_squared*neutron_mass_c2/
302 hbarc_squared*fermipot*fermipot;
307 G4double klks2=fermipot/(E-fermipot);
310 G4double ksdk=std::sqrt((E-fermipot)/E);
315 G4double k=std::sqrt(2*neutron_mass_c2*E/hbarc_squared);
320 for (thetas_o=0.*
degree; thetas_o<=90.*
degree+1e-6; thetas_o+=ang_steptheta)
322 costhetas_o_squared=std::cos(thetas_o)*std::cos(thetas_o);
324 for (phis_o=-180.*
degree; phis_o<=180.*
degree+1e-6; phis_o+=ang_stepphi)
327 if (costhetas_o_squared>=-klks2) {
332 if (std::fabs(std::sin(theta_i)/ksdk) <= 1.)
333 thetarefract = std::asin(std::sin(theta_i)/ksdk);
335 IntensS = kl4d4/costheta_i*ksdk*
S2(costheta_i_squared, klk2)*
336 SS2(costhetas_o_squared,klks2)*
337 FmuS(k,kS,theta_i,thetas_o,phis_o,b2,w2,AngCut,thetarefract)*
348 wkeit+=IntensS*ang_steptheta*ang_stepphi;
359 while (ang_stepphi>=AngCut*AngCut || ang_steptheta>=AngCut*AngCut)
361 a_max_thetas_o=max_thetas_o;
362 a_max_phis_o=max_phis_o;
367 for (thetas_o=a_max_thetas_o-ang_steptheta;
368 thetas_o<=a_max_thetas_o-ang_steptheta+1e-6;
369 thetas_o+=ang_steptheta)
371 costhetas_o_squared=std::cos(thetas_o)*std::cos(thetas_o);
372 for (phis_o=a_max_phis_o-ang_stepphi;
373 phis_o<=a_max_phis_o+ang_stepphi+1e-6;
377 if (std::fabs(std::sin(theta_i)/ksdk) <= 1.)
378 thetarefract = std::asin(std::sin(theta_i)/ksdk);
380 IntensS=kl4d4/costheta_i*ksdk*
S2(costheta_i_squared, klk2)*
381 SS2(costhetas_o_squared,klks2)*
382 FmuS(k,kS,theta_i,thetas_o,phis_o,b2,w2,AngCut,thetarefract)*
387 max_thetas_o=thetas_o;
407 G4double kl4d4=neutron_mass_c2/hbarc_squared*neutron_mass_c2/
408 hbarc_squared*fermipot*fermipot;
413 G4double costheta_i=std::cos(theta_i);
414 G4double costheta_o=std::cos(theta_o);
418 return kl4d4/costheta_i*
S2(costheta_i*costheta_i, klk2)*
419 S2(costheta_o*costheta_o,klk2)*
420 Fmu(2*neutron_mass_c2*E/hbarc_squared,theta_i,theta_o,phi_o,b*b,w*w,AngCut)*
433 G4double kl4d4=neutron_mass_c2/hbarc_squared*neutron_mass_c2/
434 hbarc_squared*fermipot*fermipot;
439 G4double klks2=fermipot/(E-fermipot);
447 G4double ksdk=std::sqrt((E-fermipot)/E);
450 G4double k=std::sqrt(2*neutron_mass_c2*E/hbarc_squared);
455 G4double costheta_i=std::cos(theta_i);
456 G4double costhetas_o=std::cos(thetas_o);
461 if(std::fabs(std::sin(theta_i)/ksdk) <= 1.)thetarefract = std::asin(std::sin(theta_i)/ksdk);
463 return kl4d4/costheta_i*ksdk*
S2(costheta_i*costheta_i, klk2)*
464 SS2(costhetas_o*costhetas_o,klks2)*
465 FmuS(k,kS,theta_i,thetas_o,phis_o,b*b,w*w,AngCut,thetarefract)*
G4double Fmu(G4double, G4double, G4double, G4double, G4double, G4double, G4double)
G4UCNMicroRoughnessHelper()
G4double IntIplus(G4double, G4double, G4double, G4int, G4int, G4double, G4double, G4double *, G4double)
G4double ProbIplus(G4double, G4double, G4double, G4double, G4double, G4double, G4double, G4double)
G4GLOB_DLL std::ostream G4cout
G4double S2(G4double, G4double)
G4double FmuS(G4double, G4double, G4double, G4double, G4double, G4double, G4double, G4double, G4double)
static G4UCNMicroRoughnessHelper * GetInstance()
~G4UCNMicroRoughnessHelper()
T max(const T t1, const T t2)
brief Return the largest of the two arguments
static const double degree
static G4UCNMicroRoughnessHelper * fpInstance
G4double SS2(G4double, G4double)
G4double ProbIminus(G4double, G4double, G4double, G4double, G4double, G4double, G4double, G4double)
G4double IntIminus(G4double, G4double, G4double, G4int, G4int, G4double, G4double, G4double *, G4double)