[Math] Calculating ray-cylinder intersection points

3dcomputational geometry

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.

Best Answer

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$.

enter image description here

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.$$