function [s,m] = pst(s,m,par,nt,nx,ny,xw,yw,output_dir,test_id);

x = 0:1:230;

h = s.zs-s.zb;
s.Cm = mean(s.ccg,2);
s.Ceqsm = mean(s.ceqsg,2);
s.Ceqbm = mean(s.ceqbg,2);
t = (0:1:nt-1)*par.tint*par.morfac/3600;
Tsim = t(end)*3600;
rhos  = 2650;

Lmix = 2*s.R.*par.Tp./max(0.01,(par.rho.*s.c.*s.H));
Lmix = min(Lmix,h+0.5*s.H);
Lmix2 = sqrt(2*s.R.*par.Tp./max(0.01,(par.rho.*s.c)));
Lmix2 = min(Lmix2,h+0.5*s.H);

% transform Su to Sz
Ssus = zeros(size(s.Susg));
Ssus(1,:) = s.Susg(1,:);
Ssus(2:end,:) = 0.5*(s.Susg(1:end-1,:)+s.Susg(2:end,:));
Sbed = zeros(size(s.Subg));
Sbed(1,:) = s.Subg(1,:);
Sbed(2:end,:) = 0.5*(s.Subg(1:end-1,:)+s.Subg(2:end,:));

dz = s.zb(:,end)-s.zb(:,1);
s.Sdzbtot = (1-par.np)*flipud(cumsum(flipud(dz)))*(xw(2)-xw(1))/Tsim;
s.Sdzavtot = (1-par.np)*flipud(cumsum(flipud(s.dzav(:,end))))*(xw(2)-xw(1))/Tsim;
s.Ssustot = mean(Ssus,2);
s.Sbedtot = mean(Sbed,2);

% sediment transports from undertow:
s.ut  = s.ue-s.u;
s.Slongwaves = s.u.*s.ccg.*h;
s.Sundertow = s.ut.*s.ccg.*h;
s.Swaveassym = s.ua.*s.ccg.*h;

% corrlation long waves and (near bed) turbulence
% for i = 1:nx
%      [R,P] = corrcoef(s.DR(i,:)',detrend(s.u(i,:))');
%      rho(i) = R(1,2); pval(i) = P(1,2);
%      [R,P] = corrcoef(s.ccg(i,:)',detrend(s.u(i,:))');
%      rho2(i) = R(1,2); pval2(i) = P(1,2);
%      [R,P] = corrcoef(s.kb(i,:)',detrend(s.u(i,:))');
%      rho3(i) = R(1,2); pval3(i) = P(1,2);
% end

% what are simulated sed.conc. at times of mobile frame measurements 
MF = 1;
hours = [0 0.1 0.3 1.0 2.04 6.0];
sec = hours*3600;
if sum(isnan(m.tMF))==0
for i = 1:length(m.tMF);    
    seccor(i) = sec(m.interval(i));
end
tMF = m.tMF+seccor;
% now link with simulations...
% windowSize = 90; % ten seconds window 
% fs = specialf('average',windowSize);
% cma = [];
% for i=1:length(xw)
%     cma(i,:) = imfilter(s.ccg(i,:),fs,'symmetric');
% end

cma = s.ccg;
for j = 1:length(nt);
    ind1 = max(1,j-45);
    ind2 = min(nt,j+45);
    cma(:,j) = mean(s.ccg(:,ind1:ind2),2);
end

ind1 = round(tMF/par.morfac/par.tint);
for i = 1:length(m.xswf)
    ind2(i) = find(m.xswf(i)==xw);
end
else
    MF = 0;
end
%% sediment concentrations
figure(11); set(gca,'FontSize',16);
plot(xw,s.Cm*par.rhos,'k-','LineWidth',1.5); hold on; plot(m.xum,m.Cm,'ks','LineWidth',1.5); 
% plot(xw,s.ceqgm*rhos,'k--','LineWidth',1.5);
xlabel('x [m]'); ylabel('c [g/l]'); axis([41 max(max(xw)) 0 45]);
% titletext = ['Simulated test and depth averaged sediment concentration (solid line) compared with the sediment concentrations obtained from suction tubes (squares).'];
% tt = title_text(titletext);
% title (tt,'FontSize',9);
pname = [output_dir,test_id '_fig8.png'];
eval(['print -dpng ' pname]);
pname = [output_dir,test_id '_fig8.eps'];
eval(['print -depsc2 ' pname]);
pname = [output_dir,test_id '_fig8.fig'];
saveas(gcf, pname, 'fig');

