Skip to content
Snippets Groups Projects
Commit 8b7c1698 authored by mey's avatar mey
Browse files

new erosanimation function for output analysis

parent 6896d6b2
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
Showing
with 235 additions and 107 deletions
mfiles/erosanimation.m
+ 235
− 107
View file @ 8b7c1698
function [B,varargout] = erosanimation(variable,varargin) function [F] = erosanimation(variable,varargin)
% Visualize output of the EROS landscape evolution model (LEM) % Visualize output of the EROS landscape evolution model as animations (LEM)
% %
% %
% The following function library is required, which can be downloaded % The following function library is required, which can be downloaded
...@@ -17,13 +17,14 @@ function [B,varargout] = erosanimation(variable,varargin) ...@@ -17,13 +17,14 @@ function [B,varargout] = erosanimation(variable,varargin)
% %
% DESCRIPTION % DESCRIPTION
% %
% erosanimation shows timeseries data either as grid average or as 2d/3d movie % erosanimation creates movie frames of landscape evolution either in
% profile view, map view or in 3d.
% %
% %
% INPUT (required) % INPUT (required)
% %
% variable variable of interest (string) % variable variable of interest (string)
% 'topo' Topographic elevation % 'topo' Topographic elevation
% 'sediment' Sediment thickness % 'sediment' Sediment thickness
% 'water' Water depth % 'water' Water depth
% 'discharge' Water discharge % 'discharge' Water discharge
...@@ -36,13 +37,15 @@ function [B,varargout] = erosanimation(variable,varargin) ...@@ -36,13 +37,15 @@ function [B,varargout] = erosanimation(variable,varargin)
% 'hum' Water discharge on the topography % 'hum' Water discharge on the topography
% 'rain' Sources (>0) and sinks (-1) of water and sediment % 'rain' Sources (>0) and sinks (-1) of water and sediment
% %
% 'profile' custom profile, second argument needs to be a DEM (GRIDobj)
% 'sprofile' stream long profile, second argument needs to be a DEM (GRIDobj)
%
% INPUT (optional) % INPUT (optional)
% %
% Parameter name/value pairs (pn,pv,...) % Parameter name/value pairs (pn,pv,...)
% %
% 'mode' visualization mode (string) (default: 'movie2') % 'mode' visualization mode (string) (default: 'movie2')
% 'average' shows the evolution of the spatial average %
% of the variable defined as required input
% 'movie2' 2d movie of variable % 'movie2' 2d movie of variable
% 'movie3' 3d movie of topographic evolution % 'movie3' 3d movie of topographic evolution
% %
...@@ -54,13 +57,7 @@ function [B,varargout] = erosanimation(variable,varargin) ...@@ -54,13 +57,7 @@ function [B,varargout] = erosanimation(variable,varargin)
% %
% OUTPUT % OUTPUT
% %
% B movie frames captured with modes 'movie2' and 'movie3' % F movie frames
% B (mode='average') 3d array of the variable of interest.
%
%
% OUTPUT (optional)
%
% meanB spatial average of variable through time
% %
% EXAMPLE % EXAMPLE
% %
...@@ -98,21 +95,28 @@ function [B,varargout] = erosanimation(variable,varargin) ...@@ -98,21 +95,28 @@ function [B,varargout] = erosanimation(variable,varargin)
% Date: 28. May, 2020 % Date: 28. May, 2020
p = inputParser; p = inputParser;
expectedInput_variable = {'topo','water','sediment','qs','flux',... expectedInput_variable = {'topo','water','sediment','flux','qs',...
'discharge','downward','stress','hum','slope','capacity','stock'}; 'discharge','downward','stress','hum','slope','capacity','stock','sprofile','profile'};
addRequired(p,'variable',@(x) any(validatestring(x,expectedInput_variable))); addRequired(p,'variable',@(x) any(validatestring(x,expectedInput_variable)));
if strcmp(variable,'sprofile') || strcmp(variable,'profile')
addRequired(p,'dem',@(x)isa(x,'GRIDobj'));
default_flowmin = 100;
addParameter(p,'flowmin',default_flowmin,@isnumeric);
parse(p,variable,varargin{:});
dem = p.Results.dem;
flowmin = p.Results.flowmin;
else
default_mode = 'movie2';
expectedInput_mode = {'movie2','movie3'};
addParameter(p,'mode',default_mode,@(x) any(validatestring(x,expectedInput_mode)));
default_viewdir = [45,45];
addParameter(p,'viewdir',default_viewdir,@isnumeric);
parse(p,variable,varargin{:});
mode = p.Results.mode;
viewdir = p.Results.viewdir;
end
default_mode = 'movie2';
expectedInput_mode = {'average','movie2','movie3'};
addParameter(p,'mode',default_mode,@(x) any(validatestring(x,expectedInput_mode)));
default_viewdir = [45,45];
addParameter(p,'viewdir',default_viewdir,@isnumeric);
parse(p,variable,varargin{:});
mode = p.Results.mode;
viewdir = p.Results.viewdir;
switch variable switch variable
case 'topo' case 'topo'
...@@ -165,86 +169,210 @@ switch variable ...@@ -165,86 +169,210 @@ switch variable
colors = 'jet'; colors = 'jet';
end end
% determine timesteps if strcmp(variable,'sprofile')
FD = FLOWobj(dem,'preprocess','c');
T = dir('*.ini'); S = STREAMobj(FD,flowacc(FD)>flowmin);
Z = dir(['*.',filetype]); S2 = modify(S,'interactive','reachselect');
[t,~] = fread_timeVec(T.name,length(Z)); marker = round(linspace(1,length(S2.distance),10));
if isempty(t) H = dir('*.alt');
t=1:length(Z); W = dir('*.water');
end S = dir('*.sed');
if isnan(t) D = dir('*.flux');
t=1:length(Z);
end [~,index] = sortrows({H.date}.');
H = H(index);
[~,index] = sortrows({Z.date}.'); W = W(index);
Z = Z(index); S = S(index);
for i = 1:length(Z) D = D(index);
[z,~] = fopengrd(Z(i).name);
B(:,:,i) = z; for i = 1:length(D)
meanB(i)=mean(z(:)); [z,~] = fopengrd(D(i).name);
end z(z==0)=NaN;
switch mode B(:,:,i) = z;
case 'average' end
plot(t,meanB)
ylabel(iylabel); sed = grd2GRIDobj(S(1).name,dem);
xlabel('time'); w = waitbar(1/length(H),['Collecting movie frames ... ']);
grid on for i=1:length(H)
B=meanB; h=figure;
case 'movie2' subplot(2,1,1)
H = dir('*.alt'); water2 = grd2GRIDobj(W(i).name,dem);
Z = dir(['*.',filetype]); dem2 = grd2GRIDobj(H(i).name,dem);
[~,index] = sortrows({H.date}.'); sed2 = grd2GRIDobj(S(i).name,dem);
H = H(index); plotdz(S2,dem,'color',[0.9 0.9 0.9]);hold on
Z = Z(index); plotdz(S2,dem-sed,'color',[0.7 0.7 0.7])
w = waitbar(1/length(H),['Collecting movie frames ... ']); plotdz(S2,dem2+water2,'color',[0 0.61 1])
for i = 1:length(H)-1
h = grd2GRIDobj(H(i+1).name); plotdz(S2,dem2,'color',[1 0.5 0.1]);hold on
z = grd2GRIDobj(Z(i+1).name); plotdz(S2,dem2-sed2,'color','k')
z.Z(z.Z==0)=NaN; xlim([0 S2.distance(end)])
% imageschs(h,z,'colormap',colors,'caxis',[min(B(:)),200],'colorbarylabel',iylabel); yl(i,1:2) = ylim;
imageschs(h,z,'colormap',colors,'caxis',[min(B(:)),max(B(:))],'colorbarylabel',iylabel); ylim([yl(1,1),yl(1,2)+20]);
title(['Time = ',num2str(t(i)),'']) title(['Time = ',num2str(i),''])
x0=10; xlabel('Distance (m)')
y0=10; ylabel('Elevation (m)');
width=2200; legend({'Initial topo','Initial bedrock','Water','Sediment','Bedrock'},'Location','northwest')
height=1900; legend('boxoff')
set(gcf,'position',[x0,y0,width,height])
set(gcf,'Visible','off')
F(i) = getframe(gcf); subplot(2,1,2)
close all flux = grd2GRIDobj(D(i).name,dem);
waitbar(i/length(H)) flux.Z(flux.Z==0)=NaN;
end imageschs(dem2,flux,'colormap','flowcolor','caxis',[nanmin(B(:)),nanmax(B(:))],'colorbarylabel','Water discharge (m^3/s)');
close(w) hold on;
B = F; plot(S2,'k--')
case 'movie3' % scatter(S2.x(marker),S2.y(marker),'k')
H = dir('*.alt'); text(S2.x(marker),S2.y(marker),num2str(round(S2.distance(marker)/1000)),'FontWeight','bold')
% Z = dir(['*.',filetype]); x0=10;
[~,index] = sortrows({H.date}.'); y0=10;
H = H(index); width=2200;
% Z = Z(index); height=1900/2;
w = waitbar(1/length(H),['Collecting movie frames ... ']); set(gcf,'position',[x0,y0,width,height])
for i = 1:length(H) set(gcf,'Visible','off')
h = grd2GRIDobj(H(i).name); F(i) = getframe(gcf);
% z = grd2GRIDobj(Z(i).name); close all
[xm,ym] = getcoordinates(h); waitbar(i/length(H))
axis off end
surface(xm,ym,h.Z,'EdgeColor','none');colorbar close(w)
view(viewdir(1),viewdir(2)) elseif strcmp(variable,'profile')
axis equal H = dir('*.alt');
c = colorbar; W = dir('*.water');
c.Label.String = 'Elevation (m)'; S = dir('*.sed');
colormap(landcolor) D = dir('*.flux');
% caxis([nanmin(B(:)),nanmax(B(:))])
title(['Time = ',num2str(t(i)),'']) [~,index] = sortrows({H.date}.');
set(gcf,'Visible','off') H = H(index);
F(i) = getframe(gcf); W = W(index);
close all S = S(index);
waitbar(i/length(H)) D = D(index);
end
close(w) for i = 1:length(D)
f = figure; [z,~] = fopengrd(D(i).name);
movie(f,F,1,10) z(z==0)=NaN;
close(f) B(:,:,i) = z;
B = F; end
sed = grd2GRIDobj(S(1).name,dem);
[d,z,x,y] = demprofile(dem);
[~,sz0] = demprofile(sed,[],x,y);
w = waitbar(1/length(H),['Collecting movie frames ... ']);
for i=1:length(H)
h=figure;
subplot(2,1,1)
water2 = grd2GRIDobj(W(i).name,dem);
dem2 = grd2GRIDobj(H(i).name,dem);
sed2 = grd2GRIDobj(S(i).name,dem);
[~,hz] = demprofile(dem2,[],x,y);
[~,wz] = demprofile(water2,[],x,y);
[~,sz] = demprofile(sed2,[],x,y);
plot(d,z,'color',[0.9 0.9 0.9]);hold on
plot(d,z-sz0,'color',[0.7 0.7 0.7])
plot(d,hz+wz,'color',[0 0.61 1])
plot(d,hz,'color',[1 0.5 0.1]);hold on
plot(d,hz-sz,'color','k')
xlim([0 d(end)])
yl(i,1:2) = ylim;
ylim([yl(1,1),yl(1,2)+20]);
title(['Time = ',num2str(i),''])
xlabel('Distance (m)')
ylabel('Elevation (m)');
subplot(2,1,2)
flux = grd2GRIDobj(D(i).name,dem);
flux.Z(flux.Z==0)=NaN;
imageschs(dem2,flux,'colormap','flowcolor','caxis',[nanmin(B(:)),nanmax(B(:))],'colorbarylabel','Water discharge (m^3/s)');
hold on;
plot(x,y,'k--')
text(x(1)+10,y(1)+10,'left');
text(x(end)+10,y(end)+10,'right');
x0=10;
y0=10;
width=2200;
height=1900/2;
set(gcf,'position',[x0,y0,width,height])
set(gcf,'Visible','off')
F(i) = getframe(gcf);
close all
waitbar(i/length(H))
end
close(w)
else
T = dir('*.ini');
Z = dir(['*.',filetype]);
[t,~] = fread_timeVec(T.name,length(Z));
if isempty(t)
t=1:length(Z);
end
if isnan(t)
t=1:length(Z);
end
[~,index] = sortrows({Z.date}.');
Z = Z(index);
for i = 1:length(Z)
[z,~] = fopengrd(Z(i).name);
z(z==0)=NaN;
B(:,:,i) = z;
end
switch mode
case 'movie2'
H = dir('*.alt');
Z = dir(['*.',filetype]);
[~,index] = sortrows({H.date}.');
H = H(index);
Z = Z(index);
w = waitbar(1/length(H),['Collecting movie frames ... ']);
for i = 1:length(H)-1
h = grd2GRIDobj(H(i+1).name);
z = grd2GRIDobj(Z(i+1).name);
z.Z(z.Z==0)=NaN;
%imageschs(h,z,'colormap',colors,'caxis',[0,1],'colorbarylabel',iylabel);
imageschs(h,z,'colormap',colors,'caxis',[nanmin(B(:)),nanmax(B(:))],'colorbarylabel',iylabel);
%set(gca,'ColorScale','log')
title(['Time = ',num2str(t(i)),''])
x0=10;
y0=10;
width=2200;
height=1900;
set(gcf,'position',[x0,y0,width,height])
set(gcf,'Visible','off')
F(i) = getframe(gcf);
close all
waitbar(i/length(H))
end
close(w)
case 'movie3'
H = dir('*.alt');
[~,index] = sortrows({H.date}.');
H = H(index);
w = waitbar(1/length(H),['Collecting movie frames ... ']);
for i = 1:length(H)
h = grd2GRIDobj(H(i).name);
[xm,ym] = getcoordinates(h);
axis off
surface(xm,ym,h.Z,'EdgeColor','none');colorbar
view(viewdir(1),viewdir(2))
axis equal
c = colorbar;
c.Label.String = 'Elevation (m)';
colormap(landcolor)
caxis([nanmin(B(:)),nanmax(B(:))])
title(['Time = ',num2str(t(i)),''])
set(gcf,'Visible','off')
F(i) = getframe(gcf);
close all
waitbar(i/length(H))
end
close(w)
end
end end
\ No newline at end of file
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment