[gdal-dev] Clipping shapefile with another produces invalid shapes

Previous Topic Next Topic
 
classic Classic list List threaded Threaded
6 messages Options
Reply | Threaded
Open this post in threaded view
|

[gdal-dev] Clipping shapefile with another produces invalid shapes

Paul Meems
Hi list,

I have a very large shapefile, with almost 700k shapes which I want to clip with a border shapefile (in red):

I'm using this command:
ogr2ogr -clipsrc "border.shp" -overwrite -explodecollections -f "ESRI Shapefile" ogr_clipped.shp Fishnet.shp
with GDAL 2.1.3, released 2017/20/01
Eventually, I'll transform the command to a call to ogr.VectorTranslate().

The result is visually correct:

But at the red asterix, I have an invalid shape. The reason is according to GEOS:
Polygon must be clockwise
I have more of these errors.
They all seem to occur when the shape due to the clipping is triangular shaped:

Is this an error/warning which can be expected or do I need to add an additional parameter?
If needed I can try to make a test set with as little data as possible.

Thanks,

Paul

_______________________________________________
gdal-dev mailing list
[hidden email]
https://lists.osgeo.org/mailman/listinfo/gdal-dev
Reply | Threaded
Open this post in threaded view
|

Re: Clipping shapefile with another produces invalid shapes

Even Rouault-2

Paul,

 

> But at the red asterix, I have an invalid shape. The reason is according to

> GEOS:

> Polygon must be clockwise

 

I can't see such error emitted by GEOS source code.

 

even@even-i700:~$ grep -r "must be clockwise" geos/

even@even-i700:~$

 

 

Are you sure this comes from it ? As far as I know, GDAL and GEOS are pretty much indifferent regarding the winding order of polygons. Other tools might be more sensitive.

 

You could possibly use the ST_ForceLHR() function of spatialite (through the SQL SQLite dialect) to post process your result:

http://www.gaia-gis.it/gaia-sins/spatialite-sql-4.2.0.html

 

Even

 

--

Spatialys - Geospatial professional services

http://www.spatialys.com


_______________________________________________
gdal-dev mailing list
[hidden email]
https://lists.osgeo.org/mailman/listinfo/gdal-dev
Reply | Threaded
Open this post in threaded view
|

Re: Clipping shapefile with another produces invalid shapes

Paul Meems
Thanks Even for your quick reply.

I was using MapWinGIS to perform the validation and I assumed it was just using GEOS.
I did a closer look at the code and it is first doing some internal checks like the minimum number of points and the order of the points.
The code is checking if the polygon shape only has 1 part. If so it should be clockwise.
After that, it is using the GEOS:IsValid method.
As I understand you correctly I should interpret this message as a warning, not as an error.

I also tried your suggestion: 
ogr2ogr -sql "select ST_ForceLHR(geometry) from ogr_clipped" -dialect sqlite ogr_clipped_lhr.shp ogr_clipped.shp
It is very fast, only a few seconds on my 700k shapes, but the result is the same. I still get the 'Polygon must be clockwise' error.
Of course, MapWinGIS also has a fix method, but it is very strict and not so fast with 700k shapes.

I'll send the resulting shapefile to my client and ask if his application can process the data.

Thanks,



Paul

Paul Meems 
Release manager, configuration manager
and forum moderator of MapWindow GIS.
www.mapwindow.org

Owner of MapWindow.nl - Support for
Dutch speaking users.
www.mapwindow.nl


The MapWindow GIS project has moved to GitHub!


Download the latest MapWinGIS mapping engine.

Download the latest MapWindow 5 open source desktop application.


2017-10-12 23:28 GMT+02:00 Even Rouault <[hidden email]>:

Paul,

 

> But at the red asterix, I have an invalid shape. The reason is according to

> GEOS:

> Polygon must be clockwise

 

I can't see such error emitted by GEOS source code.

 

even@even-i700:~$ grep -r "must be clockwise" geos/

even@even-i700:~$

 

 

Are you sure this comes from it ? As far as I know, GDAL and GEOS are pretty much indifferent regarding the winding order of polygons. Other tools might be more sensitive.

 

You could possibly use the ST_ForceLHR() function of spatialite (through the SQL SQLite dialect) to post process your result:

http://www.gaia-gis.it/gaia-sins/spatialite-sql-4.2.0.html

 

Even

 

--

Spatialys - Geospatial professional services

http://www.spatialys.com



_______________________________________________
gdal-dev mailing list
[hidden email]
https://lists.osgeo.org/mailman/listinfo/gdal-dev
Reply | Threaded
Open this post in threaded view
|

Re: Clipping shapefile with another produces invalid shapes

jratike80
In reply to this post by Even Rouault-2
Even Rouault-2 wrote
> You could possibly use the ST_ForceLHR() function of spatialite (through
> the SQL SQLite
> dialect) to post process your result:
> http://www.gaia-gis.it/gaia-sins/spatialite-sql-4.2.0.html

Interesting, that ST_ForceLHR in Spatialite gives the same result than
ST_ForceRHR in PostGIS.

Description for Spatialite:
Any Polygon Ring will be oriented accordingly to Left Hand Rule (Exterior
Ring will be clockwise oriented, and Interior Rings will be
counter-clockwise oriented).

