From 3304104910f95b357cdd4e1a021b2f5c1b3f3a58 Mon Sep 17 00:00:00 2001 From: Dongdong Tian Date: Wed, 27 Nov 2024 00:29:48 +0800 Subject: [PATCH 1/8] load_tile_map: Add the 'crs' parameter to set the CRS of the returned image --- pygmt/datasets/tile_map.py | 36 ++++++++++++++++++++++++++++++------ 1 file changed, 30 insertions(+), 6 deletions(-) diff --git a/pygmt/datasets/tile_map.py b/pygmt/datasets/tile_map.py index 0c9f1b6271d..f7bdc217b6e 100644 --- a/pygmt/datasets/tile_map.py +++ b/pygmt/datasets/tile_map.py @@ -3,7 +3,6 @@ :class:`xarray.DataArray`. """ -import contextlib from collections.abc import Sequence from typing import Literal @@ -11,17 +10,22 @@ try: import contextily + from rasterio.crs import CRS from xyzservices import TileProvider _HAS_CONTEXTILY = True except ImportError: + CRS = None TileProvider = None _HAS_CONTEXTILY = False -with contextlib.suppress(ImportError): - # rioxarray is needed to register the rio accessor +try: import rioxarray # noqa: F401 + _HAS_RIOXARRAY = True +except ImportError: + _HAS_RIOXARRAY = False + import numpy as np import xarray as xr @@ -33,6 +37,7 @@ def load_tile_map( zoom: int | Literal["auto"] = "auto", source: TileProvider | str | None = None, lonlat: bool = True, + crs: str | CRS = "EPSG:3857", wait: int = 0, max_retries: int = 2, zoom_adjust: int | None = None, @@ -42,7 +47,8 @@ def load_tile_map( The tiles that compose the map are merged and georeferenced into an :class:`xarray.DataArray` image with 3 bands (RGB). Note that the returned image is - in a Spherical Mercator (EPSG:3857) coordinate reference system. + in a Spherical Mercator (EPSG:3857) coordinate reference system (CRS) by default, + but can be customized using the ``crs`` parameter. Parameters ---------- @@ -80,6 +86,10 @@ def load_tile_map( lonlat If ``False``, coordinates in ``region`` are assumed to be Spherical Mercator as opposed to longitude/latitude. + crs + Coordinate reference system (CRS) of the returned :class:`xarray.DataArray` + image. Default is ``"EPSG:3857"`` (i.e., Spherical Mercator). The CRS can be in + either string or :class:`rasterio.crs.CRS` format. wait If the tile API is rate-limited, the number of seconds to wait between a failed request and the next try. @@ -128,6 +138,8 @@ def load_tile_map( ... raster.rio.crs.to_string() 'EPSG:3857' """ + _default_crs = "EPSG:3857" # The default CRS returned from contextily + if not _HAS_CONTEXTILY: msg = ( "Package `contextily` is required to be installed to use this function. " @@ -136,6 +148,14 @@ def load_tile_map( ) raise ImportError(msg) + if crs != _default_crs and not _HAS_RIOXARRAY: + msg = ( + f"Package `rioxarray` is required if CRS is not '{_default_crs}'. " + "Please use `python -m pip install rioxarray` or " + "`mamba install -c conda-forge rioxarray` to install the package." + ) + raise ImportError(msg) + contextily_kwargs = {} if zoom_adjust is not None: contextily_kwargs["zoom_adjust"] = zoom_adjust @@ -176,8 +196,12 @@ def load_tile_map( dims=("band", "y", "x"), ) - # If rioxarray is installed, set the coordinate reference system + # If rioxarray is installed, set the coordinate reference system. if hasattr(dataarray, "rio"): - dataarray = dataarray.rio.write_crs(input_crs="EPSG:3857") + dataarray = dataarray.rio.write_crs(input_crs=_default_crs) + + # Reproject raster image from EPSG:3857 to specified CRS, using rioxarray. + if crs != _default_crs: + dataarray = dataarray.rio.reproject(dst_crs=crs) return dataarray From cfb69d2259782b3876a3b0db0f64dba7fcb76bae Mon Sep 17 00:00:00 2001 From: Dongdong Tian Date: Wed, 27 Nov 2024 00:31:51 +0800 Subject: [PATCH 2/8] tilemap: Use the new 'crs' parameter to simplify the tilemap function --- pygmt/src/tilemap.py | 27 +++------------------------ 1 file changed, 3 insertions(+), 24 deletions(-) diff --git a/pygmt/src/tilemap.py b/pygmt/src/tilemap.py index a13c8d9c740..ab65c3ac3c9 100644 --- a/pygmt/src/tilemap.py +++ b/pygmt/src/tilemap.py @@ -9,13 +9,9 @@ from pygmt.helpers import build_arg_list, fmt_docstring, kwargs_to_strings, use_alias try: - import rioxarray # noqa: F401 from xyzservices import TileProvider - - _HAS_RIOXARRAY = True except ImportError: TileProvider = None - _HAS_RIOXARRAY = False @fmt_docstring @@ -111,38 +107,21 @@ def tilemap( kwargs : dict Extra keyword arguments to pass to :meth:`pygmt.Figure.grdimage`. - - Raises - ------ - ImportError - If ``rioxarray`` is not installed. Follow - :doc:`install instructions for rioxarray `, (e.g. via - ``python -m pip install rioxarray``) before using this function. """ kwargs = self._preprocess(**kwargs) - if not _HAS_RIOXARRAY: - raise ImportError( - "Package `rioxarray` is required to be installed to use this function. " - "Please use `python -m pip install rioxarray` or " - "`mamba install -c conda-forge rioxarray` to install the package." - ) - raster = load_tile_map( region=region, zoom=zoom, source=source, lonlat=lonlat, + crs="OGC:CRS84" if lonlat is True else "EPSG:3857", wait=wait, max_retries=max_retries, zoom_adjust=zoom_adjust, ) - - # Reproject raster from Spherical Mercator (EPSG:3857) to lonlat (OGC:CRS84) if - # bounding box region was provided in lonlat - if lonlat and raster.rio.crs == "EPSG:3857": - raster = raster.rio.reproject(dst_crs="OGC:CRS84") - raster.gmt.gtype = 1 # set to geographic type + if lonlat: + raster.gmt.gtype = 1 # Set to geographic type # Only set region if no_clip is None or False, so that plot is clipped to exact # bounding box region From 00e260580e9e8524bfe2d667882a7b2e95d91898 Mon Sep 17 00:00:00 2001 From: Dongdong Tian Date: Wed, 27 Nov 2024 14:21:15 +0800 Subject: [PATCH 3/8] Improve tests --- pygmt/tests/test_tilemap.py | 53 +++++++++++++++++++++++++++++++++++-- 1 file changed, 51 insertions(+), 2 deletions(-) diff --git a/pygmt/tests/test_tilemap.py b/pygmt/tests/test_tilemap.py index 704a96093af..6ec04decf95 100644 --- a/pygmt/tests/test_tilemap.py +++ b/pygmt/tests/test_tilemap.py @@ -2,14 +2,24 @@ Test Figure.tilemap. """ +import importlib +from unittest.mock import patch + import pytest from pygmt import Figure -contextily = pytest.importorskip("contextily") -rioxarray = pytest.importorskip("rioxarray") +try: + import contextily + + _HAS_CONTEXTILY = True +except ImportError: + _HAS_CONTEXTILY = False + +_HAS_RIOXARRAY = bool(importlib.util.find_spec("rioxarray")) @pytest.mark.mpl_image_compare +@pytest.mark.skipif(not _HAS_CONTEXTILY, reason="contextily is not installed") def test_tilemap_web_mercator(): """ Create a tilemap plot in Spherical Mercator projection (EPSG:3857). @@ -27,6 +37,10 @@ def test_tilemap_web_mercator(): @pytest.mark.benchmark @pytest.mark.mpl_image_compare +@pytest.mark.skipif( + not (_HAS_CONTEXTILY and _HAS_RIOXARRAY), + reason="contextily and rioxarray are not installed", +) def test_tilemap_ogc_wgs84(): """ Create a tilemap plot using longitude/latitude coordinates (OGC:WGS84), centred on @@ -45,6 +59,10 @@ def test_tilemap_ogc_wgs84(): @pytest.mark.mpl_image_compare @pytest.mark.parametrize("no_clip", [False, True]) +@pytest.mark.skipif( + not (_HAS_CONTEXTILY and _HAS_RIOXARRAY), + reason="contextily and rioxarray are not installed", +) def test_tilemap_no_clip(no_clip): """ Create a tilemap plot clipped to the Southern Hemisphere when no_clip is False, but @@ -60,3 +78,34 @@ def test_tilemap_no_clip(no_clip): no_clip=no_clip, ) return fig + + +@pytest.mark.skipif(_HAS_CONTEXTILY, reason="contextily is installed.") +def test_tilemap_no_contextily(): + """ + Raise an ImportError when contextily is not installed. + """ + fig = Figure() + with pytest.raises(ImportError, match="Package `contextily` is required"): + fig.tilemap( + region=[-20000000.0, 20000000.0, -20000000.0, 20000000.0], + zoom=0, + lonlat=False, + frame="afg", + ) + + +@pytest.mark.skipif(_HAS_RIOXARRAY, reason="rioxarray is installed.") +def test_tilemap_no_rioxarray(): + """ + Raise an ImportError when rioxarray is not installed and contextily is installed. + """ + fig = Figure() + # In our CI, contextily and rioxarray are installed together, so we will see the + # error about contextily, not rioxarray. Here we mock contextily as installed, to + # make sure that we see the rioxarray error message when rioxarray is not installed. + with patch("pygmt.datasets.tile_map._HAS_CONTEXTILY", True): + with pytest.raises(ImportError, match="Package `rioxarray` is required"): + fig.tilemap( + region=[-180.0, 180.0, -90, 90], zoom=0, lonlat=True, frame="afg" + ) From 3265240fda26a75d42f636b3291473dc4b634e22 Mon Sep 17 00:00:00 2001 From: Dongdong Tian Date: Mon, 30 Sep 2024 22:41:51 +0800 Subject: [PATCH 4/8] Refactor to store common contextily parameters in kwdict --- pygmt/datasets/tile_map.py | 37 +++++++++++++++++-------------------- 1 file changed, 17 insertions(+), 20 deletions(-) diff --git a/pygmt/datasets/tile_map.py b/pygmt/datasets/tile_map.py index f7bdc217b6e..a76169c1609 100644 --- a/pygmt/datasets/tile_map.py +++ b/pygmt/datasets/tile_map.py @@ -156,28 +156,25 @@ def load_tile_map( ) raise ImportError(msg) - contextily_kwargs = {} - if zoom_adjust is not None: - contextily_kwargs["zoom_adjust"] = zoom_adjust - if Version(contextily.__version__) < Version("1.5.0"): - msg = ( - "The `zoom_adjust` parameter requires `contextily>=1.5.0` to work. " - "Please upgrade contextily, or manually set the `zoom` level instead." - ) - raise TypeError(msg) - + if zoom_adjust is not None and Version(contextily.__version__) < Version("1.5.0"): + msg = ( + "The `zoom_adjust` parameter requires `contextily>=1.5.0` to work. " + "Please upgrade contextily, or manually set the `zoom` level instead." + ) + raise ValueError(msg) + + # Keyword arguments for contextily.bounds2img + contextily_kwargs = { + "zoom": zoom, + "source": source, + "ll": lonlat, + "wait": wait, + "max_retries": max_retries, + "zoom_adjust": zoom_adjust, + } west, east, south, north = region image, extent = contextily.bounds2img( - w=west, - s=south, - e=east, - n=north, - zoom=zoom, - source=source, - ll=lonlat, - wait=wait, - max_retries=max_retries, - **contextily_kwargs, + w=west, s=south, e=east, n=north, **contextily_kwargs ) # Turn RGBA img from channel-last to channel-first and get 3-band RGB only From 421c167cc81abb507a4612f89e8f5dff0434d934 Mon Sep 17 00:00:00 2001 From: Dongdong Tian Date: Tue, 3 Dec 2024 08:39:28 +0800 Subject: [PATCH 5/8] zoom_adjust only makes sense for contextily>=1.5.0 --- pygmt/datasets/tile_map.py | 17 +++++++++-------- 1 file changed, 9 insertions(+), 8 deletions(-) diff --git a/pygmt/datasets/tile_map.py b/pygmt/datasets/tile_map.py index a76169c1609..67b213ac913 100644 --- a/pygmt/datasets/tile_map.py +++ b/pygmt/datasets/tile_map.py @@ -156,13 +156,6 @@ def load_tile_map( ) raise ImportError(msg) - if zoom_adjust is not None and Version(contextily.__version__) < Version("1.5.0"): - msg = ( - "The `zoom_adjust` parameter requires `contextily>=1.5.0` to work. " - "Please upgrade contextily, or manually set the `zoom` level instead." - ) - raise ValueError(msg) - # Keyword arguments for contextily.bounds2img contextily_kwargs = { "zoom": zoom, @@ -170,8 +163,16 @@ def load_tile_map( "ll": lonlat, "wait": wait, "max_retries": max_retries, - "zoom_adjust": zoom_adjust, } + if zoom_adjust is not None: + if Version(contextily.__version__) < Version("1.5.0"): + msg = ( + "The `zoom_adjust` parameter requires `contextily>=1.5.0` to work. " + "Please upgrade contextily, or manually set the `zoom` level instead." + ) + raise ValueError(msg) + contextily_kwargs["zoom_adjust"] = zoom_adjust + west, east, south, north = region image, extent = contextily.bounds2img( w=west, s=south, e=east, n=north, **contextily_kwargs From 8346af7baf2ea3bf5f11d4c1d672b0a63fd9a91a Mon Sep 17 00:00:00 2001 From: Dongdong Tian Date: Tue, 3 Dec 2024 10:06:19 +0800 Subject: [PATCH 6/8] Call rio.reproject only if the rio accessor is available Co-authored-by: Wei Ji <23487320+weiji14@users.noreply.github.com> --- pygmt/datasets/tile_map.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/pygmt/datasets/tile_map.py b/pygmt/datasets/tile_map.py index 67b213ac913..750c59f2d3c 100644 --- a/pygmt/datasets/tile_map.py +++ b/pygmt/datasets/tile_map.py @@ -198,8 +198,8 @@ def load_tile_map( if hasattr(dataarray, "rio"): dataarray = dataarray.rio.write_crs(input_crs=_default_crs) - # Reproject raster image from EPSG:3857 to specified CRS, using rioxarray. - if crs != _default_crs: - dataarray = dataarray.rio.reproject(dst_crs=crs) + # Reproject raster image from EPSG:3857 to specified CRS, using rioxarray. + if crs != _default_crs: + dataarray = dataarray.rio.reproject(dst_crs=crs) return dataarray From 6733e196a6ec5b845e235cbe9f44581b94d2c10d Mon Sep 17 00:00:00 2001 From: Dongdong Tian Date: Tue, 3 Dec 2024 14:16:16 +0800 Subject: [PATCH 7/8] Use provider's CRS attribute if provided --- pygmt/datasets/tile_map.py | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/pygmt/datasets/tile_map.py b/pygmt/datasets/tile_map.py index 750c59f2d3c..4d305bd45ce 100644 --- a/pygmt/datasets/tile_map.py +++ b/pygmt/datasets/tile_map.py @@ -138,7 +138,9 @@ def load_tile_map( ... raster.rio.crs.to_string() 'EPSG:3857' """ - _default_crs = "EPSG:3857" # The default CRS returned from contextily + # The CRS of the returned image. Use the tile provider's CRS if available. + # Otherwise, default to EPSG:3857. + _default_crs = getattr(source, "crs", "EPSG:3857") if not _HAS_CONTEXTILY: msg = ( From c1d55b0857db751b9c126f9151496c7882f9977c Mon Sep 17 00:00:00 2001 From: Dongdong Tian Date: Tue, 3 Dec 2024 14:25:32 +0800 Subject: [PATCH 8/8] Rename _default_crs to _source_crs --- pygmt/datasets/tile_map.py | 16 ++++++++-------- 1 file changed, 8 insertions(+), 8 deletions(-) diff --git a/pygmt/datasets/tile_map.py b/pygmt/datasets/tile_map.py index 4d305bd45ce..d18bce423f8 100644 --- a/pygmt/datasets/tile_map.py +++ b/pygmt/datasets/tile_map.py @@ -138,9 +138,9 @@ def load_tile_map( ... raster.rio.crs.to_string() 'EPSG:3857' """ - # The CRS of the returned image. Use the tile provider's CRS if available. - # Otherwise, default to EPSG:3857. - _default_crs = getattr(source, "crs", "EPSG:3857") + # The CRS of the source tile provider. If the source is a TileProvider object, use + # its crs attribute if available. Otherwise, default to EPSG:3857. + _source_crs = getattr(source, "crs", "EPSG:3857") if not _HAS_CONTEXTILY: msg = ( @@ -150,9 +150,9 @@ def load_tile_map( ) raise ImportError(msg) - if crs != _default_crs and not _HAS_RIOXARRAY: + if crs != _source_crs and not _HAS_RIOXARRAY: msg = ( - f"Package `rioxarray` is required if CRS is not '{_default_crs}'. " + f"Package `rioxarray` is required if CRS is not '{_source_crs}'. " "Please use `python -m pip install rioxarray` or " "`mamba install -c conda-forge rioxarray` to install the package." ) @@ -198,10 +198,10 @@ def load_tile_map( # If rioxarray is installed, set the coordinate reference system. if hasattr(dataarray, "rio"): - dataarray = dataarray.rio.write_crs(input_crs=_default_crs) + dataarray = dataarray.rio.write_crs(input_crs=_source_crs) - # Reproject raster image from EPSG:3857 to specified CRS, using rioxarray. - if crs != _default_crs: + # Reproject raster image from the source CRS to the specified CRS. + if crs != _source_crs: dataarray = dataarray.rio.reproject(dst_crs=crs) return dataarray