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_FLP.h"
33 #include "tuning_parameters.h"
34
35 /* Compute gain to make warped filter coefficients have a zero mean log frequency response on a */
36 /* non-warped frequency scale. (So that it can be implemented with a minimum-phase monic filter.) */
37 /* Note: A monic filter is one with the first coefficient equal to 1.0. In Silk we omit the first */
38 /* coefficient in an array of coefficients, for monic filters. */
warped_gain(const silk_float * coefs,silk_float lambda,opus_int order)39 static OPUS_INLINE silk_float warped_gain(
40 const silk_float *coefs,
41 silk_float lambda,
42 opus_int order
43 ) {
44 opus_int i;
45 silk_float gain;
46
47 lambda = -lambda;
48 gain = coefs[ order - 1 ];
49 for( i = order - 2; i >= 0; i-- ) {
50 gain = lambda * gain + coefs[ i ];
51 }
52 return (silk_float)( 1.0f / ( 1.0f - lambda * gain ) );
53 }
54
55 /* Convert warped filter coefficients to monic pseudo-warped coefficients and limit maximum */
56 /* amplitude of monic warped coefficients by using bandwidth expansion on the true coefficients */
warped_true2monic_coefs(silk_float * coefs_syn,silk_float * coefs_ana,silk_float lambda,silk_float limit,opus_int order)57 static OPUS_INLINE void warped_true2monic_coefs(
58 silk_float *coefs_syn,
59 silk_float *coefs_ana,
60 silk_float lambda,
61 silk_float limit,
62 opus_int order
63 ) {
64 opus_int i, iter, ind = 0;
65 silk_float tmp, maxabs, chirp, gain_syn, gain_ana;
66
67 /* Convert to monic coefficients */
68 for( i = order - 1; i > 0; i-- ) {
69 coefs_syn[ i - 1 ] -= lambda * coefs_syn[ i ];
70 coefs_ana[ i - 1 ] -= lambda * coefs_ana[ i ];
71 }
72 gain_syn = ( 1.0f - lambda * lambda ) / ( 1.0f + lambda * coefs_syn[ 0 ] );
73 gain_ana = ( 1.0f - lambda * lambda ) / ( 1.0f + lambda * coefs_ana[ 0 ] );
74 for( i = 0; i < order; i++ ) {
75 coefs_syn[ i ] *= gain_syn;
76 coefs_ana[ i ] *= gain_ana;
77 }
78
79 /* Limit */
80 for( iter = 0; iter < 10; iter++ ) {
81 /* Find maximum absolute value */
82 maxabs = -1.0f;
83 for( i = 0; i < order; i++ ) {
84 tmp = silk_max( silk_abs_float( coefs_syn[ i ] ), silk_abs_float( coefs_ana[ i ] ) );
85 if( tmp > maxabs ) {
86 maxabs = tmp;
87 ind = i;
88 }
89 }
90 if( maxabs <= limit ) {
91 /* Coefficients are within range - done */
92 return;
93 }
94
95 /* Convert back to true warped coefficients */
96 for( i = 1; i < order; i++ ) {
97 coefs_syn[ i - 1 ] += lambda * coefs_syn[ i ];
98 coefs_ana[ i - 1 ] += lambda * coefs_ana[ i ];
99 }
100 gain_syn = 1.0f / gain_syn;
101 gain_ana = 1.0f / gain_ana;
102 for( i = 0; i < order; i++ ) {
103 coefs_syn[ i ] *= gain_syn;
104 coefs_ana[ i ] *= gain_ana;
105 }
106
107 /* Apply bandwidth expansion */
108 chirp = 0.99f - ( 0.8f + 0.1f * iter ) * ( maxabs - limit ) / ( maxabs * ( ind + 1 ) );
109 silk_bwexpander_FLP( coefs_syn, order, chirp );
110 silk_bwexpander_FLP( coefs_ana, order, chirp );
111
112 /* Convert to monic warped coefficients */
113 for( i = order - 1; i > 0; i-- ) {
114 coefs_syn[ i - 1 ] -= lambda * coefs_syn[ i ];
115 coefs_ana[ i - 1 ] -= lambda * coefs_ana[ i ];
116 }
117 gain_syn = ( 1.0f - lambda * lambda ) / ( 1.0f + lambda * coefs_syn[ 0 ] );
118 gain_ana = ( 1.0f - lambda * lambda ) / ( 1.0f + lambda * coefs_ana[ 0 ] );
119 for( i = 0; i < order; i++ ) {
120 coefs_syn[ i ] *= gain_syn;
121 coefs_ana[ i ] *= gain_ana;
122 }
123 }
124 silk_assert( 0 );
125 }
126
127 /* Compute noise shaping coefficients and initial gain values */
silk_noise_shape_analysis_FLP(silk_encoder_state_FLP * psEnc,silk_encoder_control_FLP * psEncCtrl,const silk_float * pitch_res,const silk_float * x)128 void silk_noise_shape_analysis_FLP(
129 silk_encoder_state_FLP *psEnc, /* I/O Encoder state FLP */
130 silk_encoder_control_FLP *psEncCtrl, /* I/O Encoder control FLP */
131 const silk_float *pitch_res, /* I LPC residual from pitch analysis */
132 const silk_float *x /* I Input signal [frame_length + la_shape] */
133 )
134 {
135 silk_shape_state_FLP *psShapeSt = &psEnc->sShape;
136 opus_int k, nSamples;
137 silk_float SNR_adj_dB, HarmBoost, HarmShapeGain, Tilt;
138 silk_float nrg, pre_nrg, log_energy, log_energy_prev, energy_variation;
139 silk_float delta, BWExp1, BWExp2, gain_mult, gain_add, strength, b, warping;
140 silk_float x_windowed[ SHAPE_LPC_WIN_MAX ];
141 silk_float auto_corr[ MAX_SHAPE_LPC_ORDER + 1 ];
142 const silk_float *x_ptr, *pitch_res_ptr;
143
144 /* Point to start of first LPC analysis block */
145 x_ptr = x - psEnc->sCmn.la_shape;
146
147 /****************/
148 /* GAIN CONTROL */
149 /****************/
150 SNR_adj_dB = psEnc->sCmn.SNR_dB_Q7 * ( 1 / 128.0f );
151
152 /* Input quality is the average of the quality in the lowest two VAD bands */
153 psEncCtrl->input_quality = 0.5f * ( psEnc->sCmn.input_quality_bands_Q15[ 0 ] + psEnc->sCmn.input_quality_bands_Q15[ 1 ] ) * ( 1.0f / 32768.0f );
154
155 /* Coding quality level, between 0.0 and 1.0 */
156 psEncCtrl->coding_quality = silk_sigmoid( 0.25f * ( SNR_adj_dB - 20.0f ) );
157
158 if( psEnc->sCmn.useCBR == 0 ) {
159 /* Reduce coding SNR during low speech activity */
160 b = 1.0f - psEnc->sCmn.speech_activity_Q8 * ( 1.0f / 256.0f );
161 SNR_adj_dB -= BG_SNR_DECR_dB * psEncCtrl->coding_quality * ( 0.5f + 0.5f * psEncCtrl->input_quality ) * b * b;
162 }
163
164 if( psEnc->sCmn.indices.signalType == TYPE_VOICED ) {
165 /* Reduce gains for periodic signals */
166 SNR_adj_dB += HARM_SNR_INCR_dB * psEnc->LTPCorr;
167 } else {
168 /* For unvoiced signals and low-quality input, adjust the quality slower than SNR_dB setting */
169 SNR_adj_dB += ( -0.4f * psEnc->sCmn.SNR_dB_Q7 * ( 1 / 128.0f ) + 6.0f ) * ( 1.0f - psEncCtrl->input_quality );
170 }
171
172 /*************************/
173 /* SPARSENESS PROCESSING */
174 /*************************/
175 /* Set quantizer offset */
176 if( psEnc->sCmn.indices.signalType == TYPE_VOICED ) {
177 /* Initially set to 0; may be overruled in process_gains(..) */
178 psEnc->sCmn.indices.quantOffsetType = 0;
179 psEncCtrl->sparseness = 0.0f;
180 } else {
181 /* Sparseness measure, based on relative fluctuations of energy per 2 milliseconds */
182 nSamples = 2 * psEnc->sCmn.fs_kHz;
183 energy_variation = 0.0f;
184 log_energy_prev = 0.0f;
185 pitch_res_ptr = pitch_res;
186 for( k = 0; k < silk_SMULBB( SUB_FRAME_LENGTH_MS, psEnc->sCmn.nb_subfr ) / 2; k++ ) {
187 nrg = ( silk_float )nSamples + ( silk_float )silk_energy_FLP( pitch_res_ptr, nSamples );
188 log_energy = silk_log2( nrg );
189 if( k > 0 ) {
190 energy_variation += silk_abs_float( log_energy - log_energy_prev );
191 }
192 log_energy_prev = log_energy;
193 pitch_res_ptr += nSamples;
194 }
195 psEncCtrl->sparseness = silk_sigmoid( 0.4f * ( energy_variation - 5.0f ) );
196
197 /* Set quantization offset depending on sparseness measure */
198 if( psEncCtrl->sparseness > SPARSENESS_THRESHOLD_QNT_OFFSET ) {
199 psEnc->sCmn.indices.quantOffsetType = 0;
200 } else {
201 psEnc->sCmn.indices.quantOffsetType = 1;
202 }
203
204 /* Increase coding SNR for sparse signals */
205 SNR_adj_dB += SPARSE_SNR_INCR_dB * ( psEncCtrl->sparseness - 0.5f );
206 }
207
208 /*******************************/
209 /* Control bandwidth expansion */
210 /*******************************/
211 /* More BWE for signals with high prediction gain */
212 strength = FIND_PITCH_WHITE_NOISE_FRACTION * psEncCtrl->predGain; /* between 0.0 and 1.0 */
213 BWExp1 = BWExp2 = BANDWIDTH_EXPANSION / ( 1.0f + strength * strength );
214 delta = LOW_RATE_BANDWIDTH_EXPANSION_DELTA * ( 1.0f - 0.75f * psEncCtrl->coding_quality );
215 BWExp1 -= delta;
216 BWExp2 += delta;
217 /* BWExp1 will be applied after BWExp2, so make it relative */
218 BWExp1 /= BWExp2;
219
220 if( psEnc->sCmn.warping_Q16 > 0 ) {
221 /* Slightly more warping in analysis will move quantization noise up in frequency, where it's better masked */
222 warping = (silk_float)psEnc->sCmn.warping_Q16 / 65536.0f + 0.01f * psEncCtrl->coding_quality;
223 } else {
224 warping = 0.0f;
225 }
226
227 /********************************************/
228 /* Compute noise shaping AR coefs and gains */
229 /********************************************/
230 for( k = 0; k < psEnc->sCmn.nb_subfr; k++ ) {
231 /* Apply window: sine slope followed by flat part followed by cosine slope */
232 opus_int shift, slope_part, flat_part;
233 flat_part = psEnc->sCmn.fs_kHz * 3;
234 slope_part = ( psEnc->sCmn.shapeWinLength - flat_part ) / 2;
235
236 silk_apply_sine_window_FLP( x_windowed, x_ptr, 1, slope_part );
237 shift = slope_part;
238 silk_memcpy( x_windowed + shift, x_ptr + shift, flat_part * sizeof(silk_float) );
239 shift += flat_part;
240 silk_apply_sine_window_FLP( x_windowed + shift, x_ptr + shift, 2, slope_part );
241
242 /* Update pointer: next LPC analysis block */
243 x_ptr += psEnc->sCmn.subfr_length;
244
245 if( psEnc->sCmn.warping_Q16 > 0 ) {
246 /* Calculate warped auto correlation */
247 silk_warped_autocorrelation_FLP( auto_corr, x_windowed, warping,
248 psEnc->sCmn.shapeWinLength, psEnc->sCmn.shapingLPCOrder );
249 } else {
250 /* Calculate regular auto correlation */
251 silk_autocorrelation_FLP( auto_corr, x_windowed, psEnc->sCmn.shapeWinLength, psEnc->sCmn.shapingLPCOrder + 1 );
252 }
253
254 /* Add white noise, as a fraction of energy */
255 auto_corr[ 0 ] += auto_corr[ 0 ] * SHAPE_WHITE_NOISE_FRACTION;
256
257 /* Convert correlations to prediction coefficients, and compute residual energy */
258 nrg = silk_levinsondurbin_FLP( &psEncCtrl->AR2[ k * MAX_SHAPE_LPC_ORDER ], auto_corr, psEnc->sCmn.shapingLPCOrder );
259 psEncCtrl->Gains[ k ] = ( silk_float )sqrt( nrg );
260
261 if( psEnc->sCmn.warping_Q16 > 0 ) {
262 /* Adjust gain for warping */
263 psEncCtrl->Gains[ k ] *= warped_gain( &psEncCtrl->AR2[ k * MAX_SHAPE_LPC_ORDER ], warping, psEnc->sCmn.shapingLPCOrder );
264 }
265
266 /* Bandwidth expansion for synthesis filter shaping */
267 silk_bwexpander_FLP( &psEncCtrl->AR2[ k * MAX_SHAPE_LPC_ORDER ], psEnc->sCmn.shapingLPCOrder, BWExp2 );
268
269 /* Compute noise shaping filter coefficients */
270 silk_memcpy(
271 &psEncCtrl->AR1[ k * MAX_SHAPE_LPC_ORDER ],
272 &psEncCtrl->AR2[ k * MAX_SHAPE_LPC_ORDER ],
273 psEnc->sCmn.shapingLPCOrder * sizeof( silk_float ) );
274
275 /* Bandwidth expansion for analysis filter shaping */
276 silk_bwexpander_FLP( &psEncCtrl->AR1[ k * MAX_SHAPE_LPC_ORDER ], psEnc->sCmn.shapingLPCOrder, BWExp1 );
277
278 /* Ratio of prediction gains, in energy domain */
279 pre_nrg = silk_LPC_inverse_pred_gain_FLP( &psEncCtrl->AR2[ k * MAX_SHAPE_LPC_ORDER ], psEnc->sCmn.shapingLPCOrder );
280 nrg = silk_LPC_inverse_pred_gain_FLP( &psEncCtrl->AR1[ k * MAX_SHAPE_LPC_ORDER ], psEnc->sCmn.shapingLPCOrder );
281 psEncCtrl->GainsPre[ k ] = 1.0f - 0.7f * ( 1.0f - pre_nrg / nrg );
282
283 /* Convert to monic warped prediction coefficients and limit absolute values */
284 warped_true2monic_coefs( &psEncCtrl->AR2[ k * MAX_SHAPE_LPC_ORDER ], &psEncCtrl->AR1[ k * MAX_SHAPE_LPC_ORDER ],
285 warping, 3.999f, psEnc->sCmn.shapingLPCOrder );
286 }
287
288 /*****************/
289 /* Gain tweaking */
290 /*****************/
291 /* Increase gains during low speech activity */
292 gain_mult = (silk_float)pow( 2.0f, -0.16f * SNR_adj_dB );
293 gain_add = (silk_float)pow( 2.0f, 0.16f * MIN_QGAIN_DB );
294 for( k = 0; k < psEnc->sCmn.nb_subfr; k++ ) {
295 psEncCtrl->Gains[ k ] *= gain_mult;
296 psEncCtrl->Gains[ k ] += gain_add;
297 }
298
299 gain_mult = 1.0f + INPUT_TILT + psEncCtrl->coding_quality * HIGH_RATE_INPUT_TILT;
300 for( k = 0; k < psEnc->sCmn.nb_subfr; k++ ) {
301 psEncCtrl->GainsPre[ k ] *= gain_mult;
302 }
303
304 /************************************************/
305 /* Control low-frequency shaping and noise tilt */
306 /************************************************/
307 /* Less low frequency shaping for noisy inputs */
308 strength = LOW_FREQ_SHAPING * ( 1.0f + LOW_QUALITY_LOW_FREQ_SHAPING_DECR * ( psEnc->sCmn.input_quality_bands_Q15[ 0 ] * ( 1.0f / 32768.0f ) - 1.0f ) );
309 strength *= psEnc->sCmn.speech_activity_Q8 * ( 1.0f / 256.0f );
310 if( psEnc->sCmn.indices.signalType == TYPE_VOICED ) {
311 /* Reduce low frequencies quantization noise for periodic signals, depending on pitch lag */
312 /*f = 400; freqz([1, -0.98 + 2e-4 * f], [1, -0.97 + 7e-4 * f], 2^12, Fs); axis([0, 1000, -10, 1])*/
313 for( k = 0; k < psEnc->sCmn.nb_subfr; k++ ) {
314 b = 0.2f / psEnc->sCmn.fs_kHz + 3.0f / psEncCtrl->pitchL[ k ];
315 psEncCtrl->LF_MA_shp[ k ] = -1.0f + b;
316 psEncCtrl->LF_AR_shp[ k ] = 1.0f - b - b * strength;
317 }
318 Tilt = - HP_NOISE_COEF -
319 (1 - HP_NOISE_COEF) * HARM_HP_NOISE_COEF * psEnc->sCmn.speech_activity_Q8 * ( 1.0f / 256.0f );
320 } else {
321 b = 1.3f / psEnc->sCmn.fs_kHz;
322 psEncCtrl->LF_MA_shp[ 0 ] = -1.0f + b;
323 psEncCtrl->LF_AR_shp[ 0 ] = 1.0f - b - b * strength * 0.6f;
324 for( k = 1; k < psEnc->sCmn.nb_subfr; k++ ) {
325 psEncCtrl->LF_MA_shp[ k ] = psEncCtrl->LF_MA_shp[ 0 ];
326 psEncCtrl->LF_AR_shp[ k ] = psEncCtrl->LF_AR_shp[ 0 ];
327 }
328 Tilt = -HP_NOISE_COEF;
329 }
330
331 /****************************/
332 /* HARMONIC SHAPING CONTROL */
333 /****************************/
334 /* Control boosting of harmonic frequencies */
335 HarmBoost = LOW_RATE_HARMONIC_BOOST * ( 1.0f - psEncCtrl->coding_quality ) * psEnc->LTPCorr;
336
337 /* More harmonic boost for noisy input signals */
338 HarmBoost += LOW_INPUT_QUALITY_HARMONIC_BOOST * ( 1.0f - psEncCtrl->input_quality );
339
340 if( USE_HARM_SHAPING && psEnc->sCmn.indices.signalType == TYPE_VOICED ) {
341 /* Harmonic noise shaping */
342 HarmShapeGain = HARMONIC_SHAPING;
343
344 /* More harmonic noise shaping for high bitrates or noisy input */
345 HarmShapeGain += HIGH_RATE_OR_LOW_QUALITY_HARMONIC_SHAPING *
346 ( 1.0f - ( 1.0f - psEncCtrl->coding_quality ) * psEncCtrl->input_quality );
347
348 /* Less harmonic noise shaping for less periodic signals */
349 HarmShapeGain *= ( silk_float )sqrt( psEnc->LTPCorr );
350 } else {
351 HarmShapeGain = 0.0f;
352 }
353
354 /*************************/
355 /* Smooth over subframes */
356 /*************************/
357 for( k = 0; k < psEnc->sCmn.nb_subfr; k++ ) {
358 psShapeSt->HarmBoost_smth += SUBFR_SMTH_COEF * ( HarmBoost - psShapeSt->HarmBoost_smth );
359 psEncCtrl->HarmBoost[ k ] = psShapeSt->HarmBoost_smth;
360 psShapeSt->HarmShapeGain_smth += SUBFR_SMTH_COEF * ( HarmShapeGain - psShapeSt->HarmShapeGain_smth );
361 psEncCtrl->HarmShapeGain[ k ] = psShapeSt->HarmShapeGain_smth;
362 psShapeSt->Tilt_smth += SUBFR_SMTH_COEF * ( Tilt - psShapeSt->Tilt_smth );
363 psEncCtrl->Tilt[ k ] = psShapeSt->Tilt_smth;
364 }
365 }
366