|
1 /*********************************************************************** |
|
2 Copyright (c) 2006-2011, Skype Limited. All rights reserved. |
|
3 Redistribution and use in source and binary forms, with or without |
|
4 modification, are permitted provided that the following conditions |
|
5 are met: |
|
6 - Redistributions of source code must retain the above copyright notice, |
|
7 this list of conditions and the following disclaimer. |
|
8 - Redistributions in binary form must reproduce the above copyright |
|
9 notice, this list of conditions and the following disclaimer in the |
|
10 documentation and/or other materials provided with the distribution. |
|
11 - Neither the name of Internet Society, IETF or IETF Trust, nor the |
|
12 names of specific contributors, may be used to endorse or promote |
|
13 products derived from this software without specific prior written |
|
14 permission. |
|
15 THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" |
|
16 AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE |
|
17 IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE |
|
18 ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE |
|
19 LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR |
|
20 CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF |
|
21 SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS |
|
22 INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN |
|
23 CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) |
|
24 ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE |
|
25 POSSIBILITY OF SUCH DAMAGE. |
|
26 ***********************************************************************/ |
|
27 |
|
28 #ifdef HAVE_CONFIG_H |
|
29 #include "config.h" |
|
30 #endif |
|
31 |
|
32 #include "main_FIX.h" |
|
33 |
|
34 /* Residual energy: nrg = wxx - 2 * wXx * c + c' * wXX * c */ |
|
35 opus_int32 silk_residual_energy16_covar_FIX( |
|
36 const opus_int16 *c, /* I Prediction vector */ |
|
37 const opus_int32 *wXX, /* I Correlation matrix */ |
|
38 const opus_int32 *wXx, /* I Correlation vector */ |
|
39 opus_int32 wxx, /* I Signal energy */ |
|
40 opus_int D, /* I Dimension */ |
|
41 opus_int cQ /* I Q value for c vector 0 - 15 */ |
|
42 ) |
|
43 { |
|
44 opus_int i, j, lshifts, Qxtra; |
|
45 opus_int32 c_max, w_max, tmp, tmp2, nrg; |
|
46 opus_int cn[ MAX_MATRIX_SIZE ]; |
|
47 const opus_int32 *pRow; |
|
48 |
|
49 /* Safety checks */ |
|
50 silk_assert( D >= 0 ); |
|
51 silk_assert( D <= 16 ); |
|
52 silk_assert( cQ > 0 ); |
|
53 silk_assert( cQ < 16 ); |
|
54 |
|
55 lshifts = 16 - cQ; |
|
56 Qxtra = lshifts; |
|
57 |
|
58 c_max = 0; |
|
59 for( i = 0; i < D; i++ ) { |
|
60 c_max = silk_max_32( c_max, silk_abs( (opus_int32)c[ i ] ) ); |
|
61 } |
|
62 Qxtra = silk_min_int( Qxtra, silk_CLZ32( c_max ) - 17 ); |
|
63 |
|
64 w_max = silk_max_32( wXX[ 0 ], wXX[ D * D - 1 ] ); |
|
65 Qxtra = silk_min_int( Qxtra, silk_CLZ32( silk_MUL( D, silk_RSHIFT( silk_SMULWB( w_max, c_max ), 4 ) ) ) - 5 ); |
|
66 Qxtra = silk_max_int( Qxtra, 0 ); |
|
67 for( i = 0; i < D; i++ ) { |
|
68 cn[ i ] = silk_LSHIFT( ( opus_int )c[ i ], Qxtra ); |
|
69 silk_assert( silk_abs(cn[i]) <= ( silk_int16_MAX + 1 ) ); /* Check that silk_SMLAWB can be used */ |
|
70 } |
|
71 lshifts -= Qxtra; |
|
72 |
|
73 /* Compute wxx - 2 * wXx * c */ |
|
74 tmp = 0; |
|
75 for( i = 0; i < D; i++ ) { |
|
76 tmp = silk_SMLAWB( tmp, wXx[ i ], cn[ i ] ); |
|
77 } |
|
78 nrg = silk_RSHIFT( wxx, 1 + lshifts ) - tmp; /* Q: -lshifts - 1 */ |
|
79 |
|
80 /* Add c' * wXX * c, assuming wXX is symmetric */ |
|
81 tmp2 = 0; |
|
82 for( i = 0; i < D; i++ ) { |
|
83 tmp = 0; |
|
84 pRow = &wXX[ i * D ]; |
|
85 for( j = i + 1; j < D; j++ ) { |
|
86 tmp = silk_SMLAWB( tmp, pRow[ j ], cn[ j ] ); |
|
87 } |
|
88 tmp = silk_SMLAWB( tmp, silk_RSHIFT( pRow[ i ], 1 ), cn[ i ] ); |
|
89 tmp2 = silk_SMLAWB( tmp2, tmp, cn[ i ] ); |
|
90 } |
|
91 nrg = silk_ADD_LSHIFT32( nrg, tmp2, lshifts ); /* Q: -lshifts - 1 */ |
|
92 |
|
93 /* Keep one bit free always, because we add them for LSF interpolation */ |
|
94 if( nrg < 1 ) { |
|
95 nrg = 1; |
|
96 } else if( nrg > silk_RSHIFT( silk_int32_MAX, lshifts + 2 ) ) { |
|
97 nrg = silk_int32_MAX >> 1; |
|
98 } else { |
|
99 nrg = silk_LSHIFT( nrg, lshifts + 1 ); /* Q0 */ |
|
100 } |
|
101 return nrg; |
|
102 |
|
103 } |