From: drowe67 Date: Tue, 20 Oct 2015 22:20:52 +0000 (+0000) Subject: all combos working and script to run thru them, plus comments X-Git-Url: http://git.whiteaudio.com/gitweb/?a=commitdiff_plain;h=cdae0a67fabec298505fbbc328b05630e2860bb1;p=freetel-svn-tracking.git all combos working and script to run thru them, plus comments git-svn-id: https://svn.code.sf.net/p/freetel/code@2459 01035d8c-6547-0410-b346-abe4f91aad63 --- diff --git a/codec2-dev/octave/legacyfsk.m b/codec2-dev/octave/legacyfsk.m deleted file mode 100644 index 44138e7a..00000000 --- a/codec2-dev/octave/legacyfsk.m +++ /dev/null @@ -1,390 +0,0 @@ -% legacyfsk.m -% David Rowe October 2015 -% -% Attempt at FSK waveform that can pass through legacy FM radios -% but still be optimally demodulated by SDRs. It doesn't have to be -% optimally demodulated by legacy radios. Trick is getting it to -% pass through 300-3000Hz audio filters in leagcy radios. -% -% [X] code up modulator -% [X] manchester two bit symbols -% [X] plot spectrum -% [ ] demodulate -% [ ] measure BER compared to ideal coherent FSK -% [ ] measure and adjust for tone inversion through FM demod - -1; - -fm; % analog FM library - - -function states = legacyfsk_init() - Fs = states.Fs = 96000; - Rs = states.Rs = 4800; % symbol rate over channel - Ts = states.Ts = Fs/Rs; % symbol period in samples - M = states.M = 4; % mFSK, either 2 or 4 - bpsym = state.Rb = log2(M); % bits per symbol over channel - rate = states.rate = 0.5; % Manchester code rate - nbits = 100; - nbits = states.nbits = 100; % number of payload data symbols/frame - nbits2 = states.nbits2 = nbits/rate; % number of symbols/frame over channel after manchester encoding - nsym = states.nsym = nbits2/log2(M); % number of symbols per frame - nsam = states.nsam = nsym*Ts; - - printf("Rs: %d M: %d bpsym: %d nbits: %d nbits2: %d nsym: %d nsam: %d\n", Rs, M, bpsym, nbits, nbits2, nsym, nsam); - - states.fc = states.Fs/4; - if states.M == 2 - states.f(1) = states.fc - Rs/2; - states.f(2) = states.fc + Rs/2; - else - states.f(1) = states.fc - 3*Rs/2; - states.f(2) = states.fc - Rs/2; - states.f(3) = states.fc + Rs/2; - states.f(4) = states.fc + 3*Rs/2; - end -endfunction - - -% test modulator function - -function tx = legacyfsk_mod(states, tx_bits) - Fs = states.Fs; - Ts = states.Ts; - Rs = states.Rs; - f = states.f; - M = states.M; - nsym = states.nsym; - - tx = zeros(Ts*length(tx_bits)/log2(M),1); - tx_phase = 0; - - step = log2(M); - k = 1; - for i=1:step:length(tx_bits) - if M == 2 - tone = tx_bits(i) + 1; - else - tone = (tx_bits(i:i+1) * [2 1]') + 1; - end - tx_phase_vec = tx_phase + (1:Ts)*2*pi*f(tone)/Fs; - tx((k-1)*Ts+1:k*Ts) = 2.0*cos(tx_phase_vec); k++; - tx_phase = tx_phase_vec(Ts) - floor(tx_phase_vec(Ts)/(2*pi))*2*pi; - end - -endfunction - - -function run_sim - - frames = 100; - timing_offset = 0.0; - test_frame_mode = 1; - demod = 2; - - if demod == 1 - EbNodB = 8.5+5; - else - EbNodB = 8.5; - end - EbNodB = 6; - - % init fsk modem - - more off - rand('state',1); - randn('state',1); - states = legacyfsk_init(); - Fs = states.Fs; - nbits = states.nbits; - nbits2 = states.nbits2; - Ts = states.Ts; - Rs = states.Rs; - nsam = states.nsam; - rate = states.rate; - M = states.M; - - % init analog FM modem - - fm_states.Fs = Fs; - fm_max = fm_states.fm_max = 3E3; - fd = fm_states.fd = 5E3; - fm_states.fc = states.fc; - - fm_states.pre_emp = 0; - fm_states.de_emp = 1; - fm_states.Ts = 1; - fm_states.output_filter = 1; - fm_states = analog_fm_init(fm_states); - [b, a] = cheby1(4, 1, 300/Fs, 'high'); % 300Hz HPF to simulate FM radios - - % init sim states - - rx_bits_buf = zeros(1,2*nbits2); - Terrs = Tbits = 0; - state = 0; - nerr_log = []; - - % set up the channel noise. We have log(M)*rate payload bits/symbol - % we have log2(M) bits/symbol, and rate bits per payload symbol - % TODO: explain this better as Im confused! - - EbNo = 10^(EbNodB/10); - EsNo = EbNo*rate*log2(M); - variance = states.Fs/((states.Rs)*EsNo); - printf("EbNodB: %3.1f EbNo: %3.2f EsNo: %3.2f\n", EbNodB, EbNo, EsNo); - - % set up the input bits - - if test_frame_mode == 1 - % test frame of bits, which we repeat for convenience when BER testing - test_frame = round(rand(1, nbits)); - tx_bits = []; - for i=1:frames+1 - tx_bits = [tx_bits test_frame]; - end - end - if test_frame_mode == 2 - % random bits, just to make sure sync algs work on random data - tx_bits = round(rand(1, nbits*(frames+1))); - end - if test_frame_mode == 3 - % ...10101... sequence - tx_bits = zeros(1, nbits*(frames+1)); - tx_bits(1:2:length(tx_bits)) = 1; - end - - % Manchester encoding, which removes DC term in baseband signal, - % making the waveform friendly to old-school legacy FM radios with - % voiceband filtering. The "code rate" is 0.5, which means we have - % encode one input bit into 2 output bits. The 2FSK encoder takes - % one input bit, the 4FSK encoder two input bits. - - tx_bits_encoded = zeros(1,length(tx_bits)*2); - fsk2_enc = [[1 0]; [0 1]]; - % -1.5 1.5 1.5 -1.5 -0.5 0.5 0.5 -0.5 - % 0 3 3 0 1 2 2 1 - fsk4_enc = [[0 0 1 1]; [1 1 0 0]; [0 1 1 0]; [1 0 0 1]]; - k=1; - if M == 2 - for i=1:2:length(tx_bits_encoded) - input_bit = tx_bits(k); k++; - tx_bits_encoded(i:i+1) = fsk2_enc(input_bit+1,:); - end - else - for i=1:4:length(tx_bits_encoded) - input_bits = tx_bits(k:k+1) * [2 1]'; k+=2; - tx_bits_encoded(i:i+3) = fsk4_enc(input_bits+1,:); - end - end - tx_bits_encoded(1:6*2) - % use ideal FSK modulator (note: need to try using analog FM modulator) - - tx = legacyfsk_mod(states, tx_bits_encoded); - noise = sqrt(variance)*randn(length(tx),1); - rx = tx + noise; - timing_offset_samples = round(timing_offset*Ts); - rx = [zeros(timing_offset_samples,1); rx]; - - if demod == 1 - % use analog FM demodulator, aka a $40 Baofeng - - [rx_out rx_bb] = analog_fm_demod(fm_states, rx'); - rx_out_hp = filter(b,a,rx_out); - - rx_timing = filter(ones(1,Ts),1,rx_out_hp); - - % TODO: for 4FSK determine amplitude/decn boundaries, choose closest to demod each symbol - end - - if demod == 2 - - % optimal non-coherent demod at Rs - - rx_timing_sig = zeros(1,length(rx)); - for m=1:M - phi_vec = (1:length(rx))*2*pi*states.f(m)/Fs; - dc = rx' .* exp(-j*phi_vec); - rx_filt(m,:) = abs(filter(ones(1,Ts),1,dc)); - rx_timing_sig = rx_timing_sig + rx_filt(m,1:length(rx)); - end - end - - % Estimate fine timing using line at Rs/2 that Manchester encoding provides - % We need this to sync up to Manchester codewords. TODO plot signal and - % timing "line" we extract - - Np = length(rx_timing_sig); - w = 2*pi*(Rs)/Fs; - x = (rx_timing_sig .^ 2) * exp(-j*w*(0:Np-1))'; - norm_rx_timing = angle(x)/(2*pi) - 0.42; - rx_timing = round(norm_rx_timing*Ts); - rx_timing = 0; - printf("norm_rx_timing: %4.4f rx_timing: %d\n", norm_rx_timing, rx_timing); - - % sample at optimal instant - - [tmp l] = size(rx_filt); - rx_filt_dec = rx_filt(:, Ts+rx_timing:Ts:l); - - % Max likelihood decoding of Manchester encoded symbols. Search - % through all ML possibilities to extract bits. Use energy (filter - % output sq) - - [tmp l] = size(rx_filt_dec); - if M == 2 - rx_bits = zeros(1,l/2); - k = 1; - for i=1:2:l-1 - ml = [rx_filt_dec(2,i)*rx_filt_dec(1,i+1) rx_filt_dec(1,i)*rx_filt_dec(2,i+1)]; - [mx mx_ind] = max(ml); - rx_bits(k) = mx_ind-1; k++; - end - else - rx_bits = zeros(1,l); - k = 1; - fsk4_dec = [[0 0]; [0 1]; [1 0]; [1 1]]; - for i=1:2:l-1 - %ml = [rx_filt_dec(1,i)*rx_filt_dec(4,i+1) rx_filt_dec(4,i)*rx_filt_dec(1,i+1) rx_filt_dec(2,i)*rx_filt_dec(3,i+1) rx_filt_dec(3,i)*rx_filt_dec(2,i+1)]; - ml = [(rx_filt_dec(1,i)+rx_filt_dec(4,i+1)) (rx_filt_dec(4,i)+rx_filt_dec(1,i+1)) (rx_filt_dec(2,i)+rx_filt_dec(3,i+1)) (rx_filt_dec(3,i)+rx_filt_dec(2,i+1))]; - [mx mx_ind] = max(ml); - rx_bits(k:k+1) = fsk4_dec(mx_ind,:); k+=2; - end - end - - % Run frame sync and BER logic over demodulated bits - - st = 1; - for f=1:frames - - % extract nin bits - - nin = nbits; - en = st + nin - 1; - - rx_bits_buf(1:nbits) = rx_bits_buf(nbits+1:2*nbits); - rx_bits_buf(nbits+1:2*nbits) = rx_bits(st:en); - - st += nin; - - % frame sync based on min BER - - if test_frame_mode == 1 - nerrs_min = nbits; - next_state = state; - if state == 0 - for i=1:nbits - error_positions = xor(rx_bits_buf(i:nbits+i-1), test_frame); - nerrs = sum(error_positions); - %printf("i: %d nerrs: %d nerrs_min: %d \n", i, nerrs, nerrs_min); - if nerrs < nerrs_min - nerrs_min = nerrs; - coarse_offset = i; - end - end - if nerrs_min < 3 - next_state = 1; - %printf("%d %d\n", coarse_offset, nerrs_min); - end - end - - if state == 1 - error_positions = xor(rx_bits_buf(coarse_offset:coarse_offset+nbits-1), test_frame); - nerrs = sum(error_positions); - Terrs += nerrs; - Tbits += nbits; - nerr_log = [nerr_log nerrs]; - end - - state = next_state; - - end - end - - if test_frame_mode == 1 - printf("frames: %d Tbits: %d Terrs: %d BER %4.3f\n", frames, Tbits, Terrs, Terrs/Tbits); - end - - % Bunch O'plots -------------------------------------------------------------- - - st = 1; en=10; - - Tx=fft(tx, Fs); - TxdB = 20*log10(abs(Tx(1:Fs/2))); - figure(1) - clf; - plot(TxdB) - axis([1 Fs/2 (max(TxdB)-100) max(TxdB)]) - title('Tx Spectrum'); - - figure(2) - clf - if demod == 1 - subplot(211) - plot(rx_out(st:en*states.Ts*2)); - title('After Analog FM demod'); - subplot(212) - plot(rx_out_hp(st:en*states.Ts*2)); - title('After 300Hz HPF'); - end - if demod == 2 - subplot(211); - plot(rx_filt(1,st*Ts:en*Ts)); - hold on; - plot(rx_filt(2,st*Ts:en*Ts),'g'); - if M == 4 - plot(rx_filt(3,st*Ts:en*Ts),'c'); - plot(rx_filt(4,st*Ts:en*Ts),'r'); - end - hold off; - title('Output of each filter'); - subplot(212); - plot(rx_filt_dec(1,st:en),'+'); - hold on; - plot(rx_filt_dec(2,st:en),'g+'); - if M == 4 - plot(rx_filt_dec(3,st:en),'c+'); - plot(rx_filt_dec(4,st:en),'r+'); - end - hold off; - title('Decimated output of each filter'); - end - - figure(3); - clf; - if demod == 1 - subplot(211) - h = freqz(b,a,Fs); - plot(20*log10(abs(h(1:4000)))) - title('300Hz HPF Response'); - subplot(212) - h = fft(rx_out, Fs); - plot(20*log10(abs(h(1:4000)))) - title('Spectrum of baseband modem signal after analog FM demod'); - end - -if 0 - figure(4); - clf; - subplot(211) - plot(rx_sd(st*Ts:en*Ts),'+') - title('Difference after filtering'); - subplot(212) - plot(rx_sd_sampled(st:en),'+'); - hold on; - plot(rx_sd_sampled(st:en),'g+') - hold off; - title('Both channels after sampling at ideal timing instant') -end - - figure(5) - clf; - subplot(211) - plot(rx_timing_sig(st*Ts:en*Ts).^2) - subplot(212) - F = abs(fft(rx_timing_sig(1:Fs))); - plot(F(100:8000)) - title('FFT of rx_timing_sig') -end - -run_sim; diff --git a/codec2-dev/octave/mancyfsk.m b/codec2-dev/octave/mancyfsk.m new file mode 100644 index 00000000..a5ac8ccb --- /dev/null +++ b/codec2-dev/octave/mancyfsk.m @@ -0,0 +1,500 @@ +% mancyfsk.m +% David Rowe October 2015 +% +% Manchester encoded 2FSK & 4FSK simulation. +% +% Attempt to design a FSK waveform that can pass through legacy FM +% radios but still be optimally demodulated by SDRs. It doesn't have +% to be optimally demodulated by legacy radios. Trick is getting it +% to pass through 300-3000Hz audio filters in legacy radios. +% +% [X] code up modulator +% [X] manchester two bit symbols +% [X] plot spectrum +% [X] demodulate using analog FM and ideal demods +% [X] measure BER compared to ideal coherent FSK + +1; + +fm; % analog FM library + + +function states = legacyfsk_init(M,Rs) + Fs = states.Fs = 96000; + states.Rs = Rs; % symbol rate over channel + Ts = states.Ts = Fs/Rs; % symbol period in samples + states.M = M; % mFSK, either 2 or 4 + bpsym = state.Rb = log2(M); % bits per symbol over channel + rate = states.rate = 0.5; % Manchester code rate + nbits = 100; + nbits = states.nbits = 100; % number of payload data symbols/frame + nbits2 = states.nbits2 = nbits/rate; % number of symbols/frame over channel after manchester encoding + nsym = states.nsym = nbits2/log2(M); % number of symbols per frame + nsam = states.nsam = nsym*Ts; + + %printf(" Rs: %d M: %d bpsym: %d nbits: %d nbits2: %d nsym: %d nsam: %d\n", Rs, M, bpsym, nbits, nbits2, nsym, nsam); + + states.fc = states.Fs/4; + if states.M == 2 + states.f(1) = states.fc - Rs/2; + states.f(2) = states.fc + Rs/2; + else + states.f(1) = states.fc - 3*Rs/2; + states.f(2) = states.fc - Rs/2; + states.f(3) = states.fc + Rs/2; + states.f(4) = states.fc + 3*Rs/2; + end +endfunction + + +% test modulator function + +function tx = legacyfsk_mod(states, tx_bits) + Fs = states.Fs; + Ts = states.Ts; + Rs = states.Rs; + f = states.f; + M = states.M; + nsym = states.nsym; + + tx = zeros(Ts*length(tx_bits)/log2(M),1); + tx_phase = 0; + + step = log2(M); + k = 1; + for i=1:step:length(tx_bits) + if M == 2 + tone = tx_bits(i) + 1; + else + tone = (tx_bits(i:i+1) * [2 1]') + 1; + end + tx_phase_vec = tx_phase + (1:Ts)*2*pi*f(tone)/Fs; + tx((k-1)*Ts+1:k*Ts) = 2.0*cos(tx_phase_vec); k++; + tx_phase = tx_phase_vec(Ts) - floor(tx_phase_vec(Ts)/(2*pi))*2*pi; + end + +endfunction + + +function run_sim(sim_in) + + frames = sim_in.frames; + test_frame_mode = sim_in.test_frame_mode; + M = sim_in.M; + Rs = sim_in.Rs; + demod = sim_in.demod; + EbNodB = sim_in.EbNodB; + timing_offset = sim_in.timing_offset; + + % rx timing has been adjusted experimentally + + if Rs == 4800 + if demod == 1 + rx_timing = 4; + else + rx_timing = 0; + end + end + if Rs == 2400 + if demod == 1 + rx_timing = 40; + else + rx_timing = 0; + end + end + + % init fsk modem + + more off + rand('state',1); + randn('state',1); + states = legacyfsk_init(M,Rs); + Fs = states.Fs; + nbits = states.nbits; + nbits2 = states.nbits2; + Ts = states.Ts; + nsam = states.nsam; + rate = states.rate; + + % init analog FM modem + + fm_states.Fs = Fs; + fm_max = fm_states.fm_max = 3E3; + fd = fm_states.fd = 5E3; + fm_states.fc = states.fc; + + fm_states.pre_emp = 0; + fm_states.de_emp = 1; + fm_states.Ts = 1; + fm_states.output_filter = 1; + fm_states = analog_fm_init(fm_states); + [b, a] = cheby1(4, 1, 300/Fs, 'high'); % 300Hz HPF to simulate FM radios + + % init sim states + + rx_bits_buf = zeros(1,2*nbits2); + Terrs = Tbits = 0; + state = 0; + nerr_log = []; + + % set up the channel noise. We have log(M)*rate payload bits/symbol + % we have log2(M) bits/symbol, and rate bits per payload symbol + % TODO: explain this better as Im confused! + + EbNo = 10^(EbNodB/10); + EsNo = EbNo*rate*log2(M); + variance = states.Fs/((states.Rs)*EsNo); + %printf("EbNodB: %3.1f EbNo: %3.2f EsNo: %3.2f\n", EbNodB, EbNo, EsNo); + + % set up the input bits + + if test_frame_mode == 1 + % test frame of bits, which we repeat for convenience when BER testing + test_frame = round(rand(1, nbits)); + tx_bits = []; + for i=1:frames+1 + tx_bits = [tx_bits test_frame]; + end + end + if test_frame_mode == 2 + % random bits, just to make sure sync algs work on random data + tx_bits = round(rand(1, nbits*(frames+1))); + end + if test_frame_mode == 3 + % ...10101... sequence + tx_bits = zeros(1, nbits*(frames+1)); + tx_bits(1:2:length(tx_bits)) = 1; + end + + % Manchester Encoding ----------------------------------------------------------- + + % Manchester encoding, which removes DC term in baseband signal, + % making the waveform friendly to old-school legacy FM radios with + % voiceband filtering. The "code rate" is 0.5, which means we have + % encode one input bit into 2 output bits. The 2FSK encoder takes + % one input bit, the 4FSK encoder two input bits. + + tx_bits_encoded = zeros(1,length(tx_bits)*2); + fsk2_enc = [[1 0]; [0 1]]; + % -1.5 1.5 1.5 -1.5 -0.5 0.5 0.5 -0.5 + % 0 3 3 0 1 2 2 1 + fsk4_enc = [[0 0 1 1]; [1 1 0 0]; [0 1 1 0]; [1 0 0 1]]; + k=1; + if M == 2 + for i=1:2:length(tx_bits_encoded) + input_bit = tx_bits(k); k++; + tx_bits_encoded(i:i+1) = fsk2_enc(input_bit+1,:); + end + else + for i=1:4:length(tx_bits_encoded) + input_bits = tx_bits(k:k+1) * [2 1]'; k+=2; + tx_bits_encoded(i:i+3) = fsk4_enc(input_bits+1,:); + end + end + + % FSK Modulator -------------------------------------------------------------- + + % use ideal FSK modulator (note: need to try using analog FM modulator) + + tx = legacyfsk_mod(states, tx_bits_encoded); + noise = sqrt(variance)*randn(length(tx),1); + rx = tx + noise; + timing_offset_samples = round(timing_offset*Ts); + rx = [zeros(timing_offset_samples,1); rx]; + + % Demodulator ---------------------------------------------------------------------------- + + if demod == 1 + % use analog FM demodulator, aka a $40 Baofeng + + [rx_out rx_bb] = analog_fm_demod(fm_states, rx'); + if sim_in.hpf + rx_out_hp = filter(b,a,rx_out); + else + rx_out_hp = rx_out; + end + rx_filt = filter(ones(1,Ts),1,rx_out_hp); + rx_timing_sig = rx_filt; + + % TODO: for 4FSK determine amplitude/decn boundaries, choose closest to demod each symbol + + end + + if demod == 2 + + % optimal non-coherent demod at Rs + + rx_timing_sig = zeros(1,length(rx)); + for m=1:M + phi_vec = (1:length(rx))*2*pi*states.f(m)/Fs; + dc = rx' .* exp(-j*phi_vec); + rx_filt(m,:) = abs(filter(ones(1,Ts),1,dc)); + rx_timing_sig = rx_timing_sig + rx_filt(m,1:length(rx)); + end + end + + % Fine timing estimation ------------------------------------------------------ + + % Estimate fine timing using line at Rs/2 that Manchester encoding provides + % We need this to sync up to Manchester codewords. TODO plot signal and + % timing "line" we extract + + Np = length(rx_timing_sig); + w = 2*pi*(Rs)/Fs; + x = (rx_timing_sig .^ 2) * exp(-j*w*(0:Np-1))'; + norm_rx_timing = angle(x)/(2*pi) - 0.42; + %rx_timing = round(norm_rx_timing*Ts); + %printf("norm_rx_timing: %4.4f rx_timing: %d\n", norm_rx_timing, rx_timing); + + % Max likelihood decoding of Manchester encoded symbols. Search + % through all ML possibilities to extract bits. Use energy (filter + % output sq) + + % Manchester Decoding -------------------------------------------------------- + + if M == 2 + if demod == 1 + + % sample at optimum instant + + [tmp l] = size(rx_filt); + rx_filt_dec = rx_filt(:, Ts+rx_timing:Ts:l); + + [tmp l] = size(rx_filt_dec); + rx_bits = zeros(1,l); + k = 1; + for i=1:2:l-1 + ml = [rx_filt_dec(i)-rx_filt_dec(i+1) -rx_filt_dec(i)+rx_filt_dec(i+1)]; + [mx mx_ind] = max(ml); + rx_bits(k) = mx_ind-1; k++; + end + end + + if demod == 2 + + % sample at optimum instant + + [tmp l] = size(rx_filt); + rx_filt_dec = rx_filt(:, Ts+rx_timing:Ts:l); + + [tmp l] = size(rx_filt_dec); + rx_bits = zeros(1,l); + k = 1; + for i=1:2:l-1 + %ml = [rx_filt_dec(2,i)*rx_filt_dec(1,i+1) rx_filt_dec(1,i)*rx_filt_dec(2,i+1)]; + ml = [rx_filt_dec(2,i)+rx_filt_dec(1,i+1) rx_filt_dec(1,i)+rx_filt_dec(2,i+1)]; + [mx mx_ind] = max(ml); + rx_bits(k) = mx_ind-1; k++; + end + end + else % M == 4 + if demod == 1 + % TODO: 4FSK version of demod + rx_bits=tx_bits; + end + if demod == 2 + % sample at optimal instant + + [tmp l] = size(rx_filt); + rx_filt_dec = rx_filt(:, Ts+rx_timing:Ts:l); + [tmp l] = size(rx_filt_dec); + rx_bits = zeros(1,l); + + k = 1; + fsk4_dec = [[0 0]; [0 1]; [1 0]; [1 1]]; + for i=1:2:l-1 + %ml = [rx_filt_dec(1,i)*rx_filt_dec(4,i+1) rx_filt_dec(4,i)*rx_filt_dec(1,i+1) rx_filt_dec(2,i)*rx_filt_dec(3,i+1) rx_filt_dec(3,i)*rx_filt_dec(2,i+1)]; + ml = [(rx_filt_dec(1,i)+rx_filt_dec(4,i+1)) (rx_filt_dec(4,i)+rx_filt_dec(1,i+1)) (rx_filt_dec(2,i)+rx_filt_dec(3,i+1)) (rx_filt_dec(3,i)+rx_filt_dec(2,i+1))]; + [mx mx_ind] = max(ml); + rx_bits(k:k+1) = fsk4_dec(mx_ind,:); k+=2; + end + end + end + + % useful for getting decoding right + %tx_bits(1:20) + %rx_bits(1:20) + + % Frame sync and BER logic ------------------------------------------------------------- + + st = 1; + for f=1:frames + + % extract nin bits + + nin = nbits; + en = st + nin - 1; + + rx_bits_buf(1:nbits) = rx_bits_buf(nbits+1:2*nbits); + rx_bits_buf(nbits+1:2*nbits) = rx_bits(st:en); + + st += nin; + + % frame sync based on min BER + + if test_frame_mode == 1 + nerrs_min = nbits; + next_state = state; + if state == 0 + for i=1:nbits + error_positions = xor(rx_bits_buf(i:nbits+i-1), test_frame); + nerrs = sum(error_positions); + %printf("i: %d nerrs: %d nerrs_min: %d \n", i, nerrs, nerrs_min); + if nerrs < nerrs_min + nerrs_min = nerrs; + coarse_offset = i; + end + end + if nerrs_min < 3 + next_state = 1; + %printf("%d %d\n", coarse_offset, nerrs_min); + end + end + + if state == 1 + error_positions = xor(rx_bits_buf(coarse_offset:coarse_offset+nbits-1), test_frame); + nerrs = sum(error_positions); + Terrs += nerrs; + Tbits += nbits; + nerr_log = [nerr_log nerrs]; + end + + state = next_state; + + end + end + + if test_frame_mode == 1 + if sim_in.verbose + printf(" demod: %d frames: %d EbNodB: %3.1f Tbits: %d Terrs: %d BER %4.3f\n", demod, frames, EbNodB, Tbits, Terrs, Terrs/Tbits); + else + printf(" EbNodB: %3.1f BER %4.3f\n", EbNodB, Terrs/Tbits); + end + end + + % Bunch O'plots -------------------------------------------------------------- + + close all; + + st = 1; en=20; + + Tx=fft(tx, Fs); + TxdB = 20*log10(abs(Tx(1:Fs/2))); + figure(1) + clf; + plot(TxdB) + axis([1 Fs/2 (max(TxdB)-100) max(TxdB)]) + title('Tx Spectrum'); + + figure(2) + clf + if demod == 1 + subplot(211) + plot(rx_filt(st*Ts:en*Ts)); + title('After integrator'); + subplot(212) + plot(rx_filt_dec(st:en),'+'); + title('Decimated output'); + end + if demod == 2 + subplot(211); + plot(rx_filt(1,st*Ts:en*Ts)); + hold on; + plot(rx_filt(2,st*Ts:en*Ts),'g'); + if M == 4 + plot(rx_filt(3,st*Ts:en*Ts),'c'); + plot(rx_filt(4,st*Ts:en*Ts),'r'); + end + hold off; + title('Output of each filter'); + subplot(212); + plot(rx_filt_dec(1,st:en),'+'); + hold on; + plot(rx_filt_dec(2,st:en),'g+'); + if M == 4 + plot(rx_filt_dec(3,st:en),'c+'); + plot(rx_filt_dec(4,st:en),'r+'); + end + hold off; + title('Decimated output of each filter'); + end + + figure(3) + clf; + subplot(211) + plot(rx_timing_sig(st*Ts:en*Ts).^2) + title('rx-timing-sig') + subplot(212) + F = abs(fft(rx_timing_sig(1:Fs))); + plot(F(100:8000)) + title('FFT of rx-timing-sig') + + if demod == 1 + figure(4); + clf; + h = fft(rx_out, Fs); + hdB = 20*log10(abs(h)); + plot(hdB(1:4000)) + title('Spectrum of baseband modem signal after analog FM demod'); + axis([1 4000 (max(hdB)-40) max(hdB)]) + end + + if demod == 1 + figure(5) + clf; + subplot(211) + plot(rx_out(st*Ts:en*Ts)); + title('baseband modem signal after analog FM demod'); + subplot(212) + plot(rx_out_hp(st*Ts:en*Ts)); + title('baseband modem signal after 300Hz filter'); + end +end + + +% Run various permutations of simulation here --------------------------------------- + +function run_single + + sim_in.frames = 100; + sim_in.test_frame_mode = 1; + sim_in.M = 2; + sim_in.Rs = 2400; + sim_in.demod = 1; + sim_in.EbNodB = 15; + sim_in.timing_offset = 0.0; + sim_in.hpf = 1; + sim_in.verbose = 1; + + run_sim(sim_in); +endfunction + + +function run_lots + + % adjusted a few scenarios for about 2% BER so we can compare + + sim_in.frames = 100; + sim_in.test_frame_mode = 1; + sim_in.M = 2; + sim_in.Rs = 4800; + sim_in.demod = 1; + sim_in.EbNodB = 12; + sim_in.timing_offset = 0.0; + sim_in.hpf = 1; + sim_in.verbose = 0; + + printf("Rs=4800 2FSK ideal demod\n"); + sim_in.EbNodB = 8.5; sim_in.demod = 2; run_sim(sim_in); + printf("Rs=4800 2FSK analog FM demod, not too shabby and pushes 2400bit/s thru a $40 HT!\n"); + sim_in.EbNodB = 12; sim_in.demod = 1; run_sim(sim_in); + printf("Rs=2400 2FSK analog FM demod, needs more power for same BER! Che?\n"); + sim_in.Rs = 2400; sim_in.EbNodB = 15; run_sim(sim_in); + printf("Hmm, doesnt improve with no 300Hz HPF, maybe due to less deviation?\n"); + sim_in.hpf = 0; run_sim(sim_in); + printf("Rs=2400 4FSK ideal demod, nice low Eb/No!\n"); + sim_in.demod = 2; sim_in.M = 4; sim_in.Rs = 2400; sim_in.EbNodB = 6; run_sim(sim_in); +endfunction + +%run_single; +run_lots;