views:

33

answers:

2

I need a function that returns points on a circle in three dimensions.

The circle should "cap" a line segment defined by points A and B and it's radius. each cap is perpendicular to the line segment. and centered at one of the endpoints.

Here is a shitty diagram

+1  A: 

Let N be the unit vector in the direction from A to B, i.e., N = (B-A) / length(A-B). The first step is to find two more vectors X and Y such that {N, X, Y} form a basis. That means you want two more vectors so that all pairs of {N, X, Y} are perpendicular to each other and also so that they are all unit vectors. Another way to think about this is that you want to create a new coordinate system whose x-axis lines up with the line segment. You need to find vectors pointing in the direction of the y-axis and z-axis.

Note that there are infinitely many choices for X and Y. You just need to somehow find two that work.

One way to do this is to first find vectors {N, W, V} where N is from above and W and V are two of (1,0,0), (0,1,0), and (0,0,1). Pick the two vectors for W and V that correspond to the smallest coordinates of N. So if N = (.31, .95, 0) then you pick (1,0,0) and (0,0,1) for W and V. (Math geek note: This way of picking W and V ensures that {N,W,V} spans R^3). Then you apply the Gram-Schmidt process to {N, W, V} to get vectors {N, X, Y} as above. Note that you need the vector N to be the first vector so that it doesn't get changed by the process.

So now you have two vectors that are perpendicular to the line segment and perpendicular to each other. This means the points on the circle around A are X * cos t + Y * sin t + A where 0 <= t < 2 * pi. This is exactly like the usual description of a circle in two dimensions; it is just written in the new coordinate system described above.

David Norman
Thank you. That was a very helpful answer.
Nathan
> So if N = (.31, .95, 0) then you pick (1,0,0) and (0,0,1) for W and V ------ What if one of the components of N is negative? is that smaller, or should I use the absolute value?
Nathan
The smallest in absolute value.
David Norman
+1  A: 

As David Norman noted the crux is to find two orthogonal unit vectors X,Y that are orthogonal to N. However I think a simpler way to compute these is by finding the householder reflection Q that maps N to a multiple of (1,0,0) and then to take as X the image of (0,1,0) under Q and Y as the image of (0,0,1) under Q. While this might sound complicated it comes down to:

s = (N[0] > 0.0) ? 1.0 : -1.0

t = N[0] + s; f = -1.0/(s*t);

X[0] = f*N[1]*t; X[1] = 1 + f*N[1]*N[1]; X[2] = f*N[1]*N[2];

Y[0] = f*N[2]*t; Y[1] = f*N[1]*N[2]; Y[2] = 1 + f*N[2]*N[2];

dmuir
And what do you do if N[0] == 0 or worse if N[0] is within floating point error of 0? Perhaps you first need to permute the coordinates of N so that the largest component in absolute value is first?
David Norman
N[0] near or equal to 0 is not a problem; |t| will then be close to 1 and the only division will be by something close to 1 is absolute value. X and Y will still be (nearly) orthogonal and (nearly) orthogonal to N
dmuir