Liquid-dsp based command line FM demodulator

In my last post I wrote about how I hacked Google radioreceiver into command line based FM demodulator that can be used to pipe demodulated audio into multimon-ng for AX.25 packet decoding. However it had some significant flaws. For example as far as I know it only worked reliably if input sample rate was 1024 kHz and output 48 kHz. I did not look too much into problem because I wanted to rewrite it anyway to make it more universal for supporting more modulations beside FM and AM.

I am follower of www.rtl-sdr.com and recently there was a post about new SDR software called CubicSDR. What makes it interesting for me is that it uses library called liquid-dsp for doing signal processing. I tried to find decent signal processing library that is intended for SDR, but I ended up forking Google radioreceiver out of frustation. I did not have any luck finding such library until I stumbled upon liquid-dsp.

Liquid-dsp is written in pure C and it does not have any dependencies which make it very suitable for embedded platforms. GNU Radio is very, very bloated compared to this. I am going to use it on BeagleBoneBlack or even better on the new 4 core RaspberryPi, therefore liquid-dsp embedded oriented nature fits very well for me.

FM demodulation #

Here is hypothetical demodulation loop where block of complex baseband signal is read from “source” and demodulated output is written to stdout so it can be piped to another program. Notice how to setup modulation factor kf for FSK demodulation with deviation of 3500 Hz.

#define SAMPLERATE 1024000.0f
#define FSK_DEVIATION_HZ 3500.0f
#define NUM_SAMPLES 8192

float kf = FSK_DEVIATION_HZ/SAMPLERATE; // modulation factor
float complex r[NUM_SAMPLES]; // received signal
float         y[NUM_SAMPLES]; // demodulator output

freqmod mod = freqmod_create(kf);

while (1) {
    read_complex_iq_block(r, NUM_SAMPLES); // made up function
    freqdem_demodulate_block(dem, r, NUM_SAMPLES, y);
    fwrite((uint8_t*)y, 4, NUM_SAMPLES, stdout);
}
freqmod_destroy(mod);

Lowpass filtering #

Demodulation loop has a little problem. Its input baseband signal sample rate is 1024 kHz, however we would like to demodulate FSK modulated signal which has deviation of ±3.5 kHz and data rate of 9.6 kbps. The problem is that currently unfiltered signal goes to demodulator, therefore it is unable to figure out our narrow band signal correctly, because of high frequency noise which is present in baseband signal. Lowpass filter must be added to cut off out of band high frequency signals.

#define SAMPLERATE 1024000.0f
#define CUTOFF_HZ 4200.0f
#define NUM_SAMPLES 8192

float fc = CUTOFF_HZ/SAMPLERATE; // cutoff frequency
unsigned int h_len = 64; // filter length    
float As = 70.0f; // stop-band attenuation

float complex r[NUM_SAMPLES]; // received signal
float complex rf[NUM_SAMPLES]; // filtered signal

// design filter from prototype and scale to bandwidth
firfilt_crcf q = firfilt_crcf_create_kaiser(h_len, fc, As, 0.0f);
firfilt_crcf_set_scale(q, 2.0f*fc);

while (1) {
    read_complex_iq_block(r, NUM_SAMPLES); // made up function
    for (int i=0; i<NUM_SAMPLES; i++)
        firfilt_crcf_push(q, r[i]);
        firfilt_crcf_execute(q, &rf[i]);
    }

    // here lowpass filtered signal "rf" can be given to demodulator as shown above
}
firfilt_crcf_destroy(q);

I have not yet studied how this filter response actually looks like. I just used it with cube satellite ESTCube-1 recording and those filter parameters worked very well. Notice that here fc defines filter cutoff frequency centered around 0 Hz. In this example it is ±4.2 kHz.
Here is screenshot from Baudline showing AX.25 packets after filtering step:

Resampling #

So far so good, filtering and demodulation works. Now only resampling is missing which was the main problem with Google radioreceiver. Fortunately arbitrary resampling is supported in liquid-dsp. Here is example how to use it:

#define IN_SAMPLERATE 1024000.0f
#define OUT_SAMPLERATE 48000.0f
#define NUM_SAMPLES 8192

float r= OUT_SAMPLERATE/IN_SAMPLERATE;
float As=60.0f; // resampling filter stop-band attenuation [dB]
unsigned int n= NUM_SAMPLES; // number of input samples

msresamp_crcf q = msresamp_crcf_create(r,As);
msresamp_crcf_print(q);
float delay = msresamp_crcf_get_delay(q);

// number of input samples (zero-padded)
unsigned int nx = n + (int)ceilf(delay) + 10;

// output buffer with extra padding for good measure
unsigned int ny_alloc = (unsigned int) (2*(float)nx * r);  // allocation for output

float complex r[nx]; // received signal
float complex y[ny_alloc]; // resampled signal

while (1) {
    unsigned int ny;
    read_complex_iq_block(r, NUM_SAMPLES); // made up function
    msresamp_crcf_execute(q, r, NUM_SAMPLES, y, &ny);

    // here "y" is signal with sample rate of OUT_SAMPLERATE and "ny" shows how many samples are in "y"
}
msresamp_crcf_destroy(q);

Original example uses nx in msresamp_crcf_execute, however it did not work for me. I got it to work with NUM_SAMPLES.

msresamp_crcf_execute(q, r, nx, y, &ny);

vs

msresamp_crcf_execute(q, r, NUM_SAMPLES, y, &ny);

Here is output resampled from 126 kHz to 48 kHz:

Putting it all together #

Final tool called demod can be downloaded from github.

Here sdr_fsk9600.wav is recording that contains some AX.25 packets. Those packets can be demodulated and decoded using following command:

sox -t wav sdr_fsk9600.wav -esigned-integer -b16  -r 126000 -t raw - | demod -s 126000 -r 48000 -b 4200 -m fm d=3500 | multimon-ng -t raw  -a FSK9600 /dev/stdin

Here multimon-ng fork that processes 48 kHz input instead of 22.050 kHz is used. It is because multion-ng decodes much more packets with 48 kHz input. Firstly I thought that Google radioreceiver 22.050 kHz output was not good enough for multimon-ng, however now even with liquid-dsp based approach it works better using 48 kHz sampling rate. I should really dig into multimon-ng to see why it works so badly.