0

Given: a unit sphere S with $M:=(0,0,0)$ and

a normal unit vector $n_0:=(n_x,n_y,n_z)$ to a plane E containing $M$,

further a parameter $t \in [0,2\pi)$.

The intersection of S and E defines a great circle C with center $M$.

Wanted: a parametric equation of the coordinates of points $P(t)$ in the circumference of C in polar $(\Phi (t),\Theta (t))$ or Cartesian coordinates $(x(t),y(t),z(t))$.

Had some trouble derivation, maybe there is a ready to go resolution which will be used for programming in Python language.

Anton Vrdoljak
  • 2,717
peets
  • 101

1 Answers1

2

OK. Let's start by finding two vectors, $u$ and $v$, orthogonal to $n$ and each other. My favorite approach is this: look at the two largest (in absolutely value) coordinates of $n$. Let's say they're the $x$ and $y$ coordinates. Then we build a vector $$ b = \pmatrix{-n_y\\n_x\\0} $$ which is (a) nonzero, and (b) orthogonal to $n$. [If the two largest entries had been $n_x$ and $n_z$, you'd let $b = \pmatrix{-n_z \\ 0 \\ n_x}$; I'll let you figure out the third possibility.]

Then let $u = \frac{1}{\|b\|}b$ , and $v = n \times u$.

Then define $$ x(t) = \cos(t) u_x - \sin(t) v_x\\ y(t) = \cos(t) u_y - \sin(t) v_y \\ z(t) = \cos(t) u_z - \sin(t) v_z $$ and you're done. That's perhaps better expressed as $$ P(t) = \cos(t) u + \sin (t) v . $$

Small postscript: Folks tend to not like step 1, where they have to find the two largest coordinates, arguing that there should be a simple formula for a vector $v(n)$ orthogonal to $n$. But the "can't comb the hair on a billiard ball" theorem shows that there is no such formula, at least if you want $v$ to vary continuously as a function of $n$. If you're willing to go with a discontinuous formula...well, that's what my "pick the two largest entries" actually is. There's a nifty paper somewhere that points out that by using clever tricks relating to the IEEE floating point representation of numbers between 0 and 1, there's a formula you can write in many programming languages without using anything "if"-like (including no "modulo" operations) --- I seem to recall it includes masking to select certain bits in the the FP representation, and then doing bit-wise magic --- but practically speaking, you probably want something like what I described, because a year from now it'll be a lot easier to debug!

John Hughes
  • 93,729
  • @Anton: thanks for correcting my original text; I have to learn the special formats. – peets Oct 22 '23 at 10:28
  • Thanks a lot. I use a solution with a circle in x-y plane with $n_xy$=(0,0,1) and then rotating it by an angle to the given normal vector $n_0$ with the axis of $n_0$ coordinates in x-y plane, using Quaternions. This costs a lot of computing time, but has the advantage to work in any direction of $n_0$ without exceptions. Thus I was looking for a faster way, which you described; only one cos and sin per point. The necessary exceptions for $n_0$ on an axis of x, y, z is only little overhead. The python functions can be accelerated by numba more easily. :c) – peets Oct 22 '23 at 10:34
  • v has to be normalized (length 1) - same procedure as b -> u – peets Oct 23 '23 at 06:36
  • Python program speed was about 50 times faster (without numba!) – peets Oct 23 '23 at 08:37
  • 1
    v does not have to be normalized: n and u are orthogonal unit vectors, so their cross product is also a unit vector. (It might not hurt to normalize anyhow, to suppress accumulating floating-point errors, but this isn't a procedure where there's any danger of accumulation, and all the numbers are in the sweet-spot for FP accuracy (i.e, between -1 and 1). – John Hughes Oct 23 '23 at 20:47
  • You are right. Maybe I had a software problem in my python code. I obtained smaller unit circles in some directions when plotting them. Normalizing v fixed the problem, so I did not think any further. Supposedly a wrong used python function is the reason, but I have just now no time to find it. Sorry for the confusion. – peets Oct 26 '23 at 18:32