
clear all
close all
disp(' ')
disp(date)

pathOUT='...';
fileOUT='kustlijnindicatoren_netCDF_nov2011_v2.nc';

pathIN='...';
fileIN='OUTPUT_alltransects.mat';


%% LIJKT NOODZAKELIJK OM SCHRIJVEN NAAR netCDF file MOGELIJK TE MAKEN; WAAROM?
setpref('SNCTOOLS','USE_MEXNC',true)

%% Loading inputdata

load([pathIN fileIN])
raailist=OUTPUT.raaien-(round(OUTPUT.raaien/1000000)*10000);

[BKL_xRD BKL_yRD BKL_lat BKL_lon]                   = rsp2rd2latlon(OUTPUT.BKL,OUTPUT.raaien,OUTPUT.JRKINFO);
[MKL_xRD MKL_yRD MKL_lat MKL_lon]                   = rsp2rd2latlon(OUTPUT.MKL,OUTPUT.raaien,OUTPUT.JRKINFO);
[xDV_xRD xDV_yRD xDV_lat xDV_lon]                   = rsp2rd2latlon(OUTPUT.xDV,OUTPUT.raaien,OUTPUT.JRKINFO);
[MDL_xRD MDL_yRD MDL_lat MDL_lon]                   = rsp2rd2latlon(OUTPUT.MDL,OUTPUT.raaien,OUTPUT.JRKINFO);
[Xafslag_xRD Xafslag_yRD Xafslag_lat Xafslag_lon]   = rsp2rd2latlon(OUTPUT.Xafslag,OUTPUT.raaien,OUTPUT.JRKINFO);


%% Additional parameters (from input)

    zDV         =3;
    Xaf1990     =repmat(OUTPUT.Xaf1990in(:,2),1,OUTPUT.njaar);
    DH          =(OUTPUT.MDV./(OUTPUT.xDV-Xaf1990))+zDV;
% duinratio       	  =nanmean(DH./(OUTPUT.xDV-Xaf1990),2);
% vooroeverhelling    =nanmean((OUTPUT.MKLextra.UL-OUTPUT.MKLextra.LL)./(2*(OUTPUT.MKLextra.xSea-OUTPUT.MKL)),2);
% strandhelling       =nanmean((OUTPUT.MKLextra.UL-OUTPUT.MKLextra.LL)./(2*(OUTPUT.MKL-OUTPUT.MKLextra.xLand)),2);
duinratio       	=DH./(OUTPUT.xDV-Xaf1990);
vooroeverhelling    =(OUTPUT.MKLextra.UL-OUTPUT.MKLextra.LL)./(2*(OUTPUT.MKLextra.xSea-OUTPUT.MKL));
strandhelling       =(OUTPUT.MKLextra.UL-OUTPUT.MKLextra.LL)./(2*(OUTPUT.MKL-OUTPUT.MKLextra.xLand));

%% 1. Create netCDF file

%     outputfile = [pathOUT 'netCDF_OUTPUT-test.nc'];
    outputfile = [pathOUT fileOUT];
    nc_create_empty(outputfile);


%% 2. Add overall meta info

   %%  CF convention   
   % <http://cf-pcmdi.llnl.gov/documents/cf-conventions/1.4/cf-conventions.html#description-of-file-contents>

   nc_attput(outputfile, nc_global, 'title'         , 'Coastline indicators');
   nc_attput(outputfile, nc_global, 'institution'   , 'ARCADIS - Coastal and Marine Systems');
   nc_attput(outputfile, nc_global, 'source'        , ['$HeadURL: ' mfilename '.m$']);
   nc_attput(outputfile, nc_global, 'history'       , '');
   nc_attput(outputfile, nc_global, 'references'    , 'ARCADIS reports: "Verkenning: relatie kustlijnzorg - kustveiligheid" (C04021.002681; Januari 2011) and "Nadere uitwerkingen: relatie kustlijnzorg - kustveiligheid" (C03041.002738; October 2011)');
   nc_attput(outputfile, nc_global, 'email'         , '');
   nc_attput(outputfile, nc_global, 'comment'       , 'CONCEPT RELEASE');
   nc_attput(outputfile, nc_global, 'version'       , date);
   nc_attput(outputfile, nc_global, 'Conventions'   , 'CF-1.5');
   nc_attput(outputfile, nc_global, 'CF:featureType', 'point');  % https://cf-pcmdi.llnl.gov/trac/wiki/PointObservationConventions
   
   %%  OpenEarth convention

   nc_attput(outputfile, nc_global, 'terms_for_use' , '');
   nc_attput(outputfile, nc_global, 'disclaimer'    , '');


