working on test_curves, AWGN cohpsk, dpsk working, HF results need work

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

diff --git a/codec2-dev/octave/cohpsk.m b/codec2-dev/octave/cohpsk.m
index 849efe298a23c98549910d16efc62487aaf291fb..96fbec69c6465fdd419e23bf3f24e628bb2c2370 100644 (file)
--- a/codec2-dev/octave/cohpsk.m
+++ b/codec2-dev/octave/cohpsk.m
@@ -130,7 +130,45 @@ endfunction
 
 % Symbol rate processing for tx side (modulator) -------------------------------------------------------
 
-function [tx_symb tx_bits prev_sym_tx] = bits_to_qpsk_symbols(sim_in, tx_bits, code_param, prev_sym_tx)
+function [tx_symb prev_tx_symb] = bits_to_dqpsk_symbols(sim_in, tx_bits, prev_tx_symb)
+    Nc         = sim_in.Nc;
+    Nsymbrow   = sim_in.Nsymbrow;
+
+    tx_symb = zeros(Nsymbrow,Nc);
+
+    for c=1:Nc
+      for r=1:Nsymbrow
+        i = (c-1)*Nsymbrow + r;
+        tx_symb(r,c) = qpsk_mod(tx_bits(2*(i-1)+1:2*i));  
+        tx_symb(r,c) *= prev_tx_symb(c);
+        prev_tx_symb(c) = tx_symb(r,c);
+      end
+    end
+              
+endfunction
+
+
+function [rx_symb rx_bits rx_symb_linear prev_rx_symb] = dqpsk_symbols_to_bits(sim_in, rx_symb, prev_rx_symb)
+    Nc         = sim_in.Nc;
+    Nsymbrow   = sim_in.Nsymbrow;
+
+    tx_symb = zeros(Nsymbrow,Nc);
+
+    for c=1:Nc
+      for r=1:Nsymbrow
+        tmp = rx_symb(r,c);
+        rx_symb(r,c) *= conj(prev_rx_symb(c))/abs(prev_rx_symb(c));
+        prev_rx_symb(c) = tmp;
+        i = (c-1)*Nsymbrow + r;
+        rx_symb_linear(i) = rx_symb(r,c);
+        rx_bits((2*(i-1)+1):(2*i)) = qpsk_demod(rx_symb(r,c));
+      end
+    end 
+              
+endfunction
+
+
+function [tx_symb tx_bits] = bits_to_qpsk_symbols(sim_in, tx_bits, code_param)
     ldpc_code     = sim_in.ldpc_code;
     rate          = sim_in.ldpc_code_rate;
     framesize     = sim_in.framesize;
