diff --git a/doc/examples/notebooks/jupyter_tutorial.ipynb b/doc/examples/notebooks/jupyter_tutorial.ipynb index be9d9482b95..65f38308b84 100644 --- a/doc/examples/notebooks/jupyter_tutorial.ipynb +++ b/doc/examples/notebooks/jupyter_tutorial.ipynb @@ -221,11 +221,12 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "**Only when using ipyleaflet**, the map display includes a button for querying data ( ℹ ), a button for displaying and editing computational region (□), and a button for drawing simple vector geometries (🖊) that can be saved as a GRASS native vector map.\n", + "**Only when using ipyleaflet**, the map display includes a button for querying data ( ℹ ), a button for displaying and editing computational region (□), a button for drawing simple vector geometries (🖊) that can be saved as a GRASS native vector map, and a button for plotting raster values along a line (📈).\n", "\n", "1. Query raster/vector layers: Click the **info** button ( ℹ ) to enable query mode, then select a point on the map to retrieve information. Toggle the button off when finished.\n", "2. View and edit the computational region: Click the **computational region** button (□) to display and optionally adjust the region. Update it by clicking **Update region**, then toggle the button off when done.\n", - "3. Draw and save geometries: Click the **pencil** button (🖊) to draw shapes on the map, name the vector map, and save it. The geometry will be added as a new layer and the drawing tool will close automatically.\n" + "3. Draw and save geometries: Click the **pencil** button (🖊) to draw shapes on the map, name the vector map, and save it. The geometry will be added as a new layer and the drawing tool will close automatically.\n", + "4. Profile rasters along a line: Click the **line chart** button (📈), select the rasters, and draw a line. The profile is plotted below the map and is available as the `profile` attribute of the map, e.g., `raleigh_map.profile.data` or `raleigh_map.profile.save(\"profile.png\")`.\n" ] }, { @@ -259,6 +260,25 @@ "raleigh_map.show()" ] }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "A profile can also be computed without drawing, by giving the line coordinates in the project's coordinate reference system. The figure shows the line with numbered vertices on the raster and the values along the line." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "profile = gj.Profile(\n", + " \"elevation\", [(635000, 222000), (640000, 222000), (640000, 225000)]\n", + ")\n", + "profile.show()" + ] + }, { "cell_type": "markdown", "metadata": {}, @@ -557,7 +577,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.10.12" + "version": "3.12.14" } }, "nbformat": 4, diff --git a/python/grass/jupyter/Makefile b/python/grass/jupyter/Makefile index cb39436a3e6..6d6111cac4f 100644 --- a/python/grass/jupyter/Makefile +++ b/python/grass/jupyter/Makefile @@ -9,6 +9,7 @@ MODULES = \ setup \ map \ interactivemap \ + profile \ region \ map3d \ seriesmap \ diff --git a/python/grass/jupyter/__init__.py b/python/grass/jupyter/__init__.py index 737c401da47..75150b5a9d5 100644 --- a/python/grass/jupyter/__init__.py +++ b/python/grass/jupyter/__init__.py @@ -64,6 +64,15 @@ >>> streams_map.d_vect(map="streams") >>> streams_map.show() +Values of rasters along a line can be sampled and plotted with the line +vertices marked on a map and on the chart:: + +>>> profile = gj.Profile("elevation", [(638000, 220000), (640000, 222000)]) +>>> profile.save("profile.png") + +The interactive map has a tool for drawing the line; the profile of the selected +rasters is then shown below the map and available as the *profile* attribute. + Other classes and functions are described below and in the example notebooks. Static HTML versions of the example notebooks are available on GitHub_ and interactive ones with live code are available on Binder: @@ -107,6 +116,7 @@ from .interactivemap import InteractiveMap, Raster, Vector from .map import Map from .map3d import Map3D +from .profile import Profile from .seriesmap import SeriesMap from .setup import init from .timeseriesmap import TimeSeriesMap @@ -115,6 +125,7 @@ "InteractiveMap", "Map", "Map3D", + "Profile", "Raster", "SeriesMap", "TimeSeriesMap", diff --git a/python/grass/jupyter/interactivemap.py b/python/grass/jupyter/interactivemap.py index 0813dc59646..50aa65634d3 100644 --- a/python/grass/jupyter/interactivemap.py +++ b/python/grass/jupyter/interactivemap.py @@ -16,6 +16,7 @@ import base64 import json from pathlib import Path +from .profile import Profile from .reprojection_renderer import ReprojectionRenderer from .utils import ( @@ -28,6 +29,7 @@ query_raster, query_vector, reproject_latlon, + reproject_latlon_coords, ) @@ -302,6 +304,8 @@ def _import_ipyleaflet(error): self.width = int(width) self.height = int(height) self._controllers = {} + self._profile_controller = None + self._bottom_output = None # Store vector and raster name self.raster_name = [] @@ -318,6 +322,8 @@ def _import_ipyleaflet(error): self.map = self._ipyleaflet.Map( basemap=basemap, layout=layout, scroll_wheel_zoom=True ) + # Output shown below the map, used by the profile tool for its figure. + self._bottom_output = self._ipywidgets.Output() else: self.map = self._folium.Map( @@ -419,6 +425,28 @@ def setup_query_interface(self): controller_class=InteractiveQueryController, ) + def setup_profile_interface(self): + """Sets up the raster profile interface. + + This includes creating a toggle button to activate the profile mode, + and instantiating an InteractiveProfileController which lets the user + draw a line and plots the values of the selected rasters along it. + """ + button = self._create_toggle_button( + icon="line-chart", + tooltip=_("Click to draw a line and plot raster values along it"), + controller_class=InteractiveProfileController, + ) + self._profile_controller = self._controllers[button] + return button + + @property + def profile(self): + """Profile of the last line drawn with the profile tool, None if no line""" + if self._profile_controller is None: + return None + return self._profile_controller.profile + def _create_toggle_button(self, icon, tooltip, controller_class): button = self._ipywidgets.ToggleButton( icon=icon, @@ -437,6 +465,7 @@ def _create_toggle_button(self, icon, tooltip, controller_class): rasters=self.raster_name, vectors=self.vector_name, width=self.width, + output_widget=self._bottom_output, ) self._controllers[button] = controller button.observe(self._toggle_mode, names="value") @@ -453,8 +482,9 @@ def _toggle_mode(self, change): self._controllers[change["owner"]].deactivate() def show(self): - """This function returns a folium figure or ipyleaflet map object - with a GRASS raster and/or vector overlaid on a basemap. + """This function returns a folium figure or an ipywidgets box + with the ipyleaflet map and an output area below it, with a GRASS + raster and/or vector overlaid on a basemap. If map has layer control enabled, additional layers cannot be added after calling show().""" @@ -463,9 +493,10 @@ def show(self): self.setup_query_interface(), self.setup_computational_region_interface(), self.setup_drawing_interface(), + self.setup_profile_interface(), ] button_box = self._ipywidgets.HBox( - toggle_buttons, layout=self._ipywidgets.Layout(width="150px") + toggle_buttons, layout=self._ipywidgets.Layout(width="200px") ) self.map.add( self._ipyleaflet.WidgetControl(widget=button_box, position="topright") @@ -486,7 +517,7 @@ def show(self): # ipyleaflet self.map.add(self.layer_control_object) - return self.map + return self._ipywidgets.VBox([self.map, self._bottom_output]) def save(self, filename): """Save map as an html map. @@ -846,3 +877,104 @@ def clear_popups(self): for item in reversed(list(self.map.layers)): if isinstance(item, self._ipyleaflet.Popup): self.map.remove(item) + + +class InteractiveProfileController: + """A controller for drawing a line and plotting raster values along it. + + Attributes: + map: The ipyleaflet.Map object. + _ipyleaflet: The ipyleaflet module. + _ipywidgets: The ipywidgets module. + raster_name: List of raster layers on the map. + output_widget: The ipywidgets.Output where the profile figure is shown. + draw_control: The draw control restricted to lines. + raster_select: The widget for selecting the rasters to profile. + select_control: The widget control holding raster_select. + profile: The Profile of the last drawn line, None if no line. + """ + + def __init__( + self, map_object, ipyleaflet, ipywidgets, rasters, output_widget, **kwargs + ): # pylint: disable=unused-argument + """Initializes the InteractiveProfileController. + + :param map_object: The map object. + :param ipyleaflet: The ipyleaflet module. + :param ipywidgets: The ipywidgets module. + :param rasters: List of raster layers. + :param output_widget: The ipywidgets.Output where the profile is shown. + """ + self.map = map_object + self._ipyleaflet = ipyleaflet + self._ipywidgets = ipywidgets + self.raster_name = rasters + self.output_widget = output_widget + self.profile = None + # An empty dict disables a shape; only a line makes sense for a profile. + self.draw_control = self._ipyleaflet.DrawControl( + polygon={}, + circlemarker={}, + polyline={"shapeOptions": {"color": "#d62728", "weight": 3}}, + edit=False, + remove=False, + ) + # Registered once, so toggling the tool does not stack callbacks. + self.draw_control.on_draw(self._handle_draw) + self.raster_select = self._ipywidgets.SelectMultiple( + description=_("Rasters:"), + rows=4, + style={"description_width": "initial"}, + ) + self.select_control = None + + def activate(self): + """Activates line drawing and shows the raster selection.""" + self.raster_select.options = list(self.raster_name) + self.raster_select.value = tuple(self.raster_name) + self.map.add(self.draw_control) + self.select_control = self._ipyleaflet.WidgetControl( + widget=self.raster_select, position="topright" + ) + self.map.add(self.select_control) + + def deactivate(self): + """Removes the drawn lines, the draw control, and the raster selection.""" + self.draw_control.clear() + if self.draw_control in self.map.controls: + self.map.remove(self.draw_control) + if self.select_control in self.map.controls: + self.map.remove(self.select_control) + self.select_control = None + + def _handle_draw(self, _control, action, geo_json): + """Computes and shows the profile of a newly drawn line. + + :param str action: The action type. + :param dict geo_json: The GeoJSON data. + """ + if action != "created" or geo_json["geometry"]["type"] != "LineString": + return + # Inside the output context, so errors show up in the notebook. + with self.output_widget: + self.output_widget.clear_output(wait=True) + rasters = list(self.raster_select.value) + if not rasters: + print(_("Select at least one raster to profile.")) + return + self.profile = self.create_profile( + geo_json["geometry"]["coordinates"], rasters + ) + self.profile.show() + + def create_profile(self, lonlat_coordinates, rasters): + """Creates the profile of rasters along a line. + + :param lonlat_coordinates: GeoJSON vertices as [longitude, latitude] pairs + :param rasters: List of raster names. + :return: Profile object + """ + points = reproject_latlon_coords( + [(lat, lon) for lon, lat in lonlat_coordinates] + ) + return Profile(rasters, [(east, north) for east, north, _elevation in points]) diff --git a/python/grass/jupyter/profile.py b/python/grass/jupyter/profile.py new file mode 100644 index 00000000000..55d0c41c528 --- /dev/null +++ b/python/grass/jupyter/profile.py @@ -0,0 +1,279 @@ +# MODULE: grass.jupyter.profile +# +# AUTHOR(S): Corey White +# +# PURPOSE: This module contains functions for sampling and plotting +# raster values along a line in Jupyter Notebooks +# +# SPDX-FileCopyrightText: 2026 Corey White +# SPDX-FileCopyrightText: GRASS Development Team +# SPDX-License-Identifier: GPL-2.0-or-later + +"""Raster profiles along a line""" + +import bisect +import io +import os + +from grass.tools import Tools + +from .map import Map +from .utils import get_region + + +class Profile: + """Values of one or more rasters sampled along a line. + + The line is given as coordinate pairs in the coordinate reference system + of the current project. Values are sampled with *r.profile* and the + distances of the line vertices are measured with *m.measure*. + + :Basic usage: + .. code-block:: pycon + + >>> profile = Profile("elevation", [(638000, 220000), (640000, 222000)]) + >>> profile.data["elevation"][0] + {'easting': 638000.0, 'northing': 220000.0, 'distance': 0, 'value': 109.5} + >>> profile.save("profile.png") + """ + + def __init__( + self, + rasters, + coordinates, + *, + env=None, + session=None, + resolution=None, + units=None, + ): + """Samples the rasters along the line. + + :param rasters: name of a raster or a list of raster names + :param coordinates: list of (east, north) tuples in the project CRS + :param env: runtime environment to use for execution (defaults to global) + :param session: session with environment (used if it has an env attribute) + :param resolution: sampling step along the line in map units + (defaults to the mean of the region resolutions) + :param units: distance units (meters, kilometers, feet, miles), + defaults to the project units + """ + if isinstance(rasters, str): + rasters = [rasters] + self.rasters = list(rasters) + if not self.rasters: + raise ValueError(_("At least one raster is required")) + self.coordinates = [(float(east), float(north)) for east, north in coordinates] + if len(self.coordinates) < 2: + raise ValueError(_("At least two coordinate pairs are required")) + + if env: + self._env = env.copy() + elif session and hasattr(session, "env"): + self._env = session.env.copy() + else: + self._env = os.environ.copy() + + tools = Tools(env=self._env) + units_kwargs = {"units": units} if units else {} + self.vertices, self.distance_units = self._measure_vertices(tools, units_kwargs) + + profile_kwargs = dict(units_kwargs) + if resolution is not None: + profile_kwargs["resolution"] = resolution + flat_coordinates = [value for pair in self.coordinates for value in pair] + self.data = { + raster: tools.r_profile( + input=raster, + coordinates=flat_coordinates, + flags="g", + format="json", + **profile_kwargs, + ).json + for raster in self.rasters + } + self._append_last_vertex(tools) + + def _measure_vertices(self, tools, units_kwargs): + """Returns the vertices with cumulative distances and the distance units.""" + east, north = self.coordinates[0] + vertices = [{"easting": east, "northing": north, "distance": 0}] + for east, north in self.coordinates[1:]: + previous = vertices[-1] + result = tools.m_measure( + coordinates=[previous["easting"], previous["northing"], east, north], + format="json", + **units_kwargs, + ).json + vertices.append( + { + "easting": east, + "northing": north, + "distance": previous["distance"] + result["length"], + } + ) + return vertices, result["units"]["length"] + + def _append_last_vertex(self, tools): + """Adds a sample at the last vertex, which r.profile does not output.""" + east, north = self.coordinates[-1] + values = tools.r_what( + map=self.rasters, coordinates=[east, north], format="json" + ).json[0] + distance = self.vertices[-1]["distance"] + for raster in self.rasters: + self.data[raster].append( + { + "easting": east, + "northing": north, + "distance": distance, + "value": values[raster]["value"], + } + ) + + def to_dataframe(self): + """Returns the samples as a pandas DataFrame with a raster column. + + Requires pandas. + """ + import pandas as pd # pylint: disable=import-outside-toplevel + + return pd.DataFrame( + [ + {"raster": raster, **row} + for raster, rows in self.data.items() + for row in rows + ] + ) + + def plot(self, ax=None, markers=True): + """Plots the profile of each raster as distance against value. + + :param ax: matplotlib axes to draw on, a new figure is created if None + :param bool markers: draw numbered markers at the line vertices + :return: the matplotlib axes + """ + from matplotlib.figure import Figure # pylint: disable=import-outside-toplevel + + if ax is None: + ax = Figure().subplots() + for raster, rows in self.data.items(): + distances = [row["distance"] for row in rows] + # None becomes NaN, which matplotlib leaves as a gap in the line. + values = [ + float("nan") if row["value"] is None else row["value"] for row in rows + ] + (line,) = ax.plot(distances, values, linewidth=2, label=raster) + if markers: + for number, vertex in enumerate(self.vertices, start=1): + value = self._interpolate(rows, vertex["distance"]) + if value is not None: + self._draw_marker( + ax, vertex["distance"], value, number, line.get_color() + ) + if len(self.data) > 1: + ax.legend() + ax.set_xlabel(_("Distance [{}]").format(self.distance_units)) + ax.set_ylabel(_("Value")) + return ax + + @staticmethod + def _interpolate(rows, distance): + """Returns the value at a distance by linear interpolation between samples. + + Returns None when a neighboring sample is null so that no marker + is drawn in a gap of the line. + """ + distances = [row["distance"] for row in rows] + index = bisect.bisect_right(distances, distance) - 1 + if index < 0: + return None + before = rows[index] + if index + 1 >= len(rows) or before["distance"] == distance: + return before["value"] + after = rows[index + 1] + if before["value"] is None or after["value"] is None: + return None + fraction = (distance - before["distance"]) / ( + after["distance"] - before["distance"] + ) + return before["value"] + fraction * (after["value"] - before["value"]) + + @staticmethod + def _draw_marker(ax, x, y, number, color): + """Draws a numbered marker.""" + ax.plot( + x, + y, + marker="o", + markersize=8, + markerfacecolor="white", + markeredgecolor=color, + linestyle="none", + ) + ax.annotate( + str(number), + (x, y), + textcoords="offset points", + xytext=(0, 7), + ha="center", + fontsize=8, + ) + + def _plot_map(self, ax, raster): + """Draws the raster with the line and numbered vertices.""" + import matplotlib.image as mpimg # pylint: disable=import-outside-toplevel + + # The current region is what r.profile sampled, so render exactly that. + region = get_region(env=self._env) + rendered = Map(env=self._env, use_region=True) + rendered.d_rast(map=raster) + image = mpimg.imread(rendered.filename) + ax.imshow( + image, + extent=[region["west"], region["east"], region["south"], region["north"]], + ) + eastings, northings = zip(*self.coordinates, strict=False) + ax.plot(eastings, northings, color="black", linewidth=2) + for number, (east, north) in enumerate(self.coordinates, start=1): + self._draw_marker(ax, east, north, number, "black") + ax.set_xticks([]) + ax.set_yticks([]) + ax.set_title(raster) + + def figure(self, raster=None, figsize=(10, 4)): + """Creates a figure with a map of the line and the profile chart. + + :param raster: raster shown in the map panel (defaults to the first one) + :param figsize: figure size in inches as (width, height) + :return: matplotlib Figure + """ + from matplotlib.figure import Figure # pylint: disable=import-outside-toplevel + + fig = Figure(figsize=figsize, layout="constrained") + map_ax, chart_ax = fig.subplots(1, 2, width_ratios=[1, 1.6]) + self._plot_map(map_ax, raster or self.rasters[0]) + self.plot(ax=chart_ax) + return fig + + def show(self, **kwargs): + """Displays the figure in the notebook. + + :param kwargs: arguments passed to :meth:`figure` + """ + from IPython.display import ( # pylint: disable=import-outside-toplevel + Image, + display, + ) + + buffer = io.BytesIO() + self.figure(**kwargs).savefig(buffer, format="png") + display(Image(buffer.getvalue())) + + def save(self, filename, **kwargs): + """Saves the figure to a file, the format follows the extension. + + :param filename: name of the image file + :param kwargs: arguments passed to :meth:`figure` + """ + self.figure(**kwargs).savefig(filename) diff --git a/python/grass/jupyter/tests/profile_controller_test.py b/python/grass/jupyter/tests/profile_controller_test.py new file mode 100644 index 00000000000..f56f8bcfaba --- /dev/null +++ b/python/grass/jupyter/tests/profile_controller_test.py @@ -0,0 +1,69 @@ +"""Test the profile tool of grass.jupyter.InteractiveMap""" + +from io import StringIO + +import pytest + +from grass.tools import Tools + +import grass.jupyter as gj + +ipyleaflet = pytest.importorskip( + "ipyleaflet", reason="ipyleaflet package not available" +) +ipywidgets = pytest.importorskip( + "ipywidgets", reason="ipywidgets package not available" +) + + +def lonlat_line(coordinates): + """Reproject (east, north) pairs of the current project to GeoJSON [lon, lat]""" + text = "".join(f"{east} {north}\n" for east, north in coordinates) + output = Tools().m_proj(input=StringIO(text), flags="od", separator=",").text + return [ + [float(value) for value in line.split(",")[:2]] for line in output.splitlines() + ] + + +@pytest.mark.needs_solo_run +def test_profile_tool(simple_dataset): + """Drawing a line computes the profile of the selected rasters""" + interactive_map = gj.InteractiveMap(map_backend="ipyleaflet") + interactive_map.add_raster(simple_dataset.raster_name) + assert interactive_map.profile is None + button = interactive_map.setup_profile_interface() + controller = interactive_map._controllers[button] + + button.value = True + assert controller.raster_select.value == (simple_dataset.raster_name,) + assert controller.draw_control in interactive_map.map.controls + assert controller.select_control in interactive_map.map.controls + + coordinates = [(10, 10), (110, 70)] + geo_json = { + "type": "Feature", + "geometry": {"type": "LineString", "coordinates": lonlat_line(coordinates)}, + } + controller._handle_draw(None, "created", geo_json) + profile = interactive_map.profile + assert profile is not None + assert profile.rasters == [simple_dataset.raster_name] + assert len(profile.vertices) == 2 + assert profile.vertices[0]["easting"] == pytest.approx(10, abs=0.01) + assert profile.vertices[-1]["northing"] == pytest.approx(70, abs=0.01) + assert len(profile.data[simple_dataset.raster_name]) > 2 + + button.value = False + assert controller.draw_control not in interactive_map.map.controls + assert controller.select_control is None + + +@pytest.mark.needs_solo_run +def test_show_returns_map_with_output(simple_dataset): + """The map and the output area for the profile are shown together""" + interactive_map = gj.InteractiveMap(map_backend="ipyleaflet") + interactive_map.add_raster(simple_dataset.raster_name) + box = interactive_map.show() + assert isinstance(box, ipywidgets.VBox) + assert box.children[0] is interactive_map.map + assert isinstance(box.children[1], ipywidgets.Output) diff --git a/python/grass/jupyter/tests/profile_test.py b/python/grass/jupyter/tests/profile_test.py new file mode 100644 index 00000000000..158f91db8f7 --- /dev/null +++ b/python/grass/jupyter/tests/profile_test.py @@ -0,0 +1,101 @@ +"""Test grass.jupyter.Profile""" + +import pytest + +from grass.tools import Tools + +import grass.jupyter as gj + +LINE = [(0.5, 0.5), (1.5, 0.5), (1.5, 4.5)] + + +def test_vertices(session_with_data): + """Every coordinate pair becomes a vertex with a cumulative distance""" + profile = gj.Profile("data", LINE, env=session_with_data.env) + assert len(profile.vertices) == len(LINE) + assert [vertex["distance"] for vertex in profile.vertices] == pytest.approx( + [0, 1, 5] + ) + assert (profile.vertices[-1]["easting"], profile.vertices[-1]["northing"]) == ( + 1.5, + 4.5, + ) + + +def test_data_rows(session_with_data): + """Samples run from the first to the last vertex with increasing distance""" + profile = gj.Profile("data", LINE, env=session_with_data.env) + rows = profile.data["data"] + assert set(rows[0].keys()) == {"easting", "northing", "distance", "value"} + assert rows[0] == {"easting": 0.5, "northing": 0.5, "distance": 0, "value": 5} + # One sample per cell along 5 units plus the last vertex. + assert len(rows) == 6 + distances = [row["distance"] for row in rows] + assert distances == sorted(distances) + assert rows[-1]["distance"] == pytest.approx(profile.vertices[-1]["distance"]) + # data = row() * col(), so the last vertex in row 1, column 2 has value 2. + assert rows[-1]["value"] == 2 + + +def test_multiple_rasters_aligned(session_with_data): + """Rasters are sampled at the same distances""" + tools = Tools(session=session_with_data) + tools.r_mapcalc(expression="data2 = 2 * data") + profile = gj.Profile(["data", "data2"], LINE, env=session_with_data.env) + assert profile.rasters == ["data", "data2"] + rows = profile.data["data"] + rows2 = profile.data["data2"] + assert [row["distance"] for row in rows] == [row["distance"] for row in rows2] + assert [2 * row["value"] for row in rows] == [row["value"] for row in rows2] + + +def test_null_values(session_with_data): + """Null cells become None and do not break plotting""" + pytest.importorskip("matplotlib", reason="matplotlib package not available") + tools = Tools(session=session_with_data) + tools.r_mapcalc(expression="nulls = if(col() == 2, null(), data)") + profile = gj.Profile("nulls", LINE, env=session_with_data.env) + values = [row["value"] for row in profile.data["nulls"]] + assert values[0] == 5 + assert values[1:] == [None] * 5 + ax = profile.plot() + assert len(ax.get_lines()) >= 1 + + +def test_resolution(session_with_data): + """Halving the resolution doubles the number of samples along the line""" + profile = gj.Profile("data", LINE, env=session_with_data.env) + fine = gj.Profile("data", LINE, env=session_with_data.env, resolution=0.5) + # Both include the extra sample at the last vertex. + assert len(fine.data["data"]) - 1 == 2 * (len(profile.data["data"]) - 1) + + +def test_invalid_input(session_with_data): + """Too few coordinates or no raster are rejected""" + with pytest.raises(ValueError, match="two coordinate pairs"): + gj.Profile("data", [(0.5, 0.5)], env=session_with_data.env) + with pytest.raises(ValueError, match="one raster"): + gj.Profile([], LINE, env=session_with_data.env) + + +def test_to_dataframe(session_with_data): + """DataFrame has one row per sample with the raster name""" + pd = pytest.importorskip("pandas", reason="pandas package not available") + profile = gj.Profile("data", LINE, env=session_with_data.env) + df = profile.to_dataframe() + assert isinstance(df, pd.DataFrame) + assert list(df.columns) == ["raster", "easting", "northing", "distance", "value"] + assert len(df) == len(profile.data["data"]) + assert set(df["raster"]) == {"data"} + + +def test_figure_and_save(session_with_data, tmp_path): + """Figure has a map and a chart panel and can be saved""" + pytest.importorskip("matplotlib", reason="matplotlib package not available") + profile = gj.Profile("data", LINE, env=session_with_data.env) + fig = profile.figure() + assert len(fig.axes) == 2 + filename = tmp_path / "profile.png" + profile.save(filename) + assert filename.exists() + assert filename.stat().st_size > 0 diff --git a/python/grass/jupyter/testsuite/interactivemap_test.py b/python/grass/jupyter/testsuite/interactivemap_test.py index 375ab0d2a38..adc5673906b 100644 --- a/python/grass/jupyter/testsuite/interactivemap_test.py +++ b/python/grass/jupyter/testsuite/interactivemap_test.py @@ -105,6 +105,27 @@ def test_draw(self): interactive_map._controllers[button].activate() self.assertIsNotNone(interactive_map._controllers[button].save_button_control) + @unittest.skipIf(not can_import_ipyleaflet(), "Cannot import ipyleaflet") + def test_profile(self): + """Test the profile tool with a drawn line.""" + interactive_map = gj.InteractiveMap(map_backend="ipyleaflet") + interactive_map.add_raster("elevation") + button = interactive_map.setup_profile_interface() + controller = interactive_map._controllers[button] + controller.activate() + self.assertEqual(controller.raster_select.value, ("elevation",)) + geo_json = { + "type": "Feature", + "geometry": { + "type": "LineString", + "coordinates": [[-78.68, 35.78], [-78.66, 35.79]], + }, + } + controller._handle_draw(None, "created", geo_json) + self.assertIsNotNone(interactive_map.profile) + self.assertEqual(len(interactive_map.profile.vertices), 2) + self.assertIn("elevation", interactive_map.profile.data) + @unittest.skipIf(not can_import_ipyleaflet(), "Cannot import ipyleaflet") def test_draw_computational_region(self): """Test the draw_computational_region method.""" diff --git a/python/grass/jupyter/utils.py b/python/grass/jupyter/utils.py index 210b1546db9..840b9301c4c 100644 --- a/python/grass/jupyter/utils.py +++ b/python/grass/jupyter/utils.py @@ -84,16 +84,25 @@ def reproject_region(region, from_proj, to_proj): return region +def reproject_latlon_coords(coords, env=None): + """Reproject coordinates from WGS84 to the project CRS in one m.proj call + + :param coords: list of (latitude, longitude) tuples + :param env: environment to use (defaults to global) + :return: list of (east, north, elevation) tuples in the input order + """ + text = "".join(f"{lon} {lat}\n" for lat, lon in coords) + output = Tools(env=env).m_proj(input=StringIO(text), flags="i", separator=",").text + return [tuple(map(float, line.split(","))) for line in output.splitlines()] + + def reproject_latlon(coord): """Reproject coordinates :param coord: coordinates given as tuple (latitude, longitude) :return: reprojected coordinates (returned as tuple) """ - coord_str = f"{coord[1]} {coord[0]}\n" - output = Tools().m_proj(input=StringIO(coord_str), flags="i", separator=",").text - east, north, elev = map(float, output.splitlines()[0].split(",")) - return east, north, elev + return reproject_latlon_coords([coord])[0] def _style_table(html_content):