Scilab上のFM復調はTokyoFMの実データから音声に復調できたのでOK。しかしFM変調はイケません。過去回で自前のFSK変調関数なども作ってみたのですが「定義通り」に周波数を操作すると「位相がプッツン」な波形でちょっと納得がいきません。ましてやアナログ信号のFM変調はダメね。今回は再挑戦。FM変調出来たのか?
※「手習ひデジタル信号処理」投稿順 Indexはこちら
Scilab上でFSK/FM変調の苦闘?
苦闘はかなり前の過去回で、Scilabにcomm_tbx(コミュニケーション・ツール・ボックス)をインストールしたころに遡ります。まあ、comm_tbxはデジタル変調用ツール類なので、アナログ変調向けではないのですが、ひととおりツール(関数)そろっておるし、これ使えば一般的なデジタル変復調などやり放題っと思っていた私が馬鹿でした。
変調、復調関係の関数を扱ってみるとどうも変。ツールボックスが内部で使っている関数が古くて、使用中のScilabバージョンとミスマッチなのか、警告が出たりすることも度々。動作結果が出たものも何か怪しいっす。それだけでなくHELPファイルに書かれているとおりにやっているのに動作不可のものもあります。FSK変調デス。
いたしかたなし、ということで以下の過去回ではFSK変調を自前で試みましたです。
手習ひデジタル信号処理(83) Scilab、FSK、自力更生の必要なしだけれどやってみた
信号処理素人の浅はかさ、「定義のとおりに計算すればいいのじゃろ」と、周波数パラメータそのものを変調信号で切り替えるようなコードを書いてしまいました。掲げた変調結果自体は「位相が連続」なように見えますが、滑らかな波形が得られるように引数を「選んだ」結果であります。引数によってはプッチン、と位相が飛びます。これはダメだ。
今回、上記の反省のもと、ネットを漁っていてありがたいお言葉に目が開きました。東邦大学名誉教授 佐藤先生の以下のページであります。
上記より1か所引用させていただきます。
「VCOは周波数を変化させる回路」ではなく、「位相の変化速度を制御する回路」と考えるべきです。
FM(周波数変調)だけれども、周波数を操作しようとおもっちゃいけなかったのね。。。『打たれまい、打たれまいとするから。。。いや一歩すすんで打ってもらおう』星飛雄馬開眼の一言(うろ覚えですが)であります。
しかしね、実際にScilab上でFM変調のために位相を計算するのはどうしたらよろしいの?積分するといってまともにやるのは大変そうだし。と思って検索していたら、以下のページを発見しました。いつもお世話になっておりますStack Overflow様のページです。元の題材はGnu Octave上でのFM変調、それもキャリアは方形波です。
Sinusoidal FM modulation of a sqarewave carrier in gnu-octave
Matlab、OctaveとScilabは良く似ている(といいつつScilabだけ、離れ小島だったりすること多いっす)ので、Octave用のコードをScilabに置き換えて考えるのは難しくありません。キモは下の方に小さくかかれていた cumsum()関数の利用です。幸いOctaveと同名の関数がScilabにもありました。メンドイ積分などをアカラサマにせずともこの関数一発で済みそう。
Scilab上での実装
ということでScilab上で実験してみました。まずは、サンプリング周波数Fs、テスト用のアナログ信号波形の周波数Fsig、搬送波周波数Fcを定義。そして時系列tの数列(信号波形5波長分)を作ってみました。こんな感じ。
Fs=1000; Fsig=5; Fc=100; t=[0:1/Fs:5/Fsig];
続いて、信号波形と搬送波波形を生成します。どちらも単純明快な正弦波であります。こんな感じ。
Xsig=sin(2*%pi*Fsig*t); Xc=sin(2*%pi*Fc*t);
そしてFM変調波形の計算は以下のとおり。以下は変調指数は50%という設定です。このくらいにしないとグラフ描いたときの「見栄え」が悪いので。
Fdev=Xsig * Fc/2; z=sin(2*%pi*Fc.*t + 2*%pi*cumsum(Fdev)/Fs);
しっかりcumsum関数使わせてもらってます。これでFM変調できる筈。
結果をプロット
上記で生成した波形を時間プロットしてみました。
subplot(3,1,1); plot(t, Xsig); title('Signal'); subplot(3,1,2); plot(t, Xc); title('Carrier'); subplot(3,1,3); plot(t, z); title('FM');
SignalによってCarrierがFM変調された結果がFMのところに表示されておるようであります。
念のためFFTもかけてみました。plotFFT2は過去回で作成したScilab標準のFFT関数の自前ラッパです。正規化周波数でなく、Fsから具体的な周波数に変換して表示するだけのもの。
plotFFT2(Fs, z, 1);
グラフが以下に。
搬送波周波数 100Hz を中心に変調度50%ということで+50Hz、-50Hzのところに山がある感じ。こんなんで良いのか?