1 /* Copyright (C) 2005 Analog Devices */
2 /**
3    @file filters_bfin.h
4    @brief Various analysis/synthesis filters (Blackfin version)
5 */
6 /*
7    Redistribution and use in source and binary forms, with or without
8    modification, are permitted provided that the following conditions
9    are met:
10 
11    - Redistributions of source code must retain the above copyright
12    notice, this list of conditions and the following disclaimer.
13 
14    - Redistributions in binary form must reproduce the above copyright
15    notice, this list of conditions and the following disclaimer in the
16    documentation and/or other materials provided with the distribution.
17 
18    - Neither the name of the Xiph.org Foundation nor the names of its
19    contributors may be used to endorse or promote products derived from
20    this software without specific prior written permission.
21 
22    THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
23    ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
24    LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
25    A PARTICULAR PURPOSE ARE DISCLAIMED.  IN NO EVENT SHALL THE FOUNDATION OR
26    CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
27    EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
28    PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
29    PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
30    LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
31    NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
32    SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
33 */
34 
35 #include <stdio.h>
36 
37 #define OVERRIDE_NORMALIZE16
normalize16(const spx_sig_t * x,spx_word16_t * y,spx_sig_t max_scale,int len)38 int normalize16(const spx_sig_t *x, spx_word16_t *y, spx_sig_t max_scale, int len)
39 {
40    spx_sig_t max_val=1;
41    int sig_shift;
42 
43    __asm__
44    (
45    "%0 = 0;\n\t"
46    "I0 = %1;\n\t"
47    "L0 = 0;\n\t"
48    "R1 = [I0++];\n\t"
49    "LOOP norm_max%= LC0 = %2;\n\t"
50    "LOOP_BEGIN norm_max%=;\n\t"
51       "R2 = ABS R1 || R1 = [I0++];\n\t"
52       "%0 = MAX(%0, R2);\n\t"
53    "LOOP_END norm_max%=;\n\t"
54    : "=&d" (max_val)
55    : "a" (x), "a" (len)
56    : "R1", "R2"
57    );
58 
59    sig_shift=0;
60    while (max_val>max_scale)
61    {
62       sig_shift++;
63       max_val >>= 1;
64    }
65 
66    __asm__ __volatile__
67    (
68    "I0 = %0;\n\t"
69    "L0 = 0;\n\t"
70    "I1 = %1;\n\t"
71    "L1 = 0;\n\t"
72    "R0 = [I0++];\n\t"
73    "LOOP norm_shift%= LC0 = %3 >> 1;\n\t"
74    "LOOP_BEGIN norm_shift%=;\n\t"
75       "R1 = ASHIFT R0 by %2.L || R2 = [I0++];\n\t"
76       "R3 = ASHIFT R2 by %2.L || R0 = [I0++];\n\t"
77       "R3 = PACK(R3.L, R1.L);\n\t"
78       "[I1++] = R3;\n\t"
79    "LOOP_END norm_shift%=;\n\t"
80    : : "a" (x), "a" (y), "d" (-sig_shift), "a" (len)
81    : "I0", "L0", "I1", "L1", "R0", "R1", "R2", "R3", "memory"
82    );
83    return sig_shift;
84 }
85 
86 #define OVERRIDE_FILTER_MEM2
filter_mem2(const spx_sig_t * _x,const spx_coef_t * num,const spx_coef_t * den,spx_sig_t * _y,int N,int ord,spx_mem_t * mem)87 void filter_mem2(const spx_sig_t *_x, const spx_coef_t *num, const spx_coef_t *den, spx_sig_t *_y, int N, int ord, spx_mem_t *mem)
88 {
89    spx_word32_t xy2[N+1];
90    spx_word32_t *xy = xy2+1;
91    spx_word32_t numden_a[2*ord+2];
92    spx_word16_t *numden = (spx_word16_t*) numden_a;
93    int i;
94    for (i=0;i<ord;i++)
95    {
96       numden[2*i] = num[i];
97       numden[2*i+1] = den[i];
98    }
99    __asm__ __volatile__
100    (
101    /* Register setup */
102    "R0 = %5;\n\t"      /*ord */
103 
104    "P0 = %3;\n\t"
105    "I0 = P0;\n\t"
106    "B0 = P0;\n\t"
107    "L0 = 0;\n\t"
108 
109    "P2 = %0;\n\t"
110    "I2 = P2;\n\t"
111    "L2 = 0;\n\t"
112 
113    "P4 = %6;\n\t"
114    "P0 = %1;\n\t"
115    "P1 = %2;\n\t"
116 
117    /* First sample */
118    "R1 = [P4++];\n\t"
119    "R1 <<= 1;\n\t"
120    "R2 = [P0++];\n\t"
121    "R1 = R1 + R2;\n\t"
122    "[P1++] = R1;\n\t"
123    "R1 <<= 2;\n\t"
124    "R2 <<= 2;\n\t"
125    "R2 = PACK(R1.H, R2.H);\n\t"
126    "[P2] = R2;\n\t"
127 
128    /* Samples 1 to ord-1 (using memory) */
129    "R0 += -1;\n\t"
130    "R3 = 0;\n\t"
131    "LC0 = R0;\n\t"
132    "LOOP filter_start%= LC0;\n\t"
133    "LOOP_BEGIN filter_start%=;\n\t"
134       "R3 += 1;\n\t"
135       "LC1 = R3;\n\t"
136 
137       "R1 = [P4++];\n\t"
138       "A1 = R1;\n\t"
139       "A0 = 0;\n\t"
140       "I0 = B0;\n\t"
141       "I2 = P2;\n\t"
142       "P2 += 4;\n\t"
143       "R4 = [I0++] || R5 = [I2--];\n\t"
144       "LOOP filter_start_inner%= LC1;\n\t"
145       "LOOP_BEGIN filter_start_inner%=;\n\t"
146          "A0 += R4.L*R5.L (IS), A1 -= R4.H*R5.H (IS) || R4 = [I0++] || R5 = [I2--];\n\t"
147       "LOOP_END filter_start_inner%=;\n\t"
148       "A0 += A1;\n\t"
149       "R4 = A0;\n\t"
150       "R4 <<= 1;\n\t"
151       "R2 = [P0++];\n\t"
152       "R4 = R4 + R2;\n\t"
153       "[P1++] = R4;\n\t"
154       "R4 <<= 2;\n\t"
155       "R2 <<= 2;\n\t"
156       "R2 = PACK(R4.H, R2.H);\n\t"
157       "[P2] = R2;\n\t"
158 
159    "LOOP_END filter_start%=;\n\t"
160 
161    /* Samples ord to N*/
162    "R0 = %5;\n\t"
163    "R0 <<= 1;\n\t"
164    "I0 = B0;\n\t"
165    "R0 <<= 1;\n\t"
166    "L0 = R0;\n\t"
167 
168    "R0 = %5;\n\t"
169    "R2 = %4;\n\t"
170    "R2 = R2 - R0;\n\t"
171    "R4 = [I0++];\n\t"
172    "LC0 = R2;\n\t"
173    "P3 = R0;\n\t"
174    "R0 <<= 2;\n\t"
175    "R0 += 8;\n\t"
176    "I2 = P2;\n\t"
177    "M0 = R0;\n\t"
178    "A0 = A1 = 0;\n\t"
179    "R5 = [I2--];\n\t"
180    "LOOP filter_mid%= LC0;\n\t"
181    "LOOP_BEGIN filter_mid%=;\n\t"
182       "LOOP filter_mid_inner%= LC1=P3;\n\t"
183       "LOOP_BEGIN filter_mid_inner%=;\n\t"
184          "A0 += R4.L*R5.L (IS), A1 -= R4.H*R5.H (IS) || R4 = [I0++] || R5 = [I2--];\n\t"
185       "LOOP_END filter_mid_inner%=;\n\t"
186       "R0 = (A0 += A1) || I2 += M0;\n\t"
187       "R0 = R0 << 1 || R5 = [P0++];\n\t"
188       "R0 = R0 + R5;\n\t"
189       "R0 = R0 << 2 || [P1++] = R0;\n\t"
190       "R5 = R5 << 2;\n\t"
191       "R5 = PACK(R0.H, R5.H);\n\t"
192       "A0 = A1 = 0 || [I2--] = R5\n\t"
193       "LOOP_END filter_mid%=;\n\t"
194    "I2 += 4;\n\t"
195    "P2 = I2;\n\t"
196    /* Update memory */
197    "P4 = %6;\n\t"
198    "R0 = %5;\n\t"
199    "LC0 = R0;\n\t"
200    "P0 = B0;\n\t"
201    "A0 = A1 = 0;\n\t"
202    "LOOP mem_update%= LC0;\n\t"
203    "LOOP_BEGIN mem_update%=;\n\t"
204       "I2 = P2;\n\t"
205       "I0 = P0;\n\t"
206       "P0 += 4;\n\t"
207       "R0 = LC0;\n\t"
208       "LC1 = R0;\n\t"
209       "R5 = [I2--] || R4 = [I0++];\n\t"
210       "LOOP mem_accum%= LC1;\n\t"
211       "LOOP_BEGIN mem_accum%=;\n\t"
212          "A0 += R4.L*R5.L (IS), A1 -= R4.H*R5.H (IS) || R4 = [I0++] || R5 = [I2--];\n\t"
213       "LOOP_END mem_accum%=;\n\t"
214       "R0 = (A0 += A1);\n\t"
215       "A0 = A1 = 0 || [P4++] = R0;\n\t"
216    "LOOP_END mem_update%=;\n\t"
217    "L0 = 0;\n\t"
218    : : "m" (xy), "m" (_x), "m" (_y), "m" (numden), "m" (N), "m" (ord), "m" (mem)
219    : "A0", "A1", "R0", "R1", "R2", "R3", "R4", "R5", "P0", "P1", "P2", "P3", "P4", "B0", "I0", "I2", "L0", "L2", "M0", "memory"
220    );
221 
222 }
223 
224 
225 
226 
227 #define OVERRIDE_IIR_MEM2
iir_mem2(const spx_sig_t * _x,const spx_coef_t * den,spx_sig_t * _y,int N,int ord,spx_mem_t * mem)228 void iir_mem2(const spx_sig_t *_x, const spx_coef_t *den, spx_sig_t *_y, int N, int ord, spx_mem_t *mem)
229 {
230    spx_word16_t y[N+2];
231    spx_word16_t *yy;
232    yy = y+2;
233    __asm__ __volatile__
234    (
235    /* Register setup */
236    "R0 = %5;\n\t"      /*ord */
237 
238    "P1 = %3;\n\t"
239    "I1 = P1;\n\t"
240    "B1 = P1;\n\t"
241    "L1 = 0;\n\t"
242 
243    "P3 = %0;\n\t"
244    "I3 = P3;\n\t"
245    "L3 = 0;\n\t"
246 
247    "P4 = %6;\n\t"
248    "P0 = %1;\n\t"
249    "P1 = %2;\n\t"
250 
251    /* First sample */
252    "R1 = [P4++];\n\t"
253    "R1 <<= 1;\n\t"
254    "R2 = [P0++];\n\t"
255    "R1 = R1 + R2;\n\t"
256    "[P1++] = R1;\n\t"
257    "R1 <<= 2;\n\t"
258    "W[P3] = R1.H;\n\t"
259    "R2 <<= 2;\n\t"
260 
261    /* Samples 1 to ord-1 (using memory) */
262    "R0 += -1;\n\t"
263    "R3 = 0;\n\t"
264    "LC0 = R0;\n\t"
265    "LOOP filter_start%= LC0;\n\t"
266    "LOOP_BEGIN filter_start%=;\n\t"
267       "R3 += 1;\n\t"
268       "LC1 = R3;\n\t"
269 
270       "R1 = [P4++];\n\t"
271       "A1 = R1;\n\t"
272       "I1 = B1;\n\t"
273       "I3 = P3;\n\t"
274       "P3 += 2;\n\t"
275       "LOOP filter_start_inner%= LC1;\n\t"
276       "LOOP_BEGIN filter_start_inner%=;\n\t"
277          "R4.L = W[I1++];\n\t"
278          "R5.L = W[I3--];\n\t"
279          "A1 -= R4.L*R5.L (IS);\n\t"
280       "LOOP_END filter_start_inner%=;\n\t"
281 
282       "R1 = A1;\n\t"
283       "R1 <<= 1;\n\t"
284       "R2 = [P0++];\n\t"
285       "R1 = R1 + R2;\n\t"
286       "[P1++] = R1;\n\t"
287       "R1 <<= 2;\n\t"
288       "W[P3] = R1.H;\n\t"
289       "R2 <<= 2;\n\t"
290    "LOOP_END filter_start%=;\n\t"
291 
292    /* Samples ord to N*/
293    "R0 = %5;\n\t"
294    "R0 <<= 1;\n\t"
295    "I1 = B1;\n\t"
296    "L1 = R0;\n\t"
297 
298    "R0 = %5;\n\t"
299    "R2 = %4;\n\t"
300    "R2 = R2 - R0;\n\t"
301    "R4.L = W[I1++];\n\t"
302    "LC0 = R2;\n\t"
303    "LOOP filter_mid%= LC0;\n\t"
304    "LOOP_BEGIN filter_mid%=;\n\t"
305       "LC1 = R0;\n\t"
306       "A1 = 0;\n\t"
307       "I3 = P3;\n\t"
308       "P3 += 2;\n\t"
309       "R5.L = W[I3--];\n\t"
310       "LOOP filter_mid_inner%= LC1;\n\t"
311       "LOOP_BEGIN filter_mid_inner%=;\n\t"
312          "A1 -= R4.L*R5.L (IS) || R4.L = W[I1++] || R5.L = W[I3--];\n\t"
313       "LOOP_END filter_mid_inner%=;\n\t"
314       "R1 = A1;\n\t"
315       "R1 = R1 << 1 || R2 = [P0++];\n\t"
316       "R1 = R1 + R2;\n\t"
317       "R1 = R1 << 2 || [P1++] = R1;\n\t"
318       "W[P3] = R1.H;\n\t"
319    "LOOP_END filter_mid%=;\n\t"
320 
321    /* Update memory */
322    "P4 = %6;\n\t"
323    "R0 = %5;\n\t"
324    "LC0 = R0;\n\t"
325    "P1 = B1;\n\t"
326    "LOOP mem_update%= LC0;\n\t"
327    "LOOP_BEGIN mem_update%=;\n\t"
328       "A0 = 0;\n\t"
329       "I3 = P3;\n\t"
330       "I1 = P1;\n\t"
331       "P1 += 2;\n\t"
332       "R0 = LC0;\n\t"
333       "LC1=R0;\n\t"
334       "R5.L = W[I3--] || R4.L = W[I1++];\n\t"
335       "LOOP mem_accum%= LC1;\n\t"
336       "LOOP_BEGIN mem_accum%=;\n\t"
337          "A0 -= R4.L*R5.L (IS) || R4.L = W[I1++] || R5.L = W[I3--];\n\t"
338       "LOOP_END mem_accum%=;\n\t"
339       "R0 = A0;\n\t"
340       "[P4++] = R0;\n\t"
341    "LOOP_END mem_update%=;\n\t"
342    "L1 = 0;\n\t"
343    : : "m" (yy), "m" (_x), "m" (_y), "m" (den), "m" (N), "m" (ord), "m" (mem)
344    : "A0", "A1", "R0", "R1", "R2", "R3", "R4", "R5", "P0", "P1", "P2", "P3", "P4", "B1", "I1", "I3", "L1", "L3", "memory"
345    );
346 
347 }
348 
349 #define OVERRIDE_FIR_MEM2
fir_mem2(const spx_sig_t * x,const spx_coef_t * num,spx_sig_t * y,int N,int ord,spx_mem_t * mem)350 void fir_mem2(const spx_sig_t *x, const spx_coef_t *num, spx_sig_t *y, int N, int ord, spx_mem_t *mem)
351 {
352    int i;
353    spx_coef_t den2[12];
354    spx_coef_t *den;
355    den = (spx_coef_t*)((((int)den2)+4)&0xfffffffc);
356    for (i=0;i<10;i++)
357       den[i] = 0;
358    filter_mem2(x, num, den, y, N, ord, mem);
359 }
360 
361 
362 #define OVERRIDE_COMPUTE_IMPULSE_RESPONSE
compute_impulse_response(const spx_coef_t * ak,const spx_coef_t * awk1,const spx_coef_t * awk2,spx_word16_t * y,int N,int ord,char * stack)363 void compute_impulse_response(const spx_coef_t *ak, const spx_coef_t *awk1, const spx_coef_t *awk2, spx_word16_t *y, int N, int ord, char *stack)
364 {
365    int i;
366    VARDECL(spx_word16_t *ytmp);
367    ALLOC(ytmp, N, spx_word16_t);
368    spx_word16_t *ytmp2 = ytmp;
369    y[0] = LPC_SCALING;
370    for (i=0;i<ord;i++)
371       y[i+1] = awk1[i];
372    i++;
373    for (;i<N;i++)
374       y[i] = 0;
375 
376    N-=1;
377    __asm__ __volatile__
378    (
379          "I0 = %0;\n\t"
380          "I1 = %1;\n\t"
381          "L0 = 0;\n\t"
382          "L1 = 0;\n\t"
383          "L2 = 0;\n\t"
384          "L3 = 0;\n\t"
385          "R0 = 1;\n\t"
386          "R0 <<= 13;\n\t"
387          "W[I0] = R0.L;\n\t"
388          "R0 <<= 1;\n\t"
389          "W[I1] = R0.L;\n\t"
390          "R0 = %5;\n\t"
391          "LC0 = R0;\n\t"
392          "R2 = 0;\n\t"
393          "LOOP samples%= LC0;\n\t"
394          "LOOP_BEGIN samples%=;\n\t"
395             "R2 += 1;\n\t"
396             "R2 = MIN(R2, %4);\n\t"
397             "I0 = %0;\n\t"
398             "I1 = %1;\n\t"
399             "I2 = %2;\n\t"
400             "I3 = %3;\n\t"
401             "%0 += 2;\n\t"
402             "%1 += 2;\n\t"
403             "A0 = A1 = 0;\n\t"
404             "R0.L = W[I0--] || R1.L = W[I2++];\n\t"
405             "LC1 = R2;\n\t"
406             "LOOP filter%= LC1;\n\t"
407             "LOOP_BEGIN filter%=;\n\t"
408                "A0 -= R0.L*R1.L (IS) || R0.L = W[I1--] || R1.L = W[I3++];\n\t"
409                "A1 -= R0.L*R1.L (IS) || R0.L = W[I0--] || R1.L = W[I2++];\n\t"
410             "LOOP_END filter%=;\n\t"
411             "R0 = A0, R1 = A1;\n\t"
412             "R3 = W[%1] (X);\n\t"
413             "R3 <<= 13;\n\t"
414             "R0 = R0 + R3;\n\t"
415             "R3 = R0 >>> 13;\n\t"
416             "W[%0] = R3.L;\n\t"
417             "R0 <<= 1;\n\t"
418             "R1 = R1 + R0;\n\t"
419             "R1 >>>= 13;\n\t"
420             "W[%1] = R1.L;\n\t"
421          "LOOP_END samples%=;\n\t"
422    : "=a" (ytmp2), "=a" (y)
423    : "a" (awk2), "a" (ak), "d" (ord), "m" (N), "0" (ytmp2), "1" (y)
424    : "A0", "A1", "R0", "R1", "R2", "R3", "I0", "I1", "I2", "I3", "L0", "L1", "L2", "L3", "A0", "A1"
425    );
426 }
427 
428 
429 
430 #if 0 /* Equivalent C function for filter_mem2 and compute_impulse_response */
431 #define min(a,b) ((a)<(b) ? (a):(b))
432 
433 void compute_impulse_response(const spx_coef_t *ak, const spx_coef_t *awk1, const spx_coef_t *awk2, spx_word16_t *y, int N, int ord, char *stack)
434 {
435    int i,j;
436    VARDECL(spx_word16_t *ytmp);
437    ALLOC(ytmp, N, spx_word16_t);
438 
439    y[0] = LPC_SCALING;
440    for (i=0;i<ord;i++)
441       y[i+1] = awk1[i];
442    i++;
443    for (;i<N;i++)
444       y[i] = 0;
445 
446    for (i=0;i<N;i++)
447    {
448       spx_word32_t yi = SHL32(EXTEND32(y[i]),LPC_SHIFT);
449       spx_word32_t yi2 = 0;
450       for (j=0;j<min(i,ord);j++)
451       {
452          yi = MAC16_16(yi, awk2[j], -ytmp[i-j-1]);
453          yi2 = MAC16_16(yi2, ak[j], -y[i-j-1]);
454       }
455       ytmp[i] = EXTRACT16(SHR32(yi,LPC_SHIFT));
456       yi2 = ADD32(yi2,SHL32(yi,1));
457       y[i] = EXTRACT16(SHR32(yi2,LPC_SHIFT));
458    }
459 
460 }
461 
462 
463 void filter_mem2(const spx_sig_t *_x, const spx_coef_t *num, const spx_coef_t *den, spx_sig_t *_y, int N, int ord, spx_mem_t *mem)
464 {
465    int i,j;
466    spx_word16_t xi,yi,nyi;
467    spx_word16_t x[N],y[N];
468    spx_word16_t *xx, *yy;
469    xx = x;
470    yy = y;
471    for (i=0;i<N;i++)
472    {
473       x[i] = EXTRACT16(SHR32(_x[i],SIG_SHIFT));
474    }
475 
476    for (i=0;i<ord;i++)
477    {
478       spx_word32_t yi = mem[i];
479       for (j=0;j<i;j++)
480       {
481          yi = MAC16_16(yi, num[j], x[i-j-1]);
482          yi = MAC16_16(yi, den[j], -y[i-j-1]);
483       }
484       _y[i] = ADD32(_x[i],SHL32(yi,1));
485       y[i] = EXTRACT16(SHR32(_y[i],SIG_SHIFT));
486    }
487    for (i=ord;i<N;i++)
488    {
489       spx_word32_t yi = 0;
490       for (j=0;j<ord;j++)
491       {
492          yi = MAC16_16(yi, num[j], x[i-j-1]);
493          yi = MAC16_16(yi, den[j], -y[i-j-1]);
494       }
495       _y[i] = ADD32(_x[i],SHL32(yi,1));
496       y[i] = EXTRACT16(SHR32(_y[i],SIG_SHIFT));
497    }
498 
499    for (i=0;i<ord;i++)
500    {
501       spx_mem_t m = 0;
502       for (j=0;j<ord-i;j++)
503       {
504          m = MAC16_16(m, x[N-1-j], num[j+i]);
505          m = MAC16_16(m, -y[N-1-j], den[j+i]);
506       }
507       mem[i] = m;
508    }
509 }
510 #endif
511