figure(12); set(gca,'FontSize',16);
plot(xw,s.Cm*par.rhos,'k-','LineWidth',1.5); hold on; plot(m.xum,m.Cm,'ks','LineWidth',1.5); 
% plot(xw,s.ceqgm*rhos,'k--','LineWidth',1.5);
xlabel('x [m]'); ylabel('c [g/l] & c_{eq} [g/l]'); axis([160 max(max(xw)) 0 45]);
% titletext = ['Simulated test and depth averaged sediment concentration (solid line) compared with the sediment concentrations obtained from suction tubes (squares).'];
% tt = title_text(titletext);
% title (tt,'FontSize',9);
pname = [output_dir,test_id '_fig9.png'];
eval(['print -dpng ' pname]);
pname = [output_dir,test_id '_fig9.eps'];
eval(['print -depsc2 ' pname]);
pname = [output_dir,test_id '_fig9.fig'];
saveas(gcf, pname, 'fig');

xt = [170 175 180 185 190 195 200 205];
for i = 1:length(xt)
    ind(i) = find(xw==xt(i));
end

% fact = nanmean(par.rhos*mean(s.ccg(ind,:),2)./m.Cm')

if MF==1
% scatter plot simulated and measured sediment concentrations
figure(13); set(gca,'FontSize',16);
for i = 1:length(ind1)
    loglog(m.Cmt(i),rhos*cma(ind2(i),ind1(i)),'ks','LineWidth',1.5); hold on;
end
plot([0.1 50],[0.1 50],'k','LineWidth',1.5);
plot([0.1 50],0.5*[0.1 50],'k--','LineWidth',1.5);
plot([0.1 50],2*[0.1 50],'k--','LineWidth',1.5);
xlabel('c_{measured} [g/l]'); ylabel('c_{simulated} [g/l]'); axis image;
axis([0.1 50 0.1 50]); grid on;
% titletext = ['Scatter plot of simulated time and depth averaged sediment concentrations compared with vertically integrated suction tube measurements. The solid line corresponds to a perfect match between measurements and simulations whereas simulation results between the dashed lines are within a factor two with the measurements.'];
% tt = title_text(titletext);
% title (tt,'FontSize',9);
pname = [output_dir,test_id '_fig10.png'];
eval(['print -dpng ' pname]);
pname = [output_dir,test_id '_fig10.eps'];
eval(['print -depsc2 ' pname]);
pname = [output_dir,test_id '_fig10.fig'];
saveas(gcf, pname, 'fig');
end

% wave averaged, bore averaged and near bed turbulence
figure(14); set(gca,'FontSize',16);
plot(xw,mean((s.DR/par.rho).^(2/3),2),'k:','LineWidth',1.5); hold on;
plot(xw,mean((s.DR/par.rho).^(2/3).*par.Tp./s.Tbore,2),'k--','LineWidth',1.5); hold on;
plot(xw,mean(s.kb,2),'k','LineWidth',1.5);
xlabel('x [m]'); ylabel('k [m^2/s^2]'); axis([41 max(max(xw)) 0 3.5]);
% titletext = ['Simulated wave averaged turbulence energy (dotted line), bore averaged turbulence energy (dashed line) and near-bed bore averaged turbulence energy (solid line) as function of cross-shore position.'];
% tt = title_text(titletext);
% title (tt,'FontSize',9);
pname = [output_dir,test_id '_fig11.png'];
eval(['print -dpng ' pname]);
pname = [output_dir,test_id '_fig11.eps'];
eval(['print -depsc2 ' pname]);
pname = [output_dir,test_id '_fig11.fig'];
saveas(gcf, pname, 'fig');

% mixing length
% figure(15); set(gca,'FontSize',16);
% plot(xw,mean(Lmix./h,2),'k','LineWidth',1.5); hold on; plot(xw(:,2),mean(Lmix2./h,2),'k--','LineWidth',1.5); plot(xw(:,2),mean(0.9*s.H./h,2),'k:','LineWidth',1.5);
% plot(xw,mean(min(1,(exp(h./Lmix)-1).^-1),2),'r-','LineWidth',1.5);
% xlabel('x [m]'); ylabel('L_{mix}/h [-] & (exp(h/L_{mix})-1)^-1 [-]'); axis([41 max(max(xw)) 0 1.5]);
% print mixinglength.png -dpng

% surface slope
% figure(16); set(gca,'FontSize',16);
% plot(xw,mean(s.BR,2),'k','LineWidth',1.5); hold on; 
% plot(xw,mean(tan(asin(s.BR)),2),'k--','LineWidth',1.5);
% xlabel('x [m]'); ylabel('\beta & \partial\eta/\partialx [-]'); axis([41 max(xw) 0 1.2*max(mean(tan(asin(s.BR)),2))]);
% print surfaceslope.png -dpng

% correlations:
% 1) roller dissipation and long wave flows
% 2) sediment concentration and long wave flows
% 3) near bed turbulence and long wave flows
% figure(17); set(gca,'FontSize',16);
% plot(xw,rho,'k-','LineWidth',1.5); hold on; plot(xw,rho2,'k--','LineWidth',1.5); hold on; plot(xw,rho3,'k:','LineWidth',1.5); hold on;
% grid on; xlabel('x [m]'); ylabel('\rho_{(D_R,u)} & \rho_{(c,u)} & \rho_{(k_b,u)} [-]');
% axis([41 max(max(xw)) -1 1]);

