0
$\begingroup$

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()

Catmull Rom curve

$\endgroup$
7
  • $\begingroup$ Welcome to Math.SE, @Emily. I'm not sure that here is the best place to post a request to code review. Perhaps another site around here, like Code Review, or Computer Science or one of the ones you cited posts from like Stack Overflow or Game Development. $\endgroup$ Commented May 5, 2019 at 18:53
  • $\begingroup$ Thank you. I think it's more the maths I'm struggling with though rather than the code. Would it help if I posted in plain text rather than code format? $\endgroup$
    – Emily
    Commented May 5, 2019 at 18:55
  • $\begingroup$ I think that using the code format it's the best way to check code. $\endgroup$ Commented May 5, 2019 at 18:57
  • $\begingroup$ I did not check your code, but I took a look at the comments on chf's post and one said "plasmacel: The approach from your link is simply wrong for the non-uniform case, it just happens to reduce to the correct solution if all time intervals are 1. And that's why you cannot match it with Barry-Goldman. The left and right time interval must appear separately in the equation, otherwise it can never work correctly! – Matthias Jul 5 '18 at 9:32" Is it relevant for your case? $\endgroup$ Commented May 5, 2019 at 22:33
  • $\begingroup$ Thank you for your help. I think plasmacel is not the original commentator and was providing another parameterisation in the comments so cfh, the OP's, answer still stands. I rescaled the tangent's by multiplying by (t2-t1). I'm just confused because it's managed to be done here but with the tangent being calculated as M[k] = (P[k+1] - P[k-1])/2 . I think I don't understand how a tangent can have two outputs, an x and y? $\endgroup$
    – Emily
    Commented May 6, 2019 at 7:27

1 Answer 1

-1
$\begingroup$

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!

$\endgroup$
1
  • $\begingroup$ This worked, thank you! $\endgroup$ Commented Dec 19, 2023 at 22:31

You must log in to answer this question.

Not the answer you're looking for? Browse other questions tagged .