Hi everyone,
I was working on the simulation of ultrasound propagation. I designed three homogeneous mediums. Their speed-of-sounds were 1500 m/s. Their densities were 1000 kg/m^3. Their alpha_power were set as 1.1. Their alpha_coeff were 0.002, 0.7, and 3.0, respectively.
A sensor was set to record transmitted signal.
I just found that with different alpha_coeff, the time-of-arrivals of the signals were different. The larger alpha_coeff was, the earlier the signal arrived.
The change was so significant that it was counter-intuitive. I wonder if dispersion can explain this phenomenon.
I would like to ask you, is this simulation result consistent with the physical phenomenon? If not, how should I correct my code?
The code has been attached below.
Thanks a lot.
% To compare TOFs of transmitted signals, with different [constant] att.s or wvfroms.
% * The 2-D numerical simulations are based on k-Wave toolbox
% * Emitted pulse: a 3-cycle tone-burst whose centre-frequency is 2 MHz
% * Emitter and receiver: two points about 25 cm apart, recording transmitted signal
% * Medium: homogeneous, whose attenuation power is 1.1, and attenuation coefficient
% is 0.002, 0.7 or 3.0 dB/cm/MHz^1.1
clc;clear
close all
%% Define Params
% default medium params
sos_ref = 1500; % [m/s]
density_ref = 1000; % [kg/m^3]
att_power = 1.1; % y
% define k-Wave Grid
Nx = 3072;
dx = 1e-4;
kgrid = kWaveGrid(Nx, dx, Nx, dx);
t_end = Nx * dx / sos_ref * 2;
kgrid.makeTime(sos_ref, 0.5, t_end); % CFL = 0.5
% different wvforms and att.s
small_att = 0.002; % dB/cm/MHz^y
big_att = 0.7;
bigger_att = 3.0;
wv_tb = toneBurst(1./kgrid.dt, 2e6, 3, 'SignalLength', kgrid.Nt);
%% Define numerical phantom
%------------define Medium One, no attenuation at all
medium_noAtt.sound_speed = sos_ref * ones(Nx);
medium_noAtt.density = density_ref * ones(Nx);
%------------define Medium Two, with a small atten. coeff.,
medium_smallAtt = medium_noAtt;
medium_smallAtt.alpha_power = att_power;
medium_smallAtt.alpha_coeff = small_att * ones(Nx);
%------------define Medium Three, with a big atten. coeff.,
medium_bigAtt = medium_smallAtt;
medium_bigAtt.alpha_coeff = big_att * ones(Nx);
%------------define Medium Three, with a big atten. coeff.,
medium_biggerAtt = medium_smallAtt;
medium_biggerAtt.alpha_coeff = bigger_att * ones(Nx);
%% Define source & sensor
source.p_mask = zeros(Nx, Nx);
source.p_mask(Nx / 2, 250) = 1;
sensor.mask = zeros(Nx, Nx)
sensor.mask(Nx / 2, Nx - 250) = 1;
% source.p not defined
%% Execute the simulations
input_args = {'PMLSize', 'auto',...
'PMLInside', false,...
'PlotPML', false,...
'DisplayMask', 'off',...
'DataPath', pwd,...
'DeleteData', false};
for medium_name = {'smallAtt', 'bigAtt', 'biggerAtt'}
cmd = ['medium = medium_', medium_name{1}, ';'];
eval(cmd)
source.p = wv_tb;
rf_data = kspaceFirstOrder2DC(kgrid, medium, source, sensor, input_args{:});
cmd = ['sig_', medium_name{1}, ' = rf_data;'];
eval(cmd)
end
%% Visulaization
f1 = figure('Units','centimeters', 'Position', [2, 2, 25, 40]);
plt_cnt = 0;
plt_c = 3;
plt_r = 1;
range_ = [160, 175] * 1e-6; % unit: s
for medium_name = { 'smallAtt', 'bigAtt', 'biggerAtt'}
cmd = ['sig = sig_', medium_name{1}, ';'];
eval(cmd)
env = abs(hilbert(sig));
% ----------
plt_cnt = plt_cnt + 1; subplot(plt_c, plt_r, plt_cnt)
plot(kgrid.t_array * 1e6, env, 'color','r', 'lineWidth', 1.5)
grid on
grid minor
xlabel('Time [\mus]')
xlim(range_ * 1e6)
title(medium_name{1}, 'Interpreter', 'none')
end
sgtitle('Signal Envelope')