extracting and plotting some more stats

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

diff --git a/codec2-dev/octave/cellmodem.m b/codec2-dev/octave/cellmodem.m
index ac0e0c3ba09e97c22ea83554e2186fc9c2f95a02..370053993d749c64a711ce046a46742711e01383 100644 (file)
--- a/codec2-dev/octave/cellmodem.m
+++ b/codec2-dev/octave/cellmodem.m
@@ -20,67 +20,118 @@
 %     symbol set.  Maybe 50 or 100 Hz grid for LSPs, evaluate that some how.
 
 graphics_toolkit ("gnuplot");
+rand('state',1);
 lsp;
 
-codec_cmd = "speexenc --bitrate 4000 in.raw - | speexdec - out.raw"
+% given a vector of LSP symbols, constructs a synthesised speech signals
+% anfd runs it through a codec
 
-lpc_order = 4;
-Fs = 8000;
-fo = 200;             % pitch frequency of voice 
-wo = 2*pi*fo/Fs;      % pitch in rads
-N  = Fs*0.04;         % frame length
-frames = Fs/N;        % frames to play through codec
-gain = 100;
+function [w__log mse] = run_sim(sim_in, w_log)
+  N         = sim_in.N;
+  Wo        = sim_in.Wo;
+  frames    = sim_in.frames;
+  lpc_order = sim_in.lpc_order;
 
-s = [];
-w_log = [];
-w__log = [];
+  s      = [];
+  w__log = [];
+  L      = floor(pi/Wo);
+  phi    = zeros(1,L);
 
-for f=1:frames
+  for f=1:frames
 
-  % construct LSP symbol
+    % synthesise speech signal
 
-  w = [0.1 0.4  0.5 0.8]*pi;
+    a=lsptoa(w_log(f,:));
 
-  % synthesise speech signal
+    ex = zeros(1,N);
+    for m=1:L
+      phi(m) += Wo*m*N;
+      ex += cos(phi(m) + Wo*m*(0:N-1));
+    end
 
-  a=lsptoa(w);
-  w_log = [w_log; w];
-  l=pi/wo;
-  ex = zeros(1,N);
-  for m=1:l
-    ex += cos(wo*m*(1:N));
+    s = [s filter(1, a, ex)];
   end
-  s = [s filter(1, a, ex)];
+
+  % play through codec
+
+  s *= sim_in.gain;
+  f=fopen("in.raw","wb");
+  fwrite(f,s,"short");
+  fclose(f);
+  system(sim_in.codec_cmd);
+  f=fopen("out.raw","rb");
+  s_ = fread(f,Inf,"short");
+  fclose(f);
+
+  % extract received symbols from channel and evaluate
+
+  mse = zeros(1,frames);
+  for f=1:frames
+    a_ = lpcauto(s_((f-1)*sim_in.N+1:f*N), lpc_order);
+    w_ = atolsp(a_);
+    w__log = [w__log; w_];
+    error = w__log(f,:) - w_log(f,:);
+    mse(f) = error*error';
+  end
+endfunction
+
+% constants -----------------------------------------------------
+
+sim_in.codec_cmd = "speexenc --bitrate 4000 in.raw - | speexdec - out.raw";
+
+lpc_order   = sim_in.lpc_order = 6;
+Fs          = sim_in.Fs = 8000;
+              sim_in.Fo = 100;             % pitch frequency of voice 
+              sim_in.Wo = 2*pi*Fo/Fs;      % pitch in rads
+N           = sim_in.N  = Fs*0.04;         % frame length
+frames      = sim_in.frames = 1000;        % frames to play through codec
+              sim_in.gain = 100;
+Nsym        = 8;                           % number of symbols to find
+
+% start with some LSP random vectors
+% for stable filter most be monotonically increasing on 0..pi
+
+w_log = [];
+for f=1:frames;
+  w = sort(rand(1,lpc_order)*pi);
+  w_log = [w_log; w];
 end
+[w__log mse] = run_sim(sim_in, w_log);
+
+% sort by MSE to get the best symbols
 
-% play through codec
+[sort_mse sort_ind] = sort(mse);
+symbols = w_log(sort_ind(1:Nsym),:)
+
+% Play these symbols through the codec in random order
+
+w_log = [];
+symb_ind = [];
+for f=1:frames
+  symb_ind(f) = floor(1 + rand(1,1)*Nsym);
+  w_log = [w_log; symbols(symb_ind,:)];
+end
 
-s *= gain;
-f=fopen("in.raw","wb");
-fwrite(f,s,"short");
-fclose(f);
-system(codec_cmd);
-f=fopen("out.raw","rb");
-s_ = fread(f,Inf,"short");
-fclose(f);
+[w__log mse] = run_sim(sim_in, w_log);
 
-% extract received symbols from channel and evaluate
+% now see if we can "detect" them
 
 for f=1:frames
-  a_ = lpcauto(s_((f-1)*N+1:f*N), lpc_order);
-  w_ = atolsp(a_);
-  w__log = [w__log; w_];
+
+  % check received symbol against codebook of symbols
+
+  min_e = 1E6;
+  for i=1:Nsym
+    e = w__log(f,:)*symbols(i,:)';
+    if e < min_e
+      min_e = e;
+      min_ind = i;
+    end
+  end
+
+  symb_ind_out(f) = min_ind;
 end
 
-figure(1);
-clf;
-subplot(211)
-plot(s)
-axis([1 Fs/10 -4000 4000]);
-subplot(212)
-plot(s_);
-axis([1 Fs/10 -4000 4000]);
-
-figure(2)
-plot(w_log, 'g', w__log, 'r+')
+figure(1)
+clf
+plot(symb_ind,'g',symb_ind_out,'r');

