The station concept in cmf

Meteorological data is usually measured at stations, and can be represented as such in cmf. Cells have the method get_weather(time) and get_rainfall(time) to get the current meteorological observations. If you are doing nothing, this function returns a default Weather. However, if the cell is bound to one or several measurement station (MeteoStation or RainfallStation, the functions return the condition of the bound station or an interpolated value. To simplify this connection, the stations take a position in project coordinates to find the distance of the cell to the stations of the project. The timeseries containing the forcing data are filled from any data source using standard Python.

Rainfall stations

Rainfall stations consists of a timeseries of rainfall intensities in mm/day. Note that the unit is '''always''' mm/day, regardless of the timestep of your timeseries. If you have, for example, hourly measured rainfall data in mm/h, you have to multiply the value by 24 to get the rainfall in mm/day. Since you can calculate with timeseries, you can simply move the unit conversion into your data loading script.

Create a rainfall station

In the next example, we will create rainfall station using random data. For your application you will load your data from a file, database or web service.

import cmf
import numpy as np
from datetime import datetime,timedelta
 
# Create a project
project = cmf.project()
 
# Create a timeseries, starting Jan 1st, 2010 with 1h step, with random data
raindata = cmf.timeseries.from_array(start=datetime(2010,1,1),
                                     timedelta(hours=1),
                                     numpy.random.uniform(0,30,24*365))
# Add a rainfall station to the project
rainstation = project.rainfall_stations.add(Name='Random station', 
                                            Data=raindata,
                                            Position=(0,0,0))

The last command is new: Every project has a list of rainfall stations, and you add a new station with the add method. The method expects a name for the station (very useful if you have more then one), the data as a cmf timeseries and a position in space in project coordinates. You can connect your new rainfall station with a cell c by:

rainstation.use_for_cell(c)

or you can use one of the spatial distibution methods shown at the end of this tutorial.

Meteorological stations

MeteoStation objects are much more complex than rainfall stations. The do not hold only one timeseries, but quite a number. The data is used for the various evaporation, evapotranspiration and transpiration connections available in cmf. Depending on the processes you are using in your model, not every timeseries you can populate is really used. The more physically based atmospheric connections, like PenmanMonteithET and ShuttleworthWallace need all values that are provided by the MeteoStation. However, not all values are measured at all stations. The MeteoStation object follows therefore the FAO guidelines for the Penman Monteith method, to use surrogate values, as listed below.

Timeseries provided by a meteo object:

T:: Timeseries of Temperature \( T\left[^\circ C\right] \) '''Optional''', if missing it is calculated from \( T=\frac{T_{max} + T_{min}} 2 \)

Tmax:: Timeseries of daily maximum Temperature \( T_{max} \left[^\circ C\right] \) '''Required''', but it can be generated from T, if you have a high resolution timeseries for Temperature, Tmax can be generated as follows: meteo.Tmax = meteo.T.reduce_max(meteo.T.begin, cmf.day)

Tmin:: Timeseries of daily minimum Temperature \( T_{min}\left[^\circ C\right] \) '''Required''', but it can generate from T. If you have a high resolution timeseries for Temperature, Tmin can be generated as follows: meteo.Tmin = meteo.T.reduce_min(meteo.T.begin, cmf.day)

Tground:: Timeseries of upper soil temperature \( T_{ground}\left[^\circ C\right] \) '''Optional''', if missing \( T_{ground} = T \)

Windspeed:: Timeseries of windspeed at instrument height (default 2m) in \(m/s\). '''Optional''', if not available the wind speed defaults to 2 m/s

rHmean:: Timeseries of relative mean Humidity \( rH_{mean} \left[\%\right] \) '''Optional''', not used if rHmax or Tdew is available. Use this field, if you have high resolution relative humidity data.

rHmin:: Timeseries of daily minimum relative Humidity \( rH_{min} \left[\%\right] \) '''Optional''', only used if rHmax is used

rHmax:: Timeseries of daily maximum relative Humidity \( rH_{max} \left[\%\right] \)'''Optional''', not used if Tdew is available

Tdew:: Timeseries of dew point temperature \( T_{dew}\left[^\circ C\right]\) '''Optional''', if neither Tdew, rHmax or rHmean is available then Tdew = Tmin

Sunshine:: Fractional sunshine \(\frac{n}{N}\). Sunshine duration per potential sunshine duration. '''Optional''', if not present it defaults to 0.5. If you have the total duration of sunshine, use the SetSunshineFraction method of the MeteoStation. If you have cloud coverage instead of total sunshine duration you may assume that the fractional sunshine duration equals 1-cloudfraction

Rs:: Global Radiation in \( \frac{MJ}{m^2 day} \) ''' Optional''', if not available cmf::atmosphere::global_radiation is used to calculate radiation from lat, lon, day of year and Sunshine.

T_lapse:: Temperature lapse, the slope of the temperature / height regression. Typical values are \( -0.0004 .. -0.001 \frac{^\circ C}{m} \), default is \( 0\frac{^\circ C}{m} \) (no temperature adjusting)

Creating a meteorological station

To keep the example code short, we assume, that you have your data provided in a special binary data format for timeseries. For your model, you would rather do something as it is shown in CmfTutTestData.

meteo = project.meteo_stations.add_station(name='A station name',
                                           position=(0,0,0),
                                           latitude=50.8,
                                           longitude=8.4)
# Load some data
meteo.Tmax = cmf.timeseries.from_file('Tmax.timeseries')
meteo.Tmin = cmf.timeseries.from_file('Tmin.timeseries')
meteo.Rs   = cmf.timeseries.from_file('Rs.timeseries')

Have a look at the first line. To create a meteo station, the name and position are needed as for the rainfall station, but also the position in geographic coordinates and the time zone for the calculation of radiation from sun inclination. If you want to have a daily average radiation, set meteo.daily=True. To include diurnal changes of the radiation (the night is dark) set meteo.daily=False. If you provide Rs data, the geographic position is not used and the radiation follows the given data. To assign a cell you can do meteo.use_for_cell(c) or use a concept shown in the next section.

Distributing data in space

Using the nearest station

If the forcing data for a cell should be used from one station, the station.use_for_cell idiom is fine. To assign all cells of a project to the nearest station write in your script:

project.use_nearest_meteo(z_weight=0.0)

project.use_nearest_rainfall(z_weight=0.0)

The distance is calculated using the x,y,z coordinates,of the cell \(x_c,y_c,z_c\) and the position argument of the station ''i'', \(x_i,y_i,z_i\) is calculated as:

\[ d_i(x,y,z) = w_z \left|z-z_i\right| + \sqrt{\left(x-x_i\right)^2 +\left(y-y_i\right)^2} \]

where \(w_z\) z_weight is the weight of the height difference for the distance calculation. If z_weight is not given, it defaults to 0.0 (pure euclidian hoizontal distance).

Interpolation of forcing data with IDW

Rainfall and other meterological data can be interpolated in space using the well known IDW (inverse distance weighted) interpolation technique. Since many variations of this simple interpolation technique exist, the one used in cmf sis described. Any spatial interpolation is a function of place and known data at places. The IDW function in cmf reads as:

\[ f(x,y,z,t) = \sum^N_{i=1}{f_i(t) w_i(x,y,z)} \]

\[ w_i(x,y,z) = \frac{d_i(x,y,z)^{-p}}{\sum^N_{j=0}{d_j(x,y,z)^{-p}}} \]

\[ d_i(x,y,z) = w_z \left|z-z_i\right| + \sqrt{\left(x-x_i\right)^2 + \left(y-y_i\right)^2} \]

where:

\( x,y,z \) is the actual spatial position - \( t \) is the actual time step - \( N \) is the number of stations - \( f_i(t) \) is a data value of station \(i\) and time step \(t\) - \( w_i(x,y,z) \) is a weight function of station \(i\) for position \( x,y,z \) - \( d_i(x,y,z)\) is the distance / proximity function of station \(i\) for position \(x,y,z\) - \( p \) is the weighting power. A common value is \(p=2\) - \( x_i, y_i, z_i \) is the spatial position of station \(i\) - \( w_z \) is a weight parameter. This parameter can be used change the influence of height differences on the proximity of places. \(w_z=0\) resolves in a pure horizontal distance calculation, while very heigh values (eg. 10 000) resolve in an interpolation by height. Intermediate values (e.g. 10), depending on the height differences in the study area can lead to realistic distributions of spatiotemporal values in the landscape.

Usage:

Applying IDW interpolation for all cells over all stations is used as:

project.use_IDW_meteo(z_weight=0.0, power=2.)

project.use_IDW_rainfall(z_weight=0.0, power=2.)

with z_weight \(\hat{=}\ w_z\), power \(\hat{=}\ p\)