29.1 Одномерная интерполяция

Octave поддерживает несколько методов одномерной интерполяции, большинство из которых описано в этом разделе. Дополнительные методы описаны в разделах Полиномиальная интерполяция и Интерполяция по разрозненным данным.

: yi = interp1 (x, y, xi)

: yi = interp1 (y, xi)

: yi = interp1 (…, method)

: yi = interp1 (…, extrap)

: yi = interp1 (…, "left")

: yi = interp1 (…, "right")

: pp = interp1 (…, "pp")

Одномерная интерполяция.

Интерполирует входные данные, чтобы определить значение yi в точках xi. Если x не задан, в качестве него используются индексы y (1:length (y)). Если y — матрица или N-мерный массив, интерполяция выполняется для каждого столбца y.

Метод интерполяции method может иметь одно из следующих значений:

"nearest"

Возвращает ближайший соседний элемент.

"previous"

Возвращает предыдущий соседний элемент.

"next"

Возвращает следующий соседний элемент.

"linear" (по умолчанию)

Линейная интерполяция по ближайшим соседним элементам.

"pchip"

Кусочно-кубический интерполяционный полином Эрмита — сохраняющая форму интерполяция с непрерывной первой производной.

"cubic"

Кубическая интерполяция (то же, что "pchip").

"spline"

Интерполяция кубическими сплайнами — непрерывность первой и второй производных на всей кривой.

Добавление символа ’*’ в начало любого из перечисленных выше методов заставляет interp1 считать, что точки x расположены равномерно, и использовать только x(1) и x(2). Обычно это работает быстрее и никогда не медленнее. Метод по умолчанию — "linear".

Если extrap — строка "extrap", значения за пределами конечных точек экстраполируются с использованием текущего метода method. Если extrap — число, значения за пределами конечных точек заменяются этим числом. Если extrap не задан, по умолчанию используется NA.

Если задан строковый аргумент "pp", аргумент xi указывать не следует, а interp1 возвращает кусочно-полиномиальный объект. Позже этот объект можно использовать с ppval для вычисления интерполяции. Существует эквивалентная запись: ppval (interp1 (x, y, method, "pp"), xi) == interp1 (x, y, xi, method, "extrap").

Повторяющиеся точки в x задают разрывный интерполянт. Не допускается более двух последовательных точек с одинаковым значением. Если x возрастает, разрывный интерполянт по умолчанию непрерывен справа. Если x убывает, разрывный интерполянт по умолчанию непрерывен слева. Условие непрерывности интерполянта можно задать параметрами "left" или "right", чтобы выбрать интерполянт, непрерывный слева или справа соответственно. Разрывная интерполяция разрешена только для методов "nearest" и "linear"; во всех остальных случаях значения x должны быть уникальными.

Пример использования interp1:

xf = [0:0.05:10];
yf = sin (2*pi*xf/5);
xp = [0:10];
yp = sin (2*pi*xp/5);
lin = interp1 (xp, yp, xf);
near = interp1 (xp, yp, xf, "nearest");
pch = interp1 (xp, yp, xf, "pchip");
spl = interp1 (xp, yp, xf, "spline");
plot (xf,yf,"r", xf,near,"g", xf,lin,"b", xf,pch,"c", xf,spl,"m",
      xp,yp,"r*");
legend ("original", "nearest", "linear", "pchip", "spline");

См. также: pchip, spline, interpft, interp2, interp3, interpn.

Различные методы интерполяции имеют несколько важных отличий. Метод "spline" обеспечивает непрерывность как первой, так и второй производных интерполированных значений, тогда как остальные методы этого не делают. Это означает, что результаты метода "spline" обычно более гладкие. Если интерполируемая функция действительно гладкая, "spline" даст отличные результаты. Однако если функция, значения которой нужно вычислить, в некотором смысле разрывна, интерполяция "pchip" может дать лучшие результаты.

Это можно продемонстрировать с помощью кода:

