What is the correct way to Normalize corresponding points before estimation of Fundamental matrix in OpenCV C++?

Question

I am trying to manually implement a fundamental matrix estimation function for corresponding points (based on similarities between two images). The corresponding points are obtained after performing ORB feature detection, extraction, matching and ratio test.

There is a lot of literature available on good sources about this topic. However none of them appear to give a good pseudo-code for doing the process. I went through various Chapters on Multiple View Geometry book; and also many online sources.

This source appears to give a formula for doing the normalization and I followed the formula mentioned on page 13 of this source.

Based on this formula, I created the following algorithm. I am not sure if I am doing it the right way though !

Normalization.hpp

class Normalization {
    typedef std::vector <cv::Point2f> intercepts;
    typedef std::vector<cv::Mat> matVec;
 public:
    Normalization () {}
    ~Normalization () {}

    void makeAverage(intercepts pointsVec);

    std::tuple <cv::Mat, cv::Mat> normalize(intercepts pointsVec);

    matVec getNormalizedPoints(intercepts pointsVec);

 private:
    double xAvg = 0;
    double yAvg = 0;
    double count = 0;
    matVec normalizedPts;
    double distance = 0;
    matVec matVecData;
    cv::Mat forwardTransform;
    cv::Mat reverseTransform;
};

Normalization.cpp

#include "Normalization.hpp"

typedef std::vector <cv::Point2f> intercepts;
typedef std::vector<cv::Mat> matVec;

/*******
*@brief  : The makeAverage function receives the input 2D coordinates (x, y)
*          and creates the average of x and y
*@params : The input parameter is a set of all matches (x, y pairs) in image A
************/
void Normalization::makeAverage(intercepts pointsVec) {
    count = pointsVec.size();
    for (auto& member : pointsVec) {
        xAvg = xAvg + member.x;
        yAvg = yAvg + member.y;
    }
    xAvg = xAvg / count;
    yAvg = yAvg / count;
}

/*******
*@brief  : The normalize function accesses the average distance calculated
*          in the previous step and calculates the forward and inverse transformation
*          matrices
*@params : The input to this function is a vector of corresponding points in given image
*@return : The returned data is a tuple of forward and inverse transformation matrices
*************/
std::tuple <cv::Mat, cv::Mat> Normalization::normalize(intercepts pointsVec) {
    for (auto& member : pointsVec) {
        //  Accumulate the distance for every point

        distance += ((1 / (count * std::sqrt(2))) *\
                     (std::sqrt(std::pow((member.x - xAvg), 2)\
                                + std::pow((member.y - yAvg), 2))));
    }
    forwardTransform = (cv::Mat_<double>(3, 3) << (1 / distance), \
                        0, -(xAvg / distance), 0, (1 / distance), \
                        -(yAvg / distance), 0, 0, 1);

    reverseTransform = (cv::Mat_<double>(3, 3) << distance, 0, xAvg, \
                        0, distance, yAvg, 0, 0, 1);

    return std::make_tuple(forwardTransform, reverseTransform);
}

/*******
*@brief  : The getNormalizedPoints function trannsforms the raw image coordinates into
*          transformed coordinates using the forwardTransform matrix estimated in previous step
*@params : The input to this function is a vector of corresponding points in given image
*@return : The returned data is vector of transformed coordinates
*************/
matVec Normalization::getNormalizedPoints(intercepts pointsVec) {
    cv::Mat triplet;
    for (auto& member : pointsVec) {
        triplet = (cv::Mat_<double>(3, 1) << member.x, member.y, 1);
        matVecData.emplace_back(forwardTransform * triplet);
    }
    return matVecData;
}

Is this the right way ? Are there other ways of Normalization ?

Answer 1

I think you can do it your way but in "Multiple View Geometry in Computer Vision" Hartley and Zisserman recommend isotropic scaling (p. 107):

Isotropic scaling. As a first step of normalization, the coordinates in each image are translated (by a different translation for each image) so as to bring the centroid of the set of all points to the origin. The coordinates are also scaled so that on the average a point x is of the form x = (x, y,w)T, with each of x, y and w having the same average magnitude. Rather than choose different scale factors for each coordinate direction, an isotropic scaling factor is chosen so that the x and y-coordinates of a point are scaled equally. To this end, we choose to scale the coordinates so that the average distance of a point x from the origin is equal to √ 2. This means that the “average” point is equal to (1, 1, 1)T. In summary the transformation is as follows:
(i) The points are translated so that their centroid is at the origin.
(ii) The points are then scaled so that the average distance from the origin is equal to √2.
(iii) This transformation is applied to each of the two images independently.

They state that it is important for the direct linear transformation (DLT) but even more important for the calculation of a Fundamental Matrix like you want to do. The algorithm you chose, normalized the point coordinates to (1, 1, 1) but did not apply a scaling so that the average distance from the origin is equal to √2.

Here is some code for this type of normalization. The averaging step stayed the same:

std::vector<cv::Mat> normalize(std::vector<cv::Point2d> pointsVec) {
    // Averaging
    double count = (double) pointsVec.size();
    double xAvg = 0;
    double yAvg = 0;
    for (auto& member : pointsVec) {
        xAvg = xAvg + member.x;
        yAvg = yAvg + member.y;
    }
    xAvg = xAvg / count;
    yAvg = yAvg / count;

    // Normalization
    std::vector<cv::Mat> points3d;
    std::vector<double> distances;
    for (auto& member : pointsVec) {

        double distance = (std::sqrt(std::pow((member.x - xAvg), 2) + std::pow((member.y - yAvg), 2)));
        distances.push_back(distance);
    }
    double xy_norm = std::accumulate(distances.begin(), distances.end(), 0.0) / distances.size();

    // Create a matrix transforming the points into having mean (0,0) and mean distance to the center equal to sqrt(2)
    cv::Mat_<double> Normalization_matrix(3, 3); 
    double diagonal_element = sqrt(2) / xy_norm;
    double element_13 = -sqrt(2) * xAvg / xy_norm;
    double element_23 = -sqrt(2)* yAvg/ xy_norm;

    Normalization_matrix << diagonal_element, 0, element_13,
        0, diagonal_element, element_23,
        0, 0, 1;

    // Multiply the original points with the normalization matrix
    for (auto& member : pointsVec) {
        cv::Mat triplet = (cv::Mat_<double>(3, 1) << member.x, member.y, 1);
        points3d.emplace_back(Normalization_matrix * triplet);
    }
    return points3d;
}