Skip to content

Commit

Permalink
Merge pull request #25 from drowe67/dr-bbfm
Browse files Browse the repository at this point in the history
Baseband FM PoC
  • Loading branch information
drowe67 authored Nov 7, 2024
2 parents aaa8f85 + dacd4e1 commit 490dfcc
Show file tree
Hide file tree
Showing 29 changed files with 3,155 additions and 407 deletions.
83 changes: 83 additions & 0 deletions BBFM.md
Original file line number Diff line number Diff line change
@@ -0,0 +1,83 @@
# Radio Autoencoder - Baseband FM (BBFM)

A version of the Radio Autoencoder (RADE) designed for the baseband FM channel provided by DC coupled and passband FM radios, e.g. land mobile radio (LMR) VHF/UHF use case.

# BBFM ML encoder/decoder

1. First pass training command line:
```
python3 ./train_bbfm.py --cuda-visible-devices 0 --sequence-length 400 --batch-size 512 --epochs 100 --lr 0.003 --lr-decay-factor 0.0001 --plot_loss ~/Downloads/tts_speech_16k_speexdsp.f32 model_bbfm_01 --range_EbNo --range_EbNo_start 6 --plot_loss
```
1. Inference (runs encoder and decoder, and outputs symbols `z_hat.f32`):
```
./inference_bbfm.sh model_bbfm_01/checkpoints/checkpoint_epoch_100.pth wav/brian_g8sez.wav - --write_latent z_hat.f32
```
1. Stand alone decoder, outputs speech from `z_hat.f32` to sound card:
```
./rx_bbfm.sh model_bbfm_01/checkpoints/checkpoint_epoch_100.pth z_hat.f32 -
```
1. Or save speech out to a wave file:
```
./rx_bbfm.sh model_bbfm_01/checkpoints/checkpoint_epoch_100.pth z_hat.f32 t.wav
```
1. Plot sequence of received symbols:
```
octave:4> radae_plots; do_plots_bbfm('z_hat.f32')
```
# Fading channel simulation
HF channel sim (two path Rayleigh) is pretty close to TIA-102.CAAA-E 1.6.33 Faded Channel Simulator. The measured level crossing rate (LCR) seems to meet req (f), for v=60 km/hr, f = 450 MHz, and P=1 when measured over a 10 second sample. We've used Rs=2000 symb/s here, so x-axis of plot is 1 second in time.
![LMR 60](doc/lmr_60.png)
```
octave:39> multipath_samples("lmr60",8000, 2000, 1, 10, "h_lmr60.f32")
Generating Doppler spreading samples...
fd = 25.000
path_delay_s = 2.0000e-04
Nsecplot = 1
Pav = 1.0366
P = 1
LCR_theory = 23.457
LCR_meas = 24.400
```
# Single Carrier PSK Modem
A single carrier PSK modem "back end" that connects the ML symbols to the radio. This particular modem is written in Python, and can work with DC coupled and passband BBFM radios. It uses classical DSP, rather than ML. Unlike the HF RADE waveform which used OFDM, this modem is single carrier.
1. Run a single test with some plots, Eb/No=4dB, 100ppm sample clock offset, BER should be about 0.01:
```
python3 -c "from radae import single_carrier; s=single_carrier(); s.run_test(100,sample_clock_offset_ppm=-100,plots_en=True,EbNodB=4)"
```
1. Run a suite of tests:
```
ctest -V -R bbfm_sc
```
1. Create a file of BBFM symbols, 80 symbols every 40ms, plays expected output speech:
```
./bbfm_inference.sh model_bbfm_01/checkpoints/checkpoint_epoch_100.pth wav/brian_g8sez.wav - --write_latent z.f32
```
2. Sanity check of modem, BER test using digital, BPSK symbols, the symbols in z.f32 are replaced with BPSK symbols. `t.int16` is a real valued Fs=9600Hz sample file, that could be played into a FM radio.
```
cat z.f32 | python3 sc_tx.py --ber_test > t.int16
cat t.int16 | python3 sc_rx.py --ber_test --plots > /dev/null
```
3. Send the BBFM symbols over the modem, and listen to results:
```
cat z.f32 | python3 sc_tx.py > t.int16
cat t.int16 | python3 sc_rx.py > z_hat.f32
./bbfm_rx.sh model_bbfm_01/checkpoints/checkpoint_epoch_100.pth z_hat.f32 -
```
4. Compare MSE of features passed through the system, first with z == z_hat, then with z passed through modem to get z_hat:
```
python3 loss.py features_in.f32 features_out.f32
loss: 0.033
python3 loss.py features_in.f32 features_rx_out.f32
loss: 0.035
```
This is a really good result, and likely inaudible. The `feature*.f32` files are produced as intermediate outputs form the `bbfm_inference.sh` and `bbfm_rx.sh` scripts.
42 changes: 40 additions & 2 deletions CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -168,13 +168,13 @@ add_test(NAME chirp_mpp
# Low SNR ota_test.sh, with chirp measurement, AWGN
add_test(NAME ota_test_awgn
COMMAND sh -c "cd ${CMAKE_SOURCE_DIR}; \
./test/ota_test_cal.sh ~/codec2-dev/build_linux/ -21 0.4")
./test/ota_test_cal.sh ${CODEC2_DEV_BUILD_DIR} -21 0.4")

# Low SNR ota_test.sh, with chirp measurement, MPP, high loss threshold as we only care about gross errors,
# like stuck in false sync
add_test(NAME ota_test_mpp
COMMAND sh -c "cd ${CMAKE_SOURCE_DIR}; \
./test/ota_test_cal.sh ~/codec2-dev/build_linux/ -24 0.4 --mpp --freq -25")
./test/ota_test_cal.sh ${CODEC2_DEV_BUILD_DIR} -24 0.4 --mpp --freq -25")


# Acquisition tests ------------------------------------------------------------------------------------
Expand Down Expand Up @@ -432,6 +432,7 @@ add_test(NAME radae_rx_aux_mpp
python3 loss.py features_in.f32 features_rx_out.f32 --loss 0.3 --clip_start 300")
set_tests_properties(radae_rx_aux_mpp PROPERTIES PASS_REGULAR_EXPRESSION "PASS")


# Embedding Python in C callable library ------------------------------------------------------------------------------------------------------

# Test Embedded version of streaming Tx, running with Python __main__
Expand Down Expand Up @@ -471,3 +472,40 @@ add_test(NAME radae_rx_embed_c
cat rx.f32 | PYTHONPATH='.' ${CMAKE_CURRENT_BINARY_DIR}/src/radae_rx > features_out.f32;
python3 loss.py features_in.f32 features_out.f32 --loss_test 0.15 --acq_time_test 0.5")
set_tests_properties(radae_rx_embed_c PROPERTIES PASS_REGULAR_EXPRESSION "PASS")


# BBFM -----------------------------------------------------------------------------------------------

# single carrier modem internal (inside single_carrier class) tests
add_test(NAME bbfm_sc_internal
COMMAND sh -c "cd ${CMAKE_SOURCE_DIR}; \
python3 -c 'from radae import single_carrier,single_carrier_tests; single_carrier_tests()'")
set_tests_properties(bbfm_sc_internal PROPERTIES PASS_REGULAR_EXPRESSION "ALL PASS")

# single carrier modem stand alone tx/rx, using BPSK symbol/BER test mode
add_test(NAME bbfm_sc_ber
COMMAND sh -c "cd ${CMAKE_SOURCE_DIR}; \
./bbfm_inference.sh model_bbfm_01/checkpoints/checkpoint_epoch_100.pth wav/brian_g8sez.wav /dev/null --write_latent z.f32; \
cat z.f32 | python3 sc_tx.py --ber_test > t.int16; \
cat t.int16 | python3 sc_rx.py --ber_test --target_ber 0 > /dev/null")
set_tests_properties(bbfm_sc_ber PROPERTIES PASS_REGULAR_EXPRESSION "PASS")

# single carrier modem stand alone tx/rx, using continuously valued BBFM symbols, compare loss in decoded features
add_test(NAME bbfm_sc_loss
COMMAND sh -c "cd ${CMAKE_SOURCE_DIR}; \
./bbfm_inference.sh model_bbfm_01/checkpoints/checkpoint_epoch_100.pth wav/brian_g8sez.wav /dev/null --write_latent z.f32; \
cat z.f32 | python3 sc_tx.py | python3 sc_rx.py > z_hat.f32; \
./bbfm_rx.sh model_bbfm_01/checkpoints/checkpoint_epoch_100.pth z_hat.f32 /dev/null; \
python3 loss.py features_in.f32 features_out.f32 --features_hat2 features_rx_out.f32 --compare")
set_tests_properties(bbfm_sc_loss PROPERTIES PASS_REGULAR_EXPRESSION "PASS")

# single carrier modem stand alone tx/rx with external 300-2700Hz band pass filter and no noise, measure loss. In practice, SNR will be
# quite high, so channel distortions other than noise may dominate
add_test(NAME bbfm_sc_bpf_loss
COMMAND sh -c "cd ${CMAKE_SOURCE_DIR}; \
./bbfm_inference.sh model_bbfm_01/checkpoints/checkpoint_epoch_100.pth wav/brian_g8sez.wav /dev/null --write_latent z.f32; \
cat z.f32 | python3 sc_tx.py | ${CODEC2_DEV_BUILD_DIR}/src/ch - - | python3 sc_rx.py > z_hat.f32; \
./bbfm_rx.sh model_bbfm_01/checkpoints/checkpoint_epoch_100.pth z_hat.f32 /dev/null; \
python3 loss.py features_in.f32 features_out.f32 --features_hat2 features_rx_out.f32 --compare")
set_tests_properties(bbfm_sc_bpf_loss PROPERTIES PASS_REGULAR_EXPRESSION "PASS")

7 changes: 7 additions & 0 deletions README.md
Original file line number Diff line number Diff line change
Expand Up @@ -52,6 +52,13 @@ The RDOVAE derived Python source code is released under the two-clause BSD licen
| radae_tx.[py,sh] | streaming RADAE encoder and helper script |
| radae_rx.[py,sh] | streaming RADAE decoder and helper script |
| resource_est.py | WIP estimate CPU/memory resources |
| radae_base.py | Shared ML code between models |
| radae/bbfm.py | Baseband FM PyTorch model |
| train_bbfm.py | Training for BBFM model |
| inference_bbfm.py | Baseband FM inference |
| inference_bbfm.sh | helper script for infereence_bbfm.sh |
| fm.m | Octave analog FM mod/demod simulation |
| analog_bbfm.sh | helper script for analog FM simulation |

# Installation

Expand Down
44 changes: 44 additions & 0 deletions analog_bbfm.sh
Original file line number Diff line number Diff line change
@@ -0,0 +1,44 @@
#!/bin/bash -x
#
# Analog FM simulation, for comparison to ML BBFM

CODEC2_DEV=${CODEC2_DEV:-${HOME}/codec2-dev}
OPUS=build/src
PATH=${PATH}:${OPUS}:${CODEC2_DEV}/build_linux/src
gain=6

which ch >/dev/null || { printf "\n**** Can't find ch - check CODEC2_PATH **** \n\n"; exit 1; }

source utils.sh

if [ $# -lt 3 ]; then
echo "usage (write output to file):"
echo " ./analog_bbfm.sh in.wav out.wav CNRdB"
echo "usage (play output with aplay):"
echo " ./analog_bbfm.sh in.wav - CNRdB"
exit 1
fi

if [ ! -f $1 ]; then
echo "can't find $1"
exit 1
fi

input_speech=$1
output_speech=$2
CNRdB=$3

tmp_in=$(mktemp)
tmp_out=$(mktemp)
tmp_fm=$(mktemp)

# We use hilbert clipper in ch util for speech compressor. Octave FM simulation uses 48 kHz sample rate.
# input wav -> 300-3100Hz Fs=8kHz -> ch compressor -> 300-3100Hz Fs=48kHz -> FM mod/demod
sox ${input_speech} -t .s16 -r 8000 -c 1 - sinc 0.3-3.1k | ch - - --clip 16384 --gain $gain 2>/dev/null | sox -t .s16 -r 8000 -c 1 - -t .s16 -r 48000 ${tmp_in} sinc 0.3-3.1k
echo "fm; pkg load signal; fm_mod_file('${tmp_fm}','${tmp_in}',${CNRdB}); fm_demod_file('${tmp_out}','${tmp_fm}'); quit;" | octave-cli -qf

if [ $output_speech == "-" ]; then
aplay ${tmp_out} -r 48000 -f S16_LE 2>/dev/null
elif [ $output_speech != "/dev/null" ]; then
sox -t .s16 -r 48000 -c 1 ${tmp_out} -r 8000 ${output_speech}
fi
30 changes: 30 additions & 0 deletions bbfm_bpf.m
Original file line number Diff line number Diff line change
@@ -0,0 +1,30 @@
% Octave script to explore BPF of basband symbols

function bbfm_bpf()
pkg load signal;
Fs = 8000; T = 1/Fs; Rs = 2000; M = Fs/Rs; Nsym = 6; alpha = 0.25;
rn = gen_rn_coeffs(alpha, T, Rs, Nsym, M);
bpf = fir2(100,[0 250 350 3000 3100 4000]/(Fs/2),[0.001 0.001 1 1 0.001 0.001]);

figure(1); clf;
subplot(211);
[h,w] = freqz(rn); plot(w*Fs/(2*pi), 20*log10(abs(h))); grid('minor'); ylabel('RRC');
subplot(212);
[h,w] = freqz(bpf); plot(w*Fs/(2*pi), 20*log10(abs(h))); grid; ylabel('BPF');

Nsymb = 1000;
tx_symb = 1 - 2*(rand(Nsymb,1)>0.5);
tx_pad = zeros(1,M*Nsymb);
tx_pad(1:M:end) = tx_symb;
tx = filter(rn,1,tx_pad);
tx = filter(bpf,1,tx)
rx = filter(rn,1,tx);
rx_symb = rx(1:M:end);
figure(2); clf;
subplot(211); stem(tx_symb(1:100)); ylabel('Tx symbols');
subplot(212); stem(rx_symb(1:100)); ylabel('Rx Symbols');

figure(3); clf;
plot(20*log10(abs(fft(tx)(1:length(tx)/2))))
end

142 changes: 142 additions & 0 deletions bbfm_inference.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,142 @@
"""
/* Copyright (c) 2024 modifications for radio autoencoder project
by David Rowe */
/* Copyright (c) 2022 Amazon
Written by Jan Buethe */
/*
Redistribution and use in source and binary forms, with or without
modification, are permitted provided that the following conditions
are met:
- Redistributions of source code must retain the above copyright
notice, this list of conditions and the following disclaimer.
- Redistributions in binary form must reproduce the above copyright
notice, this list of conditions and the following disclaimer in the
documentation and/or other materials provided with the distribution.
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER
OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
*/
"""

import os
import argparse

import numpy as np
import torch

from radae import BBFM, distortion_loss

parser = argparse.ArgumentParser()

parser.add_argument('model_name', type=str, help='path to model in .pth format')
parser.add_argument('features', type=str, help='path to input feature file in .f32 format')
parser.add_argument('features_hat', type=str, help='path to output feature file in .f32 format')
parser.add_argument('--latent-dim', type=int, help="number of symbols produces by encoder, default: 80", default=80)
parser.add_argument('--cuda-visible-devices', type=str, help="set to 0 to run using GPU rather than CPU", default="")
parser.add_argument('--write_latent', type=str, default="", help='path to output file of latent vectors z[latent_dim] in .f32 format')
parser.add_argument('--CNRdB', type=float, default=100, help='FM demod input CNR in dB')
parser.add_argument('--passthru', action='store_true', help='copy features in to feature out, bypassing ML network')
parser.add_argument('--h_file', type=str, default="", help='path to rate Rs fading channel magnitude samples, rate Rs time steps by Nc=1 carriers .f32 format')
parser.add_argument('--write_CNRdB', type=str, default="", help='path to output file of CNRdB per sample after fading in .f32 format')
parser.add_argument('--loss_test', type=float, default=0.0, help='compare loss to arg, print PASS/FAIL')
args = parser.parse_args()

# set visible devices
os.environ['CUDA_VISIBLE_DEVICES'] = args.cuda_visible_devices

# device
device = torch.device("cuda") if torch.cuda.is_available() else torch.device("cpu")

latent_dim = args.latent_dim

# not exposed
nb_total_features = 36
num_features = 20
num_used_features = 20

# load model from a checkpoint file
model = BBFM(num_features, latent_dim, args.CNRdB)
checkpoint = torch.load(args.model_name, map_location='cpu')
model.load_state_dict(checkpoint['state_dict'], strict=False)
checkpoint['state_dict'] = model.state_dict()

# load features from file
feature_file = args.features
features_in = np.reshape(np.fromfile(feature_file, dtype=np.float32), (1, -1, nb_total_features))
nb_features_rounded = model.num_10ms_times_steps_rounded_to_modem_frames(features_in.shape[1])
features = features_in[:,:nb_features_rounded,:]
features = features[:, :, :num_used_features]
features = torch.tensor(features)
print(f"Processing: {nb_features_rounded} feature vectors")

# default rate Rb multipath model H=1
Rb = model.Rb
Nc = 1
num_timesteps_at_rate_Rs = model.num_timesteps_at_rate_Rs(nb_features_rounded)
H = torch.ones((1,num_timesteps_at_rate_Rs,Nc))

# user supplied rate Rs multipath model, sequence of H magnitude samples
if args.h_file:
H = np.reshape(np.fromfile(args.h_file, dtype=np.float32), (1, -1, Nc))
print(H.shape, num_timesteps_at_rate_Rs)
if H.shape[1] < num_timesteps_at_rate_Rs:
print("Multipath H file too short")
quit()
H = H[:,:num_timesteps_at_rate_Rs,:]
H = torch.tensor(H)

if __name__ == '__main__':

if args.passthru:
features_hat = features_in.flatten()
features_hat.tofile(args.features_hat)
quit()

# push model to device and run test
model.to(device)
features = features.to(device)
H = H.to(device)
output = model(features,H)

# Lets check actual SNR at output of FM demod
tx_sym = output["z_hat"].cpu().detach().numpy()
S = np.mean(np.abs(tx_sym)**2)
N = np.mean(output["sigma"].cpu().detach().numpy()**2)
SNRdB_meas = 10*np.log10(S/N)
print(f"SNRdB Measured: {SNRdB_meas:6.2f}")

features_hat = output["features_hat"][:,:,:num_used_features]
features_hat = torch.cat([features_hat, torch.zeros_like(features_hat)[:,:,:16]], dim=-1)
features_hat = features_hat.cpu().detach().numpy().flatten().astype('float32')
features_hat.tofile(args.features_hat)

loss = distortion_loss(features,output['features_hat']).cpu().detach().numpy()[0]
print(f"loss: {loss:5.3f}")
if args.loss_test > 0.0:
if loss < args.loss_test:
print("PASS")
else:
print("FAIL")

# write output symbols (latent vectors)
if len(args.write_latent):
z_hat = output["z_hat"].cpu().detach().numpy().flatten().astype('float32')
z_hat.tofile(args.write_latent)

# write CNRdB after fading
if len(args.write_CNRdB):
CNRdB = output["CNRdB"].cpu().detach().numpy().flatten().astype('float32')
CNRdB.tofile(args.write_CNRdB)

Loading

0 comments on commit 490dfcc

Please sign in to comment.