FM音源プログラム (9) -- オペレータ (2)

今回は、サイン波テーブルについての話です。
「サイン波テーブル」と書きましたが、DSP など、「計算の得意な計算機」であれば、sin / cos 関数のテイラー (Taylor) 展開式
\qquad \sin(x) = x - \frac{x^3}{3!} + \frac{x^5}{5!} - \frac{x^7}{7!} + ...
\qquad \cos(x) = 1 - \frac{x^2}{2!} + \frac{x^4}{4!} - \frac{x^6}{6!} + ...
を計算して値を求めることができます。
三角関数の対称性を利用すると、x の値は 0 \le x \le \pi / 4 の範囲に還元することができます。
上の式の誤差が最も大きくなるのは x= \pi/4 ですが、その時でも、x^7 の項までで計算を打ち切った場合の誤差を調べると、2 進数で約 15 ビット精度となります。
一般の「計算の不得意な計算機」では、この程度の計算でも負担となりますから、値を事前に計算しておいた表を利用することになります。
三角関数の周期は 2\pi で、無理数ですから、そのままでは整数個のエントリを持つ表のアドレスとして使うわけにはいきません。
どっちみちスケーリングする必要があるのですから、0 〜 2\pi の範囲が、N ビットアドレスの表の 2^N 個のエントリの範囲と一致するような、
\qquad \sin(x) = {\rm SIN}[ x \cdot 2^N / 2\pi ]
というマッピングをします。
ここで、{\rm SIN}[ \cdot ] は周期 2^N のサイン波テーブルを表すものとします。
もとの三角関数の対称性
\qquad \sin(x + \pi) = -sin(x)
\qquad \sin(\pi - x) = sin(x)
より、サイン波テーブルでも
\qquad {\rm SIN}[ k + 2^{N-1} ] = -{\rm SIN}[ k ] \qquad \qquad (0 \le k < 2^{N-1})
\qquad {\rm SIN}[ 2^{N-1} - k ] = {\rm SIN}[ k ] \qquad \qquad \qquad(0 \le k < 2^{N-2})
という性質がありますから、1周期 2^N エントリの 1/4 の 2^{N-2} エントリのテーブルがあれば十分なことが分かります。
ただし、この 1/4 サイズに還元する操作のための実行時間が余計にかかることになります。
この表を補間をせずに使う、つまり、「位相」を表のアドレスとして使うときに下位ビットを切り捨ててアクセスした表のデータをそのまま使う場合には、当然、位相方向の誤差が生じます。
この誤差は、実際には振幅方向の誤差に振り替わって現れてきます。
この誤差と、表を作成する際に値のデータを量子化したために発生した誤差とが、同程度になる条件を求めてみます。
もとの三角関数で考えると、sin / cos 関数の微分係数の最大値は 1 ですから、振幅方向の量子化と、位相方向の離散化の幅が同程度になれば良いことが分かります。
もとの三角関数の振幅 ±1 を符号付き M ビットで量子化すると量子化ステップサイズは 2 / 2^M = 2^{-(M-1)} となります。
一方、位相方向に N ビットで離散化すると、周期 2\pi2^N 分割するわけですから、離散化のステップサイズは 2\pi / 2^N = 2\pi \cdot 2^{-N} となります。
離散化ステップサイズが量子化ステップサイズより小さくなる条件は、
\qquad 2\pi \cdot 2^{-N} \, < \; 2^{-(M-1)}
\qquad \pi \, < \, 2^{N-M}
\pi\, < \, 4 ですから、4 で置き換えて大きめに見積もると、
\qquad 4 \, < \, 2^{N-M}
\qquad \log_2(4) \, < \, \log_2(2^{N-M})
\qquad 2 \, < \, N - M
\qquad N \, > \, M + 2
となります。
FM音源プログラムは、最初、12 ビット DAC を持つ ADuC7026 用に作り始めたので、最終的な精度としては 12 ビットを想定していますので、 M = 12 となります。
サイン波テーブルのデータ自体としては、int16 型として作成してあります。
したがって、N = M + 2 = 12 + 2 = 14 で、サイン波1周期としては 16 K エントリのサイズとなります。
FM音源プログラムでは、当初は、この 1/4 サイズの 4 K エントリの表を使っていましたが、少しでも高速化するために、現在では 1/2 サイズの 8 K エントリの表を使っています。