added hf/multipath model, and diversity results
authordrowe67 <drowe67@01035d8c-6547-0410-b346-abe4f91aad63>
Thu, 16 Mar 2017 07:23:18 +0000 (07:23 +0000)
committerdrowe67 <drowe67@01035d8c-6547-0410-b346-abe4f91aad63>
Thu, 16 Mar 2017 07:23:18 +0000 (07:23 +0000)
git-svn-id: https://svn.code.sf.net/p/freetel/code@3066 01035d8c-6547-0410-b346-abe4f91aad63

codec2-dev/octave/ldpc_short.m

index a60924d35999ca6cdfb45a1dd38a0c1f49338699..48f62fba212ade292bef0a614114806678950347 100644 (file)
@@ -25,15 +25,21 @@ function init_cml
 end\r
 \r
 \r
-function sim_out = run_sim(HRA, Ntrials, Esvec, genie_Es, packet_size, golay = 0)\r
+function sim_out = run_sim(sim_in, HRA, Ntrials)\r
 \r
-  if golay\r
+  genie_Es    = sim_in.genie_Es;\r
+  packet_size = sim_in.packet_size;\r
+  code        = sim_in.code;\r
+  hf_en       = sim_in.hf_en;\r
+  diversity   = sim_in.diversity;\r
+  Esvec       = sim_in.Esvec;\r
+\r
+  if strcmp(code, 'golay')\r
     rate = 0.5;\r
     code_param.data_bits_per_frame = 12;\r
     framesize = 24;\r
-  else\r
-    % LDPC code\r
-\r
+  end\r
+  if strcmp(code, 'ldpc')\r
     [Nr Nc] = size(HRA);  \r
     rate = (Nc-Nr)/Nc;\r
     framesize = Nc;\r
