Hello everybody,
I am a beginner in the ultrasound field and I just started to use k-wave to modelize my experimental setup. I have some difficulties to modelize a focused single element transducer.
I created a line of source and defined a pression source (so to modelize the transducer as a piston). In front of the line I defined an acoustical lens ( so I define two different media on the grid) to create the focus of the transducer (as the hardware is).
It is working but for me it doens't reproduce the behaviour of the transducer because every point of the line generate a pressure wave in every direction and not only along one axis.
Is there a way to define the source.p as a directional source ? Is there another way to modelize it ?
Thank you in advance,
Ju
Here is my code :
%% create the computational grid
Nx = 216; % number of grid points in the x (row) direction
Ny = 216; % number of grid points in the y (column) direction
dx = 50e-3/Nx; % grid point spacing in the x direction [m]
dy = dx; % grid point spacing in the y direction [m]
kgrid = makeGrid(Nx, dx, Ny, dy);
%% define the properties of the propagation medium
% medium 1 : water
medium.sound_speed = 1500*ones(Nx, Ny); % [m/s]
medium.density = 1000*ones(Nx, Ny); % [kg/m^3]
% medium 2 : lens ( create = put 1 for the lens limit + depth)
lens_y_offset = 10; %[grid points]
lens_radius = 100; %[grid points]
lens_center_x = Nx/2; %[grid points]
lens_center_y = lens_radius+lens_y_offset; %[grid points]
lens_radians = pi/6;
Z = makeCircle(Nx, Ny, lens_center_x, lens_center_y, lens_radius, lens_radians);
% Find the limit of the half lens
[Wu_x,Wu_y] = find(Z==1);
Wu_x_sorted = sort(Wu_x);
Wu_y_sorted = sort(Wu_y,'descend');
%Loop to fill with one the half lens depth + set the density and the sound
%speed within the lens
for i = 1: length(Wu_x_sorted)
Z(Wu_x_sorted(i),min(Wu_y_sorted):Wu_y_sorted(i))=1;
medium.sound_speed(Wu_x_sorted(i),min(Wu_y_sorted):Wu_y_sorted(i)) = 1700; % [m/s]
medium.density(Wu_x_sorted(i),min(Wu_y_sorted):Wu_y_sorted(i)) = 1200; % [kg/m^3]
end
%Loop to create the symmetrical of the half lens
for i = 1:round(lens_radius*sin(lens_radians))
Z(lens_center_x+i,:) = Z(lens_center_x-i,:);
medium.sound_speed(lens_center_x+i,:) = medium.sound_speed(lens_center_x-i,:); % [m/s]
medium.density(lens_center_x+i,:) = medium.density(lens_center_x-i,:); % [kg/m^3]
end
% create the time array
[kgrid.t_array, dt] = makeTime(kgrid, medium.sound_speed);
%% define a line along X axis transducer element
transducer_line_length = lens_radius ; % [grid points]
transducer_center_x = Nx/2; %[grid points]
transducer_center_y = lens_y_offset; %[grid points]
source.p_mask = zeros(Nx,Ny);
source.p_mask(transducer_center_x-(transducer_line_length/2):transducer_center_x+(transducer_line_length/2)-1,transducer_center_y) = 1;
% define a time varying sinusoidal source
source_freq = 0.25e6; % [Hz]
source_mag = 0.5; % [Pa]
source.p = source_mag*sin(2*pi*source_freq*kgrid.t_array);
% filter the source to remove high frequencies not supported by the grid
source.p = filterTimeSeries(kgrid, medium, source.p);
% create a display mask to display the transducer
display_mask = source.p_mask;
%% create a sensor mask covering the entire computational domain using the
% opposing corners of a rectangle
sensor.mask = [1, 1, Nx, Ny].';
%% set the record mode capture the final wave-field and the statistics at
% each sensor point
sensor.record = {'p_final', 'p_max', 'p_rms'};
%% assign the input options
input_args = {'DisplayMask', display_mask, 'PMLInside', false,'PlotLayout', true, 'PlotPML', false};
% run the simulation
sensor_data = kspaceFirstOrder2D(kgrid, medium, source, sensor, input_args{:});