1176 write out daily hydrology values for calibration validation by vgro · Pull Request #1949 · ImperialCollegeLondon/virtual_ecosystem

vgro · 2026-10-22T09:57:20Z

This PR adds daily values from the hydrology model andhourly values from the abiotic model to the data object. I introduced new dimensions for these additional variables and saved then with the suffix _daily or _hourly for each relevant variable.

The hydrology csv diagnostic option is removed.

I made a small change to the diurnal temperature range, too, this needs further investigation, see #1950.

Fixes #1176

Type of change

New feature (non-breaking change which adds functionality)
Optimization (back-end change that speeds up the code)
Bug fix (non-breaking change which fixes an issue)

Key checklist

Make sure you've run the pre-commit checks: $ pre-commit run -a
All tests pass: $ poetry run pytest

Further checks

Code is commented, particularly in hard-to-understand areas
Tests added that prove fix is effective or that feature works
Relevant documentation reviewed and updated
If this PR changes the expected model inputs, are the relevant people from the
data science team included as reviewers? -> NA

…libration-validation

…alidation' of https://github.com/ImperialCollegeLondon/virtual_rainforest into 1176-write-out-daily-hydrology-values-for-calibration-validation

Codecov Comments Bot (codecov-commenter) · 2026-10-22T10:40:23Z

Codecov Report

❌ Patch coverage is 95.57522% with 5 lines in your changes missing coverage. Please review.
✅ Project coverage is 89.78%. Comparing base (f99e14e) to head (b02e669).

Files with missing lines	Patch %	Lines
virtual_ecosystem/models/abiotic/abiotic_tools.py	91.17%	3 Missing ⚠️
virtual_ecosystem/models/abiotic/microclimate.py	94.59%	2 Missing ⚠️

Additional details and impacted files

@@            Coverage Diff             @@
##             main    #1949      +/-   ##
==========================================
- Coverage   89.79%   89.78%   -0.01%     
==========================================
  Files          83       83              
  Lines        9678     9708      +30     
==========================================
+ Hits         8690     8716      +26     
- Misses        988      992       +4

☔ View full report in Codecov by Harness.
📢 Have feedback on the report? Share it here.

🚀 New features to boost your workflow:

❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

David Orme (davidorme)

vgro, some questions - some things I don't understand.

David Orme (davidorme) · 2026-10-22T13:48:28Z

+        dims: tuple[Hashable, ...]
        if values_out.ndim == 1:
-            dims = ["cell_id"]
+            dims = ("cell_id",)
        elif values_out.ndim == 2:
-            dims = ["layers", "cell_id"]
+            dims = ("layers", "cell_id")
        else:
            raise ValueError(
                f"Unsupported dimensions for variable '{var}': {values.shape}"
            )

        output[var] = DataArray(values_out, dims=dims)

+    for var in hourly_vars:
+        if var in data_record:
+            values = np.asarray(data_record[var])
+            hourly_dims: tuple[Hashable, ...]
+            if values.ndim == 2:
+                hourly_dims = ("hour", "cell_id")
+            elif values.ndim == 3:
+                hourly_dims = ("hour", "layers", "cell_id")
+            else:
+                raise ValueError(
+                    f"Unsupported hourly dimensions for variable '{var}': "
+                    f"{values.shape}"
+                )
+        elif var in output:
+            n_hours = time_dim
+            if n_hours is None:
+                if not data_record:
+                    raise ValueError(
+                        f"Cannot determine the hour dimension for static output'{var}'."
+                    )
+                n_hours = np.asarray(next(iter(data_record.values()))).shape[0]
+            values = np.broadcast_to(
+                output[var].to_numpy(),
+                (n_hours, *output[var].shape),
+            ).copy()
+            hourly_dims = ("hour", *output[var].dims)
+        else:
+            raise ValueError(
+                f"Hourly output requested for unavailable variable '{var}'."
+            )
+
+        layer_template = layer_structure.from_template()
+        coords: dict[Hashable, Any] = {
+            "hour": np.arange(values.shape[0]),
+            "cell_id": layer_template["cell_id"],
+        }
+        if "layers" in hourly_dims:
+            n_layers = values.shape[hourly_dims.index("layers")]
+            if n_layers == layer_template.sizes["layers"]:
+                coords["layers"] = layer_template["layers"]
+                coords["layer_roles"] = (
+                    "layers",
+                    layer_template["layer_roles"].to_numpy(),
+                )
+            else:
+                coords["layers"] = np.arange(n_layers)
+        output[f"{var}_hourly"] = DataArray(values, dims=hourly_dims, coords=coords)
+


I can't quite follow what the outputs are here. I think it would be faster and more efficient here to define a set of specific templates at _setup and then map variable names onto those templates.

It is hard though to see how many templates we're talking about here. It seems like four?

dims = ("cell_id",) dims = ("layers", "cell_id") hourly_dims = ("hour", "cell_id") hourly_dims = ("hour", "layers", "cell_id")

But then there are ones with fewer layers?

Probably yes, this is not a new problem but so far I have hacked around it. I'll look into it and see what makes sense. Just to confirm, do you mean the layer_structure.from_template()?

