2017-02-17 1 views
0

J'essaie de trouver des dérivés d'une spline en plusieurs points en utilisant splev en scipy. Par exemple:Dérivés d'une spline: `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] 

Je voudrais tracer les tangentes au point donné:

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

Je pense que ces tangentes sont fausses. Ma compréhension était que xp_num et yp_num sont des dérivés numériques de la spline par rapport à x et y. Donc, pour trouver dy/dx, je devrais simplement les diviser. Est-ce correct?

Finalement, je voudrais trouver tangentes d'une courbe comme this

+1

wh à fait «normalisé» faire? –

+0

@PaulPanzer: Il normalise chaque ligne. J'ai ajouté la définition de 'normaliser'. – Mahdi

+0

Etrange, votre code me semble comme s'il ne devait pas du tout créer de vecteurs dans l'intrigue. 'yp_num' et' xp_num' sont 1d, n'est-ce pas? Donc votre 'R' devrait être 1 x quelque chose, donc quand vous triez' R [1:, 0] vous devriez obtenir un tableau vide. Est-ce que j'ai râté quelque chose? –

Répondre

1

Votre problème (les dérivés évidemment faux) ne sont pas liés à la dérivée numérique puisque vous ne les utilisez pas au moins dans le code affiché . Ce qui est clairement mal à moins que votre fonction normalized fait quelque chose vraiment magique est votre divisant yp_the par xp_the puisque le premier est en effet l'augmentation, celle-ci est pas doit être constante pour obtenir

dy 
-- 
dx 

par opposition à votre

dy 
---- 
x dx 

Vous viendriez probablement de la formule pour une courbe paramétrique

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

utilisé t=x et ensuite négligé que dx/dx est constante. Genre de chose qui arrive aux meilleurs d'entre nous.

+0

Génial, c'est exactement ce que je pensais de la courbe paramétrique! Merci! – Mahdi

1

vous n'avez pas inclus votre R = normalized(yp_the/xp_the) source de

Je l'ai fait avec linalg.norm

puis j'ai changé le delta_y pour le dérivé normalisé

et donné sur carquois

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) 

enter image description here

+0

Désolé pour mon erreur. Je voulais tracer les tangentes en utilisant 'xp_num' et' yp_num'. – Mahdi