Oldukça çok ortaya çıkan bir imaj işleme problemi şudur: elimizde iki nokta grubu var, bu noktaların arasındaki eşleşmeyi biliyoruz. Öyle bir \(H\) ilişkisi bulmak istiyoruz ki verili \(x\) noktasınının (homojen) kordinatını \(x'\) noktasına taşısın, yani eldeki her veri noktasının ima ettiği eşleşmeyi bulsun.
Örnek
= [[25.8064516129,25.0],[23.87096,45.625],
x1 20.0,69.375],[28.387,92.5],
[38.709,116.875],[64.5161290323,115.0],
[64.516,89.375],[65.16,66.875],
[57.4193,45.0],[45.80645,23.75]]
[= [[93.548,66.25],[114.838,110.0],
x2 138.709,153.125],[182.580,179.375],
[241.935,204.375],[276.77,163.75],
[254.193,123.125],[212.903,73.125],
[158.064,54.375],[120.6451,40.625]]
[
= np.array(x1)
x1 = np.array(x2)
x2 0], x1[:,1], 'rd')
plt.plot(x1[:,0], x2[:,1], 'bd')
plt.plot(x2[:,0,320)
plt.xlim(0,240)
plt.ylim('vision_30vstab_02.png') plt.savefig(
Yani kırmızı noktaları mavi noktalara çeviren ilişkiyi arıyoruz. Bu transformasyonda ne var? Sağa doğru bir yer değiştirme (translation), ölçekleme (scaling), ve saat yönüne doğru bir döndürme (rotation). Bu tür 2D-2D ilişkilerine homografi adı veriliyor. Aradığımız alttaki türden bir formül [3],
\[ x' = H x\]
yani her \(x\) noktası \(H\) üzerinden \(x'\) haline gelecek. \(H\) matrisi homojen kordinatları baz alır,
\[ \left[\begin{array}{r} x' \\ y' \\ w' \end{array}\right] \left[\begin{array}{rrr} h_1 & h_2 & h_3 \\ h_3 & h_4 & h_5 \\ h_6 & h_7 & h_8 \end{array}\right] \left[\begin{array}{r} x \\ y \\ w \end{array}\right] \]
\(H\) matrisinin bazı şekilleri vardır, mesela
\[ \left[\begin{array}{r} x' \\ y' \\ 1 \end{array}\right] \left[\begin{array}{rrr} a_1 & a_2 & t_x \\ a_3 & a_4 & t_y \\ 0 & 0 & 1 \end{array}\right] \left[\begin{array}{r} x \\ y \\ 1 \end{array}\right] \]
Ya da matris içindeki bölgeleri vektör / matrisler ile özetlersek,
\[ x' = \left[\begin{array}{rr} A & t \\ 0 & 1 \end{array}\right] x \]
Üstteki transformasyona ilgin transformasyonu (affine transformation) deniyor, yamultma (warping) denen işlem budur. Bu transformasyon \(w=1\) şartını korur.
Eğer \(H\) şu türden olursa,
\[ \left[\begin{array}{r} x' \\ y' \\ 1 \end{array}\right] \left[\begin{array}{rrr} s\cos(\theta) & -s\sin(\theta) & t_x \\ s\sin(\theta) & s\cos(\theta) & t_y \\ 0 & 0 & 1 \end{array}\right] \left[\begin{array}{r} x \\ y \\ 1 \end{array}\right] \]
Ya da
\[ x' = \left[\begin{array}{rr} sR & t \\ 0 & 1 \end{array}\right] x \]
Dönüş \(R\), taşınma \(t\), dönme \(\theta\), ölçekleme \(s\). Bu transformasyona ölçeklemeye (scaling) izin veren bir katı (rigid) transformasyon deniyor. “Katı’’ demek \(s=1\), yani noktalar arası mesafeler değişmeyecek demek, daha büyük \(s\) ile tabii ölçekleme olabilir, mesafeler artabilir, ama mesafe oranları yine aynı kalır, ayrıca döndürme de -rotation- yapılabilir. Bu transformasyona yansıtsal (projective) ismi de verilir. Yansıtsal transformasyonun ilgin transformasyondan daha esnek / kuvvetli olduğu bilinir.
Not: ilgin transformasyon ve onu kestirme hesabı bazen literatürde iki boyutlu kordinat sisteminde ve \(x' = R x + t\), yani rotasyon artı yer değişimi gibi bir formda da görülebilir, biz homojen sisteme geçerek her ikisini aynı matris \(H\) içinde ve tek çarpım operasyonu ile gösterebilmiş oluyoruz. Homojen, tek matrisli formda hesap yapmak daha kolay.
Homografi hesabının kullanım alanları geniş; mesela elde olan iki imaj arasında birbirine uyan noktaları biliyorsak, \(H\)’yi hesaplayarak tüm imaj üzerinde bir değişim matrisi hesaplamış oluruz.
Yansıtsal \(H\) hesabı için direk lineer transform (direct linear transform -DLT-) tekniği var. Eldeki tüm eşleşmeler için alttaki sistemi yaratırız,
\[ \left[\begin{array}{rrrrrrrrr} -x_1 & -y_1 & -1 & 0 & 0 & 0 & x_1x_1' & y_1x_1' & x_1' \\ 0 & 0 & 0 & -x_1 & -y_1 & -1 & x_1y_1' & y_1y_1' & y_1' \\ -x_2 & -y_2 & -1 & 0 & 0 & 0 & x_2x_2' & y_2x_2' & x_2' \\ 0 & 0 & 0 & -x_2 & -y_2 & -1 & x_2y_2' & y_2y_2' & y_2' \\ & \vdots & & \vdots & & \vdots & & \vdots & \end{array}\right] \left[\begin{array}{r} h_1 \\ h_2 \\ h_3 \\ h_4 \\ h_5 \\ h_6 \\ h_7 \\ h_8 \\ h_9 \end{array}\right] = 0 \]
Bu sistem \(x' - Hx = 0\) sistemini temsil etmiş oluyor, ne kadar fazla nokta olursa üstteki matris o kadar aşağı doğru genişleyecektir (öğe ayarlaması öne göre yapılacak tabii). Mükemmel \(H\) bulunamayabilir, ama sıfıra olabildiğince yaklaşmak için üstteki problemi bir minimizasyon problemi olarak görürüz, SVD bu çözümü bize sağlar.
import numpy.linalg as lin
def H_from_points(fp,tp):
if fp.shape != tp.shape:
raise RuntimeError('number of points do not match')
= np.mean(fp[:2], axis=1)
m = np.max(np.std(fp[:2], axis=1)) + 1e-9
maxstd = np.diag([1/maxstd, 1/maxstd, 1])
C1 0][2] = -m[0]/maxstd
C1[1][2] = -m[1]/maxstd
C1[= np.dot(C1,fp)
fp
= np.mean(tp[:2], axis=1)
m = np.max(np.std(tp[:2], axis=1)) + 1e-9
maxstd = np.diag([1/maxstd, 1/maxstd, 1])
C2 0][2] = -m[0]/maxstd
C2[1][2] = -m[1]/maxstd
C2[= np.dot(C2,tp)
tp
= fp.shape[1]
nbr_correspondences = np.zeros((2*nbr_correspondences,9))
A for i in range(nbr_correspondences):
2*i] = [-fp[0][i],-fp[1][i],-1,0,0,0,
A[0][i]*fp[0][i],tp[0][i]*fp[1][i],tp[0][i]]
tp[2*i+1] = [0,0,0,-fp[0][i],-fp[1][i],-1,
A[1][i]*fp[0][i],tp[1][i]*fp[1][i],tp[1][i]]
tp[
= lin.svd(A)
U,S,V = V[8].reshape((3,3))
H
= np.dot(lin.inv(C2),np.dot(H,C1))
H
# normalize and return
return H / H[2,2]
= np.ones((len(x1),3))
x1h 2] = x1
x1h[:,:= np.ones((len(x1),3))
x2h 2] = x2
x2h[:,:= H_from_points(x1h.T,x2h.T)
A = np.dot(A, x1h.T).T
res = res.T / res[:,2]
res
0], x1[:,1], 'rd')
plt.plot(x1[:,0], x2[:,1], 'bd')
plt.plot(x2[:,0], res.T[:,1], 'bx')
plt.plot(res.T[:,0,320)
plt.xlim(0,240)
plt.ylim('vision_30vstab_03.png') plt.savefig(
Çarpı ile işaretli noktalar kestirme hesabı yapılan \(H\) ile kırmızı noktaların transform edilmesiyle elde edildi. Gerçek noktalara oldukca yakın!
İlgin transformasyon matrisinin hesabı için üstteki metotta \(h_7=0,h_8=0\) kullanmak yeterli. Alternatif bir yöntem de var, daha fazla detay için [3, sf. 76].
İmaj Bölgesi Çekip Çıkarmak
Üstteki tekniğin ilginç uygulamalarından biri; diyelim ki bir imajın belli bir bölgesindeki görüntüyü eşit kenarlı olacak şekilde çekip çıkarmak istiyorum, mesela alttaki Sudoku oyun karesini,
from scipy import ndimage
from PIL import Image
= np.array(Image.open('sudoku81.JPG').convert('L'))
im = [[257.4166, 14.9375],
corners 510.8489, 197.6145],
[59.30208, 269.65625],
[325.598958, 469.05729]]
[= np.array(corners)
corners 0], corners[:,1], 'rd')
plt.plot(corners[:,=plt.cm.Greys_r)
plt.imshow(im, cmap'vision_30vstab_04.png') plt.savefig(
Kenarları kırmızı noktalarla ben seçtim, şimdi o bölgenin alınıp eşit
kenarlı halde gösterilmesini istiyorum. Bu ne demektir? Bu seçilen her
köşe noktasının eşit kenarlı bir karenin köşelerine transform edilmesi
demektir, mesela bu köşeler \((1,300),(300,300),..\) gibi olabilir
(imajın en uç noktaları). Sonra daha önce yaptığım gibi \(H\) hesaplarım, ve o bölgedeki tüm
pikselleri alıp hesapladığım transformasyonu onlara uygularım,
scipy.ndimage.geometric_transform
bu işi yapar.
from scipy import ndimage
import scipy
import imageio
= [ [p[1],p[0],1] for p in corners]
fp = np.array(fp).T
fp = np.array([[0,0,1],[0,300,1],[300,0,1],[300,300,1]]).T
tp = H_from_points(tp,fp)
H
def warpfcn(x):
= np.array([x[0],x[1],1])
x = np.dot(H,x)
xt = xt/xt[2]
xt return xt[0],xt[1]
= ndimage.geometric_transform(im,warpfcn,(300,300))
im_g 'vision_30vstab_05.png', im_g) imageio.imwrite(
Video Stabilizasyonu
Elde tutulan kamera ile kaydedilen görüntülerde titreme çok olabilir. Mesela şurada [1] bizim cep telefonu ile kaydettiğimiz bir örnek var. Bu görüntüyü yazılım ile stabilize etmek mümkün mü? Cevap evet - ve çözüm şaşırtıcı derecede basit. [4]’ün tekniği şöyle özetlenebilir: bir video’yu baştan itibaren kare kare işlerken, her karede ilginç köşe noktaları (Harris tekniği ile) buluruz, ve bu noktaların bir sonraki resimdeki eşlerini elde ederiz, bu artık görüntü işlemede demirbaş haline gelmiş bir işlem. Sonra tüm eşlemeleri kullanarak her video karesi için bir homografi / transformasyon hesaplarız, bu transformasyon matrisi içinde \(x,y\) değişimi, yani taşınma, ve \(a\) açısı ki döndürme bilgisi vardır. Bu bilgileri \(dx,dy,da\) olarak biriktiririz.
Tüm kareler işlenince başa dönüp tüm bu değişimlerin kümülatif toplamını alarak \(x,y,a\) zaman serilerini oluştururuz. Bu zaman serileri üzerinde bir yürüyen ortalama (moving average) hesabı yaparız, bu bize pürüzüşleştirilmiş zaman serileri verir. Şimdi kümülatif serinin pürüzsüz seriden olan farklarını buluruz, ve her kare için bu farkı alıp, onunla bir \(H\) oluştururuz ve bu \(H\) ile bir önceki kare üzerinde yamultma yaparak onu “düzeltiriz’’. Bu kadar.
Bu teknik niye işliyor? İşliyor çünkü üstte gösterdiğimiz türde video’larda “beklenen” bir “akış’’, bir nokta eşleşmesi var. Düz yürüyoruz, kamera ileri dönük, ortadaki pikseller dışa doğru eşleşmeli, sağdakiler daha sağa doğru, vs. Bu beklentiyi hareketli ortalama ile hesaplamak mümkün, ve ondan olan sapmaları kameranın istenmeyen titremesi olarak algılıyoruz, ve düzeltiyoruz.
import cv2, numpy as np, pandas as pd
= "/opt/Downloads/bwalk1.mp4"
fin
= cv2.VideoCapture(fin)
cap = int(cap.get(cv2.CAP_PROP_FRAME_COUNT))
N = int(cap.get(cv2.CAP_PROP_FPS))
fps
= cap.read()
status, prev = cv2.cvtColor(prev, cv2.COLOR_BGR2GRAY)
prev_gray = prev.shape[:2]
(h,w)
= None
last_T = []
prev_to_cur_transform for k in range(N-1):
= cap.read()
status, cur = cv2.cvtColor(cur, cv2.COLOR_BGR2GRAY)
cur_gray = cv2.goodFeaturesToTrack(prev_gray,
prev_corner = 200,
maxCorners = 0.01,
qualityLevel = 30.0,
minDistance = 3)
blockSize = cv2.calcOpticalFlowPyrLK(prev_gray,
cur_corner, status, err
cur_gray,
prev_corner,None)
= []
prev_corner2 = []
cur_corner2 for i,st in enumerate(status):
if st==1:
prev_corner2.append(prev_corner[i])
cur_corner2.append(cur_corner[i])= np.array(prev_corner2)
prev_corner2 = np.array(cur_corner2)
cur_corner2 = cv2.estimateAffine2D(prev_corner2, cur_corner2)
T, inliers = T[:]
last_T
= T[0,2];
dx = T[1,2];
dy = np.arctan2(T[1,0], T[0,0])
da
prev_to_cur_transform.append([dx, dy, da])= cur[:]
prev = cur_gray[:]
prev_gray
= np.array(prev_to_cur_transform)
prev_to_cur_transform = np.cumsum(prev_to_cur_transform, axis=0)
trajectory = pd.DataFrame(trajectory)
trajectory = trajectory.rolling(window=30).mean()
smoothed_trajectory = smoothed_trajectory.bfill()
smoothed_trajectory = prev_to_cur_transform + \
new_prev_to_cur_transform - trajectory)
(smoothed_trajectory = np.array(new_prev_to_cur_transform)
new_prev_to_cur_transform
= np.zeros((2,3))
T = cv2.VideoCapture(fin)
cap = cv2.VideoWriter('/tmp/out.avi', cv2.VideoWriter_fourcc('P','I','M','1'),
out True)
fps, (w, h),
for k in range(N-1):
= cap.read()
status, cur 0,0] = np.cos(new_prev_to_cur_transform[k][2]);
T[0,1] = -np.sin(new_prev_to_cur_transform[k][2]);
T[1,0] = np.sin(new_prev_to_cur_transform[k][2]);
T[1,1] = np.cos(new_prev_to_cur_transform[k][2]);
T[0,2] = new_prev_to_cur_transform[k][0];
T[1,2] = new_prev_to_cur_transform[k][1];
T[= cv2.warpAffine(cur, T, (w,h));
cur2 ;
out.write(cur2)20); cv2.waitKey(
cv2.estimateAffine2D
çağrısı katı transformasyonu
hesaplayan bir çağrıdır, aynen H_from_points
gibi.
Üstteki kodu [1] üzerinde uygularsak stabilizasyon yapıldığını
göreceğiz. Sonuç [2]’de. C++ kodu vidstab.cpp
’de
bulunabilir. Derlemek için
g++ -c -O3 -Wall `pkg-config --cflags opencv` -D LINUX -fPIC vidstab.cpp
g++ -o vs.exe vs.o `pkg-config --libs opencv`
Canlı Zamanda (Real-Time) Stabilizasyon
[4]’ün tekniği toptan (batch) işleyen bir teknik, ortalama alınması, düzeltme yapılması için video’nun baştan sona işlenmesi, ve geriye dönülmesi gerekiyor. Düzeltme işlemini canlı olarak yapamaz mıyız?
Bu mümkün olmalı; yürüyen ortalama için [6] yazısına bakabiliriz; orada işlenen üstel ağırlıklı hareketli ortalama kullanılabilir. Bu ortalamanın özyineli (recursive) formu da vardır,
\[ z_t = \alpha g_t + (1-\alpha) z_{t-1}\]
ki \(\alpha\) kullanıcı tarafından
seçilen parametredir, en son verilere ne kadar ağırlık verileceğini
tanımlar. Algoritma şöyle olabilir: Stabilizasyon için her video karesi
işlenirken \(dx,dy,da\) farklarını
hesaplarız, bunların kümülatif toplamını da anlık hesaplarız (kolay). Bu
kümülatif \(x,y,a\)’yı üstteki tanımda
\(g_t\) olarak formüle veririz, en son
ortalama her zaman \(z_t\) içinde
olacaktır. Bu ortalamanın kumulatif olandan farkı, “sapması’’ her kare
üzerinde düzeltme amacı ile kullanılabilir. Bu kod
vsonline.py
içinde bulunabilir.
Kaynaklar
[1] Bayramlı, Veri 1 (Video), File
[2] Bayramlı, Veri 2 (Video), File
[3] Solem, Computer Vision with Python
[4] Nghia Ho, Simple Video Stabilization using OpenCV, http://nghiaho.com/?p=2093
[5] Bayramlı, OpenCV 3.0, https://burakbayramli.github.io/dersblog/sk/2017/03/opencv-30.html
[6] Bayramlı, Zaman Serileri ve Finans, ARIMA, ARCH, GARCH, Periyotlar, Yürüyen Ortalama
[7] Bayramlı, Kalman Filters and Homography: Utilizing the Matrix A https://arxiv.org/abs/1006.4910