使用图像处理计算弯管的曲率（霍夫变换抛物线检测）

Posted 2023-02-14

【中文标题】使用图像处理计算弯管的曲率（霍夫变换抛物线检测）

【英文标题】：Compute curvature of a bent pipe using image processing (Hough transform parabola detection)

【发布时间】：2016-09-28 20:45:34

【问题描述】：

我正在尝试设计一种方法来检测此管道的曲率。我尝试应用霍夫变换并发现检测到的线，但它们不位于管道表面，因此平滑它以适应贝泽曲线不起作用。请建议一些好的方法来开始这样的图像。[

hough变换检测线条得到的图像如下 [ 我正在使用标准 Matlab 代码进行概率霍夫变换线检测，该代码生成围绕结构的线段。本质上，管道的形状类似于抛物线，但对于霍夫抛物线检测，我需要在检测之前提供点的偏心率。请提出一种沿着可以拟合抛物线的曲率找到离散点的好方法。我已经给opencv和ITK打了标签，所以如果有可以在这张特定图片上实现的功能，请提出这个功能，我会尝试看看结果。

img = imread('test2.jpg');
rawimg = rgb2gray(img);
 [accum, axis_rho, axis_theta, lineprm, lineseg] = Hough_Grd(bwtu, 8, 0.01);
figure(1); imagesc(axis_theta*(180/pi), axis_rho, accum); axis xy;
 xlabel('Theta (degree)'); ylabel('Pho (pixels)');
 title('Accumulation Array from Hough Transform');
  figure(2); imagesc(bwtu); colormap('gray'); axis image;
  DrawLines_2Ends(lineseg);
  title('Raw Image with Line Segments Detected');

图像的边缘图如下，对边缘图应用霍夫变换后生成的结果也不好。我在想一个像这条曲线一样进行一般参数形状检测的解决方案可以表示为抛物线族，因此我们进行曲线拟合以估计其弯曲时的系数以分析其曲率。我需要设计一个实时程序，所以请在这个方向提出任何建议。

【问题讨论】：

嗯，只是猜测：如果您先将原件通过精巧的边缘检测器会怎样？

我尝试将霍夫变换应用于边缘图，但并没有从霍夫变换中得到更好的线检测。

我明白了。直线的结果要差得多，但抛物线的结果可能更好。也许您可以在第一种方法中使用线的交点获得可靠的点：似乎管道边缘的交叉点密度更高，但我不知道这是否能很好地推广到其他图片。我是CV新手，如果这里没有人回复，我只能建议在dsp.或cv. StackExchange站点提问。

虽然，边缘贴图实际上似乎更好的材料，但霍夫可能会失败，因为直线少得多。如果您尝试使用多个偏心率值检测抛物线，并选择最适合的那一个呢？

我也进行了过滤以消除噪音，但您的建议是正确的，即 line fit 无法正确估计 curve 。我会尝试不同偏心的建议，并使用最合适的并再次更新图片。非常感谢您的建议，它们真的很有帮助..

【参考方案1】：

我建议以下方法：

第一阶段：生成管道分段。

对图像执行阈值处理。 找到阈值图像中的连通分量。 搜索代表管道的连接组件。 代表管道的连接组件应该有一个边缘图，分为顶部和底部边缘（见附图）。 顶部和底部边缘的大小应该相似，并且它们之间的距离应该相对恒定。换句话说，它们的每像素距离的方差应该很低。

第二阶段 - 提取曲线

在这个阶段，您应该提取曲线的点以进行 Beizer 拟合。 您可以在顶部边缘或底部边缘执行此计算。 另一种选择是在管道分割的骨架上进行。

结果

管道分割。上下边缘分别用蓝色和红色标记。

代码

I = mat2gray(imread('ILwH7.jpg'));
im = rgb2gray(I);

%constant values to be used later on
BW_THRESHOLD = 0.64;
MIN_CC_SIZE = 50;
VAR_THRESHOLD = 2;
SIMILAR_SIZE_THRESHOLD = 0.85;

%stage 1 - thresholding & noise cleaning
bwIm = im>BW_THRESHOLD;
bwIm = imfill(bwIm,'holes');
bwIm = imopen(bwIm,strel('disk',1));


CC = bwconncomp(bwIm);
%iterates over the CC list, and searches for the CC which represents the
%pipe
for ii=1:length(CC.PixelIdxList)
    %ignore small CC
    if(length(CC.PixelIdxListii)<50)
        continue;
    end
    %extracts CC edges
    ccMask = zeros(size(bwIm));
    ccMask(CC.PixelIdxListii) = 1;
    ccMaskEdges = edge(ccMask);

    %finds connected components in the edges mat(there should be two).
    %these are the top and bottom parts of the pipe.
    CC2 = bwconncomp(ccMaskEdges);
    if length(CC2.PixelIdxList)~=2
        continue;
    end

    %tests that the top and bottom edges has similar sizes
    s1 = length(CC2.PixelIdxList1);
    s2 = length(CC2.PixelIdxList2);
    if(min(s1,s2)/max(s1,s2) < SIMILAR_SIZE_THRESHOLD)
        continue;
    end

    %calculate the masks of these two connected compnents
    topEdgeMask = false(size(ccMask));
    topEdgeMask(CC2.PixelIdxList1) = true;
    bottomEdgeMask = false(size(ccMask));
    bottomEdgeMask(CC2.PixelIdxList2) = true;

    %tests that the variance of the distances between the points is low
    topEdgeDists = bwdist(topEdgeMask);
    bottomEdgeDists = bwdist(bottomEdgeMask);
    var1 = std(topEdgeDists(bottomEdgeMask));
    var2 = std(bottomEdgeDists(topEdgeMask));

    %if the variances are low - we have found the CC of the pipe. break!
    if(var1<VAR_THRESHOLD && var2<VAR_THRESHOLD)
        pipeMask = ccMask;
        break;
    end


end

%performs median filtering on the top and bottom boundaries.
MEDIAN_SIZE =5;
[topCorveY, topCurveX] = find(topEdgeMask);
topCurveX = medfilt1(topCurveX);
topCurveY = medfilt1(topCurveY);
[bottomCorveY, bottomCurveX] = find(bottomEdgeMask);
bottomCurveX = medfilt1(bottomCurveX);
bottomCorveY = medfilt1(bottomCorveY);

%display results
imshow(pipeMask); hold on;
plot(topCurveX,topCorveY,'.-');
plot(bottomCurveX,bottomCorveY,'.-');

评论

在这个具体的例子中，通过阈值获取管道分割相对容易。在某些场景中，它可能更复杂。在这些情况下，您可能需要使用区域增长算法来生成管道分割。

检测表示管道的连接组件可以通过使用一些更直观的方法来完成。例如 - 它的边界的局部曲率应该很低。

【参考方案2】：

您可以找到倒置边缘图图像的连通分量 (CC)。然后，您可以以某种方式过滤这些组件，例如，基于它们的像素数，使用区域属性。这是我使用给定的 Octave 代码获得的连接组件。 现在，您可以使用 nlinfit 或任何合适的方法为每个 CC 拟合模型。

im = imread('uFBtU.png');
gr = rgb2gray(uint8(im));
er = imerode(gr, ones(3)) < .5;

[lbl, n] = bwlabel(er, 8);
imshow(label2rgb(lbl))