• 沒有找到結果。

# fundamental matrix

N/A
N/A
Protected

Share "fundamental matrix"

Copied!
20
0
0

(1)

### Structure from motion

Digital Visual Effectsg Yung-Yu Chuang

with slides by Richard Szeliski, Steve Seitz, Zhengyou Zhang and Marc Pollefyes

### Outline

• Epipolar geometry and fundamental matrix

S f i

• Structure from motion

• Factorization method

• Applications

• Applications

### The epipolar geometry

epipolar geometry demo

C C’ x x’ and X are coplanar C,C ,x,x and X are coplanar

(2)

### The epipolar geometry

What if only C C’ x are known?

What if only C,C ,x are known?

### The epipolar geometry

All points on  project on l and l’

All points on  project on l and l

### The epipolar geometry

Family of planes  and lines l and l’ intersect at e Family of planes  and lines l and l intersect at e and e’

### The epipolar geometry

epipolar pole

= intersection of baseline with image plane

epipolar geometry demo

= intersection of baseline with image plane

= projection of projection center in other image

epipolar plane = plane containing baseline epipolar plane plane containing baseline

epipolar line = intersection of epipolar plane with image

(3)

C C’

p R p’

C

T=C’-C

### TRp'p  

Two reference frames are related via the extrinsic parameters

###  

0

0 z y

T T

T T T

Multiply both sides by pT

T



 



0 0

x y

x z

T T

T T T

T

T

T

essential matrix

### Epp pp

Let M and M’ be the intrinsic matrices, then

1

1

1

1

1

### x

fundamental matrix

### The fundamental matrix F

• The fundamental matrix is the algebraic representation of epipolar geometry representation of epipolar geometry

Th f d t l t i ti fi th diti

• The fundamental matrix satisfies the condition that for any pair of corresponding points x↔x’

i th t i

in the two images

T

T

(4)

### The fundamental matrix F

F is the unique 3x3 rank 2 matrix that satisfies xTFx’=0 for all x↔x’

for all x↔x

1. Transpose: if F is fundamental matrix for (P,P’), then FT is fundamental matrix for (P’,P)

2 Epipolar lines: l F ’ & l’ FT 2. Epipolar lines: l=Fx’ & l’=FTx

3. Epipoles: on all epipolar lines, thus eTFx’=0, x’

eTF=0, similarly Fe’=0

e F 0, similarly Fe 0

4. F has 7 d.o.f. , i.e. 3x3-1(homogeneous)-1(rank2)

5. F is a correlation, projective mapping from a point x to a line l=Fx’ (not a proper correlation, i.e. not invertible)

### The fundamental matrix F

• It can be used for

– Simplifies matching

– Allows to detect wrong matchesAllows to detect wrong matches

### Estimation of F — 8-point algorithm

• The fundamental matrix F is defined by

### Fxx

for any pair of matches x and x’ in two images.

f11 f12 f13

• Let x=(u,v,1)T and x’=(u’,v’,1)T,

21 22 23

13 12 11

f f f

f f f

f f f F



f31 f32 f33

each match gives a linear equation

0 '

' '

' '

' f11uv f12uf13vu f21vv f22vf23u f31v f32 f33 uu

### 8-point algorithm

11

f f

1

´

´

´

´

´

´ 13

12 1

1 1 1 1 1 1 1 1 1 1

1

f f f v

u v v v u v u v u u u

1 0

´

´

´

´

´

´

22 21 2

2 2 2 2 2 2 2 2 2 2 2

1 1 1 1 1 1 1 1 1 1 1 1

f v f

u v v v u v u v u u u

1

´

´

´

´

´

´

31 23

f

v f u v v v u v u v u u

un n n n n n n n n n n n

33 32

 f f

33

 f

• In reality, instead of solving , we seek f to minimize subj. . Find the vector Af  0

Af j f 1

corresponding to the least singular value.

(5)

### 8-point algorithm

• To enforce that F is of rank 2, F is replaced by F’ that minimizes subject to F F' det F' 0 F that minimizes subject to . FF' det F' 0

• It is achieved by SVD. Let , where y FUΣV, let



 1

0 0

0 0

Σ 

 



0 0

0 0 Σ'

1

 , let



 



3 2

0 0

0 0

Σ



 



0 0 0

0 0

Σ' 2

then is the solution. F' UΣ'V

### 8-point algorithm

% Build the constraint matrix

