I would like to compute the gradient of a spline and visualize it in a polar plot.
I use tck, u = splprep() to obtain the spline with the cartesian coordinates and splev(u, tck, der=1) to compute its partial derivatives in x and y direction, respectively. Then, I compute the endpoints of the arrows that are supposed to visualize the gradient and convert them to polar coordinates.
The plot looks good on the first glance but if I compare the estimated direction of the gradient to the analytical solution, there are significant differences even if I increase the number of points.
MWE
from matplotlib import pyplot as plt
import numpy as np
from scipy.interpolate import splprep, splev
if __name__ == '__main__':
N = 11 # number of samples
# x = np.arange(0, N) # [Update] This did not decrease the step size for increasing N
x = np.linspace(0, 10, N)
y = np.sin(x)
tck, u = splprep([x, y], s=0) # spline
theta, r = np.arctan2(y, x), np.hypot(x, y) # convert to polar
gradient = splev(u, tck, der=1) # compute first derivative
# normalize
gradient = gradient / np.hypot(gradient[0], gradient[1])
# compare numerical and analytical solution
# direction = np.arctan2(gradient[1], gradient[0]) # [Update] this was wrong
slope = gradient[1] / gradient[0]
print(np.cos(x) - slope) # cos(x) should be the analytical solution
endpoints_x = x + gradient[0]
endpoints_y = y + gradient[1]
# convert cartesian endpoints to polar
endpoints_theta, endpoints_r = np.arctan2(endpoints_y, endpoints_x),\
np.hypot(endpoints_x, endpoints_y)
fig, ax = plt.subplots(1, 1, subplot_kw=dict(polar=True))
plt.scatter(theta, r, marker='o')
plt.plot(np.stack((theta, endpoints_theta)), np.stack((r, endpoints_r)), 'r')
plt.show()
Screenshot
Update
I found two mistakes. First, the step size in x was not decreasing by increasing N as I used np.arange(0, N). Second, I expected the numerical solution to be arctan2(gradient[1], gradient[0]) but it is simply the slope gradient[1] / gradient[0] in cartesian coordinates.
Everything works as expected now.
