% coastalIndicators2nc
%
% Generates coastal_indicators.nc (for MorphAn) which contains
%
% MKL is obtained using Winkust boundaries (Rijkswaterstaat).
%     The information are available from Deltares repository @
%     https://repos.deltares.nl/repos/OpenEarthRawData/trunk/rijkswaterstaat/kustlijnkaartenboek/raw/mkl/
%
% Dune Foot  defined as the crossing between the water level at +3.00m nap
%            and the coastline profile
%
% Indicators provided by ARCADIS and HKV: information is available in the .nc itself
%
% Data are taken from Jarkus dataset (1965 - 2010) @
% http://opendap.deltares.nl/thredds/dodsC/opendap/rijkswaterstaat/jarkus/profiles/transect.nc
%
% NB. strandlijnen dataset (1843 - 1998) is already in a .nc format @
% https://repos.deltares.nl/repos/OpenEarthRawData/trunk/rijkswaterstaat/strandlijnen/processed/strandlijnen.nc
% MHW, MLW are taken from two .mat files, obtained with MW_Jarkus.m
%
% STRANDLIJNEN dataset is NOT MERGED WITH JARKUS-STRANDLIJNEN YET! 
% - Decision on how combining them is required
%
% giorgio.santinelli@deltares.nl
% 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%% settings
% strandlijnen database - dataset 1843 - 1998
OPT.ncpath1      = 'https://repos.deltares.nl/repos/OpenEarthRawData/trunk/rijkswaterstaat/strandlijnen/processed/';
OPT.ncpath1      = 'F:\checkout\OpenEarthRawData\trunk\rijkswaterstaat\strandlijnen\processed\';
OPT.Filename1    = 'strandlijnen.nc';
% jarkus-strandlijnen database - dataset 1965 - 2010
OPT.ncpath1      = '../raw/';
OPT.Filename1_l  = 'MLW_rsp.mat';
OPT.Filename1_h  = 'MHW_rsp.mat';
% jarkus database - dataset 1965 - 2010
OPT.ncpath2      = 'http://opendap.deltares.nl/thredds/dodsC/opendap/rijkswaterstaat/jarkus/profiles/';
OPT.ncpath2      = 'F:\JARKUSdata\';
OPT.Filename2    = 'transect.nc';
% MKL from Winkust 1965 - 2010
OPT.ncpath3      = 'http://opendap.deltares.nl/thredds/dodsC/opendap/rijkswaterstaat/BKL_TKL_MKL/';
OPT.Filename3    = 'MKL.nc';
% DF
OPT.ncpath4      = 'http://opendap.deltares.nl/thredds/dodsC/opendap/rijkswaterstaat/DuneFoot/';
OPT.Filename4    = 'DF.nc';
% ARCADIS
OPT.ncpath5      = 'd:\santinel\Giorgio\Kustlijnzorg\Arcadis_kustlijnindicatoren_nc\kustlijn_indicatoren\'; % only local so far
OPT.Filename5    = 'kustlijnindicatoren_netCDF_nov2011_v2.nc';
% HKV faalkans
OPT.ncpath6      = 'http://opendap.deltares.nl/thredds/dodsC/opendap/rijkswaterstaat/faalkans_PC-Ring/';
OPT.Filename6    = 'faalkans.nc';
% HKV BeachWidth
OPT.ncpath7      = 'http://opendap.deltares.nl/thredds/dodsC/opendap/rijkswaterstaat/strandbreedte/';
OPT.Filename7    = 'strandbreedte.nc';

% general settings
OPT.outputfile   = '..\processed\coastal_indicators.nc';
OPT.disp         = true;
OPT.dump         = true;

