-
Notifications
You must be signed in to change notification settings - Fork 390
add a special projection for EPSG:27572 #2427
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
Changes from 6 commits
1a48244
aafe3f5
20f1a6a
5d38ab3
6c61a5c
950b220
bb49798
02d8842
871222c
File filter
Filter by extension
Conversations
Jump to
Diff view
Diff view
There are no files selected for viewing
| Original file line number | Diff line number | Diff line change | ||||
|---|---|---|---|---|---|---|
|
|
@@ -1815,6 +1815,59 @@ def y_limits(self): | |||||
| return self._y_limits | ||||||
|
|
||||||
|
|
||||||
| class LambertZoneII(Projection): | ||||||
| """ | ||||||
| Lambert zone II (extended) projection (https://epsg.io/27572). | ||||||
|
|
||||||
| """ | ||||||
| def __init__(self): | ||||||
| proj4_params = [('proj', 'lcc'), | ||||||
| ('lat_1', 46.8), | ||||||
| ('lat_0', 46.8), | ||||||
| ('lon_0', 0.0), | ||||||
| ('k_0', 0.99987742), | ||||||
| ('x_0', 600000.0), | ||||||
| ('y_0', 2200000.0), | ||||||
| ('pm', 'paris'), | ||||||
| ('units', 'm'), | ||||||
| ('no_defs', None)] | ||||||
|
|
||||||
| globe = Globe(semimajor_axis='6378249.2', | ||||||
| semiminor_axis='6356515.0', | ||||||
| towgs84='-168,-60,320,0,0,0,0') | ||||||
|
||||||
|
|
||||||
| super().__init__(proj4_params, globe=globe) | ||||||
|
|
||||||
| x0, x1, y0, y1 = [0, 1.2e6, 1.6e6, 2.7e6] | ||||||
|
||||||
| points = self.transform_points(self.as_geodetic(), lons, lats) |
but here it says that it is values reported using the WGS84 ellipse.
Line 678 in 694efcf
| # Geographic area of the entire dataset referenced to WGS 84 |
So, I think this may be a bug in the epsg code path and we should be doing something like:
points = self.transform_points(PlateCarree().as_geodetic(), lons, lats)
Thoughts on how that looks?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Replacing with PlateCarree does fix the initial problem mentioned in the opening message of this PR. However, the resulting bounds are not wide enough to include the southmost part of hexagonal France and Corsica.
import matplotlib.pyplot as plt
import cartopy
ax = plt.subplot(projection=cartopy.crs.epsg(27572))
ax.set_extent(cartopy.crs.epsg(27572).bounds, crs=cartopy.crs.epsg(27572))
ax.add_feature(cartopy.feature.BORDERS)
ax.add_feature(cartopy.feature.COASTLINE)
ax.scatter(564055, 2571176, transform=cartopy.crs.epsg(27572))
plt.show()There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
👍 My proposal would be to do something like this then to update the bounds manually after initialization from the EPSG crs.
class LambertZoneII(Projection):
def __init__(self):
crs = pyproj.CRS.from_epsg(27572)
super().__init__(crs.to_wkt())
# Projected bounds from https://epsg.io/27572
self.bounds = [-5242.32, 1212512.16, 1589155.51, 2706796.21]There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Ok, I've done just that. Works just fine. :-)
I have updated the test image too, as the extent is now slightly different from the one I had arbitrarily set before.
Thanks.
|
Contributor
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. This isn't a super helpful map after all with the limited boundary :) But it is small enough and shows the boundaries so perhaps still useful? Up to you whether you think it is worthwhile to keep this or not, I could be convinced either way on it.
Contributor
Author
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. Yes, unfortunately we don't see the red and blue lines as on the other test maps. But I guess this is still another useful additional unit test to have. And it does allow to see its boundaries. All in all, I would be in favour to keep it. :) |

There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Could you add some more text about this being used for France or general information about this projection so someone reading this would know more at first glance?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I have added slightly more information in the docstring.