Surface Reconstruction

Now we assume that all cameras are calibrated with respect to the common (global) coordinate system. The last phase of 3D scene reconstruction is (dense) reconstruction of surface.

The task is to reconstruct dense point cloud. We will use stereo-matching, which provides dense correspondences in a pair of images. Prior to this step, the images must be epipolarly rectified. The implementation of rectification (the rectify package) is available in the code repository.

The images of a single pair are modified by homographies Ha and Hb during rectification. Using the original camera matrices of the pair, we need to compute the camera matrices belonging to these new images, thus belonging to the set of established correspondences.

Let Pa and Pb are the two cameras (full, including K) of the image pair imga, imgb before rectification, and Fab is the fundamental matrix between the images. Then the rectified image pair and the rectifying homographies can be obtained as

import rectify

[Ha, Hb, imga_r, imgb_r] = rectify.rectify( F12, imga, imgb )

The implementation of stereo-matching (the stereo_gcs package) is available in the code repository for matlab only. When working in python, it is necessary to save the rectified data to a matlab file, run the stereo in Matlab environment, and load the results back to python.

The stereo matching works best when inicialized with correspondence seeds from sparse matching. So the inlier correspondences from a particular pair should be taken, modified by rectifying homographies to be valid for the rectified image pair, and then used for the matching. In python, the code should look as follows:

import scipy.io

task = [] # task variable for all data

# For every selected image pair with indices i_a and i_b

#   - load the image im_a, im_b, compute fundamental matrix F_ab
#   - load corresponing points u_a, u_b
#   - keep only inliers w.r.t. F_ab

    [H_a, H_b, im_a_r, im_b_r] = rectify.rectify( F_ab, im_a, im_b )

#   - modify corresponding points by H_a, H_b

#   seeds are rows of coordinates [ x_a, x_b, y ]
#      (corresponding point have the same y coordinate when rectified)
#   assuming u_a_r and u_b_r euclidean correspondences after rectification:
    seeds = np.vstack( ( u_a_r[0,:], u_b_r[0], ( u_a_r[1] + u_b_r[1] ) / 2 ) ).T

    task_i = np.array( [ im_a_r, im_b_r, seeds ], dtype=object )
    task += [ task_i ]

# now all stereo tasks are prepared, save to a matlab file

task = np.vstack( task )
scipy.io.savemat( 'stereo_in.mat', { 'task': task } )

# here run the gcs stereo in matlab, and then load the results

d = scipy.io.loadmat( 'stereo_out.mat' )
D = d['D']

# a disparity map for i-th pair is in D[i,0]
i = 0  # 1, 2, ...
Di = D[i,0]

plt.imshow( Di )

The code that runs the gcs stereo for matlab is here.

Now Di contains disparities in the first rectified image of the i-th pair (im_a_r). This means that pixel (x,y) in the first rectified image corresponds to the pixel (x+Di[y,x],) in the second rectified image (im_b_r). This correspondence should be then transformed to the original images by (inverse of) rectifying homographies, and triangulated using the cameras P_a, P_b, and gives a single point of the dense point cloud. There are pixels where disparity is not computed; the contain nan value. All other pixels give a point cloud corresponding to a part of scene surface.

This dense reconstruction should be done for different image pairs (img_i,img_j). Note that for each pair there is different fundamental matrix Fij and different homographies Ha, Hb. However, the triangulated points are all in the common coordinate system, so they can be simply merged together.

The choice of pairs to be used in stereo is not unequivocal. We recommend to (manually) select horizontal and vertical pairs constructed from neighbouring cameras. In the case of 4×3 image capturing scheme this leads to 9 horizontal and 8 vertical pairs. However, it may be needed to cover parts of surface by additional pairs.

1. Compute disparity maps between chosen image pairs. Show the results.

1. Construct the dense point cloud.
2. Export the point cloud and the cameras into ply file.