Voice Lab
mel_cepstral
void mel_cepstral(double vector_i[], double vector_o[], int order, double gain)
{
int i, k, m, lpc_index, index;
double alpha;
double previous_am[MAX_ALLOC_SIZE];
double am[MAX_ALLOC_SIZE];
double K, sum;
/*
** Information for this routine taken from:
** "Recursive Calculation of Mel-Cepstrum from LP Coefficients"
** By Keiichi Tokuda, Takao Kobayashi and Satoshi Imai
** April 1994
*/
alpha = 0.455;
K = gain;
if((order+1) >= MAX_ALLOC_SIZE) {
fprintf(stderr, "Please increase MAX_ALLOC_SIZE in window_lib.h to at least %d\n", order+2);
exit(1);
}
/* In equations 17, 18 and 19 in the article, we see the following variables:
** a() = the array of LP coefficients denoted here by "vector_i"
** a_tilde() = the current line of the matrix denoted here by "am"
** a_tilde(i-1)() = the previous line of the matrix denoted here by "previous_am"
** alpha = a constant denoted here by "alpha" and set to 0.455
** K =
** c_tilde() = output coefficient vector denoted here by "vector_o"
** M = number of LP Coefficients denoted here by "p";
*/
for(i = 0; i <= order; i++) {
/* start with a zero vector */
vector_o[i] = 0.0;
previous_am[i] = 0.0;
}
for(i = -order; i <= 0; i++) {
DEBUG_WRITE2(16, " loop i=%d\n", i);
/* added p to i into index so that index is always positive */
index = i + order;
if(index < 0) {
fprintf(stderr, "Error with index (too small)\n");
exit(1);
}
if(index > order) {
fprintf(stderr, "Error with index (too large)\n");
exit(1);
}
if(i < 0) lpc_index = -i;
if(i == 0) lpc_index = 0;
if(lpc_index < 0) {
fprintf(stderr, "Error with lpc_index (too small)\n");
exit(1);
}
if(lpc_index > order) {
fprintf(stderr, "Error with lpc_index (too large)\n");
exit(1);
}
if(index == 0) {
am[0] = vector_i[lpc_index];
for(m = 1; m <= order; m++) {
am[0] = 0;
}
} else if(index == 1) {
am[0] = vector_i[lpc_index] + alpha * previous_am[0];
am[1] = (1 - (alpha * alpha)) * previous_am[0];
for(m = 2; m <= order; m++) {
am[m] = alpha * (-am[m-1]);
}
} else {
am[0] = vector_i[lpc_index] + alpha * previous_am[0];
am[1] = (1 - (alpha * alpha)) * previous_am[0] + alpha * previous_am[1];
for(m = 2; m <= order; m++) {
am[m] = previous_am[m-1] + alpha * (previous_am[m] - am[m-1]);
}
}
for(k = 0; k <= order; k++ ) {
previous_am[k] = am[k];
}
}
vector_o[0] = log(K/ am[0]);
for(m = 1; m <= order; m++) {
sum = 0.0;
for(k = 1; k <= (m - 1); k++) {
sum += ((double) k / (double) m) * vector_o[k] * am[m-k];
}
vector_o[m] = -am[m] - sum;
}
} /* End of the mel cepstral function */
