47 typedef enum REG_TYPES {
54 static void IDW(
const std::vector<double>& vecData_in,
const std::vector<StationData>& vecStations_in,
56 static void LocalLapseIDW(
const std::vector<double>& vecData_in,
57 const std::vector<StationData>& vecStations_in,
58 const DEMObject& dem,
const size_t& nrOfNeighbors,
64 static void ODKriging(
const std::vector<double>& vecData,
65 const std::vector<StationData>& vecStations,
72 static bool allZeroes(
const std::vector<double>& vecData);
77 static double InvHorizontalDistance(
const double& X1,
const double& Y1,
const double& X2,
const double& Y2);
78 static double HorizontalDistance(
const double& X1,
const double& Y1,
const double& X2,
const double& Y2);
79 static double HorizontalDistance(
const DEMObject& dem,
const int& i,
const int& j,
80 const double& X2,
const double& Y2);
81 static void getNeighbors(
const double& x,
const double& y,
82 const std::vector<StationData>& vecStations,
83 std::vector< std::pair<double, size_t> >& list);
84 static void buildPositionsVectors(
const std::vector<StationData>& vecStations,
85 std::vector<double>& vecEastings, std::vector<double>& vecNorthings);
88 static double IDWCore(
const double& x,
const double& y,
89 const std::vector<double>& vecData_in,
90 const std::vector<double>& vecEastings,
const std::vector<double>& vecNorthings);
91 static double IDWCore(
const std::vector<double>& vecData_in,
const std::vector<double>& vecDistance_sq);
92 static double LLIDW_pixel(
const size_t& i,
const size_t& j,
93 const std::vector<double>& vecData_in,
94 const std::vector<StationData>& vecStations_in,
95 const DEMObject& dem,
const size_t& nrOfNeighbors);
97 static void steepestDescentDisplacement(
const DEMObject& dem,
const Grid2DObject& grid,
const size_t& ii,
const size_t& jj,
short &d_i_dest,
short &d_j_dest);
98 static double depositAroundCell(
const DEMObject& dem,
const size_t& ii,
const size_t& jj,
const double& precip,
Grid2DObject &grid);
100 static void WinstralSX(
const DEMObject& dem,
const double& dmax,
const double& in_bearing,
Grid2DObject& grid);
104 static double weightInvDist(
const double& d2);
105 static double weightInvDistSqrt(
const double& d2);
106 static double weightInvDist2(
const double& d2);
107 double weightInvDistN(
const double& d2);
static void LocalLapseIDW(const std::vector< double > &vecData_in, const std::vector< StationData > &vecStations_in, const DEMObject &dem, const size_t &nrOfNeighbors, Grid2DObject &grid)
Grid filling function: Similar to Interpol2D::LapseIDW but using a limited number of stations for eac...
Definition: libinterpol2D.cc:236
static void IDW(const std::vector< double > &vecData_in, const std::vector< StationData > &vecStations_in, const DEMObject &dem, Grid2DObject &grid)
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".
Definition: libinterpol2D.cc:306
static void ListonWind(const DEMObject &i_dem, Grid2DObject &VW, Grid2DObject &DW)
Grid filling function: This implementation fills a grid using a curvature and slope algorithm...
Definition: libinterpol2D.cc:344
no elevation dependence (ie: constant)
Definition: libinterpol2D.h:48
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...
Definition: libinterpol2D.cc:945
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.
Definition: libinterpol2D.cc:604
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. 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.
Definition: libinterpol2D.cc:560
linear elevation dependence
Definition: libinterpol2D.h:49
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 exposu...
Definition: libinterpol2D.cc:787
reg_types
Keywords for selecting the regression algorithm to use.
Definition: libinterpol2D.h:47
static void constant(const double &value, const DEMObject &dem, Grid2DObject &grid)
Grid filling function: This implementation fills the grid with a constant value.
Definition: libinterpol2D.cc:179
A class to represent DEMs and automatically compute some properties. This class stores elevation grid...
Definition: DEMObject.h:39
A class to represent 2D Grids. Typical application as DEM or Landuse Model.
Definition: Grid2DObject.h:37
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.
Definition: libinterpol2D.cc:416
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...
Definition: libinterpol2D.cc:39
A class to perform 1D regressions It works on a time serie and uses either ad-hoc methods or matrix a...
Definition: libfit1D.h:119
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 ...
Definition: libinterpol2D.cc:504
static void stdPressure(const DEMObject &dem, Grid2DObject &grid)
Grid filling function: This implementation builds a standard air pressure as a function of the elevat...
Definition: libinterpol2D.cc:160
double bearing(std::string bearing_str)
Converts a string bearing to a compass bearing.
Definition: IOUtils.cc:70
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 l...
Definition: libinterpol2D.cc:656
A class to perform 2D spatial interpolations. Each parameter to be interpolated declares which interp...
Definition: libinterpol2D.h:44