%% jarkus database - dataset 1965 - 2010
%% Read data
id       = nc_varget([OPT.ncpath2 OPT.Filename2],'id');
MHW      = nc_varget([OPT.ncpath2 OPT.Filename2],'mean_high_water');
MLW      = nc_varget([OPT.ncpath2 OPT.Filename2],'mean_low_water');
lat      = nc_varget([OPT.ncpath2 OPT.Filename2],'lat');
lon      = nc_varget([OPT.ncpath2 OPT.Filename2],'lon');
cross    = nc_varget([OPT.ncpath2 OPT.Filename2],'cross_shore');
along    = nc_varget([OPT.ncpath2 OPT.Filename2],'alongshore');
areaname = nc_varget([OPT.ncpath2 OPT.Filename2],'areaname');
areacode = nc_varget([OPT.ncpath2 OPT.Filename2],'areacode');
time     = nc_varget([OPT.ncpath2 OPT.Filename2],'time');

rsp_x    = nc_varget([OPT.ncpath2 OPT.Filename2],'rsp_x');
rsp_y    = nc_varget([OPT.ncpath2 OPT.Filename2],'rsp_y');
x        = nc_varget([OPT.ncpath2 OPT.Filename2],'x');
y        = nc_varget([OPT.ncpath2 OPT.Filename2],'y');
angle    = nc_varget([OPT.ncpath2 OPT.Filename2],'angle');

load([OPT.ncpath1 OPT.Filename1_l]); % MLW_rsp
load([OPT.ncpath1 OPT.Filename1_h]); % MHW_rsp
momentary_coastline = (nc_varget([OPT.ncpath3,OPT.Filename3],'momentary_coastline'))';
dune_foot           =  nc_varget([OPT.ncpath4,OPT.Filename4],'dune_foot_rsp');
MDL                 = (nc_varget([OPT.ncpath5,OPT.Filename5],'MDL'))';
Xafslag             = (nc_varget([OPT.ncpath5,OPT.Filename5],'Xafslag'))';
MDV                 = (nc_varget([OPT.ncpath5,OPT.Filename5],'MDV'))';
P                   = (nc_varget([OPT.ncpath5,OPT.Filename5],'P'))';
duneratio           = (nc_varget([OPT.ncpath5,OPT.Filename5],'duneratio'))';
nearshore_slope     = (nc_varget([OPT.ncpath5,OPT.Filename5],'nearshore_slope'))';
beach_slope         = (nc_varget([OPT.ncpath5,OPT.Filename5],'beach_slope'))';
probability_failure = (nc_varget([OPT.ncpath6,OPT.Filename6],'probability_failure'))';
beach_width_at_MHW  = (nc_varget([OPT.ncpath7,OPT.Filename7],'beach_width_at_MHW'))';
beach_width_at_MLW  = (nc_varget([OPT.ncpath7,OPT.Filename7],'beach_width_at_MLW'))';

%% create empty file
nc_create_empty(OPT.outputfile)

%% Add overall meta info
%  http://cf-pcmdi.llnl.gov/documents/cf-conventions/1.5/cf-conventions.html#description-of-file-contents
%------------------

attval =  sprintf('created on %s', datestr(now));

