Estimation of Transformations Shao-Yi Chien Department of - - PowerPoint PPT Presentation

Estimation of Transformations

簡韶逸 Shao-Yi Chien Department of Electrical Engineering National Taiwan University Fall 2018

1

Outline

• Estimation – 2D Projective Transformation

2

[Slides credit: Marc Pollefeys]

Parameter Estimation

• 2D homography

Given a set of (xi,xi’), compute H (xi’=Hxi)

• 3D to 2D camera projection

Given a set of (Xi,xi), compute P (xi=PXi)

• Fundamental matrix

Given a set of (xi,xi’), compute F (xi’TFxi=0)

• Trifocal tensor

Given a set of (xi,xi’,xi”), compute T

Number of Measurements Required

• At least as many independent equations as degrees
• f freedom required
• Example:

Hx x'

                                 1 1 λ

33 32 31 23 22 21 13 12 11

y x h h h h h h h h h y x

2 independent equations / point 8 degrees of freedom 4x2≥8

Approximate Solutions

• Minimal solution
• 4 points yield an exact solution for H
• More points
• Robust estimation algorithms, such as RANSAC
• No exact solution, because measurements are inexact

(“noise”)

• Search for “best” according to some cost function
• Algebraic or geometric/statistical cost

Gold Standard Algorithm

• Cost function that is optimal for some assumptions
• Computational algorithm that minimizes it is called

“Gold Standard” algorithm

• Other algorithms can then be compared to it

Direct Linear Transformation (DLT)

i i

Hx x   Hx x   

i i

            

i i i i

x h x h x h Hx

3 2 1 T T T

                       

i i i i i i i i i i i i i i

y x x w w y x h x h x h x h x h x h Hx x

1 2 3 1 2 3 T T T T T T

h h h x x x x x x

3 2 1

                             

T T T T T T T T T i i i i i i i i i i i i

x y x w y w

 

T i i i i

w y x      , , x

h A 

i



Direct Linear Transformation (DLT)

• Equations are linear in h

h h h x x x x x x

3 2 1

                             

T T T T T T T T T i i i i i i i i i i i i

x y x w y w h A 

i

• Only 2 out of 3 are linearly independent

(indeed, 2 eq/pt)

Direct Linear Transformation (DLT)

• Equations are linear in h

h A 

i

• Only 2 out of 3 are linearly independent

(indeed, 2 eq/pt)

h h h x x x x

3 2 1

                      

T T T T T T i i i i i i i i

x w y w

• Holds for any homogeneous

representation, e.g. (xi’,yi’,1)

Direct Linear Transformation (DLT)

• Solving for H

Ah 

h A A A A

4 3 2 1

            

size A is 8x9 or 12x9, but rank 8 Trivial solution is h=09

T is not interesting

1-D null-space yields solution of interest pick for example the one with

1 h 

Direct Linear Transformation (DLT)

• Over-determined solution

No exact solution because of inexact measurement i.e. “noise”

Ah 

h A A A

n 2 1

             

Find approximate solution

• Additional constraint needed to avoid 0, e.g.
• not possible, so minimize

1 h  Ah Ah 

DLT Algorithm

Objective Given n≥4 2D to 2D point correspondences {xi↔xi’}, determine the 2D homography matrix H such that xi’=Hxi Algorithm (i) For each correspondence xi ↔xi’ compute Ai. Usually

• nly two first rows needed.

(ii) Assemble n 2x9 matrices Ai into a single 2nx9 matrix A (iii) Obtain SVD of A. Solution for h is last column of V (iv) Determine H from h

Inhomogeneous Solution

                  ' ' h ~ ' ' ' ' ' ' ' ' ' '

i i i i i i i i i i i i i i i i i i i i i i i i

x w y w x y x x w w w y w x y y y x w w w y w x

Since h can only be computed up to scale, pick hj=1, e.g. h9=1, and solve for 8-vector

h ~

Solve using Gaussian elimination (4 points) or using linear least-squares (more than 4 points) However, if h9=0 this approach fails also poor results if h9 close to zero Therefore, not recommended Note h9=H33=0 if origin is mapped to infinity

 

1 H 1 Hx l            

 T

Degenerate Configurations

x4 x1 x3 x2 x4 x1 x3 x2 H? H’? x1 x3 x2 x4

Hx x   

i i

Constraints: i=1,2,3,4

T

l x H

4 *

 