A = [x2(1,:)‘.*x1(1,:)' x2(1,:)'.*x1(2,:)' x2(1,:)' ...[ ( , ) ( , ) ( , ) ( , ) ( , ) x2(2,:)'.*x1(1,:)' x2(2,:)'.*x1(2,:)' x2(2,:)' ...

x1(1,:)' x1(2,:)' ones(npts,1) ];

[U,D,V] = svd(A);

% Extract fundamental matrix from the column of V

% corresponding to the smallest singular value.p g g F = reshape(V(:,9),3,3)';

% E f k2 t i t

% Enforce rank2 constraint [U,D,V] = svd(F);

F = U*diag([D(1 1) D(2 2) 0])*V';

F = U diag([D(1,1) D(2,2) 0]) V ;

### 8-point algorithm

• Pros: it is linear, easy to implement and fast

C ibl i

• Cons: susceptible to noise

### Problem with 8-point algorithm

11

f f

1

´

´

´

´

´

´ 13

12 1

1 1 1 1 1 1 1 1 1 1

1

f f f v

u v v v u v u v u u u

1 0

´

´

´

´

´

´

22 21 2

2 2 2 2 2 2 2 2 2 2 2

1 1 1 1 1 1 1 1 1 1 1 1

f v f

u v v v u v u v u u u

10000 10000 100 10000 10000 100 100 ~100 1

1

´

´

´

´

´

´

31 23

f

v f u v v v u v u v u u

un n n n n n n n n n n n

~10000 ~10000 ~100 ~10000 ~10000 ~100 ~100 ~100 1

!

Orders of magnitude difference

between column of data matrix 33

32

 f f

! between column of data matrix

 least-squares yields poor results

33

 f

(6)

### Normalized 8-point algorithm

1. Transform input by , 2 C ll 8 i b i

i

i Tx

xˆ  xˆ'iTx'i ˆ'

ˆ ˆ

2. Call 8-point on to obtain 3.

i i x x ˆˆ , T

F T F'Τ ˆ

F

1

### Normalized 8-point algorithm

normalized least squares yields good results

T f i [ 1 1] [ 1 1]

(700,500)

(0,500) (-1,1) (1,1)

 2

Transform image to ~[-1,1]x[-1,1]

( , )

( , ) ( , ) ( , )

500 1

2 1 700 0

2

(0,0)

1

(0,0) (700,0) (-1,-1) (1,-1)

### Normalized 8-point algorithm

[x1, T1] = normalise2dpts(x1);

[x2, T2] = normalise2dpts(x2);

A = [x2(1,:)‘.*x1(1,:)' x2(1,:)'.*x1(2,:)' x2(1,:)' ...

x2(2,:)'.*x1(1,:)' x2(2,:)'.*x1(2,:)' x2(2,:)' ...

[ , ] p ( );

x1(1,:)' x1(2,:)' ones(npts,1) ];

[U D V] svd(A);

[U,D,V] = svd(A);

F = reshape(V(:,9),3,3)';

F reshape(V(:,9),3,3) ; [U,D,V] = svd(F);

F = U*diag([D(1,1) D(2,2) 0])*V';

% Denormalise F = T2'*F*T1;

### Normalization

function [newpts, T] = normalise2dpts(pts) c = mean(pts(1:2,:)')'; % Centroid

newp(1,:) = pts(1,:)-c(1); % Shift origin to centroid.p( , ) p ( , ) ( ); g newp(2,:) = pts(2,:)-c(2);

meandist = mean(sqrt(newp(1,:).^2 + newp(2,:).^2));

scale = sqrt(2)/meandist;

T = [scale 0 -scale*c(1) 0 scale -scale*c(2) 0 0 1 ];

t T* t newpts = T*pts;

(7)

### RANSAC

repeat

select minimal sample (8 matches) compute solution(s) for F

determine inliers

until (#inliers,#samples)>95% or too many times ( , p ) y compute F based on all inliers

(8)

### Structure from motion

Unknown Unknown Unknown Unknown camera camera viewpoints viewpoints

structure for motion: automatic recovery of camera motion and scene structure from two or more images. It is a self calibration technique and called automatic camera tracking calibration technique and called automatic camera tracking or matchmoving.

### Applications

• For computer vision, multiple-view shape reconstruction novel view synthesis and reconstruction, novel view synthesis and autonomous vehicle navigation.

F fil d ti l i ti f CGI

• For film production, seamless insertion of CGI into live-action backgrounds

### Matchmove

example #1 example #2 example #3 example #4

(9)

### Structure from motion

2D feature

tracking 3D estimation optimization (bundle adjust)

geometry fitting

SFM pipeline SFM pipeline

### Structure from motion

• Step 1: Track Features

Detect good features Shi & Tomasi SIFT – Detect good features, Shi & Tomasi, SIFT – Find correspondences between frames

• Lucas & Kanade-style motion estimationy

• window-based correlation

• SIFT matching

### KLT tracking

http://www ces clemson edu/~stb/klt/

http://www.ces.clemson.edu/ stb/klt/

### Structure from Motion

• Step 2: Estimate Motion and Structure

Si lifi d j ti d l [T i 92]

– Simplified projection model, e.g., [Tomasi 92]

– 2 or 3 views at a time [Hartley 00]

(10)

### Structure from Motion

• Step 3: Refine estimates

“B dl dj t t” i h t t – “Bundle adjustment” in photogrammetry – Other iterative methods

### Structure from Motion

• Step 4: Recover surfaces (image-based triangulation silhouettes stereo ) triangulation, silhouettes, stereo…)

Good mesh Good mesh

(11)

### Notations

• n 3D points are seen in m views

( 1) 2D i i

• q=(u,v,1): 2D image point

• p=(x,y,z,1): 3D scene point

• : projection matrix

• : projection function

• : projection function

• qij is the projection of the i-th point on image j

 j ti d th f

• ij projective depth of qij

ij

j

i

(x y z)(x/z y/z)

(x,y,z)(x/z,y/z)

ijz

### Structure from motion

• Estimate and to minimize

j

i

) );

