diff options
Diffstat (limited to 'src')
-rw-r--r-- | src/alpha.h | 89 | ||||
-rw-r--r-- | src/kernel.cu | 60 | ||||
-rw-r--r-- | src/sorting.h | 103 | ||||
-rw-r--r-- | src/transform.cu | 47 |
4 files changed, 210 insertions, 89 deletions
diff --git a/src/alpha.h b/src/alpha.h index e4a3a40..ac75834 100644 --- a/src/alpha.h +++ b/src/alpha.h @@ -4,57 +4,14 @@ #include "linalg.h" #include "intersect.h" #include "mesh.h" +#include "sorting.h" -template <class T> -__device__ void swap(T &a, T &b) -{ - T temp = a; - a = b; - b = temp; -} +#include "stdio.h" #define ALPHA_DEPTH 10 -struct HitList -{ - unsigned int size; - unsigned int indices[ALPHA_DEPTH]; - float distances[ALPHA_DEPTH]; -}; - -__device__ void add_to_hit_list(const float &distance, const int &index, HitList &h) -{ - unsigned int i; - if (h.size >= ALPHA_DEPTH) - { - if (distance > h.distances[ALPHA_DEPTH-1]) - return; - - i = ALPHA_DEPTH-1; - } - else - { - i = h.size; - h.size += 1; - } - - h.indices[i] = index; - h.distances[i] = distance; - - while (i > 0 && h.distances[i-1] > h.distances[i]) - { - swap(h.distances[i-1], h.distances[i]); - swap(h.indices[i-1], h.indices[i]); - - i -= 1; - } -} - __device__ int get_color_alpha(const float3 &origin, const float3& direction) { - HitList h; - h.size = 0; - float distance; if (!intersect_node(origin, direction, g_start_node, -1.0f)) @@ -69,6 +26,10 @@ __device__ int get_color_alpha(const float3 &origin, const float3& direction) unsigned int i; + float distances[ALPHA_DEPTH]; + unsigned int indices[ALPHA_DEPTH]; + unsigned int n=0; + do { unsigned int first_child = tex1Dfetch(node_map, *node); @@ -106,9 +67,26 @@ __device__ int get_color_alpha(const float3 &origin, const float3& direction) if (intersect_triangle(origin, direction, v0, v1, v2, distance)) { - add_to_hit_list(distance, i, h); + if (n < 1) + { + distances[0] = distance; + indices[0] = i; + } + else + { + unsigned long j = searchsorted(n, distances, distance); + + if (j <= ALPHA_DEPTH-1) + { + insert(ALPHA_DEPTH, distances, j, distance); + insert(ALPHA_DEPTH, indices, j, i); + } + } + + if (n < ALPHA_DEPTH) + n++; } - + } // triangle loop node--; @@ -118,16 +96,19 @@ __device__ int get_color_alpha(const float3 &origin, const float3& direction) } // while loop while (node != head); - if (h.size < 1) + if (n < 1) return 0; float scale = 1.0f; float fr = 0.0f; float fg = 0.0f; float fb = 0.0f; - for (i=0; i < h.size; i++) + unsigned int index; + for (i=0; i < n; i++) { - uint4 triangle_data = g_triangles[h.indices[i]]; + index = indices[i]; + + uint4 triangle_data = g_triangles[index]; float3 v0 = g_vertices[triangle_data.x]; float3 v1 = g_vertices[triangle_data.y]; @@ -138,11 +119,11 @@ __device__ int get_color_alpha(const float3 &origin, const float3& direction) if (cos_theta < 0.0f) cos_theta = dot(-normalize(cross(v1-v0,v2-v1)),-direction); - unsigned int r0 = 0xff & (g_colors[h.indices[i]] >> 16); - unsigned int g0 = 0xff & (g_colors[h.indices[i]] >> 8); - unsigned int b0 = 0xff & g_colors[h.indices[i]]; + unsigned int r0 = 0xff & (g_colors[index] >> 16); + unsigned int g0 = 0xff & (g_colors[index] >> 8); + unsigned int b0 = 0xff & g_colors[index]; - float alpha = (255 - (0xff & (g_colors[h.indices[i]] >> 24)))/255.0f; + float alpha = (255 - (0xff & (g_colors[index] >> 24)))/255.0f; fr += r0*scale*cos_theta*alpha; fg += g0*scale*cos_theta*alpha; diff --git a/src/kernel.cu b/src/kernel.cu index 307412e..4418305 100644 --- a/src/kernel.cu +++ b/src/kernel.cu @@ -62,41 +62,6 @@ __device__ void to_diffuse(Photon &p, State &s, curandState &rng, const int &max extern "C" { -/* Translate the points `a` by the vector `v` */ -__global__ void translate(int nthreads, float3 *a, float3 v) -{ - int id = blockIdx.x*blockDim.x + threadIdx.x; - - if (id >= nthreads) - return; - - a[id] += v; -} - -/* Rotate the points `a` through an angle `phi` counter-clockwise about the - axis `axis` (when looking towards +infinity). */ -__global__ void rotate(int nthreads, float3 *a, float phi, float3 axis) -{ - int id = blockIdx.x*blockDim.x + threadIdx.x; - - if (id >= nthreads) - return; - - a[id] = rotate(a[id], phi, axis); -} - -__global__ void rotate_around_point(int nthreads, float3 *a, float phi, float3 axis, float3 point) -{ - int id = blockIdx.x*blockDim.x + threadIdx.x; - - if (id >= nthreads) - return; - - a[id] -= point; - a[id] = rotate(a[id], phi, axis); - a[id] += point; -} - __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 kernel_id = blockIdx.x*blockDim.x + threadIdx.x; @@ -239,6 +204,28 @@ __global__ void process_image(int nthreads, float3 *image, int *pixels, int nima } // process_image +__global__ void distance_to_mesh(int nthreads, float3 *positions, float3 *directions, float *distances) +{ + int id = blockIdx.x*blockDim.x + threadIdx.x; + + if (id >= nthreads) + return; + + float3 position = positions[id]; + float3 direction = directions[id]; + direction /= norm(direction); + + float distance; + + int triangle_index = intersect_mesh(position, direction, distance); + + if (triangle_index == -1) + distances[id] = 1e9; + else + distances[id] = distance; + +} // distance_to_mesh + __global__ void ray_trace(int nthreads, float3 *positions, float3 *directions, int *pixels) { int id = blockIdx.x*blockDim.x + threadIdx.x; @@ -288,6 +275,9 @@ __global__ void ray_trace_alpha(int nthreads, float3 *positions, float3 *directi float3 direction = directions[id]; direction /= norm(direction); + bool hit; + float distance; + pixels[id] = get_color_alpha(position, direction); } // ray_trace diff --git a/src/sorting.h b/src/sorting.h new file mode 100644 index 0000000..ec16bbe --- /dev/null +++ b/src/sorting.h @@ -0,0 +1,103 @@ +template <class T> +__device__ void swap(T &a, T &b) +{ + T tmp = a; + a = b; + b = tmp; +} + +template <class T> +__device__ void reverse(int n, T *a) +{ + for (int i=0; i < n/2; i++) + swap(a[i],a[n-1-i]); +} + +template <class T> +__device__ void piksrt(int n, T *arr) +{ + int i,j; + T a; + + for (j=1; j < n; j++) + { + a = arr[j]; + i = j-1; + while (i >= 0 && arr[i] > a) + { + arr[i+1] = arr[i]; + i--; + } + arr[i+1] = a; + } +} + +template <class T, class U> +__device__ void piksrt2(int n, T *arr, U *brr) +{ + int i,j; + T a; + U b; + + for (j=1; j < n; j++) + { + a = arr[j]; + b = brr[j]; + i = j-1; + while (i >= 0 && arr[i] > a) + { + arr[i+1] = arr[i]; + brr[i+1] = brr[i]; + i--; + } + arr[i+1] = a; + brr[i+1] = b; + } +} + +template <class T> +__device__ unsigned long searchsorted(unsigned long n, T *arr, const T &x) +{ + unsigned long ju,jm,jl; + int ascnd; + + jl = 0; + ju = n; + + ascnd = (arr[n-1] > arr[0]); + + while (ju-jl > 1) + { + jm = (ju+jl) >> 1; + + if (x >= arr[jm] == ascnd) + jl = jm; + else + ju = jm; + } + + if (x <= arr[0]) + return 0; + else if (x == arr[n-1]) + return n-1; + else + return ju; +} + +template <class T> +__device__ void insert(unsigned long n, T *arr, unsigned long i, const T &x) +{ + unsigned long j; + for (j=n-1; j > i; j--) + arr[j] = arr[j-1]; + arr[i] = x; +} + +template <class T> +__device__ void add_sorted(unsigned long n, T *arr, const T &x) +{ + unsigned long i = searchsorted(n, arr, x); + + if (i < n) + insert(n, arr, i, x); +} diff --git a/src/transform.cu b/src/transform.cu new file mode 100644 index 0000000..57bd509 --- /dev/null +++ b/src/transform.cu @@ -0,0 +1,47 @@ +//-*-c-*- + +#include "linalg.h" +#include "rotate.h" + +extern "C" +{ + +/* Translate the points `a` by the vector `v` */ +__global__ void translate(int nthreads, float3 *a, float3 v) +{ + int id = blockIdx.x*blockDim.x + threadIdx.x; + + if (id >= nthreads) + return; + + a[id] += v; +} + +/* Rotate the points `a` through an angle `phi` counter-clockwise about the + axis `axis` (when looking towards +infinity). */ +__global__ void rotate(int nthreads, float3 *a, float phi, float3 axis) +{ + int id = blockIdx.x*blockDim.x + threadIdx.x; + + if (id >= nthreads) + return; + + a[id] = rotate(a[id], phi, axis); +} + +/* Rotate the points `a` through an angle `phi` counter-clockwise + (when looking towards +infinity along `axis`) about the axis defined + by the point `point` and the vector `axis` . */ +__global__ void rotate_around_point(int nthreads, float3 *a, float phi, float3 axis, float3 point) +{ + int id = blockIdx.x*blockDim.x + threadIdx.x; + + if (id >= nthreads) + return; + + a[id] -= point; + a[id] = rotate(a[id], phi, axis); + a[id] += point; +} + +} // extern "c" |