Remove unused tests from configure.ac
[speexdsp.git] / tmv / lsp_tm.h
1 /* Copyright (C) 2007 Hong Zhiqian */\r
2 /**\r
3    @file lsp_tm.h\r
4    @author Hong Zhiqian\r
5    @brief Various compatibility routines for Speex (TriMedia version)\r
6 */\r
7 /*\r
8    Redistribution and use in source and binary forms, with or without\r
9    modification, are permitted provided that the following conditions\r
10    are met:\r
11    \r
12    - Redistributions of source code must retain the above copyright\r
13    notice, this list of conditions and the following disclaimer.\r
14    \r
15    - Redistributions in binary form must reproduce the above copyright\r
16    notice, this list of conditions and the following disclaimer in the\r
17    documentation and/or other materials provided with the distribution.\r
18    \r
19    - Neither the name of the Xiph.org Foundation nor the names of its\r
20    contributors may be used to endorse or promote products derived from\r
21    this software without specific prior written permission.\r
22    \r
23    THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS\r
24    ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT\r
25    LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR\r
26    A PARTICULAR PURPOSE ARE DISCLAIMED.  IN NO EVENT SHALL THE FOUNDATION OR\r
27    CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,\r
28    EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,\r
29    PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR\r
30    PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF\r
31    LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING\r
32    NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS\r
33    SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.\r
34 */\r
35 \r
36 #include <ops/custom_defs.h>\r
37 #include "profile_tm.h"\r
38 \r
39 \r
40 #ifdef FIXED_POINT\r
41 \r
42 #define OVERRIDE_LSP_INTERPOLATE\r
43 void lsp_interpolate(Int16 *old_lsp, Int16 *new_lsp, Int16 *interp_lsp, int len, int subframe, int nb_subframes)\r
44 {\r
45         register int tmp  =  DIV32_16(SHL32(EXTEND32(1 + subframe),14),nb_subframes);\r
46         register int tmp2 =  16384-tmp;\r
47         register int in_0, in_1, factor, out_1, out_2, olsp, nlsp, ilsp;\r
48         register int i;\r
49 \r
50         TMDEBUG_ALIGNMEM(old_lsp);\r
51         TMDEBUG_ALIGNMEM(new_lsp);\r
52         TMDEBUG_ALIGNMEM(interp_lsp);\r
53 \r
54         LSPINTERPOLATE_START();\r
55 \r
56         factor = pack16lsb(tmp,tmp2);\r
57 \r
58         len >>= 1;\r
59         for ( i=0 ; i<len ; ++i )\r
60         {       \r
61                 olsp = ld32x(old_lsp,i);                                                        // olsp[i+1],olsp[i]\r
62                 nlsp = ld32x(new_lsp,i);                                                        // nlsp[i+1],nlsp[i]\r
63 \r
64                 in_0 = pack16lsb(nlsp,olsp);\r
65                 in_1 = pack16msb(nlsp,olsp);\r
66           \r
67                 out_1 = 8192 + ifir16(in_0,factor);\r
68                 out_2 = 8192 + ifir16(in_1,factor);\r
69 \r
70                 ilsp = pack16lsb(out_2 >> 14, out_1 >> 14);\r
71                 st32d(i << 2, interp_lsp, ilsp);\r
72 \r
73         }\r
74 \r
75         LSPINTERPOLATE_STOP();\r
76 }\r
77 \r
78 \r
79 #define OVERRIDE_CHEB_POLY_EVA\r
80 static inline Int32 cheb_poly_eva(Int16 *coef, Int16 x, int m, char *stack)\r
81 {\r
82         register int    c10, c32, c54;\r
83     register int        sum, b0, f0, f1, f2, f3;\r
84         register int    xx, f32, f10;\r
85 \r
86         CHEBPOLYEVA_START();\r
87 \r
88         xx = sex16(x);\r
89         b0 = iclipi(xx,16383);\r
90 \r
91 #if 0\r
92         c10 = ld32(coef);\r
93         c32 = ld32x(coef,1);\r
94         c54 = ld32x(coef,2);\r
95 #else   \r
96         c10 = pack16lsb(coef[1],coef[0]);\r
97         c32 = pack16lsb(coef[3],coef[2]);\r
98         c54 = pack16lsb(coef[5],coef[4]);\r
99 #endif\r
100 \r
101         f0   = ((xx * b0) >> 13) - 16384;\r
102         f1   = ((xx * f0) >> 13) - b0;\r
103         f2   = ((xx * f1) >> 13) - f0;\r
104 \r
105         if ( m == 4 )\r
106         {       sum = sex16(c54);\r
107                 f32 = pack16lsb(xx,f0);\r
108                 f10 = pack16lsb(f1,f2);\r
109 \r
110         } else\r
111         {       sum =  asri(16,c54);\r
112                 sum += ((sex16(c54) * xx) + 8192) >> 14;\r
113 \r
114                 f3   = ((xx * f2) >> 13) - f1;\r
115                 f32 = pack16lsb(f0,f1);\r
116                 f10 = pack16lsb(f2,f3);\r
117         }\r
118 \r
119         sum += (ifir16(c32,f32) + 8192) >> 14;\r
120         sum += (ifir16(c10,f10) + 8192) >> 14;\r
121 \r
122 #ifndef REMARK_ON\r
123    (void)stack;\r
124 #endif\r
125 \r
126         CHEBPOLYEVA_STOP();\r
127     return sum;\r
128 }\r
129 \r
130 \r
131 #define OVERRIDE_LSP_ENFORCE_MARGIN\r
132 void lsp_enforce_margin(Int16 *lsp, int len, Int16 margin)\r
133 {\r
134         register int i;\r
135         register int m = margin;\r
136         register int m2 = 25736-margin;\r
137         register int lsp0, lsp1, lsp2;\r
138 \r
139         TMDEBUG_ALIGNMEM(lsp);\r
140 \r
141         LSPENFORCEMARGIN_START();\r
142 \r
143         lsp0 = ld32(lsp);\r
144         lsp1 = asri(16,lsp0);\r
145         lsp0 = sex16(lsp0);\r
146         lsp2 = lsp[len-1];\r
147 \r
148         if ( lsp0 < m )\r
149         {       lsp0 = m;\r
150                 lsp[0] = m;\r
151         }\r
152 \r
153         if ( lsp2 > m2 )\r
154         {       lsp2 = m2;\r
155                 lsp[len-1] = m2;\r
156         }\r
157 \r
158         for ( i=1 ; i<len-1 ; ++i )\r
159         {       \r
160                 lsp0 += m;\r
161                 lsp2 = lsp[i+1];\r
162                 m2   = lsp2 - m;\r
163 \r
164                 if ( lsp1 < lsp0 )\r
165                 {       lsp1   = lsp0;\r
166                         lsp[i] = lsp0;\r
167                 }\r
168 \r
169                 if ( lsp1 > m2 )\r
170                 {       lsp1    = (lsp1 >> 1) + (m2 >> 1);\r
171                         lsp[i]  = lsp1;                 \r
172                 }\r
173 \r
174                 lsp0 = lsp1;\r
175                 lsp1 = lsp2;\r
176         }\r
177 \r
178         LSPENFORCEMARGIN_STOP();\r
179 }\r
180 \r
181 \r
182 #define OVERRIDE_LSP_TO_LPC\r
183 void lsp_to_lpc(Int16 *freq, Int16 *ak,int lpcrdr, char *stack)\r
184 {\r
185     VARDECL(Int16 *freqn);\r
186     VARDECL(int **xp);\r
187     VARDECL(int *xpmem);\r
188     VARDECL(int **xq);\r
189     VARDECL(int *xqmem);\r
190 \r
191     register int i, j, k;\r
192     register int xout1,xout2,xin;\r
193     register int m;\r
194     \r
195         LSPTOLPC_START();\r
196 \r
197         m = lpcrdr>>1;\r
198 \r
199         /* \r
200     \r
201        Reconstruct P(z) and Q(z) by cascading second order polynomials\r
202        in form 1 - 2cos(w)z(-1) + z(-2), where w is the LSP frequency.\r
203        In the time domain this is:\r
204 \r
205        y(n) = x(n) - 2cos(w)x(n-1) + x(n-2)\r
206     \r
207        This is what the ALLOCS below are trying to do:\r
208 \r
209          int xp[m+1][lpcrdr+1+2]; // P matrix in QIMP\r
210          int xq[m+1][lpcrdr+1+2]; // Q matrix in QIMP\r
211 \r
212        These matrices store the output of each stage on each row.  The\r
213        final (m-th) row has the output of the final (m-th) cascaded\r
214        2nd order filter.  The first row is the impulse input to the\r
215        system (not written as it is known).\r
216 \r
217        The version below takes advantage of the fact that a lot of the\r
218        outputs are zero or known, for example if we put an inpulse\r
219        into the first section the "clock" it 10 times only the first 3\r
220        outputs samples are non-zero (it's an FIR filter).\r
221     */\r
222 \r
223     ALLOC(xp, (m+1), int*);\r
224     ALLOC(xpmem, (m+1)*(lpcrdr+1+2), int);\r
225 \r
226     ALLOC(xq, (m+1), int*);\r
227     ALLOC(xqmem, (m+1)*(lpcrdr+1+2), int);\r
228     \r
229     for ( i=0; i<=m; i++ ) \r
230         {       xp[i] = xpmem + i*(lpcrdr+1+2);\r
231                 xq[i] = xqmem + i*(lpcrdr+1+2);\r
232     }\r
233 \r
234     /* work out 2cos terms in Q14 */\r
235 \r
236     ALLOC(freqn, lpcrdr, Int16);\r
237     for ( j=0,k=0 ; j<m ; ++j,k+=2 )\r
238         {       register int f;\r
239 \r
240                 f = ld32x(freq,j);\r
241 \r
242                 freqn[k] = ANGLE2X(sex16(f));\r
243                 freqn[k+1] = ANGLE2X(asri(16,f));\r
244         }\r
245     \r
246         #define QIMP  21   /* scaling for impulse */\r
247     xin = SHL32(EXTEND32(1), (QIMP-1)); /* 0.5 in QIMP format */\r
248    \r
249     /* first col and last non-zero values of each row are trivial */\r
250     \r
251     for(i=0;i<=m;i++) \r
252         {       xp[i][1] = 0;\r
253                 xp[i][2] = xin;\r
254                 xp[i][2+2*i] = xin;\r
255                 xq[i][1] = 0;\r
256                 xq[i][2] = xin;\r
257                 xq[i][2+2*i] = xin;\r
258     }\r
259 \r
260     /* 2nd row (first output row) is trivial */\r
261 \r
262     xp[1][3] = -MULT16_32_Q14(freqn[0],xp[0][2]);\r
263     xq[1][3] = -MULT16_32_Q14(freqn[1],xq[0][2]);\r
264 \r
265     xout1 = xout2 = 0;\r
266 \r
267     for( i=1 ; i<m ; ++i) \r
268         {       register int f, f0, f1, m0, m1;\r
269                 \r
270                 k  = 2*(i+1)-1;\r
271                 f  = ld32x(freqn,i);\r
272                 f0 = sex16(f);\r
273                 f1 = asri(16,f);\r
274 \r
275                 for( j=1 ; j<k ; ++j) \r
276                 {       register int _m0, _m1;  \r
277 \r
278                         _m0 = MULT16_32_Q14(f0,xp[i][j+1]);\r
279                         xp[i+1][j+2] = ADD32(SUB32(xp[i][j+2], _m0), xp[i][j]);\r
280         \r
281                         _m1 = MULT16_32_Q14(f1,xq[i][j+1]);\r
282                         xq[i+1][j+2] = ADD32(SUB32(xq[i][j+2], _m1), xq[i][j]);\r
283                 }\r
284 \r
285       \r
286                 m0 = MULT16_32_Q14(f0,xp[i][j+1]);\r
287                 xp[i+1][j+2] = SUB32(xp[i][j], m0);\r
288                 m1 = MULT16_32_Q14(f1,xq[i][j+1]);\r
289                 xq[i+1][j+2] = SUB32(xq[i][j], m1);\r
290     }\r
291 \r
292 \r
293     for( i=0,j=3 ; i<lpcrdr ; ++j,++i ) \r
294         {       register int _a0, _xp0, _xq0;\r
295 \r
296                 _xp0 = xp[m][j];\r
297                 _xq0 = xq[m][j];\r
298 \r
299         _a0 = PSHR32(_xp0 + xout1 + _xq0 - xout2, QIMP-13); \r
300                 xout1 = _xp0;\r
301                 xout2 = _xq0;\r
302      \r
303                 ak[i] = iclipi(_a0,32767); \r
304     }\r
305 \r
306         LSPTOLPC_STOP();\r
307 }\r
308 \r
309 \r
310 #endif\r