我正在尝试使用 Scipy leastsq 为二维中的一组测量点坐标找到“正方形”网格的最佳拟合(实验点大约在正方形网格上).
网格的参数是间距(x 和 y 相等)、中心位置(center_x 和 center_y)和旋转(学位)。
我定义了一个误差函数，计算每对点(实验 vs 理想网格)的欧氏距离并取平均值。我想通过 leastsq 最小化此函数，但出现错误。
函数定义如下:
import numpy as np
from scipy.optimize import leastsq
def get_spot_grid(shape, pitch, center_x, center_y, rotation=0):
x_spots, y_spots = np.meshgrid(
(np.arange(shape[1]) - (shape[1]-1)/2.)*pitch,
(np.arange(shape[0]) - (shape[0]-1)/2.)*pitch)
theta = rotation/180.*np.pi
x_spots = x_spots*np.cos(theta) - y_spots*np.sin(theta) + center_x
y_spads = x_spots*np.sin(theta) + y_spots*np.cos(theta) + center_y
return x_spots, y_spots
def get_mean_distance(x1, y1, x2, y2):
return np.sqrt((x1 - x2)**2 + (y1 - y2)**2).mean()
def err_func(params, xe, ye):
pitch, center_x, center_y, rotation = params
x_grid, y_grid = get_spot_grid(xe.shape, pitch, center_x, center_y, rotation)
return get_mean_distance(x_grid, y_grid, xe, ye)
这是实验坐标:
xe = np.array([ -23.31, -4.01, 15.44, 34.71, -23.39, -4.10, 15.28, 34.60, -23.75, -4.38, 15.07, 34.34, -23.91, -4.53, 14.82, 34.15]).reshape(4, 4)
ye = np.array([-16.00, -15.81, -15.72, -15.49, 3.29, 3.51, 3.90, 4.02, 22.75, 22.93, 23.18, 23.43, 42.19, 42.35, 42.69, 42.87]).reshape(4, 4)
我尝试以这种方式使用 leastsq:
leastsq(err_func, x0=(19, 12, 5, 0), args=(xe, ye))
但我收到以下错误:
---------------------------------------------------------------------------
TypeError Traceback (most recent call last)
<ipython-input-19-ee91cf6ce7d6> in <module>()
----> 1 leastsq(err_func, x0=(19, 12, 5, 0), args=(xe, ye))
C:\Anaconda\lib\site-packages\scipy\optimize\minpack.pyc in leastsq(func, x0, args, Dfun, full_output, col_deriv, ftol, xtol, gtol, maxfev, epsfcn, factor, diag)
369 m = shape[0]
370 if n > m:
--> 371 raise TypeError('Improper input: N=%s must not exceed M=%s' % (n, m))
372 if epsfcn is None:
373 epsfcn = finfo(dtype).eps
TypeError: Improper input: N=4 must not exceed M=1
我不知道这里有什么问题:(
请您参考如下方法:
因为 leastsq 函数假设 err_function 返回一个残差数组 docs以这种方式编写 err_function 有点困难，为什么不使用另一个 scipy 的函数 - minimize .然后你添加你的指标 - 你已经拥有的错误函数并且它有效。但是，我认为 get_spot_grid 函数中还有一处错别字(y_spots 与 y_spads)。完整代码:
import numpy as np
from scipy.optimize import leastsq, minimize
def get_spot_grid(shape, pitch, center_x, center_y, rotation=0):
x_spots, y_spots = np.meshgrid(
(np.arange(shape[1]) - (shape[1]-1)/2.)*pitch,
(np.arange(shape[0]) - (shape[0]-1)/2.)*pitch)
theta = rotation/180.*np.pi
x_spots = x_spots*np.cos(theta) - y_spots*np.sin(theta) + center_x
y_spots = x_spots*np.sin(theta) + y_spots*np.cos(theta) + center_y
return x_spots, y_spots
def get_mean_distance(x1, y1, x2, y2):
return np.sqrt((x1 - x2)**2 + (y1 - y2)**2).mean()
def err_func(params, xe, ye):
pitch, center_x, center_y, rotation = params
x_grid, y_grid = get_spot_grid(xe.shape, pitch, center_x, center_y, rotation)
return get_mean_distance(x_grid, y_grid, xe, ye)
xe = np.array([-23.31, -4.01, 15.44, 34.71, -23.39, -4.10, 15.28, 34.60, -23.75, -4.38, 15.07, 34.34, -23.91, -4.53, 14.82, 34.15]).reshape(4, 4)
ye = np.array([-16.00, -15.81, -15.72, -15.49, 3.29, 3.51, 3.90, 4.02, 22.75, 22.93, 23.18, 23.43, 42.19, 42.35, 42.69, 42.87]).reshape(4, 4)
# leastsq(err_func, x0=(19, 12, 5, 0), args=(xe, ye))
minimize(err_func, x0=(19, 12, 5, 0), args=(xe, ye))
