Matlab / Лекция 19

Лекция № 19

Matlab: Функции для решения задач интерполяции и аппроксимации

1. Функции для решения задач интерполяции

2. Функции для решения задач аппроксимации

1. Функции для решения задач интерполяции

В Matlab для решения задач интерполяции имеются следующие функции:

spline — интерполяция функции одной переменной кубическим сплайном;

interp1 — одномерная табличная интерполяция;

interp2 — двумерная табличная интерполяция;

griddata — двумерная табличная интерполяция на неравномерной сетке.

Интерполяция функции одной переменной кубическим сплайном. Функция yi = spline(x, y, xi) интерполирует значения функции y в точках xi внутри области определения функции, используя кубические сплайны.

Пример 1. Зададим синусоиду 11 точками и проведем интерполяцию кубическими сплайнами, используя мелкую сетку (рис. 1):

x = 0:10;

y = sin(x);

xi = 0:.25:10;

yi = spline(x, y, xi);

plot(x, y, 'o', xi, yi, 'g')

grid on

Рис. 1.

Большие возможности для работы со сплайнами предоставляет пользователю специализированный пакет Spline Toolbox.

Одномерная табличная интерполяция. Функция yi = interp1(x, y, xi) строит интерполирующую кривую для одномерного массива y, заданного на сетке x; выходной массив yi может быть определен на более мелкой сетке xi. Если Y — двумерный массив, то интерполирующая кривая строится для каждого столбца. По умолчанию реализована линейная интерполяция.

Функция yi = interp1(x, y, xi, ‘метод‘) позволяет задать метод интерполяции:

‘linear’	линейная
‘cubic’	кубическая
‘spline’	кубические сплайны

Принято, что аргумент x изменяется монотонно; кроме того, для кубической интерполяции предполагается, что сетка по x равномерна.

Пример 2. Зададим синусоиду 11 точками и проведем интерполяцию, используя мелкую сетку (рис. 2):

x = 0:10;

y = sin(x);

xi = 0:.25:10;

yi = interp1(x, y, xi);

plot(x, y, 'o', xi, yi, 'g')

hold on

yi = interp1(x, y, xi, 'spline' );

plot(x, y, 'ob', xi, yi, 'm')

grid on

Рис. 2.

Двумерная табличная интерполяция. Функция ZI = interp2(X, Y, Z, XI, YI) интерполирует данные, определяющие некоторую поверхность на двумерной сетке (X, Y); выходной массив ZI может быть определен на более мелкой сетке (XI, YI). По умолчанию реализована линейная интерполяция.

Функция ZI = interp2(X, Y, Z, XI, YI, ‘метод‘) позволяет задать метод интерполяции:

‘linear’	линейная
‘cubic’	кубическая

Принято, что аргументы X и Y изменяются монотонно; кроме того, для кубической интерполяции предполагается, что сетка (X, Y) равномерна.

Пример 3. Проведем интерполяцию функции peaks, используя мелкую сетку (рис. 3):

[X, Y] = meshgrid(-3:0.25:3);

Z = peaks(X, Y);

[XI, YI] = meshgrid(-3:0.125:3);

ZI = interp2(X, Y, Z, XI, YI);

mesh(X, Y, Z)

hold on

mesh(XI, YI, ZI+15)

Функция peaks осуществляет построение поверхности с нормальным законом распределения значений по координате Z.

Рис. 3.

Двумерная интерполяция на неравномерной сетке. Функция ZI = griddata(X, Y, Z, XI, YI) возвращает массив ZI, который определен на новой сетке (XI, YI) в результате интерполяции исходной функции Z, заданной на неравномерной сетке (X, Y).

Функция [XI, YI, ZI] = griddata(X, Y, Z, XI, YI) кроме массива ZI возвращает массивы XI, YI, упорядоченные по аналогии с функцией meshgrid.

Пример 4. Определим функцию на сетке, заданной 100 точками, выбранными случайно на отрезке [-2; 2]:

X = rand(100, 1) * 4 - 2;

Y = rand(100, 1) * 4 - 2;

Z = X.*exp(-X.^2 - Y.^2);

При этом функция X = rand(m, n) формирует массив размера , элементами которого являются случайные величины, распределенные по равномерному закону в интервале (0, 1).

Векторы X, Y, Z определяют 100 случайных точек на поверхности функции ZI, которую зададим на следующей равномерной сетке:

ti = -2:0.25:2;

[XI, YI] = meshgrid(ti, ti);

ZI = griddata(X, Y, Z, XI, YI);

Построим поверхность функции, полученной в результате интерполяции на неравномерной случайной сетке (рис. 4):

mesh(XI, YI, ZI)

hold on

plot3(X, Y, Z, 'r o')

Рис. 4.

2. Функции для решения задач аппроксимации

В Matlab для решения задач аппроксимации имеются следующие функции:

polyfit — аппроксимация данных полиномом;

interpft — аппроксимация периодической функции на основе быстрого преобразования Фурье.

Аппроксимация данных полиномом. Функция p = polyfit(x, y, n) находит коэффициенты полинома степени n, который аппроксимирует функцию y(x) методом наименьших квадратов. Результатом является строка p длины n +1, содержащая коэффициенты аппроксимирующего полинома.

Пример 1. Рассмотрим аппроксимацию функции :

x = (0:0.1:10)';

y = sin(x)+cos(x);

вычислим коэффициенты аппроксимирующего полинома степени 6:

p = polyfit(x, y, 6)

p = 0.0004 -0.0079 0.0323 0.2245 -1.6581 2.0673 0.8047

вычислим значения полинома в точках сетки:

f = polyval(p, x);

Построим графики функции и аппроксимирующего полинома (рис. 1):

plot(x, y, 'b o', x, f, 'r')

Рис. 1.

Рассмотрим аппроксимацию функции (рис. 2). В этом случае для качественной аппроксимации необходимо использовать полином 30 степени:

x = (0:0.1:10)';

y = sin(2*x)+cos(3*x);

p = polyfit(x, y, 30)

f = polyval(p, x);

plot(x, y, 'b o', x, f, 'r')

grid on

Рис. 2.

Аппроксимация периодической функции на основе быстрого преобразования Фурье. Функция yp = interpft(y, n) возвращает одномерный массив чисел, который является периодической функцией, определенной в n точках и аппроксимирующей одномерный массив y. Если length(x) = m, а интервал дискретности dx, то интервал дискретности для y определяется по формуле , причем n всегда превышает m.

Пример 2. Рассмотрим аппроксимацию функции , которая задана 11 точками на интервале [0; 10] (рис. 3):

x = 0:10;

y = sin(x);

xp = 0:0.25:10;

yp = interpft(y, 41);

plot(x, y, 'b o -', xp, yp, 'r')

grid on

Рис. 3.

8