MI - Fimex
Functions
interpolation.h File Reference
#include "fimex/deprecated.h"
#include "fimex/mifi_constants.h"

Go to the source code of this file.

Functions

int mifi_string_to_interpolation_method (const char *stringMethod)
 
int mifi_interpolate_f (int method, const char *proj_input, const float *infield, const double *in_x_axis, const double *in_y_axis, const int in_x_axis_type, const int in_y_axis_type, const int ix, const int iy, const int iz, const char *proj_output, float *outfield, const double *out_x_axis, const double *out_y_axis, const int out_x_axis_type, const int out_y_axis_type, const int ox, const int oy)
 
int mifi_interpolate_d (int method, char *proj_input, double *infield, double *in_x_axis, double *in_y_axis, int in_x_axis_type, int in_y_axis_type, int ix, int iy, int iz, char *proj_output, double *outfield, double *out_x_axis, double *out_y_axis, int out_x_axis_type, int out_y_axis_type, int ox, int oy)
 not implemented yet More...
 
int mifi_vector_reproject_values_f (int method, const char *proj_input, const char *proj_output, float *u_out, float *v_out, const double *out_x_axis, const double *out_y_axis, int out_x_axis_type, int out_y_axis_type, int ox, int oy, int oz)
 interpolate the vector values More...
 
int mifi_vector_reproject_values_by_matrix_f (int method, const double *matrix, float *u_out, float *v_out, int ox, int oy, int oz)
 
int mifi_vector_reproject_direction_by_matrix_f (int method, const double *matrix, float *angle_out, int ox, int oy, int oz)
 
int mifi_get_vector_reproject_matrix (const char *proj_input, const char *proj_output, const double *out_x_axis, const double *out_y_axis, int out_x_axis_type, int out_y_axis_type, int ox, int oy, double *matrix)
 
int mifi_get_vector_reproject_matrix_field (const char *proj_input, const char *proj_output, const double *in_x_field, const double *in_y_field, int ox, int oy, double *matrix)
 
int mifi_get_vector_reproject_matrix_points (const char *proj_input, const char *proj_output, int inputIsMetric, const double *out_x_points, const double *out_y_points, int on, double *matrix)
 
int mifi_get_values_f (const float *infield, float *outfield, const double x, const double y, const int ix, const int iy, const int iz)
 
int mifi_get_values_weak_extrapol_f (const float *infield, float *outfield, const double x, const double y, const int ix, const int iy, const int iz)
 
int mifi_get_values_no_extrapol_f (const float *infield, float *outfield, const double x, const double y, const int ix, const int iy, const int iz)
 
int mifi_get_values_bilinear_f (const float *infield, float *outvalues, const double x, const double y, const int ix, const int iy, const int iz)
 
int mifi_get_values_bicubic_f (const float *infield, float *outvalues, const double x, const double y, const int ix, const int iy, const int iz)
 not implemented yet More...
 
int mifi_get_values_nearest_f (const float *infieldA, const float *infieldB, float *outfield, const size_t n, const double a, const double b, const double x)
 
int mifi_get_values_linear_f (const float *infieldA, const float *infieldB, float *outfield, const size_t n, const double a, const double b, const double x)
 
int mifi_get_values_linear_weak_extrapol_f (const float *infieldA, const float *infieldB, float *outfield, const size_t n, const double a, const double b, const double x)
 
int mifi_get_values_linear_no_extrapol_f (const float *infieldA, const float *infieldB, float *outfield, const size_t n, const double a, const double b, const double x)
 
int mifi_get_values_linear_const_extrapol_f (const float *infieldA, const float *infieldB, float *outfield, const size_t n, const double a, const double b, const double x)
 
int mifi_get_values_linear_d (const double *infieldA, const double *infieldB, double *outfield, const size_t n, const double a, const double b, const double x)
 
int mifi_get_values_log_f (const float *infieldA, const float *infieldB, float *outfield, const size_t n, const double a, const double b, const double x)
 
int mifi_get_values_log_log_f (const float *infieldA, const float *infieldB, float *outfield, const size_t n, const double a, const double b, const double x)
 