nc_attput(OPT.outputfile, nc_global, 'title'         , 'coastal indicators based on Jarkus data');
nc_attput(OPT.outputfile, nc_global, 'institution'   , 'Rijkswaterstaat, Deltares, ARCADIS, HKV Lijn in Water');
nc_attput(OPT.outputfile, nc_global, 'source'        , 'measurements');
nc_attput(OPT.outputfile, nc_global, 'history'       , ['Data processed from Jarkus database and .nc files of indicators, ', attval]);
nc_attput(OPT.outputfile, nc_global, 'references'    , 'original source: http://www.rijkswaterstaat.nl , Deltares storage: http://opendap.deltares.nl/thredds/dodsC/opendap/rijkswaterstaat/ , converted with script coastalIndicators2nc.m');
nc_attput(OPT.outputfile, nc_global, 'email'         , 'giorgio.santinelli@deltares.nl');
% nc_attput(OPT.outputfile, nc_global, 'comment'       , '');
nc_attput(OPT.outputfile, nc_global, 'version'       , '0.1');
nc_attput(OPT.outputfile, nc_global, 'Conventions'   , 'CF-1.5');
% spatial reference using wkt. approach used by gdal.
nc_attput(OPT.outputfile, nc_global, 'spatial_ref', 'COMPD_CS["Amersfoort / RD New + NAP",PROJCS["Amersfoort / RD New",GEOGCS["Amersfoort",DATUM["Amersfoort",SPHEROID["Bessel 1841",6377397.155,299.1528128,AUTHORITY["EPSG","7004"]],TOWGS84[565.04,49.91,465.84,-0.40939438743923684,-0.35970519561431136,1.868491000350572,0.8409828680306614],AUTHORITY["EPSG","6289"]],PRIMEM["Greenwich",0.0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.017453292519943295],AXIS["Geodetic latitude",NORTH],AXIS["Geodetic longitude",EAST],AUTHORITY["EPSG","4289"]],PROJECTION["Oblique Stereographic",AUTHORITY["EPSG","9809"]],PARAMETER["central_meridian",5.387638888888891],PARAMETER["latitude_of_origin",52.15616055555556],PARAMETER["scale_factor",0.9999079],PARAMETER["false_easting",155000.0],PARAMETER["false_northing",463000.0],UNIT["m",1.0],AXIS["Easting",EAST],AXIS["Northing",NORTH],AUTHORITY["EPSG","28992"]],VERT_CS["Normaal Amsterdams Peil",VERT_DATUM["Normaal Amsterdams Peil",2005,AUTHORITY["EPSG","5109"]],UNIT["m",1.0],AXIS["Gravity-related height",UP],AUTHORITY["EPSG","5709"]],AUTHORITY["EPSG","7415"]]');
% nc_attput(OPT.outputfile, nc_global, 'CF:featureType', 'stationTimeSeries');
nc_attput(OPT.outputfile, nc_global, 'terms_for_use' , 'These data can be used freely for research purposes provided that the following source is acknowledged: Rijkswaterstaat, Deltares, ARCADIS, HKV Lijn in Water');
% nc_attput(OPT.outputfile, nc_global, 'disclaimer'    , 'This data is made available in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.');

%% create dimensions
% add dimensions
nc_adddim(OPT.outputfile, 'time', length(time)); % 0);
nc_adddim(OPT.outputfile, 'alongshore', length(along)); % one per transect
nc_adddim(OPT.outputfile, 'stringsize', 100);

%% add variables
clear nc
i = 0;

%% add epsg

% smaller variables first
[lon2,lat2,OPT2] = convertCoordinates(0,0,'CS1.code',28992,'CS2.code',4326);
%[CoordinateSystems, Operations , CoordSysCart ,CoordSysGeo] = GetCoordinateSystems();
epsg        = 28992;

i = i + 1;
nc(i).Name      = 'epsg';
nc(i).Nctype    = 'int';
nc(i).Dimension = {};
nc(i).Attribute = struct('Name', ...
   {'grid_mapping_name', ...
    'semi_major_axis', ...
    'semi_minor_axis', ...
    'inverse_flattening', ...
    'latitude_of_projection_origin', ...
    'longitude_of_projection_origin', ...
    'false_easting', ...
    'false_northing', ...
    'scale_factor_at_projection_origin',...
    'comment'}, ...
    'Value', ...
    {OPT2.proj_conv1.method.name,    ...
    OPT2.CS1.ellips.semi_major_axis, ...
    OPT2.CS1.ellips.semi_minor_axis, ...
    OPT2.CS1.ellips.inv_flattening,  ...
    OPT2.proj_conv1.param.value(strcmp(OPT2.proj_conv1.param.name,'Latitude of natural origin'    )), ...
    OPT2.proj_conv1.param.value(strcmp(OPT2.proj_conv1.param.name,'Longitude of natural origin'   )),...
    OPT2.proj_conv1.param.value(strcmp(OPT2.proj_conv1.param.name,'False easting'                 )),...
    OPT2.proj_conv1.param.value(strcmp(OPT2.proj_conv1.param.name,'False northing'                )),...
    OPT2.proj_conv1.param.value(strcmp(OPT2.proj_conv1.param.name,'Scale factor at natural origin')),...
    'value is equal to EPSG code'});

