5
$\begingroup$

I have data points $(x,y)$ for a plane curve, and I would like to find its curvature. While I was googling to check how could I start, I found this matlab code: 

  mx = mean(x); my = mean(y)
  X = x - mx; Y = y - my; % Get differences from means
  dx2 = mean(X.^2); dy2 = mean(Y.^2); % Get variances
  t = [X,Y]\(X.^2-dx2+Y.^2-dy2)/2; % Solve least mean squares problem
  a0 = t(1); b0 = t(2); % t is the 2 x 1 solution array [a0;b0]
  r = sqrt(dx2+dy2+a0^2+b0^2); % Calculate the radius
  a = a0 + mx; b = b0 + my; % Locate the circle's center
  curv = 1/r; % Get the curvature

No specific explanation was added to this code except:

"The circle defined by center $(a, b)$ and radius $r$ will yield the least mean square value for the expression $(x-a)^2 + (y-b)^2 - r^2$ among all possible parameters, $a$, $b$, and $r$, over the points in vectors $x$ and $y$. The circle's curvature will be $1/r$."

The code worked very well, I checked this by taking as an example a circle.

I also tried to follow what they did step by step to understand the concept used behind this code; however I was stuck on the use of the mean. I know least square is used to solve systems, but why use mean least square? Why is the mean and variance included?

Would someone help me please. Any help will be highly appreciated.

$\endgroup$
4
  • $\begingroup$ Curvature is not the same as variance. Also it is likely not very constant for the whole curve. There is one very easy way to check it anyway. Construct short pieces of circles and run the code on them. $\endgroup$
    – mathreadler
    Oct 16, 2015 at 14:23
  • $\begingroup$ Yes I do that @mathreadler, and I check the matlab code, It works nicely. However I dont understand the idea used in this code to compute the curvature. $\endgroup$
    – Nizar
    Oct 16, 2015 at 14:26
  • $\begingroup$ No one has an idea about this matlab code ?? $\endgroup$
    – Nizar
    Oct 18, 2015 at 15:52
  • $\begingroup$ Can you do least squares for a straight line fit? What difference do you find between a circle and a straight line modeling? $\endgroup$
    – Narasimham
    Oct 21, 2015 at 14:16

1 Answer 1

Reset to default
6
+25
$\begingroup$

The Matlab code may be explained by the fact that the curvature is one divided by the radius of an Osculating circle . Since all what you have is a bunch of data points, it has to be decided which data points may be relevant for the curvature at some place at the curve. At least three of them are needed to determine a circle. If more data points are taken, then a least squares procedure will be needed in order to determine a Best Fit Circle to them. For a better understanding of the latter, one might read the following document:

There is no unambiguous recipe for how many data points to take for an osculating circle at some place at the curve. It all depends upon your data and all you can have is a decent approximation. Now I don't know how Matlab works, but it seems to me that these data points form the arrays $x$ and $y$ in the code. However, speaking for myself, I'd rather ignore this particular implementation and start from scratch, with my own understanding of the matter.

Update. To demonstrate that Robert Israel's method (as explained in the reference) does indeed work:
enter image description here
The code itself is to large to fit into the margin. Serious: it really has no sense to publish all those (Delphi Pascal) programming details: I'm pretty sure that wouldn't add much to understanding.
But if someone asks me for the source code, then he or she can have it for free.

$\endgroup$
4
  • $\begingroup$ Thank you for your answer, What i did really is that I wrote another matlab code that find the center and radius of the oscillating circle ( I take the number of given points as variable greater than three), I got good results that exactly mathch the results of the above code> But really I am very curios t know how this code works?? why using the mean and the variance ??!! $\endgroup$
    – Nizar
    Oct 19, 2015 at 20:12
  • $\begingroup$ @Nizar: Without an accompanying documentation, not being an expert in Matlab, it's virtually impossible for me to do any reverse engineering about what theory is behind that code. However, I you read the PDF document mentioned in my answer, then you will see that mean $(\mu_x,\mu_y)$ and variance $(\sigma_{xx},\sigma_{xy},\sigma_{yy})$ come with the RI method in this paper. I have no idea if that may be comparable with the theory as implemented by your Matlab program. $\endgroup$
    – Han de Bruijn
    Oct 20, 2015 at 13:07
  • $\begingroup$ Thak you again for your help, but sorry I have a small question in the pdf you suggested for me : What is the role of $w_i$ in the least square formula? where it come from? what did it represent?? $\endgroup$
    – Nizar
    Oct 22, 2015 at 6:39
  • 2
    $\begingroup$ @Nizar: The $w_i$ represent weighting coefficients; they can be conveniently taken $=1$ here. (Sorry I didn't explain them in the paper; that happens if there are no referees to check content) $\endgroup$
    – Han de Bruijn
    Oct 22, 2015 at 11:20

You must log in to answer this question.

Not the answer you're looking for? Browse other questions tagged .