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 10946 : 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 10946 : set_f( a, 0.f, LFE_PLC_LPCORD + 1 );
78 10946 : set_f( a_out, 0.f, LFE_PLC_LPCORD + 1 );
79 10946 : rc = &buf[0];
80 10946 : a[0] = 1.f;
81 10946 : err = r[0];
82 10946 : a_out[0] = a[0];
83 10946 : rc[0] = ( -r[1] ) / r[0];
84 10946 : s = r[1];
85 10946 : a[1] = rc[0];
86 10946 : a_out[1] = a[1];
87 :
88 10946 : i = 1;
89 215522 : while ( i < m )
90 : {
91 1168945 : for ( j = 1; j <= i / 2; j++ )
92 : {
93 963817 : l = i - j;
94 963817 : at = a[j] + rc[i - 1] * a[l];
95 963817 : a[l] += rc[i - 1] * a[j];
96 963817 : a[j] = at;
97 : }
98 :
99 205128 : a[i] = rc[i - 1];
100 205128 : err += rc[i - 1] * s;
101 205128 : s = 0.f;
102 205128 : i++;
103 2445696 : for ( j = 0; j < i; j++ )
104 : {
105 2240568 : s += r[i - j] * a[j];
106 : }
107 205128 : rc[i - 1] = ( -s ) / err;
108 205128 : if ( fabsf( rc[i - 1] ) > 0.99945f )
109 : {
110 552 : return 1;
111 : }
112 : else
113 : {
114 4500672 : for ( j = 0; j <= m; j++ )
115 : {
116 4296096 : a_out[j] = a[j];
117 : }
118 : }
119 : }
120 :
121 10394 : return 0;
122 : }
123 :
124 :
125 : /*-----------------------------------------------------------------------------------------*
126 : * Function check_stab()
127 : *
128 : * LPC filter stability check applying given sharpening value delta
129 : *-----------------------------------------------------------------------------------------*/
130 :
131 135649 : 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 135649 : fac = 1.0f + delta;
142 135649 : fac1 = fac;
143 :
144 2848629 : for ( i = 0; i < LFE_PLC_LPCORD; i++ )
145 : {
146 2712980 : amod[i] = a[i] * fac;
147 2712980 : fac *= fac1;
148 : }
149 :
150 135649 : stable = a2rc( amod, refl, LFE_PLC_LPCORD );
151 :
152 135649 : 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 10946 : static float find_max_delta(
163 : float *a )
164 : {
165 : float delta;
166 : float eps;
167 : float fac;
168 : uint16_t stable;
169 :
170 10946 : delta = 0.0f;
171 10946 : eps = 0.01f;
172 10946 : fac = 2;
173 10946 : stable = FALSE;
174 :
175 18582 : while ( check_stab( a, eps ) )
176 : {
177 :
178 7636 : eps *= fac;
179 7636 : stable = TRUE;
180 : }
181 10946 : fac = 0.5f;
182 :
183 10946 : if ( stable )
184 : {
185 4763 : eps *= fac;
186 : }
187 :
188 23458 : while ( !stable )
189 : {
190 12512 : eps *= fac;
191 12512 : stable = check_stab( a, eps );
192 : }
193 :
194 : /* must be stable with current eps */
195 10946 : delta = eps;
196 10946 : eps *= fac;
197 :
198 : while ( 1 )
199 : {
200 104555 : delta += eps;
201 104555 : stable = check_stab( a, delta );
202 :
203 104555 : if ( !stable )
204 : {
205 51276 : if ( fabsf( eps ) > EPS_STOP )
206 : {
207 46542 : eps = -fabsf( eps ) * fac;
208 : }
209 : else
210 : {
211 4734 : eps = -fabsf( eps );
212 : }
213 : }
214 : else
215 : {
216 53279 : if ( fabsf( eps ) < EPS_STOP )
217 : {
218 10946 : break;
219 : }
220 42333 : eps = fabsf( eps ) * fac;
221 : }
222 : }
223 :
224 10946 : return delta;
225 : }
226 :
227 :
228 : /*-----------------------------------------------------------------------------------------*
229 : * Function recover_samples()
230 : *
231 : * recover lost samples by extrapolation of signal buffer
232 : *-----------------------------------------------------------------------------------------*/
233 :
234 19229 : 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 19229 : autocorr( outbuf, r, LFE_PLC_LPCORD, LFE_PLC_BUFLEN, hamm_lfe_plc, 0, 1, 1 );
245 :
246 19229 : if ( r[0] < POW_THR * LFE_PLC_BUFLEN )
247 : {
248 8283 : set_zero( rec_frame, LFE_PLC_RECLEN );
249 8283 : return;
250 : }
251 :
252 10946 : lfeplc_lev_dur( a, r, LFE_PLC_LPCORD );
253 :
254 10946 : delta = find_max_delta( a + 1 );
255 :
256 10946 : fac = 1.0f + delta;
257 10946 : att = 1.0f;
258 :
259 10946 : if ( bfi_count >= LFE_PLC_MUTE_THR )
260 : {
261 0 : att = LFE_PLC_BURST_ATT;
262 0 : fac *= att;
263 : }
264 :
265 229866 : for ( i = 1; i <= LFE_PLC_LPCORD; i++ )
266 : {
267 218920 : a[i] = a[i] * fac;
268 218920 : fac *= att * ( 1.0f + delta );
269 : }
270 :
271 10946 : set_zero( zeroes, LFE_PLC_RECLEN );
272 :
273 10946 : syn_filt( a, LFE_PLC_LPCORD, zeroes, rec_frame, LFE_PLC_RECLEN, outbuf + LFE_PLC_BUFLEN - LFE_PLC_LPCORD, 0 );
274 :
275 10946 : 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 19229 : 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 19229 : output_Fs = output_frame * FRAMES_PER_SEC;
299 19229 : fade_len = hLFE->pWindow_state->fade_len;
300 19229 : full_len = hLFE->pWindow_state->full_len;
301 19229 : dct_len = hLFE->pWindow_state->dct_len;
302 19229 : zero_pad_len = hLFE->pWindow_state->zero_pad_len;
303 19229 : pWindow_coeffs = hLFE->pWindow_state->pWindow_coeffs;
304 :
305 19229 : plc_fdel = (int16_t) ( LFE_PLC_FDEL * output_Fs / 48000 );
306 19229 : rec_frame_len = (int16_t) ( LFE_PLC_RECLEN_48K * output_Fs / 48000 );
307 19229 : recover_samples( hLFE->bfi_count, prevsynth, rec_frame );
308 :
309 19229 : set_zero( mem, 2 * LFE_PLC_FDEL / LFE_PLC_DSF );
310 :
311 19229 : 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 19229 : 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 19229 : if ( output_Fs != 48000 )
320 : {
321 9582029 : for ( i = 0; i < rec_frame_len; i++ )
322 : {
323 9571160 : rec_frame_us[i] = rec_frame_us[( i * 48000 ) / output_Fs];
324 : }
325 : }
326 :
327 57687 : for ( i = 0; i < 2; i++ )
328 : {
329 38458 : 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 38458 : ivas_tda( input_tda, ytda + i * dct_len, dct_len );
331 : }
332 :
333 19229 : return;
334 : }
|