From 4c554d621d9b66c595397a0667d212bc4ea57be1 Mon Sep 17 00:00:00 2001 From: Stan Seibert Date: Sat, 10 Sep 2011 15:05:49 -0400 Subject: Add the ability to propagate the same photons multiple times on the the GPU, and run the DAQ multiple times on the same photons in a likelihood calculation. Propagating the same photons in a warp speeds up propagation by a factor of 3 (and we could do this even better if we wanted), and this improves the statistics in a likelihood evaluation quite a bit. Running the DAQ multiple times is also an inexpensive way to improve the quality of the PDF estimates. --- src/propagate.cu | 38 ++++++++++++++++++++++++++++++++++++++ 1 file changed, 38 insertions(+) (limited to 'src') diff --git a/src/propagate.cu b/src/propagate.cu index 2d5183e..ecc32c0 100644 --- a/src/propagate.cu +++ b/src/propagate.cu @@ -8,6 +8,44 @@ extern "C" { +__global__ void +photon_duplicate(int first_photon, int nthreads, + float3 *positions, float3 *directions, + float *wavelengths, float3 *polarizations, + float *times, unsigned int *histories, + int *last_hit_triangles, + int copies, int stride) +{ + int id = blockIdx.x*blockDim.x + threadIdx.x; + + if (id >= nthreads) + return; + + int photon_id = first_photon + id; + + Photon p; + p.position = positions[photon_id]; + p.direction = directions[photon_id]; + p.polarization = polarizations[photon_id]; + p.wavelength = wavelengths[photon_id]; + p.time = times[photon_id]; + p.last_hit_triangle = last_hit_triangles[photon_id]; + p.history = histories[photon_id]; + + for (int i=1; i <= copies; i++) { + int target_photon_id = photon_id + stride * i; + + positions[target_photon_id] = p.position; + directions[target_photon_id] = p.direction; + polarizations[target_photon_id] = p.polarization; + wavelengths[target_photon_id] = p.wavelength; + times[target_photon_id] = p.time; + last_hit_triangles[target_photon_id] = p.last_hit_triangle; + histories[target_photon_id] = p.history; + } +} + + __global__ void propagate(int first_photon, int nthreads, unsigned int *input_queue, unsigned int *output_queue, curandState *rng_states, -- cgit