2014-10-02 78 views
2

嵌合我有x和y的一維陣列numpy的和我想再現y隨已知函數,以獲得「試用」。以下是我正在使用的代碼:曲線與已知功能numpy的

import numpy as np 
import matplotlib.pyplot as plt 
from scipy.optimize import curve_fit 

y = array([ 0.04022493, 0.04287536, 0.03983657, 0.0393201 , 0.03810298, 
    0.0363814 , 0.0331144 , 0.03074823, 0.02795767, 0.02413816, 
    0.02180802, 0.01861309, 0.01632699, 0.01368056, 0.01124232, 
    0.01005323, 0.00867196, 0.00940864, 0.00961282, 0.00892419, 
    0.01048963, 0.01199101, 0.01533408, 0.01855704, 0.02163586, 
    0.02630014, 0.02971127, 0.03511223, 0.03941218, 0.04280329, 
    0.04689105, 0.04960554, 0.05232003, 0.05487037, 0.05843364, 
    0.05120701]) 

x= array([ 0., 0.08975979, 0.17951958, 0.26927937, 0.35903916, 
    0.44879895, 0.53855874, 0.62831853, 0.71807832, 0.80783811, 
    0.8975979 , 0.98735769, 1.07711748, 1.16687727, 1.25663706, 
    1.34639685, 1.43615664, 1.52591643, 1.61567622, 1.70543601, 
    1.7951958 , 1.88495559, 1.97471538, 2.06447517, 2.15423496, 
    2.24399475, 2.33375454, 2.42351433, 2.51327412, 2.60303391, 
    2.6927937 , 2.78255349, 2.87231328, 2.96207307, 3.05183286, 
    3.14159265]) 

def func(x,beta): 
    return 1.0/(4.0*np.pi)*(1+beta*(3.0/2*np.cos(x)**2-1.0/2)) 

guesses = [20] 
popt,pcov = curve_fit(func,x,y,p0=guesses) 

y_fit = 1/(4.0*np.pi)*(1+popt[0]*(3.0/2*np.cos(x)**2-1.0/2)) 

plt.figure(1) 
plt.plot(x,y,'ro',x,y_fit,'k-') 
plt.show() 

代碼有效,但擬合完全關閉(請參見圖片)。任何想法爲什麼?

enter image description here 它看起來像使用式包含一個額外的參數,即,p

def func(x,beta,p): 
    return p/(4.0*np.pi)*(1+beta*(3.0/2*np.cos(x)**2-1.0/2)) 

guesses = [20,5] 
popt,pcov = curve_fit(func,x,y,p0=guesses) 

y_fit = func(angle_plot,*popt) 

plt.figure(2) 
plt.plot(x,y,'ro',x,y_fit,'k-') 
plt.show() 

print popt # [ 1.23341604 0.27362069] 

在哪一個是β-的POPT,哪一個爲p?

enter image description here

+2

它真正似乎您使用功能不符合數據以及 - 這可能是最合適的。請檢查括號。 – mdurant 2014-10-02 15:33:01

+0

@mdurant看到編輯 – diegus 2014-10-03 08:00:25

+0

第二個是'p',因爲你已經在'func' beta'後'定義它。 – 2014-10-03 15:34:01

回答

0

這將是一樣好,你可以得到(假設你得到的方程式權@mdurant建議),一個額外的截距項需要進一步完善適合:

def func(x,beta, icpt): 
    return 1.0/(4.0*np.pi)*(1+beta*(3.0/2*np.cos(x)**2-1.0/2))+icpt 

guesses = [20, 0] 
popt,pcov = curve_fit(func,x,y,p0=guesses) 

y_fit = func(x, *popt) 

plt.figure(1) 
plt.plot(x,y,'ro', x,y_fit,'k-') 
print popt #[ 0.33748816 -0.05780343] 

enter image description here

+0

它看起來我需要的第二個參數fit..In您打印是公測第一值的POPT,即0.3374? – diegus 2014-10-03 07:21:38

+0

是的,沒錯。它們的順序與'func'中的順序相同。 – 2014-10-03 15:32:13

2

這也許不是你想要的,但是,如果你只是想獲得一個不錯的選擇的數據,你可以使用np.polyfit

fit = np.polyfit(x,y,4) 
fit_fn = np.poly1d(fit) 
plt.scatter(x,y,label='data',color='r') 
plt.plot(x,fit_fn(x),color='b',label='fit') 
plt.legend(loc='upper left') 

enter image description here

注意fit給出的係數值,在這種情況下,一個四階多項式:

>>> fit 
array([-0.00877534, 0.05561778, -0.09494909, 0.02634183, 0.03936857]) 
+0

的polyfit肯定是好多了,然後我的功能,但我需要用我的功能.. – diegus 2014-10-03 08:01:44