如何拟合闭合轮廓？
我有一个表示闭合轮廓(带有噪声)的数据：
contour = [(x1, y1), (x2, y2), ...]有没有什么简单的方法来拟合轮廓？这就是numpy.polyfit。但是，如果x值重复，并且需要一些努力来确定多项式的适当次数，那么它就失败了。
回答 2
Stack Overflow用户
发布于 2012-11-29 07:11:56
从一个点到要拟合的等高线的距离是以该点为中心的极坐标中角度的周期函数。这个函数可以表示为正弦(或余弦)函数的组合，可以通过傅立叶变换精确计算。实际上，根据Parseval's theorem的说法，通过截断到前N个函数的傅里叶变换计算的线性组合最适合这N个函数。
要在实践中使用这一点，请拾取一个中心点(可能是轮廓的重心)，将轮廓转换为极坐标，并计算距中心点的距离的傅立叶变换。拟合的轮廓由最初的几个傅立叶系数给出。
剩下的一个问题是，转换为极坐标的等高线在等间距角度上没有距离值。这就是Irregular Sampling问题。由于你大概有一个相当高的样本密度，你可以很简单地通过在两个最近的点之间使用线性插值到均匀间隔的角度，或者(取决于你的数据)用一个小窗口平均来绕过这个问题。在这里，不规则采样的大多数其他解决方案都要复杂得多，而且没有必要。
编辑：示例代码，工作：
import numpy, scipy, scipy.ndimage, scipy.interpolate, numpy.fft, math
# create simple square
img = numpy.zeros( (10, 10) )
img[1:9, 1:9] = 1
img[2:8, 2:8] = 0
# find contour
x, y = numpy.nonzero(img)
# find center point and conver to polar coords
x0, y0 = numpy.mean(x), numpy.mean(y)
C = (x - x0) + 1j * (y - y0)
angles = numpy.angle(C)
distances = numpy.absolute(C)
sortidx = numpy.argsort( angles )
angles = angles[ sortidx ]
distances = distances[ sortidx ]
# copy first and last elements with angles wrapped around
# this is needed so can interpolate over full range -pi to pi
angles = numpy.hstack(([ angles[-1] - 2*math.pi ], angles, [ angles[0] + 2*math.pi ]))
distances = numpy.hstack(([distances[-1]], distances, [distances[0]]))
# interpolate to evenly spaced angles
f = scipy.interpolate.interp1d(angles, distances)
angles_uniform = scipy.linspace(-math.pi, math.pi, num=100, endpoint=False)
distances_uniform = f(angles_uniform)
# fft and inverse fft
fft_coeffs = numpy.fft.rfft(distances_uniform)
# zero out all but lowest 10 coefficients
fft_coeffs[11:] = 0
distances_fit = numpy.fft.irfft(fft_coeffs)
# plot results
import matplotlib.pyplot as plt
plt.polar(angles, distances)
plt.polar(angles_uniform, distances_uniform)
plt.polar(angles_uniform, distances_fit)
plt.show()在这种情况下，选择不同的中心点可能会有所帮助。在极端情况下，可能没有不具有此属性的中心点(如果您的轮廓看起来像like this)。在这种情况下，您仍然可以使用上面的方法来内切或限制您拥有的形状，但这本身并不是一个合适的方法来拟合它。此方法用于拟合像土豆一样的“块状”椭圆形，而不是像椒盐卷饼那样的“扭曲”椭圆形:)
Stack Overflow用户
发布于 2012-11-28 20:41:44
如果你确定了你的多项式次数，你可以简单地使用scipy.optimize中的至少函数
假设你生成了一个简单的圆。我会把它分成x和y两个部分。
data = [ [cos(t)+0.1*randn(),sin(t)+0.1*randn()] for t in rand(100)*2*np.pi ]
contour = array(data)
x,y = contour.T在给定多项式系数的情况下，编写一个简单的函数来计算每个点与0的差值。我们将曲线拟合为以原点为中心的圆。
def f(coef):
a = coef
return a*x**2+a*y**2-1我们可以简单地使用leastsq函数来找到最佳系数。
from scipy.optimize import leastsq
initial_guess = [0.1,0.1]
coef = leastsq(f,initial_guess)[0]
# coef = array([ 0.92811554])我只接受返回元组的第一个元素，因为至少返回了许多我们不需要的其他信息。
def f(coef):
a,b,cx,cy = coef
return a*(x-cx)**2+b*(y-cy)**2-1
initial_guess = [0.1,0.1,0.0,0.0]
coef = leastsq(f,initial_guess)[0]
# coef = array([ 0.92624664, 0.93672577, 0.00531 , 0.01269507])编辑：
如果出于某种原因，您需要估计拟合参数的不确定性，您可以从结果的协方差矩阵中获得此信息：
res = leastsq(f,initial_guess,full_output=True)
coef = res[0]
cov = res[1]
#cov = array([[ 0.02537329, -0.00970796, -0.00065069, 0.00045027],
# [-0.00970796, 0.03157025, 0.0006394 , 0.00207787],
# [-0.00065069, 0.0006394 , 0.00535228, -0.00053483],
# [ 0.00045027, 0.00207787, -0.00053483, 0.00618327]])
uncert = sqrt(diag(cov))
# uncert = array([ 0.15928997, 0.17768018, 0.07315927, 0.07863377])协方差矩阵的对角线是每个参数的方差，所以不确定性是它的平方根
有关拟合过程的更多信息，请访问http://www.scipy.org/Cookbook/FittingData。
我之所以使用leastsq而不是curve_fit函数，是因为curve_fit需要一个y = f(x)形式的显式函数，而不是每个隐式多项式都可以转换成这种形式(或者更好的是，几乎没有任何有趣的隐式多项式)。
https://stackoverflow.com/questions/13604611
复制相似问题
腾讯云开发者
