離散コサイン変換

二次元DCTとDFTとの比較。左はスペクトル、右はヒストグラム。低周波域での相違を示すため、スペクトルは 1/4 だけ示してある。DCTでは、パワーのほとんどが低周波領域に集中していることがわかる。

離散コサイン変換(りさんコサインへんかん、: discrete cosine transformDCT)は、離散信号周波数領域へ変換する方法の一つである。

概要[編集]

DCTは、有限数列を、余弦関数数列 cos(nk) を基底とする一次結合(つまり、適切な周波数振幅のコサインカーブの和)の係数に変換する。余弦関数は実数に対しては実数を返すので、実数列に対してはDCT係数も実数列となる。

これは、離散フーリエ変換 (DFT: discrete Fourier transform) が、実数に対しても複素数を返す exp(ink) を使うため、実数列に対しても複素数列となるのと大きな違いである。なお、DFTも偶関数数列に対しては実係数を返す、つまりコサイン成分のみとなるが、DCTはy軸で折り返して偶関数化してDFTすることと等価であり、実際にそう計算することが多い。

DCTでは、係数が実数になる上、特定の成分への集中度があがる。JPEGなどの画像圧縮、AACMP3ATRACといった音声圧縮、デジタルフィルタ等広い範囲で用いられている。

逆変換を逆離散コサイン変換: inverse discrete cosine transform (IDCT))と呼ぶ。

種類[編集]

DCTには標準的な方法が8通りあり、そのうち4つがふつうに用いられる。最も一般的な方法は type-II DCT であり、単にDCTと呼んだ場合これを指すことが多い(以下DCT-II)。同様に、DCT-IIの逆変換である type-III DCT は単に逆DCT (inverse DCT) ないしIDCTと呼ばれることが多い。

DCTに関連する変換法が二つある。一つは離散サイン変換 (DST) であり、実領域で奇関数を用いた離散フーリエ変換 (DFT) と等価である。もう一つの修正離散コサイン変換 (MDCT) は「互いに重なりのある」データのDCTに基づいている。

応用[編集]

DCT、特にDCT-IIは信号・画像処理にしばしば用いられる。特に不可逆圧縮に頻用されるが、これはDCTの持つ強力な「エネルギー圧縮」特性による。DCTで変換すると、情報が少数の低周波成分に集中する傾向が生まれ、マルコフ過程の制限に基づく信号について、非相関成分の検出に最適であるKarhunen-Loève変換に近い。

たとえばDCTはJPEGMJPEGMPEGDV等の画像圧縮に用いられる。これらの画像圧縮では、N × N のブロックに2次元DCT-IIを行い、その結果を量子化し、エントロピー圧縮する。典型的には N = 8 であり、そのブロックの行ごと、列ごとにDCT-IIの式を適用する。その結果得られる 8 × 8 行列は、要素 (0, 0) をDC(直流。周波数が 0)成分とし、行・列とも、添字が大きくなるほど垂直ないし水平方向の空間周波数が高い成分を表す。

2次元DCTとDFTとの比較。

音声圧縮に用いられるMDCTAACVorbisMP3も関連した変換法である。

DCTは、偏微分方程式をスペクトル法で解くときにも広く使われる。その場合、配列の両端での境界条件の偶奇性に対応して、異なるDCTの変種が使われる。

DCTはまた、チェビシェフ多項式とも密接に関係しており、高速DCT算法(下記)はクレーンショー・カーチス数値積分則のような、任意の関数についてのチェビシェフ近似に用いられる。

非形式的概説[編集]

フーリエ変換やそれに類似の変換(以下、類フーリエ変換とよぶ)のように、離散コサイン変換 (DCT) も関数あるいは信号を異なる周波数振幅をもつ三角関数の和として表現する。また、DCTは、離散フーリエ変換 (DFT) と同じく、離散的なデータ点からなる有限の関数に対して施される。一見してそれとわかるDCTとDFTとの違いはDCTがコサイン(余弦)関数のみを使うのに対してDFTがコサインとサイン(正弦)関数の両方を(複素数の指数関数の形式で)使うという点である。しかし、この見かけの違いはもっと本質的な違いの帰結でしかない。すなわち、DCTとDFTあるいは他の関連する変換は境界条件において異なっているということである。

有限の定義域をもつ関数に施される類フーリエ変換、すなわちDFTやDCTやフーリエ級数は、暗黙のうちにその定義域の外部に関数を「拡張」して定義しているのだと考えることができる。つまり、ある関数 f(x) を一旦三角関数の和として表現してしまうと、任意の x に対し、それがたとえ元の関数 f(x) が定義されていない x であったとしても、その x におけるその三角関数の和を計算できる。DFTやフーリエ級数では元の関数の周期的な拡張がなされていると考えることができる。DCTでは、(離散的でない)コサイン変換と同様に、元の関数を偶関数に拡張することを意味する。

しかしながら、DCTは「有限」で「離散的」な数列に対して施されるものであるから、連続なコサイン変換にはない2つの微妙な問題が引き起こされる。まず、有限の点で定義された関数は定義域に左端と右端(すなわち最小の添字と最大の添字)とをもつので、その両方それぞれで偶対称であるか奇対称であるかを指定しなければならない。次に、関数の定義域は離散的であるので、どの位置に関して関数が偶/奇対称であるのかを指定しなければならない。例えば、abcd という均等に離れた4つの点の列を考えてみよう。そして例えば、の境界で偶対称であると指定したとしよう。このときどの位置で対称なのかという微妙な相違が生ずる。すなわち、データは点 a に関して偶対称であって偶関数への拡張は dcbabcd なのだろうか、それともデータは a とその前の点との中間点に関して偶対称であって拡張は dcbaabcd であるのか(a が繰り返すのか)?

これら2重の選択が、DCTと離散サイン変換 (DST) との標準的なさまざまな変種すべてを生じさせることになる。各々の境界は偶対称であるか奇対称であるかどちらかであることができ(これは2つの境界それぞれに2つの選択肢を与える)、さらに、各々の境界であるデータ点に関して対称か、2つのデータ点の中間点に関して対称かどちらかであることができる(同様に、これは2つの境界それぞれに2つの選択肢を与える)。結局、2 × 2 × 2 × 2 = 16 種類の選択肢がある。これらの選択肢のうちの境界が偶対称であるものがDCTとよばれ、選択肢の半分の8つのタイプに対応する。残りの半分がDSTの8つのタイプとなる。

これらは境界条件が異なるだけで施される変換はすべて離散フーリエ変換であるが、これらの違いは変換を応用する際にその用途に強く影響し、さまざまなDCTの変種に対してそれぞれに有用な特性を与えている。最も直接的には、偏微分方程式スペクトル法で解くために類フーリエ変換を用いるとき、境界条件は解かせることになる問題の一部として直接指定される。あるいはまた、(DCTのタイプIVに基づいている)修正離散コサイン変換 (modified DCT, MDCT) に対しては、境界条件はMDCTの本質的な特性である時間領域エイリアシングの消去に密接に関係している。もっと微妙なあり方ながら、境界は任意の類フーリエ級数において収束の速さに影響しているので、境界条件は画像や音声圧縮に対してDCTを有用なものとしているいわゆる「エネルギー圧縮」の特性を与える原因となっている。

特に、関数に不連続性があればフーリエ級数の収束率英語版を減少させることはよく知られている。同じ原理は信号圧縮に対して類フーリエ変換の有用性を決定している。よりなめらかな関数はそれをより正確に表すために必要となるDFTやDCTの係数がより少なくてすみ、より圧縮できることになる(ここで、「なめらかさ」について語るためにDFTやDCTをそれぞれ関数のフーリエ級数とコサイン級数の近似だとみなしている)。しかし、DFTがもつ非明示的な周期性は境界において通常不連続性を作り出すことを意味する。任意に選んだ信号の断片において左と右の境界の値が共に同じ値を持つということはめったに起こることではない。対照的に、「両方」の境界が「常に」偶対称であるDCTはこれらの境界において連続した拡張を与える(ただし一般にはその傾きは不連続である)。これがなぜDCTが、とりわけ(両方の境界が偶対称である)DCTのタイプ I, II, V, VI が一般にDFTよりも信号圧縮でよい成績を収めるのかという理由である。応用上は、こうした用途には一部には計算の容易さからDCT-IIが最も好まれている。

