Extension of cubic spline interpolation
バイキュービック補間 といくつかの1次元および2次元補間 の比較。 黒 と 赤 / 黄 / 緑 / 青 の点は、それぞれ補間点と隣接サンプルに対応しています。 地面からの高さは、それぞれの値に対応しています。 数学 において 、 双三次補間は 三次スプライン補間 (データセットに三次補間を適用する手法) の拡張であり、 2次元の 規則的なグリッド 上のデータ点 を補間する 。補間された面(つまり、画像ではなくカーネルの形状)は、 双一次補間 や 最近傍補間 によって得られる対応する面よりも 滑らかになる。双三次補間は、 ラグランジュ多項式 、 三次スプライン 、または三次畳み込みアルゴリズムのいずれかを使用して実行できる 。
画像処理 において、 画像の再サンプリング では、速度が問題にならない場合、 双線形補間や最近傍補間よりも双三次補間が選ばれることが多い。4 ピクセル (2×2)のみを考慮する双線形補間とは対照的に、双三次補間は16ピクセル(4×4)を考慮します。双三次補間で再サンプリングされた画像は、 選択されたb値とc値に応じて、
異なる補間 アーティファクトが生じる可能性がある。
計算 25個の単位正方形を継ぎ合わせた 正方形上の双三次補間。Matplotlib の 実装に従った双三次補間。色は関数の値を示します。黒い点は補間対象となる所定のデータの位置です。色サンプルが放射状対称ではないことに注意してください。 [ 0 , 4 ] × [ 0 , 4 ] {\displaystyle [0,4]\times [0,4]} 上記と同じデータセットでの 双線形補間。 [ 疑わしい – 議論 ] 表面の微分は正方形の境界上で連続ではありません。 上記と同じデータセットに対する 最近傍補間。 単位正方形の4つの頂点 、、、、 における 関数値 と導関数 、 が既知であるとする 。補間された曲面は次のように表される 。 f {\displaystyle f} f x {\displaystyle f_{x}} f y {\displaystyle f_{y}} f x y {\displaystyle f_{xy}} ( 0 , 0 ) {\displaystyle (0,0)} ( 1 , 0 ) {\displaystyle (1,0)} ( 0 , 1 ) {\displaystyle (0,1)} ( 1 , 1 ) {\displaystyle (1,1)} p ( x , y ) = ∑ i = 0 3 ∑ j = 0 3 a i j x i y j . {\displaystyle p(x,y)=\sum \limits _{i=0}^{3}\sum _{j=0}^{3}a_{ij}x^{i}y^{j}.}
補間問題は16個の係数を決定することから成ります 。 関数値と対応づけることで、以下の4つの方程式が得られます。 a i j {\displaystyle a_{ij}} p ( x , y ) {\displaystyle p(x,y)}
f ( 0 , 0 ) = p ( 0 , 0 ) = a 00 , {\displaystyle f(0,0)=p(0,0)=a_{00},} f ( 1 , 0 ) = p ( 1 , 0 ) = a 00 + a 10 + a 20 + a 30 , {\displaystyle 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 , {\displaystyle f(0,1)=p(0,1)=a_{00}+a_{01}+a_{02}+a_{03},} f ( 1 , 1 ) = p ( 1 , 1 ) = ∑ i = 0 3 ∑ j = 0 3 a i j . {\displaystyle f(1,1)=p(1,1)=\textstyle \sum \limits _{i=0}^{3}\sum \limits _{j=0}^{3}a_{ij}.} 同様に、および方向の 導関数に対する 8 つの方程式 : x {\displaystyle x} y {\displaystyle y}
f x ( 0 , 0 ) = p x ( 0 , 0 ) = a 10 , {\displaystyle f_{x}(0,0)=p_{x}(0,0)=a_{10},} f x ( 1 , 0 ) = p x ( 1 , 0 ) = a 10 + 2 a 20 + 3 a 30 , {\displaystyle 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 , {\displaystyle f_{x}(0,1)=p_{x}(0,1)=a_{10}+a_{11}+a_{12}+a_{13},} f x ( 1 , 1 ) = p x ( 1 , 1 ) = ∑ i = 1 3 ∑ j = 0 3 a i j i , {\displaystyle 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 , {\displaystyle 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 , {\displaystyle 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 + 2 a 02 + 3 a 03 , {\displaystyle f_{y}(0,1)=p_{y}(0,1)=a_{01}+2a_{02}+3a_{03},} f y ( 1 , 1 ) = p y ( 1 , 1 ) = ∑ i = 0 3 ∑ j = 1 3 a i j j . {\displaystyle f_{y}(1,1)=p_{y}(1,1)=\textstyle \sum \limits _{i=0}^{3}\sum \limits _{j=1}^{3}a_{ij}j.} 混合偏微分 に関する4つの方程式 : x y {\displaystyle xy}
f x y ( 0 , 0 ) = p x y ( 0 , 0 ) = a 11 , {\displaystyle f_{xy}(0,0)=p_{xy}(0,0)=a_{11},} f x y ( 1 , 0 ) = p x y ( 1 , 0 ) = a 11 + 2 a 21 + 3 a 31 , {\displaystyle f_{xy}(1,0)=p_{xy}(1,0)=a_{11}+2a_{21}+3a_{31},} f x y ( 0 , 1 ) = p x y ( 0 , 1 ) = a 11 + 2 a 12 + 3 a 13 , {\displaystyle f_{xy}(0,1)=p_{xy}(0,1)=a_{11}+2a_{12}+3a_{13},} f x y ( 1 , 1 ) = p x y ( 1 , 1 ) = ∑ i = 1 3 ∑ j = 1 3 a i j i j . {\displaystyle f_{xy}(1,1)=p_{xy}(1,1)=\textstyle \sum \limits _{i=1}^{3}\sum \limits _{j=1}^{3}a_{ij}ij.} 上記の表現では、次の ID が使用されています。 p x ( x , y ) = ∑ i = 1 3 ∑ j = 0 3 a i j i x i − 1 y j , {\displaystyle 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 ) = ∑ i = 0 3 ∑ j = 1 3 a i j x i j y j − 1 , {\displaystyle p_{y}(x,y)=\textstyle \sum \limits _{i=0}^{3}\sum \limits _{j=1}^{3}a_{ij}x^{i}jy^{j-1},} p x y ( x , y ) = ∑ i = 1 3 ∑ j = 1 3 a i j i x i − 1 j y j − 1 . {\displaystyle 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 ) {\displaystyle p(x,y)} [ 0 , 1 ] × [ 0 , 1 ] {\displaystyle [0,1]\times [0,1]}
未知のパラメータを ベクトルに
グループ化し 、 上記の連立方程式を線形方程式の行列に再定式化することができます 。 a i j {\displaystyle a_{ij}} α = [ 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 ] T {\displaystyle \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 = [ 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 x y ( 0 , 0 ) f x y ( 1 , 0 ) f x y ( 0 , 1 ) f x y ( 1 , 1 ) ] T , {\displaystyle 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 α = x {\displaystyle A\alpha =x}
行列を反転すると、より有用な線形方程式 が得られ 、 これを 迅速かつ簡単に計算できます。 A − 1 x = α {\displaystyle A^{-1}x=\alpha } A − 1 = [ 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 ] , {\displaystyle 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],} α {\displaystyle \alpha }
16個の係数に対しては別の簡潔な行列形式も存在する。 または 、 [ 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 x y ( 0 , 0 ) f x y ( 0 , 1 ) f x ( 1 , 0 ) f x ( 1 , 1 ) f x y ( 1 , 0 ) f x y ( 1 , 1 ) ] = [ 1 0 0 0 1 1 1 1 0 1 0 0 0 1 2 3 ] [ 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 ] [ 1 1 0 0 0 1 1 1 0 1 0 2 0 1 0 3 ] , {\displaystyle {\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}},} [ 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 ] = [ 1 0 0 0 0 0 1 0 − 3 3 − 2 − 1 2 − 2 1 1 ] [ 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 x y ( 0 , 0 ) f x y ( 0 , 1 ) f x ( 1 , 0 ) f x ( 1 , 1 ) f x y ( 1 , 0 ) f x y ( 1 , 1 ) ] [ 1 0 − 3 2 0 0 3 − 2 0 1 − 2 1 0 0 − 1 1 ] , {\displaystyle {\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 ) = [ 1 x x 2 x 3 ] [ 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 ] [ 1 y y 2 y 3 ] . {\displaystyle 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 , {\displaystyle p_{x},p_{y},} p x y {\displaystyle p_{xy}} p x ( x , y ) = ∑ i = 1 3 ∑ j = 0 3 a i j i x i − 1 y j Δ x , {\displaystyle 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 ) = ∑ i = 0 3 ∑ j = 1 3 a i j x i j y j − 1 Δ y , {\displaystyle 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 x y ( x , y ) = ∑ i = 1 3 ∑ j = 1 3 a i j i x i − 1 j y j − 1 Δ x Δ y , {\displaystyle 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}},} Δ x {\displaystyle \Delta x} x {\displaystyle x} ( x , y ) {\displaystyle (x,y)} Δ y {\displaystyle \Delta y} α {\displaystyle \alpha } x = [ f ( 0 , 0 ) f ( 1 , 0 ) f ( 0 , 1 ) f ( 1 , 1 ) Δ x f x ( 0 , 0 ) Δ x f x ( 1 , 0 ) Δ x f x ( 0 , 1 ) Δ x f x ( 1 , 1 ) Δ y f y ( 0 , 0 ) Δ y f y ( 1 , 0 ) Δ y f y ( 0 , 1 ) Δ y f y ( 1 , 1 ) Δ x Δ y f x y ( 0 , 0 ) Δ x Δ y f x y ( 1 , 0 ) Δ x Δ y f x y ( 0 , 1 ) Δ x Δ y f x y ( 1 , 1 ) ] T , {\displaystyle 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},} α = A − 1 x {\displaystyle \alpha =A^{-1}x} A {\displaystyle A} x ¯ = x − x 0 x 1 − x 0 , y ¯ = y − y 0 y 1 − y 0 {\displaystyle {\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 , {\displaystyle x_{0},x_{1},y_{0},} y 1 {\displaystyle y_{1}} x {\displaystyle x} y {\displaystyle y} ( x , y ) {\displaystyle (x,y)} p ( x , y ) = ∑ i = 0 3 ∑ j = 0 3 a i j x ¯ i y ¯ j . {\displaystyle p(x,y)=\sum \limits _{i=0}^{3}\sum _{j=0}^{3}a_{ij}{\overline {x}}^{i}{\overline {y}}^{j}.}
関数値から導関数を求める 導関数が不明な場合は、通常、有限差分 などを使用して、単位正方形の角に隣接する点の関数値から近似されます 。
この方法を用いて、単導関数、 またはのいずれかを求めるには 、適切な軸における2点間の傾きを求めます 。 例えば、 いずれかの点を求めるには、 対象点の左右の点についてを求め、それらの点を結ぶ直線の傾きを計算します。同様に、 についても同様に行います 。 f x {\displaystyle f_{x}} f y {\displaystyle f_{y}} f x {\displaystyle f_{x}} f ( x , y ) {\displaystyle f(x,y)} f y {\displaystyle f_{y}}
相互導関数 を求めるには 、両軸でそれぞれ1つずつ導関数を取ります。例えば、まず の手順を使って 目標点の上と下の点の導関数を 求め、次に それらの値(通常は それらの点の の値ではなく)に対して の手順を適用して、目標点の の値を得ることができます。(あるいは、逆の方法で、まず と を計算し 、それからそれらから を 求めることもできます 。どちらの方法も結果は同じです。) f x y {\displaystyle f_{xy}} f x {\displaystyle f_{x}} x {\displaystyle x} f y {\displaystyle f_{y}} f {\displaystyle f} f x y ( x , y ) {\displaystyle f_{xy}(x,y)} f y {\displaystyle f_{y}} f x {\displaystyle f_{x}}
データセットの端で周囲の点が欠落している場合、欠落した点はいくつかの方法で近似できます。単純で一般的な方法は、既存の点から目標点までの傾きがそれ以上変化しないと仮定し、これに基づいて欠落点の仮想的な値を計算することです。
双三次畳み込みアルゴリズム 畳み込みカーネル W ( x ) {\displaystyle W(x)} 畳み込みカーネル(2Dで拡張) 双三次スプライン補間では、各グリッドセルについて、上述の線形システムの解が必要です。同様の特性を持つ補間関数は 、両次元において以下のカーネルを
用いた 畳み込み を適用することで得られます。
ここで、 は通常 -0.5 または -0.75 に設定されます。 また、 すべての非ゼロ整数 に対して となることに 注意してください 。 W ( x ) = { ( a + 2 ) | x | 3 − ( a + 3 ) | x | 2 + 1 for | x | ≤ 1 , a | x | 3 − 5 a | x | 2 + 8 a | x | − 4 a for 1 < | x | < 2 , 0 otherwise , {\displaystyle 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 {\displaystyle a} W ( 0 ) = 1 {\displaystyle W(0)=1} W ( n ) = 0 {\displaystyle W(n)=0} n {\displaystyle n}
このアプローチは、TRW Systems GroupのRifmanによって提案され、ERTS(地球資源技術衛星、後の Landsat 1 )の画像データを対象としました。 [1] [2] [3] 当初、 は-1に固定されていましたが、後に同社のSimonによってパラメータ化され、-0.5、-0.75、-1.0が意味のある値になりました。 [4] しかし、これらの提案では導出プロセスや公式が十分に示されていなかったため、後にKeysが完全な形で再提案し、 の場合 、元の関数のサンプリング間隔に関して3次収束が達成されることを示しました。 [5] a {\displaystyle a} a = − 0.5 {\displaystyle a=-0.5}
畳み込みカーネルは 次のように導出される [6] 。 W ( x ) {\displaystyle W(x)} W ( x ) = { A 1 | x | 3 + B 1 | x | 2 + C 1 | x | + D 1 for | x | ≤ 1 A 2 | x | 3 + B 2 | x | 2 + C 2 | x | + D 2 for 1 < | x | < 2 0 otherwise {\displaystyle W(x)={\begin{cases}A_{1}|x|^{3}+B_{1}|x|^{2}+C_{1}|x|+D_{1}&{\text{for }}|x|\leq 1\\A_{2}|x|^{3}+B_{2}|x|^{2}+C_{2}|x|+D_{2}&{\text{for }}1<|x|<2\\0&{\text{otherwise}}\end{cases}}}
条件: 、 、また 1 次導関数: 、の場合 、次のようになります。 W ( 0 ) = 1 {\displaystyle W(0)=1} W ( 1 ) = W ( 2 ) = 0 {\displaystyle W(1)=W(2)=0} W ′ ( 0 ) = W ′ ( 2 ) = 0 {\displaystyle W'(0)=W'(2)=0} W ′ ( 1 ) = a {\displaystyle W'(1)=a}
W ( 0 ) = D 1 = 1 , {\displaystyle W(0)=D_{1}=1,} W ( 1 − ) = A 1 + B 1 + C 1 + D 1 = 0 , {\displaystyle W(1^{-})=A_{1}+B_{1}+C_{1}+D_{1}=0,} W ( 1 + ) = A 2 + B 2 + C 2 + D 2 = 0 , {\displaystyle W(1^{+})=A_{2}+B_{2}+C_{2}+D_{2}=0,} W ( 2 − ) = 8 A 2 + 4 B 2 + 2 C 2 + D 2 = 0 , {\displaystyle W(2^{-})=8A_{2}+4B_{2}+2C_{2}+D_{2}=0,} W ′ ( 0 ) = C 1 = 0 , {\displaystyle W'(0)=C_{1}=0,} W ′ ( 1 − ) = 3 A 1 + 2 B 1 + C 1 = a , {\displaystyle W'(1^{-})=3A_{1}+2B_{1}+C_{1}=a,} W ′ ( 1 + ) = 3 A 2 + 2 B 2 + C 2 = a , {\displaystyle W'(1^{+})=3A_{2}+2B_{2}+C_{2}=a,} W ′ ( 2 − ) = 12 A 2 + 4 B 2 + C 2 = 0. {\displaystyle W'(2^{-})=12A_{2}+4B_{2}+C_{2}=0.} 上記の連立方程式を解くと、 次の式が得られます。 A 1 = a + 2 , B 1 = − ( a + 3 ) , C 1 = 0 , D 1 = 1 , {\displaystyle A_{1}=a+2,\quad B_{1}=-(a+3),\quad C_{1}=0,\quad D_{1}=1,} A 2 = a , B 2 = − 5 a , C 2 = 8 a , D 2 = − 4 a , {\displaystyle A_{2}=a,\quad B_{2}=-5a,\quad C_{2}=8a,\quad D_{2}=-4a,}
行列表記 行列表記のために畳み込みカーネルを0から1にシフト 行列表記を用いると、この式をより分かりやすく表現できます。1
次元 の場合、 0から1までの範囲で表すことができます。1次元3次畳み込み補間では、4つのサンプル点が必要であることに注意してください。各照会点に対し、左側に2つのサンプル点、右側に2つのサンプル点が配置されます。これらの点には、このテキストでは-1から2までのインデックスが付けられています。0でインデックス付けされた点から照会点までの距離は、 ここで で表されます。 p ( t ) = 1 2 [ 1 t t 2 t 3 ] [ 0 1 0 0 a 0 − a 0 − 2 a − a − 3 2 a + 3 a a a + 2 − a − 2 − a ] [ f − 1 f 0 f 1 f 2 ] {\displaystyle p(t)={\tfrac {1}{2}}{\begin{bmatrix}1&t&t^{2}&t^{3}\end{bmatrix}}{\begin{bmatrix}0&1&0&0\\a&0&-a&0\\-2a&-a-3&2a+3&a\\a&a+2&-a-2&-a\end{bmatrix}}{\begin{bmatrix}f_{-1}\\f_{0}\\f_{1}\\f_{2}\end{bmatrix}}} t {\displaystyle t} t {\displaystyle t}
であるため、これは Catmull-Romスプライン と一致する 。したがって、 を とする双三次補間は 「Catmull-Rom補間」と呼ばれることもある。 [7] a = − 0.5 {\displaystyle a=-0.5} a {\displaystyle a} − 0.5 {\displaystyle -0.5}
2次元の場合、まず に1回適用し 、次に に再度適用します 。 x {\displaystyle x} y {\displaystyle y} b − 1 = p ( t x , f ( − 1 , − 1 ) , f ( 0 , − 1 ) , f ( 1 , − 1 ) , f ( 2 , − 1 ) ) , b 0 = p ( t x , f ( − 1 , 0 ) , f ( 0 , 0 ) , f ( 1 , 0 ) , f ( 2 , 0 ) ) , b 1 = p ( t x , f ( − 1 , 1 ) , f ( 0 , 1 ) , f ( 1 , 1 ) , f ( 2 , 1 ) ) , b 2 = p ( t x , f ( − 1 , 2 ) , f ( 0 , 2 ) , f ( 1 , 2 ) , f ( 2 , 2 ) ) , {\displaystyle {\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 ) . {\displaystyle p(x,y)=p(t_{y},b_{-1},b_{0},b_{1},b_{2}).}
微分連続性 定義により 1 次導関数で連続しており、次の行列表記を使用して簡単に検証できます。
p i ′ ( 1 ) = [ 0 a 0 − a ] [ f i − 1 f i f i + 1 f i + 2 ] = p i + 1 ′ ( 0 ) = [ a 0 − a 0 ] [ f i f i + 1 f i + 2 f i + 3 ] {\displaystyle p_{i}'(1)={\begin{bmatrix}0&a&0&-a\end{bmatrix}}{\begin{bmatrix}f_{i-1}\\f_{i}\\f_{i+1}\\f_{i+2}\end{bmatrix}}=p_{i+1}'(0)={\begin{bmatrix}a&0&-a&0\end{bmatrix}}{\begin{bmatrix}f_{i}\\f_{i+1}\\f_{i+2}\\f_{i+3}\end{bmatrix}}} したがって、1 次導関数では連続です。
p i ″ ( 1 ) = [ 2 a 4 a + 6 − 2 a − 6 − 4 a ] [ f i − 1 f i f i + 1 f i + 2 ] ≠ p i + 1 ″ ( 0 ) = [ − 4 a − 2 a − 6 4 a + 6 2 a ] [ f i f i + 1 f i + 2 f i + 3 ] {\displaystyle p_{i}''(1)={\begin{bmatrix}2a&4a+6&-2a-6&-4a\end{bmatrix}}{\begin{bmatrix}f_{i-1}\\f_{i}\\f_{i+1}\\f_{i+2}\end{bmatrix}}\neq p_{i+1}''(0)={\begin{bmatrix}-4a&-2a-6&4a+6&2a\end{bmatrix}}{\begin{bmatrix}f_{i}\\f_{i+1}\\f_{i+2}\\f_{i+3}\end{bmatrix}}} したがって、2 次導関数では連続ではありません。
他の方法との比較 以下に、各方法のカーネル関数の比較を示します。
方法 カーネル関数 カーネル関数(2D拡張) 最も近い隣人 双線形補間 双三次補間
コンピュータグラフィックスでの使用 バイキュービックアルゴリズムは、画像や動画を表示用に拡大縮小する際によく使用されます( ビットマップの再サンプリングを参照)。バイキュービックアルゴリズムは、一般的な バイリニア アルゴリズムよりも優れたディテールを保持します 。
バイキュービックアルゴリズムは画像のダウンスケーリングにも使用されますが、ダウンスケーリング中の再サンプリングは実際の光学におけるズームアウトの効果とは異なるため、他の方法の方が適している場合があります。
パラメータ「a」の効果 この図の下半分は上半分を拡大したもので、左側の線のシャープネスがどのように実現されているかを示しています。バイキュービック補間によりオーバーシュートが発生し、 アキュータンス が増大します。 しかし、カーネルの負のローブにより、 オーバーシュート (ハロー)が発生します。これは クリッピングを 引き起こす可能性があり、アーティファクト( リンギングアーティファクトも参照)となりますが、 アキュータンス (見かけのシャープネス)を高めるため、望ましい効果となることもあります。この効果は、パラメータ が または の場合の方が の場合よりも 大きくなります 。 a {\displaystyle a} − 0.75 {\displaystyle -0.75} − 1.0 {\displaystyle -1.0} − 0.5 {\displaystyle -0.5}
のとき 、補間結果は元の画像(サンプリング前の連続関数)の2次導関数(剰余項を除く)までテイラー近似と一致します。 [8] [9] この効果の明確な例は、元の画像が線形関数(単純なグラデーション画像)で、補間結果が の場合にのみ線形関数になる場合です 。たとえば、入力が のとき、 補間結果は であり 、 の場合にのみ線形関数(2次項と3次項が0)になります 。 [10] a = − 0.5 {\displaystyle a=-0.5} a = − 0.5 {\displaystyle a=-0.5} [ f i − 1 f i f i + 1 f i + 2 ] = [ − 1 0 1 2 ] , {\displaystyle {\begin{bmatrix}f_{i-1}&f_{i}&f_{i+1}&f_{i+2}\end{bmatrix}}={\begin{bmatrix}-1&0&1&2\end{bmatrix}},} p ( t ) = − 2 ( 2 a + 1 ) t 3 + 3 ( 2 a + 1 ) t 2 − 2 a t , {\displaystyle p(t)=-2(2a+1)t^{3}+3(2a+1)t^{2}-2at,} a = − 0.5 {\displaystyle a=-0.5}
a = − 0.5 {\displaystyle a=-0.5} a = − 0.75 {\displaystyle a=-0.75} a = − 1.0 {\displaystyle a=-1.0}
下図は、単純なグラデーション画像を16倍に拡大したものと、標準テスト画像「カメラマン」 [11] の一部を 8倍に拡大したものを示している。前者の縦縞と後者のこめかみ周辺の人工的なまだら模様は、前述の効果によって生じたアーティファクトである。
a = − 0.5 {\displaystyle a=-0.5} a = − 0.75 {\displaystyle a=-0.75} a = − 1.0 {\displaystyle a=-1.0} 最近傍補間 (比較用)
のとき 、畳み込みカーネルの2次導関数は [10] で連続するが、補間結果は連続する2次導関数を持たないことに注意されたい。また、のとき、畳み込みカーネルの導関数は [10] における sinc関数 の導関数と一致するが、これは sinc補間 との比較以外には意味を持たない 。 a = − 0.75 {\displaystyle a=-0.75} x = 1 {\displaystyle x=1} a = − 1.0 {\displaystyle a=-1.0} x = 1 {\displaystyle x=1}
以下は、一般的に使用される実装で の設定の例です。 a {\displaystyle a}
^ Java 高度なイメージング - クラス InterpolationBicubic ^ ImageMagick v6 の例 - 再サンプリングフィルタ - キュービックフィルタファミリー ^ Stack Overflow - imresize - 双三次補間を理解しようとする ^ ab PyTorch ドキュメント > torch.nn.functional > torch.nn.functional.grid_sample ^ ab tensorflow/tensorflow/core/kernels/image/resize_bicubic_op.cc ^ jax-ml/jax - "bicubic" を使用する場合の jax.image.resize と F.interpolate の違いに注意してください #15768 ^ opencv/opencv - INTER_CUBIC は Catmull-Rom 補間をデフォルトにすべき #17720 ^ Java 高度なイメージング - クラス InterpolationBicubic2
再サンプリング位置 画像を整数倍に拡大する場合、拡大画像の重心が入力画像に対してずれないように、リサンプリング位置を入力画像のピクセル間に設定するのが一般的です(下図左参照)。 一方、リサンプリング位置を入力画像のピクセルに重なるように設定する方法もあります(下図右参照)。この方法は、変換の可逆性を確保します。
再サンプル位置の例
補間には、入力画像の左右(または上下)の2ピクセルが必要です。画像の外側の端は無効であるため、外挿が必要です。外挿の方法は実装によって異なります。 以下は、最外側のピクセルをコピーして外挿する例です。
最も外側のピクセル値をコピーして外挿する例
参照
参考文献 ^ Rifman, Samuel S. (1973). 「ERTSマルチスペクトル画像のデジタル補正」 NASA技術報告書73N28327 . ^ Bauer, Brian P. (1980). 「MSSと合成画像を用いたリサンプリング技術の調査」 (PDF) . EDC文書0044, USGS EROSデータセンター : 12. ^ ジョージ、ウォルバーグ (1994)。 デジタル画像ワーピング、第 3 版 。 130 ~ 131 ページ 。ISBN 0-8186-8944-7 。 ^ Simon, KW (1975). 「幾何学的操作のためのデジタル画像再構成と再サンプリング」 LARSシンポジウム論文67 . ^ R. Keys (1981). 「デジタル画像処理のための三次畳み込み補間」. IEEE Transactions on Acoustics, Speech, and Signal Processing . 29 (6): 1153– 1160. Bibcode :1981ITASS..29.1153K. CiteSeerX 10.1.1.320.776 . doi :10.1109/TASSP.1981.1163711. ^ キーズ(1981)、pp.1153-1154 ^ 「Catmull-Rom補間」という用語が使用されている例: INTER_CUBICはデフォルトでCatmull-Rom補間に設定される必要があります。#17720 ^ キーズ(1981)、pp.1154-1155 ^ ウォルバーグ(1994)、130ページ ^ abc サイモン(1975)、p.3A-3 ^ MIT Works がデジタル化される - カメラマン
外部リンク 標高サンプルへの補間の適用 補間理論 (双)三次補間の説明とJava/C++実装 双三次ラグランジュ補間の Excel ワークシート関数