diff --git a/.github/workflows/python-code-quality.yml b/.github/workflows/python-code-quality.yml index e38fec4d2c5..242068a9f0d 100644 --- a/.github/workflows/python-code-quality.yml +++ b/.github/workflows/python-code-quality.yml @@ -32,7 +32,7 @@ jobs: # renovate: datasource=pypi depName=flake8 FLAKE8_VERSION: "7.3.0" # renovate: datasource=pypi depName=pylint - PYLINT_VERSION: "4.0.6" + PYLINT_VERSION: "4.1.2" # renovate: datasource=pypi depName=bandit BANDIT_VERSION: "1.9.4" # renovate: datasource=pypi depName=ruff diff --git a/.github/workflows/ubuntu.yml b/.github/workflows/ubuntu.yml index a5896aedabc..b12c8ecaec9 100644 --- a/.github/workflows/ubuntu.yml +++ b/.github/workflows/ubuntu.yml @@ -107,6 +107,14 @@ jobs: env: STEPS_GET_EXCLUDE_OUTPUTS_EXTRA_EXCLUDE: ${{ steps.get-exclude.outputs.extra-exclude }} + - name: Set up Python + if: ${{ matrix.config == 'ubuntu-24.04' }} + uses: actions/setup-python@f677139bbe7f9c59b41e40162b753c062f5d49a3 # v5.2.0 + with: + python-version: 3.13 + cache: pip # zizmor: ignore[cache-poisoning] This workflow is not used for release artifacts + allow-prereleases: true + - name: Get dependencies run: | sudo apt-get update -y @@ -116,6 +124,13 @@ jobs: xargs -a <(awk '! /^ *(#|$)/' ".github/workflows/apt.txt") -r -- \ sudo apt-get install -y --no-install-recommends --no-install-suggests + - name: Install Python dependencies + if: ${{ matrix.config == 'ubuntu-24.04' }} + run: | + python -m pip install --upgrade pip + pip install -r .github/workflows/python_requirements.txt + pip install -r .github/workflows/optional_requirements.txt + - name: Create installation directory run: | mkdir "${HOME}/install" diff --git a/Dockerfile b/Dockerfile index 78f6a5f50d3..65936b9304f 100644 --- a/Dockerfile +++ b/Dockerfile @@ -4,7 +4,7 @@ # Changes to this file must be copied over to the other file. ARG GUI=without -FROM ubuntu:26.04@sha256:61ebaa5cc23ca45450db85eac015435199ec569e28ec222ea13f2aed2110b8a6 AS common_start +FROM ubuntu:26.04@sha256:88a381d5b5eeb2b35d3ad70925a362c37ce569daf43ede89ff818ec20e4d3794 AS common_start ARG BASE_NAME="ubuntu:26.04" ARG PYTHON_VERSION=3.14 diff --git a/docker/ubuntu/Dockerfile b/docker/ubuntu/Dockerfile index 78f6a5f50d3..65936b9304f 100644 --- a/docker/ubuntu/Dockerfile +++ b/docker/ubuntu/Dockerfile @@ -4,7 +4,7 @@ # Changes to this file must be copied over to the other file. ARG GUI=without -FROM ubuntu:26.04@sha256:61ebaa5cc23ca45450db85eac015435199ec569e28ec222ea13f2aed2110b8a6 AS common_start +FROM ubuntu:26.04@sha256:88a381d5b5eeb2b35d3ad70925a362c37ce569daf43ede89ff818ec20e4d3794 AS common_start ARG BASE_NAME="ubuntu:26.04" ARG PYTHON_VERSION=3.14 diff --git a/gui/wxpython/core/utils.py b/gui/wxpython/core/utils.py index 305b6c1bc3d..addda95d614 100644 --- a/gui/wxpython/core/utils.py +++ b/gui/wxpython/core/utils.py @@ -230,11 +230,7 @@ def GetValidLayerName(name): cIdx += 1 retName = "".join(retNameList) - if ( - retName - and not (retName[0] >= "A" and retName[0] <= "Z") - and not (retName[0] >= "a" and retName[0] <= "z") - ): + if retName and not ("A" <= retName[0] <= "Z") and not ("a" <= retName[0] <= "z"): retName = "x" + retName[1:] if mapset: diff --git a/gui/wxpython/gui_core/forms.py b/gui/wxpython/gui_core/forms.py index c172b54a8c9..2cdbcce334c 100644 --- a/gui/wxpython/gui_core/forms.py +++ b/gui/wxpython/gui_core/forms.py @@ -3115,7 +3115,7 @@ def OnVerbosity(self, event): def OnPageChange(self, event): sel = self.notebook.GetSelection() if not event else event.GetSelection() idx = self.notebook.GetPageIndexByName("manual") - if idx > -1 and sel == idx: + if -1 < idx == sel: # calling LoadPage() is strangely time-consuming (only first call) # FIXME: move to helpPage.__init__() if not self.manualTab.IsLoaded(): diff --git a/gui/wxpython/history/browser.py b/gui/wxpython/history/browser.py index 98db94cdfc8..6f0817accc4 100644 --- a/gui/wxpython/history/browser.py +++ b/gui/wxpython/history/browser.py @@ -16,6 +16,7 @@ @author Tomas Zigo """ +import math from datetime import datetime import wx @@ -25,7 +26,7 @@ from history.tree import HistoryBrowserTree from icons.icon import MetaIcon -import grass.script as gs +from grass.tools import Tools from grass.grassdb import history @@ -40,10 +41,10 @@ "status": _("Status:"), "mask2d": _("Mask 2D:"), "mask3d": _("Mask 3D:"), - "n": _("North:"), - "s": _("South:"), - "w": _("West:"), - "e": _("East:"), + "north": _("North:"), + "south": _("South:"), + "west": _("West:"), + "east": _("East:"), "nsres": _("North-south resolution:"), "ewres": _("East-west resolution:"), "rows": _("Number of rows:"), @@ -69,6 +70,18 @@ def make_label(key): return TRANSLATION_KEYS.get(key, "") +# Region keys in the shell format of g.region, which are also the names of the g.region parameters +REGION_KEYS = {"n": "north", "s": "south", "w": "west", "e": "east"} + + +def normalize_region(region): + """Return region settings with the keys used by g.region format=json. + + History entries recorded by older versions use the keys of the shell format. + """ + return {REGION_KEYS.get(key, key): value for key, value in region.items()} + + class HistoryInfoPanel(SP.ScrolledPanel): def __init__(self, parent, giface, title=("Command Info"), style=wx.TAB_TRAVERSAL): super().__init__(parent=parent, id=wx.ID_ANY, style=style) @@ -76,6 +89,7 @@ def __init__(self, parent, giface, title=("Command Info"), style=wx.TAB_TRAVERSA self.parent = parent self.giface = giface self.title = title + self.tools = Tools() self.region_settings = None @@ -166,7 +180,7 @@ def _general_info_filter(self, key, value): return key in filter_keys or ((key in {"mask2d", "mask3d"}) and value is True) def _region_settings_filter(self, key): - return key not in {"projection", "zone", "cells"} + return key not in {"crs", "projection", "zone", "cells"} def _updateGeneralInfoBox(self, command_info): """Update a static box for displaying general info about the command. @@ -213,7 +227,7 @@ def _updateRegionSettingsGrid(self, command_info): """ self.sizer_region_settings_grid.Clear(True) - self.region_settings = command_info["region"] + self.region_settings = normalize_region(command_info["region"]) idx = 0 for key, value in self.region_settings.items(): if self._region_settings_filter(key): @@ -248,9 +262,12 @@ def _updateRegionSettingsMatch(self): self.sizer_region_settings_match.Clear(True) # Region condition - history_region = self.region_settings + history_region = self._get_history_region() current_region = self._get_current_region() - region_matches = history_region == current_region + region_matches = all( + math.isclose(value, current_region[key], rel_tol=0, abs_tol=1e-8) + for key, value in history_region.items() + ) # Icon and button according to the condition if region_matches: @@ -320,7 +337,7 @@ def hideCommandInfo(self): def _get_current_region(self): """Get current computational region settings.""" - return gs.region() + return self.tools.g_region(flags="p", format="json").json def _get_history_region(self): """Get computational region settings of executed command.""" @@ -332,8 +349,11 @@ def _get_history_region(self): def OnUpdateRegion(self, event): """Set current region to the region of executed command.""" + parameters = {value: key for key, value in REGION_KEYS.items()} history_region = self._get_history_region() - gs.run_command("g.region", **history_region) + self.tools.g_region( + **{parameters.get(key, key): value for key, value in history_region.items()} + ) self.giface.updateMap.emit(render=False, renderVector=False) self._updateRegionSettingsMatch() self.Layout() diff --git a/gui/wxpython/tplot/frame.py b/gui/wxpython/tplot/frame.py index 9f7f025a50f..5788badd52e 100755 --- a/gui/wxpython/tplot/frame.py +++ b/gui/wxpython/tplot/frame.py @@ -493,7 +493,7 @@ def _getSTRDdata(self, timeseries): r.open() val = r.get_value(self.poi) r.close() - if val == -2147483648 and val < minmin: + if -2147483648 == val < minmin: self.timeDataR[name][row[0]]["value"] = None else: self.timeDataR[name][row[0]]["value"] = val diff --git a/python/grass/grassdb/history.py b/python/grass/grassdb/history.py index 5608f792274..3c3cba5e264 100644 --- a/python/grass/grassdb/history.py +++ b/python/grass/grassdb/history.py @@ -264,7 +264,9 @@ def get_initial_command_info(env_run): mask3d_name = f"RASTER3D_MASK@{env['MAPSET']}" # Computational region settings - region_settings = gs.region(env=env_run) + region_settings = gs.parse_command( + "g.region", flags="p", format="json", env=env_run + ) # Finalize the command info dictionary return { diff --git a/python/grass/pygrass/raster/segment.py b/python/grass/pygrass/raster/segment.py index cc301cc3cef..3f8603f761e 100644 --- a/python/grass/pygrass/raster/segment.py +++ b/python/grass/pygrass/raster/segment.py @@ -33,7 +33,7 @@ def nseg(self): ) def segments_in_mem(self): - if self.maxmem > 0 and self.maxmem < 100: + if 0 < self.maxmem < 100: seg_in_mem = (self.maxmem * self.nseg()) / 100 else: seg_in_mem = 4 * (self.rows() / self.srows + self.cols() / self.scols + 2) diff --git a/python/grass/temporal/c_libraries_interface.py b/python/grass/temporal/c_libraries_interface.py index 3804296b9d3..be9e1b2d400 100644 --- a/python/grass/temporal/c_libraries_interface.py +++ b/python/grass/temporal/c_libraries_interface.py @@ -11,6 +11,7 @@ from __future__ import annotations import logging +import math import sys from ctypes import CFUNCTYPE, POINTER, byref, c_int, c_void_p, cast from datetime import datetime @@ -825,11 +826,11 @@ def _read_raster3d_info(name, mapset): return None libraster3d.Rast3d_range_min_max(g3map, byref(min), byref(max)) - if min.value != min.value: + if math.isnan(min.value): kvp["min"] = None else: kvp["min"] = float(min.value) - if max.value != max.value: + if math.isnan(max.value): kvp["max"] = None else: kvp["max"] = float(max.value) diff --git a/python/grass/temporal/temporal_extent.py b/python/grass/temporal/temporal_extent.py index ef3ab3a8feb..c6ed6d3b5a5 100644 --- a/python/grass/temporal/temporal_extent.py +++ b/python/grass/temporal/temporal_extent.py @@ -737,8 +737,7 @@ def during(self, extent) -> bool: # Check single point of time in interval if self.D["end_time"] is None: return bool( - self.D["start_time"] >= extent.D["start_time"] - and self.D["start_time"] < extent.D["end_time"] + extent.D["start_time"] <= self.D["start_time"] < extent.D["end_time"] ) return bool( @@ -776,8 +775,7 @@ def contains(self, extent) -> bool: # Check single point of time in interval if extent.D["end_time"] is None: return bool( - self.D["start_time"] <= extent.D["start_time"] - and self.D["end_time"] > extent.D["start_time"] + self.D["start_time"] <= extent.D["start_time"] < self.D["end_time"] ) return bool( @@ -853,9 +851,10 @@ def overlaps(self, extent) -> bool: return bool( self.D["end_time"] is not None and extent.D["end_time"] is not None - and self.D["start_time"] < extent.D["start_time"] - and self.D["end_time"] < extent.D["end_time"] - and self.D["end_time"] > extent.D["start_time"] + and self.D["start_time"] + < extent.D["start_time"] + < self.D["end_time"] + < extent.D["end_time"] ) def overlapped(self, extent) -> bool: @@ -893,9 +892,10 @@ def overlapped(self, extent) -> bool: return bool( self.D["end_time"] is not None and extent.D["end_time"] is not None - and self.D["start_time"] > extent.D["start_time"] - and self.D["end_time"] > extent.D["end_time"] - and self.D["start_time"] < extent.D["end_time"] + and extent.D["start_time"] + < self.D["start_time"] + < extent.D["end_time"] + < self.D["end_time"] ) def temporal_relation(self, extent): diff --git a/raster/r.sim/r.sim.water/r.sim.water.html b/raster/r.sim/r.sim.water/r.sim.water.html index 5e280a77876..7e125f78243 100644 --- a/raster/r.sim/r.sim.water/r.sim.water.html +++ b/raster/r.sim/r.sim.water/r.sim.water.html @@ -90,7 +90,15 @@

