# forum.alglib.net

ALGLIB forum
 It is currently Fri May 24, 2024 4:18 am

 All times are UTC

### Forum rules

1. This forum can be used for discussion of both ALGLIB-related and general numerical analysis questions
2. This forum is English-only - postings in other languages will be removed.

 Page 1 of 1 [ 6 posts ]
 Print view Previous topic | Next topic
Author Message
 Post subject: Fitting circle to points with lsfitnonlinearPosted: Thu May 13, 2010 2:25 pm

Joined: Thu May 13, 2010 2:14 pm
Posts: 2
Hello,

I'm working on a project and i've to find the center of circles.
I use OpenCv to find approximatively the center, then i use lsfitnonlineariteration function on the edges points (by canny), but the program crashes.
Thanks u.

Sources:
Code:
// the result of Hough algorithm on the picture is (161, 121) with a radius of 50 px.
int centre_x = 161, centre_y = 121, r = 50;
Canny(image, edges, 50, 100);
vector<Point> pointsContour; // pointsContour contains the edges points (the superior part of the circle)
for(int i=0; i<centre_y; i++)
{
uchar *ptr = edges.ptr<uchar>(i);
for(int j=0; j<edges.cols; j++)
{
if(ptr[j])
{
pointsContour.push_back(Point(j, i));
}
}
}
int n = pointsContour.size(), m = 1, k = 2;
ap::real_1d_array y, c;
ap::real_2d_array x;
lsfitreport rep;
lsfitstate state;
x.setlength(n, m);
y.setlength(n);
c.setlength(k);
for(int i=0; i<n; i++)
{
x(i, 0) = pointsContour[i].x;
y(i) = pointsContour[i].y;
}
c(0) = centre_x;
c(1) = centre_y;

lsfitnonlinearfgh(x, y, c, n, m, k, state);
lsfitnonlinearsetcond(state, 0.1, 0.1, 0);
while(lsfitnonlineariteration(state))
{
/*
F = sqrt(r? - (x - xo)?) + yo
*/
state.f = sqrt(carre(r) - carre(state.x(0) - state.c(0))) + state.c(1);
if(state.needfg || state.needfgh)
{
/*
DF/Dxo = (-2*xo + x) / (2 * sqrt(r? - (x - xo)?))
DF/Dyo = 1;
*/
state.g(0) = (- 2.0 * state.c(0) + state.x(0)) /
(2.0 * sqrt(carre(r) - carre(state.x(0) - state.c(0))));
state.g(1) = 1.0;
}
if(state.needfgh)
{
/*
D?F/Dxo? = -(-2*xo + x) / (4 * (r? - (x - xo)?)^(3/2))
*/
state.h(0, 0) = - (-2.0 * state.c(0) + state.x(0)) /
(4.0 * (carre(r) - carre(state.x(0) - state.c(0))) * sqrt((carre(r) - carre(state.x(0) - state.c(0)))));
state.h(0, 1) = 0.0;
state.h(1, 0) = 0.0;
state.h(1, 1) = 0.0;
}
}
int info;
lsfitnonlinearresults(state, info, c, rep);
printf("%f\t%f\n", c(0), c(1));

Top

 Post subject: Re: Problem in lsfitnonlineariteration functionPosted: Thu May 13, 2010 3:56 pm

Joined: Fri May 07, 2010 7:06 am
Posts: 915
You didn't provide additional info, but I suspect that your model tries to calculate sqrt() of negative number. That's because you choose incorrect model. You try to minimize vertical distance from point to the circle. But it is much better to minimize difference between (distance(point,circle_center)^2 and circle_radius^2. It is called "total least squares". Actually, it is the only way to solve your problem.

Read http://www.math.sunysb.edu/~scott/Book3 ... ircle.html - it contains some info on your problem. You can use MinLM unit to minimize E(a,b,r) function which were mentioned there.

The only drawback is that you have to write slightly more complex code to calculate F. LSFit unit does this work for you (LSFitNonlinear is just a wrapper around MinLM unit), but it can't be used to solve total least squares problems.

Top

 Post subject: Re: Problem in lsfitnonlineariteration functionPosted: Thu May 13, 2010 6:16 pm

Joined: Thu May 13, 2010 2:14 pm
Posts: 2
Yes, the crash is due to sqrt function. I will use minlmcreatefgh function.
Thanks you very much.

Top

 Post subject: Re: Fitting circle to points with lsfitnonlinearPosted: Sun May 30, 2010 12:49 pm

Joined: Sun May 30, 2010 12:33 pm
Posts: 2
I have a very similar problem. What I have is a list of 2D points to that I want to fit a circle or ellipse with unknown center and radius. The problem is not the mathematical model for, let's say, the circle, but how to apply it in alglib.
As far as I understand lsfitnonlinear is not capabale of solving this minimization problem. So I try to use the minlm unit, but I don't see how to pass all my points to minlmcreatefgh().

Let's say we want to minimize the function E(a,b,r) = SUM[ ((xi-a)^2 + (yi-b)^2 - r^2)^2 ] according to http://www.math.sunysb.edu/~scott/Book3 ... ircle.html. Could you give me a short example on how to apply this with alglib?

Top

 Post subject: Re: Fitting circle to points with lsfitnonlinearPosted: Sun May 30, 2010 7:07 pm

Joined: Fri May 07, 2010 7:06 am
Posts: 915
I don't know what language you use, so I'll post it in pseudocode:

Code:
var
State:  MinLMState;
Rep:    MinLMReport;
S:      array of Real;
X, Y:   Real;
begin
SetLength(S, 3);
S[0]:=0; // A
S[1]:=0; // B
S[2]:=0; // R
MinLMCreateFJ(3, N, S, State);
MinLMSetCond(State, 0.0, 0.0, 0.001, 0);
while MinLMIteration(State) do
begin
if State.NeedF then
State.F:=SUM[ ((x[]-State.X[0])^2 + (y[]-State.X[1])^2 - State.X[2]^2)^2 ];
if State.NeedFiJ then
begin
State.Fi[]:=(x[]-State.X[0])^2 + (y[]-State.X[1])^2 - State.X[2]^2;
State.J[*,0]:=-2*(x[]-State.X[0]);
State.J[*,1]:=-2*(y[]-State.X[1]);
State.J[*,2]:=-2*State.X[2];
end;
end;
MinLMResults(State, S, Rep);

I hope that I typed correct code for calculation of F=SUM(Fi), Fi and J={dFi/dA, dFi/dB, dFi/dR}. Please, recheck it before actually fitting something.

Top

 Post subject: Re: Fitting circle to points with lsfitnonlinearPosted: Sun May 30, 2010 10:09 pm

Joined: Sun May 30, 2010 12:33 pm
Posts: 2
What I found is, that I have to initialize S[3] with a non-zero value, otherwise R will always be zero after the minimization, A is correct and B is wrong. If you keep that in mind and remember that r should be non-negative, everything works fine! Thank you!

Top

 Display posts from previous: All posts1 day7 days2 weeks1 month3 months6 months1 year Sort by AuthorPost timeSubject AscendingDescending
 Page 1 of 1 [ 6 posts ]

 All times are UTC

#### Who is online

Users browsing this forum: No registered users and 3 guests

 You cannot post new topics in this forumYou cannot reply to topics in this forumYou cannot edit your posts in this forumYou cannot delete your posts in this forumYou cannot post attachments in this forum

Search for: