Сплайн в 3D не может быть дифференцирован из-за ошибки AttributeError - PullRequest
Новая
       14

Сплайн в 3D не может быть дифференцирован из-за ошибки AttributeError

0 голосов
/ 29 мая 2019

Я пытаюсь приспособить сглаживающий B-сплайн к некоторым данным, и я нашел это очень полезным post здесь. Однако мне нужен не только сплайн, но и его производные, поэтому я попытался добавить в пример следующий код: 

tck_der = interpolate.splder(tck, n=1)
x_der, y_der, z_der = interpolate.splev(u_fine, tck_der)

По какой-то причине это не работает из-за проблем с типом данных. Я получаю следующую трассировку: 

Traceback (most recent call last):
  File "interpolate_point_trace.py", line 31, in spline_example
    tck_der = interpolate.splder(tck, n=1)
  File "/home/user/anaconda3/lib/python3.7/site-packages/scipy/interpolate/fitpack.py", line 657, in splder
     return _impl.splder(tck, n)
   File "/home/user/anaconda3/lib/python3.7/site-packages/scipy/interpolate/_fitpack_impl.py", line 1206, in splder
     sh = (slice(None),) + ((None,)*len(c.shape[1:]))
 AttributeError: 'list' object has no attribute 'shape'

Причина этого, по-видимому, заключается в том, что второй аргумент кортежа tck содержит список массивов numpy. Я думал, что преобразование входных данных в массив numpy также поможет, но это не меняет типы данных tck.

Отражает ли это поведение ошибку в scipy, или вход неправильно сформирован? Я попытался вручную превратить список в массив: 

tck[1] = np.array(tck[1])

но это (что меня не удивило) также выдало ошибку: 

ValueError: operands could not be broadcast together with shapes (0,8) (7,1)

Есть идеи, в чем может быть проблема? Я использовал scipy до и на 1D сплайнах, функция splder работает просто отлично, поэтому я предполагаю, что это как-то связано со сплайном, являющимся линией в 3D.

------- редактировать --------

Вот минимальный рабочий пример: 

import numpy as np
import matplotlib.pyplot as plt
from scipy import interpolate
from mpl_toolkits.mplot3d import Axes3D

total_rad = 10
z_factor = 3
noise = 0.1

num_true_pts = 200
s_true = np.linspace(0, total_rad, num_true_pts)
x_true = np.cos(s_true)
y_true = np.sin(s_true)
z_true = s_true / z_factor

num_sample_pts = 80
s_sample = np.linspace(0, total_rad, num_sample_pts)
x_sample = np.cos(s_sample) + noise * np.random.randn(num_sample_pts)
y_sample = np.sin(s_sample) + noise * np.random.randn(num_sample_pts)
z_sample = s_sample / z_factor + noise * np.random.randn(num_sample_pts)

tck, u = interpolate.splprep([x_sample, y_sample, z_sample], s=2)
x_knots, y_knots, z_knots = interpolate.splev(tck[0], tck)
u_fine = np.linspace(0, 1, num_true_pts)
x_fine, y_fine, z_fine = interpolate.splev(u_fine, tck)

# this is the part of the code I inserted: the line under this causes the crash
tck_der = interpolate.splder(tck, n=1)
x_der, y_der, z_der = interpolate.splev(u_fine, tck_der)

# end of the inserted code

fig2 = plt.figure(2)
ax3d = fig2.add_subplot(111, projection='3d')
ax3d.plot(x_true, y_true, z_true, 'b')
ax3d.plot(x_sample, y_sample, z_sample, 'r*')
ax3d.plot(x_knots, y_knots, z_knots, 'go')
ax3d.plot(x_fine, y_fine, z_fine, 'g')
fig2.show()
plt.show()

1 Ответ

2 голосов
/ 02 июля 2019

Наткнулся на ту же проблему ...

Я обошел ошибку с помощью interpolate.splder(tck, n=1) и вместо этого использовал interpolate.splev(spline_ev, tck, der=1), который возвращает производные в точках spline_ev (см. Scipy Doku ).

Если вам нужен сплайн, я думаю, что вы можете использовать interpolate.splprep() снова.

Всего что-то вроде: 

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

points = np.random.rand(10,2) * 10

(tck, u), fp, ier, msg = interpolate.splprep(points.T, s=0, k=3, full_output=True)

spline_ev = np.linspace(0.0, 1.0, 100, endpoint=True)

spline_points = interpolate.splev(spline_ev, tck)

# Calculate derivative
spline_der_points = interpolate.splev(spline_ev, tck, der=1)
spline_der = interpolate.splprep(spline_der_points.T, s=0, k=3, full_output=True)


# Plot the data and derivative
fig = plt.figure()

plt.plot(points[:,0], points[:,1], '.-', label="points")
plt.plot(spline_points[0], spline_points[1], '.-', label="tck")
plt.plot(spline_der_points[0], spline_der_points[1], '.-', label="tck_der")

#   Show tangent
plt.arrow(spline_points[0][23]-spline_der_points[0][23], spline_points[1][23]-spline_der_points[1][23], 2.0*spline_der_points[0][23], 2.0*spline_der_points[1][23])

plt.legend()

plt.show()

РЕДАКТИРОВАТЬ:

Я также открыл проблему на Github , и согласно ev-br использование interpolate.splprep не рекомендуется, и вместо него следует использовать make_interp_spline / BSpline.

...