【Python】OpenCVで特徴点の追跡【メダカの軌跡】
はじめに
タイトルのサムネイル写真は、メダカの動きをOpenCVでコーナー検出し、特徴点の追跡を行なって軌跡を描いたものです。この記事ではそのやり方を解説していきます。
Python上のOpenCVで calcOpticalFlowPyrLK
関数によって特徴点の追跡が可能です。つまりは、動画フレーム内で特徴点が移動した場所を追跡できるようになります。特徴点とは、
コーナー検出
などによって得られる画像内の物体の特徴的な座標になります。
準備
こちらのメダカビオトープの動画素材をOpenCVで加工していきます。
ここではPython3.xとOpenCVを使用しました。
$ pip list | grep opencv
opencv-contrib-python 4.6.0.66
opencv-python 4.5.5.62OpenCVがまだインストールされてない場合は、pipコマンドでインストールしましょう。numpyも使いますので合わせてインストールします。
$ pip install opencv-python
$ pip install opencv-contrib-python
$ pip install numpyHSV色空間でメダカのみを抽出する
で行なったような方法で、メダカのみをHSV色空間のフィルタリングで抽出してみます。
カラーピッカーでメダカのカラダのRGB色を取得して、HSV色空間へ変換、特定のHSV色のみを抽出した例です:
import cv2
import numpy as np
import itertools
import colorsys
def rgb2hsv(rgb):
r,g,b = rgb
h,s,v = colorsys.rgb_to_hsv(r, g, b)
return int(h*180), int(s*255), int(v)
def hsv_maks(img, rgb_list=None):
for i, hsv in enumerate([rgb2hsv(rgb) for rgb in rgb_list]):
h,s,v = hsv
lower = np.array([h-3, s-4, v-4])
upper = np.array([h+3, s+4, v+4])
src = cv2.inRange(img,lower,upper)
if i == 0:
mask_hsv = src
else:
mask_hsv = cv2.addWeighted(src1=src,alpha=1,src2=mask_hsv,beta=1,gamma=0)
return mask_hsv
def main(inpath, outpath, rgb_list):
# -------------------------------------------------
# setup
# -------------------------------------------------
vidcap = cv2.VideoCapture(inpath)
width = int(vidcap.get(cv2.CAP_PROP_FRAME_WIDTH))
height = int(vidcap.get(cv2.CAP_PROP_FRAME_HEIGHT))
fps = 30 #vidcap.get(cv2.CAP_PROP_FPS)
fmt = cv2.VideoWriter_fourcc('m', 'p', '4', 'v')
writer = cv2.VideoWriter(outpath, fmt, fps, (width, height)) # ライター作成
counter = 0
# -------------------------------------------------
# loop
# -------------------------------------------------
while True:
grabbed_frame, frame = vidcap.read()
if frame is None:
break
img_hsv = cv2.cvtColor(frame, cv2.COLOR_BGR2HSV)
mask_hsv = hsv_maks(img_hsv, rgb_list=rgb_list)
img = cv2.bitwise_and(frame, frame, mask=mask_hsv)
cv2.imshow('img', img)
img = cv2.resize(img, (width, height))
writer.write(img)
key = cv2.waitKey(5)
if key == 27: # Esc
break
counter += 1
vidcap.release()
writer.release()
if __name__ == '__main__':
inpath = '../assets/biotope_short.mov' # 1280 x 720
outpath = '../build/medaka_hsv.mp4'
rgb_list = [
[156, 134, 70], # 金色メダカ
[144,146,80], # 金色メダカのお腹
[217, 186, 137], # 金色メダカ
[122, 111, 48], # 金色メダカ
[150,120,56], # 金色メダカ
[133,106,42], # 金色メダカ
[137,122,52], # 金色メダカ
[144,145,80], # 金色メダカ
[185,149,77], # 金色メダカ
[130,117,69], # 金色メダカ
[89, 73, 47], # 黒メダカ
[115, 105, 62], # 黒メダカ
[79, 61, 29], # 黒メダカ
[100, 97, 47], # 黒メダカ
[120, 101, 46], # 黒メダカ
[108, 110, 85], # 黒メダカのお腹
[69,67,31], # 黒メダカ
[87,77,16], # 黒メダカ
[103,104,26], # 黒メダカ
[80,75,43], # 黒メダカ
[81,79,43], # 黒メダカ
[98,87,54], # 黒メダカ
[147,156,156], # 黒メダカのお腹
[133,135,139], # 黒メダカのお腹
[140,133,135], # 黒メダカのお腹
[53,46,12], # 黒メダカ
[78,69,18], # 黒メダカ
[77,85,42], # 黒メダカ
[90,82,40], # 黒メダカ
]
main(inpath, outpath, rgb_list)実行結果が次になります。 
少しわかりづらいですが、メダカのカラダの色を捉えることができてます。
hsv_maks() 関数内のlower と
upperでHSV色空間の範囲を持たせてます。範囲をゆるめればメダカのカラダの色がもう少しはっきりしますが、今度は葉っぱなど他の物体の色と被ってしまうので狭めに設定しました。
その代わり、メダカのカラダの色をカラーピッカーで取得してリスト化する作業が大変になります。
特徴点の検出（コーナー検出）
Shi-Tomasiコーナー検出のgoodFeaturesToTrack()
関数を使ってコーナー検出した例です:
import cv2
import numpy as np
import itertools
import colorsys
def rgb2hsv(rgb):
r,g,b = rgb
h,s,v = colorsys.rgb_to_hsv(r, g, b)
return int(h*180), int(s*255), int(v)
def hsv_maks(img, rgb_list=None):
for i, hsv in enumerate([rgb2hsv(rgb) for rgb in rgb_list]):
h,s,v = hsv
lower = np.array([h-3, s-4, v-4])
upper = np.array([h+3, s+4, v+4])
src = cv2.inRange(img,lower,upper)
if i == 0:
mask_hsv = src
else:
mask_hsv = cv2.addWeighted(src1=src,alpha=1,src2=mask_hsv,beta=1,gamma=0)
return mask_hsv
def fill(size, color):
w, h = size
canvas = np.zeros((h, w, 3), dtype="uint8")
cv2.rectangle(canvas, (0,0), (w,h), color, -1)
return canvas
def is_nearby(pt1, pt2, maxDistance, minDistance):
r = (pt1[0] - pt2[0])**2 + (pt1[1] - pt2[1])**2
if r < maxDistance**2 and r > minDistance**2:
return True
return False
def clamp(n, smallest, largest):
return max(smallest, min(n, largest))
def main(inpath, outpath, rgb_list):
# -------------------------------------------------
# setup
# -------------------------------------------------
vidcap = cv2.VideoCapture(inpath)
width = int(vidcap.get(cv2.CAP_PROP_FRAME_WIDTH))
height = int(vidcap.get(cv2.CAP_PROP_FRAME_HEIGHT))
fps = 30 #vidcap.get(cv2.CAP_PROP_FPS)
fmt = cv2.VideoWriter_fourcc('m', 'p', '4', 'v')
writer = cv2.VideoWriter(outpath, fmt, fps, (width, height)) # ライター作成
radius = 2
counter = 0
# -------------------------------------------------
# loop
# -------------------------------------------------
while True:
grabbed_frame, frame = vidcap.read()
if frame is None:
break
# HSVマスク作成
img_hsv = cv2.cvtColor(frame, cv2.COLOR_BGR2HSV)
mask_hsv = hsv_maks(img_hsv, rgb_list=rgb_list)
img_bgr = cv2.bitwise_and(frame, frame, mask=mask_hsv)
# グレースケール変換
gray = cv2.cvtColor(img_bgr, cv2.COLOR_BGR2GRAY)
gray = cv2.dilate(gray,None,iterations = 1)
# Shi-Tomasiコーナー検出
corners = cv2.goodFeaturesToTrack(gray, maxCorners=180, qualityLevel=0.001, minDistance=12)
corners = np.int0(corners)
# キャンバスへ描画
canvas = np.zeros_like(frame)
for i in corners:
pt = i.ravel()
cv2.circle(canvas, pt, radius, (255,255,255), -1)
for pair in itertools.combinations(corners, 2): # 組み合わせ
d1, d2 = pair
pt1 = d1.ravel()
pt2 = d2.ravel()
if is_nearby(pt1, pt2, 30, 0):
cv2.line(canvas, pt1, pt2, (255,255,255))
# 書き出し
cv2.imshow('Canvas', canvas)
canvas = cv2.resize(canvas, (width, height))
writer.write(canvas)
key = cv2.waitKey(5)
if key == 27: # Esc
break
counter += 1
vidcap.release()
writer.release()
if __name__ == '__main__':
inpath = '../assets/biotope_short.mov' # 1280 x 720
outpath = '../build/medaka_detected_corner.mp4'
rgb_list = [
[156, 134, 70], # 金色メダカ
[144,146,80], # 金色メダカのお腹
[217, 186, 137], # 金色メダカ
[122, 111, 48], # 金色メダカ
[150,120,56], # 金色メダカ
[133,106,42], # 金色メダカ
[137,122,52], # 金色メダカ
[144,145,80], # 金色メダカ
[185,149,77], # 金色メダカ
[130,117,69], # 金色メダカ
[89, 73, 47], # 黒メダカ
[115, 105, 62], # 黒メダカ
[79, 61, 29], # 黒メダカ
[100, 97, 47], # 黒メダカ
[120, 101, 46], # 黒メダカ
[108, 110, 85], # 黒メダカのお腹
[69,67,31], # 黒メダカ
[87,77,16], # 黒メダカ
[103,104,26], # 黒メダカ
[80,75,43], # 黒メダカ
[81,79,43], # 黒メダカ
[98,87,54], # 黒メダカ
[147,156,156], # 黒メダカのお腹
[133,135,139], # 黒メダカのお腹
[140,133,135], # 黒メダカのお腹
[53,46,12], # 黒メダカ
[78,69,18], # 黒メダカ
[77,85,42], # 黒メダカ
[90,82,40], # 黒メダカ
]
main(inpath, outpath, rgb_list)実行結果が次になります。 
Processingのジェネラティブアートっぽく、点と点を線で結んでみました。一定の距離以下の点と点を結ぶ方法は比較的簡単でして、やり方は の記事で解説してます。
特徴点を追跡する（モーション解析）
コーナー検出によって見つけられた特徴点を、動画フレームの前後で保持（追跡）する便利な機能がOpenCVに備えられてます。それが
calcOpticalFlowPyrLK() 関数です。
calcOpticalFlowPyrLK() は Lucas-Kanade法（英語論文）
を使ったオプティカルフローと呼ばれるものです。入力には、以前と現在の二つのグレースケール画像と、以前の特徴点の情報が必要になります。結果として移動後の特徴点の配列を得られます。
次はメダカの特徴点（軌跡）を追跡した例です:
import cv2
import numpy as np
import colorsys
def rgb2hsv(rgb):
r,g,b = rgb
h,s,v = colorsys.rgb_to_hsv(r, g, b)
return int(h*180), int(s*255), int(v)
def hsv_maks(img, rgb_list=None):
for i, hsv in enumerate([rgb2hsv(rgb) for rgb in rgb_list]):
h,s,v = hsv
lower = np.array([h-3, s-4, v-4])
upper = np.array([h+3, s+4, v+4])
src = cv2.inRange(img,lower,upper)
if i == 0:
mask_hsv = src
else:
mask_hsv = cv2.addWeighted(src1=src,alpha=1,src2=mask_hsv,beta=1,gamma=0)
return mask_hsv
def clamp(n, smallest, largest):
return max(smallest, min(n, largest))
def main(inpath, outpath, rgb_list):
# -------------------------------------------------
# setup
# -------------------------------------------------
vidcap = cv2.VideoCapture(inpath)
width = int(vidcap.get(cv2.CAP_PROP_FRAME_WIDTH))
height = int(vidcap.get(cv2.CAP_PROP_FRAME_HEIGHT))
fps = 30 #vidcap.get(cv2.CAP_PROP_FPS)
fmt = cv2.VideoWriter_fourcc('m', 'p', '4', 'v')
writer = cv2.VideoWriter(outpath, fmt, fps, (width, height)) # ライター作成
radius = 2
counter = 0
alpha = 0.0
beta = 1.0
# Parameters for lucas kanade optical flow
lk_params = dict( winSize = (15,15),
maxLevel = 2,
criteria = (cv2.TERM_CRITERIA_EPS | cv2.TERM_CRITERIA_COUNT, 10, 0.03))
grabbed_frame, frame = vidcap.read()
old_masks = []
# -------------------------------------------------
# loop
# -------------------------------------------------
while True:
grabbed_frame, frame = vidcap.read()
if frame is None:
break
mask = np.zeros_like(frame)
img_hsv = cv2.cvtColor(frame, cv2.COLOR_BGR2HSV)
mask_hsv = hsv_maks(img_hsv, rgb_list=rgb_list)
img_bgr = cv2.bitwise_and(frame, frame, mask=mask_hsv)
gray = cv2.cvtColor(img_bgr, cv2.COLOR_BGR2GRAY)
corners = cv2.goodFeaturesToTrack(gray, maxCorners=180, qualityLevel=0.001, minDistance=12)
if counter > 0:
new_pts, st, err = cv2.calcOpticalFlowPyrLK(old_gray, gray, old_pts, None, **lk_params)
good_new = new_pts[st==1]
good_old = old_pts[st==1]
for new,old in zip(good_new,good_old):
pt0 = np.int0(new.ravel())
pt1 = np.int0(old.ravel())
cv2.line(mask, pt0, pt1, (255,128,0), 2)
cv2.circle(mask, pt1, radius, (255,0,255), -1)
old_masks.insert(0, mask)
if len(old_masks) > 15:
old_masks = old_masks[:15]
a = 1
for m in old_masks:
mask = cv2.addWeighted(src1=m,alpha=a, src2=mask,beta=1,gamma=0)
a -= 0.05
alpha = clamp(alpha, 0, 1)
beta = clamp(beta, 0, 1)
cv2.imshow('img', mask)
mask = cv2.resize(mask, (width, height))
writer.write(mask)
key = cv2.waitKey(5)
if key == 27: # Esc
break
counter += 1
old_pts = corners.copy()
old_gray = gray
vidcap.release()
writer.release()
if __name__ == '__main__':
inpath = '../assets/biotope_short.mov' # 1280 x 720
outpath = '../build/medaka_flow_short.mp4'
rgb_list = [
[156, 134, 70], # 金色メダカ
[144,146,80], # 金色メダカのお腹
[217, 186, 137], # 金色メダカ
[122, 111, 48], # 金色メダカ
[150,120,56], # 金色メダカ
[133,106,42], # 金色メダカ
[137,122,52], # 金色メダカ
[144,145,80], # 金色メダカ
[185,149,77], # 金色メダカ
[130,117,69], # 金色メダカ
[89, 73, 47], # 黒メダカ
[115, 105, 62], # 黒メダカ
[79, 61, 29], # 黒メダカ
[100, 97, 47], # 黒メダカ
[120, 101, 46], # 黒メダカ
[108, 110, 85], # 黒メダカのお腹
[69,67,31], # 黒メダカ
[87,77,16], # 黒メダカ
[103,104,26], # 黒メダカ
[80,75,43], # 黒メダカ
[81,79,43], # 黒メダカ
[98,87,54], # 黒メダカ
[147,156,156], # 黒メダカのお腹
[133,135,139], # 黒メダカのお腹
[140,133,135], # 黒メダカのお腹
[53,46,12], # 黒メダカ
[78,69,18], # 黒メダカ
[77,85,42], # 黒メダカ
[90,82,40], # 黒メダカ
]
main(inpath, outpath, rgb_list)実行結果が次になります。 
▼ ぜひこちらの映像もご覧ください。
https://www.youtube.com/shorts/zQ89Vlm7oKE
calcOpticalFlowPyrLK
calcOpticalFlowPyrLK(prevImg, nextImg, prevPts, nextPts, status=..., err=..., winSize=..., maxLevel=..., criteria=..., flags: int = ..., minEigThreshold=...)| パラメータ | 意味 |
|---|---|
| prevImg | 1番目の入力画像 |
| nextImg | 2番目の入力画像 |
| prevPts | 1番目の入力画像の特徴点 |
| winSize | ピラミッドレベルにおける探索窓のサイズ |
| maxLevel | 画像ピラミッドの最大レベル数 |
| criteria | 反復探索アルゴリズムの停止基準 |
| derivLambda | 画像の空間微分の相対的な重み |
（大きな動きの場合はピラミッドのスケールをアップする）
関連記事
- 【Python】OpenCVで図形の描画からアニメーションまで【線・四角・丸・塗りつぶし】
- 【Python】OpenCVでコーナーの検出【Harris/Shi-Tomasi】
- 【Python】OpenCVで画像操作いろいろ（グレースケール・モノ・輪郭抽出・切り抜く・透過）
- 【Python】VidStabで手ぶれ補正【動画編集への道#2】
- 【Python】OpenCVで画像をアフィン変換【移動・拡大・回転・剪断】