https://qiita.com/0x20FE/items/0e229a4c31889903d019 参照

ADALM PLUTOでFM受信してみる

概要

ソフトウエア無線によるFM受信の技術メモです。RFアンプとAD変換部分は、Analog DevicesのPLUTOを使い、Windows環境のpythonにてFM受信を行ってみました。

※更新

遅延検波の解説を更新(2021/9/23)
ディエンファシスの追加(2022/6/19)　＜-忘れてた。
読みにくかったのでリライトしました(2022/8/12,2023/9/14)
[NEW] 実行環境の準備に、PLUTOの環境変数の設定を追加(2024/6/28)

関連記事
ADALM-PLUTOを有線LANに接続する (2022/8/12)

想定するアプリケーションと実行環境

想定するアプリケーション：SDR、FM受信、音声処理
実行環境：python, windows10, ADALM PLUTO
他のデバイスでも受信処理とオーディオ出力の部分はそのまま使えるかと。

実行環境の準備

PLUTOを使うために次のリンクに従ってWindows版のLibiio,python,pyadi-iioをインストールする。

初期設定では、ローカル周波数の設定範囲が(325 - 3800MHz）であるため、（teratermなど使って）シリアル通信でログインし、環境変数を設定しrebootすることで、周波数設定範囲を(70 - 6000 MHz)に拡大する。
（参考）Customizing the Pluto configuration Updating to the AD9364

シリアルコンソールから

# fw_setenv attr_name compatible
# fw_setenv attr_val ad9364
# reboot

解説

FM検波理論（遅延検波）

瞬時周波数成分を取り出す（FM検波、FM復調）方法として遅延検波を利用します。

s(t)を伝送信号としたとき、ベースバンドFM信号の瞬時周波数f(t)は
Kは定数
ベースバンドFM信号をサンプリング間隔Tでサンプルすると

（fは瞬時周波数）

この信号を１サンプル遅延させ複素共役した値と積算する。

したがって、np.angle()で複素数の偏角を求めで割ると瞬時周波数fが得られる。

遅延検波

# (i) bbIQ
# (o) dem
dem = np.angle(bbIQ[1:] * np.conj(bbIQ[:-1])) / math.pi

サンプル間隔の選定

サンプル間隔が大きくなると前記の偏角がを超え、計算結果が不連続となり異常値が発生する。逆にサンプル間隔が小さいすぎると偏角の値が小さくなり、量子化雑音の影響が大きくなる。このためサンプル間隔Tは注意深く選定する必要がある。
となるは、

ディエンファシスフィルタの導出

日本のFM放送はτ=50usのエンファシスがかかっているため、ディエンファシスフィルタをかける。
エンファシス
ディエンファシス

デジタルフィルタ構成とするためHd(s)に双一次変換を施す。

サンプリング周波数を44100Hzとすると

ディエンファシスフィルタ

# (i) aud
# (o) demp
deemp_n = np.array([1,1])*2
deemp_d = np.array([5.41,3.41])
demp = scipy.signal.lfilter(deemp_n,deemp_d,aud)

Pythonでのサウンド出力（PyAudio）

別記事を参照してください。

pyaudio

#初期化
pyAud = pyaudio.PyAudio()
#サンプリング周波数、チャンネル、フォーマット
sndStrm = pyAud.open(rate=SND_FS,channels=SND_CH_OUT,format=pyaudio.paInt16,input=False,output=True)

# 途中の繰り返し処理
# (i) pack
# (o) sndStrm
# (g) SND_CHUNK バッファ長
while pack.size > SND_CHUNK :
	# ±1に正規化された信号をSND_CHUNK毎に16bit固定少数変換して出力
    aBuf = pack[0:SND_CHUNK]*8000
    sndBuf = aBuf.astype(np.int16).tobytes()
    sndStrm.write(sndBuf)
    pack = pack[SND_CHUNK:]
    
    
#後始末    
sndStrm.stop_stream()
sndStrm.close()
pyAud.terminate()

サンプルプログラム

sdr=SdrInit()にてPlutoの設定を行い、オブジェクトを得る。
sdr.rx()でサンプルを取り出し、demLoop()にて遅延検波処理を行い音声データを得る。
遅延検波の後で、音声デバイスのサンプルレートにダウンサンプルし、デエンファシスをかける。
メインループでは、demLoopで取り出した音声データを、一旦バッファ(pack)に格納し、pyaudioのチャンクバッファサイズよりデータが溜まったところで、チャンクサイズごとに取り出し、pyaudio用のバイナリデータに変換してpyaudioのデバイスに（ブロックモードで）書き込む。

PlutoFM.py

import adi                          #pluto 
import matplotlib.pyplot as plt     #graph
import numpy as np
import scipy.signal                 # US/DS/Filter
import math                         # pi.sic.cos
import pyaudio

# System parameters
################################################
#-- application
################################################
SND_FS = 44100 # 48000
SND_CHUNK = 1024*2
SND_CH_OUT = 1

sample_rate_audio = SND_FS
################################################
#-- baseband IQ parameter
################################################
#(i) sample_rate_audio
decimate_demod = 5 
sample_rate_demod = sample_rate_audio * decimate_demod
################################################
#-- radio interface param
################################################
# PLUTO
# LO 325MHz to  3800MHz
# BW < 20MHz
# FS 521ksps to 61.4Msps # cat spec is 65.2Ksps to  61.4Msps
################################################
#(i) sample_rate_demod
decimate_capture = 4 
sdrSample_rate = sample_rate_demod * decimate_capture

sdrBufSize = 1024*20
sdrFcenter = 82500000 # 中心周波数82.5MHzの場合
#sdrFcenter = 80000000
sdrBandWidth = 400000
################################################
# Create radio object
def SdrInit():
    # (i) 
    # (o) radio Object
    sdr = adi.Pluto()
    sdr.rx_rf_bandwidth = sdrBandWidth
    sdr.rx_lo = sdrFcenter
    sdr.sample_rate = sdrSample_rate
    sdr.rx_buffer_size=sdrBufSize
    #sdr.tx_lo = 2000000000
    #sdr.tx_cyclic_buffer = True
    #sdr.tx_hardwaregain = -30
    sdr.gain_control_mode = "slow_attack"
    # Read back properties from hardware
    #print(sdr.rx_hardwaregain)
    # Get complex data back
    sdr.rx_enabled_channels = [0]
    return sdr

#=========================================================
#　FM　復調
#=========================================================
deemp_n = np.array([1,1])*2
deemp_d = np.array([5.41,3.41])
demLastIQ=0
def demLoop(demLastIQ):
    smp = sdr.rx()
    decim = scipy.signal.decimate(smp,decimate_capture)
    bbIQ = np.insert(decim,0,demLastIQ)
    #bbIQ = decim
    demLastIQ = decim[-1]
    #------------------------
    # 遅延検波
    # -+-[Delay]-- [x]---[ angle ]-->
    #  |            ^
    #  +-[Conj ]----+
    #------------------------
    dem = np.angle(bbIQ[1:] * np.conj(bbIQ[:-1])) / math.pi 
    #  decimate to audio sample rate
    aud = scipy.signal.decimate(dem,decimate_demod)
    #deEmphasis
    demp = scipy.signal.lfilter(deemp_n,deemp_d,aud)
    return demp,demLastIQ

if __name__ == "__main__":
    sdr = SdrInit()
    #====================================================
    #--- setup Audio
    #====================================================
    pyAud = pyaudio.PyAudio()
    sndStrm = pyAud.open(rate=SND_FS,channels=SND_CH_OUT,format=pyaudio.paInt16,input=False,output=True)
    #====================================================
    #--- setup FIFO
    #====================================================
    pack = np.array([])
    #====================================================
    #--- FM受信ループ
    #====================================================
    for cnt in range(2000): # while: でもよい
        #====================================================
        #--- チャンクサイズ以上たまったらサウンドデバイスに出力
        #====================================================
        while pack.size > SND_CHUNK :
            aBuf = pack[0:SND_CHUNK]*8000
            sndBuf = aBuf.astype(np.int16).tobytes()
            sndStrm.write(sndBuf)
            pack = pack[SND_CHUNK:]
            #plt.plot(aBuf,linestyle='None',marker="x")
            #plt.psd(aBuf,NFFT=1024,Fs=sample_rate_audio , Fc=0)
        #====================================================
        #--- 新たなIF信号を受信し音声データに復調しFIFOにつめる
        #====================================================
        aud,demLastIQ = demLoop(demLastIQ)
        pack = np.append(pack, aud)
    #====================================================
    # teardown 処理
    #====================================================
    sndStrm.stop_stream()
    sndStrm.close()
    pyAud.terminate()

結言

パラメータ調整は必要ですが、遅延検波で簡単にFMを復調できます。

雑多

2025年12月20日土曜日

PLUTOと復調理論とPython code 少し