Бикубическая интерполяция

В математике бикубическая интерполяция — это расширение интерполяции кубическим сплайном (метода применения кубической интерполяции к набору данных) для интерполяции точек данных на двумерной регулярной сетке . Интерполированная поверхность (имеется в виду форма ядра, а не изображение) более гладкая , чем соответствующие поверхности, полученные с помощью билинейной интерполяции или интерполяции ближайших соседей . Бикубическая интерполяция может быть выполнена с использованием полиномов Лагранжа , кубических сплайнов или кубической свертки алгоритма .

При обработке изображений бикубическая интерполяция часто выбирается вместо билинейной интерполяции или интерполяции по ближайшему соседу при повторной выборке изображения , когда скорость не является проблемой. В отличие от билинейной интерполяции, которая учитывает только 4 пикселя (2×2), бикубическая интерполяция учитывает 16 пикселей (4×4). Изображения, подвергнутые повторной выборке с помощью бикубической интерполяции, могут иметь различные артефакты интерполяции в зависимости от выбранных значений b и c.

Вычисление

Бикубическая интерполяция на квадрате $[0,4]\times [0,4]$ состоящий из 25 отдельных квадратов, соединенных вместе. Бикубическая интерполяция согласно . реализации Matplotlib Цвет указывает значение функции. Черные точки — это места интерполяции предписанных данных. Обратите внимание, что образцы цвета не являются радиально симметричными.

Билинейная интерполяция на том же наборе данных, что и выше. Производные поверхности не непрерывны по границам квадрата.

Интерполяция ближайшего соседа для того же набора данных, что и выше.

Предположим, что значения функции $f$ и производные $f_{x}$ , $f_{y}$ и $f_{xy}$ известны на четырех углах $(0,0)$ , $(1,0)$ , $(0,1)$ , и $(1,1)$ единицы площади. Тогда интерполированную поверхность можно записать как $p(x,y)=\sum \limits _{i=0}^{3}\sum _{j=0}^{3}a_{ij}x^{i}y^{j}.$

Задача интерполяции состоит в определении 16 коэффициентов $a_{ij}$ . Соответствие $p(x,y)$ со значениями функции дает четыре уравнения:

$f(0,0)=p(0,0)=a_{00},$
$f(1,0)=p(1,0)=a_{00}+a_{10}+a_{20}+a_{30},$
$f(0,1)=p(0,1)=a_{00}+a_{01}+a_{02}+a_{03},$
$f(1,1)=p(1,1)=\textstyle \sum \limits _{i=0}^{3}\sum \limits _{j=0}^{3}a_{ij}.$

Аналогично, восемь уравнений для производных в $x$ и $y$ направления:

$f_{x}(0,0)=p_{x}(0,0)=a_{10},$
$f_{x}(1,0)=p_{x}(1,0)=a_{10}+2a_{20}+3a_{30},$
$f_{x}(0,1)=p_{x}(0,1)=a_{10}+a_{11}+a_{12}+a_{13},$
$f_{x}(1,1)=p_{x}(1,1)=\textstyle \sum \limits _{i=1}^{3}\sum \limits _{j=0}^{3}a_{ij}i,$
$f_{y}(0,0)=p_{y}(0,0)=a_{01},$
$f_{y}(1,0)=p_{y}(1,0)=a_{01}+a_{11}+a_{21}+a_{31},$
$f_{y}(0,1)=p_{y}(0,1)=a_{01}+2a_{02}+3a_{03},$
$f_{y}(1,1)=p_{y}(1,1)=\textstyle \sum \limits _{i=0}^{3}\sum \limits _{j=1}^{3}a_{ij}j.$

И четыре уравнения для $xy$ смешанная частная производная :

$f_{xy}(0,0)=p_{xy}(0,0)=a_{11},$
$f_{xy}(1,0)=p_{xy}(1,0)=a_{11}+2a_{21}+3a_{31},$
$f_{xy}(0,1)=p_{xy}(0,1)=a_{11}+2a_{12}+3a_{13},$
$f_{xy}(1,1)=p_{xy}(1,1)=\textstyle \sum \limits _{i=1}^{3}\sum \limits _{j=1}^{3}a_{ij}ij.$

В приведенных выше выражениях использовались следующие тождества: $p_{x}(x,y)=\textstyle \sum \limits _{i=1}^{3}\sum \limits _{j=0}^{3}a_{ij}ix^{i-1}y^{j},$ $p_{y}(x,y)=\textstyle \sum \limits _{i=0}^{3}\sum \limits _{j=1}^{3}a_{ij}x^{i}jy^{j-1},$ $p_{xy}(x,y)=\textstyle \sum \limits _{i=1}^{3}\sum \limits _{j=1}^{3}a_{ij}ix^{i-1}jy^{j-1}.$