t = -2:2;
dt = 1;
ti =-2:0.025:2;
dti = 0.025;
y = sign (t);
ys = interp1 (t,y,ti,"spline");
yp = interp1 (t,y,ti,"pchip");
ddys = diff (diff (ys)./dti) ./ dti;
ddyp = diff (diff (yp)./dti) ./ dti;
figure (1);
plot (ti,ys,"r-", ti,yp,"g-");
legend ("spline", "pchip", 4);
figure (2);
plot (ti,ddys,"r+", ti,ddyp,"g*");
legend ("spline", "pchip");

Результат показан на рисунке 29.1 и рисунке 29.2.

Рисунок 29.1: Сравнение методов интерполяции "pchip" и "spline" для ступенчатой функции

Рисунок 29.2: Сравнение вторых производных методов интерполяции "pchip" и "spline" для ступенчатой функции

Интерполяция Фурье — это метод передискретизации, при котором сигнал преобразуется в частотную область, дополняется нулями, а затем преобразуется обратно во временную область.

: y = interpft (x, n)

: y = interpft (x, n, dim)

Интерполяция Фурье.

Если x — вектор, он передискретизируется с использованием n точек. Предполагается, что данные в x расположены равномерно. Если x — матрица или N-мерный массив, интерполяция выполняется для каждого столбца x.

Если задан dim, интерполяция выполняется вдоль измерения dim.

interpft предполагает, что интерполируемая функция периодична, поэтому для конечных точек интерполяции делаются соответствующие предположения.

См. также: interp1.

У интерполяции Фурье есть два существенных ограничения. Во-первых, предполагается, что сигнал является периодическим, поэтому непериодические сигналы будут плохо представлены на краях. Во-вторых, и сигнал, и результат его интерполяции должны быть представлены в равномерно расположенных точках. Пример использования interpft:

t = 0 : 0.3 : pi; dt = t(2)-t(1);
n = length (t); k = 100;
ti = t(1) + [0 : k-1]*dt*n/k;
y = sin (4*t + 0.3) .* cos (3*t - 0.1);
yp = sin (4*ti + 0.3) .* cos (3*ti - 0.1);
plot (ti, yp, "g", ti, interp1 (t, y, ti, "spline"), "b", ...
      ti, interpft (y, k), "c", t, y, "r+");
legend ("sin(4t+0.3)cos(3t-0.1)", "spline", "interpft", "data");

Он демонстрирует неудовлетворительную работу интерполяции Фурье для непериодических функций, что видно на рисунке 29.3.

Рисунок 29.3: Сравнение interp1 и interpft для непериодических данных

Кроме того, вспомогательные функции spline и lookup, лежащие в основе функции interp1, можно вызывать напрямую.

: pp = spline (x, y)

: yi = spline (x, y, xi)

Возвращает интерполянт в виде кубического сплайна для точек x и y.

При вызове с двумя аргументами возвращает кусочно-полиномиальный объект pp, который можно использовать с ppval для вычисления полинома в заданных точках.

При вызове с третьим входным аргументом spline вычисляет сплайн в точках xi. Третья форма вызова spline (x, y, xi) эквивалентна ppval (spline (x, y), xi).

Переменная x должна быть вектором длины n.

y может быть вектором или массивом. Если y — вектор, его длина должна быть равна либо n, либо n + 2. Если длина y равна n, используется граничное условие "not-a-knot". Если длина y равна n + 2, первое и последнее значения вектора y задают значения первой производной кубического сплайна в конечных точках.

Если y — массив, его размерность должна иметь вид [s1, s2, …, sk, n] или [s1, s2, …, sk, n + 2]. Массив внутри функции преобразуется в матрицу, ведущая размерность которой задается выражением s1 * s2 * … * sk, после чего каждая строка этой матрицы обрабатывается отдельно. Обратите внимание, что это прямо противоположно interp1, но сделано для совместимости с MATLAB.

См. также: pchip, ppval, mkpp, unmkpp.