-
Notifications
You must be signed in to change notification settings - Fork 18
Description
So far, HydPy misses a standard way to define spatial information. Here are some first ideas.
Where do we need spatial information so far?
L-Land Version 3 requires the latitutde and the longitude in decimal degrees for calculating extraterrestrial radiation. Each spatial location corresponds to the area centre of a sub-catchment. Each subcatchment is represented by at least one Element object.
The same holds for Evap Version 1.
All Conv models need information on multiple locations. Parameter InputCoordinates defines the coordinates of all relevant input points, for example of different rain gauges. Parameter OutputCoordinates defines the coordinates of all relevant output points, for example of the area centres of all relevant sub-catchments. The unit is arbitrary, but accurate interpolations usually require a projection system like UTM. Each input or output location is mapped to a Node object.
Conv Version 3 also requires additional information for each input and output location. As the names of the related parameters suggest (InputHeights and OutputHeights), this is often the terrain height. But any other type of information suitable to improve interpolation results is also fine.
Where do we need spatial information soon?
For our OpenDA adapter. Algorithms like the Ensemble Kalman Filter require a spatial noise model. And this model requires the location of certain "parts" of a catchment (e. g. the centres of sub-catchments, hydrological response units or river segments), usually defined within a separate OpenDA XML file. Handling the locations separately from HydPy's project files seems error-prone. Hence, we prefer to implement an HTTP request allowing the adapter to ask the HydPy server for the correct locations.
What do we need to discuss?
- Where is the best place to define location information?
- How much additional support is required?
A very first draw.
We could implement a Location class, allowing us to define x and y and, optionally, z and an epsg code for each Node and Element object. The following (not fully functional) example assumes such a class already exists:
from hydpy import Element, Location, Node
Node("rain_gauge",
location=Location(x=1, y=2, z=3, epsg=25832))
Node("rain_catchment",
location=Location(x=4, y=5, z=6))
Node("outlet_catchment",
location=Location(x=7, y=8))
Element("interpolator",
inlets="rain_gauge",
outlets="rain_catchment")
Element("catchment",
inputs="rain_catchment",
outlets="outlet_catchment",
location=Location(x=9, y=10, z=11))Repeating the same EPSG code multiple time is messy. One way to avoid this would be to provide a default code as an additional option:
pub.options.defaultepsgcode = 25832One could, for example, change the default code at the top of the relevant network file, when required.
In the above example, the "interpolator" element interpolates the precipitation measured at the location of node "rain_gauge" to the location of node "rain_catchment". Node "rain_catchment" passes the interpolated precipitation to element "catchment". Very likely, we want to use the same location information (the area centre of the modelled catchment) for node "rain_gauge" and element "catchment". To not repeatedly define the same location, element "catchment" could query the information from node "rain_gauge" (location="rain_gauge").
Let's assume we want to model the catchment with application model L-Land Version 3. Currently, its control files start like that:
from hydpy.models.lland_v3 import *
parameterstep("1d")
latitude(50.0)
longitude(10.0)
...L-Land Version 3 could query its latitude and its longitude from element "catchment". However, there are two points to consider:
- It would need to transform the given coordinates (in this case UTM) into decimal degrees. Realising this via a library like pyproj does not look to be complicated. However, we should first check if this breaks compatibility with the
arcpylibrary of ArcGIS Pro. Also, we need to make sure it does not complicate our builds or installing HydPy for some other reasons. - Would we convert the control parameters
LongitudeandLatitudeto derived parameters, so that the spatial location never appears within a control file?
The last point connects with the following two problems.
First, if we want OpenDA to modify, for example, the soil moisture of individual hydrological response units, we would ideally also define the location for each unit. Just to write something down:
from hydpy import Element, Location, Node
Element("catchment",
inputs="rain_catchment",
outlets="outlet_catchment",
location=Location(x=9, y=10, z=11),
sublocations=(Location(x=12, y=13), Location(x=14, y=15)))I think this information is too model-specific for a network file and would hamper its readability too much. Hence, if we actually want to support defining separate locations for individual hydrological response units (or river segments...), we should better add a new file type. Maybe, this file type should also support clarifying which variable (e.g. precipitation Nied, soil moisture BoWa, snow water equivalent WAeS) relates to which set of locations, like these two tables:
| LocationType | Variables |
|---|---|
| area centre | Nied |
| hru | BoWa, WAeS |
| Element | LocationType | Coordinates |
|---|---|---|
| catchment | area centre | x=1, y=2, z=3 |
| catchment | hru | (x=4, y=5), (x=6, y=7) |
Second, we also need the locations of selections. Here it seems natural to calculate the average of the locations of all contained elements automatically, weighted with their areas (in case of a subcatchment) or lengths (in case of a river channel). Maybe we could speak of "size" and always assume its value to be given in km² or km:
Element("catchment",
inputs="rain_catchment",
outlets="outlet_catchment",
location=Location(x=9, y=10, z=11),
size=12)| Element | LocationType | Coordinates | Size |
|---|---|---|---|
| catchment | area centre | x=1, y=2, z=3 | 20 |
| catchment | hru | (x=4, y=5), (x=6, y=7) | 5, 15 |
Or, on the Python side, combine the coordinates and the size and also allow defining its unit, for example:
Element("catchment",
inputs="rain_catchment",
outlets="outlet_catchment",
location=Location(x=9, y=10, z=11, size="12km"))