% id
i = i + 1;
nc(i).Name         = 'id';
nc(i).Nctype       = 'int';
nc(i).Dimension    = {'alongshore'};
nc(i).Attribute(1) = struct('Name', 'long_name' ,'Value', 'identifier');
nc(i).Attribute(2) = struct('Name', 'comment'   ,'Value', 'sum of area code (x1000000) and alongshore coordinate');

% areacode
i = i + 1;
nc(i).Name         = 'areacode';
nc(i).Nctype       = 'int';
nc(i).Dimension    = {'alongshore'};
nc(i).Attribute(1) = struct('Name', 'long_name' ,'Value', 'area code');
nc(i).Attribute(2) = struct('Name', 'comment'   ,'Value', 'codes for the 15 coastal areas as defined by rijkswaterstaat');

% areaname
i = i + 1;
nc(i).Name         = 'areaname';
nc(i).Nctype       = 'char';
nc(i).Dimension    = {'alongshore','stringsize'};
nc(i).Attribute(1) = struct('Name', 'long_name' ,'Value', 'area name');
nc(i).Attribute(2) = struct('Name', 'comment'   ,'Value', 'names for the 15 coastal areas as defined by rijkswaterstaat');

% alongshore
i = i + 1;
nc(i).Name         = 'alongshore';
nc(i).Nctype       = 'double'; 
nc(i).Dimension    = {'alongshore'};
nc(i).Attribute(1) = struct('Name', 'long_name' ,'Value', 'alongshore coordinate');
nc(i).Attribute(2) = struct('Name', 'units'     ,'Value', 'm');
nc(i).Attribute(3) = struct('Name', 'comment'   ,'Value', 'alongshore coordinate within the 15 coastal areas as defined by rijkswaterstaat');

% angle
i = i + 1;
nc(i).Name         = 'angle';
nc(i).Nctype       = 'double';
nc(i).Dimension    = {'alongshore'};
nc(i).Attribute(1) = struct('Name', 'long_name' ,'Value', 'angle of transect');
nc(i).Attribute(2) = struct('Name', 'units'     ,'Value', 'degrees');
nc(i).Attribute(3) = struct('Name', 'comment'   ,'Value', 'positive clockwise 0 north');

% time
i = i + 1;
nc(i).Name         = 'time';
nc(i).Nctype       = 'double';
nc(i).Dimension    = {'time'};
nc(i).Attribute(1) = struct('Name', 'standard_name' ,'Value', 'time');
nc(i).Attribute(2) = struct('Name', 'axis'          ,'Value', 'T');
nc(i).Attribute(3) = struct('Name', 'units'         ,'Value', 'days since 1970-01-01 00:00 +1:00');
nc(i).Attribute(4) = struct('Name', 'comment'       ,'Value', 'measurement year');

% indicator_value (dune_foot)
i = i + 1;
nc(i).Name         = 'dune_foot';
nc(i).Nctype       = 'double';
nc(i).Dimension    = {'alongshore','time'};
nc(i).Attribute(1) = struct('Name', 'indicator_name','Value', 'dune_foot');
nc(i).Attribute(2) = struct('Name', 'long_name'     ,'Value', 'dune foot position wrt rijkstrandpaal');
nc(i).Attribute(3) = struct('Name', 'units'         ,'Value', 'm');
nc(i).Attribute(4) = struct('Name', 'comment'       ,'Value', 'dune foot is defined by the intersection of the profile with the line +3.00m nap');
nc(i).Attribute(5) = struct('Name', '_FillValue'    ,'Value', NaN);

