1 /* Copyright (c) 2014, Cisco Systems, INC
2 Written by XiangMingZhu WeiZhou MinPeng YanWang
3
4 Redistribution and use in source and binary forms, with or without
5 modification, are permitted provided that the following conditions
6 are met:
7
8 - Redistributions of source code must retain the above copyright
9 notice, this list of conditions and the following disclaimer.
10
11 - Redistributions in binary form must reproduce the above copyright
12 notice, this list of conditions and the following disclaimer in the
13 documentation and/or other materials provided with the distribution.
14
15 THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
16 ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
17 LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
18 A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER
19 OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
20 EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
21 PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
22 PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
23 LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
24 NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
25 SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
26 */
27
28 #ifdef HAVE_CONFIG_H
29 #include "config.h"
30 #endif
31
32 #include <xmmintrin.h>
33 #include <emmintrin.h>
34 #include <smmintrin.h>
35
36 #include "main.h"
37 #include "stack_alloc.h"
38
39 /* Weighting factors for tilt measure */
40 static const opus_int32 tiltWeights[ VAD_N_BANDS ] = { 30000, 6000, -12000, -12000 };
41
42 /***************************************/
43 /* Get the speech activity level in Q8 */
44 /***************************************/
silk_VAD_GetSA_Q8_sse4_1(silk_encoder_state * psEncC,const opus_int16 pIn[])45 opus_int silk_VAD_GetSA_Q8_sse4_1( /* O Return value, 0 if success */
46 silk_encoder_state *psEncC, /* I/O Encoder state */
47 const opus_int16 pIn[] /* I PCM input */
48 )
49 {
50 opus_int SA_Q15, pSNR_dB_Q7, input_tilt;
51 opus_int decimated_framelength1, decimated_framelength2;
52 opus_int decimated_framelength;
53 opus_int dec_subframe_length, dec_subframe_offset, SNR_Q7, i, b, s;
54 opus_int32 sumSquared, smooth_coef_Q16;
55 opus_int16 HPstateTmp;
56 VARDECL( opus_int16, X );
57 opus_int32 Xnrg[ VAD_N_BANDS ];
58 opus_int32 NrgToNoiseRatio_Q8[ VAD_N_BANDS ];
59 opus_int32 speech_nrg, x_tmp;
60 opus_int X_offset[ VAD_N_BANDS ];
61 opus_int ret = 0;
62 silk_VAD_state *psSilk_VAD = &psEncC->sVAD;
63
64 SAVE_STACK;
65
66 /* Safety checks */
67 silk_assert( VAD_N_BANDS == 4 );
68 celt_assert( MAX_FRAME_LENGTH >= psEncC->frame_length );
69 celt_assert( psEncC->frame_length <= 512 );
70 celt_assert( psEncC->frame_length == 8 * silk_RSHIFT( psEncC->frame_length, 3 ) );
71
72 /***********************/
73 /* Filter and Decimate */
74 /***********************/
75 decimated_framelength1 = silk_RSHIFT( psEncC->frame_length, 1 );
76 decimated_framelength2 = silk_RSHIFT( psEncC->frame_length, 2 );
77 decimated_framelength = silk_RSHIFT( psEncC->frame_length, 3 );
78 /* Decimate into 4 bands:
79 0 L 3L L 3L 5L
80 - -- - -- --
81 8 8 2 4 4
82
83 [0-1 kHz| temp. |1-2 kHz| 2-4 kHz | 4-8 kHz |
84
85 They're arranged to allow the minimal ( frame_length / 4 ) extra
86 scratch space during the downsampling process */
87 X_offset[ 0 ] = 0;
88 X_offset[ 1 ] = decimated_framelength + decimated_framelength2;
89 X_offset[ 2 ] = X_offset[ 1 ] + decimated_framelength;
90 X_offset[ 3 ] = X_offset[ 2 ] + decimated_framelength2;
91 ALLOC( X, X_offset[ 3 ] + decimated_framelength1, opus_int16 );
92
93 /* 0-8 kHz to 0-4 kHz and 4-8 kHz */
94 silk_ana_filt_bank_1( pIn, &psSilk_VAD->AnaState[ 0 ],
95 X, &X[ X_offset[ 3 ] ], psEncC->frame_length );
96
97 /* 0-4 kHz to 0-2 kHz and 2-4 kHz */
98 silk_ana_filt_bank_1( X, &psSilk_VAD->AnaState1[ 0 ],
99 X, &X[ X_offset[ 2 ] ], decimated_framelength1 );
100
101 /* 0-2 kHz to 0-1 kHz and 1-2 kHz */
102 silk_ana_filt_bank_1( X, &psSilk_VAD->AnaState2[ 0 ],
103 X, &X[ X_offset[ 1 ] ], decimated_framelength2 );
104
105 /*********************************************/
106 /* HP filter on lowest band (differentiator) */
107 /*********************************************/
108 X[ decimated_framelength - 1 ] = silk_RSHIFT( X[ decimated_framelength - 1 ], 1 );
109 HPstateTmp = X[ decimated_framelength - 1 ];
110 for( i = decimated_framelength - 1; i > 0; i-- ) {
111 X[ i - 1 ] = silk_RSHIFT( X[ i - 1 ], 1 );
112 X[ i ] -= X[ i - 1 ];
113 }
114 X[ 0 ] -= psSilk_VAD->HPstate;
115 psSilk_VAD->HPstate = HPstateTmp;
116
117 /*************************************/
118 /* Calculate the energy in each band */
119 /*************************************/
120 for( b = 0; b < VAD_N_BANDS; b++ ) {
121 /* Find the decimated framelength in the non-uniformly divided bands */
122 decimated_framelength = silk_RSHIFT( psEncC->frame_length, silk_min_int( VAD_N_BANDS - b, VAD_N_BANDS - 1 ) );
123
124 /* Split length into subframe lengths */
125 dec_subframe_length = silk_RSHIFT( decimated_framelength, VAD_INTERNAL_SUBFRAMES_LOG2 );
126 dec_subframe_offset = 0;
127
128 /* Compute energy per sub-frame */
129 /* initialize with summed energy of last subframe */
130 Xnrg[ b ] = psSilk_VAD->XnrgSubfr[ b ];
131 for( s = 0; s < VAD_INTERNAL_SUBFRAMES; s++ ) {
132 __m128i xmm_X, xmm_acc;
133 sumSquared = 0;
134
135 xmm_acc = _mm_setzero_si128();
136
137 for( i = 0; i < dec_subframe_length - 7; i += 8 )
138 {
139 xmm_X = _mm_loadu_si128( (__m128i *)&(X[ X_offset[ b ] + i + dec_subframe_offset ] ) );
140 xmm_X = _mm_srai_epi16( xmm_X, 3 );
141 xmm_X = _mm_madd_epi16( xmm_X, xmm_X );
142 xmm_acc = _mm_add_epi32( xmm_acc, xmm_X );
143 }
144
145 xmm_acc = _mm_add_epi32( xmm_acc, _mm_unpackhi_epi64( xmm_acc, xmm_acc ) );
146 xmm_acc = _mm_add_epi32( xmm_acc, _mm_shufflelo_epi16( xmm_acc, 0x0E ) );
147
148 sumSquared += _mm_cvtsi128_si32( xmm_acc );
149
150 for( ; i < dec_subframe_length; i++ ) {
151 /* The energy will be less than dec_subframe_length * ( silk_int16_MIN / 8 ) ^ 2. */
152 /* Therefore we can accumulate with no risk of overflow (unless dec_subframe_length > 128) */
153 x_tmp = silk_RSHIFT(
154 X[ X_offset[ b ] + i + dec_subframe_offset ], 3 );
155 sumSquared = silk_SMLABB( sumSquared, x_tmp, x_tmp );
156
157 /* Safety check */
158 silk_assert( sumSquared >= 0 );
159 }
160
161 /* Add/saturate summed energy of current subframe */
162 if( s < VAD_INTERNAL_SUBFRAMES - 1 ) {
163 Xnrg[ b ] = silk_ADD_POS_SAT32( Xnrg[ b ], sumSquared );
164 } else {
165 /* Look-ahead subframe */
166 Xnrg[ b ] = silk_ADD_POS_SAT32( Xnrg[ b ], silk_RSHIFT( sumSquared, 1 ) );
167 }
168
169 dec_subframe_offset += dec_subframe_length;
170 }
171 psSilk_VAD->XnrgSubfr[ b ] = sumSquared;
172 }
173
174 /********************/
175 /* Noise estimation */
176 /********************/
177 silk_VAD_GetNoiseLevels( &Xnrg[ 0 ], psSilk_VAD );
178
179 /***********************************************/
180 /* Signal-plus-noise to noise ratio estimation */
181 /***********************************************/
182 sumSquared = 0;
183 input_tilt = 0;
184 for( b = 0; b < VAD_N_BANDS; b++ ) {
185 speech_nrg = Xnrg[ b ] - psSilk_VAD->NL[ b ];
186 if( speech_nrg > 0 ) {
187 /* Divide, with sufficient resolution */
188 if( ( Xnrg[ b ] & 0xFF800000 ) == 0 ) {
189 NrgToNoiseRatio_Q8[ b ] = silk_DIV32( silk_LSHIFT( Xnrg[ b ], 8 ), psSilk_VAD->NL[ b ] + 1 );
190 } else {
191 NrgToNoiseRatio_Q8[ b ] = silk_DIV32( Xnrg[ b ], silk_RSHIFT( psSilk_VAD->NL[ b ], 8 ) + 1 );
192 }
193
194 /* Convert to log domain */
195 SNR_Q7 = silk_lin2log( NrgToNoiseRatio_Q8[ b ] ) - 8 * 128;
196
197 /* Sum-of-squares */
198 sumSquared = silk_SMLABB( sumSquared, SNR_Q7, SNR_Q7 ); /* Q14 */
199
200 /* Tilt measure */
201 if( speech_nrg < ( (opus_int32)1 << 20 ) ) {
202 /* Scale down SNR value for small subband speech energies */
203 SNR_Q7 = silk_SMULWB( silk_LSHIFT( silk_SQRT_APPROX( speech_nrg ), 6 ), SNR_Q7 );
204 }
205 input_tilt = silk_SMLAWB( input_tilt, tiltWeights[ b ], SNR_Q7 );
206 } else {
207 NrgToNoiseRatio_Q8[ b ] = 256;
208 }
209 }
210
211 /* Mean-of-squares */
212 sumSquared = silk_DIV32_16( sumSquared, VAD_N_BANDS ); /* Q14 */
213
214 /* Root-mean-square approximation, scale to dBs, and write to output pointer */
215 pSNR_dB_Q7 = (opus_int16)( 3 * silk_SQRT_APPROX( sumSquared ) ); /* Q7 */
216
217 /*********************************/
218 /* Speech Probability Estimation */
219 /*********************************/
220 SA_Q15 = silk_sigm_Q15( silk_SMULWB( VAD_SNR_FACTOR_Q16, pSNR_dB_Q7 ) - VAD_NEGATIVE_OFFSET_Q5 );
221
222 /**************************/
223 /* Frequency Tilt Measure */
224 /**************************/
225 psEncC->input_tilt_Q15 = silk_LSHIFT( silk_sigm_Q15( input_tilt ) - 16384, 1 );
226
227 /**************************************************/
228 /* Scale the sigmoid output based on power levels */
229 /**************************************************/
230 speech_nrg = 0;
231 for( b = 0; b < VAD_N_BANDS; b++ ) {
232 /* Accumulate signal-without-noise energies, higher frequency bands have more weight */
233 speech_nrg += ( b + 1 ) * silk_RSHIFT( Xnrg[ b ] - psSilk_VAD->NL[ b ], 4 );
234 }
235
236 /* Power scaling */
237 if( speech_nrg <= 0 ) {
238 SA_Q15 = silk_RSHIFT( SA_Q15, 1 );
239 } else if( speech_nrg < 32768 ) {
240 if( psEncC->frame_length == 10 * psEncC->fs_kHz ) {
241 speech_nrg = silk_LSHIFT_SAT32( speech_nrg, 16 );
242 } else {
243 speech_nrg = silk_LSHIFT_SAT32( speech_nrg, 15 );
244 }
245
246 /* square-root */
247 speech_nrg = silk_SQRT_APPROX( speech_nrg );
248 SA_Q15 = silk_SMULWB( 32768 + speech_nrg, SA_Q15 );
249 }
250
251 /* Copy the resulting speech activity in Q8 */
252 psEncC->speech_activity_Q8 = silk_min_int( silk_RSHIFT( SA_Q15, 7 ), silk_uint8_MAX );
253
254 /***********************************/
255 /* Energy Level and SNR estimation */
256 /***********************************/
257 /* Smoothing coefficient */
258 smooth_coef_Q16 = silk_SMULWB( VAD_SNR_SMOOTH_COEF_Q18, silk_SMULWB( (opus_int32)SA_Q15, SA_Q15 ) );
259
260 if( psEncC->frame_length == 10 * psEncC->fs_kHz ) {
261 smooth_coef_Q16 >>= 1;
262 }
263
264 for( b = 0; b < VAD_N_BANDS; b++ ) {
265 /* compute smoothed energy-to-noise ratio per band */
266 psSilk_VAD->NrgRatioSmth_Q8[ b ] = silk_SMLAWB( psSilk_VAD->NrgRatioSmth_Q8[ b ],
267 NrgToNoiseRatio_Q8[ b ] - psSilk_VAD->NrgRatioSmth_Q8[ b ], smooth_coef_Q16 );
268
269 /* signal to noise ratio in dB per band */
270 SNR_Q7 = 3 * ( silk_lin2log( psSilk_VAD->NrgRatioSmth_Q8[b] ) - 8 * 128 );
271 /* quality = sigmoid( 0.25 * ( SNR_dB - 16 ) ); */
272 psEncC->input_quality_bands_Q15[ b ] = silk_sigm_Q15( silk_RSHIFT( SNR_Q7 - 16 * 128, 4 ) );
273 }
274
275 RESTORE_STACK;
276 return( ret );
277 }
278