Solid Geometry
Problems Associated
with the Perspective 
3-Point (P3P) Problem

M. Q. Rieck
                                                                                             Drake University
                                                                                                  10/29/2022

michael.rieck@drake.edu

Abstract:  A couple compelling problems in solid geometry appear to have been ignored to date in the mathematical literature. 
These problems are very easy to state, but mostly difficult to solve. However, recent advances in the analysis of the so-called 
“Perspective 3-Point (P3P) Problem” make it comparatively easy to solve these problems. After stating these problems, each is 
split into two cases, an “acute triangle case” and an “obtuse triangle case.” This is done because the analyses of these cases wind up 
being quite different. For each case of each problem, a proposed solution is given. While no proofs are presented here, the claims 
that have been proved are identified as such. Extensive experimental evidence has been collected to support all of the claims
(using Mathematica and C++).  Mathematica code for producing graphics is included.

1.  Statements of the Problems

(a)  The “Three Lines Problem”

*  Fix three directed lines, l1, l2 and l3 , in Cartesian 3-space, all passing through the origin,
           and not coplanar.
           
       *  Let α be the angle between  l2 and l3  (as directed lines), let β be the angle between  l3 and l1, 
       and let γ be the angle between  l1 and l2  . 
            
       *  Fix three positive numbers ∡A, ∡B and ∡C such that ∡A + ∡B + ∡C = π.  (These symbols 
            just represent numbers, and do not suggest the existence of points A, B and C.)
            
       *   Now ask this:  Do l1, l2, l3  contain points A, B, C, respectively, such that the triangle 
             △ABC has interior angles  ∡CAB = ∡A ,  ∡ABC = ∡B  and  ∡BCA = ∡C ?   (Radians) 
             
       *   We seek a logical combination of inequalities for answering this question, where these 
            inequalities are as simple and intuitive as possible.

(b)  The “Three Rays Problem”

*   Starting over, replace the three directed lines with three rays from the origin. 

        *   Let α, β and γ be similar to before, but now angles between rays instead of directed lines. 
        
        *   Let ∡A, ∡B and ∡C be as before.

        *   Ask the same question as before, but using the rays instead of the directed lines.

(c)  Splitting the Problems

*  Because the solution for each of these problems has a different character depending 
            on whether or not any of the angles (numbers) ∡A, ∡B , ∡C exceeds π/2, it will be handy
            to split each problem into two cases, accordingly: 
            
        *  LA refers to the Three Lines Problem when none of ∡A, ∡B , ∡C exceeds π/2. 
        
        *  LO refers to the Three Lines Problem when one of ∡A, ∡B , ∡C exceeds π/2. 
        
        *  RA refers to the Three Rays Problem when none of ∡A, ∡B , ∡C exceeds π/2. 
        
        *  RO refers to the Three Rays Problem when one of ∡A, ∡B , ∡C exceeds π/2.
        
        *  The goal now is to treat each of LA, LO, RA and RO in turn, and for each, to propose
            a system of conditions on the six quantities, ∡A, ∡B , ∡C , α, β and γ, that is capable of 
            answering the question posed by the problem.

2.  Proposed Solutions

(a)  The LA (“lines/acute”) Problem

*  This is by far the easiest of the four problems to solve. 

        *  Besides the restrictions on ∡A, ∡B and ∡C, the following evident relationships hold 
            due to the definitions:  α + β + γ < 2π , α < β + γ , β < α + γ , γ < α + β . 
                   
         *  There are no other restrictions. 

         *  This last claim has been proved; the proof proceeds by first considering the case where one
              of  l1, l2, l3 is orthogonal to the other two. 
                            
         *  The simple relationships between ∡A, ∡B , ∡C , α, β and γ, discussed thus far will be tacitly 
             understood to hold, going forward.  They are requirements for all four problems.

(b)  The LO (“lines/obtuse”) Problem

*  Besides the constraints seen in the LA problem, the LO problem has extra constraints 
             that are not at all obvious. 

         *  Evidence that extra conditions are required is easily obtained: When  l1, l2, l3  are 
              pairwise orthogonal, it is easy to see that the Lines Problem requires that △ABC be 
              an acute or right triangle. Obtuse triangles do not fit. 
             
         *  Before proposing a solution to the LO problem, it is necessary to introduce a quantity 
             D,  which is quite a complicated polynomial of  sin2∡A,  sin2∡B,  sin2∡C,  cos α ,  cos β 
             and cos γ.  It plays a critical role in the analysis, and will be defined later.  

         *  Another needed quantity q will be defined later too.

*  Henceforth, it will be convenient to let c1 = cos α , c2 = cos β , c3 = cos γ , C1 = c12 , 
             C2 = c22 , C3 = c32 ,  and C0 = c1 c2 c3 . 

         *  The claim will now be made that if ∡A > π/2, then the question in the LO problem 
             is answered in the negative if and only if the following compound condition holds: 

                   D > 0   ∧   {  q < 0
                                ∨   [ 2C2 < 1 + cos 2∡B  ∧  (1 + cos 2∡A) C3 < cos 2∡A+cos 2∡C ] 
                                ∨   [ 2C3 < 1 + cos 2∡C  ∧  (1 + cos 2∡A) C2 < cos 2∡A+cos 2∡B ]  }.
             
         *  Careful analysis and extensive testing support this claim, but as yet, it is unproven. 
         
         *  Symmetric reasoning can be applied to deal with the other two obtuse triangle cases.

(c)  The RA (“rays/acute”) Problem

*  Next, it is claimed that the question in the RA problem is answered in the affirmative 
            if and only if the following compound condition holds:
            
                                  [ α < ∡A   →   c2 cos∡C + c3 cos∡B > 0 ] 
                          ∧     [ β < ∡B   →   c3 cos∡A + c1 cos∡C > 0 ] 
                          ∧     [ γ < ∡C    →   c1 cos∡B + c2 cos∡A > 0 ] 
                          ∧     [ (α < ∡A  ∧  β > ∡B  ∧  γ > ∡C )  →  D > 0 ] 
                          ∧     [ (α > ∡A  ∧  β < ∡B  ∧  γ > ∡C)   →  D > 0 ] 
                          ∧     [ (α > ∡A  ∧  β > ∡B  ∧  γ < ∡C)   →  D > 0 ] . 
               
                (Here “→” is logical implication, and of course p → q is equivalent to ⌉ p ⋁ q.)
            
        *  This condition has been proved to be necessary in order that the RA question
            be answered in the affirmative. 
            
       *   While extensive tests suggest that it is also sufficient, this has not yet been proved.

(d)  The RO (“rays/obtuse”) Problem

*  This problem is the most challenging, and while a proposed solution exists, and has
             been extensively tested, it is somewhat complicated, and will be omitted here.

         *  The proposed solution can be found in the “Four Problems” document posted at 
         
                                  https://github.com/mqrieck/tetrahedron_test.cpp

(e)  The Quantities D and q

*  Though the claims made here concerning the LO and the RA problems have not yet been
             completely proved, they have emerged from a careful analysis involving recently published
             aspects of the P3P problem, in particular, the 2021 paper by Rieck and Wang. 
             
        *   The quantity D is discussed in the paper, where it is referred to as the “discriminant of
             Grunert’s system,” and where it is demonstrated that D = 0 if and only if the P3P problem 
             has a repeated solution. 
             
        *   Details will now be given for efficient ways to compute D and q.

*   Fix a triangle △ A’ B’  C’  on the unit circle in the Cartesian plane (xy-plane), similar to
             △AB C.  Further, do this so that there are angles (Subscript[ϕ, 1]) , (Subscript[ϕ, 2]) and (Subscript[ϕ, 3]) such that (Subscript[ϕ, 1])+(Subscript[ϕ, 2])+(Subscript[ϕ, 3]) = 0, 
             and the coordinates of A’ , B’ and C’ are, respectively, ((Subscript[x, 1]), (Subscript[y, 1])) = (cos (Subscript[ϕ, 1]), sin (Subscript[ϕ, 1])),   ((Subscript[x, 2]), (Subscript[y, 2]))
             = (cos (Subscript[ϕ, 2]), sin (Subscript[ϕ, 2])), and ((Subscript[x, 3]), (Subscript[y, 3])) = (cos (Subscript[ϕ, 3]), sin (Subscript[ϕ, 3])).   (This is always possible.) 

        *   Define these quantities: 
        
                  H  =  1 + 2 C0 - C1 - C2 - C3 ,  
                  L   =  2 [ ((Subscript[y, 1])+(Subscript[y, 2])+(Subscript[y, 3]))(1 - C0) + (1-(Subscript[x, 1]))(Subscript[y, 1])(C1 - 1) + (1-(Subscript[x, 2]))(Subscript[y, 2])(C2 - 1) + (1-(Subscript[x, 3]))(Subscript[y, 3])(C3 - 1) ] ,
                  R   =  2 [ (1+(Subscript[x, 1])+(Subscript[x, 2])+(Subscript[x, 3]))(1 - C0) + (1+(Subscript[x, 1]))(Subscript[x, 1])(C1 - 1) + (1+(Subscript[x, 2]))(Subscript[x, 2])(C2 - 1) + (1+(Subscript[x, 3]))(Subscript[x, 3])(C3 - 1) ] ,
                  K  =  - H - R , 
                  D  =  4 H ( 2 K + 3 H )3 - ( K2 + L2 + 12 H K + 9 H2 )2   (sign different elsewhere).

         *   Further define these quantities: 
                  (Subscript[j, 2])  =  (1 - C2)(1 - C3) (Subscript[x, 2]) ,  (Subscript[k, 2])  =  (1 - C2)(1 - C3) (Subscript[y, 2]) ,  (Subscript[m, 2]) =  (1 - C3)(C2 - cos 2∡B) ,
                  (Subscript[j, 3])  =  (1 - C2)(1 - C3) (Subscript[x, 3]) ,  (Subscript[k, 3])  =  (1 - C2)(1 - C3) (Subscript[y, 3]) ,  (Subscript[m, 3]) =  (1 - C2)(C3 - cos 2∡C) ,
                   q   =  ((Subscript[j, 2]) (Subscript[m, 3]) - (Subscript[j, 3]) (Subscript[m, 2]))2 + ((Subscript[k, 2]) (Subscript[m, 3]) - (Subscript[k, 3]) (Subscript[m, 2]))2 - ((Subscript[j, 2]) (Subscript[k, 3]) - (Subscript[j, 3]) (Subscript[k, 2]))2.

3.  Demonstrations

(a)  The Three Lines Problem

*   It is straightforward to check that 

                   H = 4 sin[(α + β + γ)/2] sin[(-α + β + γ)/2] sin[(α - β + γ)/2] sin[(α + β - γ)/2], 

           and that the requirements imposed on α, β, γ mean that H > 0, which will be assumed. 
           
      *   The following code demonstrates the claims made concerning the Three Lines Problem.

      *   The triangle △ABC is set by specifying two side lengths (b and c) and an interior angle (∡A), 
           but only the shape of the triangle matters, not the scale.  A value for c3 is also fixed.     

      *   For these fixed parameters, actual pairs (c1, c2) that can occur (for lines l1, l2, l3) are 
           generated, and displayed as points in the (c1, c2)-plane.  This is the blue region.

*   These points are generated by extending the method of John Sullivan that is discussed in 
           John Wetzel’s paper “An Ancient Elliptic Locus.”
           
      *   Another region is computed based on the solution proposed in these slides.  It shows the
           points (c1, c2) that are alleged to not be allowable.  This is the yellow region. 
            
      *   These two regions always appear to partition the ellipse interior given by H > 0, as needed.
      
      *   The red curve is where D = 0, and generally, part of this curve serves as a boundary between 
           the blue region and the yellow region.
           
      *   Further discussion of all this, and some images produced from the following demos can 
           be found in my “Four Problems” document.