%% 2. Create dimensions

   nc_add_dimension(outputfile, 'Time'   , length(OUTPUT.jaren))
   nc_add_dimension(outputfile, 'ID'	 , length(OUTPUT.raaien))


%% 3. Create variables

   clear nc
   ifld = 0;
   
            ifld = ifld + 1;
   nc(ifld).Name         = 'Year';
   nc(ifld).Nctype       = 'int32';
   nc(ifld).Dimension    = {'Time'};
   nc(ifld).Attribute(1) = struct('Name', 'long_name'      ,'Value', 'Time');
   nc(ifld).Attribute(2) = struct('Name', 'units'          ,'Value', 'years');
   nc(ifld).Attribute(3) = struct('Name', 'axis'           ,'Value', 'T');
   nc(ifld).Attribute(4) = struct('Name', 'comment'        ,'Value', 'Year of JarKus measurement');

            ifld = ifld + 1;
   nc(ifld).Name         = 'ID';
   nc(ifld).Nctype       = 'int32';
   nc(ifld).Dimension    = {'ID'};
   nc(ifld).Attribute(1) = struct('Name', 'long_name'      ,'Value', 'Transect_id');
   nc(ifld).Attribute(2) = struct('Name', 'units'          ,'Value', '-');
   nc(ifld).Attribute(3) = struct('Name', 'axis'           ,'Value', 'X');
   nc(ifld).Attribute(4) = struct('Name', 'comment'        ,'Value', '= (areacode * 1000000) + transectnumber');

   
            ifld = ifld + 1;
   nc(ifld).Name         = 'BKL';
   nc(ifld).Nctype       = 'double';
   nc(ifld).Dimension    = {'ID'};
   nc(ifld).Attribute(1) = struct('Name', 'long_name'      ,'Value', 'BKL position');
   nc(ifld).Attribute(2) = struct('Name', 'units'          ,'Value', 'm');
   nc(ifld).Attribute(3) = struct('Name', 'axis'           ,'Value', 'Y');
   nc(ifld).Attribute(4) = struct('Name', 'comment'        ,'Value', '`Basis KustLijn` (BKL) position w.r.t. RSP');
   nc(ifld).Attribute(5) = struct('Name', '_FillValue'     ,'Value', NaN);
   
   nc=nc_add_newcoordsys(nc,'BKL','w.r.t. RSP');
   ifld=length(nc);
   

   ifld = ifld + 1;
   nc(ifld).Name         = 'MKL';
   nc(ifld).Nctype       = 'double';
   nc(ifld).Dimension    = {'Time','ID'};
   nc(ifld).Attribute(1) = struct('Name', 'long_name'      ,'Value', 'MKL position');
   nc(ifld).Attribute(2) = struct('Name', 'units'          ,'Value', 'm');
   nc(ifld).Attribute(3) = struct('Name', 'axis'           ,'Value', 'Y');
   nc(ifld).Attribute(4) = struct('Name', 'comment'        ,'Value', '`Momentane KustLijn` (MKL) position w.r.t. RSP');
   nc(ifld).Attribute(5) = struct('Name', '_FillValue'     ,'Value', NaN);
   
   nc=nc_add_newcoordsys(nc,'MKL','w.r.t. RSP');
   ifld=length(nc);
   
      
            ifld = ifld + 1;
   nc(ifld).Name         = 'xDV';
   nc(ifld).Nctype       = 'double';
   nc(ifld).Dimension    = {'Time','ID'};
   nc(ifld).Attribute(1) = struct('Name', 'long_name'      ,'Value', 'Duinvoet position');
   nc(ifld).Attribute(2) = struct('Name', 'units'          ,'Value', 'm');
   nc(ifld).Attribute(3) = struct('Name', 'axis'           ,'Value', 'Y');
   nc(ifld).Attribute(4) = struct('Name', 'comment'        ,'Value', '`Duinvoet` (xDV) position w.r.t. RSP (choice: dunefootlevel at NAP +3m)');
   nc(ifld).Attribute(5) = struct('Name', '_FillValue'     ,'Value', NaN);
   
   nc=nc_add_newcoordsys(nc,'xDV','w.r.t. RSP');
   ifld=length(nc);
   

            ifld = ifld + 1;
   nc(ifld).Name         = 'MDL';
   nc(ifld).Nctype       = 'double';
   nc(ifld).Dimension    = {'Time','ID'};
   nc(ifld).Attribute(1) = struct('Name', 'long_name'      ,'Value', 'MDL position');
   nc(ifld).Attribute(2) = struct('Name', 'units'          ,'Value', 'm');
   nc(ifld).Attribute(3) = struct('Name', 'axis'           ,'Value', 'Y');
   nc(ifld).Attribute(4) = struct('Name', 'comment'        ,'Value', '`Momentane DuinLijn` (MDL) position w.r.t. RSP');
   nc(ifld).Attribute(5) = struct('Name', '_FillValue'     ,'Value', NaN);
   
   nc=nc_add_newcoordsys(nc,'MDL','w.r.t. RSP');
   ifld=length(nc);
   
   
            ifld = ifld + 1;
   nc(ifld).Name         = 'Xafslag';
   nc(ifld).Nctype       = 'double';
   nc(ifld).Dimension    = {'Time','ID'};
   nc(ifld).Attribute(1) = struct('Name', 'long_name'      ,'Value', 'Xafslag');
   nc(ifld).Attribute(2) = struct('Name', 'units'          ,'Value', 'm');
   nc(ifld).Attribute(3) = struct('Name', 'axis'           ,'Value', 'Y');
   nc(ifld).Attribute(4) = struct('Name', 'comment'        ,'Value', '`Positie Maatgevend Afslagpunt` (Xafslag) w.r.t. RSP');
   nc(ifld).Attribute(5) = struct('Name', '_FillValue'     ,'Value', NaN);
   
   nc=nc_add_newcoordsys(nc,'Xafslag','w.r.t. RSP');
   ifld=length(nc);
   
   
               ifld = ifld + 1;
   nc(ifld).Name         = 'MDV';
   nc(ifld).Nctype       = 'double';
   nc(ifld).Dimension    = {'Time','ID'};
   nc(ifld).Attribute(1) = struct('Name', 'long_name'      ,'Value', 'MDV');
   nc(ifld).Attribute(2) = struct('Name', 'units'          ,'Value', 'm^3/m');
   nc(ifld).Attribute(3) = struct('Name', 'axis'           ,'Value', 'Y');
   nc(ifld).Attribute(4) = struct('Name', 'comment'        ,'Value', '`Momentaan Duin Volume` (MDV)'); 
   nc(ifld).Attribute(5) = struct('Name', '_FillValue'     ,'Value', NaN);
   
   
            ifld = ifld + 1;
   nc(ifld).Name         = 'P';
   nc(ifld).Nctype       = 'double';
   nc(ifld).Dimension    = {'Time','ID'};
   nc(ifld).Attribute(1) = struct('Name', 'long_name'      ,'Value', 'Probability of exceeding Xaf1990');
   nc(ifld).Attribute(2) = struct('Name', 'units'          ,'Value', '1/year');
   nc(ifld).Attribute(3) = struct('Name', 'axis'           ,'Value', 'Y');
   nc(ifld).Attribute(4) = struct('Name', 'comment'        ,'Value', 'Probability of exceeding (P) the `Maatgevend Afslagpunt van 1990` (Xaf1990)'); 
   nc(ifld).Attribute(5) = struct('Name', '_FillValue'     ,'Value', NaN);
   
   
               ifld = ifld + 1;
   nc(ifld).Name         = 'duneratio';
   nc(ifld).Nctype       = 'double';
   nc(ifld).Dimension    = {'Time','ID'};
   nc(ifld).Attribute(1) = struct('Name', 'long_name'      ,'Value', 'Characteristic dune ratio');
   nc(ifld).Attribute(2) = struct('Name', 'units'          ,'Value', '-');
   nc(ifld).Attribute(3) = struct('Name', 'axis'           ,'Value', 'Y');
   nc(ifld).Attribute(4) = struct('Name', 'comment'        ,'Value', '= mean( DH/(xDV-Xaf1990) )'); 
   nc(ifld).Attribute(5) = struct('Name', '_FillValue'     ,'Value', NaN);
   
   
               ifld = ifld + 1;
   nc(ifld).Name         = 'nearshore_slope';
   nc(ifld).Nctype       = 'double';
   nc(ifld).Dimension    = {'Time','ID'};
   nc(ifld).Attribute(1) = struct('Name', 'long_name'      ,'Value', 'Characteristic nearshore_slope');
   nc(ifld).Attribute(2) = struct('Name', 'units'          ,'Value', '-');
   nc(ifld).Attribute(3) = struct('Name', 'axis'           ,'Value', 'Y');
   nc(ifld).Attribute(4) = struct('Name', 'comment'        ,'Value', '= mean( (UL-LL)/(2*(xSea-MKL)) )'); 
   nc(ifld).Attribute(5) = struct('Name', '_FillValue'     ,'Value', NaN);
   
   
               ifld = ifld + 1;
   nc(ifld).Name         = 'beach_slope';
   nc(ifld).Nctype       = 'double';
   nc(ifld).Dimension    = {'Time','ID'};
   nc(ifld).Attribute(1) = struct('Name', 'long_name'      ,'Value', 'Characteristic beach_slope');
   nc(ifld).Attribute(2) = struct('Name', 'units'          ,'Value', '-');
   nc(ifld).Attribute(3) = struct('Name', 'axis'           ,'Value', 'Y');
   nc(ifld).Attribute(4) = struct('Name', 'comment'        ,'Value', '= mean( (UL-LL)/(2*(MKL-xLand)) )'); 
   nc(ifld).Attribute(5) = struct('Name', '_FillValue'     ,'Value', NaN);