形式的定義[編集]

形式的には、1次元のDCT F: RNRN は、ある可逆な線形写像(ただし、R実数の集合)、または同じことであるが、ある正則な N × N 正方行列であって、以下に示された式で表される。ただし、これらの式では、N 個の実数列 x0, ..., xN−1N 個の実数列 X0, ..., XN−1 に変換される。

DCT-I[編集]

フーリエ変換や他の類似の変換と同じように式全体にかかる定数係数にはばらつきがあり、文献やライブラリによっては、DFTとの対応からこの式を 2 倍したものや、逆変換との対称性から {2/(N−1)}1/2 倍したものによって定義している場合もあるので注意を要する。また、x0xN−1 の項を 21/2 倍し、対応して X0XN−1 を 2−1/2 倍していることもある。後者の変更によって、ある定数倍を除いて変換は直交変換となるが、このときには実偶関数に対するDFTとは直接の関連を失うことになる。

DCT-Iでは、境界条件から 2(N − 1) 周期に拡張された関数を考えていることに対応するので、N ≥ 2 でないと定義できないことに注意されたい。他のタイプのDCTはすべて、N ≥ 1 であればよい。なお、N = 2 のときは上式の総和の項は消え、X0 = 1/2 (x0 + x1), X1 = 1/2 (x0x1) となる。

DCT-Iは(全体が2倍になる違いを除いて)、2(N − 1) 個の実数をもつ偶対称関数のDFTとまったく同じものである。たとえば、DCT-Iで N = 5 とし、5個の実数を abcde とすると、これは8個の実数 abcdedcb(偶対称)に対するDFTを 2 で割ったものになる。

DCT-Iは次の境界条件の場合に対応している:

  • xnn = 0 に関して偶対称、n = N − 1 に関して偶対称。
  • Xk についても同様。

DCT-II[編集]

DCT-IIは信号の圧縮分野などの応用では最も広く用いられている方法で、単にDCT (the DCT) と呼ばれることもある。DCT-Iと同様の理由により、これを2倍したものや、(2/N)1/2 倍したものとして定められている場合もあり、また直交化のために X0 の項のみ 2−1/2 倍されている場合もある。最後の場合DFTとの直接の対応は失われる。

このタイプは境界の両側に要素の間隔の半分のシフトを含んだ偶対称への拡張を考える。例えば N = 5 のときの実数列を abcde とすれば、2N = 10 個の実数列 abcdeedcba となる。両端の要素 a, e が繰り返される点がDCT-Iとは異なっている。ただし半分のシフトを行っているため、DFTとの対応を考える場合にはさらに倍にして偶数の添字の要素を 0 とした 4N 個の実数列をとる。すなわち、4N 個の実数列 y0, ..., y4N−1 を、

  • y2n = 0   (0 ≤ n < 2N である n について),
  • y2n+1 = y4N−2n−1 = xn   (0 ≤ n < N である n について)

を満たすものとすると、DCT-IIはこの実数列 yn をDFTで変換し 2 で割ったものと一致する。

DCT-IIは次の境界条件に対応する:

  • xnn = −1/2 に関して偶対称、n = N − 1/2 に関して偶対称。
  • Xkk = 0 に関して偶対称、k = N に関して奇対称。

DCT-III[編集]

DCT-IIIは(ある定数倍を無視すれば)DCT-IIの逆変換である。そのため、単に「逆DCT」(the inverse DCT, IDCT) と呼ばれることがある。

x0 の項を 倍することもある(対応する変形は上記DCT-II参照)。そうすると、DCT-IIとDCT-IIIとは互いに転置になる。DCT-IIIの行列は直交になるが、DFTとの直接の対応関係は失われる。

DCT-IIIは次の境界条件にあたる:

  • xnn = 0 で偶対称かつ n = N で奇対称。
  • Xkk = −1/2 で偶対称かつ k = N − 1/2 で偶対称。

DCT-IV[編集]

DCT-IVの行列は(定数倍を無視すれば)直交である。