% indicator_value (MHW_rsp)
i = i + 1;
nc(i).Name         = 'mean_high_water';
nc(i).Nctype       = 'double';
nc(i).Dimension    = {'alongshore','time'};
nc(i).Attribute(1) = struct('Name', 'indicator_name','Value', 'mean_high_water');
nc(i).Attribute(2) = struct('Name', 'long_name'     ,'Value', 'mean high water position wrt rijkstrandpaal');
nc(i).Attribute(3) = struct('Name', 'units'         ,'Value', 'm');
nc(i).Attribute(4) = struct('Name', 'comment'       ,'Value', 'mean high water position is defined by the intersection of the profile with mean high water line');
nc(i).Attribute(5) = struct('Name', '_FillValue'    ,'Value', NaN);

% indicator_value (MLW_rsp)
i = i + 1;
nc(i).Name         = 'mean_low_water';
nc(i).Nctype       = 'double';
nc(i).Dimension    = {'alongshore','time'};
nc(i).Attribute(1) = struct('Name', 'indicator_name','Value', 'mean_low_water');
nc(i).Attribute(2) = struct('Name', 'long_name'     ,'Value', 'mean low water position wrt rijkstrandpaal');
nc(i).Attribute(3) = struct('Name', 'units'         ,'Value', 'm');
nc(i).Attribute(4) = struct('Name', 'comment'       ,'Value', 'mean low water position is defined by the intersection of the profile with mean low water line');
nc(i).Attribute(5) = struct('Name', '_FillValue'    ,'Value', NaN);

% indicator_value (momentary_coastline)
i = i + 1;
nc(i).Name         = 'momentary_coastline';
nc(i).Nctype       = 'double';
nc(i).Dimension    = {'alongshore','time'};
nc(i).Attribute(1) = struct('Name', 'indicator_name','Value', 'momentary_coastline');
nc(i).Attribute(2) = struct('Name', 'long_name'     ,'Value', 'MKL position');
nc(i).Attribute(3) = struct('Name', 'units'         ,'Value', 'm');
nc(i).Attribute(4) = struct('Name', 'comment'       ,'Value', 'MKL position defined by Rijkswaterstaat');
nc(i).Attribute(5) = struct('Name', '_FillValue'    ,'Value', NaN);

% indicator_value (MDL)
i = i + 1;
nc(i).Name         = 'MDL';
nc(i).Nctype       = 'double';
nc(i).Dimension    = {'alongshore','time'};
nc(i).Attribute(1) = struct('Name', 'indicator_name','Value', 'MDL');
nc(i).Attribute(2) = struct('Name', 'long_name'     ,'Value', 'MDL position');
nc(i).Attribute(3) = struct('Name', 'units'         ,'Value', 'm');
nc(i).Attribute(4) = struct('Name', 'comment'       ,'Value', '`Momentane DuinLijn` (MDL) position wrt RSP');
nc(i).Attribute(5) = struct('Name', '_FillValue'    ,'Value', NaN);

% indicator_value (Xafslag)
i = i + 1;
nc(i).Name         = 'Xafslag';
nc(i).Nctype       = 'double';
nc(i).Dimension    = {'alongshore','time'};
nc(i).Attribute(1) = struct('Name', 'indicator_name','Value', 'Xafslag');
nc(i).Attribute(2) = struct('Name', 'long_name'     ,'Value', 'Xafslag');
nc(i).Attribute(3) = struct('Name', 'units'         ,'Value', 'm');
nc(i).Attribute(4) = struct('Name', 'comment'       ,'Value', '`Positie Maatgevend Afslagpunt` (Xafslag) wrt RSP');
nc(i).Attribute(5) = struct('Name', '_FillValue'    ,'Value', NaN);

