Calculating ray-cylinder intersection points How is it done? I couldn't find a general formula/algorithm for this problem. I've read this article, I've managed to calculate the intersection points around the $z$-axis like they explained, but I didn't understand how to transform, rotate, and scale such that I can calculate the intersection points in an arbitrary location.
I'm trying to build a simple 3D engine, so detailed explanation/general formula would be great.


*

*My ray is defined with: origin point, normalized direction vector.

*My (infinite) cylinder is defined with: ray, radius.
If it's possible to extend the answers also for a finite cylinder-ray intersection, that would be more than great. Thanks in advance!

P.S. I've read many articles and comments (in Stack Overflow and elsewhere). Like I said, couldn't find a helpful solution for my problem.
P.P.S. I know how to calculate ray-triangle and ray-sphere intersections with arbitrary locations, if that helps.
 A: Consider the shortest line segment between the cylinder axis and the ray and form a reference frame such that the cylindre axis is $z$ and the line segment is $x$. This frame is built by using the cross product of the two direction vectors (assumed to be unit), $\vec x=\vec c\times\vec r$, then $\vec z=\vec c,\vec y=\vec z\times\vec x$. Normalize these vectors. Then project the vector that joins a point of both lines onto $x$ to get the shortest segment. This projections is the vector $\vec d=(\vec{o_co_r}\cdot \vec x)\,\vec x$.

The origin point $\vec o=\vec o_c+s\,\vec c=\vec o_r+t\,\vec r-\vec d$ is found by solving
$$\begin{cases}
(\vec o_c+s\,\vec c)\cdot\vec y=(\vec o_r+t\,\vec r)\cdot\vec y,\\
(\vec o_c+s\,\vec c)\cdot\vec z=(\vec o_r+t\,\vec r)\cdot\vec z
\end{cases}$$ for $s$ or $t$.
In the new frame, the cylindre has the implicit equation
$$x^2+y^2=r^2$$
and the ray the parametric equations
$$\begin{cases}x=d,\\y=\beta t,\\z=\gamma t\end{cases}$$ where $d=\vec d\cdot\vec x,\beta=\vec r\cdot \vec y,\gamma=\vec r\cdot\vec z$.
Now
$$d^2+\beta^2t^2=r^2$$ gives the two intersections of the infinite cylindre with the infinite ray. You can restrict to $t\ge0$ for a half-ray.
If the cylindre has finite extent, the two basis will have the equations $z=z_{min}$ and $z=z_{max}$ (wrt to the origin $o$), giving the intersections by
$$z=\gamma t.$$
When you have the range of $t$ inside the infinite cylindre, and the range of $t$ between the two basis, it is an easy matter to find the common range. Then you can compute the two intersection points in the auxiliary frame and back-transform to the initial frame. Notice that the matrix of the back transform is the transpose of the direct transform, and that you invert
$$\vec q=R(\vec p-\vec o)$$ by
$$\vec p=R^{-1}\vec q+\vec o.$$
A: The points at which the ray intersects the cylinder are the only ones on the line that are at a distance equal to the radius from the cylinder’s axis. Since you’re starting from a description of the cylinder as axis and radius, you can use a standard formula for the distance from a point to a line to find these points instead of trying to come up with an equation for the cylinder or trying to come up with a transformation into some standard configuration.  
Let $\mathbf x(t) = \mathbf p_0+t \mathbf v$ be the parameterization of the ray with the given starting point and direction vector. Choose two points $\mathbf x_1$ and $\mathbf x_2$ on the cylinder’s axis: since that’s also defined by a ray (line?) you can choose the origin point of that line for $\mathbf x_1$ and add any convenient multiple of the axis direction vector to it for the other. Letting $r$ be the cylinder’s radius, the point-line distance formula gives following the quadratic equation in $t$: $${\lvert (\mathbf x(t)-\mathbf x_1)\times(\mathbf x(t)-\mathbf x_2)\rvert^2 \over \lvert \mathbf x_1-\mathbf x_2 \rvert^2}=r^2.$$ Expand and solve for $t$, rejecting any negative solutions, then compute $\mathbf x(t)$ for each resulting value of $t$. The one with the lesser $t$-value is the nearer to the origin of the ray.  
For a finite cylinder, you can then project these points onto the cylinder’s axis and perform a range check. If you choose for $\mathbf x_1$ and $\mathbf x_2$ above the two points on the cylinder’s axis that bound the cylinder, then if $\mathbf p$ is a solution to the infinite intersection, it lies on the bounded cylinder iff $0\le(\mathbf p-\mathbf x_1)\cdot(\mathbf x_2-\mathbf x_1)\le(\mathbf x_2-\mathbf x_1)\cdot(\mathbf x_2-\mathbf x_1)$.