Эта процедура дает поверхность $p(x,y)$ на единичной площади $[0,1]\times [0,1]$ которая непрерывна и имеет непрерывные производные. Затем можно выполнить бикубическую интерполяцию на регулярной сетке произвольного размера путем объединения таких бикубических поверхностей, гарантируя совпадение производных на границах.

Группировка неизвестных параметров $a_{ij}$ в векторе $\alpha =\left[{\begin{smallmatrix}a_{00}&a_{10}&a_{20}&a_{30}&a_{01}&a_{11}&a_{21}&a_{31}&a_{02}&a_{12}&a_{22}&a_{32}&a_{03}&a_{13}&a_{23}&a_{33}\end{smallmatrix}}\right]^{T}$ и позволяя $x=\left[{\begin{smallmatrix}f(0,0)&f(1,0)&f(0,1)&f(1,1)&f_{x}(0,0)&f_{x}(1,0)&f_{x}(0,1)&f_{x}(1,1)&f_{y}(0,0)&f_{y}(1,0)&f_{y}(0,1)&f_{y}(1,1)&f_{xy}(0,0)&f_{xy}(1,0)&f_{xy}(0,1)&f_{xy}(1,1)\end{smallmatrix}}\right]^{T},$ приведенную выше систему уравнений можно переформулировать в матрицу линейного уравнения $A\alpha =x$ .

Обращение матрицы дает более полезное линейное уравнение $A^{-1}x=\alpha$ , где $A^{-1}=\left[{\begin{smallmatrix}{\begin{array}{rrrrrrrrrrrrrrrr}1&0&0&0&0&0&0&0&0&0&0&0&0&0&0&0\\0&0&0&0&1&0&0&0&0&0&0&0&0&0&0&0\\-3&3&0&0&-2&-1&0&0&0&0&0&0&0&0&0&0\\2&-2&0&0&1&1&0&0&0&0&0&0&0&0&0&0\\0&0&0&0&0&0&0&0&1&0&0&0&0&0&0&0\\0&0&0&0&0&0&0&0&0&0&0&0&1&0&0&0\\0&0&0&0&0&0&0&0&-3&3&0&0&-2&-1&0&0\\0&0&0&0&0&0&0&0&2&-2&0&0&1&1&0&0\\-3&0&3&0&0&0&0&0&-2&0&-1&0&0&0&0&0\\0&0&0&0&-3&0&3&0&0&0&0&0&-2&0&-1&0\\9&-9&-9&9&6&3&-6&-3&6&-6&3&-3&4&2&2&1\\-6&6&6&-6&-3&-3&3&3&-4&4&-2&2&-2&-2&-1&-1\\2&0&-2&0&0&0&0&0&1&0&1&0&0&0&0&0\\0&0&0&0&2&0&-2&0&0&0&0&0&1&0&1&0\\-6&6&6&-6&-4&-2&4&2&-3&3&-3&3&-2&-1&-2&-1\\4&-4&-4&4&2&2&-2&-2&2&-2&2&-2&1&1&1&1\end{array}}\end{smallmatrix}}\right],$ что позволяет $\alpha$ рассчитать быстро и легко.

Возможна и другая краткая матричная форма для 16 коэффициентов: ${\begin{bmatrix}f(0,0)&f(0,1)&f_{y}(0,0)&f_{y}(0,1)\\f(1,0)&f(1,1)&f_{y}(1,0)&f_{y}(1,1)\\f_{x}(0,0)&f_{x}(0,1)&f_{xy}(0,0)&f_{xy}(0,1)\\f_{x}(1,0)&f_{x}(1,1)&f_{xy}(1,0)&f_{xy}(1,1)\end{bmatrix}}={\begin{bmatrix}1&0&0&0\\1&1&1&1\\0&1&0&0\\0&1&2&3\end{bmatrix}}{\begin{bmatrix}a_{00}&a_{01}&a_{02}&a_{03}\\a_{10}&a_{11}&a_{12}&a_{13}\\a_{20}&a_{21}&a_{22}&a_{23}\\a_{30}&a_{31}&a_{32}&a_{33}\end{bmatrix}}{\begin{bmatrix}1&1&0&0\\0&1&1&1\\0&1&0&2\\0&1&0&3\end{bmatrix}},$ или ${\begin{bmatrix}a_{00}&a_{01}&a_{02}&a_{03}\\a_{10}&a_{11}&a_{12}&a_{13}\\a_{20}&a_{21}&a_{22}&a_{23}\\a_{30}&a_{31}&a_{32}&a_{33}\end{bmatrix}}={\begin{bmatrix}1&0&0&0\\0&0&1&0\\-3&3&-2&-1\\2&-2&1&1\end{bmatrix}}{\begin{bmatrix}f(0,0)&f(0,1)&f_{y}(0,0)&f_{y}(0,1)\\f(1,0)&f(1,1)&f_{y}(1,0)&f_{y}(1,1)\\f_{x}(0,0)&f_{x}(0,1)&f_{xy}(0,0)&f_{xy}(0,1)\\f_{x}(1,0)&f_{x}(1,1)&f_{xy}(1,0)&f_{xy}(1,1)\end{bmatrix}}{\begin{bmatrix}1&0&-3&2\\0&0&3&-2\\0&1&-2&1\\0&0&-1&1\end{bmatrix}},$ где $p(x,y)={\begin{bmatrix}1&x&x^{2}&x^{3}\end{bmatrix}}{\begin{bmatrix}a_{00}&a_{01}&a_{02}&a_{03}\\a_{10}&a_{11}&a_{12}&a_{13}\\a_{20}&a_{21}&a_{22}&a_{23}\\a_{30}&a_{31}&a_{32}&a_{33}\end{bmatrix}}{\begin{bmatrix}1\\y\\y^{2}\\y^{3}\end{bmatrix}}.$

