如何在MATLAB中可视化球体的相交? [英] How do I visualize the intersection of spheres in MATLAB?

问题描述

似乎在几个地方都提出了这个问题(包括在SO上).在可视化三边测量问题的结果时，我最近遇到了对此的需求.

It seems this question has been asked in a few places (including on SO). I recently came across the need for this when visualizing results of a trilateration problem.

几乎在每种情况下，答案都会引导查询着眼于 Wolfram进行数学学习，但不包括任何代码.数学确实是一个很好的参考，但是如果我问一个关于编程的问题，那么一些代码也可能会有所帮助. (对代码问题的回答避免诸如编写代码很简单"之类的尖刻注释时，当然也很感激.)

In almost every case, the answer directs the inquiry to look at Wolfram for the math but excludes any code. The math really is a great reference, but if I'm asking a question on programming, some code might help as well. (It certainly is also appreciated when answers to a code question avoid pithy comments like "writing the code is trivial").

那么如何才能在MATLAB中可视化球体的相交呢?我在下面有一个简单的解决方案.

So how can one visualize the intersection of spheres in MATLAB? I have a simple solution below.

推荐答案

我写了一个小脚本来做到这一点.随时提出建议和编辑.通过检查每个球体的表面是否落在所有其他球体的体积之内来进行工作.

I wrote a small script to do just this. Feel free to make suggestions and edits. It works by checking if the surface of each sphere falls within the volume of all of the other spheres.

对于球面相交，在sphere()函数调用中使用更多的面更好(但更慢).这应该在可视化中给出更密集的结果.对于仅球形的可视化，较小的数目(约50)就足够了.请参阅注释，以了解如何可视化每个注释.

For sphere intersection, it's better (but slower) to use a larger number of faces in the sphere() function call. This should give denser results in the visualization. For the sphere-alone visualization, a smaller number (~50) should suffice. See the comments for how to visualize each.

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秒在笔记本电脑上运行组合的可视化效果.

For the 6 spheres used in these examples, it took an average of 1.934 seconds to run the combined visualization on my laptop.

6个球体的交点:

Intersection of 6 spheres:

实际6个领域:

下面，我将两者结合在一起，以便您可以在球体视图中看到相交处.

Below, I've combined the two so you can see the intersection in the view of the spheres.

对于这些示例:

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

我希望这对希望将这种效果可视化的其他人有所帮助.

I hope this helps anyone else who may desire to visualize this effect.

这篇关于如何在MATLAB中可视化球体的相交?的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！