title: fit a ellipse in Python given a set of points xi=(xi,yi)
body: <p>I am computing a series of index from a 2D points (x,y). One index is the ratio between minor and major axis. To fit the ellipse i am using the following <a href="http://nicky.vanforeest.com/misc/fitEllipse/fitEllipse.html#the-approach" rel="nofollow noreferrer">post</a></p>
  
  <p>when i run these function the final results looks strange because the center and the axis length are not in scale with the 2D points</p>
  
  <pre><code>center =  [  560415.53298363+0.j  6368878.84576771+0.j]
  angle of rotation =  (-0.0528033467597-5.55111512313e-17j)
  axes =  [0.00000000-557.21553487j  6817.76933256  +0.j]
  </code></pre>
  
  <p>thanks in advance for help</p>
  
  <pre><code>import numpy as np
  from numpy.linalg import eig, inv
  
  def fitEllipse(x,y):
      x = x[:,np.newaxis]
      y = y[:,np.newaxis]
      D =  np.hstack((x*x, x*y, y*y, x, y, np.ones_like(x)))
      S = np.dot(D.T,D)
      C = np.zeros([6,6])
      C[0,2] = C[2,0] = 2; C[1,1] = -1
      E, V =  eig(np.dot(inv(S), C))
      n = np.argmax(np.abs(E))
      a = V[:,n]
      return a
  
  def ellipse_center(a):
      b,c,d,f,g,a = a[1]/2, a[2], a[3]/2, a[4]/2, a[5], a[0]
      num = b*b-a*c
      x0=(c*d-b*f)/num
      y0=(a*f-b*d)/num
      return np.array([x0,y0])
  
  def ellipse_angle_of_rotation( a ):
      b,c,d,f,g,a = a[1]/2, a[2], a[3]/2, a[4]/2, a[5], a[0]
      return 0.5*np.arctan(2*b/(a-c))
  
  def ellipse_axis_length( a ):
      b,c,d,f,g,a = a[1]/2, a[2], a[3]/2, a[4]/2, a[5], a[0]
      up = 2*(a*f*f+c*d*d+g*b*b-2*b*d*f-a*c*g)
      down1=(b*b-a*c)*( (c-a)*np.sqrt(1+4*b*b/((a-c)*(a-c)))-(c+a))
      down2=(b*b-a*c)*( (a-c)*np.sqrt(1+4*b*b/((a-c)*(a-c)))-(c+a))
      res1=np.sqrt(up/down1)
      res2=np.sqrt(up/down2)
      return np.array([res1, res2])
  
  if __name__ == '__main__':
  
      points = [(560036.4495758876, 6362071.890493258),
       (560036.4495758876, 6362070.890493258),
       (560036.9495758876, 6362070.890493258),
       (560036.9495758876, 6362070.390493258),
       (560037.4495758876, 6362070.390493258),
       (560037.4495758876, 6362064.890493258),
       (560036.4495758876, 6362064.890493258),
       (560036.4495758876, 6362063.390493258),
       (560035.4495758876, 6362063.390493258),
       (560035.4495758876, 6362062.390493258),
       (560034.9495758876, 6362062.390493258),
       (560034.9495758876, 6362061.390493258),
       (560032.9495758876, 6362061.390493258),
       (560032.9495758876, 6362061.890493258),
       (560030.4495758876, 6362061.890493258),
       (560030.4495758876, 6362061.390493258),
       (560029.9495758876, 6362061.390493258),
       (560029.9495758876, 6362060.390493258),
       (560029.4495758876, 6362060.390493258),
       (560029.4495758876, 6362059.890493258),
       (560028.9495758876, 6362059.890493258),
       (560028.9495758876, 6362059.390493258),
       (560028.4495758876, 6362059.390493258),
       (560028.4495758876, 6362058.890493258),
       (560027.4495758876, 6362058.890493258),
       (560027.4495758876, 6362058.390493258),
       (560026.9495758876, 6362058.390493258),
       (560026.9495758876, 6362057.890493258),
       (560025.4495758876, 6362057.890493258),
       (560025.4495758876, 6362057.390493258),
       (560023.4495758876, 6362057.390493258),
       (560023.4495758876, 6362060.390493258),
       (560023.9495758876, 6362060.390493258),
       (560023.9495758876, 6362061.890493258),
       (560024.4495758876, 6362061.890493258),
       (560024.4495758876, 6362063.390493258),
       (560024.9495758876, 6362063.390493258),
       (560024.9495758876, 6362064.390493258),
       (560025.4495758876, 6362064.390493258),
       (560025.4495758876, 6362065.390493258),
       (560025.9495758876, 6362065.390493258),
       (560025.9495758876, 6362065.890493258),
       (560026.4495758876, 6362065.890493258),
       (560026.4495758876, 6362066.890493258),
       (560026.9495758876, 6362066.890493258),
       (560026.9495758876, 6362068.390493258),
       (560027.4495758876, 6362068.390493258),
       (560027.4495758876, 6362068.890493258),
       (560027.9495758876, 6362068.890493258),
       (560027.9495758876, 6362069.390493258),
       (560028.4495758876, 6362069.390493258),
       (560028.4495758876, 6362069.890493258),
       (560033.4495758876, 6362069.890493258),
       (560033.4495758876, 6362070.390493258),
       (560033.9495758876, 6362070.390493258),
       (560033.9495758876, 6362070.890493258),
       (560034.4495758876, 6362070.890493258),
       (560034.4495758876, 6362071.390493258),
       (560034.9495758876, 6362071.390493258),
       (560034.9495758876, 6362071.890493258),
       (560036.4495758876, 6362071.890493258)]
  
  
      a_points = np.array(points)
      x = a_points[:, 0]
      y = a_points[:, 1]
      from pylab import *
      plot(x,y)
      show()
      a = fitEllipse(x,y)
      center = ellipse_center(a)
      phi = ellipse_angle_of_rotation(a)
      axes = ellipse_axis_length(a)
  
      print "center = ",  center
      print "angle of rotation = ",  phi
      print "axes = ", axes
  
      from pylab import *
      plot(x,y)
      plot(center[0:1],center[1:], color = 'red')
      show()
  </code></pre>
  
  <p>each vertex is a xi,y,i point
  <img src="https://i.stack.imgur.com/GxBZO.png" alt="enter image description here"></p>
  
  <p>plot of 2D point and center of fit ellipse
  <img src="https://i.stack.imgur.com/NiRu9.png" alt="enter image description here"></p>
  
  <p>using OpenCV i have the following result:</p>
  
  <pre><code>import cv
  
  PointArray2D32f = cv.CreateMat(1, len(points), cv.CV_32FC2)
  for (i, (x, y)) in enumerate(points):
      PointArray2D32f[0, i] = (x, y)
      # Fits ellipse to current contour.
     (center, size, angle) = cv.FitEllipse2(PointArray2D32f)
  
  (center, size, angle) 
  ((560030.625, 6362066.5),(10.480490684509277, 17.20206642150879),144.34889221191406)
  </code></pre>