function mkl2netcdf(filename, m, varargin)
% 
% mkl2netcdf
%
%	Creates mkl.nc file. 
%	'm' struct is called from convertMKL2nc.m
%	Uses matlab nc functions
%
%	---------------------------------------
% 
%       Giorgio Santinelli
%
%       Giorgio.Santinelli@deltares.nl
%
%       Deltares
%       P.O. Box 177
%       2600 MH Delft
%       The Netherlands
%
%	---------------------------------------
%
% $Id: mkl2netcdf.m 8819 2015-11-25 15:16:00Z santinel $
% $Date: 2015-11-25 16:16:00 +0100 (Wed, 25 Nov 2015) $
% $Author: santinel $
% $Revision: 8819 $
% $HeadURL: https://svn.oss.deltares.nl/repos/openearthrawdata/trunk/rijkswaterstaat/kustlijnkaartenboek/scripts/mkl2netcdf.m $


OPT = struct(...
    'publisher_name', getenv('USERNAME'),...
    'publisher_email', 'info@deltares.nl',...
    'historyatt', '$HeadURL $ $Id $',...
    'processing_level', 'preliminary');

OPT = setproperty(OPT, varargin);

datefmt = 'yyyy-mm-ddTHH:MMZ'; % date format
STRINGSIZE = 100;

    Attributes = [...
    struct( 'Name', 'naming_authority', 'Value', 'deltares.nl' ),... 
    struct( 'Name', 'id', 'Value', sprintf('MKL_release%s', datestr(nowutc, 'yyyymmdd'))),...
    struct( 'Name', 'Metadata_Conventions', 'Value', 'Unidata Dataset Discovery v1.0'),...
    struct( 'Name', 'title', 'Value', 'momentary coastline'),...
    struct( 'Name', 'summary', 'Value', 'Coastline indicator MKL along the Dutch coast since 1965, as described in the '' Kustlijnkaarten rapport'' http://www.rijkswaterstaat.nl/water/veiligheid/bescherming_tegen_het_water/veiligheidsmaatregelen/kustlijnzorg/kustlijnkaarten/'),...
    struct( 'Name', 'keywords', 'Value', 'Coastline, MKL, Dutch coast'),...
    struct( 'Name', 'keywords_vocabulary', 'Value', 'http://www.eionet.europa.eu/gemet'),...
    struct( 'Name', 'standard_name_vocabulary', 'Value', 'http://cf-convention.github.io/'),...
    struct( 'Name', 'history', 'Value', 'Data processed from xls/csv data - Winkust/Morphan output'),...
    struct( 'Name', 'HEADurl', 'Value', 'https://svn.oss.deltares.nl/repos/openearthrawdata/trunk/rijkswaterstaat/kustlijnkaartenboek/'),...
    struct( 'Name', 'Years covered', 'Value', ['1965-' num2str(m.jarkusyear(end))]),...
    struct( 'Name', 'institution', 'Value', 'Rijkswaterstaat'),...
    struct( 'Name', 'source'     , 'Value', 'on shore and off shore measurements'),...
    struct( 'Name', 'references' , 'Value', 'Original source: Rijkswaterstaat (Morphan), Deltares storage: http://opendap.deltares.nl/thredds/dodsC/opendap/rijkswaterstaat/ converted with script convertMKL2nc.m'),...
    struct( 'Name', 'Conventions', 'Value', 'CF-1.6'),...
    ... % Creator Search attributes
    struct( 'Name', 'creator_name', 'Value', 'Rijkswaterstaat'),...
    struct( 'Name', 'creator_url', 'Value', 'http://www.rijkswaterstaat.nl'),...
    struct( 'Name', 'creator_email', 'Value', 'info@rijkswaterstaat.nl'),...
    struct( 'Name', 'date_created', 'Value', datestr(nowutc, datefmt)),...
    struct( 'Name', 'date_modified', 'Value', datestr(nowutc, datefmt)),...
    struct( 'Name', 'date_issued', 'Value', datestr(nowutc, datefmt)),...
    ... % Publisher Search attributes
    struct( 'Name', 'publisher_name', 'Value', OPT.publisher_name),...
    struct( 'Name', 'publisher_url', 'Value', 'http://www.deltares.nl'),...
    struct( 'Name', 'publisher_email', 'Value', OPT.publisher_email),...
    ... % Extent Search attributes
    ... % Other Extent Information attributes
    struct( 'Name', 'geospatial_vertical_units', 'Value', 'm'),...
    struct( 'Name', 'geospatial_vertical_positive', 'Value', 'up'),...
    ...% Other attributes
	struct( 'Name', 'processing_level', 'Value', OPT.processing_level),...
	struct( 'Name', 'license', 'Value', [sprintf('These data can be used freely for research purposes provided that the following source is acknowledged: %s. ', 'RIJKSWATERSTAAT')...
                '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.']),...
    struct( 'Name', 'time_coverage_start', 'Value', '1965'),...
    struct( 'Name', 'time_coverage_end', 'Value', num2str(m.jarkusyear(end))),...
    struct( 'Name', 'time_coverage_resolution', 'Value', 'year')...
    ];

