Gage adjustment#
Concept#
The objective of this module is the adjustment of radar-based rainfall estimates by rain gage observations. However, this module could also be applied to adjust satellite rainfall by rain gage observations, remotely sensed soil moisture patterns by ground truthing moisture sensors, or any dense spatial point pattern which could be adjusted by sparse point measurements (ground truth).
Basically, we only need two data sources:
point observations (e.g. rain gage observations)
set of (potentially irregular) unadjusted point values (e.g. remotely sensed rainfall)
[Goudenhoofdt et al., 2009] provide an excellent overview of adjustment procedures. The general idea is that we quantify the error of the remotely sensed rainfall at the rain gage locations, assuming the rain gage observation to be accurate.
The error can be assumed to be purely additive
(AdjustAdd
), purely multiplicative
(AdjustMultiply
, AdjustMFB
)
or a mixture of both (AdjustMixed
).
If the error is assumed to be heterogeneous in space
(AdjustAdd
, AdjustMultiply
,
AdjustMixed
), the error at the rain gage locations is
interpolated to the radar bin locations and then used to adjust (correct)
the raw radar rainfall estimates. In case of the AdjustMFB approach, though,
the multiplicative error is assumed to be homogeneous in space.
Quick start#
The basic procedure consists of creating an adjustment object from the class
you want to use for adjustment. After that, you can call the object with the
actual data that is to be adjusted. The following example is using the
additive error model with default settings. obs_coords
and
raw_coords
represent arrays with coordinate pairs for the gage
observations and the radar bins, respectively. obs
and raw
are
arrays containing the actual data:
adjuster = AdjustAdd(obs_coords, raw_coords)
adjusted = adjuster(obs, raw)
Both obs
and raw
need to be flat (1-dimensional) arrays of shape (n, )
that have the same length as the obs_coords
and raw_coords
arrays,
respectively.
The user can specify the approach that should be used to interpolate the error
in space, as well as the keyword arguments which control the behaviour of the
interpolation approach. For this purpose, all interpolation classes from the
wradlib.ipol
module are available and can be passed by using the
ipclass
argument. The default interpolation class is
Inverse Distance Weighting (Idw
). If you want to use
e.g. linear barycentric interpolation:
import wradlib.ipol as ipol
adjuster = AdjustAdd(obs_coords, raw_coords, ipclass=ipol.Linear)
adjusted = adjuster(obs, raw)
Warning
Be aware that there are a lot of control parameters that can dramatically influence the behaviour of the adjustment (which gauges are considered, how is an error interpolation carried out, …). Read the docs carefully and try to experiment with the effects of the different control parameters. There might be situations in which the algorithms decide - based on the control parameter - not to do an adjustment and just return the unadjusted values.
Cross validation#
Another helpful feature is an easy-to-use method for leave-one-out
cross-validation [B4]. Cross validation is a standard
procedure for verifying rain gage adjustment or interpolation procedures. You
can start the cross validation in the same way as you start the actual
adjustment, however, you call the xvalidate
method instead. The result of the cross validation are pairs of observation
and the corresponding adjustment result at the observation location. Using the
wradlib.verify
module, you can compute error metrics for the cross
validation results:
adjuster = AdjustAdd(obs_coords, raw_coords)
observed, estimated = adjuster.xvalidate(obs, raw)
from wradlib.verify import ErrorMetrics
metrics = ErrorMetrics(observed, estimated)
metrics.report()
The basic adjustment class that inherits to all other classes. |
|
Multiplicative gage adjustment using one correction factor for the entire domain. |
|
Gage adjustment using a multiplicative error model |
|
Gage adjustment using an additive error model. |
|
Gage adjustment using a mixed error model (additive and multiplicative). |
|
Get the raw values in the neighbourhood of the observation points |
|
Same behaviour as the other adjustment classes, but returns an interpolation of rain gage observations |
|
Same behaviour as the other adjustment classes, but simply returns the unadjusted data. |
- class wradlib.adjust.AdjustBase(obs_coords, raw_coords, *, nnear_raws=9, stat='median', mingages=5, minval=0.0, mfb_args=None, ipclass=<class 'wradlib.ipol.Idw'>, **ipargs)[source]#
Bases:
IpolBase
The basic adjustment class that inherits to all other classes.
All methods except the
__call__
method are inherited to the following adjustment classes.- Parameters:
obs_coords (
numpy.ndarray
) – array of floats of shape (number of points, 2) x and y coordinate pairs of observation locations (e.g. rain gauges).raw_coords (
numpy.ndarray
) – array of floats of shape (number of points, 2) x and y coordinate pairs of raw (unadjusted) radar fieldnnear_raws (
int
) – Defaults to 9. This parameter controls the number of radar bins or grid cells (in the neighbourhood of a rain gauge) which is used to compute the value of the radar observation AT a rain gauge.stat (
str
) – Defaults to ‘median’. Must be either ‘mean’, ‘median’, or ‘best’. This parameter controls the statistic that is used to compute the value of the radar observation AT a rain gauge based on the neighbourhood specified by parameternnear_raws
.mingages (
int
) – Defaults to 5. Minimum number of valid gages required for an adjustment. If less valid gauges are available, the adjustment procedure will return unadjusted raw values. If you do not want to use this feature, you need to setmingages=0
.minval (
float
) – If the gage or radar observation is below this threshold, the location will not be used for adjustment. For additive adjustment, this value should be set to zero (default value). For multiplicative adjustment, values larger than zero might be chosen in order to minimize artifacts.mfb_args (
dict
) – Only used for AdjustMFB - This set of parameters controls how the mean field bias is computed. Items of the dictionary are:method: string defaults to ‘linregr’ which fits a regression line through observed and estimated values and then gets the bias from the inverse of the slope. Other values: ‘mean’ or ‘median’ compute the mean or the median of the ratios between gauge and radar observations.
minslope, minr, maxp: When using method=’linregr’, these parameters control whether a linear regression turned out to be robust (minimum allowable slope, minimum allowable correlation, maximim allowable p-value). If the regression result is not considered robust, no adjustment will take place.
Ipclass (
wradlib.ipol.IpolBase
) – an interpolation class fromwradlib.ipol
Not used for AdjustMFB - default value isIdw
(Inverse Distance Weighting).ipargs (
dict
) – keyword arguments to create an instance of ipclass Not used for AdjustMFB - forIdw
, these keyword arguments would e.g. bennear
orp
.
Examples
See Adjusting radar-base rainfall estimates by rain gauge observations.
- xvalidate(obs, raw)[source]#
Leave-One-Out Cross Validation, applicable to all gage adjustment classes.
This method will be inherited to other Adjust classes. It should thus be applicable to all adjustment procedures without any modification. This way, the actual adjustment procedure has only to be defined once in the
__call__
method.The output of this method can be evaluated by using the verify.ErrorMetrics class.
- Parameters:
obs (
numpy.ndarray
) – array of floatsraw (
numpy.ndarray
) – array of floats
- Returns:
obs (
numpy.ndarray
) – array of floats valid observations at those locations which have a valid radar observationestatobs (
numpy.ndarray
) – array of floats estimated values at the valid observation locations
- class wradlib.adjust.AdjustMFB(obs_coords, raw_coords, *, nnear_raws=9, stat='median', mingages=5, minval=0.0, mfb_args=None, ipclass=<class 'wradlib.ipol.Idw'>, **ipargs)[source]#
Bases:
AdjustBase
Multiplicative gage adjustment using one correction factor for the entire domain.
This method is also known as the Mean Field Bias correction.
Note
Inherits from
wradlib.adjust.AdjustBase
For a complete overview of parameters for the initialisation of adjustment objects, as well as an extensive example, please see
wradlib.adjust.AdjustBase
.- Returns:
output (
numpy.ndarray
) – array of adjusted radar values
- class wradlib.adjust.AdjustMultiply(obs_coords, raw_coords, *, nnear_raws=9, stat='median', mingages=5, minval=0.0, mfb_args=None, ipclass=<class 'wradlib.ipol.Idw'>, **ipargs)[source]#
Bases:
AdjustBase
Gage adjustment using a multiplicative error model
First, an instance of AdjustMultiply has to be created. Calling this instance then does the actual adjustment. The motivation behind this performance. In case the observation points are always the same for different time steps, the computation of neighbours and inverse distance weights only needs to be performed once during initialisation.
AdjustMultiply automatically takes care of invalid gage or radar observations (e.g. NaN, Inf or other typical missing data flags such as -9999). However, in case e.g. the observation data contain missing values, the computation of the inverse distance weights needs to be repeated in
__call__
which is at the expense of performance.Note
Inherits from
wradlib.adjust.AdjustBase
For a complete overview of parameters for the initialisation of adjustment objects, as well as an extensive example, please see
wradlib.adjust.AdjustBase
.- Returns:
output (
numpy.ndarray
) – array of adjusted radar values
- class wradlib.adjust.AdjustAdd(obs_coords, raw_coords, *, nnear_raws=9, stat='median', mingages=5, minval=0.0, mfb_args=None, ipclass=<class 'wradlib.ipol.Idw'>, **ipargs)[source]#
Bases:
AdjustBase
Gage adjustment using an additive error model.
First, an instance of AdjustAdd has to be created. Calling this instance then does the actual adjustment. The motivation behind this performance. In case the observation points are always the same for different time steps, the computation of neighbours and inverse distance weights only needs to be performed once.
AdjustAdd automatically takes care of invalid gage or radar observations (e.g. NaN, Inf or other typical missing data flags such as -9999). However, in case e.g. the observation data contains missing values, the computation of the inverse distance weights needs to be repeated in
__call__
which is at the expense of performance.Note
Inherits from
wradlib.adjust.AdjustBase
For a complete overview of parameters for the initialisation of adjustment objects, as well as an extensive example, please see
wradlib.adjust.AdjustBase
.- Returns:
output (
numpy.ndarray
) – array of adjusted radar values
- class wradlib.adjust.AdjustMixed(obs_coords, raw_coords, *, nnear_raws=9, stat='median', mingages=5, minval=0.0, mfb_args=None, ipclass=<class 'wradlib.ipol.Idw'>, **ipargs)[source]#
Bases:
AdjustBase
Gage adjustment using a mixed error model (additive and multiplicative).
The mixed error model assumes that you have both a multiplicative and an additive error term. The intention is to overcome the drawbacks of the purely additive and multiplicative approaches (see
AdjustAdd
andAdjustMultiply
). The formal representation of the error model according to [Pfaff, 2010] is:\[R_{gage} = R_{radar} \cdot (1 + \delta) +0 \epsilon\]\(\delta\) and \(\epsilon\) have to be assumed to be independent and normally distributed. The present implementation is based on teh Least Squares estimation of \(\delta\) and \(\epsilon\) for each rain gage location. \(\delta\) and \(\epsilon\) are then interpolated and used to correct the radar rainfall field.
The least squares implementation uses the equation for the error model plus the condition to minimize (\(\delta^2 + \epsilon^2\)) for each gage location. The idea behind this is that \(\epsilon\) dominates the adjustment for small deviations between radar and gage while \(\delta\) dominates in case of large deviations.
Usage: First, an instance of AdjustMixed has to be created. Calling this instance then does the actual adjustment. The motivation behind this is performance. In case the observation points are always the same for different time steps, the computation of neighbours and inverse distance weights only needs to be performed once during initialisation.
AdjustMixed automatically takes care of invalid gage or radar observations (e.g. NaN, Inf or other typical missing data flags such as -9999). However, in case e.g. the observation data contain missing values, the computation of the inverse distance weights needs to be repeated in
__call__
which is at the expense of performance.Note
Inherits from
wradlib.adjust.AdjustBase
For a complete overview of parameters for the initialisation of adjustment objects, as well as an extensive example, please see
wradlib.adjust.AdjustBase
.- Returns:
output (
numpy.ndarray
) – array of adjusted radar values
- class wradlib.adjust.RawAtObs(obs_coords, raw_coords, *, nnear=9, stat='median')[source]#
Bases:
object
Get the raw values in the neighbourhood of the observation points
- Parameters:
obs_coords (
numpy.ndarray
) – array of float coordinate pairs of observations pointsraw_coords (
numpy.ndarray
) – array of float coordinate pairs of raw (unadjusted) fieldnnear (
int
) – number of neighbours which should be considered in the vicinity of each point in obsstat (
str
) – function name
- class wradlib.adjust.GageOnly(obs_coords, raw_coords, *, nnear_raws=9, stat='median', mingages=5, minval=0.0, mfb_args=None, ipclass=<class 'wradlib.ipol.Idw'>, **ipargs)[source]#
Bases:
AdjustBase
Same behaviour as the other adjustment classes, but returns an interpolation of rain gage observations
First, an instance of GageOnly has to be created. Calling this instance then does the actual adjustment. The motivation behind this performance. In case the observation points are always the same for different time steps, the computation of neighbours and inverse distance weights only needs to be performed once during initialisation.
GageOnly automatically takes care of invalid gage or radar observations (e.g. NaN, Inf or other typical missing data flags such as -9999). However, in case e.g. the observation data contain missing values, the computation of the inverse distance weights needs to be repeated in
__call__
which is at the expense of performance.Note
Inherits from
wradlib.adjust.AdjustBase
For a complete overview of parameters for the initialisation of adjustment objects, as well as an extensive example, please see
wradlib.adjust.AdjustBase
.- Returns:
output (
numpy.ndarray
) – array of adjusted radar values
- class wradlib.adjust.AdjustNone(obs_coords, raw_coords, *, nnear_raws=9, stat='median', mingages=5, minval=0.0, mfb_args=None, ipclass=<class 'wradlib.ipol.Idw'>, **ipargs)[source]#
Bases:
AdjustBase
Same behaviour as the other adjustment classes, but simply returns the unadjusted data.
This class can be used for benchmark verification experiments as a control for unadjusted data.
Note
Inherits from
wradlib.adjust.AdjustBase
For a complete overview of parameters for the initialisation of adjustment objects, as well as an extensive example, please see
wradlib.adjust.AdjustBase
.- Returns:
output (
numpy.ndarray
) – array of unadjusted radar values