5
$\begingroup$

I'm currently trying to use the Mathematica's FindFit command to fit a curve of some data to the bottom half of a circle. The data is of the form

Data = {{a1, a2}, {b1, b2}, ...}

My code looks like as follows

s = FindFit[Data, y0 - Sqrt[R^2 - (x - x0)^2], {x0, y0, R}, x]

The window is printing the below error,

FindFit::nrlnum: The function value {74724.4 + (-1. + y)^2, 75066.6 + (-1. + y)^2, 75409.6 + (-1. + y)^2, 75753.3 + (-1. + y)^2, 76097.9 + (-1. + y)^2, 76443.2 + (-1. + y)^2, <<40>>, 91274.4 + (-1. + y)^2, 91652.5 + (-1. + y)^2, 92031.4 + (1. + y)^2, 92411.1 + (-1. + y)^2, <<71>>} is not a list of real numbers with dimensions {121} at {x0, y0, R} = {1., 1., 1.}. >>

Any help would be appreciated.

Improve this question
$\endgroup$
6
  • $\begingroup$ Welcome to the site, John! I have formatted your post for you. Your problem is a common syntax error. See mathematica.stackexchange.com/a/18395/8 - underscore is a reserved character. For this reason, I have marked your post as a duplicate of that question. $\endgroup$ Commented Nov 14, 2014 at 6:25
  • $\begingroup$ Well my actual script contains no underscores in the nomenclature. I am new to both mathematica as well as stackexchange so I guess I should have stayed as close to my script as possible. My data name is simple "Data" in my script so there is another underlying error present. I have removed this underscore in my question above to avoid this confusion. $\endgroup$ Commented Nov 14, 2014 at 15:12
  • 1
    $\begingroup$ There is an undefined variable y in either your data or in the actual code that you used, but did not show here. Please inspect your data, and/or quit the kernel and try again. $\endgroup$ Commented Nov 14, 2014 at 15:52
  • 2
    $\begingroup$ Here are some relevant posts from MSE and other forums.1 2 3 4 $\endgroup$ Commented Nov 14, 2014 at 17:24
  • $\begingroup$ There does appear to be an undefined variable within the code. Quitting the kernel appears to have fixed this problem however the window is now printing "FindFit::nrlnum: "The function value {-7.44087-273.373\ I,-7.37627-273.998\ I,-7.29269-274.623\ I,-7.21934-275.248\ I,<<44>>,-2.97508-303.373\ I,-2.93673-303.998\ I,<<71>>} is not a list of real numbers with dimensions {121} at {x0,y0,R} = {1.,1.,1.}."". $\endgroup$ Commented Nov 14, 2014 at 17:26

3 Answers 3

Reset to default
14
$\begingroup$

Here is an approach using linear regression:

