Interpolation of scattered data

Scattered data interpolation enables approximation of functions from unstructured point sets where data points are not aligned to a regular grid. The Interpolation class provides factory methods for creating scattered interpolators in 2D and N dimensions, including nearest-neighbor, radial basis function (RBF), and triangulated linear methods.

What is Scattered Data?

Scattered data consists of function samples at irregularly distributed points that do not lie on a structured grid. Unlike grid-based data, scattered points can be located anywhere in the domain, making them ideal for:

  • Sensor measurements at arbitrary locations

  • GPS coordinates with elevation data

  • Weather station observations

  • Laboratory experiments with non-uniform sampling

  • Finite element analysis node data

Scattered interpolation methods construct a continuous function that passes through (or approximates) the given data points, allowing evaluation at arbitrary query locations within or outside the sampled region.

2D vs. N-Dimensional Methods: Some methods are specialized for 2D data and exploit 2D geometry (e.g., Delaunay triangulation), while others work in arbitrary dimensions (e.g., nearest-neighbor and RBF interpolation).

Supported Scattered Interpolation Methods

The scattered interpolation API provides three main methods:

Scattered interpolation method comparison

Method

Dimensions

Smoothness

Domain

Use Cases

Nearest-neighbor

Any

C0 (piecewise constant)

Global (all space)

Fast lookups, categorical data

Triangulated linear (2D)

2D only

C0 (continuous)

Convex hull of points

General 2D interpolation, visualization

Radial Basis Function (RBF)

Any

C∞ (smooth)

Global (all space)

Smooth approximations, geospatial modeling

Nearest-Neighbor Interpolation:

The NearestScatterSurface returns the value of the closest data point to the query location. It is piecewise constant, creating Voronoi-like regions around each data point. This method is extremely fast and works in any number of dimensions, but produces discontinuous results. Use for categorical data or when computational speed is paramount.

Triangulated Linear Interpolation (2D):

The TriangulatedSurface2D performs Delaunay triangulation of the 2D point set and uses piecewise linear interpolation within each triangle. The result is continuous but has discontinuous gradients at triangle edges. Triangulated interpolation is well-suited for general 2D scattered data and is commonly used in computer graphics and geographic information systems (GIS).

The interpolant is only defined over the triangulation (typically the convex hull of the data points). Queries outside this region can be handled using ExtrapolationMode.

Radial Basis Function (RBF) Interpolation:

The RbfSurface constructs a smooth global interpolant using radial basis functions. RBF methods represent the interpolant as a weighted sum of radial functions centered at each data point. Common radial basis functions include Gaussian, multiquadric (MQ), inverse multiquadric (IMQ), and thin-plate splines.

RBF interpolation produces smooth, infinitely differentiable surfaces and works in arbitrary dimensions. However, it requires solving a linear system of size equal to the number of data points, which can be computationally expensive for large datasets. RBF interpolation is global: the interpolant is defined everywhere in space, and there is no concept of extrapolation.

Creating Scattered Interpolators

Scattered interpolators are created using factory methods in the Interpolation class. The canonical representation for N-dimensional data is an Array2D<T> of shape (dimensions, nPoints) along with a 1D array of values. Convenience overloads are provided for common cases (2D and 3D with separate coordinate arrays).

Nearest-Neighbor Interpolation:

To create a nearest-neighbor interpolator, pass the scattered point coordinates and corresponding function values:

C# 
// Define scattered 2D points using Vectors
Vector<double> x = Vector.Create(0.0, 1.0, 2.0, 0.5, 1.5);
Vector<double> y = Vector.Create(0.0, 0.0, 0.0, 1.0, 1.0);
Vector<double> values = Vector.Create(1.0, 2.0, 3.0, 1.5, 2.5);

// Create nearest-neighbor interpolator
var surface = Numerics.NET.Interpolation.ScatteredNearestInterpolator(
    x, y, values);

// Evaluate at query points
double z1 = surface.Evaluate(0.6, 0.2);
double z2 = surface.Evaluate(1.2, 0.8);
Console.WriteLine($"Nearest neighbor at (0.6, 0.2): {z1}");
Console.WriteLine($"Nearest neighbor at (1.2, 0.8): {z2}");

Radial Basis Function (RBF) Interpolation:

RBF interpolators require the same inputs as nearest-neighbor, but also allow specification of the radial basis function type and shape parameter:

C# 
// Define scattered 2D points using Vectors
Vector<double> x = Vector.Create(0.0, 1.0, 2.0, 0.5, 1.5, 1.0);
Vector<double> y = Vector.Create(0.0, 0.0, 0.0, 1.0, 1.0, 0.5);
Vector<double> values = Vector.Create(1.0, 2.0, 3.0, 1.5, 2.5, 2.0);