Description for PostGIS:
Forces the orientation of the vertices in a polygon to follow a
Right-Hand-Rule, in which the area that is bounded by the polygon is to the
right of the boundary. In particular, the exterior ring is orientated in a
clockwise direction and the interior rings in a counter-clockwise direction.
This function is a synonym for ST_ForcePolygonCW

Perhaps the rule of thumb would be not to talk just about right hand/left
hand rules in GIS because they are ambiguous, but mention also what it
really means like in the descriptions above, or in the GeoJSON
specification:

A linear ring MUST follow the right-hand rule with respect to the
area it bounds, i.e., exterior rings are counterclockwise, and
holes are clockwise.

-Jukka Rahkonen-




--
Sent from: http://osgeo-org.1560.x6.nabble.com/GDAL-Dev-f3742093.html
_______________________________________________
gdal-dev mailing list
[hidden email]
https://lists.osgeo.org/mailman/listinfo/gdal-dev
Reply | Threaded
Open this post in threaded view
|

Re: Clipping shapefile with another produces invalid shapes

Even Rouault-2
In reply to this post by Paul Meems

On vendredi 13 octobre 2017 09:06:46 CEST Paul Meems wrote:

> Thanks Even for your quick reply.

>

> I was using MapWinGIS to perform the validation and I assumed it was just

> using GEOS.

> I did a closer look at the code and it is first doing some internal checks

> like the minimum number of points and the order of the points.

> The code is checking if the polygon shape only has 1 part. If so it should

> be clockwise.

> After that, it is using the GEOS:IsValid method.

> As I understand you correctly I should interpret this message as a warning,

> not as an error.

 

Different software, standard, formats have differen rules on if there must be a winding order and which one it is...

 

A return of false by GEOS:IsValid() should be examined a bit more closely (but apparently you get this error before calling IsValid())

 

>

> I also tried your suggestion:

> ogr2ogr -sql "select ST_ForceLHR(geometry) from ogr_clipped" -dialect

> sqlite ogr_clipped_lhr.shp ogr_clipped.shp

 

Jukka answered on the LHR vs RHR confusion, but I completely missed the fact that you were outputting to shapefile. In that case, whatever the original orientation, or if you modified it with a SQL statement, the OGR shapefile writer will force the polygons to follow the shapefile specification, which requires the outer ring to be in clockwise order. So you shouldn't get this error... Which means that either there's an issue in some case in the shapefile driver (numerical precision issue), or in the check in MapWinGIS.

 

Even

 

--

Spatialys - Geospatial professional services

http://www.spatialys.com


_______________________________________________
gdal-dev mailing list
[hidden email]
https://lists.osgeo.org/mailman/listinfo/gdal-dev
Reply | Threaded
Open this post in threaded view
|

Re: Clipping shapefile with another produces invalid shapes

Paul Meems
Thanks Even and Jukka for the clarification.

Due to historical reasons, MapWinGIS is not using GDAL/OGR to read/write shapefiles. It is using its own implementation.
For now, I will clip with OGR first and then fix with MapWinGIS, this will produce valid shapes.
It does take a bit longer.

Should I create a test dataset to reproduce this issue. I'm not capable of debugging GDAL myself.




Paul

Paul Meems 
Release manager, configuration manager
and forum moderator of MapWindow GIS.
www.mapwindow.org

Owner of MapWindow.nl - Support for
Dutch speaking users.
www.mapwindow.nl


The MapWindow GIS project has moved to GitHub!


Download the latest MapWinGIS mapping engine.

Download the latest MapWindow 5 open source desktop application.


2017-10-13 13:01 GMT+02:00 Even Rouault <[hidden email]>:

On vendredi 13 octobre 2017 09:06:46 CEST Paul Meems wrote:

> Thanks Even for your quick reply.

>

> I was using MapWinGIS to perform the validation and I assumed it was just

> using GEOS.

> I did a closer look at the code and it is first doing some internal checks

> like the minimum number of points and the order of the points.

> The code is checking if the polygon shape only has 1 part. If so it should

> be clockwise.

> After that, it is using the GEOS:IsValid method.

> As I understand you correctly I should interpret this message as a warning,

> not as an error.

 

Different software, standard, formats have differen rules on if there must be a winding order and which one it is...

 

A return of false by GEOS:IsValid() should be examined a bit more closely (but apparently you get this error before calling IsValid())

 

>

> I also tried your suggestion:

> ogr2ogr -sql "select ST_ForceLHR(geometry) from ogr_clipped" -dialect

> sqlite ogr_clipped_lhr.shp ogr_clipped.shp

 

Jukka answered on the LHR vs RHR confusion, but I completely missed the fact that you were outputting to shapefile. In that case, whatever the original orientation, or if you modified it with a SQL statement, the OGR shapefile writer will force the polygons to follow the shapefile specification, which requires the outer ring to be in clockwise order. So you shouldn't get this error... Which means that either there's an issue in some case in the shapefile driver (numerical precision issue), or in the check in MapWinGIS.

 

Even

 

--

Spatialys - Geospatial professional services

http://www.spatialys.com



_______________________________________________
gdal-dev mailing list
[hidden email]
https://lists.osgeo.org/mailman/listinfo/gdal-dev