circfit[pts_] := Module[{reg, lm, bf, exp, center, rad}, reg = {2 #1, 2 #2, #2^2 + #1^2} & @@@ pts; lm = LinearModelFit[reg, {1, x, y}, {x, y}]; bf = lm["BestFitParameters"]; exp = (x - #2)^2 + (y - #3)^2 - #1 - #2^2 - #3^2 & @@ bf; {center, rad} = {{#2, #3}, Sqrt[#2^2 + #3^2 + #1]} & @@ bf; circlefit[{"expression" -> exp, "center" -> center, "radius" -> rad}]]; circlefit[list_][field_] := field /. list; circlefit[list_]["Properties"] := list /. Rule[field_, _] :> field; circlefit /: ReplaceAll[fields_, circlefit[list_]] := fields /. list; Format[circlefit[list_], StandardForm] := HoldForm[circlefit]["<" <> ToString@Length@list <> ">"]

This assumes that underlying data is a circle and aim is to find center and radius. There are doubtless much better ways.

For fun (and done quickly, so apologies for some inefficiencies and ugliness):

Manipulate[ pt1 = Table[{x, cy + Sqrt[rad^2 - (x - cx)^2]}, {x, -4, 4, 0.1}]; pt2 = Table[{x, cy - Sqrt[rad^2 - (x - cx)^2]}, {x, -4, 4, 0.1}]; pt = Cases[Join[pt1, pt2], {_Real, _Real}]; pts = pt + RandomReal[{-p, p}, {Length[pt], 2}]; With[{fit = circfit[pts]}, Column[{Show[ ContourPlot[Evaluate@fit["expression"], {x, -1, 8}, {y, -1, 8}, Contours -> {0}, ContourShading -> None, ContourStyle -> Directive[Red, Thick]], ListPlot[pts], AspectRatio -> Automatic, PlotLabel -> fit["expression"] == 0, ImageSize -> 400, PlotRange -> {{-1, 8}, {-1, 8}}], Grid[{{"variable", "value", "model"}, {"center", {cx, cy}, fit["center"]}, {"radius", rad, fit["radius"]}}] }]], {p, 0.05, 1}, {cx, 0, 2}, {cy, 0, 2}, {rad, 1, 5} ]

Animated Gif showing the result of the manipulate

Improve this answer
$\endgroup$
5
  • $\begingroup$ I have to wonder what those much better ways might be? I'd be surprised if there is something better than a least squares solution, which is what you have provided. $\endgroup$ Commented Nov 16, 2014 at 1:24
  • $\begingroup$ @DanielLichtblau...I guess I mean better ways of implementing the idea. I'm still not sure what the underlying problem is,e.g. whether estimation of known circle or circular approximation of closed figure...where perhaps interpolation better for latter... $\endgroup$ Commented Nov 16, 2014 at 3:02
  • $\begingroup$ Thank you. I was able to get it working by simply Taylor expanding the lower half of the circle and fitting the curve to the produced polynomial. Your method appears to provide more correct values however. $\endgroup$ Commented Nov 17, 2014 at 17:37
  • 1
    $\begingroup$ Knowing the underlying problem as ubpdqn points out is critical. If all one needs is a circle that looks good visually, then the above method works fine (and that's not to take away from the above method). But if statistical inferences need to be made about the parameters, then a model that mimics the data generation process is needed. Then maximum likelihood or a Bayesian approach could be used to estimate the model parameters and associated measures of precision. $\endgroup$ Commented Jun 12, 2015 at 16:43
  • 1
    $\begingroup$ @DanielLichtblau you could use MLESAC or RANSAC if there is extreme noise present. That's still least squares/MLE but more robust to outliers, also once appeared in this image processing answer $\endgroup$ Commented Jun 12, 2015 at 18:10
6
$\begingroup$

Not exactly the problem, but I have this code languishing on my hard drive to compute a circle center and radius given a bunch of pts. From there it is straight-forward to get the polynomial.

First make some data...

pts = Table[RandomReal[{.9, 1.1}] {Sin[#], Cos[#]} & @ RandomReal[{0, 2 π}], {10}] (* {{0.347046, -0.930661}, {0.218395, -1.07027}, {0.890144, -0.179376}, \ {-0.828852, -0.367061}, {0.770728, 0.607237}, {-0.118834, 0.934153}, {0.21855, -0.949232}, {0.136237, 1.02154}, {-0.881469, -0.285087}, {-0.878687, -0.377516}}*)

Define a penalty function

penalty[c_, radius_, pts_] := (Map[radius - Norm[# - c] &, pts])^2 // Total

Use Minimize[ ] to find the center of the circle and radius

res = Minimize[penalty[{x, y}, rr, pts], {x, y, rr}] (* {0.0268756, {x -> 0.0247297, y -> -0.0119147, rr -> 0.969958}} *) circ = Circle[{x, y}, rr] /. (res // Last) Show[circ // Graphics, ListPlot[pts]]

enter image description here

Improve this answer
$\endgroup$
2
  • $\begingroup$ I made a small adjustment to your code to eliminate the non-localized Symbol theta and replace it with an anonymous function. I hope you do not mind. $\endgroup$ Commented Jun 6, 2017 at 0:38
  • $\begingroup$ Nope, thanks for cleaning it up. $\endgroup$ Commented Jun 6, 2017 at 14:28
1
$\begingroup$

I added a third element to the data list to give a "z" value that is a radius of 1.

Data = {{a1, a2, 1}, {b1, b2, 1}, ...}

Then, the FindFit function works well.

FindFit[data, ((x - a)^2 + (y - b)^2)/r^2, {a, b, r}, {x, y}]
Improve this answer
$\endgroup$
1
  • $\begingroup$ for this to work reliably you need to supply a good initial guess for a,b (eg use the mean value of the data ) $\endgroup$ Commented Jun 5, 2017 at 16:54

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.