我试图计算一个由点数组(x,y,z)给出的曲面的曲率。最初我试图拟合一个多项式方程z=a+bx+cx^2+dy+exy+fy^2) 然后计算高斯曲率
$ K = \frac{F_{xx}\cdot F_{yy}-{F_{xy}}^2}{(1+{F_x}^2+{F_y}^2)^2} $
但是,如果曲面是复杂的,则问题是拟合的。我找到了这个Matlab程序来数值计算曲率。我想知道如何在Python中做同样的事情。
function [K,H,Pmax,Pmin] = surfature(X,Y,Z),
% SURFATURE - COMPUTE GAUSSIAN AND MEAN CURVATURES OF A SURFACE
% [K,H] = SURFATURE(X,Y,Z), WHERE X,Y,Z ARE 2D ARRAYS OF POINTS ON THE
% SURFACE. K AND H ARE THE GAUSSIAN AND MEAN CURVATURES, RESPECTIVELY.
% SURFATURE RETURNS 2 ADDITIONAL ARGUEMENTS,
% [K,H,Pmax,Pmin] = SURFATURE(...), WHERE Pmax AND Pmin ARE THE MINIMUM
% AND MAXIMUM CURVATURES AT EACH POINT, RESPECTIVELY.
% First Derivatives
[Xu,Xv] = gradient(X);
[Yu,Yv] = gradient(Y);
[Zu,Zv] = gradient(Z);
% Second Derivatives
[Xuu,Xuv] = gradient(Xu);
[Yuu,Yuv] = gradient(Yu);
[Zuu,Zuv] = gradient(Zu);
[Xuv,Xvv] = gradient(Xv);
[Yuv,Yvv] = gradient(Yv);
[Zuv,Zvv] = gradient(Zv);
% Reshape 2D Arrays into Vectors
Xu = Xu(:); Yu = Yu(:); Zu = Zu(:);
Xv = Xv(:); Yv = Yv(:); Zv = Zv(:);
Xuu = Xuu(:); Yuu = Yuu(:); Zuu = Zuu(:);
Xuv = Xuv(:); Yuv = Yuv(:); Zuv = Zuv(:);
Xvv = Xvv(:); Yvv = Yvv(:); Zvv = Zvv(:);
Xu = [Xu Yu Zu];
Xv = [Xv Yv Zv];
Xuu = [Xuu Yuu Zuu];
Xuv = [Xuv Yuv Zuv];
Xvv = [Xvv Yvv Zvv];
% First fundamental Coeffecients of the surface (E,F,G)
E = dot(Xu,Xu,2);
F = dot(Xu,Xv,2);
G = dot(Xv,Xv,2);
m = cross(Xu,Xv,2);
p = sqrt(dot(m,m,2));
n = m./[p p p];
% Second fundamental Coeffecients of the surface (L,M,N)
L = dot(Xuu,n,2);
M = dot(Xuv,n,2);
N = dot(Xvv,n,2);
[s,t] = size(Z);
% Gaussian Curvature
K = (L.*N - M.^2)./(E.*G - F.^2);
K = reshape(K,s,t);
% Mean Curvature
H = (E.*N + G.*L - 2.*F.*M)./(2*(E.*G - F.^2));
H = reshape(H,s,t);
% Principal Curvatures
Pmax = H + sqrt(H.^2 - K);
Pmin = H - sqrt(H.^2 - K);
我希望我在这里不会太晚。我处理的问题完全相同(我工作的公司的产品)。
首先必须考虑的是,这些点必须表示矩形网格。X是二维数组,Y是二维数组,Z是二维数组。如果您有一个非结构化的浊点,只有一个矩阵形状的Nx3(第一列是X,第二列是Y,第三列是Z),那么您不能应用这个matlab函数。
我已经开发了一个与这个Matlab脚本相当的Python,在这里我只计算Z矩阵的平均曲率(我想你可以从这个脚本中得到灵感,并调整它以获得你想要的所有曲率),忽略X和Y,假设网格是正方形的。我认为你可以“理解”我在做什么和如何做,并根据你的需要调整它:
如果其他人偶然发现这个问题,为了完整性,我提供了以下代码,灵感来自heltonbiker。
下面是一些python代码,用于计算高斯曲率,如“使用几何固有权重计算距离图像的表面曲率”中的方程式(3)所述*,T.Kurita和p.Boulanger,1992。
注:
虽然很晚了,但是没有什么坏处。我修改了“surference”函数以用于Python。 免责声明:我不是原始代码的作者。信用卡到期的地方。
相关问题 更多 >
编程相关推荐