diff options
Diffstat (limited to 'src/kernel.cu')
-rw-r--r-- | src/kernel.cu | 238 |
1 files changed, 98 insertions, 140 deletions
diff --git a/src/kernel.cu b/src/kernel.cu index a36ff91..be0087d 100644 --- a/src/kernel.cu +++ b/src/kernel.cu @@ -9,14 +9,12 @@ #include "photon.h" #include "alpha.h" -#define RED_WAVELENGTH 685 -#define BLUE_WAVELENGTH 445 -#define GREEN_WAVELENGTH 545 - __device__ void fAtomicAdd(float *addr, float data) { while (data) + { data = atomicExch(addr, data+atomicExch(addr, 0.0f)); + } } __device__ __noinline__ void to_diffuse(Photon &p, State &s, curandState &rng, const int &max_steps) @@ -28,9 +26,9 @@ __device__ __noinline__ void to_diffuse(Photon &p, State &s, curandState &rng, c int command; - command = fill_state(s, p, rng); + fill_state(s, p); - if (command == BREAK) + if (p.last_hit_triangle == -1) break; command = propagate_to_boundary(p, s, rng); @@ -87,104 +85,76 @@ __global__ void rotate(int nthreads, float3 *points, float phi, float3 axis) points[id] = rotate(points[id], phi, axis); } - __global__ void build_rgb_lookup(int nthreads, float3 *positions, float3 *directions, float3 *rgb_lookup1, float3 *rgb_lookup2, unsigned long long seed, int runs, int max_steps) +__global__ void update_xyz_lookup(int nthreads, int total_threads, int offset, float3 position, curandState *rng_states, float wavelength, float3 xyz, float3 *xyz_lookup1, float3 *xyz_lookup2, int max_steps) { - int id = blockIdx.x*blockDim.x + threadIdx.x; - - if (id >= nthreads) - return; - - curandState rng; - curand_init(seed, id, 0, &rng); - - Photon p0; - p0.position = positions[id]; - p0.direction = directions[id]; - p0.direction /= norm(p0.direction); - p0.polarization = uniform_sphere(&rng); - p0.time = 0.0f; - p0.history = 0x0; + int kernel_id = blockIdx.x*blockDim.x + threadIdx.x; + int id = kernel_id + offset; - float distance; - - int hit_triangle = intersect_mesh(p0.position, p0.direction, distance); - - if (hit_triangle != id) + if (kernel_id >= nthreads || id >= total_threads) return; - // note triangles from built geometry mesh are in order + curandState rng = rng_states[kernel_id]; - uint4 triangle_data = g_triangles[hit_triangle]; + uint4 triangle_data = g_triangles[id]; float3 v0 = g_vertices[triangle_data.x]; float3 v1 = g_vertices[triangle_data.y]; float3 v2 = g_vertices[triangle_data.z]; - float cos_theta = dot(normalize(cross(v1-v0, v2-v1)), -p0.direction); + float a = curand_uniform(&rng); + float b = uniform(&rng, 0.0f, (1.0f - a)); + float c = 1.0f - a - b; - if (cos_theta < 0.0f) - cos_theta = dot(-normalize(cross(v1-v0, v2-v1)), -p0.direction); + float3 direction = a*v0 + b*v1 + c*v2 - position; + direction /= norm(direction); - Photon p; - State s; + float distance; + int triangle_index = intersect_mesh(position, direction, distance); - for (int i=0; i < runs; i++) + if (triangle_index != id) { - p = p0; - p.wavelength = RED_WAVELENGTH; + rng_states[kernel_id] = rng; + return; + } - to_diffuse(p, s, rng, max_steps); + float cos_theta = dot(normalize(cross(v1-v0,v2-v1)),-direction); - if (p.history & REFLECT_DIFFUSE) - { - if (s.inside_to_outside) - { - fAtomicAdd(&rgb_lookup1[p.last_hit_triangle].x, cos_theta); - } - else - { - fAtomicAdd(&rgb_lookup2[p.last_hit_triangle].x, cos_theta); - } - } + if (cos_theta < 0.0f) + cos_theta = dot(-normalize(cross(v1-v0,v2-v1)),-direction); - p = p0; - p.wavelength = GREEN_WAVELENGTH; + Photon p; + p.position = position; + p.direction = direction; + p.wavelength = wavelength; + p.polarization = uniform_sphere(&rng); + p.last_hit_triangle = -1; + p.time = 0; + p.history = 0; - to_diffuse(p, s, rng, max_steps); + State s; + to_diffuse(p, s, rng, max_steps); - if (p.history & REFLECT_DIFFUSE) + if (p.history & REFLECT_DIFFUSE) + { + if (s.inside_to_outside) { - if (s.inside_to_outside) - { - fAtomicAdd(&rgb_lookup1[p.last_hit_triangle].y, cos_theta); - } - else - { - fAtomicAdd(&rgb_lookup2[p.last_hit_triangle].y, cos_theta); - } + fAtomicAdd(&xyz_lookup1[p.last_hit_triangle].x, cos_theta*xyz.x); + fAtomicAdd(&xyz_lookup1[p.last_hit_triangle].y, cos_theta*xyz.y); + fAtomicAdd(&xyz_lookup1[p.last_hit_triangle].z, cos_theta*xyz.z); } - - p = p0; - p.wavelength = BLUE_WAVELENGTH; - - to_diffuse(p, s, rng, max_steps); - - if (p.history & REFLECT_DIFFUSE) + else { - if (s.inside_to_outside) - { - fAtomicAdd(&rgb_lookup1[p.last_hit_triangle].z, cos_theta); - } - else - { - fAtomicAdd(&rgb_lookup2[p.last_hit_triangle].z, cos_theta); - } + fAtomicAdd(&xyz_lookup2[p.last_hit_triangle].x, cos_theta*xyz.x); + fAtomicAdd(&xyz_lookup2[p.last_hit_triangle].y, cos_theta*xyz.y); + fAtomicAdd(&xyz_lookup2[p.last_hit_triangle].z, cos_theta*xyz.z); } } -} // build_rgb_lookup + rng_states[kernel_id] = rng; -__global__ void render(int nthreads, curandState *rng_states, float3 *positions, float3 *directions, float *red, float *green, float *blue, float3 *rgb_lookup1, float3 *rgb_lookup2, int max_steps) +} // update_xyz_lookup + +__global__ void update_xyz_image(int nthreads, curandState *rng_states, float3 *positions, float3 *directions, float wavelength, float3 xyz, float3 *xyz_lookup1, float3 *xyz_lookup2, float3 *image, int nlookup_calls, int max_steps) { int id = blockIdx.x*blockDim.x + threadIdx.x; @@ -193,81 +163,69 @@ __global__ void render(int nthreads, curandState *rng_states, float3 *positions, curandState rng = rng_states[id]; - Photon p0; - p0.position = positions[id]; - p0.direction = directions[id]; - p0.direction /= norm(p0.direction); - p0.polarization = uniform_sphere(&rng); - p0.time = 0.0f; - p0.history = 0x0; + Photon p; + p.position = positions[id]; + p.direction = directions[id]; + p.direction /= norm(p.direction); + p.wavelength = wavelength; + p.polarization = uniform_sphere(&rng); + p.last_hit_triangle = -1; + p.time = 0; + p.history = 0; State s; - Photon p = p0; - p.wavelength = RED_WAVELENGTH; - to_diffuse(p, s, rng, max_steps); if (p.history & REFLECT_DIFFUSE) { if (s.inside_to_outside) { - red[id] = rgb_lookup1[p.last_hit_triangle].x; - } - else - { - red[id] = rgb_lookup2[p.last_hit_triangle].x; - } - } - else - { - red[id] = 0.0; - } - - p = p0; - p.wavelength = GREEN_WAVELENGTH; - - to_diffuse(p, s, rng, max_steps); - - if (p.history & REFLECT_DIFFUSE) - { - if (s.inside_to_outside) - { - green[id] = rgb_lookup1[p.last_hit_triangle].y; + image[id].x += xyz.x*xyz_lookup1[p.last_hit_triangle].x/nlookup_calls; + image[id].y += xyz.y*xyz_lookup1[p.last_hit_triangle].y/nlookup_calls; + image[id].z += xyz.z*xyz_lookup1[p.last_hit_triangle].z/nlookup_calls; } else { - green[id] = rgb_lookup2[p.last_hit_triangle].y; + image[id].x += xyz.x*xyz_lookup2[p.last_hit_triangle].x/nlookup_calls; + image[id].y += xyz.y*xyz_lookup2[p.last_hit_triangle].y/nlookup_calls; + image[id].z += xyz.z*xyz_lookup2[p.last_hit_triangle].z/nlookup_calls; } } - else - { - green[id] = 0.0; - } - - p = p0; - p.wavelength = BLUE_WAVELENGTH; - - to_diffuse(p, s, rng, max_steps); - - if (p.history & REFLECT_DIFFUSE) - { - if (s.inside_to_outside) - { - blue[id] = rgb_lookup1[p.last_hit_triangle].z; - } - else - { - blue[id] = rgb_lookup2[p.last_hit_triangle].z; - } - } - else - { - blue[id] = 0.0; - } - + rng_states[id] = rng; -} // render +} // update_xyz_image + +__global__ void process_image(int nthreads, float3 *image, int *pixels, int nimages) +{ + int id = blockIdx.x*blockDim.x + threadIdx.x; + + if (id >= nthreads) + return; + + float3 rgb = image[id]/nimages; + + if (rgb.x < 0.0f) + rgb.x = 0.0f; + if (rgb.y < 0.0f) + rgb.y = 0.0f; + if (rgb.z < 0.0f) + rgb.z = 0.0f; + + if (rgb.x > 1.0f) + rgb.x = 1.0f; + if (rgb.y > 1.0f) + rgb.y = 1.0f; + if (rgb.z > 1.0f) + rgb.z = 1.0f; + + unsigned int r = floorf(rgb.x*255.0f); + unsigned int g = floorf(rgb.y*255.0f); + unsigned int b = floorf(rgb.z*255.0f); + + pixels[id] = r << 16 | g << 8 | b; + +} // process_image /* Trace the rays starting at `positions` traveling in the direction `directions` to their intersection with the global mesh. If the ray @@ -322,9 +280,9 @@ __global__ void propagate(int nthreads, curandState *rng_states, float3 *positio int command; - command = fill_state(s, p, rng); + fill_state(s, p); - if (command == BREAK) + if (p.last_hit_triangle == -1) break; command = propagate_to_boundary(p, s, rng); |