Merge pull request #156 from raphaelquast/dev

Merge for v6.2.1

docs/api.rst

@@ -1734,8 +1734,11 @@ Pre-defined WebMap services:
     S1GBM
     S2_cloudless
     GEBCO
+    GMRT
+    GLAD
     CAMS
     DLR_basemaps
+    OpenPlanetary
 
 **Services specific for Austria (Europe)**
 
@@ -2343,8 +2346,17 @@ Logos can be re-positioned and re-sized with the :ref:`layout_editor`!
 🌈 Colorbars (with a histogram)
 -------------------------------
 
-| Before adding a colorbar, you must plot the data using ``m.plot_map()``.
-| A colorbar with a colored histogram on top can then be added to the map via ``m.add_colorbar``.
+Before adding a colorbar, you must plot the data using ``m.plot_map(vmin=..., vmax=...)``.
+
+- ``vmin`` and ``vmax`` hereby specify the value-range used for assigning colors (e.g. the limits of the colorbar).
+- If no explicit limits are provided, the min/max values of the data are used.
+
+Once a dataset has been plotted, a colorbar with a colored histogram on top can be added to the map by calling ``m.add_colorbar()``.
+
+
+.. note::
+    | The colorbar always represents the dataset that was used in the last call to ``m.plot_map()``.
+    | If you need multiple colorbars, use an individual ``Maps`` object for each dataset! (e.g. via ``m2  = m.new_layer()``)
 
 
 .. note::
@@ -2392,11 +2404,6 @@ Logos can be re-positioned and re-sized with the :ref:`layout_editor`!
     |                                                                    |                                          |
     +--------------------------------------------------------------------+------------------------------------------+
 
-.. note::
-    | You must plot a dataset first! (e.g. by calling ``m.plot_map()``)
-    | The colorbar always represents the dataset that was used in the last call to ``m.plot_map()``.
-    | If you need multiple colorbars, use an individual ``Maps`` object for each dataset! (e.g. via ``m2  = m.new_layer()``)
-
 
 
 Once the colorbar has been created, the colorbar-object can be accessed via ``m.colorbar``.

eomaps/_data_manager.py

