📄 npp.c
字号:
/* ================================================================== */
/* */
/* Microsoft Speech coder ANSI-C Source Code */
/* SC1200 1200 bps speech coder */
/* Fixed Point Implementation Version 7.0 */
/* Copyright (C) 2000, Microsoft Corp. */
/* All rights reserved. */
/* */
/* ================================================================== */
/*---------------------------------------------------------------------
* enh_fun.c - Speech Enhancement Functions
*
* Author: Rainer Martin, AT&T Labs-Research
*
* Last Update: $Id: enh_fun.c,v 1.2 1999/02/19 17:55:12 martinr Exp martinr $
*
*---------------------------------------------------------------------
*/
#include "npp.h"
#include "fs_lib.h"
#include "mat_lib.h"
#include "constant.h"
#include "dsp_sub.h"
#include "fft_lib.h"
#include "global.h"
#define X003_Q15 983 /* 0.03 * (1 << 15) */
#define X005_Q15 1638 /* 0.05 * (1 << 15) */
#define X006_Q15 1966 /* 0.06 * (1 << 15) */
#define X018_Q15 5898 /* 0.18 * (1 << 15) */
#define X065_Q15 21299 /* 0.65 * (1 << 15) */
#define X132_Q11 2703 /* 1.32 * (1 << 11) */
#define X15_Q13 12288 /* 1.5 * (1 << 13) */
#define X22_Q11 4506 /* 2.2 * (1 << 11) */
#define X44_Q11 9011 /* 4.4 * (1 << 11) */
#define X88_Q11 18022 /* 8.8 * (1 << 11) */
const static Shortword sqrt_tukey_256_180[ENH_WINLEN] = { /* Q15 */
677, 1354, 2030, 2705, 3380, 4053, 4724, 5393,
6060, 6724, 7385, 8044, 8698, 9349, 9996, 10639,
11277, 11911, 12539, 13162, 13780, 14391, 14996, 15595,
16188, 16773, 17351, 17922, 18485, 19040, 19587, 20126,
20656, 21177, 21690, 22193, 22686, 23170, 23644, 24108,
24561, 25004, 25437, 25858, 26268, 26668, 27056, 27432,
27796, 28149, 28490, 28818, 29134, 29438, 29729, 30008,
30273, 30526, 30766, 30992, 31205, 31405, 31592, 31765,
31924, 32070, 32202, 32321, 32425, 32516, 32593, 32656,
32705, 32740, 32761, 32767, 32767, 32767, 32767, 32767,
32767, 32767, 32767, 32767, 32767, 32767, 32767, 32767,
32767, 32767, 32767, 32767, 32767, 32767, 32767, 32767,
32767, 32767, 32767, 32767, 32767, 32767, 32767, 32767,
32767, 32767, 32767, 32767, 32767, 32767, 32767, 32767,
32767, 32767, 32767, 32767, 32767, 32767, 32767, 32767,
32767, 32767, 32767, 32767, 32767, 32767, 32767, 32767,
32767, 32767, 32767, 32767, 32767, 32767, 32767, 32767,
32767, 32767, 32767, 32767, 32767, 32767, 32767, 32767,
32767, 32767, 32767, 32767, 32767, 32767, 32767, 32767,
32767, 32767, 32767, 32767, 32767, 32767, 32767, 32767,
32767, 32767, 32767, 32767, 32767, 32767, 32767, 32767,
32767, 32767, 32767, 32767, 32767, 32767, 32767, 32767,
32767, 32767, 32767, 32767, 32761, 32740, 32705, 32656,
32593, 32516, 32425, 32321, 32202, 32070, 31924, 31765,
31592, 31405, 31205, 30992, 30766, 30526, 30273, 30008,
29729, 29438, 29134, 28818, 28490, 28149, 27796, 27432,
27056, 26668, 26268, 25858, 25437, 25004, 24561, 24108,
23644, 23170, 22686, 22193, 21690, 21177, 20656, 20126,
19587, 19040, 18485, 17922, 17351, 16773, 16188, 15595,
14996, 14391, 13780, 13162, 12539, 11911, 11277, 10639,
9996, 9349, 8698, 8044, 7385, 6724, 6060, 5393,
4724, 4053, 3380, 2705, 2030, 1354, 677, 0
};
/* ====== Entities from Enhanced_Data ====== */
static Shortword enh_i = 0; /* formerly D->I */
static Shortword lambdaD[ENH_VEC_LENF]; /* overestimated noise */
/* psd(noise_bias * noisespect)*/
static Shortword lambdaD_shift[ENH_VEC_LENF];
static Shortword SN_LT; /* long term SNR */
static Shortword SN_LT_shift;
static Shortword n_pwr;
static Shortword n_pwr_shift;
static Shortword YY[ENH_VEC_LENF]; /* signal periodogram of current frame */
static Shortword YY_shift[ENH_VEC_LENF];
static Shortword sm_shift[ENH_VEC_LENF];
static Shortword noise_shift[ENH_VEC_LENF];
static Shortword noise2_shift[ENH_VEC_LENF];
static Shortword av_shift[ENH_VEC_LENF];
static Shortword av2_shift[ENH_VEC_LENF];
static Shortword act_min[ENH_VEC_LENF]; /* current minimum of long window */
static Shortword act_min_shift[ENH_VEC_LENF];
static Shortword act_min_sub[ENH_VEC_LENF];
/* current minimum of sub-window */
static Shortword act_min_sub_shift[ENH_VEC_LENF];
static Shortword vk[ENH_VEC_LENF];
static Shortword vk_shift[ENH_VEC_LENF];
static Shortword ksi[ENH_VEC_LENF]; /* a priori SNR */
static Shortword ksi_shift[ENH_VEC_LENF];
static Shortword var_rel[ENH_VEC_LENF];
/* est. relative var. of smoothedspect */
static Shortword var_rel_av; /* average relative var. of smoothedspect */
/* only for minimum statistics */
static Shortword smoothedspect[ENH_VEC_LENF]; /* smoothed signal spectrum */
static Shortword var_sp_av[ENH_VEC_LENF];
/* estimated mean of smoothedspect */
static Shortword var_sp_2[ENH_VEC_LENF];
/* esitmated 2nd moment of smoothedspect */
static Shortword noisespect[ENH_VEC_LENF]; /* estimated noise psd */
static Shortword noisespect2[ENH_VEC_LENF];
static Shortword alphacorr; /* correction factor for alpha_var, Q15 */
static Shortword alpha_var[ENH_VEC_LENF]; /* variable smoothing */
/* parameter for psd */
static Shortword circb[NUM_MINWIN][ENH_VEC_LENF]; /* ring buffer */
static Shortword circb_shift[NUM_MINWIN][ENH_VEC_LENF];
static Shortword initial_noise[ENH_WINLEN]; /* look ahead noise */
/* samples (Q0) */
static Shortword speech_in_npp[ENH_WINLEN]; /* input of one frame */
static Shortword ybuf[2*ENH_WINLEN + 2];
/* buffer for FFT, this can be eliminated if */
/* we can write a better real-FFT program for DSP */
static Longword temp_yy[ENH_WINLEN + 2];
/* ====== Prototypes ====== */
static void smoothing_win(Shortword initial_noise[]);
static void compute_qk(Shortword qk[], Shortword gamaK[],
Shortword gamaK_shift[], Shortword GammaQ_TH);
static void gain_log_mmse(Shortword qk[], Shortword Gain[],
Shortword gamaK[], Shortword gamaK_shift[],
Shortword m);
static Shortword ksi_min_adapt(BOOLEAN n_flag, Shortword ksi_min,
Shortword sn_lt, Shortword sn_lt_shift);
static void smoothed_periodogram(Shortword YY_av, Shortword yy_shift);
static void bias_compensation(Shortword biased_spect[],
Shortword bias_shift[],
Shortword biased_spect_sub[],
Shortword bias_sub_shift[]);
static Shortword noise_slope();
static Shortword comp_data_shift(Shortword num1, Shortword shift1,
Shortword num2, Shortword shift2);
static void min_search(Shortword biased_spect[], Shortword bias_shift[],
Shortword biased_spect_sub[],
Shortword bias_sub_shift[]);
void enh_init();
static void minstat_init();
static void process_frame(Shortword inspeech[], Shortword outspeech[]);
static void gain_mod(Shortword qk[], Shortword GainD[], Shortword m);
/****************************************************************************
**
** Function: npp
**
** Description: The noise pre-processor
**
** Arguments:
**
** Shortword sp_in[] ---- input speech data buffer (Q0)
** Shortword sp_out[] ---- output speech data buffer (Q0)
**
** Return value: None
**
*****************************************************************************/
void npp(Shortword sp_in[], Shortword sp_out[])
{
static Shortword speech_overlap[ENH_OVERLAP_LEN];
static BOOLEAN first_time = TRUE;
Shortword speech_out_npp[ENH_WINLEN]; /* output of one frame */
if (first_time){
if (rate == RATE1200)
v_equ(initial_noise, sp_in, ENH_WINLEN);
else {
v_zap(initial_noise, ENH_OVERLAP_LEN);
v_equ(&initial_noise[ENH_OVERLAP_LEN], sp_in, ENH_WIN_SHIFT);
}
enh_init(); /* Initialize enhancement routine */
v_zap(speech_in_npp, ENH_WINLEN);
first_time = FALSE;
}
/* Shift input buffer from the previous frame */
v_equ(speech_in_npp, &(speech_in_npp[ENH_WIN_SHIFT]), ENH_OVERLAP_LEN);
v_equ(&(speech_in_npp[ENH_OVERLAP_LEN]), sp_in, ENH_WIN_SHIFT);
/* ======== Process one frame ======== */
process_frame(speech_in_npp, speech_out_npp);
/* Overlap-add the output buffer. */
v_add(speech_out_npp, speech_overlap, ENH_OVERLAP_LEN);
v_equ(speech_overlap, &(speech_out_npp[ENH_WIN_SHIFT]), ENH_OVERLAP_LEN);
/* ======== Output enhanced speech ======== */
v_equ(sp_out, speech_out_npp, ENH_WIN_SHIFT);
}
/*****************************************************************************/
/* Subroutine gain_mod: compute gain modification factor based on */
/* signal absence probabilities qk */
/*****************************************************************************/
static void gain_mod(Shortword qk[], Shortword GainD[], Shortword m)
{
register Shortword i;
Shortword temp, temp1, temp2, temp3, temp4, shift;
Shortword temp_shift, temp2_shift;
Longword L_sum, L_tmp;
for (i = 0; i < m; i++){
/* temp = 1.0 - qk[i]; */
temp = sub(ONE_Q15, qk[i]); /* Q15 */
if (temp == 0)
temp = 1;
shift = norm_s(temp);
temp = shl(temp, shift);
temp_shift = negate(shift);
/* temp2 = temp*temp; */
L_sum = L_mult(temp, temp);
shift = norm_l(L_sum);
temp2 = extract_h(L_shl(L_sum, shift));
temp2_shift = sub(shl(temp_shift, 1), shift);
/* GM[i] = temp2 / (temp2 + (temp + ksi[i])*qk[i]*exp(-vk[i])) */
/* exp(-vk) = 2^(2*vk*(-0.5*log2(e))), and -23637 is -0.5*log2(e). */
L_sum = L_mult(vk[i], -23637); /* vk*(-0.5*log2(e)), Q15 */
shift = add(vk_shift[i], 1);
L_sum = L_shr(L_sum, sub(15, shift)); /* 2^x x Q16 */
/* temp + ksi[i] */
shift = sub(ksi_shift[i], temp_shift);
if (shift > 0){
temp4 = add(ksi_shift[i], 1);
temp3 = add(shr(temp, add(shift, 1)), shr(ksi[i], 1));
} else {
temp4 = add(temp_shift, 1);
temp3 = add(shr(temp, 1), shl(ksi[i], sub(shift, 1)));
}
/* (temp + ksi[i])*qk[i] */
L_tmp = L_mult(temp3, qk[i]);
/* temp + ksi[i])*qk[i]*exp(-vk[i] */
L_sum = L_add(L_sum, L_deposit_h(temp4)); /* add exp */
shift = extract_h(L_sum); /* this is exp part */
temp4 = (Shortword) (extract_l(L_shr(L_sum, 1)) & 0x7fff);
temp1 = shr(mult(temp4, 9864), 3); /* change to base 10 Q12 */
temp1 = pow10_fxp(temp1, 14); /* out Q14 */
L_tmp = L_mpy_ls(L_tmp, temp1); /* Q30 */
temp3 = norm_l(L_tmp);
temp1 = extract_h(L_shl(L_tmp, temp3));
shift = add(shift, sub(1, temp3)); /* make L_tmp Q31 */
/* temp2 + (temp + ksi[i])*qk[i]*exp(-vk[i]) */
temp1 = shr(temp1, 1);
temp2 = shr(temp2, 1);
temp = sub(shift, temp2_shift);
if (temp > 0){
temp3 = add(temp1, shr(temp2, temp));
temp4 = shr(temp2, temp);
} else {
temp3 = add(shl(temp1, temp), temp2);
temp4 = temp2;
}
temp = divide_s(temp4, temp3); /* temp is previously known as GM[]. */
/* limit lowest GM value */
if (temp < GM_MIN)
temp = GM_MIN; /* Q15 */
GainD[i] = mult(GainD[i], temp); /* modified gain */
}
}
/*****************************************************************************/
/* Subroutine compute_qk: compute the probability of speech absence */
/* This uses an harddecision approach due to Malah (ICASSP 1999). */
/*****************************************************************************/
static void compute_qk(Shortword qk[], Shortword gamaK[],
Shortword gamaK_shift[], Shortword GammaQ_TH)
{
register Shortword i;
static BOOLEAN first_time = TRUE;
static Shortword qla[ENH_VEC_LENF];
/* Initializing qla[] */
if (first_time){
fill(qla, X05_Q15, ENH_VEC_LENF);
first_time = FALSE;
}
/* qla[] = alphaq * qla[]; */
v_scale(qla, ENH_ALPHAQ, ENH_VEC_LENF);
for (i = 0; i < ENH_VEC_LENF; i++){
/* if (gamaK[i] < GammaQ_TH) */
if (comp_data_shift(gamaK[i], gamaK_shift[i], GammaQ_TH, 0) < 0)
/* qla[] += betaq; */
qla[i] = add(qla[i], ENH_BETAQ);
}
v_equ(qk, qla, ENH_VEC_LENF);
}
/*****************************************************************************/
/* Subroutine gain_log_mmse: compute the gain factor and the auxiliary */
/* variable vk for the Ephraim&Malah 1985 log spectral estimator. */
/* Approximation of the exponential integral due to Malah, 1996 */
/*****************************************************************************/
static void gain_log_mmse(Shortword qk[], Shortword Gain[],
Shortword gamaK[], Shortword gamaK_shift[],
Shortword m)
{
register Shortword i;
Shortword ksi_vq, temp1, temp2, shift;
Longword L_sum;
for (i = 0; i < m; i++){
/* 1.0 - qk[] */
temp1 = sub(ONE_Q15, qk[i]);
shift = norm_s(temp1);
temp1 = shl(temp1, shift);
temp2 = sub(ksi_shift[i], (Shortword) (-shift));
/* ksi[] + 1.0 - qk[] */
if (temp2 > 0){
temp1 = shr(temp1, add(temp2, 1));
temp1 = add(temp1, shr(ksi[i], 1));
temp2 = shr(ksi[i], 1);
} else {
temp1 = add(shr(temp1,1), shl(ksi[i], sub(temp2, 1)));
temp2 = shl(ksi[i], sub(temp2, 1));
}
/* ksi_vq = ksi[i] / (ksi[i] + 1.0 - qk[i]); */
ksi_vq = divide_s(temp2, temp1); /* Q15 */
/* vk[i] = ksi_vq * gamaK[i]; */
L_sum = L_mult(ksi_vq, gamaK[i]);
shift = norm_l(L_sum);
vk[i] = extract_h(L_shl(L_sum, shift));
vk_shift[i] = sub(gamaK_shift[i], shift);
/* The original floating version compares vk[] against 2^{-52} == */
/* 2.220446049250313e-16. Tian uses 32767*2^{-52} instead. */
if (comp_data_shift(vk[i], vk_shift[i], 32767, -52) < 0){
vk[i] = 32767; /* MATLAB eps */
vk_shift[i] = -52;
}
if (comp_data_shift(vk[i], vk_shift[i], 26214, -3) < 0){
/* eiv = -2.31 * log10(vk[i]) - 0.6; */
temp1 = log10_fxp(vk[i], 15); /* Q12 */
L_sum = L_shl(L_deposit_l(temp1), 14); /* Q26 */
L_sum = L_add(L_sum, L_shl(L_mult(vk_shift[i], 9864), 10));
L_sum = L_mpy_ls(L_sum, -18923);
/* -18923 = -2.31 (Q13). out Q24 */
L_sum = L_sub(L_sum, 10066330L); /* 10066330L = 0.6 (Q24), Q24 */
} else if (comp_data_shift(vk[i], vk_shift[i], 25600, 8) > 0){
/* eiv = pow(10, -0.52 * 200 - 0.26); */
L_sum = 1; /* Q24 */
/* vk[] = 200; */
vk[i] = 25600;
vk_shift[i] = 8;
} else if (comp_data_shift(vk[i], vk_shift[i], 32767, 0) > 0){
/* eiv = pow(10, -0.52 * vk[i] - 0.26); */
L_sum = L_mult(vk[i], -17039); /* -17039 == -0.52 Q15 */
L_sum = L_sub(L_sum, L_shr(L_deposit_h(8520), vk_shift[i]));
/* 8520 == 0.26 Q15 */
L_sum = L_shr(L_sum, sub(14, vk_shift[i]));
L_sum = L_mpy_ls(L_sum, 27213); /* Q15 to base 2 */
shift = extract_h(L_shl(L_sum, 1)); /* integer part */
temp1 = (Shortword) (extract_l(L_sum) & 0x7fff);
temp1 = shr(mult(temp1, 9864), 3); /* Q12 to base 10 */
temp1 = pow10_fxp(temp1, 14); /* output Q14 */
L_sum = L_shl(L_deposit_l(temp1), 10); /* Q24 now */
L_sum = L_shl(L_sum, shift); /* shift must < 0 */
} else {
/* eiv = -1.544 * log10(vk[i]) + 0.166; */
temp1 = vk[i];
if (vk_shift[i] != 0)
temp1 = shl(temp1, vk_shift[i]);
temp1 = log10_fxp(temp1, 15); /* Q12 */
L_sum = L_shl(L_deposit_l(temp1), 13); /* Q25 */
L_sum = L_mpy_ls(L_sum, -25297); /* -25297 = -1.544, Q14. Q24 */
L_sum = L_add(L_sum, 2785018L); /* 2785018L == 0.166, Q24 */
}
/* Now "eiv" is kept in L_sum. */
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -