added some comments and a PER demo
authordrowe67 <drowe67@01035d8c-6547-0410-b346-abe4f91aad63>
Wed, 12 Oct 2016 00:30:28 +0000 (00:30 +0000)
committerdrowe67 <drowe67@01035d8c-6547-0410-b346-abe4f91aad63>
Wed, 12 Oct 2016 00:30:28 +0000 (00:30 +0000)
git-svn-id: https://svn.code.sf.net/p/freetel/code@2891 01035d8c-6547-0410-b346-abe4f91aad63

codec2-dev/octave/test_ldpc_fsk_lib.m

index af158a58298c462f022c746ab3527e2bb904abec..8f9ddb7e6004bc26e2dc27f8e3f84aed286ea399 100644 (file)
@@ -1,7 +1,9 @@
 % test_ldpc_fsk_lib\r
 % David Rowe 16 April 2016\r
 %\r
-% Some tests for ldpc_fsk_lib, and C versions ldpc_enc and ldpc_dec\r
+% A series of tests for ldpc_fsk_lib, and C versions ldpc_enc and ldpc_dec.\r
+% Gradually builds up complete C command line for SSTV balloon system,\r
+% using Octave versions of LDPC and FSK modem as reference points.\r
 \r
 1;\r
 \r
@@ -57,7 +59,7 @@ function detected_data = decode_from_file(filename)
 end\r
 \r
 \r
-% plots a BER curve for the decoder.  Takes a while to run, uses parallel cores\r
+% plots a BER curve for the LDPC decoder.  Takes a while to run, uses parallel cores\r
 \r
 function plot_curve\r
   num_cores = 4;              % set this to the number of cores you have\r
@@ -243,7 +245,7 @@ function save_rtlsdr(filename, s)
   %iq(2:2:2*l) = 127 + im*(127/mx); \r
   iq(1:2:2*l) = 127 + 32*re; \r
   iq(2:2:2*l) = 127 + 32*im; \r
-  figure(3); clf; plot(iq);\r
+  figure(3); clf; plot(iq); title('simulated IQ signal from RTL SDR');\r
   fs = fopen(filename,"wb");\r
   fwrite(fs,iq,"uchar");\r
   fclose(fs);\r
@@ -273,16 +275,20 @@ endfunction
 \r
 \r
 % Using simulated SSTV packet, generate complex fsk mod signals, 8-bit\r
-% unsigned IQ for feeding into demod chain\r
+% unsigned IQ for feeding into C demod chain.  Can also be used to\r
+% generate BER curves.  Found bugs in UW size and our use of csdr\r
+% re-sampler using this function, and by gradually and carefully\r
+% building up the C command line.\r
+\r
 #{\r
 todo: [X] uncoded BER\r
           [X] octave fsk demod\r
           [X] use C demod\r
-      [ ] compare uncoded BER to unsigned 8 bit IQ to regular 16-bit\r
-          [ ] generate complex rx signal with noise\r
-          [ ] used cmd line utils to drive demod\r
-      [ ] test with resampler\r
-      [ ] measure effect on PER with coding\r
+      [X] compare uncoded BER to unsigned 8 bit IQ to regular 16-bit\r
+          [X] generate complex rx signal with noise\r
+          [X] used cmd line utils to drive demod\r
+      [X] test with resampler\r
+      [X] measure effect on PER with coding\r
 #}\r
 \r
 function [n_uncoded_errs n_uncoded_bits] = run_sstv_sim(sim_in, EbNodB)\r
@@ -303,7 +309,7 @@ function [n_uncoded_errs n_uncoded_bits] = run_sstv_sim(sim_in, EbNodB)
   % note fixed frame of bits used for BER testing\r
 \r
   tx_codeword = gen_sstv_frame;\r
-  length(tx_codeword)\r
+\r
   % init FSK modem\r
 \r
   fsk_horus_as_a_lib = 1;\r
@@ -313,13 +319,13 @@ function [n_uncoded_errs n_uncoded_bits] = run_sstv_sim(sim_in, EbNodB)
   states.dA(1:states.M) = 1;\r
   states.tx_real = 0;  % Octave fsk_mod generates complex valued output\r
 \r
-  % Set up simulated tx tones to sit in teh middle of cdsr passband\r
+  % Set up simulated tx tones to sit in the middle of cdsr passband\r
 \r
   filt_low_norm = 0.1; filt_high_norm = 0.4;\r
   fc = states.Fs*(filt_low_norm + filt_high_norm)/2;\r
   %fc = 1800;\r
-  f1 = fc - states.Rs/2\r
-  f2 = fc + states.Rs/2\r
+  f1 = fc - states.Rs/2;\r
+  f2 = fc + states.Rs/2;\r
   states.ftx = [f1 f2];\r
 \r
   % set up AWGN channel \r
@@ -340,7 +346,7 @@ function [n_uncoded_errs n_uncoded_bits] = run_sstv_sim(sim_in, EbNodB)
   noise_real = sqrt(variance)*randn(length(tx),1);\r
   noise_complex = sqrt(variance/2)*(randn(length(tx),1) + j*randn(length(tx),1));\r
 \r
-  % demodulate\r
+  % demodulate -----------------------------------------------------\r
 \r
   if demod_type == 1\r
 \r
@@ -396,7 +402,7 @@ function [n_uncoded_errs n_uncoded_bits] = run_sstv_sim(sim_in, EbNodB)
   end\r
 \r
   if demod_type == 3\r
-    % C demod driven by command line kung fu\r
+    % C demod driven by csdr command line kung fu\r
 \r
     assert(states.tx_real == 0, "need complex signal for this test");\r
     rx = tx + noise_complex;\r
@@ -413,8 +419,6 @@ function [n_uncoded_errs n_uncoded_bits] = run_sstv_sim(sim_in, EbNodB)
     rx = tx + noise_complex;\r
     SNRdB = 10*log10(var(tx)/var(noise_real));\r
     \r
-    % interpolate by factor or 2, then use straight line interpolation\r
-\r
     printf("resampling ...\n");\r
     rx_resample_fract = fractional_resample(rx, 1.08331);\r
     %rx_resample_fract = fractional_resample(rx_resample_fract, 1/1.08331);\r
@@ -423,71 +427,96 @@ function [n_uncoded_errs n_uncoded_bits] = run_sstv_sim(sim_in, EbNodB)
     printf("run C cmd line chain ...\n");\r
 %    system("cat fsk_demod_resample.iq | csdr convert_u8_f | csdr bandpass_fir_fft_cc 0.1 0.4 0.05 | csdr realpart_cf | csdr convert_f_s16 | ../build_linux/src/fsk_demod 2X 8 9600 1200 - fsk_demod.bin");\r
     system("cat fsk_demod_resample.iq | csdr convert_u8_f | csdr bandpass_fir_fft_cc 0.1 0.4 0.05 | csdr realpart_cf | csdr convert_f_s16 | ../unittest/tsrc - - 0.9230968 | ../build_linux/src/fsk_demod 2X 8 9600 1200 - fsk_demod.bin");\r
-%    system("cat fsk_demod_resample.iq | csdr convert_u8_f | csdr bandpass_fir_fft_cc 0.1 0.4 0.05 | csdr fractional_decimator_ff 1.08331 | csdr realpart_cf | csdr convert_f_s16 | ../build_linux/src/fsk_demod 2X 8 9600 1200 - fsk_demod.bin");\r
+%    system("cat fsk_demod_resample.iq | csdr convert_u8_f | csdr bandpass_fir_fft_cc 0.1 0.4 0.05 | csdr realpart_cf | csdr fractional_decimator_ff 1.08331 | csdr convert_f_s16 | ../build_linux/src/fsk_demod 2X 8 9600 1200 - fsk_demod.bin");\r
     f = fopen("fsk_demod.bin","rb"); rx_bit_stream = fread(f, "uint8")'; fclose(f);\r
   end\r
 \r
- % state machine. Look for SSTV UW.  When found count bit errors over one frame of bits\r
 \r
-  state = "wait for uw";\r
-  start_uw_ind = 16*10+1; end_uw_ind = start_uw_ind + 5*10 - 1;\r
-  uw_rs232 = tx_codeword(start_uw_ind:end_uw_ind); luw = length(uw_rs232);\r
-  start_frame_ind =  end_uw_ind + 1;\r
-  nbits = length(rx_bit_stream);\r
-  uw_thresh = 5;\r
-  n_uncoded_errs = 0;\r
-  n_uncoded_bits = 0;\r
-  n_packets_rx = 0;\r
-  last_i = 0;\r
+  if demod_type == 5\r
+\r
+    % C demod with resampler and use C code to measure PER, in this\r
+    % test we don't need to run state machine below as C code gives us\r
+    % the ouputs we need\r
 \r
-  % might as well include RS232 framing bits in uncoded error count\r
+    assert(states.tx_real == 0, "need complex signal for this test");\r
+    rx = tx + noise_complex;\r
+    SNRdB = 10*log10(var(tx)/var(noise_real));\r
+    \r
+    printf("resampling ...\n");\r
+    rx_resample_fract = fractional_resample(rx, 1.08331);\r
+    save_rtlsdr("fsk_demod_resample.iq", rx_resample_fract);\r
 \r
