134680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)/* 234680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) Copyright 1992, 1993, 1994 by Jutta Degener and Carsten Bormann, 334680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) Technische Universitaet Berlin 434680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 534680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) Any use of this software is permitted provided that this notice is not 634680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) removed and that neither the authors nor the Technische Universitaet Berlin 734680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) are deemed to have made any representations as to the suitability of this 834680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) software for any purpose nor are held responsible for any defects of 934680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) this software. THERE IS ABSOLUTELY NO WARRANTY FOR THIS SOFTWARE. 1034680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 1134680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) As a matter of courtesy, the authors request to be informed about uses 1234680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) this software has found, about bugs in this software, and about any 1334680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) improvements that may be of general interest. 1434680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 1534680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) Berlin, 28.11.1994 1634680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) Jutta Degener 1734680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) Carsten Bormann 1834680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 1934680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 2034680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) Code modified by Jean-Marc Valin 2134680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 2234680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) Speex License: 2334680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 2434680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) Redistribution and use in source and binary forms, with or without 2534680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) modification, are permitted provided that the following conditions 2634680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) are met: 2734680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 2834680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) - Redistributions of source code must retain the above copyright 2934680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) notice, this list of conditions and the following disclaimer. 3034680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 3134680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) - Redistributions in binary form must reproduce the above copyright 3234680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) notice, this list of conditions and the following disclaimer in the 3334680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) documentation and/or other materials provided with the distribution. 3434680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 3534680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) - Neither the name of the Xiph.org Foundation nor the names of its 3634680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) contributors may be used to endorse or promote products derived from 3734680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) this software without specific prior written permission. 3834680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 3934680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS 4034680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT 4134680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR 4234680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE FOUNDATION OR 4334680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, 4434680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, 4534680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR 4634680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF 4734680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING 4834680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS 4934680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. 5034680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)*/ 5134680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 5234680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)#ifdef HAVE_CONFIG_H 5334680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)#include "config.h" 5434680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)#endif 5534680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 5634680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)#include "lpc.h" 5734680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 5834680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)#ifdef BFIN_ASM 5934680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)#include "lpc_bfin.h" 6034680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)#endif 6134680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 6234680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)/* LPC analysis 6334680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) * 6434680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) * The next two functions calculate linear prediction coefficients 6534680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) * and/or the related reflection coefficients from the first P_MAX+1 6634680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) * values of the autocorrelation function. 6734680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) */ 6834680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 6934680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)/* Invented by N. Levinson in 1947, modified by J. Durbin in 1959. 7034680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) */ 7134680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 7234680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)/* returns minimum mean square error */ 7334680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)spx_word32_t _spx_lpc( 7434680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)spx_coef_t *lpc, /* out: [0...p-1] LPC coefficients */ 7534680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)const spx_word16_t *ac, /* in: [0...p] autocorrelation values */ 7634680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)int p 7734680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)) 7834680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles){ 7934680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) int i, j; 8034680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) spx_word16_t r; 8134680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) spx_word16_t error = ac[0]; 8234680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 8334680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) if (ac[0] == 0) 8434680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) { 8534680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) for (i = 0; i < p; i++) 8634680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) lpc[i] = 0; 8734680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) return 0; 8834680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) } 8934680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 9034680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) for (i = 0; i < p; i++) { 9134680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 9234680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) /* Sum up this iteration's reflection coefficient */ 9334680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) spx_word32_t rr = NEG32(SHL32(EXTEND32(ac[i + 1]),13)); 9434680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) for (j = 0; j < i; j++) 9534680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) rr = SUB32(rr,MULT16_16(lpc[j],ac[i - j])); 9634680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)#ifdef FIXED_POINT 9734680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) r = DIV32_16(rr+PSHR32(error,1),ADD16(error,8)); 9834680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)#else 9934680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) r = rr/(error+.003*ac[0]); 10034680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)#endif 10134680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) /* Update LPC coefficients and total error */ 10234680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) lpc[i] = r; 10334680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) for (j = 0; j < i>>1; j++) 10434680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) { 10534680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) spx_word16_t tmp = lpc[j]; 10634680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) lpc[j] = MAC16_16_P13(lpc[j],r,lpc[i-1-j]); 10734680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) lpc[i-1-j] = MAC16_16_P13(lpc[i-1-j],r,tmp); 10834680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) } 10934680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) if (i & 1) 11034680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) lpc[j] = MAC16_16_P13(lpc[j],lpc[j],r); 11134680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 11234680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) error = SUB16(error,MULT16_16_Q13(r,MULT16_16_Q13(error,r))); 11334680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) } 11434680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) return error; 11534680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)} 11634680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 11734680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 11834680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)#ifdef FIXED_POINT 11934680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 12034680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)/* Compute the autocorrelation 12134680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) * ,--, 12234680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) * ac(i) = > x(n) * x(n-i) for all n 12334680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) * `--' 12434680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) * for lags between 0 and lag-1, and x == 0 outside 0...n-1 12534680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) */ 12634680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 12734680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)#ifndef OVERRIDE_SPEEX_AUTOCORR 12834680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)void _spx_autocorr( 12934680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)const spx_word16_t *x, /* in: [0...n-1] samples x */ 13034680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)spx_word16_t *ac, /* out: [0...lag-1] ac values */ 13134680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)int lag, 13234680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)int n 13334680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)) 13434680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles){ 13534680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) spx_word32_t d; 13634680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) int i, j; 13734680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) spx_word32_t ac0=1; 13834680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) int shift, ac_shift; 13934680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 14034680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) for (j=0;j<n;j++) 14134680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) ac0 = ADD32(ac0,SHR32(MULT16_16(x[j],x[j]),8)); 14234680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) ac0 = ADD32(ac0,n); 14334680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) shift = 8; 14434680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) while (shift && ac0<0x40000000) 14534680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) { 14634680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) shift--; 14734680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) ac0 <<= 1; 14834680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) } 14934680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) ac_shift = 18; 15034680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) while (ac_shift && ac0<0x40000000) 15134680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) { 15234680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) ac_shift--; 15334680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) ac0 <<= 1; 15434680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) } 15534680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 15634680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 15734680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) for (i=0;i<lag;i++) 15834680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) { 15934680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) d=0; 16034680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) for (j=i;j<n;j++) 16134680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) { 16234680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) d = ADD32(d,SHR32(MULT16_16(x[j],x[j-i]), shift)); 16334680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) } 16434680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 16534680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) ac[i] = SHR32(d, ac_shift); 16634680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) } 16734680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)} 16834680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)#endif 16934680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 17034680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 17134680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)#else 17234680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 17334680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 17434680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 17534680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)/* Compute the autocorrelation 17634680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) * ,--, 17734680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) * ac(i) = > x(n) * x(n-i) for all n 17834680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) * `--' 17934680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) * for lags between 0 and lag-1, and x == 0 outside 0...n-1 18034680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) */ 18134680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)void _spx_autocorr( 18234680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)const spx_word16_t *x, /* in: [0...n-1] samples x */ 18334680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)float *ac, /* out: [0...lag-1] ac values */ 18434680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)int lag, 18534680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)int n 18634680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)) 18734680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles){ 18834680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) float d; 18934680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) int i; 19034680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) while (lag--) 19134680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) { 19234680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) for (i = lag, d = 0; i < n; i++) 19334680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) d += x[i] * x[i-lag]; 19434680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) ac[lag] = d; 19534680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) } 19634680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) ac[0] += 10; 19734680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)} 19834680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 19934680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles)#endif 20034680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 20134680572440d7894ef8dafce81d8039ed80726a2Torne (Richard Coles) 202