simple unquantised three resonator model
authordrowe67 <drowe67@01035d8c-6547-0410-b346-abe4f91aad63>
Mon, 5 Dec 2016 04:52:04 +0000 (04:52 +0000)
committerdrowe67 <drowe67@01035d8c-6547-0410-b346-abe4f91aad63>
Mon, 5 Dec 2016 04:52:04 +0000 (04:52 +0000)
git-svn-id: https://svn.code.sf.net/p/freetel/code@2925 01035d8c-6547-0410-b346-abe4f91aad63

codec2-dev/octave/newamp.m
codec2-dev/octave/newamp1_batch.m [new file with mode: 0644]
codec2-dev/octave/newamp1_fbf.m [new file with mode: 0644]

index 1a4cd9a41c57f1d0b0def8b3b7cb6f1ce7f395ea..54b66c23f15db07e72d2415e046735db26616e59 100644 (file)
@@ -861,31 +861,50 @@ function [v] = est_voicing_bits(v1, v5)
 endfunction
 
 
-function AmdB_ = piecewise_model(AmdB, Wo) 
+function [AmdB_ residual] = piecewise_model(AmdB, Wo) 
     L = length(AmdB);
     l1000 = floor(L/4);     
-    l2000 = floor(2*L/4);     
-    l3000 = floor(3*L/4);     
     AmdB_ = ones(1,L);
+    mask_sample_freqs_kHz = (1:L)*Wo*4/pi;
+
+    % fit a resonator to max of first 300 - 1000 Hz
 
-    % fit a resonator to max of first 1000 Hz
+    fmin = 0.150;
+    lmin = floor(L*fmin/4);
+    [mx mx_ind] = max(AmdB(lmin+1:l1000));
+    mx_ind += lmin;
+    AmdB_ = parabolic_resonator(mx_ind*Wo*4/pi, mask_sample_freqs_kHz) + mx;  
 
-    [mx mx_ind] = max(AmdB(1:l1000));
-    mask_sample_freqs_kHz = (1:l1000)*Wo*4/pi;
-    AmdB_(1:l1000) = parabolic_resonator(mx_ind*Wo*4/pi, mask_sample_freqs_kHz) + mx;  
+    % fit a 2nd resonator, must be above 1000Hz
 
-    % fit a resonator to max of second 1000 Hz
+    fr1 = mx_ind*Wo*4/pi;
+    fmin = 1;
+    lmin = round(L*fmin/4);
 
-    [mx mx_ind] = max(AmdB(l1000+1:l2000));
-    mx_ind += l1000;
-    mask_sample_freqs_kHz = (l1000+1:l2000)*Wo*4/pi;
-    AmdB_(l1000+1:l2000) = parabolic_resonator(mx_ind*Wo*4/pi, mask_sample_freqs_kHz) + mx;  
+    [mx mx_ind] = max(AmdB(lmin+1:L));
+    mx_ind += lmin;
+    AmdB_ = max(AmdB_, parabolic_resonator(mx_ind*Wo*4/pi, mask_sample_freqs_kHz) + mx);  
+    fr2 = mx_ind*Wo*4/pi;
 
-    % fit a resonator to max of third 1000 Hz
+    % fit a third resonator, must be +/- 300 Hz after 2nd resonator
 
-    [mx mx_ind] = max(AmdB(l2000+1:L));
-    mx_ind += l2000;
-    mask_sample_freqs_kHz = (l2000+1:L)*Wo*4/pi;
-    AmdB_(l2000+1:L) = parabolic_resonator(mx_ind*Wo*4/pi, mask_sample_freqs_kHz) + mx;  
+    residual = AmdB;
+    residual(1:lmin) = 0;
 
+    fr2 = mx_ind*Wo*4/pi;
+    fmin = fr2 - 0.300;
+    fmax = fr2 + 0.300;
+    
+    lmin = max(1, round(L*fmin/4));
+    lmax = min(L, round(L*fmax/4));
+    residual(lmin:lmax) = 0;
+
+    [mx mx_ind] = max(residual);
+    AmdB_ = max(AmdB_, parabolic_resonator(mx_ind*Wo*4/pi, mask_sample_freqs_kHz) + mx);  
+    fr3 = mx_ind*Wo*4/pi;
+   
+    %figure(3);
+    %plot(mask_sample_freqs_kHz, residual);
+
+    printf("\nfr1: %f fr2: %f fr3: %f\n", fr1, fr2, fr3);
 endfunction
diff --git a/codec2-dev/octave/newamp1_batch.m b/codec2-dev/octave/newamp1_batch.m
new file mode 100644 (file)
index 0000000..e58879a
--- /dev/null
@@ -0,0 +1,77 @@
+% newamp1_batch.m
+%
+% Copyright David Rowe 2016
+% This program is distributed under the terms of the GNU General Public License 
+% Version 2
+%
+% Octave script to batch process model parameters using the new
+% amplitude model.  Used for generating samples we can listen to.
+%
+% Usage:
+%   ~/codec2-dev/build_linux/src$ ./c2sim ../../raw/hts1a.raw --dump hts1a
+%   $ cd ~/codec2-dev/octave
+%   octave:14> newamp1_batch("../build_linux/src/hts1a")
+%   ~/codec2-dev/build_linux/src$ ./c2sim ../../raw/hts1a.raw --amread hts1a_am.out -o - | play -t raw -r 8000 -s -2 -
+% Or with a little more processing:
+%   codec2-dev/build_linux/src$ ./c2sim ../../raw/hts2a.raw --amread hts2a_am.out --awread hts2a_aw.out --phase0 --postfilter --Woread hts2a_Wo.out -o - | play -q -t raw -r 8000 -s -2 -
+
+
+% process a whole file and write results
+
+function  newamp1_batch(samname, optional_Am_out_name, optional_Aw_out_name)
+  newamp;
+  more off;
+
+  max_amp = 80;
+
+  model_name = strcat(samname,"_model.txt");
+  model = load(model_name);
+  [frames nc] = size(model);
+
+  if nargin == 2
+    Am_out_name = optional_Am_out_name;
+  else
+    Am_out_name = sprintf("%s_am.out", samname);
+  end
+
+  fam  = fopen(Am_out_name,"wb"); 
+
+  % encoder loop ------------------------------------------------------
+
+  for f=1:frames
+    printf("%d ", f);   
+
+    Wo = model(f,1);
+    L = min([model(f,2) max_amp-1]);
+    Am = model(f,3:(L+2));
+    AmdB = 20*log10(Am);
+
+    % find mask
+
+    #{
+    maskdB = mask_model(AmdB, Wo, L);
+    AmdB_ = maskdB;
+    [mx mx_ind] = max(AmdB_);
+    AmdB_(mx_ind) += 6;
+    #}
+
+    AmdB_ = piecewise_model(AmdB, Wo);
+    #{
+    l1000 = floor(L/4);     
+    AmdB_ = AmdB;
+    [mx mx_ind] = max(AmdB_(1:l1000));
+    mask_sample_freqs_kHz = (1:l1000)*Wo*4/pi;
+    AmdB_(1:l1000) = parabolic_resonator(mx_ind*Wo*4/pi, mask_sample_freqs_kHz) + mx;  
+    #}
+
+    Am_ = zeros(1,max_amp);
+    Am_(2:L) = 10 .^ (AmdB_(1:L-1)/20);  % C array doesnt use A[0]
+    fwrite(fam, Am_, "float32");
+    Am_ = zeros(1,max_amp);
+  end
+
+  fclose(fam);
+  printf("\n")
+
+endfunction
+  
diff --git a/codec2-dev/octave/newamp1_fbf.m b/codec2-dev/octave/newamp1_fbf.m
new file mode 100644 (file)
index 0000000..e170721
--- /dev/null
@@ -0,0 +1,98 @@
+% newamp1_fbf.m
+%
+% Copyright David Rowe 2016
+% This program is distributed under the terms of the GNU General Public License 
+% Version 2
+%
+% Interactive Octave script to explore frame by frame operation of new amplitude
+% modelling model.
+%
+% Usage:
+%   Make sure codec2-dev is compiled with the -DDUMP option - see README for
+%    instructions.
+%   ~/codec2-dev/build_linux/src$ ./c2sim ../../raw/hts1a.raw --dump hts1a
+%   $ cd ~/codec2-dev/octave
+%   octave:14> newamp1_fbf("../build_linux/src/hts1a",50)
+
+
+
+function newamp1_fbf(samname, f=10)
+  newamp;
+  more off;
+  plot_spectrum = 1;
+  mask_en = 1;
+
+  % load up text files dumped from c2sim ---------------------------------------
+
+  sn_name = strcat(samname,"_sn.txt");
+  Sn = load(sn_name);
+  sw_name = strcat(samname,"_sw.txt");
+  Sw = load(sw_name);
+  model_name = strcat(samname,"_model.txt");
+  model = load(model_name);
+  [frames tmp] = size(model);
+
+  % Keyboard loop --------------------------------------------------------------
+
+  k = ' ';
+  do 
+    figure(1);
+    clf;
+    s = [ Sn(2*f-1,:) Sn(2*f,:) ];
+    size(s);
+    plot(s);
+    axis([1 length(s) -20000 20000]);
+    title('Time Domain Speech');
+
+    Wo = model(f,1);
+    L = model(f,2);
+    Am = model(f,3:(L+2));
+    AmdB = 20*log10(Am);
+    Am_freqs_kHz = (1:L)*Wo*4/pi;
+
+    #{   
+    [maskdB Am_freqs_kHz] = mask_model(AmdB, Wo, L);
+    AmdB_ = maskdB;
+    [mx mx_ind] = max(AmdB_);
+    AmdB_(mx_ind) += 6;
+    #}
+
+    [AmdB_ residual] = piecewise_model(AmdB, Wo);
+
+    figure(2);
+    clf;
+    title('Frequency Domain');
+
+    axis([1 4000 -20 80]);
+    hold on;
+    plot((1:L)*Wo*4000/pi, AmdB,";Am;r+-");
+    plot(Am_freqs_kHz*1000, AmdB_, ';model;c');
+    plot(Am_freqs_kHz*1000, residual, ';model;g');
+
+    hold off;
+
+    % interactive menu ------------------------------------------
+
+    printf("\rframe: %d  menu: n-next  b-back  q-quit  m-mask_en", f);
+    fflush(stdout);
+    k = kbhit();
+
+    if (k == 'm')
+      if mask_en
+        mask_en = 0;
+      else
+        mask_en = 1; 
+      end
+    endif
+    if (k == 'n')
+      f = f + 1;
+    endif
+    if (k == 'b')
+      f = f - 1;
+    endif
+  until (k == 'q')
+  printf("\n");
+
+endfunction
+