Fit data into a circle

Hi, I’m new to Julia. I’ve been looking for a way in Julia to fit my data into a circle.
My data is a set of (x, y) coordinates, and I want to fit them into a circle to get the centre and radius etc. I know some libraries that do this in Matlab/python but haven’t found similar things in Julia. Is there anything you would recommend? or any reading material that I might have missed?
thank you all very much!

The https://github.com/JuliaNLSolvers/LsqFit.jl package will do what you want.

1 Like

Minimize \Sigma_i ((x_i - c_x)^2 + (y_i - c_y)^2 - r^2)^2 by choosing c_x, c_y and r. Use Optim.jl.

4 Likes

Couldn’t you just compute the mean to find the center and the rms distance from the mean to find the radius? i.e.

using Statistics
center = (mean(x), mean(y))
radius = hypot(stdm(x, center[1]), stdm(y, center[2]))
6 Likes

This requires the points to be evenly distributed around the circle.

Not really. If you analytically derive the optimal c_x, c_y and r from the objective in my comment above, you will arrive at the neat recipe suggested by Steven.

Edit: this is wrong. See below.

2 Likes

I’m sorry, I thought @stevengj was suggesting the mean, stdm approach as an alternative to running an iterative solver or optimizer.

This is not the cleanest code in the world (some allocations that could be trimmed at least), and seems a bit fragile with respect to initial conditions, but here’s a working example using LsqFit.

using LsqFit
	
data = [(-2.0,1.0), (0.0, 3.0), (2.0, 1.0)]
	
ptrs = 1:length(data)
	
y = zeros(length(data))
	
function circle_model(t,p)
	out = []
	x0,y0,r = p
	for ptr ∈ ptrs
		x,y = data[ptr]
		push!(out,(x-x0)^2+(y-y0)^2-r^2)
	end
	out
end
	
p0 = [0.0, 0.0, 1.0] # using 0.0 as initial radius fails
	
fit = curve_fit(circle_model, ptrs, y, p0)
	
fit.param
1 Like

There is no need for an iterative solver if there is an analytical solution :slight_smile:

Thank you all very much!
I tried both approaches and the results are slightly different (see pic).
I was wondering if they are essentially the same then why still slightly different?
(I used a particularly bad/uneven dataset to test)
left: mean stdm
right: LsqFit

compare

If the data is uneven, you likely want a fitting approach, not mean/stdm. (The two are not equivalent.)

1 Like

You can’t just average all the x’s and all the y’s to get the center. Because you have more points in the lower right, the center point is biased towards that direction. Looking at this visually, I think the least squares fit got the answer right.

[EDIT: an interesting test would be to filter out that clump of points with x<30400 and y>17000 and repeat the two techniques. The mean will find the center of the remaining arc, clearly not the center. The LsqFit should still find the center although it may get lost more easily.]

I am curious if some other solver would be more efficient but I don’t know if speed matters to you.

Yes agreed (see pic)compare
For this specific purpose, the speed might not be a big issue.
next step I’ll try to find if there’s a way to fit data into a 10by10-circle array. I know the diameter and distance between centres in the ideal situation, but really data are noisy and some of them may not be a circle but a dot or an arc instead.
do you think I can use LsqFit to achieve this?
Thanks!

Thanks for running my suggested test case, It looks like the LsqFit is hanging in there pretty good.

I don’t understand what you mean by “10by10-circle array”. Do you have one set of data with 100 of these patterns that you’ll need to separate somehow?

https://ieeexplore.ieee.org/document/1246564
Here is an overview of a few methods, if you plan to write a circle fit package :wink:

A. Full Least-Squares (FLS) Method
B. Average of Intersections Method
C. Reduced Least-Squares (RLS) Method
D. Modified Least-Squares (MLS) Methods
E. Ka˙sa’s Method

1 Like

If I have time I may look into point E. I think this was also extended for hyperspheres.

I wrote a small package using the loss function of Kasa’s method in combination with LsqFit.lmfit.

This is just a very simple implementation and I did not check if the standard error is correct at all.
PRs welcome :slight_smile:

3 Likes

Do you have one set of data with 100 of these patterns that you’ll need to separate somehow?

Yes that’s the case. see the illustration.


there are 100 (or more) possible positions but not all of them are filled. I need to get each of them separately and run the data analysis pipeline for each of them. I know the distance between the centres and the radius in theory.
previously I was using Matlab’s computer vision toolbox to cluster them based on distance, without using the prior knowledge of the “theoretical distribution”. But I kind of feel maybe there is a way I can fit the entire dataset into a model and get these 100+ circles in one go.
This is already beyond what I was originally asking, and I’m still learning and exploring different approaches. it would be great if you have any idea or suggestion about this.
Thanks a lot for all the help!

1 Like

Do you know the positions of the circles? Otherwise you could try a circle Hough transform to get an initial guess.
I think Images.jl has an implementation.

1 Like

I know the relative positions, like the pitch, but not the absolute position in each data set. (i.e. I don’t know the precise coordinates of the any of these circles, and the array may be tilted or have some other slight deformation

Otherwise you could try a circle Hough transform to get an initial guess.
I think Images.jl has an implementation.

Thanks a lot! I’ll check it out.