Wrong coordinates returned by `AreaDefintion.get_lonlats` for some projections in out-of-Earth locations

In the context of https://github.com/pytroll/pyresample/pull/546, I was testing that the new boundary extraction methodology works for all projections available in cartopy, and I stumbled in wrong values returned by AreaDefinition.get_lonlats().
After a quick investigation, the error seems to arise from wrong results by pyproj transformations !

Here below I provide a code that enable to visualize the wrong coordinates as well as snapshot of what is going on for some of the problematic projections.

Problem description

The wrong results occurs only for some CRS in out-of-Earth locations.
Usually pyproj returns np.inf values in out-of-Earth locations, while with the projections listed below, AreaDefinition.get_lonlats() returns wrong/bad longitude/latitude values. In a couple of cases, even returns latitudes outside the (-90, 90) bounds.

NOTE1: These errors causes wrong resampling with pyresample !

NOTE2: Note that for all cartopy projections not listed below, pyproj returns the correct results !

Expected Output

Usually pyproj returns np.inf values in out-of-Earth locations.

Code Sample, a minimal, complete, and verifiable piece of code

import numpy as np 
import pyproj 
import cartopy
from pyresample.future.geometry import AreaDefinition 
import matplotlib.pyplot as plt 


list_bugs = [ 
    cartopy.crs.EquidistantConic(central_longitude=0.0, central_latitude=0.0, false_easting=0.0, false_northing=0.0, standard_parallels=(20.0, 50.0), globe=None),
    cartopy.crs.AzimuthalEquidistant(central_longitude=0.0, central_latitude=0.0, false_easting=0.0, false_northing=0.0, globe=None),
    cartopy.crs.AlbersEqualArea(central_longitude=0.0, central_latitude=0.0, false_easting=0.0, false_northing=0.0, standard_parallels=(20.0, 50.0), globe=None),
    cartopy.crs.LambertConformal(central_longitude=-96.0, central_latitude=39.0, false_easting=0.0, false_northing=0.0, standard_parallels=(33, 45), globe=None, cutoff=-30),
    cartopy.crs.EqualEarth(central_longitude=0, false_easting=None, false_northing=None, globe=None),
    cartopy.crs.Hammer(central_longitude=0, false_easting=None, false_northing=None, globe=None),
    cartopy.crs.Sinusoidal(central_longitude=0.0, false_easting=0.0, false_northing=0.0, globe=None),
    cartopy.crs.EckertI(central_longitude=0, false_easting=None, false_northing=None, globe=None),
    cartopy.crs.EckertII(central_longitude=0, false_easting=None, false_northing=None, globe=None),
    cartopy.crs.EckertIII(central_longitude=0, false_easting=None, false_northing=None, globe=None),
    cartopy.crs.EckertIV(central_longitude=0, false_easting=None, false_northing=None, globe=None),
    cartopy.crs.EckertV(central_longitude=0, false_easting=None, false_northing=None, globe=None),
    cartopy.crs.EckertVI(central_longitude=0, false_easting=None, false_northing=None, globe=None),
]

cartopy_crs = list_bugs[1]
for cartopy_crs in list_bugs:
    name = type(cartopy_crs).__name__
    x_min, x_max = cartopy_crs.x_limits
    y_min, y_max = cartopy_crs.y_limits
    area_extent = [x_min, y_min, x_max, y_max]
    proj_dict = cartopy_crs.to_dict()
    shape = (200, 300)
    area_def = AreaDefinition(
        crs=proj_dict,
        shape=shape,
        area_extent=area_extent
    )
    lons, lats = area_def.get_lonlats()
    
    plt.imshow(lons)
    plt.title(name+"Longitude"); plt.colorbar()
    plt.show()
    
    plt.imshow(lats)
    plt.title(name+" Latitude"); plt.colorbar()
    plt.show() 

To further investigate the erroneous results, just take values in the lower left corner (which are out of Earth in the specified projections) and look at the returned values:

from pyresample.geometry import get_geodetic_crs_with_no_datum_shift, TransformDirection, CRS
self = area_def
proj_wkt = self.crs_wkt

x, y = area_def.get_proj_coords(data_slice=(-1, 0))
crs = CRS.from_wkt(proj_wkt)
gcrs = get_geodetic_crs_with_no_datum_shift(crs)
transformer = pyproj.Transformer.from_crs(gcrs, crs, always_xy=True)
lon, lat = transformer.transform(x, y, direction=TransformDirection.INVERSE)
print(lon, lat)   # should be inf, inf 

Examples

Equidistant Conic

Azimuthal Equidistant

Albers Equal Area

Equal Earth

Lamber Conformal

Sinusoidal

Eckert

Start by running the provided cartopy projection example and inspect AreaDefinition.get_lonlats(), get_proj_coords(), and the inverse pyproj Transformer call shown in the issue. Compare out-of-Earth results with the expected inf values for the listed projections, then add regression coverage demonstrating correct invalid-coordinate handling and verify that resampling no longer uses the erroneous values.