CAIRS - Continuous Assimilation of Integrating Rain Sensors
CAIRS is a framework to reconstruct rain fields by assimilating signals of fundamentally different rain sensors.
In particular, the integration characteristics of sensors are explicitly considered. For example, non-recording standard rain gauges integrate over time and deliver information such as the daily rainfall sums. The rain-induced attenuation of micro wave links (MWL) can be used to measure the path-integrated intensities - an example of a sensor with spatial integration.
Sensor signals with different scales (e.g. continuous, binary) can be assimilated. Furthermore, CAIRS is formulated continuously in time and space. This is helpful because it enables a natural consideration of signals with irregular time-intervals.
Note, CAIRS is still under development and the interface may change.
CAIRS can then be installed with the Julia command
After that, CAIRS behaves like a normal package. For example, it can
be updated with
First, the package CAIRS must be loaded. For convinience, it is also
recommended to load the packages
using CAIRS using Distributions using Base.Dates
Every sensor must be characterized. In the simplest case a sensor measures the rain intensity at a point. In this case (the logarithm of) the signal distribution must be defined conditioned on the intesity at this coordinate:
function log_p_gauge(S::Float64, R::Vector) # non-linear continuous rain gauge mu = 0.1+R^2.0 # Note, the signal and can be non-linearly # related to the rain intensity. sigma = 0.005 ## log of normal density, p(S|R) logpdf(Normal(mu, sigma), S) # doesn't have to be normal end sensor_gauge = Sensor(log_p_gauge)
For integrating sensors, also the integration domain must be
specified. For example, a micro wave link (MWL) with length
6 may be
function log_p_MWL(S::Float64, I::Float64) R_mean = I/6.0 sigma = 0.1 ## log of normal density, p(S|I) logpdf(Normal(R_mean, sigma), S) end sensor_MWL = Sensor(log_p_MWL, Coor(6, 0, 0)) # integrates along a path of length 6
The prior of the rain field is modeled as Gaussian process (GP). A GP is described by a mean and a covariance function.
This functions can be specified by the user. The mean function returns
the prior mean of the rain intensity at a given coordinate. It must
take a single argument of type
Coor. The covariance function must
return the covariance of the rain intensities at two given point, given
by two arguments of type
Coor. Note, it is not checked if the
provided function is a valid covariance function!
However, helpers to construct valid functions are provided. The functions
cov_exponential() create a simple constant
mean, and a separable gamma-exponential covariance function. Only the
parameters must be provided:
mean_GP = mean_constant(mean=2.0) cov_GP = cov_exponential(sigma=10.0, # standard deviation of GP l_spatial=1.5, # spatial correlation length l_temporal=Minute(1), # temporal correlation length gamma=1.0) # exponent for smoothness in [0, 2]
Other types of covariance functions will be added in future.
The next step is to import the signals. Every signal must have an
attached sensor. Signals can be constructed with the function
Signal or more conveniently with
add_signal() expected that the signals of every sensor are
stored in a separate file. The file must contain two columns:
- Column 1: date and time
- Column 2: signal values
## path to example data that come with the CAIRS package path1 = joinpath(Pkg.dir("CAIRS"), "example", "data", "Sensor1.csv") path2 = joinpath(Pkg.dir("CAIRS"), "example", "data", "Sensor2.csv") sig = Signal # create an empty array for Signals add_signal!(sig, # add signal to vector 'sig' path1, # file name sensor_gauge, # sensor Coor(5, 6), # coordinate of the sensor date_format="d.m.yyyy HH:MM:SS", delim=',') # delimitation character add_signal!(sig, path2, sensor_MWL, # MWL link Coor(4.2, 2), # coordinate of one end point of the sensor 0.9, # rotation around the point defined above in [rad] date_format="d.m.yyyy HH:MM:SS", delim=',')
Information about a signal can be printed with
Writing the sensor positions in a file is useful for plotting:
Definition of prediction points
The location for which a prediction of the rain intesity is desired must be
defined as an
Vector of coordinates. Coordinates are
Coor(x, y, time). Time can be a number or a
### create a simple grid nn = 20 loc_pred = [Coor(i, j, time) for i=linspace(0, 10, nn), j=linspace(0, 10, nn), time=DateTime(2013, 11, 22, 13, 15, 00) : Minute(1): DateTime(2013, 11, 22, 13, 20, 00) ]
This produced a regular grid, but the point could also be irregularly distributed. Also, not only predictions for coordinates but also for intesities integrated over a domain can be made. Domains are defined by the function
The assimilation of the signals and the computation of the predictions are done with
R_pred = predict(loc_pred, # vector or array with locations for predictions sig, # vector of signals mean_GP, # mean function of prior cov_GP, # covariance function of prior n_sample_calib = 20000, # number of iterations of the Gibbs sampler burn_in = 5000, # number of removed samples (and length of adaptation) n_sample_pred = 6000, # number of samples for predictions delta = Second(90)) # consider all signals within time 'delta' # from prediction points
Write a summary of the samples in a file that is used for visualization:
Visualization with R
One possibility to visualize the result is to use R. A simple
R-script to produce rain maps comes with CAIRE. It requires that R and
tripack are installed.
pathRscript = joinpath(Pkg.dir("CAIRS"), "R", "compute_rain_map.r") run(`Rscript $pathRscript rain_field.csv sensor_coor.csv out.pdf`)
Note, here it is assumed that
Rscript is in PATH.
Scheidegger, A. and Rieckermann, J. (2014) "Bayesian assimilation of rainfall sensors with fundamentally different integration characteristics" in WRaH Proceedings, Washington, DC.