function Process_hydrodynamic_runs(info, dirs)
% --- location of x = 0
Tstart = 10*60; %3600*2;
Tend   = 1800; %3600*5;

boers_id = {'1A','1B','1C'};

for ll=1:3

L1 = load([dirs.data filesep boers_id{ll} filesep 'Boers_data_processed.mat']);

if strcmp(boers_id{ll},'1A')
    Tp =  2.1;
    name = '1A_SB_BOI'
elseif strcmp(boers_id{ll},'1B')
    Tp = 2.1;
    name = '1B_SB_BOI'
elseif strcmp(boers_id{ll},'1C')
    Tp = 3.4;
    name = '1C_SB_BOI'
end



fname            = [dirs.output filesep name  filesep 'xboutput.nc'];

point_zs         = nc_varget(fname, 'point_zs'); point_zs = point_zs(:,:);
point_u          = nc_varget(fname, 'point_u');  point_u = point_u(:,:);
point_zb         = nc_varget(fname, 'point_zb'); point_zb = squeeze(point_zb(1,:));
globalx          = nc_varget(fname, 'globalx');
t                = nc_varget(fname, 'globaltime'); t = t(:);
zb               = nc_varget(fname, 'zb'); 
zs               = nc_varget(fname, 'zs'); %M.zs = squeeze(M.zs(1,:,:)); M.zs = M.zs(:,:);
point_time       = nc_varget(fname, 'pointtime'); point_time = point_time(:,:);
pointx           = nc_varget(fname, 'pointx');
point_H          = nc_varget(fname, 'point_H'); point_H = point_H(:,:);
globaltime       = nc_varget([fname],'globaltime');


id_wanted           = point_time>=Tstart & point_time<= Tend;

% seperation between IG and HF
fband   = [(1/Tp/20) (1/Tp/2) 10];
Fs      = 1/(point_time(3)-point_time(2));

% --- initialise arrays

LF_eta              = NaN(length(point_time(id_wanted)),length(pointx)); 
VLF_eta             = NaN(length(point_time(id_wanted)),length(pointx)); 
LF_tot_eta          = NaN(length(point_time(id_wanted)),length(pointx)); 

Hm0_HF      = NaN(1,length(pointx)); 
Hm0_LF      = NaN(1,length(pointx)); 
Hm0_VLF     = NaN(1,length(pointx)); 
Hm0_LF_tot  = NaN(1,length(pointx)); 
Groupiness  = NaN(1,length(pointx)); 
setup       = NaN(1,length(pointx)); 

for jj = 1:length(pointx)
        % --- Filtering LF freq
        LF_eta(:,jj)    = bandpassfilter(point_zs(id_wanted,jj),Fs,[fband(1) fband(2)]);

        % --- VLF
        VLF_eta(:,jj)   = bandpassfilter(point_zs(id_wanted,jj),Fs,[0 fband(1)]);

        LF_tot_eta(:,jj)   = bandpassfilter(point_zs(id_wanted,jj),Fs,[0 fband(2)]);

        % Calculating wave heights + wave groups
        Hm0_HF(1,jj)                = rms(point_H(id_wanted,jj))*sqrt(2); 
        Hm0_LF(1,jj)                = 4 * std(LF_eta(:,jj));
        Hm0_VLF(1,jj)               = 4 * std(VLF_eta(:,jj));
        Hm0_LF_tot(1,jj)            = 4 * std(LF_tot_eta(:,jj));

        Groupiness(1,jj)     = std( point_H(id_wanted,jj) )/mean( point_H(id_wanted,jj) );

        % setup
        setup(1,jj) = nanmean(point_zs(id_wanted,jj));  
end

[C,ia,ic] = unique(pointx);
figure(1)
set(gcf, 'Position', [10 10 600 600]);
subplot('position',[0.1 0.7 0.7 0.2]);
plot(L1.x,L1.Hm0_HF,'o','markersize',5, 'Color', 'k',  'DisplayName', ['H_{m0} HF data' ]); hold on;
plot(L1.x,L1.Hm0_LF,'<','markersize',5, 'Color', 'k',  'DisplayName', ['H_{m0} LF data' ]); hold on;

