OPT = struct(...
    ...'maindir', 'F:\Sources\OpenEarthModels\tudelft\dune_erosion_coastline_curvature\simulations\',...
    'maindir', 'd:\heijer\sources\dune_erosion_coastline_curvature\simulations\00000000\xbeach\',...
    'destdir', fullfile(fileparts(which(mfilename)), '..', '..', 'dune_erosion_bayesian_network', 'netica'),...
    'regexp_select', 'run(36).*?0\\xbsum\.nc',...
    'regexp_legend', 'run\d{2}.*?radius\d{4,5}(.*?rev\d{4}\D?)?',...
    'strrep_legend', {{{'\' '; '} {'s0' 's '} {'00' '00 m'}}},...
    'save', false ...
    );

[dirs dircont files] = dirlisting(OPT.maindir);

ncfiles = files(~cellfun(@isempty, regexp(files, OPT.regexp_select)));

for iradii = 1:length(ncfiles)
    ncfile = ncfiles{iradii};
    
    result(iradii) = decc_analyse(ncfile,...
        's', false,...
        'zb', true,...
        'Bal', false,...
        'Fy_mean', false,...
        'v_mean', false,...
        'u_mean', false,...
        'Svsg_cs', false,...
        ...'ee_mean', true,...
        'contour', -2.2);
end

m = regexp(ncfiles, '(?<=radius)\d+', 'match');
radii = cellfun(@(x,y) repmat(str2double(x), 1, y), [m{:}], num2cell(cellfun(@length, {result.Ero})), 'UniformOutput', false);

erofac_straight_coast = @(a) 1+.08*sind(a*2);
a = rand(1,1000)*90;
netica_write_case_file(fullfile(OPT.destdir, 'curvature_wave_angle.cas'),...
    'phi', [abs([result.phi]) a],...
    'erosionfactor', [[result.Ero]/170 erofac_straight_coast(a)],...
    'radius', [[radii{:}] 2e4*ones(size(a))]);

%% add jarkus info
tr = jarkus_transects(...
    'year', 2010,...
...    'areacode', 2:12,...
    'output', {'areacode', 'angle', 'id', 'time', 'x', 'y', 'cross_shore', 'altitude'});
[curvatures radii relativeAngle distances] = jarkus_curvature(tr, 'method', 'circle');
id = 1:length(radii);
radii(~isfinite(radii) | radii<=0) = 2e4;
% radii(radii<0) = 0;
netica_write_case_file(fullfile(OPT.destdir, 'jarkus_curvature_ori.cas'),...
    'radius', radii,...
    'wavedir', rand(length(tr.areacode(id)),1)*50+270,...
    'areacode', tr.areacode(id),...
    'theta', tr.angle(id));