% netcdfFaalkans
%
% Generates faalkans.nc which contains the following output from HKV: 
% faalkans    = probability of failure for dune retreat (first dune row)
%
% faalkans is obtained by HKV.
%     The information are available from Deltares repository @
%     https://repos.deltares.nl/repos/OpenEarthRawData/trunk/rijkswaterstaat/kustlijnkaartenboek/raw/mkl/
%
% Vincent Vuik
% v.vuik@hkv.nl
% 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%


%% settings

% jarkus database for transect info
%OPT.ncpath3      = 'http://opendap.deltares.nl/thredds/dodsC/opendap/rijkswaterstaat/jarkus/profiles/';
OPT.ncpath       = 'http://opendap.tudelft.nl/thredds/dodsC/data2/deltares/rijkswaterstaat/jarkus/profiles/';
OPT.Filename     = 'transect.nc';
OPT.fullncpath   = [OPT.ncpath OPT.Filename];

% general settings
OPT.fillvalue    = NaN;
OPT.refdatenum   = datenum(1970,1,1);
OPT.outputfile   = 'faalkans.nc';
OPT.disp         = true;
OPT.dump         = true;

%% Read data form jarkus transect.nc
id         = nc_varget([OPT.fullncpath],'id');
along      = nc_varget([OPT.fullncpath],'alongshore');
areaname   = nc_varget([OPT.fullncpath],'areaname');
areacode   = nc_varget([OPT.fullncpath],'areacode');
time       = nc_varget([OPT.fullncpath],'time');

jarkusyear = str2num(datestr(nc_cf_time([OPT.fullncpath]),'yyyy'));

%% create the matrix for each parameter, with a structure similar to the jarkus file

bestand               = ['faalkans_2014.xlsx'];
raaien                = xlsread(bestand,'totaal','A2:A1983');
jaren                 = xlsread(bestand,'totaal','B1:AU1');
faalkans              = xlsread(bestand,'totaal','B2:AU1983');
kwaliteit             = xlsread(bestand,'totaal','AW2:AX1983'); % quality code for trend & magnitude

for i=1:length(raaien);
    lokid(i)          = find(id==raaien(i)); 
end
Nj                    = length(jaren);
Nr                    = length(id);
Pf                    = NaN.*ones(Nj,Nr);
Pf( 1:Nj,lokid)       = faalkans';
Qtrend                = NaN.*ones(1,Nr);
Qmgntd                = NaN.*ones(1,Nr);
Qtrend(1,lokid)       = kwaliteit(:,1)';
Qmgntd(1,lokid)       = kwaliteit(:,2)';

%% 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'         , 'probability of failure of dunes');
nc_attput(OPT.outputfile, nc_global, 'institution'   , 'HKV Lijn in Water');
nc_attput(OPT.outputfile, nc_global, 'email'         , 'v.vuik@hkv.nl');
nc_attput(OPT.outputfile, nc_global, 'version'       , '4.0');
nc_attput(OPT.outputfile, nc_global, 'report'        , 'HKV (2012), Indicatoren voor kustlijnzorg, report pr2063.30');

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

%% add variables
clear nc;
i = 0;

%% add epsg

% smaller variables first
[lon_dummy,lat_dummy,OPT2] = convertCoordinates(0,0,'CS1.code',28992,'CS2.code',4326);
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');

% 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');

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

% quality code trend
i = i + 1;
nc(i).Name         = 'reliability_trend_probability_failure';
nc(i).Nctype       = 'double';
nc(i).Dimension    = {'alongshore'};
nc(i).Attribute(1) = struct('Name', 'long_name' ,'Value', 'Reliability of computed trend 1965-2010 in probability of failure');
nc(i).Attribute(2) = struct('Name', 'comment'   ,'Value', '0=bad, 1=poor, 2=good, -1=not assigned');

% quality code magnitude
i = i + 1;
nc(i).Name         = 'reliability_magnitude_probability_failure';
nc(i).Nctype       = 'double';
nc(i).Dimension    = {'alongshore'};
nc(i).Attribute(1) = struct('Name', 'long_name' ,'Value', 'Reliability of computed magnitude of probability of failure');
nc(i).Attribute(2) = struct('Name', 'comment'   ,'Value', '0=bad, 1=poor, 2=good, -1=not assigned');

%% 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, 'time'                    , time);

nc_varput(OPT.outputfile, 'probability_failure'     , Pf);
nc_varput(OPT.outputfile, 'reliability_trend_probability_failure'     , Qtrend);
nc_varput(OPT.outputfile, 'reliability_magnitude_probability_failure' , Qmgntd);

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

%% end