Project 3 : Camera Calibration and Fundamental Matrix Estimation with RANSAC

In this project, camera and scene geometry related theories and techniques are implemented. In the first part, the camera projection matrix which maps 3D world coordinates to image coordinates is estimated. In the second part, the fundamental matrix which relates points in one scene to epipolar lines in another image is also estimated by point correspondence method. In the third part, the estimation of fundamental matrix using RANSAC is implemented.

  1. Camera Projection Matrix
  2. Fundamental Matrix Estimation
  3. Fundamental Matrix with RANSAC

Figure 1 shows the implementation results performing by RANSAC.

Example of implementation of RANSAC.

Code and Algorithms

Part 1: Camera Projection Matrix

In this part, the projection matrix that goes from world 3D coordinates to 2D image coordinates is computed. Linear regression is set up to find the elements of matrix M. The function calculate_projection_matrix.m is shown below.


function M = calculate_projection_matrix( Points_2D, Points_3D )
[row_3d,~]=size(Points_3D);

x=1;
for y=1:2:2*row_3d
    A(y,:)=[Points_3D(x,:),1,0,0,0,0,(-1*Points_2D(x,1))*Points_3D(x,:),-1*Points_2D(x,1)];
    A(y+1,:)=[0,0,0,0,Points_3D(x,:),1,(-1*Points_2D(x,2))*Points_3D(x,:),-1*Points_2D(x,2)];
    x=x+1;
end

[~,~,V]=svd(A);
M=V(:,end);
M=reshape(M,[],3)';

fprintf('Randomly setting matrix entries as a placeholder\n')
end

After the calculation of matrix M, the camera C can be found by equation C=-Q^-1*m4. The code for finding camera C is:


function  Center  = compute_camera_center( M )

Q=M(1:3,1:3);
m4=M(:,4);
Center=-inv(Q)*m4;

end

Part 2: Fundamental Matrix Estimation

This part estimates the mapping of points in one image to lines in another image by means of the fundamental matrix. For this purpose, the decomposition of matrix F using singular value decomposition is used. In addition, a rank 2 matrix is estimated by setting the smallest singular value in S to zero thus generating S2. Finally, the fundamental matrix is calculated as F = U S2 V'. To test the result, the drawing of the epipolar lines on one image which correspond to a point in the other point is presented in later of this report.

The function estimate_fundamental_matrix is shown below:


function  F_matrix  = estimate_fundamental_matrix(Points_a,Points_b)
[row_a,~]=size(Points_a);
[row_b,~]=size(Points_b);

mean_a=mean(Points_a);
mean_b=mean(Points_b);

T_a=[1,0,-1*mean_a(1);% normalize points (extra)
   0,1,-1*mean_a(2);
   0,0,1];
T_b=[1,0,-1*mean_b(1);
   0,1,-1*mean_b(2);
   0,0,1];

a_calc=Points_a;
b_calc=Points_b;

a_calc(:,1)=a_calc(:,1)-mean_a(1);
a_calc(:,2)=a_calc(:,2)-mean_a(2);
s_a=1/std2(a_calc);

b_calc(:,1)=b_calc(:,1)-mean_b(1);
b_calc(:,2)=b_calc(:,2)-mean_b(2);
s_b=1/std2(b_calc);

S_a=[s_a,0,0;
   0,s_a,0;
   0,0,1
    ];
S_b=[s_b,0,0;
   0,s_b,0;
   0,0,1
    ];

a_calc=zeros(row_a,3);
b_calc=zeros(row_b,3);

for i=1:row_a
   a_calc(i,:)=S_a*T_a*[Points_a(i,1);Points_a(i,2);1];
end
Points_a=a_calc(:,1:2);

for i=1:row_b
   b_calc(i,:)=S_b*T_b*[Points_b(i,1);Points_b(i,2);1];
end
Points_b=b_calc(:,1:2);

for i=1:row_a
  A(i,:)=[Points_a(i,1)*Points_b(i,1),Points_b(i,1)*Points_a(i,2),Points_b(i,1),Points_b(i,2)*Points_a(i,1),Points_a(i,2)*Points_b(i,2),Points_b(i,2),Points_a(i,1),Points_a(i,2),1];                                                     
end

