#ifndef ICLMEX_H #define ICLMEX_H #include #include /** \mainpage ICLMEx (ICL Model Extraction package) The ICLMex library currently provides functions to extract parameterized models from given sets of x and y coordinates. The algorithms base on the paper Direct Least Square Fitting of Ellipses of Andrew W. Fitzgibbon, Maurizio Milu and Robert B. Fischer. The basic structure is the GeneralModel class. It formalizes all information needed to extract a model from a given set of x and y positions. NOTE: A further todo will be to generalize the algorithms to an abitrary input feature dimension. Yet, no generalization is implemented, so the model extraction functions and data structures are restricted to 2D problems. The following will give a small insight of the used algorithm:

The direct least square fitting algorithm

In the following, models are formalized by \f$ F=AX=0 \f$, where the vector \f$\vec{X}=(f_1(x,y),f_2(x,y),...)^T\f$ contains combined features of the input data \f$x \;and \;y\f$, and the components \f$a_i\f$ of \f$\vec{A}=(a_1,a_2,a_3,...)^T\f$ are coefficient to the features. So models are expressed as the intersection of polynoms over x and y and the \f$F=0\f$ plane. The following examples shall illustrate this:

model examples

In the following five examples are given for the primitives: ellipse, so called restricted ellipse, circle, straight line and straight lines visiting the origin

ellipse models

an elliptic model for example can be parameterized as follows: \f[ ax^2 + bxy + cy^2 + dx + ey + f = 0; \f] this can be rewritten using the above notation: \f[ F(A,X)=0 \; with \; A=(a,b,c,d,e,f)^T \; and \; X=(x^2,xy,y^2,x,y,1)^T \f]

restricted ellipse models

If we leave out the ellipse models mixed term (b), the major axis of the ellipse remain one of the x- and the y-axis. We restate the formula as follows: \f[ F(A,X)=0 \; with \; A=(a,b,c,d,e)^T \; and \; X=(x^2,y^2,x,y,1)^T \f]

circle models

An additional restriction of the ellipse model will forbid the ellipse to have different stretch factors along the two major axis, which complies a circle model: \f[ F(A,X)=0 \; with \; A=(a,b,c,d)^T \; and \; X=(x^2+y^2,x,y,1)^T \f]

straight line models

If we leave out the quadratic Terms of the above model, we get a model for straight lines: \f[ F(A,X)=0 \; with \; A=(a,b,c)^T \; and \; X=(x,y,1)^T \f] in this case the vector \f$(a,b)^T\f$ complies the normal vector of the straight line, and the distance to the origin is defined by c

straight line through the origin models

The above models coefficient c can be left out, to get a model which restricts the lines to visit the origin. \f[ F(A,X)=0 \; with \; A=(a,b)^T \; and \; X=(x,y)^T \f]

Direct least square fitting

The following shall give a coarse insight in the implemented algorithms; further information can be found in the above mentioned paper. We have a set of data points coordinates \f$\vec{x}\f$ and \f$\vec{y}\f$, and want to minimize the function: \f[ \|Da\|^2 \f] where D is the so called design matrix, where row i complies \f$X(x_i,y_i)\f$. So D is in case of searching ellipses: \f[ D=\left(\begin{array}{cccccc} x_1^2 & x_1y_1 & y_1^2 & x_1 & y_1 & 1\\ x_2^2 & x_2y_2 & y_2^2 & x_2 & y_2 & 2\\ x_3^2 & x_3y_3 & y_3^2 & x_3 & y_3 & 3\\ & & ... & & & \\ x_N^2 & x_Ny_N & y_N^2 & x_N & y_N & N\\ \end{array}\right) \f] To forbid the minimization procedure to find the trivial solution \f$A=0\f$ we have to introduce a so called constraint matrix \f$C\f$, which is in the simplest case just the \f$dim(A)\times dim(A)\f$ identity matrix. The minimization can be tackled as a constraint problem introducing Lagrange multipliers: Minimize \f$E=\|Da\|^2\f$ subject to \f$a^TCa=1\f$. In case of \f$C=Id_{dim(A)}\f$ this constraint is equivalent to the constraint \f$\|a\|^2=1\f$. The above mentioned paper discusses also other choices for \f$C\f$. Introducing Lagrange multipliers \f$\lambda\f$ we arrive the System: \f[ \begin{array}{r} 2D^TDa-2\lambda Ca = 0 \\ a^TCa = 1 \end{array} \f] By further introducing the so called scatter matrix \f$S=D^TD\f$ this can be restated as follows: \f[ \begin{array}{r} Sa = \lambda Ca = 0 \\ a^TCa = 1 \end{array} \f] So we have to solve a generalized eigenvector problem. More details should not be discussed here.

Implemented algorithm

The paper gives a short summary of the algorithm to simplify implementation in terms of a 6-line MatLab implementation, which is sketched in the following: -# create the Design matrix \f$D\f$ -# calculate \f$S=D^TD\f$ -# create the constraint matrix \f$C\f$ (here \f$Id\f$) -# find eigenvectors and eigenvalues of \f$S^{-1}C\f$ -# find the largest valid eigenvalue \f$\lambda_{max}\f$ -# optimal parameters are the components of the eigenvector \f$\hat{\nu}_{max}\f$ */ namespace icl{ /// fits a model to given xs and ys coordinate vector /** after the call, the given model reference has the found parameters. This parameters can easily read out with the models "[]" operator. @see GeneralizedModel */ template void fitModel(T *xs, T *ys, unsigned int nPoints, GeneralModel &model); /// this fitModel function implements a metha algorithm to get a better robustness against non gaussian outliers /** if subSetCoutn == -1 it is set to the dim of the model A simple idea to become more robust against non gaussian outliers, this implementation takes subSetCount subsets of the given xs, and ys sets, of size subSetSize, and calculates the seached model params by using the above function for this subset. Single outliers will effect only a small percentage of all of the evaluated subsets. The rule of thumb is to estimate the median of each parameter for the result. (@see Paper: Further Five Point Fit Ellipse Fitting by Paul L. Rosin) */ template void fitModel(T *xs, T *ys,unsigned int nPoints, GeneralModel &model, int subSetSize, int subSetCount=-1); /// draws a model with calculated params /** The model must be fitted with the above function before it can be drawed into an image. */ template void drawModel(GeneralModel &model, Img *image, X *color); /// draws a model with calculated params /** The model must be fitted with the above function before it can be drawed into an image. */ template void drawQuadraticModel(QuadraticModel &model, Img *image, X *color); } #endif