Define:

 

4 4 4 4 *

x x l x x H     k

T

 

3 , 2 , 1 , x l x x H

4 *

    i

i i T

Then, H* is rank-1 matrix and thus not a homography (case A) (case B) If H* is unique solution, then no homography mapping xi→xi’(case B) If further solution H exist, then also αH*+βH (case A) (2-D null-space in stead of 1-D null-space)

Solutions from Lines

i i

l H l

T

  Ah 

2D homographies from 2D lines Minimum of 4 lines Minimum of 5 points or 5 planes 3D Homographies (15 dof) 2D affinities (6 dof) Minimum of 3 points or lines Conic provides 5 constraints

Solutions from Mixed Type

• 2D homography
• cannot be determined uniquely from the

correspondence of 2 points and 2 line

• can from 3 points and 1 line or 1 point and 3 lines

16

Cost Functions

• Algebraic distance
• Geometric distance
• Reprojection error
• Comparison
• Geometric interpretation
• Sampson error

Algebraic Distance

Ah

DLT minimizes

Ah  e

residual vector

i

e

partial vector for each (xi↔xi’) algebraic error vector

 

2 2 2 alg

h x x x x Hx , x                 

T T T T T T i i i i i i i i i i i

x w y w e d algebraic distance

 

2 2 2 1 2 2 1 alg

x , x a a d   where

 

2 1 T 3 2 1

x x , , a    a a a

 

2 2 2 2 alg

Ah Hx , x e e d

i i i i i

   

 

Not geometrically/statistically meaningfull, but given good normalization it works fine and is very fast (use for initialization for non-linear minimization)

Geometric Distance

measured coordinates estimated coordinates true coordinates

x

x ˆ

x

 

2 H

x H , x argmin H ˆ

i i i

d  



Error in one image

e.g. calibration pattern

 

 

2 2 1

• H

Hx , x x H , x argmin H ˆ

i i i i i

d d    



Symmetric transfer error d(.,.) Euclidean distance (in image)

 

   

i i i i i i i i i

d d

i i

x ˆ H ˆ x ˆ subject to x ˆ , x x ˆ , x argmin x ˆ , x ˆ , H ˆ

2 2 x ˆ , x ˆ H,

      





Reprojection error

Symmetric Transfer Error v.s. Reprojection Error

 

 

2 2 1

• Hx

, x x H x,    d d

   

2 2

x ˆ , x x ˆ x,    d d

Symmetric Transfer Error Reprojection Error

Comparison of Geometric and Algebraic Distances

                   

i i i i i i i i i i

w x x w y w w y e ˆ ˆ ˆ ˆ h A Error in one image

 

T i i i i

w y x      , , x

 

x H ˆ , ˆ , ˆ x ˆ      

T i i i i

w y x

                     

3 2 1

h h h x x x x

T T T T T T i i i i i i i i

x w y w

     

2 2 2 i i alg

ˆ ˆ ˆ ˆ x ˆ , x

i i i i i i i i

w x x w y w w y d              

     

 

 

i i i i i i i i i i

w w d w x w x w y w y d                    ˆ / x ˆ , x / ˆ / ˆ ˆ / ˆ / x ˆ , x

i i alg 2 / 1 2 2 2 i i

typical, but not , except for affinities

1  

i

w

i 3x

h ˆ  

i

w

 For affinities, DLT can minimize geometric distance

these two distance metrics are related, but not identical

represents 2 quadrics in 4 (quadratic in X)

Geometric Interpretation of Reprojection Error

Estimating homography ~ fit surface to points X=(x,y,x’,y’)T in 4

Hx x   

i i

H

ν            

2 2 2 2 2 2 2 i i

x ˆ , x x ˆ , x ˆ ˆ ˆ ˆ X ˆ X

i i i i i i i i i i i i

d d y y x x y y x x                 

     

2 H 2 2

, X x ˆ , x x ˆ , x



i i i i i

d d d



   

Analog to conic fitting

 

Cx x C , x

T 2 alg

 d

 

2

C , x



d

Sampson Error

2

X ˆ X 

H

ν

Vector that minimizes the geometric error is the closest point on the variety to the measurement

X ˆ

X

between algebraic and geometric error

X ˆ

Sampson error: 1st order approximation of

 

X Ah

H

  C

   

X H H X H

δ X X δ X      C C C X X ˆ δX  

 

X ˆ

