matlab kwave模拟超声换能器学习笔记

持续更新中。。。 

d​​​​​​kwave官方示例http://www.k-wave.org/documentation/k-wave.php

kwave官方手册下载http://www.k-wave.org/manual/k-wave_user_manual_1.1.pdf

一、均匀介质传播的模拟

建立物理空间网格

% create the computational grid
Nx = 128;           % number of grid points in the x (row) direction
Ny = 128;           % number of grid points in the y (column) direction
dx = 0.1e-3;        % grid point spacing in the x direction [m]
dy = 0.1e-3;        % grid point spacing in the y direction [m]
kgrid = kWaveGrid(Nx, dx, Ny, dy);

官方手册中提到，可以使用一些计算来定义网格

在定义的grid最外层默认添加完美匹配层，该层是声波的吸收层，二维情况下占用定义格点的最外围20个单元，注意不要在此区域添加换能器。

定义传播介质

% define the properties of the propagation medium
medium.sound_speed = 1500;  % [m/s]
medium.alpha_coeff = 0.75;  % [dB/(MHz^y cm)]
medium.alpha_power = 1.5;

定义初始压力分布

% create initial pressure distribution using makeDisc
disc_magnitude = 5; % [Pa]
disc_x_pos = 50;    % [grid points]
disc_y_pos = 50;    % [grid points]
disc_radius = 8;    % [grid points]
disc_1 = disc_magnitude * makeDisc(Nx, Ny, disc_x_pos, disc_y_pos, disc_radius);

disc_magnitude = 3; % [Pa]
disc_x_pos = 80;    % [grid points]
disc_y_pos = 60;    % [grid points]
disc_radius = 5;    % [grid points]
disc_2 = disc_magnitude * makeDisc(Nx, Ny, disc_x_pos, disc_y_pos, disc_radius);

source.p0 = disc_1 + disc_2;

定义sensor mask

sensor mask是声压被记录的区域，可以用以下三种方式给出

1. 一个0-1二值矩阵来给出所有需要记录声压的格点（需要与上一步定义的grid大小一致）
2. 区域的两个对角坐标（仅对于线、矩形或cubic）
3. 坐标集合（的矩阵，D是物理空间的维数，坐标原点在定义的格点中心，也可以使用makeCartCircle等函数来创建坐标集合）

% define a centered circular sensor
sensor_radius = 4e-3;   % [m]
num_sensor_points = 50;
sensor.mask = makeCartCircle(sensor_radius, num_sensor_points);

运行模拟

仿真函数kspaceFirstOrder{1,2,3}D，需要定义完以上四个部分后可以调用。

图1: kwave仿真pipeline[1]

 

调用函数kspaceFirstOrder2D来进行仿真

% run the simulation
sensor_data = kspaceFirstOrder2D(kgrid, medium, source, sensor);

函数返回值sensor_data记录了sensor.mask中的时间序列，可以用索引sensor_data(sensor_point_index, time_index)取值。

sensor_data可视化代码

% plot the simulated sensor data
figure;
imagesc(sensor_data, [-1, 1]);
colormap(getColorMap);
ylabel('Sensor Position');
xlabel('Time Step');
colorbar;

Reference

[1] Bradley Treeby, Ben Cox, and Jiri Jaros, k-wave_user_manual_1.1, 2016