Convert all CRLF in the SILK code, tabs to spaces, and trailing
[opus.git] / silk / silk_VAD.c
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, (subject to the limitations in the disclaimer below)
5 are permitted provided that the following conditions 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 Skype Limited, nor the names of specific
12 contributors, may be used to endorse or promote products derived from
13 this software without specific prior written permission.
14 NO EXPRESS OR IMPLIED LICENSES TO ANY PARTY'S PATENT RIGHTS ARE GRANTED
15 BY THIS LICENSE. THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND
16 CONTRIBUTORS ''AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING,
17 BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND
18 FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE
19 COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
20 INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT
21 NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF
22 USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON
23 ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
24 (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
25 OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
26 ***********************************************************************/
27
28 #include <stdlib.h>
29 #include "silk_main.h"
30
31 /**********************************/
32 /* Initialization of the Silk VAD */
33 /**********************************/
34 opus_int silk_VAD_Init(                              /* O    Return value, 0 if success                  */
35     silk_VAD_state              *psSilk_VAD         /* I/O  Pointer to Silk VAD state                   */
36 )
37 {
38     opus_int b, ret = 0;
39
40     /* reset state memory */
41     SKP_memset( psSilk_VAD, 0, sizeof( silk_VAD_state ) );
42
43     /* init noise levels */
44     /* Initialize array with approx pink noise levels (psd proportional to inverse of frequency) */
45     for( b = 0; b < VAD_N_BANDS; b++ ) {
46         psSilk_VAD->NoiseLevelBias[ b ] = SKP_max_32( SKP_DIV32_16( VAD_NOISE_LEVELS_BIAS, b + 1 ), 1 );
47     }
48
49     /* Initialize state */
50     for( b = 0; b < VAD_N_BANDS; b++ ) {
51         psSilk_VAD->NL[ b ]     = SKP_MUL( 100, psSilk_VAD->NoiseLevelBias[ b ] );
52         psSilk_VAD->inv_NL[ b ] = SKP_DIV32( SKP_int32_MAX, psSilk_VAD->NL[ b ] );
53     }
54     psSilk_VAD->counter = 15;
55
56     /* init smoothed energy-to-noise ratio*/
57     for( b = 0; b < VAD_N_BANDS; b++ ) {
58         psSilk_VAD->NrgRatioSmth_Q8[ b ] = 100 * 256;       /* 100 * 256 --> 20 dB SNR */
59     }
60
61     return( ret );
62 }
63
64 /* Weighting factors for tilt measure */
65 static const opus_int32 tiltWeights[ VAD_N_BANDS ] = { 30000, 6000, -12000, -12000 };
66
67 /***************************************/
68 /* Get the speech activity level in Q8 */
69 /***************************************/
70 opus_int silk_VAD_GetSA_Q8(                          /* O    Return value, 0 if success                  */
71     silk_encoder_state          *psEncC,            /* I/O  Encoder state                               */
72     const opus_int16             pIn[]               /* I    PCM input                                   */
73 )
74 {
75     opus_int   SA_Q15, pSNR_dB_Q7, input_tilt;
76     opus_int   decimated_framelength, dec_subframe_length, dec_subframe_offset, SNR_Q7, i, b, s;
77     opus_int32 sumSquared, smooth_coef_Q16;
78     opus_int16 HPstateTmp;
79     opus_int16 X[ VAD_N_BANDS ][ MAX_FRAME_LENGTH / 2 ];
80     opus_int32 Xnrg[ VAD_N_BANDS ];
81     opus_int32 NrgToNoiseRatio_Q8[ VAD_N_BANDS ];
82     opus_int32 speech_nrg, x_tmp;
83     opus_int   ret = 0;
84     silk_VAD_state *psSilk_VAD = &psEncC->sVAD;
85
86     /* Safety checks */
87     SKP_assert( VAD_N_BANDS == 4 );
88     SKP_assert( MAX_FRAME_LENGTH >= psEncC->frame_length );
89     SKP_assert( psEncC->frame_length <= 512 );
90     SKP_assert( psEncC->frame_length == 8 * SKP_RSHIFT( psEncC->frame_length, 3 ) );
91
92     /***********************/
93     /* Filter and Decimate */
94     /***********************/
95     /* 0-8 kHz to 0-4 kHz and 4-8 kHz */
96     silk_ana_filt_bank_1( pIn,          &psSilk_VAD->AnaState[  0 ], &X[ 0 ][ 0 ], &X[ 3 ][ 0 ], psEncC->frame_length );
97
98     /* 0-4 kHz to 0-2 kHz and 2-4 kHz */
99     silk_ana_filt_bank_1( &X[ 0 ][ 0 ], &psSilk_VAD->AnaState1[ 0 ], &X[ 0 ][ 0 ], &X[ 2 ][ 0 ], SKP_RSHIFT( psEncC->frame_length, 1 ) );
100
101     /* 0-2 kHz to 0-1 kHz and 1-2 kHz */
102     silk_ana_filt_bank_1( &X[ 0 ][ 0 ], &psSilk_VAD->AnaState2[ 0 ], &X[ 0 ][ 0 ], &X[ 1 ][ 0 ], SKP_RSHIFT( psEncC->frame_length, 2 ) );
103
104     /*********************************************/
105     /* HP filter on lowest band (differentiator) */
106     /*********************************************/
107     decimated_framelength = SKP_RSHIFT( psEncC->frame_length, 3 );
108     X[ 0 ][ decimated_framelength - 1 ] = SKP_RSHIFT( X[ 0 ][ decimated_framelength - 1 ], 1 );
109     HPstateTmp = X[ 0 ][ decimated_framelength - 1 ];
110     for( i = decimated_framelength - 1; i > 0; i-- ) {
111         X[ 0 ][ i - 1 ]  = SKP_RSHIFT( X[ 0 ][ i - 1 ], 1 );
112         X[ 0 ][ i ]     -= X[ 0 ][ i - 1 ];
113     }
114     X[ 0 ][ 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 = SKP_RSHIFT( psEncC->frame_length, SKP_min_int( VAD_N_BANDS - b, VAD_N_BANDS - 1 ) );
123
124         /* Split length into subframe lengths */
125         dec_subframe_length = SKP_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             sumSquared = 0;
133             for( i = 0; i < dec_subframe_length; i++ ) {
134                 /* The energy will be less than dec_subframe_length * ( SKP_int16_MIN / 8 ) ^ 2.            */
135                 /* Therefore we can accumulate with no risk of overflow (unless dec_subframe_length > 128)  */
136                 x_tmp = SKP_RSHIFT( X[ b ][ i + dec_subframe_offset ], 3 );
137                 sumSquared = SKP_SMLABB( sumSquared, x_tmp, x_tmp );
138
139                 /* Safety check */
140                 SKP_assert( sumSquared >= 0 );
141             }
142
143             /* Add/saturate summed energy of current subframe */
144             if( s < VAD_INTERNAL_SUBFRAMES - 1 ) {
145                 Xnrg[ b ] = SKP_ADD_POS_SAT32( Xnrg[ b ], sumSquared );
146             } else {
147                 /* Look-ahead subframe */
148                 Xnrg[ b ] = SKP_ADD_POS_SAT32( Xnrg[ b ], SKP_RSHIFT( sumSquared, 1 ) );
149             }
150
151             dec_subframe_offset += dec_subframe_length;
152         }
153         psSilk_VAD->XnrgSubfr[ b ] = sumSquared;
154     }
155
156     /********************/
157     /* Noise estimation */
158     /********************/
159     silk_VAD_GetNoiseLevels( &Xnrg[ 0 ], psSilk_VAD );
160
161     /***********************************************/
162     /* Signal-plus-noise to noise ratio estimation */
163     /***********************************************/
164     sumSquared = 0;
165     input_tilt = 0;
166     for( b = 0; b < VAD_N_BANDS; b++ ) {
167         speech_nrg = Xnrg[ b ] - psSilk_VAD->NL[ b ];
168         if( speech_nrg > 0 ) {
169             /* Divide, with sufficient resolution */
170             if( ( Xnrg[ b ] & 0xFF800000 ) == 0 ) {
171                 NrgToNoiseRatio_Q8[ b ] = SKP_DIV32( SKP_LSHIFT( Xnrg[ b ], 8 ), psSilk_VAD->NL[ b ] + 1 );
172             } else {
173                 NrgToNoiseRatio_Q8[ b ] = SKP_DIV32( Xnrg[ b ], SKP_RSHIFT( psSilk_VAD->NL[ b ], 8 ) + 1 );
174             }
175
176             /* Convert to log domain */
177             SNR_Q7 = silk_lin2log( NrgToNoiseRatio_Q8[ b ] ) - 8 * 128;
178
179             /* Sum-of-squares */
180             sumSquared = SKP_SMLABB( sumSquared, SNR_Q7, SNR_Q7 );          /* Q14 */
181
182             /* Tilt measure */
183             if( speech_nrg < ( 1 << 20 ) ) {
184                 /* Scale down SNR value for small subband speech energies */
185                 SNR_Q7 = SKP_SMULWB( SKP_LSHIFT( silk_SQRT_APPROX( speech_nrg ), 6 ), SNR_Q7 );
186             }
187             input_tilt = SKP_SMLAWB( input_tilt, tiltWeights[ b ], SNR_Q7 );
188         } else {
189             NrgToNoiseRatio_Q8[ b ] = 256;
190         }
191     }
192
193     /* Mean-of-squares */
194     sumSquared = SKP_DIV32_16( sumSquared, VAD_N_BANDS ); /* Q14 */
195
196     /* Root-mean-square approximation, scale to dBs, and write to output pointer */
197     pSNR_dB_Q7 = ( opus_int16 )( 3 * silk_SQRT_APPROX( sumSquared ) ); /* Q7 */
198
199     /*********************************/
200     /* Speech Probability Estimation */
201     /*********************************/
202     SA_Q15 = silk_sigm_Q15( SKP_SMULWB( VAD_SNR_FACTOR_Q16, pSNR_dB_Q7 ) - VAD_NEGATIVE_OFFSET_Q5 );
203
204     /**************************/
205     /* Frequency Tilt Measure */
206     /**************************/
207     psEncC->input_tilt_Q15 = SKP_LSHIFT( silk_sigm_Q15( input_tilt ) - 16384, 1 );
208
209     /**************************************************/
210     /* Scale the sigmoid output based on power levels */
211     /**************************************************/
212     speech_nrg = 0;
213     for( b = 0; b < VAD_N_BANDS; b++ ) {
214         /* Accumulate signal-without-noise energies, higher frequency bands have more weight */
215         speech_nrg += ( b + 1 ) * SKP_RSHIFT( Xnrg[ b ] - psSilk_VAD->NL[ b ], 4 );
216     }
217
218     /* Power scaling */
219     if( speech_nrg <= 0 ) {
220         SA_Q15 = SKP_RSHIFT( SA_Q15, 1 );
221     } else if( speech_nrg < 32768 ) {
222         if( psEncC->frame_length == 10 * psEncC->fs_kHz ) {
223             speech_nrg = SKP_LSHIFT_SAT32( speech_nrg, 16 );
224         } else {
225             speech_nrg = SKP_LSHIFT_SAT32( speech_nrg, 15 );
226         }
227
228         /* square-root */
229         speech_nrg = silk_SQRT_APPROX( speech_nrg );
230         SA_Q15 = SKP_SMULWB( 32768 + speech_nrg, SA_Q15 );
231     }
232
233     /* Copy the resulting speech activity in Q8 */
234     psEncC->speech_activity_Q8 = SKP_min_int( SKP_RSHIFT( SA_Q15, 7 ), SKP_uint8_MAX );
235
236     /***********************************/
237     /* Energy Level and SNR estimation */
238     /***********************************/
239     /* Smoothing coefficient */
240     smooth_coef_Q16 = SKP_SMULWB( VAD_SNR_SMOOTH_COEF_Q18, SKP_SMULWB( SA_Q15, SA_Q15 ) );
241
242     if( psEncC->frame_length == 10 * psEncC->fs_kHz ) {
243         smooth_coef_Q16 >>= 1;
244     }
245
246     for( b = 0; b < VAD_N_BANDS; b++ ) {
247         /* compute smoothed energy-to-noise ratio per band */
248         psSilk_VAD->NrgRatioSmth_Q8[ b ] = SKP_SMLAWB( psSilk_VAD->NrgRatioSmth_Q8[ b ],
249             NrgToNoiseRatio_Q8[ b ] - psSilk_VAD->NrgRatioSmth_Q8[ b ], smooth_coef_Q16 );
250
251         /* signal to noise ratio in dB per band */
252         SNR_Q7 = 3 * ( silk_lin2log( psSilk_VAD->NrgRatioSmth_Q8[b] ) - 8 * 128 );
253         /* quality = sigmoid( 0.25 * ( SNR_dB - 16 ) ); */
254         psEncC->input_quality_bands_Q15[ b ] = silk_sigm_Q15( SKP_RSHIFT( SNR_Q7 - 16 * 128, 4 ) );
255     }
256
257     return( ret );
258 }
259
260 /**************************/
261 /* Noise level estimation */
262 /**************************/
263 void silk_VAD_GetNoiseLevels(
264     const opus_int32                 pX[ VAD_N_BANDS ],  /* I    subband energies                            */
265     silk_VAD_state              *psSilk_VAD         /* I/O  Pointer to Silk VAD state                   */
266 )
267 {
268     opus_int   k;
269     opus_int32 nl, nrg, inv_nrg;
270     opus_int   coef, min_coef;
271
272     /* Initially faster smoothing */
273     if( psSilk_VAD->counter < 1000 ) { /* 1000 = 20 sec */
274         min_coef = SKP_DIV32_16( SKP_int16_MAX, SKP_RSHIFT( psSilk_VAD->counter, 4 ) + 1 );
275     } else {
276         min_coef = 0;
277     }
278
279     for( k = 0; k < VAD_N_BANDS; k++ ) {
280         /* Get old noise level estimate for current band */
281         nl = psSilk_VAD->NL[ k ];
282         SKP_assert( nl >= 0 );
283
284         /* Add bias */
285         nrg = SKP_ADD_POS_SAT32( pX[ k ], psSilk_VAD->NoiseLevelBias[ k ] );
286         SKP_assert( nrg > 0 );
287
288         /* Invert energies */
289         inv_nrg = SKP_DIV32( SKP_int32_MAX, nrg );
290         SKP_assert( inv_nrg >= 0 );
291
292         /* Less update when subband energy is high */
293         if( nrg > SKP_LSHIFT( nl, 3 ) ) {
294             coef = VAD_NOISE_LEVEL_SMOOTH_COEF_Q16 >> 3;
295         } else if( nrg < nl ) {
296             coef = VAD_NOISE_LEVEL_SMOOTH_COEF_Q16;
297         } else {
298             coef = SKP_SMULWB( SKP_SMULWW( inv_nrg, nl ), VAD_NOISE_LEVEL_SMOOTH_COEF_Q16 << 1 );
299         }
300
301         /* Initially faster smoothing */
302         coef = SKP_max_int( coef, min_coef );
303
304         /* Smooth inverse energies */
305         psSilk_VAD->inv_NL[ k ] = SKP_SMLAWB( psSilk_VAD->inv_NL[ k ], inv_nrg - psSilk_VAD->inv_NL[ k ], coef );
306         SKP_assert( psSilk_VAD->inv_NL[ k ] >= 0 );
307
308         /* Compute noise level by inverting again */
309         nl = SKP_DIV32( SKP_int32_MAX, psSilk_VAD->inv_NL[ k ] );
310         SKP_assert( nl >= 0 );
311
312         /* Limit noise levels (guarantee 7 bits of head room) */
313         nl = SKP_min( nl, 0x00FFFFFF );
314
315         /* Store as part of state */
316         psSilk_VAD->NL[ k ] = nl;
317     }
318
319     /* Increment frame counter */
320     psSilk_VAD->counter++;
321 }