Showing posts with label GEOFRAME. Show all posts
Showing posts with label GEOFRAME. Show all posts

Friday, August 21, 2026

Six Notebooks to represent water budgets from compartment models (with EPNs inside)

Since the paper on GEOtop (Rigon et al., 2006), which claimed that a hydrological model should look at the water and energy budgets, rather than at single fluxes, in order to gain real knowledge of the system, I have been struggling to build models that respond to this effort. Estimating a discharge, an evaporation, a snow water equivalent one at a time is not the same thing as closing a budget: it is the budget — the simultaneous accounting of inputs, outputs and storage variations — that constrains the pieces to be mutually coherent, and it is from that coherence that understanding comes.

Detail of a Giacomo Balla artwork

The effort has grown in various directions since then. On the conceptual and strategic side, it evolved into the vision of Digital eARth Twin Hydrology systems, the DARTHs (Rigon et al., 2022), discussed in several posts on this blog (see the posts on DARTHs). On the software side, it required a complete rebuilding of the informatics of these systems, deployed in the GEOframe system: a component-based infrastructure built on the OMS3 framework (David et al., 2013; Formetta et al., 2014), whose NewAGE branch (Bancheri et al., 2020) is the one we routinely use for operational and research modelling (see the GEOframe posts). On the formal side, it produced the Extended Petri Nets (EPN) formalism (Bancheri, Serafin and Rigon, 2019), which gives a graphical and mathematical grammar for writing budget-based models: places for storages, transitions for fluxes, controllers for the variables that regulate rates without exchanging mass, splitters for the partition of fluxes.

In all of this effort, one persistent issue has been how to represent the water budget in a compact and coherent way (on the energy budget we are still working). Compact, because a budget involves many time series at once — inputs, outputs, several storages — and they must be readable on a single temporal abscissa. Coherent, because the representation should not merely juxtapose the series but make their mutual constraint visible: at every time step, inputs minus outputs must equal the sum of storage variations, and a figure of the budget should let the reader see whether this closure holds, column by column.

