grass.jupyter: Use InteractiveMap with custom projections

python/grass/jupyter/interactivemap.py

@@ -16,7 +16,7 @@
 import base64
 import json
 from pathlib import Path
-from .reprojection_renderer import ReprojectionRenderer
+from .reprojection_renderer import ReprojectionRenderer, CustomRenderer
 
 from .utils import (
     get_region_bounds_latlon,
@@ -213,6 +213,7 @@
         use_region=False,
         saved_region=None,
         map_backend=None,
+        crs=None,
     ):
         """Creates a blank folium/ipyleaflet map centered on g.region.
 
@@ -248,6 +249,7 @@
         :param bool use_region: use computational region of current mapset
         :param str saved_region: name of saved computation region
         :param str map_backend: "ipyleaflet" or "folium" or None
+        :param str crs: custom CRS to use for rendering (default None)
         """
         self._ipyleaflet = None
         self._folium = None
@@ -305,14 +307,41 @@
         # Store Region
         self.region = None
 
+        if crs:
+            if crs.lower().startswith("epsg"):
+                crs_code = crs[4:]
+                crs_dict = {"name": crs_code, "custom": False}
+            else:
+                # Assume custom CRS, get current CRS from g.proj
+                import subprocess
+
+                result = subprocess.run(
+                    ["g.proj", "-jf"],
+                    capture_output=True,
+                    text=True,
+                    check=True,
+                )
+                crs_proj4 = result.stdout.strip()
+                crs_dict = {"name": "Custom", "custom": True, "proj4def": crs_proj4}
+            self._renderer = CustomRenderer(crs=crs_dict)
+        else:
+            self._renderer = ReprojectionRenderer(
+                use_region=use_region, saved_region=saved_region
+            )
+
         if self._ipyleaflet:
             basemap = xyzservices.providers.query_name(tiles)
             if API_key and basemap.get("accessToken"):
                 basemap["accessToken"] = API_key
             layout = self._ipywidgets.Layout(width=f"{width}px", height=f"{height}px")
