Carrier frequency offset in Bell 103 FSK system

Question

This question is a sequel of the discussion that started in this question. Dan Boschen suggested some very nice solutions to decode a Bell 103 FSK signal. I am especially interested in decoding the 300bps signals with a C program running on a microcontroller. Unfortunately, the real input shows some deviations I cannot yet explain, which complicate everything. Andy Walls has indicated that it is Carrier Frequency Offset (CFO). Hence this post.

The real input signal

signal = [-2654, -7612, -3886, 7878, 11615, 3576, -7827, -8718, 2024, 10608, 7037, -4689, -10767, -3449, 7954, 9506, -268, -10213, -7943, 3860, 10379, 4329, -7315, -10597, -993, 9291, 8319, -1073, -8453, -9421, -3838, 4587, 9746, 8721, 2079, -6025, -10183, -7285, 543, 7892, 9928, 5577, -2463, -9105, -9553, -3460, 4873, 9778, 8289, 1465, -6545, -10267, -6913, 1409, 9233, 9670, 1549, -7870, -10014, -2476, 7280, 10301, 4253, -5740, -10644, -5449, 4735, 10329, 6648, -3054, -10283, -7899, 1853, 9598, 8580, -116, -9005, -9653, -1175, 8088, 9893, 4392, -3651, -9721, -9107, -2347, 5821, 9970, 7589, 329, -7404, -10201, -5967, 2302, 8859, 9547, 4147, -4083, -9749, -8730, -1712, 6344, 10047, 7267, -227, -7833, -10173, -4173, 5971, 10306, 5370, -4504, -10636, -6707, 3355, 10120, 7727, -1563, -9704, -8781, 437, 9021, 9372, 1417, -8025, -10149, -2603, 7177, 10249, 4220, -5744, -10668, -5861, 3143, 9068, 8965, 2886, -5287, -10183, -8154, -627, 7100, 9972, 6363, -1451, -8575, -9938, -4485, 3901, 9493, 8821, 2499, -5667, -10182, -7705, 6, 7571, 9986, 5825, -2953, -10279, -8524, 1050, 9229, 8920, 568, -8609, -9954, -1889, 7670, 10083, 3544, -6378, -10624, -4769, 5398, 10408, 6134, -3699, -10422, -7309, 2637, 9894, 8222, -642, -8674, -9449, -3731, 4583, 9831, 8655, 2073, -6095, -10210, -7331, 538, 7869, 9865, 5497, -2561, -9185, -9618, -3499, 4834, 9739, 8214, 1379, -6616, -10322, -6953, 1080, 8409, 9903, 3054, -6925, -10372, -4007, 6091, 10442, 5543, -4415, -10569, -6726, 3321, 10072, 7715, -1583, -9742, -8847, 336, 8928, 9332, 1383, -8081, -10213, -2681, 7110, 10289, 4978, -3519, -9480, -9000, -2202, 5933, 10045, 7614, 315, -7409, -10204, -5942, 2302, 8857, 9509, 4089, -4137, -9793, -8746, -1724, 6328, 10021, 7222, -275, -7838, -10123, -4857, 4825, 10468, 6696, -2987, -10231, -7801, 1940, 9661, 8653, -69, -8957, -9569, -1109, 8162, 9911, 2851, -6971, -10531, -4089, 5999, 10393, 5513, -4426, -10578, -6787, 2807, 9047, 8738, 2828, -5427, -10221, -8193, -622, 7121, 9974, 6362, -1471, -8559, -9908, -4432, 3953, 9537, 8864, 2520, -5641, -10149, -7648, 56, 7589, 9994, 5989, -2088, -9526, -9544, -573, 8418, 9480, 2021, -7641, -10408, -3409, 6581, 10342, 4891, -5077, -10617, -6097, 4070, 10274, 7276, -2246, -9966, -8345, 1145, 9352, 9059, 713, -8506, -10084, -3341, 4988, 9829, 8253, 1308, -6650, -10324, -6863, 1158, 8236, 9730, 4961, -3183, -9470, -9348, -2856, 5407, 9862, 7865, 757, -7108, -10304, -6431, 1743, 8569, 9581, 3392, -6440, -10648, -4769, 5371, 10406, 6074, -3750, -10452, -7322, 2604, 9869, 8181, -838, -9383, -9254, -389, 8576, 9648, 2108, -7540, -10389, -3403, 6571, 10352, 5141, -3853, -9537, -8504, -1471, 6542, 10119, 7243, -323, -7846, -10123, -5404, 2918, 9108, 9246, 3457, -4768, -10013, -8423, -1103, 6785, 10003, 6737, -959, -8269, -10056, -4821, 4225, 10333, 7837, -1482, -9649];

The signal consists of a repeating 101010 pattern, with a total of 427 samples. The samples are obtained at PCM Fs = 8 kHz, which has been verified by use of a logic analyzer. The resulting frequency spectrum is shown below. It deviates from the expected 1070 and 1270 Hz. In particular, the big lobe at 1170 Hz exactly in between the two expected frequencies, and the two smaller side lobes.

Unfortunately I do not have an oscilloscope available to verify the 1070 and 1270 Hz signals are actually there. I am hoping to have one tomorrow, or Monday. However, the hardware generating the real signal is Bell 103 compliant and has not had problems in previous tests (although those were some time ago).

Answer 1

You actually haven't written a question here.

However, processing the signal samples you provided, I don't see CFO

I do see that the first symbol is way off frequency.

