splprepからスプライン誘導体を評価するには？

Question

私はBスプラインを使ってnumpyで閉じた曲線の曲率を計算したいと思っています。私はスムーズな結果を得るためにデータよりもスプライン表現の導関数を評価したいと思います。しかし、以下のコードはエラーを返します：私が正しくsplder機能を使用していた場合

Traceback (most recent call last): 
     File "<stdin>", line 16, in <module> 
     File "/Users/jfl/anaconda/lib/python2.7/site-packages/scipy/interpolate/fitpack.py", line 1212, in splder 
     c = (c[1:-1-k] - c[:-2-k]) * k/dt 
    TypeError: unsupported operand type(s) for -: 'list' and 'list' 

shell returned 1 

は、だから私は疑問に思って...

import numpy as np 
import scipy.interpolate as intplt 
import matplotlib.pyplot as plt 


t = np.linspace(0,2*np.pi,100) 
x = np.cos(t) 
y = np.sin(t) 


pts = np.vstack((x,y)) 
tck, u = intplt.splprep(pts, u=None,k=3, s=0.0, per=1) 
u_new = np.linspace(u.min(), u.max(), 1000) 
x_new, y_new = intplt.splev(u_new, tck, der=0) 

tck_der1 = intplt.splder(tck) 
tck_der2 = intplt.splder(tck_der1) 


xp, yp = intplt.splev(u_new, tck_der1, der=0) 
xpp, ypp = intplt.splev(u_new, tck_der2, der=0) 



plt.figure() 
plt.plot(x,y,".") 
plt.plot(x_new,y_new) 

plt.figure() 
curvature = np.abs(xp* ypp - yp* xpp)/np.power(xp** 2 + yp** 2, 3/2) 
plt.plot(u_new,curvature) 

plt.show() 

Answer 1

splderは（splrepと機能y = f(x)を補間）スカラースプラインをサポートしていますsplprepで得られたベクトルスプラインではありません。いずれにせよ、あなたはそれを必要としない：ちょうどsplevのderパラメータを使用します。

xp, yp = intplt.splev(u_new, tck, der=1) 
xpp, ypp = intplt.splev(u_new, tck, der=2) 

これは誘導体の一般的な有限差分評価ではありません。 splevはスプライン構造を実際に使用して計算します。これはあなたがしたいことです。

私はあなたが上記を試してみましたが、出力に不満だったと思います：

しかし、これは面白いものばかりmatplotlibのです。プロットウィンドウは、+ 9.998e-1から+ 9.998e-1 + 0.0006までです。言い換えれば、曲率はほぼ一定（1）なので、matplotlibは補間からくる必然的なノイズを増幅します。合理的なプロットウィンドウを設定するだけで、問題は消えます。

（または、素敵な非定数の曲率を得るために、あなたの例としてx = 2*np.cos(t)のようなものを使用。）