two modems running, 4.7dB between them

git-svn-id: https://svn.code.sf.net/p/freetel/code@2452 01035d8c-6547-0410-b346-abe4f91aad63

diff --git a/codec2-dev/octave/legacyfsk.m b/codec2-dev/octave/legacyfsk.m
index a74caedb837f89bdcc282e5b2fb7736816e9d108..150e2a12380bfb59e93dabfa1c3243e7c964b779 100644 (file)
--- a/codec2-dev/octave/legacyfsk.m
+++ b/codec2-dev/octave/legacyfsk.m
@@ -11,6 +11,7 @@
 %     [X] plot spectrum
 % [ ] demodulate
 % [ ] measure BER compared to ideal coherent FSK
+% [ ] measure and adjust for tone inversion through FM demod
 
 1;
 
@@ -23,6 +24,10 @@ function states = legacyfsk_init()
   Ts = states.Ts = Fs/Rs;
   nbits = states.nbits = 100;                  % number of payload data symbols/frame
   nbits2 = states.nbits2 = states.nbits*2;     % number of bits/frame over channel after manchester encoding
+
+  states.fc = states.Fs/4;
+  states.f1 = states.fc - Rs/2;
+  states.f2 = states.fc + Rs/2;
 endfunction
 
 
@@ -34,7 +39,7 @@ function tx = legacyfsk_mod(states, tx_bits)
     Ts = states.Ts;
     Fs = states.Fs;
     Rs = states.Rs;
-    f1 = 24E3-Rs/2; f2 = 24E3+Rs/2;
+    f1 = states.f1; f2 = states.f2;
 
     for i=1:length(tx_bits)
       if tx_bits(i) == 0
@@ -50,9 +55,16 @@ endfunction
 
 function run_sim
 
-  frames = 10;
-  EbNodB = 12.5;
+  frames = 100;
+  timing_offset = 0.0;
   test_frame_mode = 1;
+  demod = 1;
+
+  if demod == 1
+    EbNodB = 8.5+4.7;
+  else
+    EbNodB = 8.5;
+  end
 
   % init fsk modem
 
@@ -64,13 +76,14 @@ function run_sim
   nbits = states.nbits;
   nbits2 = states.nbits2;
   Ts = states.Ts;
+  Rs = states.Rs;
 
   % init analog FM modem
 
   fm_states.Fs = Fs;  
   fm_max = fm_states.fm_max = 3E3;
   fd = fm_states.fd = 5E3;
-  fm_states.fc = 24E3;
+  fm_states.fc = states.fc;
 
   fm_states.pre_emp = 0;
   fm_states.de_emp  = 1;
@@ -79,7 +92,7 @@ function run_sim
   fm_states = analog_fm_init(fm_states);
   [b, a] = cheby1(4, 1, 300/Fs, 'high');   % 300Hz HPF to simulate FM radios
 
-  rx_bits_buf = zeros(1,2*nbits);
+  rx_bits_buf = zeros(1,2*nbits2);
   Terrs = Tbits = 0;
   state = 0;
   nerr_log = [];
@@ -117,16 +130,16 @@ function run_sim
   % old-school legacy FM radios.
 
   tx_bits_mapped = zeros(1,length(tx_bits)*2);
-  j = 1;
+  k= 1;
   for i=1:2:length(tx_bits_mapped)
-    if tx_bits(j)
+    if tx_bits(k)
       tx_bits_mapped(i) = 1;
       tx_bits_mapped(i+1) = 0;
     else
       tx_bits_mapped(i) = 0;
       tx_bits_mapped(i+1) = 1;
     end
-    j++;
+    k++;
   end
 
   % use ideal FSK modulator (note: need to try using analog FM modulator)
@@ -134,20 +147,58 @@ function run_sim
   tx = legacyfsk_mod(states, tx_bits_mapped);
   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_sd = filter(ones(1,Ts),1,rx_out_hp);
+  end
+
+  if demod == 2
 
-  % use analog FM demodulator
+    % optimal non-coherent demod at Rs
 
