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.

1 files changed, 38 insertions, 0 deletions

diff --git a/src/propagate.cu b/src/propagate.cu
index 2d5183e..ecc32c0 100644
--- a/src/propagate.cu
+++ b/src/propagate.cu
@@ -9,6 +9,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,
 	  float3 *positions, float3 *directions,