如何找到椭圆的方程

Question

我希望能够使用一般圆锥曲线的方程，在给定五个或六个点的情况下找到椭圆的方程：

A x2 + B xy + C y2 + D x + E y + F = 0。

起初我尝试使用六个点。这是我的Python代码：

    import numpy as np
    def conic_section(p1, p2, p3, p4, p5, p6):
        def row(point):
            return [point[0]*point[0], point[0]*point[1], point[1]*point[1],                         
            point[0], point[1], 1]
        matrix=np.matrix([row(p1),row(p2),row(p3),row(p4),row(p5), row(p6)])
        b=[0,0,0,0,0,0]
        return np.linalg.solve(matrix,b)

    print conic_section(np.array([6,5]), np.array([2,9]), np.array([0,0]),         
    np.array([11, 5.5]), np.array([6, 7]), np.array([-1,-1]))

问题在于这会返回解[0,0,0,0,0,0]，因为等式右边是零向量。

然后我试图通过减去F并将其除以来更改圆锥曲线：

A x2 + B xy + C y2 + D x + E y + F = 0

A x2 + B xy + C y2 + D x + E y = -F

A' x2 + B xy + C' y2 + D' x + E' y = -1.

但是这样做不起作用的原因是，如果我的一个点是（0,0），那么我最终会得到一个有一行零的矩阵，但是等式的右边在向量中的条目中具有-1。换句话说，如果我的一个点是（0,0）-那么“F”应该是0，所以我不能将它除掉。

任何帮助都将不胜感激。 谢谢。

- Jbag1212

---

你说你不能除掉它是什么意思？线性代数完全可以处理这个问题。假设这六个点是独立的，并且在椭圆上，那么系统可以解出方程。 - Willem Van Onsem

请提供六个可能会导致问题的点。 - Willem Van Onsem

如果一个椭圆有任何一个点与原点相接触，那么它方程中的F值等于0。:) 在这种情况下，如果你有一个没有常数的椭圆，我认为你不能使用这种方法。如果我错了，请有人纠正我。 - cs95

如果F=0，那么您不能将方程式除以F，因为你无法除以零。其中五个不起作用的点是：np.array([6,5]), np.array([2,9]), np.array([0,0]), np.array([11, 5.5]), np.array([6, 7])。 - Jbag1212

如果你知道你有一个椭圆（而不是更一般的圆锥曲线），A 必须是非零的。由于你可以通过任意因子缩放 A 到 F，你可以在你的线性方程组中添加一个额外的约束条件 A = 1（如果你有六个点，则删除其中一个；五个点足以确定圆锥曲线）。 - Mark Dickinson

Answer 1

看起来你有椭圆上的确切点，不需要近似，并通过某种奇怪的方式使用 Braikenridge-Maclaurin 构造来得到圆锥曲线。

五个点(x[i],y[i])确定了这个显式方程的椭圆(mathworld page, eq. 8)

因此，要找到椭圆方程，您可以通过对第一行的子式进行余子式展开（代数余子式）来构建行列式。例如，系数A是从x1y1到右下角的子矩阵的行列式值，系数B是没有xiyi列的子矩阵的行列式的否定值，依此类推。

- MBo

Answer 2

椭圆方程（不考虑平移和旋转）：

目标是解决变量A到F的线性方程:

使用：

from math import sin, cos, pi, sqrt

import matplotlib.pyplot as plt
import numpy as np
from numpy.linalg import eig, inv

# basis parameters of the ellipse
a = 7
b = 4

def ellipse(t, a, b):
    return a*cos(t), b*sin(t)

points = [ellipse(t, a, b) for t in np.linspace(0, 2*pi, 100)]
x, y = [np.array(v) for v in list(zip(*points))]

fig = plt.figure()
plt.scatter(x, y)
plt.show()

def fit_ellipse(x, y):
    x = x[:, np.newaxis]
    y = y[:, np.newaxis]
    D =  np.column_stack((x**2, x*y, y**2, 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))
    return V[:, n]

A, B, C, D, E, F = fit_ellipse(x, y)
K = D**2/(4*A) + E**2/(4*C) - F

# a, b
print('a:', sqrt(K/A), 'b:', sqrt(K/C))

输出：

a: 6.999999999999998 b: 4.0

参见：

http://mathworld.wolfram.com/ConicSection.html https://fr.wikipedia.org/wiki/Ellipse_(math%C3%A9matiques)#Forme_matricielle http://nicky.vanforeest.com/misc/fitEllipse/fitEllipse.html

- glegoux

---

这些方程式中的未知数是 A 到 F，且它们是线性的。问题在于这些方程式是齐次线性方程式。 - Mark Dickinson