(* Be sure to activate this first to compute discriminant (D) and also compute q *)
disc[x1_,x2_,x3_,y1_,y2_,y3_,c1_, c2_, c3_] :=  Module[ {L,R,H},
4H(((H-2R))^3) -((((L^2)+(R^2)-10H R -2(H^2)))^2)  /. {
L -> 2 ((-1+(c1^2)) (1-x1) y1+(-1+(c2^2)) (1-x2) y2+(-1+(c3^2)) (1-x3) y3+(1-c1 c2 c3) (y1+y2+y3)), 
R -> 2 ((-1+(c1^2)) x1 (1+x1)+(-1+(c2^2)) x2 (1+x2)+(-1+(c3^2)) x3 (1+x3)+(1-c1 c2 c3) (1+x1+x2+x3)), 
H ->  1-(c1^2)-(c2^2)+2 c1 c2 c3-(c3^2) }  ];
q[x2_,y2_,x3_,y3_,B_,C_,c2_,c3_] := Module[
{ j2,k2,m2,j3,k3,m3 },  
(j2 m3 - j3 m2)^2 +  (k2 m3 - k3 m2)^2 - (j2 k3 - j3 k2)^2  /. {
j2 -> (1-c2^2)(1-c3^2) x2 , k2 -> (1-c2^2)(1-c3^2) y2 , m2 ->(1-c3^2) (c2^2-Cos[2 B]), 
j3 -> (1-c2^2)(1-c3^2) x3 , k3 -> (1-c2^2)(1-c3^2) y3 , m3 ->(1-c2^2) (c3^2-Cos[2 C]) } ];

(* Demonstration of the Three Lines Problem *) 
Manipulate[ Module[
{a,B,C,costurn,sinturn,x10,x20,x30,y10,y20,y30,x1,x2,x3,y1,y2,y3,t1,t2,t3,DD,HH,LL,RR,EE, region, plot, ellipse,contour, points, discrCurve, lines, L, R, c3, type, ϕ1,ϕ2,ϕ3,showList}, 
c3 = Cos[ϕ]; 
a = Sqrt[b^2+c^2-2b c Cos[A]];
B = ArcCos[(a^2+c^2-b^2)/(2a c)];
C = ArcCos[(a^2+b^2-c^2)/(2a b)];
If[ A > B+C , type = "obtuseA", 
If[ B > C+A, type = "obtuseB", 
If[ C > A+B, type = "obtuseC", 
type = "acute" ]]];
(*
costurn = Cos[2*(B-C)/3];
sinturn = Sin[2*(B-C)/3];
x10 = 1; y10 = 0;
x20 = Cos[2*C]; y20 =  Sin[2*C];
x30 = Cos[2*B]; y30 = -Sin[2*B];
x1 = x10 * costurn - y10 * sinturn;
y1 = x10 * sinturn + y10 * costurn;
x2 = x20 * costurn - y20 * sinturn;
y2 = x20 * sinturn + y20 * costurn;
x3 = x30 * costurn - y30 * sinturn;
y3 = x30 * sinturn + y30 * costurn;
t1 = y1/(1+x1); t2 = y2/(1+x2); t3 = y3/(1+x3); 
ϕ1 = ArcTan[x1,y1]; ϕ2 = ArcTan[x2,y2];  ϕ3 = ArcTan[x3,y3];
*)
ϕ1 = 2(A+2B)/3; ϕ2 = -2(2A+B)/3; ϕ3 = 2(A-B)/3; 
x1 = Cos[ϕ1]; x2 = Cos[ϕ2]; x3 = Cos[ϕ3]; 
y1 = Sin[ϕ1]; y2 = Sin[ϕ2]; y3 = Sin[ϕ3]; 
t1 = Tan[ϕ1/2]; t2 = Tan[ϕ2/2]; t3 = Tan[ϕ3/2]; 

If[ generate,  
plot = ParametricPlot[ {
(c Cos[ϕ] Sin[θ]+b (Cos[A] Cos[θ]-Cos[ψ] Sin[A] Sin[θ]) Sin[ϕ])/Sqrt[(c^2) ((Sin[θ])^2)+(b^2) ((Sin[ϕ])^2)+2 b c Sin[θ] Sin[ϕ] (Cos[A] Cos[θ+ϕ]-Cos[ψ] Sin[A] Sin[θ+ϕ])],   
(c Sin[θ]+b Sin[ϕ] (Cos[A] Cos[θ+ϕ]-Cos[ψ] Sin[A] Sin[θ+ϕ])) /Sqrt[((c^2) ((Sin[θ])^2)+(b^2) ((Sin[ϕ])^2)+2 b c Sin[θ] Sin[ϕ] (Cos[A] Cos[θ+ϕ]-Cos[ψ] Sin[A] Sin[θ+ϕ]))]
}, {ψ,-N[π]+.0001,N[π]-.0001}, {θ, .0001, N[2π]-.0001},PlotRange->{{-1.1,1.1},{-1.1,1.1}}, PlotStyle->{RGBColor[.15,.35,.6],Opacity[op1]}, PlotPoints->pp1 ] ]; 

If[ forbidden, 
If[ type == "obtuseA", 
region = RegionPlot[ 
   1 - c1^2 - c2^2 - c3^2 + 2 c1 c2 c3 >= 0   &&   
disc[x1,x2,x3,y1,y2,y3,c1,c2,c3] >= 0   &&  ( 
q[x2,y2,x3,y3,B,C,c2,c3] <= 0   ||   
( 2 (c2^2) ≤ 1+Cos[2 B] && (c3^2) (1+Cos[2 A]) ≤ Cos[2 A]+Cos[2 C] ) ||
((c2^2) (1+Cos[2 A]) ≤ Cos[2 A]+Cos[2 B] && 2 (c3^2) ≤ 1+Cos[2 C] )  )  , 
 {c1,-1.1,1.1},{c2,-1.1,1.1}  , PlotStyle->{Yellow, Opacity[op2]}, BoundaryStyle->None, PlotPoints->pp2 ]]; 

If[ type == "obtuseB", 
region = RegionPlot[ 
   1 - c1^2 - c2^2 - c3^2 + 2 c1 c2 c3 >= 0   &&   
disc[x1,x2,x3,y1,y2,y3,c1,c2,c3] >= 0   &&   (
q[x3,y3,x1,y1,C,A,c3,c1] <= 0     ||  
( 2 (c3^2) ≤ 1+Cos[2 C] && (c1^2) (1+Cos[2 B]) ≤ Cos[2 A]+Cos[2 B] ) ||
( (c3^2) (1+Cos[2 B]) ≤ Cos[2 B]+Cos[2 C] && 2 (c1^2)≤1+Cos[2 A] ) )  , 
 {c1,-1.1,1.1},{c2,-1.1,1.1}  , PlotStyle->{Yellow, Opacity[op2]}, BoundaryStyle->None, PlotPoints->pp2 ]]; 

If[ type == "obtuseC", 
region = RegionPlot[ 
   1 - c1^2 - c2^2 - c3^2 + 2 c1 c2 c3 >= 0   &&   
disc[x1,x2,x3,y1,y2,y3,c1,c2,c3] >= 0   &&   (
q[x1,y1,x2,y2,A,B,c1,c2] <= 0     || 
( 2 (c1^2) ≤ 1+Cos[2 A] && (c2^2) (1+Cos[2 C]) ≤ Cos[2 B]+Cos[2 C] )  ||
((c1^2) (1+Cos[2 C]) ≤ Cos[2 A]+Cos[2 C] && 2 (c2^2) ≤ 1+Cos[2 B] ) )  ,  
{c1,-1.1,1.1},{c2,-1.1,1.1}  , PlotStyle->{Yellow, Opacity[op2]}, BoundaryStyle->None, PlotPoints->pp2 ]]; 

If[ type == "acute", 
region = RegionPlot[ False, {c1,-1.1,1.1},{c2,-1.1,1.1}] ]
];

ellipse = ContourPlot[ 1-c1^2-c2^2-c3^2+2c1 c2 c3==0, {c1,-1,1}, {c2,-1,1},PlotRange->{{-1.1,1.1},{-1.1,1.1}},ContourStyle->{Black,Thickness[.005]}, PlotPoints->pp3 ];

discrCurve = ContourPlot[  0 ==( 
-4 (((1-2 R))^3)+(((-2+(L^2)-10 R+(R^2)))^2) //. {
L -> (2 ((-1+(c1^2)) (1-x1) y1+(-1+(c2^2)) (1-x2) y2+(-1+(c3^2)) (1-x3) y3+(1-c1 c2 c3) (y1+y2+y3)))/(1-(c1^2)-(c2^2)+2 c1 c2 c3-(c3^2)),R -> ((2 ((-1+(c1^2)) x1 (1+x1)+(-1+(c2^2)) x2 (1+x2)+(-1+(c3^2)) x3 (1+x3)+(1-c1 c2 c3) (1+x1+x2+x3)))/(1-(c1^2)-(c2^2)+2 c1 c2 c3-(c3^2))),
c3 -> Cos[ϕ] }) ,
 {c1,-1.1,1.1},{c2,-1.1,1.1},PlotRange->{{-1.1,1.1},{-1.1,1.1}}, RegionFunction->Function[{c1,c2}, 1-c1^2-c2^2-c3^2+2c1 c2 c3 >= 0],
ContourStyle->{Red,Thickness[.003]}, PlotPoints->pp3 ];

showList = { ellipse }; 
If[ generate, AppendTo[showList,plot]]; 
If[ forbidden, AppendTo[showList,region]]; 
If[ showCurve, showList = Append[ showList,discrCurve ]]; 
Show[ showList, ImageSize->600,Axes->False,Frame->False,PlotLabel-> Style[ 
(* "phis: " <> ToString[SetPrecision[ϕ1/π,3]] <> "π " <> ToString[SetPrecision[ϕ2/π,3]] <> "π " <> ToString[SetPrecision[ϕ3/π,3]] <> "π\n" <> *)
"side lengths: " <> ToString[SetPrecision[a,3]] <> "  " <> ToString[SetPrecision[b,3]]<> "  " <> ToString[SetPrecision[c,3]] <> "\ninterior angles: "<>ToString[SetPrecision[A/π,3]] <> "π  " <> ToString[SetPrecision[B/π,3]]<> "π  " <> ToString[SetPrecision[C/π,3]]<>"π" <> "\n" <> type <> " ,   ϕ = " <> ToString[SetPrecision[ϕ/π,3]]<> "π,  cos ϕ = " <> ToString[Chop[SetPrecision[Cos[ϕ],3]]] , 20 ]] ] , 
{{b,1.4}, .1, 5},{{c,1}, .1, 5},{{A,.3π}, 0, N[π]},{{ϕ,.575π}, 0.0001, N[π]-.0001}, 
(*{{b,2},.1,5},{{c,3},.1,5},{{A,.2π},0,N[π]},{{ϕ,.55π},0.0001,N[π]-.0001},*)
(*{{b,2},.1,5},{{c,3},.1,5},{{A,.2π},0,N[π]},{{ϕ,.425π},0.0001,N[π]-.0001},*)
{{op1,.4},0,1},{{op2,.4},0,1},{pp1,50},{pp2,50},{pp3,50},
{generate,{True,False}},{forbidden,{True,False}},{showCurve,{False,True}},ControlPlacement->Right ]

(b)  The Three Rays Problem

*   The story and the code here are similar to the Three Lines Problem.

*   However, this code is only designed to work when △ABC is an acute triangle, or when
     ∡C > π/2 , but not when ∡A > π/2 nor when ∡B > π/2 .

