% Simulations In One Dimension Example
%
% This example provides a demonstration of using k-Wave for the simulation
% and detection of photoacoustic waves within a one-dimensional
% heterogeneous propagation medium. It builds on the Homogeneous
% Propagation Medium and Heterogeneous Propagation Medium examples.
%
% author: Bradley Treeby
% date: 14th July 2009
% last update: 16th July 2009
%  
% This example is part of the k-Wave Toolbox (http://www.k-wave.org)

clear all;

% =========================================================================
% SIMULATION
% =========================================================================

% create the computational grid
Nx = 512;       % number of pixels in the x (column) direction
dx = 0.05e-3;   % pixel width [m]
kgrid = makeGrid(Nx, dx);

% define the properties of the propagation medium
c = 1500*ones(Nx, 1);           % [m/s]
c(1:round(Nx/3)) = 2000;        % [m/s]
rho = 1000*ones(Nx, 1);         % [kg/m^3]
rho(round(4*Nx/5):end) = 1500;  % [kg/m^3]

% create initial pressure distribution using a smoothly shaped sinusoid
x_pos = 280;    % pixels
width = 100;    % pixels
height = 1;
in = 0:pi/(width/2):2*pi;
p0 = [zeros(x_pos,1); ((height/2)*sin(in-pi/2)+(height/2)).'; zeros(Nx - x_pos  - width - 1, 1)];

% create a Cartesian sensor mask
sensor_mask = [-10e-3, 10e-3];  % [mm]

% use makeTime to determine a suitable dt
[t_array dt] = makeTime(kgrid, max(c(:)));

% replace t_array with a longer array to capture the reflections
t_end = 2.5*kgrid.x_size/max(c(:));
t_array = 0:dt:t_end;

% run the simulation
sensor_data = kspaceFirstOrder1D(p0, kgrid, c, rho, t_array, sensor_mask, 'PlotLayout', true);

% =========================================================================
% VISUALISATION
% =========================================================================

% plot the recorded time signals
figure;
plot(sensor_data(1, :), 'b-');
hold on;
plot(sensor_data(2, :), 'r-');
axis tight;
set(gca, 'YLim', [-0.1 0.7]);
ylabel('Pressure');
xlabel('Time Step');
legend('Sensor Position 1', 'Sensor Position 2');