2014-10-22 172 views
2

我試圖擬合二維高斯到一些灰度圖像數據,這是由一個二維數組給出的。 lmfit庫實現了一個易於使用的Model類,它應該能夠做到這一點。 不幸的是,文檔(http://lmfit.github.io/lmfit-py/model.html)只提供了一維擬合的例子。對於我的情況,我簡單地用兩個獨立變量構造lmfit模型。Python lmfit:擬合2D模型

下面的代碼似乎對我有效,但會導致scipy拋出「minpack.error:函數調用的結果不是一個適當的浮點數組。

Tom總結:如何將2D(x1,x2) - >(y)數據輸入到lmfit模型。

以下是我的方法: 所有內容都包含在GaussianFit2D類中,但以下是重要部分: 這就是高斯函數。文檔說有關用戶定義的函數

Of course, the model function will have to return an array that will be the same size as the data being modeled. Generally this is handled by also specifying one or more independent variables.

我真的不明白這是什麼意思應該,因爲對於給定的值X1,X2的唯一合理的結果是標量值。

def _function(self, x1, x2, amp, wid, cen1, cen2): 
    val = (amp/(np.sqrt(2*np.pi)*wid)) * np.exp(-((x1-cen1)**2+(x2-cen2)**2)/(2*wid**2)) 
    return val 

這裏生成的模型:

def _buildModel(self, **kwargs): 
    model = lmfit.Model(self._function, independent_vars=["x1", "x2"], 
         param_names=["amp", "wid", "cen1", "cen2"]) 
    return model 

這是獲取數據的功能,構建模型,而params,並呼籲lmfit符合():終於在這裏

def fit(self, data, freeX, **kwargs): 
    freeX = np.asarray(freeX, float) 
    model = self._buildModel(**kwargs) 
    params = self._generateModelParams(model, **kwargs) 

    model.fit(data, x1=freeX[0], x2=freeX[1], params=params) 

ANF這個適合的功能被稱爲:

data = np.asarray(img, float) 
    gaussFit = GaussianFit2D() 
    x1 = np.arange(len(img[0, :])) 
    x2 = np.arange(len(img[:, 0])) 
    fit = gaussFit.fit(data, [x1, x2]) 
+0

你需要'lmfit',或其他工具('curve_fit'或SciPy的'leastsq')也未嘗不可? – Evert 2014-10-22 07:57:01

+0

lmfit很迷人,因爲它很容易爲參數,約束等添加初始值......如果我不管理它,那麼scipy可能是第二次嘗試。 – m0e 2014-10-23 06:52:19

回答

0

好吧,與開發人員一起寫下來,並從他們那裏得到答案(在此感謝馬特)。

其基本思想是將所有輸入平坦化爲1D數據,從lmfit隱藏> 1維輸入。 以下是你如何做到的。 修改你的函數:

def function(self, x1, x2): 
     return (x1+x2).flatten() 

拼合你的二維輸入數組要適合:

... 
data = data.flatten() 
... 

修改兩個1D X變量,例如,你有他們的任意組合:

... 
x1n = [] 
x2n = [] 
    for i in x1: 
     for j in x2: 
       x1n.append(i) 
       x2n.append(j) 
x1n = np.asarray(x1n) 
x2n = np.asarray(x2n) 
... 

並將任何東西扔進鉗工:

model.fit(data, x1=x1n, x2=x2n, params=params) 
0

這是一個供您參考的例子,希望它能幫助您。

import numpy 
from lmfit import Model 

def gaussian(x, cenu, cenv, wid): 

    u = x[:, 0] 
    v = x[:, 1] 
    return (1/(2*numpy.pi*wid**2)) * numpy.exp(-(u-cenu)**2/(2*wid**2)-(v-cenv)**2/(2*wid**2)) 

data = numpy.empty((25,3)) 
x = numpy.arange(-2,3,1) 
y = numpy.arange(-2,3,1) 
xx, yy = numpy.meshgrid(x, y) 
data[:,0] = xx.flatten() 
data[:,1] = yy.flatten() 


data[:, 2]= gaussian(data[:,0:2],0,0,0.5) 

print 'xx\n', xx 
print 'yy\n',yy 
print 'data to be fit\n', data[:, 2] 

cu = 0.9 
cv = 0.5 
wid = 1 
gmod = Model(gaussian) 

gmod.set_param_hint('cenu', value=cu, min=cu-2, max=cu+2) 
gmod.set_param_hint('cenv', value=cv, min=cv -2, max=cv+2) 
gmod.set_param_hint('wid', value=wid, min=0.1, max=5) 
params = gmod.make_params() 
result = gmod.fit(data[:, 2], x=data[:, 0:2], params=params) 
print result.fit_report(min_correl=0.25) 
print result.best_values 
print result.best_fit 
+0

解釋你的答案也會有幫助。 – Fred 2015-04-16 16:11:35