One step in this direction is the production of six Jupyter Notebooks, obtained with the help of Claude (Anthropic's AI assistant) and certainly improvable, which you find on OSF. As usual, this is a seed for improvements.

What the notebooks contain

The six parts share a single Python library (hydrobudget.py, accompanied by a test suite) — with a second module, epn_registry.py, for the model dictionary of Part 6 — and exchange data through feather files and JSON declarations, so that each notebook can be read — and modified — on its own.

Part 1 — Synthetic data. A minimal three-bucket model (snowpack, root zone, groundwater, with air temperature as forcing) generates the datasets used throughout: a one-year demonstration run, a four-year analysis run with an imposed dry year, and — this is a methodological point I care about — an independent thirty-year baseline from which all climatologies and drought thresholds are computed, so that the reference statistics are never contaminated by the period being analyzed. Everything is saved to feather, and every generated dataset is closure-checked at birth — including, now, a small hillslope-and-riparian system, which is the net whose drawing in Part 4 has a crossing: it could previously be drawn but never verified, for want of a table. This is shortcut since I was not having real data. 

Part 2 — Representing the budget. The central figure of the series: inputs drawn as bars hanging from the top of the panel (the hyetograph convention), outputs rising from the bottom, storage variations as signed stacked bars below, optionally the total storage in a third panel — all on one shared temporal abscissa, at any time step from hourly to yearly, with fluxes summed and states sampled at period end under aggregation. The budget check is built into the figure: a line marks in − out at every step, and each column of storage variations must reach it exactly. On the synthetic data the maximum column-wise closure error is of the order of 10⁻¹⁴ mm; on real data, the columns where the bars miss the line become a diagnostic of where measurements or model do not close. The bookkeeping is declared, not hard-coded — one learns, for instance, that snowfall is not a flux input to the soil but an input to the snowpack storage, and the figure enforces the distinction. What makes the residual informative, rather than merely small, is how it behaves when the control volume is declared wrongly: the notebook now tabulates four candidate volumes — the soil view as drawn, the same view with Δ(snowpack) added, the whole system, and the whole system with Δ(snowpack) removed. The two consistent pairings close to about 10⁻¹³ mm; both inconsistent ones fail by the same 116.32 mm, which is the largest monthly change in snowpack storage. Equal failure at the scale of exactly one term is the signature of a term placed on the wrong side of a boundary, not of an accumulating error.

Part 3 — Droughts. Once the budget is an object, droughts become excursions of that object: persistent negative anomalies of storages and accumulated fluxes with respect to their (baseline) climatology, identified by the classical threshold-level method of run theory (Yevjevich, 1967), with pooling, severities and intensities in explicit units. Because each drought type lives in a different variable of the budget, the notebook separates them — meteorological, snow, soil moisture, streamflow, groundwater — and displays their propagation chain on a timeline (in the spirit of Van Loon, 2015). Snow droughts are further classified as dry, warm, warm-and-dry or other, using precipitation ratios and temperature anomalies (Harpold et al., 2017): this is why the temperature is generated and stored in Part 1. The reference is now computed by day of year rather than by calendar month, and the notebook shows both: the step discontinuities of a monthly reference are visible in the anomaly it produces, and therefore in the spells the threshold method extracts from it — which is a way of seeing how much the choice of reference, and not only the choice of threshold, is doing. There is certainly more to work out, in this part.

Part 4 — Representing the EPNs. The same systems drawn as Extended Petri Nets, following the graphical conventions of Bancheri, Serafin and Rigon (2019): colored circles for places, squares for transitions inheriting the color of the place they exit, framed triangles for controllers with dashed information arcs, diamonds for splitters with their partition fractions, everything laid on an integer grid with mathematical symbols as labels. Two details I find pleasing: arc routing can circumnavigate nodes through waypoints, and when a crossing between arcs is genuinely unavoidable — the notebook distinguishes by computation the avoidable from the necessary ones — it is denoted by a small open circle, the old circuit-diagram convention for "crossing without exchange". The rules that place the nodes are worth stating, because they are computable and not aesthetic preferences: places and transitions stay on grid nodes, arrows are straight lines or compositions of straight lines except for the controllers, a multi-part arrow is permitted precisely in order to avoid a crossing, and among the admissible configurations the more compact is chosen. The justification for not trading one against the other is that a net which obeys the no-crossing rule in some layout can always be deformed onto the grid, so compactness and crossing-freedom are not in competition. Part 6 gives these rules a test they cannot have on a hand-made example: applied blind to 47 published nets, none of which declares a single coordinate, 34 come out with no crossing at all. The basic idea is that giving as input the EPN representation of the model, the program knows which data to expect and consequently which graphs to produce, see Part 5 below. There are various improvements that can be envisioned also here. 

Part 5 — The EPN as budget analysis. The point of the whole exercise: the EPN declaration is not only a drawing, it is the bookkeeping. From the same JSON that draws the net, the code derives the whole-system budget and one budget per place, verifies the closure node by node, and plots the two-panel budget figure of every place of the net. The same JSON then acts as a data contract for simulation outputs: we propose a format that GEOframe-NewAGE results should obey — one tidy table per hydrological response unit, one column per EPN transition (internal fluxes included, which is an argument for making them first-class outputs of any simulator), one column per place state — so that if simulator and analysis share the same EPN, every node must close to numerical precision, and any residual is either a format violation or a bug. Coherence between the GEOframe outputs and the EPN scheme of the simulator, in other words, becomes something one can check mechanically. One lesson from making the contract strict: the closure check used to skip the missing steps, so a table that was 95 % gaps reported a better residual than an honest one — precisely backwards. It now reports the sample size beside the error (on a deliberately broken demo table, n_steps = 4 of_possible = 364), so a good number obtained from four steps cannot pass for a good number, and a declared column that is present but mostly empty now fails the contract outright.

Part 6 — The EPN dictionary. If an EPN is a grammar, then it should be possible to write a vocabulary in it, and the vocabulary is worth more than any single sentence. Part 6 transcribes the models of the MARRMoT collection (Knoben et al., 2019) as EPNs — 47 nets over 46 model identifiers, since one model is drawn twice in the source with two different structures — and stores each as a JSON entry carrying the net, the governing equations as they are printed beside it, the constitutive closures, the parameters, and the page it was read from. The collection then becomes something one can query rather than leaf through: find the models with a snowpack, the ones with an IUH transition, the ones with at least four storages. Sizes run from one storage to nine, with a median of four; 43 of the nets have a junction, 33 a splitter, one an IUH.

The reason for storing the equations next to the drawing is that it makes the transcription falsifiable. The deck prints the ODE of every store beside its net, so dS/dt = Ps − Es − Perc is an independent statement of which fluxes touch that store and with which sign, and the code compares the two readings: 173 of the 187 places carry a printed equation, and of those 161 agree with the arcs drawn around them. Twelve places in ten models do not. I have not guessed at those: each is recorded in a table with the arc missing on one side or the other, because deciding which of the two readings is right means going back to the original paper, not to my transcription. Fourteen places print no equation at all in the deck, so they cannot be checked — a gap in the source, and it seemed better to count it than to hide it.

The second half of the notebook is about adding an entry, since a dictionary nobody can extend is a catalogue. A small two-bucket model is walked through the three gates an entry must pass — the schema, the equation cross-check, and a redraw — and each gate is broken first, so that the reader sees the error message rather than a description of it.

Two things I did not expect to learn from doing this. The first is that transcribing a drawing is the weak link in the whole exercise: the errors the cross-check caught in my own transcriptions were exactly the kind a careful reader makes — a flux attached to the wrong store, a splitter branch that should have been an evaporative loss, a melt term turned into a self-loop. The second is that the disagreements which survive are interesting in themselves. A published net and a published equation that do not say the same thing are a fact about the literature, not about my code.

A seed

The notebooks run top to bottom without errors, the figures are checked geometrically before they ship, and the library carries 112 tests — the dictionary's own invariants among them: that an entry's stored counts match its net, that its status flag is the computed verdict of the cross-check rather than an opinion, that the search index is recomputable from the entries it summarises — but none of this makes them finished. The synthetic model is deliberately minimal; the drought thresholds are heuristic; the energy budget is absent; the connection to real GEOframe runs is, for now, a contract waiting for its first signatures. As usual, this is a seed for improvements: take them, break them, and tell me where.

One practical warning for anyone who does. The .ipynb files are generated, by the scripts build_01.pybuild_06.py: an edit made in Jupyter is discarded by the next build. This is deliberate — it is what makes the series reproducible, and the build is now byte-identical, so two builds of unchanged sources give identical files and a real change shows up as a real difference — but it means that the natural thing to do, which is to open a notebook and fix something in place, is the one thing that will quietly lose the fix. Edit the builder.

Files are here: https://osf.io/jt7z8/files/osfstorage

References

  • Bancheri, M., Serafin, F., & Rigon, R. (2019). The representation of hydrological dynamical systems using Extended Petri Nets (EPN). Water Resources Research, 55(11), 8895–8921. https://doi.org/10.1029/2019WR025099
  • Bancheri, M., Rigon, R., & Manfreda, S. (2020). The GEOframe-NewAge modelling system applied in a data-scarce environment. Water, 12(1), 86. https://doi.org/10.3390/w12010086
  • David, O., Ascough II, J. C., Lloyd, W., Green, T. R., Rojas, K. W., Leavesley, G. H., & Ahuja, L. R. (2013). A software engineering perspective on environmental modeling framework design: The Object Modeling System. Environmental Modelling & Software, 39, 201–213. https://doi.org/10.1016/j.envsoft.2012.03.006
  • Formetta, G., Antonello, A., Franceschi, S., David, O., & Rigon, R. (2014). Hydrological modelling with components: A GIS-based open-source framework. Environmental Modelling & Software, 55, 190–200. https://doi.org/10.1016/j.envsoft.2014.01.019
  • Harpold, A. A., Dettinger, M., & Rajagopal, S. (2017). Defining snow drought and why it matters. Eos, 98. https://doi.org/10.1029/2017EO068775
  • Knoben, W. J. M., Freer, J. E., Fowler, K. J. A., Peel, M. C., & Woods, R. A. (2019). Modular Assessment of Rainfall–Runoff Models Toolbox (MARRMoT) v1.2: an open-source, extendable framework providing implementations of 46 conceptual hydrologic models as continuous state-space formulations. Geoscientific Model Development, 12(6), 2463–2480. https://doi.org/10.5194/gmd-12-2463-2019
  • Rigon, R., Bertoldi, G., & Over, T. M. (2006). GEOtop: A distributed hydrological model with coupled water and energy budgets. Journal of Hydrometeorology, 7(3), 371–388. https://doi.org/10.1175/JHM497.1
  • Rigon, R., Formetta, G., Bancheri, M., Tubini, N., D'Amato, C., David, O., & Massari, C. (2022). HESS Opinions: Participatory Digital eARth Twin Hydrology systems (DARTHs) for everyone. Hydrology and Earth System Sciences, 26, 4773–4800. https://doi.org/10.5194/hess-26-4773-2022
  • Van Loon, A. F. (2015). Hydrological drought explained. WIREs Water, 2(4), 359–392. https://doi.org/10.1002/wat2.1085
  • Yevjevich, V. (1967). An objective approach to definitions and investigations of continental hydrologic droughts. Hydrology Papers 23, Colorado State University.

Saturday, June 13, 2026

For future working: A note for a comprehensive object oriented implementation of hydrological dynamical systems

Here is an ambition worth stating plainly: a single engine that can solve any hydrological dynamical system, where building a new model means instantiating a few new classes rather than touching the engine at all. Not a model — a kit. The bricks are the storages and the fluxes; the engine assembles them, differentiates them, and solves them; and the same solvers are reused no matter what physics you bolt on. You add a new process by extending the kit, never by editing the machinery — open to extension, closed to modification — and you do it by writing to interfaces rather than to concrete classes, so the solver never needs to know, or care, what it is solving. The reward is a system that is at once generic, efficient, and expandable, with strikingly little code to maintain.
NOt forgetting that this is very preliminary material, food for thinking, probably bugged material, please find:

This is not a new creed. It is the generic-programming philosophy that Berti set out for scientific computing two decades ago — efficient, reusable components built on the right abstractions, so that code stops being rewritten for every variant — and it is exactly the spirit in which Niccolò and I built WHETGEO-1D (Tubini and Rigon, 2022), where a ClosureEquation interface and a factory let you swap van Genuchten for Brooks–Corey without disturbing a line of the solver. What I want to do here is take that same philosophy and lift it from the one-dimensional soil column up to the whole topology of a dynamical system, and the thing that makes it possible is an old friend.

The friend is the Extended Petri Net. A hydrological dynamical system is one: storages are places, fluxes are transitions, and a third set of objects, the controllers, complete the causal wiring. Follow that picture all the way down into the software, and the thing every modeller secretly dreads — assembling, by hand, the system of equations to hand to a solver — almost disappears.

The grammar, and why "universal" is not a boast

Before the software, the picture. An EPN is drawn with a small, fixed vocabulary, and once you have it the claim of universality stops being rhetorical. A place is a circle, a transition a square; a square carrying a black dot is a driven input, the rain falling into the net; a dashed square is the outlet to the world. A diamond is a splitter, where one quantity divides into branches whose fractions sum to one — precipitation parting into snow and rain. A dotted circle is a collector, a summing junction with no storage of its own. And a little histogram marks a flux computed by convolution, a unit hydrograph. That is very nearly the whole alphabet.

Why believe it is enough? Because Marialaura, Anna De Nardi and I drew all forty-six models of the MARRMoT collection (Knoben et al., 2019) in exactly this vocabulary — Collie, GR4J, TOPMODEL, HBV, Sacramento, VIC, the Tank cascades, the whole zoo — and each one turns out to be just a different choice of places, fluxes and wiring. If a single alphabet spells every word in the dictionary, it is the right alphabet. The classes that follow are nothing more than this grammar given types.

One class of storages, many kinds of flux

The first thing the Petri-net view tells you is an asymmetry. A place is universal: it is a state variable whose time variation we study, and one storage is, as an object, exactly like another. A bucket of soil water, a snowpack, a channel reach — all of them are just \(S_i\) with a balance to satisfy. So in an object-oriented design, storages are a single class.

Fluxes are the opposite. A flux can be a known time series (a forcing), or a constitutive law that declares a mathematical form and a dependence on the state. So a transition wants to be an interface with concrete implementations — an external flux here, a Darcy-Buckingham law there, a power-law discharge somewhere else. Each flux knows two things that matter for what follows: where it moves mass (its source and target places) and which state variables it reads. Those two are not the same, and keeping them apart is the whole trick.

Once you take that seriously, the "many kinds of flux" become a small, nameable family. There is the external flux (a forcing) and the constitutive flux (a law of the state); the controller, a quantity derived from the state that gates a flux without carrying any mass; the splitter and the collector of the grammar above; the stoichiometric flux, one rate metered into several currencies, which we will need the moment evapotranspiration appears; and the convolution flux, the one transition that carries a memory. They all implement the same interface. The only place the software must be careful is with the threshold laws that conceptual models love — the \(\text{if } S>S_{\max}\), the \(\max(\cdot,0)\) — which are not differentiable at the kink; those enter in regularised form, the same medicine the soil-water retention curve already takes.

The two graphs

An EPN actually carries two graphs over the same nodes, and conflating them is, I think, why generic model engines so often turn into a tangle.

The first is the incidence graph — the Petri net proper. It records which flux moves water between which storages, and it is summarised by the incidence matrix \(\mathbf{C}\), with \(+1\) where a flux enters a storage and \(-1\) where it leaves. From it, the entire model collapses into one line:

\[ \frac{\mathrm{d}\mathbf{S}}{\mathrm{d}t} \;=\; \mathbf{C}\,\mathbf{Q}(\mathbf{S},t). \]

That is conservation, and nothing more. The second graph is the causal one, directed on the storages alone: there is an arc from \(S_j\) to \(S_i\) whenever the equation for \(S_i\) contains a flux that reads \(S_j\). This is the graph that decides how the equations are coupled.

And here is where the controllers finally make sense. A controller is an arc that lives in the second graph but not in the first — a flux that is gated or modulated by some storage without any mass passing through that storage. The read arcs and inhibitor arcs of Petri-net theory. They carry no water, so they are absent from \(\mathbf{C}\); but they carry causation, so they are present in the dependency graph and therefore in the Jacobian. That is exactly what we mean when we say controllers "complete the causal structure of the model".

Assembling the system the solver can eat

If you want an implicit step — and for cyclically coupled, stiff budgets you almost always do — backward Euler turns the master equation into a residual,

\[ \mathbf{F}(\mathbf{S}^{\,n+1}) \;:=\; \bigl(\mathbf{S}^{\,n+1}-\mathbf{S}^{\,n}\bigr) -\Delta t\,\mathbf{C}\,\mathbf{Q}(\mathbf{S}^{\,n+1},t^{\,n+1}) \;=\; \mathbf{0}, \]

and a Newton method needs its Jacobian,

\[ \mathbf{J} \;=\; \mathbf{I} \;-\; \Delta t\,\mathbf{C}\,\frac{\partial \mathbf{Q}}{\partial \mathbf{S}}. \]

This is the passage people worry about: how do you build this nonlinear system automatically for an arbitrary model? The answer is that the structure is already in the objects. The sparsity of \(\partial\mathbf{Q}/\partial\mathbf{S}\) is just the list of dependencies each flux declared. So the assembler becomes a single scatter loop over the fluxes — exactly like assembling element matrices in a finite-element code. You walk the fluxes once; each one drops its value into the balance of its two endpoints and its derivative into the matching rows of the Jacobian. You never build \(\mathbf{C}\) or the flux Jacobian as dense matrices; they assemble themselves, entry by entry.

The derivatives I would get from automatic differentiation. If each flux is written over a differentiable scalar (in Java, Hipparchus' DerivativeStructure), one forward evaluation returns both the flux value and its partials. Concretely, a linear reservoir is just

DerivativeStructure evaluate(AdContext ctx, double t) {
    return ctx.value(from).divide(tau);   // Q = S_from / tau
}

and the engine differentiates it for you. The consequence is the thing I actually care about: you write the physics once, per flux, and never again touch the solver. Adding a new constitutive law costs you its forward formula and nothing else.

Loops, sequences, and where the parallelism hides

There is one more gift in the causal graph. Run Tarjan's algorithm on it and you get its strongly connected components. A component with more than one storage is a loop — a knot of mutually dependent equations that has to be solved simultaneously. Every other component is a single equation that can be solved in turn. The graph of components is a DAG, and its topological order is simply the order in which to solve. In the algebra of the Jacobian this is a block-triangular structure; it is the same decomposition that equation-based modelling languages call BLT. And components sitting at the same level of that DAG do not depend on one another, so they can go to different threads. The parallelism was never something we had to impose — it was sitting in the dependency graph all along, waiting to be read off.

The solver, finally, gets to be ignorant. It sees only a residual, a Jacobian, and a starting guess. That is enough to let an off-the-shelf Newton–Raphson handle the easy blocks and a Nested Newton (Casulli–Zanolli) handle the monotone Richards-type ones — without either of them ever knowing what they are solving.

Is it efficient? Yes — but structurally

The natural worry is that all this generality must cost something at runtime. It does not, or rather, the cost lands exactly where it should. The only expensive operation — the implicit Newton solve — happens only inside the loops; the rest of the network, which is most of it, is cheap explicit updates. And the whole structural analysis — finding the loops, the solution order, the sparsity of the Jacobian — depends only on the graph, which never changes in time. So you do it once, at the start, and pay nothing more for it on the millions of timesteps that follow. The efficiency is not a clever trick in the inner loop; it is a consequence of letting the structure tell you where the hard work actually is. The one honest caveat is that you must respect that structure: solve everything monolithically and you throw the gift away. Couple only what is genuinely stiff, and let the rest run loose and parallel.

Net3, upgraded

This is where Francesco Serafin's Net3 comes in, and the fit is almost too neat — Net3 grew out of the same EPN picture. Net3 takes a model as a directed acyclic graph and runs it in parallel, which is wonderful for a river network (a tree is a DAG) but awkward for coupled dynamics, because coupling makes loops, and a DAG cannot hold a loop. The repair is the one move we already have: take each loop, each strongly connected component, and crush it into a single super-node whose insides are the simultaneous solve. The graph of super-nodes is a DAG again, and Net3 schedules it exactly as before. The coupled solving lives inside the node; the parallel routing lives between nodes. Net3 keeps doing what it is good at, and inherits the one thing it could not do.

Three budgets, one net — and the trouble with ET

The original EPN paper already hints that we might want to solve more than water: the energy budget, the carbon budget, all at once. The beauty is that the equation is the same, \(\dot{\mathbf{S}}_b = \mathbf{C}_b\,\mathbf{Q}_b\); only the meaning of the parameters changes. Water places hold storages, energy places hold internal energy, carbon places hold pools. What ties them together are the controllers — quantities derived from one budget's state that reach into another. Temperature, born of the energy budget, governs evapotranspiration in the water budget; soil moisture, born of the water budget, governs the thermal conductivity in the energy budget and the stomata in the carbon budget. Each of these is, once again, an arc that lives in the causal graph but in nobody's incidence matrix.

And then there is the genuinely awkward case, the one that makes the whole thing interesting: a single quantity that belongs to two budgets at once. Evapotranspiration is the textbook offender. It is a loss of water, at rate \(E\); it is also a loss of energy, at rate \(\lambda E\), the latent heat carried away. How do you stop the two budgets from quarrelling about how much of it happened?

The clean answer is to stop thinking of separate budgets and to write one net over all the currencies at once, letting the incidence matrix carry not just \(\pm 1\) but real stoichiometric coefficients. Then ET is a single column with an entry of \(-1\) in the water rows and \(-\lambda\) in the energy rows. You evaluate it once — it has one rate — and drop that one rate into both budgets, scaled appropriately. The consequence is the thing I find quietly satisfying: the latent heat is exactly \(\lambda\) times the water loss at every step of the iteration, not just at the end. The two budgets cannot disagree, because there is only one number. Conservation across currencies stops being something you check and becomes something the structure guarantees.

If you have ever wondered what Penman–Monteith really is, this is it. Penman solves for the surface temperature and the evaporation rate together, between the energy and the mass budgets — which is precisely solving one of these little coupled super-nodes with a shared ET column and a shared temperature controller. The multi-currency net is just Penman–Monteith let off its leash: the same closure, now for any number of budgets, solved numerically rather than by hand. And when all the controllers descend from a single free energy, the cross-couplings ought to be symmetric — Onsager reciprocity — which gives a quiet, structural way to check that the coupled model is thermodynamically honest.

Routing, travel times, and the unit hydrograph

One flux refuses to be memoryless, and it rewards a closer look, because its kernel is not arbitrary — it is a travel-time distribution. The instantaneous unit hydrograph carries the effective rainfall to the outlet by convolution, and the cleanest way to hold it is not the IUH \(f\) itself but its integral, the S-function

\[ s_f(T) \;:=\; \int_0^T f(y)\,\mathrm{d}y, \]

which is, up to the catchment area, the cumulative distribution of travel times — the very residence-time object that runs through the age-ranked budget story (Rigon, Bancheri and Green, 2016). Write the discharge against travel time and a single rain record contributes differences of \(s_f\) over its own length; the records then simply superpose, because the routing is linear and time-invariant. That is the entire numerical recipe, worked out impulse by impulse in the illustrated guide (Rigon et al., 2022b; Rigon et al., 2016).

The interesting part is what the kernel's origin decides. A parametric hydrograph — an exponential, which is just a single linear reservoir; or a Nash cascade — reduces exactly to a little chain of reservoir places, so it folds back into the basic vocabulary and asks for no new type. A geomorphological hydrograph does not: the width function, read off a digital elevation model as the area of the catchment at each flow distance from the outlet and mapped to time by a velocity, gives \(s_f\) straight from the shape of the basin and is no finite chain of reservoirs. That is exactly why the convolution flux has to be a first-class citizen rather than sugar over a cascade. And the linearity is a genuine assumption: when the hydrograph changes shape with the storm — an event-specific GIUH — the kernel becomes a controller of the forcing, and the clean convolution gives way to the general, time-varying case.

There is a quiet bonus for anyone thinking of running the model in real time. Because each incoming rain record commits the discharge for the next several steps — water already fallen, travel times already fixed — the convolution hands you a short forecast for nothing, unmodifiable by rain that has not yet arrived. Feeding the engine live data is then only a matter of swapping the file behind a forcing for a streaming source; the engine never notices whether the rain fell last year or a minute ago.

Why I like this

What pleases me here is how much the architecture buys by simply refusing to mix two things up. Conservation lives in the incidence matrix. Causation lives in the dependency graph. The physics lives in the fluxes, written one at a time. And the assembler — the small piece of code that compiles all of it into a system of equations — is one loop. It is the kind of separation that, once you see it, makes you wonder why the equations ever felt like the hard part. They were never the hard part. The hard part was deciding what was a place, what was a flux, and what was only a controller.

And that, in the end, is the whole point of the kit. Every new model — a different catchment, a snow scheme, a coupled carbon budget — is a handful of new flux classes implementing the same interface, dropped into the same engine, solved by the same Newton. Nothing in the machinery changes; the machinery was closed to modification from the start. The code stays small not because we were clever line by line, but because we let the abstraction carry the weight. That is what generic programming promised for scientific computing, and it is satisfying to watch hydrology turn out to be such a natural place to collect on the promise.

References

Bancheri, M., Serafin, F., and Rigon, R. (2019). The Representation of Hydrological Dynamical Systems Using Extended Petri Nets (EPN). Water Resources Research, 55(11), 8895–8921. doi:10.1029/2019WR025099.

Berti, G. (2000). Generic Software Components for Scientific Computing. PhD thesis, BTU Cottbus. See also Berti, G. (2006), GrAL — the grid algorithms library, Future Generation Computer Systems, 22(1–2), 110–122.

Knoben, W. J. M., Freer, J. E., Fowler, K. J. A., Peel, M. C., and Woods, R. A. (2019). Modular Assessment of Rainfall–Runoff Models Toolbox (MARRMoT) v1.2: an open-source, extendable framework providing implementations of 46 conceptual hydrologic models as continuous state-space formulations. Geoscientific Model Development, 12, 2463–2480. doi:10.5194/gmd-12-2463-2019.

Rigon, R., Bancheri, M., Formetta, G., and de Lavenne, A. (2016a). The geomorphological unit hydrograph from a historical-critical perspective. Earth Surface Processes and Landforms, 41(1), 27–37. doi:10.1002/esp.3855.

Rigon, R., Bancheri, M., and Green, T. R. (2016b). Age-ranked hydrological budgets and a travel time description of catchment hydrology. Hydrology and Earth System Sciences, 20(12), 4929–4947. doi:10.5194/hess-20-4929-2016.

Rigon, R., Formetta, G., Bancheri, M., Tubini, N., D'Amato, C., David, O., and Massari, C. (2022a). HESS Opinions: Participatory Digital eARth Twin Hydrology systems (DARTHs) for everyone — a blueprint for hydrologists. Hydrology and Earth System Sciences, 26, 4773–4800. doi:10.5194/hess-26-4773-2022.

Rigon, R., Franceschi, S., Formetta, G., Bancheri, M., and Tubini, N. (2022b). An illustrated guide to IUH/GIUH estimation. Authorea preprint. doi:10.22541/au.164192110.08629205/v1.

Serafin, F. (2019). Enabling Modeling Framework with Surrogate Modeling Capabilities and Complex Networks (the Net3 subsystem). PhD thesis, University of Trento. See also Serafin, F., David, O., Carlson, J. R., Green, T. R., and Rigon, R. (2021), Environmental Modelling & Software, 146, 105231.

Tubini, N. and Rigon, R. (2022). Implementing the Water, HEat and Transport model in GEOframe (WHETGEO-1D v.1.0): algorithms, informatics, design patterns, open science features, and 1D deployment. Geoscientific Model Development, 15, 75–104. doi:10.5194/gmd-15-75-2022.

Friday, July 18, 2025

Integrating GLEAM Earth Observation data strategy usage within GEOframe

To better understand the previous discussions, let's examine the specific case of GLEAM 3.0 (Martens et al., 2017) and 4.0. GLEAM (Miralles et al., 2025) is a global evapotranspiration product built on multiple Earth Observation (EO) resources. Since evapotranspiration (ET) cannot be measured directly, it must be inferred or modeled from available data. AS you know, GEOframe is our system for doing hydrology. For GEOframe methodologies in catchments application, please see this previous post.

Actually GLEAM operates at 0.1-degree spatial resolution (approximately 11 km grid cells), which is adequate for global analyses but insufficient for our purposes. Our objective requires information at 1km (approximateli 0.01-degree resolution), particularly for applications in complex terrain such as the Alps, where significant topographic variation occurs within very small areas.

From Miralles et al. 2025. References can be recovered there


Earth Observation Resources in GLEAM

According to Miralles et al. (2025), GLEAM utilizes the EO and reanalysis resources listed in Table 1. Below is an analysis of how each resource could be integrated with GEOframe:

Radiation

Current GLEAM approach: 0.1-0.5 degree resolution GEOframe implementation: Point-wise calculations using local solar radiation, filtered through cloud interception and atmospheric scattering models, then topographically corrected using digital elevation models.
Limitations: The coarse satellite resolution is inadequate for rugged terrain. Additionally, GEOframe's current empirical methods lack reliability, often requiring radiation estimates from randomly selected points within catchments, reducing representativeness.
Potential improvements: Satellite data could be fused with ground measurements to enhance overall accuracy.

Air Temperature

GEOframe implementation: Kriging interpolation with drift using ground station data.
Assessment: Generally reliable since temperature varies gradually across space, though canopy effects may introduce complications.

Precipitation

GEOframe implementation: Ground station measurements interpolated using kriging (typically without drift, as drift was found insignificant). Event-specific variograms are employed.
Key challenge: Determining whether kriging interpolation accurately captures storm spatial patterns. Satellite and radar data could provide valuable validation and improvement opportunities. See also the last post here.

Wind Speed

Application: Required for Penman-Monteith formulations.
Potential integration: ERA5 reanalysis data could supplement ground station measurements through data fusion/assimilation approaches, pending reliability validation (see Azimi et al., 2025)

Vapor Pressure Deficit (VPD)

Current use: Input for Penman-Monteith solutions (D'Amato and Rigon, 2025).
Technical note: VPD represents the temperature difference between emitting surfaces and air. Understanding EO estimation methods could enable valuable comparative analyses with the Prospero model (Bottazzi et al., 2021; D'Amato et al, 2025), where VPD emerges from energy budget calculations.

Carbon Dioxide Concentration

Application: Controls transpiration conductance in both Jarvis and Ball-Berry-Leuning parameterizations.
Current status: Available as input parameter in GEOframe for Penman-Monteith and Prospero models but not utilized in Priestley-Taylor formulations.

Snow Water Equivalent (SWE)

GEOframe approach: Calculated from precipitation, temperature, and snowpack evolution models. It should not be confused with Snow Covered Area (SCA).
Improvement opportunities: MODIS snow products offer superior resolution compared to GLEAM's 25 km resolution. For mountainous terrain with 5000 m elevation changes within 25 km, higher-resolution products could provide significant improvements.

Surface Soil Moisture

GEOframe implementation: Prognostic variable within root zone compartment.
Integration potential: Could enable GEOframe calibration if EO resolution and reliability improve. GEOframe soil moisture could be upscaled to match EO data resolution for comparative analysis.

Vegetation Optical Depth (VOD)

Definition: Proxy for vegetation biomass and cumulative transpiration (assuming linear correlation).
Current status: Not implemented in GEOframe.
Integration potential: Could be connected to Leaf Area Index (LAI), which is used in both interception and transpiration calculations.

Fraction of Absorbed Photosynthetic Radiation (fPAR)

Current status: Not used in GEOframe, which employs total radiation instead.
Advantage: Available at appropriate spatial scales for potential integration.

Leaf Area Index (LAI)

Applications: Useful for both interception and transpiration calculations (when using Penman-Monteith or Prospero models).
Integration potential: High, given its direct relevance to existing GEOframe processes.

Vegetation Height

Value: Excellent spatial resolution and direct application in aerodynamic resistance calculations.
Integration status: Not currently used but could be easily incorporated into GEOframe.

Land Cover Fraction

Current status: Not implemented in GEOframe.
Potential application: Could enhance transpiration estimations.

Soil Properties

Current status: Not utilized, as no GEOframe parameters currently depend on soil characteristics.
Future applications: Could become relevant if replacing reservoir-based root zone approaches with simplified versions of WHETGEO (Tubini and Rigon, 2022) or GEOSPACE (D'Amato and Rigon, 2025b).

GLEAM Methodological Components and GEOframe Integration

Rainfall Interception

GLEAM approach: Utilizes the van Dijk-Bruijnzeel model (van Dijk et al., 2001), developed from global experimental datasets.
GEOframe current implementation: Uses the Gash model.
Integration opportunity: Adding a van Dijk-Bruijnzeel component to GEOframe could enhance model compatibility and performance. Notably, Zhong et al. (2022) successfully constrained interception estimates using fPAR, providing valuable insights for EO integration.

Potential Evapotranspiration

GLEAM approach: Employs the Penman equation for potential ET estimation.
GEOframe compatibility: This methodology is already available in GEOframe, enabling direct reproduction of GLEAM's approach. However, D'Amato et al. (2025) implement a more sophisticated Penman-Monteith formulation than GLEAM, allowing for comparative analyses.
Aerodynamic conductance: GLEAM uses Thom's equation, which differs from GEOframe's current formulation but could be easily implemented. Both approaches require roughness length and zero displacement height parameters that can be derived from EO vegetation retrievals.

Soil Moisture Integration

GLEAM4 advancement: Incorporates data assimilation using European Space Agency (ESA) Climate Change Initiative (CCI) surface soil moisture data through a Newtonian Nudging scheme. The method decomposes soil moisture into anomalies and computes uncertainties using triple collocation (Miralles et al., 2025).
GEOframe current approach: Relies solely on root zone reservoirs for ET sources.
Enhancement opportunities:
  • Adding ET sources from groundwater reservoirs through minor modifications to GEOframe's groundwater component
  • Implementing GLEAM4's multi-layer running water balance approach, which considers constant root depth per land cover fraction
Physical realism considerations: While GLEAM4 moved from reservoir models (similar to GEOframe) to multi-layer approaches, GEOframe could implement GEOSPACE (D'Amato et al., 2025) to achieve superior physical realism compared to GLEAM4. However, the reliability of satellite-derived soil moisture estimates requires careful validation.

Stress Function Formulations

GEOframe current approach: GEOET (GEOframe's transpiration component) employs empirical schemes following Jarvis or Ball-Berry-Leuning (BBL) parameterizations.
GLEAM4 innovation: Introduces an innovative deep neural network approach replacing traditional semi-empirical stress computations. As described by the authors: "GLEAM4 replaces the original semi-empirical computation based on soil moisture and vegetation optical depth (VOD) with the deep neural network approach presented in Koppa et al., 2022."
Neural network advantages: The approach recognizes that actual-to-potential transpiration ratios are controlled by numerous environmental variables with non-linear interactions, including:
  • Soil moisture and VOD
  • Vapor pressure deficit (VPD)
  • Incoming solar radiation (SWi)
  • Air temperature (Ta)
  • CO2 concentration
  • Wind speed (u)
  • Leaf Area Index (LAI)
Training methodology: The neural network learns universal transpiration stress functions using global eddy-covariance and sap flow data, with separate parameterizations for tall and short vegetation.
Implementation potential: Incorporating this neural network approach could represent a significant alternative for GEOframe, moving beyond traditional empirical formulations to data-driven, physically-informed methods.

Implementation Considerations

Full understanding and implementation of these methodological improvements requires careful examination of Miralles et al. (2025) and its supporting literature. The integration of these approaches could substantially enhance GEOframe's capabilities while maintaining compatibility with global EO datasets.

Conclusions

While GLEAM provides a comprehensive framework using multiple EO resources, significant opportunities exist for improving spatial resolution and integrating these datasets with process-based models like GEOframe. The main challenges involve resolution limitations and the need for validation of empirical methods against ground-truth data. The methodological advances in GLEAM4, particularly the neural network-based stress functions and improved data assimilation schemes, offer promising directions for enhancing GEOframe's predictive capabilities. Overall GLEAM should not be considered as a EO dataset but a modeling product. Describing it as an EO product makes it more objective that it is actually. Other global products on ET are available. Ecostress (Pierrat et al., 2025) is a recent notable example. A scrutiny similat to the one applied to GLEAM can be made with that platform too, but I let you as an exercise.

This post is part of the dissemination material of the Space It Up project funded by the Italian Space Agency, ASI, and the Ministry of University and Research, MUR, under contract n. 2024-5-E.0 - CUP n. I53D24000060005.


References


  • Azimi, Shima, Christian Massari, Gaia Roati, Silvia Barbetta, and Riccardo Rigon. 2025. “A New Tool for Correcting the Spatial and Temporal Pattern of Global Precipitation Products across Mountainous Terrain: Precipitation and Hydrological Analysis.” Journal of Hydrology 660 (133530): 133530. https://doi.org/10.1016/j.jhydrol.2025.133530.
  • Bottazzi, M., M. Bancheri, M. Mobilia, and G. Bertoldi. 2021. “Comparing Evapotranspiration Estimates from the Geoframe-Prospero Model with Penman–Monteith and Priestley-Taylor Approaches under Different Climate Conditions.” WATER. https://www.mdpi.com/2073-4441/13/9/1221.
  • D’Amato, Concetta, and Riccardo Rigon. 2025. “Elementary Mathematics Helps to Shed Light on the Transpiration Budget under Water Stress.” Ecohydrology: Ecosystems, Land and Water Process Interactions, Ecohydrogeomorphology 18 (2). https://doi.org/10.1002/eco.70009.
  • D’Amato, Concetta, Niccolò Tubini, and Riccardo Rigon. 2025. “A Component Based Modular Treatment of the Soil-Plant-Atmosphere Continuum: The GEOSPACE Framework (v.1.2.9).” https://doi.org/10.5194/egusphere-2024-4128.
  • D’Amato, Concetta, Niccolò Tubini, and Riccardo Rigon. 2025. “A Component Based Modular Treatment of the Soil-Plant-Atmosphere Continuum: The GEOSPACE Framework (v.1.2.9).” https://doi.org/10.5194/egusphere-2024-4128.
  • Dijk, A. I. J. M. van, and L. A. Bruijnzeel. 2001. “Modelling Rainfall Interception by Vegetation of Variable Density Using an Adapted Analytical Model. Part 1. Model Description.” Journal of Hydrology 247 (3–4): 230–38. https://doi.org/10.1016/s0022-1694(01)00392-4.
  • Koppa, Akash, Dominik Rains, Petra Hulsman, Rafael Poyatos, and Diego G. Miralles. 2022. “A Deep Learning-Based Hybrid Model of Global Terrestrial Evaporation.” Nature Communications 13 (1): 1912. https://doi.org/10.1038/s41467-022-29543-7.
  • Martens, Brecht, Diego G. Miralles, Hans Lievens, Robin Van Der Schalie, Richard A. M. De Jeu, Diego Fernández-Prieto, Hylke E. Beck, Wouter A. Dorigo, and Niko E. C. Verhoest. 2017. “GLEAM v3: Satellite-Based Land Evaporation and Root-Zone Soil Moisture.” Geoscientific Model Development 10 (5): 1903–25. http://www.geosci-model-dev-discuss.net/gmd-2016-162/.
  • Pierrat, Zoe Amie, Adam J. Purdy, Gregory Halverson, Joshua B. Fisher, Kanishka Mallick, Madeleine Pascolini-Campbell, Youngryel Ryu, et al. 2025. “Evaluation of ECOSTRESS Collection 2 Evapotranspiration Products: Strengths and Uncertainties for Evapotranspiration Modeling.” Water Resources Research 61 (6). https://doi.org/10.1029/2024wr039404.
  • Miralles, Diego G., Olivier Bonte, Akash Koppa, Oscar M. Baez-Villanueva, Emma Tronquo, Feng Zhong, Hylke E. Beck, et al. 2025. “GLEAM4: Global Land Evaporation and Soil Moisture Dataset at 0.1 Resolution from 1980 to near Present.” Scientific Data 12 (1): 416. https://doi.org/10.1038/s41597-025-04610-y.
  • Zhong, Feng, Shanhu Jiang, Albert I. J. M. van Dijk, Liliang Ren, Jaap Schellekens, and Diego G. Miralles. 2022. “Revisiting Large-Scale Interception Patterns Constrained by a Synthesis of Global Experimental Data.” https://doi.org/10.5194/hess-2022-155.

Tuesday, June 3, 2025

OMS Runner Library: Streamlining Hydrological Model Execution

 The OMS Runner Library v1.2.2 represents a significant advancement in hydrological modeling workflow automation, specifically designed to simplify the execution of OMS3 (Object Modeling System) simulations. For hydrologists and water resources engineers working with GEOframe and OMS3, this Python library addresses the seamless integration and execution of simulation models across different computing platforms. What follows assume a lot of knowlege that you can get by looking to some of our Winter Schools or some of our lab classes as  Physical Hydrology (in Italian) or  Biosphere Atmosphere and Climate Interactions. 

What is OMS3?

The Object Modeling System (OMS3) is a Java-based framework widely used in environmental and hydrological modeling. It provides a robust platform for developing, coupling, and executing complex simulation models. However, working with OMS3 often requires dealing with Java classpaths, configuration files, and platform-specific execution commands – tasks that can be time-consuming and error-prone, especially for researchers focused on scientific analysis rather than software engineering.

The Solution: Python Integration

The OMS Runner Library bridges this gap by providing a comprehensive Python interface for OMS3 operations. This is particularly valuable because Python has become the lingua franca of scientific computing, with most hydrologists already familiar with its ecosystem of tools like pandas, matplotlib, and Jupyter notebooks.

The library automatically handles the complexities of Java environment detection, ensuring that Java JDK 11 is properly configured across Windows, macOS, and Linux systems. This cross-platform compatibility is crucial for research teams working in diverse computing environments, from field laptops running Windows to high-performance computing clusters running Linux.

Please find:

Version 1.2.4

Version 1.2.2

Key Capabilities

One of the library's standout features is its intelligent simulation management. It can automatically discover simulation files within a project, maintain configuration databases, and execute models either individually or in sophisticated batch processing workflows. For hydrologists working with multiple scenarios – such as climate change impact assessments or calibration procedures – the parallel execution capabilities can reduce computational time.

The library supports various execution patterns: sequential processing for dependent simulations, parallel execution for independent model runs, and asynchronous background processing for long-running computations. This flexibility allows researchers to optimize their workflows based on available computational resources and modeling requirements.

Practical Applications

In practical hydrological applications, this translates to significant productivity gains. A researcher studying watershed responses to different precipitation scenarios can now set up dozens of model runs with just a few lines of Python code, monitor their progress through Jupyter notebooks, and automatically collect results for analysis. The library's integration with popular Python data analysis tools means results can be immediately processed, visualized, and shared.

Users can explore more about GEOframe's capabilities and latest developments at the GEOframe blog, where detailed tutorials and case studies demonstrate advanced hydrological modeling workflows.

The comprehensive logging and error handling features are particularly valuable in operational hydrology contexts, where model reliability and traceability are paramount. The library maintains detailed execution histories, facilitates debugging, and provides clear diagnostic information when issues arise.


Saturday, January 25, 2025

GEOSPACE or Soil-Plants-Atmosphere-Continuum Estimator in GEOframe first paper

The soil-plant-atmosphere continuum (SPAC) system is a complex and interconnected network of physical phenomena, encompassing heat transfer, evapotranspiration, precipitation, water absorption, soil water flow, substance transport, and gas exchange. These processes govern the exchange of energy, matter, and water within the SPAC system. To better understand and model SPAC interactions, interdisciplinary approaches are essential due to the inherent complexity of the system. Instead of relying on a single monolithic model, we propose a component-based modeling approach, where each component addresses a specific aspect of the system. Object-oriented programming (OOP) is adopted as the foundational framework for this approach, providing flexibility and adaptability to accommodate the ever-changing nature of the SPAC system.

Please find the paper by clicking on the Figure

The Soil Plant Atmosphere Continuum Estimator in GEOframe (GEOSPACE) is presented in this paper, in particular the one-dimensional development, GEOSPACE-1D. The framework is a tool designed to facilitate robust, reliable and transparent simulations of SPAC interactions. It embraces the principles of open-source software and modular design, aiming to promote open, reusable, and reproducible research practices. By implementing the OOP, GEOSPACE-1D breaks down the complexity of SPAC modeling into smaller, self-contained structures, each responsible for a specific scientific or mathematical concept. This modular architecture adheres to the "open to extensions, closed to modifications" philosophy, enabling easy model extension without disrupting existing components. Equations are implemented in an abstract manner, emphasizing the use of common interfaces over concrete classes, a hallmark of contemporary OOP. GEOSPACE-1D adopts a generic programming framework, where distinct classes adhere to a common interface. This compartmentalization serves two critical purposes: validating individual processes against analytical solutions and facilitating the integration of novel processes into the system.

The paper emphasizes the significance of modeling the coupling between infiltration and evapotranspiration for accurate hydrological simulations. It explores the interplay between plant transpiration, soil evaporation, and soil moisture dynamics, highlighting the need to account for these interactions in SPAC models. The paper concludes by underlining the importance of modularity, transparency, and openness in SPAC modeling, principles that underlie the development of GEOSPACE-1D and its components. Overall, GEOSPACE-1D represents a promising approach to SPAC modeling, providing a flexible and extensible framework for studying complex interactions within the Earth's Critical Zone. It is worth recalling that the fundamental premise of GEOSPACE-1D is not to create a single soil-plant-atmosphere model, but to establish a system that allows the creation of a series of soil-plant-atmosphere models, adapted to the specific needs of the user's case study.

Friday, September 13, 2024

A new tool for correcting the spatial and temporal pattern of global precipitation products across mountainous terrain: precipitation and hydrological analysis


This study primarily aims to integrate global precipitation data into hydrological models at the catchment scale, a common practice in hydrological research. Specifically, the study investigates how biased spatial and temporal patterns in precipitation data affectmodel performance and uncertainty. The European Meteorological Observations(EMO) and Climate Hazards Group InfraRed Precipitation with Station data (CHIRPS) global datasets are utilized as inputs for the GEOframe-NewAGE hydrological model to simulate the hydrological processes of the mountainous Aosta Valley catchment in northwestern Italy. The uncertainty of the hydrological model forced with global precipitation data is assessed using a proposed method called Empirical Conditional Probability (EcoProb). The results show that, although traditional performance metrics suggest similar outcomes for the model forced with EMO and CHIRPS, the proposed uncertainty analysis reveals higher uncertainty when CHIRPS is used as the precipitation input. To leverage all useful information in the global precipitation data, the spatial correlation of CHIRPS was combined with a subset of raingauges using the EcoProb method to modify the EMO precipitation data. This approach enabled the integration of the advantages of EMO and CHIRPS, which offer higher temporal and spatial correlation with ground observation, respectively, into a unified precipitation product. The combined dataset, referred to as the EcoProbSet product in this study, outperformed both the CHIRPS and EMO products, reducing the uncertainty introduced into hydrological models compared to the original global datasets.

You can find the paper preprint by clicking on the Figure above. 

Monday, August 5, 2024

Mumbai GEOframe School !

 We have just completed our effort with the GEOframe Mumbai Monsoon School, inserted in a larger initiative, of the GISE HUB which included one day long SCPP workshop on "Recent Advances in Hydrological Modelling" on 31st July. Besides being trained on GEOframe, hands on training on Dynamic Budyko model was provided by prof. Basudev Biswal (GS) and his postdoc  Prashant Istalkar. Lectures on the 31st July covered a wide range of topics including flood inundation modelling, socio-hydrology, land-surface modeling, climate-change impact assessment, machine learning models, and complex networks.


Great thanks to
Sumit Sen and Basudev Biswal for organizing the School. Hospitality was superb, discussions enriching and seeing the dedication and smartness of students an encouraging academic experience. We hope that the School will have follows up both at IIT and UniTrento and exchanges could continue in the future. For further information, see also the Linkedin post by Basudev here
The GEOframe material of the School is available to anyone and the slides and videos (when uploaded) will be available at the GEOframe blog page dedicated to the School.
The success of the School, from our side, is the outcome of many that are listed in this "people of GEOframe" presentation available here
For students who want to complete a personal exercise with GEOframe, the GEOframe team is available to assist. Upon completion, each student will receive a University of Trento T-shirt.

Sunday, June 9, 2024

On catchment analysis (modeling)

In a series of papers (Abera et al., 2016, Abera et al. 2017a, Abera et al., 2017b, Azimi et al., 2023), we have outlined a methodology for studying basins, focusing on specific locations BUT looking especially to the methodologies. They are also summarized in slides that I typically use in my hydrological modeling classes. These slides summarize the analysis requirements in seven key steps, supported by various notebooks that implement the methodologies.


Each time we begin a new catchment analysis, please ensure these methodological suggestions are considered. Overlooking them can be quite frustrating. Consistently revisit and apply the reference material to build upon previous work and past achievements. Criticize previous methods if necessary, but do not disregard them.
There are two critical steps that are often neglected. The first is data analysis—specifically, examining data coherence and comparing multiple data sources. This preliminary analysis can provide significant insights before any modeling begins, but it is rarely pursued. Instead, input data are directly used in the model, leading to issues later because something seems off.
The second neglected step is validation. There is a tendency to be satisfied with performance metrics like KGE or Nash-Sutcliffe, but these should be starting points, not final assessments. Other benchmarks, such as those proposed by Addor et al., (2018) should be used to critically evaluate the results, not just applied mechanically.
Recently, Azimi et al., 2023 introduced a more refined analysis method (called in future papers EcoProb), which allows for finer discrimination of model behavior by separating the ranges of response. This method should be considered for more precise analysis.
Additionally, since Abera (and likely earlier), we have tried to refocus our analysis on not just discharges but also on budgets. Understanding budget behavior can reveal significant insights and prevent errors, but it is often sidelined. We need to improve in this area.
Mapping is another crucial aspect. While we often rely on time series plots, spatial representation is essential to show the irreducible spatial heterogeneity. This is evident in soil moisture studies, such as the recent work by Andreis et al., and should also apply to other quantities like snow cover and depth.
I have worked with many of you to create effective graphs and maps. Using a full range of colors is beneficial, but please remember that some journals (AGU and EGU) require color-blind friendly plots. Address this requirement from the beginning to avoid last-minute modifications.

P.S. I - One distinguishing feature of GEOframe compared to other systems is its ability to explore multiple working hypotheses. Although this capability exists, it has not been utilized so far. Let's make full use of it moving forward.

P.S. II - When I read a paper from collaborators, I assume that all materials, including data, software, notebooks, and .sim files, are organized and shared as supplemental material for reproducibility. To achieve this, it is crucial to maintain order and keep the material up-to-date from the beginning. Otherwise, it becomes a nightmare.

References

Abera, Wuletawu, Luca Brocca, and Riccardo Rigon. 2016. “Comparative Evaluation of Different Satellite Rainfall Estimation Products and Bias Correction in the Upper Blue Nile (UBN) Basin.” Atmospheric Research 178-179 (September): 471–83. https://doi.org/10.1016/j.atmosres.2016.04.017.

Abera, Wuletawu, Giuseppe Formetta, Luca Brocca, and Riccardo Rigon. 2017. “Modeling the Water Budget of the Upper Blue Nile Basin Using the JGrass-NewAge Model System and Satellite Data.” Hydrology and Earth System Sciences 21 (6): 3145–65. https://doi.org/10.5194/hess-21-3145-2017.

Abera, Wuletawu, Giuseppe Formetta, Marco Borga, and Riccardo Rigon. 2017. “Estimating the Water Budget Components and Their Variability in a Pre-Alpine Basin with JGrass-NewAGE.” Advances in Water Resources 104 (June): 37–54. https://doi.org/10.1016/j.advwatres.2017.03.010.

Addor, N., G. Nearing, C. Prieto, A. J. Newman, N. Le Vine, and M. P. Clark. 2018. “A Ranking of Hydrological Signatures Based on Their Predictability in Space.” Water Resources Research 54 (11): 8792–8812. https://doi.org/10.1029/2018wr022606.

Azimi, Shima, Christian Massari, Giuseppe Formetta, Silvia Barbetta, Alberto Tazioli, Davide Fronzi, Sara Modanesi, Angelica Tarpanelli, and Riccardo Rigon. 2023. “On Understanding Mountainous Carbonate Basins of the Mediterranean Using Parsimonious Modeling Solutions.” Hydrology and Earth System Sciences 27 (24): 4485–4503. https://doi.org/10.5194/hess-27-4485-2023.

Thursday, May 16, 2024

GEOframe-New AGE material for beginners

Dear User or Dear Explorer,

Here we aim to summarize some of the material related to GEOframe-NewAGE. The main source is certainly the GEOframe blog:

However, for a logical introduction, it may be useful to start here:



The most recent material on GEOframe-NewAGE is from the latest school, accessible from this point:

By following the links for each day, you can download the slides and watch the lesson videos. 

Another useful resource is provided through hydrological modeling tutorials:

GEOframe's infrastructure is based on the Object Modelling System v3:

Various papers and applications related to GEOframe have been written and developed; you can find them here:

For any further assistance, the GEOframe crew can be reached at geoframe-schools@googlegroups.com.

Please feel free to reach out  us if you have any questions. Next Winter School  on GEOframe-NewAGE will be in January 2025 from 7 to11 in Trento University.  Next Summer School (on Land-Atmosphere interactions will be June 2-6 2025.  This Summer will be holding one Summer School on GEOframe-NewAGE in El Cairo and one in Mumbay (both the last weeks of July).

Tuesday, December 26, 2023

Code Washing

 This time, I want to address the concerning issue of students inappropriately reusing open-source code without a clear understanding of open-source licenses.



It's crucial for students to grasp the essence of open source licenses, understanding that they are not just permissions to copy but guidelines for responsible use. Engaging with open-source code should involve a genuine learning process, encouraging students to comprehend and apply the principles embedded in the code they explore.
Merely having access to code doesn't grant the right to take it, make superficial changes, or translating from a programming language to another, remove original authors, and claim the altered code as their own. While open source encourages learning through code exposure, wholesale copying with only minor alterations, especially without restructuring for object-oriented code, doesn't constitute "creating a new code base."
In such instances, phrases like 'I looked at Mickey Mouse code, but I am using my own code' are, at the very least, misleading and likely a form of plagiarism. I term this practice "code washing." My plea: steer clear of it and adhere to ethical behavior.
The notion of "code washing" not only undermines the integrity of individual work but also compromises the collaborative spirit of open source. It's essential to emphasize that acknowledging and respecting the original authors not only aligns with ethical standards but also fosters a culture of transparency and collaboration in the coding community.

Friday, November 17, 2023

Some pills on what we do for agriculture droughts

 Just to introduce the debate about droughts simulation, agriculture, new technologies that can be used for improving agriculture. Below the presentation.



Just click on the Figure to see the presentation given at the Festival della Meteorologia 2023


Tuesday, October 10, 2023

Notes about the dynamic nature of the GEOframe-Po Project

Here below you can find some provisional notes, to be improved in the next days about our Deployment of the GEOframe system to the river Po for the basin Authority of the river Po.  

Basin extraction

it's not a straightforward operation. In fact, it has never been done systematically all over Italy. It serves two opposing needs: to be objective and to align with the official grid provided by basin Authorities and Regions. The initial phase relies mainly on slope analysis and requires processing digital terrain data, which have become available only in recent years, especially if we refer to data produced with laser altimetry. The starting point is the Digital Elevation Models (DEMs) provided by the regions, which have been reprojected and standardized to correct reference systems. The initiation of the hydrographic networks is determined by an area threshold, while sub-basins, for the Po river, are delineated to have an average area of 10 km2. Procedures have been standardized in geographic information systems (GIS) over the last twenty years, but for this specific task, the Horton Machine library developed by Univrsity of Trento and HydroloGIS was used (Abera et al., 2016, serving as reference), incorporating some innovative elements: a parser to aggregate smaller basins into adjacent larger ones and addressing certain topological situations, especially those in flat areas for the subsequent use with GEOframe.
The tools was named GEOframeInputBuilder.

The extraction of lakes, particularly the large Lombard lakes and Lake Garda, required special attention and made the process less automated. Visual analysis reveals a differentiated geometry between mountain basins and lowland inter-basins, since the early years of fluvial geomorphology, but now objectively observed. The database, now available, enables statistical analysis of their geometry and topology, which previously relied on more qualitative cartographic analysis. The basin initiation with an area threshold is functional to the hydrological modelling but the reader should be aware that this topic is a very alive hydrological research topic, especially along with the work by Gianluca Botter and coworkers [insert CITATION].

The grid, as currently constructed, will be distributed for free use and will serve as a fundamental standard for further cartographic-digital and hydrological analyses and developments.

Photo by Luigi Ghirri



Interpolation

Interpolation techniques have seen significant development between the 1980s and 90s [insert citation], but especially geostatistical methods have slowly made their way into the practice of digital analysis of meteorological forcings in the hydrological cycle. These require the definition of an estimation model of the correlation between measurements, known as a variogram, the robustness of which is fundamental to the reliability of the result.
The starting database is made up of measurements collected by ground stations from regional entities operating on the Po basin. These data have been analyzed, cleaned, and subsequently interpolated, currently on each centroid of the sub-basins identified in the first phase of the work. The interpolation was carried out for precipitation and temperatures on a daily scale, as a first step to produce hourly or sub-hourly interpolation at any point of a suitable one-kilometer grid.
The interpolation technique used was kriging with drift to account for orographic effects, especially on temperature. For the interpolation of the experimental variogram, a ?linear? Exponential? What else? model was used using the interpolators implemented in GEOframe.
The interpolation covered the entire period from 1990 to today, and the data are stored in CSV files in folders containing the data for each individual sub-basin.

It is clear that the procedure is a first approximation that will serve as the basis for future improvements. For example, the extension of the interpolation on the one-kilometer grid is one aspect. The next improvement could be to introduce high-resolution reanalysis data, combining geostatistical techniques with simulations of atmospheric circulation and any data coming from radar and satellite. Convergent research come from atmospheric physics and meteorology whose resolution is arrived at the scales useful for hydrology. Some work should be done for connecting better the two communities.

Setup:

GEOframe-NewAGE allows numerous configurations, as various components are available for the same phenomenon. For the basic configuration of each single Hydrologic Response Unit (HRU), the one already partially tested in [insert citation] called Embedded Reservoir Model (ERM) was chosen, the description of which can be found in the cited bibliography or in the linked videos. In summary, the ERM model is composed of a component for interception, one for snow, when present, a fast surface runoff separator based on the Hymod model, a nonlinear reservoir for the description of the root zone, and a second nonlinear reservoir for groundwater. Structurally, it is not much different from the HBV Model [insert citation]. In the basic configuration, flood propagation is neglected.
For the part of evapotranspiration, a simple Priestley-Taylor model was used, where however the radiation is provided through a rather accurate model [insert citations].
Each of these ERM models was then connected to the others through the Net3 infrastructure [insert citation] to form a directed acyclic graph in which each node represents an HRU. Potentially, each HRU can be characterized not only by its own topographic and meteorological data, but also by its own models.
In the basic configuration, however, the same model structure is usually used for all HRUs while the values of the model parameters are obtained by subsequent calibration with spatially differentiated parameters, if the available data allow it.
The potential setup variants are numerous, encompassing at least three options for snow modeling, three for evapotranspiration modeling, and an array of choices for reservoir modeling. The inclusion or exclusion of flow propagation modules, as well as the potential elimination or addition of compartments to be modeled and their diverse connections, further expand the possibilities. An overview of potential topological configurations is presented, for instance, in [insert MaRmot citation]. As even a novice reader can comprehend, the possible combinations multiply far beyond exponentially with the number of connected Hydrological Response Units (HRUs), which can, in turn, be linked in various manners. This complexity underscores why our comprehensive study on the Po River necessitates distribution and further refinement by others to enhance the precision of the results and better align them with local needs which cannot be gained by a single yet very productive team of people. In turn this open the question on how the re-analysis performed by external researchers or teams can be accepted and inserted back into the main project.

The analysis of multiple configurations is therefore entrusted to later phases of the project.

Calibration

Among the phases of a simulation, the calibration phase is the most time-consuming. It essentially consists of a large number of attempts to combine the model parameters to reproduce the measured data as faithfully as possible. The space of possible parameters is generally very large, even for a single simulation HRU. Therefore, the tools for calibration try to use intelligent strategies (including ML) to quickly guess which are the best parameter configurations.

The goodness of the simulated values' fit to the measured ones is usually quantified through some goodness of fit (GOF) algorithms. In our case, these are generally the KGE [insert citation] or the NS [insert citation]. An analysis of the various GOFs can be found in [insert citation], whose result can be further detailed, in the validation phase (see below), with additional indicators such as those presented, for example, in Addor et al., 2017. Another method of analysis, post-hoc of the goodness of the simulations, much more refined, is that presented in [insert Shima work citation]. The latter can also serve as a Bias corrector of the final result and it is going to be systematically applied to the results of the Po project.

From an algorithmic point of view, the calibration carried out in the project is based on the LUCA model [insert citation], which is a sophisticated implementation of SCEM-UA [insert citation], but a particle swarm [insert citation] could also be used. The calibration procedure follows some standards. Having a set of data to base the calibration on, the data are usually divided into two subsets, one used for calibration and another for the so-called validation phase. In the former, the problem of having available input and output data is solved, determining the parameters (or models) in a way similar to what is done in normal ML techniques (which, for this purpose, could probably be used profitably). In the latter, the performance of the model solution on data not used for parameter determination (and should be "independent" of the former) is evaluated. As already mentioned, in the validation phase, additional GOF indicators can be used to better discern the performance of the adopted solution.

A note concerns the word "validation". This is the term used but does not imply any ontological meaning about the nature of the truth described by the model, but only a practical meaning related to the reliability of the model in predicting a certain sequence of numbers.
The calibration/validation procedure can be implemented for a single variable, in the specific case, usually the flow in a section of the hydrographic network, or for more variables, for example, snow cover, soil water content, evapotranspiration, if these measurements are available. These latter possible measures, however, have a different character from the discharge as, while discharge is an aggregate variable, resulting from the concentration of the fallen water on the watershed area in a single point, the others remain variables distributed spatially, before being aggregated for the purposes of the watersheds budget, and therefore the methods of determining the goodness of reproduction of the measured data follow more articulated paths, if not more complex. The good thing is that GEOframe allows you to calibrate the various quantities separately, as each of them is modeled by "different components" that can be used separately from the overall model. The use case is performed throufh quite a lot of manual intervention so far and could be made more automatic.

In any case, if the target variables are more than one, we speak of multi-objective calibration, while if there are variables measured at multiple sites, we speak of multi-site calibration [insert citation].

I would like further to suggest an enhancement to our analysis and move from the daily to hourly time scale. This is particularly crucial for understanding processes within smaller watersheds, approximately on a 1km^2 scale, where many significant phenomena demonstrate sub-daily dynamics.


Simulation/ Analysis/ECP

The validation phase is already a simulation stage (with predetermined parameters) and represents the normal completion of operations in a production phase. This production phase is usually understood in the hydrological literature as hindcasting, that is, as functional to the understanding of past events for which an explanation is sought in a quantitative framework. This involves the use of more accurate analysis and indicators than those used in the calibration/validation phase which require a certain speed. One of these is the analysis through empirical conditional distributions, as illustrated in Azimi et al., 2023. These analyses can eventually lead to a rethinking of the setup and calibration/validation phases in order to obtain more accurate results. As shown in Azimi et al (2023, 2024), ECPs can also be used as bias correctors and improve the overall statistical performance of the model's results, at least if it shows a certain stationarity of temporal behavior, that is, if, for example, the effects attributable to global warming do not significantly impact the structure of the model (including its parameters). The determination of the "reliability" of the models is then a key concept in the development of digital twins of the hydrological system (Rigon et al, 2022).

Another matter, and much less frequented by hydrologists, is that of forecasting future events. These future events, obviously, have to do with the time series input to hydrological models and therefore require forecasts of precipitation, temperature, wind speed, and air humidity. It is known that the meteorological system (global and local) is affected by a lack of predictability (predictability) due to deterministic chaos effects [insert citation]. To date, weather predictions have reliability, with respect to weather categories, of a few days, they have the ability to predict atmospheric temperatures, but they are still very imprecise in determining the amount of precipitation, in essence, they can be used to predict the hydrological future but with questionable quantitative value. The theoretical reason for this debacle has been somewhat said, but there are also others, for example, the heterogeneity of ground conditions and the absence of a description of the soil-atmosphere feedbacks, both conditions not described in meteorological models. Hydrological forecasts can therefore only be of a statistical nature and produce scenarios [insert citation], which are not devoid of meaning and practical value. In this area between Hydrology and meteorology the search for a common ground is mandatory for any evolution. In GEOframe, however, the input data treatment/modelling is quite well separated from the hydrological computation and any new source of data can be easily (but not without person/months work) included.

Distribution of results and participatory science

A fundamental aspect, already widely discussed in Rigon et al., 2022, is to understand how the results of a model can be shared with various users, but also how the model, its setup (including a very expensive phase of DEM analysis, weather data interpolation, and calibration/validation) can be shared, saving other researchers time. GEOframe is built in such a way that this is possible (share ability is by design of the informatics) and some experiences have already been made in this sense. Some within the Trento working group, others with external research groups from the University of Milan (whose work is to be incorporated) and the Polytechnic of Turin, where the basic data and models already pre-digested by the University of Trento served for further developments and analysis on some parts of the Po basin already processed.
The question on how to preserve, make use of multiple contributions to code, data, simulation configurations and simulations, is still open though.
It should be clarified that the GEOframe system is not only a set of data and models, but also a library of analysis tools, especially developed through Python Notebooks and often documented through a series of slides and video lessons [add the links here] and Schools [https://abouthydrology.blogspot.com/2021/10/the-geoframe-schools-index.html]. Although this system can be improved and automated, it has allowed the group from the Polytechnic of Turin to dramatically shorten the modeling times of a series of basins in Piedmont and will allow, for the moment in the planning stage, the sharing of the setup and analysis of the Lombard area of the large Alpine lakes. Other analyses, developed in parallel on areas such as Friuli by the University of Udine, can easily be inserted into a possible national system that covers all of Italy, even though they were developed separately.
From the informatics point of view organizing all of this information through appropriate repositories would be mandatory in the future for an effcient use of the resources.

Conclusions

The GEOframe-Po project is more than just a collection of models; it envisions a comprehensive system that encompasses a variety of input and output datasets, model configurations, and the flexibility to operate on diverse platforms such as laptops, servers, and the cloud (leveraging the OMS/CSIP platform). The interfaces, as evidenced by the available self-instruction materials, can range from simple text-based designs to more sophisticated visual tools, including augmented reality devices.
The system is designed for continuous improvement and customization, with the ability to implement changes with minimal overhead. This was a strategic requirement pursued at various levels of the information technology aspect of the project [insert citations]. The current models can be broadly categorized as physically based, with the majority of the implementation comprising what is referred to in literature as "lumped" models. However, the system is designed to accommodate extensions to more distributed models, a possibility that has already been partially realized in some research lines of the Trento group.
The integration of machine learning techniques into the system is also possible [insert citation], even though they have not been utilized to date. The design of the GEOframe-Po project, therefore, represents a flexible, adaptable, and forward-thinking approach to modeling and data analysis.