add a special projection for EPSG:27572 (#2427)

* add a projection for L2E

this projection is relevant for mainland France, and is necessary
because `cartopy.crs.epsg(27572)` does not work when used as a
``projection`` argument for geoaxes due to the bounds not being projected wide enough.

lib/cartopy/crs.py

@@ -18,6 +18,7 @@
 import warnings
 
 import numpy as np
+import pyproj
 from pyproj import Transformer
 from pyproj.exceptions import ProjError
 import shapely.geometry as sgeom
@@ -685,7 +686,9 @@ def __init__(self, *args, **kwargs):
             y1 = self.area_of_use.north
             lons = np.array([x0, x0, x1, x1])
             lats = np.array([y0, y1, y1, y0])
-            points = self.transform_points(self.as_geodetic(), lons, lats)
+            points = self.transform_points(
+                PlateCarree().as_geodetic(), lons, lats
+            )
             x = points[:, 0]
             y = points[:, 1]
             self.bounds = (x.min(), x.max(), y.min(), y.max())
@@ -1815,6 +1818,20 @@ def y_limits(self):
         return self._y_limits
 
 
+class LambertZoneII(Projection):
+    """
+    Lambert zone II (extended) projection (https://epsg.io/27572), a
+    legacy projection that covers hexagonal France and Corsica.
+
+    """
+    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]
+
+
 class LambertAzimuthalEqualArea(Projection):
     """
     A Lambert Azimuthal Equal-Area projection.

lib/cartopy/tests/crs/test_lambert_conformal.py

@@ -114,3 +114,29 @@ def test_single_npole(self):
         assert_array_almost_equal(n_pole_crs.y_limits,
                                   expected_y,
                                   decimal=0)
+
+
+class TestLambertZoneII:
+    def setup_class(self):
+        self.point_a = (1.4868268900254693, 48.13277955695077)
+        self.point_b = (-2.3188020040300126, 48.68412929316207)
+        self.src_crs = ccrs.PlateCarree()
+        self.nan = float('nan')
+
+    def test_default(self):
+        proj = ccrs.LambertZoneII()
+        res = proj.transform_point(*self.point_a, src_crs=self.src_crs)
+        np.testing.assert_array_almost_equal(res,
+                                             (536690.18620, 2348515.62248),
+                                             decimal=5)
+        res = proj.transform_point(*self.point_b, src_crs=self.src_crs)
+        np.testing.assert_array_almost_equal(res,
+                                             (257199.57387, 2419655.71471),
+                                             decimal=5)
+
+    def test_nan(self):
+        proj = ccrs.LambertZoneII()
+        res = proj.transform_point(0.0, float('nan'), src_crs=self.src_crs)
+        assert np.all(np.isnan(res))
+        res = proj.transform_point(float('nan'), 0.0, src_crs=self.src_crs)
+        assert np.all(np.isnan(res))

lib/cartopy/tests/mpl/test_mpl_integration.py

@@ -172,6 +172,7 @@ def test_simple_global():
     pytest.param((ccrs.InterruptedGoodeHomolosine, dict(emphasis='land')),
                  id='InterruptedGoodeHomolosine'),
     ccrs.LambertCylindrical,
+    ccrs.LambertZoneII,
     pytest.param((ccrs.Mercator, dict(min_latitude=-85, max_latitude=85)),
                  id='Mercator'),
     ccrs.Miller,