function [ density, sound_speed ] = makeMediumLayers3D(kgrid, thck_wall, thck_abs, rho_abs)
%Create 3 different medium including rigid walls and absorbent material on the walls
%Author: Pavel Horak
%       makeMediumLayers3D( Nx, Ny, Nz,thck_wall, thck_abs, rho_abs )
%       Nx,Ny,Nz,dx,dy,dz   -      [grid points]
%       thck_wall           -      [m]
%       thck_abs            -      [m]
%       rho_abs             -      [kg/m^3]
Nx=kgrid.Nx;
Ny=kgrid.Ny;
Nz=kgrid.Nz;
dx=kgrid.dx;
dy=kgrid.dy;
dz=kgrid.dz;
%fixed constants:
rho_air = 1.225;    % [kg/m^3]
c_air=343;          % [m/s]
rho_wall=1000;      % [kg/m^3] 2300
c_wall=1700;        % [kg/m^3] 1700
c_abs=c_air;

if thck_wall==0
sound_speed = c_air*ones(Nx, Ny, Nz);
density = rho_air*ones(Nx, Ny, Nz);
else
%convert metres to gridpoints
thck_wallx=round(thck_wall/dx);
thck_wally=round(thck_wall/dy);
thck_wallz=round(thck_wall/dz);
thck_absx=round(thck_abs/dx);
thck_absy=round(thck_abs/dy);
thck_absz=round(thck_abs/dz);


%air
sound_speed = c_air*ones(Nx, Ny, Nz);
density = rho_air*ones(Nx, Ny, Nz);
%walls
wall = zeros(Nx, Ny,Nz);
wall(1:thck_wallx, :,:) = 1;
wall(:, 1:thck_wally,:) = 1;
wall(:,:, 1:thck_wallz) = 1;
wall(end - thck_wallx + 1:end, :,:) = 1;
wall(:, end - thck_wally + 1:end,:) = 1;
wall(:,:, end - thck_wallz + 1:end) = 1;

%absorber
absorber = zeros(Nx, Ny,Nz);
absorber(thck_wallx:thck_wallx+thck_absx, :,:) = 1;
absorber(:, thck_wally:thck_wally+thck_absy,:) = 1;
absorber(:,:, thck_wallz:thck_wallz+thck_absz) = 1;
absorber(end -thck_wallx-thck_absx + 1:end-thck_wallx+1, :,:) = 1;
absorber(:, end -thck_wallx-thck_absx + 1:end-thck_wallx+1,:) = 1;
absorber(:,:, end -thck_wallx-thck_absx + 1:end-thck_wallx+1) = 1;

%define medium
sound_speed(wall == 1) = c_wall;
density(wall == 1) = rho_wall;
sound_speed(absorber == 1) = c_abs;
density(absorber == 1) = rho_abs;
end
end


