MeteoIODoc 20240328.18c31bd1
mio::Interpol2D Class Reference

Detailed Description

A class to perform 2D spatial interpolations. Each parameter to be interpolated declares which interpolation method to use. Then the class computes the interpolation for each 2D grid point, combining the inputs provided by the available data sources.

Author
Mathias Bavay

#include <libinterpol2D.h>

Public Types

enum  reg_types { R_CST , R_LIN }
 Keywords for selecting the regression algorithm to use. More...
 

Static Public Member Functions

static void stdPressure (const DEMObject &dem, Grid2DObject &grid)
 Grid filling function: This implementation builds a standard air pressure as a function of the elevation. More...
 
static void constant (const double &value, const DEMObject &dem, Grid2DObject &grid)
 Grid filling function: This implementation fills the grid with a constant value. More...
 
static void IDW (const std::vector< double > &vecData_in, const std::vector< StationData > &vecStations_in, const DEMObject &dem, Grid2DObject &grid, const double &scale, const double &alpha=1.)
 Grid filling function: This implementation fills a grid using Inverse Distance Weighting. for example, the air temperatures measured at several stations would be given as values, the stations positions as positions and projected to a grid. No elevation detrending is performed, the DEM is only used for checking if a grid point is "nodata". More...
 
static void LocalLapseIDW (const std::vector< double > &vecData_in, const std::vector< StationData > &vecStations_in, const DEMObject &dem, const size_t &nrOfNeighbors, const double &MaxDistance, Grid2DObject &grid, const double &scale, const double &alpha=1.)
 Grid filling function: Similar to Interpol2D::LapseIDW but using a limited number of stations for each cell. More...
 
static void ListonWind (const DEMObject &i_dem, Grid2DObject &VW, Grid2DObject &DW, const double &eta)
 Grid filling function: This implementation fills a grid using a curvature and slope algorithm, as described in G. E. Liston and K. Elder, "A meteorological distribution system for high-resolution terrestrial modeling (MicroMet)", Journal of Hydrometeorology, 7.2, 2006. More...
 
static void CurvatureCorrection (DEMObject &dem, const Grid2DObject &ta, Grid2DObject &grid)
 Distribute precipitation in a way that reflects snow redistribution on the ground, according to (Huss, 2008) This method modifies the solid precipitation distribution according to the local slope and curvature. See "Quantitative evaluation of different hydrological modelling approaches in a partly glacierized Swiss watershed", Magnusson et All., Hydrological Processes, 2010, under review. and "Modelling runoff from highly glacierized alpine catchments in a changing climate", Huss et All., Hydrological Processes, 22, 3888-3902, 2008. More...
 
static void SteepSlopeRedistribution (const DEMObject &dem, const Grid2DObject &ta, Grid2DObject &grid)
 redistribute precip from steeper slopes to gentler slopes by following the steepest path from top to bottom and gradually depositing precip during descent More...
 
static void PrecipSnow (const DEMObject &dem, const Grid2DObject &ta, Grid2DObject &grid)
 Distribute precipitation in a way that reflects snow redistribution on the ground, according to (Huss, 2008) This method modifies the solid precipitation distribution according to the local slope and curvature: all pixels whose slope is greater than 60° will not receive any snow at all. All pixels whose slope is less than 40° will receive full snow and any pixel between 40° and 60° sees a linear correction between 100% and 0% snow. After this step, a curvature correction is applied: pixels having the minimu curvature see 50% snow more, pixels having the maximum curvature see 50% snow less and pixels ate the middle of the curvature range are unaffected. More...
 
static void ODKriging (const std::vector< double > &vecData, const std::vector< StationData > &vecStations, const DEMObject &dem, const Fit1D &variogram, Grid2DObject &grid)
 Ordinary Kriging matrix formulation This implements the matrix formulation of Ordinary Kriging, as shown (for example) in "Statistics for spatial data", Noel A. C. Cressie, John Wiley & Sons, revised edition, 1993, pp122. More...
 
static void RyanWind (const DEMObject &dem, Grid2DObject &VW, Grid2DObject &DW)
 compute the change of wind direction by the local terrain This is according to Ryan, "a mathematical model for diagnosis and prediction of surface winds in mountainous terrain", 1977, journal of applied meteorology, 16, 6. More...
 