RealQ[x_]:=Element[x,Reals]===True ; 
(* Demonstration of the Three Rays Problem *) 
(* Only good for testing acute and obtuseC cases, not obtuseA nor obtuseB cases *)
Manipulate[ Module[
{a,B,C,costurn,sinturn,x10,x20,x30,y10,y20,y30,x1,x2,x3,y1,y2,y3,t1,t2,t3,DD,HH,LL,RR,EE, region, plot, ellipse,contour, points, discrCurve, lines, L, R, c3, type, ϕ1,ϕ2,ϕ3,showList, pt1, pt2}, 

c3 = Cos[ϕ]; 
a = Sqrt[b^2+c^2-2b c Cos[A]];
B = ArcCos[(a^2+c^2-b^2)/(2a c)];
C = ArcCos[(a^2+b^2-c^2)/(2a b)];
If[ A > B+C , type = "obtuseA", 
If[ B > C+A, type = "obtuseB", 
If[ C > A+B, type = "obtuseC", 
type = "acute" ]]];
(*
costurn = Cos[2*(B-C)/3];
sinturn = Sin[2*(B-C)/3];
x10 = 1; y10 = 0;
x20 = Cos[2*C]; y20 =  Sin[2*C];
x30 = Cos[2*B]; y30 = -Sin[2*B];
x1 = x10 * costurn - y10 * sinturn; y1 = x10 * sinturn + y10 * costurn;
x2 = x20 * costurn - y20 * sinturn; y2 = x20 * sinturn + y20 * costurn;
x3 = x30 * costurn - y30 * sinturn; y3 = x30 * sinturn + y30 * costurn;
ϕ1 = ArcTan[x1,y1]; ϕ2 = ArcTan[x2,y2];  ϕ3 = ArcTan[x3,y3] ;
*)
ϕ1 = 2(A+2B)/3; ϕ2 = -2(2A+B)/3; ϕ3 = 2(A-B)/3; 
x1 = Cos[ϕ1]; x2 = Cos[ϕ2]; x3 = Cos[ϕ3]; 
y1 = Sin[ϕ1]; y2 = Sin[ϕ2]; y3 = Sin[ϕ3]; 
t1 = Tan[ϕ1/2]; t2 = Tan[ϕ2/2]; t3 = Tan[ϕ3/2]; 

If[ generate, plot = ParametricPlot[ {
(c Cos[ϕ] Sin[θ]+b (Cos[A] Cos[θ]-Cos[ψ] Sin[A] Sin[θ]) Sin[ϕ])/Sqrt[(c^2) ((Sin[θ])^2)+(b^2) ((Sin[ϕ])^2)+2 b c Sin[θ] Sin[ϕ] (Cos[A] Cos[θ+ϕ]-Cos[ψ] Sin[A] Sin[θ+ϕ])],   
(c Sin[θ]+b Sin[ϕ] (Cos[A] Cos[θ+ϕ]-Cos[ψ] Sin[A] Sin[θ+ϕ])) /Sqrt[((c^2) ((Sin[θ])^2)+(b^2) ((Sin[ϕ])^2)+2 b c Sin[θ] Sin[ϕ] (Cos[A] Cos[θ+ϕ]-Cos[ψ] Sin[A] Sin[θ+ϕ]))]
}, {ψ,-N[π]+.0001,N[π]-.0001}, {θ,.0001,N[π]-ϕ-.0001},PlotRange->{{-1.1,1.1},{-1.1,1.1}}, PlotStyle->{RGBColor[.15,.35,.6],Opacity[op1]}, PlotPoints->pp1 ]];

If[ forbidden, 
If[ type == "obtuseA",   (* wrong *) 
region = RegionPlot[    1 - c1^2 - c2^2 - c3^2 + 2 c1 c2 c3 >= 0  , 
 {c1,-1.1,1.1},{c2,-1.1,1.1}  , PlotStyle->{Yellow, Opacity[op2]},PlotPoints->pp2, BoundaryStyle->None ]]; 

If[ type == "obtuseB",   (* wrong *) 
region = RegionPlot[ 
   1 - c1^2 - c2^2 - c3^2 + 2 c1 c2 c3 >= 0   , 
 {c1,-1.1,1.1},{c2,-1.1,1.1}  , PlotStyle->{Yellow, Opacity[op2]},PlotPoints->pp2, BoundaryStyle->None ]]; 

If[ type == "obtuseC",  
If[ Cos[C] < c3 < -Cos[C],
pt1 = {
Sign[c3]Sqrt[ ((c3^2) Cos [2 ϕ1+ϕ2]-Cos[ϕ1-ϕ2] Cos[2 ϕ1+ϕ2]+Cos[ϕ1+2 ϕ2]-(c3^2) Cos[ϕ1+2 ϕ2]-(Sqrt[2]) (Sqrt[1-2 (c3^2)+Cos[ϕ1-ϕ2]]) Abs[Sin[((ϕ1-ϕ2)/2)]] Sin[2 ϕ1+ϕ2]) /(1-(c3^2)+(c3^2) Cos[2 ϕ1+ϕ2]-Cos[ϕ1-ϕ2] Cos[2 ϕ1+ϕ2]-(Sqrt[2]) (Sqrt[1-2 (c3^2)+Cos[ϕ1-ϕ2]]) Abs[Sin[((ϕ1-ϕ2)/2)]] Sin[2 ϕ1+ϕ2])]
,
Sqrt[(Cos[2 ϕ1+ϕ2]-(c3^2) Cos[2 ϕ1+ϕ2]+(c3^2) Cos[ϕ1+2 ϕ2]-Cos[ϕ1-ϕ2] Cos[ϕ1+2 ϕ2]-(Sqrt[2]) (Sqrt[1-2 (c3^2)+Cos[ϕ1-ϕ2]]) Abs[Sin[((ϕ1-ϕ2)/2)]] Sin[ϕ1+2 ϕ2]) /(1-(c3^2)+(c3^2) Cos[ϕ1+2 ϕ2]-Cos[ϕ1-ϕ2] Cos[ϕ1+2 ϕ2]-(Sqrt[2]) (Sqrt[1-2 (c3^2)+Cos[ϕ1-ϕ2]]) Abs[Sin[((ϕ1-ϕ2)/2)]] Sin[ϕ1+2 ϕ2])]
}; 
pt2 = {
Sqrt[ ((c3^2) Cos[2 ϕ1+ϕ2]-Cos[ϕ1-ϕ2] Cos[2 ϕ1+ϕ2]+Cos[ϕ1+2 ϕ2]-(c3^2) Cos[ϕ1+2 ϕ2]+(Sqrt[2]) (Sqrt[1-2 (c3^2)+Cos[ϕ1-ϕ2]]) Abs[Sin[((ϕ1-ϕ2)/2)]] Sin[2 ϕ1+ϕ2]) /(1-(c3^2)+(c3^2) Cos[2 ϕ1+ϕ2]-Cos[ϕ1-ϕ2] Cos[2 ϕ1+ϕ2]+(Sqrt[2])(Sqrt[1-2 (c3^2)+Cos[ϕ1-ϕ2]]) Abs[Sin[((ϕ1-ϕ2)/2)]] Sin[2 ϕ1+ϕ2])]
,
Sign[c3]Sqrt[(Cos[2 ϕ1+ϕ2]-(c3^2) Cos[2 ϕ1+ϕ2]+(c3^2) Cos[ϕ1+2 ϕ2]-Cos[ϕ1-ϕ2] Cos[ϕ1+2 ϕ2]+(Sqrt[2]) (Sqrt[1-2 (c3^2)+Cos[ϕ1-ϕ2]]) Abs[Sin[((ϕ1-ϕ2)/2)]] Sin[ϕ1+2 ϕ2]) /(1-(c3^2)+(c3^2) Cos[ϕ1+2 ϕ2]-Cos[ϕ1-ϕ2] Cos[ϕ1+2 ϕ2]+(Sqrt[2]) (Sqrt[1-2 (c3^2)+Cos[ϕ1-ϕ2]]) Abs[Sin[((ϕ1-ϕ2)/2)]] Sin[ϕ1+2 ϕ2])]
}; 
If[ RealQ[pt1[[1]]] && RealQ[pt1[[2]]] && RealQ[pt2[[1]]] && RealQ[pt1[[2]]] , 
   region = RegionPlot[ 
   1 - c1^2 - c2^2 - c3^2 + 2 c1 c2 c3 >= 0   &&  (
    (c1 < pt1[[1]] && c2 < pt1[[2]]) ||
 (c1 < pt2[[1]] && c2 < pt2[[2]]) ||
 (c1 < pt1[[1]] && c2 < Cos[B] && disc[x1,x2,x3,y1,y2,y3,c1,c2,c3] < 0) || 
 (c2 < pt2[[2]] && c1 < Cos[A] && disc[x1,x2,x3,y1,y2,y3,c1,c2,c3] < 0) || 
  Cos[B]c1 + Cos[A] c2 <= 0   || 
   (   (* the forbidden region for the lines problem *)
disc[x1,x2,x3,y1,y2,y3,c1,c2,c3] >= 0   &&   (
q[x1,y1,x2,y2,A,B,c1,c2] <= 0   || 
( 2 (c1^2) ≤ 1+Cos[2 A] && (c2^2) (1+Cos[2 C]) ≤ Cos[2 B]+Cos[2 C] )  ||
((c1^2) (1+Cos[2 C]) ≤ Cos[2 A]+Cos[2 C] && 2 (c2^2) ≤ 1+Cos[2 B] )  
)
 ) )  ,  {c1,-1.1,1.1},{c2,-1.1,1.1}  , PlotStyle->{Yellow, Opacity[op2]}, PlotPoints->pp2,BoundaryStyle->None ]  ,
region = RegionPlot[ 
   1 - c1^2 - c2^2 - c3^2 + 2 c1 c2 c3 >= 0   , 
 {c1,-1.1,1.1},{c2,-1.1,1.1}  , PlotStyle->{Yellow, Opacity[op]}, PlotPoints->pp1,BoundaryStyle->None ]] ]; 
If[c3 >= -Cos[C], 
region = RegionPlot[ 
   1 - c1^2 - c2^2 - c3^2 + 2 c1 c2 c3 >= 0   && (
   (c1 < 0 && c2 < Cos[B])|| 
 (c2 < 0 && c1 < Cos[A])||
(c1 < Cos[A] && c2 < Cos[B] &&  disc[x1,x2,x3,y1,y2,y3,c1,c2,c3] < 0) 
),  {c1,-1.1,1.1},{c2,-1.1,1.1}  , PlotStyle->{Yellow, Opacity[op2]}, PlotPoints->pp2,BoundaryStyle->None ] ];
If[c3 <= Cos[C], 
region = RegionPlot[ 
   1 - c1^2 - c2^2 - c3^2 + 2 c1 c2 c3 >= 0   && (
   c1 > Cos[A]||  c2 > Cos[B]
),  {c1,-1.1,1.1},{c2,-1.1,1.1}  , PlotStyle->{Yellow, Opacity[op2]}, PlotPoints->pp2,BoundaryStyle->None ] ]
]; 

If[ type == "acute",
region = RegionPlot[ 
(   1 - c1^2 - c2^2 - c3^2 + 2 c1 c2 c3 >= 0   &&   (
(* sphere-based rules" *)
      A + β + γ > 2π ||
    α + B + γ > 2π  ||
    α + β + C > 2π  || 
(*
     β + γ - α  > 2(B+C)  || 
    γ + α - β  > 2(C+A)  ||
    α + β - γ  > 2(A+B)  ||
*)
(* toroid-based rules *)
(*
    (α < A    &&    β > Max[B,C+α] )   || 
   (α < A    &&    γ > Max[C,B+α] )   ||
   (β < B    &&    γ > Max[C,A+β] )   ||
   (β < B    &&    α > Max[A,C+β] )   ||
   (γ < C    &&    α > Max[A,B+γ] )   ||   
   (γ < C    &&    β > Max[B,A+γ] )   ||
*)
   ( α < A  &&  Cos[C]c2+Cos[B]c3 < 0 ) ||
   ( β < B  &&  Cos[A]c3+Cos[C]c1 < 0 )  ||
   ( γ < C  &&  Cos[B]c1+Cos[A]c2 < 0 )  || 
(* deltoid-based rules *) 
Not[ disc[x1,x2,x3,y1,y2,y3,c1,c2,c3] > 0  || 
     ( ( α >= A || β < B || γ < C ) && ( α < A || β >= B || γ < C ) && ( α < A || β < B || γ >= C ) ) ]
) /. {α->ArcCos[c1],β->ArcCos[c2],γ->ArcCos[c3]}  )  , 
{c1,-1.1,1.1},{c2,-1.1,1.1} , PlotStyle->{Yellow, Opacity[op2]},PlotPoints->pp2, BoundaryStyle->None ] ]
];

ellipse = ContourPlot[ 1-c1^2-c2^2-c3^2+2c1 c2 c3==0, {c1,-1,1}, {c2,-1,1},PlotRange->{{-1.1,1.1},{-1.1,1.1}},ContourStyle->{Black,Thickness[.005]},PlotPoints->pp3 ];

discrCurve = ContourPlot[  0 ==( 
-4 (((1-2 R))^3)+(((-2+(L^2)-10 R+(R^2)))^2) //. {
L -> (2 ((-1+(c1^2)) (1-x1) y1+(-1+(c2^2)) (1-x2) y2+(-1+(c3^2)) (1-x3) y3+(1-c1 c2 c3) (y1+y2+y3)))/(1-(c1^2)-(c2^2)+2 c1 c2 c3-(c3^2)),R -> ((2 ((-1+(c1^2)) x1 (1+x1)+(-1+(c2^2)) x2 (1+x2)+(-1+(c3^2)) x3 (1+x3)+(1-c1 c2 c3) (1+x1+x2+x3)))/(1-(c1^2)-(c2^2)+2 c1 c2 c3-(c3^2))),
c3 -> Cos[ϕ] }) ,
 {c1,-1.1,1.1},{c2,-1.1,1.1},PlotRange->{{-1.1,1.1},{-1.1,1.1}}, RegionFunction->Function[{c1,c2}, 1-c1^2-c2^2-c3^2+2c1 c2 c3 >= 0],
ContourStyle->{Red,Thickness[.003]}, PlotPoints->pp3 ];

showList = { ellipse }; 
If[ generate, AppendTo[showList,plot]]; 
If[ forbidden, AppendTo[showList,region]]; 
If[ showCurve, showList = Append[ showList,discrCurve ]]; 
Show[ showList, ImageSize->600,Axes->False,Frame->False,PlotLabel-> Style[ 
"side lengths: " <> ToString[SetPrecision[a,3]] <> "  " <> ToString[SetPrecision[b,3]]<> "  " <> ToString[SetPrecision[c,3]] <> "\ninterior angles: "<>ToString[SetPrecision[A/π,3]] <> "π  " <> ToString[SetPrecision[B/π,3]]<> "π  " <> ToString[SetPrecision[C/π,3]]<>"π" <> "\n" <> type <> " ,   ϕ = " <> ToString[SetPrecision[ϕ/π,3]]<> "π,  cos ϕ = " <> ToString[Chop[SetPrecision[Cos[ϕ],3]]] , 20 ]] ] , 
    {{b,1.4}, .1, 5},{{c,1}, .1, 5},{{A,.3π}, 0, N[π]},{{ϕ,.575π}, 0.0001, N[π]-.0001}, 
(* {{b,1.9},.1,5},{{c,3},.1,5},{{A,.2π},0,N[π]},{{ϕ,.575π},0.0001,N[π]-.0001}, *)
{{op1,.4},0,1},{{op2,.4},0,1},{pp1,50},{pp2,50},{pp3,50},
{generate,{True,False}},{forbidden,{True,False}},{showCurve,{False,True}},ControlPlacement->Right ]

