diff options
author | Dan Zhu <zxdan@google.com> | 2019-07-23 18:16:37 +0000 |
---|---|---|
committer | Gerrit Code Review <noreply-gerritcodereview@google.com> | 2019-07-23 18:16:37 +0000 |
commit | e7548da489bed380fa6f32308bc9531f9dc7d917 (patch) | |
tree | 0dbf835c4c3c6954d9cdd9ab60d550745533ca62 | |
parent | 5db8ff42db397f864f7707d7a1890fd70cd8fff6 (diff) | |
parent | c7093c18b067e8690bde47558cefe32284c84db3 (diff) | |
download | libvpx-e7548da489bed380fa6f32308bc9531f9dc7d917.tar libvpx-e7548da489bed380fa6f32308bc9531f9dc7d917.tar.gz libvpx-e7548da489bed380fa6f32308bc9531f9dc7d917.tar.bz2 libvpx-e7548da489bed380fa6f32308bc9531f9dc7d917.zip |
Merge "Add Horn & Schunck Estimator"
-rw-r--r-- | tools/3D-Reconstruction/MotionEST/HornSchunck.py | 205 | ||||
-rw-r--r-- | tools/3D-Reconstruction/MotionEST/MotionEST.py | 16 |
2 files changed, 215 insertions, 6 deletions
diff --git a/tools/3D-Reconstruction/MotionEST/HornSchunck.py b/tools/3D-Reconstruction/MotionEST/HornSchunck.py new file mode 100644 index 000000000..0bf431cf6 --- /dev/null +++ b/tools/3D-Reconstruction/MotionEST/HornSchunck.py @@ -0,0 +1,205 @@ +#!/usr/bin/env python +# coding: utf-8 +import numpy as np +import numpy.linalg as LA +from scipy.ndimage.filters import gaussian_filter +from scipy.sparse import csc_matrix +from scipy.sparse.linalg import inv +from MotionEST import MotionEST +"""Horn & Schunck Model""" + + +class HornSchunck(MotionEST): + """ + constructor: + cur_f: current frame + ref_f: reference frame + blk_sz: block size + alpha: smooth constrain weight + sigma: gaussian blur parameter + """ + + def __init__(self, cur_f, ref_f, blk_sz, alpha, sigma, max_iter=100): + super(HornSchunck, self).__init__(cur_f, ref_f, blk_sz) + self.cur_I, self.ref_I = self.getIntensity() + #perform gaussian blur to smooth the intensity + self.cur_I = gaussian_filter(self.cur_I, sigma=sigma) + self.ref_I = gaussian_filter(self.ref_I, sigma=sigma) + self.alpha = alpha + self.max_iter = max_iter + self.Ix, self.Iy, self.It = self.intensityDiff() + + """ + Build Frame Intensity + """ + + def getIntensity(self): + cur_I = np.zeros((self.num_row, self.num_col)) + ref_I = np.zeros((self.num_row, self.num_col)) + #use average intensity as block's intensity + for i in xrange(self.num_row): + for j in xrange(self.num_col): + r = i * self.blk_sz + c = j * self.blk_sz + cur_I[i, j] = np.mean(self.cur_yuv[r:r + self.blk_sz, c:c + self.blk_sz, + 0]) + ref_I[i, j] = np.mean(self.ref_yuv[r:r + self.blk_sz, c:c + self.blk_sz, + 0]) + return cur_I, ref_I + + """ + Get First Order Derivative + """ + + def intensityDiff(self): + Ix = np.zeros((self.num_row, self.num_col)) + Iy = np.zeros((self.num_row, self.num_col)) + It = np.zeros((self.num_row, self.num_col)) + sz = self.blk_sz + for i in xrange(self.num_row - 1): + for j in xrange(self.num_col - 1): + """ + Ix: + (i ,j) <--- (i ,j+1) + (i+1,j) <--- (i+1,j+1) + """ + count = 0 + for r, c in {(i, j + 1), (i + 1, j + 1)}: + if 0 <= r < self.num_row and 0 < c < self.num_col: + Ix[i, j] += ( + self.cur_I[r, c] - self.cur_I[r, c - 1] + self.ref_I[r, c] - + self.ref_I[r, c - 1]) + count += 2 + Ix[i, j] /= count + """ + Iy: + (i ,j) (i ,j+1) + ^ ^ + | | + (i+1,j) (i+1,j+1) + """ + count = 0 + for r, c in {(i + 1, j), (i + 1, j + 1)}: + if 0 < r < self.num_row and 0 <= c < self.num_col: + Iy[i, j] += ( + self.cur_I[r, c] - self.cur_I[r - 1, c] + self.ref_I[r, c] - + self.ref_I[r - 1, c]) + count += 2 + Iy[i, j] /= count + count = 0 + #It: + for r in xrange(i, i + 2): + for c in xrange(j, j + 2): + if 0 <= r < self.num_row and 0 <= c < self.num_col: + It[i, j] += (self.ref_I[r, c] - self.cur_I[r, c]) + count += 1 + It[i, j] /= count + return Ix, Iy, It + + """ + Get weighted average of neighbor motion vectors + for evaluation of laplacian + """ + + def averageMV(self): + avg = np.zeros((self.num_row, self.num_col, 2)) + """ + 1/12 --- 1/6 --- 1/12 + | | | + 1/6 --- -1/8 --- 1/6 + | | | + 1/12 --- 1/6 --- 1/12 + """ + for i in xrange(self.num_row): + for j in xrange(self.num_col): + for r, c in {(-1, 0), (1, 0), (0, -1), (0, 1)}: + if 0 <= i + r < self.num_row and 0 <= j + c < self.num_col: + avg[i, j] += self.mf[i + r, j + c] / 6.0 + for r, c in {(-1, -1), (-1, 1), (1, -1), (1, 1)}: + if 0 <= i + r < self.num_row and 0 <= j + c < self.num_col: + avg[i, j] += self.mf[i + r, j + c] / 12.0 + return avg + + def est(self): + count = 0 + """ + u_{n+1} = ~u_n - Ix(Ix.~u_n+Iy.~v+It)/(IxIx+IyIy+alpha^2) + v_{n+1} = ~v_n - Iy(Ix.~u_n+Iy.~v+It)/(IxIx+IyIy+alpha^2) + """ + denom = self.alpha**2 + np.power(self.Ix, 2) + np.power(self.Iy, 2) + while count < self.max_iter: + avg = self.averageMV() + self.mf[:, :, 1] = avg[:, :, 1] - self.Ix * ( + self.Ix * avg[:, :, 1] + self.Iy * avg[:, :, 0] + self.It) / denom + self.mf[:, :, 0] = avg[:, :, 0] - self.Iy * ( + self.Ix * avg[:, :, 1] + self.Iy * avg[:, :, 0] + self.It) / denom + count += 1 + self.mf *= self.blk_sz + + def est_mat(self): + row_idx = [] + col_idx = [] + data = [] + + N = 2 * self.num_row * self.num_col + b = np.zeros((N, 1)) + for i in xrange(self.num_row): + for j in xrange(self.num_col): + """(IxIx+alpha^2)u+IxIy.v-alpha^2~u IxIy.u+(IyIy+alpha^2)v-alpha^2~v + """ + u_idx = i * 2 * self.num_col + 2 * j + v_idx = u_idx + 1 + b[u_idx, 0] = -self.Ix[i, j] * self.It[i, j] + b[v_idx, 0] = -self.Iy[i, j] * self.It[i, j] + #u: (IxIx+alpha^2)u + row_idx.append(u_idx) + col_idx.append(u_idx) + data.append(self.Ix[i, j] * self.Ix[i, j] + self.alpha**2) + #IxIy.v + row_idx.append(u_idx) + col_idx.append(v_idx) + data.append(self.Ix[i, j] * self.Iy[i, j]) + + #v: IxIy.u + row_idx.append(v_idx) + col_idx.append(u_idx) + data.append(self.Ix[i, j] * self.Iy[i, j]) + #(IyIy+alpha^2)v + row_idx.append(v_idx) + col_idx.append(v_idx) + data.append(self.Iy[i, j] * self.Iy[i, j] + self.alpha**2) + + #-alpha^2~u + #-alpha^2~v + for r, c in {(-1, 0), (1, 0), (0, -1), (0, 1)}: + if 0 <= i + r < self.num_row and 0 <= j + c < self.num_col: + u_nb = (i + r) * 2 * self.num_col + 2 * (j + c) + v_nb = u_nb + 1 + + row_idx.append(u_idx) + col_idx.append(u_nb) + data.append(-1 * self.alpha**2 / 6.0) + + row_idx.append(v_idx) + col_idx.append(v_nb) + data.append(-1 * self.alpha**2 / 6.0) + for r, c in {(-1, -1), (-1, 1), (1, -1), (1, 1)}: + if 0 <= i + r < self.num_row and 0 <= j + c < self.num_col: + u_nb = (i + r) * 2 * self.num_col + 2 * (j + c) + v_nb = u_nb + 1 + + row_idx.append(u_idx) + col_idx.append(u_nb) + data.append(-1 * self.alpha**2 / 12.0) + + row_idx.append(v_idx) + col_idx.append(v_nb) + data.append(-1 * self.alpha**2 / 12.0) + M = csc_matrix((data, (row_idx, col_idx)), shape=(N, N)) + M_inv = inv(M) + uv = M_inv.dot(b) + + for i in xrange(self.num_row): + for j in xrange(self.num_col): + self.mf[i, j, 0] = uv[i * 2 * self.num_col + 2 * j + 1, 0] * self.blk_sz + self.mf[i, j, 1] = uv[i * 2 * self.num_col + 2 * j, 0] * self.blk_sz diff --git a/tools/3D-Reconstruction/MotionEST/MotionEST.py b/tools/3D-Reconstruction/MotionEST/MotionEST.py index 0e04bdd40..68cf7e743 100644 --- a/tools/3D-Reconstruction/MotionEST/MotionEST.py +++ b/tools/3D-Reconstruction/MotionEST/MotionEST.py @@ -53,24 +53,28 @@ class MotionEST(object): h = min(self.blk_sz, self.height - cur_y) w = min(self.blk_sz, self.width - cur_x) cur_blk = self.cur_yuv[cur_y:cur_y + h, cur_x:cur_x + w, :] - ref_x = cur_x + mv[1] - ref_y = cur_y + mv[0] - if 0 <= ref_x < self.width and 0 <= ref_y < self.height: + ref_x = int(cur_x + mv[1]) + ref_y = int(cur_y + mv[0]) + if 0 <= ref_x < self.width - w and 0 <= ref_y < self.height - h: ref_blk = self.ref_yuv[ref_y:ref_y + h, ref_x:ref_x + w, :] else: ref_blk = np.zeros((h, w, 3)) - return self.metric(cur_blk, ref_blk) + return metric(cur_blk, ref_blk) """ distortion of motion field """ - def distortion(self, metric=MSE): + def distortion(self, mask=None, metric=MSE): loss = 0 + count = 0 for i in xrange(self.num_row): for j in xrange(self.num_col): + if not mask is None and mask[i, j]: + continue loss += self.dist(i, j, self.mf[i, j], metric) - return loss / self.num_row / self.num_col + count += 1 + return loss / count """ evaluation |