%% Simple example to plot evolution of a profile with forcing

iline=11;  %% choose the transect you want to see



% locate vlugtenburg directory on this computer
curdir=pwd;
ind_v=findstr(curdir, 'vlugtenburg');
pathvlugtenburg=(curdir(1:ind_v+10));  % gives back cometing like F:\check_outs\tudelft\vlugtenburg\

%%  load survey data 
close all
NetCDFfilename=[pathvlugtenburg '\Vlugtenburg.nc'];
%nc_dump(NetCDFfilename);  % give a quick overview of the contents

dates=nc_varget(NetCDFfilename,'time');     % index -1 because netcdf starts counting at zero.. Weird Yanks.
alt=nc_varget(NetCDFfilename,'altitude');  % Format is surveydate,linenr,xloc
surveytype=nc_varget(NetCDFfilename,'survey_type');
x=nc_varget(NetCDFfilename,'cross_shore');

Surveys=list_surveys;
indjet=Surveys.indjet; indwalk=Surveys.indwalk; pairs=Surveys.Survey_pairs;  % rewrite to local variables... bit easier later on


%% load wave data
wavefilename=[pathvlugtenburg '\MeteoHydro\Wavedata\EURHs.mat'];
load (wavefilename)
%% load waterlevel
wlfilename=[pathvlugtenburg '\MeteoHydro\Waterlevel\HvHh.mat'];
load (wlfilename)

% Filter signal waterlevel in a daiy signal
opth=2;  %% choose 1 for mean daily waterlevel, 2 for maximum daily waterlevel
dailyh(:,1)=(h(1,1):1:h(end,1))';  % t vector with timesteps of a day

dt=h(2,1)-h(1,1);
dt_new=dailyh(2,1)-dailyh(1,1);
for i=1:length(dailyh)-1
    date_i=dailyh(i,1);
    indh=find(h(:,1) > date_i  & h(:,1) <= date_i  + dt_new);
    if opth==1
    dailyh(i,2)=sum(h(indh,2))/length(indh);
    end
    if opth==2
    dailyh(i,2)=max(h(indh,2));
    end

end

%% figure all profiles

iline=7;
figure
subplot(211)
hold all
for idate=1:length(pairs)
    
subplot(211)
    plot(x,squeeze(alt(pairs(idate,1),iline,:)))
plot(x,squeeze(alt(pairs(idate,2),iline,:)))

subplot(212)
  plot(x,squeeze(alt(pairs(idate,1),iline,:))*-1)
%pause
end
title(['Profile ' num2str(iline)])
%legend ([dates(indjet,:)])
xlim([0 1300])
xlabel('cross shore distance [m]')





%%
tmin=datenum('01.05.2009', 'dd.mm.yyyy');
tmax=datenum('01.08.2011', 'dd.mm.yyyy');

figure
% subplot top: wave conditions
subplot(6,1,[1:2])
plot(Hs(:,1),Hs(:,2),'r');
hold on
bar(dailyh(:,1),dailyh(:,2))
legend('Hs @ Eur','h @ HvH')
for idate=1:length(dates)
    if surveytype(idate,1) == 'J'
        plot([ datenum(dates(idate,:)) datenum(dates(idate,:))], [0 500],'-.k')
    end
end
xlim([tmin tmax])
ylabel('[cm]')
datetick('x',20, 'keeplimits')

title ('Onshore migration nourishment Vlugtenburg')


% subplot bottom : timestack surveyline
opt2=1;  %% choose between a timestack of the depths (opt2=1) or a timestack of the sed erosion (opt2=2)

subplot(6,1,[3:6])

if opt2==1
pcolor(datenum(dates(indjet,:)),x,squeeze(alt(indjet,iline,:))'); 
hold on;
pcolor(datenum(dates(indwalk,:)),x,squeeze(alt(indwalk,iline,:))');
shading flat
caxis([-8 3]); %colorbar
load('ColormapTimestack.mat')
colormap(cmap);

% Show line top nourishment
%plot(datenum(dates(indjet,:)),x(ind_nour),'w','linewidth',2)
end

if opt2==2
SedEr_j=diff(squeeze(alt(indjet,iline,:)));      
SedEr_w=diff(squeeze(alt(indwalk,iline,:)));    
pcolor(datenum(dates(indjet(2:end),:)),x,SedEr_j'); 
hold on;
pcolor(datenum(dates(indwalk(2:end),:)),x,SedEr_w');
shading flat
caxis([-1 1]);
cmap=colormap;
cmap(28:35,:)=1;
colormap(cmap)
%colorbar
end


% plot lines in the timestack where the surveys happened
for idate=1:length(dates)
    if surveytype(idate,1) == 'J'
        plot([ datenum(dates(idate,:)) datenum(dates(idate,:))], [0 800],'-.k')
    end
end

ylim([100 800])
xlim([tmin tmax])
datetick('x',20, 'keeplimits')
colorbar('location','West')
xlabel('date')
ylabel('crossshore dist')
title(['Timestack profile evolution line ' num2str(iline)])
finame=(['images\Onshoremigrationnourishment' ]);
%print (finame,'-dmeta')

%% figure all profiles
if 0
iline=19;
figure
hold all
for idate=1:length(indjet)
    plot(x,mean(squeeze(alt(indjet(idate),3:19,:)))-3*(idate-1))
end
title(['Mean Profile '])
legend ([dates(indjet,:)])
xlim([0 1300])
xlabel('cross shore distance [m]')

end

%% figure all profiles
iline=7;
figure
hold all
for idate=19:length(indjet)
    plot(x,squeeze(alt(indjet(idate),iline,:))-0*(idate-1))

end
title(['Profile ' num2str(iline)])
legend ([dates(indjet,:)])
xlim([0 1300])
xlabel('cross shore distance [m]')
grid on