
用Field II和MATLAB实现超声相控阵聚焦成像从零开始的工程实践指南超声相控阵技术作为现代医学影像领域的核心技术之一其聚焦成像能力直接决定了图像分辨率和诊断准确性。本文将带您从零开始通过Field II仿真平台和MATLAB编程环境完整实现超声相控阵的聚焦成像流程。不同于理论讲解我们聚焦于可执行的代码实践每个步骤都配有物理意义解释和常见错误排查方法确保您不仅能运行代码更能理解背后的工程逻辑。1. 环境准备与基础概念在开始编码前我们需要明确几个关键概念。相控阵超声成像通过精确控制各个阵元的发射时序相位延迟实现声波的定向发射和接收。与平面波成像相比聚焦成像在焦点区域能获得更高的信噪比和分辨率这是以帧率为代价换取的。1.1 Field II安装与配置Field II是丹麦理工大学开发的超声仿真工具箱MATLAB环境下运行。安装步骤从官网下载最新版Field II包将文件夹添加到MATLAB路径addpath(genpath(您的FieldII路径)); field_init; % 初始化工具箱验证安装field_info; % 应显示版本信息注意Field II使用前必须初始化且不同版本API可能有差异。本文基于2023版编写。1.2 核心参数物理意义超声成像质量取决于多个关键参数以下是需要重点理解的参数对照表参数符号物理意义典型值影响维度trans.fc中心频率5-10MHz分辨率/穿透深度trans.numele阵元数量64-256横向分辨率trans.pitch阵元间距0.1-0.3mm栅瓣伪影trans.c声速1540m/s深度校准userset.txfocus发射焦点深度30-60mm焦区分辨率2. 换能器建模与参数设置准确的换能器模型是仿真的基础。我们将构建一个64阵元的线性相控阵这是临床常用的配置。2.1 换能器几何定义trans.fc 6e6; % 6MHz中心频率 trans.numele 64; % 64个阵元 trans.width 136.9e-6; % 阵元宽度(米) trans.pitch 171.1e-6; % 阵元中心间距 trans.kerf trans.pitch - trans.width; % 阵元间间隙 trans.height 14e-3; % 阵元高度( elevation方向) trans.elevationFocus 60e-3; % 高程方向焦点深度 trans.c 1540; % 组织声速(m/s) % 计算阵元物理位置 trans.ElementPos trans.pitch * (-((trans.numele-1)/2):((trans.numele-1)/2));这段代码定义了换能器的基本几何特性。常见错误包括单位不一致米/毫米混用阵元位置计算错误导致不对称声速设置不当影响深度校准2.2 成像区域参数userset.theta -pi/4; % 扫描角度范围(-45度) userset.fs 100e6; % 采样频率100MHz userset.dep1 0e-3; % 成像起始深度 userset.dep2 80e-3; % 成像结束深度 userset.lat1 trans.ElementPos(1); % 横向起始位置 userset.lat2 trans.ElementPos(end); % 横向结束位置 userset.txfocus 30e-3; % 发射焦点深度30mm userset.Raynum 128; % 扫描线数量 userset.rayDelta 2*(-userset.theta)/(userset.Raynum-1); % 角度增量 userset.aperture 64*trans.pitch; % 有效孔径 trans.radius -(userset.aperture/2)/tan(-userset.theta); % 曲率半径3. 发射与接收波束形成聚焦成像的核心在于精确计算各阵元的延时模式。我们采用延时叠加(DAS)算法实现波束形成。3.1 延时计算原理对于每个扫描线需要计算从发射焦点到成像点的距离从成像点到各接收阵元的距离总传播时间决定采样时刻% 计算单条扫描线的延时模式 function delays calc_delay(xT, x_focus, z_focus, c) distances sqrt((xT - x_focus).^2 z_focus.^2); delays (distances - min(distances)) / c; end3.2 完整波束形成流程% 初始化发射和接收孔径 emit xdc_linear_array(trans.numele, trans.width, trans.height, ... trans.kerf, 1, 1, [0 0 0]); rcv xdc_linear_array(trans.numele, trans.width, trans.height, ... trans.kerf, 1, 1, [0 0 0]); % 计算各扫描线参数 Angles userset.theta:userset.rayDelta:(userset.theta (userset.Raynum-1)*userset.rayDelta); for i 1:userset.Raynum txorigin -(trans.radius * tan(Angles(i))); txsteer Angles(i); [x_focus, z_focus] calc_focal(txorigin, txsteer); % 设置发射延时 tx_dd calc_delay(trans.ElementPos, x_focus, z_focus, trans.c); xdc_apodization(emit, 0, ones(1, trans.numele)); xdc_center_focus(emit, [0 0 0]); xdc_focus_times(emit, 0, tx_dd); % 设置接收参数 xdc_apodization(rcv, 0, ones(1, trans.numele)); xdc_focus_times(rcv, 0, zeros(1, trans.numele)); % 仿真散射点回波 point_pos [0 0 10; 0 0 20; 0 0 30; 0 0 40; 0 0 50; 0 0 60; 0 0 70]/1000; point_amp 20*ones(1, size(point_pos,1)); [v_temp, tstart(i)] calc_scat_multi(emit, rcv, point_pos, point_amp); rf_data(1:size(v_temp,1), :, i) v_temp; end4. 图像重建与显示获得射频数据后需要通过波束形成和后期处理生成最终的B超图像。4.1 延时叠加算法实现function dasdata das(rf_an, tau, trans) [nsamples, nchannels] size(rf_an); dasdata zeros(size(xx,1), size(xx,2)); for k 1:length(x) % 计算发射距离 dTX sqrt(x(k)^2 (z(k)-trans.radius)^2) trans.radius; % 计算各通道接收距离 dRX sqrt((trans.ElementPos-x(k)).^2 z(k).^2); % 总延时 tau (dTX dRX)/trans.c max(tx_dd(ii)); % 插值获取各通道数据 for ch 1:nchannels t_idx round(tau(ch)*userset.fs); if t_idx 0 t_idx nsamples dasdata(k) dasdata(k) rf_an(t_idx, ch); end end end end4.2 图像后处理与显示% 希尔伯特变换获取包络 rf hilbert(rf_data); % 定义成像网格 [xx, zz] meshgrid(linspace(userset.lat1, userset.lat2, 256), ... linspace(userset.dep1, userset.dep2, 512)); x xx(:); z zz(:); % 执行DAS波束形成 dasdata zeros(size(xx)); for ii 1:userset.Raynum rf_an rf(:,:,ii); dasdata dasdata das(rf_an, tau, trans); end % 对数压缩和显示 FrameData abs(dasdata); maxd max(FrameData(:)); dBRange 60; % 动态范围60dB figure; imagesc(xx*1000, zz*1000, 20*log10(FrameData/maxd), [-dBRange, 0]); colormap(gray); colorbar; xlabel(Lateral Position (mm)); ylabel(Depth (mm)); title(超声相控阵聚焦成像结果);5. 性能优化与调试技巧实际工程实现中计算效率和图像质量需要平衡。以下是几个关键优化点5.1 计算加速策略并行计算将不同扫描线的处理分配到多个CPU核心parfor i 1:userset.Raynum % 各扫描线独立处理 endGPU加速将波束形成移植到GPUgpu_rf gpuArray(rf_data); % 在GPU上执行DAS5.2 常见问题排查图像出现伪影检查阵元位置定义是否正确验证延时计算是否考虑了声束偏转确认采样率是否满足Nyquist定理分辨率不足提高中心频率权衡穿透深度增加扫描线数量权衡帧率优化孔径大小和焦点位置深度标定错误确认声速参数与实际介质匹配检查时间基准校正验证采样时钟精度在完成基础实现后可以尝试以下进阶改进动态接收聚焦提升分辨率相干复合成像降低噪声自适应波束形成算法