H

 C

 

δ X X

X H H

    C C e  

X

Jδ

Find the vector that minimizes subject to

X

δ e  

X

Jδ

X

δ X J

H

with

   C

Find the vector that minimizes subject to

X

δ e  

X

Jδ

X

δ

Use Lagrange multipliers: minimize derivatives

 

Jδ 2λ

• δ

δ

X X X

 e

T T T

J 2λ

• δ

2

X



 

Jδ 2

X

  e

λ J δ X

T

  λ JJ

T

   e

  e

1 T

JJ λ



  

  e

1 T T X

JJ J δ



  

X

δ X X ˆ  

  e

e

1 T T X T X 2 X

JJ δ δ δ



 

Sampson Error

Sampson Error

2

X ˆ X 

H

ν

Vector that minimizes the geometric error is the closest point on the variety to the measurement

X ˆ

X

between algebraic and geometric error

X ˆ

Sampson error: 1st order approximation of

 

X Ah

H

  C

   

X H H X H

δ X X δ X      C C C X X ˆ δX  

 

X ˆ

H

 C

 

δ X X

X H H

    C C e  

X

Jδ

Find the vector that minimizes subject to

X

δ e  

X

Jδ

X

δ X J

H

with

   C

  e

e

1 T T X T X 2 X

JJ δ δ δ



 

(Sampson error)

Sampson Approximation

A few points (i) For a 2D homography X=(x,y,x’,y’) (ii) is the algebraic error vector (iii) is a 2x4 matrix, e.g. (iv) Similar to algebraic error in fact, same as Mahalanobis distance (v) Sampson error independent of linear reparametrization (cancels out in between e and J) (vi) Must be summed for all points (vii) Close to geometric error, but much fewer unknowns

  e

e

1 T T 2 X

JJ δ





e e e

T 2



2 JJT

e

 

31 21 3 T 2 T 11

/ h x h x h y h w x y w J

i i i i i i

           

X J

H

   C

 

X

H

C e 

 



 e

e

1 T T JJ

Statistical Cost Function and Maximum Likelihood Estimation

• Optimal cost function related to noise model of measurement
• Assume zero-mean isotropic Gaussian noise (assume outliers

removed)

 

  



2 2

2 / x x, 2

πσ 2 1

x Pr

 d

e 

   

   

2 2 i

2 x H , x

/ 2 i

πσ 2 1

H | x Pr



i

d

e

i







 

     



     constant x H , x H | x Pr log

2 i 2 i

σ 2 1

i

d

Error in one image Maximum Likelihood Estimate

 





2 i

x H , x

i

d

Equivalent to minimizing the geometric error function

Statistical Cost Function and Maximum Likelihood Estimation

• Optimal cost function related to noise model of measurement
• Assume zero-mean isotropic Gaussian noise (assume outliers

removed)

 

  



2 2

2 / x x, 2

πσ 2 1

x Pr

 d

e 

   

 

   

2 2 i 2 i

2 x H , x x , x

/ 2 i

πσ 2 1

H | x Pr



       





 

i i

d d

e

i

Error in both images Maximum Likelihood Estimate

   

2 i 2 i

x ˆ , x x ˆ , x

i i

d d   



Equivalent to minimizing the reprojection error function

Mahalanobis Distance

• General Gaussian case

Measurement X with covariance matrix Σ

   

X X X X X X

1 T 2

    

  2 2

X X X X

 

    

Error in two images (independent)

2 2

X X X X

i i

i i i i i  

    



Varying covariances

Invariance to Transforms ?

x T x ~     Tx x ~  Hx x  

x ~ H ~ x ~   T H ~ T H

1

• ?

  Tx H ~ x T    Tx H ~ T x

• 1

  

will result change? for which algorithms? for which transformations?

Non-invariance of DLT

Given and H computed by DLT, and Does the DLT algorithm applied to yield ?

i i i i

x T x ~ , Tx x ~     

i i

x x  

i i

x ~ x ~  

• 1

HT T H ~  

 

i i i i i

e Tx HT T x T x ~ H ~ x ~ ~

• 1

       

 

i i i

e

* *

T Hx x T      

   

h A , R ~ , ~ h ~ A ~

i 2 1 2 1 i

s e e s e e   

T T

Effect of change of coordinates on algebraic error         1 t sR T

        s s R t

• R

T

T *

for similarities so

 

 

i i i i

