assert( (length(xp) >=3) && (length(yp) >= 3) );
y = zeros(1,length(x));
+ k = 1;
for i=1:length(x)
xi = x(i);
% k is index into xp of where we start 3 points used to form parabola
- k = 1;
while ((xp(k+1) < xi) && (k < (length(xp)-2)))
k++;
end
x1 = xp(k); y1 = yp(k); x2 = xp(k+1); y2 = yp(k+1); x3 = xp(k+2); y3 = yp(k+2);
+ %printf("k: %d i: %d xi: %f x1: %f y1: %f\n", k, i, xi, x1, y1);
a = ((y3-y2)/(x3-x2)-(y2-y1)/(x2-x1))/(x3-x1);
b = ((y3-y2)/(x3-x2)*(x2-x1)+(y2-y1)/(x2-x1)*(x3-x2))/(x3-x1);
rate_K_sample_freqs_kHz = rate_K_sample_freqs_Hz/1000;
endfunction
+
function [rate_K_surface rate_K_sample_freqs_kHz] = resample_const_rate_f_mel(model, K)
rate_K_sample_freqs_kHz = mel_sample_freqs_kHz(K);
rate_K_surface = resample_const_rate_f(model, rate_K_sample_freqs_kHz);
AmdB_peak = max(AmdB);
AmdB(find(AmdB < (AmdB_peak-50))) = AmdB_peak-50;
-
+
rate_L_sample_freqs_kHz = (1:L)*Wo*4/pi;
%rate_K_surface(f,:) = interp1(rate_L_sample_freqs_kHz, AmdB, rate_K_sample_freqs_kHz, "spline", "extrap");
endfunction
+% construct energy quantiser table, and save to text file to include in C
+
+function energy_q = create_energy_q
+ energy_q = 10 + 40/16*(0:15);
+endfunction
+
+function save_energy_q(fn)
+ energy_q = create_energy_q;
+ f = fopen(fn, "wt");
+ fprintf(f, "1 %d\n", length(energy_q));
+ for n=1:length(energy_q)
+ fprintf(f, "%f\n", energy_q(n));
+ end
+ fclose(f);
+endfunction
+
+
+% save's VQ in format that can be compiled by Codec 2 build system
+
+function save_vq(vqset, filenameprefix)
+ [Nvec order stages] = size(vqset);
+ for s=1:stages
+ fn = sprintf("%s_%d.txt", filenameprefix, s);
+ f = fopen(fn, "wt");
+ fprintf(f, "%d %d\n", order, Nvec);
+ for n=1:Nvec
+ for k=1:order
+ fprintf(f, "% 8.4f ", vqset(n,k,s));
+ end
+ fprintf(f, "\n");
+ end
+ fclose(f);
+ end
+endfunction
+
+
% --------------------------------------------------------------------------------
-% Experimental functions used for masking, piecewise models
+% Experimental functions used for masking, piecewise models, not part of newamp1
% --------------------------------------------------------------------------------
K = 20;
[frames tmp] = size(rate_K_surface_c);
[rate_K_surface sample_freqs_kHz] = resample_const_rate_f_mel(model(1:frames,:), K);
-rate_K_surface
+
+ melvq;
+ load train_120_vq; m=5;
+
+ for f=1:frames
+ mean_f(f) = mean(rate_K_surface(f,:));
+ rate_K_surface_no_mean(f,:) = rate_K_surface(f,:) - mean_f(f);
+ end
+
+ [res rate_K_surface_no_mean_ ind] = mbest(train_120_vq, rate_K_surface_no_mean, m);
+
+ for f=1:frames
+ rate_K_surface_no_mean_(f,:) = post_filter(rate_K_surface_no_mean_(f,:), sample_freqs_kHz, 1.5);
+ end
+
+ rate_K_surface_ = zeros(frames, K);
+ energy_q = create_energy_q;
+ for f=1:frames
+ [mean_f_ indx] = quantise(energy_q, mean_f(f));
+ indexes(f,3) = indx - 1;
+ rate_K_surface_(f,:) = rate_K_surface_no_mean_(f,:) + mean_f_;
+ end
+
figure(1);
- mesh(rate_K_surface);
+ mesh(rate_K_surface_);
figure(2);
- mesh(rate_K_surface_c);
+ mesh(rate_K_surface__c);
figure(3);
- mesh(rate_K_surface - rate_K_surface_c);
-
+ mesh(rate_K_surface_ - rate_K_surface__c);
+ axis([1 K 1 frames -1 1])
#{
for f=1:frames
${D}/lspmelvq3.txt
)
+set(CODEBOOKSGE ${D}/gecb.txt)
+
set(CODEBOOKSNEWAMP1
${D}/train_120_1.txt
${D}/train_120_2.txt
)
-set(CODEBOOKSGE ${D}/gecb.txt)
+set(CODEBOOKSNEWAMP1_ENERGY
+ ${D}/newamp1_energy_q.txt
+)
# when crosscompiling we need a native executable
if(CMAKE_CROSSCOMPILING)
DEPENDS generate_codebook ${CODEBOOKSNEWAMP1}
)
+# codebooknewamp1_energy.c
+add_custom_command(
+ OUTPUT ${CMAKE_CURRENT_BINARY_DIR}/codebooknewamp1_energy.c
+ COMMAND generate_codebook newamp1_energy_cb ${CODEBOOKSNEWAMP1_ENERGY} > ${CMAKE_CURRENT_BINARY_DIR}/codebooknewamp1_energy.c
+ DEPENDS generate_codebook ${CODEBOOKSNEWAMP1_ENERGY}
+)
+
#
# codec2 library sources
#
codebookjvm.c
codebookmel.c
codebooklspmelvq.c
- codebooknewamp1.c
codebookge.c
+ codebooknewamp1.c
+ codebooknewamp1_energy.c
golay23.c
freedv_api.c
freedv_vhf_framing.c
--- /dev/null
+1 16
+10.000000
+12.500000
+15.000000
+17.500000
+20.000000
+22.500000
+25.000000
+27.500000
+30.000000
+32.500000
+35.000000
+37.500000
+40.000000
+42.500000
+45.000000
+47.500000
extern const struct lsp_codebook ge_cb[];
extern const struct lsp_codebook lspmelvq_cb[];
extern const struct lsp_codebook newamp1vq_cb[];
+extern const struct lsp_codebook newamp1_energy_cb[];
#endif
x1 = xp[k]; y1 = yp[k]; x2 = xp[k+1]; y2 = yp[k+1]; x3 = xp[k+2]; y3 = yp[k+2];
- printf("k: %d np: %d i: %d xi: %f x1: %f y1: %f\n", k, np, i, xi, x1, y1);
+ //printf("k: %d np: %d i: %d xi: %f x1: %f y1: %f\n", k, np, i, xi, x1, y1);
a = ((y3-y2)/(x3-x2)-(y2-y1)/(x2-x1))/(x3-x1);
b = ((y3-y2)/(x3-x2)*(x2-x1)+(y2-y1)/(x2-x1)*(x3-x2))/(x3-x1);
AmdB_peak = AmdB[m];
}
rate_L_sample_freqs_kHz[m] = m*model->Wo*4.0/M_PI;
- printf("m: %d AmdB: %f AmdB_peak: %f sf: %f\n", m, AmdB[m], AmdB_peak, rate_L_sample_freqs_kHz[m]);
+ //printf("m: %d AmdB: %f AmdB_peak: %f sf: %f\n", m, AmdB[m], AmdB_peak, rate_L_sample_freqs_kHz[m]);
}
/* clip between peak and peak -50dB, to reduce dynamic range */
}
+/*---------------------------------------------------------------------------*\
+
+ FUNCTION....: rate_K_mbest_encode
+ AUTHOR......: David Rowe
+ DATE CREATED: Jan 2017
+
+ Two stage rate K newamp1 VQ quantiser using mbest search.
+
+\*---------------------------------------------------------------------------*/
+
+float rate_K_mbest_encode(int *indexes, float *x, float *xq, int ndim, int mbest_entries)
+{
+ int i, j, n1, n2;
+ const float *codebook1 = newamp1vq_cb[0].cb;
+ const float *codebook2 = newamp1vq_cb[1].cb;
+ struct MBEST *mbest_stage1, *mbest_stage2;
+ float target[ndim];
+ float w[ndim];
+ int index[MBEST_STAGES];
+ float mse, tmp;
+
+ for(i=0; i<ndim; i++)
+ w[i] = 1.0;
+
+ mbest_stage1 = mbest_create(mbest_entries);
+ mbest_stage2 = mbest_create(mbest_entries);
+ for(i=0; i<MBEST_STAGES; i++)
+ index[i] = 0;
+
+ /* Stage 1 */
+
+ mbest_search(codebook1, x, w, ndim, newamp1vq_cb[0].m, mbest_stage1, index);
+ MBEST_PRINT("Stage 1:", mbest_stage1);
+
+ /* Stage 2 */
+
+ for (j=0; j<mbest_entries; j++) {
+ index[1] = n1 = mbest_stage1->list[j].index[0];
+ for(i=0; i<ndim; i++)
+ target[i] = x[i] - codebook1[ndim*n1+i];
+ mbest_search(codebook2, target, w, ndim, newamp1vq_cb[1].m, mbest_stage2, index);
+ }
+ MBEST_PRINT("Stage 2:", mbest_stage2);
+
+ n1 = mbest_stage2->list[0].index[1];
+ n2 = mbest_stage2->list[0].index[0];
+ mse = 0.0;
+ for (i=0;i<ndim;i++) {
+ tmp = codebook1[ndim*n1+i] + codebook2[ndim*n2+i];
+ mse += (x[i]-tmp)*(x[i]-tmp);
+ xq[i] = tmp;
+ }
+
+ mbest_destroy(mbest_stage1);
+ mbest_destroy(mbest_stage2);
+
+ indexes[0] = n1; indexes[1] = n2;;
+
+ return mse;
+}
+
+
+/*---------------------------------------------------------------------------*\
+ FUNCTION....: post_filter
+ AUTHOR......: David Rowe
+ DATE CREATED: Jan 2017
+
+ Post Filter, has a big impact on speech quality after VQ. When used
+ on a mean removed rate K vector, it raises formants, and supresses
+ anti-formants. As it manipulates amplitudes, we normalise energy to
+ prevent clipping or large level variations. pf_gain of 1.2 to 1.5
+ (dB) seems to work OK. Good area for further investigations and
+ improvements in speech quality.
+
+\*---------------------------------------------------------------------------*/
+
+void post_filter_newamp1(float vec[], float sample_freq_kHz[], int K, float pf_gain)
+{
+ int k;
+
+ /*
+ vec is rate K vector describing spectrum of current frame lets
+ pre-emp before applying PF. 20dB/dec over 300Hz. Postfilter
+ affects energy of frame so we measure energy before and after
+ and normalise. Plenty of room for experiment here as well.
+ */
+
+ float pre[K];
+ float e_before = 0.0;
+ float e_after = 0.0;
+ for(k=0; k<K; k++) {
+ pre[k] = 20.0*log10f(sample_freq_kHz[k]/0.3);
+ vec[k] += pre[k];
+ e_before += powf(10.0, 2.0*vec[k]/20.0);
+ vec[k] *= pf_gain;
+ e_after += powf(10.0, 2.0*vec[k]/20.0);
+ }
+
+ float gain = e_after/e_before;
+ float gaindB = 10*log10f(gain);
+
+ for(k=0; k<K; k++) {
+ vec[k] -= gaindB;
+ vec[k] -= pre[k];
+ }
+}
float ftomel(float fHz);
void mel_sample_freqs_kHz(float rate_K_sample_freqs_kHz[], int K);
void resample_const_rate_f(MODEL *model, float rate_K_vec[], float rate_K_sample_freqs_kHz[], int K);
+float rate_K_mbest_encode(int *indexes, float *x, float *xq, int ndim, int mbest_entries);
+void post_filter_newamp1(float vec[], float sample_freq_kHz[], int K, float pf_gain);
#endif
MODEL model;
void *nlp_states;
float pitch, prev_uq_Wo;
- int i,m,f;
+ int i,m,f,k;
if (argc != 2) {
printf("usage: ./tnewamp1 RawFile\n");
int K = 20;
float rate_K_sample_freqs_kHz[K];
- float model_octave[FRAMES][MAX_AMP+2]; // model params in matrix format, useful for C <-> Octave
- float rate_K_surface[FRAMES][K]; // rate K vecs for each frame, form a surface that makes pretty graphs
+ float model_octave[FRAMES][MAX_AMP+2]; // model params in matrix format, useful for C <-> Octave
+ float rate_K_surface[FRAMES][K]; // rate K vecs for each frame, form a surface that makes pretty graphs
+ float rate_K_surface_no_mean[FRAMES][K]; // mean removed surface
+ float rate_K_surface_no_mean_[FRAMES][K]; // quantised mean removed surface
+ float mean[FRAMES];
+ float mean_[FRAMES];
+ float rate_K_surface_[FRAMES][K]; // quantised rate K vecs for each frame
for(f=0; f<FRAMES; f++)
for(m=0; m<MAX_AMP+2; m++)
model_octave[f][m] = 0.0;
mel_sample_freqs_kHz(rate_K_sample_freqs_kHz, K);
+
//for(int k=0; k<K; k++)
// printf("k: %d sf: %f\n", k, rate_K_sample_freqs_kHz[k]);
two_stage_pitch_refinement(&model, Sw);
estimate_amplitudes(&model, Sw, W, 1);
- /* Resample at rate K ----------------------------------------*/
+ /* newamp1 processing ----------------------------------------*/
resample_const_rate_f(&model, &rate_K_surface[f][0], rate_K_sample_freqs_kHz, K);
- for(int k=0; k<K; k++)
- printf("k: %d sf: %f sv: %f\n", k, rate_K_sample_freqs_kHz[k], rate_K_surface[f][k]);
- printf("\n");
+ float sum = 0.0;
+ for(k=0; k<K; k++)
+ sum += rate_K_surface[f][k];
+ mean[f] = sum/K;
+ for(k=0; k<K; k++)
+ rate_K_surface_no_mean[f][k] = rate_K_surface[f][k] - mean[f];
+
+ int vq_indexes[2];
+ rate_K_mbest_encode(vq_indexes, &rate_K_surface_no_mean[f][0], &rate_K_surface_no_mean_[f][0], K, 5);
+
+ post_filter_newamp1(&rate_K_surface_no_mean_[f][0], rate_K_sample_freqs_kHz, K, 1.5);
+
+ int energy_index;
+ float w[1] = {1.0};
+ float se;
+ energy_index = quantise(newamp1_energy_cb[0].cb, &mean[f], w, newamp1_energy_cb[0].k, newamp1_energy_cb[0].m, &se);
+ mean_[f] = newamp1_energy_cb[0].cb[energy_index];
+
+ for(k=0; k<K; k++)
+ rate_K_surface_[f][k] = rate_K_surface_no_mean_[f][k] + mean_[f];
+
+ //for(int k=0; k<K; k++)
+ // printf("k: %d sf: %f sv: %f\n", k, rate_K_sample_freqs_kHz[k], rate_K_surface[f][k]);
+ //printf("\n");
/* log vectors */
assert(fout != NULL);
fprintf(fout, "# Created by tnewamp1.c\n");
octave_save_float(fout, "rate_K_surface_c", (float*)rate_K_surface, FRAMES, K, K);
+ octave_save_float(fout, "mean_c", (float*)mean, FRAMES, 1, 1);
+ octave_save_float(fout, "rate_K_surface_no_mean_c", (float*)rate_K_surface_no_mean, FRAMES, K, K);
+ octave_save_float(fout, "rate_K_surface_no_mean__c", (float*)rate_K_surface_no_mean_, FRAMES, K, K);
+ octave_save_float(fout, "mean__c", (float*)mean_, FRAMES, 1, 1);
+ octave_save_float(fout, "rate_K_surface__c", (float*)rate_K_surface_, FRAMES, K, K);
octave_save_float(fout, "model_c", (float*)model_octave, FRAMES, MAX_AMP+2, MAX_AMP+2);
fclose(fout);
projects / freetel-svn-tracking.git / commitdiff