DCT-IVの変種のひとつで、各変換のデータが「重なり合っている」変形を、修正離散コサイン変換 (modified DCT, MDCT) と呼ぶ。

DCT-IVは次の境界条件に対応する:

  • xnn = −1/2 で偶対称、n = N − 1/2 で奇対称。
  • Xk についても同様。

DCT-V – VIII[編集]

DCTタイプ I – IV は、実偶関数への偶数次DFTと等価である。原理的には、実際にはさらに、実偶関数への奇数次DFTに対応する4タイプのDCTタイプ V – VIII が存在する (Martucci, 1994)。タイプ V – VIII は、cos 関数の引数の分母に N + 1/2 の係数がある。ただし、タイプ V – VIII は実際にはほとんど使われない。

(自明な実偶配列である、1つの数 a への長さ 1(奇数)のDFTは、N = 1 のDCT-Vである)

逆変換[編集]

DCT-Iの逆変換は、DCT-Iの 2/(N − 1) 倍である。DCT-IVの逆変換は、DCT-IVの 2/N 倍である。DCT-IIの逆変換はDCT-IIIの 2/N 倍で、DCT-IIIの逆変換はDCT-IIの 2/N 倍である。

DFT同様、これらの変換公式の最前部にある標準化係数は便宜的なもので、扱いによって異なる。たとえば変換式を 倍(DCT-Iでは {2/(N − 1)}1/2 倍)する著者もおり、その場合何も乗算しなくても逆変換になる。

計算法[編集]

上記の公式を直接使うと、計算量は O(N2) となるが、高速フーリエ変換 (FFT) と同様の技法を使って、計算量を O(N logN) に減らせる。また、計算量 O(N) の事前処理と事後処理を加えることで、FFTそのものを使ってもDCTを計算できる。

当然ながら、ふつうは、最も効率がいいのはDCT専用のアルゴリズムであり、FFTはそれに及ばない(例外については後述する)。とはいうものの、DCTに特化したアルゴリズム(少なくとも2の冪乗個のデータに関しては、現在知られている中で最も計算量の少ないものを含め)は通常、FFTのアルゴリズムと密接に関連している。というのも、DCTは本質的には偶である実数データに対するDFTであるから、高速DCTのアルゴリズムの設計には、FFTを元にして、データの対称性に基づき冗長な計算を減らすことができる。この設計は自動化もできる (Frigo & Johnson, 2005)。クーリー・テューキーのアルゴリズム英語版に基づくものが最も一般的だが、FFTのアルゴリズムならこれに限らず何でも用いることができる。たとえば、ウィノグラードのアルゴリズム (Winograd algorithm) を用いると、加算の回数が増える代わりに乗算の回数を最小化することができ、一般的には効率があがる。同様なアルゴリズムは Feig & Winograd (1992) によってもDCT向きに提唱されている。DFT、DCT、および類似の変換法のアルゴリズムは互いに密接に関連しているため、どれかの変換法で改善が行われると、理論的には他の変換法にも即座に応用することができる (Duhamel & Vetterli, 1990)。

理論的には、FFTそのものを変更なしで用いた場合、DCT専用のアルゴリズムに比べいくらかのオーバーヘッドを伴うことになるが、この方法には明瞭な利点がある。高度に最適化されたFFTプログラムが広く出回っていることである。かくして実際には、一般的な N 長のデータを扱う場合、FFTを元にしたアルゴリズムの方が容易に性能を出せることが多い(現代の主なハードウェアの速度は、単純な計算量で決まるようなものではなく、プログラムの最適化によって、それに応じたハードウェアの改良も行われる)。一方、DCT専用アルゴリズムは、少量かつ固定長のデータ(たとえば、JPEGで用いられる 8 × 8 のDCT-II)向けや、音声圧縮用途の小規模なDCT(ないしMDCT)向けに広く用いられている。このような組み込みシステム用には、プログラムコードが短くて済むことも重要だからである。