sd d x ~ H ~ , x ~ Hx , x

alg alg

  

Non-invariance of DLT

(T*: cofactor matrix)

Non-invariance of DLT

Given and H computed by DLT, and Does the DLT algorithm applied to yield ?

i i i i

x T x ~ , Tx x ~     

i i

x x  

i i

x ~ x ~  

• 1

HT T H ~  

 

   

1 H ~ subject to x ~ H ~ , x ~ minimize 1 H subject to x ~ H ~ , x ~ minimize 1 H subject to Hx , x minimize

2 alg 2 alg 2 alg

       

  

i i i i i i i i i

d d d

Invariance of Geometric Error

 

 

   

i i i i i i i i

sd d d d Hx , x Hx T , x T Tx HT T , x T x ~ H ~ , x ~

• 1

          

Given and H, and Assume T’ is a similarity transformations

, x T x ~ , Tx x ~

i i i i

    

i i

x x   , x ~ x ~

i i

 

• 1

HT T H ~  

Normalizing Transformations

• Since DLT is not invariant,

what is a good choice of coordinates?

e.g. Isotropic scaling

• Translate centroid to origin
• Scale to a average distance to the origin
• Independently on both images

2

1 norm

1 2 / 2 / T



             h h w w h w

Or

Importance of Normalization

h h h 1 1

3 2 1

                              

i i i i i i i i i i i i i i

x y x x x y x y y y x y y x

~102 ~102 ~102 ~102 ~104 ~104 ~102 1 1

• rders of magnitude difference!

Without normalization With normalization

Normalized DLT Algorithm

Objective Given n≥4 2D to 2D point correspondences {xi↔xi’}, determine the 2D homography matrix H such that xi’=Hxi Algorithm (i) Normalize points (ii) Apply DLT algorithm to (iii) Denormalize solution

, x ~ x ~

i i

 

i norm i i norm i

x T x ~ , x T x ~     

norm

• 1

norm

T H ~ T H  

Employ this algorithm instead of the original DLT algorithm!

• More accurate
• Invariant to arbitrary choices of the scale and coordinate origin

Normalization is also called pre-conditioning

Iterative Minimization Methods

Required to minimize geometric error (i) Often slower than DLT (ii) Require initialization (iii) No guaranteed convergence, local minima (iv) Stopping criterion required Therefore, careful implementation required: (i) Cost function (ii) Parameterization (minimal or not) (iii) Cost function ( parameters ) (iv) Initialization (v) Iterations

Parameterization

Parameters should cover complete space and allow efficient estimation of cost

• Minimal or over-parameterized? e.g. 8 or 9

(minimal often more complex, also cost surface) (good algorithms can deal with over-parameterization) (sometimes also local parameterization)

• Parametrization can also be used to restrict

transformation to particular class, e.g. affine

Function Specifications

(i) Measurement vector XN with covariance Σ (ii) Set of parameters represented by vector P N (iii) Mapping f : M →N. Range of mapping is surface S representing allowable measurements (iv) Cost function: squared Mahalanobis distance Goal is to achieve , or get as close as possible in terms of Mahalanobis distance

         

P X P X P X

1 T 2

f f f     

 

 

X P  f

Error in one image

 





2 i

x H , x

i

d

 

n

f Hx ,..., Hx , Hx h :

2 1



 

h X f 

 

 

2 2 1

• Hx

, x x H , x

i i i i i

d d   



Symmetric transfer error

 

n n

f Hx ,..., Hx , Hx , x H ,..., x H , x H h :

2 1

• 1

2

• 1

1

• 1

   

 

h X f 

Reprojection error

   

2 i 2 i

x ˆ , x x ˆ , x

i i

d d   



 

h X f 

X composed of 2n inhomogeneous coordinates of the points 𝑦𝑗

′

X composed of 4n-vector inhomogeneous coordinates of the points𝑦𝑗 and 𝑦𝑗

′

X composed of 4n-vector

Initialization

• Typically, use linear solution
• If outliers, use robust algorithm
• Alternative, sample parameter space

Iteration Methods

Many algorithms exist

• Newton’s method
• Levenberg-Marquardt
• Powell’s method
• Simplex method

Levenberg-Marquardt Algorithm

For a mapping function f with parameter vector To an estimated measurement vector We want to find p that can minimize , where f(p) can be approximated as with small and J= Find to minimize

44

Levenberg-Marquardt Algorithm

