矩阵中从一个点开始的螺旋循环

Question

我有一个2D网格，想从X、Y开始并保存窗口(W)和重叠(OP)的角落。我尝试了这些代码（链接），但没有一个适合我的目的。

如图所示，我想从一个随机点（黑色单元格）开始，并在螺旋循环中保存每个新窗口的角落位置（由黑色圆圈表示）。该算法可用于任何网格大小（不仅限于正方形）和任何起始点位置。

Matlab也有一个类似于我想要的函数（spiral），但它不接受网格、窗口大小和重叠(OP)。

我期望这个图的输出如下：(8,12) (11,12) (11,9) (8,9) (4,9) (4,12) (4,15) ...

我正在使用以下代码，从一个角落开始，使用定义的W、OP和矩阵大小逐步填充矩阵：

W = [10 12];
OP = [4 3];

M = zeros(100,110);

for i=[1:W(1)-OP(1):size(M,1)-W(1), size(M,1)-W(1)+1]
  for j=[1:W(2)-OP(2):size(M,2)-W(2), size(M,2)-W(2)+1]
      block = rand(W(1),W(2));
      M(i:i+W(1)-1, j:j+W(2)-1) = block;
      imagesc(M); axis equal tight xy
      pause(.1)
  end;
end;

更明确地说，我应该如何更改“上面”的代码，以便从位置（x，y）开始，并根据W、OP和size（M）螺旋填充整个矩阵。

谢谢！

- Sam

---

这个图很不清晰，我无法确定你想要保存哪些位置。 - rayryeng

你能否使用最小的样本输入数据并告诉我们期望的输出结果？此外，重叠区域可能超过一个元素宽，对吗？ - Divakar

我编辑了这个问题。是的，OP可以是多个元素。 - Sam

我在文档中没有看到任何“spiral”函数。你有参考资料吗？ - Ratbert

@Sam：我猜你的例子应该是(8,12) (11,12) (11,9) (8,9) (5,9) (5,12) (5,15)，所以步长始终为3？ - Daniel

Answer 1

基本问题

让数据定义为：

step = 3;  %// step size
x0 = 8;    %// x coordinate of origin
y0 = 12;   %// y coordinate of origin
N = 32;    %// number of steps

那么螺旋的坐标可以通过以下方式作为复平面中的值来获得^†:

z = x0+1j*y0 + step*cumsum([0 -1j.^(-floor(sqrt(4*(0:N)+1))-1)]);

当然，x和y坐标随后是：

x = real(z);
y = imag(z);

使用上述给定的示例值，plot(z,'o-')（或plot(x,y,'o-')）会产生以下图形

† 关键是生成序列1,2,3,3,4,4,5,5,5,6,6,6,7,7,7,7,8,8,8,8...，我感激OEIS 解决了那部分。该序列结果为4n+1的平方根的整数部分，其中n=1,2,3，...。

如何包含重叠和窗口大小

为考虑重叠，请按丹尼尔的建议从步骤中减去其值。

为考虑窗口大小，N应足够大，以便螺旋达到窗口边界之外的某个点；然后仅保留前面的点。

由于事先计算出需要多大的N很困难，因此一种可能的方法是在循环中以指数方式增加N，直到足够大。指数增加确保循环迭代次数较小。下面的代码使用2的幂次方作为N。

%// Data
step = 3;     %// step size
overlap = 1;  %// overlap
x0 = 20;      %// x coordinate of origin
y0 = 15;      %// y coordinate of origin
xmin = 0;     %// window boundary: min x
xmax = 40;    %// window boundary: max x
ymax = 30;    %// window boundary: min y
ymin = 0;     %// window boundary: max y

%// Computations
stepov = step-overlap;
N = 8; %// Initial value. Will be increased as needed
done = false;
while ~done
    z = x0+1j*y0 + stepov*cumsum([0 -1j.^(-floor(sqrt(4*(0:N)+1))-1)]);
        %// compute coordinates of N points
    ind = find(real(z)<xmin | real(z)>xmax | imag(z)<ymin | imag(z)>ymax, 1);
        %// find index of first z out of boundary, if any
    done = ~isempty(ind); %// exit if we have reached outside window boundary
    N = N*2; %// increase number of steps for next try