(c)  Singular Surface in (c1, c2, c3) - Space

*   The equation D = 0 defines an interesting surface in (c1, c2, c3) - space.

       *   We continue to restrict attention to the region where H > 0.
       
       *   Here one fixes the angles B and C (meaning  ∡B and ∡C) for △ABC.
       
       *   If △ABC is acute, the D = 0 surface is homeomorphic to a sphere, and separates an 
            inside region, where D < 0, from an outside region, where D > 0.
             
      *   But if △ABC is obtuse, the D = 0 surface is self-intersecting, and separates the H > 0 
            region into three subregions, two where D > 0, and one where D < 0.
            
       *   This is essentially why the solution to the LO problem is non-trivial.
       
       *   The following code provides a nice rendering of this surface.

       *   A good understanding of the P3P problem is needed to understand how this code 
            generates surfaces and curves, which will not be explained here.

Manipulate[ Module[
{AA,BB,CC,costurn,sinturn,x10,x20,x30,y10,y20,y30,x1,x2,x3,y1,y2,y3,t1,t2,t3,ϕ3,t,Z,points,lines, plot0,plot1,plot2,plot3,infTet,contours, sqrt,curve0,curve1,curve2,curve3,edge01,edge02,edge03,edge11,edge12,edge13,edge21,edge22, edge23,edge31,edge32,edge33,showlist,showlabel,gold,type,beltcolor=White,edgecolor, tetcolor = Black},

edgecolor = gold = RGBColor[0.9961,0.7734,0.1523];
ϕ3 = - ϕ1 - ϕ2;
BB = (ϕ1-ϕ3)/2;
CC = (ϕ2-ϕ1)/2;
AA = π - BB - CC;
If [ ϕ1 >= ϕ2 || ϕ3 <= -Pi || AA ≤ 0 || BB < 0 || CC < 0 , type = "ERROR!!!", 
If [ AA > π/2 ||  BB > π/2 ||  CC > π/2, type="obtuse", type = "acute"]];
x1 = Cos[ϕ1]; y1 = Sin[ϕ1];
x2 = Cos[ϕ2]; y2 = Sin[ϕ2];
x3 = Cos[ϕ3]; y3 = Sin[ϕ3];
t1 = y1/(1+x1); t2 = y2/(1+x2); t3 = y3/(1+x3);

plot0 = ParametricPlot3D[  
{   (1+Z-Cos[θ-ϕ2]-Cos[θ-ϕ3]+Cos[ϕ2-ϕ3])/
((Sqrt[2+Z-2 Cos[θ-ϕ2]]) (Sqrt[2+Z-2 Cos[θ-ϕ3]])),
(1+Z-Cos[θ-ϕ3]-Cos[θ+ϕ2+ϕ3]+Cos[ϕ2+2 ϕ3])/((Sqrt[2+Z-2 Cos[θ-ϕ3]]) (Sqrt[2+Z-2 Cos[θ+ϕ2+ϕ3]])), 
(1+Z-Cos[θ-ϕ2]-Cos[θ+ϕ2+ϕ3]+Cos[2 ϕ2+ϕ3])/((Sqrt[2+Z-2 Cos[θ-ϕ2]]) (Sqrt[2+Z-2 Cos[θ+ϕ2+ϕ3]]))   } , 
{θ,θL,θH}, {Z,ZL,ZH}, PlotRange->{{c1L,c1H},{c2L,c2H},{c3L,c3H}}, PlotPoints->pps, PlotStyle->{Yellow,Opacity[op]},Mesh->mesh ];

plot1 = ParametricPlot3D[  
{   (1+Z-Cos[θ-ϕ2]-Cos[θ-ϕ3]+Cos[ϕ2-ϕ3])/
((Sqrt[2+Z-2 Cos[θ-ϕ2]]) (Sqrt[2+Z-2 Cos[θ-ϕ3]])),
- ((1+Z-Cos[θ-ϕ3]-Cos[θ+ϕ2+ϕ3]+Cos[ϕ2+2 ϕ3])/((Sqrt[2+Z-2 Cos[θ-ϕ3]]) (Sqrt[2+Z-2 Cos[θ+ϕ2+ϕ3]]))), 
- ((1+Z-Cos[θ-ϕ2]-Cos[θ+ϕ2+ϕ3]+Cos[2 ϕ2+ϕ3])/((Sqrt[2+Z-2 Cos[θ-ϕ2]]) (Sqrt[2+Z-2 Cos[θ+ϕ2+ϕ3]])))   } , 
{θ,θL,θH}, {Z,ZL,ZH}, PlotRange->{{c1L,c1H},{c2L,c2H},{c3L,c3H}}, PlotPoints->pps, PlotStyle->{Red,Opacity[op]},Mesh->mesh ];
 
plot2 = ParametricPlot3D[  
{   - ((1+Z-Cos[θ-ϕ2]-Cos[θ-ϕ3]+Cos[ϕ2-ϕ3])/((Sqrt[2+Z-2 Cos[θ-ϕ2]]) (Sqrt[2+Z-2 Cos[θ-ϕ3]]))),
(1+Z-Cos[θ-ϕ3]-Cos[θ+ϕ2+ϕ3]+Cos[ϕ2+2 ϕ3])/((Sqrt[2+Z-2 Cos[θ-ϕ3]]) (Sqrt[2+Z-2 Cos[θ+ϕ2+ϕ3]])),
- ((1+Z-Cos[θ-ϕ2]-Cos[θ+ϕ2+ϕ3]+Cos[2 ϕ2+ϕ3])/((Sqrt[2+Z-2 Cos[θ-ϕ2]]) (Sqrt[2+Z-2 Cos[θ+ϕ2+ϕ3]])))   } , 
{θ,θL,θH}, {Z,ZL,ZH}, PlotRange->{{c1L,c1H},{c2L,c2H},{c3L,c3H}}, PlotPoints->pps, PlotStyle->{Darker[Green],Opacity[op]},Mesh->mesh ];

plot3 = ParametricPlot3D[  
{   - ((1+Z-Cos[θ-ϕ2]-Cos[θ-ϕ3]+Cos[ϕ2-ϕ3])/((Sqrt[2+Z-2 Cos[θ-ϕ2]]) (Sqrt[2+Z-2 Cos[θ-ϕ3]]))),
- ((1+Z-Cos[θ-ϕ3]-Cos[θ+ϕ2+ϕ3]+Cos[ϕ2+2 ϕ3])/((Sqrt[2+Z-2 Cos[θ-ϕ3]]) (Sqrt[2+Z-2 Cos[θ+ϕ2+ϕ3]]))), 
(1+Z-Cos[θ-ϕ2]-Cos[θ+ϕ2+ϕ3]+Cos[2 ϕ2+ϕ3])/((Sqrt[2+Z-2 Cos[θ-ϕ2]]) (Sqrt[2+Z-2 Cos[θ+ϕ2+ϕ3]]))   }, 
{θ,θL,θH}, {Z,ZL,ZH}, PlotRange->{{c1L,c1H},{c2L,c2H},{c3L,c3H}}, PlotPoints->pps, PlotStyle->{Blue,Opacity[op]},Mesh->mesh ];

curve0 = ParametricPlot3D[  
{ Sqrt[ (Cos[θ-ϕ1]+ Cos[ϕ2-ϕ3]) / (Cos[θ-ϕ1]+ 1)]  ,
 Sqrt[ (Cos[θ-ϕ2]+ Cos[ϕ3-ϕ1]) / (Cos[θ-ϕ2]+ 1)]  ,
 Sqrt[ (Cos[θ-ϕ3]+ Cos[ϕ1-ϕ2]) / (Cos[θ-ϕ3]+ 1)]   } 
 , {θ,θL,θH}, PlotRange->{{c1L,c1H},{c2L,c2H},{c3L,c3H}}, PlotPoints->ppc, PlotStyle->{beltcolor,Thickness[bt]} ];

curve1 = ParametricPlot3D[  
{ Sqrt[ (Cos[θ-ϕ1]+ Cos[ϕ2-ϕ3]) / (Cos[θ-ϕ1]+ 1)]  ,
 -Sqrt[ (Cos[θ-ϕ2]+ Cos[ϕ3-ϕ1]) / (Cos[θ-ϕ2]+ 1)]  ,
 -Sqrt[ (Cos[θ-ϕ3]+ Cos[ϕ1-ϕ2]) / (Cos[θ-ϕ3]+ 1)]   } 
 , {θ,θL,θH}, PlotRange->{{c1L,c1H},{c2L,c2H},{c3L,c3H}}, PlotPoints->ppc, PlotStyle->{beltcolor,Thickness[bt]} ];

curve2 = ParametricPlot3D[  
{ -Sqrt[ (Cos[θ-ϕ1]+ Cos[ϕ2-ϕ3]) / (Cos[θ-ϕ1]+ 1)]  ,
 Sqrt[ (Cos[θ-ϕ2]+ Cos[ϕ3-ϕ1]) / (Cos[θ-ϕ2]+ 1)]  ,
 -Sqrt[ (Cos[θ-ϕ3]+ Cos[ϕ1-ϕ2]) / (Cos[θ-ϕ3]+ 1)]   } 
 , {θ,θL,θH}, PlotRange->{{c1L,c1H},{c2L,c2H},{c3L,c3H}}, PlotPoints->ppc, PlotStyle->{beltcolor,Thickness[bt]} ];

curve3 = ParametricPlot3D[  
{ -Sqrt[ (Cos[θ-ϕ1]+ Cos[ϕ2-ϕ3]) / (Cos[θ-ϕ1]+ 1)]  ,
 -Sqrt[ (Cos[θ-ϕ2]+ Cos[ϕ3-ϕ1]) / (Cos[θ-ϕ2]+ 1)]  ,
 Sqrt[ (Cos[θ-ϕ3]+ Cos[ϕ1-ϕ2]) / (Cos[θ-ϕ3]+ 1)]   } 
 , {θ,θL,θH}, PlotRange->{{c1L,c1H},{c2L,c2H},{c3L,c3H}}, PlotPoints->ppc, PlotStyle->{beltcolor,Thickness[bt]} ];

 edge01 = ParametricPlot3D[  {
((1+Z+Cos[ϕ2]+Cos[ϕ2-ϕ3]+Cos[ϕ3])/((Sqrt[2+Z+2 Cos[ϕ2]]) (Sqrt[2+Z+2 Cos[ϕ3]]))),
(1+Z+Cos[ϕ3]+Cos[ϕ2+ϕ3]+Cos[ϕ2+2 ϕ3])/((Sqrt[2+Z+2 Cos[ϕ3]]) (Sqrt[2+Z+2 Cos[ϕ2+ϕ3]])),
(1+Z+Cos[ϕ2]+Cos[ϕ2+ϕ3]+Cos[2 ϕ2+ϕ3])/((Sqrt[2+Z+2 Cos[ϕ2]]) (Sqrt[2+Z+2 Cos[ϕ2+ϕ3]]))
},{Z,ZL, ZH}, PlotRange->{{c1L,c1H},{c2L,c2H},{c3L,c3H}}, PlotPoints->ppc, PlotStyle->  {edgecolor, Thickness[et]} ];

edge02 = ParametricPlot3D[  {
(1+Z+Cos[ϕ2-ϕ3]-Sin[(π/6)+ϕ2]-Sin[(π/6)+ϕ3])/((Sqrt[2+Z-2 Sin[(π/6)+ϕ2]]) (Sqrt[2+Z-2 Sin[(π/6)+ϕ3]])),
(1+Z+Cos[ϕ2+2 ϕ3]-Sin[(π/6)-ϕ2-ϕ3]-Sin[(π/6)+ϕ3])/((Sqrt[2+Z-2 Sin[(π/6)-ϕ2-ϕ3]]) (Sqrt[2+Z-2 Sin[(π/6)+ϕ3]])),
(1+Z+Cos[2 ϕ2+ϕ3]-Sin[(π/6)+ϕ2]-Sin[(π/6)-ϕ2-ϕ3])/((Sqrt[2+Z-2 Sin[(π/6)+ϕ2]]) (Sqrt[2+Z-2 Sin[(π/6)-ϕ2-ϕ3]]))
  }, {Z,ZL, ZH}, PlotRange->{{c1L,c1H},{c2L,c2H},{c3L,c3H}}, PlotPoints->ppc, PlotStyle->  {edgecolor, Thickness[et]} ];

edge03 = ParametricPlot3D[  {
(1+Z+Cos[ϕ2-ϕ3]-Sin[(π/6)-ϕ2]-Sin[(π/6)-ϕ3])/((Sqrt[2+Z-2 Sin[(π/6)-ϕ2]]) (Sqrt[2+Z-2 Sin[(π/6)-ϕ3]])),
(1+Z+Cos[ϕ2+2 ϕ3]-Sin[(π/6)-ϕ3]-Sin[(π/6)+ϕ2+ϕ3])/((Sqrt[2+Z-2 Sin[(π/6)-ϕ3]]) (Sqrt[2+Z-2 Sin[(π/6)+ϕ2+ϕ3]])),
(1+Z+Cos[2 ϕ2+ϕ3]-Sin[(π/6)-ϕ2]-Sin[(π/6)+ϕ2+ϕ3])/((Sqrt[2+Z-2 Sin[(π/6)-ϕ2]]) (Sqrt[2+Z-2 Sin[(π/6)+ϕ2+ϕ3]]))
}, {Z,ZL, ZH}, PlotRange->{{c1L,c1H},{c2L,c2H},{c3L,c3H}}, PlotPoints->ppc, PlotStyle->  {edgecolor, Thickness[et]} ];

edge11 = ParametricPlot3D[  {
((1+Z+Cos[ϕ2]+Cos[ϕ2-ϕ3]+Cos[ϕ3])/((Sqrt[2+Z+2 Cos[ϕ2]]) (Sqrt[2+Z+2 Cos[ϕ3]]))),
- ((1+Z+Cos[ϕ3]+Cos[ϕ2+ϕ3]+Cos[ϕ2+2 ϕ3])/((Sqrt[2+Z+2 Cos[ϕ3]]) (Sqrt[2+Z+2 Cos[ϕ2+ϕ3]]))),
- ((1+Z+Cos[ϕ2]+Cos[ϕ2+ϕ3]+Cos[2 ϕ2+ϕ3])/((Sqrt[2+Z+2 Cos[ϕ2]]) (Sqrt[2+Z+2 Cos[ϕ2+ϕ3]])))
},{Z,ZL, ZH}, PlotRange->{{c1L,c1H},{c2L,c2H},{c3L,c3H}}, PlotPoints->ppc, PlotStyle->  {edgecolor, Thickness[et]} ];

edge12 = ParametricPlot3D[  {
(1+Z+Cos[ϕ2-ϕ3]-Sin[(π/6)+ϕ2]-Sin[(π/6)+ϕ3])/((Sqrt[2+Z-2 Sin[(π/6)+ϕ2]]) (Sqrt[2+Z-2 Sin[(π/6)+ϕ3]])),
- ((1+Z+Cos[ϕ2+2 ϕ3]-Sin[(π/6)-ϕ2-ϕ3]-Sin[(π/6)+ϕ3])/((Sqrt[2+Z-2 Sin[(π/6)-ϕ2-ϕ3]]) (Sqrt[2+Z-2 Sin[(π/6)+ϕ3]]))),
- ((1+Z+Cos[2 ϕ2+ϕ3]-Sin[(π/6)+ϕ2]-Sin[(π/6)-ϕ2-ϕ3])/((Sqrt[2+Z-2 Sin[(π/6)+ϕ2]]) (Sqrt[2+Z-2 Sin[(π/6)-ϕ2-ϕ3]])))
}, {Z,ZL, ZH}, PlotRange->{{c1L,c1H},{c2L,c2H},{c3L,c3H}}, PlotPoints->ppc, PlotStyle->  {edgecolor, Thickness[et]} ];

edge13 = ParametricPlot3D[  {
(1+Z+Cos[ϕ2-ϕ3]-Sin[(π/6)-ϕ2]-Sin[(π/6)-ϕ3])/((Sqrt[2+Z-2 Sin[(π/6)-ϕ2]]) (Sqrt[2+Z-2 Sin[(π/6)-ϕ3]])),
- ((1+Z+Cos[ϕ2+2 ϕ3]-Sin[(π/6)-ϕ3]-Sin[(π/6)+ϕ2+ϕ3])/((Sqrt[2+Z-2 Sin[(π/6)-ϕ3]]) (Sqrt[2+Z-2 Sin[(π/6)+ϕ2+ϕ3]]))),
- ((1+Z+Cos[2 ϕ2+ϕ3]-Sin[(π/6)-ϕ2]-Sin[(π/6)+ϕ2+ϕ3])/((Sqrt[2+Z-2 Sin[(π/6)-ϕ2]]) (Sqrt[2+Z-2 Sin[(π/6)+ϕ2+ϕ3]])))
}, {Z,ZL, ZH}, PlotRange->{{c1L,c1H},{c2L,c2H},{c3L,c3H}}, PlotPoints->ppc, PlotStyle->  {edgecolor, Thickness[et]} ];

edge21 = ParametricPlot3D[  {
- ((1+Z+Cos[ϕ2]+Cos[ϕ2-ϕ3]+Cos[ϕ3])/((Sqrt[2+Z+2 Cos[ϕ2]]) (Sqrt[2+Z+2 Cos[ϕ3]]))),
(1+Z+Cos[ϕ3]+Cos[ϕ2+ϕ3]+Cos[ϕ2+2 ϕ3])/((Sqrt[2+Z+2 Cos[ϕ3]]) (Sqrt[2+Z+2 Cos[ϕ2+ϕ3]])),
- ((1+Z+Cos[ϕ2]+Cos[ϕ2+ϕ3]+Cos[2 ϕ2+ϕ3])/((Sqrt[2+Z+2 Cos[ϕ2]]) (Sqrt[2+Z+2 Cos[ϕ2+ϕ3]])))
},{Z,ZL, ZH}, PlotRange->{{c1L,c1H},{c2L,c2H},{c3L,c3H}}, PlotPoints->ppc, PlotStyle->  {edgecolor, Thickness[et]} ];

edge22 = ParametricPlot3D[  {
- ((1+Z+Cos[ϕ2-ϕ3]-Sin[(π/6)+ϕ2]-Sin[(π/6)+ϕ3])/((Sqrt[2+Z-2 Sin[(π/6)+ϕ2]]) (Sqrt[2+Z-2 Sin[(π/6)+ϕ3]]))),
(1+Z+Cos[ϕ2+2 ϕ3]-Sin[(π/6)-ϕ2-ϕ3]-Sin[(π/6)+ϕ3])/((Sqrt[2+Z-2 Sin[(π/6)-ϕ2-ϕ3]]) (Sqrt[2+Z-2 Sin[(π/6)+ϕ3]])),
- ((1+Z+Cos[2 ϕ2+ϕ3]-Sin[(π/6)+ϕ2]-Sin[(π/6)-ϕ2-ϕ3])/((Sqrt[2+Z-2 Sin[(π/6)+ϕ2]]) (Sqrt[2+Z-2 Sin[(π/6)-ϕ2-ϕ3]])))
}, {Z,ZL, ZH}, PlotRange->{{c1L,c1H},{c2L,c2H},{c3L,c3H}}, PlotPoints->ppc, PlotStyle->  {edgecolor, Thickness[et]} ];

edge23 = ParametricPlot3D[  {
- ((1+Z+Cos[ϕ2-ϕ3]-Sin[(π/6)-ϕ2]-Sin[(π/6)-ϕ3])/((Sqrt[2+Z-2 Sin[(π/6)-ϕ2]]) (Sqrt[2+Z-2 Sin[(π/6)-ϕ3]]))),
(1+Z+Cos[ϕ2+2 ϕ3]-Sin[(π/6)-ϕ3]-Sin[(π/6)+ϕ2+ϕ3])/((Sqrt[2+Z-2 Sin[(π/6)-ϕ3]]) (Sqrt[2+Z-2 Sin[(π/6)+ϕ2+ϕ3]])),
- ((1+Z+Cos[2 ϕ2+ϕ3]-Sin[(π/6)-ϕ2]-Sin[(π/6)+ϕ2+ϕ3])/((Sqrt[2+Z-2 Sin[(π/6)-ϕ2]]) (Sqrt[2+Z-2 Sin[(π/6)+ϕ2+ϕ3]])))
}, {Z,ZL, ZH}, PlotRange->{{c1L,c1H},{c2L,c2H},{c3L,c3H}}, PlotPoints->ppc, PlotStyle->  {edgecolor, Thickness[et]} ];

edge31 = ParametricPlot3D[  {
- ((1+Z+Cos[ϕ2]+Cos[ϕ2-ϕ3]+Cos[ϕ3])/((Sqrt[2+Z+2 Cos[ϕ2]]) (Sqrt[2+Z+2 Cos[ϕ3]]))),
- ((1+Z+Cos[ϕ3]+Cos[ϕ2+ϕ3]+Cos[ϕ2+2 ϕ3])/((Sqrt[2+Z+2 Cos[ϕ3]]) (Sqrt[2+Z+2 Cos[ϕ2+ϕ3]]))),
(1+Z+Cos[ϕ2]+Cos[ϕ2+ϕ3]+Cos[2 ϕ2+ϕ3])/((Sqrt[2+Z+2 Cos[ϕ2]]) (Sqrt[2+Z+2 Cos[ϕ2+ϕ3]]))
},{Z,ZL, ZH}, PlotRange->{{c1L,c1H},{c2L,c2H},{c3L,c3H}}, PlotPoints->ppc, PlotStyle->  {edgecolor, Thickness[et]} ];

edge32 = ParametricPlot3D[  {
- ((1+Z+Cos[ϕ2-ϕ3]-Sin[(π/6)+ϕ2]-Sin[(π/6)+ϕ3])/((Sqrt[2+Z-2 Sin[(π/6)+ϕ2]]) (Sqrt[2+Z-2 Sin[(π/6)+ϕ3]]))),
- ((1+Z+Cos[ϕ2+2 ϕ3]-Sin[(π/6)-ϕ2-ϕ3]-Sin[(π/6)+ϕ3])/((Sqrt[2+Z-2 Sin[(π/6)-ϕ2-ϕ3]]) (Sqrt[2+Z-2 Sin[(π/6)+ϕ3]]))),
(1+Z+Cos[2 ϕ2+ϕ3]-Sin[(π/6)+ϕ2]-Sin[(π/6)-ϕ2-ϕ3])/((Sqrt[2+Z-2 Sin[(π/6)+ϕ2]]) (Sqrt[2+Z-2 Sin[(π/6)-ϕ2-ϕ3]]))
}, {Z,ZL, ZH}, PlotRange->{{c1L,c1H},{c2L,c2H},{c3L,c3H}}, PlotPoints->ppc, PlotStyle->  {edgecolor, Thickness[et]} ];

edge33 = ParametricPlot3D[  {
- ((1+Z+Cos[ϕ2-ϕ3]-Sin[(π/6)-ϕ2]-Sin[(π/6)-ϕ3])/((Sqrt[2+Z-2 Sin[(π/6)-ϕ2]]) (Sqrt[2+Z-2 Sin[(π/6)-ϕ3]]))),
- ((1+Z+Cos[ϕ2+2 ϕ3]-Sin[(π/6)-ϕ3]-Sin[(π/6)+ϕ2+ϕ3])/((Sqrt[2+Z-2 Sin[(π/6)-ϕ3]]) (Sqrt[2+Z-2 Sin[(π/6)+ϕ2+ϕ3]]))),
(1+Z+Cos[2 ϕ2+ϕ3]-Sin[(π/6)-ϕ2]-Sin[(π/6)+ϕ2+ϕ3])/((Sqrt[2+Z-2 Sin[(π/6)-ϕ2]]) (Sqrt[2+Z-2 Sin[(π/6)+ϕ2+ϕ3]]))
}, {Z,ZL, ZH}, PlotRange->{{c1L,c1H},{c2L,c2H},{c3L,c3H}}, PlotPoints->ppc, PlotStyle->  {edgecolor, Thickness[et]} ];

points = ListPointPlot3D[{{-Cos[AA],Cos[BB],Cos[CC]},{Cos[AA],-Cos[BB],Cos[CC]},{Cos[AA],Cos[BB],-Cos[CC]},{-Cos[AA],-Cos[BB],-Cos[CC]}} , PlotStyle-> {Black,PointSize[.02]}];

lines = ListLinePlot3D[{
{{-Cos[AA],-Cos[BB],-Cos[CC]},{-Cos[AA],Cos[BB],Cos[CC]}},{{-Cos[AA],-Cos[BB],-Cos[CC]},{Cos[AA],-Cos[BB],Cos[CC]}},{{-Cos[AA],-Cos[BB],-Cos[CC]},{Cos[AA],Cos[BB],-Cos[CC]}},
{{-Cos[AA],Cos[BB],Cos[CC]},{Cos[AA],-Cos[BB],Cos[CC]},{Cos[AA],Cos[BB],-Cos[CC]},{-Cos[AA],Cos[BB],Cos[CC]}}
} , PlotStyle->{{tetcolor,Thickness[.006]}} ];

infTet = ContourPlot3D[ 1-c1^2-c2^2-c3^2+2c1 c2 c3 == 0 , {c1,c1L,c1H}, {c2,c2L,c2H}, {c3,c3L,c3H}, ContourStyle->{Gray,Opacity[.05]}, Mesh->None]; 

showlist = {} ;
If[ surface0, AppendTo[showlist, plot0] ];
If[ surface1, AppendTo[showlist, plot1] ];
If[ surface2, AppendTo[showlist, plot2] ];
If[ surface3, AppendTo[showlist, plot3] ];
If[ inflTetr, AppendTo[showlist, infTet] ];
If[ belt, showlist = Join[showlist, {curve0  ,curve1,curve2,curve3  }] ];
If[ edges, showlist = Join[showlist, {edge01, edge02,edge03,edge11, edge12,edge13,edge21, edge22,edge23,edge31, edge32,edge33}] ]; 
If[ tetrahedron,showlist = Join[showlist, {lines (*, points*)}] ]; 
If[ label, 
showlabel = Style[ ToString[SetPrecision[ϕ1/π,3]] <> "π  " <> ToString[SetPrecision[ϕ2/π,3]] <> "π  " <> ToString[SetPrecision[ϕ3/π,3]] <> "π\n" <> 
ToString[SetPrecision[ϕ1,3]] <> "  " <> ToString[SetPrecision[ϕ2,3]] <> "  " <> ToString[SetPrecision[ϕ3,3]] <> "\n" <> 
ToString[SetPrecision[AA/π,3]] <> "π  " <> ToString[SetPrecision[BB/π,3]] <> "π  " <> ToString[SetPrecision[CC/π,3]] <> "π\n" <> 
type <>"\n" <> 
ToString[SetPrecision[θL/π,3]] <> "π  " <> 
ToString[SetPrecision[θH/π,3]] <> "π\n" <> 
ToString[SetPrecision[θL,3]] <> "  " <> 
ToString[SetPrecision[θH,3]] <> "\n" <> 
ToString[SetPrecision[c1L,3]] <> "  " <> 
ToString[SetPrecision[c1H,3]] <> "  " <> 
ToString[SetPrecision[c2L,3]] <> "  " <> 
ToString[SetPrecision[c2H,3]] <> "  " <> 
ToString[SetPrecision[c3L,3]] <> "  " <> 
ToString[SetPrecision[c3H,3]] , 16 ], 
showlabel = "" ];
Show[ showlist , Axes->axes, AxesLabel->{c1,c2,c3}, Boxed->box, BoxRatios->{br,1,1}, PlotLabel-> showlabel, ImageSize->600 ]
] ,    {{ϕ1,.20 Pi},0,N[Pi]},{{ϕ2,.35 Pi},0,N[Pi]}, 
{{c1L,-1},-1,1}, {{c1H,1},-1,1}, {{c2L,-1},-1,1}, {{c2H,1},-1,1}, {{c3L,-1},-1,1}, {{c3H,1},-1,1}, {{θL,-N[Pi]+.001},-N[Pi],N[Pi]},{{θH,N[Pi]+.001},-N[Pi],N[Pi]+.001},{{ZL,0},0,10},{{ZH,10},0,10},{br,1},{ppc,200},{pps,100}, {op,.5},{bt,.005}, {et,.003}, {mesh, 0}, {surface0, {True,False}},{surface1, {True,False}},{surface2, {True,False}},{surface3, {True,False}},{inflTetr, {False,True}},{tetrahedron, {False,True}},{belt, {True,False}},{edges, {True,False}},{axes, {False,True}},{box, {False, True}},{label, {True,False}},ControlPlacement->Right ]

