00001 /* 00002 Copyright 1992, 1993, 1994 by Jutta Degener and Carsten Bormann, 00003 Technische Universitaet Berlin 00004 00005 Any use of this software is permitted provided that this notice is not 00006 removed and that neither the authors nor the Technische Universitaet Berlin 00007 are deemed to have made any representations as to the suitability of this 00008 software for any purpose nor are held responsible for any defects of 00009 this software. THERE IS ABSOLUTELY NO WARRANTY FOR THIS SOFTWARE. 00010 00011 As a matter of courtesy, the authors request to be informed about uses 00012 this software has found, about bugs in this software, and about any 00013 improvements that may be of general interest. 00014 00015 Berlin, 28.11.1994 00016 Jutta Degener 00017 Carsten Bormann 00018 */ 00019 00020 00021 /* LPC- and Reflection Coefficients 00022 * 00023 * The next two functions calculate linear prediction coefficients 00024 * and/or the related reflection coefficients from the first P_MAX+1 00025 * values of the autocorrelation function. 00026 */ 00027 00028 /* Invented by N. Levinson in 1947, modified by J. Durbin in 1959. 00029 */ 00030 inline float /* returns minimum mean square error */ 00031 wld( 00032 float * lpc, /* [0...p-1] LPC coefficients */ 00033 const float * ac, /* in: [0...p] autocorrelation values */ 00034 float * ref, /* out: [0...p-1] reflection coef's */ 00035 int p 00036 ) 00037 { 00038 int i, j; float r, error = ac[0]; 00039 00040 lpc[0]=1; //for compatibility 00041 lpc++; //for compatibility 00042 00043 if (ac[0] == 0) { 00044 for (i = 0; i < p; i++) ref[i] = 0; return 0; } 00045 00046 for (i = 0; i < p; i++) { 00047 00048 /* Sum up this iteration's reflection coefficient. 00049 */ 00050 r = -ac[i + 1]; 00051 for (j = 0; j < i; j++) r -= lpc[j] * ac[i - j]; 00052 ref[i] = r /= error; 00053 00054 /* Update LPC coefficients and total error. 00055 */ 00056 lpc[i] = r; 00057 for (j = 0; j < i/2; j++) { 00058 float tmp = lpc[j]; 00059 lpc[j] += r * lpc[i-1-j]; 00060 lpc[i-1-j] += r * tmp; 00061 } 00062 if (i % 2) lpc[j] += lpc[j] * r; 00063 00064 error *= 1.0 - r * r; 00065 } 00066 return error; 00067 } 00068 00069 00070 /* Compute the autocorrelation 00071 * ,--, 00072 * ac(i) = > x(n) * x(n-i) for all n 00073 * `--' 00074 * for lags between 0 and lag-1, and x == 0 outside 0...n-1 00075 */ 00076 inline void autocorr( 00077 const float * x, /* in: [0...n-1] samples x */ 00078 float *ac, /* out: [0...lag-1] ac values */ 00079 int lag, int n) 00080 { 00081 float d; int i; 00082 lag++; //Conpensate for the old routine 00083 while (lag--) { 00084 for (i = lag, d = 0; i < n; i++) d += x[i] * x[i-lag]; 00085 ac[lag] = d; 00086 } 00087 }