endfunction
+% choose largest sample in band, idea is we care more about finding
+% peaks, can handle some error in frequency. x are non linear
+% (arbitrary) sampling points in kHz
+
+function y = interp_largest(f0_Hz, AmdB, x_kHz)
+ L = length(AmdB);
+ x = x_kHz*1000;
+ y = zeros(1,length(x));
+ bw = x(2) - x(1);
+ k = 1;
+
+ for i=1:length(x)
+
+ % determine limits of this band
+
+ if i>1
+ bw = x(i) - x(i-1);
+ end
+ band_low = x(i) - bw/2; band_high = x(i) + bw/2;
+
+ % map band limits to harmonics
+
+ if x(i) < f0_Hz
+ m_low = m_high = 1;
+ else
+ m_low = round(band_low/f0_Hz); m_high = round(band_high/f0_Hz)-1;
+ m_low = max(1, m_low); m_high = min(L, m_high); m_high = max(m_low, m_high);
+ end
+
+ printf("L: %d f0: %f i: %d band_low: %f band_high: %f m_low: %d m_high: %d\n",L, f0_Hz, i, band_low, band_high, m_low, m_high);
+ % find max in band
+
+ y(i) = max(AmdB(m_low:m_high));
+ end
+
+endfunction
+
% simple linear interpolator
function y = interp_linear(xp, yp, x)
function rate_K_sample_freqs_kHz = mel_sample_freqs_kHz(K)
- mel_start = ftomel(200); mel_end = ftomel(3700);
+ mel_start = ftomel(100); mel_end = ftomel(0.95*4000);
step = (mel_end-mel_start)/(K-1);
mel = mel_start:step:mel_end;
rate_K_sample_freqs_Hz = 700*((10 .^ (mel/2595)) - 1);
endfunction
+function [rate_K_vec_corrected orig_error error nasty_error_log nasty_error_m_log] = correct_rate_K_vec(rate_K_vec, rate_K_sample_freqs_kHz, AmdB, AmdB_, K, Wo, L, Fs)
+
+ % aliasing correction --------------------------------------
+
+ % The mel sample rate decreases as frequency increases. Look for
+ % any regions above 1000Hz where we have missed definition of a
+ % spectral peak (formant) due to aliasing. Adjust the rate K
+ % sample levels to restore peaks. Theory is that correct
+ % definition of a formant is less important than the frequency of
+ % the formant. As long as we define a formant in that general
+ % frequency area it will sound OK.
+
+ Am_freqs_kHz = (1:L)*Wo*4/pi;
+
+ % Lets see where we have made an error
+
+ error = orig_error = AmdB - AmdB_;
+
+ Ncorrections = 3; % maximum number of rate K samples to correct
+ error_thresh = 3; % only worry about errors larger than thresh
+
+ start_m = floor(L*1000/(Fs/2));
+ error(1:start_m) = 0; % first 1000Hz is densly sampled so ignore
+ nasty_error_m_log = []; nasty_error_log = [];
+
+
+ rate_K_vec_corrected = rate_K_vec;
+ for i=1:Ncorrections
+ [mx mx_m] = max(error);
+
+ if mx > error_thresh
+ nasty_error_log = [nasty_error_log mx];
+ nasty_error_m_log = [nasty_error_m_log mx_m];
+
+ % find closest rate K sample to nasty error
+
+ nasty_error_freq = mx_m*Wo*Fs/(2*pi*1000)
+ [tmp closest_k] = min(abs(rate_K_sample_freqs_kHz - nasty_error_freq));
+ rate_K_vec_corrected(closest_k) = AmdB(mx_m);
+
+ % zero out error in this region and look for another large error region
+
+ k = max(1, closest_k-1);
+ rate_K_prev_sample_kHz = rate_K_sample_freqs_kHz(k)
+ k = min(K, closest_k+1);
+ rate_K_next_sample_kHz = rate_K_sample_freqs_kHz(k)
+
+ [tmp st_m] = min(abs(Am_freqs_kHz - rate_K_prev_sample_kHz));
+ [tmp en_m] = min(abs(Am_freqs_kHz - rate_K_next_sample_kHz));
+ if closest_k == K
+ en_m = L;
+ end
+ error(st_m:en_m) = 0;
+ end
+ end
+endfunction
+
+
% Take a rate K surface and convert back to time varying rate L
function [model_ AmdB_] = resample_rate_L(model, rate_K_surface, rate_K_sample_freqs_kHz)
endfunction
+function rate_K_surface_no_slope = remove_slope(rate_K_surface)
+ [frames K] = size(rate_K_surface);
+ rate_K_surface_no_slope = zeros(frames,K);
+ for f=1:frames
+ [gradient intercept] = linreg(1:K, rate_K_surface(f,:), K);
+ printf("f: %d gradient: %f intercept: %f\n", f, gradient, intercept);
+ rate_K_surface_no_slope(f,:) = rate_K_surface(f,:) - (intercept + gradient*(1:K));
+ end
+endfunction
#{
Octave script to batch process model parameters using the new
- amplitude model, version 1. Outputs anotehr set of model parameters
+ amplitude model, version 1. Outputs another set of model parameters
that can be fed to c2sim for listening tests. The companion
newamp1_fbf.m script is used to visualise the processing frame by frame
newamp;
more off;
- max_amp = 80;
+ max_amp = 160;
postfilter = 0; % optional postfiler that runs on Am, not used atm
synth_phase = 1;
%model_ = experiment_filter(model);
%model_ = experiment_filter_dec_filter(model);
- %[model_ surface] = experiment_mel_freq(model, 1, 1, voicing);
+ [model_ surface] = experiment_mel_freq(model, 0, 1, voicing);
+
%model_ = experiment_dec_abys(model, 8, 1, 1, 1, voicing);
+#{
[model_ voicing_ indexes] = experiment_rate_K_dec(model, voicing); % encoder/decoder, lets toss away results except for indexes
- %[model_ voicing_] = model_from_indexes(indexes); % decoder uses just indexes, outputs vecs for synthesis
[model_ voicing_] = model_from_indexes_fbf(indexes); % decoder uses just indexes, outputs vecs for synthesis
+#}
+
+ %[model_ voicing_] = model_from_indexes(indexes); % decoder uses just indexes, outputs vecs for synthesis
%model_ = experiment_dec_linear(model_);
%model_ = experiment_energy_rate_linear(model, 1, 0);
% synthesis phase spectra from magnitiude spectra using minimum phase techniques
- fft_enc = 128;
+ fft_enc = 256;
phase = determine_phase(model_, f, fft_enc);
assert(length(phase) == fft_enc);
%Aw = zeros(1, fft_enc*2);
function [model_ rate_K_surface] = experiment_mel_freq(model, vq_en=0, plots=1, voicing)
[frames nc] = size(model);
- K = 20;
- [rate_K_surface rate_K_sample_freqs_kHz] = resample_const_rate_f_mel(model, K);
-
+ K = 20; Fs = 8000;
+
+ for f=1:frames
+ 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;
+ [rate_K_vec rate_K_sample_freqs_kHz] = resample_const_rate_f_mel(model(f,:), K);
+ [tmp_ AmdB_] = resample_rate_L(model(f,:), rate_K_vec, rate_K_sample_freqs_kHz);
+ [rate_K_vec_corrected orig_error error nasty_error_log nasty_error_m_log] = correct_rate_K_vec(rate_K_vec, rate_K_sample_freqs_kHz, AmdB, AmdB_, K, Wo, L, Fs);
+ %rate_K_surface(f,:) = rate_K_vec;
+ rate_K_surface(f,:) = rate_K_vec_corrected;
+ end
+
if plots
figure(1); clf; mesh(rate_K_surface);
figure(2); clf; plot_dft_surface(rate_K_surface)
Fs = 100; dF = Fs/frames;
mesh(1:K, (1:frames/2)*dF, 20*log10(abs(dft_surface(1:frames/2,:,:))));
endfunction
+
function newamp1_fbf(samname, f=10)
newamp;
more off;
- quant_en = 0; pf_en = 0; plot_phase = 1;
+ quant_en = 0; pf_en = 0; plot_phase = 0;
melvq;
-
- K=20; load train_120_vq; m=5;
+
+ Fs = 8000; K=20; load train_120_vq; m=5;
% load up text files dumped from c2sim ---------------------------------------
AmdB = 20*log10(Am);
Am_freqs_kHz = (1:L)*Wo*4/pi;
- % plots for mel sampling
+ % resample at rate K
[rate_K_vec rate_K_sample_freqs_kHz] = resample_const_rate_f_mel(model(f,:), K);
+ [tmp_ AmdB_] = resample_rate_L(model(f,:), rate_K_vec, rate_K_sample_freqs_kHz);
+ [rate_K_vec_corrected orig_error error nasty_error_log nasty_error_m_log] = correct_rate_K_vec(rate_K_vec, rate_K_sample_freqs_kHz, AmdB, AmdB_, K, Wo, L, Fs);
+ [tmp_ AmdB_corrected] = resample_rate_L(model(f,:), rate_K_vec, rate_K_sample_freqs_kHz);
+
+ % quantisation simulation ---------------------------------
mean_f = mean(rate_K_vec);
rate_K_vec_no_mean = rate_K_vec - mean_f;
rate_K_vec_ = rate_K_vec_no_mean_ + mean_f;
[model_ AmdB_] = resample_rate_L(model(f,:), rate_K_vec_, rate_K_sample_freqs_kHz);
+ %[model_ AmdB_largest] = resample_rate_L(model(f,:), rate_K_largest, rate_K_sample_freqs_kHz);
% plots ----------------------------------
axis([1 4000 -20 80]);
hold on;
+
plot((1:L)*Wo*4000/pi, AmdB,";Am;b+-");
- plot(rate_K_sample_freqs_kHz*1000, rate_K_vec, ';rate K mel;g+-');
- if quant_en >= 1
- plot((1:L)*Wo*4000/pi, AmdB_,";Am quant;k+-");
- end
- if quant_en == 2
- plot(rate_K_sample_freqs_kHz*1000, rate_K_vec_, ';rate K mel quant;r+-');
- end
+ stem(rate_K_sample_freqs_kHz*1000, rate_K_vec, 'g');
+ plot(rate_K_sample_freqs_kHz*1000, rate_K_vec, ';rate K;g:');
+ plot(rate_K_sample_freqs_kHz*1000, rate_K_vec_corrected, ';rate K;g-', 'linewidth', 2);
+
+ plot((1:L)*Wo*4000/pi, -10 + orig_error ,";AmdB - AmdB_ error;k+-");
+ plot((1:L)*Wo*4000/pi, -10 + error ,";AmdB - AmdB_ error;r-");
+ plot(nasty_error_m_log*Wo*4000/pi, AmdB(nasty_error_m_log), 'ro', 'markersize', 10);
+
hold off;
+ #{
figure(3);
clf;
title('Frequency Domain 2');
plot(rate_K_sample_freqs_kHz*1000, rate_K_vec_no_mean_, ';rate K mel no mean quant;r+-');
end
hold off;
+ #}
if plot_phase
phase_512 = determine_phase(model_, 1, 512);
if k == 'm'
quant_en++;
- if quant_en > 2
+ if quant_en > 1
quant_en = 0;
end
endif