我如何拟合3D数据 [英] how do I fit 3D data

问题描述

我有一个3D点列表，我想将其适合一个球体

I have a list of 3D points, and I want to fit the to a sphere:

R^2 = (x-x0)^2 + (y-y0)^2 + (z-z0)^2

所以我想，我要表达z并使用4个参数(x0，y0，z0和R)拟合2D数据:

So I thought, I'd express z, and fit 2D data with 4 parameters (x0, y0, z0 and R):

z = sqrt(R^2 - (x-x0)^2 - (y-y0)^2) + z0

这是一个代码(它是较大项目的一部分):

Here's a code (it is a part of larger project):

#!/usr/bin/python

from scipy import *
from scipy.optimize import leastsq

Coordinates = load("data.npy")

xyCoords = Coordinates[:, [0, 1]]
zCoords  = Coordinates[:, 2]

p0 = [149.33499, 148.95999, -218.84893225608857, 285.72893713890107]

fitfunc = lambda p, x: sqrt(p[3]**2 - (x[0] - p[0])**2 - (x[1] - p[1])**2) + x[2]
errfunc = lambda p, x, z: fitfunc(p, x) - z
p1, flag = leastsq(errfunc, p0, args=(xyCoords, zCoords))

print p1

我得到了错误:

ValueError: operands could not be broadcast together with shapes (2) (1404)

这里是 data.npy 的链接.

推荐答案

您需要正确定义fitfunc:

fitfunc = lambda p, x: sqrt(p[3]**2 - (x[:, 0] - p[0])**2 - (x[:, 1] - p[1])**2) + p[2]

我认为您的方法不是很健壮，因为当您采用sqrt时，有两种解决方案，一种是肯定的，一种是否定的，而您只是在考虑肯定的.因此，除非您的所有观点都在一个领域的上半部，否则您的方法将行不通.最好将r设置为fitfunc:

I don't think that your approach is very robust, because when you take the sqrt there are two solutions, one positive, one negative, and you are only considering the positive. So unless all your points are on the top half of a sphere, your approach will not work. It is probably better to make r your fitfunc:

import numpy as np
from scipy.optimize import leastsq

# test data: sphere centered at 'center' of radius 'R'
center = (np.random.rand(3) - 0.5) * 200
R = np.random.rand(1) * 100
coords = np.random.rand(100, 3) - 0.5
coords /= np.sqrt(np.sum(coords**2, axis=1))[:, None]
coords *= R
coords += center

p0 = [0, 0, 0, 1]

def fitfunc(p, coords):
    x0, y0, z0, R = p
    x, y, z = coords.T
    return np.sqrt((x-x0)**2 + (y-y0)**2 + (z-z0)**2)

errfunc = lambda p, x: fitfunc(p, x) - p[3]

p1, flag = leastsq(errfunc, p0, args=(coords,))

>>> center
array([-39.77447344, -69.89096249,  44.46437355])
>>> R
array([ 69.87797469])
>>> p1
array([-39.77447344, -69.89096249,  44.46437355,  69.87797469])

这篇关于我如何拟合3D数据的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！