%% Define and create dimensions
% %  unlimited "time" dimension in Matlab seems to give variable "time" a 
% %  size equal to the the maximum dimension among all variables [2268]. 
% %  This leads to storing a useless amount of data. 
%    time_dim = struct('Name', 'time', 'Length', 0, 'Unlimited', true);
    time_dim = struct('Name', 'time', 'Length', length(m.time), 'Unlimited', false);
    bounds2_dim = struct('Name', 'bounds2', 'Length', 2, 'Unlimited', false);
    alongshore_dim = struct('Name', 'alongshore', 'Length', length(m.id), 'Unlimited', false);
    cross_shore_dim = struct('Name', 'cross_shore', 'Length', length(m.cross), 'Unlimited', false);
    stringsize_dim = struct('Name', 'stringsize', 'Length', STRINGSIZE, 'Unlimited', false);
    Dimensions = [time_dim, bounds2_dim, alongshore_dim, cross_shore_dim, stringsize_dim];
    
%% add epsg
    % smaller variables first
    epsg        = 28992;
    [~,~,OPT]=convertCoordinates(0,0,'CS1.code',epsg,'CS2.code',4326);
    %[CoordinateSystems, Operations , CoordSysCart ,CoordSysGeo] = GetCoordinateSystems();
    s.Name      = 'epsg';
    s.Datatype    = 'int32';
    s.Dimensions = [];
    s.Attributes = 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', ...
        {OPT.proj_conv1.method.name,    ...
        OPT.CS1.ellips.semi_major_axis, ...
        OPT.CS1.ellips.semi_minor_axis, ...
        OPT.CS1.ellips.inv_flattening,  ...
        OPT.proj_conv1.param.value(strcmp(OPT.proj_conv1.param.name,'Latitude of natural origin'    )), ...
        OPT.proj_conv1.param.value(strcmp(OPT.proj_conv1.param.name,'Longitude of natural origin'   )),...
        OPT.proj_conv1.param.value(strcmp(OPT.proj_conv1.param.name,'False easting'                 )),...
        OPT.proj_conv1.param.value(strcmp(OPT.proj_conv1.param.name,'False northing'                )),...
        OPT.proj_conv1.param.value(strcmp(OPT.proj_conv1.param.name,'Scale factor at natural origin')),...
        'value is equal to EPSG code'});
    Variables = s;

