如何在MATLAB中可视化多个球体的交集？

Question

如何在MATLAB中可视化多个球体的交集？

3

似乎这个问题已经在几个地方（包括SO）被问过了。当我需要可视化三边测量问题的结果时，我最近遇到了这个需求。

几乎每种情况下，答案都会引导查询者去看Wolfram的数学内容，但不包括任何代码。这些数学确实是一个很好的参考，但如果我在编程方面有问题，一些代码也可能有所帮助。（当然，如果回答一个代码问题时要避免简短的评论，比如“编写代码很容易”也会受到赞赏）。

那么如何在MATLAB中可视化球体的交点呢？我以下提供了一个简单的解决方案。

- marcman

1个回答

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- marcman · Accepted Answer

我写了一个小脚本来实现这个功能。欢迎提出建议和修改意见。它的工作原理是检查每个球体表面是否落在所有其他球体的体积之内。

对于球体交集，最好（但速度较慢）在sphere()函数调用中使用更多的面数。这应该会在可视化中给出更密集的结果。对于独立的球体可视化，较小的数字（~50）就足够了。请参阅注释以了解如何可视化每个内容。

close all
clear
clc

% centers   : 3 x N matrix of [X;Y;Z] coordinates
% dist      : 1 x N vector of sphere radii

%% Plot spheres (fewer faces)
figure, hold on % One figure to rule them all
[x,y,z] = sphere(50); % 50x50-face sphere
for i = 1 : size(centers,2)
    h = surfl(dist(i) * x + centers(1,i), dist(i) * y + centers(2,i), dist(i) * z + centers(3,i));
    set(h, 'FaceAlpha', 0.15)
    shading interp
end

%% Plot intersection (more faces)
% Create a 1000x1000-face sphere (bigger number = better visualization)
[x,y,z] = sphere(1000);

% Allocate space
xt = zeros([size(x), size(centers,2)]);
yt = zeros([size(y), size(centers,2)]);
zt = zeros([size(z), size(centers,2)]);
xm = zeros([size(x), size(centers,2), size(centers,2)]);
ym = zeros([size(y), size(centers,2), size(centers,2)]);
zm = zeros([size(z), size(centers,2), size(centers,2)]);

% Calculate each sphere
for i = 1 : size(centers, 2)
    xt(:,:,i) = dist(i) * x + centers(1,i);
    yt(:,:,i) = dist(i) * y + centers(2,i);
    zt(:,:,i) = dist(i) * z + centers(3,i);
end

% Determine whether the points of each sphere fall within another sphere
% Returns booleans
for i = 1 : size(centers, 2)
    [xm(:,:,:,i), ym(:,:,:,i), zm(:,:,:,i)] = insphere(xt, yt, zt, centers(1,i), centers(2,i), centers(3,i), dist(i)+0.001);
end

% Exclude values of x,y,z that don't fall in every sphere
xmsum = sum(xm,4);
ymsum = sum(ym,4);
zmsum = sum(zm,4);
xt(xmsum < size(centers,2)) = 0;
yt(ymsum < size(centers,2)) = 0;
zt(zmsum < size(centers,2)) = 0;

% Plot intersection
for i = 1 : size(centers,2)
    xp = xt(:,:,i);
    yp = yt(:,:,i);
    zp = zt(:,:,i);
    zp(~(xp & yp & zp)) = NaN;
    surf(xt(:,:,i), yt(:,:,i), zp, 'EdgeColor', 'none');
end

以下是insphere函数的代码：

function [x_new,y_new,z_new] = insphere(x,y,z, x0, y0, z0, r)
    x_new = (x - x0).^2 + (y - y0).^2 + (z - z0).^2 <= r^2;
    y_new = (x - x0).^2 + (y - y0).^2 + (z - z0).^2 <= r^2;
    z_new = (x - x0).^2 + (y - y0).^2 + (z - z0).^2 <= r^2;
end

可视化示例

这些示例中使用的 6 个球体，平均需要 1.934 秒在我的笔记本电脑上运行组合可视化。

6 个球体的交集：

实际的 6 个球体：

下面，我将两者结合起来，让您可以看到球体的交集。

对于这些示例:

centers =

   -0.0065   -0.3383   -0.1738   -0.2513   -0.2268   -0.3115
    1.6521   -5.7721   -1.7783   -3.5578   -2.9894   -5.1412
    1.2947   -0.2749    0.6781    0.2438    0.4235   -0.1483

dist =

    5.8871    2.5280    2.7109    1.6833    1.9164    2.1231

我希望这可以帮助其他想要实现这种效果的人。