オプティカルフローの実装(Python)
オプティカルフローを求めるための方法の1つであるLucas Kanade法を実装してみました.
勾配法
Lucal Kanade法を実装する前にまず勾配法について学ぶ必要があります.勾配法とは,連続する2枚の画像において追跡対象の物体の移動量は微小であると仮定し,オプティカルフローを求める方法です.まず時刻tにおける画像上の輝度値をI(x,y,t)としましょう.そして時間Δtの間にその画素が(Δx,Δy)動いたとしましょう.すると,
そしてIの微小変化を表現するするために,(1式)をテイラー展開で一次近似します.x,y,tが変化した分だけIも変化しその増分が加えられることを意味しています
ここで(2)について両辺のI(x,y,t)はほぼ同じとみなせるので移行して打ち消せます.そしてΔtで割ると,
そしてΔtが限りなく0に近いと仮定すると微分の形にもっていけます.つまり,
とが時間変化に対する画素の移動方向を表しており,これがオプティカルフローとなります.また,とが時間変化に対する輝度値の変化となります.そして最終的には,
と表され,時間変化に対するIの変化は,画素値の勾配と画素のオプティカルフローとの内積で表されることになります.
ここでオプティカルフロー(u,v)を求めたいわけですが,未知変数2つに対して式が一つしかありません.なのでこのままでは解けません.しかし,Lucal Kanade法ではある仮定を置くことで式の数を増やして(u,v)を求めることができます.
Lucal Kanade法
ではその仮定とは何かというと,注目画素の周辺の局所領域では(u,v)の値は同じだという仮定です.
こっからの数式は気力がなくなったのでそのうち書きます.今のところは他のサイトを参照・・
ソースコード
汚いコードですんません.なんか間違ってる気がしますが大体こんな感じです.
from scipy import ndimage import cv2 import numpy as np def opticalFlow(pre, next): image1 = cv2.cvtColor(pre, cv2.COLOR_RGB2GRAY).astype(np.float) image2 = cv2.cvtColor(next, cv2.COLOR_RGB2GRAY).astype(np.float) #輝度勾配 dx = ndimage.gaussian_filter(image1, sigma=1.0, order=(0,1)) dy = ndimage.gaussian_filter(image1, sigma=1.0, order=(1,0)) #Iの時間変化 つまり次の時間の画像と前時間の画像との差分が時間変化 de = image2 - image1 height, width = image1.shape for i in range(1, height-1, 20): for j in range(1, width-1, 20): #局所領域(3*3) patchDx = dx[i-1:i+1+1,j-1:j+1+1] patchDy = dy[i-1:i+1+1,j-1:j+1+1] patchDe = de[i-1:i+1+1,j-1:j+1+1] A = np.zeros((2,2)) A[0,0] = np.sum(patchDx**2) A[0,1] = np.sum(patchDx * patchDy) A[1,0] = A[0,1] A[1,1] = np.sum(patchDy**2) B = np.zeros((2)) B[0] = -np.sum(patchDx * patchDe) B[1] = -np.sum(patchDy * patchDe) #ランクが2じゃないと逆行列もとめられん if np.linalg.matrix_rank(A) == 2: u, v = np.linalg.solve(A, B) #フローを描画 next = cv2.line(next, (j,i), (int(j+u), int(i+v)), (255,0,0), 1) return next def main(): cap = cv2.VideoCapture(0) cv2.namedWindow("frame") preFrame = None while(True): ret, frame = cap.read() if preFrame != None: flow = opticalFlow(preFrame, frame) cv2.imshow("frame", flow) if cv2.waitKey(20) & 0xFF == 27: break preFrame = frame cap.release() cv2.destroyAllWindows() main()