From 51b4168f8083a9c7aaf61dfc1f82c6cfc3efb4e5 Mon Sep 17 00:00:00 2001 From: drowe67 Date: Mon, 5 Dec 2016 04:52:04 +0000 Subject: [PATCH] simple unquantised three resonator model git-svn-id: https://svn.code.sf.net/p/freetel/code@2925 01035d8c-6547-0410-b346-abe4f91aad63 --- codec2-dev/octave/newamp.m | 53 +++++++++++------ codec2-dev/octave/newamp1_batch.m | 77 ++++++++++++++++++++++++ codec2-dev/octave/newamp1_fbf.m | 98 +++++++++++++++++++++++++++++++ 3 files changed, 211 insertions(+), 17 deletions(-) create mode 100644 codec2-dev/octave/newamp1_batch.m create mode 100644 codec2-dev/octave/newamp1_fbf.m diff --git a/codec2-dev/octave/newamp.m b/codec2-dev/octave/newamp.m index 1a4cd9a4..54b66c23 100644 --- a/codec2-dev/octave/newamp.m +++ b/codec2-dev/octave/newamp.m @@ -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 index 00000000..e58879ae --- /dev/null +++ b/codec2-dev/octave/newamp1_batch.m @@ -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 index 00000000..e170721d --- /dev/null +++ b/codec2-dev/octave/newamp1_fbf.m @@ -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 + + -- 2.25.1