0

我有两个矩阵。每个矩阵的维度为 3 * k，表示 k 个圆，每列以 [xyr] 的形式表示，其中 (x,y) 是圆的中心，r 是它的半径。所以一个有 3 个圆的矩阵是 [x1 x2 x3; y1 y2 y3; r1 r2 r3]

我需要找到两个矩阵之间重叠的圆圈和重叠区域。假设，考虑第一个矩阵中的一个圆。现在考虑第二个矩阵中的每个其他圆圈。现在，我需要找出第二个矩阵中的哪些圆圈与所考虑的圆圈重叠。我需要为第一个矩阵中的每个圆圈执行此操作。类似地，对于第二个矩阵中的每个圆相对于第一个矩阵。

所以我需要两个矩阵中的每个圆（其中有 k1+k2 个），与另一个矩阵重叠的面积有多少，另一个矩阵中重叠的圆是什么。

显然，可能有多个重叠的圆圈。两个矩阵的维度 k 可能不同。我有每个矩阵的 2 个额外矩阵，一个按 x 坐标排序，另一个按 y 坐标排序，如果这有助于计算。问题是矩阵中有很多圆圈，我正在寻找一种有效的方法来做。此外，我想将此扩展到两个以上的矩阵，然后一个有效的算法将大大提高执行时间。

此链接中给出了这两个图像的预览（我拥有相应的圆圈矩阵）： https ://www.dropbox.com/s/om5has5uw91dj6p/overlap.jpg

4

2 回答 2

2

你需要两件事：

  1. 找出哪些圆对重叠（距离 < 半径总和）
  2. 计算这些对的重叠区域

已编辑我修改了代码以添加一些矢量化并捕获一个圆圈完全在另一个圆圈内的情况——Wolfram 公式被打破的情况。与 1000 x 2000 圈相比，展示了速度。

这是一些执行这些操作的代码（带有简单示例）。它似乎相当有效。用 1.4 秒来比较 1000 和 2000 圈（166 万重叠）。

% circles [x; y; r]
m1 = rand(3, 1000);
m2 = rand(3, 2000);

tic
% distances between circles:
dx = bsxfun(@minus, m1(1,:), m2(1,:)');
dy = bsxfun(@minus, m1(2,:), m2(2,:)');
sumr = bsxfun(@plus, m1(3,:), m2(3,:)');

% distance between centers:
dist = sqrt(dx.^2 + dy.^2);

% pairs that have overlap:
pairs = find(dist < sumr);

s1 = size(m1,2);
s2 = size(m2,2);

% calculate overlaps
c1 = repmat(1:s1, [s2 1]); % circle from m1 corresponding to value in pairs()
c2 = repmat(1:s2, [s1 1])'; % ditto for m2

o = overlap(m1(3, c1(pairs)), m2(3, c2(pairs)), dist(pairs)');
toc

将以下内容放在路径中自己的文件overlap.m中：这将计算重叠

function o = overlap(r1, r2, d)
r12 = r1.^2;
r22 = r2.^2;
d2 = d.^2;
a1 = acos((d2 + r12 - r22)./(2.*d.*r1));
a2 = acos((d2 + r22 - r12)./(2.*d.*r2));
o = r12.*a1 + r22.*a2-.5*sqrt((-d+r1+r2)*(d+r1-r2)*(d-r1+r2)*(d+r1+r2));

% fix situation where r1 - r2 > d:
% i.e. one circle fully inside the other
f = find(abs(imag(a1)) + abs(imag(a2)) > 0);
o(f) = pi * min(r1(f),r2(f)).^2;
于 2013-09-12T15:46:50.477 回答
0

我想过使用光栅化来近似圆之间的重叠。这个想法是创建N点来近似圆形路径，然后将它们传递给poly2mask函数。这将返回一个在定义圆的位置填充的二进制掩码图像（该算法应以亚像素精度工作）。对所有圆对执行此操作，然后我们可以通过对两个掩码应用逻辑 AND 来计算交集，并计算结果中剩余的“on”像素的数量（或者更好地bwarea用于稍微更好的估计）。同样，这只是一个近似值，因为我们使用的是离散的像素网格......

从好的方面来说，我们不必担心数学方程。事实上，这种方法适用于任何多边形形状，在这些多边形形状中，很难为交叉点提出封闭形式的解决方案。

不幸的是，这种蛮力方法被证明比我预期的要慢（没有@Floris的代码那么快），无论如何我都会发布我的尝试。我从另一个代码中借用了检查两个圆圈是否感兴趣的想法，如果您拥有的大多数圆圈相距甚远，这应该会大大减少时间。

只是为了好玩，我添加了一些代码来可视化这个过程！

ANIMATION = true;

% size of image we are working in
%sz = [180 360];
sz = [100 100];

% generated two sets of random circles (specified as columns [x;y;r])
% (I am just trying to create circles wholly visible within the box)
k1 = 100; k2 = 80;
M1 = bsxfun(@times, rand(3,k1), [sz(:)*0.6+sz(:)*0.2;0.1*min(sz)+0.1*min(sz)]);
M2 = bsxfun(@times, rand(3,k2), [sz(:)*0.6+sz(:)*0.2;0.1*min(sz)+0.1*min(sz)]);

% animation
if ANIMATION
    clf
    hImg = imshow(zeros(sz), 'InitialMag','fit');
    hLine(1) = line(NaN, NaN, 'Color','r', 'LineWidth',2);
    hLine(2) = line(NaN, NaN, 'Color','b', 'LineWidth',2);
    axis on
end

% used to approximate circles
num = 50;
t = linspace(0, 2*pi, num);
ct = cos(t); st = sin(t);

% test to find which circles intersect
dist = pdist2(M1(1:2,:)', M2(1:2,:)');
sumr = bsxfun(@plus, M1(3,:)', M2(3,:));
skipIdx = (sumr < dist);

% compute overlap between all pairs
overlap = zeros(k1,k2);
for i=1:k1
    for j=1:k2
        % early skip if circles dont interset
        if skipIdx(i,j), continue, end

        % compute each circle points
        x1 = M1(3,i)*ct + M1(1,i); y1 = M1(3,i)*st + M1(2,i);
        x2 = M2(3,j)*ct + M2(1,j); y2 = M2(3,j)*st + M2(2,j);

        % rasterize circles
        BW1 = poly2mask(x1, y1, sz(1), sz(2));
        BW2 = poly2mask(x2, y2, sz(1), sz(2));

        % compute area of intersection of the two masks
        %overlap(i,j) = sum(BW1(:)&BW2(:));
        overlap(i,j) = bwarea(BW1&BW2);

        if ANIMATION
            % update animation
            set(hImg, 'CData',BW1&BW2)
            set(hLine(1), 'XData',x1, 'YData',y1)
            set(hLine(2), 'XData',x2, 'YData',y2)
            title(sprintf('%g',overlap(i,j)))
            drawnow
        end
    end
end

为了了解近似值有多准确，我捕获了一个迭代，其中一个圆圈完全在另一个圆圈内。这样我们就可以将实际的相交区域与我们的近似值进行比较。

对于下面显示的图像，我们有：

K>> M1(:,i)
ans =
   58.2106    % x
   24.0996    % y
    8.9387    % r

K>> pi*M1(3,i)^2
ans =
  251.0122

K>> overlap(i,j)
ans =
  252.5000

我认为这已经足够接近了 :) 使用更高分辨率的像素网格会提高近似值，但肯定会进一步减慢它的速度。

路口

于 2013-09-13T22:06:58.373 回答