当前位置： 首页 > news >正文

超声原始数据重构成B扫成像的MATLAB实现

news 2025/7/26 15:20:10

基本原理

  1. 原始RF数据：超声探头接收的原始射频信号
  2. 包络检测：提取RF信号的幅度信息
  3. 对数压缩：压缩动态范围以适应显示
  4. 时间增益补偿(TGC)：补偿深度相关的信号衰减
  5. 扫描转换：将数据映射到图像坐标系

MATLAB代码

function ultrasound_bscan_reconstruction()% 超声B扫成像重建close all; clc;% 生成模拟超声RF数据[rf_data, fs, c] = simulate_ultrasound_data();% 显示原始RF数据display_raw_rf(rf_data, fs, c);% 重建B超图像b_mode_image = reconstruct_bscan(rf_data, fs, c);% 显示最终B超图像display_bscan_image(b_mode_image);
endfunction [rf_data, fs, c] = simulate_ultrasound_data()% 模拟超声RF数据生成c = 1540;               % 声速 (m/s)fs = 100e6;             % 采样频率 (Hz)num_samples = 1500;     % 每条扫描线采样点数num_lines = 128;        % 扫描线数量depth = 0.04;           % 成像深度 (m)% 创建时间和空间轴t = (0:num_samples-1)/fs;   % 时间轴z = t * c / 2;              % 深度轴% 生成扫描线位置 (线性阵列)x_pos = linspace(-0.02, 0.02, num_lines);% 初始化RF数据矩阵rf_data = zeros(num_samples, num_lines);% 添加组织模型中的散射点scatterers = [% x(m), z(m), 反射强度0.000, 0.010, 1.0;0.005, 0.015, 0.8;-0.004, 0.020, 0.9;0.000, 0.030, 1.2;0.008, 0.025, 0.7;-0.007, 0.018, 0.6;];% 添加高反射界面 (模拟器官边界)z_boundary = 0.012:0.002:0.038;for i = 1:length(z_boundary)scatterers = [scatterers; linspace(-0.015, 0.015, 20)', repmat(z_boundary(i), 20, 1), repmat(1.5, 20, 1)];end% 为每条扫描线生成RF信号for line = 1:num_lines% 探头脉冲响应 (高斯调制正弦波)t_pulse = (-10:10)/fs * 5;pulse = exp(-(t_pulse*1e6).^2) .* sin(2*pi*5e6*t_pulse);% 初始化扫描线信号line_signal = zeros(num_samples, 1);% 添加散射点响应for s = 1:size(scatterers, 1)dx = scatterers(s, 1) - x_pos(line);dz = scatterers(s, 2);distance = sqrt(dx^2 + dz^2);  % 到散射点的距离% 计算往返时间time_of_flight = 2 * distance / c;% 找到最近的采样点sample_idx = round(time_of_flight * fs) + 1;% 确保在有效范围内if sample_idx > 10 && sample_idx < num_samples - 10% 添加脉冲响应start_idx = sample_idx - 10;end_idx = sample_idx + 10;line_signal(start_idx:end_idx) = line_signal(start_idx:end_idx) + ...scatterers(s, 3) * pulse';endend% 添加噪声和系统增益变化line_signal = line_signal + 0.02 * randn(size(line_signal));rf_data(:, line) = line_signal;end% 添加系统增益随深度变化 (模拟TGC效果)tgc_gain = linspace(1, 5, num_samples)';rf_data = rf_data .* repmat(tgc_gain, 1, num_lines);
endfunction display_raw_rf(rf_data, fs, c)% 显示原始RF数据figure('Name', '原始RF数据', 'Position', [100, 100, 900, 600]);% 计算深度轴num_samples = size(rf_data, 1);t = (0:num_samples-1)/fs;depth = t * c / 2;% 显示单条RF信号subplot(2, 2, [1, 3]);plot(t*1e6, rf_data(:, round(size(rf_data, 2)/2)));xlabel('时间 (\mus)');ylabel('幅度');title('单条扫描线RF信号');grid on;% 显示RF数据图像subplot(2, 2, 2);imagesc([], depth*1000, rf_data);colormap(gray);colorbar;title('RF数据矩阵');ylabel('深度 (mm)');xlabel('扫描线索引');axis tight;% 显示RF数据幅度谱subplot(2, 2, 4);[pxx, f] = pwelch(rf_data(:, round(size(rf_data, 2)/2)), [], [], [], fs);plot(f/1e6, 10*log10(pxx));xlabel('频率 (MHz)');ylabel('功率谱密度 (dB/Hz)');title('RF信号频谱');grid on;xlim([0, 20]);
endfunction b_mode_image = reconstruct_bscan(rf_data, fs, c)% 重建B超图像% 参数设置f0 = 5e6;       % 中心频率 (Hz)bw = 0.7;       % 带宽系数% 获取数据尺寸[num_samples, num_lines] = size(rf_data);% 步骤1: 带通滤波rf_filtered = bandpass_filter(rf_data, fs, f0, bw);% 步骤2: 包络检测 (希尔伯特变换)analytic_signal = hilbert(rf_filtered);envelope = abs(analytic_signal);% 步骤3: 对数压缩log_compressed = 20 * log10(envelope + eps);  % 加eps避免log(0)% 步骤4: 时间增益补偿 (TGC)depth = ((0:num_samples-1)/fs) * c / 2;tgc = 1.0 + 1.5 * depth / max(depth);  % 线性增益补偿tgc_corrected = log_compressed .* tgc';% 步骤5: 降采样 (可选)dec_factor = 2;b_mode_downsampled = tgc_corrected(1:dec_factor:end, :);% 步骤6: 扫描转换 (此处为线性阵列，直接显示)b_mode_image = b_mode_downsampled;
endfunction filtered_data = bandpass_filter(data, fs, f0, bw)% 带通滤波器设计nyquist = fs/2;f_low = max(0, f0*(1 - bw/2));    % 低频截止f_high = min(nyquist, f0*(1 + bw/2)); % 高频截止% 设计滤波器[b, a] = butter(4, [f_low, f_high]/nyquist, 'bandpass');% 应用滤波器filtered_data = filtfilt(b, a, data);
endfunction display_bscan_image(b_mode_image)% 显示B超图像figure('Name', 'B超图像重建', 'Position', [100, 100, 1000, 800]);% 归一化图像img = b_mode_image - min(b_mode_image(:));img = img / max(img(:));% 显示图像imagesc(img);colormap(gray);title('重建的B超图像');xlabel('扫描线索引');ylabel('深度采样点');% 添加深度刻度 (示例值)depth_ticks = linspace(0, size(b_mode_image, 1), 6);depth_labels = round(linspace(0, 40, 6)); % 假设最大深度40mmset(gca, 'YTick', depth_ticks, 'YTickLabel', depth_labels);ylabel('深度 (mm)');% 添加颜色条colorbar;% 添加图像处理选项uicontrol('Style', 'pushbutton', 'String', '增强对比度', ...'Position', [20 20 100 30], 'Callback', @enhance_contrast);uicontrol('Style', 'pushbutton', 'String', '平滑处理', ...'Position', [140 20 100 30], 'Callback', @apply_smoothing);uicontrol('Style', 'pushbutton', 'String', '边缘增强', ...'Position', [260 20 100 30], 'Callback', @edge_enhancement);% 回调函数function enhance_contrast(~, ~)img_enhanced = imadjust(img);imagesc(img_enhanced);colormap(gray);title('对比度增强的B超图像');endfunction apply_smoothing(~, ~)h = fspecial('gaussian', [5 5], 1.5);img_smoothed = imfilter(img, h);imagesc(img_smoothed);colormap(gray);title('平滑处理后的B超图像');endfunction edge_enhancement(~, ~)img_edges = img - imgaussfilt(img, 2) + 0.5;img_edges(img_edges < 0) = 0;img_edges(img_edges > 1) = 1;imagesc(img_edges);colormap(gray);title('边缘增强的B超图像');end
end

1. 模拟超声原始数据生成

  • 创建模拟组织环境（散射点模型）
  • 生成超声脉冲响应（高斯调制正弦波）
  • 模拟扫描线采集过程
  • 添加时间增益补偿(TGC)和随机噪声

2. 原始RF数据显示

  • 单条扫描线RF信号波形
  • RF数据矩阵图像
  • RF信号频谱分析

3. B扫成像重建处理

  • 带通滤波：提取探头频带内的信号
  • 包络检测：使用希尔伯特变换提取信号幅度
  • 对数压缩：压缩动态范围（40-60dB）
  • 时间增益补偿：补偿深度相关的信号衰减
  • 降采样：减少数据量（可选步骤）

4. B超图像显示

  • 显示重建的B超图像
  • 添加深度刻度
  • 交互式图像处理按钮：1对比度增强 2平滑处理 3边缘增强

参考 超声原始数据的重构成B扫成像的方法 youwenfan.com/contentcsa/65820.html

重建过程关键步骤详解

  1. 带通滤波：

    • 目的：去除低频噪声和高频干扰
    • 方法：设计Butterworth带通滤波器
    • 实现：bandpass_filter函数

  2. 包络检测：

    • 目的：提取RF信号的幅度信息

    • 方法：希尔伯特变换获取解析信号

    • 数学表达：

      analytic_signal = hilbert(rf_signal);
envelope = abs(analytic_signal);

  3. 对数压缩：

    • 目的：将大动态范围压缩到适合显示的0-255范围

    • 方法：取对数并缩放

    • 公式：

      log_compressed = 20 * log10(envelope + eps);

  4. 时间增益补偿(TGC)：

    • 目的：补偿深度增加导致的信号衰减

    • 方法：应用与深度成正比的增益曲线

    • 实现：

      tgc = 1.0 + 1.5 * depth / max_depth;
tgc_corrected = log_compressed .* tgc';

  5. 扫描转换：

    • 目的：将数据映射到图像坐标系
    • 注意：本代码模拟线性阵列探头，直接显示数据矩阵
    • 对于扇形扫描，需要插值到笛卡尔坐标系

运行结果

运行此代码将生成四个图形窗口：

  1. 原始RF数据显示（包含RF信号波形、数据矩阵和频谱）
  2. 重建的B超图像（显示组织结构和散射点）
  3. 交互式界面（提供图像增强选项）

B超图像将显示模拟的器官边界和散射点，您可以使用界面上的按钮进行图像增强处理。

这个实现提供了超声成像的基本原理演示，实际临床应用中的系统会更加复杂，包含更多优化和校正步骤。

http://www.dtcms.com/a/294388.html

相关文章：

  • 【AI News | 20250722】每日AI进展
  • now能减少mysql的压力吗
  • 【Android】用 ViewPager2 + Fragment + TabLayout 实现标签页切换
  • linux性能调整和故障排查
  • LeetCode热题100--24. 两两交换链表中的节点--中等
  • Linux文件——Ext2文件系统(3)_软硬链接
  • Ubuntu 1804 编译ffmpeg qsv MediaSDK libva 遇到的问题记录
  • #Linux内存管理# 详细介绍madvise函数的工作原理
  • Elasticsearch(ES)安装
  • 分布式电商系统：缓存策略、负载均衡与容灾方案
  • 解决 Electron 中 window.open 打开新窗口的各种“坑”
  • Python 程序设计讲义（6）：Python 的基本用法——运算符与表达式
  • API 汇总：ONLYOFFICE 文档最近更新
  • 背包DP之0/1背包
  • 11-1 浅层神经网络及计算前向传播
  • 局部重要性注意力LIA,通过区域重要性图与门控机制实现高阶信息交互，自适应增强有用特征、抑制冗余信息，平衡模型性能与效率。
  • VR-Doh: 革新3D建模的虚拟现实体验
  • DPVR亮相青岛品牌日，崂山科创力量引领AI眼镜新浪潮
  • 基于PLC的轨检小车控制器设计
  • .NET-键控服务依赖注入
  • 【实战】Dify从0到100进阶--文档解读（13）API前端再开发
  • 苍穹外卖DAY11
  • 【LeetCode数据结构】栈和队列的应用——设计循环队列问题详解
  • 【后端】FastAPI的Pydantic 模型
  • Excel 将数据导入到SQLServer数据库
  • Java TCP 通信详解：从基础到实战，彻底掌握面向连接的网络编程
  • 通用表格识别技术的应用，深刻改变人们处理表格数据的方式
  • 如何最简单、通俗地理解Python的numpy库？
  • Ubuntu22.04.5 LTS安装与使用Docker
  • 【优选算法-多源 BFS】多源 BFS：解决多个起点的广度优先搜索