static void Winstral (const DEMObject &dem, const Grid2DObject &TA, const double &dmax, const double &in_bearing, Grid2DObject &grid)
 Alter a precipitation field with the Winstral Sx exposure coefficient This implements the wind exposure coefficient (Sx) for one bearing as in "Simulating wind fields and snow redistribution using terrain‐based parameters to model snow accumulation and melt over a semi‐arid mountain catchment.", Winstral, Adam, and Danny Marks, Hydrological Processes 16.18 (2002), pp3585-3603. More...
 
static void Winstral (const DEMObject &dem, const Grid2DObject &TA, const Grid2DObject &DW, const Grid2DObject &VW, const double &dmax, Grid2DObject &grid)
 
static bool allZeroes (const std::vector< double > &vecData)
 check if the points measurements are all at zero This check can be performed to trigger optimizations: it is quicker to fill the grid directly with zeroes instead of running a complicated algorithm. More...
 
static double getTanMaxSlope (const Grid2DObject &dem, const double &dmin, const double &dmax, const double &bearing, const size_t &ii, const size_t &jj)
 compute the max slope angle looking toward the horizon in a given direction. The search distance is limited between dmin and dmax from the starting point (i,j) and the elvation difference must be at least 2 meters (otherwise it is considered flat). If a slope start to be positive/negative before turning negative/positive, the first one would be returned (so a pixel just behind a ridge would still see an uphill slope to the ridge even if the slope behind the ridge would be greater). More...
 

Member Enumeration Documentation

◆ reg_types

Keywords for selecting the regression algorithm to use.

Enumerator
R_CST 

no elevation dependence (ie: constant)

R_LIN 

linear elevation dependence

Member Function Documentation

◆ allZeroes()

bool mio::Interpol2D::allZeroes ( const std::vector< double > &  vecData)
static

check if the points measurements are all at zero This check can be performed to trigger optimizations: it is quicker to fill the grid directly with zeroes instead of running a complicated algorithm.

Returns
true if all data is set to zero

◆ constant()

void mio::Interpol2D::constant ( const double &  value,
const DEMObject dem,
Grid2DObject grid 
)
static

Grid filling function: This implementation fills the grid with a constant value.

Parameters
valuevalue to put in the grid
demarray of elevations (dem). This is needed in order to know if a point is "nodata"
grid2D array to fill

◆ CurvatureCorrection()

void mio::Interpol2D::CurvatureCorrection ( DEMObject dem,
const Grid2DObject ta,
Grid2DObject grid 
)
static

Distribute precipitation in a way that reflects snow redistribution on the ground, according to (Huss, 2008) This method modifies the solid precipitation distribution according to the local slope and curvature. See "Quantitative evaluation of different hydrological modelling approaches in a partly glacierized Swiss watershed", Magnusson et All., Hydrological Processes, 2010, under review. and "Modelling runoff from highly glacierized alpine catchments in a changing climate", Huss et All., Hydrological Processes, 22, 3888-3902, 2008.

Parameters
demarray of elevations (dem). The slope must have been updated as it is required for the DEM analysis.
taarray of air temperatures used to determine if precipitation is rain or snow
grid2D array of precipitation to fill
Author
Florian Kobierska, Jan Magnusson, Rob Spence and Mathias Bavay

◆ getTanMaxSlope()

double mio::Interpol2D::getTanMaxSlope ( const Grid2DObject dem,
const double &  dmin,
const double &  dmax,
const double &  bearing,
const size_t &  i,
const size_t &  j 
)
static

compute the max slope angle looking toward the horizon in a given direction. The search distance is limited between dmin and dmax from the starting point (i,j) and the elvation difference must be at least 2 meters (otherwise it is considered flat). If a slope start to be positive/negative before turning negative/positive, the first one would be returned (so a pixel just behind a ridge would still see an uphill slope to the ridge even if the slope behind the ridge would be greater).

This is exactly identical with the Winstral Sx factor for a single direction. Or the upwind slope for Ryan.

Parameters
[in]demDEM to work with
[in]dminminimum search distance (ie all points at less than dmin are skipped)
[in]dmaxmaximum search distance
[in]bearingdirection of the search
[in]ix index of the cell to start the search from
[in]jy index of the cell to start the search from
Returns
tan of the maximum slope angle from the (i,j) cell in the given direction

◆ IDW()

void mio::Interpol2D::IDW ( const std::vector< double > &  vecData_in,
const std::vector< StationData > &  vecStations_in,
const DEMObject dem,
Grid2DObject grid,
const double &  scale,
const double &  alpha = 1. 
)
static