% indicator_value (MDV)
i = i + 1;
nc(i).Name         = 'MDV';
nc(i).Nctype       = 'double';
nc(i).Dimension    = {'alongshore','time'};
nc(i).Attribute(1) = struct('Name', 'indicator_name','Value', 'MDV');
nc(i).Attribute(2) = struct('Name', 'long_name'     ,'Value', 'MDV');
nc(i).Attribute(3) = struct('Name', 'units'         ,'Value', 'm^3/m');
nc(i).Attribute(4) = struct('Name', 'comment'       ,'Value', '`Momentaan Duin Volume` (MDV)');
nc(i).Attribute(5) = struct('Name', '_FillValue'    ,'Value', NaN);

% indicator_value (P)
i = i + 1;
nc(i).Name         = 'P';
nc(i).Nctype       = 'double';
nc(i).Dimension    = {'alongshore','time'};
nc(i).Attribute(1) = struct('Name', 'indicator_name','Value', 'P');
nc(i).Attribute(2) = struct('Name', 'long_name'     ,'Value', 'Probability of exceeding Xaf1990');
nc(i).Attribute(3) = struct('Name', 'units'         ,'Value', '1/year');
nc(i).Attribute(4) = struct('Name', 'comment'       ,'Value', 'Probability of exceeding (P) the `Maatgevend Afslagpunt van 1990` (Xaf1990)');
nc(i).Attribute(5) = struct('Name', '_FillValue'    ,'Value', NaN);

% indicator_value (duneratio)
i = i + 1;
nc(i).Name         = 'duneratio';
nc(i).Nctype       = 'double';
nc(i).Dimension    = {'alongshore','time'};
nc(i).Attribute(1) = struct('Name', 'indicator_name','Value', 'duneratio');
nc(i).Attribute(2) = struct('Name', 'long_name'     ,'Value', 'Characteristic dune ratio');
nc(i).Attribute(3) = struct('Name', 'units'         ,'Value', '-');
nc(i).Attribute(4) = struct('Name', 'comment'       ,'Value', ' = mean( DH/(xDV-Xaf1990) )');
nc(i).Attribute(5) = struct('Name', '_FillValue'    ,'Value', NaN);

% indicator_value (nearshore_slope)
i = i + 1;
nc(i).Name         = 'nearshore_slope';
nc(i).Nctype       = 'double';
nc(i).Dimension    = {'alongshore','time'};
nc(i).Attribute(1) = struct('Name', 'indicator_name','Value', 'nearshore_slope');
nc(i).Attribute(2) = struct('Name', 'long_name'     ,'Value', 'Characteristic nearshore_slope');
nc(i).Attribute(3) = struct('Name', 'units'         ,'Value', 'm');
nc(i).Attribute(4) = struct('Name', 'comment'       ,'Value', ' = mean( (UL-LL)/(2*(xSea-MKL)) )');
nc(i).Attribute(5) = struct('Name', '_FillValue'    ,'Value', NaN);

% indicator_value (beach_slope)
i = i + 1;
nc(i).Name         = 'beach_slope';
nc(i).Nctype       = 'double';
nc(i).Dimension    = {'alongshore','time'};
nc(i).Attribute(1) = struct('Name', 'indicator_name','Value', 'beach_slope');
nc(i).Attribute(2) = struct('Name', 'long_name'     ,'Value', 'Characteristic beach_slope');
nc(i).Attribute(3) = struct('Name', 'units'         ,'Value', 'm');
nc(i).Attribute(4) = struct('Name', 'comment'       ,'Value', ' = mean( (UL-LL)/(2*(MKL-xLand)) )');
nc(i).Attribute(5) = struct('Name', '_FillValue'    ,'Value', NaN);