DESCRIPTION

the time series flag -t and output_step parameter defining the time step in minutes for writing output files. Files are saved with a suffix representing time since the start of simulation in minutes -(e.g. wdepth.05, wdepth.10). +(e.g. wdepth.05, wdepth.10) and are timestamped with that time. +The simulation advances in time steps which usually do not fall exactly on the output times. +A map holds the state at the time step closest to the time in its name, +so at most half a time step earlier or later. +When the time step is longer than output_step, there are fewer time steps than +output times, and a time step writes only the maps for the output time closest to it. +The series always ends with maps named by the duration which hold the state +at the end of the run, also when the duration is not a multiple of output_step +or when the simulation stopped early. Monitoring of water depth at specific points is supported. A vector map with observation points and a path to a logfile must be provided. For each point in the vector map which is located in @@ -209,12 +217,12 @@

Run summary

mean_source_rateMean rainfall excessm/s mean_infiltrationMean infiltration rate, 0 without infiltration inputm/s threadsThreads used for the computationcount - outputsOne entry per set of written maps: one per output_step with -t, otherwise a single entry + outputsOne entry per set of written maps: with -t, one per written output step, the last one named by duration, otherwise a single entry

-Each entry of outputs contains the simulated_time (s) and timestamp -of the written maps, the number of walkers_remaining at that time, and +Each entry of outputs contains the simulated_time (s) when the maps +were written, their timestamp, the number of walkers_remaining at that time, and the names of the depth, discharge, error and walkers maps, or null for maps which were not requested. diff --git a/raster/r.sim/r.sim.water/r.sim.water.md b/raster/r.sim/r.sim.water/r.sim.water.md index 692f0d06c90..ec2fded33af 100644 --- a/raster/r.sim/r.sim.water/r.sim.water.md +++ b/raster/r.sim/r.sim.water/r.sim.water.md @@ -80,7 +80,17 @@ water depth and discharge maps can be saved during simulation using the time series flag **-t** and **output_step** parameter defining the time step in minutes for writing output files. Files are saved with a suffix representing time since the start of simulation in minutes (e.g. -wdepth.05, wdepth.10). Monitoring of water depth at specific points is +wdepth.05, wdepth.10) and are timestamped with that time. The +simulation advances in time steps which usually do not fall exactly on +the output times. A map holds the state at the time step closest to the +time in its name, so at most half a time step earlier or later. When the +time step is longer than **output_step**, there are fewer time steps +than output times, and a time step writes only the maps for the output +time closest to it. The series always ends with maps named by the +**duration** which hold the state at the end of the run, also when the +duration is not a multiple of **output_step** or when the simulation +stopped early. +Monitoring of water depth at specific points is supported. A vector map with observation points and a path to a logfile must be provided. For each point in the vector map which is located in the computational region the water depth is logged each time step in the @@ -197,12 +207,12 @@ stored in the history of the output raster maps under the same keys (see | `mean_source_rate` | Mean rainfall excess | m/s | | `mean_infiltration` | Mean infiltration rate, 0 without infiltration input | m/s | | `threads` | Threads used for the computation | count | -| `outputs` | One entry per set of written maps: one per **output_step** with **-t**, otherwise a single entry | | +| `outputs` | One entry per set of written maps: with **-t**, one per written output step, the last one named by **duration**, otherwise a single entry | | -Each entry of `outputs` contains the `simulated_time` (s) and `timestamp` -of the written maps, the number of `walkers_remaining` at that time, and -the names of the `depth`, `discharge`, `error` and `walkers` maps, or -`null` for maps which were not requested. +Each entry of `outputs` contains the `simulated_time` (s) when the maps +were written, their `timestamp`, the number of `walkers_remaining` at +that time, and the names of the `depth`, `discharge`, `error` and +`walkers` maps, or `null` for maps which were not requested. Summary of a time series run with two output steps in JSON: diff --git a/raster/r.sim/r.sim.water/tests/r_sim_water_nprocs_test.py b/raster/r.sim/r.sim.water/tests/r_sim_water_nprocs_test.py index 8510050fa6c..73f4fe40412 100644 --- a/raster/r.sim/r.sim.water/tests/r_sim_water_nprocs_test.py +++ b/raster/r.sim/r.sim.water/tests/r_sim_water_nprocs_test.py @@ -18,6 +18,8 @@ RAIN = np.where( (abs(row - SIZE // 2) <= SIZE // 4) & (abs(col - SIZE // 2) <= SIZE // 4), 50, 0 ) +# Infiltration only where it rains, at a lower rate than the rain +INFILTRATION = np.where(RAIN > 0, 20, 0) @pytest.fixture(scope="module") @@ -76,3 +78,22 @@ def test_all_walkers_move_with_nprocs(session): single = total_water(depth(session, nprocs=1, **options)) multiple = total_water(depth(session, nprocs=4, **options)) assert multiple == pytest.approx(single, rel=1e-6) + + +def test_infiltration_takes_exactly_the_capacity(session): + """Infiltration removes exactly the capacity of the cells, also in parallel. + + The walkers created in a cell outweigh its infiltration capacity, so the + capacity is used up in the first step, before any walker moves, and it + is not restored. The water left for the rest of the simulation is then + the rain minus the infiltration, so the total is that fraction of the + total without infiltration. + """ + expected = total_water( + depth(session, nprocs=1, infil_value=0, diffusion_coeff=0) + ) * (1 - INFILTRATION.sum() / RAIN.sum()) + for nprocs in (1, 4): + actual = total_water( + depth(session, nprocs=nprocs, infil=INFILTRATION, diffusion_coeff=0) + ) + assert actual == pytest.approx(expected, rel=1e-6) diff --git a/raster/r.sim/r.sim.water/tests/r_sim_water_summary_test.py b/raster/r.sim/r.sim.water/tests/r_sim_water_summary_test.py index a7c3f0fd74f..d62e3b5abc6 100644 --- a/raster/r.sim/r.sim.water/tests/r_sim_water_summary_test.py +++ b/raster/r.sim/r.sim.water/tests/r_sim_water_summary_test.py @@ -120,8 +120,8 @@ def test_json_summary_time_series(session_tools): assert [output["depth"] for output in outputs] == ["depth.01", "depth.02"] assert outputs[0]["simulated_time"] < outputs[1]["simulated_time"] assert summary["simulated_time"] >= outputs[1]["simulated_time"] - for output in outputs: - assert output["timestamp"] == timestamp_for(output["simulated_time"]) + for minute, output in enumerate(outputs, start=1): + assert output["timestamp"] == f"{minute} minutes" assert output["discharge"] is None assert session_tools.r_info(map=output["depth"], format="json")["rows"] == 5 diff --git a/raster/r.sim/r.sim.water/tests/r_sim_water_test.py b/raster/r.sim/r.sim.water/tests/r_sim_water_test.py index 8dcf2f02736..0e6f67a61f3 100644 --- a/raster/r.sim/r.sim.water/tests/r_sim_water_test.py +++ b/raster/r.sim/r.sim.water/tests/r_sim_water_test.py @@ -529,12 +529,12 @@ def test_duration_affects_time_series_progression(long_slope_session): """A longer duration must create more time-series output maps. With output_step=5, duration=10 produces maps at t=5,10 while - duration=20 produces maps at t=5,10,15 and, if walkers survive to the - end, also t=20. The exact count of the longer run is not pinned: whether - the final t=20 map appears depends on whether the last walkers leave the - domain before the final step, which varies with seed and platform. We - assert the robust lower bound (at least three maps) and that the longer - simulation yields strictly more maps than the shorter one. + duration=20 produces maps at t=5,10,15,20. The exact count of the longer + run is not pinned: when the last walkers leave the domain early, which + varies with seed and platform, the steps after that are not written, + only the final t=20 map. We assert the robust lower bound (at least + three maps) and that the longer simulation yields strictly more maps + than the shorter one. """ tools = Tools(session=long_slope_session) @@ -573,8 +573,8 @@ def test_duration_affects_time_series_progression(long_slope_session): f"10-min simulation with output_step=5 should produce 2 maps, got {len(maps_10)}" ) - # 20-min run should produce at least 3 time-series maps (t=5, t=10, t=15; - # t=20 as well when the last walkers survive to the final step). + # 20-min run should produce at least 3 time-series maps (t=5, t=10 and + # t=20, and t=15 as well when the last walkers stay past it). maps_20 = list(tools.g_list(type="raster", pattern="depth_20min*", format="json")) assert len(maps_20) >= 3, ( f"20-min simulation with output_step=5 should produce at least 3 maps, " diff --git a/raster/r.sim/r.sim.water/tests/r_sim_water_time_series_test.py b/raster/r.sim/r.sim.water/tests/r_sim_water_time_series_test.py new file mode 100644 index 00000000000..dbcdd628559 --- /dev/null +++ b/raster/r.sim/r.sim.water/tests/r_sim_water_time_series_test.py @@ -0,0 +1,169 @@ +# SPDX-License-Identifier: GPL-2.0-or-later + +"""Tests of when and under which names r.sim.water writes time series maps""" + +import os + +import pytest + +import grass.script as gs +from grass.experimental.mapset import TemporaryMapsetSession +from grass.tools import Tools + +# On this slope, the time step is about 10.7 s, so a 1-minute output step +# does not fall on an iteration and the old schedule, which wrote every +# (int)(60 / 10.7) = 5 iterations and named maps by the rounded simulated +# time, wrote dep.04 twice. +SLOPE = { + "elevation": "elevation", + "man_value": 0.1, + "duration": 10, + "random_seed": 1, + "nprocs": 1, +} + + +@pytest.fixture(scope="module") +def slope_project(tmp_path_factory): + """Project with a 50 m by 500 m plane at 10 m rising 0.1 m per row""" + project = tmp_path_factory.mktemp("simwe_time_series") / "simwe" + gs.create_project(project) + with gs.setup.init(project, env=os.environ.copy()) as session: + tools = Tools(session=session) + tools.g_region(s=0, n=500, w=0, e=50, res=10, flags="s") + tools.r_mapcalc(expression="elevation = row() * 0.1") + yield session + + +@pytest.fixture +def slope_tools(slope_project): + with ( + TemporaryMapsetSession(env=slope_project.env) as session, + Tools(session=session, capture_stderr=True) as tools, + ): + yield tools + + +def minutes_of_timestamp(tools, name): + """Minutes of a raster timestamp such as '1 minute' or '10 minutes'""" + number, unit = tools.r_timestamp(map=name).text.split() + assert unit in {"minute", "minutes"} + return int(number) + + +def check_series(tools, summary, minutes): + """Series maps are named and timestamped by the given minutes. + + Each map is written at the iteration closest to its time, except the + last one which is written at the end of the run. + """ + outputs = summary["outputs"] + names = [f"dep.{minute:02d}" for minute in minutes] + assert [output["depth"] for output in outputs] == names + listed = [item["name"] for item in tools.g_list(type="raster", format="json")] + assert sorted(listed) == sorted([*names, "elevation"]) + for minute, output in zip(minutes[:-1], outputs[:-1], strict=True): + assert output["timestamp"] == f"{minute} minutes" + assert minutes_of_timestamp(tools, output["depth"]) == minute + assert abs(output["simulated_time"] - minute * 60) <= summary["time_step"] / 2 + assert outputs[-1]["timestamp"] == f"{minutes[-1]} minutes" + assert minutes_of_timestamp(tools, names[-1]) == minutes[-1] + assert outputs[-1]["simulated_time"] == summary["simulated_time"] + + +def test_every_output_step_written_once(slope_tools): + """Each output step gets one map named and timestamped by its time""" + summary = slope_tools.r_sim_water( + **SLOPE, + depth="dep", + walkers_output="walkers", + output_step=1, + flags="tp", + format="json", + ).json + check_series(slope_tools, summary, list(range(1, 11))) + # Walker maps follow the same schedule with an underscore separator. + walkers = [output["walkers"] for output in summary["outputs"]] + assert walkers == [f"walkers_{minute:02d}" for minute in range(1, 11)] + listed = slope_tools.g_list(type="vector", format="json") + assert sorted(item["name"] for item in listed) == walkers + + +def test_final_step_named_by_duration(slope_tools): + """Duration which is not a multiple of output_step still ends the series""" + summary = slope_tools.r_sim_water( + **SLOPE, depth="dep", output_step=3, flags="tp", format="json" + ).json + check_series(slope_tools, summary, [3, 6, 9, 10]) + + +@pytest.mark.parametrize( + ("mintimestep", "time_step", "minutes"), + [ + (17.5, 70, [1, 2, 4, 5, 6, 7, 8, 9, 10]), + (30, 120, [2, 4, 6, 8, 10]), + (45, 180, [3, 6, 9, 10]), + ], +) +def test_time_step_longer_than_output_step( + slope_tools, mintimestep, time_step, minutes +): + """An iteration writes only the output step closest to it + + The simulated time per iteration is the minimum time step times the + time coefficient 4. With 70 s, the iteration ending at 210 s is equally + close to 3 and 4 minutes and writes the later one. With 120 s and 180 s, + the iterations fall exactly on every second and third output step. The + old schedule wrote every (int)(60 / time step) = 0 iterations, that is, + no maps at all. + """ + result = slope_tools.r_sim_water( + **SLOPE, + depth="dep", + output_step=1, + mintimestep=mintimestep, + flags="tp", + format="json", + ) + summary = result.json + assert summary["time_step"] == time_step + check_series(slope_tools, summary, minutes) + assert "longer than output_step" in result.stderr + + +def test_final_step_after_early_stop(slope_tools): + """The series ends with the duration also when all walkers left early + + With a very low Manning's n, all walkers leave the slope before the + first output step. + """ + summary = slope_tools.r_sim_water( + **(SLOPE | {"man_value": 0.002, "duration": 2}), + depth="dep", + output_step=1, + flags="tp", + format="json", + ).json + assert summary["stopped_early"] is True + assert summary["simulated_time"] < 60 + assert [output["depth"] for output in summary["outputs"]] == ["dep.02"] + assert summary["outputs"][0]["timestamp"] == "2 minutes" + + +def test_no_output_step_past_duration(slope_tools): + """An output step past the duration is not written + + With a time step of 160 s and a duration of 11 minutes, the last + iteration ends at 640 s. The output step closest to it is 12 minutes, + which is past the duration. + """ + summary = slope_tools.r_sim_water( + **(SLOPE | {"duration": 11}), + depth="dep", + output_step=3, + mintimestep=40, + flags="tp", + format="json", + ).json + assert summary["time_step"] == 160 + check_series(slope_tools, summary, [3, 6, 9, 11]) diff --git a/raster/r.sim/simlib/hydro.c b/raster/r.sim/simlib/hydro.c index 2bcc53bcc88..457872cf865 100644 --- a/raster/r.sim/simlib/hydro.c +++ b/raster/r.sim/simlib/hydro.c @@ -44,6 +44,8 @@ void main_loop(const Setup *setup, const Geometry *geometry, int iblock; double conn = 1.0; double addac; + // Time of the last time series output step written [s] + double series_time = 0.; // nblock is reserved for Monte Carlo replicas. A future // change will allow nblock > 1, give each replica an @@ -124,12 +126,6 @@ void main_loop(const Setup *setup, const Geometry *geometry, for (i = 1; i <= setup->miter; i++) { /* iteration loop depending on simulation time and deltap */ G_percent(i, setup->miter, 1); - if (setup->iterout > 0 && i % setup->iterout == 0) { - /* nfiterw = i / iterout + 10; - nfiterh = i / iterout + 40; */ - G_debug(2, "iblock=%d i=%d miter=%d nwalk=%d nwalka=%d", iblock, - i, setup->miter, sim->nwalk, sim->nwalka); - } if (sim->nwalka == 0 && i > 1) { summary->stopped_early = true; @@ -178,39 +174,54 @@ void main_loop(const Setup *setup, const Geometry *geometry, } if (grids->zz[k][l] != UNDEF) { - if (grids->inf[k][l] != UNDEF && - grids->inf[k][l] > 0) { - // Walker's contribution to water depth in this - // cell for this timestep [m] - double decr = factor * sim->w[lw].m; - // Compare with the depth the cell can absorb - // this timestep [m] - if (grids->inf[k][l] * setup->deltap > decr) { - // The cell can absorb the full walker. - // Reduce infiltration rate [m/s]. - grids->inf[k][l] -= decr / setup->deltap; - // Eliminate the walker - sim->w[lw].m = 0.; - continue; - } - else { - // The cell can't absorb the full walker. - // Reduce the walker mass by the equivalent - // of what an infiltration-rate source would - // generate as walker weight. - sim->w[lw].m -= sim->rwalk * - grids->inf[k][l] / - setup->sisum; - // Cell's infiltration capacity is fully - // exhausted - grids->inf[k][l] = 0.; - // Eliminate walker if needed - if (sim->w[lw].m < 0.) { + /* Lock only cells which may have capacity left. */ + double inf; +#pragma omp atomic read + inf = grids->inf[k][l]; + bool eliminated = false; + if (inf > 0) { +#pragma omp critical(infiltration) + if (grids->inf[k][l] != UNDEF && + grids->inf[k][l] > 0) { + // Walker's contribution to water depth in + // this cell for this timestep [m] + double decr = factor * sim->w[lw].m; + // Compare with the depth the cell can + // absorb this timestep [m] + if (grids->inf[k][l] * setup->deltap > + decr) { + // The cell can absorb the full walker. + // Reduce infiltration rate [m/s]. +#pragma omp atomic + grids->inf[k][l] -= + decr / setup->deltap; + // Eliminate the walker sim->w[lw].m = 0.; - continue; + eliminated = true; + } + else { + // The cell can't absorb the full + // walker. Reduce the walker mass by the + // equivalent of what an + // infiltration-rate source would + // generate as walker weight. + sim->w[lw].m -= sim->rwalk * + grids->inf[k][l] / + setup->sisum; + // Cell's infiltration capacity is fully + // exhausted +#pragma omp atomic write + grids->inf[k][l] = 0.; + // Eliminate walker if needed + if (sim->w[lw].m < 0.) { + sim->w[lw].m = 0.; + eliminated = true; + } } } } + if (eliminated) + continue; /* Add walker weight to water depth or * concentration. The captured sum includes the @@ -304,13 +315,33 @@ void main_loop(const Setup *setup, const Geometry *geometry, /* Total remaining walkers for this iteration */ sim->nwalka = nwalka; + // Output step j is written at the iteration closest to + // j * output_step and named by that time. With a time step + // longer than output_step, an iteration is the closest one to + // several steps and writes the step closest to it. + bool write_series = false; + if (settings->ts && settings->iterout > 0) { + double iteration_time = simulated_seconds(setup, i); + double step_time = + settings->iterout * + floor(iteration_time / settings->iterout + 0.5); + // A step more than half a time step ahead is closer to the + // next iteration. No step is past the duration. + if (step_time > series_time && + step_time <= + iteration_time + 0.5 * time_step_seconds(setup) && + step_time <= settings->timesec) { + series_time = step_time; + write_series = true; + } + } + /* Changes made by Soeren 8. Mar 2011 to replace the site walker * output implementation */ /* Save all walkers located within the computational region and with valid z coordinates */ if (outputs->outwalk != NULL && - (i == setup->miter || - (setup->iterout > 0 && i % setup->iterout == 0))) { + (i == setup->miter || write_series)) { sim->nstack = 0; for (lw = 0; lw < sim->nwalk; lw++) { @@ -341,15 +372,18 @@ void main_loop(const Setup *setup, const Geometry *geometry, } /* lw loop */ } - if (settings->ts && setup->iterout > 0 && i % setup->iterout == 0) { + if (write_series) { + G_debug(2, "iblock=%d i=%d miter=%d nwalk=%d nwalka=%d", iblock, + i, setup->miter, sim->nwalk, sim->nwalka); /* call output for iteration output */ if (outputs->erdep != NULL) erod(grids->gama, setup, geometry, grids); /* divergence of gama field */ double itime = simulated_seconds(setup, i); - int ii = output_data(itime, conn, setup, geometry, settings, - sim, inputs, outputs, grids, summary); + int ii = + output_data(itime, series_time, conn, setup, geometry, + settings, sim, inputs, outputs, grids, summary); if (ii != 1) G_fatal_error(_("Unable to write raster maps")); } @@ -442,13 +476,16 @@ void main_loop(const Setup *setup, const Geometry *geometry, } /* Write final maps here because we know the last time stamp here */ - if (!settings->ts) { + // A time series always ends with a step named by the duration, also when + // the duration is not a multiple of output_step or the run stopped early. + if (!settings->ts || series_time < settings->timesec) { // All blocks have completed; gama is the eventual cumulative total, // so no extrapolation is needed. conn = 1.0; double itime = simulated_seconds(setup, summary->iterations_completed); - int ii = output_data(itime, conn, setup, geometry, settings, sim, - inputs, outputs, grids, summary); + double name_time = settings->ts ? settings->timesec : itime; + int ii = output_data(itime, name_time, conn, setup, geometry, settings, + sim, inputs, outputs, grids, summary); if (ii != 1) G_fatal_error(_("Cannot write raster maps")); } diff --git a/raster/r.sim/simlib/input.c b/raster/r.sim/simlib/input.c index 9987711e71d..ddb6298fbe8 100644 --- a/raster/r.sim/simlib/input.c +++ b/raster/r.sim/simlib/input.c @@ -339,10 +339,11 @@ int grad_check(Setup *setup, const Geometry *geometry, const Settings *settings, (int)(settings->timesec / (setup->deltap * setup->timec)); /* number of iterations = number of cells to pass */ - setup->iterout = - (int)(settings->iterout / - (setup->deltap * setup->timec)); /* number of cells to pass for - time series output */ + if (settings->ts && settings->iterout > 0 && + time_step_seconds(setup) > settings->iterout) + G_warning(_("Time step of %.2f s is longer than output_step of %d s, " + "some output steps may be skipped"), + time_step_seconds(setup), settings->iterout); fprintf(stderr, "\n"); G_message(_("Min elevation \t= %.2f m\nMax elevation \t= %.2f m\n"), zmin, diff --git a/raster/r.sim/simlib/output.c b/raster/r.sim/simlib/output.c index f322ac65221..393e2b12eb7 100644 --- a/raster/r.sim/simlib/output.c +++ b/raster/r.sim/simlib/output.c @@ -143,12 +143,14 @@ static void write_history(const char *name, double tt, const Setup *setup, G_write_raster_timestamp(name, timestamp); } -/* conn is the sequential-block extrapolation factor nblock/iblock: scales +/* tt is the simulated time reached [s] and name_time the time [s] the maps + * are named and timestamped by, rounded to minutes. + * conn is the sequential-block extrapolation factor nblock/iblock: scales * the cumulative partial sum in gama into an estimator of the eventual total * for snapshots taken before all blocks have run. With nblock = 1 (or after * the final block) conn = 1.0 and the output formulas are unchanged. err is * written from gammas, which already has conn baked into its accumulator. */ -int output_data(double tt, double conn, const Setup *setup, +int output_data(double tt, double name_time, double conn, const Setup *setup, const Geometry *geometry, const Settings *settings, const Simulation *sim, const Inputs *inputs, const Outputs *outputs, const Grids *grids, Summary *summary) @@ -187,7 +189,7 @@ int output_data(double tt, double conn, const Setup *setup, ndigit = 5; /* Convert to minutes */ - tt_minutes = (int)(tt / 60. + 0.5); + tt_minutes = (int)(name_time / 60. + 0.5); /* Create timestamp */ snprintf(timestamp_buf, sizeof(timestamp_buf), "%d minutes", tt_minutes); diff --git a/raster/r.sim/simlib/simlib.h b/raster/r.sim/simlib/simlib.h index be11d77350b..fe71bea6012 100644 --- a/raster/r.sim/simlib/simlib.h +++ b/raster/r.sim/simlib/simlib.h @@ -33,14 +33,13 @@ typedef struct { double hbeta; // Weighting factor for water flow velocity vector double hhmax; // Threshold water depth [m] double frac; // Water diffusion constant - int iterout; // Time interval for creating output maps [minutes] - int timesec; // Time how long the simulation runs [minutes] + int iterout; // Time interval for creating output maps [seconds] + int timesec; // Time how long the simulation runs [seconds] bool ts; // Time series output double mintimestep; // Minimum time step for the simulation [seconds] } Settings; typedef struct { - int iterout; // Number of iterations for creating output maps int miter; // Total number of iterations double chmean; // Mean Manning's n double si0; // Mean rainfall excess (or sediment concentration?) @@ -177,7 +176,7 @@ void main_loop(const Setup *setup, const Geometry *geometry, const Settings *settings, Simulation *sim, ObservationPoints *points, const Inputs *inputs, const Outputs *outputs, Grids *grids, Summary *summary); -int output_data(double tt, double conn, const Setup *setup, +int output_data(double tt, double name_time, double conn, const Setup *setup, const Geometry *geometry, const Settings *settings, const Simulation *sim, const Inputs *inputs, const Outputs *outputs, const Grids *grids, Summary *summary); diff --git a/scripts/i.pansharpen/i.pansharpen.py b/scripts/i.pansharpen/i.pansharpen.py index 535fd032532..703f876c9e0 100755 --- a/scripts/i.pansharpen/i.pansharpen.py +++ b/scripts/i.pansharpen/i.pansharpen.py @@ -776,7 +776,7 @@ def matchhist(original, target, matched): for j in arrays[target]: # find the grey value in target that corresponds to the cdf # closest to the original cdf - if j[1] <= i[1] + min_difference and j[1] >= i[1] - min_difference: + if i[1] - min_difference <= j[1] <= i[1] + min_difference: # build a reclass rules file from the original grey value and # corresponding grey value from target out_line = "%d = %d\n" % (i[0], j[0]) diff --git a/scripts/r.in.wms/wms_drv.py b/scripts/r.in.wms/wms_drv.py index 754c6acb1f7..be171f97865 100644 --- a/scripts/r.in.wms/wms_drv.py +++ b/scripts/r.in.wms/wms_drv.py @@ -348,11 +348,13 @@ def _computeRequestData(self, bbox, tl_corner, tile_span, tile_size, mat_num_bbo for col in ["min_col", "max_col"]: for row in ["min_row", "max_row"]: if ( - self.t_num_bbox["min_row"] <= self.t_num_bbox[row] - and self.t_num_bbox[row] <= mat_num_bbox["max_row"] + self.t_num_bbox["min_row"] + <= self.t_num_bbox[row] + <= mat_num_bbox["max_row"] ) and ( - self.t_num_bbox["min_col"] <= self.t_num_bbox[col] - and self.t_num_bbox[col] <= mat_num_bbox["max_col"] + self.t_num_bbox["min_col"] + <= self.t_num_bbox[col] + <= mat_num_bbox["max_col"] ): self.intersects = True diff --git a/scripts/r.tileset/r.tileset.py b/scripts/r.tileset/r.tileset.py index 5397a850e3b..fd71601bf80 100644 --- a/scripts/r.tileset/r.tileset.py +++ b/scripts/r.tileset/r.tileset.py @@ -225,9 +225,9 @@ def bboxesIntersect(bbox_1, bbox_2) -> bool: cin = [False, False] for i in (0, 1): if ( - (bi_a1[i] <= bi_b1[i] and bi_a2[i] >= bi_b1[i]) + (bi_a1[i] <= bi_b1[i] <= bi_a2[i]) or (bi_a1[i] <= bi_b1[i] and bi_a2[i] >= bi_b2[i]) - or (bi_b1[i] <= bi_a1[i] and bi_b2[i] >= bi_a1[i]) + or (bi_b1[i] <= bi_a1[i] <= bi_b2[i]) or (bi_b1[i] <= bi_a1[i] and bi_b2[i] >= bi_a2[i]) ): cin[i] = True