Multi-view Scene Flow Estimation: A View Centered Variational Approach

Abstract

We present a novel method for recovering the 3D structure and scene flow from calibrated multi-view sequences. We propose a 3D point cloud parametrization of the 3D structure and scene flow that allows us to directly estimate the desired unknowns. A unified global energy functional is proposed to incorporate the information from the available sequences and simultaneously recover both depth and scene flow. The functional enforces multi-view geometric consistency and imposes brightness constancy and piecewise smoothness assumptions directly on the 3D unknowns. It inherently handles the challenges of discontinuities, occlusions, and large displacements. The main contribution of this work is the fusion of a 3D representation and an advanced variational framework that directly uses the available multi-view information. This formulation allows us to advantageously bind the 3D unknowns in time and space. Different from optical flow and disparity, the proposed method results in a nonlinear mapping between the images’ coordinates, thus giving rise to additional challenges in the optimization process. Our experiments on real and synthetic data demonstrate that the proposed method successfully recovers the 3D structure and scene flow despite the complicated nonconvex optimization problem.

Appendices

Appendix A: Mapping Between Images

Our 3D parameterization in the presented framework introduces a nonlinear transformation of the 3D unknowns, Z and V, to each of the image’s plane. A notable challenge in the minimization of the proposed functional arises from the nontrivial mapping of the images’ coordinates to the reference camera coordinate system.

Using our parametrization, each pixel in the reference camera, (x,y), and its corresponding depth, Z(x,y), specify a 3D point, P (see Eq. (5)). It follows that projecting P onto the ith camera maps (x,y,Z(x,y)) to the point p _i =(x _i ,y _i )^T. That is,

$$ \textbf {p}_i = \mathit{Proj}\bigl(\textbf {P},M^i\bigr) = f^i\bigl(x,y,Z(x,y)\bigr), $$

(17)

where f ⁱ is the mapping to the corresponding ith image. More precisely, f ⁱ is given by substituting Eq. (5) into Eq. (1). For example, the component x _i is given by:

$$ x_i = {a\cdot Z + b \over c\cdot Z + d}~. $$

(18)

The coefficients a,b,c and d depend on the reference camera coordinates, (x,y):

(19)

where M ⁱ is the 3×4 projection matrix of the ith camera (subscripts denote the row and column indices). The expression for y _i is equivalently computed.

Similarly, at time step t+1, projecting $\widehat{\textbf {P}}= \textbf {P}+ \textbf {V}$ maps (x,y,Z(x,y),V(x,y)) to $\widehat{\textbf {p}}_{i}$, denoted by a mapping, $\widehat{f}^{i}$:

$$ \widehat{\textbf {p}}_i = \mathit{Proj}\bigl(\,\widehat{\textbf {P}},M^i\bigr) = \widehat{f}^i\bigl(x,y,Z(x,y),\textbf {V}(x,y) \bigr). $$

(20)

Analogously to Eq. (18), the component $\widehat{x}_{i}$ is given by:

$$ \widehat{x}_i = {a\cdot Z + M^i_{{11}}\cdot u + M^i_{{12}}\cdot v + M^i_{{13}}\cdot w+ b \over c\cdot Z + M^i_{{31}} \cdot u + M^i_{{32}}\cdot v + M^i_{{33}} \cdot w + d}, $$

(21)

where the coefficients a,b,c and d are defined in Eq. (19).

Appendix B: Image Derivatives with Respect to the 3D Unknowns

A first step toward the numerical solution of the resulting Euler-Lagrange equations (Eq. (12) or Eq. (13)) requires computing the derivatives of the intensity functions with respect to the 3D unknowns. To produce the final expressions for these derivatives, the nonlinear relation between the 3D unknowns and the image plane has to be carefully considered (see Appendix A). This appendix shows how these computations are performed. The mathematical analysis is preformed in the continuous domain. Thus, the frames as well as the 3D unknowns are regarded as continuous functions. Finally, the resulting equations are discretized by using standard approximations for the derivatives.

For simplicity, given a time step, t, we use the intensity functions $I_{i}^{t}$ and $I_{i}^{t+1}$ to abbreviate I _i (p _i ,t) and $I_{i}(\widehat{\textbf {p}}_{i}, t+1)$, respectively. We next elaborate on the computation of derivatives of $I_{i}^{t}$ and $I_{i}^{t+1}$ with respect to Z and u, denoted by $\partial_{Z}I_{i}^{t}$, $\partial_{Z}I_{i}^{t+1}$ and $\partial_{u}I_{i}^{t+1}$ (the other derivatives with respect to v and w are similarly computed).

$I_{i}^{t}$ can be regarded as a function of the reference image coordinates, (x,y), and the corresponding depth, Z(x,y), by considering a composition of two functions: the ith intensity function and the mapping transformation, defined in Appendix A. That is,

