如何提高图像希尔伯特扫描的性能？

Question

这种基于希尔伯特曲线的图像扫描方法。曲线看起来像（从 1 到 6 阶）：

它可以用于图像扫描。因此，例如，我的三阶曲线代码是：

Hilbert=[C(1,1) C(1,2) C(2,2) C(2,1) C(3,1) C(4,1) C(4,2) C(3,2) C(3,3) C(4,3) C(4,4) C(3,4)...
      C(2,4) C(2,3) C(1,3) C(1,4) C(1,5) C(2,5) C(2,6) C(1,6) C(1,7) C(1,8) C(2,8) C(2,7)...
      C(3,7) C(3,8) C(4,8) C(4,7) C(4,6) C(3,6) C(3,5) C(4,5) C(5,5) C(6,5) C(6,6) C(5,6)...
      C(5,7) C(5,8) C(6,8) C(6,7) C(7,7) C(7,8) C(8,8) C(8,7) C(8,6) C(7,6) C(7,5) C(8,5)...
      C(8,4) C(8,3) C(7,3) C(7,4) C(6,4) C(5,4) C(5,3) C(6,3) C(6,2) C(5,2) C(5,1) C(6,1)...
      C(7,1) C(7,2) C(8,2) C(8,1)];

而且它工作得很快。我为 8 阶和 9 阶曲线制作了相同的函数，但它的工作速度非常非常慢。 9阶，也许，永远不会结束。至少，我没有耐心等待结束 - 2小时后我刚刚关闭了程序。但 7 阶曲线运行 15 秒。怎么了？我可以做同样的事情，但更快吗？是的，程序需要读取 512 * 512 个数组元素，但让它更快也不是不可能。

那么，我到底需要什么 - 我有数组元素的坐标，它们按照应该读取的顺序排列。我需要可接受的时间来读取它们并写入新数组。怎么做？

附言英语对我来说仍然很难，如果有不清楚的地方 - 请问我。

Answer 1

在线快速搜索，您可以在 Steve Eddins 博客上找到 a post about Hilbert curves。这是他生成曲线的实现：

function [x,y] = hilbert_curve(order)
    A = zeros(0,2);
    B = zeros(0,2);
    C = zeros(0,2);
    D = zeros(0,2);

    north = [ 0  1];
    east  = [ 1  0];
    south = [ 0 -1];
    west  = [-1  0];

    for i=1:order
        AA = [B ; north ; A ; east  ; A ; south ; C];
        BB = [A ; east  ; B ; north ; B ; west  ; D];
        CC = [D ; west  ; C ; south ; C ; east  ; A];
        DD = [C ; south ; D ; west  ; D ; north ; B];

        A = AA;
        B = BB;
        C = CC;
        D = DD;
    end

    subs = [0 0; cumsum(A)] + 1;
    x = subs(:,1);
    y = subs(:,2);
end

返回的 xy 坐标是[1,2^order] 范围内的整数。如下所示，该函数足够快：

>> for order=1:10, tic, [x,y] = hilbert_curve(order); toc; end
Elapsed time is 0.001478 seconds.
Elapsed time is 0.000603 seconds.
Elapsed time is 0.000228 seconds.
Elapsed time is 0.000251 seconds.
Elapsed time is 0.000361 seconds.
Elapsed time is 0.000623 seconds.
Elapsed time is 0.001288 seconds.
Elapsed time is 0.007269 seconds.
Elapsed time is 0.029440 seconds.
Elapsed time is 0.117728 seconds.

现在让我们用覆盖了曲线的图像来测试它。我们将图像大小调整为 128x128，以便我们可以看到图案而不会过度拥挤，但您绝对可以为您的情况使用 512x512：

%// some grayscale square image
img = imread('cameraman.tif');

%// scale it down for better visualization
N = 128
assert(N > 0 && mod(N,2)==0);
img = imresize(img, [N N]);

%// space-filling Hilbert curve
order = log2(N)
[x,y] = hilbert_curve(order);

%// show image with curve overlayed
imshow(img, 'InitialMagnification',400)
h = line(x, y);

让我们放大一点，以更好地了解曲线是如何覆盖所有像素的：

>> zoom(10)
>> set(h, 'Marker','.')

最后你可以使用下标索引到图像矩阵中：

>> ind = sub2ind([N N], x, y);
>> pix = img(ind);    %// linear indexing

地点：

>> whos ind
  Name          Size             Bytes  Class     Attributes

  ind       16384x1             131072  double