diff options
Diffstat (limited to 'src/daq.cu')
| -rw-r--r-- | src/daq.cu | 198 |
1 files changed, 0 insertions, 198 deletions
diff --git a/src/daq.cu b/src/daq.cu deleted file mode 100644 index 5a68846..0000000 --- a/src/daq.cu +++ /dev/null @@ -1,198 +0,0 @@ -// -*-c++-*- -#include <curand_kernel.h> - -__device__ unsigned int -float_to_sortable_int(float f) -{ - return __float_as_int(f); - //int i = __float_as_int(f); - //unsigned int mask = -(int)(i >> 31) | 0x80000000; - //return i ^ mask; -} - -__device__ float -sortable_int_to_float(unsigned int i) -{ - return __int_as_float(i); - //unsigned int mask = ((i >> 31) - 1) | 0x80000000; - //return __int_as_float(i ^ mask); -} - -extern "C" -{ - -__global__ void -reset_earliest_time_int(float maxtime, int ntime_ints, unsigned int *time_ints) -{ - int id = threadIdx.x + blockDim.x * blockIdx.x; - if (id < ntime_ints) { - unsigned int maxtime_int = float_to_sortable_int(maxtime); - time_ints[id] = maxtime_int; - } -} - -__global__ void -run_daq(curandState *s, unsigned int detection_state, float time_rms, - int first_photon, int nphotons, float *photon_times, - unsigned int *photon_histories, int *last_hit_triangles, - int *solid_map, int nsolids, unsigned int *earliest_time_int, - unsigned int *channel_q, unsigned int *channel_histories) -{ - - int id = threadIdx.x + blockDim.x * blockIdx.x; - - if (id < nphotons) { - curandState rng = s[id]; - int photon_id = id + first_photon; - int triangle_id = last_hit_triangles[photon_id]; - - if (triangle_id > -1) { - int solid_id = solid_map[triangle_id]; - unsigned int history = photon_histories[photon_id]; - - if (solid_id < nsolids && (history & detection_state)) { - float time = photon_times[photon_id] + curand_normal(&rng) * time_rms; - unsigned int time_int = float_to_sortable_int(time); - - atomicMin(earliest_time_int + solid_id, time_int); - atomicAdd(channel_q + solid_id, 1); - atomicOr(channel_histories + solid_id, history); - } - - } - - s[id] = rng; - - } - -} - -__global__ void -convert_sortable_int_to_float(int n, unsigned int *sortable_ints, - float *float_output) -{ - int id = threadIdx.x + blockDim.x * blockIdx.x; - - if (id < n) - float_output[id] = sortable_int_to_float(sortable_ints[id]); -} - -__global__ void -bin_hits(int nchannels, unsigned int *channel_q, float *channel_time, - unsigned int *hitcount, int tbins, float tmin, float tmax, int qbins, - float qmin, float qmax, unsigned int *pdf) -{ - int id = threadIdx.x + blockDim.x * blockIdx.x; - - if (id >= nchannels) - return; - - unsigned int q = channel_q[id]; - float t = channel_time[id]; - - if (t < 1e8 && t >= tmin && t < tmax && q >= qmin && q < qmax) { - hitcount[id] += 1; - - int tbin = (t - tmin) / (tmax - tmin) * tbins; - int qbin = (q - qmin) / (qmax - qmin) * qbins; - - // row major order (channel, t, q) - int bin = id * (tbins * qbins) + tbin * qbins + qbin; - pdf[bin] += 1; - } -} - -__global__ void -accumulate_pdf_eval(int time_only, int nchannels, unsigned int *event_hit, - float *event_time, float *event_charge, float *mc_time, - unsigned int *mc_charge, // quirk of DAQ! - unsigned int *hitcount, unsigned int *bincount, - float min_twidth, float tmin, float tmax, - float min_qwidth, float qmin, float qmax, - float *nearest_mc, int min_bin_content) -{ - int id = threadIdx.x + blockDim.x * blockIdx.x; - - if (id >= nchannels) - return; - - // Was this channel hit in the Monte Carlo? - float channel_mc_time = mc_time[id]; - if (channel_mc_time >= 1e8) - return; // Nothing else to do - - // Is this channel inside the range of the PDF? - float distance; - int channel_bincount = bincount[id]; - if (time_only) { - if (channel_mc_time < tmin || channel_mc_time > tmax) - return; // Nothing else to do - - // This MC information is contained in the PDF - hitcount[id] += 1; - - // Was this channel hit in the event-of-interest? - int channel_event_hit = event_hit[id]; - if (!channel_event_hit) - return; // No need to update PDF value for unhit channel - - // Are we inside the minimum size bin? - float channel_event_time = event_time[id]; - distance = fabsf(channel_mc_time - channel_event_time); - if (distance < min_twidth/2.0) - bincount[id] = channel_bincount + 1; - - } - else { // time and charge PDF - float channel_mc_charge = mc_charge[id]; // int->float conversion because DAQ just returns an integer - - if (channel_mc_time < tmin || channel_mc_time > tmax || - channel_mc_charge < qmin || channel_mc_charge > qmax) - return; // Nothing else to do - - // This MC information is contained in the PDF - hitcount[id] += 1; - - // Was this channel hit in the event-of-interest? - int channel_event_hit = event_hit[id]; - if (!channel_event_hit) - return; // No need to update PDF value for unhit channel - - // Are we inside the minimum size bin? - float channel_event_time = event_time[id]; - float channel_event_charge = event_charge[id]; - float normalized_time_distance = fabsf(channel_event_time - channel_mc_time)/min_twidth/2.0; - float normalized_charge_distance = fabsf(channel_event_charge - channel_mc_charge)/min_qwidth/2.0; - distance = sqrt(normalized_time_distance*normalized_time_distance + normalized_charge_distance*normalized_charge_distance); - - if (distance < 1.0f) - bincount[id] = channel_bincount + 1; - } - - // Do we need to keep updating the nearest_mc list? - if (channel_bincount + 1 >= min_bin_content) - return; // No need to perform insertion into nearest_mc because bincount is a better estimate of the PDF value - - // insertion sort the distance into the array of nearest MC points - int offset = min_bin_content * id; - int entry = min_bin_content - 1; - - // If last entry less than new entry, nothing to update - if (distance > nearest_mc[offset + entry]) - return; - - // Find where to insert the new entry while sliding the rest - // to the right - entry--; - while (entry >= 0) { - if (nearest_mc[offset+entry] >= distance) - nearest_mc[offset+entry+1] = nearest_mc[offset+entry]; - else - break; - entry--; - } - - nearest_mc[offset+entry+1] = distance; -} - -} // extern "C" |