I would suggest you look at some of the validation tests in the org.hipparchus.optim.nonlinear.vector.leastsquare
package. The first step would be to look at one of the implementations of test problems (say CircleProblem
for example) and how it is used in LevenbergMarquardtOptimizerTest
to create a LeastSquaresBuilder
.
The two main things to implement are the model function and its Jacobian whith respect to the parameters:
f = \left\{
\begin{gather}
A x_1 + B y_1 + C\\
D x_1 + E y_1 + F\\
\ldots\\
A x_n + B y_n + C\\
D x_n + E y_n + F\\
\end{gather}
\right.
J = \begin{bmatrix}
\frac{\partial x'_1}{\partial A} & \frac{\partial x'_1}{\partial B} & \frac{\partial x'_1}{\partial C} & \frac{\partial x'_1}{\partial D} & \frac{\partial x'_1}{\partial E} & \frac{\partial x'_1}{\partial F}\\
\frac{\partial y'_1}{\partial A} & \frac{\partial y'_1}{\partial B} & \frac{\partial y'_1}{\partial C} & \frac{\partial y'_1}{\partial D} & \frac{\partial y'_1}{\partial E} & \frac{\partial y'_1}{\partial F}\\
\ldots&\ldots&\ldots&\ldots&\ldots&\ldots\\
\frac{\partial x'_n}{\partial A} & \frac{\partial x'_n}{\partial B} & \frac{\partial x'_n}{\partial C} & \frac{\partial x'_n}{\partial D} & \frac{\partial x'_n}{\partial E} & \frac{\partial x'_n}{\partial F}\\
\frac{\partial y'_n}{\partial A} & \frac{\partial y'_n}{\partial B} & \frac{\partial y'_n}{\partial C} & \frac{\partial y'_n}{\partial D} & \frac{\partial y'_n}{\partial E} & \frac{\partial y'_n}{\partial F}\\
\end{bmatrix}
= \begin{bmatrix}
x_1 & y_1 & 1 & 0 & 0 & 0\\
0 & 0 & 0 & x_1 & y_1 & 1\\
\ldots&\ldots&\ldots&\ldots&\ldots&\ldots\\
x_n & y_n & 1 & 0 & 0 & 0\\
0 & 0 & 0 & x_n & y_n & 1
\end{bmatrix}
When you have implemented these two mathematical functions (for example in a class that provides two methods getModelFunction()
and getModelFunctionJacobian()
), you can put them in a LeastSquaresProblem
thanks to a LeastSquaresBuilder
, using something along these lines:
LeastSquaresProblem myProblem =
new LeastSquaresBuilder()
.checkerPair(new SimpleVectorValueChecker(1e-12, 1e-12))
.maxEvaluations(100)
.maxIterations(100)
.model(imageTransform.getModelFunction(),
imageTransform.getModelFunctionJacobian())
.target(targetArrayContainingXprimesYprimes)
.start(startArrayContainingInitialABCDEFWhichMayBeFullOfZeros)
.build();
The last step is to solve this least squares problem. This is done by calling:
Optimum optimum =
new LevenbergMarquardtOptimizer().optimize(myProblem);
The optimum
instance will contain a lot of things related to the optimum solution of your problem. This includes the number of iterations, the value of the cost function, covariances, residuals for all the test points, RMS and of course using getPoint()
the vector containing the parameters (in your case A, B, C, D, E, and F).
There is no need to bother about Householder transformations. Everything is taken care of in the Levenberg-Marquardt algorithm. The Householder transformations are really low level parts and there is a lot of stuff on top of that to build a least squares solver. Levenberg-Marquardt is a good solver. It is a kind of hybrid method between a gradient method and a Gauss-Newton, its most renown property is that it is able to start even far from the solution (which a Gauss-Newton alone can’t do). This is why I suggest to just use an array full of zeros for the start point. If your two images are only slightly different, you could also start with a point where only A and E are set to 1.0, so your initial model corresponds to the identity transform.