画像データの勾配を可視化

画像のデータ構造

私たちが画像を認識するときは２次元としてみていますが、基本的な画像のデータ構造は座標情報、256段階の明るさと3チャンネルの色情報で構成されています。

今回は画像データを1チャンネルの行列として読み込み、三次元空間にプロットすることで勾配を可視化してみようと思います。

画像の行列化

numpyを使ってjpg画像を行列に変換していきます。

例えば500×500ピクセルの画像の場合、行列のサイズも同様になります。

from PIL import Image
import numpy as np
img = np.array(Image.open('sample_1.jpg').convert('L'))
print(im.shape)

画像の断面図

画像を行列に変換できたので、特定の行や列を指定するとその断面図を見ることができます。

試しに、いつも使っている酵母の画像を行列化して、その断面図を見てみます。

https://en.wikipedia.org/wiki/Saccharomyces_cerevisiae

440×440のデータなので、323行目の断面図を見てみます。

img = np.array(Image.open('sample_1.jpg').convert('L'))
y = img[322]
fig = plt.figure()
plt.plot([i for i in range(len(y))],y)
fig.savefig('fig1-a.jpg',dpi = 500)

縦軸が明度、横軸が323行目における座標です。

以下が出力結果です。

わかりやすくするために、323行目に線を入れてみます。

img[322] = [0 for i in range(len(img[0]))]
cv2.imwrite('result.jpg',img)

以下が出力結果です。

連続的な断面図の出力

上記で作成した断面図はそのほかに画像配列の行の数だけ存在することになります。

これを、連続的に出力することでCTスキャンのような連続的な出力結果を見ることができます。

for i in range(img.shape[0]):
    fig = plt.figure()
    plt.plot([i for i in range(len(img[i]))],img[i])
    fig.savefig('fig.jpg')
    time.sleep(1)
    print(i)

ここで、matplotlibによるグラフの出力に少し時間がかかることを考慮してパッチ処理を適用します。

三次元空間にプロット

上記の連続断面図は、三次元空間を利用することで一度に表すことができます。

プロット数が膨大になるので、tqdmで処理状況を可視化しておきます。

from mpl_toolkits.mplot3d import Axes3D
import matplotlib.pyplot as plt
from tqdm import tqdm
import numpy as np

fig = plt.figure(figsize = (8, 8))
img = np.array(Image.open('sample_1.jpg').convert('L'))
ax = fig.add_subplot(111, projection='3d')
x = [i for i in range(img.shape[0])]
y = [i for i in range(img.shape[1])]
for i in tqdm(x):
    for j in y:
        ax.scatter(i, j, img[i][j], s = 1, c = "blue")
fig.savefig('3d.jpg',dpi=500)

以下が出力結果です。

プロット数が多すぎてちょっとみにくいですね。

ラプラシアンカーネルによる平滑化処理後のマッピング

そのままの画像ではプロットが散在しすぎて見づらいため、ラプラシアンカーネルによる平滑化処理を行った画像を使用して見たいと思います。

ラプラシアンカーネル（二次微分）の計算と実装

畳み込み演算まで一気に実装していきます。

import cv2
import numpy as np
import matplotlib.pyplot as plt

kernel = np.array([[1, 1, 1],
                   [1, -8, 1],
                   [1, 1, 1]])
                   
def filter2d(src, kernel):
    m, n = kernel.shape
    h = src.shape[0]
    w = src.shape[1]
    dst = np.zeros((h, w))
    for y in range(int((m-1)/2), h - int((m-1)/2)):
        for x in range(int((m-1)/2), w - int((m-1)/2)):
            dst[y][x] = np.sum(src[y-int((m-1)/2):y+int((m-1)/2)+1, x-int((m-1)/2):x+int((m-1)/2)+1]*kernel)
    return dst

img_raw = cv2.imread('sample_1.jpg')
img_gray = cv2.cvtColor(img_raw,cv2.COLOR_RGB2GRAY)
dst = filter2d(img_gray, kernel)
cv2.imwrite('result_Laplacian.jpg',dst)

上記で得られた画像を行列化します。この時、明度が50以下のピクセルは0とみなす処理を加えました。

from mpl_toolkits.mplot3d import Axes3D
import matplotlib.pyplot as plt
from tqdm import tqdm
import numpy as np
fig = plt.figure(figsize = (8, 8))
img = np.array(Image.open('result_Laplacian.jpg').convert('L'))
ax = fig.add_subplot(111, projection='3d')
x = [i for i in range(img.shape[0])]
y = [i for i in range(img.shape[1])]
for i in tqdm(x):
    for j in y:
        k = img[i][j]
        if k<50:
            ax.scatter(i, j,0, s = 1, c='white')
        else:
            ax.scatter(i, j,k, s = 1, c='black')
plt.show()
fig.savefig('3d.jpg',dpi=500)

以下が出力結果です。