Find 4 points equidistant from a centroid in 3D space where all points lie on a predefined plane

Question

I have a plane defined as an xyz vector and a point that lies on the plane.

I would like to generate xyz coordinates for 4 points (N_points) on the plane surrounding the defined point (centroid) at a defined distance/radius (r).

My current solution only works in 2D. I would like to expand this to work in 3D but my knowledge of geometry is failing me. Any ideas would be much appreciated.

def circlePoints(r, N_points, plane=(1,1,1), centroid=(0,0,0), rotation=0):
    (plane_x, plane_y, plane_z) = plane
    (centroid_x, centroid_y, centroid_z) = centroid

    step = (np.pi*2) / N_points
    rot=rotation
    i=0
    points=[]
    for i in xrange(N_points):
        x = round(centroid_x + ((np.sin(rot)*r) / plane_x), 2)
        y = round(centroid_y + ((np.cos(rot)*r) / plane_y), 2)
        z=0 #?
        points.append((x,y,z))
        rot+=step
    return points

print circlePoints(1, 4, [1,2,0], [2,3,1])
print circlePoints(1, 4)

Answer 1

We need to find two vectors perpendicular to plane (the normal). We can do so by the following procedure:

• Normalize plane
• Set a vector k = (1, 0, 0)
• Calculate math.abs(np.dot(k, plane))
• If > 0.9 then set k = (0, 1, 0)
• Calculate a = np.cross(k, plane)) and b = np.cross(plane, a)
• You now have two vectors in the plane. You can get any points in the plane by adding some number times these two vectors and adding to centeroid
• If you want specific distances, you need to normalize a and b

Code:

import numpy as np
import math

def normalize(a):
    b = 1.0 / math.sqrt(np.sum(a ** 2))
    return a * b

def circlePoints(r, N_points, plane=(1,1,1), centroid=(0,0,0)):
    p = normalize(np.array(plane))
    k = (1, 0, 0)
    if math.fabs(np.dot(k, p)) > 0.9:
        k = (0, 1, 0)
    a = normalize(np.cross(k, p))
    b = normalize(np.cross(p, a))
    step = (np.pi * 2) / N_points
    ang = [step * i for i in xrange(N_points)]
    return [(np.array(centroid) + \
            r * (math.cos(rot) * a + math.sin(rot) * b)) \
            for rot in ang]

print circlePoints(10, 5, (1, 1, 1), (0, 0, 0))