% sediment transports
figure(18); set(gca,'FontSize',16);
plot(xw,s.Sdzbtot,'k-','LineWidth',2.5); hold on; plot(x,m.Sdzb,'k--','LineWidth',2.5);
plot(xw,s.Sdzavtot,'r-','LineWidth',1.5); plot(xw,s.Ssustot,'r--','LineWidth',1.5); plot(xw,s.Sbedtot,'r-.','LineWidth',1.5); 
xlabel('x [m]'); ylabel('S [m^3/m/s]'); axis([41 max(max(xw)) 1.2*min(min(s.Sdzbtot),min(m.Sdzb)) 2*max(max(s.Sdzbtot),max(m.Sdzb))]);
% titletext = ['Measured (thick dashed line) and simulated (thick solid line) test av-eraged sediment transport from bed level changes. The simulated transport is separated in a transport due to avalanching (dashed-dotted line) and a transport related to the hydrody-namics (dotted line).'];
% tt = title_text(titletext);
% title (tt,'FontSize',9);
pname = [output_dir,test_id '_fig12.png'];
eval(['print -dpng ' pname]);
pname = [output_dir,test_id '_fig12.eps'];
eval(['print -depsc2 ' pname]);
pname = [output_dir,test_id '_fig12.fig'];
saveas(gcf, pname, 'fig');

% decompose flow-based sediment transports
figure(19); set(gca,'FontSize',16);
plot(xw,s.Ssustot,'k-','LineWidth',2.5); hold on;
plot(xw,mean(s.Slongwaves,2),'k--','LineWidth',1.5); plot(xw,mean(s.Sundertow,2),'k-.','LineWidth',1.5); plot(xw,mean(s.Swaveassym,2),'k:','LineWidth',1.5); 
%plot(xw(:,2),mean(Slongwaves+Sundertow+Swaveassym,2),'k--','LineWidth',2.5); 
xlabel('x [m]'); ylabel('S [m^3/m/s]'); axis([41 max(xw) 1.2*min(min(s.Sdzbtot),min(m.Sdzb)) 1.2*max(max(mean(s.Swaveassym,2)),max(m.Sdzb))]);
% titletext = ['Simulated test averaged sediment transport related to the hydrodynamics (solid line) divided into wave asymmetry related sediment transport (dotted line), long wave related sediment transport (dashed line) and sediment transport associated with the short wave undertow (dashed-dotted line).'];
% tt = title_text(titletext);
% title (tt,'FontSize',9);
pname = [output_dir,test_id '_fig13.png'];
eval(['print -dpng ' pname]);
pname = [output_dir,test_id '_fig13.eps'];
eval(['print -depsc2 ' pname]);
pname = [output_dir,test_id '_fig13.fig'];
saveas(gcf, pname, 'fig');


