前回の記事に続いて、本記事では二次微分によるエッジ検出器について記していきます。

二次微分を利用したエッジ検出の理論

前回と同様に、二次微分の数学的表現からフィルタの形を求めてみます。

二次微分では、前回の記事にも記したように、色の急激な変化点にゼロ交差点が発生します。
そのゼロ交差点を用いてエッジ検出を行っていきます。

連続系の$x$方向、$y$方向の二次微分の和をラプラシアンといい、$\nabla^2$で表されます。

\nabla^2 = \frac{\partial}{\partial x^2} + \frac{\partial}{\partial y^2}

また、離散系の二次微分は、前回記事の$I(x, y)$を用いて以下のように計算ができます。

\begin{align}
I_{xx}(x, y) &= {I(x + 1, y) - I(x, y)} - {I(x, y) - I(x - 1, y)}\\
&= I(x + 1, y) + I(x - 1, y) - 2I(x, y)
\end{align}

\begin{align}
I_{yy}(x, y) &= {I(x, y + 1) - I(x, y)} - {I(x, y) - I(x, y - 1)}\\
&= I(x, y + 1) + I(x, y - 1) - 2I(x, y)
\end{align}

以上のことから、ラプラシアンのカーネルは以下のように表されます。

F = 
\begin{pmatrix}
0 & 1 & 0 \\
1 & -4 & 1 \\
0 & 1 & 0
\end{pmatrix}

斜め方向を考慮した場合は以下のようになります。

F = 
\begin{pmatrix}
1 & 1 & 1 \\
1 & -8 & 1 \\
1 & 1 & 1
\end{pmatrix}

ただ、このままだと、ノイズの影響を避けられないため、通常はガウシアンフィルタと併用されます。
そのため、ガウシアンフィルタで平滑化を行い、ラプラシアンフィルタでエッジ検出する流れをラプラシアン・ガウシアンフィルタ(Laplacian of Gasussian filter, LoG)と言います。
ラプラシアン、ガウシアンフィルタは通常以下の式で表されます。

\begin{align}
LoG(x, y, \sigma) &= \frac{\partial}{\partial x^2} G + \frac{\partial}{\partial y^2} \\
&= - \frac{1}{\pi\sigma^2}\biggl(1 - \frac{x^2 + y^2}{2\sigma^2}\biggr) \exp\biggl(- \frac{x^2 + y^2}{2\sigma^2}\biggr)
\end{align}

エッジ検出のPythonでの実装

ラプラシアン・ガウシアンフィルタの実装

前回のコードとあまり変わりません。後半を修正したという感じになります。

Pythonでの実装

import cv2
import numpy as np

# カーネル関数
def kernel_function(x, y, sigma):
    z = (1 / (2 * np.pi * sigma ** 2)) * np.exp(-(x ** 2 + y ** 2) / (2 * sigma ** 2))
    return z

# ガウシアンフィルタ用のカーネルを作成する関数
def create_kernel(size=(3, 3), sigma=1):
    kernel = np.zeros(size)
    kerner_x, kernel_y = int((size[1] - 1) / 2), int((size[0] - 1) / 2)
    for x in range(size[1]):
        for y in range(size[0]):
            kernel[y][x] = kernel_function(x - kerner_x, y - kernel_y, sigma)
    # カーネルの規格化
    kernel_norm = kernel / np.sum(kernel)
    return kernel_norm

# ガウシアンフィルタの計算
def gaussian_filter(img, kernel):
    # カーネルサイズを取得
    line, column = kernel.shape
    # フィルタをかける画像の高さと幅を取得
    height, width = img.shape

    # 畳み込み演算をしない領域の幅を指定
    height_ignore = int((line - 1) / 2)
    width_ignore = int((column - 1) / 2)

    # 出力画像用の配列
    img_filter = img.copy()

    # フィルタリングの計算
    for y in range(height_ignore, height - height_ignore):
        for x in range(width_ignore, width - width_ignore):
            img_filter[y][x] = np.sum(img[y - height_ignore: y + height_ignore + 1, x - width_ignore: x + width_ignore + 1] * kernel)
    return img_filter

# フィルタ計算する関数
def laplacian_filter(src, kernel):
    # カーネルサイズ
    m, n = kernel.shape
    # 畳み込み演算をしない領域の幅
    ignore = int((m - 1) / 2)
    height, width = src.shape
    
    img_edge = np.zeros((height, width))
    
    for y in range(ignore, height - ignore):
        for x in range(ignore, weight - ignore):
            # 畳み込み演算
            img_edge[y][x] = np.sum(src[y - ignore:y + ignore + 1, x - ignore: x + ignore + 1] * kernel)
    
    return img_edge

# 画像を白黒画像としてインポート
img_gray = cv2.imread('rankuru_noise.jpg', cv2.IMREAD_GRAYSCALE)

# カーネルを定義(行数と列数は奇数である必要あり)
gaussian_kernel = create_kernel(size=(3, 3), sigma=5)

# ガウシアンフィルタを適用
img_filter = gaussian_filter(img_gray, gaussian_kernel)

# ラプラシアンフィルタのカーネルの定義
laplacian_kernel = np.array([
    [1, 1,  1],
    [1, -8, 1],
    [1, 1,  1],
])

# ラプラシアンフィルタを適用
img_edge = laplacian_filter(img_filter, laplacian_kernel)

上記でラプラシアン・ガウシアンフィルタが適用できます。
計算結果を表示してみます。

import matplotlib.pyplot as plt
fig, ax = plt.subplots(nrows=1, ncols=2, figsize=(14, 6))
ax = ax.flatten()
ax[0].imshow(gray_img, cmap='gray')
ax[1].imshow(dst, cmap='gray_r')
ax[0].set_xticks([]);ax[0].set_yticks([])
ax[1].set_xticks([]);ax[1].set_yticks([])
plt.show()

元画像と得られた画像の比較は以下のようになります。

この画像では少し分かりにくいので、結果をもっと分かりやすくするために更なる手順が必要だと考えられます。

まとめ

二次微分を用いたエッジ検出である「ラプラシアン・ガウシアンフィルタ」について実装を行いました。
問題点などあれば指摘していただけると幸いです。

以上になります。