Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

pyshp-2.2.0 test failure #2012

Closed
manisandro opened this issue Mar 9, 2022 · 7 comments · Fixed by #2028
Closed

pyshp-2.2.0 test failure #2012

manisandro opened this issue Mar 9, 2022 · 7 comments · Fixed by #2028
Milestone

Comments

@manisandro
Copy link

I'm building proj9 for fedora, and rebuilding cartopy against it yielded the following test failure:

=================================== FAILURES ===================================
___________________________ TestLakes.test_geometry ____________________________
[gw0] linux -- Python 3.10.2 /usr/bin/python3

self = <cartopy.tests.test_shapereader.TestLakes object at 0x7f7df66ec190>

    def test_geometry(self):
        lake_geometry = self.test_lake_geometry
        assert lake_geometry.type == 'Polygon'
    
        polygon = lake_geometry
    
        expected = np.array([(-84.85548682324658, 11.147898667846633),
                            (-85.29013729525353, 11.176165676310276),
                            (-85.79132117383625, 11.509737046754324),
                            (-85.8851655748783, 11.900100816287136),
                            (-85.5653401354239, 11.940330918826362),
                            (-85.03684526237491, 11.5216484643976),
                            (-84.85548682324658, 11.147898667846633),
                            (-84.85548682324658, 11.147898667846633)])
    
>       assert_array_almost_equal(expected, polygon.exterior.coords)
E       AssertionError: 
E       Arrays are not almost equal to 6 decimals
E       
E       Mismatched elements: 12 / 16 (75%)
E       Max absolute difference: 0.75447591
E       Max relative difference: 0.00887234
E        x: array([[-84.855487,  11.147899],
E              [-85.290137,  11.176166],
E              [-85.791321,  11.509737],...
E        y: array([[-84.855487,  11.147899],
E              [-84.855487,  11.147899],
E              [-85.036845,  11.521648],...

../../BUILDROOT/python-cartopy-0.20.1-1.fc37.x86_64/usr/lib64/python3.10/site-packages/cartopy/tests/test_shapereader.py:45: AssertionError
---------------------------- Captured stdout setup -----------------------------
@manisandro
Copy link
Author

Actually this looks like a bug in pyshp-2.2.0. With pyshp-2.1.3 the test passes.

@manisandro manisandro changed the title Proj9 test failure pyshp-2.2.0 test failure Mar 9, 2022
@manisandro
Copy link
Author

So upon further investigation, pyshp-2.2.0 enforces the positive orientation of the outer polygon ring, while pyshp-2.1.3 does not. The lake_geometry in the test is negatively oriented, and the ring of expected points is also negatively orientated. Since pyshp-2.2.0 reverses the orientation of the ring, assert_array_almost_equal(expected, polygon.exterior.coords) will fail. See also GeospatialPython/pyshp@298d981

@manisandro manisandro reopened this Mar 9, 2022
@manisandro
Copy link
Author

I just noticed there's also the issue that depending on whether pyshp or fiona are used to read the shapefiles, the returned ring orientations will differ, as fiona won't reverse the orientation. Not sure what convention cartopy expects, but it might also something that needs to be solved by shapely (which deals with transforming the geometry from pyshp into a shapely.geometry which is what cartopy ends up using AFAICS.

CC @QuLogic as Fedora cartopy maintainer.

@QuLogic
Copy link
Member

QuLogic commented Mar 9, 2022

Yes, I commented on that commit; I'm still confused as to whether this change is intentional.

@karimbahgat
Copy link
Contributor

The switch to enforcing ring orientation in pyshp 2.2 was indeed intentional and in an effort to conform to the new GeoJSON RFC7946 specification, where it is one of the requirements. In theory the switch should not have caused any types of issues, since the original 2008 GeoJSON specification did not say anything about ring orientation. The exception as evidenced in this issue is for test cases where expected coordinate sequences have to be made explicit.

However, after reading more up on the GeoJSON RFC7946 spec, I have come to conclude that the next version pyshp will probably revert back to not enforcing ring orientation as before (instead keeping the rings in whichever orientation they were stored in the shapefile). I had not realized this before, but RFC7946 also requires that all coordinates must be in WGS84 coordinates, and some other rules about geometry splitting around the antimeridian. Supporting these requirements would then mean that pyshp would have to be dependent on pyproj in order to do projection transformations for non-WGS84 shapefiles and possibly shapely to do clipping where necessary. This would also make things much slower.

Indeed, the reason fiona does not enforce ring orientation is that it relies on GDAL to generate the geojson, and the GDAL GeoJSON driver still uses the 2008 specification by default. Given the backwards incompatible changes, the consensus in the GDAL community appears to have been that the switch to RFC7946 should be manually set and should be "all-or-nothing". Since pyshp cannot reasonably be expected to support "all" of RFC7946, I think it's better to only support the 2008 GeoJSON spec and just make this clear in the docs. After all, the purpose of pyshp is not to generate perfectly valid geojson, but rather to have quick access to the geometry contents in a human-readable standard - for these purposes 2008 GeoJSON is sufficient enough. If anyone wants to generate RFC7946-valid GeoJSON for file-based storage or sending across the internet, the burden should then fall on the user or some third-party library to do the necessary coordinate projection, clipping, etc.

While this means that cartopy's test error will disappear when I revert this in the next pyshp version, I do suspect this issue might come up again soon. As the only official standard, it's likely that RFC7946 will only become more widely used, and if GDAL makes RFC7946 the default in future versions then fiona will produce geojson which enforces ring orientation as pyshp 2.2 does, and the cartopy test will fail once again. Given this, maybe a more robust fix here is to enforce a particular ring orientation in the test itself? This might be as simple as sorting the rings before comparing them: assert_array_almost_equal(sorted(expected), sorted(polygon.exterior.coords))?

@QuLogic
Copy link
Member

QuLogic commented Mar 26, 2022

Thanks for clarifying things @karimbahgat.

@QuLogic QuLogic added this to the 0.21 milestone Apr 26, 2022
karimbahgat added a commit to GeospatialPython/pyshp that referenced this issue Apr 27, 2022
This is only required in the RFC7946 GeoJSON specification, which it won't be possible to fully support in pyshp. Instead sticking to the original 2008 GeoJSON, keeping all rings in the same orientation as stored in the shapefile. More details at SciTools/cartopy#2012.
@karimbahgat
Copy link
Contributor

For info, this test should no longer fail for pyshp v2.3.0. Geojson ring orientation has now been reverted back to preserving the original shapefile ring orientation (not enforcing ring orientation based on exterior/interior). See a9ba72da52743850164185dd8db5fd821494fd11 for details.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging a pull request may close this issue.

3 participants