[r2 sci_HF(ll) relbias_HF(ll) bss] = xb_skill([L1.x', L1.Hm0_HF'], [pointx(ia),Hm0_HF(ia)']);
plot(pointx,Hm0_HF,'b-','linewidth',2, 'DisplayName', 'H_{m0 HF}');hold on;
hold on
[r2 sci_LF(ll) relbias_LF(ll) bss] = xb_skill([L1.x', L1.Hm0_LF'], [pointx(ia),Hm0_LF(ia)']);
plot(pointx,Hm0_LF,'r-','linewidth',2,'DisplayName', 'H_{m0} LF'); hold on;

leg=legend('Location',[.85 .7 .1 .2]);
set(leg);
grid on
ylabel('Wave height (m)','FontSize',12)
xlim([0 30])

subplot('position',[0.1 0.4 0.7 0.2]);
[r2 sci_eta(ll) relbias_eta(ll) bss] = xb_skill([L1.x', L1.setup'], [pointx(ia),setup(ia)']);
plot(L1.x, L1.setup, 'o', 'markersize',5,'Color', 'k', 'DisplayName', 'Setup data'); hold on;
plot(pointx,setup,'b-','linewidth',2, 'DisplayName', 'setup'); hold on;
leg=legend('Location',[.85 .4 .1 .2]);
set(leg)
grid on
ylabel('Setup (m)','FontSize',12)
xlim([0 30])

subplot('position',[0.1 0.1 0.7 0.2]);
plot(globalx,zb(1,:),'k','linewidth',2, 'DisplayName', 'Bed');hold on;
plot(globalx,squeeze( zs(1,1,:)),'b','linewidth',2, 'DisplayName', 'mwl')
leg=legend('Location',[.85 .1 .1 .2]);
set(leg)
xlabel('x (m)','FontSize',12)
ylabel('depth (m)','FontSize',12)
grid on
xlim([0 30])

print('-dpng','-r300',[dirs.output,filesep 'Overview',boers_id{ll}])
close all

end

figure();
subplot('position',[0.1 0.7 0.7 0.2]);
bar(categorical(boers_id),  [relbias_HF', sci_HF'], 'grouped')
grid on
title('H_{m0} HF','FontSize',8)
legend({'rel. bias', 'SCI'},'FontSize',8,'Location',[.85 .7 .1 .2])
ylim([-0.05 0.1])
subplot('position',[0.1 0.4 0.7 0.2]);
bar(categorical(boers_id),  [relbias_LF', sci_LF'], 'grouped')
grid on
title('H_{m0} LF','FontSize',8)
ylim([-0.4 0.4])

subplot('position',[0.1 0.1 0.7 0.2]);
bar(categorical(boers_id),  [relbias_eta', sci_eta'], 'grouped')
grid on
title('setup','FontSize',8)
ylim([-1 1])

print('-dpng','-r300',[dirs.output,filesep 'Overview'])



dlm = {'$H_{m0,HF}$ Rel. bias','$H_{m0,HF}$ SCI','$H_{m0,LF}$ Rel. bias','$H_{m0,LF}$ SCI','setup Rel. bias','setup SCI'};
matrix2latex([relbias_HF; sci_HF; relbias_LF; sci_LF; relbias_eta; sci_eta], 'filename', 'skill.tex', ...
    'rowjustification','p{1.5cm}',...
    'caption',  'The statistical scores for the Boers experiments. The scatter index (SCI) and relative bias (rel. bias) are shown for the short-wave height, infragravity wave height and setup.', ...
    'rowlabel', dlm, ...
    'collabel', boers_id, ...
    'format',	{'%s', '%4.2f', '%4.2f','%4.2f'});

end