NVIDIA Jetson Nano 開発者キットでソフトウェア無線（動作確認編）

# はじめに

[**セットアップ編**](https://qiita.com/tsutof/items/75a8b8d9ac9f039b4016) で、Jetson Nano でソフトウェア無線（SDR: Software Defined Radio）を利用するための環境（cuSignal + SoapySDR）構築について紹介しました。受信機には [RTL-SDR BLOG V.3](https://www.rtl-sdr.com/rtl-sdr-blog-v-3-dongles-user-guide/) を使います。今回の動作確認編では、構築した環境上で、FM 放送の受信を例に、動作確認してみます。

# cuSignal に付属する Jupyter ノートブック

[cuSignal](https://github.com/rapidsai/cusignal) には、その使い方を示したいくつかの Jupyter ノートブックが付属し、ソフトウェア無線関連のノートブックもあります。

| ノートブック | 内容 | Jetson Nano で動作 |
|:-:|:-:|:-:|
| [online_signal_processing_tools.ipynb](https://github.com/rapidsai/cusignal/blob/branch-0.20/notebooks/sdr/online_signal_processing_tools.ipynb) | cuSignal をソフトウェア無線で利用する際の基本的なサンプル | :heavy_check_mark: |
| [sdr_integration.ipynb](https://github.com/rapidsai/cusignal/blob/branch-0.20/notebooks/sdr/sdr_integration.ipynb) | rtlsdr を利用した FM 放送受信サンプル | 部分的（設定を変更すれば動作） |
| [sdr_wfm_demod.ipynb](https://github.com/rapidsai/cusignal/blob/branch-0.20/notebooks/sdr/sdr_wfm_demod.ipynb) | SoapySDR を利用した FM 放送受信サンプル | :heavy_check_mark:  |

# 動作確認を解説

cuSignal 付属のノートブックもとても役に立ちましたが、自分自身の理解のためにノートブックを作成したので、それを以下に解説します。なお、作成したノートブックは [Gist で公開](https://gist.github.com/tsutof/6daf9e8221428d0e43b649de981d5498) しています。

:point_down: ライブラリのインポート

```python
import SoapySDR
from SoapySDR import * #SOAPY_SDR_ constants
import matplotlib.pyplot as plt
import numpy
import importlib
```

:point_down: **sample_rate** は RTL-SDR がデータをサンプルするサンプリングレートです。FM 放送音声のサンプリングレートである **audio_fs** との違いに注意。復調後にダウンサンプリングするための比率を **resample_factor** で定義しています。

```python
sample_rate = int(2.4e6)
audio_fs = int(48e3)
resample_factor = sample_rate // audio_fs
read_elements = 1024
channel = 0
```

:point_down: **fm_freq** は FM 放送局の周波数です。ここでは、80.0MHz の TOKYO FM を設定しています。

```python
fm_freq = 80.0e6 # FM Station Frequency
gain = 40
sample_time = 1 # Sampling time in sec
```

:point_down: メモリ消費は大きいですが、動作が分かりやすいので、まず、受信データを **sample_time** 秒間蓄えてから、変調します。**total_elements** はそのサンプル数です。SoapySDR の readStream 関数は（私のやり方が間違っている可能性もありますが）一度に多くのサンプルを読み出すと動作が不安定だったので、複数回に分けて読み出し、**num_reads** がその回数。

```python
total_elements = int(sample_rate * sample_time)
num_reads = total_elements // read_elements
```

:point_down: 受信デバイス（ここでは RTL-SDR）の初期設定を行います。[このページ:page_facing_up:](https://github.com/pothosware/SoapySDR/wiki/PythonSupport#basic-example) を参考にしました。

```python
#enumerate devices
results = SoapySDR.Device.enumerate()
for result in results: print(result)

#create device instance
#args can be user defined or from the enumeration result
args = dict(driver="rtlsdr")
sdr = SoapySDR.Device(args)

#query device info
print(sdr.listAntennas(SOAPY_SDR_RX, channel))
print(sdr.listGains(SOAPY_SDR_RX, channel))
freqs = sdr.getFrequencyRange(SOAPY_SDR_RX, channel)
for freqRange in freqs: print(freqRange)

#apply settings
sdr.setSampleRate(SOAPY_SDR_RX, channel, sample_rate)
sdr.setFrequency(SOAPY_SDR_RX, channel, fm_freq)
sdr.setGain(SOAPY_SDR_RX, channel, gain)
```

    {driver=rtlsdr, label=Generic RTL2832U OEM :: 00000001, manufacturer=Realtek, product=RTL2838UHIDIR, serial=00000001, tuner=Rafael Micro R820T}
    ('RX',)
    ('TUNER',)
    2.3999e+07, 1.764e+09

:point_down: 受信データ読み込み用のバッファ・メモリを用意します。

```python
buff_len = total_elements
cpu_buff = numpy.array([0] * buff_len, numpy.complex64)
print('buff_len={}'.format(buff_len))
```

    buff_len=2400000

:point_down: 受信機からデータを読み込みます。先ほど述べたように、readStream 関数で一度に多くのデータを読み出すと動作が不安定だったため、複数回に分けて読み出しています。最初の読み出しに失敗することもあったので、空読みしてから、本当の読み出しを行っています。

```python
#setup a stream (complex floats)
rxStream = sdr.setupStream(SOAPY_SDR_RX, SOAPY_SDR_CF32)
sdr.activateStream(rxStream) #start streaming

while(True):
    sr = sdr.readStream(rxStream, [cpu_buff], 1024)
    if sr.ret == 1024:
        break

offset = 0
for i in range(num_reads):
    sr = sdr.readStream(
        rxStream, [cpu_buff[offset:offset + read_elements]], read_elements, timeoutUs=int(8e12))
    if sr.ret != read_elements:
        print('readStream error sr.ret={}'.format(sr.ret))
        break
    offset += read_elements

#shutdown the stream
sdr.deactivateStream(rxStream) #stop streaming
sdr.closeStream(rxStream)
```

:point_down: パワースペクトラム密度を求める関数。**sig** を **cusignal** の別名と、**scipy.signal** の別名に切り替えることで、両方のライブラリでの実行に対応します。

```python
def welch(buff, sample_rate):
    return sig.welch(buff, sample_rate, nfft=1024, scaling='density', return_onesided=False)
```

:point_down: FM 復調を行う関数。こちらも、**cupy** と **numpy** に対応。受信機から読み出せるデータは I-Q 表現と呼ばれる複素数形式で、[np.angle](https://numpy.org/doc/stable/reference/generated/numpy.angle.html) 関数で位相角が求まり、[np.unwrap](https://numpy.org/doc/stable/reference/generated/numpy.unwrap.html) 関数で位相角のシフトを求めます。[np.diff](https://numpy.org/doc/stable/reference/generated/numpy.diff.html) 関数で微分します。次に、ダウンサンプリングで、受信機のサンプリングレートから、FM 放送の音声サンプリングレートへ変換。最後に正規化を行って完了です。

```python
def demodulate(buff, resample_factor):
    b = buff
    b = np.diff(np.unwrap(np.angle(b)))
    b = sig.resample_poly(b, 1, resample_factor, window='flattop')
    b /= np.pi
    return b
```

:sweat_smile: 私は、説明できる程の知見がないので、FM 変調の仕組みはこの辺でご勘弁を。もう少し詳しく知りたい方は、[Interface 2021年5月号](https://interface.cqpub.co.jp/magazine/202105/) の特集「Pythonで無線信号処理」をお勧めします。

## CuPy + cuSignal

まずは、CuPy と cuSignal で、FM 復調を実行します。

:point_down: **cupy** の別名を **np** へ、**cusignal** の別名を **sig** に設定します。

```python
np = importlib.import_module('cupy')
print(np.__name__)
sig = importlib.import_module('cusignal')
print(sig.__name__)
```

    cupy
    cusignal

:point_down: GPU が利用できるメモリを割り当て、そこへ、先程、受信機から読み出したデータをコピーします。

```python
gpu_buff = sig.get_shared_mem(buff_len, dtype=np.complex64)
gpu_buff[:] = cpu_buff
```

:point_down: パワースペクトラム密度を計算。

```python
f, Pxx_den = welch(gpu_buff, sample_rate)
```

:point_down: 計算結果をグラフ表示

```python
plt.semilogy(np.asnumpy(np.fft.fftshift(f/1e4)), np.asnumpy(np.fft.fftshift(Pxx_den)))
plt.show()
```

![output_13_0.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/384208/5b80c7a8-7545-c4cd-a57f-06d944520d6d.png)

なぜ中心周波数に凹みがあるのか？不明。私のコードに問題があるのかも知れません。   

:point_down: FM 復調を行います。 

```python
%%time
b = demodulate(gpu_buff, resample_factor)
```

処理時間は後の方でまとめてあります。

:point_down: 音声波形を表示。

```python
b = np.asnumpy(b).astype(np.float32)
```

```python
x = [i for i in range(len(b))]
plt.plot(x, b, label="test")
```

    [<matplotlib.lines.Line2D at 0x7f311a5550>]

![output_17_1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/384208/f6fad242-2388-8765-eba8-11213ea1dc5c.png)

## NumPy + SciPy Signal

次に、NumPy と SciPy Signal で、FM 復調を実行します。

:point_down: **numpy** の別名を **np** へ、**scipy.signal** の別名を **sig** に設定します。

```python
np = importlib.import_module('numpy')
print(np.__name__)
sig = importlib.import_module('scipy.signal')
print(sig.__name__)
```

    numpy
    scipy.signal

:point_down: パワースペクトラム密度を計算。

```python
f, Pxx_den = welch(cpu_buff, sample_rate)
```

:point_down: 計算結果をグラフ表示

```python
plt.semilogy(np.fft.fftshift(f/1e4), np.fft.fftshift(Pxx_den))
plt.show()
```

![output_21_0.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/384208/3f87fd91-5e6d-2938-2fe2-42f9337b1426.png)
    
データが同じなので、CuPy + cuSignal で計算した場合と同じ結果です。

:point_down: FM 復調を行います。

```python
%%time
b = demodulate(cpu_buff, resample_factor)
```

処理時間は後の方でまとめてあります。

:point_down: 音声波形を表示。

```python
x = [i for i in range(len(b))]
plt.plot(x, b, label="test")
```

    [<matplotlib.lines.Line2D at 0x7f311da5e0>]

![output_23_1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/384208/9b971087-4de5-7bc6-dcff-6373696cf71f.png)

こちらも、データが同じなので、CuPy + cuSignal で復調した場合と同じ結果です。
    
:point_down: 最後に音声データをファイルに保存します。

```python
from scipy.io import wavfile
wavfile.write('demod.wav', rate=audio_fs, data=1e2*b)
```

# 処理時間

音声1秒間分のFM復調処理（10回計測した平均値）
動作電力モード：MAXN（sudo jetson_clocksも実行）

|   |CPU Time <br> User [msec] | CPU Time <br> Sys [msec] | CPU Time <br> Total [msec] | Wall Time [msec] |
|---|---|---|---|---|
| CuPy + cuSignal <br>（GPU利用） | 13.6 | 2.4 | 16.0 | **126.0** |
| NumPy + SciPy Signal <br>（CPUのみ） | 530.0 | 86.8 | 616.0 | **619.0** |

# 最後に

前回の [セットアップ編](https://qiita.com/tsutof/items/75a8b8d9ac9f039b4016) と、今回の動作確認編で、Jetson Nano でソフトウェア無線を試す環境が整いました。しかし、まだ、FM 復調を試しただけで、ソフトウェア無線の目的である、ハードウェアを変更せずに、いろいろな通信方式に対応することには、ほぼ遠いです。Jetson Nano を使うからには、ディープラーニング推論も組み合わせて、面白い使い方を探ってみたいです。例えば、AI ノイズ除去など。私の勉強不足で、そこまで行くのは大変ですが、もし、何かできたら、続きの Qiita 記事を作成したいと思います。