int mifi_points2position (double *points, const int n, const double *axis, const int num, const int axis_type)
 find position in array of position in projection More...
 
static int mifi_3d_array_position (int x, int y, int z, int ix, int iy, int iz)
 
int mifi_project_values (const char *proj_input, const char *proj_output, double *in_out_x_vals, double *in_out_y_vals, const int num)
 project values so that the projetion (x,y) => (x_proj), (y_proj) can be expressed as x_proj(x,y), y_proj(x,y) More...
 
int mifi_project_axes (const char *proj_input, const char *proj_output, const double *in_x_axis, const double *in_y_axis, const int ix, const int iy, double *out_xproj_axis, double *out_yproj_axis)
 project axes so that the projetion (x,y) => (x_proj), (y_proj) can be expressed as x_proj(x,y), y_proj(x,y) More...
 
int mifi_fill2d_f (size_t nx, size_t ny, float *field, float relaxCrit, float corrEff, size_t maxLoop, size_t *nChanged)
 Method to fill undefined values in a 2d field. More...
 
int mifi_creepfill2d_f (size_t nx, size_t ny, float *field, unsigned short repeat, char setWeight, size_t *nChanged)
 Method to fill undefined values in a 2d field in stable time. More...
 
int mifi_creepfillval2d_f (size_t nx, size_t ny, float *field, float defaultVal, unsigned short repeat, char setWeight, size_t *nChanged)
 Method to fill undefined values in a 2d field in stable time. More...
 
int mifi_griddistance (size_t nx, size_t ny, const double *lonVals, const double *latVals, float *gridDistX, float *gridDistY)
 
size_t mifi_compute_vertical_velocity (size_t nx, size_t ny, size_t nz, double dx, double dy, const float *gridDistX, const float *gridDistY, const double *ap, const double *b, const float *zs, const float *ps, const float *u, const float *v, const float *t, float *w)
 
size_t mifi_bad2nanf (float *posPtr, float *endPtr, float badVal)
 
size_t mifi_nanf2bad (float *posPtr, float *endPtr, float badVal)
 
 MIFI_DEPRECATED (int mifi_isnanf(float val))
 
 MIFI_DEPRECATED (int mifi_isnand(double val))
 

Function Documentation

static int mifi_3d_array_position ( int  x,
int  y,
int  z,
int  ix,
int  iy,
int  iz 
)
inlinestatic
size_t mifi_bad2nanf ( float *  posPtr,
float *  endPtr,
float  badVal 
)

Convert bad-values to nan. The mifi_ functions don't handle bad values generally, but forward this work to the floating-point IEEE NaN's. This function converts a general bad value to a nan in a float array.

Parameters
posPtrstart pointer of the float array
endPtrend-pointer of the float array (excluded from conversion)
badValbad value to be converted to nan
Returns
0 (in fimex < 0.61, this was possibly the number of conversions, but it was never strict)

Referenced by mifi_3d_array_position().

size_t mifi_compute_vertical_velocity ( size_t  nx,
size_t  ny,
size_t  nz,
double  dx,
double  dy,
const float *  gridDistX,
const float *  gridDistY,
const double *  ap,
const double *  b,
const float *  zs,
const float *  ps,
const float *  u,
const float *  v,
const float *  t,
float *  w 
)

Compute vertical velocity from continuity equation, integrating over all model-levels. Derived from compw.f: J.E. Haugen, (C) 1995 DNMI

Parameters
nxsize of grid in x direction
nysize of grid in y direction
nzsize of levels
dxx-grid-distance in projection-plane (m) (use radian*R for spherical coordinates)
dyy-grid-distance in projection-plan (m) (use radian*R for spherical coordinates)
gridDistXdistance in m on surface between two grid points in x-direction (nx*ny)
gridDistYdistance in m on surface between two grid points in y-direction (nx*ny)
apap parameter of full hybrid-sigma coordinates in Pa
bb parameter of full hybrid-sigma coordinates (dimensionless)
zssurface-geopotential in m (nx*ny)
pssurface-pressure in Pa (nx*ny)
ux-velocity in m/s (nx*ny*nz)
vy-velocity in m/s (nx*ny*nz)
tabs. temperature in K (nx*ny*nz)
woutput, vertical velocity in m/s, must be preallocated (nx*ny*nz)
Returns
MIFI_OK/MIFI_ERROR

