I fitted a plane into a set of points in R3 using minimization of L2 error (taken from this post: https://gist.github.com/amroamroamro/1db8d69b4b65e8bc66a6).

data = cluster_data[:ALLOWED_PEAKS, :]

grid_steps = 20

# regular grid covering the domain of the data
mn = np.min(data, axis=0)
mx = np.max(data, axis=0)

X, Y = np.meshgrid(np.linspace(mn[0], mx[0], grid_steps),
                   np.linspace(mn[1], mx[1], grid_steps))
XX = X.flatten()
YY = Y.flatten()

# best-fit linear plane (1st-order)
A = np.c_[data[:, 0], data[:, 1], np.ones(data.shape[0])]
C, _, _, _ = scipy.linalg.lstsq(A, data[:, 2])    # coefficients

# evaluate it on grid
Z = C[0] * X + C[1] * Y + C[2]

gci = int((grid_steps / 2) - 1 ) # grid_center_index

centroid = np.array([X[gci, gci], Y[gci, gci], Z[gci, gci]]).reshape(1,-1)  # take center point of grid steps

surface_norm = np.array([-C[0], -C[1], 1]).reshape(1,-1)

surface_norm_unit = normalize(surface_norm, norm='l2')

centroid_norm = centroid + surface_norm_unit

The plane looks like a great fit by visual inspection.

Now I want to find the normal vector of the plane. When looking at the code, the plane was determined by using an equation of the form

cz = ax + by + d

When assuming c = 1, and using Wikipedia's knowledge about normal vectors on planes

For a plane given by the equation ax+by+cz+d=0, the vector (a,b,c) is a normal.

my normal vector should be [ C[0], C[1], -1] or [ -C[0], -C[1], 1].

When plotting this normal vector, I receive a line that does not appear to be perpendicular on the plane.

enter image description here

Especially when looking from on side and making the plane appear to be a line (very thin, I hope you can see it)

enter image description here

Am I missing something fundamental about normal vectors? Is my assumption incorrect, that on the second image it should be orthogonal to the thin line? Here is my code for plotting:

# plot points and fitted surface using Matplotlib
fig = plt.figure(figsize=(10, 10))
ax = fig.gca(projection='3d')
ax.plot_surface(X, Y, Z, rstride=1, cstride=1, alpha=0.2)
ax.scatter(data[:, 0], data[:, 1], data[:, 2], c='r', s=50)
ax.scatter(centroid[0,0], centroid[0,1], centroid[0,2], c='b', s=100)
ax.scatter(centroid_norm[0,0], centroid_norm[0,1], centroid_norm[0,2], c='b', s=100)

When multiplying the calculated normal vector with a vector on the plane, the result is 2.77555756e-17, so nearly 0, which should be an indicator of both vectors being orthogonal.

  • 1
    $\begingroup$ Your plot is deceiving you because of the gross coordinate distortions. You can evaluate angles visually only when all three coordinates are drawn on identical scales. Try that and see whether the plane looks perpendicular to the normal vector. $\endgroup$
    – whuber
    Feb 1 '18 at 14:53
  • 1
    $\begingroup$ @whuber , thank you very much, I feel like a fool. If I want to visualize I need the same scale. Just to be safe, I tested it with scaled axes on the graph an it works great. If you want to add it as a solution, I would mark it as accepted. $\endgroup$ Feb 1 '18 at 16:43

I was experiencing the same issue.

whuber is right about this. I used plt3d.autoscale(False),it worked for me. This answer is related to this question: https://stackoverflow.com/q/3461869

For you this might work ax = fig.gca(projection='3d',autoscale_on=False)

ref: https://matplotlib.org/3.1.1/api/axes_api.html#matplotlib.axes.Axes

If that did not work, find an appropriate method to turn off auto-scaling from matplotlib site.

It messed up my plot but I think it could be fixed by limiting the number of points to be drawn or by increasing the size of the canvas. I did not bother about fixing it.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.