# How to find intersections of a curve and tetrahedra?

I have a curve going through tetrahedron elements. How to find the intersections of this curve with the tetrahedrons?

The curve is constructed by a few points and represents the axis of the artery going through tissue and I have tissue elements which are tetrahedron elements. I need to find those tetrahedron elements and the intersection points. I attached the image of my problem below. Pink elements are the ones that I have found manually and the blue elements show the elements constructing the tissue. The curve is white and the yellow colour shows an artery with the axis.

• Welcome to Math.SE. To make it easier for us to see the image you may put an exclamation mark (!) before the initial square bracket. – Ertxiem Apr 4 at 22:44
• I do not have 10 reputation score to do that. I get the required reputation score. – Elias Apr 4 at 22:59
• Use 3D barycentric coordinates en.wikipedia.org/wiki/Barycentric_coordinate_system (but you have first to understand how they work in 2D). Besides, couldn't your issue be treated in 2D by projecting on the 3 cordinate planes (xOy, yOz, zOx) ? – Jean Marie Apr 4 at 23:03
• I really liked the idea of simplifying it to 2D on each coordinate planes. It seems it can be treated that way. I need to find the most efficient way to do this. I have more than 200000 tets. – Elias Apr 4 at 23:18
• @JeanMarie : Could you please give me a hint how to make use of barycentric coordinate system? Does it have an advantage over the other approaches? – Elias Apr 5 at 1:29

A brute force approach:

Your problem isn't that easy and will require some sophistication to achieve efficiency. Anyway, a first solution can obtained by brute-forcing.

In any case, you have to discretize the curve as a contiguous sequence of line segments, i.e. a polyline, otherwise the equations are too complicated.

The brute-force approach consists in:

• for every tetrahedron:

• for every segment:

• determine if the segment crosses the tetrahedron.

You can also perform the search with the two loops swapped.

For the intersection operation, preprocess the tetrahedra to get the equations of the planes of the faces, oriented in such a way that the internal points make the equation positive, $$ax+by+cz+d\ge0$$. Now take the line of support of the segment and intersect the four faces. Use the parametric equation of the segment, $$p=p_0+t\cdot p_0p_1$$, and solve for $$t$$. You will get four intersections but at most two of them are valid: check which intersections belong to the inside of the tetrahedron by checking the sign when plugged in the equations of the three other planes.

Finally, consider the range of $$t$$ values that cross the tetrahedron, and take the intersection with the interval $$[0,1]$$ to restrict to the line segment. Notice that several segments may intersect a given tetrahedron.

To accelerate the process, several techniques can be used.

• First it is advisable to compute the axis-aligned bounding boxes of the tetrahedra and also of the segments. This will allow quick rejection of non-intersecting elements.

• Much more efficiency can be gained if you arrange the bounding boxes in a hierarchy of bounding boxes. When a segment does not interfere with a high-level box, all lower-level tests are avoided.

• If your tetrahedrization allows it, you can also work with the topology of the mesh: when you exit a tetrahedron, by a face, you perforce enter the neighboring tetrahedron that shares that face. Then by following the polyline, you will enter and leave tetrahedra at the same time that you follow the trajectory. This way, you only visit the relevant tetrahedra. But degeneracies are t be feared (segment meeting an edge or a vertex).

• Nice techniques! simple and efficient. Thanks! – Elias Apr 4 at 23:50
• @Elias: there is much more to be said about this problem. I hope you can get started with my hints. – Yves Daoust Apr 4 at 23:53
• I am going to start with the last technique. I have the first element I find the point that the curve leaves the tet and the neighboring tet will be the the element that I want. I implement this and hopefully I can do what I want. – Elias Apr 5 at 0:04

I'll start by some notation.

Let $$P$$ be a point in the curve.

The tetrahedron has faces numbered by $$k$$, with $$k \in \{1,2,3,4\}$$. The vector normal to each face can be denoted by $$\vec{n_k}$$. Choose them all pointing outward (or all pointing inward) to make our life easier. Let $$C_k$$ be the centre of face $$k$$.

You can construct 4 vectors $$\vec{u_k} = \vec{PC_k} = C_k - P$$, starting at point $$P$$ and ending at each center $$C_k$$. If $$P$$ is interior to the tetrahedron, then the dot product $$\vec{u_k} \cdot \vec{n_k} > 0$$ for all $$k$$, meaning that the vectors $$\vec{u_k}$$ and $$\vec{n_k}$$ point in (roughly) the same direction. A negative dot product would mean that that pair of vectors points in (roughly) opposite directions, so the point would be on the outside of that face.

By the way, to make the computations faster, $$C_k$$ can be any vertex that belongs to that face instead of the centre.

To find the interception point, we will need to look at the point in the curve after $$P$$. Let $$Q$$ be the point immediately after $$P$$. There is an interception if $$P$$ is interior to a tetrahedron and $$Q$$ is exterior to the same tetrahedron, or the other way around.

The interception is at $$R = P + x \vec{PQ}$$, with $$x \in (0;1)$$ such that $$\vec{RC_k} \cdot n_k = 0$$ for one of the faces. Or, equivalently, if $$x \in (0;1)$$, with $$x = \frac{\vec{PC_k} \cdot \vec{n_k}}{\vec{PQ} \cdot \vec{n_k}}$$

My approach has the limitation of assuming that each point in the curve will be inside of one tetrahedron.

• I am not the downvoter, but I am sorry, you don't answer the question : what you say hasn't any operational character. – Jean Marie Apr 4 at 23:07
To do so, I find barrycentric coordinates for point P and point Q which makes it easier to find if a point is inside the element or not. If all components are positive then the point is inside otherwise outside. If for example 2nd component is negative then the point Q is outside with respect to face $$(V_1V_3V_4)$$ where Vs are vertices of the element. Then I use parametric equation of the segment to find where the segment crossed the element as in Yves' answer. Also, I can find the face it crosses and based on that I find the adjacent element and do the same thing to get the next intersection and element. There are some special cases like when the intersection is close to the edge or vertex that should be considered though.