-  nbits_frame = code_param.data_bits_per_frame*10/8;  \r
+    printf("run C cmd line chain - uncoded PER\n");\r
+    system("cat fsk_demod_resample.iq | csdr convert_u8_f | csdr bandpass_fir_fft_cc 0.1 0.4 0.05 | csdr realpart_cf | csdr convert_f_s16 | ../unittest/tsrc - - 0.9230968 | ../build_linux/src/fsk_demod 2X 8 9600 1200 - - | ../src/drs232 - /dev/null -v");\r
 \r
-  uw_errs = zeros(1, nbits);\r
-  for i=luw:nbits\r
-    uw_errs(i) = sum(xor(rx_bit_stream(i-luw+1:i), uw_rs232));\r
+    printf("run C cmd line chain - LDPC coded PER\n");\r
+    system("cat fsk_demod_resample.iq | csdr convert_u8_f | csdr bandpass_fir_fft_cc 0.1 0.4 0.05 | csdr realpart_cf | csdr convert_f_s16 | ../unittest/tsrc - - 0.9230968 | ../build_linux/src/fsk_demod 2XS 8 9600 1200 - - | ../src/drs232_ldpc - /dev/null -v");\r
   end\r
 \r
-  for i=luw:nbits\r
-    next_state = state;\r
-    if strcmp(state, 'wait for uw')\r
-      if uw_errs(i) <= uw_thresh\r
-        next_state = 'count errors';\r
-        tx_frame_ind = start_frame_ind;\r
-        rx_frame_ind = i + 1;\r
-        n_uncoded_errs_this_frame = 0;\r
-        %printf("%d %s %s\n", i, state, next_state);\r
-        if last_i\r
-          printf("i: %d i-last_i: %d ", i, i-last_i);\r
+  if demod_type != 5\r
+    % state machine. Look for SSTV UW.  When found count bit errors over one frame of bits\r
+\r
+    state = "wait for uw";\r
+    start_uw_ind = 16*10+1; end_uw_ind = start_uw_ind + 5*10 - 1;\r
+    uw_rs232 = tx_codeword(start_uw_ind:end_uw_ind); luw = length(uw_rs232);\r
+    start_frame_ind =  end_uw_ind + 1;\r
+    nbits = length(rx_bit_stream);\r
+    uw_thresh = 5;\r
+    n_uncoded_errs = 0;\r
+    n_uncoded_bits = 0;\r
+    n_packets_rx = 0;\r
+    last_i = 0;\r
+\r
+    % might as well include RS232 framing bits in uncoded error count\r
+\r
+    nbits_frame = code_param.data_bits_per_frame*10/8;  \r
+\r
+    uw_errs = zeros(1, nbits);\r
+    for i=luw:nbits\r
+      uw_errs(i) = sum(xor(rx_bit_stream(i-luw+1:i), uw_rs232));\r
+    end\r
+\r
+    for i=luw:nbits\r
+      next_state = state;\r
+      if strcmp(state, 'wait for uw')\r
+        if uw_errs(i) <= uw_thresh\r
+          next_state = 'count errors';\r
+          tx_frame_ind = start_frame_ind;\r
+          rx_frame_ind = i + 1;\r
+          n_uncoded_errs_this_frame = 0;\r
+          %printf("%d %s %s\n", i, state, next_state);\r
+          if last_i\r
+            printf("i: %d i-last_i: %d ", i, i-last_i);\r
+          end\r
         end\r
       end\r
-    end\r
-    if strcmp(state, 'count errors')\r
-      n_uncoded_errs_this_frame += xor(rx_bit_stream(i), tx_codeword(tx_frame_ind));\r
-      n_uncoded_bits++;\r
-      tx_frame_ind++;\r
-      if tx_frame_ind == (start_frame_ind+nbits_frame)\r
-        n_uncoded_errs += n_uncoded_errs_this_frame;\r
-        printf("n_uncoded_errs_this_frame: %d\n", n_uncoded_errs_this_frame);\r
-        frame_rx232_rx = rx_bit_stream(rx_frame_ind:rx_frame_ind+nbits_frame-1);\r
-        %tx_codeword(start_frame_ind+1:start_frame_ind+10)\r
-        %frame_rx232_rx(1:10)\r
-        sstv_checksum(frame_rx232_rx);\r
-        last_i = i;\r
-        n_packets_rx++;\r
-        next_state = 'wait for uw';\r
+      if strcmp(state, 'count errors')\r
+        n_uncoded_errs_this_frame += xor(rx_bit_stream(i), tx_codeword(tx_frame_ind));\r
+        n_uncoded_bits++;\r
+        tx_frame_ind++;\r
+        if tx_frame_ind == (start_frame_ind+nbits_frame)\r
+          n_uncoded_errs += n_uncoded_errs_this_frame;\r
+          printf("n_uncoded_errs_this_frame: %d\n", n_uncoded_errs_this_frame);\r
+          frame_rx232_rx = rx_bit_stream(rx_frame_ind:rx_frame_ind+nbits_frame-1);\r
+          %tx_codeword(start_frame_ind+1:start_frame_ind+10)\r
+          %frame_rx232_rx(1:10)\r
+          sstv_checksum(frame_rx232_rx);\r
+          last_i = i;\r
+          n_packets_rx++;\r
+          next_state = 'wait for uw';\r
+        end\r
       end\r