Referenced by mifi_3d_array_position().

int mifi_creepfill2d_f ( size_t  nx,
size_t  ny,
float *  field,
unsigned short  repeat,
char  setWeight,
size_t *  nChanged 
)

Method to fill undefined values in a 2d field in stable time.

This method will fill undefined values by interpolation of neighboring defined values + the average. A value is assumed to be defined if it is defined in the input field, or if it has been defined through the interpolation method (the defined fields will 'creep' into the undefined area).

The results are very similar to mifi_fill2d_f, but the time will vary only with the size of the undefined area, not with the smoothness of the defined values.

Parameters
nxsize of field in x-direction
nysize of field in x-direction
fieldthe data-field to be filled (input/output)
repeatnumber of times values should be re-smoothed (depending on grid-size, 20-100 (linear with time used)).
setWeightdefault weight of original values (versus derived values with weight = 1). Must be >= 1, e.g. 2 the higher the value, the smoother the approxamation from the undefined border to average.
nChangednumber of changed values (output)
Returns
error-code, usually MIFI_OK

Referenced by mifi_3d_array_position(), and MetNoFimex::InterpolatorCreepFill2d::operator()().

int mifi_creepfillval2d_f ( size_t  nx,
size_t  ny,
float *  field,
float  defaultVal,
unsigned short  repeat,
char  setWeight,
size_t *  nChanged 
)

Method to fill undefined values in a 2d field in stable time.

This method will fill undefined values by interpolation of neighboring defined values + the defaultValue. A value is assumed to be defined if it is defined in the input field, or if it has been defined through the interpolation method (the defined fields will 'creep' into the undefined area).

The results are very similar to mifi_fill2d_f, but the time will vary only with the size of the undefined area, not with the smoothness of the defined values.

Parameters
nxsize of field in x-direction
nysize of field in x-direction
fieldthe data-field to be filled (input/output)
defaultValthe defaultValue to be used
repeatnumber of times values should be re-smoothed (depending on grid-size, 20-100 (linear with time used)).
setWeightdefault weight of original values (versus derived values with weight = 1). Must be >= 1, e.g. 2 the higher the value, the smoother the approxamation from the undefined border to average.
nChangednumber of changed values (output)
Returns
error-code, usually MIFI_OK

Referenced by mifi_3d_array_position(), and MetNoFimex::InterpolatorCreepFillVal2d::operator()().

MIFI_DEPRECATED ( int   mifi_isnanffloat val)

check if the value is a nan

Parameters
valthe value to test
Returns
0 on false, otherwise true
Warning
this function should only be used in C++, which doesn't define the isnan macro defined in C99
Deprecated:
use mifi_isnan() template from fimex/Utils.h

Referenced by mifi_3d_array_position().

MIFI_DEPRECATED ( int   mifi_isnanddouble val)

check if the value is a nan

Parameters
valthe value to test
Returns
0 on false, otherwise true
Warning
this function should only be used in C++, which doesn't define the isnan macro defined in C99
Deprecated:
use mifi_isnan() template from fimex/Utils.h
int mifi_fill2d_f ( size_t  nx,
size_t  ny,
float *  field,
float  relaxCrit,
float  corrEff,
size_t  maxLoop,
size_t *  nChanged 
)

Method to fill undefined values in a 2d field.

Solves Laplace's equation with Neumann boundary conditions (dA/dn = 0) in rectangular coordinates by an iterative method to fill-in reasonable values at gridpoints containing values with MIFI_UNDEFINED_F or NaNs

Translated to C from Fortran code by H.Engedahl and A.Foss (1990-93).

Parameters
nxsize of field in x-direction
nysize of field in x-direction
fieldthe data-field to be filled (input/output)
relaxCritrelaxation criteria. Usually 4 orders of magnitude lower than data in field.
corrEffCoef. of overrelaxation, between +1.2 and +2.0
maxLoopMax. allowed no. of scans in relaxation procedure.
nChangednumber of changed values (output)
Returns
error-code, usually MIFI_OK