4.  Afterthoughts

*   The proposed solutions to LO and RA can written be in terms of polynomial inequalities
          in c1 , c2 and c3 , without explicit mention of α , β nor γ . 

     *   More precisely, each of the proposed condition for can be written as a disjunction of 
          conjunctions of such polynomial inequalities.

     *   Perhaps something like the Wu-Ritt decomposition of a system of polynomial equations 
          could be used to establish the correctness of the proposed solutions.

*   The Perspective 3-Point Problem is similar to the Three Rays Problem, but its focus is 
          essentially on finding all the ways to choose the triangle △ABC, and solving it requires 
          extracting the roots of a polynomial of degree greater than two.

*   As I pointed out in an unpublished paper, presented at a regional MAA meeting, the 
          P3P Problem can be reformulated as a problem on a sphere involving four great circles. 

     *   The results here concerning the Three Lines Problem can be interpreted in this context. 

     *   This in turn produces results in the “elliptic plane,” i.e., in the projective plane with a metric
          inherited from the sphere. 

     *   The spherical interpretation of the triviality of the LA Problem can be used to make a
          somewhat interesting puzzle, involving a sphere with three great circles drawn on it, 
          along with a movable ring (same size as the circles) with markings on it. The goal is to place 
          the ring on the sphere so that it becomes a great circle, with its markings on the other circles.

