I've implemented Catmull-Rom in python and now I want to be able to get the normal of points so I'm first finding the tangent, to look something like this. For Catmull Rom I have 4 points P0, P1, P2 and P3. Using cfh's answer I've found the tangents of the two end points, P1 and P2
#m1 represents tangent at starting point P1
m1 = (t2-t1)*(((P1-P0)/(t1-t0))-((P2-P0)/(t2-t0))+((P2-P1)/(t2-t1)))
#m2 epresents tangent at ending point P2
m1 = (t2-t1)*(((P2-P1)/(t2-t1))-((P3-P1)/(t3-t1))+((P3-P2)/(t3-t2)))
I've then put these into the standard formula for a cubic spline following this post
cubic_spline = (2*t**3 - 3*t**2 + 1) * P1 + (t**3 - 2*t**2 + t) * m1 + (-2*t**3 + 3*t**2) * P2 + (t**3 - t**2) * m2
cubic_spline_deriv = (6*t**2 - 6*t)*P1 + (3*t**2 - 4*t + 1)*m1 + (-6*t**2 + 6*t)*P2 + (3*t**2 - 2*t)*m2
I think I'm massively misunderstanding the maths though.
I think that the output of cubic_spline_deriv is a tangent vector. However, when I try to plot it with plt.quiver I'm way off finding a tangent. I'd also expect cubic_spline to give the same output as C, but it doesn't.
Would it be at all possible for anyone to please take a look and see where I'm massively misunderstanding please?
My whole code is below:
#https://stackoverflow.com/questions/34894837/how-to-get-connect-two-part-of-curve-and-get-the-points-position-of-connecting-c
##https://stackoverflow.com/questions/37214786/emulating-excels-scatter-with-smooth-curve-spline-function-in-matplotlib-for
def CatmullRomSpline(P0, P1, P2, P3, nPoints=100):
"""
P0, P1, P2, and P3 should be (x,y) point pairs that define the Catmull-Rom spline.
nPoints is the number of points to include in this curve segment.
"""
# Convert the points to numpy so that we can do array multiplication
P0, P1, P2, P3 = map(np.array, [P0, P1, P2, P3])
# Calculate t0 to t4
alpha = 0.5
def tj(ti, Pi, Pj):
xi, yi = Pi
xj, yj = Pj
return ( ( (xj-xi)**2 + (yj-yi)**2 )**0.5 )**alpha + ti
t0 = 0
t1 = tj(t0, P0, P1)
t2 = tj(t1, P1, P2)
t3 = tj(t2, P2, P3)
# Only calculate points between P1 and P2
t = np.linspace(t1,t2,nPoints)
# Reshape so that we can multiply by the points P0 to P3
# and get a point for each value of t.
t = t.reshape(len(t),1)
#print(t)
A1 = (t1-t)/(t1-t0)*P0 + (t-t0)/(t1-t0)*P1
A2 = (t2-t)/(t2-t1)*P1 + (t-t1)/(t2-t1)*P2
A3 = (t3-t)/(t3-t2)*P2 + (t-t2)/(t3-t2)*P3
#print(A1)
#print(A2)
#print(A3)
B1 = (t2-t)/(t2-t0)*A1 + (t-t0)/(t2-t0)*A2
B2 = (t3-t)/(t3-t1)*A2 + (t-t1)/(t3-t1)*A3
C = (t2-t)/(t2-t1)*B1 + (t-t1)/(t2-t1)*B2
#m1 represents tangent at starting point P1
m1 = (t2-t1)*(((P1-P0)/(t1-t0))-((P2-P0)/(t2-t0))+((P2-P1)/(t2-t1)))
#m2 epresents tangent at ending point P2
m2 = (t2-t1)*(((P2-P1)/(t2-t1))-((P3-P1)/(t3-t1))+((P3-P2)/(t3-t2)))
cubic_spline = (2*t**3 - 3*t**2 + 1) * P1 + (t**3 - 2*t**2 + t) * m1 + (-2*t**3 + 3*t**2) * P2 + (t**3 - t**2) * m2
cubic_spline_deriv = (6*t**2 - 6*t)*P1 + (3*t**2 - 4*t + 1)*m1 + (-6*t**2 + 6*t)*P2 + (3*t**2 - 2*t)*m2
return C, cubic_spline, cubic_spline_deriv
def CatmullRomChain(P):
"""
Calculate Catmull Rom for a chain of points and return the combined curve.
"""
sz = len(P)
# The curve C will contain an array of (x,y) points.
#C is a list of x and y coordinates
C = []
Q = []
Qd = []
for i in range(sz-3):
c, cubic_spline, cubic_spline_deriv = CatmullRomSpline(P[i], P[i+1], P[i+2], P[i+3])
C.extend(c)
Q.extend(cubic_spline)
Qd.extend(cubic_spline_deriv)
return C, Q, Qd
# Define a set of points for curve to go through
Points = [[0,1.5],[2,2],[3,1],[4,0.5],[5,1],[6,2],[7,3]]
# Calculate the Catmull-Rom splines through the points
c, q, qd = CatmullRomChain(Points)
u = qd[0][0]
v = qd[0][1]
print(c)
print(q)
print(len(qd))
# Convert the Catmull-Rom curve points into x and y arrays and plot
x,y = zip(*c)
plt.plot(x,y)
# Plot the control points
px, py = zip(*Points)
plt.plot(px,py,'or')
plt.quiver(2, 2, u, v)
plt.show()
if catmull-rom have v0 v1 v2 v3..
float t2 = t * t;
float m0 = (v2 - v0) * .5f;
float m1 = (v3 - v1) * .5f;
tangent = (t2 - t) * 6f * v1 + (3f * t2 - 4f * t + 1f) * m0 + (-6f * t2 + 6f * t) * v2 + (3f * t2 - 2f * t) * m1;
if you want normalized >> tangent.normalized
if you want normal vector >> rotate right 90 degree
if tangent is Vector2, Very Simple >> normal = (tangent.y, -tangent.x)
Have a nice day!