%% Define and create variables
    s.Name       = 'id';
    s.Datatype   = 'int32';
    s.Dimensions = alongshore_dim;
    s.Attributes = struct('Name' ,{'long_name' , 'comment'},...
                         'Value',{'identifier' , 'sum of area code (*1e6) and alongshore coordinate'});
    Variables = [Variables, s];
    
    flag_values = 1:17;%unique(m.areaCode);
    flag_meanings = rws_kustvak(flag_values);
    flag_meanings = str2line(flag_meanings,'s',', ');
    
    s.Name       = 'areacode';
    s.Datatype   = 'int32';
    s.Dimensions = alongshore_dim;
    s.Attributes = struct('Name' ,{'long_name','flag_values','flag_meanings','flag_comment'       ,'comment'},...
                          'Value',{'area code', flag_values , flag_meanings ,'points to: areaname','codes for the 17 coastal areas (kustvakken) as defined by Rijkswaterstaat'});
    Variables = [Variables, s];

    s.Name      = 'areaname';
    s.Datatype  = 'char';
    s.Dimensions = [stringsize_dim, alongshore_dim];
    s.Attributes = struct('Name',{'long_name','flag_comment'        ,'comment'},...
                         'Value',{'area name','indexed in: areacode','names for the 17 coastal areas (kustvakken) as defined by Rijkswaterstaat'});
    Variables = [Variables, s];

    s.Name       = 'alongshore';
    s.Datatype   = 'double';
    s.Dimensions = alongshore_dim;
    s.Attributes = struct('Name' ,{'long_name'             , 'units', 'comment'},...
                         'Value',{'alongshore coordinate', 'm'     , 'alongshore coordinate within the 17 coastal areas (kustvakken) as defined by Rijkswaterstaat'});
    Variables = [Variables, s];

    s.Name       = 'cross_shore';
    s.Datatype   = 'double';
    s.Dimensions = cross_shore_dim;
    s.Attributes = struct('Name' ,{'long_name'             , 'units', 'comment'},...
                         'Value',{'cross-shore coordinate', 'm'    , 'cross-shore coordinate relative to the rsp (rijks strand paal)'});
    Variables = [Variables, s];

    s.Name       = 'time';
    s.Datatype   = 'double';
    s.Dimensions = time_dim;
    s.Attributes = struct('Name' ,{'standard_name'          ,'axis' ,'units'                 ,'cell_methods','bounds'     ,'comment'         },...
                         'Value',{'time'                   ,'T'    ,'days since 1970-01-01' ,'mean'        ,'time_bounds','measurement year (date is artificially set to July, 1st); see time_MKL'});
    Variables = [Variables, s];

    s.Name      = 'time_bounds';
    s.Datatype    = 'double';
    s.Dimensions = [bounds2_dim, time_dim];
    s.Attributes = struct('Name' ,{'standard_name'          ,'units'                },...
                         'Value',{'time'                   ,'days since 1970-01-01'});
    Variables = [Variables, s];

    % ! important: time has to be the last! 
    % http://stackoverflow.com/questions/16938779/matlab-nc-unlimited-in-the-wrong-index-when-writing-unlimited-time-dimension-i
    
    % MKL_lat
    s.Name       = 'momentary_coastline_lat';
    s.Datatype   = 'double';
    s.Dimensions = [alongshore_dim, time_dim];
    s.Attributes = struct('Name' ,{'standard_name'         , 'long_name',             'units',        'axis', 'comment'},...
                         'Value',{'latitude'               , 'MKL latitude position', 'degree_north', 'X',    'MKL latitude position as defined by Rijkswaterstaat'});
    Variables = [Variables, s];
    
    % MKL_lon
    s.Name       = 'momentary_coastline_lon';
    s.Datatype   = 'double';
    s.Dimensions = [alongshore_dim, time_dim];
    s.Attributes = struct('Name' ,{'standard_name', 'long_name',              'units',        'axis', 'comment'},...
                          'Value',{'longitude'    , 'MKL longitude position', 'degrees_east', 'Y',    'MKL longitude position as defined by Rijkswaterstaat'});
    Variables = [Variables, s];
    
    % MKL_x
    s.Name       = 'momentary_coastline_x';
    s.Datatype   = 'double';
    s.Dimensions = [alongshore_dim, time_dim];
    s.Attributes = struct('Name' ,{'standard_name'          , 'long_name',      'units', 'axis', 'comment'},...
                          'Value',{'projection_x_coordinate', 'MKL x position', 'm',     'X',    'MKL x position as defined by Rijkswaterstaat'});
    Variables = [Variables, s];
    
    % MKL_y
    s.Name       = 'momentary_coastline_y';
    s.Datatype   = 'double';
    s.Dimensions = [alongshore_dim, time_dim];
    s.Attributes = struct('Name' ,{'standard_name'          , 'long_name',      'units', 'axis', 'comment'},...
                          'Value',{'projection_y_coordinate', 'MKL y position', 'm',     'Y',    'MKL y position as defined by Rijkswaterstaat'});
    Variables = [Variables, s];
    
    % MKL_rsp
    s.Name       = 'momentary_coastline';
    s.Datatype   = 'double';
    s.Dimensions = [alongshore_dim, time_dim];
    s.Attributes = struct('Name' ,{'long_name'   , 'units', 'comment'},...
                          'Value',{'MKL position', 'm'    , 'MKL position with respect to the beach pole (rijks strand paal) as defined by Rijkswaterstaat'});
    Variables = [Variables, s];
 
    % time_MKL
    s.Name       = 'time_MKL';
    s.Datatype   = 'double';
    s.Dimensions = [alongshore_dim, time_dim];
    s.Attributes = struct('Name' ,{'long_name', 'units'                    , 'comment'},...
                          'Value',{'MKL date',  'days since 1970-01-01'    , 'MKL date as defined by Rijkswaterstaat'});
    Variables = [Variables, s];

    % hoog
    s.Name         = 'high_boundary_MKL';
    s.Datatype     = 'double';
    s.Dimensions   = [alongshore_dim, time_dim];
    s.Attributes   = struct('Name' ,{'long_name',         'units', 'comment'},...
                          'Value',{'MKL high boundary', 'm'    , 'MKL high boundary with respect to NAP, as defined by Rijkswaterstaat'});
    Variables = [Variables, s];

    % laag
    s.Name       = 'low_boundary_MKL';
    s.Datatype   = 'double';
    s.Dimensions = [alongshore_dim, time_dim];
    s.Attributes = struct('Name' ,{'long_name',        'units', 'comment'},...
                          'Value',{'MKL low boundary', 'm'    , 'MKL low boundary with respect to NAP, as defined by Rijkswaterstaat'});
    Variables = [Variables, s];

    % XMaxSea
    s.Name       = 'seaward_boundary_MKL';
    s.Datatype   = 'double';
    s.Dimensions = [alongshore_dim, time_dim];
    s.Attributes = struct('Name' ,{'long_name',            'units', 'comment'},...
                          'Value',{'MKL seaward boundary', 'm'    , 'MKL seaward boundary, as defined by Rijkswaterstaat'});
    Variables = [Variables, s];
               
    % XMinLand
    s.Name       = 'landward_boundary_MKL';
    s.Datatype   = 'double';
    s.Dimensions = [alongshore_dim, time_dim];
    s.Attributes = struct('Name' ,{'long_name',             'units', 'comment'},...
                          'Value',{'MKL landward boundary', 'm'    , 'MKL landward boundary, as defined by Rijkswaterstaat'});
    Variables = [Variables, s];

    % Volume
    s.Name       = 'volume_MKL';
    s.Datatype   = 'double';
    s.Dimensions = [alongshore_dim, time_dim];
    s.Attributes = struct('Name' ,{'long_name',  'units', 'comment'},...
                          'Value',{'MKL volume', 'm'    , 'MKL volume in m^2/m, as defined by Rijkswaterstaat'});
    Variables = [Variables, s];


