Hello, I noticed that the reflectivity is not correctly calculated when the complex refractive index is specified. When a photon propagate from material 1 to material 2, reflectivity must be calculated from the complex refractive index of two materials. However, in G4OpBoundaryProcess::GetReflectivity(), only the complex refractive index of material 2 is used for the calculation of reflectivity. It seems that the calculation assumes refractive index of material 1 to be always 1 (i.e. vacuum). Therefore, when I set the complex refractive index for the optical surface by using AddProperty("REALRINDEX", ...) and AddProperty("IMAGINARYRINDEX", ...), reflectivity is not calculated properly. In fact, in my simulation, I tried to change the refractive index of material 1, but the reflectivity did not change. Is this designed to be like this?
This issue has been brought up before. Please, see: http://hypernews.slac.stanford.edu/HyperNews/geant4/get/opticalphotons/472.html and the lengthy discussion following this posting. In particular, the last entry: http://hypernews.slac.stanford.edu/HyperNews/geant4/get/opticalphotons/472/2/1/1/1/1.html I'd be grateful if someone would tell me what exactly should be changed in the code; e.g. how the code should read correctly.
Thank you for the quick response. > http://hypernews.slac.stanford.edu/HyperNews/geant4/get/opticalphotons/472/2/1/1/1/1.html Yes, this post seems to suggest exactly what we have to do. Here I copied those formulas in the post: n1 * sin(theta_1) = n2 * sin (theta_2), where r_s = (n1*cos(theta_1)-n2*cos(theta_2))/(n1*cos(theta_1)+n2*cos(theta_2)) r_t = (n1*cos(theta_2)-n2*cos(theta_1))/(n1*cos(theta_2)+n2*cos(theta_1)) In those formulas, if you substitute 1 to n1, it is equivalent to the formulas in line 1228-1236 of G4OpBoundaryProcess.cc: CosPhi=std::sqrt(u-((std::sin(incidentangle)*std::sin(incidentangle))/(N*N))); numeratorTE = std::cos(incidentangle) - N*CosPhi; denominatorTE = std::cos(incidentangle) + N*CosPhi; rTE = numeratorTE/denominatorTE; numeratorTM = N*std::cos(incidentangle) - CosPhi; denominatorTM = N*std::cos(incidentangle) + CosPhi; rTM = numeratorTM/denominatorTM; Therefore, for example, I think you can modify the codes above like this: G4complex N1 = complex refractive index of pre-steppoint; CosPhi=std::sqrt(u-((std::sin(incidentangle)*std::sin(incidentangle))*(N1*N1)/(N*N))); numeratorTE = N1*std::cos(incidentangle) - N*CosPhi; denominatorTE = N1*std::cos(incidentangle) + N*CosPhi; rTE = numeratorTE/denominatorTE; numeratorTM = N*std::cos(incidentangle) - N1*CosPhi; denominatorTM = N*std::cos(incidentangle) + N1*CosPhi; rTM = numeratorTM/denominatorTM; (please also check the code by yourself!)
Created attachment 351 [details] G4OpBoundaryProcess.cc
Thanks Kei Ieki, G4OpBoundaryProcess::GetReflectivity has been changed accordingly and uses the PreStepPoint material's complex index of refraction if provided by the user, else: N1(Rindex1,0) Please see attached. The code is in tag op-V10-01-03 and will be in G4.10.2 Peter