// Create RBF interpolator with default radial basis function
var surface = Numerics.NET.Interpolation.ScatteredRbfInterpolator(
    x, y, values);

// Evaluate at query points - RBF provides smooth interpolation
double z1 = surface.Evaluate(0.5, 0.5);
double z2 = surface.Evaluate(1.5, 0.5);
Console.WriteLine($"RBF interpolation at (0.5, 0.5): {z1:F3}");
Console.WriteLine($"RBF interpolation at (1.5, 0.5): {z2:F3}");

RBF interpolation also supports N-dimensional data:

C# 
// RBF also works in N dimensions
// Define 3D scattered points using Matrix (each row is a point)
Matrix<double> points = Matrix.Create(new double[,] {
    { 0.0, 0.0, 0.0 },
    { 1.0, 0.0, 0.0 },
    { 0.0, 1.0, 0.0 },
    { 0.0, 0.0, 1.0 },
    { 0.5, 0.5, 0.5 }
});
Vector<double> values3D = Vector.Create(1.0, 2.0, 3.0, 4.0, 2.5);

var surface3D = Numerics.NET.Interpolation.ScatteredRbfInterpolator(
    points, values3D);

Vector<double> queryPoint = Vector.Create(0.3, 0.3, 0.3);
double znd = surface3D.Evaluate(queryPoint);
Console.WriteLine($"3D RBF interpolation: {znd:F3}");

Triangulated Linear Interpolation (2D):

For 2D scattered data, triangulated linear interpolation performs Delaunay triangulation and uses piecewise linear interpolation within triangles:

C# 
// Define scattered 2D points using Vectors
Vector<double> x = Vector.Create(0.0, 2.0, 1.0, 0.5, 1.5);
Vector<double> y = Vector.Create(0.0, 0.0, 1.5, 0.8, 0.8);
Vector<double> values = Vector.Create(1.0, 3.0, 4.0, 2.0, 3.5);

// Create triangulated linear interpolator
// This performs Delaunay triangulation and linear interpolation within triangles
var surface = Numerics.NET.Interpolation.ScatteredLinear2DInterpolator(
    x, y, values);

// Evaluate at query points inside the convex hull
double z1 = surface.Evaluate(1.0, 0.5);
double z2 = surface.Evaluate(1.2, 0.9);
Console.WriteLine($"Triangulated interpolation at (1.0, 0.5): {z1:F3}");
Console.WriteLine($"Triangulated interpolation at (1.2, 0.9): {z2:F3}");

Evaluating Scattered Interpolators

Evaluation of scattered interpolators follows the same pattern as grid interpolators. For 2D surfaces, use surface.Evaluate(x, y). For N-dimensional surfaces, pass a coordinate array or span.

Performance characteristics differ by method:

  • Nearest-neighbor: Currently uses brute-force search with O(n) evaluation cost per query, where n is the number of data points. Spatial indexing acceleration is planned for future versions.

  • Triangulated linear: Uses point location in the triangulation, typically O(log n) for a well-structured triangulation.

  • RBF: Evaluates a weighted sum of radial basis functions, which is O(n) per query, where n is the number of data points.

Derivatives and Gradients

Derivative support depends on the interpolation method:

  • Nearest-neighbor: Gradient and Hessian are effectively zero (piecewise constant function).

  • Triangulated linear: Gradient is constant within each triangle and discontinuous at triangle edges. Hessian is zero within triangles.

  • RBF: Gradient and Hessian are available analytically for most radial basis functions and are smooth everywhere.

The following example demonstrates gradient computation for an RBF surface:

C# 
// Define scattered 2D points using Vectors
Vector<double> x = Vector.Create(0.0, 1.0, 2.0, 0.5, 1.5, 1.0);
Vector<double> y = Vector.Create(0.0, 0.0, 0.0, 1.0, 1.0, 0.5);
Vector<double> values = Vector.Create(1.0, 2.0, 3.0, 1.5, 2.5, 2.0);

// RBF surfaces support gradients
var surface = Numerics.NET.Interpolation.ScatteredRbfInterpolator(
    x, y, values);

// Compute gradient at a point
Vector<double> gradient = Vector.Create(0.0, 0.0);
surface.Gradient(1.0, 0.5, gradient);
Console.WriteLine($"RBF gradient at (1.0, 0.5): [{gradient[0]:F3}, {gradient[1]:F3}]");

Domain and Extrapolation

