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 задают разрывную интерполянт. Может быть не более 2 последовательных точек с одинаковым значением. Если 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.

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

: y = interpft (x, n) ¶

: y = interpft (x, n, dim) ¶

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

Если x — вектор, то 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.

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

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