実際のところ、通常のFFTを用いたDCTアルゴリズムといっても、それはしばしば、実数の偶関数データに対するより大規模なFFTから冗長な処理をそぎ落としたものと等価であり、計算量から見ても最適でありうる。たとえば、DCT-IIは 4N の偶対称な実数データ(偶数番目の要素が 0)に対するDFTと等価である。FFTを用いた一般的な計算法(たとえばFFTPACKFFTWに用いられている)の一つは Makhoul (1980) による。この手法は、実で偶なDFT(DCT-IIに対応する)における radix-4 時間間引きFFTの1ステップと見ることもできる(基数を 4 にする radix-4 ステップによって 4N 個のデータに対するDFTが4つのDFTに分解されることになり、それぞれのDFTは N 個の実数データに対するものとなる。4つのDFTのうち2つは 0 で、データが偶対称であることから、残りの2つは互いに等しくなる。かくして N 個の実数データに対するFFT 1回と、O(N) のバタフライ演算で計算できることとなる)。偶数の添字を持つ要素が 0 であるから、radix-4 ステップは split-radix ステップと正確に同じものである。続いて N 個の実数データに対するFFTを実データ split-radix FFT(Sorensen et al., 1987等)を用いて行えば、最終的な算法全体は、すでに述べた 2 の冪乗データに対するDCT-IIアルゴリズムのうち、最も計算量が少ないものに匹敵する(実数演算の回数が 2N log2NN + 2 のオーダーである[1])。したがって、計算量の点ではDCTをFFTで計算することが本質的に悪であるというわけではなく、単に使おうとしているFFTアルゴリズムの最適化の問題であることがある。アルゴリズムではなく実装上の問題であるが、データ量 N が小さい場合は、独立したFFTルーチンを呼び出すための関数呼び出しに伴うオーバーヘッドの方が問題になりうるほどである。

注釈[編集]

  1. ^ 正確には、実数演算の回数、特に実数乗算の回数は、変換式のスケーリングに幾分依存する。計算量 2N log2NN + 2 はDCT-IIについて前述の定義を用いた場合で、式全体が でスケーリングされていれば、乗算を2回節約できる。出力を個別にスケーリングすることが許されるなら、さらに乗算を減らせる。size-8 であるJPEGに関する結果を参照されたい (Arai et al. , 1988)。

参考文献[編集]

  • K.R. Rao and P. Yip, Discrete Cosine Transform: Algorithms, Advantages, Applications, 1990, Academic Press:Boston; K.R. Rao, P. Yip, and V. Britanak, 2006, 2 sub ed., ISBN 0-12-580251-X; 安田浩, 藤原洋訳 『画像符号化技術 — DCTとその国際標準』, 1992, オーム社, ISBN 4-274-03401-1
  • A. V. Oppenheim, R. W. Schafer, and J. R. Buck, Discrete-Time Signal Processing, second edition (Prentice-Hall, New Jersey, 1999).
  • S. A. Martucci, "Symmetric convolution and the discrete sine and cosine transforms," IEEE Trans. Sig. Processing SP-42, 1038–1051 (1994).
  • Matteo Frigo and Steven G. Johnson: FFTW, http://www.fftw.org/. フリー(GPL ライセンス)の C ライブラリで、任意の大きさの 1 次元と多次元の DCT(タイプ I–IV)を高速に計算できる。 M. Frigo and S. G. Johnson, "The Design and Implementation of FFTW3," Proceedings of the IEEE 93 (2), 216–231 (2005) も参照。
  • E. Feig, S. Winograd. "Fast algorithms for the discrete cosine transform," IEEE Transactions on Signal Processing 40 (9), 2174–2193 (1992).
  • P. Duhamel and M. Vetterli, "Fast Fourier transforms: a tutorial review and a state of the art," Signal Processing 19, 259–299 (1990).
  • John Makhoul, "A fast cosine transform in one and two dimensions," IEEE Trans. Acoust. Speech Sig. Proc. 28 (1), 27–34 (1980).
  • H. V. Sorensen, D. L. Jones, M. T. Heideman, and C. S. Burrus, "Real-valued fast Fourier transform algorithms," IEEE Trans. Acoust. Speech Sig. Processing ASSP-35, 849–863 (1987).
  • Y. Arai, T. Agui, and M. Nakajima, "A fast DCT-SQ scheme for images," Trans. IEICE 71 (11), 1095–1097 (1988).

関連項目[編集]