fitting an image with 2D equation in python
I have an image and I want to fit it to 2D equation in order to extract nx and ny parameters. First I defined 2D function and residuals from fit then I read the image file and then I tried to fit it using leastsq method, this is my code:
#!/usr/bin/python import pyfits import numpy as np import numpy.random as npr import matplotlib.pyplot as plt import scipy.optimize nx=870 ny=901 # define 2D function def fun(nx,ny): n=(1+((nx**2+ny**2)**0.5/150)**2)**-3.7 return n vfun=np.vectorize(fun) nxlist=np.linspace(-nx,nx,870) nylist=np.linspace(-ny,ny,901) X,Y=np.meshgrid(nxlist,nylist) Z=vfun(X,Y) def residuals(p,y,nx,ny): nx,ny = p err = y-fun(nx,ny) return err def peval (nx,ny,p): nx,ny=p return fun(nx,ny) # read image file def image(): h = pyfits.open('image.fits') IM = h.data # copy the image data into a numpy (numerical python) array return IM y_true = image() y_meas = y_true+0.1*np.random.randn(ny,nx) # add noise colmap = plt.get_cmap('CMRmap') # load CMRmap colormap plt.imshow(y_meas, cmap=colmap, origin='lower') # plot image using gray colorbar plt.show() # initial values p0=[300,500] plsq = scipy.optimize.leastsq(residuals,p0,args=(y_meas,nx,ny)) print plsq
However, I got this error message
File "image_fit_test.py", line 51, in <module> plsq = scipy.optimize.leastsq(residuals,p0,args=(y_meas,nx,ny)) File "/.../anaconda/lib/python2.7/site-packages/scipy/optimize/minpack.py", line 364, in leastsq gtol, maxfev, epsfcn, factor, diag) minpack.error: Result from function call is not a proper array of floats.
Please can someone suggest any solution and where is this going wrong?
thank you in advance.
Simply replace residuals() with the following should solve your problem:
def residuals(p,y,nx,ny): nx,ny = p err = y-fun(nx,ny) return err.flatten()
Basically I suspect the function call of residuals(p0, meas, nx, ny) would return a 2d array of the shape of (nx, ny), which results in the minpack.error exception. You need to pass a 1d array (or a float) to leastsq().