?? openlpc.fixed.c
字號:
temp >>= PRECISION;
temp = (temp - ftofix32(.0170881256f)) * x;
temp >>= PRECISION;
temp = (temp + ftofix32(.0308918810f)) * x;
temp >>= PRECISION;
temp = (temp - ftofix32(.0501743046f)) * x;
temp >>= PRECISION;
temp = (temp + ftofix32(.0889789874f)) * x;
temp >>= PRECISION;
temp = (temp - ftofix32(.2145988016f)) * x;
temp >>= PRECISION;
temp = (temp + ftofix32(1.570796305f)) * fixsqrt32(itofix32(1) - x);
temp >>= PRECISION;
return sign * (ftofix32(M_PI/2) - (fixed32)temp);
}
#define PREEMPH
#define bcopy(a, b, n) memmove(b, a, n)
#define LPC_FILTORDER 10
#define FS 8000.0 /* Sampling rate */
#define MAXWINDOW 1000 /* Max analysis window length */
typedef struct openlpc_e_state{
int framelen, buflen;
fixed32 s[MAXWINDOW], y[MAXWINDOW], h[MAXWINDOW];
fixed32 xv1[3], yv1[3],
xv2[2], yv2[2],
xv3[1], yv3[3],
xv4[2], yv4[2];
fixed32 w[MAXWINDOW], r[LPC_FILTORDER+1];
} openlpc_e_state_t;
typedef struct openlpc_d_state{
fixed32 Oldper, OldG, Oldk[LPC_FILTORDER + 1];
fixed32 bp[LPC_FILTORDER+1];
fixed32 exc;
fixed32 gainadj;
int pitchctr, framelen, buflen;
} openlpc_d_state_t;
#define FC 200.0 /* Pitch analyzer filter cutoff */
#define DOWN 5 /* Decimation for pitch analyzer */
#define MINPIT 40.0 /* Minimum pitch (observed: 74) */
#define MAXPIT 320.0 /* Maximum pitch (observed: 250) */
#define MINPER (int)(FS / (DOWN * MAXPIT) + .5) /* Minimum period */
#define MAXPER (int)(FS / (DOWN * MINPIT) + .5) /* Maximum period */
#define REAL_MINPER (DOWN * MINPER) /* converted to samples units */
#define WSCALE 1.5863 /* Energy loss due to windowing */
#define BITS_FOR_LPC 38
#define ARCSIN_Q /* provides better quantization of first two k[] at low bitrates */
#if BITS_FOR_LPC == 38
/* (38 bit LPC-10, 2.7 Kbit/s @ 20ms, 2.4 Kbit/s @ 22.5 ms */
static int parambits[LPC_FILTORDER] = {6,5,5,4,4,3,3,3,3,2};
#elif BITS_FOR_LPC == 32
/* (32 bit LPC-10, 2.4 Kbit/s, not so good */
static int parambits[LPC_FILTORDER] = {5,5,5,4,3,3,2,2,2,1};
#else /* BITS_FOR_LPC == 80 */
/* 80-bit LPC10, 4.8 Kbit/s */
static int parambits[LPC_FILTORDER] = {8,8,8,8,8,8,8,8,8,8};
#endif
static fixed32 logmaxminper;
static int sizeofparm; /* computed by lpc_init */
static void auto_correl1(fixed32 *w, int n, fixed32 *r)
{
int i, k;
fixed64 temp, temp2;
for (k=0; k <= MAXPER; k++, n--) {
temp = 0;
for (i=0; i < n; i++) {
temp2 = w[i];
temp += temp2 * w[i+k];
}
r[k] = (fixed32)(temp >> PRECISION);
}
}
static void auto_correl2(fixed32 *w, int n, fixed32 *r)
{
int i, k;
fixed64 temp, temp2;
for (k=0; k <= LPC_FILTORDER; k++, n--) {
temp = 0;
for (i=0; i < n; i++) {
temp2 = w[i];
temp += temp2 * w[i+k];
}
r[k] = (fixed32)(temp >> PRECISION);
}
}
static void durbin(fixed32 r[], int p, fixed32 k[], fixed32 *g)
{
int i, j;
fixed32 a[LPC_FILTORDER+1], at[LPC_FILTORDER+1], e;
for (i=0; i <= p; i++)
a[i] = at[i] = 0;
e = r[0];
for (i=1; i <= p; i++) {
k[i] = -r[i];
for (j=1; j < i; j++) {
at[j] = a[j];
k[i] -= fixmul32(a[j], r[i-j]);
}
if (e == 0) { /* fix by John Walker */
*g = 0;
return;
}
k[i] = fixdiv32(k[i], e);
a[i] = k[i];
for (j=1; j < i; j++)
a[j] = at[j] + fixmul32(k[i], at[i-j]);
e = fixmul32(e, (itofix32(1) - fixmul32(k[i], k[i])));
}
if (e < 0) {
e = 0; /* fix by John Walker */
}
*g = fixsqrt32(e);
}
static void calc_pitch(fixed32 w[], int len, fixed32 *per)
{
int i, j, rpos;
fixed32 d[MAXWINDOW / DOWN], r[MAXPER + 1], rmax;
fixed32 rval, rm, rp;
fixed32 x, y;
fixed32 vthresh;
/* decimation */
for (i=0, j=0; i < len; i+=DOWN)
d[j++] = w[i];
auto_correl1(d, len / DOWN, r);
/* find peak between MINPER and MAXPER */
x = itofix32(1);
rpos = 0;
rmax = 0;
for (i = 1; i <= MAXPER; i++) {
rm = r[i-1];
rp = r[i+1];
y = rm+r[i]+rp; /* find max of integral from i-1 to i+1 */
if (y > rmax && r[i] > rm && r[i] > rp && i > MINPER) {
rmax = y;
rpos = i;
}
}
/* consider adjacent values */
rm = r[rpos-1];
rp = r[rpos+1];
if(rpos > 0) {
x = fixdiv32(((rpos-1) * rm + rpos * r[rpos] + (rpos+1) * rp), (rm+r[rpos]+rp));
}
/* normalize, so that 0. < rval < 1. */
rval = (r[0] == 0 ? 0 : fixdiv32(r[rpos], r[0]));
/* periods near the low boundary and at low volumes
are usually spurious and
manifest themselves as annoying mosquito buzzes */
*per = 0; /* default: unvoiced */
if ( x > itofix32(MINPER) && /* x could be < MINPER or even < 0 if rpos == MINPER */
x < itofix32(MAXPER + 1) /* same story */
) {
vthresh = ftofix32(0.6);
if(r[0] > ftofix32(0.002)) /* at low volumes (< 0.002), prefer unvoiced */
vthresh = ftofix32(0.25); /* drop threshold at high volumes */
if(rval > vthresh)
*per = x * DOWN;
}
}
/* Initialization of various parameters */
openlpc_encoder_state *create_openlpc_encoder_state(void)
{
openlpc_encoder_state *state;
state = (openlpc_encoder_state *)malloc(sizeof(openlpc_encoder_state));
return state;
}
void init_openlpc_encoder_state(openlpc_encoder_state *st, int framelen)
{
int i, j;
st->framelen = framelen;
memset(st->y, 0, sizeof(st->y));
st->buflen = framelen * 3 / 2;
/* (st->buflen > MAXWINDOW) return -1;*/
for(i=0, j=0; i<sizeof(parambits)/sizeof(parambits[0]); i++) {
j += parambits[i];
}
sizeofparm = (j + 7) / 8 + 2;
for (i = 0; i < st->buflen; i++) {
st->s[i] = 0;
/* this is only calculated once, but used each frame, */
/* so we will use floating point for accuracy */
st->h[i] = ftofix32(WSCALE*(0.54 - 0.46 * cos(2 * M_PI * i / (st->buflen-1.0))));
}
/* init the filters */
st->xv1[0] = st->xv1[1] = st->xv1[2] = st->yv1[0] = st->yv1[1] = st->yv1[2] = 0;
st->xv2[0] = st->xv2[1] = st->yv2[0] = st->yv2[1] = 0;
st->xv3[0] = st->yv3[0] = st->yv3[1] = st->yv3[2] = 0;
st->xv4[0] = st->xv4[1] = st->yv4[0] = st->yv4[1] = 0;
logmaxminper = fixlog32(fixdiv32(itofix32(MAXPER), itofix32(MINPER)));
}
void destroy_openlpc_encoder_state(openlpc_encoder_state *st)
{
if(st != NULL)
{
free(st);
st = NULL;
}
}
/* LPC Analysis (compression) */
int openlpc_encode(const short *buf, unsigned char *parm, openlpc_encoder_state *st)
{
int i, j;
fixed32 per, gain, k[LPC_FILTORDER+1];
fixed32 per1, per2;
fixed32 xv10, xv11, xv12, yv10, yv11, yv12, xv30, yv30, yv31, yv32;
#ifdef PREEMPH
fixed32 xv20, xv21, yv20, yv21, xv40, xv41, yv40, yv41;
#endif
xv10 = st->xv1[0];
xv11 = st->xv1[1];
xv12 = st->xv1[2];
yv10 = st->yv1[0];
yv11 = st->yv1[1];
yv12 = st->yv1[2];
xv30 = st->xv3[0];
yv30 = st->yv3[0];
yv31 = st->yv3[1];
yv32 = st->yv3[2];
/* convert short data in buf[] to signed lin. data in s[] and prefilter */
for (i=0, j=st->buflen - st->framelen; i < st->framelen; i++, j++) {
/* special handling here for the intitial conversion */
fixed32 u = (fixed32)(buf[i] << (PRECISION - 15));
/* Anti-hum 2nd order Butterworth high-pass, 100 Hz corner frequency */
/* Digital filter designed by mkfilter/mkshape/gencode A.J. Fisher
mkfilter -Bu -Hp -o 2 -a 0.0125 -l -z */
xv10 = xv11;
xv11 = xv12;
#ifdef FAST_FILTERS
xv12 = ((u * 15) >> 4) + (u >> 7) + ((u * 11) >> 14); /* /GAIN */
yv10 = yv11;
yv11 = yv12;
yv12 = (fixed32)((xv10 + xv12) - (xv11 + xv11)
- ((yv10 * 7) >> 3) - ((yv10 * 5) >> 8)
+ ((yv11 * 15) >> 3) + (yv11 >> 6) );
#else
xv12 = fixmul32(u, ftofix32(0.94597831)); /* /GAIN */
yv10 = yv11;
yv11 = yv12;
yv12 = (fixed32)((xv10 + xv12) - (xv11 + xv11)
+ fixmul32(ftofix32(-0.8948742499), yv10) + fixmul32(ftofix32(1.8890389823), yv11));
#endif
u = st->s[j] = yv12; /* also affects input of next stage, to the LPC filter synth */
/* low-pass filter s[] -> y[] before computing pitch */
/* second-order Butterworth low-pass filter, corner at 300 Hz */
/* Digital filter designed by mkfilter/mkshape/gencode A.J. Fisher
MKFILTER.EXE -Bu -Lp -o 2 -a 0.0375 -l -z */
#ifdef FAST_FILTERS
xv30 = ((u * 3) >> 6) + (u >> 13); /* GAIN */
yv30 = yv31;
yv31 = yv32;
yv32 = xv30 - ((yv30 * 23) >> 5) + (yv30 >> 9)
+ ((yv31 * 107) >> 6) - (yv31 >> 9);
#else
xv30 = fixmul32(u, ftofix32(0.04699658)); /* GAIN */
yv30 = yv31;
yv31 = yv32;
yv32 = xv30 + fixmul32(ftofix32(-0.7166152306), yv30) + fixmul32(ftofix32(1.6696186545), yv31);
#endif
st->y[j] = yv32;
}
st->xv1[0] = xv10;
st->xv1[1] = xv11;
st->xv1[2] = xv12;
st->yv1[0] = yv10;
st->yv1[1] = yv11;
st->yv1[2] = yv12;
st->xv3[0] = xv30;
st->yv3[0] = yv30;
st->yv3[1] = yv31;
st->yv3[2] = yv32;
#ifdef PREEMPH
/* operate optional preemphasis s[] -> s[] on the newly arrived frame */
xv20 = st->xv2[0];
xv21 = st->xv2[1];
yv20 = st->yv2[0];
yv21 = st->yv2[1];
xv40 = st->xv4[0];
xv41 = st->xv4[1];
yv40 = st->yv4[0];
yv41 = st->yv4[1];
for (j=st->buflen - st->framelen; j < st->buflen; j++) {
fixed32 u = st->s[j];
/* handcoded filter: 1 zero at 640 Hz, 1 pole at 3200 */
#define TAU (FS / 3200.f)
#define RHO (0.1f)
xv20 = xv21; /* e(n-1) */
#ifdef FAST_FILTERS
xv21 = ((u * 3) >> 1) +((u * 43) >> 9); /* e(n) , add 4 dB to compensate attenuation */
yv20 = yv21;
yv21 = ((yv20 * 11) >> 4) + ((yv20 * 7) >> 10) /* u(n) */
+ ((xv21 * 23) >> 5) + ((xv21 * 7) >> 11)
- ((xv20 * 11) >> 4) - ((xv20 * 7) >> 10);
#else
xv21 = fixmul32(u, ftofix32(1.584)); /* e(n) , add 4 dB to compensate attenuation */
yv20 = yv21;
yv21 = fixmul32(ftofix32(TAU/(1.0f+RHO+TAU)), yv20) /* u(n) */
+ fixmul32(ftofix32((RHO+TAU)/(1.0f+RHO+TAU)), xv21)
- fixmul32(ftofix32(TAU/(1.0f+RHO+TAU)), xv20);
#endif
u = yv21;
/* cascaded copy of handcoded filter: 1 zero at 640 Hz, 1 pole at 3200 */
?? 快捷鍵說明
復制代碼
Ctrl + C
搜索代碼
Ctrl + F
全屏模式
F11
切換主題
Ctrl + Shift + D
顯示快捷鍵
?
增大字號
Ctrl + =
減小字號
Ctrl + -