将平面适配到点云中(3D)

2022-05-23 00:00:00 python arrays slice curve-fitting tuples

问题描述

我有一个包含多列数据的数据文件,我想从这个数据文件中提取3列(表示坐标),并将它们放到另一个文件中,然后使用新创建的文件,我想使用scipy.Optimize.curve_fit来拟合一个平面或曲面(或您想怎么称呼它)。以下是我的代码:

# -*- coding: utf-8 -*-

from pylab import *
import matplotlib.pyplot as plt 
import numpy as np
from scipy import stats
from scipy.optimize import curve_fit



### processing function    
def store(var,textfile):
    data=loadtxt(textfile,skiprows=1)
    p0=[]
    p1=[]
    p2=[]
    for i in range(0,len(data)):
        p0.append(float(data[i,2]))
        p1.append(float(data[i,3]))
        p2.append(float(data[i,4]))
    var.append(p0)
    var.append(p1)
    var.append(p2)

#extracting the data from a textfile
datafile1='cracks_0101005_5k_tensionTestCentreCrack_l0.001a0_r0.01.txt'
a1=[]
store(a1, datafile1)


rcParams.update({'legend.numpoints':1,'font.size': 20,'axes.labelsize':25,'xtick.major.pad':10,'ytick.major.pad':10,'legend.fontsize':14})
lw=2
ms=10

#fitting a surface(curve) into the data

def func(data, m, n, o): 
    return m*data[:,0] + n*data[:,2] + o 

guess=(1,1,1)

params, pcov = curve_fit(func, a1[::2, :2], a1[:,1], guess)
print (params)

我收到以下错误消息:

Traceback (most recent call last):
  File "fitcurve.py", line 41, in <module>
    params, pcov = curve_fit(func, a1[::2, :2], a1[:,1], guess)
TypeError: list indices must be integers, not tuple

您能告诉我我做错了什么吗?

只是为了更清楚地说明: 我试着把Y作为我的相依函数,所以它是X和Z的函数。 显然a1[]是列表而不是数组,对吗? 但是,即使我将其更改为数组Myarray=np.asarray(a1),我也会收到其他一些奇怪的消息。 如果有人能帮我理解一下这个问题,我将不胜感激。

干杯


解决方案

以下是我在代码中注意到的可能错误:

您希望将y作为x,z的函数,因此要发送的X数组可能是a1[:, ::2]。但这意味着func已经获得了m,2数组,所以这里它必须是return m*data[:,0] + n*data[:,1] + o 尽管如此,我认为它应该是两个参数的拟合,而不是三个。您可以根据相应的结果计算可能的m, n, o

import matplotlib 
matplotlib.use('Qt5Agg')
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
import numpy as np
from matplotlib import style
from random import random
style.use('ggplot')

from scipy.optimize import curve_fit

""" make some data and save to file """
data = []
a, b, s = 1.233, -2.17, 5.2
for i in range(88):
    x = 10 * (2 * random() - 1)
    y = 10 * (2 * random() - 1)
    z = a * x +  b * y+ s * (2 * random() - 1) * 0.5
    data += [[x, y, z]]
data = np.array(data)
np.savetxt("data.txt", data)

""" get the data and use unpack to directly write into x,y,z variables"""
xData, yData, zData = np.loadtxt("data.txt", unpack=True)
"""...while I actally need the packed version as well, so I could load again"""
#allData = np.loadtxt("data.txt")
"""...or..."""
allData = np.array( zip(xData, yData, zData) )

def func(data, m, o): 
    return m * data[:,0] + o * data[:, 1] 
    
guess = (1, 1)
params, pcov = curve_fit(func, allData[:, ::2], allData[:,1], guess)

""" showing data and fit result"""
x = np.linspace(-10, 10, 10)
y = np.linspace(-10, 10, 10)
X, Y = np.meshgrid(x, y)
Z = -params[0] / params[1] * X + 1 / params[1] * Y

fig1 = plt.figure(1)
ax = fig1.add_subplot( 1, 1, 1, projection='3d')
ax.scatter(xData, yData, zData)
ax.plot_wireframe(X, Y, Z, color='#4060a6', alpha=0.6 )
ax.set_title(
    "({:1.2f},{:1.2f})".format(
        -params[0] / params[1], 1 / params[1]
    )
)
plt.show()

请注意,当您安装y = m * x + o * z时,我绘制了z = a * x + b * yMITb = 1 / oa = -m / o的线框,即n = 1。您可以相应地重新调整m, n, o的比例。

相关文章