ogl_beamforming

Ultrasound Beamforming Implemented with OpenGL
git clone anongit@rnpnr.xyz:ogl_beamforming.git
Log | Files | Refs | Feed | Submodules | LICENSE

hadamard.glsl (2472B)


      1 /* See LICENSE for license details. */
      2 layout(local_size_x = 32, local_size_y = 32, local_size_z = 1) in;
      3 
      4 layout(std430, binding = 1) readonly restrict buffer buffer_1 {
      5 	int rf_data[];
      6 };
      7 
      8 layout(std430, binding = 2) writeonly restrict buffer buffer_2 {
      9 	vec2 out_data[];
     10 };
     11 
     12 layout(r8i, binding = 0) readonly restrict uniform iimage2D hadamard;
     13 
     14 void main()
     15 {
     16 	/* NOTE: each invocation takes a time sample and a receive channel.
     17 	 * It first maps the column to the correct column in the rf data then
     18 	 * does the dot product with the equivalent row of the hadamard matrix.
     19 	 * The result is stored to the equivalent row, column index of the output.
     20 	 */
     21 	uint time_sample = gl_GlobalInvocationID.x;
     22 	uint channel     = gl_GlobalInvocationID.y;
     23 	uint acq         = gl_GlobalInvocationID.z;
     24 
     25 	/* NOTE: offsets for storing the results in the output data */
     26 	uint out_off = dec_data_dim.x * dec_data_dim.y * acq + dec_data_dim.x * channel + time_sample;
     27 
     28 	/* NOTE: channel mapping is stored as u16s so we must do this to extract the final value */
     29 	uint ch_array_idx = (channel / 8);
     30 	uint ch_vec_idx   = (channel % 8) / 2;
     31 	uint ch_elem_lfs  = ((~channel) & 1u) * 16;
     32 	uint rf_channel   = (channel_mapping[ch_array_idx][ch_vec_idx] << ch_elem_lfs) >> 16;
     33 
     34 	/* NOTE: stride is the number of samples between acquistions; off is the
     35 	 * index of the first acquisition for this channel and time sample  */
     36 	uint rf_stride = dec_data_dim.x;
     37 	uint rf_off    = rf_raw_dim.x * rf_channel + time_sample;
     38 
     39 	/* NOTE: rf_data index and stride considering the data is i16 not i32 */
     40 	uint ridx       = rf_off / 2;
     41 	uint ridx_delta = rf_stride / 2;
     42 
     43 	/* NOTE: rf_data is i16 so each access grabs two time samples at time.
     44 	 * We need to shift arithmetically (maintaining the sign) to get the
     45 	 * desired element. If the time sample is even we take the upper half
     46 	 * and if its odd we take the lower half. */
     47 	uint lfs = ((~time_sample) & 1u) * 16;
     48 
     49 	/* NOTE: Compute N-D dot product */
     50 	int sum = 0;
     51 	switch (decode) {
     52 	case DECODE_MODE_NONE: {
     53 		ridx += ridx_delta * acq;
     54 		sum   = (rf_data[ridx] << lfs) >> 16;
     55 	} break;
     56 	case DECODE_MODE_HADAMARD: {
     57 		/* NOTE: offset to get the correct column in hadamard matrix */
     58 		uint hoff = dec_data_dim.z * acq;
     59 
     60 		for (int i = 0; i < dec_data_dim.z; i++) {
     61 			int data = (rf_data[ridx] << lfs) >> 16;
     62 			sum  += imageLoad(hadamard, ivec2(i, acq)).x * data;
     63 			ridx += ridx_delta;
     64 		}
     65 	} break;
     66 	}
     67 	out_data[out_off] = vec2(float(sum), 0);
     68 }