ogl_beamforming

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

Commit: a4c7561750ffdfaf52ef28da28e44f3756d8b254
Parent: 2cd2cab85c3516a6f6ec071a6043d3684f602b3b
Author: Tyler Henry
Date:   Thu, 30 Jul 2026 12:24:05 -0600

core: add basic readi support

Diffstat:
Mbeamformer.meta | 24++++++++++++++++--------
Mbeamformer_core.c | 22++++++++++++++++------
Mbeamformer_internal.h | 5++++-
Mgenerated/beamformer.c | 35+++++++++++++++++++++++++----------
Mshaders/das.glsl | 48+++++++++++++++++++++++++++++++++++++++++++++++-
5 files changed, 108 insertions(+), 26 deletions(-)

diff --git a/beamformer.meta b/beamformer.meta @@ -1,11 +1,12 @@ -@Constant(16) ChunkChannelCount -@Constant(4) FilterSlots -@Constant(4096) MaxBacklogFrames -@Constant(256) MaxChannelCount -@Constant(256) MaxEmissionsCount -@Constant(16) MaxComputeShaderStages -@Constant(16) MaxParameterBlocks -@Constant(3) MaxRawDataFramesInFlight +@Constant(16) ChunkChannelCount +@Constant(4) FilterSlots +@Constant(4096) MaxBacklogFrames +@Constant(256) MaxChannelCount +@Constant(256) MaxEmissionsCount +@Constant(16) MaxComputeShaderStages +@Constant(16) MaxParameterBlocks +@Constant(3) MaxRawDataFramesInFlight +@Constant(65536) MaxHadamardElements @Enumeration ShaderResourceKind { @@ -203,6 +204,8 @@ { [contrast_mode ContrastMode ] [emission_parameters EmissionParameters] + [readi_group_count U32] + [readi_group U32] } @Struct Parameters @@ -269,6 +272,7 @@ [focal_vectors V2 MaxChannelCount] [sparse_elements S16 MaxChannelCount] [transmit_receive_orientations U16 MaxChannelCount] + [hadamard_matrix F16 MaxHadamardElements] } @Emit @@ -404,6 +408,7 @@ @Shader(das.glsl) DAS { @Constant MaxChannelCount + @Constant MaxHadamardElements @Enumeration AcquisitionKind @Enumeration InterpolationMode @@ -439,6 +444,8 @@ [SingleFocus B32] [FocusDepth F32] [TransmitAngle F32] + + [ReadiGroupCount U32] } @PushConstants @@ -455,6 +462,7 @@ [output_size_z U32] [cycle_t U32] [channel_offset S32] + [readi_group U32] } } diff --git a/beamformer_core.c b/beamformer_core.c @@ -271,14 +271,17 @@ das_valid_points(iv3 points) } function void -beamformer_update_hadamard(BeamformerComputePlan *cp, i32 order, b32 row_major, Arena arena) +beamformer_update_hadamard(BeamformerComputePlan *cp, i32 order, b32 row_major, Arena arena, b32 das_matrix) { f16 *hadamard = make_hadamard_transpose(&arena, order, row_major); if (hadamard) { - u64 offset = offsetof(BeamformerComputeArrayParameters, Hadamard); - u64 size = sizeof(*((BeamformerComputeArrayParameters *)0)->Hadamard) * order * order; + u64 offset = das_matrix ? offsetof(BeamformerComputeArrayParameters, DasHadamard) + : offsetof(BeamformerComputeArrayParameters, DecodeHadamard); + u64 size = das_matrix ? sizeof(*((BeamformerComputeArrayParameters *)0)->DasHadamard) + : sizeof(*((BeamformerComputeArrayParameters *)0)->DecodeHadamard); + size *= order * order; vk_buffer_range_upload(&cp->array_parameters, hadamard, offset, size, 0); - cp->hadamard_order = order; + if(!das_matrix) cp->hadamard_order = order; } } @@ -737,8 +740,11 @@ plan_compute_pipeline(BeamformerComputePlan *cp, BeamformerParameterBlock *pb, A db->InterpolationMode = pb->parameters.interpolation_mode; db->TransmitAngle = pb->parameters.focal_vector.E[0]; db->FocusDepth = pb->parameters.focal_vector.E[1]; + db->ReadiGroupCount = pb->parameters.readi_group_count; db->TransmitReceiveOrientation = pb->parameters.transmit_receive_orientation; + cp->readi_group = pb->parameters.readi_group; + // NOTE(rnp): old gcc will miscompile an assignment memory_copy(cp->xdc_transform.E, pb->parameters.xdc_transform.E, sizeof(cp->xdc_transform)); @@ -1060,7 +1066,10 @@ beamformer_commit_parameter_block(BeamformerCtx *ctx, BeamformerComputePlan *cp, if (pb->parameters.decode_mode != BeamformerDecodeMode_None && cp->hadamard_order != (i32)cp->acquisition_count) { - beamformer_update_hadamard(cp, (i32)cp->acquisition_count, vk_gpu_info()->cooperative_matrix, arena); + beamformer_update_hadamard(cp, (i32)cp->acquisition_count, vk_gpu_info()->cooperative_matrix, arena, false); + + if (pb->parameters.readi_group_count > 1) + beamformer_update_hadamard(cp, (i32)pb->parameters.readi_group_count, false, arena, true); } }break; @@ -1131,7 +1140,7 @@ do_compute_shader(BeamformerCtx *ctx, VulkanHandle cmd, BeamformerComputePlan *c case BeamformerShaderKind_Decode:{ BeamformerDecodePushConstants pc = { - .hadamard_buffer = cp->array_parameters.gpu_pointer + offsetof(BeamformerComputeArrayParameters, Hadamard), + .hadamard_buffer = cp->array_parameters.gpu_pointer + offsetof(BeamformerComputeArrayParameters, DecodeHadamard), .rf_buffer = pp_input_pointer, }; @@ -1238,6 +1247,7 @@ do_compute_shader(BeamformerCtx *ctx, VulkanHandle cmd, BeamformerComputePlan *c .output_size_z = cp->output_points.z, .cycle_t = das_cycle_t++, .channel_offset = channel_offset, + .readi_group = cp->readi_group, .array_parameters = cp->array_parameters.gpu_pointer + offsetof(BeamformerComputeArrayParameters, FocalVectors), }; memory_copy(pc.voxel_transform.E, cp->das_voxel_transform.E, sizeof(pc.voxel_transform)); diff --git a/beamformer_internal.h b/beamformer_internal.h @@ -266,10 +266,11 @@ typedef struct { // X(kind, format, elements) #define BEAMFORMER_COMPUTE_ARRAY_PARAMETERS_LIST \ - X(Hadamard, f16, BeamformerMaxChannelCount * BeamformerMaxChannelCount) \ + X(DecodeHadamard, f16, BeamformerMaxChannelCount * BeamformerMaxChannelCount) \ X(FocalVectors, v2, BeamformerMaxChannelCount) \ X(SparseElements, i16, BeamformerMaxChannelCount) \ X(TransmitReceiveOrientations, u16, BeamformerMaxChannelCount) \ + X(DasHadamard, f16, BeamformerMaxChannelCount * BeamformerMaxChannelCount) \ typedef enum { #define X(k, ...) BeamformerComputeArrayParameterKind_##k, @@ -325,6 +326,8 @@ struct BeamformerComputePlan { v2 xdc_element_pitch; m4 xdc_transform; + u32 readi_group; + GPUBuffer array_parameters; BeamformerFilter filters[BeamformerFilterSlots]; diff --git a/generated/beamformer.c b/generated/beamformer.c @@ -11,6 +11,7 @@ #define BeamformerMaxComputeShaderStages (16) #define BeamformerMaxParameterBlocks (16) #define BeamformerMaxRawDataFramesInFlight (3) +#define BeamformerMaxHadamardElements (65536) typedef enum { BeamformerShaderResourceKind_Buffer = 0, @@ -212,6 +213,7 @@ typedef struct { b32 SingleFocus; f32 FocusDepth; f32 TransmitAngle; + u32 ReadiGroupCount; } BeamformerDASBakeParameters; typedef struct { @@ -251,6 +253,7 @@ typedef struct { u32 output_size_z; u32 cycle_t; i32 channel_offset; + u32 readi_group; } BeamformerDASPushConstants; typedef struct { @@ -364,6 +367,8 @@ typedef struct { typedef struct { BeamformerContrastMode contrast_mode; BeamformerEmissionParameters emission_parameters; + u32 readi_group_count; + u32 readi_group; } BeamformerExtraParameters; typedef struct { @@ -392,6 +397,8 @@ typedef struct { u32 decimation_rate; BeamformerContrastMode contrast_mode; BeamformerEmissionParameters emission_parameters; + u32 readi_group_count; + u32 readi_group; } BeamformerParameters; typedef struct { @@ -420,6 +427,8 @@ typedef struct { u32 decimation_rate; BeamformerContrastMode contrast_mode; BeamformerEmissionParameters emission_parameters; + u32 readi_group_count; + u32 readi_group; i16 channel_mapping[BeamformerMaxChannelCount]; i16 sparse_elements[BeamformerMaxEmissionsCount]; u8 transmit_receive_orientations[BeamformerMaxEmissionsCount]; @@ -448,6 +457,7 @@ typedef struct { v2 focal_vectors[BeamformerMaxChannelCount]; i16 sparse_elements[BeamformerMaxChannelCount]; u16 transmit_receive_orientations[BeamformerMaxChannelCount]; + f16 hadamard_matrix[BeamformerMaxHadamardElements]; } BeamformerDASArrayParameters; typedef union { @@ -629,6 +639,7 @@ read_only global MetaStructMember *meta_struct_members_by_id[] = { {14, 56, 1, 0}, {8, 60, 1, 0}, {8, 64, 1, 0}, + {18, 68, 1, 0}, }, (MetaStructMember []){ {18, 0, 1, 0}, @@ -689,6 +700,7 @@ read_only global str8 *meta_struct_member_names_by_id[] = { str8_comp("SingleFocus"), str8_comp("FocusDepth"), str8_comp("TransmitAngle"), + str8_comp("ReadiGroupCount"), }, (str8 []){ str8_comp("SizeX"), @@ -706,7 +718,7 @@ read_only global str8 *meta_struct_member_names_by_id[] = { read_only global MetaStructInfo meta_struct_info_by_id[] = { {str8_comp("DecodeBakeParameters"), 11, 44, 0}, {str8_comp("FilterBakeParameters"), 12, 48, 0}, - {str8_comp("DASBakeParameters"), 17, 68, 0}, + {str8_comp("DASBakeParameters"), 18, 72, 0}, {str8_comp("ReshapeBakeParameters"), 9, 36, 0}, }; @@ -809,6 +821,7 @@ read_only global str8 beamformer_shader_global_header_strings[] = { "};\n" "\n"), str8_comp("#define MaxChannelCount (256)\n\n"), + str8_comp("#define MaxHadamardElements (65536)\n\n"), str8_comp("" "#define AcquisitionKind_FORCES 0\n" "#define AcquisitionKind_UFORCES 1\n" @@ -836,9 +849,10 @@ read_only global str8 beamformer_shader_global_header_strings[] = { "\n"), str8_comp("" "struct DASArrayParameters {\n" - " f32vec2 focal_vectors[MaxChannelCount];\n" - " int16_t sparse_elements[MaxChannelCount];\n" - " uint16_t transmit_receive_orientations[MaxChannelCount];\n" + " f32vec2 focal_vectors[MaxChannelCount];\n" + " int16_t sparse_elements[MaxChannelCount];\n" + " uint16_t transmit_receive_orientations[MaxChannelCount];\n" + " float16_t hadamard_matrix[MaxHadamardElements];\n" "};\n" "\n"), str8_comp("" @@ -858,6 +872,7 @@ read_only global str8 beamformer_shader_global_header_strings[] = { " uint32_t output_size_z;\n" " uint32_t cycle_t;\n" " int32_t channel_offset;\n" + " uint32_t readi_group;\n" "};\n" "\n"), str8_comp("" @@ -933,18 +948,18 @@ read_only global b8 beamformer_shader_primitive_is_vertex[] = { read_only global i32 *beamformer_shader_header_vectors[] = { (i32 []){0, 1, 2}, (i32 []){3, 4, 5, 6}, - (i32 []){7, 8, 9, 10, 3, 4, 11, 12, 13}, - (i32 []){14}, - 0, + (i32 []){7, 8, 9, 10, 11, 3, 4, 12, 13, 14}, (i32 []){15}, - (i32 []){16, 17}, - (i32 []){18}, + 0, + (i32 []){16}, + (i32 []){17, 18}, + (i32 []){19}, }; read_only global i32 beamformer_shader_header_vector_lengths[] = { 3, 4, - 9, + 10, 1, 0, 1, diff --git a/shaders/das.glsl b/shaders/das.glsl @@ -320,6 +320,51 @@ RESULT_TYPE FORCES(const vec3 xdc_world_point) return result; } +RESULT_TYPE READI_FORCES(const vec3 xdc_world_point) +{ + RESULT_TYPE result = RESULT_TYPE(0); + + float z_delta_squared = xdc_world_point.z * xdc_world_point.z; + float transmit_y_delta = xdc_world_point.y - xdc_element_pitch.y * ChannelCount / 2; + float transmit_yz_squared = transmit_y_delta * transmit_y_delta + z_delta_squared; + + // NOTE(tkh): The row we use matches the acquisition group, the column is the element group we are beamforming. + s32 hadamard_offset = s32(readi_group) * s32(ReadiGroupCount); + + for (f32 chunk_channel = 0; chunk_channel < f32(ChunkChannelCount); chunk_channel += 1.0f) { + float rx_channel = float(channel_offset) + chunk_channel; + float receive_x_delta = xdc_world_point.x - rx_channel * xdc_element_pitch.x; + float a_arg = abs(FNumber * receive_x_delta / xdc_world_point.z); + + if (a_arg < 0.5f) { + s32 channel_rf_offset = s32(rf_element_offset) + s32(chunk_channel) * SampleCount * AcquisitionCount; + channel_rf_offset -= s32(InterpolationMode == InterpolationMode_Cubic); + + float receive_index = sample_index(sqrt(receive_x_delta * receive_x_delta + z_delta_squared)); + float apodization = apodize(a_arg); + + // NOTE(tkh): Iterating over groups of tx elements, each group is AcquisitionCount sequential elements. + // The first element in each group is beamformed using the first acquisition, the second element in each group is beamformed using the second acquisition, etc. + for (s32 tx_group = 0; tx_group < s32(ReadiGroupCount); tx_group++) { + s32 rf_offset = channel_rf_offset; + f16 hadamard_value = ArrayParameters(array_parameters).data.hadamard_matrix[hadamard_offset + tx_group]; + float group_apodization = apodization * f32(hadamard_value); + + for (s32 tx_event = 0; tx_event < AcquisitionCount; tx_event++) { + s32 tx_element = tx_group * AcquisitionCount + tx_event; + float transmit_x_delta = xdc_world_point.x - xdc_element_pitch.x * tx_element; + float transmit_index = sqrt(transmit_yz_squared + transmit_x_delta * transmit_x_delta) * SamplingFrequency / SpeedOfSound; + + SAMPLE_TYPE value = group_apodization * sample_rf(rf_offset, receive_index + transmit_index); + result += RESULT_STORE(value); + rf_offset += SampleCount; + } + } + } + } + return result; +} + void main() { uvec3 out_voxel = gl_GlobalInvocationID; @@ -337,7 +382,8 @@ void main() case AcquisitionKind_FORCES: case AcquisitionKind_UFORCES: { - sum = FORCES(world_point); + sum = ReadiGroupCount > 1 ? READI_FORCES(world_point) + : FORCES(world_point); }break; case AcquisitionKind_HERCULES: case AcquisitionKind_UHERCULES: