Fix bug that caused photons to NAN_ABORT if they hit a triangle

at exactly normal incidence.

The plane of incidence was undefined in that case, but should have
been the plane normal to polarization vector.

1 files changed, 9 insertions, 1 deletions

diff --git a/src/photon.h b/src/photon.h
index f76ca54..45ad79b 100644
--- a/src/photon.h
+++ b/src/photon.h
@@ -179,7 +179,15 @@ __device__ void propagate_at_boundary(Photon &p, State &s, curandState &rng)
 	float refracted_angle = asinf(sinf(incident_angle)*s.refractive_index1/s.refractive_index2);
 
 	float3 incident_plane_normal = cross(p.direction, s.surface_normal);
-	incident_plane_normal /= norm(incident_plane_normal);
+	float incident_plane_normal_length = norm(incident_plane_normal);
+
+	// Photons at normal incidence do not have a unique plane of incidence,
+	// so we have to pick the plane normal to be the polarization vector
+	// to get the correct logic below
+	if (incident_plane_normal_length < 1e-6f)
+	  incident_plane_normal = p.polarization;
+	else
+	  incident_plane_normal /= incident_plane_normal_length;
 
 	float normal_coefficient = dot(p.polarization, incident_plane_normal);
 	float normal_probability = normal_coefficient*normal_coefficient;