python – 通过轨道数据拟合椭圆
我为它绕太阳运行的行星的(x,y,z)坐标生成了一堆数据.现在我想通过这些数据拟合椭圆.
我试图做的: 我基于五个参数创建了一个虚拟椭圆:半长轴&定义尺寸和偏差的偏心率形状和三个旋转椭圆的欧拉角.由于我的数据并不总是以原点为中心,因此我还需要翻译需要额外三个变量(dx,dy,dz)的椭圆. 我的问题在于最后一部分:最小化偏差并找到变量的值. 为了最大限度地减少偏差,我使用scipy.optimize.minimize来尝试近似最佳拟合变量,但它只是做得不够好: 我最好的一个看起来像是Here is an image,这是一个非常慷慨准确的初始猜测. (蓝色=数据,红色=适合) Here is the entire code.(无需数据,它会生成自己的假数据) 简而言之,我使用这个scipy函数: initial_guess = [0.3,0.2,0.1,0.7,3,0.0,-0.1,0.0] bnds = ((0.2,0.5),(0.1,0.3),(0,2*np.pi),(-0.5,(-0.3,0.3)) #reasonable bounds for the variables result = optimize.minimize(deviation,initial_guess,args=(data,),method='L-BFGS-B',bounds=bnds,tol=1e-8) #perform minimalisation semi_major,eccentricity,inclination,periapsis,longitude,dx,dz = result["x"] 要最小化此错误(或偏差)功能: def deviation(variables,data): """ This function calculates the cumulative seperation between the ellipse fit points and data points and returns it """ num_pts = len(data[:,0]) semi_major,dz = variables dummy_ellipse = generate_ellipse(num_pts,semi_major,dz,dz) deviations = np.zeros(len(data[:,0])) pair_deviations = np.zeros(len(data[:,0])) # Calculate separation between each pair of points for j in range(len(data[:,0])): for i in range(len(data[:,0])): pair_deviations[i] = np.sqrt((data[j,0]-dummy_ellipse[i,0])**2 + (data[j,1]-dummy_ellipse[i,1])**2 + (data[j,2]-dummy_ellipse[i,2])**2) deviations[j] = min(pair_deviations) # only pick the closest point to the data point j. total_deviation = sum(deviations) return total_deviation (我的代码可能有点乱,而且效率低,我是新手) 我可能在编码中犯了一些逻辑错误,但我认为它归结为scipy.minimize.optimize函数.我不知道它是如何工作的以及对它的期望.我还建议在处理这么多变量时尝试马尔可夫链蒙特卡罗.我确实看了一下司仪,但现在它已经超出了我的脑海. 解决方法
首先,你的目标函数中有一个拼写错误,它会阻止其中一个变量的优化:
dummy_ellipse = generate_ellipse(...,dz) 应该 dummy_ellipse = generate_ellipse(...,dz) 另外,将sqrt取出并最小化欧氏距离的平方和使得优化器在数值上更容易一些. 由于BFGS求解器假设的min(),你的目标函数也无处可辨,因此它的性能不是最理想的. 此外,从分析几何角度来解决问题可能有所帮助:3d中的椭圆被定义为两个方程的解 f1(x,z,p) = 0 f2(x,p) = 0 其中p是椭圆的参数.现在,为了使参数适合数据集,您可以尝试最小化 F(p) = sum_{j=1}^N [f1(x_j,y_j,z_j,p)**2 + f2(x_j,p)**2] 总和超过数据点. 更好的是,在这个问题公式中你可以使用optimize.leastsq,这在最小二乘问题中可能更有效. (编辑:李大同) 【声明】本站内容均来自网络,其相关言论仅代表作者个人观点,不代表本站立场。若无意侵犯到您的权利,请及时与联系站长删除相关内容! |