# 3D Points within 3D Polygon in PostGIS

I've been working with PostGIS for a while now, and I'm able to run queries over massive datasets reasonably easy and fast. However, I do have an issue now that I don't know how to solve. I did find some workarounds, although they mean lower performance, higher execution times, and more complex queries.

The objective is straightforward: to get all the points within a volume. But, the volume is not a simple extrusion of a polygon, it is instead of the difference of a higher volume minus some smaller volumes here and there.

My first attempt: to tackle the problem was inspired by this question, although it didn't work:

1. Create the master polygon, using a third axis (Z) as the base for the volume:

ST_Polygon(ST_GeomFromText('LINESTRING(Longs, Lats, Zs...)'), 4326);

1. Extrude the polygon using the height of the volume, getting a POLYHEDRALSURFACE as a resulting geometry:

ST_Extrude(Previous_polygon, x=0, y=0, z=height);

1. Having a table with the information regarding each master volume and its small volumes to be removed from. I created a new geometry removing all these small volumes from the master volume:

ST_3DDifference(Master volume, combination of smaller volumes to be removed);

1. Query the table with the positions and joining the resulting master volume on

3dpoints.geometry &&& Master_Volume.geometry.

As I said, it didn't work because the points were also found outside the master volume. Not sure why! I guess I'll try to come back to that option should the following solution be too much of a hassle.

The second attempt I tried was:

1. Query the table with the positions straightaway, merging the volume as a 2D polygon, and including the filter on the Z-axis, something like:

ST_Contains(Master_Polygon.geometry, 2dpoints.geometry) AND 2dpoints.Z BETWEEN Master_Polygon.MinZ AND Master_Polygon.MaxZ

1. Of course, this included all the point within the smaller volumes I didn't want to consider, so then I added the condition to remove those points:

ST_Contains(Master_Polygon.geometry, 2dpoints.geometry) AND 2dpoints.Z BETWEEN Master_Polygon.MinZ AND Master_Polygon.MaxZ AND NOT (ST_Contains(Smaller_Polygon.geometry, 2dpoints.geometry) AND 2dpoints.Z BETWEEN Smaller_Polygon.MinZ AND Smaller_Polygon.MaxZ)

And it does work, although only when the smaller volume (polygon) is only one. Should the lower volumes be more than one, then the condition is not enough as it would create a different set of points based on each condition. This gives all the points within the master volume but not in the lower volume A would be returned or all the points within the master volume but not in the lower volume B.

Question: Is there a way to get all the points in the master volume by removing all those points with a condition which are in any of the smaller volumes?

I apologize for the long post, but its due to the complicated explanation of the issue. I did find a solution, and it's by splitting the master volume into smaller volumes and merging it after the query. But it is not quite right, especially when I will have to update the geometries and make the master volume to change accordingly automatically.

• a very naive option would be to count the geometry occurances and select only those that were foumd more than once... Commented May 18, 2018 at 14:48
• I have done similar things. I built a line of sight application recently, with extruded trees and buildings. The 3D functions were horrifically slow for the amount of data I was processing, so I went back to 2D to get candidates for intersection, before doing actual 3D calcs on the dumped polyhedral surfaces from the 3D volumes. Commented May 18, 2018 at 14:51
• any progress yet? Commented May 22, 2018 at 14:13

I played around with some test setups a while ago; it seems ST_3DUnion, ST_3DDifference and ST_3DIntersection work on Volumes; you could try sth. like this:

WITH
volume AS (
SELECT ST_3DDifference(
<master_volume>,
ST_3DUnion(<smaller_volumes>)
) AS geom
)
SELECT ST_3DIntersection(vol.geom, pts.geom) AS geom
FROM volume AS vol,
<points> AS pts

Since the spatial relation filters won´t work, a selective JOIN is impossible; thus this will take damn long and returns GEOMETRYCOLLECTION EMPTY for all non-intersecting points that you then need to filter out.

Just wanted to share; I´m not sure if this is actually going to work, maybe someone would know exactly if and why...

It seems to return correct results with simple test setups like

SELECT ST_AsText(
ST_3DIntersection(
ST_3DDifference(
ST_MakeSolid(ST_Extrude(ST_Buffer('POINTZ(0 0 0)'::geometry, 20), 0, 0, 20)),
ST_3DUnion(ST_MakeSolid(ST_Extrude(ST_Buffer('POINTZ(0 0 0)'::geometry, 1), 0, 0, 10)),
ST_MakeSolid(ST_Extrude(ST_Buffer(ST_GeomFromText('POINTZ(0 0 5)'), 1), 0, 0, 10)))
),
'POINTZ(0 0 5)'::geometry     --'POINTZ(2 0 5)'::geometry
)
) AS geom

but that could also just be completely random...

• Correct, they work on volumes. The ST_MakeSolid function is significantly underdocumented, imho, in all of this. Sometimes they work on polyhedral surfaces also, it kind of depends what you are expecting back. Commented May 18, 2018 at 14:47
• @JohnPowellakaBarça ah nice, thanks, I never made it to larger real life setups with my testing. do you have any more experience with them, or other functions that also work? some drawbacks (apart from being reeeaally slow...)? Commented May 18, 2018 at 14:52
• @JohnPowellakaBarça could you speed up the cross join with a &&& filter in WHERE? Commented May 18, 2018 at 14:54
• The speed is the main issue. I am still learning how it all fits together, but I have settled on a 2D/3D combo approach. Use 2D to get candidates, and then actual 3D only for those results. Commented May 18, 2018 at 14:55
• No idea about &&&, but seems possible. Test away :D. Commented May 18, 2018 at 14:55