end
z = z(1:ind-1); %// only keep values that are within the boundary
x = real(z);
y = imag(z);

使用代码中指定的数据，得到的图形如下。请注意，最后一个点是（38,0）。下一个点将是（38，-2），超出了窗口边界。

- Luis Mendo

---

@Benoit_11 谢谢！我喜欢用复数自然地解决问题 :-) - Luis Mendo

@Sam：你必须将 step 设置为 gridsize-op。在生成左下角时，只有到下一个左下角的距离才是相关的。 - Daniel

谢谢Luis。这非常接近我所需要的，除了边界！我该如何考虑它？因为当它到达右侧的末尾时，我想将代码回滚到左侧部分。我的意思是，我想覆盖我的矩阵的所有部分。 - Sam

是的，你说得对，我很抱歉，因为我认为我的问题已经足够清楚了。 - Sam

@Sam 但是“重叠”（OP）是输入还是结果？您不能指定任意重叠并确保窗口被完全填充。 - Luis Mendo

Answer 2

这里有一段代码，它可以产生期望的输出。只需要对spiral_generic进行小的修改就可以满足您的要求：

function demo()
spiral_generic([10,11],[3,4])
W = [10 12];
OP = [4 3];
%make sure your start point is really on the grid of r and c, this is not checked!
start = [19,28];
M = zeros(100,110);
r=[1:W(1)-OP(1):size(M,1)-W(1), size(M,1)-W(1)+1];
c=[1:W(2)-OP(2):size(M,2)-W(2), size(M,2)-W(2)+1];
startindex=[find(r==start(1),1,'first'),find(c==start(2),1,'first')];
A=spiral_generic([numel(r),numel(c)],startindex);
[~,idx]=sort(A(:));
[ridx,cidx]=ind2sub(size(A),idx);
%blocks contains the lower left corners in order of processing.
blocks=[r(ridx);c(cidx)];
for blockindex=blocks
       block = rand(W(1),W(2));
       M(blockindex(1):blockindex(1)+W(1)-1, blockindex(2):blockindex(2)+W(2)-1) = block;
       imagesc(M);
       pause(.1)
end

end
function A = spiral_generic(n, P)
% Makes NxN matrix filled up spirally starting with point P
  r = max([P - 1, n - P]);              % Radius of the bigger matrix
  M = spiral(2 * r + 1);                % Bigger matrix itself
  M = permute(M,[2,1]);                 % changing start direction of the spiral
  M = M(:,end:-1:1);                    % chaning spin orientation
  C = r + 1 - (P - 1);                  % Top-left corner of A in M
  A = M(C(1):C(1)+n(1)-1, C(2):C(2)+n(2)-1);  % Get the submatrix
  [~, order] = sort(A(:));              % Get elements' order
  A(order) = 1:(n(1)*n(2));                     % Fill with continous values
end

- Daniel

---

我能定义起始点吗？（在你的代码中是“center”），因为我的起始点不在矩阵的中心。如果您查看上面的图像，您会发现黑色单元格不在中心。 - Sam

是的，中心可以改变。您没有描述边界到达时它应该如何行为，所以floor(min(center./window_size))*2+1这个术语可能需要更改。 - Daniel

我编辑了问题以精确展示我的需求。 - Sam

我得到了这个（https://www.dropbox.com/s/dh44712m9azctav/spiral_out.png?dl=0），但它与我需要的不同！ - Sam

我使用以下代码进行了更正：for blockindex=1:size(blocks,2) block = rand(W(1),W(2)); M(blocks(1,blockindex):blocks(1,blockindex)+W(1)-1, blocks(2,blockindex):blocks(2,blockindex)+W(2)-1) = block; imagesc(M); axis equal tight xy pause(.1) end - Sam