Расширение на прямолинейные сетки.

Часто приложения требуют бикубической интерполяции с использованием данных в прямолинейной сетке, а не в единичном квадрате. В этом случае тождества для $p_{x},p_{y},$ и $p_{xy}$ становиться $p_{x}(x,y)=\textstyle \sum \limits _{i=1}^{3}\sum \limits _{j=0}^{3}{\frac {a_{ij}ix^{i-1}y^{j}}{\Delta x}},$ $p_{y}(x,y)=\textstyle \sum \limits _{i=0}^{3}\sum \limits _{j=1}^{3}{\frac {a_{ij}x^{i}jy^{j-1}}{\Delta y}},$ $p_{xy}(x,y)=\textstyle \sum \limits _{i=1}^{3}\sum \limits _{j=1}^{3}{\frac {a_{ij}ix^{i-1}jy^{j-1}}{\Delta x\Delta y}},$ где $\Delta x$ это $x$ расстояние между ячейкой, содержащей точку $(x,y)$ и подобное для $\Delta y$ . В этом случае наиболее практичный подход к вычислению коэффициентов $\alpha$ это позволить $x=\left[{\begin{smallmatrix}f(0,0)&f(1,0)&f(0,1)&f(1,1)&\Delta xf_{x}(0,0)&\Delta xf_{x}(1,0)&\Delta xf_{x}(0,1)&\Delta xf_{x}(1,1)&\Delta yf_{y}(0,0)&\Delta yf_{y}(1,0)&\Delta yf_{y}(0,1)&\Delta yf_{y}(1,1)&\Delta x\Delta yf_{xy}(0,0)&\Delta x\Delta yf_{xy}(1,0)&\Delta x\Delta yf_{xy}(0,1)&\Delta x\Delta yf_{xy}(1,1)\end{smallmatrix}}\right]^{T},$ затем решить $\alpha =A^{-1}x$ с $A$ как раньше. Далее нормализованные интерполяционные переменные вычисляются как ${\begin{aligned}{\overline {x}}&={\frac {x-x_{0}}{x_{1}-x_{0}}},\\{\overline {y}}&={\frac {y-y_{0}}{y_{1}-y_{0}}}\end{aligned}}$ где $x_{0},x_{1},y_{0},$ и $y_{1}$ являются $x$ и $y$ координаты точек сетки, окружающих точку $(x,y)$ . Тогда интерполирующая поверхность становится $p(x,y)=\sum \limits _{i=0}^{3}\sum _{j=0}^{3}a_{ij}{\overline {x}}^{i}{\overline {y}}^{j}.$

Нахождение производных от значений функции

Если производные неизвестны, они обычно аппроксимируются по значениям функции в точках, соседних с углами единичного квадрата, например, с использованием конечных разностей .

Чтобы найти любую из отдельных производных, $f_{x}$ или $f_{y}$ , используя этот метод, найдите наклон между двумя окружающими точками на соответствующей оси. Например, чтобы вычислить $f_{x}$ для одной из точек найдите $f(x,y)$ для точек слева и справа от целевой точки и вычислите их наклон, и аналогично для $f_{y}$ .

Чтобы найти перекрестную производную $f_{xy}$ , возьмите производную по обеим осям, по одной. Например, можно сначала использовать $f_{x}$ процедура поиска $x$ производные точек выше и ниже целевой точки, затем используйте $f_{y}$ процедуру над этими значениями (а не, как обычно, над значениями $f$ для этих точек), чтобы получить значение $f_{xy}(x,y)$ для целевой точки. (Или можно сделать это в обратном направлении, предварительно вычислив $f_{y}$ а потом $f_{x}$ из тех. Оба дают эквивалентные результаты.)

