possible fix for always assuming WGS84

.gitignore

@@ -130,3 +130,6 @@ dmypy.json
 
 # PyCharm:
 .idea
+
+# vscode:
+.vscode

morecantile/models.py

@@ -40,7 +40,6 @@
 NumType = Union[float, int]
 BoundsType = Tuple[NumType, NumType]
 LL_EPSILON = 1e-11
-WGS84_CRS = pyproj.CRS.from_epsg(4326)
 axesInfo = Annotated[List[str], Field(min_length=2, max_length=2)]
 
 
@@ -486,7 +485,7 @@ class TileMatrixSet(BaseModel, arbitrary_types_allowed=True):
     ]
 
     # Private attributes
-    _geographic_crs: pyproj.CRS = PrivateAttr(default=WGS84_CRS)
+    _geographic_crs: pyproj.CRS = PrivateAttr()
     _to_geographic: pyproj.Transformer = PrivateAttr()
     _from_geographic: pyproj.Transformer = PrivateAttr()
 
@@ -495,7 +494,7 @@ def __init__(self, **data):
         super().__init__(**data)
 
         self._geographic_crs = pyproj.CRS.from_user_input(
-            data.get("_geographic_crs", WGS84_CRS)
+            data.get("_geographic_crs", self.crs._pyproj_crs.geodetic_crs)
         )
 
         try:
@@ -656,7 +655,7 @@ def custom(
         title: Optional[str] = None,
         id: Optional[str] = None,
         ordered_axes: Optional[List[str]] = None,
-        geographic_crs: pyproj.CRS = WGS84_CRS,
+        geographic_crs: pyproj.CRS = None,
         screen_pixel_size: float = 0.28e-3,
         decimation_base: int = 2,
         **kwargs: Any,
@@ -689,8 +688,8 @@ def custom(
             Tile Matrix Set title
         id: str, optional
             Tile Matrix Set identifier
-        geographic_crs: pyproj.CRS
-            Geographic (lat,lon) coordinate reference system (default is EPSG:4326)
+        geographic_crs: pyproj.CRS, optional
+            Geographic (lat,lon) coordinate reference system (default is CRS's geodetic CRS)
         ordered_axes: list of str, optional
             Override Axis order (e.g `["N", "S"]`) else default to CRS's metadata
         screen_pixel_size: float, optional
@@ -705,6 +704,9 @@ def custom(
         TileMatrixSet
 
         """
+        # TODO is this correct to always assume?
+        geographic_crs = geographic_crs or crs.geodetic_crs
+
         matrix_scale = matrix_scale or [1, 1]
 
         if ordered_axes:

tests/test_models.py

@@ -6,6 +6,7 @@
 from collections.abc import Iterable
 
 import pyproj
+import pyproj.crs
 import pytest
 from pydantic import ValidationError
 from rasterio.crs import CRS as rioCRS
@@ -289,6 +290,22 @@ def test_zoom_for_res():
     assert tms.zoom_for_res(5000) == 6
 
 
+def test_custom_not_earth():
+    crs = pyproj.CRS.from_user_input("IAU_2015:49900")
+    extent = [-90, -180, 90, 180]
+    mars_tms = morecantile.TileMatrixSet.custom(extent, crs, id="MarsGeographicCRS")
+    assert "4326" not in mars_tms.geographic_crs.to_wkt()
+    assert "4326" not in mars_tms.rasterio_geographic_crs.to_wkt()
+    wkt_mars_web_mercator = 'PROJCRS["Mars (2015) - Sphere XY / Pseudo-Mercator",BASEGEOGCRS["Mars (2015) - Sphere",DATUM["Mars (2015) - Sphere",ELLIPSOID["Mars (2015) - Sphere",3396190,0,LENGTHUNIT["metre",1,ID["EPSG",9001]]],ANCHOR["Viking 1 lander : 47.95137 W"]],PRIMEM["Reference Meridian",0,ANGLEUNIT["degree",0.0174532925199433,ID["EPSG",9122]]]],CONVERSION["Popular Visualisation Pseudo-Mercator",METHOD["Popular Visualisation Pseudo Mercator",ID["EPSG",1024]],PARAMETER["Latitude of natural origin",0,ANGLEUNIT["degree",0.0174532925199433],ID["EPSG",8801]],PARAMETER["Longitude of natural origin",0,ANGLEUNIT["degree",0.0174532925199433],ID["EPSG",8802]],PARAMETER["False easting",0,LENGTHUNIT["metre",1],ID["EPSG",8806]],PARAMETER["False northing",0,LENGTHUNIT["metre",1],ID["EPSG",8807]]],CS[Cartesian,2],AXIS["easting (X)",east,ORDER[1],LENGTHUNIT["metre",1,ID["EPSG",9001]]],AXIS["northing (Y)",north,ORDER[2],LENGTHUNIT["metre",1,ID["EPSG",9001]]],USAGE[SCOPE["Web mapping and visualisation."],AREA["World between 85.06 S and 85.06 N."],BBOX[-85.850511287,-180,85.0511287,180]],REMARK["Use semi-major radius as sphere radius for interoperability. Source of IAU Coordinate systems: doi:10.1007/s10569-017-9805-5"]]'
+    crs_mars_web_mercator = pyproj.CRS.from_wkt(wkt_mars_web_mercator)
+    extent_wm = [-85.850511287, -180, 85.0511287, 180]
+    mars_tms_wm = morecantile.TileMatrixSet.custom(
+        extent_wm, crs_mars_web_mercator, id="MarsWebMercator"
+    )
+    assert "4326" not in mars_tms_wm.geographic_crs.to_wkt()
+    assert "4326" not in mars_tms_wm.rasterio_geographic_crs.to_wkt()
+
+
 def test_schema():
     """Translate Model to Schema."""
     tms = morecantile.tms.get("WebMercatorQuad")
@@ -300,13 +317,15 @@ def test_schema():
         "+proj=stere +lat_0=90 +lon_0=0 +k=2 +x_0=0 +y_0=0 +R=3396190 +units=m +no_defs"
     )
     extent = [-13584760.000, -13585240.000, 13585240.000, 13584760.000]
-    with pytest.warns(UserWarning):
-        tms = morecantile.TileMatrixSet.custom(extent, crs, id="MarsNPolek2MOLA5k")
+    tms = morecantile.TileMatrixSet.custom(extent, crs, id="MarsNPolek2MOLA5k")
     assert tms.model_json_schema()
     assert tms.model_dump(exclude_none=True)
     json_doc = json.loads(tms.model_dump_json(exclude_none=True))
     assert json_doc["crs"] == "http://www.opengis.net/def/crs/IAU/2015/49930"
-
+    # test to ensure a correct geographic crs
+    mars_tms = TileMatrixSet.model_validate(json_doc)
+    assert "4326" not in mars_tms.geographic_crs.to_wkt()
+    assert "4326" not in mars_tms.rasterio_geographic_crs.to_wkt()
     crs = pyproj.CRS.from_epsg(3031)
     extent = [-948.75, -543592.47, 5817.41, -3333128.95]  # From https:///epsg.io/3031
     tms = morecantile.TileMatrixSet.custom(extent, crs, id="MyCustomTmsEPSG3031")