Skip to content

Commit

Permalink
Merge pull request #20 from hawahee/master
Browse files Browse the repository at this point in the history
Thanks @hawahee 
and sorry for the late merge!
curvature: Fix inaccuracies in the implemented Evans method
Cheers, Wolfgang
  • Loading branch information
wschwanghart committed Oct 4, 2021
2 parents 333a8f7 + 8269eab commit 042db7f
Showing 1 changed file with 11 additions and 5 deletions.
16 changes: 11 additions & 5 deletions @GRIDobj/curvature.m
Original file line number Diff line number Diff line change
Expand Up @@ -124,13 +124,19 @@
correctedges = true;
shape = 'valid';
end

% First-order partial derivatives:
[fx,fy] = gradient(dem,cs);

if correctedges
dem = padarray(dem,[1 1],'symmetric');
end

% First-order partial derivatives:
% kernel for dz/dx
kernel = [-1 0 1; -1 0 1; -1 0 1]./(6*cs);
fx = conv2(dem,kernel,shape);
% kernel for dz/dy
kernel = [1 1 1; 0 0 0; -1 -1 -1]./(6*cs);
fy = conv2(dem,kernel,shape);

% Second order derivatives according to Evans method (see Olaya 2009)
%
% z1 z2 z3
Expand Down Expand Up @@ -159,9 +165,9 @@
case 'profc'
curv = - (fx.^2 .* fxx + 2*fx.*fy.*fxy + fy.^2.*fyy)./((fx.^2 + fy.^2).*(1 + fx.^2 + fy.^2).^(3/2));
case 'tangc'
curv = - (fy.^2 .* fxx + 2*fx.*fy.*fxy + fx.^2.*fyy)./((fx.^2 + fy.^2).*(1 + fx.^2 + fy.^2).^(1/2));
curv = - (fy.^2 .* fxx - 2*fx.*fy.*fxy + fx.^2.*fyy)./((fx.^2 + fy.^2).*(1 + fx.^2 + fy.^2).^(1/2));
case 'planc'
curv = - (fy.^2 .* fxx + 2*fx.*fy.*fxy + fx.^2.*fyy)./((fx.^2 + fy.^2).^(3/2));
curv = - (fy.^2 .* fxx - 2*fx.*fy.*fxy + fx.^2.*fyy)./((fx.^2 + fy.^2).^(3/2));
case 'meanc'
curv = - ((1+fy.^2).*fxx - 2.*fxy.*fx.*fy + (1+fx.^2).*fyy)./ ...
(2.* (fx.^2+fy.^2+1).^(3/2));
Expand Down

0 comments on commit 042db7f

Please sign in to comment.