The Octave code that generated this picture is

% Octave requires an explcit load of the signal processing packasge, once.
pkg load signal

Fs = 8000; % sample rate
f1 = 1070;
f2 = 1270;

% Given input signal
x = [-2654, -7612, -3886, 7878, 11615, 3576, -7827, -8718, 2024, 10608, 7037, -4689, -10767, -3449, 7954, 9506, -268, -10213, -7943, 3860, 10379, 4329, -7315, -10597, -993, 9291, 8319, -1073, -8453, -9421, -3838, 4587, 9746, 8721, 2079, -6025, -10183, -7285, 543, 7892, 9928, 5577, -2463, -9105, -9553, -3460, 4873, 9778, 8289, 1465, -6545, -10267, -6913, 1409, 9233, 9670, 1549, -7870, -10014, -2476, 7280, 10301, 4253, -5740, -10644, -5449, 4735, 10329, 6648, -3054, -10283, -7899, 1853, 9598, 8580, -116, -9005, -9653, -1175, 8088, 9893, 4392, -3651, -9721, -9107, -2347, 5821, 9970, 7589, 329, -7404, -10201, -5967, 2302, 8859, 9547, 4147, -4083, -9749, -8730, -1712, 6344, 10047, 7267, -227, -7833, -10173, -4173, 5971, 10306, 5370, -4504, -10636, -6707, 3355, 10120, 7727, -1563, -9704, -8781, 437, 9021, 9372, 1417, -8025, -10149, -2603, 7177, 10249, 4220, -5744, -10668, -5861, 3143, 9068, 8965, 2886, -5287, -10183, -8154, -627, 7100, 9972, 6363, -1451, -8575, -9938, -4485, 3901, 9493, 8821, 2499, -5667, -10182, -7705, 6, 7571, 9986, 5825, -2953, -10279, -8524, 1050, 9229, 8920, 568, -8609, -9954, -1889, 7670, 10083, 3544, -6378, -10624, -4769, 5398, 10408, 6134, -3699, -10422, -7309, 2637, 9894, 8222, -642, -8674, -9449, -3731, 4583, 9831, 8655, 2073, -6095, -10210, -7331, 538, 7869, 9865, 5497, -2561, -9185, -9618, -3499, 4834, 9739, 8214, 1379, -6616, -10322, -6953, 1080, 8409, 9903, 3054, -6925, -10372, -4007, 6091, 10442, 5543, -4415, -10569, -6726, 3321, 10072, 7715, -1583, -9742, -8847, 336, 8928, 9332, 1383, -8081, -10213, -2681, 7110, 10289, 4978, -3519, -9480, -9000, -2202, 5933, 10045, 7614, 315, -7409, -10204, -5942, 2302, 8857, 9509, 4089, -4137, -9793, -8746, -1724, 6328, 10021, 7222, -275, -7838, -10123, -4857, 4825, 10468, 6696, -2987, -10231, -7801, 1940, 9661, 8653, -69, -8957, -9569, -1109, 8162, 9911, 2851, -6971, -10531, -4089, 5999, 10393, 5513, -4426, -10578, -6787, 2807, 9047, 8738, 2828, -5427, -10221, -8193, -622, 7121, 9974, 6362, -1471, -8559, -9908, -4432, 3953, 9537, 8864, 2520, -5641, -10149, -7648, 56, 7589, 9994, 5989, -2088, -9526, -9544, -573, 8418, 9480, 2021, -7641, -10408, -3409, 6581, 10342, 4891, -5077, -10617, -6097, 4070, 10274, 7276, -2246, -9966, -8345, 1145, 9352, 9059, 713, -8506, -10084, -3341, 4988, 9829, 8253, 1308, -6650, -10324, -6863, 1158, 8236, 9730, 4961, -3183, -9470, -9348, -2856, 5407, 9862, 7865, 757, -7108, -10304, -6431, 1743, 8569, 9581, 3392, -6440, -10648, -4769, 5371, 10406, 6074, -3750, -10452, -7322, 2604, 9869, 8181, -838, -9383, -9254, -389, 8576, 9648, 2108, -7540, -10389, -3403, 6571, 10352, 5141, -3853, -9537, -8504, -1471, 6542, 10119, 7243, -323, -7846, -10123, -5404, 2918, 9108, 9246, 3457, -4768, -10013, -8423, -1103, 6785, 10003, 6737, -959, -8269, -10056, -4821, 4225, 10333, 7837, -1482, -9649];

% Scale to between +/- 1.0 to avoid ridiculously large numbers on plots
xs = x/max(abs(x));

% Run through a Hilbert filter to make the signal analytic
xsc = hilbert(xs);

% FM demodulate using an approximate d(phi)/dt
y = arg(xsc(2:end) .* conj(xsc(1:(end-1))));
ys = y * (Fs/2)/pi;

% Plot the output of the FM discrimator
t = [1:length(ys)]/Fs;
f1line = ones(1, length(t))*f1;
f2line = ones(1, length(t))*f2;
plot(t, ys, t, f1line, t, f2line);
ylim([f1-300 f2+300]);
title('FM Discriminator Output');
xlabel('Time (seconds)');
ylabel('Frequency (Hz)');

By the way, the Octave/MatLab 'hilbert' implementation does a terrible thing - it zeros FFT bins to perform filtering - which should never be done in a real life, streaming Hilbert filter implementation. Oppenheim and Schafer have a section in their textbook regarding practical Hilbert filters.