手習ひデジタル信号処理(92) Scilab、デシメーション処理、サンプルデータで動作確認

Joseph Halfmoon

前回まで「お求めやすい」ソフトウエア無線受信USBドングルRTL-SDRから取得した信号をScilabへ引き込むべくPythonプログラムなど作成。一応の完成後、間が空いてしまいました。忘却力の年寄はほぼ完全に忘れております。今回からデータを受け取る側のScilabでの処理をやっていこうと思います。また直ぐに忘れる?

※「手習ひデジタル信号処理」投稿順 Indexはこちら

FMラジオ放送の復調処理のおさらい

前回作成のRTL-SDRで受信した信号の「輸出」用のプログラム、短い時間範囲ですが確認用に音声波形も出力できるようになってました。それと同じことをScilab上でも処理できれば音声波形への復調(リアルタイムじゃないけれど)ができるというものです。

Python上でコードした大まかな処理ステップは以下のようです。

    1. オーディオ信号のサンプリング周波数Faは、44100Hz
    2. RTL-SDRで取り出すIQ信号のサンプリング周波数Fsは、Faの48倍の2.1168MHzに設定
    3. Fsでサンプリングされた信号を8:1のダウンサンプリング
    4. 上記信号の位相の変化を検出して復調処理
    5. 復調処理した信号をさらに6:1のダウンサンプリングしてFa周波数の音声信号とする

ここで2回登場するダウンサンプリング処理ですが、Pythonプログラム上では、以下のScipy所蔵の関数にお願いしてました。

scipy.signal.resample()

サンプル周波数の変更(ダウンサンプリング、アップサンプリング)ができる関数です。ダウンサンプリングの場合、自動で元信号にローパスフィルタ(アンチ・エイリアス)を掛けてからサンプリングしなおしてくれる優れものです。フィルタ特性などは「よきにはからって」くれます。

一方、お世話になっておりますScilabのツールボックス comm_tbx には、

downsample

という関数がありますが、「素のダウンサンプル」でアンチエイリアスフィルタなどは別途自分でやれ、ということみたいです。フィルタ(ポリフェーズ化された)を適用してからダウンサンプルしてくれるデシメーション用の以下の関数もありますが、

polyphase_decimation

一度にやってくれるだけで、こちらもフィルタの設定などは自分でやれということみたいっす。どうもScipyのように「お任せ」ではできないようです。あたりまえか。

Scilab上でのデシメーション処理の準備

そこで、自力でデシメーション処理(アンチエイリアスフィルタ付)をすることといたしました。まず、過去回で作成済の小ネタの自前関数共を召喚いたします。こんな感じ。

exec('genSin.sci', -1)
exec('genFIRLP2.sci', -1)
exec('plotFFT2.sci', -1)

また、具体的な周波数やデータ点数なども値を決めときます。

Fa = 44100;
Fs = Fa * 48;
Fcutoff = Fs / 8 / 2;
nData=1024*48;
被テスト波形の生成

まずは、正弦波をいくつか組み合わせて被テスト波形を作り、それをデシメーション処理してみて、所望の周波数が得られているか確認したいと思います。

サンプル波形の生成は自前関数を使って以下で。

[t5k, y5k]=genSin(Fs, nData, 5000, 0, 1);
[t10k, y10k]=genSin(Fs, nData, 10000, 0, 1);
[t50k, y50k]=genSin(Fs, nData, 50000, 0, 1);
[t300k, y300k]=genSin(Fs, nData, 300000, 0, 1);
sample=y5k+y10k+y50k+y300k;

8:1のデシメーション処理をするので、Fs=2.1168MHzに対して、デシメーション処理後のサンプリング周波数は264600Hzということになります。上記にはナイキスト周波数132300Hz以上の300kHz波形が含まれておりますが、こいつは抑制しないとエイリアシングが起こるということです。

まずは生のサンプル波形のFFTプロット。自前関数でプロットっす。

plotFFT2(Fs, sample, 1);

こんな感じ。sample

下の周波数が見ずらいので、下の方を拡大したグラフが以下に。sampleZoom

5k、10k、50k、そして300kHzにピークがならんでおります。

アンチエイリアス無でダウンサンプリング

さて上記のサンプル波形を「素でダウンサンプル」してみるとどうなるか。当然エイリアスが折り返されてくる筈。

sampleN8 = downsample(sample, 8);
plotFFT2(Fs/8, sampleN8, 1);

プロット結果が以下に。黄色のマーカー部分が折り返されてきている「いらない」やつデス。aliasAfterDownSampleNoFilter

 

アンチエイリアス・フィルタ

どんなフィルタでLPFすべきか、3個くらい作成して特性を比較してみました。

coeff10=genFIRLP2(Fs, 10, Fcutoff, 1);
coeff100=genFIRLP2(Fs, 100, Fcutoff, 1);
coeff200=genFIRLP2(Fs, 200, Fcutoff, 1);

最初は、フィルタ長10点。短いですけど、「急峻」とはお世辞にも言えず。

coeff10magPlot

 

フィルタ長100点の場合が以下に。頭の平なところが微妙に波打っているけどまあ急峻?coeff100magPlot

フィルタ長200点。フィルタの規模的には上記の2倍に巨大化しているけれどもそれほど特性が違うわけでもない?coeff200magPlot

今回は真ん中の奴で実験してみることに。

アンチエイリアス・フィルタ適用後ダウンサンプリング

さて、サンプル波形に上記のアンチエイリアス・フィルタを適用してから、ダウンサンプルするコードが以下に。

af100=filter(coeff100, 1, sample);
sampleAF8=downsample(af100, 8);
plotFFT2(Fs/8, sampleAF8, 1);

FFTかけたグラフが以下に。sampleAF8

黄色のマーカ部分、-40dBくらいになっているので、まあ「許せる」のではないかと。知らんけど。

手習ひデジタル信号処理(91) PythonでFMラジオデータを音声再生、RTL-SDR に戻る

手習ひデジタル信号処理(93) Scilab、デシメーション処理、実データで動作確認 へ進む