@@ -47,10 +53,30 @@ function sim_out = run_sim(HRA, Ntrials, Esvec, genie_Es, packet_size, golay = 0
     code_param.P_matrix = [];\r
     code_param.data_bits_per_frame = length(code_param.H_cols) - length( code_param.P_matrix ); \r
   end\r
+  if strcmp(code, 'diversity')\r
+    rate = 0.5;\r
+    code_param.data_bits_per_frame = 7;\r
+    framesize = 14;\r
+  end\r
 \r
   mod_order = 2; \r
   bps = code_param.bits_per_symbol = log2(mod_order);\r
 \r
+  % init HF model\r
+\r
+  if hf_en\r
+\r
+    % some typical values\r
+\r
+    Rs = 50; dopplerSpreadHz = 1.0; path_delay = 1E-3*Rs;\r
+\r
+    nsymb = Ntrials*framesize*bps;\r
+    spread1 = doppler_spread(dopplerSpreadHz, Rs, nsymb);\r
+    spread2 = doppler_spread(dopplerSpreadHz, Rs, nsymb);\r
+    hf_gain = 1.0/sqrt(var(spread1)+var(spread2));\r
+    % printf("nsymb: %d lspread1: %d\n", nsymb, length(spread1));\r
+  end\r
+\r
   % loop around each Esvec point\r
   \r
   for ne = 1:length(Esvec)\r
@@ -60,32 +86,59 @@ function sim_out = run_sim(HRA, Ntrials, Esvec, genie_Es, packet_size, golay = 0
 \r
     Terrs = 0;  Tbits = 0;  Ferrs = 0; \r
     tx_bits = rx_bits = [];\r
+    hfi = 1;\r
 \r
     for nn = 1: Ntrials        \r
       data = round( rand( 1, code_param.data_bits_per_frame ) );\r
       tx_bits = [tx_bits data];\r
-      if golay\r
+      if strcmp(code, 'golay')\r
         codeword = egolayenc(data);\r
-      else\r
+      end\r
+      if strcmp(code, 'ldpc')\r
         codeword = LdpcEncode( data, code_param.H_rows, code_param.P_matrix );\r
       end\r
+      if strcmp(code, 'diversity')\r
+        codeword = [data data];\r
+      end\r
+\r
       code_param.code_bits_per_frame = length( codeword );\r
       Nsymb = code_param.code_bits_per_frame/bps;      \r
        \r
       % modulate\r
+\r
       s = 1 - 2 * codeword;   \r
       code_param.symbols_per_frame = length( s );\r
               \r
+      if hf_en\r
+\r
+        % simplified rate Rs simulation model that doesn't include\r
+        % ISI, just freq filtering.  We assume perfect phase estimation\r
+        % so it's just amplitude distortion.\r
+\r
+        for i=1:length(s)\r
+\r
+          % 1.5*Rs carrier spacing, symbols mapped to 14 carriers\r
+          % OK we'd probably use QPSK in practice but meh a few approximations....\r
+\r
+          w = 1.5*mod(i,14)*2*pi;  \r
+          hf_model(i) = hf_gain*(spread1(hfi) + exp(-j*w*path_delay)*spread2(hfi));\r
+          s(i) *= abs(hf_model(i));\r
+          hfi++;\r
+        end\r
+      end\r
+\r
       variance = 1/(2*EsNo);\r
       noise = sqrt(variance)* randn(1,code_param.symbols_per_frame); \r
       r = s + noise;\r
 \r
       Nr = length(r);  \r
 \r
-      if golay\r
+      if strcmp(code, 'golay')\r
         detected_data = egolaydec(r < 0);\r
         detected_data = detected_data(code_param.data_bits_per_frame+1:framesize);\r
-      else\r
+      end\r
+\r
+      if strcmp(code, 'ldpc')\r
 \r
         % in the binary case the LLRs are just a scaled version of the rx samples ...\r
         % if the EsNo is known -- use the following \r
@@ -106,6 +159,10 @@ function sim_out = run_sim(HRA, Ntrials, Esvec, genie_Es, packet_size, golay = 0
         detected_data = detected_data(1:code_param.data_bits_per_frame);\r
       end\r
 \r
+      if strcmp(code, 'diversity')\r
+        detected_data = (r(1:7) + r(8:14)) < 0;\r
+      end\r
+\r
       rx_bits = [rx_bits detected_data];\r
 \r
       error_positions = xor( detected_data, data );\r
@@ -121,10 +178,10 @@ function sim_out = run_sim(HRA, Ntrials, Esvec, genie_Es, packet_size, golay = 0
     % one big string of bits lets us account for cases were FEC framesize\r
     % is less than packet size\r
 \r
+    error_positions = xor(tx_bits, rx_bits);\r
     Perrs = 0; Tpackets = 0;\r
     for i=1:packet_size:length(tx_bits)-packet_size\r
-      error_positions = xor( tx_bits(i:i+packet_size-1), rx_bits(i:i+packet_size-1));\r
-      Nerrs = sum(error_positions);\r
+      Nerrs = sum(error_positions(i:i+packet_size-1));\r
       if Nerrs>0,  Perrs = Perrs +1;  end\r
       Tpackets++;\r
     end\r
@@ -143,7 +200,88 @@ function sim_out = run_sim(HRA, Ntrials, Esvec, genie_Es, packet_size, golay = 0
     sim_out.FERvec = FERvec;\r
     sim_out.TERvec  = TERvec;\r
     sim_out.framesize = framesize;\r
+    sim_out.error_positions = error_positions;\r
+  end\r
+endfunction\r
+\r
+\r
+function plot_curves(hf_en)\r
+\r
+  Ntrials = 500;\r
+\r
+  sim_in.genie_Es    = 1;\r
+  sim_in.packet_size = 28;\r
+  sim_in.code        = 'ldpc';\r
+  sim_in.hf_en       = hf_en;\r
+  sim_in.diversity   = 0;\r
+\r
+  if hf_en\r
+    Esvec = -3:0.5:6; \r
+  else\r
+    Esvec = 0:0.5:6; \r
+  end\r
+  sim_in.Esvec = Esvec;\r
+\r
+  load HRA_112_112.txt\r
+  load HRA_112_56.txt\r
+  load HRA_56_56.txt\r
+  load HRA_56_28.txt\r
+\r
+  sim_out1 = run_sim(sim_in, HRA_112_112, Ntrials);\r
+  sim_out2 = run_sim(sim_in, HRA_112_56 , Ntrials);\r
+  sim_out3 = run_sim(sim_in, HRA_56_56  , Ntrials*2);\r
+  sim_out4 = run_sim(sim_in, HRA_56_28  , Ntrials*2);\r
+  sim_in.code = 'golay';\r
+  sim_out5 = run_sim(sim_in, [], Ntrials*10);\r
+  sim_in.code = 'diversity';\r
+  sim_out6 = run_sim(sim_in, [], Ntrials*10);\r
+\r
+  if hf_en\r
+    Ebvec_theory = 0.5:9;\r
+    EbNoLin = 10.^(Ebvec_theory/10);\r
+    uncoded_BER_theory = 0.5.*(1-sqrt(EbNoLin./(EbNoLin+1)));\r
+  else\r
+    Ebvec_theory = -2:0.5:6;\r
+    uncoded_BER_theory = 0.5*erfc(sqrt(10.^(Ebvec_theory/10)));\r
   end\r
+\r
+  % need standard packet size to compare\r
+  % packet error if bit 0, or bit 1, or bit 2 .....\r
+  %              or bit 0 and bit 1\r
+  % no packet error if all bits ok (1-p(0))*(1-p(1))\r
+  % P(packet error) = p(0)+p(1)+....\r
+\r
+  uncoded_PER_theory = 1 - (1-uncoded_BER_theory).^sim_in.packet_size;\r
+\r
+  figure(1); clf;\r
+  semilogy(Ebvec_theory,  uncoded_BER_theory, 'b+-;BPSK theory;','markersize', 10, 'linewidth', 2)\r
+  hold on;\r
+  semilogy(sim_out1.Ebvec, sim_out1.BERvec, 'g+-;rate 1/2 HRA 112 112;','markersize', 10, 'linewidth', 2)\r
+  semilogy(sim_out2.Ebvec, sim_out2.BERvec, 'r+-;rate 2/3 HRA 112 56;','markersize', 10, 'linewidth', 2)\r
+  semilogy(sim_out3.Ebvec, sim_out3.BERvec, 'c+-;rate 1/2 HRA 56 56;','markersize', 10, 'linewidth', 2)\r
+  semilogy(sim_out4.Ebvec, sim_out4.BERvec, 'k+-;rate 2/3 HRA 56 28;','markersize', 10, 'linewidth', 2)\r
+  semilogy(sim_out5.Ebvec, sim_out5.BERvec, 'm+-;rate 1/2 Golay (24,12);','markersize', 10, 'linewidth', 2)\r
+  semilogy(sim_out6.Ebvec, sim_out6.BERvec, 'bo-;rate 1/2 Diversity;','markersize', 10, 'linewidth', 2)\r
+  hold off;\r
+  xlabel('Eb/No')\r
+  ylabel('BER')\r
+  grid\r
+  legend("boxoff");\r
+\r
+  figure(2); clf;\r
+  semilogy(Ebvec_theory,  uncoded_PER_theory, 'b+-;BPSK theory;','markersize', 10, 'linewidth', 2)\r
+  hold on;\r
+  semilogy(sim_out1.Ebvec, sim_out1.PERvec, 'g+-;rate 1/2 HRA 112 112;','markersize', 10, 'linewidth', 2)\r
+  semilogy(sim_out2.Ebvec, sim_out2.PERvec, 'r+-;rate 2/3 HRA 112 56;','markersize', 10, 'linewidth', 2)\r
+  semilogy(sim_out3.Ebvec, sim_out3.PERvec, 'c+-;rate 1/2 HRA 56 56;','markersize', 10, 'linewidth', 2)\r
+  semilogy(sim_out4.Ebvec, sim_out4.PERvec, 'k+-;rate 2/3 HRA 56 28;','markersize', 10, 'linewidth', 2)\r
+  semilogy(sim_out5.Ebvec, sim_out5.PERvec, 'm+-;rate 1/2 Golay (24,12);','markersize', 10, 'linewidth', 2)\r
+  semilogy(sim_out6.Ebvec, sim_out6.PERvec, 'bo-;rate 1/2 Diversity;','markersize', 10, 'linewidth', 2)\r
+  hold off;\r
+  xlabel('Eb/No')\r
+  ylabel('PER')\r
+  grid\r
+  legend("boxoff");\r
 endfunction\r
 \r
 % Start simulation here ----------------------------------------------\r
@@ -154,61 +292,6 @@ more off;
 format;\r
 init_cml;\r
 \r
-Ntrials = 500;\r
-Esvec = -3:0.5:3; \r
-packet_size = 28;\r
-\r
-\r
-load HRA_112_112.txt\r
-load HRA_112_56.txt\r
-load HRA_56_56.txt\r
-load HRA_56_28.txt\r
-\r
-sim_out1 = run_sim(HRA_112_112, Ntrials, Esvec, 1, packet_size);\r
-sim_out2 = run_sim(HRA_112_56 , Ntrials, Esvec, 1, packet_size);\r
-sim_out3 = run_sim(HRA_56_56  , Ntrials*2, Esvec, 1, packet_size);\r
-sim_out4 = run_sim(HRA_56_28  , Ntrials*2, Esvec, 1, packet_size);\r
-sim_out5 = run_sim([], Ntrials*10, Esvec, 1, packet_size, 1);\r
-\r
-Ebvec_theory = -2:0.5:6;\r
-uncoded_BER_theory = 0.5*erfc(sqrt(10.^(Ebvec_theory/10)));\r
-\r
-% need standard packet size to compare\r
-% packet error if bit 0, or bit 1, or bit 2 .....\r
-%              or bit 0 and bit 1\r
-% no packet error if all bits ok (1-p(0))*(1-p(1))\r
-% P(packet error) = p(0)+p(1)+....\r
-\r
-uncoded_PER_theory = 1 - (1-uncoded_BER_theory).^packet_size;\r
-\r
-figure(1); clf;\r
-semilogy(Ebvec_theory,  uncoded_BER_theory, 'b+-;BPSK theory;','markersize', 10, 'linewidth', 2)\r
-hold on;\r
-semilogy(sim_out1.Ebvec, sim_out1.BERvec, 'g+-;rate 1/2 HRA 112 112;','markersize', 10, 'linewidth', 2)\r
-semilogy(sim_out2.Ebvec, sim_out2.BERvec, 'r+-;rate 2/3 HRA 112 56;','markersize', 10, 'linewidth', 2)\r
-semilogy(sim_out3.Ebvec, sim_out3.BERvec, 'c+-;rate 1/2 HRA 56 56;','markersize', 10, 'linewidth', 2)\r
-semilogy(sim_out4.Ebvec, sim_out4.BERvec, 'k+-;rate 2/3 HRA 56 28;','markersize', 10, 'linewidth', 2)\r
-semilogy(sim_out5.Ebvec, sim_out5.BERvec, 'm+-;rate 1/2 Golay (24,12);','markersize', 10, 'linewidth', 2)\r
-hold off;\r
-xlabel('Eb/No')\r
-ylabel('BER')\r
-grid\r
-legend("boxoff");\r
-\r
-\r
-figure(2); clf;\r
-semilogy(Ebvec_theory,  uncoded_PER_theory, 'b+-;BPSK theory;','markersize', 10, 'linewidth', 2)\r
-hold on;\r
-semilogy(sim_out1.Ebvec, sim_out1.PERvec, 'g+-;rate 1/2 HRA 112 112;','markersize', 10, 'linewidth', 2)\r
-semilogy(sim_out2.Ebvec, sim_out2.PERvec, 'r+-;rate 2/3 HRA 112 56;','markersize', 10, 'linewidth', 2)\r
-semilogy(sim_out3.Ebvec, sim_out3.PERvec, 'c+-;rate 1/2 HRA 56 56;','markersize', 10, 'linewidth', 2)\r
-semilogy(sim_out4.Ebvec, sim_out4.PERvec, 'k+-;rate 2/3 HRA 56 28;','markersize', 10, 'linewidth', 2)\r
-semilogy(sim_out5.Ebvec, sim_out5.PERvec, 'm+-;rate 1/2 Golay (24,12);','markersize', 10, 'linewidth', 2)\r
-hold off;\r
-xlabel('Eb/No')\r
-ylabel('PER')\r
-grid\r
-legend("boxoff");\r
-\r
+plot_curves('hf');\r
 \r
 \r