I have a simulation that imports SEM (https://imgur.com/8mFivA5) images of aerogels with different pore sizes, the pores may be filled with fluid or air. I want to define the region inside the pores with different values of density and sound speed. Is there a efficient way to do this?
Here is my code so far!
clearvars;
% create the computational grid
Nx = 200 ; % number of grid points in the x (row) direction
Ny = 200; % number of grid points in the y (column) direction
dx = 5e-6; % grid point spacing in the x direction [m]
dy = 5e-6; % grid point spacing in the y direction [m]
kgrid = kWaveGrid(Nx, dx, Ny, dy);
%load image
img1 = loadImage('x-ca-1.bmp');
%resize
img1 = resize(img1, [Nx,Ny]);
img_indices1 = find(img1 ==1);
% define the medium properties
%Medium%
medium.sound_speed = 1500*ones(Nx, Ny); % [m/s]
medium.density = 1000 * ones(Nx, Ny);
%aerogel
medium.sound_speed(img_indices1) = 300;
medium.density(img_indices1) = 500;
%Define source
mag = 10;
source.p0 = mag * makeLine(Nx, Ny,[10,99],[10,101]);
% define a 2D binary sensor mask in the shape of a line
x_offset1 = 25; % [grid points]
x_offset2 = 175;
width = 1; % [grid points]
sensor.mask = zeros(Nx, Ny);
sensor.mask(x_offset1, Ny/2 - width/2 + 1:Ny/2 + width/2) = 1;
sensor.mask(x_offset2, Ny/2 - width/2 + 1:Ny/2 + width/2) = 1;
% define the frequency response of the sensor elements
center_freq = 40e6; % [Hz]
bandwidth = 80; % [%]
sensor.frequency_response = [center_freq, bandwidth];
% set the input options
input_args = {'DataCast', 'gpuArray-single'};
%run simulations
[sensor_data] = kspaceFirstOrder2D(kgrid, medium, source, sensor, input_args{:});
% compute the amplitude spectra of the recorded time series
[f, sensor_data1] = spect(sensor_data(1,:), 1/kgrid.dt);
[f, sensor_data2] = spect(sensor_data(2,:), 1/kgrid.dt);
% plot the amplitude spectra
[f_sc, scale, prefix] = scaleSI(max(f));
plot(f * scale, sensor_data1, 'k-',...
f * scale, sensor_data2, 'r-');
xlabel(['Frequency [' prefix 'Hz]']);
ylabel('Amplitude [au]');
legend('sensor-1', 'sensor-2');