GNU Octave：霍夫变换

Question

我正在尝试使用 hough 变换，不幸的是它似乎没有输出与 Drawn 线对应的 r 和 theta。我一直试图在这个网站和其他网站上找到答案，但到目前为止我尝试的一切都失败了。

I=zeros(80, 80);
for n=1:25;
  I(n+20, n+2)=1;
  I(n+1, n*2+17)=1;
  I(n+1, n*2+16)=1;
  end

hough = houghtf(I,"line", pi*[0:360]/180);
threshHough = hough>.9*max(hough(:));
[r, theta] = find(threshHough>0)
%theta = (theta-91)*pi/180
%r=r-size(hough,1)/2

imshow(I)

Answer 1

Octave中的houghtf函数将一行参数化为

r = x*cos(theta) + y*sin(theta)

输出是一个 NxM 矩阵，其中第一维表示 r，第二维表示 theta。

documentation关于r参数不是很清楚。它只表示 N = 2*diag_length - 1，图像的对角线长度为 diag_length。但它甚至没有告诉你原点在哪里。经过一些实验，我发现原点在 (N-1)/2。因此，在指数 (ii,jj),

找到峰值后

r = ii - (size(hough,1)-1)/2;

文档中对 theta 参数的描述要好得多：返回的 bin 对应于输入中给定的值。所以 theta 是您传递给 houghft 的数组 pi*[0:360]/180 的索引。你可以这样写：

angles = pi*[0:360]/180;
theta = angles(jj);

请注意霍夫变换处理的是线的方向，而不是线的方向。您可以使用 theta 或 theta+180 参数化同一行。因此，上面代码生成的 hough 图像是多余的，所有方向都重复了两次（实际上，0 度方向的结果在 180 度和 360 度重复了两次）。相反，使用

angles = pi*[0:179]/180;

接下来，简单地对 Hough 变换进行阈值处理并获取所有像素并不是最好的方法：一些峰值会更强，导致领结形检测，您将只想使用其中最高的观点。 bwmorph 函数的 'shrink' 方法是一种将阈值后的检测减少到单个点的快捷方式，但不能保证在最高点。另请注意，对于 0 度的线条，领结在图像边缘被切成两半，另一半以 179 度角向后出现。这需要额外的努力来修复。下面，我通过扩展 angles 一点，然后删除重复的点来修复它。

最后，当对线进行参数化时，似乎 Octave 选择 x 作为图像的第一个索引，y 作为第二个索引，这与 imshow 相反

将它们放在一起，我们得到：

% Create test image with some lines
I = zeros(80, 80);
for n=1:25;
  I(n+20,n+2) = 1;
  I(n+1,n*2+17) = 1;
  I(n+1,n*2+16) = 1;
end
I(34:73,55) = 1;
I(60,15:64) = 1;

% Compute the Hough transform
angles = pi*[-10:189]/180;
hough = houghtf(I,"line",angles);

% Detect maxima in the Hough transform -- this is a bit of a hack
detect = hough>.5*max(hough(:));
detect = bwmorph(detect,'shrink',inf);
[ii, jj] = find(detect);
r = ii - (size(hough,1)-1)/2;
theta = angles(jj);

% Remove duplicate points by deleting any point with an angle 
% outside of the interval [0,180).
dup = theta<-1e-6 | theta>=pi-1e-6;
r(dup) = [];
theta(dup) = [];

% Compute line parameters (using Octave's implicit singleton expansion)
r = r(:)';
theta = theta(:)';
x = repmat([1;80],1,length(r)); % 2xN matrix, N==length(r)
y = (r - x.*cos(theta))./sin(theta); % solve line equation for y
% The above goes wrong when theta==0, fix that:
horizontal = theta < 1e-6;
x(:,horizontal) = r(horizontal);
y(:,horizontal) = [1;80];

% Plot
figure
imshow(I)
hold on
plot(y,x,'r-','linewidth',2) % note the reversed x and y!

由于我们的检测方式，对角线偏离了一个像素。我们从不寻找局部最大值的位置，我们将所有像素都放在阈值以上，然后选择阈值中间的点。

GNU Octave：霍夫变换

GNU Octave: Hough Transform

octave

hough-transform