0

I'm trying to plot a Bezier surface in matplotlib. I have my functions for x,y, and z (I didnt type them out, I had sympy automatically generate them) but whenever I run my code, I just get a blank plot. By the way, px, py, and pz refer to matrix of control points. Can anyone tell me what I did wrong and what I need to correct? Or maybe suggest another way to plot a bezier surface using python?

fig = plt.figure()
ax = fig.gca(projection='3d')

u = np.linspace(0, 1, 10) 
v = np.linspace(0, 1, 10)

x = px[0][0]*u**3*v**3 + px[0][1]*3*u**3*v**2*(-v + 1) + px[0][2]*3*u**3*v*   (-v + 1)**2 + px[0][3]*u**3*(-v + 1)**3 + \
px[1][0]*3*u**2*v**3*(-u + 1) + px[1][1]*9*u**2*v**2*(-u + 1)*(-v + 1) + px[1][2]*9*u**2*v*(-u + 1)*(-v + 1)**2 + \
px[1][3]*3*u**2*(-u + 1)*(-v + 1)**3 + px[2][0]*3*u*v**3*(-u + 1)**2 + px[2][1]*9*u*v**2*(-u + 1)**2*(-v + 1) + \
px[2][2]*9*u*v*(-u + 1)**2*(-v + 1)**2 + px[2][3]*3*u*(-u + 1)**2*(-v + 1)**3 + px[3][0]*v**3*(-u + 1)**3 + \
px[3][1]*3*v**2*(-u + 1)**3*(-v + 1) + px[3][2]*3*v*(-u + 1)**3*(-v + 1)**2 + px[3][3]*(-u + 1)**3*(-v + 1)**3

y = py[0][0]*u**3*v**3 + py[0][1]*3*u**3*v**2*(-v + 1) + py[0][2]*3*u**3*v*(-v + 1)**2 + py[0][3]*u**3*(-v + 1)**3 +\
py[1][0]*3*u**2*v**3*(-u + 1) + py[1][1]*9*u**2*v**2*(-u + 1)*(-v + 1) + py[1][2]*9*u**2*v*(-u + 1)*(-v + 1)**2 + \
py[1][3]*3*u**2*(-u + 1)*(-v + 1)**3 + py[2][0]*3*u*v**3*(-u + 1)**2 + py[2][1]*9*u*v**2*(-u + 1)**2*(-v + 1) + \
py[2][2]*9*u*v*(-u + 1)**2*(-v + 1)**2 + py[2][3]*3*u*(-u + 1)**2*(-v + 1)**3 + py[3][0]*v**3*(-u + 1)**3 + \
py[3][1]*3*v**2*(-u + 1)**3*(-v + 1) + py[3][2]*3*v*(-u + 1)**3*(-v + 1)**2 + py[3][3]*(-u + 1)**3*(-v + 1)**3

z = pz[0][0]*u**3*v**3 + pz[0][1]*3*u**3*v**2*(-v + 1) + pz[0][2]*3*u**3*v*(-v + 1)**2 + pz[0][3]*u**3*(-v + 1)**3 + \
pz[1][0]*3*u**2*v**3*(-u + 1) + pz[1][1]*9*u**2*v**2*(-u + 1)*(-v + 1) + pz[1][2]*9*u**2*v*(-u + 1)*(-v + 1)**2 +\
pz[1][3]*3*u**2*(-u + 1)*(-v + 1)**3 + pz[2][0]*3*u*v**3*(-u + 1)**2 + pz[2][1]*9*u*v**2*(-u + 1)**2*(-v + 1) + \
pz[2][2]*9*u*v*(-u + 1)**2*(-v + 1)**2 + pz[2][3]*3*u*(-u + 1)**2*(-v + 1)**3 + pz[3][0]*v**3*(-u + 1)**3 + \
pz[3][1]*3*v**2*(-u + 1)**3*(-v + 1) + pz[3][2]*3*v*(-u + 1)**3*(-v + 1)**2 + pz[3][3]*(-u + 1)**3*(-v + 1)**3

ax.plot_surface(x, y, z, rstride=4, cstride=4, color='b')

plt.show()
Pablo
  • 4,821
  • 12
  • 52
  • 82

1 Answers1

1

Here's a couple of examples to plot bezier curves with matplotlib:

EXAMPLE1

From stackoverflow

import numpy as np
from scipy.misc import comb
from matplotlib import pyplot as plt


def bernstein_poly(i, n, t):
    return comb(n, i) * (t**(n - i)) * (1 - t)**i


