function [] = multipleLayers(varargin)
OPT = struct(...
    'x', 1,...
    'dir', cd,...
    'nd',20,...
    'dzbed',0.1,...
    'pbbed',1,...
    'zb',1,...
    'XB',1);

OPT = setproperty(OPT, varargin{:});
currentdir = cd;
%initial time
i=1;
% ngd = 2; 
% nd = 20; 
% dzbed = 0.1;

% x = OPT.x(1,:)';

for i = 1:OPT.XB.nt 
     for m = 1:OPT.nd 
        zc(:,:,m,i) = zb(:,:,i)-(m-1)*OPT.dzbed; 
     end 
end 

% zb = squeeze(zb(:,2,:));
% zc = squeeze(zc(:,2,:,:));
% pbbed = squeeze(pbbed(:,2,:,:,:));
% % pcolor(repmat(x,[1 20]),zc(:,:,i),pbbed(:,:,1,i))
% shading flat

% pbbed = xb_readvar('pbbed.dat',OPT.XB);
sed1 = squeeze(pbbed(:,:,:,1,end));
sed2 = squeeze(pbbed(:,:,:,2,end));

temp1 = cat(1,sed1(:,:,1),sed1(:,:,2));
temp2 = cat(1,sed2(:,:,1),sed2(:,:,2));
for i = 3:OPT.nd
    temp1 = cat(1,temp1,sed1(:,:,i));
    temp2 = cat(1,temp2,sed2(:,:,i));
end
temp1 = temp1(:);
temp2 = temp2(:);

cd(dir)
save gdist1.inp temp1 -ascii
save gdist2.inp temp2 -ascii
cd(currentdir)
end

clear all; clc

x = load(fullfile('Z:','Checkouts','OpenEarthModels','tudelft','Chandeleurs_Islands','InputData','Grids','x.grd'));
x = x(1,:)';
XB = xb_getdimensions;
zb = xb_readvar('zb.dat',XB);
pbbed = xb_readvar('pbbed.dat',XB);

ngd = 2; 
nd = 20; 
dzbed = 0.1;
for i = 1:XB.nt 
     for m = 1:nd 
        zc(:,:,m,i) = zb(:,:,i)-(m-1)*0.1; 
     end 
end 

zb = squeeze(zb(:,2,:));
zc = squeeze(zc(:,2,:,:));
pbbed = squeeze(pbbed(:,2,:,:,:));
figure();
pcolor(repmat(x,[1 20]),zc(:,:,i),pbbed(:,:,1,i))
shading flat