[~,~,Q]=svd(A);% perform singular value decomposition
X=reshape(Q(:,end),[3,3])';
[U,S2,V]=svd(X);
[~,index]=min(S2(:));% smallest singular value
[m_r,m_c]=ind2sub(size(S2),index);
S2(m_r,m_c)=0;
F_matrix=U*S2*V';
F_matrix=transpose(S_b*T_b)*F_matrix*(S_a*T_a);%use the scaling matrices to adjust fundamental matrix so that it can operate on the original pixel coordiantes
end

Part 3: Fundamental Matrix with RANSAC

In this part, the fundamental matrix with unreliable point correspondences computed with SIFT (in project 2) is computed. Since last squares regression is not appropriate in this scenario because of multiple outliers, RANSAC is needed in conjunction with fundamental matrix estimation. In order to do this, 9 of point coreespondeces are iteratively chose to solve for the fundamental matrix using the function, and later the number of inliers are counted.


function  [Best_Fmatrix, inliers_a, inliers_b] = ransac_fundamental_matrix(matches_a, matches_b)
[row,~]=size(matches_a);
maximum=0;
threshold=0.002;%threshold between inlier and outlier
inliers_a=[];
inliers_b=[];

for iteration_num=1:700
    a_calc=[];
    b_calc=[];
    rand_pk=randperm(row,9);%implement 9-point correspondence
    
    for i=1:9
        points_a(i,:)=matches_a(rand_pk(i),:);
        points_b(i,:)=matches_b(rand_pk(i),:);
    end
    F=estimate_fundamental_matrix(points_a,points_b);%calculate Fundamental matrix
    
    r=0;
    for i=1:row
        if (abs([matches_b(i,:),1]*F*transpose([matches_a(i,:),1]))<=threshold) %threshold between inlier and outlier
            r=r+1;
            a_calc=vertcat(a_calc,matches_a(i,:));
            b_calc=vertcat(b_calc,matches_b(i,:));
        end
    end
    
    if(r>=maximum)
        maximum=r;
        inliers_a=a_calc;
        inliers_b=b_calc;
        Best_Fmatrix=F;
    end
    
end

if size(inliers_a,1)>100  % number of inliers
    inliers_a=inliers_a(1:100,:);
    inliers_b=inliers_b(1:100,:);
end

end

Some Extra Works

For extra credit, the fundamental matrix is improved by normalizing the coordinates before computing the fundamental matrix. The normalizing part is implemented in the function estimate_fundamental_matrix.m. As shown in estimate_fundamental_matrix.m above, the transformation matrix T is the product of the scale and offset matrix. And finally, the fundamental matrix is adjusted by the using of scaling matrices since the goal is make it operates on the original pixel coordinates. The results in later part of this report will show the improvement.

Results and Discussion

Part 1: Results



Part 1 Results

For the results in Part 1, the Projection Matrix is:

0.4583 -0.2947 -0.0140 0.0040

-0.0509 -0.0546 -0.5411 -0.0524

0.1090 0.1783 -0.0443 0.5968

The total residual is: 0.0445. The estimated location of camera is: <-1.5127, -2.3517, 0.2826>, which is basically the same as in project description.

Part 2: Results



Part 2 Results

For the results in Part 2, the Fundamental Matrix after normalization is:

-0.0000 0.0000 -0.0004

0.0000 -0.0000 0.0032

-0.0000 -0.0044 0.1063

As the images shown above, it is obvious that the epipolar lines on the left image correspond to points in the other image. All of the epipolar lines crossing through the corresponding point in the other image can be observed perfectly.

Part 3: Results

Mount Rushmore Results




Mount Rushmore Results: 700 iteration, normalized

Notre Dame Results




Notre Dame Results: 700 iteration, normalized

Gaudi Results




Gaudi Results: 700 iteration, normalized

Woodruff Dorm Results




Woodruff Dorm Results: 700 iteration, normalized

It can be observed from above image pairs that the fundamental matrix estimation works fine. Most of the spurious feature matches are rejected. Since there were too many inliers during the debug process, the number of inliers shown above is limitted to 100. By this setting, the image pairs can be better observed. In addition, when compared with the image pairs without normalization, it is obvious that the epipolar lines in normalized image are much more consistent and precise. However, the epilolar lines in un-normalized image pairs sometimes make no sense. (For your reference, the un-normalization image pair is also shown below).




Gaudi Results: Un-normalization