c++opencvlevenberg-marquardt

OpenCV::LMSolver getting a simple example to run

PROBLEM: The documentation for cv::LMSolver at opencv.org is very thin, to say the least. Finding some useful examples on the internet, also, was not possible.

APPROACH: So, I did write some simple code:

#include <opencv2/calib3d.hpp>
#include <iostream>

using namespace cv;
using namespace std;

struct Easy : public LMSolver::Callback {
   Easy() = default;

   virtual bool compute(InputArray f_param, OutputArray f_error, OutputArray f_jacobian) const override
   { 
        Mat  param = f_param.getMat();

        if( f_error.empty() ) f_error.create(1, 1, CV_64F);             // dim(error) = 1
        Mat  error = f_error.getMat();                          

        vector<double> x{param.at<double>(0,0), param.at<double>(1,0)}; // dim(param) = 2
        double         error0 = calc(x);
        error.at<double>(0,0) = error0;

        if( ! f_jacobian.needed() ) return true; 
        else if( f_jacobian.empty() ) f_jacobian.create(1, 2, CV_64F); 
        Mat  jacobian = f_jacobian.getMat();

        double e = 1e-10; // estimate derivatives in epsilon environment       
        jacobian.at<double>(0, 0) = (calc({x[0] + e, x[1]    }) - error0) / e; // d/dx0 (error)
        jacobian.at<double>(0, 1) = (calc({x[0],     x[1] + e}) - error0) / e; // d/dx1 (error)
        return true;
   }

   double calc(const vector<double> x) const { return x[0]*x[0] + x[1]*x[1]; }
};

int main(int argc, char** argv) 
{
    Ptr<Easy>      callback   = makePtr<Easy>();
    Ptr<LMSolver>  solver     = LMSolver::create(callback, 100000, 1e-37);
    Mat            parameters = (Mat_<double>(2,1) << 5, 100);
    solver->run(parameters);
    cout << parameters << endl;
}

QUESTIONS:

Solution

  • What does the return value of LMSolver::Callback::compute() report to the caller?

    Thankfully, opencv is opensource, so we might be able to figure this out simply by checking out the code.

    Looking at the source code on Github, I found that all of the calls to compute() look like:

    if( !cb->compute(x, r, J) )
    return -1;

    Returning false simply causes the solver to bail out. So it seems that the return value of the callback's compute() is simply whether the generation of the jacobian was successful or not.

    Currently, it finds the minimum at (-9e-07,4e-5). How can the precision be improved?

    If anything, you should at least compare the return value of run() against your maximum iteration count to make sure that it did, in fact, converge as much as it could.