如何提高Sobel边缘检测器的效率

Question

如何提高Sobel边缘检测器的效率

pythonperformancenumpyedge-detectionsobel

3

我正在使用Python从头开始编写一个计算机视觉库，以便与rpi相机配合使用。目前，我已经实现了将图像转换为灰度和其他一些基本的图像操作，这些操作在我的model B rpi3上都运行得比较快。

然而，我的边缘检测函数使用Sobel算子（维基百科描述）比其他函数慢得多，尽管它确实起作用。以下是代码：

def sobel(img):
    xKernel = np.array([[-1,0,1],[-2,0,2],[-1,0,1]])
    yKernel = np.array([[-1,-2,-1],[0,0,0],[1,2,1]])
    sobelled = np.zeros((img.shape[0]-2, img.shape[1]-2, 3), dtype="uint8")
    for y in range(1, img.shape[0]-1):
        for x in range(1, img.shape[1]-1):
            gx = np.sum(np.multiply(img[y-1:y+2, x-1:x+2], xKernel))
            gy = np.sum(np.multiply(img[y-1:y+2, x-1:x+2], yKernel))
            g = abs(gx) + abs(gy) #math.sqrt(gx ** 2 + gy ** 2) (Slower)
            g = g if g > 0 and g < 255 else (0 if g < 0 else 255)
            sobelled[y-1][x-2] = g
    return sobelled

并使用这张猫的灰度图像运行它：

我得到了这个看起来正确的回复：

该库的应用，尤其是这个功能，是在一个下棋机器人上，其中边缘检测将有助于识别棋子的位置。问题在于它需要超过15秒才能运行，这是一个重大问题，因为它会显著增加机器人下棋所需的时间。

我的问题是：如何加速它？

到目前为止，我已经尝试了一些方法：

而不是对gx和gy值进行平方、相加，然后开方得到总梯度，我只需要对绝对值进行求和。这样可以显著提高速度。
使用来自rpi摄像头的较低分辨率图像。这显然是一种简单的方法，可以使这些操作运行得更快，但它并不是非常可行，因为即使在最小可用分辨率480x360下，速度仍然大大降低，远低于相机的最大分辨率3280x2464。
编写嵌套的for循环来执行矩阵卷积，以代替np.sum(np.multiply(...))。这最终略微变慢，这让我感到惊讶，因为由于np.multiply返回一个新的数组，我认为使用循环应该更快。不过我认为这可能是因为numpy主要是用C编写的，或者新数组实际上没有被存储，所以不需要太长时间，但我不太确定。

任何帮助都将不胜感激 - 我认为需要改进的主要是第三点，即矩阵乘法和求和。

- Joe Iddon

你尝试过OpenCV的Sobel吗？还有，你尝试过2D卷积吗？ - Divakar

@Divakar 是的，我已经使用 OpenCV 实现了整个棋子检测，但我正在尝试用 Python 从头开始编写它。2D 卷积非常广泛，我以为我已经实现了它... - Joe Iddon

我不是很清楚 - 你是在说你不能使用Scipy中的2D卷积吗？还是你已经尝试过它，结果变得更慢了？ - Divakar

3个回答

1

我曾遇到同样的问题，并且使用Numba库中的@jit将我的代码加速了约600倍（请参见链接：https://numba.pydata.org/numba-doc/latest/user/5minguide.html）。在我的函数上方添加@jit（nopython=True）就足以完成任务。

- Matthijs Geerlings

1

就它的价值而言，这是一个补充：

Sobel x kernel:
-1 0 1
-2 0 2
-1 0 1

你不需要一个单独的内核。三分之一的操作总是得到零。只需不计算它们。其余部分可以简化：

sum = -inputvalue[y-1][x-1] - 2 * inputvalue[y][x-1] - inputvalue[y+1][x-1]
+ inputvalue[y-1][x+1] + 2 * inputvalue[y][x+1] + inputvalue[y+1][x+1]

相比于朴素的方法通过核函数循环9次乘法和9次加法，这里只需要2次乘法、3次减法和3次加法，而且不需要循环即可完成计算，这将显著减少计算时间。

我对上面提到的numpy示例中实现的1000倍速度提升感到惊讶。但是这种方法帮助我大大提高了速度：）

- scaszoo

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- alkasm · Accepted Answer

