0
$\begingroup$

I have 2 points on a sphere, given by $(\theta, \phi)$. How can I find $n$ equally spaced points on the smaller arc of the great circle between them, without converting the coordinates to another system?

I am sure that the smaller arc does not cross the singularity line $\phi=0,2\pi$.

coordinates = {{θ1, ϕ1}, {θ2, ϕ2}}

Thanks in advance.

Improve this question
$\endgroup$
4
  • $\begingroup$ "without converting the coordinates to another system" - so you don't allow conversion to Cartesian coordinates, even if it can make things a bit easier? Anyway, you might want to see this. $\endgroup$
    – J. M.'s lack of A.I. ♦
    Apr 6, 2017 at 16:59
  • $\begingroup$ No, because I have done some after-processing of the points in the case when they lie on the arc phi = 0, 2*pi. If I convert them, and then again, I will lose this processing. $\endgroup$
    – Morgan
    Apr 6, 2017 at 17:03
  • $\begingroup$ This math.stackexchange.com/questions/1476391/… gives you the cross product in spherical. That will be normal to the plane holding the great circle through your two points and seems to be a large step towards what you want. $\endgroup$
    – Bill
    Apr 6, 2017 at 23:08
  • $\begingroup$ I provide a formula in mathematica.stackexchange.com/questions/180422/… but it's hideously ugly and needs to be simplified. $\endgroup$
    – user1722
    Aug 22, 2018 at 15:24

2 Answers 2

Reset to default
2
$\begingroup$

I remain mystified by the OP's insistence against a Cartesian conversion, since it eases the task of subdivision in this case. Nevertheless, here is a simple demonstration of "slerping" (previously shown here):

sphericalDistance[{u1_, v1_}, {u2_, v2_}] := 
         InverseHaversine[Haversine[v1 - v2] + Sin[v1] Sin[v2] Haversine[u1 - u2]]

With[{n = 10, p1 = N[{π/6, π/4}], p2 = N[{-π/2, π/2}]},
     Block[{d, h},
           d = sphericalDistance[p1, p2]; h = Range[0, 1, 1/(n - 1)];
           arcDots = (Sin[Transpose[{1 - h, h}] d]/Sin[d]).
                     (Append[Sin[#2] Through[{Cos, Sin}[#1]], Cos[#2]] & @@@ {p1, p2});
           Graphics3D[{Sphere[], Red, Sphere[arcDots, 0.02]}, Boxed -> False]]]

spherical linear interpolation

If the longitude and colatitude angles are needed, here is some code for obtaining them:

{ArcTan @@ Most[#], ArcCos[Last[#]/Norm[#]]} & /@ arcDots
Improve this answer
$\endgroup$
2
  • $\begingroup$ Good, only in my question I have theta, then phi, so your #1 and #2 are the other way around. $\endgroup$
    – Morgan
    Apr 25, 2017 at 19:17
  • $\begingroup$ @Morgan, I'll leave the conversion to your preferred convention up to you; I'm accustomed to putting the longitude before the colatitude. $\endgroup$
    – J. M.'s lack of A.I. ♦
    Apr 26, 2017 at 2:39
1
$\begingroup$

A first thought:

f[s_, t_] := {Sin[s] Cos[t], Sin[s] Sin[t], Cos[s]}
fun[{a_, b_}, {c_, d_}, n_] := 
 Module[{p1 = f[a, b], p2 = f[c, d], ax, angle, pts},
  ax = Cross[p1, p2];
  angle = VectorAngle[p1, p2];
  pts = Table[RotationMatrix[j angle/n, ax].p1, {j, 0, n, 1}];
  Graphics3D[{Opacity[0.5], Sphere[], Opacity[1], Red, 
    PointSize[0.02], Point[{p1, p2}], PointSize[0.01], Black, 
    Point[pts]}]]

For example:

Manipulate[
 fun[u, v, 
  n], {u, {0, 0}, {2 Pi, 2 Pi}}, {{v, {1, 1}}, {0, 0}, {2 Pi, 
   2 Pi}}, {n, {5, 10, 20}}]

enter image description here

Improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

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