Referenced by mifi_3d_array_position(), and MetNoFimex::InterpolatorFill2d::operator()().

int mifi_get_values_bicubic_f ( const float *  infield,
float *  outvalues,
const double  x,
const double  y,
const int  ix,
const int  iy,
const int  iz 
)

not implemented yet

The bicubic convolution algorithm assigns a value f(x,y) = X * M * F * Mt * Yt with x, y between (0 <= x < 1), X = (1,x, x2, x3), Y = (1, y, y^2, y^3) and F a 4*4 matrix consisting of the original values of f(-1,-1) to f(2,2).

M is the convolution matrix with a = -0.5 as described by wikipedia (or Catmull-Rom for a = 1, not used here)

Mt and Yt are the transposed matrices/vector.

See also
http://en.wikipedia.org/wiki/Bicubic_interpolation
http://java.sun.com/products/java-media/jai/forDevelopers/jai-apidocs/javax/media/jai/InterpolationBicubic.html
int mifi_get_values_bilinear_f ( const float *  infield,
float *  outvalues,
const double  x,
const double  y,
const int  ix,
const int  iy,
const int  iz 
)

Bilinear interpolation requires a neighborhood extending one pixel to the right and below the central sample. If the fractional subsample position is given by (xfrac, yfrac), the resampled pixel value will be:

  (1 - yfrac) * [(1 - xfrac)*s00 + xfrac*s01] +
  yfrac       * [(1 - xfrac)*s10 + xfrac*s11]

This is documented by the following diagram:

                      s00    s01

                          .      < yfrac

                      s10    s11
                          ^
                         xfrac
See also
http://java.sun.com/products/java-media/jai/forDevelopers/jai-apidocs/javax/media/jai/InterpolationBilinear.html
Warning
if any of the 4 used values of infield is undefined or outside of infield, the return value will be undefined
int mifi_get_values_f ( const float *  infield,
float *  outfield,
const double  x,
const double  y,
const int  ix,
const int  iy,
const int  iz 
)

Get the nearest neighbor of a value. Values are rounded to array-position.

Parameters
infield3d fortran array of size ix,iy,iz
outfield1d array of size iz containing the values
x,y
ix,iy,iz
int mifi_get_values_linear_const_extrapol_f ( const float *  infieldA,
const float *  infieldB,
float *  outfield,
const size_t  n,
const double  a,
const double  b,
const double  x 
)

Same as mifi_get_values_nearest_f() but with constant extrapolation

Parameters
infieldA
infieldB
outfield
n
a
b
x
Returns
int mifi_get_values_linear_d ( const double *  infieldA,
const double *  infieldB,
double *  outfield,
const size_t  n,
const double  a,
const double  b,
const double  x 
)

This is the same as mifi_get_values_linear_f() for double input/output values.

int mifi_get_values_linear_f ( const float *  infieldA,
const float *  infieldB,
float *  outfield,
const size_t  n,
const double  a,
const double  b,
const double  x 
)

Linear interpolation/extrapolation of values in the arrays infieldA and infieldB at position a and b to a field at outfield at position x with o(x) = in(a) + (x - a) * (in(b) - in(a)) / (b - a) (that describes a linear function o(x) = m*x + c )

This interpolation can be used for linear time-interpolation and uses full extrapolation.

Parameters
infieldAarray of size n with values of input at position a
infieldBarray of size n with values of input at position b
outfieldarray of size n with values of input at position x, output
nsize of arrays
aposition of infieldA
bposition of infieldB
xposition of outfield
Returns
MIFI_OK return-value set for compatibility with mifi_get_values_log_f()
int mifi_get_values_linear_no_extrapol_f ( const float *  infieldA,
const float *  infieldB,
float *  outfield,
const size_t  n,
const double  a,
const double  b,
const double  x 
)

Same as mifi_get_values_nearest_f() but without any extrapolation