@@ -170,17 +208,6 @@ function [tx_symb tx_bits prev_sym_tx] = bits_to_qpsk_symbols(sim_in, tx_bits, c
     for c=Nc+1:Nc:Nc*Nchip
       tx_symb(:,c:c+Nc-1) = tx_symb(:,1:Nc);
     end
-            
-    % Optionally DQPSK encode
- 
-    if strcmp(modulation,'dqpsk')
-      for c=1:Nc*Nchip
-        for r=1:Nsymbrowpilot
-          tx_symb(r,c) *= prev_sym_tx(c);
-          prev_sym_tx(c) = tx_symb(r,c);
-        end
-      end               
-    end
 
     % ensures energy/symbol is normalised when spreading
 
@@ -658,6 +685,10 @@ function sim_out = ber_test(sim_in)
 
     [spread spread_2ms hf_gain] = init_hf_model(Fs, Rs, Nsymbrowpilot*Ntrials);
 
+    if strcmp(modulation,'dqpsk')
+      Nsymbrowpilot = Nsymbrow;
+    end
+
     % Start Simulation ----------------------------------------------------------------
 
     for ne = 1:length(Esvec)
@@ -692,6 +723,7 @@ function sim_out = ber_test(sim_in)
         w_offset     = pi/16;
 
         ct_symb_buf = zeros(2*Nsymbrowpilot, Nc);
+        prev_tx_symb = prev_rx_symb = ones(1,Nc);
 
         % simulation starts here-----------------------------------
  
@@ -703,10 +735,21 @@ function sim_out = ber_test(sim_in)
               tx_bits = round(rand(1,framesize));                       
             end
 
-            [tx_symb tx_bits prev_sym_tx] = bits_to_qpsk_symbols(sim_in, tx_bits, code_param, prev_sym_tx);
-   
-            tx_bits_buf(1:framesize) = tx_bits_buf(framesize+1:2*framesize);
-            tx_bits_buf(framesize+1:2*framesize) = tx_bits;
+            if strcmp(modulation,'qpsk')
+
+              [tx_symb tx_bits] = bits_to_qpsk_symbols(sim_in, tx_bits, code_param);
+
+              % one frame delay on bits for qpsk
+
+              tx_bits_buf(1:framesize) = tx_bits_buf(framesize+1:2*framesize);
+              tx_bits_buf(framesize+1:2*framesize) = tx_bits;
+
+            end
+            if strcmp(modulation, 'dqpsk')
+              [tx_symb prev_tx_symb] = bits_to_dqpsk_symbols(sim_in, tx_bits, prev_tx_symb);
+              tx_bits_buf(1:framesize) = tx_bits;
+            end
+
 
             s_ch = tx_symb;
 
@@ -759,16 +802,20 @@ function sim_out = ber_test(sim_in)
             ct_symb_buf(1:Nsymbrowpilot,:) = ct_symb_buf(Nsymbrowpilot+1:2*Nsymbrowpilot,:);
             ct_symb_buf(Nsymbrowpilot+1:2*Nsymbrowpilot,:) = s_ch;
 
-            [rx_symb rx_bits rx_symb_linear amp_ phi_ EsNo_ sim_in] = qpsk_symbols_to_bits(sim_in, ct_symb_buf(1:Nsymbrowpilot+Npilotsframe,:));                                 
-
-            phi_log = [phi_log; phi_];
-            amp_log = [amp_log; amp_];
-
+            if strcmp(modulation,'qpsk')
+              [rx_symb rx_bits rx_symb_linear amp_ phi_ EsNo_ sim_in] = qpsk_symbols_to_bits(sim_in, ct_symb_buf(1:Nsymbrowpilot+Npilotsframe,:));                                 
+              phi_log = [phi_log; phi_];
+              amp_log = [amp_log; amp_];
+            end
+            if strcmp(modulation,'dqpsk')
+              [rx_symb rx_bits rx_symb_linear prev_rx_symb] = dqpsk_symbols_to_bits(sim_in, s_ch, prev_rx_symb);                                 
+            end
+                        
             % Wait until we have enough frames to do pilot assisted phase estimation
 
             if nn > 1
               rx_symb_log = [rx_symb_log rx_symb_linear];
-              EsNo__log = [EsNo__log EsNo_];
+              %EsNo__log = [EsNo__log EsNo_];
 
               % Measure BER
 
@@ -851,55 +898,57 @@ function sim_out = ber_test(sim_in)
 
        end
 
-        % set up time axis to include gaps for pilots
-
-        [m1 n1] = size(phi_log);
-        phi_x = [];
-        phi_x_counter = 1;
-        p = Ns;
-        for r=1:m1
-          if p == Ns
-            phi_x_counter++;
-            p = 0;
-          end
-          p++;
-          phi_x = [phi_x phi_x_counter++];        
-        end
-
-        phi_x -= Nsymbrowpilot; % account for delay in pilot buffer
-
-        figure(5);
-        clf
-        subplot(211)
-        plot(phi_x, phi_log(:,2),'r+;Estimated HF channel phase;')
-        if hf_sim
-          hold on;
-          [m n] = size(hf_model);
-          plot(angle(hf_model(1:m,2)),'g;HF channel phase;')
-          hold off;
-        end
-        ylabel('Phase (rads)');
-        legend('boxoff');
+       if strcmp(modulation,'qpsk')
+          % set up time axis to include gaps for pilots
+
+         [m1 n1] = size(phi_log);
+         phi_x = [];
+         phi_x_counter = 1;
+         p = Ns;
+         for r=1:m1
+           if p == Ns
+             phi_x_counter++;
+             p = 0;
+           end
+           p++;
+           phi_x = [phi_x phi_x_counter++];        
+         end
+
+         phi_x -= Nsymbrowpilot; % account for delay in pilot buffer
+
+         figure(5);
+         clf
+         subplot(211)
+         plot(phi_x, phi_log(:,2),'r+;Estimated HF channel phase;')
+         if hf_sim
+           hold on;
+           [m n] = size(hf_model);
+           plot(angle(hf_model(1:m,2)),'g;HF channel phase;')
+           hold off;
+         end
+         ylabel('Phase (rads)');
+         legend('boxoff');
+
+         subplot(212)
+         plot(phi_x, amp_log(:,2),'r+;Estimated HF channel amp;')
+         if hf_sim
+           hold on;
+           plot(abs(hf_model(1:m,2)))
+           hold off;
+         end
+         ylabel('Amplitude');
+         xlabel('Time (symbols)');
+         legend('boxoff');
+       end
 
-        subplot(212)
-        plot(phi_x, amp_log(:,2),'r+;Estimated HF channel amp;')
-        if hf_sim
-          hold on;
-          plot(abs(hf_model(1:m,2)))
-          hold off;
-        end
-        ylabel('Amplitude');
-        xlabel('Time (symbols)');
-        legend('boxoff');
-
-        figure(4)
-        clf
-        subplot(211)
-        stem(Nerrs_log)
-        subplot(212)
-        if ldpc_code
-          stem(ldpc_Nerrs_log)
-        end
+       figure(4)
+       clf
+       subplot(211)
+       stem(Nerrs_log)
+       subplot(212)
+       if ldpc_code
+         stem(ldpc_Nerrs_log)
+       end
 
    end
 

diff --git a/codec2-dev/octave/tcohpsk.m b/codec2-dev/octave/tcohpsk.m
index 58e805c23db87ef2ce81dab4540d56313f43ba47..109f125d1d5fd7c051c2d142443e076a40934834 100644 (file)
--- a/codec2-dev/octave/tcohpsk.m
+++ b/codec2-dev/octave/tcohpsk.m
@@ -134,7 +134,7 @@ for i=1:frames
 
   tx_bits_log = [tx_bits_log tx_bits];
 
-  [tx_symb tx_bits prev_tx_sym] = bits_to_qpsk_symbols(acohpsk, tx_bits, [], []);
+  [tx_symb tx_bits] = bits_to_qpsk_symbols(acohpsk, tx_bits, [], []);
   tx_symb_log = [tx_symb_log; tx_symb];
 
   tx_fdm_frame = [];

diff --git a/codec2-dev/octave/test_cohpsk.m b/codec2-dev/octave/test_cohpsk.m
index bdfced8bc8e4967c557bfe38ccff1714ceed2c70..50832bceefa6f9ebdf1cfb454511265e325830ae 100644 (file)
--- a/codec2-dev/octave/test_cohpsk.m
+++ b/codec2-dev/octave/test_cohpsk.m
@@ -27,63 +27,87 @@ graphics_toolkit ("gnuplot");
 cohpsk;
 
 function test_curves
-
+     
   sim_in = standard_init();
+  sim_in.do_write_pilot_file = 0;
+
+  % single test point ---------------------------------------
 
   sim_in.verbose          = 1;
   sim_in.plot_scatter     = 1;
 
   sim_in.Esvec            = 10; 
+  sim_in.framesize        = 32;
   sim_in.hf_sim           = 1;
-  sim_in.Ntrials          = 1000;
+  sim_in.Ntrials          = 100;
   sim_in.Rs               = 50;
-  sim_in.Nc               = 9;
-  sim_in.Np               = 4;
-  sim_in.Ns               = 8;
+  sim_in.Nc               = 4;
+  sim_in.Np               = 2;
+  sim_in.Ns               = 4;
   sim_in.Nchip            = 1;
   sim_in.modulation       = 'qpsk';
-  sim_in.ldpc_code_rate   = 0.5;
+  sim_in.ldpc_code_rate   = 1;
   sim_in.ldpc_code        = 0;
 
   sim_qpsk                = ber_test(sim_in);
 
+  % AWGN curves ----------------------------------------------------
+
+  sim_in.Ntrials          = 500;
   sim_in.hf_sim           = 0;
   sim_in.plot_scatter     = 0;
-  sim_in.Esvec            = 10:20; 
+  sim_in.Esvec            = 5:10; 
+
   Ebvec = sim_in.Esvec - 10*log10(2);
   BER_theory = 0.5*erfc(sqrt(10.^(Ebvec/10)));
-
-  sim_in.Np               = 0;
-  sim_in.Nchip            = 1;
-
+  
   sim_in.modulation       = 'dqpsk';
   sim_dqpsk               = ber_test(sim_in, 'dqpsk');
-  sim_in.hf_sim           = 1;
-  sim_in.hf_mag_only      = 1;
+
   sim_in.modulation       = 'qpsk';
-  sim_qpsk_hf_ideal       = ber_test(sim_in, 'qpsk');
+  sim_in.Ns               = 4;
+  sim_in.Np               = 2;
+  sim_qpsk_pilot          = ber_test(sim_in, 'qpsk');
+
+  % HF curves ----------------------------------------------------
+
+  sim_in.Ntrials          = 200;
+  sim_in.hf_sim           = 1;
+  sim_in.plot_scatter     = 0;
+  sim_in.Esvec            = 5:20; 
+
+  Ebvec = sim_in.Esvec - 10*log10(2);
+  BER_theory = 0.5*erfc(sqrt(10.^(Ebvec/10)));
+  
   sim_in.modulation       = 'dqpsk';
-  sim_in.hf_mag_only      = 0;
   sim_dqpsk_hf            = ber_test(sim_in, 'dqpsk');
+
   sim_in.modulation       = 'qpsk';
-  sim_in.Ns               = 4
+  sim_in.Ns               = 4;
   sim_in.Np               = 2;
-  sim_qpsk_hf_pilot       = ber_test(sim_in, 'qpsk');
+  sim_qpsk_pilot_hf       = ber_test(sim_in, 'qpsk');
+
+  % plot results ---------------------------------------------------
 
   figure(1); 
   clf;
   semilogy(Ebvec, BER_theory,'r;QPSK theory;')
   hold on;
+
   semilogy(sim_dqpsk.Ebvec, sim_dqpsk.BERvec,'c;DQPSK AWGN;')
-  semilogy(sim_qpsk_hf_ideal.Ebvec, sim_qpsk_hf_ideal.BERvec,'b;QPSK HF ideal;')
-  semilogy(sim_dqpsk_hf.Ebvec, sim_dqpsk_hf.BERvec,'k;DQPSK HF;')
-  semilogy(sim_qpsk_hf_pilot.Ebvec, sim_qpsk_hf_pilot.BERvec,'r;QPSK Np=2 Ns=4 HF;')
+  semilogy(sim_qpsk_pilot.Ebvec, sim_qpsk_pilot.BERvec,'b;QPSK pilot AWGN;')
+
+  %semilogy(sim_qpsk_hf_ideal.Ebvec, sim_qpsk_hf_ideal.BERvec,'b;QPSK HF ideal;')
+  semilogy(sim_dqpsk_hf.Ebvec, sim_dqpsk_hf.BERvec,'c;DQPSK HF;')
+  semilogy(sim_qpsk_pilot_hf.Ebvec, sim_qpsk_pilot_hf.BERvec,'b;QPSK pilot HF;')
+
   hold off;
 
   xlabel('Eb/N0')
   ylabel('BER')
   grid("minor")
   axis([min(Ebvec) max(Ebvec) 1E-3 1])
+  legend("boxoff");
 endfunction
 
 
@@ -100,22 +124,22 @@ function test_single
   sim_in.Ns               = 4;
   sim_in.Np               = 2;
   sim_in.Nchip            = 1;
-% sim_in.ldpc_code_rate   = 0.5;
-% sim_in.ldpc_code        = 1;
   sim_in.ldpc_code_rate   = 1;
   sim_in.ldpc_code        = 0;
 
-  sim_in.Ntrials          = 500;
-  sim_in.Esvec            = 8; 
+  sim_in.Ntrials          = 100;
+  sim_in.Esvec            = 10; 
   sim_in.hf_sim           = 0;
   sim_in.hf_mag_only      = 0;
   sim_in.modulation       = 'qpsk';
 
-  sim_in.do_write_pilot_file = 1;
+  sim_in.modulation       = 'dqpsk';
+
+  sim_in.do_write_pilot_file = 0;
 
   sim_qpsk_hf             = ber_test(sim_in);
 
-  fep=fopen("errors_450.bin","wb"); fwrite(fep, sim_qpsk_hf.ldpc_errors_log, "short"); fclose(fep);
+  %fep=fopen("errors_450.bin","wb"); fwrite(fep, sim_qpsk_hf.ldpc_errors_log, "short"); fclose(fep);
 endfunction
 
 
@@ -569,8 +593,8 @@ endfunction
 % Start simulations ---------------------------------------
 
 more off;
-%test_curves();
-test_single();
+test_curves();
+%test_single();
 %rate_Fs_tx("tx_zero.raw");
 %rate_Fs_tx("tx.raw");
 %rate_Fs_rx("tx_-4dB.wav")