1

I have two Geometry shapes of types Polygon or MultiPolygon [with Lat & Lng points] that I want to get their intersected shape. I do this as follows:

land_shape.STIntersection(grid_shape) as shape

Then, I check to see when

shape.STIsValid() = 1 and geography::STGeomFromText(shape.STAsText(), 4326).STIsValid() <> 1

and I get many results back. There are instances where the intersected shape as a geometry is valid, but as a geography isn't. This is an issue because I build on the geography shape to get the area.

I tested to see if converting both inputs (land_shape & grid_shape) to geography would also be invalid, but they are not invalid. So, inputs are not corrupt.

I'm not sure what could cause the intersected shape to be invalid as a geography. I tried using IsValidDetailed() and the reasons were:

24409: Not valid because some portion of polygon ring (1) lies in the interior of a polygon.
24413: Not valid because of two overlapping edges in curve (1).

But, why were they valid as a geometry then? I know that MakeValid() can fix it, but I just need to understand why that behavior happens to begin with, and if I can eliminate this issue without using SQL Server's MakeValid() as I read it isn't entirely accurate.

Update: Even after applying MakeValid() I still get very large areas returned How I calculate the area is as follows:

geography :: STGeomFromText(land_shape.STAsText(), 4326).STIntersection( geography :: STGeomFromText(grid_shape.STAsText(), 4326) ).STArea()

I thought initially that the problem was the Geography shapes being invalid, but looks like it might not just be that.

I tried applying ReorientObject() but still same issue

geography :: STGeomFromText(land_shape.STAsText(), 4326).ReorientObject().STIntersection( geography :: STGeomFromText(grid_shape.STAsText(), 4326).ReorientObject() ).STArea()
HR1
  • 487
  • 4
  • 14
  • Can't you first convert them to geography and then do the intersections? – siggemannen Apr 05 '23 at 08:41
  • @siggemannen I already do this || geography :: STGeomFromText(land_shape.STAsText(), 4326).STIntersection( geography :: STGeomFromText(grid_shape.STAsText(), 4326) ).STArea() || But I was still getting very large areas, so checking if the shape geography that calculated the area was valid was just an approach to look into the problem – HR1 Apr 05 '23 at 11:13
  • Strange. Could you paste your areas somewhere? Or it is the shape of the secret area51 compounds / uncovered gold mines / treasures – siggemannen Apr 05 '23 at 13:03

1 Answers1

1

Saving Geometry as WKT and reading it as Geography does not properly converts from geometry to geography. The two types have different interpretations of that WKT, even though they both use similar WKT and they use same coordinates. Edges between two vertices in Geometry are the straight lines on the flat map. Edges between the same two vertices in Geography follows the geodesic line - so it is a very different line.

Consider e.g. two loops:

  'linestring(30 40, 80 40, 80 60, 30 60, 30 40)'

and

  'linestring(50 42, 60 42, 60 50, 50 50, 50 42)'

On a flat map the second loop is strictly inside the first one - so you can make a polygon with first loop as shell, and second loop as hole - and it will be a valid polygon. But if you use the same WKT and try to interpret it as Geography - the second loop intersects the first one, and if you construct a polygon the same way - it is no longer valid.

enter image description here

Thus a polygon

 'polygon((30 40, 80 40, 80 60, 30 60, 30 40),(50 42, 60 42, 60 50, 50 50, 50 42))'

is valid as Geometry, but not valid as Geography. There are also issues around anti-meridian and poles. The proper conversion means approximating "flat" edges of Geometry with geodesic edges of Geography and vice versa - usually involving densification of the edge, i.e. adding some points in the middle. E.g. this is what Google BigQuery does in this case:

select st_geogfromtext(
    'polygon((30 40, 80 40, 80 60, 30 60, 30 40),(50 42, 60 42, 60 50, 50 50, 50 42))', 
    planar => TRUE)

note the planar argument, telling BigQuery the input comes from a flat geometry. The result is pretty large, but it follows the original "flat" shape with about 10 meter precision:

POLYGON((30 40, 30.1953125 40, 30.390625 40, 30.5859375 40, 30.78125 40, 30.9765625 40, 31.171875 40, 31.3671875 40, 31.5625 40, 31.7578125 40, 31.953125 40, 32.1484375 40, 32.34375 40, 32.5390625 40, 32.734375 40, 32.9296875 40, 33.125 40, 33.3203125 40, ....

One more issue is ring orientation - which does not matter for Geometry, but is very important for Geography, see https://learn.microsoft.com/en-us/sql/relational-databases/spatial/polygon?view=sql-server-ver16#orientation-of-spatial-data

Michael Entin
  • 7,189
  • 3
  • 21
  • 26
  • But how would you suggest I accurately get the intersected area 'STArea()' for the shapes I've converted to Geography but were not valid using sqlserver? – HR1 Apr 05 '23 at 06:04
  • I realized that using MakeValid() doesn't solve this issue because I still get very large areas for shapes I know don't intersect – HR1 Apr 05 '23 at 06:05
  • I know and gave the answer to the original question - what is causing it. I don't know the answer of how to solve it within SqlServer, sorry. You might need 3rd party tools. But you mentioning large area reminded me there is another issue - ring orientation, see the last paragraph in updated answer. The incorrect orientation is likely the source of huge area (near the area of whole globe). – Michael Entin Apr 05 '23 at 06:14
  • One suggestion: instead of intersecting Geometries and then converting to Geography, convert to Geography first, then intersect Geographies. It might help with some shapes, but is not guaranteed to work in all cases. – Michael Entin Apr 05 '23 at 06:18
  • I also tried loading this into ArcMap exluding the geography intersection shapes that are not valid, and I still get "The number of points is less than required for feature." So I'm guessing the issue you described will not necessarily be detected by STIsValid(). Not sure how to else detect it – HR1 Apr 05 '23 at 06:45
  • What I mean is there are still valid geographies that yield very large areas – HR1 Apr 05 '23 at 06:48
  • Valid Geographies with very large areas are sign of "inverted" polygons due to wrong loop orientation (see last link in my answer). Reorienting, when the result is still valid, usually gives correct area. – Michael Entin Apr 05 '23 at 07:00
  • SQL server offers a ReorientObject() function, I tried using it and I still get same big area. – HR1 Apr 05 '23 at 11:14