Matlab - 用 surf 和 min 函数减去两个 3D 图形

Question

我正在尝试制作一个看起来像这样的冲浪图：

到目前为止我有：

x = [-1:1/100:1];
y = [-1:1/100:1];

[X,Y] = meshgrid(x,y);

Triangle1 = -abs(X) + 1.5;
Triangle2 = -abs(Y) + 1.5;

Z = min(Triangle1, Triangle2);

surf(X,Y,Z);
shading flat
colormap winter;
hold on;

[X,Y,Z] = sphere();
Sphere = surf(X, Y, Z + 1.5 );% sphere with radius 1 centred at (0,0,1.5)
hold off;

此代码生成的图形如下所示：

• 竖立一个底面为方形 ([-1,1]x[-1,1]) 且顶点位于原点 (0,0) 上方 c = 1.5 处的金字塔。

  李>

• 金字塔的顶部被挖空，方法是移除落在以顶点为中心的半径为 r=1 的球体内的部分。

所以我需要保留金字塔内部的球体表面部分并删除其余部分。请注意，每个图中的 y 轴是不同的，这就是为什么第二个图看起来有点浓缩的原因。是的，有一个金字塔进入球体，从那个角度很难看到。

我将使用 70（方位角）和 35（仰角）的视角。并确保轴正确缩放（如图所示）。在移除球体的适当表面后，我将使用 AXIS TIGHT 选项来获得适当的尺寸。

Answer 1

这是我的谦虚建议：

N = 400; % resolution
x = linspace(-1,1,N);
y = linspace(-1,1,N);
[X,Y] = meshgrid(x,y);
Triangle1 = -abs(X)+1.5 ;
Triangle2 = -abs(Y)+1.5 ;
Z = min(Triangle1, Triangle2);
Trig = alphaShape(X(:),Y(:),Z(:),2);
[Xs,Ys,Zs] = sphere(N-1);
Sphere = alphaShape(Xs(:),Ys(:),Zs(:)+2,2);
% get all the points from the pyramid that are within the sphere:
inSphere = inShape(Sphere,X(:),Y(:),Z(:));
Zt = Z;
Zt(inSphere) = nan; % remove the points in the sphere
surf(X,Y,Zt)
shading interp
view(70,35)
axis tight

我使用alphaShape 对象从金字塔中删除所有不需要的点，然后在没有它们的情况下绘制它：

我知道，它并不完美，因为您看不到金字塔内圆圈的底部，但我为实现这一目标所做的所有尝试都失败了。我的基本想法是像这样将它们绘制在一起：

hold on;
Zc = Zs;
inTrig = inShape(Trig,Xs(:),Ys(:),Zs(:)+1.5);
Zc(~inTrig) = nan;
surf(Xs,Ys,Zc+1.5)
hold off

但结果不是那么好，因为你无法真正看到金字塔内的圆圈。

不管怎样，我把这个贴在这里，因为它可能会给你一个工作方向。

Answer 2

EBH 方法的替代方案。

在 MATLAB 中很难使用 3d 中两个形状相减的通用算法。相反，如果您还记得半径为r 以(x0,y0,z0) 为中心的球体的方程式是

r^2 = (x-x0)^2 + (y-y0)^2 + (z-z0)^2

然后求解z 得到z = z0 +/- sqrt(r^2-(x-x0)^2-(y-y0)^2)，其中在平方根前面使用+ 得到球体的顶部，- 得到底部。在这种情况下，我们只对球体的底部感兴趣。为了得到最终的表面，我们只需在金字塔和半球之间取最小z。

请注意，半球的域由实心圆r^2-(x-x0)^2-(y-y0)^2 >= 0 定义。我们将域外的任何项定义为无穷大，以便在取最小值时忽略它们。

N = 400; % resolution
z0 = 1.5;  % sphere z offset
r = 1;     % sphere radius
x = linspace(-1,1,N);
y = linspace(-1,1,N);
[X,Y] = meshgrid(x,y);
% pyramid
Triangle1 = -abs(X)+1.5 ;
Triangle2 = -abs(Y)+1.5 ;
Pyramid = min(Triangle1, Triangle2);    
% half-sphere (hemisphere)
sqrt_term = r^2 - X.^2 - Y.^2;
HalfSphere = -sqrt(sqrt_term) + z0;
HalfSphere(sqrt_term < 0) = inf;
Z = min(HalfSphere, Pyramid);
surf(X,Y,Z)
shading interp
view(70,35)
axis tight