diff options
Diffstat (limited to 'src/photon.h')
-rw-r--r-- | src/photon.h | 15 |
1 files changed, 13 insertions, 2 deletions
diff --git a/src/photon.h b/src/photon.h index 4203a50..fa3b85b 100644 --- a/src/photon.h +++ b/src/photon.h @@ -1,6 +1,7 @@ #ifndef __PHOTON_H__ #define __PHOTON_H__ +#include "stdio.h" #include "linalg.h" #include "materials.h" #include "rotate.h" @@ -50,6 +51,11 @@ enum enum {BREAK, CONTINUE, PASS}; // return value from propagate_to_boundary +__device__ float get_theta(const float3 &a, const float3 &b) +{ + return acosf(fmaxf(-1.0f,fminf(1.0f,dot(a,b)))); +} + __device__ void fill_state(State &s, Photon &p) { p.last_hit_triangle = intersect_mesh(p.position, p.direction, s.distance_to_boundary, p.last_hit_triangle); @@ -97,6 +103,9 @@ __device__ void fill_state(State &s, Photon &p) s.absorption_length = interp_property(p.wavelength, material1.absorption_length); s.scattering_length = interp_property(p.wavelength, material1.scattering_length); + printf("wavelength = %f", p.wavelength); + printf("scattering length = %f\n", s.scattering_length); + } // fill_state __device__ void rayleigh_scatter(Photon &p, curandState &rng) @@ -147,6 +156,8 @@ __device__ int propagate_to_boundary(Photon &p, State &s, curandState &rng) { if (scattering_distance <= s.distance_to_boundary) { + printf("scattering distance = %f\n", scattering_distance); + p.time += scattering_distance/(SPEED_OF_LIGHT/s.refractive_index1); p.position += scattering_distance*p.direction; @@ -169,7 +180,7 @@ __device__ int propagate_to_boundary(Photon &p, State &s, curandState &rng) __device__ void propagate_at_boundary(Photon &p, State &s, curandState &rng) { - float incident_angle = acosf(dot(s.surface_normal, -p.direction)); + float incident_angle = get_theta(s.surface_normal,-p.direction); float refracted_angle = asinf(sinf(incident_angle)*s.refractive_index1/s.refractive_index2); float3 incident_plane_normal = cross(p.direction, s.surface_normal); @@ -260,7 +271,7 @@ __device__ int propagate_at_surface(Photon &p, State &s, curandState &rng) else { // specularly reflect - float incident_angle = acosf(dot(s.surface_normal, -p.direction)); + float incident_angle = get_theta(s.surface_normal,-p.direction); float3 incident_plane_normal = cross(p.direction, s.surface_normal); incident_plane_normal /= norm(incident_plane_normal); |