%% Runup location alongshore

for line = 1:22
    year = 2009;
    month = 7:12;
    for i = 1: length(month)
        [cross_shore,altitude(:,i)] = combine_surveys(year,month(i),line);
    end
    year = 2010;
    month = 1:5;
    for i = 1: length(month)
        [cross_shore,altitude(:,i+6)] = combine_surveys(year,month(i),line);
    end
    pro_var(:,line) =  nanvar(altitude')';
    pro_mean(:,line) = nanmean(altitude')';
    if line ==11 %remove outlier
        index_runup(1,line) = find(cross_shore > 77 & cross_shore < 150 & pro_var(:,line) > 0.01 ,1,'first');
    else
        index_runup(1,line) = find(cross_shore > 50 & cross_shore < 150 & pro_var(:,line) > 0.01 ,1,'first');
    end
%     index_closure(1,line) = find(pro_var);
    cross_runup_lim(:,line) = cross_shore(index_runup(1,line));
    Z_runup_lim(:,line) = pro_mean(index_runup(1,line),line);

end
dist_alongshore = nc_varget('vlugtenburg.nc','dist_alongshore');




figure
subplot(2,1,1)
plot(dist_alongshore,Z_runup_lim,'o:r','linewidth',2); figure(gcf)
ylabel({'Vertical coordiante of runup line'; '[m NAP]'})
subplot(2,1,2)
plot(dist_alongshore,cross_runup_lim,'s:b','linewidth',2); figure(gcf)
xlabel('Along shore distance [m]')
ylabel({'Cross shore location runup line';'[m]'})

%% Slope in intertidal area
%%% we hebben de cross shore locatie van de Runup
%%% nu eerst zeewaardse referentie kiezen.
%%% Ik probeer nu eerst z = 0 m NAP
if 0
    zee_ref = 4;

    [value0 index1] = min(abs(pro_mean+zee_ref));
    cross_zero = cross_shore(index1);

    for i=1:1:22

        p(i,:) = polyfit(cross_shore(index_runup(i):index1(i)),pro_mean(index_runup(i):index1(i),i),1);
        if 0
            figure
            plot(cross_shore,pro_mean(:,i),'b','linewidth',2)
            hold all
            plot(cross_shore(index_runup(i):index1(i)),polyval(p(i,:),cross_shore(index_runup(i):index1(i))),'r','linewidth',2)
            hline(0)
            title(['Line ' num2str(i) ' -- P = ' num2str(p(i,:))])
            ylim([-10 10])
        end
    end

    figure
    plot(dist_alongshore,p(1:22,1), 'DisplayName', 'p(1:22,1)', 'YDataSource', 'p(1:22,1)'); figure(gcf)
    title(['sea reference = ' num2str(-zee_ref) ' m NAP'])
    ylabel('profile gradient')
    xlabel('Along shore distance [m]')
end
%% make summarizing plots
figure(56)
close(56)
zee_ref = 4;
[value0 index1] = min(abs(pro_mean+zee_ref));
cross_zero = cross_shore(index1);
% index_runup = 80*ones(1,22)
for i=1:1:22
    p1(i,:) = polyfit(cross_shore(index_runup(i):index1(i)),pro_mean(index_runup(i):index1(i),i),1);
end

for i=3
    figure(55)
    subplot(2,1,1)
    plot(cross_shore,pro_mean(:,i),'b','linewidth',1)
    hold all
    plot(cross_shore(index_runup(i):index1(i)),polyval(p1(i,:),cross_shore(index_runup(i):index1(i))),'k','linewidth',2)
    %         hline(0)
    %         title(['Line ' num2str(i) ' -- P = ' num2str(p1(i,:))])
    ylim([-10 10])
    set(gca,'XDir','reverse')
end
for i=17
    figure(55)
    subplot(2,1,2)
    plot(cross_shore,pro_mean(:,i),'b','linewidth',1)
    hold all
    plot(cross_shore(index_runup(i):index1(i)),polyval(p1(i,:),cross_shore(index_runup(i):index1(i))),'k','linewidth',2)
    %         hline(0)
    %         title(['Line ' num2str(i) ' -- P = ' num2str(p1(i,:))])
    ylim([-10 10])
    set(gca,'XDir','reverse')
end


figure(56)
% subplot(2,1,1)
plot(dist_alongshore,p1(1:22,1),'o:k','linewidth',1); figure(gcf)
% title(['sea reference = ' num2str(-zee_ref) ' m NAP'])
ylabel('profile gradient')
ylim([-0.0261 -0.0080])
text(dist_alongshore(22)+50,p1(22)-0.003,['Foreshore until ' num2str(-zee_ref) ' m NAP'] ,'HorizontalAlignment','left')%,'Color',[0 1 0])
% xlabel('Along shore distance [m]')

zee_ref = 8;
[value0 index0] = min(abs(pro_mean+zee_ref));
cross_zero = cross_shore(index0);
for i=1:1:22
    p2(i,:) = polyfit(cross_shore(index_runup(i):index0(i)),pro_mean(index_runup(i):index0(i),i),1);
end
for i=3
    figure(55)
    subplot(2,1,1)
    %         plot(cross_shore,pro_mean(:,i))
    %         hold all
    plot(cross_shore(index_runup(i):index0(i)),polyval(p2(i,:),cross_shore(index_runup(i):index0(i))),'r','linewidth',2)
    hline(0)
    %         title(['Line ' num2str(i) ' -- P = ' num2str(p2(i,:))])
    title(['Line ' num2str(i) ' (= ' num2str(floor(dist_alongshore(i))) ' m alongshore distance)'])
    %     xlabel('Cross shore distance [m]')
    ylabel('Height [m] NAP')
    ylim([-10 10])
end
for i=17
    figure(55)
    subplot(2,1,2)
    plot(cross_shore,pro_mean(:,i),'b','linewidth',1)
    hold all
    plot(cross_shore(index_runup(i):index0(i)),polyval(p2(i,:),cross_shore(index_runup(i):index0(i))),'r','linewidth',2)
    hline(0)
    %         title(['Line ' num2str(i) ' -- P = ' num2str(p2(i,:))])
    title(['Line ' num2str(i) ' (= ' num2str(floor(dist_alongshore(i))) ' m alongshore distance)'])
    ylim([-10 10])
    xlabel('Cross shore distance [m]')
    ylabel('Height [m] NAP')

end


%
figure(56)
% subplot(2,1,2)
hold all
plot(dist_alongshore,p2(1:22,1),'s:r','linewidth',1); figure(gcf)
% title(['sea reference = ' num2str(-zee_ref) ' m NAP'])
% ylabel('profile gradient')
xlabel('Along shore distance [m]')
text(dist_alongshore(22)+50,p2(22)-0.003,['Active profile until ' num2str(-zee_ref) ' m NAP'],'HorizontalAlignment','left','Color',[1 0 0])
% figure
% plot(p(1:22,2), 'DisplayName', 'p(1:22,1)', 'YDataSource', 'p(1:22,1)'); figure(gcf)

%% correlate gradients with runup heigth
% get rid of unrealistic value first

val_index = Z_runup_lim<3;
figure(65);

subplot(1,2,2)
plot(p2(val_index,1),Z_runup_lim(val_index),'*k');
hold all
set(gca,'XDir','reverse')
xlabel('profile gradient')
ylabel({'Vertical coordiante of runup line'; '[m NAP]'})

% fit line !!
runupVSgradient = polyfit(p2(val_index,1)',Z_runup_lim(val_index),1);

%now plot result in figure
plot(p2(val_index,1)',polyval(runupVSgradient,p2(val_index,1)'),'r')
R=corrcoef(polyval(runupVSgradient,p2(val_index,1)'),Z_runup_lim(val_index))
ylim([1.4 3.4])
title ('off-shore slope')




figure(65);
subplot(1,2,1)
plot(p1(val_index,1),Z_runup_lim(val_index),'xk');
hold all
set(gca,'XDir','reverse')
xlabel('profile gradient')
ylabel({'Vertical coordiante of runup line'; '[m NAP]'})
title ('in-shore slope')

% fit line !!
runupVSgradient = polyfit(p1(val_index,1)',Z_runup_lim(val_index),1);

%now plot result in figure
plot(p1(val_index,1)',polyval(runupVSgradient,p1(val_index,1)'),'b')
R=corrcoef(polyval(runupVSgradient,p1(val_index,1)'),Z_runup_lim(val_index))
ylim([1.4 3.4])

set(gcf,'Position',[360   692   653   230])


%% plot mean contour plot
[along_i cross_i] = meshgrid(dist_alongshore,cross_shore);

figure
[xi,yi] = meshgrid(dist_alongshore,cross_shore);
[C,h]=contour(yi,xi,pro_mean,[-8:3]);
hold on
contour(yi,xi,pro_mean,[0 0],'linewidth',3)
plot(cross_i, along_i, '-.','Color',[0.9 0.9 0.9])%,'linewidth',5)
for k = 1:22;text(cross_shore(1)+200,dist_alongshore(k),num2str(k),...
        'HorizontalAlignment','right');end
title('Mean topography in local axis system')
ylabel('along shore [m]')
xlabel('cross shore [m]')
axis equal
xlim([0 1050])
set(h,'ShowText','on')
set(gca,'XDir','reverse')


%
% line=10
% figure;
% plot(cross_shore,pro_var(:,line));
% figure;
% plot(cross_shore,pro_mean(:,line));

% determine runup limit
% cross_runup_lim = cross_shore(find(cross_shore > 0 & cross_shore < 200 & pro_var > 0.005 ,1,'first'));
% Z_runup_lim = pro_mean(find(cross_shore == cross_runup_lim))

%% icce plots
if 1
%%
    
    %%
    figure(77)
    % xlabel('Along shore distance [m]')
    plot(dist_alongshore,p1(1:22,1),'o:g','linewidth',1); figure(gcf)
    % title(['sea reference = ' num2str(-zee_ref) ' m NAP'])
    ylabel('profile gradient')
    xlabel('Along shore distance [m]')
    ylim([-0.025 -0.019])
    title(['Nearshore until -4 m NAP'])%,'Color',[0 1 0])
    % text(dist_alongshore(22)+50,p1(22)-0.003,['Nearshore until -4 m NAP'] ,'HorizontalAlignment','left')%,'Color',[0 1 0])ylim([-0.025 -0.019])
    set(gcf,'position',[360 502 400 160])
        xlim([0 1800])



    figure(78)
    plot(dist_alongshore,p2(1:22,1),'s:m','linewidth',1); figure(gcf)
    % title(['sea reference = ' num2str(-zee_ref) ' m NAP'])
    % ylabel('profile gradient')
    xlabel('Along shore distance [m]')
    title(['Active profile until -8 m NAP'])
    % text(dist_alongshore(22)+50,p2(22)-0.003,['Active profile until -8 m
    % NAP'],'HorizontalAlignment','left','Color',[1 0 0])
    ylabel('profile gradient')
    ylim([-0.015 -0.009])
    set(gcf,'position',[360 502 400 160])
        xlim([0 1800])

    figure(111)
    plot(dist_alongshore,Z_runup_lim,'o:r','linewidth',2); figure(gcf)
    ylabel({'Vertical level of B_{ma}'; '[m NAP]'})
    xlabel('Along shore distance [m]')
    set(gcf,'position',[360 502 400 160])
    xlim([0 1800])
    ylim([1.8 2.8])

    %%
    figure(551)
    for i=3
        %     subplot(2,1,1)
        plot(cross_shore,pro_mean(:,i),'b','linewidth',2)
        hold all
        plot(cross_shore(index_runup(i):index1(i)),polyval(p1(i,:),cross_shore(index_runup(i):index1(i))),'g','linewidth',2)
        hline(0)
        %         title(['Line ' num2str(i) ' -- P = ' num2str(p1(i,:))])
        ylim([-10 10])
        set(gca,'XDir','reverse')
        set(gcf,'position',[360 502 400 160])
        title(['Line ' num2str(i) ' (= ' num2str(floor(dist_alongshore(i))) ' m alongshore distance)'])
        ylim([-10 10])
        xlim([0 1050])
        xlabel('Cross shore distance [m]')
        ylabel('Height [m] NAP')
    end
%%
        figure(561)
        hold off
        plot(p1(val_index,1),Z_runup_lim(val_index),'xk');
hold all
set(gca,'XDir','reverse')
xlabel('profile gradient')
ylabel({'Level of B_{ma}'; '[m NAP]'})
title ('Nearshore corr')

% fit line !!
runupVSgradient = polyfit(p1(val_index,1)',Z_runup_lim(val_index),1);

%now plot result in figure
plot(p1(val_index,1)',polyval(runupVSgradient,p1(val_index,1)'),'b','linewidth',2)
R=corrcoef(polyval(runupVSgradient,p1(val_index,1)'),Z_runup_lim(val_index))
ylim([1.4 3.4])
xlim([-0.0238 -0.021])

        set(gcf,'position',[360 502 200 160])        
        
%%

    figure(552)
    hold off
    for i=3
        %     subplot(2,1,1)
        plot(cross_shore,pro_mean(:,i),'b','linewidth',2)
        hold all
    plot(cross_shore(index_runup(i):index0(i)),polyval(p2(i,:),cross_shore(index_runup(i):index0(i))),'m','linewidth',2)
          plot(cross_shore(index_runup(i):index1(i)),polyval(p1(i,:),cross_shore(index_runup(i):index1(i))),'g','linewidth',2)
    hline(0)
        %         title(['Line ' num2str(i) ' -- P = ' num2str(p1(i,:))])
        ylim([-10 10])
        set(gca,'XDir','reverse')
        set(gcf,'position',[360 502 400 160])
        title(['Line ' num2str(i) ' (= ' num2str(floor(dist_alongshore(i))) ' m alongshore distance)'])
        ylim([-10 10])
        xlim([0 1050])
        xlabel('Cross shore distance [m]')
        ylabel('Height [m] NAP')
    end
%%
     figure(562)
        hold off
plot(p2(val_index,1),Z_runup_lim(val_index),'*k');
hold all
set(gca,'XDir','reverse')
xlabel('profile gradient')
ylabel({'Level of B_{ma}'; '[m NAP]'})

% fit line !!
runupVSgradient = polyfit(p2(val_index,1)',Z_runup_lim(val_index),1);

%now plot result in figure
plot(p2(val_index,1)',polyval(runupVSgradient,p2(val_index,1)'),'r','linewidth',2)
R=corrcoef(polyval(runupVSgradient,p2(val_index,1)'),Z_runup_lim(val_index))
ylim([1.4 3.4])
xlim([-0.015 -0.009])

title ('Active profile corr')


        set(gcf,'position',[360 502 200 160])        
        
%%
end