Find to minimize The least-square solution: Augmented normal equation (with damping term ):

45

Hessian

Gold Standard Algorithm

Objective Given n≥4 2D to 2D point correspondences {xi↔xi’}, determine the Maximum Likelyhood Estimation of H (this also implies computing optimal xi’=Hxi) Algorithm (i) Initialization: compute an initial estimate using normalized DLT

• r RANSAC

(ii) Geometric minimization of -Either Sampson error:

• Minimize the Sampson error
• Minimize using Levenberg-Marquardt over 9 entries of h
• r Gold Standard error:
• compute initial estimate for optimal {xi}
• minimize cost
• ver {H,x1,x2,…,xn}
• if many points, use sparse method

   

2 i 2 i

x ˆ , x x ˆ , x

i i

d d   



Robust Estimation

• What if set of matches contains gross outliers?

RANSAC: RANdom SAmple Consensus

Objective Robust fit of model to data set S which contains outliers Algorithm (i) Randomly select a sample of s data points from S and instantiate the model from this subset. (ii) Determine the set of data points Si which are within a distance threshold t of the model. The set Si is the consensus set of samples and defines the inliers of S. (iii) If the subset of Si is greater than some threshold T, re- estimate the model using all the points in Si and terminate (iv) If the size of Si is less than T, select a new subset and repeat the above. (v) After N trials the largest consensus set Si is selected, and the model is re-estimated using all the points in the subset Si

Distance Threshold

Choose t so probability for inlier is α (e.g. 0.95)

• Often empirically
• Zero-mean Gaussian noise σ then follows

distribution with m=codimension of model

2 

d

2 m



(dimension+codimension=dimension space)

Codimension Model t 2 1 Line (l), Fundamental matrix (F) 3.84σ2 2 Homography (H), Camera Matrix (P) 5.99σ2 3 Trifocal tensor (T) 7.81σ2

How Many Samples?

Choose N so that, with probability p, at least one random sample is free from outliers. e.g. p=0.99

   

 

s

e p N     1 1 log / 1 log

 

 

p e

N s

    1 1 1

proportion of outliers e

s 5% 10% 20% 25% 30% 40% 50% 2 2 3 5 6 7 11 17 3 3 4 7 9 11 19 35 4 3 5 9 13 17 34 72 5 4 6 12 17 26 57 146 6 4 7 16 24 37 97 293 7 4 8 20 33 54 163 588 8 5 9 26 44 78 272 1177

Acceptable Consensus Set?

• Typically, terminate when inlier ratio reaches expected ratio
• f inliers

 n

e T   1

Adaptively Determining the Number

• f Samples

e is often unknown a priori, so pick worst case, e.g. 50%,

and adapt if more inliers are found, e.g. 80% would yield

e=0.2

• N=∞, sample_count =0
• While N >sample_count repeat
• Choose a sample and count the number of inliers
• Set e=1-(number of inliers)/(total number of points)
• Recompute N from e
• Increment the sample_count by 1
• Terminate

   

   

s

e p N     1 1 log / 1 log

Robust Maximum Likelyhood Estimation

Previous MLE algorithm considers fixed set of inliers Better, robust cost function (reclassifies)

   

      





• utlier

inlier ρ with ρ

2 2 2 2 2 2 i

t e t t e e e d i R

Other Robust Algorithms

• RANSAC maximizes number of inliers
• LMedS minimizes median error

Automatic Computation of H

Objective Compute homography between two images Algorithm (i) Interest points: Compute interest points in each image (ii) Putative correspondences: Compute a set of interest point matches based on some similarity measure (iii) RANSAC robust estimation: Repeat for N samples (a) Select 4 correspondences and compute H (b) Calculate the distance d for each putative match (c) Compute the number of inliers consistent with H (d<t) Choose H with most inliers (iv) Optimal estimation: re-estimate H from all inliers by minimizing ML cost function with Levenberg-Marquardt (v) Guided matching: Determine more matches using prediction by computed H Optionally iterate last two steps until convergence

Determine Putative Correspondences

• Compare interest points

Similarity measure:

• SAD, SSD, ZNCC on small neighborhood
• If motion is limited, only consider interest points with

similar coordinates

• More advanced approaches exist, based on invariance…
• Such as SIFT

Example: robust computation

Interest points (500/image) Putative correspondences (268) Outliers (117) Inliers (151) Final inliers (262)