References

-  R. M. Haralick, C.-N. Lee, K. Ottenberg and N. Nölle, “Review and Analysis of Solutions of the 
    Three Point Perspective Pose Estimation Problem” (1997), J. Computer Vision, v. 13, n. 3, pp. 
    331-356.
    
-  M. Q. Rieck, “On the Discriminant of Grunert’s System of Algebraic Equations and Related
   Topics” (2018), J. Math. Imaging and Vision, v. 60, n. 5, pp. 737-762.
   
-  M. Q. Rieck, B. Wang, “Locating Perspective Three-Point Problem Solutions in Spatial Regions”
    (2021)  J. Math. Imaging and Vision, v. 63, n. 8, pp. 953-973.
    
-  J. E. Wetzel, “An Ancient Elliptic Locus” (2010), Amer. Math. Monthly, v. 117, n. 2, pp. 161-167.

-   https://github.com/mqrieck/tetrahedron_test.cpp/

(*** NEW CODE SECTION ***)

(* Multiple slices for the Lines and Rays Problems *) 
Manipulate[ Module[
{a,BB,CC,costurn,sinturn,x10,x20,x30,y10,y20,i,plots,y30,infTet,x1,x2,x3,y1,y2,y3,t1,t2,t3,DD,HH,LL,RR,EE, region, plot, L, R, c3, type, ϕ1,ϕ2,ϕ3}, 
c3 = Cos[ϕ]; 
a = Sqrt[b^2+c^2-2b c Cos[AA]];
BB = ArcCos[(a^2+c^2-b^2)/(2a c)];
CC = ArcCos[(a^2+b^2-c^2)/(2a b)];
If[ AA > BB+CC , type = "obtuseA", 
If[ BB > CC+AA, type = "obtuseB", 
If[ CC > AA+BB, type = "obtuseC", 
type = "acute" ]]];
ϕ1 = 2(AA+2BB)/3; ϕ2 = -2(2AA+BB)/3; ϕ3 = 2(AA-BB)/3; 
x1 = Cos[ϕ1]; x2 = Cos[ϕ2]; x3 = Cos[ϕ3]; 
y1 = Sin[ϕ1]; y2 = Sin[ϕ2]; y3 = Sin[ϕ3]; 
t1 = Tan[ϕ1/2]; t2 = Tan[ϕ2/2]; t3 = Tan[ϕ3/2]; 
infTet = ContourPlot3D[ 1-c1^2-c2^2-c3^2+2c1 c2 c3 == 0 , {c1,-1,1}, {c2,-1,1}, {c3,-1,1}, ContourStyle->{Gray,Opacity[.05]}, Mesh->None]; 
plots = {infTet}; 
For[i=0, i<=n, i++, ϕ = N[π i/n]; 
If[ useRays,
plot = ParametricPlot3D[ {
(c Cos[ϕ] Sin[θ]+b (Cos[AA] Cos[θ]-Cos[ψ] Sin[AA] Sin[θ]) Sin[ϕ])/Sqrt[(c^2) ((Sin[θ])^2)+(b^2) ((Sin[ϕ])^2)+2 b c Sin[θ] Sin[ϕ] (Cos[AA] Cos[θ+ϕ]-Cos[ψ] Sin[AA] Sin[θ+ϕ])],   
(c Sin[θ]+b Sin[ϕ] (Cos[AA] Cos[θ+ϕ]-Cos[ψ] Sin[AA] Sin[θ+ϕ])) /Sqrt[((c^2) ((Sin[θ])^2)+(b^2) ((Sin[ϕ])^2)+2 b c Sin[θ] Sin[ϕ] (Cos[AA] Cos[θ+ϕ]-Cos[ψ] Sin[AA] Sin[θ+ϕ]))], 
Cos[ϕ]
}, {ψ,-N[π]+.0001,N[π]-.0001}, {θ,.0001,N[π]-ϕ-.0001},PlotRange->{{-1.1,1.1},{-1.1,1.1},{-1.1,1.1}}, PlotStyle->{Black,Opacity[op]}, PlotPoints->pp ] , 
 plot = ParametricPlot3D[ {
(c Cos[ϕ] Sin[θ]+b (Cos[AA] Cos[θ]-Cos[ψ] Sin[AA] Sin[θ]) Sin[ϕ])/Sqrt[(c^2) ((Sin[θ])^2)+(b^2) ((Sin[ϕ])^2)+2 b c Sin[θ] Sin[ϕ] (Cos[AA] Cos[θ+ϕ]-Cos[ψ] Sin[AA] Sin[θ+ϕ])],
(c Sin[θ]+b Sin[ϕ] (Cos[AA] Cos[θ+ϕ]-Cos[ψ] Sin[AA] Sin[θ+ϕ])) /Sqrt[((c^2) ((Sin[θ])^2)+(b^2) ((Sin[ϕ])^2)+2 b c Sin[θ] Sin[ϕ] (Cos[AA] Cos[θ+ϕ]-Cos[ψ] Sin[AA] Sin[θ+ϕ]))],
Cos[ϕ]
}, {ψ,-N[π]+.0001,N[π]-.0001}, {θ,.0001,N[2π]-.0001}, PlotRange->{{-1.1,1.1},{-1.1,1.1},{-1.1,1.1}}, PlotStyle->{Black,Opacity[op]}, PlotPoints->pp ]
];
plots = Append[plots,plot]
];
Show[ plots, ImageSize->800,Axes->False,Boxed->False,Frame->False,MeshStyle->None,PlotLabel-> Style[ 
(* "phis: " <> ToString[SetPrecision[ϕ1/π,3]] <> "π " <> ToString[SetPrecision[ϕ2/π,3]] <> "π " <> ToString[SetPrecision[ϕ3/π,3]] <> "π\n" <> *)
"side lengths: " <> ToString[SetPrecision[a,3]] <> "  " <> ToString[SetPrecision[b,3]]<> "  " <> ToString[SetPrecision[c,3]] <> "\ninterior angles: "<>ToString[SetPrecision[AA/π,3]] <> "π  " <> ToString[SetPrecision[BB/π,3]]<> "π  " <> ToString[SetPrecision[CC/π,3]]<>"π" <> "\n" <> type  , 20 ]] 
] , {{b,2},.1,5},{{c,3},.1,5},{{AA,.2π},0,N[π]},{{op,.5},0,1},{n,16}, 
{pp,50},{useRays,{True,False}}, ControlPlacement->Right ]