def bezier_curve(points, nTimes=1000):
    nPoints = len(points)
    xPoints = np.array([p[0] for p in points])
    yPoints = np.array([p[1] for p in points])

    t = np.linspace(0.0, 1.0, nTimes)

    polynomial_array = np.array(
        [bernstein_poly(i, nPoints - 1, t) for i in range(0, nPoints)])

    xvals = np.dot(xPoints, polynomial_array)
    yvals = np.dot(yPoints, polynomial_array)

    return xvals, yvals


if __name__ == "__main__":
    nPoints = 4
    points = np.random.rand(nPoints, 2) * 200
    xpoints = [p[0] for p in points]
    ypoints = [p[1] for p in points]

    xvals, yvals = bezier_curve(points, nTimes=1000)
    plt.plot(xvals, yvals)
    plt.plot(xpoints, ypoints, "ro")
    for nr in range(len(points)):
        plt.text(points[nr][0], points[nr][1], nr)

    plt.show()

EXAMPLE2

From matplotlib tutorial

import matplotlib.pyplot as plt
from matplotlib.path import Path
import matplotlib.patches as patches

verts = [
    (0., 0.),  # P0
    (0.2, 1.), # P1
    (1., 0.8), # P2
    (0.8, 0.), # P3
    ]

codes = [Path.MOVETO,
         Path.CURVE4,
         Path.CURVE4,
         Path.CURVE4,
         ]

path = Path(verts, codes)

fig = plt.figure()
ax = fig.add_subplot(111)
patch = patches.PathPatch(path, facecolor='none', lw=2)
ax.add_patch(patch)

xs, ys = zip(*verts)
ax.plot(xs, ys, 'x--', lw=2, color='black', ms=10)

ax.text(-0.05, -0.05, 'P0')
ax.text(0.15, 1.05, 'P1')
ax.text(1.05, 0.85, 'P2')
ax.text(0.85, -0.05, 'P3')

ax.set_xlim(-0.1, 1.1)
ax.set_ylim(-0.1, 1.1)
plt.show()

If you ask me... my favourite way to think about parametric curves is using their matrix form, like showed in this website

EDIT:

If you want to extend it to the 3d case, here's a working example:

import matplotlib as mpl
import numpy as np
from scipy.misc import comb
from matplotlib import pyplot as plt
from mpl_toolkits.mplot3d import Axes3D


def bernstein_poly(i, n, t):
    return comb(n, i) * (t**(n - i)) * (1 - t)**i


def bezier_curve(points, nTimes=1000):
    nPoints = len(points)
    xPoints = np.array([p[0] for p in points])
    yPoints = np.array([p[1] for p in points])
    zPoints = np.array([p[2] for p in points])

    t = np.linspace(0.0, 1.0, nTimes)

    polynomial_array = np.array(
        [bernstein_poly(i, nPoints - 1, t) for i in range(0, nPoints)])

    xvals = np.dot(xPoints, polynomial_array)
    yvals = np.dot(yPoints, polynomial_array)
    zvals = np.dot(zPoints, polynomial_array)

    return xvals, yvals, zvals


if __name__ == "__main__":
    nPoints = 4
    points = np.random.rand(nPoints, 3) * 200
    xpoints = [p[0] for p in points]
    ypoints = [p[1] for p in points]
    zpoints = [p[2] for p in points]

    xvals, yvals, zvals = bezier_curve(points, nTimes=1000)

    fig = plt.figure()
    ax = fig.gca(projection='3d')
    ax.plot(xvals, yvals, zvals, label='bezier')
    ax.plot(xpoints, ypoints, zpoints, "ro")
    for nr in range(len(points)):
        ax.text(points[nr][0], points[nr][1], points[nr][2], nr)

    plt.show()
Community
  • 1
  • 1
BPL
  • 9,632
  • 9
  • 59
  • 117
  • Hey, thanks this is really helpful. I actually thought about using path after seeing it on matplotlib, but I wasn't sure how to extend it to 3 dimension surfaces. Can you explain how to extend either example to 3 dimensional surfaces? Is it as simple as just adding a "zvals" for the first example? – Jason Irukulapati Aug 15 '16 at 19:58
  • @JasonIrukulapati Sure, I've edited my answer and provided you a simple example at the end of it. I always encourage learning the 2d case and once you feel comfortable with it extending to 3d it becomes trivial – BPL Aug 15 '16 at 20:15
  • @JasonIrukulapati You're welcome :) , just be aware the right way to thanks somebody in StackOverflow is marking the answer it helped you out as valid. – BPL Aug 15 '16 at 20:40