Subcatchment masking by JoostBuitink · Pull Request #885 · Deltares/Wflow.jl

JoostBuitink · 2026-04-22T15:21:18Z

Issue addressed

Fixes ENH: run a subset of a wflow model #270.
Adds a logging message on which profile type is being used
Fixes an issue where the field capacity variable was not yet included in the standard_names dict
Fixes the maximum subsurface flow not being adjusted for subdaily timesteps
(only when using layered or layered_exponential profiles)

Explanation

Explain how you addressed the bug/feature request, what choices you made and why.

Checklist

Updated tests or added new tests
Branch is up to date with master
Tests & pre-commit hooks pass
Updated documentation if needed
Updated changelog.qmd if needed

Additional Notes (optional)

Note that for the testing, I decided to use the landuse map, as no subcatchment map was present in the example models.

SouthEndMusic

Some general comments. I don't understand everything well enough to approve the PR

SouthEndMusic · 2026-04-23T07:08:07Z

+represents your entire catchment, and all cells that have a value are interpreted as active
+cells. You can update this to another layer which contains a subset of the original data.
+
+Alternatively, you can also point this to a layer with subcatchment information (e.g. based


Suggested change

Alternatively, you can also point this to a layer with subcatchment information (e.g. based

Alternatively, you can point this to a layer with subcatchment information (e.g. based

SouthEndMusic · 2026-04-23T07:09:12Z

+When using this functionality, please take the following points into account:
+
+- Ensure that all upstream regions are included (unless you make use of external inflows)
+- When using water demand and allocation: ensure that there are no allocations area's with


Suggested change

- When using water demand and allocation: ensure that there are no allocations area's with

- When using water demand and allocation: ensure that there are no allocation areas with

SouthEndMusic · 2026-04-23T07:09:52Z

+
+- Ensure that all upstream regions are included (unless you make use of external inflows)
+- When using water demand and allocation: ensure that there are no allocations area's with
+  the same IDs outside your active domain (a warning is shown in the logging)


Suggested change

the same IDs outside your active domain (a warning is shown in the logging)

the same ID outside your active domain (a warning is shown in the logging)

SouthEndMusic · 2026-04-23T07:11:30Z

+- When using local_inertial for routing: a ghost node will be added at the end of each river
+  section (with a set distance and water level). This might cause some difference in the
+  results, as it no longer contains the "real" downstream boundary condition. If this is
+  important, ensure that the gauges for evaluation have sufficient downstream pixels in your


Does pixel here refer to a cell? I feel that's a better word to refer to the concept

Thanks, yes cell is indeed a better word!

SouthEndMusic · 2026-04-23T07:18:11Z

+        config,
+        "land_water_allocation_area__count",
+        LandHydrologySBM;
+        metadata = ParameterMetadata(; type = Int, allow_missing = true),


Why not add this to a standard name map?

It is already part of the standard names, but I needed to overwrite the metadata in order to get the full array (so not only for the active cells, as is done in the lines above this change).

vers-w

LGTM and good to have this functionality!
See my comments below.

vers-w · 2026-05-11T13:49:14Z

+represents your entire catchment, and all cells that have a value are interpreted as active
+cells. You can update this to another layer which contains a subset of the original data.
+
+Alternatively, you can point this to a layer with subcatchment information (e.g. based on


Suggested change

Alternatively, you can point this to a layer with subcatchment information (e.g. based on

Alternatively, you can point `subbasin_location__count` to a layer with subcatchment information (e.g. based on

vers-w · 2026-05-11T13:50:12Z

+```toml
+[input]
+subbasin_location__count = "subcatchment_gauges_obs"
+subbasin_active_area__count = [1016, 1013, 1011, 206, 203, 503]


This is more in line with the already used subbasin_location__count entry:

Suggested change

subbasin_active_area__count = [1016, 1013, 1011, 206, 203, 503]

subbasin_active_location__count = [1016, 1013, 1011, 206, 203, 503]

vers-w · 2026-05-11T13:55:50Z

+- When writing output at specific cells (`index`) and/or x-y positions (`coordinate`): this
+  can cause unreliable results and/or give errors. If you need output at specific locations,


I guess this means a check is not performed for index and coordinate (removing entries that are not part of active indices)? If that's the case maybe add this at the start. I think the "this can cause unreliable results" is a bit unclear.

This was indeed not checked. I adjusted the code to show a warning rather than an error, which means that this message is no longer relevant, and users are a bit more flexible with using the same toml file

vers-w · 2026-05-11T14:00:33Z

+- Fixed a bug where the maximum subsurface flow was not adjusted for subdaily timesteps
+  (only when using `layered` or `layered_exponential` profiles)


Nice catch!
Note that this also needs fixing in the release v1 branch and these lines need then to be moved below v1.0.3 in the changelog.

vers-w · 2026-05-11T16:53:02Z

        # find the node id of the downstream cell
        to_node = searchsortedfirst(indices, to_index)
+        if to_node > length(indices) || indices[to_node] != to_index
+            @warn "Invalid drainage direction value at node `$from_node` (LDD=`$ldd_val`), assuming pit"


Suggested change

@warn "Invalid drainage direction value at node `$from_node` (LDD=`$ldd_val`), assuming pit"

@warn "Invalid drainage direction value at node `$from_node` (LDD=`$ldd_val`), assign pit value at node"

vers-w · 2026-05-11T18:18:04Z

+            config;
+            indices = indices,


Suggested change

config;

indices = indices,

config,

indices,

vers-w · 2026-05-11T18:20:41Z

+        # update 2d map to only include active cells, and convert to missing for inactive cells
+        map_2d = mask_to_indices(map_2d, indices; fill_value = missing)


Suggest to remove this part and add indices as argument to ncread function.

vers-w · 2026-05-11T19:35:23Z

 end
+
+"Mask a 2D array, keeping only values at the given `indices` and setting the rest to `fill_value`."
+function mask_to_indices(


As suggested for the io.jl part to remove this function (and use indices directly in ncread), I think this is now only used for the river location. It is not clear to me why a Union type is required? I would suggest to either not make fill_value optional and use the correct type in the function call (based on array_2d), so you only need eltype(array_2d) or to make use of the metadata of the parameter that has the required type and fill information.

the function is removed and incorporated with ncread

vers-w · 2026-05-11T19:43:08Z

            sel=network.indices,
        )
-        replace!(x -> ismissing(x) ? 0 : x, reservoirs)
+        # check if any reservoirs are found in the active model domain, if not disable reservoir model component


Nice to have this check! For consistency would be good if reservoir network info is also adapted if reservoirs are not found in the active model domain.

And probably good to have a similar check for drains?

I have added a similar check for drains!

vers-w · 2026-05-11T19:55:44Z

+```
+
+::: {.callout-important}
+When using this functionality, please take the following points into account:


Just wondering: any considerations needed when running with multiple threads/this works fine with this functionality?

Not in particular; since the new outlets are assumed pits, the multithreading recognizes them as individual basins

vers-w · 2026-05-12T05:57:12Z

+  results, as it no longer contains the "real" downstream boundary condition. If this is
+  important, ensure that the gauges for evaluation have sufficient downstream cells in your
+  selected model domain.


Is it not better (easier) to recommend to add boundary conditions at these points?

yes true, added that as well!

vers-w · 2026-05-12T06:02:42Z

+  can cause unreliable results and/or give errors. If you need output at specific locations,
+  use the `map` functionality instead, as this will still write output at the targeted
+  locations when running only a subsection of your model.
+:::


Probably good to add what a user needs to take into account for the Darcy based groundwater flow concept?

also made the output writing when using coordinates or index values a bit more relaxed, to throw a warning rather than an error when invalid cells are selected

vers-w

LGTM, except for the failing tests ;-)!

JoostBuitink added 13 commits April 21, 2026 13:29

support running with a mask

ef3c7cb

ensure output maps are correct

c9443b8

specify basin IDs to run

f745ea9

add check for allocation areas

e192e44

redundant line

97d2e24

update docstring

a0d8be5

add test

ec65896

update docs

f577655

add logging message

e9eb1b9

fix test

f0d9c34

update test

9bbd503

fix docs formatting

2de479a

Update changelog.qmd

8fba1fa

JoostBuitink marked this pull request as ready for review April 23, 2026 06:25

SouthEndMusic reviewed Apr 23, 2026

View reviewed changes

JoostBuitink added 5 commits April 24, 2026 10:35

address comments

e142603

add logging message for the kh_profile type

0e799d1

fix theta_fc standardname

defc31f

bugfix ssf

b02a6a4

add note on z_layered

84dd1d1

vers-w requested changes May 11, 2026

View reviewed changes

vers-w reviewed May 12, 2026

View reviewed changes

address review comments part 1

c95b9c3

also made the output writing when using coordinates or index values a bit more relaxed, to throw a warning rather than an error when invalid cells are selected

JoostBuitink mentioned this pull request May 27, 2026

Detecting cells which are causing is_cyclic #666

Open

JoostBuitink added 5 commits May 29, 2026 16:43

add some more sel=indices statements

f1144f5

remove mask_to_indices function

a1ebe55

Merge branch 'master' into subcatchment_masking

9a93121

small fixes

4eec54b

Merge branch 'master' into subcatchment_masking

08c6d7e

JoostBuitink requested a review from vers-w June 1, 2026 12:59

JoostBuitink added 5 commits June 1, 2026 15:34

Merge branch 'master' into subcatchment_masking

bda0cbb

fix some merge conflicts

370c90d

update test

f011b38

Merge branch 'master' into subcatchment_masking

495957b

Merge branch 'master' into subcatchment_masking

51a24ee

vers-w approved these changes Jun 10, 2026

View reviewed changes

JoostBuitink added 6 commits June 10, 2026 09:59

some bugfixes

a49240a

Merge branch 'master' into subcatchment_masking

1239287

fix tests

8cf16d1

Merge branch 'master' into subcatchment_masking

7b7f64f

Merge branch 'master' into subcatchment_masking

6fff784

Update sbm_metadata.json

bc82fcd

JoostBuitink merged commit aa99309 into master Jun 10, 2026
11 of 14 checks passed

JoostBuitink deleted the subcatchment_masking branch June 10, 2026 13:46

SouthEndMusic mentioned this pull request Jun 23, 2026

Increasing efficiency during global optimisation of Wflow model #951

Closed

	Alternatively, you can also point this to a layer with subcatchment information (e.g. based
	Alternatively, you can point this to a layer with subcatchment information (e.g. based

	- When using water demand and allocation: ensure that there are no allocations area's with
	- When using water demand and allocation: ensure that there are no allocation areas with

	the same IDs outside your active domain (a warning is shown in the logging)
	the same ID outside your active domain (a warning is shown in the logging)

	Alternatively, you can point this to a layer with subcatchment information (e.g. based on
	Alternatively, you can point `subbasin_location__count` to a layer with subcatchment information (e.g. based on

	subbasin_active_area__count = [1016, 1013, 1011, 206, 203, 503]
	subbasin_active_location__count = [1016, 1013, 1011, 206, 203, 503]

		- When writing output at specific cells (`index`) and/or x-y positions (`coordinate`): this
		can cause unreliable results and/or give errors. If you need output at specific locations,

		- Fixed a bug where the maximum subsurface flow was not adjusted for subdaily timesteps
		(only when using `layered` or `layered_exponential` profiles)

	@warn "Invalid drainage direction value at node `$from_node` (LDD=`$ldd_val`), assuming pit"
	@warn "Invalid drainage direction value at node `$from_node` (LDD=`$ldd_val`), assign pit value at node"

		# update 2d map to only include active cells, and convert to missing for inactive cells
		map_2d = mask_to_indices(map_2d, indices; fill_value = missing)

Uh oh!

Conversation

JoostBuitink commented Apr 22, 2026 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Issue addressed

Explanation

Checklist

Additional Notes (optional)

Uh oh!

SouthEndMusic left a comment

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

vers-w left a comment

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

vers-w left a comment

Choose a reason for hiding this comment

Uh oh!

Uh oh!

Reviewers

Assignees

Labels

Projects

Milestone

Development

Uh oh!

3 participants

JoostBuitink commented Apr 22, 2026 •

edited

Loading