(* Multiple slices in three directions for the Lines and Rays Problems *) 
Manipulate[ Module[
{a,BB,CC,costurn,sinturn,x10,x20,x30,y10,y20,i,plots,y30,infTet,e,x1,x2,x3,y1,y2,y3,t1,t2,t3,DD,HH,LL,RR,EE, region, plot, L, R, c3, type, ϕ1,ϕ2,ϕ3, 
cA,cB,cC,Lims}, 
c3 = Cos[ϕ]; 
a = Sqrt[b^2+c^2-2b c Cos[AA]];
BB = ArcCos[(a^2+c^2-b^2)/(2a c)];
CC = ArcCos[(a^2+b^2-c^2)/(2a b)];
If[ AA > BB+CC , type = "obtuseA", 
If[ BB > CC+AA, type = "obtuseB", 
If[ CC > AA+BB, type = "obtuseC", 
type = "acute" ]]];
ϕ1 = 2(AA+2BB)/3; ϕ2 = -2(2AA+BB)/3; ϕ3 = 2(AA-BB)/3; 
x1 = Cos[ϕ1]; x2 = Cos[ϕ2]; x3 = Cos[ϕ3]; 
y1 = Sin[ϕ1]; y2 = Sin[ϕ2]; y3 = Sin[ϕ3]; 
t1 = Tan[ϕ1/2]; t2 = Tan[ϕ2/2]; t3 = Tan[ϕ3/2]; 
cA = Cos[AA]; cB = Cos[BB]; cC = Cos[CC]; 
plots = {};
If[showInfTet, 
infTet = ContourPlot3D[ 1-c1^2-c2^2-c3^2+2c1 c2 c3 == 0 , {c1,c1L,c1H}, {c2,c2L,c2H}, {c3,c3L,c3H}, ContourStyle->{{Gray,Opacity[op2]}}, Mesh->None,
PlotPoints->pp]; 
plots = Append[plots,infTet]
];
If[showLims, 
Lims = ContourPlot3D[ {c1 == cA,c2 == cB,c3 == cC} , {c1,c1L,c1H}, {c2,c2L,c2H}, {c3,c3L,c3H}, 
RegionFunction->Function[{x,y,z},1+2 x y z - x^2 -y^2-z^2>=0 ], 
RegionBoundaryStyle->None, 
ContourStyle->{{Yellow,Opacity[op3]},{Cyan,Opacity[op3]},{Magenta,Opacity[op3]}}, Mesh->None,
PlotPoints->pp ]; 
plots = Append[plots,Lims]
];
If[ useEnds, e=0,e=1 ];
For[i=e, i<=n-e, i++, ϕ = N[π i/n]; 
If[ useRays,
plot = ParametricPlot3D[ {
(c Cos[ϕ] Sin[θ]+b (Cos[AA] Cos[θ]-Cos[ψ] Sin[AA] Sin[θ]) Sin[ϕ])/Sqrt[(c^2) ((Sin[θ])^2)+(b^2) ((Sin[ϕ])^2)+2 b c Sin[θ] Sin[ϕ] (Cos[AA] Cos[θ+ϕ]-Cos[ψ] Sin[AA] Sin[θ+ϕ])],   
(c Sin[θ]+b Sin[ϕ] (Cos[AA] Cos[θ+ϕ]-Cos[ψ] Sin[AA] Sin[θ+ϕ])) /Sqrt[((c^2) ((Sin[θ])^2)+(b^2) ((Sin[ϕ])^2)+2 b c Sin[θ] Sin[ϕ] (Cos[AA] Cos[θ+ϕ]-Cos[ψ] Sin[AA] Sin[θ+ϕ]))], 
Cos[ϕ]
}, {ψ,-N[π]+.0001,N[π]-.0001}, {θ,.0001,N[π]-ϕ-.0001},PlotRange->{{c1L,c1H},{c2L,c2H},{c3L,c3H}}, PlotStyle->{Red,Opacity[op1]}, PlotPoints->pp ] , 
 plot = ParametricPlot3D[ {
(c Cos[ϕ] Sin[θ]+b (Cos[AA] Cos[θ]-Cos[ψ] Sin[AA] Sin[θ]) Sin[ϕ])/Sqrt[(c^2) ((Sin[θ])^2)+(b^2) ((Sin[ϕ])^2)+2 b c Sin[θ] Sin[ϕ] (Cos[AA] Cos[θ+ϕ]-Cos[ψ] Sin[AA] Sin[θ+ϕ])],
(c Sin[θ]+b Sin[ϕ] (Cos[AA] Cos[θ+ϕ]-Cos[ψ] Sin[AA] Sin[θ+ϕ])) /Sqrt[((c^2) ((Sin[θ])^2)+(b^2) ((Sin[ϕ])^2)+2 b c Sin[θ] Sin[ϕ] (Cos[AA] Cos[θ+ϕ]-Cos[ψ] Sin[AA] Sin[θ+ϕ]))],
Cos[ϕ]
}, {ψ,-N[π]+.0001,N[π]-.0001}, {θ,.0001,N[2π]-.0001}, PlotRange->{{c1L,c1H},{c2L,c2H},{c3L,c3H}}, PlotStyle->{Red,Opacity[op1]}, PlotPoints->pp ]
];
plots = Append[plots,plot]
];
For[i=e, i<=n-e, i++, ϕ = N[π i/n]; 
If[ useRays,
plot = ParametricPlot3D[ {
Cos[ϕ], 
(a Cos[ϕ] Sin[θ]+c (Cos[BB] Cos[θ]-Cos[ψ] Sin[BB] Sin[θ]) Sin[ϕ])/Sqrt[(a^2) ((Sin[θ])^2)+(c^2) ((Sin[ϕ])^2)+2 c a Sin[θ] Sin[ϕ] (Cos[BB] Cos[θ+ϕ]-Cos[ψ] Sin[BB] Sin[θ+ϕ])],   
(a Sin[θ]+c Sin[ϕ] (Cos[BB] Cos[θ+ϕ]-Cos[ψ] Sin[BB] Sin[θ+ϕ])) /Sqrt[((a^2) ((Sin[θ])^2)+(c^2) ((Sin[ϕ])^2)+2 c a Sin[θ] Sin[ϕ] (Cos[BB] Cos[θ+ϕ]-Cos[ψ] Sin[BB] Sin[θ+ϕ]))]
}, {ψ,-N[π]+.0001,N[π]-.0001}, {θ,.0001,N[π]-ϕ-.0001},PlotRange->{{c1L,c1H},{c2L,c2H},{c3L,c3H}}, PlotStyle->{Darker[Green],Opacity[op1]}, PlotPoints->pp ] , 
 plot = ParametricPlot3D[ {
Cos[ϕ],
(a Cos[ϕ] Sin[θ]+c (Cos[BB] Cos[θ]-Cos[ψ] Sin[BB] Sin[θ]) Sin[ϕ])/Sqrt[(a^2) ((Sin[θ])^2)+(c^2) ((Sin[ϕ])^2)+2 c a Sin[θ] Sin[ϕ] (Cos[BB] Cos[θ+ϕ]-Cos[ψ] Sin[BB] Sin[θ+ϕ])],
(a Sin[θ]+c Sin[ϕ] (Cos[BB] Cos[θ+ϕ]-Cos[ψ] Sin[BB] Sin[θ+ϕ])) /Sqrt[((a^2) ((Sin[θ])^2)+(c^2) ((Sin[ϕ])^2)+2 c a Sin[θ] Sin[ϕ] (Cos[BB] Cos[θ+ϕ]-Cos[ψ] Sin[BB] Sin[θ+ϕ]))]
}, {ψ,-N[π]+.0001,N[π]-.0001}, {θ,.0001,N[2π]-.0001}, PlotRange->{{c1L,c1H},{c2L,c2H},{c3L,c3H}}, PlotStyle->{Darker[Green],Opacity[op1]}, PlotPoints->pp ]
];
plots = Append[plots,plot]
];
For[i=e, i<=n-e, i++, ϕ = N[π i/n]; 
If[ useRays,
plot = ParametricPlot3D[ {
(b Sin[θ]+a Sin[ϕ] (Cos[CC] Cos[θ+ϕ]-Cos[ψ] Sin[CC] Sin[θ+ϕ])) /Sqrt[((b^2) ((Sin[θ])^2)+(a^2) ((Sin[ϕ])^2)+2 a b Sin[θ] Sin[ϕ] (Cos[CC] Cos[θ+ϕ]-Cos[ψ] Sin[CC] Sin[θ+ϕ]))],Cos[ϕ], 
(b Cos[ϕ] Sin[θ]+a (Cos[CC] Cos[θ]-Cos[ψ] Sin[CC] Sin[θ]) Sin[ϕ])/Sqrt[(b^2) ((Sin[θ])^2)+(a^2) ((Sin[ϕ])^2)+2 a b Sin[θ] Sin[ϕ] (Cos[CC] Cos[θ+ϕ]-Cos[ψ] Sin[CC] Sin[θ+ϕ])]
}, {ψ,-N[π]+.0001,N[π]-.0001}, {θ,.0001,N[π]-ϕ-.0001},PlotRange->{{c1L,c1H},{c2L,c2H},{c3L,c3H}}, PlotStyle->{Blue,Opacity[op1]}, PlotPoints->pp ] , 
 plot = ParametricPlot3D[ {
(b Sin[θ]+a Sin[ϕ] (Cos[CC] Cos[θ+ϕ]-Cos[ψ] Sin[CC] Sin[θ+ϕ])) /Sqrt[((b^2) ((Sin[θ])^2)+(a^2) ((Sin[ϕ])^2)+2 a b Sin[θ] Sin[ϕ] (Cos[CC] Cos[θ+ϕ]-Cos[ψ] Sin[CC] Sin[θ+ϕ]))],
Cos[ϕ],
(b Cos[ϕ] Sin[θ]+a (Cos[CC] Cos[θ]-Cos[ψ] Sin[CC] Sin[θ]) Sin[ϕ])/Sqrt[(b^2) ((Sin[θ])^2)+(a^2) ((Sin[ϕ])^2)+2 a b Sin[θ] Sin[ϕ] (Cos[CC] Cos[θ+ϕ]-Cos[ψ] Sin[CC] Sin[θ+ϕ])]
}, {ψ,-N[π]+.0001,N[π]-.0001}, {θ,.0001,N[2π]-.0001}, PlotRange->{{c1L,c1H},{c2L,c2H},{c3L,c3H}}, PlotStyle->{Blue,Opacity[op1]}, PlotPoints->pp ]
];
plots = Append[plots,plot]
];
Show[ plots, ImageSize->800,Axes->False,Boxed->False,Frame->False,MeshStyle->None,BoxRatios->{1,1,1},PlotLabel-> Style[ 
(* "phis: " <> ToString[SetPrecision[ϕ1/π,3]] <> "π " <> ToString[SetPrecision[ϕ2/π,3]] <> "π " <> ToString[SetPrecision[ϕ3/π,3]] <> "π\n" <> *)
"side lengths: " <> ToString[SetPrecision[a,3]] <> "  " <> ToString[SetPrecision[b,3]]<> "  " <> ToString[SetPrecision[c,3]] <> 
"\nvertices: ("<>ToString[SetPrecision[x1,3]] <>"," <> ToString[SetPrecision[y1,3]] <> ") (" <> 
ToString[SetPrecision[x2,3]] <>"," <> ToString[SetPrecision[y2,3]] <> ") (" <> 
ToString[SetPrecision[x3,3]] <>"," <> ToString[SetPrecision[y3,3]] <>
 ")\ninterior angles: "<>ToString[SetPrecision[AA/π,3]] <> "π  " <> ToString[SetPrecision[BB/π,3]]<> "π  " <> ToString[SetPrecision[CC/π,3]]<>"π" <> "\ntheir cosines: "<>ToString[SetPrecision[cA,3]] <> "  " <> ToString[SetPrecision[cB,3]]<> "  " <> ToString[SetPrecision[cC,3]] <> "\n" <> type  , 20 ]] 
] , {{b,3},.1,5},{{c,3},.1,5},{{AA,π/3.},0,N[π]},
{{c1L,-1},-1,1},{{c1H,1},-1,1}, {{c2L,-1},-1,1},{{c2H,1},-1,1}, 
{{c3L,-1},-1,1},{{c3H,1},-1,1},{{op1,.5},0,1},{{op2,.1},0,1},{{op3,.5},0,1},
{n,30}, {pp,50},{useRays,{True,False}},{useEnds,{True,False}},{showInfTet,{True,False}},{showLims,{True,False}}, ControlPlacement->Right ]

