% Showing the results (positon of the interface) of example 1
% TO 091224

load name            % get basename of this model
load(basename);      % load data for this model (generaged in mf_adapt

stitle=sprintf('%s: Interf. at T=',basename);  % initial title text

daysperyear=365.25;

%% Grid  xGr,yGr and zGr, ZETA and IBOUND are loaded

[xGr,yGr,xm,ym]=modelsize(xGr,yGr);

[NROW,NCOL,NLAY]=size(IBOUND);

H  =readDat([basename,'.HDS']);  % use only H.totim
ZTA=readBud([basename,'.ZTA']);  % get interface positions

NT   =length(ZTA);      % total number of time steps in ZTA file
NSURF=length(ZTA(1).label);

%% Figure of rotating interface
figure; hold on

xlabel('x [m]'); ylabel('y [m]'); zlabel('z [m]'); grid on;

set(gca,'xlim',xGr([1 end]),'ylim',sort(yGr([1 end])),'zlim',sort(zGr([1 end])));

for it=NT          % for all records in bud file (and H-file)
    for iSurf=1:length(ZTA(1).label)
        idx=strmatch(sprintf('ZETAPLANE%d',iSurf'),ZTA(it).label);
        values=ZTA(it).term{idx}(:,:);
        values(ym<ym(8),:)=NaN;
        hdl=surf(xm,ym,values);
   
    end
   title(sprintf('%s, T=%d d',basename,H(it).totim));
    view(10,10);
end

%% Figure of rotating interface
figure; hold on

stl='-:;.';

xlabel('x [m]'); ylabel('z [m]'); grid on;

set(gca,'xlim',xGr([1 end]),'ylim',sort(zGr([1 end])));

leg=[];
for it=NT          % for all records in bud file (and H-file)
    for iSurf=1:length(ZTA(1).label)
        idx=strmatch(sprintf('ZETAPLANE%d',iSurf'),ZTA(it).label);
        plot(xm,ZTA(it).term{idx}(ym==ym(9),:),['k' stl(iSurf)]);
        leg{end+1}=sprintf('it=%d, iSurf=%d',it,iSurf);
    end
    title(sprintf('%s, T=%d d',basename,H(it).totim));
end
legend(leg);

%% Computed positions using ZTA
figure; hold on

xlabel('x [m]'); ylabel('y [m]'); grid on;

set(gca,'xlim',xGr([1 end]),'ylim',sort(yGr([1 end])),'zlim',sort(zGr([1 end])));

conts=-(5:5:50);
for it=NT          % for all records in bud file (and H-file)
    for iSurf=1:length(ZTA(1).label)
        idx=strmatch(sprintf('ZETAPLANE%d',iSurf'),ZTA(it).label);
        contour3(xm,ym,ZTA(it).term{idx}(:,:),conts,['k' stl(iSurf)]);
    end
    title(sprintf('%s, T=%d d',basename,H(it).totim));
    view(2);
end
