Hi everyone,
I am doing well in photoacoustic image reconstruction by using time reversal and interpolating data in fluid model based on "2D Time Reversal Reconstruction For A Circular Sensor Example". Hence, I am trying to do the same in the elastic model because the medium includes the skull. I encountered a problem in that the circular sensors emit waves in an axisymmetric manner and the reconstructed image is incorrect and axisymmetric. But everything went well by interpolating data and creating an equivalent continuous circle. The codes of sensor, source, time reversal are as follows.
Yours,
Zak
------------------------------------------------------------
% assign the source
sampling_freq = 1/kgrid.dt; % [Hz]
source.s_mask = makeDisc(Nx, Ny, disc_x_pos, disc_y_pos, disc_radius);
source.sxx = -disc_magnitude.* toneBurst(sampling_freq, tone_burst_freq, tone_burst_cycles,'Plot',false);
source.syy = source.sxx;
% define sensor
sensor_radius = 12e-2; % [m]
sensor_angle = 4*pi/2; % [rad]
sensor_pos = [0, 0]; % [m]
num_sensor_points = 40;
sensor.mask = makeCartCircle(sensor_radius, num_sensor_points, sensor_pos, sensor_angle);
% sensor.mask = makeCircle(Nx, Ny, sensor_pos(1,1), sensor_pos(1,2), round(sensor_radius/dx), sensor_angle);
% set the input options
input_args = {'Smooth', false, 'PMLInside', false, 'PlotPML', false,...
'PlotSim',true};
% run the elastic simulation
sensor_data_elastic = pstdElastic2D(kgrid, medium, source, sensor, input_args{:});
% =========================================================================
% ELASTIC RECONSTRUCTION
% =========================================================================
% redefine kgrid_recon
kgrid_recon = kWaveGrid(Nx, dx, Ny, dy);
kgrid_recon.makeTime(cp1,cfl);
% add noise to the recorded sensor data
signal_to_noise_ratio = 40; % [dB]
sensor_data_elastic = addNoise(sensor_data_elastic, signal_to_noise_ratio, 'peak');
% reset the initial pressure
source.sxx = 0;
source.syy = source.sxx;
% assign the time reversal data(elastic)
sensor.time_reversal_boundary_data = sensor_data_elastic;
% run the time-reversal reconstruction
p0_recon_elastic = pstdElastic2D(kgrid_recon, medium, source, sensor, input_args{:});
% create a binary sensor mask of an equivalent continuous circle
sensor_radius_grid_points = round(sensor_radius / kgrid_recon.dx);
binary_sensor_mask = makeCircle(kgrid_recon.Nx, kgrid_recon.Ny, ...
kgrid_recon.Nx/2 + 1, kgrid_recon.Ny/2 + 1, sensor_radius_grid_points, sensor_angle);
% assign to sensor structure
sensor.mask = binary_sensor_mask;
cart_sensor_mask = makeCartCircle(sensor_radius, num_sensor_points, sensor_pos, sensor_angle);
% interpolate data to remove the gaps and assign to sensor structure
sensor.time_reversal_boundary_data = interpCartData(kgrid_recon, sensor_data_elastic, cart_sensor_mask, binary_sensor_mask);
% run the time-reversal reconstruction
p0_recon_interp_elastic = pstdElastic2D(kgrid_recon, medium, source, sensor, input_args{:});------------------------------------------------------------------