%% 4 Create variables with attibutes

   for ifld=1:length(nc)
      disp(['Adding :',num2str(ifld),' ',nc(ifld).Name])
      nc_addvar(outputfile, nc(ifld));   
   end
   clear ifld
   
%% 5 Fill variables

   nc_varput(outputfile, 'Year' ,       OUTPUT.jaren);
   nc_varput(outputfile, 'ID' ,         OUTPUT.raaien);
   nc_varput(outputfile, 'BKL' ,        OUTPUT.BKL');
       nc_varput(outputfile, 'BKL_xRD' ,    BKL_xRD');
       nc_varput(outputfile, 'BKL_yRD' ,    BKL_yRD');
       nc_varput(outputfile, 'BKL_lat' ,    BKL_lat');
       nc_varput(outputfile, 'BKL_lon' ,    BKL_lon');
   nc_varput(outputfile, 'MKL' ,        OUTPUT.MKL');
       nc_varput(outputfile, 'MKL_xRD' ,    MKL_xRD');
       nc_varput(outputfile, 'MKL_yRD' ,    MKL_yRD');
       nc_varput(outputfile, 'MKL_lat' ,    MKL_lat');
       nc_varput(outputfile, 'MKL_lon' ,    MKL_lon');
   nc_varput(outputfile, 'xDV' ,        OUTPUT.xDV');
       nc_varput(outputfile, 'xDV_xRD' ,    xDV_xRD');
       nc_varput(outputfile, 'xDV_yRD' ,    xDV_yRD');
       nc_varput(outputfile, 'xDV_lat' ,    xDV_lat');
       nc_varput(outputfile, 'xDV_lon' ,    xDV_lon');
   nc_varput(outputfile, 'MDV' ,        OUTPUT.MDV');
   nc_varput(outputfile, 'MDL' ,        OUTPUT.MDL');
       nc_varput(outputfile, 'MDL_xRD' ,    MDL_xRD');
       nc_varput(outputfile, 'MDL_yRD' ,    MDL_yRD');
       nc_varput(outputfile, 'MDL_lat' ,    MDL_lat');
       nc_varput(outputfile, 'MDL_lon' ,    MDL_lon');
   nc_varput(outputfile, 'Xafslag' ,	OUTPUT.Xafslag');
       nc_varput(outputfile, 'Xafslag_xRD' ,Xafslag_xRD');
       nc_varput(outputfile, 'Xafslag_yRD' ,Xafslag_yRD');
       nc_varput(outputfile, 'Xafslag_lat' ,Xafslag_lat');
       nc_varput(outputfile, 'Xafslag_lon' ,Xafslag_lon');
   nc_varput(outputfile, 'P' ,          OUTPUT.P_f');
   nc_varput(outputfile, 'duneratio' ,  duinratio');
   nc_varput(outputfile, 'nearshore_slope' ,vooroeverhelling');
   nc_varput(outputfile, 'beach_slope' ,strandhelling');


%% 6 Check

   nc_dump(outputfile);


%%
return
























%%
%%
%% EXAMPLE CODE:
%%
%%


%% Set default options, and overrule with user defined keyword,value pairs
   OPT.refdatenum = datenum(1970,1,1);
   OPT.fillvalue  = nan;
   OPT.filename   = 'nc_cf_example2nc.inp';
   
   OPT = setproperty(OPT,varargin);

%% 0. Read raw data
% Make a function that returns all data + meta-data in one struct, e.g.:
%
%  * D             = knmi_potwind(OPT.filename);
%  * D             = knmi_etmgeg (OPT.filename);
%  * D             = donar_read  (OPT.filename);
%  
% Below is just a example taht creates soem ranodm ata:

   D.datenum     = floor(now) + [0:2:24]./24;
   D.version     = 0;
   D.lat         = [ 4  5  6];
   D.lon         = [52 53 54 55];
   D.temperature = repmat(nan,[length(D.datenum) length(D.lat) length(D.lon)]);
   for i=1:length(D.datenum)
   for j=1:length(D.lat)
   for k=1:length(D.lon)
      D.temperature(i,j,k) = i.^2 + (D.lat(j).^2 + D.lon(k));
   end
   end
   end
   D.timezone    = '+00:00';

%% 1. Create file
   outputfile = [filename(OPT.filename),'.nc'];
   nc_create_empty (outputfile); % only change extension with respect to input file

%% 2. Add overall meta info

   %%  CF convention   
   % <http://cf-pcmdi.llnl.gov/documents/cf-conventions/1.4/cf-conventions.html#description-of-file-contents>

   nc_attput(outputfile, nc_global, 'title'         , '');
   nc_attput(outputfile, nc_global, 'institution'   , 'Deltares');
   nc_attput(outputfile, nc_global, 'source'        , '');
   % Insert SVN keyword $HeadURL$ that will list the url of this mfile. This
   % approach automatically ensures that the name of the script that made
   % the netCDF file is included in the netCDF file itself.
   nc_attput(outputfile, nc_global, 'history'       , ['tranformation to NetCDF: $HeadURL$']);
   % Provide web-link to originator of dataset (can be url in OpenEarthRawData).
   nc_attput(outputfile, nc_global, 'references'    , '<http://openearth.deltares.nl>');
   nc_attput(outputfile, nc_global, 'email'         , '');
   
   nc_attput(outputfile, nc_global, 'comment'       , '');
   nc_attput(outputfile, nc_global, 'version'       , D.version);
						    
   nc_attput(outputfile, nc_global, 'Conventions'   , 'CF-1.4');
   nc_attput(outputfile, nc_global, 'CF:featureType', 'stationTimeSeries');  % https://cf-pcmdi.llnl.gov/trac/wiki/PointObservationConventions
   
   %%  OpenEarth convention

   nc_attput(outputfile, nc_global, 'terms_for_use' , 'These data can be used freely for research purposes provided that the following source is acknowledged: KNMI.');
   nc_attput(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.');

%% 2. Create dimensions

   nc_add_dimension(outputfile, 'time'  , length(D.datenum))
   nc_add_dimension(outputfile, 'lat'   , length(D.lat))
   nc_add_dimension(outputfile, 'lon'   , length(D.lon))

%% 3. Create variables

   clear nc
   ifld = 0;
   
   %%% Define dimensions in this order: [time,z,y,x]
   %
   % * For standard names vocabulary by CF group see:
   %   <http://cf-pcmdi.llnl.gov/documents/cf-standard-names/standard-name-table/current/standard-name-table>
   %   From this list the following quantities are used in this mfule:
   %   contains: latitude, longitude
   % * For standard units vocabulary UDUNITS by UNIDATA see:
   %   <http://www.unidata.ucar.edu/software/udunits/>

   %%% Latitude
   % Prescribed as dimension associated with variable by CF convention in:
   % <http://cf-pcmdi.llnl.gov/documents/cf-conventions/1.4/cf-conventions.html#latitude-coordinate>
   
      ifld = ifld + 1;
   nc(ifld).Name         = 'lat';   % This name is required when extracing the data with nc_varget(,'lat').
   nc(ifld).Nctype       = 'float'; % no double needed
   nc(ifld).Dimension    = {'lat'}; % Should conform with dimension lat defined above
   nc(ifld).Attribute(1) = struct('Name', 'long_name'      ,'Value', 'latitude');      % Name free of choice, will appear in plots
   nc(ifld).Attribute(2) = struct('Name', 'units'          ,'Value', 'degrees_north'); % Note: 1st type of degrees, chosen from UDUNITS list.
   nc(ifld).Attribute(3) = struct('Name', 'standard_name'  ,'Value', 'latitude');      % Prescribed by CF convention and CF standard name table.

   %%% Longitude
   % Prescribed as dimension associated with variable by CF convention in:
   % <http://cf-pcmdi.llnl.gov/documents/cf-conventions/1.4/cf-conventions.html#longitude-coordinate>
   
      ifld = ifld + 1;
   nc(ifld).Name         = 'lon';   % This name is required when extracing the data with nc_varget(,'lat').
   nc(ifld).Nctype       = 'float'; % no double needed
   nc(ifld).Dimension    = {'lon'}; % Should conform with dimension lat defined above
   nc(ifld).Attribute(1) = struct('Name', 'long_name'      ,'Value', 'longitude');    % Name free of choice, will appear in plots
   nc(ifld).Attribute(2) = struct('Name', 'units'          ,'Value', 'degrees_east'); % Note: 2nd type of degrees, chosen from UDUNITS list.
   nc(ifld).Attribute(3) = struct('Name', 'standard_name'  ,'Value', 'longitude');    % Prescribed by CF convention and CF standard name table.

   %%% local x
   % Should be associated with a coordinate system.
  
      ifld = ifld + 1;
   nc(ifld).Name         = 'x';
   nc(ifld).Nctype       = 'float'; % no double needed
   nc(ifld).Dimension    = {'lat'};
   nc(ifld).Attribute(1) = struct('Name', 'long_name'      ,'Value', 'x');                       % Name free of choice, will appear in plots
   nc(ifld).Attribute(2) = struct('Name', 'units'          ,'Value', 'meter');                   % Chosen from UDUNITS list.
   nc(ifld).Attribute(3) = struct('Name', 'standard_name'  ,'Value', 'projection_x_coordinate'); % Prescribed by CF standard name table.


   %%% local y
   % Should be associated with a coordinate system.
   
      ifld = ifld + 1;
   nc(ifld).Name         = 'y';
   nc(ifld).Nctype       = 'float'; % no double needed
   nc(ifld).Dimension    = {'lon'};
   nc(ifld).Attribute(1) = struct('Name', 'long_name'      ,'Value', 'y');                       % Name free of choice, will appear in plots
   nc(ifld).Attribute(2) = struct('Name', 'units'          ,'Value', 'meter');                   % Chosen from UDUNITS list.
   nc(ifld).Attribute(3) = struct('Name', 'standard_name'  ,'Value', 'projection_y_coordinate'); % Prescribed by CF standard name table.

   %%% Time
   % <http://cf-pcmdi.llnl.gov/documents/cf-conventions/1.4/cf-conventions.html#time-coordinate>
   % 
   % time is a dimension, so there are two options:
   %
   % * The variable name needs the same as the dimension
   %   <http://cf-pcmdi.llnl.gov/documents/cf-conventions/1.4/cf-conventions.html#id2984551>
   %
   % * There needs to be an indirect mapping through the coordinates attribute
   %   <http://cf-pcmdi.llnl.gov/documents/cf-conventions/1.4/cf-conventions.html#id2984605>
   
   OPT.timezone = timezone_code2iso(D.timezone);

      ifld = ifld + 1;
   nc(ifld).Name         = 'time';
   nc(ifld).Nctype       = 'double'; % float not sufficient as datenums are big: doubble
   nc(ifld).Dimension    = {'time'}; % {'locations','time'} % does not work in ncBrowse, nor in Quickplot (is indirect time mapping)
   nc(ifld).Attribute(1) = struct('Name', 'long_name'      ,'Value', 'time');
   nc(ifld).Attribute(2) = struct('Name', 'units'          ,'Value',['days since ',datestr(OPT.refdatenum,'yyyy-mm-dd'),' 00:00:00 ',OPT.timezone]);
   nc(ifld).Attribute(3) = struct('Name', 'standard_name'  ,'Value', 'time');
   nc(ifld).Attribute(4) = struct('Name', '_FillValue'     ,'Value', OPT.fillvalue);
  %nc(ifld).Attribute(5) = struct('Name', 'bounds'         ,'Value', '');
   
  %%% Parameters with standard names
  % * <http://cf-pcmdi.llnl.gov/documents/cf-standard-names/standard-name-table/current/>

      ifld = ifld + 1;
   nc(ifld).Name         = 'T';
   nc(ifld).Nctype       = 'float';
   nc(ifld).Dimension    = {'time','lat','lon'};
   nc(ifld).Attribute(1) = struct('Name', 'long_name'      ,'Value', 'air temperature');
   nc(ifld).Attribute(2) = struct('Name', 'units'          ,'Value', 'degree_Celsius'); % Note: 3rd type of degrees
   nc(ifld).Attribute(3) = struct('Name', 'standard_name'  ,'Value', 'air_temperature');
   nc(ifld).Attribute(4) = struct('Name', '_FillValue'     ,'Value', OPT.fillvalue);

%% 4 Create variables with attibutes
% When variable definitons are created before actually writing the
% data in the next cell, netCDF can nicely fit all data into the
% file without the need to relocate any info. So there is no 
% need to use NC_PADHEADER.

   for ifld=1:length(nc)
      disp(['Adding :',num2str(ifld),' ',nc(ifld).Name])
      nc_addvar(outputfile, nc(ifld));   
   end
   
%% 5 Fill variables

   nc_varput(outputfile, 'time' , D.datenum-OPT.refdatenum);
   nc_varput(outputfile, 'lat'  , D.lat);
   nc_varput(outputfile, 'lon'  , D.lon);
   nc_varput(outputfile, 'T'    , D.temperature);

%% 6 Check

   nc_dump(outputfile);
   
%% For more information see: <OpenEarth.Deltares.nl>