-            self.map = self._ipyleaflet.Map(
-                basemap=basemap, layout=layout, scroll_wheel_zoom=True
-            )
+            if crs:
+                self.map = self._ipyleaflet.Map(
+                    basemap=basemap, layout=layout, scroll_wheel_zoom=True, crs=crs_dict
+                )
+            else:
+                self.map = self._ipyleaflet.Map(
+                    basemap=basemap, layout=layout, scroll_wheel_zoom=True
+                )
 
         else:
             self.map = self._folium.Map(
@@ -325,10 +354,6 @@
         self.layer_control_object = None
         self.region_rectangle = None
 
-        self._renderer = ReprojectionRenderer(
-            use_region=use_region, saved_region=saved_region
-        )
-
     def add_vector(self, name, title=None, **kwargs):
         """Imports vector into temporary WGS84 location, re-formats to a GeoJSON and
         adds to map.

python/grass/jupyter/reprojection_renderer.py

@@ -18,6 +18,7 @@
 import tempfile
 import weakref
 from pathlib import Path
+
 import grass.script as gs
 from .map import Map
 from .utils import (
@@ -168,3 +169,145 @@ def render_vector(self, name):
         )
 
         return json_file
+
+
+class CustomRenderer:
+    """Class for handling layers with a custom CRS without reprojection."""
+
+    def __init__(self, crs, use_region=False, saved_region=None, work_dir=None):
+        """
+        Initialize CustomRenderer.
+
+        param str crs: CRS string (EPSG code or Proj4 definition)
+        param str work_dir: Path to directory where files should be written
+        """
+        self.crs = crs
+
+        # Temporary folder for all our files
+        if not work_dir:
+            # Resource managed by weakref.finalize.
+            self._tmp_dir = (
+                # pylint: disable=consider-using-with
+                tempfile.TemporaryDirectory()
+            )
+
+            def cleanup(tmpdir):
+                tmpdir.cleanup()
+
+            weakref.finalize(self, cleanup, self._tmp_dir)
+        else:
+            self._tmp_dir = work_dir
+
+        # Set up the environment and location
+        self._src_env = os.environ.copy()
+        if not crs.get("custom"):
+            self._rcfile_custom, self._custom_env = setup_location(
+                "custom_location", self._tmp_dir.name, self._src_env, epsg=crs["name"]
+            )
+        else:
+            # If custom, pass the Proj4 string
+            self._rcfile_custom, self._custom_env = setup_location(
+                "custom_location",
+                self._tmp_dir.name,
+                src_env=self._src_env,
+                proj4=crs["proj4def"],
+            )
+
+        self._region_manager = RegionManagerForInteractiveMap(
+            use_region, saved_region, self._src_env, self._custom_env
+        )
+
+    def get_bbox(self):
+        """Return bounding box of computation region in WGS84"""
+        return self._region_manager.bbox
+
+    def render_raster(self, name):
+        """
+        Render a raster layer in its native CRS and save it as a PNG.
+
+        :param name: Name of the raster
+        :type name: str
+        :return: Path to the PNG file and the reprojected bounds
+        :rtype: tuple
+        """
+        file_info = gs.find_file(name, element="cell", env=self._src_env)
+        full_name = file_info["fullname"]
+        name = file_info["name"]
+        mapset = file_info["mapset"]
+
+        # Set region based on the raster's native bounds
+        self._region_manager.set_region_from_raster(full_name)
+
+        # Reproject raster into the custom CRS location
+        env_info = gs.gisenv(env=self._src_env)
+        tgt_name = full_name.replace("@", "_")
+        gs.run_command(
+            "r.proj",
+            input=full_name,
+            output=tgt_name,
+            mapset=mapset,
+            location=env_info["LOCATION_NAME"],
+            dbase=env_info["GISDBASE"],
+            resolution=self._region_manager.resolution,
+            env=self._custom_env,
+        )
+
+        # Wite the raster to a PNG file using Map
+        raster_info = gs.raster_info(tgt_name, env=self._custom_env)
+        filename = os.path.join(self._tmp_dir.name, f"{tgt_name}.png")
+        img = Map(
+            width=int(raster_info["cols"]),
+            height=int(raster_info["rows"]),
+            env=self._custom_env,
+            filename=filename,
+            use_region=True,
+        )
+        img.run("d.rast", map=tgt_name)
+
+        # Reproject the raster's bounds to WGS84 for overlaying the PNG
+        old_bounds = get_region(self._src_env)
+        from_proj = get_location_proj_string(env=self._src_env)
+        to_proj = get_location_proj_string(env=self._custom_env)
+        bounds = reproject_region(old_bounds, from_proj, to_proj)
+        new_bounds = [
+            [bounds["north"], bounds["west"]],
+            [bounds["south"], bounds["east"]],
+        ]
+
+        return filename, new_bounds
+
+    def render_vector(self, name):
+        """Reproject vector to WGS84 and save geoJSON in working directory. Return
+        geoJSON filename.
+
+        param str name: name of vector"""
+        # Find full name of vector
+        file_info = gs.find_file(name, element="vector")
+        full_name = file_info["fullname"]
+        name = file_info["name"]
+        mapset = file_info["mapset"]
+        new_name = full_name.replace("@", "_")
+        # set bbox
+        self._region_manager.set_bbox_vector(full_name)
+        # Reproject vector into WGS84 Location
+        env_info = gs.gisenv(env=self._src_env)
+        gs.run_command(
+            "v.proj",
+            input=name,
+            output=new_name,
+            mapset=mapset,
+            project=env_info["LOCATION_NAME"],
+            dbase=env_info["GISDBASE"],
+            env=self._custom_env,
+        )
+        # Convert to GeoJSON
+        json_file = Path(self._tmp_dir.name) / f"{new_name}.json"
+        gs.run_command(
+            "v.out.ogr",
+            input=new_name,
+            output=json_file,
+            format="GeoJSON",
+            env=self._custom_env,
+        )
+
+        return json_file

python/grass/jupyter/utils.py

@@ -331,7 +331,7 @@ def estimate_resolution(raster, mapset, location, dbase, env):
     return (cell_ew + cell_ns) / 2.0
 
 
-def setup_location(name, path, epsg, src_env):
+def setup_location(name, path, src_env, epsg=None, proj4=None):
     """Setup temporary location with different projection but
     same computational region as source location
 
@@ -346,7 +346,10 @@ def setup_location(name, path, epsg, src_env):
     # Create new environment
     rcfile, new_env = gs.create_environment(path, name, "PERMANENT")
     # Location and mapset
-    gs.create_location(path, name, epsg=epsg, overwrite=True)
+    if proj4:
+        gs.create_location(path, name, proj4=proj4, overwrite=True)
+    elif epsg:
+        gs.create_location(path, name, epsg=epsg, overwrite=True)
     # Reproject region
     set_target_region(src_env, new_env)
     return rcfile, new_env