diff --git a/codec2-dev/octave/cohpsk_demod_plot.m b/codec2-dev/octave/cohpsk_demod_plot.m
index fe0c78a2a6ffa0d230e1c7f4ff7d5fd66305076a..cbe5ca55094c9a777ef9845e49d7aeb4d34f8c2f 100644 (file)
--- a/codec2-dev/octave/cohpsk_demod_plot.m
+++ b/codec2-dev/octave/cohpsk_demod_plot.m
@@ -47,3 +47,10 @@ figure(5);
 clf;
 plot(error_positions_hist_c);
 title('Error Position Histogram');
+
+figure(6)
+y = 1:length(rx_amp_log_c);
+x = 1:Nc*Nd;
+mesh(x,y,20*log10(rx_amp_log_c));
+grid
+title('Channel Amplitude dB');

diff --git a/codec2-dev/octave/tcohpsk.m b/codec2-dev/octave/tcohpsk.m
index 9f565a5652fe43af4b71699b2423840e0f2422f7..8ffeb86685168f8dcea7f42fd6919964648a58f2 100644 (file)
--- a/codec2-dev/octave/tcohpsk.m
+++ b/codec2-dev/octave/tcohpsk.m
@@ -44,8 +44,9 @@
 %  [ ] pilot based EsNo estimation
 %  [ ] filter reqd with compression?
 %      + make sure not too much noise passed into noise floor
-%  [ ] different diversity combination
-%  [ ] histogram of bit errors
+%  [X] different diversity combination
+%      + taking largest symbol didn't help
+%  [X] histogram of bit errors
 %      + lot of data
 %      + ssb filter
 %      + compression

diff --git a/codec2-dev/src/cohpsk.c b/codec2-dev/src/cohpsk.c
index 2c880f81a8b76a4c0602cea2b928302de47b87f1..0802de0df989463ad94265c46f094b892d965d66 100644 (file)
--- a/codec2-dev/src/cohpsk.c
+++ b/codec2-dev/src/cohpsk.c
@@ -581,6 +581,8 @@ int sync_state_machine(struct COHPSK *coh, int sync, int next_sync)
         /* check that sync is still good, fall out of sync on consecutive bad frames */
 
         corr_with_pilots(&corr, &mag, coh, coh->ct, coh->f_fine_est);
+        coh->ratio = fabsf(corr)/mag;        
+
         // printf("%f\n", cabsolute(corr)/mag);
 
         if (fabsf(corr)/mag < 0.8) 

diff --git a/codec2-dev/src/cohpsk_ch.c b/codec2-dev/src/cohpsk_ch.c
index 16c3791bf3dc7fcf13014b0bfef360aa683a8f9e..539de49864820cefcdcc2ff29471f28fed7258cc 100644 (file)
--- a/codec2-dev/src/cohpsk_ch.c
+++ b/codec2-dev/src/cohpsk_ch.c
@@ -107,7 +107,7 @@ int main(int argc, char *argv[])
         inclip = atof(argv[6]);
     }
     else {
-        fprintf(stderr, "usage: %s InputRealModemRawFileFs7500Hz OutputRealModemRawFileFs7500Hz SNR3000Hz FoffHz FadingEn InputClip0to1\n", argv[0]);
+        fprintf(stderr, "usage: %s InputRealModemRawFileFs7500Hz OutputRealModemRawFileFs7500Hz No(dB/Hz) FoffHz FadingEn InputClip0to1\n", argv[0]);
         exit(1);
     }
     fprintf(stderr, "NodB: %4.2f foff: %4.2f Hz fading: %d inclip: %4.2f\n", NodB, foff_hz, fading_en, inclip);

diff --git a/codec2-dev/src/cohpsk_demod.c b/codec2-dev/src/cohpsk_demod.c
index 47fc42f1ba90f82b4e58daf413636caf3ae06bee..8a483049c8b5c33ae1daa781bc70e2e06693dff5 100644 (file)
--- a/codec2-dev/src/cohpsk_demod.c
+++ b/codec2-dev/src/cohpsk_demod.c
@@ -130,6 +130,7 @@ int main(int argc, char *argv[])
 
                 f_est_log[frames-1] = cohpsk->f_est;
                 ratio_log[frames-1] = cohpsk->ratio;
+                //fprintf(stderr,"ratio: %f\n", cohpsk->ratio);
 
                 //printf("frames: %d log_data_r: %d\n", frames, log_data_r);
                 if (frames == logframes)

diff --git a/codec2-dev/src/cohpsk_put_test_bits.c b/codec2-dev/src/cohpsk_put_test_bits.c
index 086c85cf61353f45f9468d67a7f584e8e3eba6bf..21e8a6ee78a169da15edb0d190826ae82b0485f7 100644 (file)
--- a/codec2-dev/src/cohpsk_put_test_bits.c
+++ b/codec2-dev/src/cohpsk_put_test_bits.c
@@ -126,7 +126,7 @@ int main(int argc, char *argv[])
 
     if (foct != NULL) {
         octave_save_int(foct, "nerr_log_c", nerr_log, 1, logframes);  
-        octave_save_int(foct, "error_positions_hist_c", error_positions_hist, 1, logframes);  
+        octave_save_int(foct, "error_positions_hist_c", error_positions_hist, 1, COHPSK_BITS_PER_FRAME);  
         fclose(foct);
     }