From 8a45fdf4d72ebfca4e52b8f4d629f60f9861633b Mon Sep 17 00:00:00 2001 From: drowe67 Date: Thu, 5 Apr 2012 09:47:54 +0000 Subject: [PATCH] modified freq offset estimation to fix bug that occurs with 1 symbol delay. Improves freq est quite a lot git-svn-id: https://svn.code.sf.net/p/freetel/code@363 01035d8c-6547-0410-b346-abe4f91aad63 --- codec2-dev/octave/README_fdmdv.txt | 16 ++++++- codec2-dev/octave/fdmdv.m | 72 +++++++++++++++++------------- codec2-dev/octave/fdmdv_demod.m | 44 +++++++++++++++--- codec2-dev/octave/fdmdv_ut.m | 54 +++++++++++++++++----- 4 files changed, 137 insertions(+), 49 deletions(-) diff --git a/codec2-dev/octave/README_fdmdv.txt b/codec2-dev/octave/README_fdmdv.txt index ecf22445..d6f415ce 100644 --- a/codec2-dev/octave/README_fdmdv.txt +++ b/codec2-dev/octave/README_fdmdv.txt @@ -6,11 +6,23 @@ sox -r 8000 -s -2 mod_dqpsk.raw -s -2 mod_dqpsk_8008hz.raw rate -h 8008 TODO -[ ] Get file based mod and demod working again -[ ] write file based ch simulator +[X] Get file based mod and demod working again +[ ] timing wraps around + + what is affect of bouncing back and forth over boundary? + + could mean we are going back and forth a symbol + + do we need to logic to lose or gain a frame? + + think so, add or remove samples, or a whole frame [ ] demod outputs ber (maybe after settling time) [ ] try to run from shell script [ ] run a few tests [ ] start coding in C and repeat tests [ ] arb bit stream input to Octave mod and demod +Tests + +[ ] fdmdv_demod('mod_dqpsk_8008hz.raw',1400*60); +[ ] fdmdv_demod('mod_dqpsk_7992hz.raw',1400*60); +[ ] mod_dqpsk_awgn_4dB_8008hz.raw +[ ] mod_dqpsk_good_4dB_8008hz.raw +[ ] mod_dqpsk_moderate_4dB_8008hz.raw +[ ] mod_dqpsk_moderate_4dB_7992hz.raw diff --git a/codec2-dev/octave/fdmdv.m b/codec2-dev/octave/fdmdv.m index af8b8ee1..35017d5d 100644 --- a/codec2-dev/octave/fdmdv.m +++ b/codec2-dev/octave/fdmdv.m @@ -1,7 +1,7 @@ % fdmdv.m % % Functions that implement a Frequency Divison Multiplexed Modem for -% Digital Voice (FMDV)over HF channels. +% Digital Voice (FDMDV) over HF channels. % % Copyright David Rowe 2012 % This program is distributed under the terms of the GNU General Public License @@ -284,30 +284,45 @@ endfunction % Estimate frequency offset of FDM signal using BPSK pilot. This is quite % sensitive to pilot tone level wrt other carriers -function foff = rx_est_freq_offset(rx_fdm, pilot) +function foff = rx_est_freq_offset(rx_fdm, pilot, pilot_prev) global M; - global Nc; - global Fs; - global Rs; - global Fcentre; - global freq; - global freq_rx_pilot; global Npilotbaseband; - global Npilotlpf; - global Npilotcoeff; - global pilot_baseband; - global pilot_lpf; - global pilot_coeff; + global pilot_baseband1; + global pilot_baseband2; + global pilot_lpf1; + global pilot_lpf2; % down convert latest M samples of pilot by multiplying by - % ideal two-tone BPSK pilot signal + % ideal BPSK pilot signal we have generated locally - pilot_baseband(1:Npilotbaseband-M) = pilot_baseband(M+1:Npilotbaseband); - c = Nc+1; + pilot_baseband1(1:Npilotbaseband-M) = pilot_baseband1(M+1:Npilotbaseband); + pilot_baseband2(1:Npilotbaseband-M) = pilot_baseband2(M+1:Npilotbaseband); for i=1:M - pilot_baseband(Npilotbaseband-M+i) = rx_fdm(i)*conj(pilot(i)); + pilot_baseband1(Npilotbaseband-M+i) = rx_fdm(i) * conj(pilot(i)); + pilot_baseband2(Npilotbaseband-M+i) = rx_fdm(i) * conj(pilot_prev(i)); end + [foff1 max1 pilot_lpf1] = lpf_peak_pick(pilot_baseband1, pilot_lpf1); + [foff2 max2 pilot_lpf2] = lpf_peak_pick(pilot_baseband2, pilot_lpf2); + + if max1 > max2 + foff = foff1; + else + foff = foff2; + end +endfunction + + +% LPF and peak pick part of freq est, put in a function as we call it twice + +function [foff imax pilot_lpf] = lpf_peak_pick(pilot_baseband, pilot_lpf) + global M; + global Npilotlpf; + global Npilotcoeff; + global Fs; + global Mpilotfft; + global pilot_coeff; + % LPF cutoff 200Hz, so we can handle max +/- 200 Hz freq offset pilot_lpf(1:Npilotlpf-M) = pilot_lpf(M+1:Npilotlpf); @@ -320,20 +335,12 @@ function foff = rx_est_freq_offset(rx_fdm, pilot) % decimate to improve DFT resolution, window and DFT Mpilot = Fs/(2*200); % calc decimation rate given new sample rate is twice LPF freq - Mpilotfft = 256; s = pilot_lpf(1:Mpilot:Npilotlpf) .* hanning(Npilotlpf/Mpilot)'; S = abs(fft(s, Mpilotfft)); - %figure(3) - %plot(real(pilot_baseband)) - %plot(real(s)) - %figure(4) - %plot(abs(fft(pilot_baseband))) - %plot(S) - % peak pick and convert to Hz - [x ix] = max(S); + [imax ix] = max(S); r = 2*200/Mpilotfft; % maps FFT bin to frequency in Hz if ix > Mpilotfft/2 @@ -552,6 +559,7 @@ phase_rx = ones(Nc+1,1); % Freq offset estimator constants +global Mpilotfft = 256; global Npilotcoeff = 30; % number of pilot LPF coeffs global pilot_coeff = fir1(Npilotcoeff, 200/(Fs/2))'; % 200Hz LPF global Npilotbaseband = Npilotcoeff + 4*M; % number of pilot baseband samples reqd for pilot LPF @@ -559,10 +567,14 @@ global Npilotlpf = 4*M; % number of samples we D % Freq offset estimator states constants -global pilot_baseband; -pilot_baseband = zeros(1, Npilotbaseband); % pilot baseband samples -global pilot_lpf -pilot_lpf = zeros(1, Npilotlpf); % LPF pilot samples +global pilot_baseband1; +global pilot_baseband2; +pilot_baseband1 = zeros(1, Npilotbaseband); % pilot baseband samples +pilot_baseband2 = zeros(1, Npilotbaseband); % pilot baseband samples +global pilot_lpf1 +global pilot_lpf2 +pilot_lpf1 = zeros(1, Npilotlpf); % LPF pilot samples +pilot_lpf2 = zeros(1, Npilotlpf); % LPF pilot samples % Timing estimator states diff --git a/codec2-dev/octave/fdmdv_demod.m b/codec2-dev/octave/fdmdv_demod.m index 29cd6dc6..ed26e7a2 100644 --- a/codec2-dev/octave/fdmdv_demod.m +++ b/codec2-dev/octave/fdmdv_demod.m @@ -33,6 +33,7 @@ function fdmdv_demod(rawfilename, nbits) pilot_freq = freq(Nc+1); pilot_phase = 1; pilot_filter_mem = zeros(1, Nfilter); + prev_pilot = zeros(M,1); % BER stats @@ -40,6 +41,8 @@ function fdmdv_demod(rawfilename, nbits) total_bits = 0; bit_errors_log = []; sync_log = []; + test_frame_sync_log = []; + test_frame_sync_state = 0; rx_symbols_log = []; rx_timing_log = []; @@ -53,8 +56,10 @@ function fdmdv_demod(rawfilename, nbits) % frequency offset estimation and correction [pilot pilot_rx_bit pilot_symbol pilot_filter_mem pilot_phase] = generate_pilot_fdm(pilot_rx_bit, pilot_symbol, pilot_filter_mem, pilot_phase, pilot_freq); - foff = rx_est_freq_offset(rx_fdm, pilot); + foff = rx_est_freq_offset(rx_fdm, pilot, prev_pilot); + prev_pilot = pilot; foff_log = [ foff_log foff ]; + foff = 0; foff_rect = exp(j*2*pi*foff/Fs); for i=1:M @@ -80,16 +85,38 @@ function fdmdv_demod(rawfilename, nbits) sync_log = [sync_log sync]; % count bit errors if we find a test frame + % Allow 15 frames for filter memories to fill and time est to settle [test_frame_sync bit_errors] = put_test_bits(rx_bits); if (test_frame_sync == 1 && f > 15) total_bit_errors = total_bit_errors + bit_errors; total_bits = total_bits + Ntest_bits; bit_errors_log = [bit_errors_log bit_errors]; - else - bit_errors_log = [bit_errors_log -1]; end + % test frame sync state machine, just for more informative plots + + next_test_frame_sync_state = test_frame_sync_state; + if (test_frame_sync_state == 0) + if (test_frame_sync == 1) + next_test_frame_sync_state = 1; + test_frame_count = 0; + end + end + + if (test_frame_sync_state == 1) + % we only expect another test_frame_sync pulse every 4 symbols + test_frame_count++; + if (test_frame_count == 4) + test_frame_count = 0; + if ((test_frame_sync == 0)) + next_test_frame_sync_state = 0; + end + end + end + test_frame_sync_state = next_test_frame_sync_state; + test_frame_sync_log = [test_frame_sync_log test_frame_sync_state]; + end % --------------------------------------------------------------------- @@ -134,11 +161,16 @@ function fdmdv_demod(rawfilename, nbits) figure(4) clf; - subplot(211) + subplot(311) stem(sync_log) + axis([0 frames 0 1.5]); title('BPSK Sync') - subplot(212) + subplot(312) stem(bit_errors_log); - title('Bit Errors for test data') + title('Bit Errors for test frames') + subplot(313) + plot(test_frame_sync_log); + axis([0 frames 0 1.5]); + title('Test Frame Sync') endfunction diff --git a/codec2-dev/octave/fdmdv_ut.m b/codec2-dev/octave/fdmdv_ut.m index b9de56bd..40557397 100644 --- a/codec2-dev/octave/fdmdv_ut.m +++ b/codec2-dev/octave/fdmdv_ut.m @@ -12,7 +12,7 @@ fdmdv; % load modem code % Simulation Parameters -------------------------------------- frames = 50; -EbNo_dB = 73; +EbNo_dB = 7.3; Foff_hz = 0; modulation = 'dqpsk'; hpa_clip = 150; @@ -25,8 +25,6 @@ rx_phase_log = 0; rx_timing_log = 0; tx_pwr = 0; noise_pwr = 0; -total_bit_errors = 0; -total_bits = 0; rx_fdm_log = []; rx_baseband_log = []; rx_bits_offset = zeros(Nc*Nb*2); @@ -35,8 +33,15 @@ prev_rx_symbols = ones(Nc+1,1); foff_log = []; tx_baseband_log = []; tx_fdm_log = []; -sync_log = []; + +% BER stats + +total_bit_errors = 0; +total_bits = 0; bit_errors_log = []; +sync_log = []; +test_frame_sync_log = []; +test_frame_sync_state = 0; % pilot states, used for copy of pilot at rx @@ -45,6 +50,7 @@ pilot_symbol = sqrt(2); pilot_freq = freq(Nc+1); pilot_phase = 1; pilot_filter_mem = zeros(1, Nfilter); +prev_pilot = zeros(M,1); % fixed delay simuation @@ -148,7 +154,8 @@ for f=1:frames % frequency offset estimation and correction [pilot pilot_rx_bit pilot_symbol pilot_filter_mem pilot_phase] = generate_pilot_fdm(pilot_rx_bit, pilot_symbol, pilot_filter_mem, pilot_phase, pilot_freq); - foff = rx_est_freq_offset(rx_fdm_delay, pilot); + foff = rx_est_freq_offset(rx_fdm_delay, pilot, prev_pilot); + prev_pilot = pilot; foff_log = [ foff_log foff ]; %foff = 0; foff_rect = exp(j*2*pi*foff/Fs); @@ -188,10 +195,30 @@ for f=1:frames total_bit_errors = total_bit_errors + bit_errors; total_bits = total_bits + Ntest_bits; bit_errors_log = [bit_errors_log bit_errors]; - else - bit_errors_log = [bit_errors_log -1]; + end + + % test frame sync state machine, just for more informative plots + + next_test_frame_sync_state = test_frame_sync_state; + if (test_frame_sync_state == 0) + if (test_frame_sync == 1) + next_test_frame_sync_state = 1; + test_frame_count = 0; + end end + if (test_frame_sync_state == 1) + % we only expect another test_frame_sync pulse every 4 symbols + test_frame_count++; + if (test_frame_count == 4) + test_frame_count = 0; + if ((test_frame_sync == 0)) + next_test_frame_sync_state = 0; + end + end + end + test_frame_sync_state = next_test_frame_sync_state; + test_frame_sync_log = [test_frame_sync_log test_frame_sync_state]; end % --------------------------------------------------------------------- @@ -216,7 +243,7 @@ printf("Eb/No (meas): %2.2f (%2.2f) dB %d bits %d errors BER: (%1.4f) PAPR: % figure(1) clf; [n m] = size(rx_symbols_log); -plot(real(rx_symbols_log(1:Nc+1,10:m)),imag(rx_symbols_log(1:Nc+1,10:m)),'+') +plot(real(rx_symbols_log(1:Nc+1,15:m)),imag(rx_symbols_log(1:Nc+1,15:m)),'+') axis([-2 2 -2 2]); title('Scatter Diagram'); @@ -243,12 +270,17 @@ title('FDM Rx Spectrum'); figure(4) clf; -subplot(211) +subplot(311) stem(sync_log) +axis([0 frames 0 1.5]); title('BPSK Sync') -subplot(212) +subplot(312) stem(bit_errors_log); -title('Bit Errors for test data') +title('Bit Errors for test frames') +subplot(313) +plot(test_frame_sync_log); +axis([0 frames 0 1.5]); +title('Test Frame Sync') % TODO % + handling sample slips, extra plus/minus samples -- 2.25.1