2017-02-17 1 views
0

splevscipyで使用して、いくつかの点でスプラインの派生を見つけようとしています。スプラインの派生物:scipy splev

import numpy as np 
from scipy.interpolate import splprep, splev 
import matplotlib.pyplot as plt 

# function to normalize each row 
def normalized(a, axis=-1, order=2): 
    l2 = np.atleast_1d(np.linalg.norm(a, order, axis)) 
    l2[l2==0] = 1 
    return a/np.expand_dims(l2, axis) 

# points and spline 
pts = np.array([[0,0],[1,1],[2,np.sqrt(2)],[4,2],[9,3]]) 
tck, u = splprep(pts.T, u=None, k=3, s=0.0, per=0) 

# compute new points and derivatives 
u_new = np.linspace(u.min(), u.max(), 5*u.shape[0]) 
x_new, y_new = splev(u_new, tck, der=0) 
xp_num, yp_num = splev(pts, tck, der=1) # numerical derivatices 
xp_the, yp_the= pts[1:,0], 0.5/np.sqrt(pts[1:,0]) # analytical derivatives 
R = normalized(yp_num/xp_num) 
X,Y = pts[1:,0], pts[1:,1] 
U,V = X + R[1:,0], Y+R[1:,1] 

私は与えられた点で接線をプロットしたいと思います:たとえば

plt.plot(x_new,y_new,'-b') 
plt.plot(pts[:,0],pts[:,1],'--ro') 
plt.quiver(X,Y,U,V, angles='xy', scale_units='xy') 

enter image description here

私は、これらの接線が間違っていると思います。私の理解はxp_numyp_numxyに関するスプラインの数値派生です。だからdy/dxを見つけるには、私は単にそれらを分けてください。これは正しいです?

結局、私はあなたが、少なくともあなたが投稿コードでそれらを使用していないので、(明らかに間違っ誘導体は)数値誘導体に関連していないthis

+1

wh正規化されていますか? –

+0

@PaulPanzer:各行を正規化します。私は 'normalize'の定義を追加しました。 – Mahdi

+0

奇妙なことに、あなたのコードは、プロット内に全くベクトルを作成しないように見えます。 'yp_num'と' xp_num'は1dですね。したがって、あなたの 'R'は1Xなので、' R [1:、0]をスライスすると、空の配列が得られるはずです。私はここに何かを逃していますか –

答えて

1

あなたの問題のような曲線の接線を見つけたいのですが。前者は確かに増分であるため、どのようなあなたのnormalized機能が本当に魔法何かをしない限り、明らかに間違っていることはあなたがxp_theによってyp_theを分割され、後者は、あなたの

とは対照的に、

dy 
-- 
dx 

を取得するために一定でなければなりません

dy 
---- 
x dx 

あなたはおそらく使用パラメトリック曲線

.  dy 
     -- 
dy dt 
-- = ---- 
dx dx 
     -- 
     dt 

ための式から来ていましたとなり、dx/dxが一定であることが見落とされる。私たちの最善を尽くすことの種類。

+0

驚くべきことに、私はパラメトリック曲線を考えました。ありがとう! – Mahdi

1

あなたR = normalized(yp_the/xp_the)ソース

が含まれていませんでした私はその後、私は正規化された派生

ためdelta_Yを変更し、矢筒に

import numpy as np 
from scipy.interpolate import splprep, splev 
import matplotlib.pyplot as plt 

# points and spline 
pts = np.array([[0,0],[1,1],[2,np.sqrt(2)],[4,2],[9,3]]) 
tck, u = splprep(pts.T, u=None, k=3, s=0.0, per=0) 

# compute new points and derivatives 
u_new = np.linspace(u.min(), u.max(), 5*u.shape[0]) 
x_new, y_new = splev(u_new, tck, der=0) 
xp_num, yp_num = splev(pts, tck, der=1) # numerical derivatices 
xp_the, yp_the= pts[1:,0], 0.5/np.sqrt(pts[1:,0]) # analytical derivatives 
#R = normalized(yp_the/xp_the) 
N = np.linalg.norm(np.array([xp_the, yp_the]), axis=0) 

X,Y = pts[1:,0], pts[1:,1] 
#U,V = X + R[1:,0], Y+R[1:,1] 
U,V = X + xp_the/N, Y + X*yp_the/N # delta Y = dy/dx * x 

plt.axes().set_aspect('equal', 'datalim') 

plt.plot(x_new,y_new,'-b') 
plt.plot(pts[:,0],pts[:,1],'--ro') 
#plt.quiver(X,Y,U,V, scale=10, pivot='mid')# angles='xy', scale_units=’xy’, scale=1 

plt.plot((X, U), (Y, V), '-k', lw=3) 
をあきらめたlinalg.norm

でそれをやりました

enter image description here

+0

私のミスについては申し訳ありません。私は 'xp_num'と' yp_num'を使って接線をプロットするつもりでした。 – Mahdi