@@ -117,10 +117,14 @@ def set_props(
             self._y0min, self._y0max = np.nanmin(self.y0), np.nanmax(self.y0)
 
         # estimate the radius (used as margin on data selection)
-        self._r = self.m._shapes._estimate_radius(self.m, "out")
-        if self._r is not None and all(np.isfinite(i) for i in self._r):
-            self._radius_margin = [i * self._radius_margin_factor for i in self._r]
-        else:
+        try:
+            self._r = self.m._shapes._estimate_radius(self.m, "out")
+            if self._r is not None and all(np.isfinite(i) for i in self._r):
+                self._radius_margin = [i * self._radius_margin_factor for i in self._r]
+            else:
+                self._radius_margin = None
+        except Exception:
+            self._r = None
             self._radius_margin = None
 
         if update_coll_on_fetch:
@@ -294,10 +298,13 @@ def _set_lims(self):
         # (1 x radius) to the map
         if self._r is not None and all(np.isfinite(self._r)):
             rx, ry = self._r
-            x0min = self._x0min - rx
-            y0min = self._y0min - ry
-            x0max = self._x0max + rx
-            y0max = self._y0max + ry
+        else:
+            rx, ry = (0, 0)
+
+        x0min = self._x0min - rx
+        y0min = self._y0min - ry
+        x0max = self._x0max + rx
+        y0max = self._y0max + ry
 
         ymin, ymax = self.m.ax.projection.y_limits
         xmin, xmax = self.m.ax.projection.x_limits
@@ -746,10 +753,10 @@ def _get_val_from_index(self, ind):
         """
 
         # TODO
-        # Note datashader transposes the data by itself!
+        # Note datashader transposes the data by itself if 1D coords are provided!
         # (to pick the correct value, we need to pick the transposed one!)
 
-        if self.m.shape.name == "shade_raster":
+        if self.m.shape.name == "shade_raster" and self.x0_1D is not None:
             val = self.z_data.T.flat[ind]
         else:
             val = self.z_data.flat[ind]

eomaps/_shapes.py

@@ -407,11 +407,13 @@ def _get_geod_circle_points(self, x, y, crs, radius, n=20):
 
             # transform from in-crs to lon/lat
             radius_t = shapes._get_transformer(
-                self._m.get_crs(crs), CRS.from_epsg(4326)
+                self._m.get_crs(crs),
+                self._m.CRS.PlateCarree(globe=self._m.crs_plot.globe),
             )
             # transform from lon/lat to the plot_crs
             plot_t = shapes._get_transformer(
-                CRS.from_epsg(4326), CRS.from_user_input(self._m.crs_plot)
+                self._m.CRS.PlateCarree(globe=self._m.crs_plot.globe),
+                CRS.from_user_input(self._m.crs_plot),
             )
 
             lon, lat = radius_t.transform(x, y)

eomaps/_version.py

@@ -1 +1 @@
-__version__ = "6.2"
+__version__ = "6.2.1"

eomaps/_webmap.py

@@ -651,7 +651,7 @@ def __init__(
         name : str
             The name of the API.
         service_type : str, optional
-            the service-type to use ("wms" or "wmts"). The default is "wmts".
+            the service-type to use ("wms", "wmts" or "xyz"). The default is "wmts".
         layers : set, optional
             A set of default layers used for delayed fetching (and autocompletion...)
             As soon as one of the layers is accessed, the API is fetched and the
@@ -668,6 +668,8 @@ def __init__(
         self._service_type = service_type
         self._params = _params
         self._fetched = False
+        self._REST_API = None
+
         if layers is None:
             layers = set()
         self._layers = layers
@@ -695,35 +697,35 @@ def _fetch_services(self):
         # set _fetched to True immediately to avoid issues in __getattribute__
         self._fetched = True
 
-        print(f"EOmaps: ... fetching services for '{self._name}'")
-
-        self._REST_API = _REST_API(self._REST_url, _params=self._params)
-
-        found_folders = set()
-        for foldername, services in self._REST_API._structure.items():
-            setattr(
-                self,
-                foldername,
-                _multi_REST_WMSservice(
-                    m=self._m,
-                    services=services,
-                    service_type=self._service_type,
-                    url=self._REST_url,
-                ),
-            )
-            found_folders.add(foldername)
+        if self._REST_API is None:
+            print(f"EOmaps: ... fetching services for '{self._name}'")
+            self._REST_API = _REST_API(self._REST_url, _params=self._params)
+
+            found_folders = set()
+            for foldername, services in self._REST_API._structure.items():
+                setattr(
+                    self,
+                    foldername,
+                    _multi_REST_WMSservice(
+                        m=self._m,
+                        services=services,
+                        service_type=self._service_type,
+                        url=self._REST_url,
+                    ),
+                )
+                found_folders.add(foldername)
 
-        new_layers = found_folders - self._layers
-        if len(new_layers) > 0:
-            print(f"EOmaps: ... found some new folders: {new_layers}")
+            new_layers = found_folders - self._layers
+            if len(new_layers) > 0:
+                print(f"EOmaps: ... found some new folders: {new_layers}")
 
-        invalid_layers = self._layers - found_folders
-        if len(invalid_layers) > 0:
-            print(f"EOmaps: ... could not find the folders: {invalid_layers}")
-        for i in invalid_layers:
-            delattr(self, i)
+            invalid_layers = self._layers - found_folders
+            if len(invalid_layers) > 0:
+                print(f"EOmaps: ... could not find the folders: {invalid_layers}")
+            for i in invalid_layers:
+                delattr(self, i)
 
-        print("EOmaps: done!")
+            print("EOmaps: done!")
 
 
 class _REST_WMSservice(_WebServiec_collection):
@@ -754,6 +756,13 @@ def _url(self):
                 url = WMSurl
             else:
                 url = None
+        elif self._service_type == "xyz":
+            suffix = "/tile/{z}/{y}/{x}"
+            WMSurl = url + suffix
+            if requests.get(url + "/tile/{2}/{1}/{1}").status_code == 200:
+                url = WMSurl
+            else:
+                url = None
         return url
 
     def _fetch_layers(self):
@@ -774,6 +783,10 @@ def _fetch_layers(self):
                     self._layers["layer_" + _sanitize(lname)] = _wmts_layer(
                         self._m, wmts, lname
                     )
+            elif self._service_type == "xyz":
+                self._layers["xyz_layer"] = _xyz_tile_service(
+                    self._m, url, 19, "xyz_layer"
+                )
 
     @property
     @lru_cache()
@@ -917,7 +930,10 @@ def _getz(d, zmax):
         # see https://wiki.openstreetmap.org/wiki/Zoom_levels
         # see https://stackoverflow.com/a/75251360/9703451
 
-        z = int(np.clip(np.ceil(np.log2(1 / d * 40075016.68557849)), 0, zmax))
+        equatorial_circumfence = 40075016.68557849
+        # (e.g. self._crs.globe.semiminor_axis * np.pi * 2)
+
+        z = int(np.clip(np.ceil(np.log2(1 / d * equatorial_circumfence)), 0, zmax))
         return z
 
     def getz(self, extent, target_resolution, zmax):
@@ -996,7 +1012,7 @@ def _image_and_extent(
 
             # reproject the extent to the output-crs
 
-            target_extent = ogc_clients._target_extents(extent, self._crs, output_proj)
+            target_extent = ogc_clients._target_extents(extent, wms_proj, output_proj)
             if len(target_extent) > 0:
                 target_extent = target_extent[0]
             else:
@@ -1015,7 +1031,7 @@ def _image_and_extent(
             original_extent = extent
             img, extent = warp_array(
                 img,
-                source_proj=self._crs,
+                source_proj=wms_proj,
                 source_extent=original_extent,
                 target_proj=output_proj,
                 target_res=regrid_shape,
@@ -1175,6 +1191,16 @@ def _do_add_layer(self, m, layer, **kwargs):
         m.BM.add_bg_artist(self._artist, layer=layer)
 
 
+class _xyz_tile_service_nonearth(_xyz_tile_service):
+    def __call__(self, *args, **kwargs):
+        print(
+            f"EOmaps WARNING: The WebMap service '{self.name}' shows images from a "
+            "different celestrial body projected to an earth-based crs! "
+            "Units used in scalebars, geod_crices etc. represent earth-based units!"
+        )
+        super().__call__(*args, **kwargs)
+
+
 # ------------------------------------------------------------------------------
 # The following is very much copied from cartopy.mpl.slippy_image_artist.py
 # https://github.com/SciTools/cartopy/blob/main/lib/cartopy/mpl/slippy_image_artist.py