Parameters
infieldA
infieldB
outfield
n
a
b
x
Returns
int mifi_get_values_linear_weak_extrapol_f ( const float *  infieldA,
const float *  infieldB,
float *  outfield,
const size_t  n,
const double  a,
const double  b,
const double  x 
)

Same as mifi_get_values_nearest_f() but with extrapolation limited to a - (a-b), or b + (a-b).

Parameters
infieldA
infieldB
outfield
n
a
b
x
Returns
int mifi_get_values_log_f ( const float *  infieldA,
const float *  infieldB,
float *  outfield,
const size_t  n,
const double  a,
const double  b,
const double  x 
)

Logarithmic interpolation/extrapolation of values in the arrays infieldA and infieldB at position a and b to a field at outfield at position x with o(x) = m*log(x) + c

This interpolation can be used for i.e. log(p)-interpolation. It is tested against results from ncl int2p and vintp2p_ecmwf log(p) interpolation.

Parameters
infieldAarray of size n with values of input at position a
infieldBarray of size n with values of input at position b
outfieldarray of size n with values of input at position x, output
nsize of arrays
aposition of infieldA
bposition of infieldB
xposition of outfield
Returns
MIFI_OK on success, MIFI_ERROR if log of a, b or x undefined
int mifi_get_values_log_log_f ( const float *  infieldA,
const float *  infieldB,
float *  outfield,
const size_t  n,
const double  a,
const double  b,
const double  x 
)

