ogl_beamforming

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

Commit: b149668559f6045fa4f87bcafa266f64cecdbf00
Parent: 2ee3075a97ba51a835629182d8964da1afaa9b16
Author: Randy Palamar
Date:   Fri, 21 Jun 2024 21:46:19 -0600

do cubic interp between time samples

Diffstat:
Mmain.c | 2--
Mshaders/render.glsl | 2+-
Mshaders/uforces.glsl | 41+++++++++++++++++++++++++----------------
3 files changed, 26 insertions(+), 19 deletions(-)

diff --git a/main.c b/main.c @@ -157,7 +157,6 @@ init_fragment_shader_ctx(FragmentShaderCtx *ctx, uv3 out_data_dim) ctx->shader = LoadShader(NULL, "shaders/render.glsl"); ctx->out_data_dim_id = glGetUniformLocation(ctx->shader.id, "u_out_data_dim"); - glUniform3uiv(ctx->out_data_dim_id, 1, out_data_dim.E); /* TODO: add min max uniform */ /* output texture for image blitting */ @@ -227,7 +226,6 @@ reload_shaders(BeamformerCtx *ctx, Arena a) UnloadShader(ctx->fsctx.shader); ctx->fsctx.shader = updated_fs; ctx->fsctx.out_data_dim_id = GetShaderLocation(updated_fs, "u_out_data_dim"); - glUniform3uiv(ctx->fsctx.out_data_dim_id, 1, ctx->out_data_dim.E); } } diff --git a/shaders/render.glsl b/shaders/render.glsl @@ -22,7 +22,7 @@ vec3 hsv2rgb(vec3 hsv) void main() { ivec2 coord = ivec2(fragTexCoord * u_out_data_dim.xy); - float smp = out_data[coord.y * u_out_data_dim.x + coord.x]; + float smp = out_data[coord.y * u_out_data_dim.x + coord.x]; smp = 20 * log(abs(smp)) + 50; v_out_colour = vec4(hsv2rgb(vec3(smp, 0.8, 0.95)), 1); diff --git a/shaders/uforces.glsl b/shaders/uforces.glsl @@ -14,51 +14,60 @@ layout(location = 4) uniform uvec3 u_out_data_dim; layout(location = 5) uniform float u_sound_speed = 1452; layout(location = 6) uniform float u_sampling_frequency = 2.0833e7; layout(location = 7) uniform float u_focal_depth = 0.07; -//layout(location = 10) uniform sampler2D u_element_positions; +//layout(location = 9) uniform sampler2D u_element_positions; + +float cubic(float a, float b, float t) +{ + float param = -2 * t * t * t + 3 * t * t; + return (1 - param) * a + param * b; +} void main() { - vec3 scale = vec3(u_rf_data_dim) / vec3(u_out_data_dim); vec2 pixel = vec2(gl_GlobalInvocationID.xy); - ivec3 rf_coord = ivec3(gl_GlobalInvocationID.xyz * scale); - ivec3 out_coord = ivec3(gl_GlobalInvocationID.xyz); + ivec2 out_coord = ivec2(gl_GlobalInvocationID.xy); /* NOTE: Convert pixel to physical coordinates */ /* TODO: Send these in like the 3D program */ - - vec2 xdc_upper_left = vec2(-0.0096, -0.0096); - vec2 xdc_bottom_right = vec2( 0.0096, 0.0096); //vec2 xdc_upper_left = texture(u_element_positions, ivec2(0, 0)).xy; //vec2 xdc_bottom_right = texture(u_element_positions, ivec2(1, 1)).xy; + vec2 xdc_upper_left = vec2(-0.0096, -0.0096); + vec2 xdc_bottom_right = vec2( 0.0096, 0.0096); vec2 xdc_size = abs(xdc_upper_left - xdc_bottom_right); + /* TODO: image extent can be different than xdc_size */ /* TODO: for now assume y-dimension is along transducer center */ vec3 image_point = vec3( xdc_upper_left.x + pixel.x * xdc_size.x / u_out_data_dim.x, 0, - pixel.y * 80e-3 / u_out_data_dim.y + pixel.y * 60e-3 / u_out_data_dim.y + 10e-3 ); - float dx = xdc_size.x / float(u_rf_data_dim.y); - float dzsign = sign(image_point.z - u_focal_depth); - /* TODO: Send this into the GPU */ float sparse_elems[] = {17, 33, 49, 65, 80, 96, 112}; + float x = image_point.x - xdc_upper_left.x; + float dx = xdc_size.x / float(u_rf_data_dim.y); + float dzsign = sign(image_point.z - u_focal_depth); + float sum = 0; - float x = image_point.x - xdc_upper_left.x; /* NOTE: skip first acquisition since its garbage */ uint ridx = u_rf_data_dim.y * u_rf_data_dim.x; for (uint i = 1; i < u_rf_data_dim.z; i++) { - vec3 focal_point = vec3(sparse_elems[i - 1] * dx, 0, u_focal_depth); + vec3 focal_point = vec3(sparse_elems[i - 1] * dx, 0, u_focal_depth); float transmit_dist = u_focal_depth + dzsign * distance(image_point, focal_point); vec2 rdist = vec2(x, image_point.z); for (uint j = 0; j < u_rf_data_dim.y; j++) { - float dist = transmit_dist + length(rdist); - uint rx = uint(dist * u_sampling_frequency / u_sound_speed); + float dist = transmit_dist + length(rdist); + float rsample = dist * u_sampling_frequency / u_sound_speed; + + /* NOTE: do cubic interp between adjacent time samples */ + uint rx1 = uint(floor(rsample)); + uint rx2 = uint(ceil(rsample)); + float param = rsample - float(rx1); + sum += cubic(rf_data[ridx + rx1], rf_data[ridx + rx2], param); - sum += rf_data[ridx + rx]; rdist.x -= dx; ridx += u_rf_data_dim.x; }