如何解决ValueError：形状2,210和2,210不对齐：210dim 1！= 2dim 0
我想使用scipy.optimize.curve_fit将强度分布函数拟合到2D图像数据，并且无法在我的代码中找到错误：
# Define doughnut beam intensity distribution function
def doughnut(x,y,x0,y0,A,FWHM):
'''2D intensity distribution function of doughnut beams (DOI: 10.1126/science.aak9913,https://science.sciencemag.org/content/sci/suppl/2016/12/21/science.aak9913.DC1/Balzarotti_SM.pdf).
Parameters
----------
x,y : float
X and Y coordinates,orthogonal to beam axis
x0 : float
X offset
y0 : float
Y offset
A : float
Peak intensity
FWHM : float
Full width at half maximum
'''
return A*np.exp(1)*4*np.log(2)*(np.dot(x+x0,x+x0) + np.dot(y+y0,y+y0))/FWHM**2*np.exp(-4*np.log(2)*(np.dot(x+x0,y+y0))/FWHM**2)
# Read image file names
pathname = '/home/user/doughnut_beam/'
filenameList = [filename for filename in os.listdir(pathname)
if filename.endswith('.tif')]
# Open image files,fit doughnut beam intensity distribution function
for filename in filenameList:
img = Image.open(pathname + filename)
X,Y = img.size
xRange = np.arange(1,X+1)
yRange = np.arange(1,Y+1)
xGrid,yGrid = np.meshgrid(xRange,yRange)
xyGrid = np.vstack((xGrid.ravel(),yGrid.ravel())) # scipy.optimize.curve_fit requires 2xN-array
imgArray = np.array(img)
imgArrayFlat = imgArray.ravel() # Flatten 2D pixel data into 1D array for scipy.optimize.curve_fit
params_opt,params_cov = curve_fit(doughnut,xyGrid,imgArrayFlat)
这是Jupyter Notebook的输出：
---------------------------------------------------------------------------
ValueError Traceback (most recent call last)
<ipython-input-44-eaa3ebdb6469> in <module>()
17 imgArrayFlat = imgArray.ravel() # Flatten 2D pixel data into 1D array for scipy.optimize.curve_fit
18
---> 19 params_opt,imgArrayFlat)
/usr/lib/python3/dist-packages/scipy/optimize/minpack.py in curve_fit(f,xdata,ydata,p0,sigma,absolute_sigma,check_finite,bounds,method,jac,**kwargs)
749 # Remove full_output from kwargs,otherwise we're passing it in twice.
750 return_full = kwargs.pop('full_output',False)
--> 751 res = leastsq(func,Dfun=jac,full_output=1,**kwargs)
752 popt,pcov,infodict,errmsg,ier = res
753 cost = np.sum(infodict['fvec'] ** 2)
/usr/lib/python3/dist-packages/scipy/optimize/minpack.py in leastsq(func,args,Dfun,full_output,col_deriv,ftol,xtol,gtol,maxfev,epsfcn,factor,diag)
381 if not isinstance(args,tuple):
382 args = (args,)
--> 383 shape,dtype = _check_func('leastsq','func',func,n)
384 m = shape[0]
385 if n > m:
/usr/lib/python3/dist-packages/scipy/optimize/minpack.py in _check_func(checker,argname,thefunc,numinputs,output_shape)
25 def _check_func(checker,26 output_shape=None):
---> 27 res = atleast_1d(thefunc(*((x0[:numinputs],) + args)))
28 if (output_shape is not None) and (shape(res) != output_shape):
29 if (output_shape[0] != 1):
/usr/lib/python3/dist-packages/scipy/optimize/minpack.py in func_wrapped(params)
461 if transform is None:
462 def func_wrapped(params):
--> 463 return func(xdata,*params) - ydata
464 elif transform.ndim == 1:
465 def func_wrapped(params):
<ipython-input-43-3e0adae6fbe0> in doughnut(x,FWHM)
17 Full width at half maximum
18 '''
---> 19 return A*np.exp(1)*4*np.log(2)*(np.dot(x+x0,y+y0))/FWHM**2)
ValueError: shapes (2,210) and (2,210) not aligned: 210 (dim 1) != 2 (dim 0)
更新：由于某些原因，使用numpy.dot对函数定义中的（偏移）变量x+x0和y+y0进行平方运算不起作用。只需更改为**运算符即可得出正确的图：
# UPDATED: Define doughnut beam intensity distribution function
def doughnut(x,orthogonal to beam axis
x0 : float
X offset
y0 : float
Y offset
A : float
Peak intensity
FWHM : float
Full width at half maximum
'''
return A*np.exp(1)*4*np.log(2)*((x+x0)**2 + (y+y0)**2)/FWHM**2*np.exp(-4*np.log(2)*((x+x0)**2 + (y+y0)**2)/FWHM**2)
fig = plt.figure()
ax = fig.gca(projection='3d')
# Make data
X = np.arange(-10,10,0.25)
Y = np.arange(-10,0.25)
X,Y = np.meshgrid(X,Y)
Z = doughnut(X,Y,x0=0,y0=0,A=1.5,FWHM=7)
# Plot the surface
surf = ax.plot_surface(X,Z)
plt.show()
=> Plot
但是：现在，在尝试拟合数据时出现了一个新错误：
---------------------------------------------------------------------------
TypeError Traceback (most recent call last)
<ipython-input-61-eaa3ebdb6469> in <module>()
17 imgArrayFlat = imgArray.ravel() # Flatten 2D pixel data into 1D array for scipy.optimize.curve_fit
18
---> 19 params_opt,diag)
384 m = shape[0]
385 if n > m:
--> 386 raise TypeError('Improper input: N=%s must not exceed M=%s' % (n,m))
387 if epsfcn is None:
388 epsfcn = finfo(dtype).eps
TypeError: Improper input: N=5 must not exceed M=2
解决方法
这应该可以解决问题。看一下for_fitting函数，看看如何以curve_fit接受的方式打包所有内容。
import matplotlib.pyplot as plt
import numpy as np
from scipy.optimize import curve_fit
def doughnut(y,x,y0,x0,A,FWHM):
"""2D intensity distribution function of doughnut beams (DOI: 10.1126/science.aak9913,https://science.sciencemag.org/content/sci/suppl/2016/12/21/science.aak9913.DC1/Balzarotti_SM.pdf).
Parameters
----------
y,x : float
X and Y coordinates,orthogonal to beam axis
y0 : float
Y offset
x0 : float
X offset
A : float
Peak intensity
FWHM : float
Full width at half maximum
"""
return (
A
* np.e
* 4
* np.log(2)
* ((x + x0) ** 2 + (y + y0) ** 2)
/ FWHM ** 2
* np.exp(-4 * np.log(2) * ((x + x0) ** 2 + (y + y0) ** 2) / FWHM ** 2)
)
fig0,(ax0,ax1,ax2) = plt.subplots(1,3,sharex=True,sharey=True)
# Make data
X = np.arange(-10,10,0.25)
Y = np.arange(-10,0.25)
X,Y = np.meshgrid(X,Y)
true_params = (0,100,7)
Z = doughnut(Y,X,*true_params)
# Plot the surface
ax0.matshow(Z,extent=(-10,-10))
ax0.set_title("Ground Truth")
def for_fitting(xdata,FWHM):
yy,xx = xdata
return doughnut(yy,xx,FWHM).ravel()
noisy_data = np.random.poisson(Z) + np.random.randn(*Z.shape)
ax1.matshow(noisy_data,-10))
ax1.set_title("Noisy Data")
opt_params,cov = curve_fit(for_fitting,(Y,X),noisy_data.ravel(),p0=(0,1))
print(opt_params)
fit_Z = doughnut(Y,*opt_params)
ax2.matshow(fit_Z,-10))
ax2.set_title("Fit")
fig0.tight_layout()
plt.show()
版权声明：本文内容由互联网用户自发贡献，该文观点与技术仅代表作者本人。本站仅提供信息存储空间服务，不拥有所有权，不承担相关法律责任。如发现本站有涉嫌侵权/违法违规的内容， 请发送邮件至 dio@foxmail.com 举报，一经查实，本站将立刻删除。