-  [rx_out rx_bb] = analog_fm_demod(fm_states, rx');
-  rx_out_hp = filter(b,a,rx_out);
+    phi1_vec = (1:length(rx))*2*pi*states.f1/Fs;
+    phi2_vec = (1:length(rx))*2*pi*states.f2/Fs;
 
-  % filter using manchester code templates, aka integration or matched filter
+    f1_dc = rx' .* exp(-j*phi1_vec);
+    f2_dc = rx' .* exp(-j*phi2_vec);
+    
+    rx_filt_one = abs(filter(ones(1,Ts),1,f1_dc));
+    rx_filt_zero = abs(filter(ones(1,Ts),1,f2_dc));
 
-  rx_filt_one = filter(manchester_one,1,rx_out_hp);
-  rx_filt_zero = filter(manchester_zero,1,rx_out_hp);
+    rx_sd = rx_filt_one - rx_filt_zero;
+  end
+
+  % From here on code is common to both demods....
+
+  % Estimate fine timing using line at Rs/2 that Manchester encoding provides
+  
+  Np = length(rx_sd);
+  w = 2*pi*(Rs/2)/Fs;
+  x = rx_sd * exp(-j*w*(0:Np-1))';
+  norm_rx_timing = angle(x)/(2*pi)
   
-  rx_timing = 5;
-  rx_filt_one_dec = rx_filt_one(rx_timing:2*Ts:length(rx_filt_one));
-  rx_filt_zero_dec = rx_filt_zero(rx_timing:2*Ts:length(rx_filt_zero));
+  % Sample at ideal timing estimate (rate Rs)
+
+  %norm_rx_timing=0.25
+  rx_timing = round(norm_rx_timing*Ts)
+  rx_sd_sampled = rx_sd(Ts+rx_timing:Ts:Np);
+
+  % Manchester decoding, we combine energy of two bits at Rs to provide
+  % one bit at Rs/2
+
+  rx_dec = rx_sd_sampled(1:2:frames*nbits2) - rx_sd_sampled(2:2:frames*nbits2);
+  rx_bits = rx_dec < 0;
+  %tx_bits(1:20)
+  %rx_bits(1:20)
+
+  % Run frame sync and BER logic over demodulated bits
 
   st = 1;
   for f=1:frames
@@ -156,11 +207,11 @@ function run_sim
 
     nin = nbits;
     en = st + nin - 1;
-    rx_bits = rx_filt_one_dec(st:en) < rx_filt_zero_dec(st:en);
-    st += nin;
 
     rx_bits_buf(1:nbits) = rx_bits_buf(nbits+1:2*nbits);
-    rx_bits_buf(nbits+1:2*nbits) = rx_bits;
+    rx_bits_buf(nbits+1:2*nbits) = rx_bits(st:en);
+
+    st += nin;
 
     % frame sync based on min BER
 
@@ -171,7 +222,7 @@ function run_sim
         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);
+          %printf("i: %d nerrs: %d nerrs_min: %d \n", i, nerrs, nerrs_min);
           if nerrs < nerrs_min
             nerrs_min = nerrs;
             coarse_offset = i;
@@ -202,6 +253,8 @@ function run_sim
 
   % Bunch O'plots --------------------------------------------------------------
 
+  st = 1; en=20;
+
   Tx=fft(tx, Fs);
   TxdB = 20*log10(abs(Tx(1:Fs/2)));
   figure(1)
@@ -212,34 +265,43 @@ function run_sim
 
   figure(2)
   clf
-  subplot(211)
-  st = 1; en=20;
-  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');
+  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
+    plot(rx_filt_one(st:en*states.Ts*2));
+    hold on;
+    plot(rx_filt_zero(st:en*states.Ts*2),'g');
+    hold off;
+  end
 
   figure(3);
   clf;
-  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');
+  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
 
   figure(4);
   clf;
   subplot(211)
-  plot(rx_filt_one(st:en*states.Ts*2) - rx_filt_zero(st:en*states.Ts*2),'+')
-  title('Difference after matched filtering');
+  plot(rx_sd(st*Ts:en*Ts),'+')
+  title('Difference after filtering');
   subplot(212)
-  plot(rx_filt_one_dec(st:en),'+');
+  plot(rx_sd_sampled(st:en),'+');
   hold on;
-  plot(rx_filt_zero_dec(st:en),'g+')
+  plot(rx_sd_sampled(st:en),'g+')
   hold off;
   title('Both channels after sampling at ideal timing instant')