VirtualPhotonics / VTS

Virtual Tissue Simulator
https://virtualphotonics.org
Other
34 stars 9 forks source link

Bugfix/136 montecarlo tissues getdistancetoboundary logic in error #137

Closed hayakawa16 closed 6 months ago

hayakawa16 commented 6 months ago

Draft PR for now.

hayakawa16 commented 6 months ago

Thank you @lmalenfant!!

janakarana commented 6 months ago

Is it possible to include unit tests to check the outputs of GetRefractedDirection and GetAngleRelativeToBoundaryNormal in MultiConcentricInclusionTissue.cs for different input settings ?

  1. SurfaceNormal and photon.DP.Direction < 90
  2. SurfaceNormal and photon.DP.Direction > 90
  3. SurfaceNormal and photon.DP.Direction = 90
hayakawa16 commented 6 months ago

Hi @janakarana, thanks for your comments.

For GetRefractedDirection, I have multiple checks (15 in all) including the Cases 1-4 documented in the header. I also have tests when going from 1.4 to 1.0 and outside critical angle. Case 1: is exiting the top at 45 degree angle. This is at 45 degree angle with surface normal (top normal is pointed up). Case 2: is exiting bottom at 45 degree angle. This is at 45 degree angle with surface normal (bottom normal is pointed down). Case 3: is entering top at 30 degree angle. This is in opposite direction from normal at ~30 degrees (or ~120 degrees). Case 4: is entering bottom at 30 degree angle. This is in opposite direction from normal at ~30 degrees (or ~120 degrees). Do you have other suggestions?

For GetAngleRelativeToBoundaryNormal, I will definitely fill out the tests more thoroughly. This gets called by Photon anytime photon crosses tissue region so didn't cover it as fully. I will now.

hayakawa16 commented 6 months ago

In forming the following conclusion, I've just checked sites that describe Fresnel and also the Monte Carlo C code I inherited from Dunn/Smithies. Please push back if you think I'm in error.

GetAngleRelativeToBoundaryNormal is solely called by Photon/CrossRegionOrReflect. In that method it is used to determine the probability of reflecting using Fresnel's law. In Fresnel's law (for particles, not waves), the cos(theta_i) is dot product of the incoming ray with the normal "interface" (reference wikipedia section "Configuration"), it is not with the "outward" directional normal. So the cos(theta_i) is determined by taking the dot product of the incoming ray with the normal interface (assumed to be in same direction as the incoming ray, so dot product is not negative).

I just checked my old C code and what gets sent into Fresnel. This code only handles layer regions. When the photon is crossing down into a layer, cos(theta_i) is the photon direction cosine Uz (Uz is positive into the tissue so Uz is positive). When the photon is crossing up into a layer, the negative of the direction cosine -Uz is used. In other words, the cos(theta_i) is always positive. This got translated to our C# code in MultiLayerTissue/GetAngleRelativeToBoundaryNormal as Math.Abs(photon.DP.Direction.Uz).

Our tissue regions' SurfaceNormal is directional, it is the outward normal for curved surfaces. So if it is used to determine the cos(theta_i), the absolute value needs to be taken. Therefore, our tissues GetAngleRelativeToBoundaryNormal should take the absolute value of the dot product with the incoming ray and the region's SurfaceNormal.

I am currently writing unit tests based on code that is assuming this conclusion.

janakarana commented 6 months ago

Thanks hayakawa16 or this clarification. I agree with your logic. It is correct. Reflection at a cylinder (line 246-253 in MultiConcentricInclusionTissue.cs) is correct and provides the correct direction. I have a question about Line 298: if (currentN > nextN && sinTheta2Squared > 1.0) return GetReflectedDirection(currentPosition, currentDirection);

Why do you use sinTheta2Squared > 1.0 condition? Sine angle is never greater than 1.

janakarana commented 6 months ago

Sorry. I understood why that condition (sinTheta2Squared > 1.0) is for. It is for the total internal reflection. Everything looks good to me.

janakarana commented 6 months ago

Everything looks good

sonarcloud[bot] commented 6 months ago

Quality Gate Passed Quality Gate passed

Issues
0 New issues
0 Accepted issues

Measures
0 Security Hotspots
100.0% Coverage on New Code
0.2% Duplication on New Code

See analysis details on SonarCloud

hayakawa16 commented 6 months ago

Hi @janakarana, thank you very much for your thorough review! In case you'd like to see the reference I used for this code, I have attached it. I also referred this website in section "Fresnel": https://blog.demofox.org/2017/01/09/raytracing-reflection-refraction-fresnel-total-internal-reflection-and-beers-law/

DeGreve2006.pdf