From 079e8bf269a124b4b8ffc90ccc02de0381d04eb3 Mon Sep 17 00:00:00 2001 From: drowe67 Date: Tue, 16 May 2017 03:03:30 +0000 Subject: [PATCH] added analog signal simulation and interleaving of (224,112) code, and a basic sync system for rx git-svn-id: https://svn.code.sf.net/p/freetel/code@3137 01035d8c-6547-0410-b346-abe4f91aad63 --- codec2-dev/octave/ofdm_rx.m | 170 ++++++++++++++++++++++++++++-------- codec2-dev/octave/ofdm_tx.m | 57 +++++++++--- 2 files changed, 182 insertions(+), 45 deletions(-) diff --git a/codec2-dev/octave/ofdm_rx.m b/codec2-dev/octave/ofdm_rx.m index ba3b1425..c99b8b0c 100644 --- a/codec2-dev/octave/ofdm_rx.m +++ b/codec2-dev/octave/ofdm_rx.m @@ -11,7 +11,7 @@ [ ] way to fall out of sync #} -function ofdm_rx(filename) +function ofdm_rx(filename, interleave_frames = 1) ofdm_lib; ldpc; more off; @@ -35,18 +35,6 @@ function ofdm_rx(filename) [code_param framesize rate] = ldpc_init_user(HRA_112_112, modulation, mod_order, mapping); assert(Nbitsperframe == code_param.code_bits_per_frame); - % fixed test frame of tx bits - - rand('seed', 100); - tx_bits = round(rand(1,code_param.data_bits_per_frame)); - tx_codeword = LdpcEncode(tx_bits, code_param.H_rows, code_param.P_matrix); - - % init logs and BER stats - - rx_bits = []; rx_np = []; timing_est_log = []; delta_t_log = []; foff_est_hz_log = []; - phase_est_pilot_log = []; - Terrs = Tbits = 0; - % load real samples from file Ascale= 2E5; @@ -54,6 +42,46 @@ function ofdm_rx(filename) Nsam = length(rx); Nframes = floor(Nsam/Nsamperframe); prx = 1; + % buffers for interleaved frames + + Nsymbolsperframe = code_param.code_bits_per_frame/bps + Nsymbolsperinterleavedframe = interleave_frames*Nsymbolsperframe + rx_np = zeros(1, Nsymbolsperinterleavedframe); + rx_amp = zeros(1, Nsymbolsperinterleavedframe); + + % OK generate tx frame for BER calcs + + rand('seed', 100); + tx_bits = []; tx_codewords = []; + for f=1:interleave_frames + atx_bits = round(rand(1,code_param.data_bits_per_frame)); + tx_bits = [tx_bits atx_bits]; + acodeword = LdpcEncode(atx_bits, code_param.H_rows, code_param.P_matrix); + tx_codewords = [tx_codewords acodeword]; + end + + % used for rx frame sync on interleaved symbols - we demod the + % entire interleaved frame to raw bits + + tx_symbols = []; + for s=1:Nsymbolsperinterleavedframe + tx_symbols = [tx_symbols qpsk_mod( tx_codewords(2*(s-1)+1:2*s) )]; + end + + tx_symbols = gp_interleave(tx_symbols); + + tx_bits_raw = []; + for s=1:Nsymbolsperinterleavedframe + tx_bits_raw = [tx_bits_raw qpsk_demod(tx_symbols(s))]; + end + + % init logs and BER stats + + rx_bits = []; rx_np_log = []; timing_est_log = []; delta_t_log = []; foff_est_hz_log = []; + phase_est_pilot_log = []; + Terrs = Tbits = Terrs_coded = Tbits_coded = Tpackets = Tpacketerrs = 0; + Nbitspervocframe = 28; + % 'prime' rx buf to get correct coarse timing (for now) prx = 1; @@ -61,7 +89,7 @@ function ofdm_rx(filename) states.rxbuf(Nrxbuf-nin+1:Nrxbuf) = rx(prx:nin); prx += nin; - state = 'searching'; + state = 'searching'; frame_count = 0; % main loop ---------------------------------------------------------------- @@ -80,30 +108,53 @@ function ofdm_rx(filename) prx += states.nin; [rx_bits_raw states aphase_est_pilot_log arx_np arx_amp] = ofdm_demod(states, rxbuf_in); - rx_codeword = ldpc_dec(code_param, max_iterations, demod_type, decoder_type, arx_np, min(EsNo,30), arx_amp); - rx_bits = rx_codeword(1:code_param.data_bits_per_frame); + frame_count++; - % measure errors and iterate start machine + % update sliding windows of rx-ed symbols and symbol amplitudes - errors = xor(tx_bits, rx_bits); - Nerrs = sum(errors); - next_state = state; + rx_np(1:Nsymbolsperinterleavedframe-Nsymbolsperframe) = rx_np(Nsymbolsperframe+1:Nsymbolsperinterleavedframe); + rx_np(Nsymbolsperinterleavedframe-Nsymbolsperframe+1:Nsymbolsperinterleavedframe) = arx_np; + rx_amp(1:Nsymbolsperinterleavedframe-Nsymbolsperframe) = rx_amp(Nsymbolsperframe+1:Nsymbolsperinterleavedframe); + rx_amp(Nsymbolsperinterleavedframe-Nsymbolsperframe+1:Nsymbolsperinterleavedframe) = arx_amp; + + printf("f: %d state: %s frame_count: %d\n", f, state, frame_count); + + % If looking for sync: check raw BER on frame just received + % against all possible positions in the interleaver frame. + % iterate state machine ------------------------------------ + + next_state = state; if strcmp(state,'searching') - if Nerrs == 0 - next_state = 'synced'; + + % If looking for sync: check raw BER on frame just received + % against all possible positions in the interleaver frame. + + for ff=1:interleave_frames + st = (ff-1)*Nbitsperframe+1; en = st + Nbitsperframe - 1; + errors = xor(tx_bits_raw(st:en), rx_bits_raw); + Nerrs = sum(errors); + printf(" ff: %d Nerrs: %d\n", ff, Nerrs); + if Nerrs/Nbitsperframe < 0.1 + next_state = 'synced'; + frame_count = ff; + if f < interleave_frames + % point trying a LDPC decode if we don't have a full frame! + frame_count -= interleave_frames; + end + end end end state = next_state; if strcmp(state,'searching') - % attempt coarse timing estimate (i.e. detect start of frame) + % still searching? Attempt coarse timing estimate (i.e. detect start of frame) st = M+Ncp + Nsamperframe + 1; en = st + 2*Nsamperframe; [ct_est foff_est] = coarse_sync(states, states.rxbuf(st:en), states.rate_fs_pilot_samples); if states.verbose - printf("Nerrs: %d ct_est: %4d foff_est: %3.1f\n", Nerrs, ct_est, foff_est); + printf(" Nerrs: %d ct_est: %4d foff_est: %3.1f\n", Nerrs, ct_est, foff_est); end % calculate number of samples we need on next buffer to get into sync @@ -114,29 +165,75 @@ function ofdm_rx(filename) states.sample_point = states.timing_est = 1; states.foff_est_hz = foff_est; + end + + if strcmp(state,'synced') && (frame_count == interleave_frames) + + % de-interleave QPSK symbols - else + arx_np = gp_deinterleave(rx_np); + arx_amp = gp_deinterleave(rx_amp); - % we are in sync so log states and bit errors + % LDPC decode - rx_np = [rx_np arx_np]; + rx_bits = []; + for ff=1:interleave_frames + st = (ff-1)*Nbitsperframe/bps+1; en = st + Nbitsperframe/bps - 1; + rx_codeword = ldpc_dec(code_param, max_iterations, demod_type, decoder_type, arx_np(st:en), min(EsNo,30), arx_amp(st:en)); + rx_bits = [rx_bits rx_codeword(1:code_param.data_bits_per_frame)]; + end + + errors_coded = xor(tx_bits, rx_bits); + Nerrs_coded = sum(errors_coded); + Terrs_coded += Nerrs_coded; + Tbits_coded += code_param.data_bits_per_frame*interleave_frames; + Nerrs_coded_log(f) = Nerrs_coded; + + printf(" Nerrs_coded: %d\n", Nerrs_coded); + + % we are in sync so log states and bit/packet error stats + + rx_np_log = [rx_np_log arx_np]; timing_est_log = [timing_est_log states.timing_est]; delta_t_log = [delta_t_log states.delta_t]; foff_est_hz_log = [foff_est_hz_log states.foff_est_hz]; phase_est_pilot_log = [phase_est_pilot_log; aphase_est_pilot_log]; - % measure bit errors + % measure uncoded bit errors + rx_bits_raw = []; + for s=1:Nsymbolsperinterleavedframe + rx_bits_raw = [rx_bits_raw qpsk_demod(rx_np(s))]; + end + errors = xor(tx_bits_raw, rx_bits_raw); + Nerrs = sum(errors); Terrs += Nerrs; Nerrs_log(f) = Nerrs; - Tbits += code_param.data_bits_per_frame; + Tbits += code_param.code_bits_per_frame*interleave_frames; + + % measure packet errors based on Codec 2 packet size + + Nvocframes = floor(code_param.data_bits_per_frame*interleave_frames/Nbitspervocframe); + for fv=1:Nvocframes + st = (fv-1)*Nbitspervocframe + 1; + en = st + Nbitspervocframe - 1; + Nvocpacketerrs = sum(xor(tx_bits(st:en), rx_bits(st:en))); + if Nvocpacketerrs + Tpacketerrs++; + end + Tpackets++; + end + + frame_count = 0; end end - printf("Coded BER: %5.4f Tbits: %d Terrs: %d\n", Terrs/Tbits, Tbits, Terrs); + printf("Coded BER: %5.4f Tbits: %5d Terrs: %5d PER: %5.4f Tpacketerrs: %5d Tpackets: %5d\n", + Terrs_coded/Tbits_coded, Tbits_coded, Terrs_coded, Tpacketerrs/Tpackets, Tpacketerrs, Tpackets); + printf("Raw BER..: %5.4f Tbits: %5d Terrs: %5d\n", Terrs/Tbits, Tbits, Terrs); figure(1); clf; - plot(rx_np,'+'); + plot(rx_np_log,'+'); axis([-2 2 -2 2]); title('Scatter'); @@ -155,11 +252,14 @@ function ofdm_rx(filename) figure(4); clf; plot(foff_est_hz_log) - mx = max(abs(foff_est_hz_log)); - axis([1 max(Nframes,2) -mx mx]); + axis([1 max(Nframes,2) -2 2]); title('Fine Freq'); + ylabel('Hz') figure(5); clf; - %plot(Nerrs_log); - + subplot(211) + title('Nerrs Log') + stem(Nerrs_log); + subplot(212) + stem(Nerrs_coded_log); endfunction diff --git a/codec2-dev/octave/ofdm_tx.m b/codec2-dev/octave/ofdm_tx.m index a414fc4a..d706a491 100644 --- a/codec2-dev/octave/ofdm_tx.m +++ b/codec2-dev/octave/ofdm_tx.m @@ -19,9 +19,10 @@ #} -function ofdm_tx(filename, Nsec, EbNodB=100, channel='awgn', freq_offset_Hz=0) +function ofdm_tx(filename, Nsec, interleave_frames = 1, EbNodB=100, channel='awgn', freq_offset_Hz=0) ofdm_lib; ldpc; + gp_interleaver; % init modem @@ -38,20 +39,56 @@ function ofdm_tx(filename, Nsec, EbNodB=100, channel='awgn', freq_offset_Hz=0) [code_param framesize rate] = ldpc_init_user(HRA_112_112, modulation, mod_order, mapping); assert(Nbitsperframe == code_param.code_bits_per_frame); - % generate fixed test frame of tx bits and run OFDM modulator - % todo: maybe extend this to 4 or 8 frames, one is a bit short + % Generate fixed test frame of tx bits and run OFDM modulator Nrows = Nsec*Rs; Nframes = floor((Nrows-1)/Ns); - rand('seed', 100); - tx_bits = round(rand(1,code_param.data_bits_per_frame)); - tx_codeword = LdpcEncode(tx_bits, code_param.H_rows, code_param.P_matrix); + % Adjust Nframes so we have an integer number of interleaver frames + % in simulation + + Nframes = interleave_frames*round(Nframes/interleave_frames); + % OK generate an interleaver frame of codewords from random data bits + + rand('seed', 100); + tx_bits = tx_symbols = []; + for f=1:interleave_frames + atx_bits = round(rand(1,code_param.data_bits_per_frame)); + tx_bits = [tx_bits atx_bits]; + codeword = LdpcEncode(atx_bits, code_param.H_rows, code_param.P_matrix); + for b=1:2:Nbitsperframe + tx_symbols = [tx_symbols qpsk_mod(codeword(b:b+1))]; + end + end + + tx_symbols = gp_interleave(tx_symbols); + + atx = []; + for f=1:interleave_frames + st = (f-1)*Nbitsperframe/bps+1; en = st + Nbitsperframe/bps-1; + atx = [atx ofdm_txframe(states, tx_symbols(st:en))]; + end + tx = []; - for f=1:Nframes - tx = [tx ofdm_mod(states, tx_codeword)]; + for f=1:Nframes/interleave_frames + tx = [tx atx]; end + + % not very pretty way to process analog signals with exactly the same channel + + analog_hack = 0; + if analog_hack + s = load_raw('../raw/ve9qrp_10s.raw')'; + tx = hilbert(s); + + % normalise power to same as ofdm tx + + nom_tx_pwr = 2/(Ns*(M*M)) + Nc/(M*M); + tx_pwr = var(tx); + tx *= sqrt(nom_tx_pwr/tx_pwr); + end + Nsam = length(tx); % channel simulation @@ -61,7 +98,7 @@ function ofdm_tx(filename, Nsec, EbNodB=100, channel='awgn', freq_offset_Hz=0) woffset = 2*pi*freq_offset_Hz/Fs; SNRdB = EbNodB + 10*log10(700/3000); - printf("EbNo: %3.1f dB SNR(3k) est: %3.1f dB foff: %3.1fHz\n", EbNodB, SNRdB, freq_offset_Hz); + printf("EbNo: %3.1f dB SNR(3k) est: %3.1f dB foff: %3.1fHz ", EbNodB, SNRdB, freq_offset_Hz); % set up HF model --------------------------------------------------------------- @@ -108,7 +145,7 @@ function ofdm_tx(filename, Nsec, EbNodB=100, channel='awgn', freq_offset_Hz=0) noise = sqrt(variance/2)*0.5*randn(1,Nsam); rx = real(rx) + noise; - printf("measured SNR: %3.2f dB\n", 10*log10(var(real(tx))/var(noise))); + printf("measured SNR: %3.2f dB\n", 10*log10(var(real(tx))/var(noise))+10*log10(4000) - 10*log10(3000)); Ascale = 2E5; frx=fopen(filename,"wb"); fwrite(frx, Ascale*rx, "short"); fclose(frx); -- 2.25.1