Suppose I have a mesh of the sphere (points on the sphere and a triangulation).

What is a good and efficient way to find the triangle which contains a point on the sphere (the point does not need to be a triangulation node)? It would be nice if the operation could be vectorized, to improve Matlab speed.

The purpose of this is to construct a simple spherical (linear) interpolation operator on the sphere.

  • $\begingroup$ I don't see the difference between this problem and finding a point in a triangle in the plane. Can you elaborate? $\endgroup$ Commented Oct 9, 2015 at 11:36
  • $\begingroup$ I don't have Matlab right now, but have you tried using triangulation.pointlocation directly? Not sure if this works for 3D points but "2D" simplices (triangles). I believe triangulation.cartesiantobarycentric does work in this case (the 3D point is 1st projected onto the triangle's plane), so the pointlocation() function may also ... $\endgroup$
    – GeoMatt22
    Commented Oct 9, 2015 at 16:51
  • $\begingroup$ @WolfgangBangerth: I thought that the surface problem was more complicated than the plane one and needed other methods. $\endgroup$ Commented Oct 9, 2015 at 20:45
  • $\begingroup$ @GeoMatt22: I'll try to see if I can use the functionalities of the functions you mention. $\endgroup$ Commented Oct 9, 2015 at 20:45
  • 1
    $\begingroup$ The sphere problem is as easy as the 2d problem because projecting an arbitrary point onto the sphere is so simple. In essence, you only have to work in polar coordinates, which by and large are two dimensional if you ensure that no vertex lies at either pole. $\endgroup$ Commented Oct 11, 2015 at 1:47

2 Answers 2


A generic way of looking up the element in which a point of given coordinates may be is to sort them into a quadtree (octree in 3D). The leaves of the tree will contain only element having an intersection with a given range of coordinates (in your case this could be ranges in $(\theta, \phi)$), the number of elements per leaf is kept smaller than a small enough constant $n$ (around 10 to 20), and the test of whether the point belongs to an element is performed on all of these. This gives a cost of order $n + log (N/n)$ for an $N$-element triangulation.

If your question is about how to test optimally that the point belongs to a triangle, please state it explicitly.


In the case of the sphere there is a simple and computationally efficient method in the case the triangles are not far from being equilateral (thanks for the correction in the comments). Suppose we are working on the unit sphere, so the center of the sphere is at the origin. Then in order to find the triangle which contains a certain point $P$ we construct the exterior normals to each triangle and search the one which makes the minimal angle with the vector $\vec{OP}$. The triangle corresponding to the closest normal will contain the point.

  • $\begingroup$ Can you prove it? I believe this may be wrong if triangles are far enough from equilateral. $\endgroup$
    – Joce
    Commented Feb 8, 2016 at 16:24
  • $\begingroup$ Besides, the cost here is in $N$, number of elements in the mesh, so you should locate only a constant number of point to ensure that your algorithm has linear cost -- so unsuitable for interpolation, which requires $N$ look-ups. $\endgroup$
    – Joce
    Commented Feb 8, 2016 at 16:34
  • $\begingroup$ Yes, you are right. It seems that the triangles need to be almost equilateral. In what I needed, that was the case. $\endgroup$ Commented Feb 8, 2016 at 21:49

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.