The concept of "domain" and "extrapolation" differs between scattered interpolation methods:

  • Nearest-neighbor and RBF: These methods are global, meaning the interpolant is defined everywhere in space. There is no concept of "out of bounds" or extrapolation—every query point has a defined value.

  • Triangulated linear (2D): The interpolant is only defined over the convex hull (or triangulation) of the data points. Queries outside this region can be handled using ExtrapolationMode:

    • Throw: Throws an exception for out-of-triangulation queries.

    • Constant: Returns a fill value (typically NaN) for out-of-triangulation queries.

    • Clamp: Clamps query coordinates to the triangulation boundary and interpolates at the boundary.

The triangulation boundary is determined geometrically. For a convex hull, a point is "inside" if it lies within the convex region spanned by the data points.

  Note

For RBF interpolation, the interpolant smoothly extrapolates to infinity based on the radial basis function's asymptotic behavior. For Gaussian RBF, values approach zero at large distances. For multiquadric, values grow unboundedly.

Examples

Example 1: Interpolating scattered elevation data

Suppose you have GPS coordinates with elevation measurements at scattered locations. You can use triangulated linear interpolation for a quick, continuous surface:

C# 
// Define scattered 2D points using Vectors
Vector<double> x = Vector.Create(0.0, 2.0, 1.0, 0.5, 1.5);
Vector<double> y = Vector.Create(0.0, 0.0, 1.5, 0.8, 0.8);
Vector<double> values = Vector.Create(1.0, 3.0, 4.0, 2.0, 3.5);

// Create triangulated linear interpolator
// This performs Delaunay triangulation and linear interpolation within triangles
var surface = Numerics.NET.Interpolation.ScatteredLinear2DInterpolator(
    x, y, values);

// Evaluate at query points inside the convex hull
double z1 = surface.Evaluate(1.0, 0.5);
double z2 = surface.Evaluate(1.2, 0.9);
Console.WriteLine($"Triangulated interpolation at (1.0, 0.5): {z1:F3}");
Console.WriteLine($"Triangulated interpolation at (1.2, 0.9): {z2:F3}");

For smoother results, use RBF interpolation:

C# 
// Define scattered 2D points using Vectors
Vector<double> x = Vector.Create(0.0, 1.0, 2.0, 0.5, 1.5, 1.0);
Vector<double> y = Vector.Create(0.0, 0.0, 0.0, 1.0, 1.0, 0.5);
Vector<double> values = Vector.Create(1.0, 2.0, 3.0, 1.5, 2.5, 2.0);

// Create RBF interpolator with default radial basis function
var surface = Numerics.NET.Interpolation.ScatteredRbfInterpolator(
    x, y, values);

// Evaluate at query points - RBF provides smooth interpolation
double z1 = surface.Evaluate(0.5, 0.5);
double z2 = surface.Evaluate(1.5, 0.5);
Console.WriteLine($"RBF interpolation at (0.5, 0.5): {z1:F3}");
Console.WriteLine($"RBF interpolation at (1.5, 0.5): {z2:F3}");

Example 2: N-dimensional sensor data

For scattered sensor measurements in 3D or higher dimensions, RBF interpolation provides a global smooth approximation:

C# 
// RBF also works in N dimensions
// Define 3D scattered points using Matrix (each row is a point)
Matrix<double> points = Matrix.Create(new double[,] {
    { 0.0, 0.0, 0.0 },
    { 1.0, 0.0, 0.0 },
    { 0.0, 1.0, 0.0 },
    { 0.0, 0.0, 1.0 },
    { 0.5, 0.5, 0.5 }
});
Vector<double> values3D = Vector.Create(1.0, 2.0, 3.0, 4.0, 2.5);

var surface3D = Numerics.NET.Interpolation.ScatteredRbfInterpolator(
    points, values3D);

Vector<double> queryPoint = Vector.Create(0.3, 0.3, 0.3);
double znd = surface3D.Evaluate(queryPoint);
Console.WriteLine($"3D RBF interpolation: {znd:F3}");

This approach generalizes to any number of dimensions, making it suitable for multidimensional parameter spaces and high-dimensional data analysis.

Example 3: Comparing scattered methods

For a given 2D dataset, you can compare the three methods:

  • Nearest-neighbor: Fast, discontinuous.

  • Triangulated linear: Continuous, fast, linear within triangles.

  • RBF: Smooth, global, slower for large datasets.

The choice depends on your data characteristics, desired smoothness, and performance requirements. For visualization and smooth derivative calculations, prefer RBF or cubic splines (if data can be gridded). For general-purpose 2D interpolation with reasonable performance, use triangulated linear.

See Also