(* 
  Taking ArcCos's, we obtain the "spikes" previously generated using very
 different methods
*) 
Manipulate[ Module[
{a,BB,CC,costurn,sinturn,x10,x20,x30,y10,y20,i,plots,y30,lines,e,x1,x2,x3,y1,y2,y3,t1,t2,t3,DD,HH,LL,RR,EE, region, plot, L, R, c3, type, ϕ1,ϕ2,ϕ3}, 
c3 = Cos[ϕ]; 
a = Sqrt[b^2+c^2-2b c Cos[AA]];
BB = ArcCos[(a^2+c^2-b^2)/(2a c)];
CC = ArcCos[(a^2+b^2-c^2)/(2a b)];
If[ AA > BB+CC , type = "obtuseA", 
If[ BB > CC+AA, type = "obtuseB", 
If[ CC > AA+BB, type = "obtuseC", 
type = "acute" ]]];
ϕ1 = 2(AA+2BB)/3; ϕ2 = -2(2AA+BB)/3; ϕ3 = 2(AA-BB)/3; 
x1 = Cos[ϕ1]; x2 = Cos[ϕ2]; x3 = Cos[ϕ3]; 
y1 = Sin[ϕ1]; y2 = Sin[ϕ2]; y3 = Sin[ϕ3]; 
t1 = Tan[ϕ1/2]; t2 = Tan[ϕ2/2]; t3 = Tan[ϕ3/2]; 
lines = ListLinePlot3D[{{{0,0,0},{0,π,π}},{{0,0,0},{π,0,π}},{{0,0,0},{π,π,0}},{{0,π,π},{π,0,π}},{{π,0,π},{π,π,0}},{{π,π,0},{0,π,π}}},
PlotRange->{{c1L,c1H},{c2L,c2H},{c3L,c3H}}, PlotStyle->{LightGray}];
plots = {lines};
If[ useEnds, e=0,e=1 ];
For[i=e, i<=n-e, i++, ϕ = N[π i/n]; 
If[ useRays,
plot = ParametricPlot3D[ {
ArcCos[(c Cos[ϕ] Sin[θ]+b (Cos[AA] Cos[θ]-Cos[ψ] Sin[AA] Sin[θ]) Sin[ϕ])/Sqrt[(c^2) ((Sin[θ])^2)+(b^2) ((Sin[ϕ])^2)+2 b c Sin[θ] Sin[ϕ] (Cos[AA] Cos[θ+ϕ]-Cos[ψ] Sin[AA] Sin[θ+ϕ])]],   
ArcCos[(c Sin[θ]+b Sin[ϕ] (Cos[AA] Cos[θ+ϕ]-Cos[ψ] Sin[AA] Sin[θ+ϕ])) /Sqrt[((c^2) ((Sin[θ])^2)+(b^2) ((Sin[ϕ])^2)+2 b c Sin[θ] Sin[ϕ] (Cos[AA] Cos[θ+ϕ]-Cos[ψ] Sin[AA] Sin[θ+ϕ]))]], 
ArcCos[Cos[ϕ]]
}, {ψ,-N[π]+.0001,N[π]-.0001}, {θ,.0001,N[π]-ϕ-.0001},PlotRange->{{c1L,c1H},{c2L,c2H},{c3L,c3H}}, PlotStyle->{Red,Opacity[op]}, PlotPoints->pp ] , 
 plot = ParametricPlot3D[ {
ArcCos[(c Cos[ϕ] Sin[θ]+b (Cos[AA] Cos[θ]-Cos[ψ] Sin[AA] Sin[θ]) Sin[ϕ])/Sqrt[(c^2) ((Sin[θ])^2)+(b^2) ((Sin[ϕ])^2)+2 b c Sin[θ] Sin[ϕ] (Cos[AA] Cos[θ+ϕ]-Cos[ψ] Sin[AA] Sin[θ+ϕ])]],
ArcCos[(c Sin[θ]+b Sin[ϕ] (Cos[AA] Cos[θ+ϕ]-Cos[ψ] Sin[AA] Sin[θ+ϕ])) /Sqrt[((c^2) ((Sin[θ])^2)+(b^2) ((Sin[ϕ])^2)+2 b c Sin[θ] Sin[ϕ] (Cos[AA] Cos[θ+ϕ]-Cos[ψ] Sin[AA] Sin[θ+ϕ]))]],
ArcCos[Cos[ϕ]]
}, {ψ,-N[π]+.0001,N[π]-.0001}, {θ,.0001,N[2π]-.0001}, PlotRange->{{c1L,c1H},{c2L,c2H},{c3L,c3H}}, PlotStyle->{Red,Opacity[op]}, PlotPoints->pp ]
];
plots = Append[plots,plot]
];
For[i=e, i<=n-e, i++, ϕ = N[π i/n]; 
If[ useRays,
plot = ParametricPlot3D[ {
ArcCos[Cos[ϕ]], 
ArcCos[(a Cos[ϕ] Sin[θ]+c (Cos[BB] Cos[θ]-Cos[ψ] Sin[BB] Sin[θ]) Sin[ϕ])/Sqrt[(a^2) ((Sin[θ])^2)+(c^2) ((Sin[ϕ])^2)+2 c a Sin[θ] Sin[ϕ] (Cos[BB] Cos[θ+ϕ]-Cos[ψ] Sin[BB] Sin[θ+ϕ])]],   
ArcCos[(a Sin[θ]+c Sin[ϕ] (Cos[BB] Cos[θ+ϕ]-Cos[ψ] Sin[BB] Sin[θ+ϕ])) /Sqrt[((a^2) ((Sin[θ])^2)+(c^2) ((Sin[ϕ])^2)+2 c a Sin[θ] Sin[ϕ] (Cos[BB] Cos[θ+ϕ]-Cos[ψ] Sin[BB] Sin[θ+ϕ]))]]
}, {ψ,-N[π]+.0001,N[π]-.0001}, {θ,.0001,N[π]-ϕ-.0001},PlotRange->{{c1L,c1H},{c2L,c2H},{c3L,c3H}}, PlotStyle->{Darker[Green],Opacity[op]}, PlotPoints->pp ] , 
 plot = ParametricPlot3D[ {
ArcCos[Cos[ϕ]],
ArcCos[(a Cos[ϕ] Sin[θ]+c (Cos[BB] Cos[θ]-Cos[ψ] Sin[BB] Sin[θ]) Sin[ϕ])/Sqrt[(a^2) ((Sin[θ])^2)+(c^2) ((Sin[ϕ])^2)+2 c a Sin[θ] Sin[ϕ] (Cos[BB] Cos[θ+ϕ]-Cos[ψ] Sin[BB] Sin[θ+ϕ])]],
ArcCos[(a Sin[θ]+c Sin[ϕ] (Cos[BB] Cos[θ+ϕ]-Cos[ψ] Sin[BB] Sin[θ+ϕ])) /Sqrt[((a^2) ((Sin[θ])^2)+(c^2) ((Sin[ϕ])^2)+2 c a Sin[θ] Sin[ϕ] (Cos[BB] Cos[θ+ϕ]-Cos[ψ] Sin[BB] Sin[θ+ϕ]))]]
}, {ψ,-N[π]+.0001,N[π]-.0001}, {θ,.0001,N[2π]-.0001}, PlotRange->{{c1L,c1H},{c2L,c2H},{c3L,c3H}}, PlotStyle->{Darker[Green],Opacity[op]}, PlotPoints->pp ]
];
plots = Append[plots,plot]
];
For[i=e, i<=n-e, i++, ϕ = N[π i/n]; 
If[ useRays,
plot = ParametricPlot3D[ {
ArcCos[(b Sin[θ]+a Sin[ϕ] (Cos[CC] Cos[θ+ϕ]-Cos[ψ] Sin[CC] Sin[θ+ϕ])) /Sqrt[((b^2) ((Sin[θ])^2)+(a^2) ((Sin[ϕ])^2)+2 a b Sin[θ] Sin[ϕ] (Cos[CC] Cos[θ+ϕ]-Cos[ψ] Sin[CC] Sin[θ+ϕ]))]],
ArcCos[Cos[ϕ]], 
ArcCos[(b Cos[ϕ] Sin[θ]+a (Cos[CC] Cos[θ]-Cos[ψ] Sin[CC] Sin[θ]) Sin[ϕ])/Sqrt[(b^2) ((Sin[θ])^2)+(a^2) ((Sin[ϕ])^2)+2 a b Sin[θ] Sin[ϕ] (Cos[CC] Cos[θ+ϕ]-Cos[ψ] Sin[CC] Sin[θ+ϕ])]]
}, {ψ,-N[π]+.0001,N[π]-.0001}, {θ,.0001,N[π]-ϕ-.0001},PlotRange->{{c1L,c1H},{c2L,c2H},{c3L,c3H}}, PlotStyle->{Blue,Opacity[op]}, PlotPoints->pp ] , 
 plot = ParametricPlot3D[ {
ArcCos[(b Sin[θ]+a Sin[ϕ] (Cos[CC] Cos[θ+ϕ]-Cos[ψ] Sin[CC] Sin[θ+ϕ])) /Sqrt[((b^2) ((Sin[θ])^2)+(a^2) ((Sin[ϕ])^2)+2 a b Sin[θ] Sin[ϕ] (Cos[CC] Cos[θ+ϕ]-Cos[ψ] Sin[CC] Sin[θ+ϕ]))]],
ArcCos[Cos[ϕ]],
ArcCos[(b Cos[ϕ] Sin[θ]+a (Cos[CC] Cos[θ]-Cos[ψ] Sin[CC] Sin[θ]) Sin[ϕ])/Sqrt[(b^2) ((Sin[θ])^2)+(a^2) ((Sin[ϕ])^2)+2 a b Sin[θ] Sin[ϕ] (Cos[CC] Cos[θ+ϕ]-Cos[ψ] Sin[CC] Sin[θ+ϕ])]]
}, {ψ,-N[π]+.0001,N[π]-.0001}, {θ,.0001,N[2π]-.0001}, PlotRange->{{c1L,c1H},{c2L,c2H},{c3L,c3H}}, PlotStyle->{Blue,Opacity[op]}, PlotPoints->pp ]
];
plots = Append[plots,plot]
];
Show[ plots, ImageSize->800,Axes->False,Boxed->False,Frame->False,MeshStyle->None,BoxRatios->{1,1,1},PlotLabel-> Style[ 
(* "phis: " <> ToString[SetPrecision[ϕ1/π,3]] <> "π " <> ToString[SetPrecision[ϕ2/π,3]] <> "π " <> ToString[SetPrecision[ϕ3/π,3]] <> "π\n" <> *)
"side lengths: " <> ToString[SetPrecision[a,3]] <> "  " <> ToString[SetPrecision[b,3]]<> "  " <> ToString[SetPrecision[c,3]] <> "\ninterior angles: "<>ToString[SetPrecision[AA/π,3]] <> "π  " <> ToString[SetPrecision[BB/π,3]]<> "π  " <> ToString[SetPrecision[CC/π,3]]<>"π" <> "\n" <> type  , 20 ]] 
] , {{b,2.7},.1,5},{{c,3},.1,5},{{AA,.7π},0,N[π]},
{{c1L,0},0,π},{{c1H,π},0,π}, {{c2L,0},0,π},{{c2H,π},0,π}, 
{{c3L,0},0,π},{{c3H,π},0,π},{{op,.5},0,1},{n,40}, {pp,50},
{useRays,{True,False}},{useEnds,{True,False}}, ControlPlacement->Right ]