Line data Source code
1 : /******************************************************************************************************
2 :
3 : (C) 2022-2025 IVAS codec Public Collaboration with portions copyright Dolby International AB, Ericsson AB,
4 : Fraunhofer-Gesellschaft zur Foerderung der angewandten Forschung e.V., Huawei Technologies Co. LTD.,
5 : Koninklijke Philips N.V., Nippon Telegraph and Telephone Corporation, Nokia Technologies Oy, Orange,
6 : Panasonic Holdings Corporation, Qualcomm Technologies, Inc., VoiceAge Corporation, and other
7 : contributors to this repository. All Rights Reserved.
8 :
9 : This software is protected by copyright law and by international treaties.
10 : The IVAS codec Public Collaboration consisting of Dolby International AB, Ericsson AB,
11 : Fraunhofer-Gesellschaft zur Foerderung der angewandten Forschung e.V., Huawei Technologies Co. LTD.,
12 : Koninklijke Philips N.V., Nippon Telegraph and Telephone Corporation, Nokia Technologies Oy, Orange,
13 : Panasonic Holdings Corporation, Qualcomm Technologies, Inc., VoiceAge Corporation, and other
14 : contributors to this repository retain full ownership rights in their respective contributions in
15 : the software. This notice grants no license of any kind, including but not limited to patent
16 : license, nor is any license granted by implication, estoppel or otherwise.
17 :
18 : Contributors are required to enter into the IVAS codec Public Collaboration agreement before making
19 : contributions.
20 :
21 : This software is provided "AS IS", without any express or implied warranties. The software is in the
22 : development stage. It is intended exclusively for experts who have experience with such software and
23 : solely for the purpose of inspection. All implied warranties of non-infringement, merchantability
24 : and fitness for a particular purpose are hereby disclaimed and excluded.
25 :
26 : Any dispute, controversy or claim arising under or in relation to providing this software shall be
27 : submitted to and settled by the final, binding jurisdiction of the courts of Munich, Germany in
28 : accordance with the laws of the Federal Republic of Germany excluding its conflict of law rules and
29 : the United Nations Convention on Contracts on the International Sales of Goods.
30 :
31 : *******************************************************************************************************/
32 :
33 : #include <stdint.h>
34 : #include "options.h"
35 : #include "prot.h"
36 : #include "ivas_prot.h"
37 : #include "ivas_rom_com.h"
38 : #include <math.h>
39 : #ifdef DEBUGGING
40 : #include "debug.h"
41 : #endif
42 : #include "wmc_auto.h"
43 :
44 : /*------------------------------------------------------------------------------------------*
45 : * Local constants
46 : *------------------------------------------------------------------------------------------*/
47 :
48 : #define LFE_PLC_DSF ( 48000 / LFE_PLC_FS )
49 : #define LFE_PLC_LPCORD ( MAX_LP_FILTER_ORDER )
50 : #define POW_THR ( 1.0e-8f )
51 : #define LFE_PLC_RECLEN_48K ( ( IVAS_LFE_NUM_COEFFS_IN_SUBGRP + 1 ) * L_FRAME48k / IVAS_LFE_NUM_COEFFS_IN_SUBGRP + LFE_PLC_FDEL )
52 : #define LFE_PLC_RECLEN ( ( LFE_PLC_RECLEN_48K / LFE_PLC_DSF ) )
53 : #define LFE_PLC_MUTE_THR ( 10 )
54 : #define LFE_PLC_BURST_ATT ( powf( powf( 10.0f, -3.0f / 20.0f ), 1.0f / ( LFE_PLC_FS * 20 / 1000 ) ) ) /* attenuate 3dB per frame starting with 10th consecutive loss */
55 : #define EPS_STOP 1e-5f
56 :
57 : /*---------------------------------------------------------------------*
58 : * lfeplc_lev_dur()
59 : *
60 : * Wiener-Levinson-Durbin algorithm to compute LP parameters from the autocorrelations
61 : * of input signal
62 : *---------------------------------------------------------------------*/
63 :
64 : /*! r: stability flag */
65 90 : static int16_t lfeplc_lev_dur(
66 : float *a_out, /* o : LP coefficients (a[0] = 1.0) */
67 : const float *r, /* i : vector of autocorrelations */
68 : const int16_t m /* i : order of LP filter */
69 : )
70 : {
71 : int16_t i, j, l;
72 : float buf[TCXLTP_LTP_ORDER];
73 : float *rc; /* reflection coefficients 0,...,m-1 */
74 : float s, at, err;
75 : float a[LFE_PLC_LPCORD + 1];
76 :
77 90 : set_f( a, 0.f, LFE_PLC_LPCORD + 1 );
78 90 : set_f( a_out, 0.f, LFE_PLC_LPCORD + 1 );
79 90 : rc = &buf[0];
80 90 : a[0] = 1.f;
81 90 : err = r[0];
82 90 : a_out[0] = a[0];
83 90 : rc[0] = ( -r[1] ) / r[0];
84 90 : s = r[1];
85 90 : a[1] = rc[0];
86 90 : a_out[1] = a[1];
87 :
88 90 : i = 1;
89 1800 : while ( i < m )
90 : {
91 9810 : for ( j = 1; j <= i / 2; j++ )
92 : {
93 8100 : l = i - j;
94 8100 : at = a[j] + rc[i - 1] * a[l];
95 8100 : a[l] += rc[i - 1] * a[j];
96 8100 : a[j] = at;
97 : }
98 :
99 1710 : a[i] = rc[i - 1];
100 1710 : err += rc[i - 1] * s;
101 1710 : s = 0.f;
102 1710 : i++;
103 20520 : for ( j = 0; j < i; j++ )
104 : {
105 18810 : s += r[i - j] * a[j];
106 : }
107 1710 : rc[i - 1] = ( -s ) / err;
108 1710 : if ( fabsf( rc[i - 1] ) > 0.99945f )
109 : {
110 0 : return 1;
111 : }
112 : else
113 : {
114 37620 : for ( j = 0; j <= m; j++ )
115 : {
116 35910 : a_out[j] = a[j];
117 : }
118 : }
119 : }
120 :
121 90 : return 0;
122 : }
123 :
124 :
125 : /*-----------------------------------------------------------------------------------------*
126 : * Function check_stab()
127 : *
128 : * LPC filter stability check applying given sharpening value delta
129 : *-----------------------------------------------------------------------------------------*/
130 :
131 1146 : static uint16_t check_stab(
132 : const float *a,
133 : const float delta )
134 : {
135 : float amod[LFE_PLC_LPCORD], refl[LFE_PLC_LPCORD];
136 : float fac;
137 : float fac1;
138 : int16_t i;
139 : uint16_t stable;
140 :
141 1146 : fac = 1.0f + delta;
142 1146 : fac1 = fac;
143 :
144 24066 : for ( i = 0; i < LFE_PLC_LPCORD; i++ )
145 : {
146 22920 : amod[i] = a[i] * fac;
147 22920 : fac *= fac1;
148 : }
149 :
150 1146 : stable = a2rc( amod, refl, LFE_PLC_LPCORD );
151 :
152 1146 : return stable;
153 : }
154 :
155 :
156 : /*-----------------------------------------------------------------------------------------*
157 : * Function find_max_delta()
158 : *
159 : * Find maximum LPC filter sharpening by iteration to get a filter that is almost instable
160 : *-----------------------------------------------------------------------------------------*/
161 :
162 90 : static float find_max_delta(
163 : float *a )
164 : {
165 : float delta;
166 : float eps;
167 : float fac;
168 : uint16_t stable;
169 :
170 90 : delta = 0.0f;
171 90 : eps = 0.01f;
172 90 : fac = 2;
173 90 : stable = FALSE;
174 :
175 171 : while ( check_stab( a, eps ) )
176 : {
177 :
178 81 : eps *= fac;
179 81 : stable = TRUE;
180 : }
181 90 : fac = 0.5f;
182 :
183 90 : if ( stable )
184 : {
185 48 : eps *= fac;
186 : }
187 :
188 162 : while ( !stable )
189 : {
190 72 : eps *= fac;
191 72 : stable = check_stab( a, eps );
192 : }
193 :
194 : /* must be stable with current eps */
195 90 : delta = eps;
196 90 : eps *= fac;
197 :
198 : while ( 1 )
199 : {
200 903 : delta += eps;
201 903 : stable = check_stab( a, delta );
202 :
203 903 : if ( !stable )
204 : {
205 444 : if ( fabsf( eps ) > EPS_STOP )
206 : {
207 402 : eps = -fabsf( eps ) * fac;
208 : }
209 : else
210 : {
211 42 : eps = -fabsf( eps );
212 : }
213 : }
214 : else
215 : {
216 459 : if ( fabsf( eps ) < EPS_STOP )
217 : {
218 90 : break;
219 : }
220 369 : eps = fabsf( eps ) * fac;
221 : }
222 : }
223 :
224 90 : return delta;
225 : }
226 :
227 :
228 : /*-----------------------------------------------------------------------------------------*
229 : * Function recover_samples()
230 : *
231 : * recover lost samples by extrapolation of signal buffer
232 : *-----------------------------------------------------------------------------------------*/
233 :
234 660 : static void recover_samples(
235 : const int16_t bfi_count,
236 : float *outbuf,
237 : float *rec_frame )
238 : {
239 : int16_t i;
240 : float zeroes[LFE_PLC_RECLEN];
241 : float delta, fac, att;
242 : float r[LFE_PLC_LPCORD + 1], a[LFE_PLC_LPCORD + 1];
243 :
244 660 : autocorr( outbuf, r, LFE_PLC_LPCORD, LFE_PLC_BUFLEN, hamm_lfe_plc, 0, 1, 1 );
245 :
246 660 : if ( r[0] < POW_THR * LFE_PLC_BUFLEN )
247 : {
248 570 : set_zero( rec_frame, LFE_PLC_RECLEN );
249 570 : return;
250 : }
251 :
252 90 : lfeplc_lev_dur( a, r, LFE_PLC_LPCORD );
253 :
254 90 : delta = find_max_delta( a + 1 );
255 :
256 90 : fac = 1.0f + delta;
257 90 : att = 1.0f;
258 :
259 90 : if ( bfi_count >= LFE_PLC_MUTE_THR )
260 : {
261 0 : att = LFE_PLC_BURST_ATT;
262 0 : fac *= att;
263 : }
264 :
265 1890 : for ( i = 1; i <= LFE_PLC_LPCORD; i++ )
266 : {
267 1800 : a[i] = a[i] * fac;
268 1800 : fac *= att * ( 1.0f + delta );
269 : }
270 :
271 90 : set_zero( zeroes, LFE_PLC_RECLEN );
272 :
273 90 : syn_filt( a, LFE_PLC_LPCORD, zeroes, rec_frame, LFE_PLC_RECLEN, outbuf + LFE_PLC_BUFLEN - LFE_PLC_LPCORD, 0 );
274 :
275 90 : return;
276 : }
277 :
278 :
279 : /*-----------------------------------------------------------------------------------------*
280 : * Function ivas_lfe_tdplc()
281 : *
282 : * MDCT interface recover lost samples by extrapolation of signal buffer
283 : *-----------------------------------------------------------------------------------------*/
284 :
285 660 : void ivas_lfe_tdplc(
286 : LFE_DEC_HANDLE hLFE, /* i/o: LFE decoder handle */
287 : float *prevsynth, /* i : previous frame synthesis */
288 : float *ytda, /* o : output time-domain buffer */
289 : const int16_t output_frame /* i : output frame length */
290 : )
291 : {
292 : float rec_frame_us[LFE_PLC_RECLEN_48K], input_tda[L_FRAME48k];
293 : float rec_frame[LFE_PLC_RECLEN], mem[2 * LFE_PLC_FDEL / LFE_PLC_DSF];
294 : int32_t output_Fs;
295 : int16_t i, fade_len, full_len, dct_len, zero_pad_len, plc_fdel, rec_frame_len;
296 : const float *pWindow_coeffs;
297 :
298 660 : output_Fs = output_frame * FRAMES_PER_SEC;
299 660 : fade_len = hLFE->pWindow_state->fade_len;
300 660 : full_len = hLFE->pWindow_state->full_len;
301 660 : dct_len = hLFE->pWindow_state->dct_len;
302 660 : zero_pad_len = hLFE->pWindow_state->zero_pad_len;
303 660 : pWindow_coeffs = hLFE->pWindow_state->pWindow_coeffs;
304 :
305 660 : plc_fdel = (int16_t) ( LFE_PLC_FDEL * output_Fs / 48000 );
306 660 : rec_frame_len = (int16_t) ( LFE_PLC_RECLEN_48K * output_Fs / 48000 );
307 660 : recover_samples( hLFE->bfi_count, prevsynth, rec_frame );
308 :
309 660 : set_zero( mem, 2 * LFE_PLC_FDEL / LFE_PLC_DSF );
310 :
311 660 : modify_Fs( prevsynth + LFE_PLC_BUFLEN - LFE_PLC_FDEL / LFE_PLC_DSF, LFE_PLC_FDEL / LFE_PLC_DSF, LFE_PLC_FS, rec_frame_us, 48000, mem, 0 );
312 660 : modify_Fs( rec_frame, LFE_PLC_RECLEN, LFE_PLC_FS, rec_frame_us, 48000, mem, 0 );
313 :
314 : /*samples are generated with 48k sampling rate
315 : and then converted to required sampling rate by simple decimation
316 : as signal is already bandlimited*/
317 :
318 : /*decimation to correct sampling rate*/
319 660 : if ( output_Fs != 48000 )
320 : {
321 13932 : for ( i = 0; i < rec_frame_len; i++ )
322 : {
323 13920 : rec_frame_us[i] = rec_frame_us[( i * 48000 ) / output_Fs];
324 : }
325 : }
326 :
327 1980 : for ( i = 0; i < 2; i++ )
328 : {
329 1320 : ivas_dct_windowing( fade_len, full_len, dct_len, zero_pad_len, pWindow_coeffs, output_frame, input_tda, rec_frame_us + plc_fdel, rec_frame_us + plc_fdel + fade_len + i * dct_len );
330 1320 : ivas_tda( input_tda, ytda + i * dct_len, dct_len );
331 : }
332 :
333 660 : return;
334 : }
|