%% write schema to file        

ncSchema = struct(... 'Filename', strrep(filename, '.nc', 'test.nc'),...
        'Name', '/',...
        'Dimensions', Dimensions,...
        'Variables', Variables,...
        'Attributes', Attributes,...
        'Groups', [],...
        'Format', 'classic');

    ncwriteschema(filename, ncSchema)

%% Store index variables
    ncwrite(filename, 'time'        , m.time    );
    ncwrite(filename, 'time_bounds' , m.time_bounds);

    ncwrite(filename, 'id'          , m.id);
    ncwrite(filename, 'areacode'    , m.areacode);
    areanames = m.areaname;
    areanames(:, size(areanames,2)+1:STRINGSIZE) = ' ';
    ncwrite(filename, 'areaname', areanames);

    ncwrite(filename, 'alongshore', m.along);
    ncwrite(filename, 'cross_shore', m.cross);
    
    ncwrite(filename, 'momentary_coastline_lat'        , m.MKL_lat');
    ncwrite(filename, 'momentary_coastline_lon'        , m.MKL_lon');
    ncwrite(filename, 'momentary_coastline_x'          , m.MKL_x');
    ncwrite(filename, 'momentary_coastline_y'          , m.MKL_y');
    
    ncwrite(filename, 'momentary_coastline' , m.Xmkl');
    
    ncwrite(filename, 'time_MKL'             , m.timeMKL');
    ncwrite(filename, 'high_boundary_MKL'    , m.hoog');
    ncwrite(filename, 'low_boundary_MKL'     , m.laag');
    ncwrite(filename, 'seaward_boundary_MKL' , m.XzeeMKL');
    ncwrite(filename, 'landward_boundary_MKL', m.XMinLand');
    ncwrite(filename, 'volume_MKL'           , m.Volume');
    
%% check
ncinfo(filename);

%%%%%%


end % function mkl2netcdf.m