So I think you can start with a maximal set of coordinates. And then there is a long but limited set of possible outputs for the model where we know what coordinates are used. So, you could use the maximal set of coords to define a smaller number of specific templates, or you could define which coords are required (or easier, not required) for each variable. The reduced coordinates object then defines the expected shape of the data. The code below is a sketch that defines the dims per variable, removes most of the heavy lifting of coordinate juggling - still not sure about the reduced layers though.

# Define the most complex shape all_coords = layer_template.coords.assign(hours=np.arange(24)) # Define unrequired coords unused_dims = { "var_name_1": ("hours", "layers",), "var_name_1_hourly": ("layers",), } # Build outputs outputs = {} for var, data in vars: try: outputs[var] = DataArray(data, coords=all_coords.drop_dims(unused_dims[var])) except xr.CoordinateValidationError: raise ValueError(f"Unexpected data shape with {var}: {data.shape} not {tuple(c.dims.values())}")

I'm not sure about the best specific implementation, but I agree with David Orme (@davidorme) that everything that could be offloaded to the setup because it is information already known, so it is run only once, will be a gain. In the end, all these nested conditionals and loops will slow down the execution unnecessary. Besides, the code will be shorter and simpler (this method is already in the long side).

Diego Alonso Álvarez (dalonsoa)

I mostly agree with David's comments about trying to do some templating upfront. Also, some functions are becoming a bit large, so when possible I would split them into smaller chunks.

Diego Alonso Álvarez (dalonsoa) · 2026-10-23T05:02:16Z

+        dims: tuple[Hashable, ...]
        if values_out.ndim == 1:
-            dims = ["cell_id"]
+            dims = ("cell_id",)
        elif values_out.ndim == 2:
-            dims = ["layers", "cell_id"]
+            dims = ("layers", "cell_id")
        else:
            raise ValueError(
                f"Unsupported dimensions for variable '{var}': {values.shape}"
            )

        output[var] = DataArray(values_out, dims=dims)

+    for var in hourly_vars:
+        if var in data_record:
+            values = np.asarray(data_record[var])
+            hourly_dims: tuple[Hashable, ...]
+            if values.ndim == 2:
+                hourly_dims = ("hour", "cell_id")
+            elif values.ndim == 3:
+                hourly_dims = ("hour", "layers", "cell_id")
+            else:
+                raise ValueError(
+                    f"Unsupported hourly dimensions for variable '{var}': "
+                    f"{values.shape}"
+                )
+        elif var in output:
+            n_hours = time_dim
+            if n_hours is None:
+                if not data_record:
+                    raise ValueError(
+                        f"Cannot determine the hour dimension for static output'{var}'."
+                    )
+                n_hours = np.asarray(next(iter(data_record.values()))).shape[0]
+            values = np.broadcast_to(
+                output[var].to_numpy(),
+                (n_hours, *output[var].shape),
+            ).copy()
+            hourly_dims = ("hour", *output[var].dims)
+        else:
+            raise ValueError(
+                f"Hourly output requested for unavailable variable '{var}'."
+            )
+
+        layer_template = layer_structure.from_template()
+        coords: dict[Hashable, Any] = {
+            "hour": np.arange(values.shape[0]),
+            "cell_id": layer_template["cell_id"],
+        }
+        if "layers" in hourly_dims:
+            n_layers = values.shape[hourly_dims.index("layers")]
+            if n_layers == layer_template.sizes["layers"]:
+                coords["layers"] = layer_template["layers"]
+                coords["layer_roles"] = (
+                    "layers",
+                    layer_template["layer_roles"].to_numpy(),
+                )
+            else:
+                coords["layers"] = np.arange(n_layers)
+        output[f"{var}_hourly"] = DataArray(values, dims=hourly_dims, coords=coords)
+


I'm not sure about the best specific implementation, but I agree with David Orme (@davidorme) that everything that could be offloaded to the setup because it is information already known, so it is run only once, will be a gain. In the end, all these nested conditionals and loops will slow down the execution unnecessary. Besides, the code will be shorter and simpler (this method is already in the long side).

Diego Alonso Álvarez (dalonsoa) · 2026-10-23T05:17:33Z

-                daily_lists=daily_lists,
-                water_balance_fluxes=water_balance_fluxes,
-                daily_diagnostics=daily_diagnostics,
+        for name, values in daily_lists.items():


Do the dimensions and coordinates of these change every iteration? Because if they do not, you could defined them upfront (in the setup or in the first iteration, and then cache the result), and use them in every subsequent calls. It might not save much time, but will make the code much simpler, just populating the templates for the given variable. Much like pointed out for the hourly case.

Yup - I agree. Creating template arrays in setup along with the mapping of var to template should make it compact (apart from the mapping) and clean.

The plants model has a similar thing in _setup - some of these could even be moved into core if enough models use them, but the hour and day dims are going to be local templates.

virtual_ecosystem/virtual_ecosystem/models/plants/plants_model.py

Lines 426 to 449 in f99e14e

# Define xarray templates