即使你正在建立自己的库，你也应该绝对使用卷积库，它们将在后端使用C或Fortran执行结果操作，这将快得多。

但是，如果您愿意自己实现卷积，请使用线性可分离滤波器。这是一个想法:

图像：

1 2 3 4 5
2 3 4 5 1
3 4 5 1 2

Sobel x 核：

-1 0 1
-2 0 2
-1 0 1

结果：

8, 3, -7

在卷积的第一个位置，你将计算9个值。首先，为什么？你永远不会加中间那一列，不要浪费时间去乘它。但是这不是线性可分滤波器的重点。这个想法很简单。当你将核置于第一个位置时，你将用[1, 2, 1]乘以第三列。但是两步后，你将用[-1,-2,-1]乘以第三列。真浪费！你已经计算出了它，现在只需要否定它。这就是线性可分离滤波器的思想。请注意，您可以将过滤器分解为两个向量的矩阵外积:

[1]
[2]  *  [-1, 0, 1]
[1]

在这里进行外积运算会得到相同的矩阵。因此，这里的想法是将操作分成两个部分。首先用行向量乘以整个图像，然后再用列向量。取行向量。

-1 0 1

横跨整个图像，我们最终得到

2  2  2
2  2 -3
2 -3 -3

然后将列向量传递进去进行乘法和求和运算，我们再次得到

8, 3, -7

还有一个可能有用的巧妙技巧（取决于您在内存和效率之间的权衡）：

请注意，在单行乘法中，您会忽略中间值，并仅从左值减去右值。这意味着实际上您正在执行以下操作：从这两个图像中减去了另一个图像：

3 4 5     1 2 3
4 5 1  -  2 3 4
5 1 2     3 4 5

如果您从图像中剪切掉前两列，则得到左矩阵，如果您剪切掉后两列，则得到右矩阵。因此，您可以简单地计算卷积的第一部分，如下所示：

result_h = img[:,2:] - img[:,:-2]

然后，您可以循环遍历Sobel算子的其余列。或者，您甚至可以进一步进行相同的操作。这一次是针对垂直情况，您只需要添加第一行和第三行以及两倍的第二行；或者使用numpy加法：

result_v = result_h[:-2] + result_h[2:] + 2*result_h[1:-1]

完成了！我可能会在不久的将来添加一些时间计算。对于一些草率的Jupyter笔记本定时（即1000x1000图像），进行简略计算：

新方法（图像总和）：每个循环8.18毫秒±399微秒（7次运行的平均值±标准偏差，每次循环100次）

旧方法（双重循环）：每个循环7.32秒±207毫秒（7次运行的平均值±标准偏差，每个循环1次）

没错，你没看错：速度提升了1000倍。

这里是比较这两种方法的一些代码：

import numpy as np

def sobel_x_orig(img):
    xKernel = np.array([[-1,0,1],[-2,0,2],[-1,0,1]])
    sobelled = np.zeros((img.shape[0]-2, img.shape[1]-2))
    for y in range(1, img.shape[0]-1):
        for x in range(1, img.shape[1]-1):
            sobelled[y-1, x-1] = np.sum(np.multiply(img[y-1:y+2, x-1:x+2], xKernel))
    return sobelled

def sobel_x_new(img):
    result_h = img[:,2:] - img[:,:-2]
    result_v = result_h[:-2] + result_h[2:] + 2*result_h[1:-1]
    return result_v

img = np.random.rand(1000, 1000)
sobel_new = sobel_x_new(img)
sobel_orig = sobel_x_orig(img)

assert (np.abs(sobel_new-sobel_orig) < 1e-12).all()

当然，1e-12 是一个相当严格的容差，但是这是每个元素的容差，所以应该没问题。但我有一张 float 图像，对于 uint8 图像来说，你会发现差异更大。

请注意，您可以对 任何线性可分离滤波器 进行此操作！这包括高斯滤波器。还要注意的是，总体而言，这需要大量的运算。在 C、Fortran 或其他语言中，通常只需对单行/列向量进行两次卷积实现，因为最终需要循环遍历每个矩阵的每个元素；无论是只添加还是相乘，以这种方式添加图像值在 C 中并不比执行卷积快。但通过 numpy 数组进行循环遍历非常缓慢，因此这种方法在 Python 中速度要快得多。