+      state = next_state;\r
     end\r
-    state = next_state;\r
-  end\r
 \r
-  uncoded_ber = n_uncoded_errs/n_uncoded_bits;\r
-  printf("EbNodB: %4.1f SNRdB: %4.1f pkts: %d bits: %d errs: %d BER: %4.3f\n", \r
-          EbNodB, SNRdB, n_packets_rx, n_uncoded_bits, n_uncoded_errs, uncoded_ber);  \r
+    uncoded_ber = n_uncoded_errs/n_uncoded_bits;\r
+    printf("EbNodB: %4.1f SNRdB: %4.1f pkts: %d bits: %d errs: %d BER: %4.3f\n", \r
+            EbNodB, SNRdB, n_packets_rx, n_uncoded_bits, n_uncoded_errs, uncoded_ber);  \r
 \r
-  figure(2);\r
-  plot(uw_errs);\r
+    figure(2);\r
+    plot(uw_errs);\r
+    title('Unique Word Hamming Distance')\r
+  end\r
 \r
 endfunction\r
 \r
@@ -524,7 +553,7 @@ rand('state',1);
 \r
 % ------------------ select which demo/test to run here ---------------\r
 \r
-demo = 10;\r
+demo = 11;\r
 \r
 if demo == 1\r
   printf("simple_ut....\n");\r
@@ -577,6 +606,8 @@ if demo == 7
   test_c_decoder;\r
 end\r
 \r
+% generates simulated demod soft decision symbols to drive C ldpc decoder with\r
+\r
 if demo == 8\r
   frames = 100;\r
   EsNodB = 3;\r
@@ -631,35 +662,49 @@ if demo == 9
 end\r
 \r
 \r
+% Plots uncoded BER curves for two different SSTV simulations.  Used\r
+% to compare results with different processing steps as we build up C\r
+% command line.  BER curves are powerful ways to confirm system is\r
+% operating as expected, several bugs were found using this system.\r
+\r
 if demo == 10\r
-  sim_in.frames = 100;\r
-  EbNodBvec = 7;\r
+  sim_in.frames = 10;\r
+  EbNodBvec = 7:10;\r
 \r
-  sim_in.demod_type = 4;\r
-  ber_octave = [];\r
+  sim_in.demod_type = 3;\r
+  ber_test1 = [];\r
   for i = 1:length(EbNodBvec)\r
     [n_uncoded_errs n_uncoded_bits] = run_sstv_sim(sim_in, EbNodBvec(i));\r
-    ber_octave(i) = n_uncoded_errs/n_uncoded_bits;\r
+    ber_test1(i) = n_uncoded_errs/n_uncoded_bits;\r
   end\r
-  #{\r
+  \r
   sim_in.demod_type = 4;\r
   ber_c = [];\r
   for i = 1:length(EbNodBvec)\r
     [n_uncoded_errs n_uncoded_bits] = run_sstv_sim(sim_in, EbNodBvec(i));\r
-    ber_c(i) = n_uncoded_errs/n_uncoded_bits;\r
+    ber_test2(i) = n_uncoded_errs/n_uncoded_bits;\r
   end\r
 \r
   figure(1);\r
   clf;\r
-  semilogy(EbNodBvec,  ber_octave, '+-;first test;')\r
+  semilogy(EbNodBvec,  ber_test1, '+-;first test;')\r
   grid;\r
   xlabel('Eb/No (dB)')\r
   ylabel('BER')\r
 \r
   hold on;\r
-  semilogy(EbNodBvec,  ber_c, 'g+-;second test;')\r
+  semilogy(EbNodBvec,  ber_test2, 'g+-;second test;')\r
   legend("boxoff");\r
   hold off;\r
-  #}\r
+  \r
 end\r
 \r
+\r
+% Measure PER of complete coded and uncoded system\r
+\r
+if demo == 11\r
+  sim_in.frames = 100;\r
+  EbNodB = 13;\r
+  sim_in.demod_type = 5;\r
+  run_sstv_sim(sim_in, EbNodB);\r
+end\r