我有以下三个点在三维空间与他们各自的法向量。
A、B、C是位置，而N_A、N_B和N_C是它们各自的法向矢量。

A = np.array([  348.92065834, -1402.3305998,     32.69313966])
N_A = np.array([-0.86925426,  0.02836434, -0.49355091])

B =  np.array([282.19332067,  82.52027998,  -5.92595371])
N_B = np.array([-0.82339849,  0.43041935,  0.3698028])

C = np.array([247.37475615,  -3.70129865, -22.10494737])
N_C = np.array([-0.83989222, 0.23796899, 0.48780305])

三个点几乎在一个平面上，但两个最近的点B和C之间的方向略有变化。从那里到点A，我假设在X-Y坐标中可能存在曲率，如x1c 0d1x所示，但在X-Z坐标中可能存在抛物面，如

所示
假设的条件是，如果曲率非常明显，圆柱体可以拟合这三个点。考虑到它们各自在X和Z坐标上的法向量，B和C的法向量朝下，而A的法向量朝上。所以，总的来说，它可能是一个抛物面。问题是如何拟合它们，考虑它们的法向量。如果不可能，然后如何从X和Y方向用一定曲率来拟合它们。
这是一个情节的代码

import numpy as np
import matplotlib.pyplot as plt
def set_axes_radius(ax, origin, radius):

    ax.set_xlim3d([origin[0] - radius, origin[0] + radius])
    ax.set_ylim3d([origin[1] - radius, origin[1] + radius])
    ax.set_zlim3d([origin[2] - radius, origin[2] + radius])

def set_axes_equal(ax, zoom=1.):
    '''
        Make axes of 3D plot have equal scale so that spheres appear as spheres,
        cubes as cubes, etc..  This is one possible solution to Matplotlib's
        ax.set_aspect("equal") and ax.axis("equal") not working for 3D.
        input:
          ax:   a matplotlib axis, e.g., as output from plt.gca().

    '''

    limits = np.array([
        ax.get_xlim3d(),
        ax.get_ylim3d(),
        ax.get_zlim3d(),
    ])

    origin = np.mean(limits, axis=1)
    radius = 0.5 * np.max(np.abs(limits[:, 1] - limits[:, 0])) / zoom
    set_axes_radius(ax, origin, radius)
%matplotlib qt 

# positions and their respective normal vectors

A = np.array([  348.92065834, -1402.3305998,     32.69313966])
N_A = np.array([-0.86925426,  0.02836434, -0.49355091])

B =  np.array([282.19332067,  82.52027998,  -5.92595371])
N_B = np.array([-0.82339849,  0.43041935,  0.3698028])

C = np.array([247.37475615,  -3.70129865, -22.10494737])
N_C = np.array([-0.83989222, 0.23796899, 0.48780305])

# A plane is given by

# a*x + b*y + c*z + d = 0

# where (a, b, c) is the normal.

# If the point (x, y, z) lies on the plane, then solving for d yield:

# d = -(a*x + b*y + c*z)

d_A = -np.sum(N_A * A)
d_B = -np.sum(N_B * B)
d_C = -np.sum(N_C * C)

# Create a meshgrid:

delta = 200

xlim_A = A[0] - delta, A[0] + delta
ylim_A = A[1] - delta, A[1] + delta
xx_A, yy_A = np.meshgrid(np.arange(*xlim_A), np.arange(*ylim_A))

xlim_B = B[0] - delta, B[0] + delta
ylim_B = B[1] - delta, B[1] + delta
xx_B, yy_B = np.meshgrid(np.arange(*xlim_B), np.arange(*ylim_B))

xlim_C = C[0] - delta, C[0] + delta
ylim_C = C[1] - delta, C[1] + delta
xx_C, yy_C = np.meshgrid(np.arange(*xlim_C), np.arange(*ylim_C))

# Solving the equation above for z:

# z = -(a*x + b*y +d) / c

zz_A = -(N_A[0] * xx_A + N_A[1] * yy_A + d_A) / N_A[2]
zz_B = -(N_B[0] * xx_B + N_B[1] * yy_B + d_B) / N_B[2]
zz_C = -(N_C[0] * xx_C + N_C[1] * yy_C + d_C) / N_C[2]

fig = plt.figure()
ax = fig.add_subplot(projection='3d')
ax.plot_surface(xx_A, yy_A, zz_A, alpha=0.5, color='g')
ax.plot_surface(xx_B, yy_B, zz_B, alpha=0.5, color='cyan')
ax.plot_surface(xx_C, yy_C, zz_C, alpha=0.5, color='crimson')

# Plot point.

x_A, y_A, z_A = A
x_B, y_B, z_B = B
x_C, y_C, z_C = C

Plane_A, = ax.plot(x_A, y_A, z_A, marker='o', markersize=5, color='g')
Plane_A.set_label('A position')
ax.legend()

Plane_B, = ax.plot(x_B, y_B, z_B, marker='o', markersize=5, color='cyan')
Plane_B.set_label('B position')
ax.legend()

Plane_C, = ax.plot(x_C, y_C, z_C, marker='o', markersize=5, color='crimson')
Plane_C.set_label('C position')
ax.legend()

# Plot normal.

dx_A, dy_A, dz_A = delta * N_A
ax.quiver(x_A, y_A, z_A, dx_A, dy_A, dz_A, arrow_length_ratio=0.15, linewidth=3, color='g')

dx_B, dy_B, dz_B = delta * N_B
ax.quiver(x_B, y_B, z_B, dx_B, dy_B, dz_B, arrow_length_ratio=0.15, linewidth=3, color='cyan')

dx_C, dy_C, dz_C = delta * N_C
ax.quiver(x_B, y_C, z_C, dx_C, dy_C, dz_C, arrow_length_ratio=0.15, linewidth=3, color='crimson')

# Enforce equal axis aspects so that the normal also appears to be normal.

ax.set_xlim(xmax=1500,xmin=-1500)
ax.set_ylim(ymax=400, ymin=-400)

zlim = max(A[2], B[2], C[2]) - delta, max(A[2], B[2], C[2]) + delta
ax.set_zlim(*zlim)

ax = plt.gca()

# ax.set_box_aspect([1,1,1])

set_axes_equal(ax)
ax.set_xlabel('X', fontsize=20)
ax.set_ylabel('Y', fontsize=20)
ax.set_zlabel('Z', fontsize=20)

plt.show()