3

I want to compute surface areas of a large CityGML data set. However, I ran into some issues, which I could not resolve or fully understand.

The error message is the following:

Polygon is invalid : points don't lie in the same plane

Minimal examples:

  1. SELECT st_3darea( ST_SetSrid(ST_GeometryFromText( 'POLYGON Z ((0 0 0, 0 1 0, 1 1 0, 1 0 0, 0 0 0))' ),4258)) /*<- This works fine */

  2. SELECT st_3darea( ST_SetSrid(ST_GeometryFromText( 'POLYGON Z ((0 0 0, 0 11, 1 1 0, 1 0 0, 0 0 0))' ),4258)) /*<- This does not work: Error Polygon is invalid : points don't lie in the same plane */

  3. SELECT st_3darea( ST_SetSrid(ST_GeometryFromText( 'POLYGON Z ((0 0 0, 0 11, 1 11, 1 0 0, 0 0 0))' ),4258)) /*<- This works fine again */

  4. SELECT st_3darea( ST_SetSrid(ST_GeometryFromText( 'POLYGON Z ((0 0 0, 0 1 1, 1 0 0, 0 0 0))' ),4258)) /*<- This works fine, removed one point */

From my understanding in 2), the function is not able to transform the coordinates into a flat plane to proceed with its calculations. Somehow, I am not sure if it is a bug or a feature (is the polygons in 2) invalid, why?)

As answered by@geozelot, its not a bug, but PostGIS does not support curved geometries.

My idea for a solution:

Is there a generic workflow/function to split polygons (e.g. 2)) in a form like in 4) triangles and then add the sums of the areas?

I would need it automated for the CityGML data set, which I do not know how to do e.g. non-automated

SELECT st_3darea( ST_SetSrid(ST_GeometryFromText( 'POLYGON Z ((0 0 0, 0 1 1, 1 0 0, 0 0 0))' ),4258)) + SELECT st_3darea( ST_SetSrid(ST_GeometryFromText( 'POLYGON Z ((0 1 1, 1 1 1, 1 0 0, 0 1 1))' ),4258))

1 Answer 1

3

The vertices given in 2) cannot be elements of the same 2-dimensional plane; i.e. vertex (0 1 1) cannot be contained in the plane defined by the other vertices. Simple POLYGON geometries (or faces), no matter the Z component, must not be curved.

PostGIS does have (limited) support for CURVED geometries that could be used for non-linear faces, e.g.

SELECT ST_AsText(ST_SetSrid(ST_GeometryFromText( 'CURVEPOLYGON(CIRCULARSTRING(0 0 0, 0 1 1, 1 1 0, 1 0 0, 0 0 0))' ),4258));

however, they won't solve your issue.


Note that ST_IsValid won't detect this as it does not consider z values; ST_IsPlanar, however, should detect this.

3
  • thank you for the insight and good to know that the support is limited. This is what i suspected in "From my understanding". In the real life data set, this occurs often. What could be a work around? I suggested to split the vertices to form triangles, but lack the experience to generalize it. Any thoughts on that? Commented Jun 9, 2021 at 10:39
  • 1
    @fx-9750GPLUS I don't have SFCGAL available at the moment, but you could try ST_Tesselate and dump the resulting TINs. ST_DelaunayTriangles(GEOS) or ST_ConstrainedDelaunayTriangles (SFCGAL) will create triangular Polyons, but have a deterministic order of vertex usage that cannot be influenced.
    – geozelot
    Commented Jun 9, 2021 at 11:06
  • I tried this on my problem.ST_Tesselate would not work with the same error: "Polygon is invalid : points don't lie in the same plane". For ST_DelaunayTriangles(GEOS) an empty geometerycollection is returned and ST_ConstrainedDelaunayTriangles (SFCGAL) returns a wrong and very small polygon. Commented Jul 7, 2021 at 8:34

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.