Grid filling function: This implementation fills a grid using Inverse Distance Weighting. for example, the air temperatures measured at several stations would be given as values, the stations positions as positions and projected to a grid. No elevation detrending is performed, the DEM is only used for checking if a grid point is "nodata".

Parameters
vecData_ininput values to use for the IDW
vecStations_inposition of the "values" (altitude and coordinates)
demarray of elevations (dem). This is needed in order to know if a point is "nodata"
grid2D array to fill
scaleThe scale factor is used to smooth the grid. It is added to the distance before applying the weights in order to come into the tail of "1/d".
alphaThe weights are computed as 1/dist^alpha, so give alpha=1 for standards 1/dist weights.

◆ ListonWind()

void mio::Interpol2D::ListonWind ( const DEMObject i_dem,
Grid2DObject VW,
Grid2DObject DW,
const double &  eta 
)
static

Grid filling function: This implementation fills a grid using a curvature and slope algorithm, as described in G. E. Liston and K. Elder, "A meteorological distribution system for high-resolution terrestrial modeling (MicroMet)", Journal of Hydrometeorology, 7.2, 2006.

Parameters
i_demarray of elevations (dem). The slope must have been updated as it is required for the DEM analysis.
VW2D array of Wind Velocity to fill
DW2D array of Wind Direction to fill
eta(curvature length scale)

◆ LocalLapseIDW()

void mio::Interpol2D::LocalLapseIDW ( const std::vector< double > &  vecData_in,
const std::vector< StationData > &  vecStations_in,
const DEMObject dem,
const size_t &  nrOfNeighbors,
const double &  MaxDistance,
Grid2DObject grid,
const double &  scale,
const double &  alpha = 1. 
)
static

Grid filling function: Similar to Interpol2D::LapseIDW but using a limited number of stations for each cell.

Parameters
vecData_ininput values to use for the IDW
vecStations_inposition of the "values" (altitude and coordinates)
demarray of elevations (dem)
nrOfNeighborsnumber of neighboring stations to use for each pixel
MaxDistancemaximum allowed distance of neighboring stations to the target grid point
grid2D array to fill
scaleThe scale factor is used to smooth the grid. It is added to the distance before applying the weights in order to come into the tail of "1/d".
alphaThe weights are computed as 1/dist^alpha, so give alpha=1 for standards 1/dist weights.

◆ ODKriging()

void mio::Interpol2D::ODKriging ( const std::vector< double > &  vecData,
const std::vector< StationData > &  vecStations,
const DEMObject dem,
const Fit1D variogram,
Grid2DObject grid 
)
static

Ordinary Kriging matrix formulation This implements the matrix formulation of Ordinary Kriging, as shown (for example) in "Statistics for spatial data", Noel A. C. Cressie, John Wiley & Sons, revised edition, 1993, pp122.

First, Ordinary kriging assumes stationarity of the mean of all random variables. We start by solving the following system:

\begin{eqnarray*} \mathbf{\lambda} & = & \mathbf{\Gamma_0^{-1}} \cdot \mathbf{\gamma^*} \\ \left[ \begin{array}{c} \lambda_1 \\ \vdots \\ \lambda_i \\ \mu \end{array} \right] & = & { \left[ \begin{array}{cccc} \Gamma_{1,1} & \cdots & \Gamma_{1,i} & 1 \\ \vdots & \ddots & \vdots & \vdots \\ \Gamma_{i,1} & \cdots & \Gamma_{i,i} & 1 \\ 1 & \cdots & 1 & 0 \end{array} \right] }^{-1} \cdot \left[ \begin{array}{c} \gamma^*_1 \\ \vdots \\ \gamma^*_i \\ 1 \end{array} \right] \end{eqnarray*}

where the \(\lambda_i\) are the interpolation weights (at each station i), \(\mu\) is the Lagrange multiplier (used to minimize the error), \(\Gamma_{i,j}\) is the covariance between the stations i and j and \(\gamma^*_i\) the covariances between the station i and the local position where the interpolation has to be computed. This covariance is computed based on distance, using the variogram that gives covariance = f(distance). The variogram is established by fitting a statistical model to all the (distance, covariance) points originating from the station measurements. The statistical model of the variogram enables computing the covariance for any distance, therefore it is possible to compute the \(\gamma^*_i\).

Once the \(\lambda_i\) have been computed, the locally interpolated value is computed as

\[ \mathbf{X^*} = \sum \mathbf{\lambda_i} * \mathbf{X_i} \]

