Finding and plotting parabola in matlab

Question

I have been trying to find a parabola in an image. For the starting purposes I took an image with a black parabola on a white background. Then I found the black pixels on the image using the find command by

[yi xi] = find(im<10); % im is the image with black parabola and white background

After that I randomly took 3 points from the collection and solved the equation for the parabola using the symbolic toolbox using

syms x y;
%solve them for the parabola equation
A  = [ x^2 x y 1 ;x0^2 x0 y0 1; x1^2 x1 y1 1; x2^2 x2 y2 1];

where

%(x0,y0) = (104,137)

%(x1,y1) = (244,161)

%(x2,y2) = (300,229)

S = solve(det(A),y);

Then I get the coeffcients a,b,c as

a = 0.0100

b = -1.6800

c = 254.1900

where a, b and c are

a*x^2 + b*x + c = y;

Since now I have got the eqn I plot the parabola by putting the coefficient values and taking

 xx = 1:300;

 yy = a*xx.^2 + b*xx +c ;

then I plot the parabola on the image as

 plot(xx,yy,'-');

For the confirmation that I have taken correct points I also plot the selected points on the image and they lie exactly on the parabola in the image. So that is not the problem.

The problems are :

1. The parabola that I plot (blue) doesn't lie on the parabola of the image(black).

2. When I put the value of x co-ordinates in the above equation. The value of y is not the same as of the y co-ordinate.

   for eg: (104,137)

   0.0100*104*104 -1.68*104 + 254.19 = 108.16 - 174.72 + 254.19 = 187.63 whereas it should be 137

My parabola is wrong. Any help will be appreciated. The images are

Answer 1

I think you must be rounding somewhere in your calculation of a, b, and c.

I had a go with the three points you mentioned using the function fit (with fit type poly2) and my a,b,c were 0.0053, -1.6802 and 254.1895 (or add more decimal places when using format long).

Similarly, when using same the code you give, the output of solve is:

S = (73*x^2)/13720 - (5763*x)/3430 + 87187/343
S2 = double(coeffs(S))

S2 =

  254.1895   -1.6802    0.0053

This gives me the same a, b, and c as with fit/poly2 (just from looking at it, the output of 73/13720 cannot be 0.01). Also, if I plot this curve over the original points using the same code you give, it works fine. So the only remaining source of error that I can see is some sort of rounding in whatever code you use to extract values of a, b, and c from the output of solve.

Answer 2

The following works. The habit of imshow to swap the intuitive meaning of ordinate and abscissa as used in plot can make overlays problematic at times. Your problem might stem from how you define your coordinate system when plotting with different function (imshow or related image display routines versus plot) in particular the position of the origin in the image. Here's my code:

% create image of parabola
k =[-0.3 10 10];
x = [1:0.1:50];
y = round(k(1)*x.^2 + k(2)*x + k(3));
crd = unique(round([x(:) y(:)]),'rows');
Nx = 100;
Ny = 100;
img = ones(Nx,Ny)*255;
iok = find((crd(:,1)>0) & (crd(:,1)<=Nx) & (crd(:,2)>0) & (crd(:,2)<=Ny));
indx = sub2ind([Nx Ny],crd(iok,1),crd(iok,2));
img(indx) = 0;
imshow(img)

Next we fit the parabola

% pick three points:
x0 = crd(1,1);
y0 = crd(1,2);
x1 = crd(80,1);
y1 = crd(80,2);
x2 = crd(160,1);
y2 = crd(160,2);
% plot these on the original image
hold on
plot(y0,x0,y1,x1,y2,x2,'Markersize',20,'Linestyle','none','Marker','x','Color','r')

Solved the equation precisely as you showed, then the result is

S = -1094/3627*x^2+12073/1209*x+37415/3627

Overlay the fitted equation

a= -1094/3627;
b = 12073/1209;
c = 37415/3627;
xx = 1:100;
yy = a*xx.^2 + b*xx +c ;
% then I plot the parabola on the image as
plot(yy,xx,'g--');

A not so pretty plot (green= fit, red crosses=picked points, blue=parabola plotted without swapping order of x and y):