1

閉じた輪郭（ノイズあり）を表すデータがあります：

contour = [(x1, y1), (x2, y2), ...]

輪郭に合わせる簡単な方法はありますか？numpy.polyfitがあります。ただし、x値が繰り返されると失敗し、多項式の適切な次数を決定するためにある程度の努力が必要になります。

4

2 に答える 2

3

フィットしようとしている点から輪郭までの距離は、その点を中心とする極座標の角度の周期関数です。この関数は、フーリエ変換によって正確に計算できる正弦 (または余弦) 関数の組み合わせとして表すことができます。実際には、最初の N 関数に切り捨てられたフーリエ変換によって計算された線形結合は、 Parseval の定理によると、それらの N 関数に最もよく適合します。

これを実際に使用するには、中心点 (おそらく輪郭の重心) を選択し、輪郭を極座標に変換して、中心点からの距離のフーリエ変換を計算します。適合した輪郭は、最初のいくつかのフーリエ係数によって与えられます。

残っている 1 つの問題は、極座標に変換された等高線が等間隔の角度で距離値を持たないことです。これが不規則サンプリング問題です。おそらくかなり高密度のサンプルがあるため、等間隔の角度に最も近い 2 つのポイント間の線形補間を使用するか、(データによっては) 小さなウィンドウで平均化することで、これを非常に簡単に回避できます。不規則なサンプリングに対する他のほとんどのソリューションは、ここでは非常に複雑で不要です。

編集:サンプルコード、作品:

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()

PS 選択した中心点を通る角度に沿ったいくつかの光線が 2 回交差するほど輪郭が非凸 (再入可能) である場合、注意が必要な特殊なケースが 1 つあります。この場合、別の中心点を選択すると役立つ場合があります。極端な場合、このプロパティを持たない中心点がない可能性があります (輪郭がこのように見える場合)。その場合でも、上記の方法を使用して形状に内接または外接することができますが、これはそれ自体をフィッティングするための適切な方法ではありません。このメソッドは、プレッツェルのような「ねじれた」楕円形ではなく、ジャガイモのような「ゴツゴツした」楕円形をフィッティングすることを目的としています:)

于 2012-11-28T23:11:55.647 に答える
1

多項式の次数を修正する場合は、scipy.optimize の leastsq 関数を使用するだけです。

単純な円を生成するとします。それを 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])

返されたタプルの最初の要素のみを取ります。これは、leastsq が必要のない他の多くの情報を返すためです。

一般的な中心を持つ楕円など、より複雑な多項式を当てはめる必要がある場合は、より複雑な関数を使用できます。

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 をご覧ください。

使用しやすいcurve_fit関数ではなく、leastsqを使用した理由は、curve_fitには形式の明示的な関数が必要であり、y = f(x)すべての暗黙的な多項式がその形式に変換できるわけではないためです（または、興味深い暗黙的な多項式はほとんどありません）まったく）

于 2012-11-28T12:41:43.980 に答える