2014-10-26 14 views
5

Bir dairenin merkezini bulmak için bazı veri noktalarını sığdırmaya çalışıyorum. Aşağıdaki noktalarının tamamı dairenin çevresi gürültülü veri noktalar şunlardır: Ben scipy http://wiki.scipy.org/Cookbook/Least_Squares_Circle gibi bazı kütüphane kullanmaya çalışıyordu ama sorun mevcut işlevleri kullanarak yaşıyorumPython'da en küçük kare kullanılarak dairenin merkezi nasıl bulunur?

data = [(2.2176383052987667, 4.218574252410221), 
(3.3041214516913033, 5.223500807396272), 
(4.280815855023374, 6.461487709813785), 
(4.946375258539319, 7.606952538212697), 
(5.382428804463699, 9.045717060494576), 
(5.752578028217334, 10.613667377465823), 
(5.547729017414035, 11.92662513852466), 
(5.260208374620305, 13.57722448066025), 
(4.642126672822957, 14.88238955729078), 
(3.820310290976751, 16.10605425390148), 
(2.8099420132544024, 17.22588), 
(1.5731539516426183, 18.17052077121059), 
(0.31752822350872545, 18.75261434891438), 
(-1.2408437559671106, 19.119355580780265), 
(-2.680901948575409, 19.15018791257732), 
(-4.190406775175328, 19.001321726517297), 
(-5.533990404926917, 18.64857428377178), 
(-6.903383826792998, 17.730112542165955), 
(-8.082883753215347, 16.928080323602334), 
(-9.138397388219254, 15.84088004983959), 
(-9.92610373064812, 14.380575762984085), 
(-10.358670204629814, 13.018017342781242), 
(-10.600053524240247, 11.387283417089911), 
(-10.463673966507077, 10.107554951600699), 
(-10.179820255235496, 8.429558128401448), 
(-9.572153386953028, 7.1976672709797676), 
(-8.641475289758178, 5.8312286526738175), 
(-7.665976739804268, 4.782663065707469), 
(-6.493033077746997, 3.8549965442534684), 
(-5.092340806635571, 3.384419909199452), 
(-3.6530364510489073, 2.992272643733981), 
(-2.1522365767310796, 3.020780664301393), 
(-0.6855406924835704, 3.0767643753777447), 
(0.7848958776292426, 3.6196842530995332), 
(2.0614188482646947, 4.32795711960546), 
(3.2705467984691508, 5.295836809444288), 
(4.359297538484424, 6.378324784240816), 
(4.981264502955681, 7.823851404553242)] 

.

Orada örneğin geçerli:

# == METHOD 2 == 
from scipy  import optimize 

method_2 = "leastsq" 

def calc_R(xc, yc): 
    """ calculate the distance of each 2D points from the center (xc, yc) """ 
    return sqrt((x-xc)**2 + (y-yc)**2) 

def f_2(c): 
    """ calculate the algebraic distance between the data points and the mean circle centered at c=(xc, yc) """ 
    Ri = calc_R(*c) 
    return Ri - Ri.mean() 

center_estimate = x_m, y_m 
center_2, ier = optimize.leastsq(f_2, center_estimate) 

xc_2, yc_2 = center_2 
Ri_2  = calc_R(*center_2) 
R_2  = Ri_2.mean() 
residu_2 = sum((Ri_2 - R_2)**2) 

Ama bu tek xy kullandığı görünen? Bu işlevi veri örneğime nasıl bağlayacağınıza dair herhangi bir fikir var mı?

cevap

2

Veri noktalarınız oldukça temiz görünüyor ve aykırı değerler görmüyorum, pek çok çember uydurma algoritması çalışacaktır. doğrusal en küçük kareler çözüldü, sonra

AX+BY+C=X²+Y²

(X-Xc)^2+(Y-Yc)^2=R²

2XcX+2YcY+R²-Xc²-Yc²=X²+Y² olarak yeniden yazılır:

Sana sihirli sorunu doğrusallaştırılmasıyla çalışır Coope yöntemiyle, başlangıç ​​için tavsiye .

2

Çevreye uyum sağlama konusunda tecrübem yok, ancak genel olarak elipslerin daha genel durumlarıyla çalıştım. Bunu gürültülü verilerle doğru bir şekilde yapmak önemsiz değildir. Bu problem için, Halir ve Flusser tarafından Numerically stable direct least squares fitting of ellipses'da açıklanan algoritma oldukça iyi çalışır. Kağıt, Numpy'ye çevirmek için kolay olması gereken Matlab kodunu içerir. Belki de bu algoritmayı bir elipse sığdırmak için kullanabilir ve sonra iki eksenin ortalamasını yarıçap olarak alabilirsin. Makaledeki bazı referanslar da uygun dairelerden bahsediyor, bunlara bakmak isteyebilirsiniz.

+0

Gerçekten de öyle, önemsiz değil. Önerinize bakacak. – mimoralea

3

bir yukarı Bas Swinckels yazılan izleyin Bana yardım merkezi bulabilirsiniz Yukarıdaki kodu kullanarak bir elips

https://github.com/bdhammel/least-squares-ellipse-fitting

uydurma Halíř ve Flusser yöntemini uygulayan kodumu sonrası düşündü gibi aşağıdaki yöntem.

from ellipses import LSqEllipse 
import numpy as np 
import matplotlib.pyplot as plt 
from matplotlib.patches import Ellipse 

lsqe = LSqEllipse() 
lsqe.fit(data) 
center, width, height, phi = lsqe.parameters() 

plt.close('all') 
fig = plt.figure(figsize=(6,6)) 
ax = fig.add_subplot(111) 
ax.axis('equal') 
ax.plot(data[0], data[1], 'ro', label='test data', zorder=1) 

ellipse = Ellipse(xy=center, width=2*width, height=2*height, angle=np.rad2deg(phi), 
       edgecolor='b', fc='None', lw=2, label='Fit', zorder = 2) 
ax.add_patch(ellipse) 

plt.legend() 
plt.show() 

enter image description here

İlgili konular