homography_error.h

Go to the documentation of this file.

// Copyright (c) 2011 libmv authors.
//
// Permission is hereby granted, free of charge, to any person obtaining a copy
// of this software and associated documentation files (the "Software"), to
// deal in the Software without restriction, including without limitation the
// rights to use, copy, modify, merge, publish, distribute, sublicense, and/or
// sell copies of the Software, and to permit persons to whom the Software is
// furnished to do so, subject to the following conditions:
//
// The above copyright notice and this permission notice shall be included in
// all copies or substantial portions of the Software.
//
// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
// IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
// AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS
// IN THE SOFTWARE.
 
#ifndef LIBMV_MULTIVIEW_HOMOGRAPHY_ERRORS_H_
#define LIBMV_MULTIVIEW_HOMOGRAPHY_ERRORS_H_
 
#include "libmv/multiview/projection.h"
 
namespace libmv {
namespace homography {
namespace homography2D {
 
struct AsymmetricError {
  static void Residuals(const Mat& H, const Mat& x1, const Mat& x2, Mat2X* dx) {
    dx->resize(2, x1.cols());
    Mat3X x2h_est;
    if (x1.rows() == 2) {
      x2h_est = H * EuclideanToHomogeneous(static_cast<Mat2X>(x1));
    } else {
      x2h_est = H * x1;
    }
    dx->row(0) = x2h_est.row(0).array() / x2h_est.row(2).array();
    dx->row(1) = x2h_est.row(1).array() / x2h_est.row(2).array();
    if (x2.rows() == 2) {
      *dx = x2 - *dx;
    } else {
      *dx = HomogeneousToEuclidean(static_cast<Mat3X>(x2)) - *dx;
    }
  }
  static void Residuals(const Mat& H, const Vec& x1, const Vec& x2, Vec2* dx) {
    Vec3 x2h_est;
    if (x1.rows() == 2) {
      x2h_est = H * EuclideanToHomogeneous(static_cast<Vec2>(x1));
    } else {
      x2h_est = H * x1;
    }
    if (x2.rows() == 2) {
      *dx = x2 - x2h_est.head<2>() / x2h_est[2];
    } else {
      *dx = HomogeneousToEuclidean(static_cast<Vec3>(x2)) -
            x2h_est.head<2>() / x2h_est[2];
    }
  }
  static double Error(const Mat& H, const Mat& x1, const Mat& x2) {
    Mat2X dx;
    Residuals(H, x1, x2, &dx);
    return dx.squaredNorm();
  }
  static double Error(const Mat& H, const Vec& x1, const Vec& x2) {
    Vec2 dx;
    Residuals(H, x1, x2, &dx);
    return dx.squaredNorm();
  }
};
 
struct SymmetricError {
  static double Error(const Mat& H, const Vec& x1, const Vec& x2) {
    // TODO(keir): This is awesomely inefficient because it does a 3x3
    // inversion for each evaluation.
    Mat3 Hinv = H.inverse();
    return AsymmetricError::Error(H, x1, x2) +
           AsymmetricError::Error(Hinv, x2, x1);
  }
  // TODO(julien) Add residuals function \see AsymmetricError
};
struct AlgebraicError {
  // TODO(julien) Make an AlgebraicError2Rows and AlgebraicError3Rows
 
  static void Residuals(const Mat& H, const Mat& x1, const Mat& x2, Mat3X* dx) {
    dx->resize(3, x1.cols());
    Vec3 col;
    for (int i = 0; i < x1.cols(); ++i) {
      Residuals(H, x1.col(i), x2.col(i), &col);
      dx->col(i) = col;
    }
  }
  static void Residuals(const Mat& H, const Vec& x1, const Vec& x2, Vec3* dx) {
    Vec3 x2h_est;
    if (x1.rows() == 2) {
      x2h_est = H * EuclideanToHomogeneous(static_cast<Vec2>(x1));
    } else {
      x2h_est = H * x1;
    }
    if (x2.rows() == 2) {
      *dx = SkewMat(EuclideanToHomogeneous(static_cast<Vec2>(x2))) * x2h_est;
    } else {
      *dx = SkewMat(x2) * x2h_est;
    }
    // TODO(julien) This is inefficient since it creates an
    // identical 3x3 skew matrix for each evaluation.
  }
  static double Error(const Mat& H, const Mat& x1, const Mat& x2) {
    Mat3X dx;
    Residuals(H, x1, x2, &dx);
    return dx.squaredNorm();
  }
  static double Error(const Mat& H, const Vec& x1, const Vec& x2) {
    Vec3 dx;
    Residuals(H, x1, x2, &dx);
    return dx.squaredNorm();
  }
};
// TODO(keir): Add error based on ideal points.
 
}  // namespace homography2D
// TODO(julien) add homography3D errors
}  // namespace homography
}  // namespace libmv
 
#endif  // LIBMV_MULTIVIEW_HOMOGRAPHY_ERRORS_H_