self.data_object_templates = {

"cnp_pft": xr.DataArray(

data=np.zeros((self.grid.n_cells, len(self.flora.pft_name), 3)),

coords={

"cell_id": self.data["cell_id"],

"pft": list(self.flora.pft_name),

"element": ["C", "N", "P"],

},

),

"cnp": xr.DataArray(

data=np.zeros((self.grid.n_cells, 3)),

coords={"cell_id": self.data["cell_id"], "element": ["C", "N", "P"]},

),

"cell": xr.zeros_like(self.data["elevation"]),

"pft": xr.DataArray(

data=np.zeros((self.grid.n_cells, len(self.flora.pft_name))),

coords={

"cell_id": self.data["cell_id"],

"pft": list(self.flora.pft_name),

},

),

}

Diego Alonso Álvarez (dalonsoa) · 2026-10-23T05:19:18Z

 variable_type = "float"

+[[variable]]
+axis = ["spatial"]


I'm I bit unsure about the approach here, but these variables will have another dimension, even if it is not validated anywhere. I think it should be disclosed here, if just as a sort of documentation.

Yeah - this is a longstanding lag in the documentation/internal. The axis data is currently woefully incomplete and we currently don't have validation on the non-spatial axes.

I agree this should be declared here, but it does not work with the current setup, I am getting pydantic ValueErrors if I add a new dimension. I am not very familiar with the BaseModel and the axis validation, is this something to hold or can I include this in the core easily?

This is definitely a different PR. Convert them to comments in the TOML if you want to capture the data, but otherwise this is a much bigger fix than this PR.

vgro added 6 commits October 21, 2026 14:26

remove daily diagnostics to csv

38021db

daily hydrology values added

5744fe0

test fix

4ab5a6f

added hourly data to data object

ae241bf

start day diemnsion with 0

9493ad3

fix diurnal temperature range

97f32b2

vgro linked an issue Oct 22, 2026 that may be closed by this pull request

Write out daily hydrology values for calibration validation #1176

Open

vgro and others added 3 commits October 22, 2026 10:57

Merge branch 'main' into 1176-write-out-daily-hydrology-values-for-ca…

d6b76d1

…libration-validation

Merge branch 'main' into 1176-write-out-daily-hydrology-values-for-ca…

d17dabf

…libration-validation

Merge branch '1176-write-out-daily-hydrology-values-for-calibration-v…

ca56020

…alidation' of https://github.com/ImperialCollegeLondon/virtual_rainforest into 1176-write-out-daily-hydrology-values-for-calibration-validation

vgro requested review from Diego Alonso Álvarez (dalonsoa) and David Orme (davidorme) and removed request for David Orme (davidorme) October 22, 2026 10:47

David Orme (davidorme) reviewed Oct 22, 2026

View reviewed changes

Diego Alonso Álvarez (dalonsoa) requested changes Oct 23, 2026

View reviewed changes

vgro added 11 commits October 26, 2026 09:08

removed _hourly variables for static variables

96da694

added output templates to abiotic model setup

a810380

added output templates to hydrology model

b530b71

replaced soil_moisture with soil_temperature in tests

26fd756

small test fixes

9004e6f

renamed AbioticOutputTemplate

c3a0f81

added HydrologyOutputTemplate class

318a7e8

move output template class to _tools.py

6caf66e

remove unnecessary copies

0a57b15

moved to suggested plant implementation

da04fa5

docs fix

b02e669

	# Define xarray templates
	self.data_object_templates = {
	"cnp_pft": xr.DataArray(
	data=np.zeros((self.grid.n_cells, len(self.flora.pft_name), 3)),
	coords={
	"cell_id": self.data["cell_id"],
	"pft": list(self.flora.pft_name),
	"element": ["C", "N", "P"],
	},
	),
	"cnp": xr.DataArray(
	data=np.zeros((self.grid.n_cells, 3)),
	coords={"cell_id": self.data["cell_id"], "element": ["C", "N", "P"]},
	),
	"cell": xr.zeros_like(self.data["elevation"]),
	"pft": xr.DataArray(
	data=np.zeros((self.grid.n_cells, len(self.flora.pft_name))),
	coords={
	"cell_id": self.data["cell_id"],
	"pft": list(self.flora.pft_name),
	},
	),
	}

Conversation

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

Uh oh!

Type of change

Key checklist

Further checks

Uh oh!

Codecov Comments Bot (codecov-commenter) commented Oct 22, 2026 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Codecov Report

Uh oh!

David Orme (davidorme) left a comment

Choose a reason for hiding this comment

Uh oh!

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

David Orme (davidorme) Oct 22, 2026 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Diego Alonso Álvarez (dalonsoa) 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!

David Orme (davidorme) Oct 26, 2026 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Reviewers

Uh oh!

Assignees

Labels

Projects

Milestone

Development

Uh oh!

4 participants

vgro commented Oct 22, 2026 •

edited

Loading

Codecov Comments Bot (codecov-commenter) commented Oct 22, 2026 •

edited

Loading

David Orme (davidorme) Oct 22, 2026 •

edited

Loading

David Orme (davidorme) Oct 26, 2026 •

edited

Loading