$$ I_i^t\bigl(x,y,Z(x,y)\bigr) = I_i(f^i\bigl(x,y,Z(x,y),t\bigr). $$

(22)

Similarly, $I_{i}^{t+1}$ can be regarded as a function of (x,y,Z(x,y)) and V. That is,

$$ I_i^{t+1}(x,y,Z,\textbf {V})= I_i \bigl(\widehat{f}^i(x,y,Z,\textbf {V}),t+1\bigr). $$

(23)

Considering Eqs. (17)–(20), the chain rule is applied for computing the partial derivatives:

(24)

(25)

(26)

The derivatives $\partial_{Z}\textbf {p}_{i}^{T} = ( \partial_{Z}x_{i}, \partial_{Z}y_{i})^{T}$ are directly computed from Eqs. (18)–(19).

To compute the derivative of $I_{i}^{t}$ with respect to p _i , $(\nabla I_{i}^{t})^{T}$, we use a warping approach. As discussed in Appendix A, a nonlinear mapping relates each of the image’s plane to the reference camera. By warping $I_{i}^{t}$ toward the reference image using the estimated Z, the values of $I_{i}^{t}$ can be directly related to the reference image values, $I_{0}^{t}$. Specifically, the required derivatives, $\nabla I_{i}^{t}$ are then computed using the warped image. Let $I_{i,w}^{t}$ be the warped image of $I_{i}^{t}$. That is,

(27)

The warped image gradient is related to the original image by:

(28)

where J is the Jacobian matrix of the change of coordinates, (x _i ,y _i )→(x,y). Therefore, the original image derivatives are obtained by multiplying Eq. (28) by J ⁻¹, leading to:

(29)

The Jacobian matrix, J, is obtained by computing the derivatives of p _i with respect to x and y. In particular, J involves the derivatives of Z(x,y), namely ∂ _x Z and ∂ _y Z. Following the explanation above, $\nabla I_{i,w}^{t+1}$ is similarly computed. In this case, the Jacobian matrix, J, additionally involves the derivatives of u,v and w with respect to the reference camera coordinates.

Appendix C: Linearizion

This appendix describes the linearizion process of the resulting Euler-Lagrange equations and the numerical approximations used. At each pyramid level, a linear system of equations is obtained and small increments in the 3D unknowns, dZ, and d V, are estimated. The total solution, Z+dZ, and V+d V, is then used to initialize the next finer level (see Sect. 2.3.2).

Considering equations (12)–(13), there are two sources of nonlinearity:

• 1. nonlinearized data term;

• 2. nonquadratic cost function Ψ.

Following the numerical approach suggested by Brox et al. (2004), two nested fixed point iterations are used at each pyramid level to remove the nonlinearity.

The outer iteration is responsible for removing the nonlinearity resulting from the nonlinear data term, using fixed point iteration on Z and V. Let k be the outer index iteration. The solution at the (k+1)th iteration is decomposed of the previous solution and small, unknown increments. That is, Z ^k+1=Z ^k+dZ ^k and V ^k+1=V ^k+d V ^k, where d V ^k=(du ^k,dv ^k,dw ^k)^T.

The first step toward linearizion is approximating the nonlinear expression given in Eq. (11) using first order Taylor expansion. We use $\Delta_{i}^{k},~\widehat{\Delta}_{i}^{k}$ and $\Delta_{i}^{t, k}$ to denote the expressions given in Eq. (11) using the fixed values Z ^k and V ^k. That is,

(30)

where $\textbf {p}_{i}^{k} = \mathit{Proj}(\textbf {P}^{k}, M^{i})$ and P ^k is given by placing Z ^k in Eq. (5). The expressions for $\widehat{\textbf {p}}_{i}^{k}$ and $\widehat{\textbf {P}}^{k}$ are analogously given. Using these notations, the first order Taylor expansions for these expressions are given by:

(31)

Equation (31) is computed by using the first order Taylor expansion for the following expressions:

(32)

(33)

where P ^k+1=P ^k+d P ^k is given by placing Z ^k+dZ ^k (Eq. (5)). Similarly, $\widehat{\textbf {P}}^{k+1}=\widehat{\textbf {P}}^{k} + d\widehat{\textbf {P}}^{k}$ where $d\widehat{\textbf {P}}^{k}= dZ^{k} + d\textbf {V}^{k}$. The computation of the image derivatives with respect to the 3D unknowns is detailed in Appendix A.

Therefore, deriving the associated Euler-Lagrange equations with respect to the unknown increments dZ ^k and du ^k results in:

(34)

(35)

The dependency of the above two equations in the increments, dZ ^k and du ^k, is obtained by substituting Eq. (31) into $\Delta_{i}^{k+1},\Delta_{i}^{t, k+1}$, and, $\widehat{\Delta}_{i}^{k+1}$. The equations for dv ^k and dw ^k are similar to Eq. (35).

Applying the above approximations (Eq. (31)), the resulting Euler-Lagrange equations are a nonlinear system of equations in the unknowns dZ ^k and d V ^k. The remaining nonlinearity is originated by Ψ′. Therefore, an additional fixed point iterations loop for Ψ′ expressions is preformed. Finally, after standard discretization of the derivatives, a linear system of equations is introduced. The solution is obtained by applying the successive overrelaxation (SOR) method.