% indicator_value (probability_failure)
i = i + 1;
nc(i).Name         = 'probability_failure';
nc(i).Nctype       = 'double';
nc(i).Dimension    = {'alongshore','time'};
nc(i).Attribute(1) = struct('Name', 'indicator_name','Value', 'probability_failure');
nc(i).Attribute(2) = struct('Name', 'long_name'     ,'Value', 'Probability of failure of the first dunerow');
nc(i).Attribute(3) = struct('Name', 'units'         ,'Value', '1/year');
nc(i).Attribute(4) = struct('Name', 'comment'       ,'Value', 'Computed by PC-Ring, based on the DUROS+ module');
nc(i).Attribute(5) = struct('Name', '_FillValue'    ,'Value', NaN);

% indicator_value (beach_width_at_MHW)
i = i + 1;
nc(i).Name         = 'beach_width_at_MHW';
nc(i).Nctype       = 'double';
nc(i).Dimension    = {'alongshore','time'};
nc(i).Attribute(1) = struct('Name', 'indicator_name','Value', 'beach_width_at_MHW');
nc(i).Attribute(2) = struct('Name', 'long_name'     ,'Value', 'Beach width at MHW');
nc(i).Attribute(3) = struct('Name', 'units'         ,'Value', 'm');
nc(i).Attribute(4) = struct('Name', 'comment'       ,'Value', 'Distance between dune boundary and sea boundary at MHW');
nc(i).Attribute(5) = struct('Name', '_FillValue'    ,'Value', NaN);

% indicator_value (beach_width_at_MLW)
i = i + 1;
nc(i).Name         = 'beach_width_at_MLW';
nc(i).Nctype       = 'double';
nc(i).Dimension    = {'alongshore','time'};
nc(i).Attribute(1) = struct('Name', 'indicator_name','Value', 'beach_width_at_MLW');
nc(i).Attribute(2) = struct('Name', 'long_name'     ,'Value', 'Beach width at MLW');
nc(i).Attribute(3) = struct('Name', 'units'         ,'Value', 'm');
nc(i).Attribute(4) = struct('Name', 'comment'       ,'Value', 'Distance between dune boundary and sea boundary at MLW');
nc(i).Attribute(5) = struct('Name', '_FillValue'    ,'Value', NaN);

%% create variables with attributes
for i = 1:length(nc)
    if OPT.disp;disp(['adding ',num2str(i),' ',nc(i).Name]);end
    nc_addvar(OPT.outputfile, nc(i));
end

%% store variables
nc_varput(OPT.outputfile, 'id'              , id);
nc_varput(OPT.outputfile, 'areacode'        , areacode);
nc_varput(OPT.outputfile, 'areaname'        , areaname);
nc_varput(OPT.outputfile, 'alongshore'      , along);
nc_varput(OPT.outputfile, 'angle'           , angle);
nc_varput(OPT.outputfile, 'time'            , time);

nc_varput(OPT.outputfile, 'dune_foot'           , dune_foot);
nc_varput(OPT.outputfile, 'mean_high_water'     , MHW_rsp);
nc_varput(OPT.outputfile, 'mean_low_water'      , MLW_rsp);
nc_varput(OPT.outputfile, 'momentary_coastline' , momentary_coastline);
nc_varput(OPT.outputfile, 'MDL'                 , MDL);
nc_varput(OPT.outputfile, 'Xafslag'             , Xafslag);
nc_varput(OPT.outputfile, 'MDV'                 , MDV);
nc_varput(OPT.outputfile, 'P'                   , P);
nc_varput(OPT.outputfile, 'duneratio'           , duneratio);
nc_varput(OPT.outputfile, 'nearshore_slope'     , nearshore_slope);
nc_varput(OPT.outputfile, 'beach_slope'         , beach_slope);
nc_varput(OPT.outputfile, 'probability_failure' , probability_failure);
nc_varput(OPT.outputfile, 'beach_width_at_MHW'  , beach_width_at_MHW);
nc_varput(OPT.outputfile, 'beach_width_at_MLW'  , beach_width_at_MLW);

%% check
if OPT.dump
    nc_dump(OPT.outputfile);
end

%% end