where \(X_i\) is the value measured at station i.

Parameters
vecDatavector containing the values as measured at the stations
vecStationsvector of stations
demdigital elevation model
variogramvariogram regression model
grid2D array of precipitation to fill
Author
Mathias Bavay

◆ PrecipSnow()

void mio::Interpol2D::PrecipSnow ( const DEMObject dem,
const Grid2DObject ta,
Grid2DObject grid 
)
static

Distribute precipitation in a way that reflects snow redistribution on the ground, according to (Huss, 2008) This method modifies the solid precipitation distribution according to the local slope and curvature: all pixels whose slope is greater than 60° will not receive any snow at all. All pixels whose slope is less than 40° will receive full snow and any pixel between 40° and 60° sees a linear correction between 100% and 0% snow. After this step, a curvature correction is applied: pixels having the minimu curvature see 50% snow more, pixels having the maximum curvature see 50% snow less and pixels ate the middle of the curvature range are unaffected.

For more, see "Quantitative evaluation of different hydrological modelling approaches in a partly glacierized Swiss watershed", Magnusson et All., Hydrological Processes, 2010, under review. and "Modelling runoff from highly glacierized alpine catchments in a changing climate", Huss et All., Hydrological Processes, 22, 3888-3902, 2008.

Parameters
demarray of elevations (dem). The slope must have been updated as it is required for the DEM analysis.
taarray of air temperatures used to determine if precipitation is rain or snow
grid2D array of precipitation to fill
Author
Florian Kobierska, Jan Magnusson and Mathias Bavay

◆ RyanWind()

void mio::Interpol2D::RyanWind ( const DEMObject dem,
Grid2DObject VW,
Grid2DObject DW 
)
static

compute the change of wind direction by the local terrain This is according to Ryan, "a mathematical model for diagnosis and prediction of surface winds in mountainous terrain", 1977, journal of applied meteorology, 16, 6.

Parameters
demarray of elevations (dem). The slope and azimuth must have been updated as they are required for the DEM analysis.
VW2D array of wind speed to fill
DW2D array of wind direction to fill
Author
Mathias Bavay

◆ stdPressure()

void mio::Interpol2D::stdPressure ( const DEMObject dem,
Grid2DObject grid 
)
static

Grid filling function: This implementation builds a standard air pressure as a function of the elevation.

Parameters
demarray of elevations (dem)
grid2D array to fill

◆ SteepSlopeRedistribution()

void mio::Interpol2D::SteepSlopeRedistribution ( const DEMObject dem,
const Grid2DObject ta,
Grid2DObject grid 
)
static

redistribute precip from steeper slopes to gentler slopes by following the steepest path from top to bottom and gradually depositing precip during descent

Parameters
demarray of elevations (dem). The slope must have been updated as it is required for the DEM analysis.
taarray of air temperatures used to determine if precipitation is rain or snow
grid2D array of precipitation to fill
Author
Rob Spence and Mathias Bavay

◆ Winstral() [1/2]

void mio::Interpol2D::Winstral ( const DEMObject dem,
const Grid2DObject TA,
const double &  dmax,
const double &  in_bearing,
Grid2DObject grid 
)
static

Alter a precipitation field with the Winstral Sx exposure coefficient This implements the wind exposure coefficient (Sx) for one bearing as in "Simulating wind fields and snow redistribution using terrain‐based parameters to model snow accumulation and melt over a semi‐arid mountain catchment.", Winstral, Adam, and Danny Marks, Hydrological Processes 16.18 (2002), pp3585-3603.

A linear correlation between erosion coefficients and eroded mass is assumed, that is that the points with maximum erosion get all their precipitation removed. The eroded mass is then distributed on the cells with positive Sx (with a linear correlation between positive Sx and deposited mass) and enforcing mass conservation within the domain.

Remarks
Only cells with an air temperature below freezing participate in the redistribution
Parameters
demdigital elevation model
TAair temperature grid (in order to discriminate between solid and liquid precipitation)
dmaxsearch radius
in_bearingwind direction to consider
grid2D array of precipitation to fill
Author
Mathias Bavay

◆ Winstral() [2/2]

void mio::Interpol2D::Winstral ( const DEMObject dem,
const Grid2DObject TA,
const Grid2DObject DW,
const Grid2DObject VW,
const double &  dmax,
Grid2DObject grid 
)
static

The documentation for this class was generated from the following files: