From b409a4b0c6131e01fc5a03c0fc31caa4829b0dec Mon Sep 17 00:00:00 2001
From: Johnathan Corgan <jcorgan@corganenterprises.com>
Date: Mon, 18 Jul 2011 15:54:43 -0700
Subject: gr-vocoder: re-implemented gr-codec2-vocoder inside gr-vocoder
---
gr-vocoder/lib/codec2/quantise.c | 851 +++++++++++++++++++++++++++++++++++++++
1 file changed, 851 insertions(+)
create mode 100644 gr-vocoder/lib/codec2/quantise.c
(limited to 'gr-vocoder/lib/codec2/quantise.c')
diff --git a/gr-vocoder/lib/codec2/quantise.c b/gr-vocoder/lib/codec2/quantise.c
new file mode 100644
index 0000000000..ff8d156b5b
--- /dev/null
+++ b/gr-vocoder/lib/codec2/quantise.c
@@ -0,0 +1,851 @@
+/*---------------------------------------------------------------------------*\
+
+ FILE........: quantise.c
+ AUTHOR......: David Rowe
+ DATE CREATED: 31/5/92
+
+ Quantisation functions for the sinusoidal coder.
+
+\*---------------------------------------------------------------------------*/
+
+/*
+ All rights reserved.
+
+ This program is free software; you can redistribute it and/or modify
+ it under the terms of the GNU Lesser General Public License version 2.1, as
+ published by the Free Software Foundation. This program is
+ distributed in the hope that it will be useful, but WITHOUT ANY
+ WARRANTY; without even the implied warranty of MERCHANTABILITY or
+ FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public
+ License for more details.
+
+ You should have received a copy of the GNU Lesser General Public License
+ along with this program; if not, see <http://www.gnu.org/licenses/>.
+
+*/
+
+#include <assert.h>
+#include <ctype.h>
+#include <stdio.h>
+#include <stdlib.h>
+#include <string.h>
+#include <math.h>
+
+#include "defines.h"
+#include "dump.h"
+#include "quantise.h"
+#include "lpc.h"
+#include "lsp.h"
+#include "fft.h"
+
+#define LSP_DELTA1 0.01 /* grid spacing for LSP root searches */
+
+/*---------------------------------------------------------------------------*\
+
+ FUNCTION HEADERS
+
+\*---------------------------------------------------------------------------*/
+
+float speech_to_uq_lsps(float lsp[], float ak[], float Sn[], float w[],
+ int order);
+
+/*---------------------------------------------------------------------------*\
+
+ FUNCTIONS
+
+\*---------------------------------------------------------------------------*/
+
+int lsp_bits(int i) {
+ return lsp_cb[i].log2m;
+}
+
+#if VECTOR_QUANTISATION
+/*---------------------------------------------------------------------------*\
+
+ quantise_uniform
+
+ Simulates uniform quantising of a float.
+
+\*---------------------------------------------------------------------------*/
+
+void quantise_uniform(float *val, float min, float max, int bits)
+{
+ int levels = 1 << (bits-1);
+ float norm;
+ int index;
+
+ /* hard limit to quantiser range */
+
+ printf("min: %f max: %f val: %f ", min, max, val[0]);
+ if (val[0] < min) val[0] = min;
+ if (val[0] > max) val[0] = max;
+
+ norm = (*val - min)/(max-min);
+ printf("%f norm: %f ", val[0], norm);
+ index = fabs(levels*norm + 0.5);
+
+ *val = min + index*(max-min)/levels;
+
+ printf("index %d val_: %f\n", index, val[0]);
+}
+
+#endif
+
+/*---------------------------------------------------------------------------*\
+
+ quantise_init
+
+ Loads the entire LSP quantiser comprised of several vector quantisers
+ (codebooks).
+
+\*---------------------------------------------------------------------------*/
+
+void quantise_init()
+{
+}
+
+/*---------------------------------------------------------------------------*\
+
+ quantise
+
+ Quantises vec by choosing the nearest vector in codebook cb, and
+ returns the vector index. The squared error of the quantised vector
+ is added to se.
+
+\*---------------------------------------------------------------------------*/
+
+long quantise(const float * cb, float vec[], float w[], int k, int m, float *se)
+/* float cb[][K]; current VQ codebook */
+/* float vec[]; vector to quantise */
+/* float w[]; weighting vector */
+/* int k; dimension of vectors */
+/* int m; size of codebook */
+/* float *se; accumulated squared error */
+{
+ float e; /* current error */
+ long besti; /* best index so far */
+ float beste; /* best error so far */
+ long j;
+ int i;
+
+ besti = 0;
+ beste = 1E32;
+ for(j=0; j<m; j++) {
+ e = 0.0;
+ for(i=0; i<k; i++)
+ e += pow((cb[j*k+i]-vec[i])*w[i],2.0);
+ if (e < beste) {
+ beste = e;
+ besti = j;
+ }
+ }
+
+ *se += beste;
+
+ return(besti);
+}
+
+/*---------------------------------------------------------------------------*\
+
+ lspd_quantise
+
+ Scalar lsp difference quantiser.
+
+\*---------------------------------------------------------------------------*/
+
+void lspd_quantise(
+ float lsp[],
+ float lsp_[],
+ int order
+)
+{
+ int i,k,m;
+ float lsp_hz[LPC_MAX];
+ float lsp__hz[LPC_MAX];
+ float dlsp[LPC_MAX];
+ float dlsp_[LPC_MAX];
+ float wt[1];
+ const float *cb;
+ float se;
+ int indexes[LPC_MAX];
+
+ /* convert from radians to Hz so we can use human readable
+ frequencies */
+
+ for(i=0; i<order; i++)
+ lsp_hz[i] = (4000.0/PI)*lsp[i];
+
+ dlsp[0] = lsp_hz[0];
+ for(i=1; i<order; i++)
+ dlsp[i] = lsp_hz[i] - lsp_hz[i-1];
+
+ /* simple uniform scalar quantisers */
+
+ wt[0] = 1.0;
+ for(i=0; i<order; i++) {
+ if (i)
+ dlsp[i] = lsp_hz[i] - lsp__hz[i-1];
+ else
+ dlsp[0] = lsp_hz[0];
+
+ k = lsp_cbd[i].k;
+ m = lsp_cbd[i].m;
+ cb = lsp_cbd[i].cb;
+ indexes[i] = quantise(cb, &dlsp[i], wt, k, m, &se);
+ dlsp_[i] = cb[indexes[i]*k];
+
+ if (i)
+ lsp__hz[i] = lsp__hz[i-1] + dlsp_[i];
+ else
+ lsp__hz[0] = dlsp_[0];
+ }
+ for(; i<order; i++)
+ lsp__hz[i] = lsp__hz[i-1] + dlsp[i];
+
+ /* convert back to radians */
+
+ for(i=0; i<order; i++)
+ lsp_[i] = (PI/4000.0)*lsp__hz[i];
+}
+
+/*---------------------------------------------------------------------------*\
+
+ lspd_vq_quantise
+
+ Vector lsp difference quantiser.
+
+\*---------------------------------------------------------------------------*/
+
+void lspdvq_quantise(
+ float lsp[],
+ float lsp_[],
+ int order
+)
+{
+ int i,k,m,ncb, nlsp;
+ float dlsp[LPC_MAX];
+ float dlsp_[LPC_MAX];
+ float wt[LPC_ORD];
+ const float *cb;
+ float se;
+ int index;
+
+ dlsp[0] = lsp[0];
+ for(i=1; i<order; i++)
+ dlsp[i] = lsp[i] - lsp[i-1];
+
+ for(i=0; i<order; i++)
+ dlsp_[i] = dlsp[i];
+
+ for(i=0; i<order; i++)
+ wt[i] = 1.0;
+
+ /* scalar quantise dLSPs 1,2,3,4,5 */
+
+ for(i=0; i<5; i++) {
+ if (i)
+ dlsp[i] = (lsp[i] - lsp_[i-1])*4000.0/PI;
+ else
+ dlsp[0] = lsp[0]*4000.0/PI;
+
+ k = lsp_cbdvq[i].k;
+ m = lsp_cbdvq[i].m;
+ cb = lsp_cbdvq[i].cb;
+ index = quantise(cb, &dlsp[i], wt, k, m, &se);
+ dlsp_[i] = cb[index*k]*PI/4000.0;
+
+ if (i)
+ lsp_[i] = lsp_[i-1] + dlsp_[i];
+ else
+ lsp_[0] = dlsp_[0];
+ }
+ dlsp[i] = lsp[i] - lsp_[i-1];
+ dlsp_[i] = dlsp[i];
+
+ //printf("lsp[0] %f lsp_[0] %f\n", lsp[0], lsp_[0]);
+ //printf("lsp[1] %f lsp_[1] %f\n", lsp[1], lsp_[1]);
+
+#ifdef TT
+ /* VQ dLSPs 3,4,5 */
+
+ ncb = 2;
+ nlsp = 2;
+ k = lsp_cbdvq[ncb].k;
+ m = lsp_cbdvq[ncb].m;
+ cb = lsp_cbdvq[ncb].cb;
+ index = quantise(cb, &dlsp[nlsp], wt, k, m, &se);
+ dlsp_[nlsp] = cb[index*k];
+ dlsp_[nlsp+1] = cb[index*k+1];
+ dlsp_[nlsp+2] = cb[index*k+2];
+
+ lsp_[0] = dlsp_[0];
+ for(i=1; i<5; i++)
+ lsp_[i] = lsp_[i-1] + dlsp_[i];
+ dlsp[i] = lsp[i] - lsp_[i-1];
+ dlsp_[i] = dlsp[i];
+#endif
+ /* VQ dLSPs 6,7,8,9,10 */
+
+ ncb = 5;
+ nlsp = 5;
+ k = lsp_cbdvq[ncb].k;
+ m = lsp_cbdvq[ncb].m;
+ cb = lsp_cbdvq[ncb].cb;
+ index = quantise(cb, &dlsp[nlsp], wt, k, m, &se);
+ dlsp_[nlsp] = cb[index*k];
+ dlsp_[nlsp+1] = cb[index*k+1];
+ dlsp_[nlsp+2] = cb[index*k+2];
+ dlsp_[nlsp+3] = cb[index*k+3];
+ dlsp_[nlsp+4] = cb[index*k+4];
+
+ /* rebuild LSPs for dLSPs */
+
+ lsp_[0] = dlsp_[0];
+ for(i=1; i<order; i++)
+ lsp_[i] = lsp_[i-1] + dlsp_[i];
+}
+
+void check_lsp_order(float lsp[], int lpc_order)
+{
+ int i;
+ float tmp;
+
+ for(i=1; i<lpc_order; i++)
+ if (lsp[i] < lsp[i-1]) {
+ printf("swap %d\n",i);
+ tmp = lsp[i-1];
+ lsp[i-1] = lsp[i]-0.05;
+ lsp[i] = tmp+0.05;
+ }
+}
+
+void force_min_lsp_dist(float lsp[], int lpc_order)
+{
+ int i;
+
+ for(i=1; i<lpc_order; i++)
+ if ((lsp[i]-lsp[i-1]) < 0.01) {
+ lsp[i] += 0.01;
+ }
+}
+
+/*---------------------------------------------------------------------------*\
+
+ lpc_model_amplitudes
+
+ Derive a LPC model for amplitude samples then estimate amplitude samples
+ from this model with optional LSP quantisation.
+
+ Returns the spectral distortion for this frame.
+
+\*---------------------------------------------------------------------------*/
+
+float lpc_model_amplitudes(
+ float Sn[], /* Input frame of speech samples */
+ float w[],
+ MODEL *model, /* sinusoidal model parameters */
+ int order, /* LPC model order */
+ int lsp_quant, /* optional LSP quantisation if non-zero */
+ float ak[] /* output aks */
+)
+{
+ float Wn[M];
+ float R[LPC_MAX+1];
+ float E;
+ int i,j;
+ float snr;
+ float lsp[LPC_MAX];
+ float lsp_hz[LPC_MAX];
+ float lsp_[LPC_MAX];
+ int roots; /* number of LSP roots found */
+ int index;
+ float se;
+ int k,m;
+ const float * cb;
+ float wt[LPC_MAX];
+
+ for(i=0; i<M; i++)
+ Wn[i] = Sn[i]*w[i];
+ autocorrelate(Wn,R,M,order);
+ levinson_durbin(R,ak,order);
+
+ E = 0.0;
+ for(i=0; i<=order; i++)
+ E += ak[i]*R[i];
+
+ for(i=0; i<order; i++)
+ wt[i] = 1.0;
+
+ if (lsp_quant) {
+ roots = lpc_to_lsp(ak, order, lsp, 5, LSP_DELTA1);
+ if (roots != order)
+ printf("LSP roots not found\n");
+
+ /* convert from radians to Hz to make quantisers more
+ human readable */
+
+ for(i=0; i<order; i++)
+ lsp_hz[i] = (4000.0/PI)*lsp[i];
+
+ /* simple uniform scalar quantisers */
+
+ for(i=0; i<10; i++) {
+ k = lsp_cb[i].k;
+ m = lsp_cb[i].m;
+ cb = lsp_cb[i].cb;
+ index = quantise(cb, &lsp_hz[i], wt, k, m, &se);
+ lsp_hz[i] = cb[index*k];
+ }
+
+ /* experiment: simulating uniform quantisation error
+ for(i=0; i<order; i++)
+ lsp[i] += PI*(12.5/4000.0)*(1.0 - 2.0*(float)rand()/RAND_MAX);
+ */
+
+ for(i=0; i<order; i++)
+ lsp[i] = (PI/4000.0)*lsp_hz[i];
+
+ /* Bandwidth Expansion (BW). Prevents any two LSPs getting too
+ close together after quantisation. We know from experiment
+ that LSP quantisation errors < 12.5Hz (25Hz setp size) are
+ inaudible so we use that as the minimum LSP separation.
+ */
+
+ for(i=1; i<5; i++) {
+ if (lsp[i] - lsp[i-1] < PI*(12.5/4000.0))
+ lsp[i] = lsp[i-1] + PI*(12.5/4000.0);
+ }
+
+ /* as quantiser gaps increased, larger BW expansion was required
+ to prevent twinkly noises */
+
+ for(i=5; i<8; i++) {
+ if (lsp[i] - lsp[i-1] < PI*(25.0/4000.0))
+ lsp[i] = lsp[i-1] + PI*(25.0/4000.0);
+ }
+ for(i=8; i<order; i++) {
+ if (lsp[i] - lsp[i-1] < PI*(75.0/4000.0))
+ lsp[i] = lsp[i-1] + PI*(75.0/4000.0);
+ }
+
+ for(j=0; j<order; j++)
+ lsp_[j] = lsp[j];
+
+ lsp_to_lpc(lsp_, ak, order);
+#ifdef DUMP
+ dump_lsp(lsp);
+#endif
+ }
+
+#ifdef DUMP
+ dump_E(E);
+#endif
+ #ifdef SIM_QUANT
+ /* simulated LPC energy quantisation */
+ {
+ float e = 10.0*log10(E);
+ e += 2.0*(1.0 - 2.0*(float)rand()/RAND_MAX);
+ E = pow(10.0,e/10.0);
+ }
+ #endif
+
+ aks_to_M2(ak,order,model,E,&snr, 1); /* {ak} -> {Am} LPC decode */
+
+ return snr;
+}
+
+/*---------------------------------------------------------------------------*\
+
+ aks_to_M2()
+
+ Transforms the linear prediction coefficients to spectral amplitude
+ samples. This function determines A(m) from the average energy per
+ band using an FFT.
+
+\*---------------------------------------------------------------------------*/
+
+void aks_to_M2(
+ float ak[], /* LPC's */
+ int order,
+ MODEL *model, /* sinusoidal model parameters for this frame */
+ float E, /* energy term */
+ float *snr, /* signal to noise ratio for this frame in dB */
+ int dump /* true to dump sample to dump file */
+)
+{
+ COMP Pw[FFT_DEC]; /* power spectrum */
+ int i,m; /* loop variables */
+ int am,bm; /* limits of current band */
+ float r; /* no. rads/bin */
+ float Em; /* energy in band */
+ float Am; /* spectral amplitude sample */
+ float signal, noise;
+
+ r = TWO_PI/(FFT_DEC);
+
+ /* Determine DFT of A(exp(jw)) --------------------------------------------*/
+
+ for(i=0; i<FFT_DEC; i++) {
+ Pw[i].real = 0.0;
+ Pw[i].imag = 0.0;
+ }
+
+ for(i=0; i<=order; i++)
+ Pw[i].real = ak[i];
+ fft(&Pw[0].real,FFT_DEC,1);
+
+ /* Determine power spectrum P(w) = E/(A(exp(jw))^2 ------------------------*/
+
+ for(i=0; i<FFT_DEC/2; i++)
+ Pw[i].real = E/(Pw[i].real*Pw[i].real + Pw[i].imag*Pw[i].imag);
+#ifdef DUMP
+ if (dump)
+ dump_Pw(Pw);
+#endif
+
+ /* Determine magnitudes by linear interpolation of P(w) -------------------*/
+
+ signal = noise = 0.0;
+ for(m=1; m<=model->L; m++) {
+ am = floor((m - 0.5)*model->Wo/r + 0.5);
+ bm = floor((m + 0.5)*model->Wo/r + 0.5);
+ Em = 0.0;
+
+ for(i=am; i<bm; i++)
+ Em += Pw[i].real;
+ Am = sqrt(Em);
+
+ signal += pow(model->A[m],2.0);
+ noise += pow(model->A[m] - Am,2.0);
+ model->A[m] = Am;
+ }
+ *snr = 10.0*log10(signal/noise);
+}
+
+/*---------------------------------------------------------------------------*\
+
+ FUNCTION....: encode_Wo()
+ AUTHOR......: David Rowe
+ DATE CREATED: 22/8/2010
+
+ Encodes Wo using a WO_LEVELS quantiser.
+
+\*---------------------------------------------------------------------------*/
+
+int encode_Wo(float Wo)
+{
+ int index;
+ float Wo_min = TWO_PI/P_MAX;
+ float Wo_max = TWO_PI/P_MIN;
+ float norm;
+
+ norm = (Wo - Wo_min)/(Wo_max - Wo_min);
+ index = floor(WO_LEVELS * norm + 0.5);
+ if (index < 0 ) index = 0;
+ if (index > (WO_LEVELS-1)) index = WO_LEVELS-1;
+
+ return index;
+}
+
+/*---------------------------------------------------------------------------*\
+
+ FUNCTION....: decode_Wo()
+ AUTHOR......: David Rowe
+ DATE CREATED: 22/8/2010
+
+ Decodes Wo using a WO_LEVELS quantiser.
+
+\*---------------------------------------------------------------------------*/
+
+float decode_Wo(int index)
+{
+ float Wo_min = TWO_PI/P_MAX;
+ float Wo_max = TWO_PI/P_MIN;
+ float step;
+ float Wo;
+
+ step = (Wo_max - Wo_min)/WO_LEVELS;
+ Wo = Wo_min + step*(index);
+
+ return Wo;
+}
+
+/*---------------------------------------------------------------------------*\
+
+ FUNCTION....: speech_to_uq_lsps()
+ AUTHOR......: David Rowe
+ DATE CREATED: 22/8/2010
+
+ Analyse a windowed frame of time domain speech to determine LPCs
+ which are the converted to LSPs for quantisation and transmission
+ over the channel.
+
+\*---------------------------------------------------------------------------*/
+
+float speech_to_uq_lsps(float lsp[],
+ float ak[],
+ float Sn[],
+ float w[],
+ int order
+)
+{
+ int i, roots;
+ float Wn[M];
+ float R[LPC_MAX+1];
+ float E;
+
+ for(i=0; i<M; i++)
+ Wn[i] = Sn[i]*w[i];
+ autocorrelate(Wn, R, M, order);
+ levinson_durbin(R, ak, order);
+
+ E = 0.0;
+ for(i=0; i<=order; i++)
+ E += ak[i]*R[i];
+
+ roots = lpc_to_lsp(ak, order, lsp, 5, LSP_DELTA1);
+ if (roots != order) {
+ /* for some reason LSP roots could not be found */
+ /* some alpha testers are reporting this condition */
+ fprintf(stderr, "LSP roots not found!\nroots = %d\n", roots);
+ for(i=0; i<=order; i++)
+ fprintf(stderr, "a[%d] = %f\n", i, ak[i]);
+
+ /* some benign LSP values we can use instead */
+ for(i=0; i<order; i++)
+ lsp[i] = (PI/order)*(float)i;
+ }
+
+ return E;
+}
+
+/*---------------------------------------------------------------------------*\
+
+ FUNCTION....: encode_lsps()
+ AUTHOR......: David Rowe
+ DATE CREATED: 22/8/2010
+
+ From a vector of unquantised (floating point) LSPs finds the quantised
+ LSP indexes.
+
+\*---------------------------------------------------------------------------*/
+
+void encode_lsps(int indexes[], float lsp[], int order)
+{
+ int i,k,m;
+ float wt[1];
+ float lsp_hz[LPC_MAX];
+ const float * cb;
+ float se;
+
+ /* convert from radians to Hz so we can use human readable
+ frequencies */
+
+ for(i=0; i<order; i++)
+ lsp_hz[i] = (4000.0/PI)*lsp[i];
+
+ /* simple uniform scalar quantisers */
+
+ wt[0] = 1.0;
+ for(i=0; i<order; i++) {
+ k = lsp_cb[i].k;
+ m = lsp_cb[i].m;
+ cb = lsp_cb[i].cb;
+ indexes[i] = quantise(cb, &lsp_hz[i], wt, k, m, &se);
+ }
+}
+
+/*---------------------------------------------------------------------------*\
+
+ FUNCTION....: decode_lsps()
+ AUTHOR......: David Rowe
+ DATE CREATED: 22/8/2010
+
+ From a vector of quantised LSP indexes, returns the quantised
+ (floating point) LSPs.
+
+\*---------------------------------------------------------------------------*/
+
+void decode_lsps(float lsp[], int indexes[], int order)
+{
+ int i,k;
+ float lsp_hz[LPC_MAX];
+ const float * cb;
+
+ for(i=0; i<order; i++) {
+ k = lsp_cb[i].k;
+ cb = lsp_cb[i].cb;
+ lsp_hz[i] = cb[indexes[i]*k];
+ }
+
+ /* convert back to radians */
+
+ for(i=0; i<order; i++)
+ lsp[i] = (PI/4000.0)*lsp_hz[i];
+}
+
+/*---------------------------------------------------------------------------*\
+
+ FUNCTION....: bw_expand_lsps()
+ AUTHOR......: David Rowe
+ DATE CREATED: 22/8/2010
+
+ Applies Bandwidth Expansion (BW) to a vector of LSPs. Prevents any
+ two LSPs getting too close together after quantisation. We know
+ from experiment that LSP quantisation errors < 12.5Hz (25Hz setp
+ size) are inaudible so we use that as the minimum LSP separation.
+
+\*---------------------------------------------------------------------------*/
+
+void bw_expand_lsps(float lsp[],
+ int order
+)
+{
+ int i;
+
+ for(i=1; i<5; i++) {
+ if (lsp[i] - lsp[i-1] < PI*(12.5/4000.0))
+ lsp[i] = lsp[i-1] + PI*(12.5/4000.0);
+ }
+
+ /* As quantiser gaps increased, larger BW expansion was required
+ to prevent twinkly noises. This may need more experiment for
+ different quanstisers.
+ */
+
+ for(i=5; i<8; i++) {
+ if (lsp[i] - lsp[i-1] < PI*(25.0/4000.0))
+ lsp[i] = lsp[i-1] + PI*(25.0/4000.0);
+ }
+ for(i=8; i<order; i++) {
+ if (lsp[i] - lsp[i-1] < PI*(75.0/4000.0))
+ lsp[i] = lsp[i-1] + PI*(75.0/4000.0);
+ }
+}
+
+/*---------------------------------------------------------------------------*\
+
+ FUNCTION....: apply_lpc_correction()
+ AUTHOR......: David Rowe
+ DATE CREATED: 22/8/2010
+
+ Apply first harmonic LPC correction at decoder. This helps improve
+ low pitch males after LPC modelling, like hts1a and morig.
+
+\*---------------------------------------------------------------------------*/
+
+void apply_lpc_correction(MODEL *model)
+{
+ if (model->Wo < (PI*150.0/4000)) {
+ model->A[1] *= 0.032;
+ }
+}
+
+/*---------------------------------------------------------------------------*\
+
+ FUNCTION....: encode_energy()
+ AUTHOR......: David Rowe
+ DATE CREATED: 22/8/2010
+
+ Encodes LPC energy using an E_LEVELS quantiser.
+
+\*---------------------------------------------------------------------------*/
+
+int encode_energy(float e)
+{
+ int index;
+ float e_min = E_MIN_DB;
+ float e_max = E_MAX_DB;
+ float norm;
+
+ e = 10.0*log10(e);
+ norm = (e - e_min)/(e_max - e_min);
+ index = floor(E_LEVELS * norm + 0.5);
+ if (index < 0 ) index = 0;
+ if (index > (E_LEVELS-1)) index = E_LEVELS-1;
+
+ return index;
+}
+
+/*---------------------------------------------------------------------------*\
+
+ FUNCTION....: decode_energy()
+ AUTHOR......: David Rowe
+ DATE CREATED: 22/8/2010
+
+ Decodes energy using a WO_BITS quantiser.
+
+\*---------------------------------------------------------------------------*/
+
+float decode_energy(int index)
+{
+ float e_min = E_MIN_DB;
+ float e_max = E_MAX_DB;
+ float step;
+ float e;
+
+ step = (e_max - e_min)/E_LEVELS;
+ e = e_min + step*(index);
+ e = pow(10.0,e/10.0);
+
+ return e;
+}
+
+/*---------------------------------------------------------------------------*\
+
+ FUNCTION....: encode_amplitudes()
+ AUTHOR......: David Rowe
+ DATE CREATED: 22/8/2010
+
+ Time domain LPC is used model the amplitudes which are then
+ converted to LSPs and quantised. So we don't actually encode the
+ amplitudes directly, rather we derive an equivalent representation
+ from the time domain speech.
+
+\*---------------------------------------------------------------------------*/
+
+void encode_amplitudes(int lsp_indexes[],
+ int *energy_index,
+ MODEL *model,
+ float Sn[],
+ float w[])
+{
+ float lsps[LPC_ORD];
+ float ak[LPC_ORD+1];
+ float e;
+
+ e = speech_to_uq_lsps(lsps, ak, Sn, w, LPC_ORD);
+ encode_lsps(lsp_indexes, lsps, LPC_ORD);
+ *energy_index = encode_energy(e);
+}
+
+/*---------------------------------------------------------------------------*\
+
+ FUNCTION....: decode_amplitudes()
+ AUTHOR......: David Rowe
+ DATE CREATED: 22/8/2010
+
+ Given the amplitude quantiser indexes recovers the harmonic
+ amplitudes.
+
+\*---------------------------------------------------------------------------*/
+
+float decode_amplitudes(MODEL *model,
+ float ak[],
+ int lsp_indexes[],
+ int energy_index,
+ float lsps[],
+ float *e
+)
+{
+ float snr;
+
+ decode_lsps(lsps, lsp_indexes, LPC_ORD);
+ bw_expand_lsps(lsps, LPC_ORD);
+ lsp_to_lpc(lsps, ak, LPC_ORD);
+ *e = decode_energy(energy_index);
+ aks_to_M2(ak, LPC_ORD, model, *e, &snr, 1);
+ apply_lpc_correction(model);
+
+ return snr;
+}
--
cgit v1.2.3