| Summary: | Reflectivity calculation with complex refractive index is incorrect | ||
|---|---|---|---|
| Product: | Geant4 | Reporter: | iekikei |
| Component: | processes/optical | Assignee: | gum |
| Status: | RESOLVED FIXED | ||
| Severity: | major | ||
| Priority: | P5 | ||
| Version: | 10.1 | ||
| Hardware: | PC | ||
| OS: | Linux | ||
| Attachments: | G4OpBoundaryProcess.cc | ||
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 |
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?