对 rgb 颜色分量执行 sobel 过滤器

Question

我正在尝试将图像分离为其 rgb 颜色分量，然后尝试对每个颜色分量图像执行 sobel h 和 v 过滤器。我不太确定问题是什么，但我得到了以下除以 0 的错误。我认为这可能是因为我的 u 和 v 原来是完全相同的数组。

错误

/Users/Sam/PycharmProjects/temp/A2.py:51: RuntimeWarning: divide by zero encountered in true_divide
  theta = 0.5 * np.arctan(2 * gxy / (gxx - gyy))
/Users/Sam/PycharmProjects/temp/A2.py:51: RuntimeWarning: invalid value encountered in true_divide
  theta = 0.5 * np.arctan(2 * gxy / (gxx - gyy))
/Users/Sam/PycharmProjects/temp/A2.py:54: RuntimeWarning: invalid value encountered in sqrt
  fTheta = np.sqrt(0.5 * ((gxx + gyy) + (gxx - gyy) * np.cos(2 * theta) + (2 * gxy * np.sin(2 * theta))))
Traceback (most recent call last):
  File "/Users/Sam/PycharmProjects/A1/venv/lib/python3.6/site-packages/numpy/core/fromnumeric.py", line 51, in _wrapfunc
    return getattr(obj, method)(*args, **kwds)

代码

import skimage.filters as filt
import numpy as np
def color_dot_product(A, B):
    return np.sum(A.conj() * B, axis=2)


mushroom = io.imread('mushroom.jpg')
I = util.img_as_float(mushroom)
  red = I[:, :, 0]  # Zero out contribution from green
  blue = I[:,:,1]
  green = I[:,:,2]

  # Compute horizontal and vertical derivatives of red, blue and green channels through applying sobel filters
  u = I.copy()
  u[:, :, 0] = filt.sobel_h(red)
  u[:, :, 1] = filt.sobel_h(green)
  u[:, :, 2] = filt.sobel_h(blue)
  v = I.copy()
  v[:, :, 0] = filt.sobel_v(red)
  v[:, :, 1] = filt.sobel_v(green)
  v[:, :, 2] = filt.sobel_v(blue)
  gxx = color_dot_product(u, u)
  gyy = color_dot_product(v, v)
  gxy = color_dot_product(u, v)

  # Calculate gradient direction (direction of the largest colour change)
  theta = 0.5 * np.arctan(2 * gxy / (gxx - gyy))

  # Calculate the magnitude of the rate of change
  fTheta = np.sqrt(0.5 * ((gxx + gyy) + (gxx - gyy) * np.cos(2 * theta) + (2 * gxy * np.sin(2 * theta))))

Answer 1

当除以这样的图像时，很可能某个像素的值为零，因此您会得到除以零。

这通常可以通过在除数之前显式测试每个像素是否为零来避免，或者如果除数为非负数，则添加一个非常小的值。

但在这种情况下，您根本不需要划分。函数atan2 接受两个输入参数，因此atan2(y,x) 等价于atan(y/x)。除了它返回范围内的角度 (-π,π] （即，它为您提供完整的 360 度角度范围）。我喜欢上面的维基百科，因为atan2 是一个通用函数，存在于所有语言中。除了Numpy 它被称为arctan2（我猜想与众不同的奇怪愿望？）。所以替换：

theta = 0.5 * np.arctan(2 * gxy / (gxx - gyy))

与：

theta = 0.5 * np.arctan2(2 * gxy, gxx - gyy)

---

您计算的gxx、gyy、gxy 三元组看起来很像结构张量。如果这是您打算计算的内容，则需要为这三个组件中的每一个添加额外的模糊。这会导致梯度信息在局部进行平均，从而产生更少的零梯度位置。