( ( log )

, , , , , (

1 1 1

1 j i ij

m

j n

i ij n

m p p w P Π p q

Π

Π

### 

 

j

otherwise

j in view visible

is if 0

1 i

ij

w p





otherwise

j 0

• Assume isotropic Gaussian noise, it is reduced top ,

2 1

1, , , , , ) ( )

( j i ij

m n

wij Π p q p

p Π

Π

1 1 1

1, , , , , ) ( )

( j i ij

j i

ij n

m p p w Π p q

Π

Π

### 

• Start from a simpler projection model

• Start from a simpler projection model

### Orthographic projection

• Special case of perspective projection

Di t f th COP t th PP i i fi it – Distance from the COP to the PP is infinite

Image World

– Also called “parallel projection”: (x, y, z) → (x, y)

### SFM under orthographic projection

2D image Orthographic projection

incorporating 3D rotation 3D scene

image offset point incorporating 3D rotation 3D scene

point

offset

### tΠtΠpq  

211 23 31 21 2 23 31 21

• Trick

– Choose scene origin to be centroid of 3D pointsg p – Choose image origins to be centroid of 2D points – Allows us to drop the camera translation:Allows us to drop the camera translation:

(12)

###    

projection of n features in one image:

n 3 3

n 2

2

n

### 

1 2 n

2

1 q q p p p

q

1 1 12

11

q q qn Π

projection of n features in m images

32

1 2 2

22 21

1 12

11

n n

n

p p

Π p q

q

q

n 3 3

n 2m 2m

2 1

qm qm qmn Πm

3 n 2m

2m

measurement

motion

### S

shape

Key Observation: rank(W) <= 3

3 3 2

2

### W  MS

known solve for

• Factorization Technique

n 3 3 m 2 n

2m

Factorization Technique

– W is at most rank 3 (assuming no noise)

– We can use singular value decomposition to factor W:

33 2

2

### M ' S '

We can use singular value decomposition to factor W:

n 33 m n 2

2m

– S’ differs from S by a linear transformation A:

### W  

1

– Solve for A by enforcing metric constraints on M

### Metric constraints

• Orthographic Camera

R f  th l

  T 01 10

– Rows of  are orthonormal:

• Enforcing “Metric” Constraints

  0 1

– Compute A such that rows of M have these properties

### M 'AMM ' 

Trick (not in original Tomasi/Kanade paper, but in followup work)

• Constraints are linear in AAT:

T T T

T A A G where G AA

' ' ' '

1 0

0 1

• Solve for G first by writing equations for every iin M

• Then G = AATby SVD (since U = V)

0 1

Then G AA by SVD (since U V)

### W

n 2m 33 n 2m n

2m

• SVD gives this solutiong

– Provides optimal rank 3 approximation W’ of W

n 2m n 2m n 2m

###  E

• Approach

– Estimate W’, then use noise-free factorization of W’

as before

– Result minimizes the SSD between positions of image features and projection of the reconstruction

(13)

### ResultsExtensions to factorization methods

• Projective projection Wi h i i d

• With missing data

• Projective projection with missing data

### Levenberg-Marquardt method

• LM can be thought of as a combination of steepest descent and the Newton method steepest descent and the Newton method.

When the current solution is far from the correct one the algorithm behaves like a correct one, the algorithm behaves like a

steepest descent method: slow, but guaranteed to converge When the current solution is close to converge. When the current solution is close to the correct solution, it becomes a Newton’s method

method.

(14)

find try to ,

ts measuremen of

set a

Given x

Here minimal is

distance squared

that the so

vector parameter

best

the p

T

).

ˆ ( with ˆ,

Here, minimal.

is distance

squared

p x

x

x  f

### Levenberg-Marquardt method

• μ=0 → Newton’s method

d h d

• μ→∞ → steepest descent method

• Strategy for choosing μ

– If error is not reduced, keep trying larger μ until it does

– If error is reduced, accept it and reduce μ for the next iteration

• Bundle adjustment (BA) is a technique for simultaneously refining the 3D structure and simultaneously refining the 3D structure and camera parameters

It i bl f bt i i ti l

• It is capable of obtaining an optimal

reconstruction under certain assumptions on

i d l F G i

image error models. For zero-mean Gaussian image errors, BA is the maximum likelihood

ti t estimator.

(15)

• n 3D points are seen in m views

i h j i f h i h i i j

• xij is the projection of the i-th point on image j

• ajj is the parameters for the j-th camera

• bi is the parameters for the i-th point

• BA attempts to minimize the projection error

• BA attempts to minimize the projection error

Euclidean distance

predicted projection Euclidean distance

3 views and 4 points

(16)

Multiplied by

### Issues in SFM

N li l di i

• Nonlinear lens distortion

• Degeneracy and critical surfaces

• Prior knowledge and scene constraints

• Multiple motions

• Multiple motions

(17)

every 50th frame of a 800-frame sequencey q

lifetime of 3192 tracks from the previous sequencep q

track length histogramg g

(18)

### Nonlinear lens distortion

effect of lens distortion

### Prior knowledge and scene constraints

add a constraint that several lines are parallelp

### Prior knowledge and scene constraints

add a constraint that it is a turntable sequenceq

(19)

### Jurassic park2d3 boujou

Enemy at the Gate, Double Negative

### 2d3 boujou

Enemy at the Gate, Double Negative

(20)

### VideoTrace

http://www.acvt.com.au/research/videotrace/

### Project #3 MatchMove

• It is more about using tools in this project

Y h i h lib i

• You can choose either calibration or structure from motion to achieve the goal

• Calibration

• Voodoo/Icarus

• Examples from previous classes #1 #2

• Examples from previous classes, #1, #2

### References

• Richard Hartley, In Defense of the 8-point Algorithm, ICCV, 1995.

• Carlo Tomasi and Takeo Kanade,,Shape and Motion from Image p g Streams: A Factorization Method, Proceedings of Natl. Acad. Sci., 1993.

Manolis Lourakis and Antonis Argyros The Design and

• Manolis Lourakis and Antonis Argyros, The Design and

Implementation of a Generic Sparse Bundle Adjustment Software Package Based on the Levenberg-Marquardt Algorithm, FORTH- ICS/TR 320 2004

ICS/TR-320 2004.

• N. Snavely, S. Seitz, R. Szeliski,Photo Tourism: Exploring Photo Collections in 3D, SIGGRAPH 2006.,

• A. Hengel et. al.,VideoTrace: Rapid Interactive Scene Modelling from Video, SIGGRAPH 2007.

This paper presents a calibration technique which estimates the mismatch, and adjusts the equivalent gain of the phase noise canceling circuitry to improve the degree of phase

fundamental theorem for line integrals.) A force field that is a gradient field is called a conservative field.. Since the line integral over a closed path is zero, the work done by

integrals given by the Fundamental Theorem, the notation is traditionally used for an antiderivative of f and is called an indefinite integral....

• Goal is to construct a no-arbitrage interest rate tree consistent with the yields and/or yield volatilities of zero-coupon bonds of all maturities.. – This procedure is

• The scene with depth variations and the camera has movement... Planar scene (or a

• The scene with depth variations and the camera has movement... Planar scene (or a

(It is also acceptable to have either just an image region or just a text region.) The layout and ordering of the slides is specified in a language called SMIL.. SMIL is covered in

Because simultaneous localization, mapping and moving object tracking is a more general process based on the integration of SLAM and moving object tracking, it inherits the