# How is ST_PointOnSurface calculated?

The PostGIS documentation states that `ST_PointOnSurface` returns "a POINT guaranteed to lie on the surface". It seems like this function could be trivially implemented to give results that satisfy the documentation but provide little real-world utility, though I'm certain that PostGIS provides a non-trivial implementation.

This introduction to PostGIS provides a nice comparison and contrast of `ST_Centroid` with `ST_PointOnSurface` and says that "[ST_PointOnSurface] is substantially more computationally expensive than the centroid operation".

Is there a more thorough explanation of how `ST_PointOnSurface` is calculated? I've been using `ST_Centroid`, but have encountered some edge cases in my data where the centroid is outside the geometry. I believe that `ST_PointOnSurface` is the correct substitute, but the function name and documentation leave room for uncertainty.

Further, is the computational expense of `ST_PointOnSurface` incurred even if the centroid does lie within the geometry already?

• It exists precisely since the centroid of non-convex polygons is not always included in it. It has nothing to do with heights and DEMs if that's the confusing part of the name. Implementation details are best checked in the code, but I believe you'd get a better answer on GIS.se. Commented Nov 4, 2013 at 18:16
• Good point on GIS.se. Is there a way to migrate this question there? I understand why both functions exist. I find the name confusing because there are infinitely many points on the surface of the polygon goemetries I'm working with. However, only a small subset of those points serves my purpose. I want to know I'm getting a point that makes sense for how I want to use it.
– mjobrien
Commented Nov 4, 2013 at 20:16
• I'm curious about cases like an (rough) oval or a rectangle, where both ST_PointOnSurface and ST_Centroid would return a point within the object, albeit at slightly different locations. In this case is there any value left to using ST_PointOnSurface, given that its big selling point no longer is a factor? Commented May 21 at 3:49

Based on a few experiments, I think `ST_PointOnSurface()` works roughly like this, if the geometry is a polygon:

1. Trace an east-west ray, lying half-way between the northern and southern extents of the polygon.
2. Find the longest segment of the ray that intersects the polygon.
3. Return the point that is half-way along said segment.

That may not make sense, so here's a sketch of a polygon with a ray dividing it into a northern and southern parts:

``````             _
/ \             <-- northern extent
/   \
/     \
/       \
/         \      __
/           \    /  \
/_ _ _ P _ _ _\  / _ _\  P = point-on-surface
/               \/      \
/                         \
/            C              \   C = centroid
/                             \
/                              /
/______________________________/  <-- southern extent
``````

Thus, `ST_PointOnSurface()` and `ST_Centroid()` are usually different points, even on convex polygons.

The only reason for the "surface" in the name, I think, is that if the geometry has 3D lines then the result will be simply be one of the vertexes.

I would agree that more explanation (and better naming) would have been useful and hope a GEOS programmer might shed some more light on the matter.

• Looking at the libgeos code, I believe you are right. The horizontal bisector is found, then the midpoint of the widest intersection is used. Commented Nov 6, 2013 at 6:41
• Does this mean that this call always produce the same point always? As in, same if called again with the same geometry? Commented Oct 21, 2020 at 11:11
• @iamthadiyan - I would assume so. Commented Jan 9, 2021 at 0:30
• I think the orginator of this function was a French since the word `surface` in the name of the function in French means "Area" Commented Aug 22, 2023 at 6:23

The Algorithm mentioned by Martin F is not guaranteed to work for every polygon. I am not sure if this is the algorithm used by said function or not, however consider the example below:

In this example, the longest line segment that the bisector crosses is between the points (-0.147,2.5) and (3.207,2.5), and the midpoint of this segment is the point (1.53,2.5), which is certainly not on the surface of the polygon.

Edit: I do realise this should've been a reply to Mark F's answer, however I didn't know how to include or reference an image in a reply.