从零构建CPWC超声成像仿真：Field II实战与模块化工作流解析

1. CPWC超声成像仿真入门指南第一次接触CPWC超声成像仿真时我被各种专业术语和复杂的数学公式搞得晕头转向。经过几个月的实战摸索终于总结出一套小白也能快速上手的方法。CPWC相干平面波复合是近年来超声成像领域的热门技术它通过发射多个角度的平面波再对回波信号进行相干叠加能够显著提升图像分辨率和信噪比。Field II作为超声仿真领域的瑞士军刀是进行CPWC仿真的不二之选。这个由丹麦理工大学开发的MATLAB工具箱已经成为了行业标准。不过对于新手来说Field II的文档确实不太友好很多关键细节需要自己摸索。下面我就把踩过的坑和实战经验分享给大家。整个仿真流程可以分为三个主要阶段前期准备环境配置和探头建模、核心计算平面波发射接收循环、后期处理数据封装和图像重建。每个阶段都有需要注意的细节比如在探头建模时阵元间距的设置会直接影响成像质量在核心计算阶段延迟时间的计算精度决定了最终的图像分辨率。2. 环境准备与探头建模2.1 Field II环境初始化在开始任何仿真工作前首先要确保Field II正确安装并配置好MATLAB环境。我建议创建一个专门的工作目录把所有相关脚本和数据都放在这里。初始化时先用field_init(0)命令清空之前的设置这是个好习惯能避免很多奇怪的bug。声速和采样频率是两个最关键的参数。在人体组织中声速通常设为1540 m/s这个值会影响所有与时间相关的计算。采样频率则要根据探头中心频率来选择一般至少是中心频率的4倍。比如对于7.5MHz的探头采样频率建议设置在30MHz以上。% Field II初始化示例代码 field_init(0); set_field(c, 1540); % 设置声速[m/s] set_field(fs, 100e6); % 设置采样频率[Hz] set_field(use_rectangles, 1); % 使用矩形阵元2.2 探头参数设置探头是超声成像系统的眼睛它的参数设置直接影响成像质量。以常用的L11-4v线性阵列探头为例它有128个阵元中心频率7.5MHz。每个阵元的宽度和高度需要精确设置通常宽度为0.2mm高度为5mm。阵元间距pitch也是一个关键参数一般等于或略大于阵元宽度以避免栅瓣效应。% 探头建模示例代码 probe uff.linear_array(); probe.N 128; % 阵元数量 probe.pitch 0.3e-3; % 阵元间距[m] probe.element_width 0.2e-3; % 阵元宽度[m] probe.element_height 5e-3; % 阵元高度[m] probe.f0 7.5e6; % 中心频率[Hz]在实际项目中我发现探头几何参数的微小误差都会导致仿真结果出现明显偏差。建议先用简单的点散射体测试探头设置是否正确确认无误后再进行复杂仿真。3. 脉冲建模与孔径定义3.1 发射脉冲设计超声成像的质量很大程度上取决于发射脉冲的特性。一个理想的脉冲应该具有较宽的带宽和良好的对称性。我通常使用高斯脉冲叠加方波激励的方式来生成发射脉冲这种方法简单有效参数调整也很直观。脉冲的分数带宽fractional bandwidth是个重要参数一般设置在0.6-0.7之间。太小的带宽会影响轴向分辨率太大的带宽则可能导致信号失真。脉冲持续时间通常设为2-3个周期既能保证足够的能量又不会导致分辨率下降。% 脉冲建模示例代码 pulse uff.pulse(); pulse.fractional_bandwidth 0.65; % 分数带宽 pulse.center_frequency probe.f0; % 中心频率 % 计算冲激响应 t0 (-1/pulse.fractional_bandwidth/probe.f0):dt:(1/pulse.fractional_bandwidth/probe.f0); impulse_response gauspuls(t0, probe.f0, pulse.fractional_bandwidth); impulse_response impulse_response - mean(impulse_response); % 去除直流分量3.2 发射与接收孔径在Field II中我们需要分别定义发射和接收孔径。发射孔径通常使用全阵列发射而接收孔径则可以灵活设置不同的接收模式。xdc_linear_array函数用于创建线性阵列孔径创建后需要设置激励信号、冲激响应等参数。一个容易忽略的细节是挡板baffle设置它模拟了探头外壳对声场的影响。在实际应用中我建议将挡板设为软挡板soft baffle这样更接近真实探头的声学特性。% 孔径定义示例代码 Th xdc_linear_array(probe.N, probe.element_width, probe.element_height, ... probe.kerf, 1, 1, [0 0 Inf]); % 发射孔径 Rh xdc_linear_array(probe.N, probe.element_width, probe.element_height, ... probe.kerf, 1, 1, [0 0 Inf]); % 接收孔径 % 设置发射孔径参数 xdc_excitation(Th, excitation); xdc_impulse(Th, impulse_response); xdc_baffle(Th, 0); % 软挡板 % 设置接收孔径参数 xdc_impulse(Rh, impulse_response); xdc_baffle(Rh, 0);4. 平面波序列与体模设置4.1 平面波角度计算CPWC技术的核心在于使用多个角度的平面波进行照射。角度范围的计算很关键通常由F数F-number决定。F数越小角度范围越大但可能导致边缘区域的信号质量下降。经过多次测试我发现F数在1.5-2.0之间通常能取得较好的平衡。角度数量也是一个需要权衡的参数。角度太少会导致复合效果不佳太多则会增加计算量。对于大多数应用场景11-21个角度是个不错的选择。在MATLAB中可以使用linspace函数均匀分布这些角度。% 平面波序列定义示例代码 F_number 1.7; % F数 alpha_max atan(1/(2*F_number)); % 最大角度[rad] Na 15; % 平面波数量 alpha linspace(-alpha_max, alpha_max, Na); % 角度向量4.2 虚拟体模创建体模是验证成像算法的重要工具。在仿真中我们通常使用点散射体来模拟组织。点散射体的位置和幅度可以根据需要灵活设置。为了测试系统的分辨率我建议设置一组间距逐渐减小的点散射体。Field II提供了calc_scat_multi函数来计算散射体产生的回波信号。在实际应用中要注意散射体的数量不宜过多否则会显著增加计算时间。通常100-200个散射体就能很好地验证系统性能。% 虚拟体模示例代码 point_position [linspace(-5e-3,5e-3,20); zeros(1,20); linspace(10e-3,15e-3,20)]; point_amplitudes ones(size(point_position,1),1); % 散射体幅度5. 核心计算与数据采集5.1 平面波发射接收循环这是整个仿真最耗时的部分。对于每个发射角度我们需要设置相应的发射延迟然后计算接收信号。发射延迟的计算公式为x*sin(α)/c其中x是阵元位置α是发射角度c是声速。这个延迟实现了平面波的倾斜发射。在接收端我们通常不设置动态聚焦而是保留原始的通道数据RF data以便后续处理。Field II的calc_scat_multi函数会返回一个矩阵包含所有接收阵元的时域信号。% 核心循环示例代码 CPW zeros(cropat, probe.N, Na); % 预分配内存 for n 1:Na % 设置发射延迟 delays probe.geometry(:,1) * sin(alpha(n)) / c0; xdc_focus_times(Th, 0, delays); % 计算回波信号 [v,t] calc_scat_multi(Th, Rh, point_position, point_amplitudes); % 存储数据 CPW(1:size(v,1),:,n) v; end5.2 数据封装与存储采集到的原始数据需要妥善封装以便后续处理。我推荐使用UFFUltrasound File Format格式存储数据这是一种开放的超声数据格式支持多种超声数据处理工具。在封装数据时要确保包含所有必要的元数据如探头参数、采样频率、发射序列等。这些信息对后续的图像重建和数据分析至关重要。一个好的做法是为每个数据集添加详细的注释说明采集条件和处理历史。% 数据封装示例代码 channel_data uff.channel_data(); channel_data.sampling_frequency fs; channel_data.sound_speed c0; channel_data.sequence seq; channel_data.data CPW;6. 图像重建与结果分析6.1 延迟叠加算法实现图像重建的核心是延迟叠加Delay and SumDAS算法。这个算法需要对每个图像点计算来自各个阵元的信号的正确延迟然后进行叠加。在CPWC中还需要对不同角度获取的数据进行相干复合。重建过程首先需要定义扫描区域scan通常是一个二维网格。然后设置适当的接收孔径窗函数如Tukey窗来控制旁瓣水平。F数在这里再次出现它决定了接收孔径的有效大小。% 图像重建示例代码 scan uff.linear_scan(); scan.x_axis linspace(-10e-3,10e-3,256); scan.z_axis linspace(5e-3,20e-3,256); pipe pipeline(); pipe.channel_data channel_data; pipe.scan scan; pipe.receive_apodization.window uff.window.tukey25; pipe.receive_apodization.f_number F_number; b_data pipe.go({midprocess.das() postprocess.coherent_compounding()});6.2 结果可视化与性能评估重建后的图像可以使用MATLAB的绘图功能进行可视化。为了客观评估图像质量我通常会计算以下几个指标分辨率通过点散射体的扩散程度评估、对比度通过均匀区域和病灶区域的信号差异评估和信噪比。在实际项目中我发现CPWC相比传统的聚焦成像确实能显著提升图像质量特别是在深部区域。但也要注意CPWC对运动非常敏感在存在组织运动的情况下可能需要额外的运动补偿算法。% 结果可视化示例代码 figure; b_data.plot(); title(CPWC重建图像); colorbar;经过完整的仿真流程后你会得到一组可以用于算法开发和性能评估的超声图像数据。这个过程虽然复杂但掌握了Field II的使用方法后你就能快速验证各种新的成像算法和探头设计了。