Tali Basha (now Dekel) published a Multi-View Scene Flow (MVSF)
at CVPR 2010, and IJCV 2013 [home]
(with code) [author]
[tr].
We search for possibilites for user interaction.
Given n cameras Ci with frames Ii at t=0 and t=1, and an initial estimate of Z0 in camera C0, the objective is to dermine the true depth Z0 and the 3D-motion U,V,W. A coarse-to-fine pyramid with η=.9 is used to remove the nonlinearity in the data term, an outer loop to remove the nonlinearity from the perspective division, and an inner loop to remove the nonlinearity from the nonquadratic penalizer ψ(s2)=[s2+ε2].5. In the algorithm sketch below, only U is used for illustration, V/W/Z0 go analoguously:
U0 | for pyramid level "L"=14..0: | unknowns | helpers | images | maps | coordinates | |||||
U1 | U=0 or U=UL+1 | U | I_0t, I_0tXY | MapX_Id | |||||||
U1 | for outer iteration "O"=0..3: | ||||||||||
U2 | keep U fixed; rewarp Ii,0 and Ii,1 with current U; dU=0 | dU, UXY | occl0t | I_0tw, I_0twXY | MapX_I_0t | X, Xt, Xxyz, Xtxy, projxy0t,xyz, projxyt,uvw | |||||
U2 | for inner iteration "I"=0..10: | ||||||||||
U3 | keep dU fixed; build Mx=B with current dU; x=0 | dUXY | occl0t, div0wz | MapX_I_0tin | Zin, Xin, Xtin | ||||||
U3 | for SOR iteration "S"=0..40: | ||||||||||
U4 | keep ψ{flow,stereo1,stereo2} fixed; solve Mx=B | x | |||||||||
U3 | dU = x | ||||||||||
U2 | U += dU | ||||||||||
U1 | UL-1 = upsample(U) |
Due to the 3D to 2D projection nature of the computation, a lot of coordinates are juggled around, as are their derivatives. The list below shows the variables in the MVSF algorithm in a sequenced overview. All equations refer to the technical report [tr] which has the most details. There, BC{m,s1,s2} → BC{flow,st.init,stereo} refers to brightness constancy and S{m,s} → S{flow,stereo} to smoothness.
U1 | MapX_Id | = | Id(Ny,Nx) | |||||||
U2 | X | = | ComputeXY(Z0, MapX_Id) | (eq.2) | ||||||
Xt | = | X+U | , | Zt | = | Z0+W | (eq.3,4) | |||
MapX_I_0 | = | ComputeProjMaps(X, Z0) | (eq.5,16-17) | |||||||
MapX_I_t | = | ComputeProjMaps(Xt, Zt) | ||||||||
occl0 | = | ComputeOcclusion(X, Z0, MapX_I_0) | (eq.7, alg.1) | |||||||
occlt | = | ComputeOcclusion(Xt, Zt, MapX_I_t) | (eq.7) | |||||||
I_0w | = | Smooth(Remap(I_0, MapX_I_0)) | (eq.15) | |||||||
I_tw | = | Smooth(Remap(I_t, MapX_I_t)) | ||||||||
I_0wX | = | Derivative(I_0w) | ||||||||
I_twX | = | Derivative(I_tw) | ||||||||
Z0X | = | GradientFast(Z0) | ||||||||
UX | = | GradientFast(U) | ||||||||
Xxyz | = | Gradient3DPoints(Z0, Z0X, MapX_Id) | , | Xtxy | = | Xxy+UXY | ||||
projx0,xyz | = | GradientProjCoor(X, Z0, Xxyz, Z0X) | ||||||||
projxt,xyz | = | GradientProjCoor(Xt, Zt, Xtxyz, ZtX) | ||||||||
projxt,uvw | = | GradientProjCoor(Xt, Zt, 01, 01) | ||||||||
U3 | Zin | = | Z0+dZ0 | , | Ztin | = | Zin+W | |||
Xin | = | ComputeXY(Zin, MapX_Id) | , | Xtin | = | Xin+U+dU | ||||
MapX_I_0in | = | ComputeProjMaps(Xin, Zin) | ||||||||
MapX_I_tin | = | ComputeProjMaps(Xtin, Ztin) | ||||||||
dZ0X | = | GradientFast(dZ0) | ||||||||
dUX | = | GradientFast(dU) | ||||||||
div0wz | = | divergence(Z0X, dZ0X, UX, dUX) | (eq.12,13) | |||||||
M, B, occl0t | = | MB(Z0, dZ0, U, dU, div0wz, occl0t, MapX_I_0tin, I_0(X), I_0tw(X), projx0t,xy, projxt,uvw) | ||||||||
U4 | x | = | SOR(M, B) |
This is the system of linear equations to be solved and the best
chance to insert user input into the optimization. Green
entries mark optional parts (when solving for UVWZ
instead of UVW). Note that the full M is actually a
sparse (Nx*Ny*|uvwz|)2 matrix,
and the M here represents the diagonal and neighbour
locations. In the following, one 4-part for UVWZ is shown.
The SOR solver in the section after this shows each of the UVWZ
uses one row of M: U will use line0, V line1,
W line2, Z line3.
M | ⋅ | x | = | B | |||||||||||||
⌜ | au1 | av1 | aw1 | az1 | bu↑ | bu↓ | bu← | bu→ | ⌝ | ⌜ | dU | ⌝ | ⌜ | bu | ⌝ | ||
au2 | av2 | aw2 | az2 | bv↑ | bv↓ | bv← | bv→ | . | dV | = | bv | ||||||
au3 | av3 | aw3 | az3 | bw↑ | bw↓ | bw← | bw→ | dW | bw | ||||||||
⌞ | au4 | av4 | aw4 | az4 | bz↑ | bz↓ | bz← | bz→ | ⌟ | ⌞ | dZ | ⌟ | ⌞ | bz | ⌟ |
where:
M0..3 | |||||||||
au1 | = | α*μu | * [ ∑i=1..4(bfi) ] | , | bf2 | = | .5 * (div0[j,i] + div0[j,i-1]) | , where bf1,bf2,bf3,bf4 are right(1),left(2),down(3),up(4) (Ss or Sm, eq.12,13) | |
+ | ψflow | * [ itwU ]2 | , | ψflow | = | (12 + [it+(itwU*dU)+(itwV*dV)+(itwW*dW) - i0 + (itwZ-i0wZ)*dZ]2)-.5 (BCm, Δti=1+, eq.11) | |||
+ | ψstereo | * [ itwU - iltwU ]2 | , | ψstereo | = | (12 + [it+(itwU*dU)+(itwV*dV)+(itwW*dW)+(itwZ*dZ) - ilt-(iltwU*dU)-(iltwV*dV)-(iltwW*dW)-(iltwZ*dZ)]2)-.5 (BCs2, Δ^i, eq.11) | |||
+ | ψflowL | * [ iltwU ]2 | , | ψflowL | = | (12 + [ilt+(iltwU*dU)+(iltwV*dV)+(iltwW*dW)+(iltwZ*dZ) - il0]2)-.5 (BCm, Δti=0, eq.11) | |||
aw2 | = | ψflow | * [ itwW * itwV ] | , | itwU | = | (itX*projxt,u)+(itY*projyt,u) | , itX=(itwX*projyt,y - itwY*projyt,x)/det | , det=projxt,x*projyt,y - projxt,y*projyt,x (J, eq.22) |
+ | ψstereo | * [ (itwW - iltwW) * (itwV - iltwV) ] | , | iltwU | = | (iltX*projxlt,u)+(iltY*projylt,u) | , iltX=(iltwX*projylt,y - iltwY*projylt,x)/det | , det=projxlt,x*projylt,y - projxlt,y*projylt,x (J, eq.22) | |
+ | ψflowL | * [ iltwW * iltwV ] | , | it | = | Images_tw[1+], ilt=Images_tw[0] | , i0=Images_0w[1+], il0=Images_0[0] | , itwX=Images_tw_X[1+], iltwX=Images_tw_X[0] | |
az4 | = | αz*μz | * [ ∑i=1..4(bfzi) ] | , | i0wZ | = | (i0X*projx0,z)+(i0Y*projy0,z) (eq.18) | , i0X=(i0wX*projy0,y - i0wY*projy0,x)/det | , det=projx0,x*projy0,y - projx0,y*projy0,x (J, eq.22) |
+ | ψflow | * [ itwZ - i0wZ ]2 | , | itwZ | = | (itX*projxt,z)+(itY*projyt,z) (eq.20) | , itX=(itwX*projyt,y - itwY*projyt,x)/det | , det=projxt,x*projyt,y - projxt,y*projyt,x (J, eq.22) | |
+ | ψst.init | * [ i0wZ ]2 | , | ψst.init | = | (12 + [i0+(i0wZ*dZ) - il0]2)-.5 (BCs1, Δi=1+, eq.11) | |||
+ | ψstereo | * [ itwZ - iltwZ ]2 | |||||||
M4..7 | |||||||||
bu↑ | = | -α*μu | * bf4 | , | bf4 | = | .5 * (div0[j,i] + div0[j-1,i]) | , div0=(.0012 + [(UX+dUX)2+(UY+dUY)2+(VX+dVX)2+(VY+dVY)2])-.5 (Sm(uv), eq.13) | |
bw← | = | -αw*μw | * bfw2 | , | bfw2 | = | .5 * (divw[j,i] + divw[j,i-1]) | , divw=(.0012 + [μw * ((WX+dWX)2+(WY+dWY)2)])-.5 (Sm(w), eq.13) | |
bz→ | = | -αz*μz | * bfz1 | , | bfz1 | = | .5 * (divz[j,i] + divz[j,i+1]) | , divz=(.0012 + [μz * ((Z0X+dZ0X)2+(Z0Y+dZ0Y)2)])-.5 (Ss, eq.12) | |
B | |||||||||
bu | = | -α*μu | * [ ∑i=1..4(bfi*bui) - ∑i=1..4(bfi) * bu0 ] | , | bu0 | = | U[j,i], bu1=U[j,i+1], ... | ||
- | ψflow | * [ itwU * (it-i0) ] | |||||||
- | ψstereo | * [ (it-ilt) * (itwU-iltwU) ] | |||||||
- | ψflowL | * [ iltwU * (ilt-il0) ] | |||||||
bw | = | -αw*μw | * [ ∑i=1..4(bfwi*bwi) - ∑i=1..4(bfwi) * bw0 ] | , | bw0 | = | W[j,i], bw1=W[j,i+1], bw2=W[j,i-1], ... | ||
- | ψflow | * [ itwW * (it-i0) ] | |||||||
- | ψstereo | * [ (it-ilt) * (itwW-iltwW) ] | |||||||
- | ψflowL | * [ iltwW * (ilt-il0) ] | |||||||
bz |
= |
-αz*μz |
* [ ∑i=1..4(bfzi*bzi) - ∑i=1..4(bfzi) * bz0 ] | , |
bz0 | = |
Z0[j,i],
..., bz4=Z0[j-1,i] |
||
- |
ψflow | * [ (itwZ-i0wZ)
* (i0-it) ] |
|||||||
- | ψst.init | * [ (i0-il0) * i0wZ ] | |||||||
- | ψstereo | * [ (it-ilt) * (itwZ-iltwZ) ] | |||||||
- | ψflowL | * [ iltwZ * (ilt-il0) ] |
Neighbours: The order 1..4 is always right(1),left(2),down(3),up(4),
and the variables bf[i] and bu[i]
always run concurrent, e.g. where bf2 asks left, bu2
asks left too. In M4..8, that order is
reversed.
Observations:
Successive over-relaxation is a faster variant of the Jacobi
iteration, using the diagonal and upper+lower triangular
submatrixes. The row index i ∈ Nx*Ny*|uvwz|
has separate entries for uvwz.
rnew | = | (1.-ΩSOR)*rold + ΩSOR * (b - sum)/diag | = | (1.-ΩSOR)*rold + ΩSOR * (bu - sumU)/au1 |
diag | = | Mi,0..3 where uvwz=0123 | = | au1, av2, aw3, or az4 |
b | = | Bi | = | bu, bv, bw, or bz |
sumU | = | Mi,1*SORi+1 + Mi,2*SORi+2 + Mi,3*SORi+3 | = | av1*SORv + aw1*SORw + az1*SORz |
+ | Mi,4*SORi,1↑ + Mi,5*SORi,1↓ + Mi,6*SORi,1← + Mi,7*SORi,1→ | + | bu↑*SORu,1↑ + bu↓*SORu,1↓ + bu←*SORu,1← + bu→*SORu,1→ | |
sumZ | = | Mi,0*SORi-3 + Mi,1*SORi-2 + Mi,2*SORi-1 | = | au4*SORu + av4*SORv + aw4*SORw |
+ | Mi,4*SORi,1↑ + Mi,5*SORi,1↓ + Mi,6*SORi,1← + Mi,7*SORi,1→ | + | bz↑*SORz,1↑ + bz↓*SORz,1↓ + bz←*SORz,1← + bz→*SORz,1→ |
For Jacobi iterations, set ΩSOR=1 so that rnew is not dependent on rold.
EOF (Sep:2014)