Log-log interpolation/extrapolation of values in the arrays infieldA and infieldB at position a and b to a field at outfield at position x that describes a function: o(x) = m*log(log(x)+c

This interpolation can be used for i.e. log(log(p))-interpolation.

Warning
It is tested against results from ncl vintp2p_ecmwf log(log(p)) interpolation, but results vary slightly (~1%) for unknown reason.
Parameters
infieldAarray of size n with values of input at position a
infieldBarray of size n with values of input at position b
outfieldarray of size n with values of input at position x, output
nsize of arrays
aposition of infieldA
bposition of infieldB
xposition of outfield
Returns
MIFI_OK on success, MIFI_ERROR if log of a, b or x undefined
int mifi_get_values_nearest_f ( const float *  infieldA,
const float *  infieldB,
float *  outfield,
const size_t  n,
const double  a,
const double  b,
const double  x 
)

NearestNeighbor "interpolation"/extrapolation of values in the arrays infieldA and infieldB at position a and b to a field at outfield at position x

Parameters
infieldAarray of size n with values of input at position a
infieldBarray of size n with values of input at position b
outfieldarray of size n with values of input at position x, output
nsize of arrays
aposition of infieldA
bposition of infieldB
xposition of outfield
Returns
MIFI_OK return-value set for compatibility with mifi_get_values_log_f()
int mifi_get_values_no_extrapol_f ( const float *  infield,
float *  outfield,
const double  x,
const double  y,
const int  ix,
const int  iy,
const int  iz 
)

Get the nearest neighbor of a value. Values are rounded to array-position. No extrapolation will occur.

Parameters
infield3d fortran array of size ix,iy,iz
outfield1d array of size iz containing the values
x,y
ix,iy,iz
int mifi_get_values_weak_extrapol_f ( const float *  infield,
float *  outfield,
const double  x,
const double  y,
const int  ix,
const int  iy,
const int  iz 
)

Get the nearest neighbor of a value. Values are rounded to array-position. Extrapolates to once the distance of the two leftmost values (or rightmost).

Parameters
infield3d fortran array of size ix,iy,iz
outfield1d array of size iz containing the values
x,y
ix,iy,iz
int mifi_get_vector_reproject_matrix ( const char *  proj_input,
const char *  proj_output,
const double *  out_x_axis,
const double *  out_y_axis,
int  out_x_axis_type,
int  out_y_axis_type,
int  ox,
int  oy,
double *  matrix 
)

calculate the vector reprojection matrix used in mifi_vector_reproject_values_f

Parameters
proj_inputproj4-string of projection of infield
proj_outputproj4-string of projection of outfield
out_x_axisfield of size ox. Axis needs to be strong monotonous and if longitude/latitude in degree
out_y_axisfield of size oy. Axis needs to be strong monotonous and if longitude/latitude in degree
out_x_axis_typeone of MIFI_LATITUDE, MIFI_LONGITUDE, MIFI_PROJ_AXIS
out_y_axis_typeone of MIFI_LATITUDE, MIFI_LONGITUDE, MIFI_PROJ_AXIS
oxx-dimension of outfield
oyy-dimension of outfield
matrixmatrix of size (4*ox*oy)
Returns
MIFI_OK or error value
int mifi_get_vector_reproject_matrix_field ( const char *  proj_input,
const char *  proj_output,
const double *  in_x_field,
const double *  in_y_field,
int  ox,
int  oy,
double *  matrix 
)

calculate the vector reprojection matrix used in mifi_vector_reproject_values_f without changing the axes, just the direction of the vectors.

Parameters
proj_inputproj4-string of projection of infield
proj_outputproj4-string of projection of outfield
in_x_fieldfield of size ox*oy with the values in the input-projection in_x_field[x+oy*y] = inXAxis[x]. The values must be in radian or m.
in_y_fieldfield of size ox*oy with the values in the input-projection in_y_field[x+oy*y] = inYAxis[y]. The values must be in radian or m.
oxx-dimension of in/outfield
oyy-dimension of in/outfield
matrixmatrix of size (4*ox*oy)
Returns
MIFI_OK or error value
int mifi_get_vector_reproject_matrix_points ( const char *  proj_input,
const char *  proj_output,
int  inputIsMetric,
const double *  out_x_points,
const double *  out_y_points,
int  on,
double *  matrix 
)

calculate the vector reprojection matrix when projecting to a list of n-points

Parameters
proj_inputproj4-string of projection of infield
proj_outputproj4-string of projection of outfield
inputIsMetric1 if projection, 0 if (rotated) latlon
out_x_pointsoutput-points values in the output-projection out_x_points[on]. The values must be in radian or m.
out_y_pointsoutput-points values in the output-projection out_y_points[on]. The values must be in radian or m.
onnumber of output-points
matrixpreallocated matrix of size (4*on)
Returns
MIFI_OK or error value
int mifi_griddistance ( size_t  nx,
size_t  ny,
const double *  lonVals,
const double *  latVals,
float *  gridDistX,
float *  gridDistY 
)

Calculate the real distance in m between neigboring grid-cells (center to center). Distances are calculated using the great-circle distance. The size of a gridcell is gridDistX*gridDistY.

Parameters
nxpoints in x direction
nypoints in y direction
lonValslongitude values in degree (j*ny + i)
latValslatitude values in degree (j*ny + i)
gridDistXoutput x-distance between neighbors in m (j*ny +i)
gridDistYoutput y-distance between neighbors in m (j*ny +i)
Returns
error-code or MIFI_OK

Referenced by mifi_3d_array_position().

int mifi_interpolate_d ( int  method,
char *  proj_input,
double *  infield,
double *  in_x_axis,
double *  in_y_axis,
int  in_x_axis_type,
int  in_y_axis_type,
int  ix,
int  iy,
int  iz,
char *  proj_output,
double *  outfield,
double *  out_x_axis,
double *  out_y_axis,
int  out_x_axis_type,
int  out_y_axis_type,
int  ox,
int  oy 
)

not implemented yet

double version of mifi_interpolate_f

See also
mifi_interpolate_f
int mifi_interpolate_f ( int  method,
const char *  proj_input,
const float *  infield,
const double *  in_x_axis,
const double *  in_y_axis,
const int  in_x_axis_type,
const int  in_y_axis_type,
const int  ix,
const int  iy,
const int  iz,
const char *  proj_output,
float *  outfield,
const double *  out_x_axis,
const double *  out_y_axis,
const int  out_x_axis_type,
const int  out_y_axis_type,
const int  ox,
const int  oy 
)

Interpolation between two projections. Missing values are set to MIFI_UNDEFINED_F which is implemented as C99 nanf. The coordinates of a cell give the midpoint of a cell, i.e. cell (10,20) spans ([9.5..10.5[,[19.5-20.5[)

Parameters
methodone of MIFI_INTERPOL_NEAREST_NEIGHBOR MIFI_INTERPOL_BILINEAR MIFI_INTERPOL_BICUBIC
proj_inputproj4-string of projection of infield
infieldreal rectangular array of dimension infield[iz,iy,ix]
in_x_axisfield of size ix. Axis needs to be strong monotonous and if longitude/latitude in degree
in_y_axisfield of size iy. Axis needs to be strong monotonous and if longitude/latitude in degree
in_x_axis_typeone of MIFI_LATITUDE, MIFI_LONGITUDE, MIFI_PROJ_AXIS
in_y_axis_typeone of MIFI_LATITUDE, MIFI_LONGITUDE, MIFI_PROJ_AXIS
ixx-dimension of infield
iyy-dimension of infield
izz-dimension of infield and outfield. The z-dim allows you to convert several fields at once without calculating the projection again and again.
proj_outputproj4-string of projection of outfield
outfieldreal rectangular array of dimension outfield[iz,oy,ox]
out_x_axisfield of size ox. Axis needs to be strong monotonous and if longitude/latitude in degree
out_y_axisfield of size oy. Axis needs to be strong monotonous and if longitude/latitude in degree
out_x_axis_typeone of MIFI_LATITUDE, MIFI_LONGITUDE, MIFI_PROJ_AXIS
out_y_axis_typeone of MIFI_LATITUDE, MIFI_LONGITUDE, MIFI_PROJ_AXIS
oxx-dimension of outfield
oyy-dimension of outfield
size_t mifi_nanf2bad ( float *  posPtr,
float *  endPtr,
float  badVal 
)

Convert nan back to bad-values. See mifi_bad2nanf

Parameters
posPtrstart pointer of the float array
endPtrend-pointer of the float array (excluded from conversion)
badValvalue NaNs will be converted to
Returns
0 (in fimex < 0.61, this was possibly the number of conversions, but it was never strict)

Referenced by mifi_3d_array_position().

int mifi_points2position ( double *  points,
const int  n,
const double *  axis,
const int  num,
const int  axis_type 
)

find position in array of position in projection

Points2position uses linear splines to find the array-position of points in the given axis. Undefined (NaN/Inf) values will be translated to -999 to be far outside the position-range (0-(num-1)).

Parameters
pointsthe values will get changed from points in axis coordinates to array coordinates
nnumber of values in points
axiscoordinate axis
numnumber of elements in coordinate axis
axis_typetype of axis, one of MIFI_LONGITUDE, MIFI_LATITUDE, MIFI_PROJ_AXIS
int mifi_project_axes ( const char *  proj_input,
const char *  proj_output,
const double *  in_x_axis,
const double *  in_y_axis,
const int  ix,
const int  iy,
double *  out_xproj_axis,
double *  out_yproj_axis 
)

project axes so that the projetion (x,y) => (x_proj), (y_proj) can be expressed as x_proj(x,y), y_proj(x,y)

all axes must be given or will be returned in radians when converted from/to latlon

Parameters
proj_inputinput projection proj string
proj_outputoutput projection proj string
in_x_axisx-axis in input-projection
in_y_axisy-axis in input-projection
ixsize of x-axis
iysize of y-axis
out_xproj_axisoutput-values of x_proj(x,y), field needs to be allocated in at least ix*iy size
out_yproj_axisoutput-values of y_proj(x,y), field needs to be allocated in at least ix*iy size
Returns
error-code

Referenced by mifi_3d_array_position().

int mifi_project_values ( const char *  proj_input,
const char *  proj_output,
double *  in_out_x_vals,
double *  in_out_y_vals,
const int  num 
)

project values so that the projetion (x,y) => (x_proj), (y_proj) can be expressed as x_proj(x,y), y_proj(x,y)

all values must be given or will be returned in radians when converted from/to latlon

Parameters
proj_inputinput projection proj string
proj_outputoutput projection proj string
in_out_x_valsx-values, will be input and output
in_out_y_valsy-values, will be input and output
numsize of arrays
Returns
error-code

Referenced by mifi_3d_array_position().

int mifi_string_to_interpolation_method ( const char *  stringMethod)

Convert interpolation methods in string-format to mifi_interpol_method (see mifi_constants.h)

Parameters
stringMethodone of nearestneighbor, bilinear,bicubic, forward_max, forward_min, forward_mean, forward_median, forward_sum, coord_nearestneighbor, coord_kdtree
Returns
mifi_interpol_method enum, or MIFI_INTERPOL_UNKNOWN on error
int mifi_vector_reproject_direction_by_matrix_f ( int  method,
const double *  matrix,
float *  angle_out,
int  ox,
int  oy,
int  oz 
)

Calculate the reprojected directions with a known matrix. Directions are the angle between the projections y-Axis (0degree) clockwise to 360degree.

Parameters
method(one of MIFI_VECTOR_KEEP_SIZE) (MIFI_VECTOR_RESIZE no longer available)
matrixreprojection matrix of size (4,ox,oy)
angle_outangles in degrees at position in the output-projection (i.e. by prevously applying mifi_interpolate_f). The values here will be changed!
oxx-dimension of outfield
oyy-dimension of outfield
ozz-dimension of the outfield
Returns
MIFI_OK or error value
int mifi_vector_reproject_values_by_matrix_f ( int  method,
const double *  matrix,
float *  u_out,
float *  v_out,
int  ox,
int  oy,
int  oz 
)

calculate the reprojected vectors with a known matrix for mifi_vector_reproject_values_f

Parameters
method(one of MIFI_VECTOR_KEEP_SIZE) (MIFI_VECTOR_RESIZE no longer available)
matrixreprojection matrix of size (4,ox,oy)
u_outvalues of u, with position in the output-projection (i.e. by prevously applying mifi_interpolate_f). The values here will be changed!
v_outvalues of v, with position in the output-projection (i.e. by prevously applying mifi_interpolate_f). The values here will be changed!
oxx-dimension of outfield
oyy-dimension of outfield
ozz-dimension of the outfield
Returns
MIFI_OK or error value
int mifi_vector_reproject_values_f ( int  method,
const char *  proj_input,
const char *  proj_output,
float *  u_out,
float *  v_out,
const double *  out_x_axis,
const double *  out_y_axis,
int  out_x_axis_type,
int  out_y_axis_type,
int  ox,
int  oy,
int  oz 
)

interpolate the vector values

When reprojecting a vector (i.e. wind (u, v)) from one projection to another, not only the base-position of the vector will change, but also the angle of the vector might change due to rotation and streching within the projection. Thus, the values of (u,v) have to be changed accordingly to projection.

This function allows to only rotate the vector values (MIFI_VECTOR_KEEP_SIZE) which is useful to keep the windspeed constant, even if the projected plane has a different scale, or to completely reproject the vector (MIFI_VECTOR_RESIZE).

This function is implemented by using a first order tailor expansion of the projection: (u', v') = A (u,v) with A a matrix defined at each point (x,y) through

1 proj(x,y)_x' = a11*x+a21*y
2 proj(x,y)_y' = a12*x+a22*y

and the same formulas for (x+delta, y) and (x, y+delta) (with delta a small value against the x or y)

Parameters
method(one of MIFI_VECTOR_KEEP_SIZE, MIFI_VECTOR_RESIZE)
proj_inputproj4-string of projection of infield
proj_outputproj4-string of projection of outfield
u_outvalues of u, with position in the output-projection (i.e. by prevously applying mifi_interpolate_f). The values here will be changed!
v_outvalues of v, with position in the output-projection (i.e. by prevously applying mifi_interpolate_f). The values here will be changed!
out_x_axisfield of size ox. Axis needs to be strong monotonous and if longitude/latitude in degree
out_y_axisfield of size oy. Axis needs to be strong monotonous and if longitude/latitude in degree
out_x_axis_typeone of MIFI_LATITUDE, MIFI_LONGITUDE, MIFI_PROJ_AXIS
out_y_axis_typeone of MIFI_LATITUDE, MIFI_LONGITUDE, MIFI_PROJ_AXIS
oxx-dimension of outfield
oyy-dimension of outfield
ozz-dimension of the outfield
Returns
MIFI_OK or error value