%%
ncfile = fullfile(fileparts(which(mfilename)), '..', 'data', 'jarkus_volumes.nc');

id = desb_illustration_id;
st = find(id == nc_varget(ncfile, 'id')) - 1 ;
volumes = nc_varget(ncfile, 'volume_sediment', [0 st 0], [-1 1 -1]);
csb = nc_varget(ncfile, 'cross_shore_boundary', [st 0], [1 -1]);

tvolumes = sum(volumes,2);

durosnc = fullfile(fileparts(which(mfilename)), '..', 'simulations', sprintf('%08i', id), 'durosplus', 'run00', sprintf('%08i.nc', id));
erosion_volume = nc_varget(durosnc, 'erosion_volume');

ax1 = subplot(121);
plot(tvolumes/diff(csb([1 end])), erosion_volume, 'o')
xlabel({'Total'; 'transect volume [m^3/m/m]'})
ylabel('Erosion volume [m^3/m]')
v1 = axis;
text(v1(2), v1(4), sprintf('\\rho = %.2f', corr(tvolumes, erosion_volume)),...
    'HorizontalAlignment', 'right',...
    'VerticalAlignment', 'top')
text(v1(1), v1(4), '(a)',...
    'HorizontalAlignment', 'left',...
    'VerticalAlignment', 'top')

ax2 = subplot(122);
plot(volumes(:,3)/diff(csb(3:4)), erosion_volume, 'o')
xlabel({'Sedimentation zone (2)'; ' transect volume [m^3/m/m]'})
v2 = axis;
text(v2(2), v2(4), sprintf('\\rho = %.2f', corr(volumes(:,3), erosion_volume)),...
    'HorizontalAlignment', 'right',...
    'VerticalAlignment', 'top')
text(v2(1), v2(4), '(b)',...
    'HorizontalAlignment', 'left',...
    'VerticalAlignment', 'top')

axes(...
    'Position', [0 0 1 1],...
    'Visible', 'off');
text(0.5, 1,...
    sprintf('transect %i', id),...
    'HorizontalAlignment', 'center',...
    'VerticalAlignment', 'top')

%%
set(ax1, 'XLimMode', 'manual',...
    'YLimMode', 'manual')
set(ax2, 'YLimMode', 'manual')
figureStyle(@styleThesis)
tightenfigure
[fpath fname] = fileparts(which(mfilename));
figurepath = fullfile(fpath, '..', 'figures');
set(gcf, 'tag', '$HeadURL$ $Revision$')
print(gcf, '-depsc2', fullfile(figurepath, [fname '.eps']))
if ispc
    saveas(gcf, fullfile(figurepath, [fname '.emf']))
end
saveas(gcf, fullfile(figurepath, [fname '.fig']))