GNURadio QuadratureDemodブロックのアルゴリズム

TL;DR

FSKをGNURadioで復調する時に使う、QuadratureDemodブロックは、$n$番目のサンプリングで得られたI/Q値を \[ \begin{aligned} x_n &= A\mathrm{e}^{j\omega t + \phi}\\ &= \mathrm{exp}(2 \pi j\frac{F_n}{F_s}n + \phi) \end{aligned} \] としたとき、

\[ \mathrm{arg}(x_n \times x_{n-1}^*) \] を計算している。これはI/Q信号を極座標表示した時の偏角の時間微分であり、まさに周波数である。 ただし、$F_s$はサンプリング周波数、$F_n$は信号波の周波数とする。

---

GNURadioのWikiにある、Quadrature Demodを見ていたら、以下の記述があった。

Mathematically, this block calculates the product of the one-sample delayed input and the conjugate undelayed signal, and then calculates the argument of the resulting complex number… (数学的に、このブロックは、1サンプル遅れた入力と、遅れていない信号の複素共役の積を求め、出てきた複素数の偏角を計算する)

とある。この文章の後、数式がかかれているようだが、なぜか数式がレンダリングされない。というか、$\TeX$の構文が画像になっている。せっかくなので、ちょっと式をフォローしてみよう。

任意のサンプリングにおける周波数の算出

SDRで受信して出てきたI/Q信号は、複素数で書ける。任意の時刻$t$における信号$x$は、

\[ x = A\mathrm{e}^{j\omega t} \] と書ける。ここで、$\omega = 2 \pi f$は角振動数であり、$A$は振幅である。

SDRでは、サンプリング周波数$F_s$の周期$T_s = \frac{1}{F_s}$毎に取得されたI/Q信号が得られる。よって、$n$番目のサンプリングを得た時刻$t$は、$t = T_s n$であるから、 $n$番目のサンプリングで得たI/Q信号は

\[ x_n = A\mathrm{exp} [j\omega_n T_s n + \phi] \] となる。ここで、$\omega_n$はn番目のサンプリングで得られたI/Q信号の角周波数、また$\mathrm{exp}(a) = \mathrm{e}^a$である。

求める値は、

1サンプル遅れた入力と、遅れていない信号の複素共役の積を求め、出てきた複素数の偏角

であるから、以下の式になる。

\[ Y_n = \mathrm{arg} [x_{n-1} \times x_n^* ] \]

まず、右辺括弧内を計算しよう。計算は簡単。手の運動。

\[ \begin{aligned} y_n &= x_{n-1} \times x_n^* \\ &= A\mathrm{exp} [j\omega_{n-1} T_s (n-1) + \phi] \times A\mathrm{exp} [-(j\omega_{n} T_s n + \phi)] \\ &= A^2\mathrm{exp} [j\omega_{n-1} T_s (n-1) + \phi -(j\omega_{n} T_s n + \phi)] \\ &= A^2\mathrm{exp} [jT_s \{ n(\omega_{n-1} - \omega_n) - \omega_{n-1} \}] \\ &= A^2\mathrm{exp} [jTs \{ \frac{\omega_{n-1} - \omega_n}{\omega_n}n - \omega_{n-1} \}] \end{aligned} \]

ここで、FSKの周波数偏移を$\delta f$とすると、 $|\omega_{n-1} - \omega_n| \le 2\delta\omega$である。FSKの場合、キャリア周波数に比べて周波数偏差は充分小さい。たとえば、920MHz帯で50kbpsであれば、周波数偏移25kHzである。よって、

\[ |\omega_{n-1} - \omega_n| \le 2\delta\omega \ll \omega_n \]

となり、 \[ \frac{\omega_{n-1} - \omega_n}{\omega_n} \approx 0 \] と近似することができる。 前述の例に適用して確認してみると、

\[ \frac{2\times 50 \times 10^3}{920 \times 10^6} \approx \times 10^{-4} \] となる。あまり良い近似ではないが、実用上問題ないだろう。

以上のことから、 \[ \begin{aligned} y_n &\approx A^2\mathrm{exp}(-jT_s \omega_{n-1}) \\ Y_n &\approx \mathrm{arg}(y_n) \\ &= -T_s \omega_{n-1} = -2\pi \frac{F_{n-1}}{F_s} \end{aligned} \] となる。

結果の評価

げげげっ。$F_{n-1}$も$F_s$も周波数なので0より大きいから、$Y_n$が負の値になってしまう。 物理的に考えると、$Y_n$は、複素数の偏角の時間微分、すなわち周波数を計算しているのであるから$Y_n$は正の実数でなければならない。

考えてみると、$Y_n$は、$n-1$番目の偏角から$n$番目の偏角を引いている。これは引き算の順序が逆である。そういえば、件の文章の後の数式もちゃんと$\TeX$を読んでみると、

\[ y[n] = …. = \mathrm{arg} (A\mathrm{e}^{j 2\pi \frac{f}{f_s}n + \phi_0} \mathrm{e}^{-j 2\pi \frac{f}{f_s} (n-1)}) \]

みたいなことが書いてある。これって、$n$番目から$n-1$番目の偏角を引いている… だまされた。

ということで、GNURadioのQuadrature DemodのWikiおよびドキュメントは間違っていて、 たとえば、

Mathematically, this block calculates the product of the undelayed input and the complex conjugate of the one-sample delayed signal, and then calculates the argument of the resulting complex number

と修正すべきである。

一応、ソースを見て私が間違ってないか確認しなきゃ。