====== 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.

===== Epipolar Rectification and Stereo-matching =====

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 [[http://cmp.felk.cvut.cz/cmp/courses/TDV/code|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 )


===== Stereo-matching =====

The implementation of stereo-matching (the ''stereo_gcs'' package) is available in the [[http://cmp.felk.cvut.cz/cmp/courses/TDV/code|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 [[http://cmp.felk.cvut.cz/cmp/courses/TDV/code/run_stereo_for_py.m|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 4x3 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.


=== Task 4-1 ===

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

=== Task 4-2 ===

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