From 3f8b4720172cfe8f67a535b2317be61eeddc33ab Mon Sep 17 00:00:00 2001
From: drowe67 <drowe67@01035d8c-6547-0410-b346-abe4f91aad63>
Date: Wed, 14 Jun 2017 01:15:46 +0000
Subject: [PATCH] first pass at quantising dct coeffs using uniform quantisers,
 sounds reasonable ... something like 170 bits which is a bit off the pace but
 a gd start.  Prob worth trying kmeans again just in case there was a bug

git-svn-id: https://svn.code.sf.net/p/freetel/code@3189 01035d8c-6547-0410-b346-abe4f91aad63
---
 codec2-dev/octave/newamp.m        | 49 +++++++++++++++++++++++++++++--
 codec2-dev/octave/newamp1_batch.m | 13 ++++----
 2 files changed, 52 insertions(+), 10 deletions(-)

diff --git a/codec2-dev/octave/newamp.m b/codec2-dev/octave/newamp.m
index 1af8a4ac..a4580d8a 100644
--- a/codec2-dev/octave/newamp.m
+++ b/codec2-dev/octave/newamp.m
@@ -281,7 +281,7 @@ endfunction
 % (ie std dev) of each coeff.  Those coeffs with the greatest change
 % need the most bits to quantise
 
-function [map rms_map mx unwrapped_dcts] = create_map_rms(rate_K_surface, nr, nc)
+function [map rms_map mx mx_ind unwrapped_dcts] = create_map_rms(rate_K_surface, nr, nc)
   unwrapped_dcts = dct_blocks(rate_K_surface);
   [mx mx_ind] = sort(std(unwrapped_dcts));
   mx_ind = fliplr(mx_ind); mx = fliplr(mx);
@@ -296,9 +296,33 @@ function [map rms_map mx unwrapped_dcts] = create_map_rms(rate_K_surface, nr, nc
 endfunction
 
 
-% Design quantisers for each DCT coeff
+% plot histogram of each 2D DCT coeff, so we can get a feel for
+% quantiser design
+
+function plot_dct2_hists(rate_K_surface, nr, nc)
+  [map rms_map mx mx_ind unwrapped_dcts] = create_map_rms(rate_K_surface, nr, nc);
+  Ncoeff = nr*nc;
+  fign = 1; subplotn = 1;
+  close all; figure(fign); clf;
+  Nplot = 60;
+  for i=1:Nplot   
+    subplot(5,4,subplotn);
+    d = unwrapped_dcts(:,mx_ind(i));
+    d = round(d);
+    hist(d);
+    subplotn++;
+    if (subplotn > 20) && (i != Nplot)
+      subplotn = 1;
+      fign++;
+      figure(fign); clf;
+    end
+  end
+endfunction
+
 
-function [quantisers nbits] = design_quantisters(rate_K_surface, nr, nc, nbits_max)
+% Design kmeans quantisers for each DCT coeff.  This didn't work very well.
+
+function [quantisers nbits] = design_quantisters_kmeans(rate_K_surface, nr, nc, nbits_max)
   [map rms_map mx unwrapped_dcts] = create_map_rms(rate_K_surface, nr, nc);
   nq = nr*nc;
   quantisers = zeros(nq, 2^nbits_max); nbits = zeros(nq,1);
@@ -320,6 +344,25 @@ function [quantisers nbits] = design_quantisters(rate_K_surface, nr, nc, nbits_m
 endfunction
 
 
+% Uniform quantisers designed to fit limits of each DCT coeff
+
+function [quantisers nlevels] = design_quantisters_uniform(rate_K_surface, nr, nc, nlevels_max)
+  [map rms_map mx mx_ind unwrapped_dcts] = create_map_rms(rate_K_surface, nr, nc);
+  
+  nq = nr*nc;
+  quantisers = zeros(nq, nlevels_max); nlevels = zeros(nq, 1);
+
+  for i=1:nq
+    d = unwrapped_dcts(:,mx_ind(i));
+    d = floor(d/16);
+    q_min = floor(min(d));
+    q_max = ceil(max(d));
+    nlevels(i) = q_max-q_min+1;
+    quantisers(i,1:nlevels(i)) = 16*(q_min:q_max);
+  end
+endfunction
+
+
 % Determine a phase spectra from a magnitude spectra
 % from http://www.dsprelated.com/showcode/20.php
 % Haven't _quite_ figured out how this works but have to start somewhere ....
diff --git a/codec2-dev/octave/newamp1_batch.m b/codec2-dev/octave/newamp1_batch.m
index 9901f7df..82d25c8d 100644
--- a/codec2-dev/octave/newamp1_batch.m
+++ b/codec2-dev/octave/newamp1_batch.m
@@ -500,7 +500,7 @@ function [model_ rate_K_surface] = experiment_mel_freq(model, vq_en=0, plots=1,
   K = 20; Fs = 8000; correct_rate_K_en = 1;
 
   quantisers = load("dct2quant.txt"); 
-  nbits = load("dct2quantnbits.txt");
+  nlevels = load("dct2quantnlevels.txt");
   map = load("dct2map.txt");
 
   for f=1:frames
@@ -560,16 +560,15 @@ function [model_ rate_K_surface] = experiment_mel_freq(model, vq_en=0, plots=1,
     % quantise, replace on map
 
     E = zeros(Nt,K);
-    E = D;
     for r=1:Nt
       for c=1:K
         quantiser_num = map(r,c);
-        nlevels = 2^nbits(quantiser_num);
-        if quantiser_num < 10
-          printf("r %d c %d quantiser_num %d nlevels %d\n", r,c,quantiser_num,nlevels);
-          levels = quantisers(quantiser_num,1:nlevels);
+        if quantiser_num < 100
+          printf("r %d c %d quantiser_num %d nlevels %d ", r, c, quantiser_num, nlevels(quantiser_num));
+          levels = quantisers(quantiser_num, 1:nlevels(quantiser_num));
           quant_out = quantise(levels, D(r,c));
-          printf("%f %f\n", D(r,c), quant_out );
+          %quant_out = 8*round(D(r,c)/8);
+          printf("in %f out %f\n", D(r,c), quant_out );
           E(r,c) = quant_out;
         end
       end
-- 
2.25.1