На краях набора данных, когда отсутствуют некоторые окружающие точки, недостающие точки можно аппроксимировать несколькими методами. Простой и распространенный метод состоит в том, чтобы предположить, что наклон от существующей точки до целевой точки продолжается без дальнейших изменений, и использовать это для расчета гипотетического значения для недостающей точки.

Алгоритм бикубической свертки

Интерполяция бикубическим сплайном требует решения описанной выше линейной системы для каждой ячейки сетки. Интерполятор с аналогичными свойствами можно получить, применив свертку со следующим ядром в обоих измерениях: $W(x)={\begin{cases}(a+2)|x|^{3}-(a+3)|x|^{2}+1&{\text{for }}|x|\leq 1,\\a|x|^{3}-5a|x|^{2}+8a|x|-4a&{\text{for }}1<|x|<2,\\0&{\text{otherwise}},\end{cases}}$ где $a$ обычно устанавливается равным -0,5 или -0,75. Обратите внимание, что $W(0)=1$ и $W(n)=0$ для всех ненулевых целых чисел $n$ .

Этот подход был предложен Кизом, который показал, что $a=-0.5$ производит сходимость третьего порядка относительно интервала дискретизации исходной функции. ^[1]

Если мы воспользуемся матричным обозначением для общего случая $a=-0.5$ , мы можем выразить уравнение более удобным способом: $p(t)={\tfrac {1}{2}}{\begin{bmatrix}1&t&t^{2}&t^{3}\end{bmatrix}}{\begin{bmatrix}0&2&0&0\\-1&0&1&0\\2&-5&4&-1\\-1&3&-3&1\end{bmatrix}}{\begin{bmatrix}f_{-1}\\f_{0}\\f_{1}\\f_{2}\end{bmatrix}}$ для $t$ от 0 до 1 для одного измерения. Обратите внимание, что для одномерной интерполяции кубической свертки требуются 4 точки выборки. Для каждого запроса слева расположены два образца, а справа — два образца. В этом тексте эти точки пронумерованы от −1 до 2. Расстояние от точки с индексом 0 до точки запроса обозначается $t$ здесь.

Для двух измерений, впервые примененных один раз в $x$ и снова в $y$ : ${\begin{aligned}b_{-1}&=p(t_{x},f_{(-1,-1)},f_{(0,-1)},f_{(1,-1)},f_{(2,-1)}),\\[1ex]b_{0}&=p(t_{x},f_{(-1,0)},f_{(0,0)},f_{(1,0)},f_{(2,0)}),\\[1ex]b_{1}&=p(t_{x},f_{(-1,1)},f_{(0,1)},f_{(1,1)},f_{(2,1)}),\\[1ex]b_{2}&=p(t_{x},f_{(-1,2)},f_{(0,2)},f_{(1,2)},f_{(2,2)}),\end{aligned}}$ $p(x,y)=p(t_{y},b_{-1},b_{0},b_{1},b_{2}).$

Использование в компьютерной графике

Нижняя половина этого рисунка представляет собой увеличение верхней половины, показывающее, как создается кажущаяся резкость левой линии. Бикубическая интерполяция вызывает перерегулирование, что увеличивает резкость .

Бикубический алгоритм часто используется для масштабирования изображений и видео для отображения (см. передискретизацию растрового изображения ). Он сохраняет мелкие детали лучше, чем обычный билинейный алгоритм.

Однако из-за отрицательных лепестков ядра это вызывает перерегулирование (ореолы). Это может вызвать отсечение и является артефактом (см. также артефакты звона ), но увеличивает остроту (кажущуюся резкость) и может быть желательным.

См. также

Пространственное сглаживание
Поверхность Безье
Билинейная интерполяция
Кубический сплайн Эрмита , одномерный аналог бикубического сплайна.
Передискретизация Ланцоша
Интерполяция естественных соседей
Синк-фильтр
Сплайн-интерполяция
Трикубическая интерполяция
Интерполяция направленной кубической свертки

Ссылки

^ Р. Киз (1981). «Интерполяция кубической свертки для цифровой обработки изображений». Транзакции IEEE по акустике, речи и обработке сигналов . 29 (6): 1153–1160. Бибкод : 1981ITASS..29.1153K . CiteSeerX 10.1.1.320.776 . дои : 10.1109/ТАССП.1981.1163711 .

Внешние ссылки

[Keys-1] Р. Киз (1981). «Интерполяция кубической свертки для цифровой обработки изображений». Транзакции IEEE по акустике, речи и обработке сигналов . 29 (6): 1153–1160. Бибкод : 1981ITASS..29.1153K . CiteSeerX 10.1.1.320.776 . дои : 10.1109/ТАССП.1981.1163711 .

[1]