Lines rotation - problem Given $3$ infinite lines $L_1$, $L_2$, $L_3$ in $R^3$ ($3D$ space), each line is represented by two $3D$ points, the line $L_1$ is rotating around the axis line $L_3$ until it intersects with the line $L_2$ (or in other words, until there's a plane that $L_1$ and $L_2$ span).
Rotation of a line around axis line means that every point on the line is always in the same distance from the axis line.

I need to find the angle that $L_1$ needs to rotate (around $L_3$) in which it intersects with $L_2$ (or none if the lines don't intersect).
  Anyone knows how to find it? I've been trying to solve it for quite a while.

A few things I already tried:
The surface of revolution created by rotating the line $L_1$ around $L_3$ is hyperboloid (rotated in space). I was trying to find the equation of the hyperboloid (given $L_1$ and $L_3$) and then to find the intersection point(s) between the hyperboloid and the line $L_2$. But I didn't succeed to find the hyperboloid equation.
Also, I wrote code (in matlab) that rotates a line around another line. so I can have the equation of $L_1$ rotated by the angle theta around $L_3$ (equation depending by one variable - theta ($\theta$)). and I need to find for which $\theta$ there's a point intersection between $L_1$ and $L_2$. The problem is that the equation is extremely long so I couldn't really continue and find the theta in which there's an intersection point.
I'd be happy for any suggestion that can help.
Thank you very much
David
 A: Sketch of a solution: Call the rotation $R$. The point $P$ at which $RL_1$ and $L_2$ intersect (if at all), is the image $RQ$ of a point $Q$ on $L_1$ that is the same distance from the axis of rotation $L_3$ as $P$. This problem is thus equivalent to finding $P\in L_2$ and $Q\in L_1$ that are equidistant from $L_3$ and lie in a plane perpendicular to it. One way to find these two points is to parameterize $L_1$ and $L_2$ and then use well-known formulas for the distance between a point and line in 3-D. Once you have these two points, computing the angle formed by them and the rotation axis is straightforward.  
Update: A better approach to finding $P$ and $Q$ is to parameterize $L_3$ as $S(\lambda)=(1-\lambda)A+\lambda B$, where $A$ and $B$ are any two convenient points on this line. You can then find $P(\lambda)$ and $Q(\lambda)$, the intersections of $L_1$ and $L_2$, respectively, with the plane perpendicular to $L_3$ through $S$. Solve $\|P-S\|=\|Q-S\|$, which you can square to get a simple quadratic equation, for $\lambda$, and you have the three points you need to compute the required rotation angle. There may, of course be zero, one or two solutions.

Update (expanding on my comment below): Once you have the points $P$, $Q$ and $S$, an appropriate transformation matrix can be constructed without explicitly computing any angles.  
Let $\mathbf u$ be a unit direction vector for $L_3$ ($B-A$ normalized will do), $\mathbf v=(P-S)/\|P-S\|$ and $\mathbf w=(Q-S)/\|Q-S\|$. Since you haven’t specified an orientation for $L_3$, it’s not clear what you mean by a counterclockwise rotation, but let’s say it means that you’d like the rotation angle to be in the range $[0,\pi]$. Examine the sign of the scalar triple product $\mathbf u\cdot\mathbf v\times\mathbf w=\det\begin{bmatrix}\mathbf u&\mathbf v&\mathbf w\end{bmatrix}$. If this is negative, then reverse $\mathbf u$. If it’s zero, then either $L_1$ and $L_2$ already intersect at $P$ or $\mathbf w=-\mathbf v$, in which case the rotation angle will be $\pi$ and there’s no need to make any adjustments. In fact, there’s no need to make any adjustment at all: the orientation changes cancel when the change-of-basis matrices are combined below, so you end up with the same transformation matrix either way.  
Now, both $\mathscr B_1=(\mathbf u,\mathbf v,\mathbf u\times\mathbf v)$ and $\mathscr B_2=(\mathbf u,\mathbf w,\mathbf u\times\mathbf w)$ are right-handed orthonormal bases for $\mathbb R^3$. This means that the matrices $$\begin{align}R_1&=\begin{bmatrix}\mathbf u&\mathbf v&\mathbf u\times\mathbf v\end{bmatrix} \\ R_2&=\begin{bmatrix}\mathbf u&\mathbf w&\mathbf u\times\mathbf w\end{bmatrix}\end{align}$$ represent rotations that map the standard basis onto $\mathscr B_1$ and $\mathscr B_2$, respectively. Therefore, the homogeneous transformation matrix that represents the rotation you’re looking for is $$R=\left[\begin{array}{c|c}R_2&S\\\hline0&1\end{array}\right]\left[\begin{array}{c|c}R_1&S\\\hline0&1\end{array}\right]^{-1}.$$ You can avoid the matrix inversion by taking advantage of block structure and the fact that $R_1$ is orthogonal: $$\left[\begin{array}{c|c}R_1&S\\\hline0&1\end{array}\right]^{-1}=\left[\begin{array}{c|c}R_1^{-1}&-R_1^{-1}S\\\hline0&1\end{array}\right]=\left[\begin{array}{c|c}R_1^T&-R_1^TS\\\hline0&1\end{array}\right]$$ and multiplying it all out yields $$R=\left[\begin{array}{c|c}R_2R_1^T&(I-R_2R_1^T)S\\\hline0&1\end{array}\right].$$ Note that you may use any convenient point on the line instead of $S$ for this part of the construction.
A: I like to use Plucker coordinates for lines. 


*

*A line through two points A and B has coordinates $$ \ell_{AB} = \begin{pmatrix} {\bf r}_B - {\bf r}_A \\ {\bf r}_B \times {\bf r}_A \end{pmatrix} $$

*A line through a point A along a direction z has coordinates $$ \ell_{Az} = \begin{pmatrix} {\bf z} \\ {\bf r}_A \times {\bf z} \end{pmatrix} $$

*Two lines $\ell_1 = ({\bf e}_1, \;{\bf m}_1)$ and $\ell_2 = ({\bf e}_2, \;{\bf m}_2)$ have minimum distance between them $$ d = \frac{ | {\bf e}_1 \cdot {\bf m}_2 +{\bf e}_2 \cdot {\bf m}_1 |}{\| {\bf e}_1 \times {\bf e}_2 \|} $$

*A line $\ell_1 = ({\bf e}_1, \;{\bf m}_1)$ is rotated about a point A about the z direction has coordinates $$ \ell_1' = \begin{pmatrix} {\bf e}_1' \\ {\bf m}_1' \end{pmatrix} = \begin{pmatrix} {\rm R}_z {\bf e}_1 \\ \left({\bf r}_A - {\rm R}_z {\bf r}_A\right) \times {\rm R}_z {\bf e}_1+{\rm R}_z{\bf m}_1 \end{pmatrix}$$  where ${\rm R}_z$ is the 3×3 rotation matrix about axis z with angle $\theta$. Can be defined with the Rodrigues' rotation formula. $$\begin{aligned} {\rm R}_z & = 1+\sin \theta {\bf K} + \left(1-\cos \theta\right) {\bf K}^2   \end{aligned}$$ and K the Cross Product 3×3 Matrix of the rotation axis.
Your task is to find the angle $\theta$ that satisfies coincidence (from the distance formula)
$$ {\bf e}_1' \cdot {\bf m}_2 + {\bf m}_1' \cdot {\bf e}_2 =0 $$
The solution of the above is equivalent to a 4-th order polynomial (use $\theta = 2 \tan^{-1}(t)$ substitution, or $\cos \theta = \frac{1-t^2}{1+t^2}$ and $\sin \theta = \frac{2 t}{1+t^2})$
A: I finaly managed to write code that solves the problem. 
I didn't write the code for all the possibilites, but it works for most of them.
Here's a sketch of the solution (My solution is based on @amd suggestion).


*

*part1 - rotate all lines (L1, L2, L3) so that line L3 (axis line) will be the same as z axis.

*part2 - change representation of lines Li (i=1 or 2) so that the point of Li is the intersection between the line and Z=0. and the vector is between the points intersection of Z=0 and Z=1 (special cases are for lines with constant z. I didn't take care of those cases).
that way we can use the same 't' in the following lines representations:
L1: (p1x,p1y,0)+t*(v1x,v1y,1)
L2: (p2x,p2y,0)+t*(v2x,v2y,1)
for every value of t in both lines we get two points in the same height (height is the z value).
so we want to find the t values in which the distances from the z axis are the same. for these t values we get points that will intersect.
we can get 0,1 or 2 solutions here.

*part 3 - calculate the angle between the points calculated before (in the same z level).
Also, since the rotation is always counter clockwise I sometimes needed to take the complementary angle to 360 degrees.
For anyone who wants to try the code (matlab), I published it here (including simulation of the rotations): 
https://github.com/David1324/lines_rotations
A: EDIT1:
Interesting problem. I am not sure if I am making things more complicated than necessary in the hint. 
Let $L_4$ be perpendicular to $L_3$. Generate (find equations of) two one-sheeted hyperbolas swept around their axes of symmetry $L_3,L_4$ made by sweeping straight generators $L_1,L_2$. The intersection of these surfaces produces 3D hyperbolic space curves.
Out of the intersection we can pick up required point of the $L2$ generator at a particular rotation.
wlog $L_3,L_4$ can be the $ z-,y-$ axes. 

Others may solve it perhaps more directly using quaternions.
