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_fx.h"
36 : #include "ivas_cnst.h"
37 : #include <math.h>
38 : #include <assert.h>
39 : #include "wmc_auto.h"
40 :
41 : #include "ivas_prot_fx.h"
42 :
43 :
44 : /*-----------------------------------------------------------------------*
45 : * Local constants
46 : *-----------------------------------------------------------------------*/
47 :
48 : #define IVAS_PCA_SM_FAC 0.0234375f
49 : #define IVAS_PCA_N_ITER_QR 10
50 :
51 : #define IVAS_PCA_SM_FAC_FX 50331648 // Q31
52 : #define ONE_MINUS_IVAS_PCA_SM_FAC_FX 2097152000 // MAX_32 - IVAS_PCA_SM_FAC_FX -> Q31
53 : /*-----------------------------------------------------------------------*
54 : * Local function definitions
55 : *-----------------------------------------------------------------------*/
56 :
57 0 : static void ivas_bitstream_write_int32(
58 : BSTR_ENC_HANDLE hMetaData,
59 : const Word32 val,
60 : const Word16 bits )
61 : {
62 : /* MSBs */
63 0 : push_next_indice( hMetaData, (UWord16) L_shr( val, 16 ), sub( bits, 16 ) );
64 :
65 : /* LSBs */
66 0 : push_next_indice( hMetaData, (UWord16) L_and( val, 0x0000FFFF ), 16 );
67 :
68 0 : return;
69 : }
70 :
71 :
72 4 : static void pca_enc_reset_fx(
73 : PCA_ENC_STATE *hPCA )
74 : {
75 : Word16 i;
76 :
77 : /* reset states for interpolation and multiplexing */
78 4 : eye_matrix_fx( hPCA->prev_eigVec_fx, FOA_CHANNELS, MAX_WORD16 );
79 4 : set16_fx( hPCA->prev_ql_fx, 0, IVAS_PCA_INTERP );
80 4 : hPCA->prev_ql_fx[0] = MAX_WORD16;
81 4 : move16();
82 4 : set16_fx( hPCA->prev_qr_fx, 0, IVAS_PCA_INTERP );
83 4 : hPCA->prev_qr_fx[0] = MAX_WORD16;
84 4 : move16();
85 20 : FOR( i = 0; i < FOA_CHANNELS; i++ )
86 : {
87 16 : hPCA->prev_D_fx[i] = ONE_IN_Q13; // 0.25 in Q15
88 16 : move16();
89 : }
90 4 : eye_matrix_fx( hPCA->mem_eigVec_interp_fx, FOA_CHANNELS, MAX_WORD16 );
91 4 : set32_fx( hPCA->old_r_sm_fx, 0, FOA_CHANNELS * FOA_CHANNELS );
92 4 : hPCA->old_r_sm_q = 31;
93 4 : move16();
94 :
95 4 : return;
96 : }
97 :
98 0 : static void pca_transform_sub_fx(
99 : Word16 *eigVec, // Q15
100 : Word32 *transformed_data[8], /* i : input/transformed audio channels Q11 */
101 : const Word16 start,
102 : const Word16 len,
103 : const Word16 n_channels )
104 : {
105 : Word16 i, j, k;
106 : Word32 temp;
107 : Word32 buffer_data[FOA_CHANNELS];
108 : Word32 temp2;
109 0 : FOR( j = 0; j < len; j++ )
110 : {
111 0 : FOR( k = 0; k < n_channels; k++ )
112 : {
113 0 : buffer_data[k] = transformed_data[k][j + start]; // Q11
114 0 : move32();
115 : }
116 0 : FOR( k = 0; k < n_channels; k++ )
117 : {
118 0 : temp = 0;
119 0 : move32();
120 0 : FOR( i = 0; i < n_channels; i++ )
121 : {
122 0 : temp2 = Mpy_32_16_1( buffer_data[i], eigVec[k * IVAS_PCA_INTERP + i] ); // Q11
123 0 : temp = L_add( temp, temp2 );
124 : }
125 0 : transformed_data[k][( j + start )] = temp; // Q11
126 0 : move32();
127 : }
128 : }
129 :
130 0 : return;
131 : }
132 :
133 0 : static void pca_enc_transform_fx(
134 : PCA_ENC_STATE *hPCA,
135 : Word16 *ql_fx, // Q15
136 : Word16 *qr_fx, // Q15
137 : Word32 *transformed_data_fx[8], // Q11
138 : const Word16 input_frame,
139 : const Word16 n_channels )
140 : {
141 : Word16 time_slot;
142 : Word16 slot_len;
143 : Word16 eigVec_interp_fx[FOA_CHANNELS * FOA_CHANNELS]; /* eigenvectors in current frame Q15*/
144 : Word16 ql_interp_fx[IVAS_PCA_LEN_INTERP_Q], qr_interp_fx[IVAS_PCA_LEN_INTERP_Q];
145 :
146 0 : quat_shortestpath_fx( hPCA->prev_ql_fx, ql_fx, hPCA->prev_qr_fx, qr_fx );
147 0 : pca_interp_preproc_fx( hPCA->prev_ql_fx, hPCA->prev_qr_fx, ql_fx, qr_fx, IVAS_PCA_N_SLOTS, ql_interp_fx, qr_interp_fx );
148 :
149 0 : slot_len = idiv1616( input_frame, IVAS_PCA_N_SLOTS );
150 :
151 0 : FOR( time_slot = 0; time_slot < IVAS_PCA_N_SLOTS; time_slot++ )
152 : {
153 : /* convert from double quaternion to 4D matrix */
154 0 : dquat2mat_fx( &ql_interp_fx[IVAS_PCA_INTERP * time_slot], &qr_interp_fx[IVAS_PCA_INTERP * time_slot], eigVec_interp_fx );
155 0 : pca_transform_sub_fx( eigVec_interp_fx, transformed_data_fx, imult1616( slot_len, time_slot ), slot_len, n_channels );
156 : }
157 :
158 0 : return;
159 : }
160 :
161 1750 : static void pca_update_state_fx(
162 : PCA_ENC_STATE *hPCA,
163 : Word16 *ql, // Q15
164 : Word16 *qr, // Q15
165 : Word16 *eigVec, // Q15
166 : const Word16 n_channels )
167 : {
168 1750 : Copy( qr, hPCA->prev_qr_fx, IVAS_PCA_INTERP ); // Q15
169 1750 : Copy( ql, hPCA->prev_ql_fx, IVAS_PCA_INTERP ); // Q15
170 1750 : Copy( eigVec, hPCA->prev_eigVec_fx, imult1616( n_channels, n_channels ) ); // Q15
171 :
172 1750 : return;
173 : }
174 :
175 1750 : static void swap_eigvec_fx(
176 : Word32 *eigVec, // Q31
177 : const Word32 i,
178 : const Word32 j )
179 : {
180 : Word32 eigVec_tmp[FOA_CHANNELS];
181 : Word32 k;
182 :
183 8750 : FOR( k = 0; k < FOA_CHANNELS; k++ )
184 : {
185 7000 : eigVec_tmp[k] = eigVec[k * FOA_CHANNELS + i];
186 7000 : move32();
187 : }
188 8750 : FOR( k = 0; k < FOA_CHANNELS; k++ )
189 : {
190 7000 : eigVec[k * FOA_CHANNELS + i] = eigVec[k * FOA_CHANNELS + j];
191 7000 : move32();
192 : }
193 8750 : FOR( k = 0; k < FOA_CHANNELS; k++ )
194 : {
195 7000 : eigVec[k * FOA_CHANNELS + j] = eigVec_tmp[k];
196 7000 : move32();
197 : }
198 :
199 1750 : return;
200 : }
201 :
202 1750 : static void sort4_D_eigVec_fx(
203 : Word32 *D, // Q
204 : Word32 *eigVec // Q31
205 : )
206 : {
207 : Word32 tempr;
208 :
209 1750 : IF( LT_32( D[0], D[1] ) )
210 : {
211 0 : SWAP( D[0], D[1] );
212 0 : swap_eigvec_fx( eigVec, 0, 1 );
213 : }
214 :
215 1750 : IF( LT_32( D[2], D[3] ) )
216 : {
217 0 : SWAP( D[2], D[3] );
218 0 : swap_eigvec_fx( eigVec, 2, 3 );
219 : }
220 :
221 1750 : IF( LT_32( D[0], D[2] ) )
222 : {
223 0 : SWAP( D[0], D[2] );
224 0 : swap_eigvec_fx( eigVec, 0, 2 );
225 : }
226 :
227 1750 : IF( LT_32( D[1], D[3] ) )
228 : {
229 0 : SWAP( D[1], D[3] );
230 0 : swap_eigvec_fx( eigVec, 1, 3 );
231 : }
232 :
233 1750 : IF( LT_32( D[1], D[2] ) )
234 : {
235 0 : SWAP( D[1], D[2] );
236 0 : swap_eigvec_fx( eigVec, 1, 2 );
237 : }
238 :
239 : /* swap last 2 values */
240 1750 : SWAP( D[2], D[3] );
241 1750 : swap_eigvec_fx( eigVec, 2, 3 );
242 :
243 1750 : return;
244 : }
245 :
246 : /*-------------------------------------------------------------------------
247 : * ivas_pca_enc_init()
248 : *
249 : * Initialize PCA encoder
250 : *------------------------------------------------------------------------*/
251 :
252 4 : void ivas_pca_enc_init_fx(
253 : PCA_ENC_STATE *hPCA /* i/o: PCA encoder structure */
254 : )
255 : {
256 4 : hPCA->prev_bypass_decision = PCA_MODE_INACTIVE;
257 4 : move32();
258 4 : pca_enc_reset_fx( hPCA );
259 :
260 4 : return;
261 : }
262 :
263 : /*-------------------------------------------------------------------------
264 : * ivas_pca_enc()
265 : *
266 : * PCA encoder
267 : *------------------------------------------------------------------------*/
268 :
269 1750 : void ivas_pca_enc_fx(
270 : const ENCODER_CONFIG_HANDLE hEncoderConfig, /* i : configuration structure */
271 : PCA_ENC_STATE *hPCA, /* i : PCA encoder structure */
272 : BSTR_ENC_HANDLE hMetaData, /* i/o: MetaData handle */
273 : Word32 *data_fx[8], // Q11
274 : const Word16 input_frame, /* i : input frame length */
275 : const Word16 n_channels /* i : number of channels */
276 : )
277 : {
278 : Word32 eigVec_tmp_fx[FOA_CHANNELS * FOA_CHANNELS]; /* transpose / tmp matrix for permutation */
279 : Word32 min_dot_fx, min_dot2_fx;
280 : Word16 i, j, k, l;
281 : Word16 fac_cost_fx;
282 : Word16 len_subfr;
283 : Word16 path[FOA_CHANNELS];
284 : Word32 r_fx[FOA_CHANNELS * FOA_CHANNELS]; /* covariance matrix */
285 : Word64 r_fx_64[FOA_CHANNELS * FOA_CHANNELS]; /* covariance matrix */
286 : Word32 D_fx[FOA_CHANNELS];
287 : Word32 eigVec_fx[FOA_CHANNELS * FOA_CHANNELS]; /* eigenvectors in current frame */
288 : Word32 index[2];
289 : Word32 *ptr_sig_fx[FOA_CHANNELS];
290 : Word16 temp_fx;
291 : Word32 temp_fx32;
292 : Word16 D_tmp_fx[FOA_CHANNELS];
293 : Word32 fac_D_fx, sum_D_fx;
294 : Word32 ivas_total_brate;
295 : Word32 alpha_fx;
296 : Word32 one_minus_alpha_fx;
297 : Word32 r_sm_fx[16];
298 : Word16 bypass_decision;
299 : Word16 det_fx, temp_e;
300 : Word64 W_tmp;
301 : Word16 cost_mtx_fx[FOA_CHANNELS * FOA_CHANNELS]; /* cost matrix */
302 : Word16 ql_fx[IVAS_PCA_INTERP], qr_fx[IVAS_PCA_INTERP];
303 : Word32 L_tmp, L_tmp1;
304 : Word16 eigVec_fx16[FOA_CHANNELS * FOA_CHANNELS];
305 : Word32 dist_alt_fx32, dotl_fx, dotr_fx;
306 : Word16 q;
307 :
308 1750 : ivas_total_brate = hEncoderConfig->ivas_total_brate;
309 1750 : move32();
310 : /* if PCA is disabled, just pass-through */
311 1750 : IF( hEncoderConfig->Opt_PCA_ON == 0 )
312 : {
313 : /* write by-pass indicator */
314 0 : push_next_indice( hMetaData, PCA_MODE_INACTIVE, 1 );
315 :
316 0 : return;
317 : }
318 :
319 : /* handle bitrate switching */
320 1750 : IF( NE_32( ivas_total_brate, PCA_BRATE ) )
321 : {
322 0 : pca_enc_reset_fx( hPCA );
323 :
324 0 : IF( NE_32( hEncoderConfig->last_ivas_total_brate, PCA_BRATE ) )
325 : {
326 0 : eye_matrix_fx( hPCA->mem_eigVec_interp_fx, FOA_CHANNELS, MAX_16 );
327 :
328 : /* copy input data into output directly as previous frame was already in by-pass mode */
329 0 : FOR( k = 0; k < n_channels; k++ )
330 : {
331 : // ToDo: TBV
332 : }
333 : }
334 : ELSE
335 : {
336 : /* set PCA by-pass mode in current frame and interpolate transform as previous frame used PCA */
337 0 : pca_enc_transform_fx( hPCA, ql_fx, qr_fx, data_fx, input_frame, n_channels ); // output Q9
338 : }
339 :
340 0 : pca_update_state_fx( hPCA, ql_fx, qr_fx, eigVec_fx16, n_channels );
341 0 : hPCA->prev_bypass_decision = PCA_MODE_INACTIVE;
342 0 : move16();
343 :
344 0 : return;
345 : }
346 :
347 : /*-----------------------------------------------------------------*
348 : * Covariance
349 : *-----------------------------------------------------------------*/
350 :
351 1750 : len_subfr = shr( input_frame, 1 );
352 :
353 5250 : FOR( i = 0; i < input_frame; i += len_subfr )
354 : {
355 17500 : FOR( k = 0; k < FOA_CHANNELS; k++ )
356 : {
357 14000 : ptr_sig_fx[k] = &data_fx[k][i]; // Q11
358 : }
359 :
360 :
361 3500 : cov_subfr_fx( ptr_sig_fx, r_fx_64, n_channels, shr( input_frame, 1 ) );
362 : Word16 tmp_q[FOA_CHANNELS * FOA_CHANNELS];
363 59500 : FOR( Word16 idx = 0; idx < FOA_CHANNELS * FOA_CHANNELS; idx++ )
364 : {
365 56000 : Word16 w_norm = W_norm( r_fx_64[idx] );
366 56000 : W_tmp = W_shl( r_fx_64[idx], w_norm );
367 56000 : L_tmp = W_extract_h( W_tmp );
368 56000 : r_fx[idx] = L_tmp;
369 56000 : tmp_q[idx] = add( sub( w_norm, 32 ), Q11 );
370 56000 : move32();
371 56000 : move16();
372 : }
373 :
374 :
375 3500 : q = tmp_q[0];
376 3500 : move16();
377 56000 : FOR( Word16 idx = 1; idx < FOA_CHANNELS * FOA_CHANNELS; idx++ )
378 : {
379 52500 : if ( GT_16( q, tmp_q[idx] ) )
380 : {
381 184 : q = tmp_q[idx];
382 184 : move16();
383 : }
384 : }
385 :
386 :
387 59500 : FOR( Word16 idx = 0; idx < FOA_CHANNELS * FOA_CHANNELS; idx++ )
388 : {
389 56000 : r_fx[idx] = L_shr( r_fx[idx], ( sub( tmp_q[idx], q ) ) ); // q
390 56000 : move32();
391 : }
392 :
393 3500 : alpha_fx = IVAS_PCA_SM_FAC_FX; // Q31
394 3500 : one_minus_alpha_fx = ONE_MINUS_IVAS_PCA_SM_FAC_FX; // Q31
395 3500 : move32();
396 3500 : move32();
397 :
398 3500 : Word16 min_q = s_min( q, hPCA->old_r_sm_q );
399 :
400 59500 : FOR( k = 0; k < FOA_CHANNELS * FOA_CHANNELS; k++ )
401 : {
402 56000 : L_tmp = L_shr( Mpy_32_32( alpha_fx, r_fx[k] ), sub( q, min_q ) ); // min_q
403 56000 : L_tmp1 = L_shr( Mpy_32_32( one_minus_alpha_fx, hPCA->old_r_sm_fx[k] ), sub( hPCA->old_r_sm_q, min_q ) ); // min_q
404 56000 : r_sm_fx[k] = L_add( L_tmp, L_tmp1 ); // min_q
405 56000 : hPCA->old_r_sm_fx[k] = r_sm_fx[k]; // min_q
406 56000 : move32();
407 56000 : move32();
408 : }
409 3500 : hPCA->old_r_sm_q = min_q;
410 3500 : move16();
411 : }
412 :
413 : /* conditioning */
414 8750 : FOR( k = 0; k < FOA_CHANNELS; k++ )
415 : {
416 7000 : temp_fx32 = r_sm_fx[k * FOA_CHANNELS + k]; // min_q
417 7000 : move32();
418 7000 : IF( LT_32( temp_fx32, L_shr( 64424, sub( 31, hPCA->old_r_sm_q ) ) ) ) // IVAS_PCA_COV_THRES in Q31
419 : {
420 0 : temp_fx32 = L_shr( 64424, sub( 31, hPCA->old_r_sm_q ) ); // IVAS_PCA_COV_THRES in Q31 /*hPCA->old_r_sm_q */
421 : }
422 7000 : r_sm_fx[k * FOA_CHANNELS + k] = temp_fx32; /* pointer reuse */ // hPCA->old_r_sm_q
423 7000 : move32();
424 : }
425 :
426 :
427 : /*-----------------------------------------------------------------*
428 : * Eigenvalue decomposition
429 : *-----------------------------------------------------------------*/
430 1750 : eig_qr_fx( r_sm_fx, hPCA->old_r_sm_q, IVAS_PCA_N_ITER_QR, eigVec_fx, D_fx, n_channels );
431 :
432 : /* force positive eigenvalues */
433 1750 : sum_D_fx = 0;
434 1750 : Word16 sum_e = 0;
435 1750 : move32();
436 1750 : move16();
437 8750 : FOR( k = 0; k < n_channels; k++ )
438 : {
439 7000 : IF( D_fx[k] < 0 )
440 : {
441 0 : D_fx[k] = L_negate( D_fx[k] );
442 0 : move32();
443 :
444 0 : FOR( l = 0; l < n_channels; l++ )
445 : {
446 0 : eigVec_fx[l * n_channels + k] = L_negate( eigVec_fx[l * n_channels + k] );
447 0 : move32();
448 : }
449 : }
450 7000 : sum_D_fx = BASOP_Util_Add_Mant32Exp( sum_D_fx, sum_e, D_fx[k], sub( 31, hPCA->old_r_sm_q ), &sum_e );
451 : }
452 :
453 : /* normalize */
454 : Word16 tmp, exp;
455 1750 : tmp = BASOP_Util_Divide3232_Scale( 1, L_add( sum_D_fx, EPSILON_FX ), &exp );
456 1750 : exp = add( exp, sub( 0, sum_e ) );
457 1750 : IF( sum_D_fx == 0 )
458 : {
459 0 : fac_D_fx = MAX_32;
460 0 : move32();
461 : }
462 : ELSE
463 : {
464 1750 : fac_D_fx = L_shl( L_deposit_h( tmp ), exp ); // Q31
465 : }
466 :
467 :
468 8750 : FOR( k = 0; k < n_channels; k++ )
469 : {
470 7000 : D_fx[k] = Mpy_32_32( fac_D_fx, D_fx[k] );
471 7000 : move32();
472 : }
473 :
474 :
475 : /* sorting (sanity check) with SPAR inversion of last two channels */
476 1750 : sort4_D_eigVec_fx( D_fx, eigVec_fx );
477 :
478 : /* normalize and compute amount of decorrelation */
479 1750 : dist_alt_fx32 = 0;
480 1750 : Word16 dist_alt_e16 = 0;
481 1750 : move32();
482 1750 : move16();
483 :
484 8750 : FOR( k = 0; k < FOA_CHANNELS; k++ )
485 : {
486 7000 : dist_alt_fx32 = BASOP_Util_Add_Mant32Exp( dist_alt_fx32, dist_alt_e16, BASOP_Util_Loge( Mpy_32_32( r_fx[k * FOA_CHANNELS + k], fac_D_fx ), sub( 31, hPCA->old_r_sm_q ) ), 6, &dist_alt_e16 ); // Q25
487 7000 : dist_alt_fx32 = BASOP_Util_Add_Mant32Exp( dist_alt_fx32, dist_alt_e16, L_negate( BASOP_Util_Loge( D_fx[k], sub( 31, hPCA->old_r_sm_q ) ) ), 6, &dist_alt_e16 ); // Q25
488 : }
489 :
490 : /*-----------------------------------------------------------------*
491 : * Eigenvector alignment
492 : *-----------------------------------------------------------------*/
493 :
494 1750 : IF( EQ_16( hPCA->prev_bypass_decision, PCA_MODE_ACTIVE ) )
495 : {
496 : /* compute absolute cost matrix */
497 0 : FOR( k = 0; k < n_channels; k++ ) /* column */
498 : {
499 0 : FOR( l = 0; l < n_channels; l++ ) /* row */
500 : {
501 0 : fac_cost_fx = mult( hPCA->prev_D_fx[l], extract_l( D_fx[k] ) ); // hPCA->old_r_sm_q
502 0 : temp_fx = 0;
503 0 : temp_e = 0;
504 0 : move16();
505 0 : move16();
506 0 : FOR( i = 0; i < n_channels; i++ ) /* row */
507 : {
508 0 : temp_e = BASOP_Util_Add_MantExp( temp_fx, temp_e, mult( hPCA->prev_eigVec_fx[i * n_channels + l], extract_l( eigVec_fx[i * n_channels + k] ) ), 0, &temp_fx );
509 : }
510 0 : IF( temp_fx < 0 )
511 : {
512 0 : temp_fx = negate( temp_fx );
513 : }
514 0 : cost_mtx_fx[k * n_channels + l] = shr( mult( temp_fx, fac_cost_fx ), temp_e ); // hPCA->old_r_sm_q
515 0 : move16();
516 : }
517 : }
518 :
519 : /* find optimal permutation */
520 0 : exhst_4x4_fx( cost_mtx_fx, path, 1 );
521 : }
522 : ELSE
523 : {
524 : /* no alignment needed if previous PCA is inactive */
525 8750 : FOR( i = 0; i < 4; i++ )
526 : {
527 7000 : path[i] = i;
528 7000 : move16();
529 : }
530 : }
531 :
532 : /* permute eigenvectors */
533 8750 : FOR( k = 0; k < n_channels; k++ )
534 : {
535 7000 : j = path[k];
536 7000 : move16();
537 : /* copy j-th column to column k */
538 35000 : FOR( l = 0; l < n_channels; l++ )
539 : {
540 28000 : eigVec_tmp_fx[l * n_channels + k] = eigVec_fx[l * n_channels + j]; // Q31
541 28000 : move16();
542 : }
543 7000 : D_tmp_fx[k] = extract_l( D_fx[j] );
544 7000 : move16();
545 : }
546 :
547 29750 : FOR( k = 0; k < n_channels * n_channels; k++ )
548 : {
549 28000 : eigVec_fx[k] = eigVec_tmp_fx[k]; // Q31
550 28000 : move32();
551 : }
552 :
553 : /* check for sign inversions */
554 8750 : FOR( k = 0; k < n_channels; k++ ) /* column */
555 : {
556 7000 : temp_fx = 0;
557 7000 : temp_e = 0;
558 : // Word64 W_add_fx = 1;
559 7000 : move16();
560 7000 : move16();
561 35000 : FOR( i = 0; i < n_channels; i++ ) /* row */
562 : {
563 28000 : temp_e = BASOP_Util_Add_MantExp( temp_fx, temp_e, mult( hPCA->prev_eigVec_fx[i * n_channels + k], extract_l( eigVec_fx[i * n_channels + k] ) ), 0, &temp_fx );
564 : }
565 7000 : move32();
566 7000 : move16();
567 :
568 7000 : IF( temp_fx < 0 )
569 : {
570 17905 : FOR( i = 0; i < n_channels; i++ )
571 : {
572 14324 : eigVec_fx[i * n_channels + k] = L_negate( eigVec_fx[i * n_channels + k] ); // Q31
573 14324 : move32();
574 : }
575 : }
576 : }
577 :
578 : /* force rotation matrix(det = +1) */
579 1750 : Copy_Scale_sig32_16( eigVec_fx, eigVec_fx16, FOA_CHANNELS * FOA_CHANNELS, -16 ); // Q15
580 1750 : det_fx = mat_det4_fx( eigVec_fx16 );
581 1750 : IF( det_fx < 0 )
582 : {
583 0 : swap_eigvec_fx( eigVec_fx, 2, 3 );
584 0 : temp_fx = D_tmp_fx[3];
585 0 : D_tmp_fx[3] = D_tmp_fx[2];
586 0 : D_tmp_fx[2] = temp_fx;
587 0 : move16();
588 0 : move16();
589 0 : move16();
590 : }
591 :
592 : /* update state */
593 8750 : FOR( k = 0; k < n_channels; k++ )
594 : {
595 7000 : hPCA->prev_D_fx[k] = D_tmp_fx[k]; // Q31
596 7000 : move16();
597 : }
598 :
599 : /*-----------------------------------------------------------------*
600 : * Rotation matrix parametrization and quantization
601 : *-----------------------------------------------------------------*/
602 :
603 : /* convert frrm rotation matrix to double quaternion */
604 1750 : Copy_Scale_sig_32_16( eigVec_fx, eigVec_fx16, FOA_CHANNELS * FOA_CHANNELS, -16 ); // Q15
605 1750 : mat2dquat_fx( eigVec_fx16, ql_fx, qr_fx );
606 :
607 1750 : dotl_fx = Dot_product( hPCA->prev_ql_fx, ql_fx, 4 ); // Q31
608 1750 : dotr_fx = Dot_product( hPCA->prev_qr_fx, qr_fx, 4 ); // Q31
609 1750 : IF( LT_32( dotl_fx, dotr_fx ) )
610 : {
611 540 : min_dot_fx = dotl_fx; // Q31
612 540 : move32();
613 : }
614 : ELSE
615 : {
616 1210 : min_dot_fx = dotr_fx; // Q31
617 1210 : move32();
618 : }
619 :
620 1750 : IF( LT_16( ql_fx[0], qr_fx[0] ) )
621 : {
622 540 : min_dot2_fx = L_deposit_h( ql_fx[0] ); // Q31
623 : }
624 : ELSE
625 : {
626 1210 : min_dot2_fx = L_deposit_h( qr_fx[0] ); // Q31
627 : }
628 :
629 1750 : bypass_decision = PCA_MODE_ACTIVE;
630 1750 : move16();
631 1750 : if ( LT_32( L_shr( dist_alt_fx32, sub( 31, dist_alt_e16 ) ), IVAS_PCA_THRES_DIST_ALT_FX ) )
632 : {
633 1750 : bypass_decision = PCA_MODE_INACTIVE;
634 1750 : move16();
635 : }
636 1750 : if ( LT_32( min_dot_fx, IVAS_PCA_THRES_MIN_DOT_FX ) )
637 : {
638 1720 : bypass_decision = PCA_MODE_INACTIVE;
639 1720 : move16();
640 : }
641 1750 : if ( LT_32( min_dot2_fx, IVAS_PCA_THRES_MIN_DOT2_FX ) )
642 : {
643 1133 : bypass_decision = PCA_MODE_INACTIVE;
644 1133 : move16();
645 : }
646 :
647 :
648 : /* if PCA is inactive */
649 1750 : IF( EQ_16( bypass_decision, PCA_MODE_INACTIVE ) )
650 : {
651 1750 : eye_matrix_fx( eigVec_fx16, 4, MAX_16 );
652 1750 : set16_fx( ql_fx, 0, 4 );
653 1750 : set16_fx( qr_fx, 0, 4 );
654 1750 : ql_fx[0] = MAX_16;
655 1750 : qr_fx[0] = MAX_16;
656 1750 : move16();
657 1750 : move16();
658 :
659 : /* write by-pass indicator */
660 1750 : push_next_indice( hMetaData, PCA_MODE_INACTIVE, 1 );
661 1750 : IF( EQ_16( hPCA->prev_bypass_decision, PCA_MODE_INACTIVE ) )
662 : {
663 1750 : eye_matrix_fx( hPCA->mem_eigVec_interp_fx, FOA_CHANNELS, MAX_16 );
664 : }
665 : ELSE
666 : {
667 : /* set PCA by-pass mode in current frame and interpolate transform as previous frame used PCA */
668 0 : pca_enc_transform_fx( hPCA, ql_fx, qr_fx, data_fx, input_frame, n_channels );
669 : }
670 :
671 1750 : pca_update_state_fx( hPCA, ql_fx, qr_fx, eigVec_fx16, n_channels );
672 1750 : hPCA->prev_bypass_decision = PCA_MODE_INACTIVE;
673 1750 : move16();
674 1750 : return;
675 : }
676 :
677 : /* force ql to have first component with positive sign */
678 :
679 0 : IF( ql_fx[0] < 0 )
680 : {
681 0 : FOR( i = 0; i < 4; i++ )
682 : {
683 0 : ql_fx[i] = negate( ql_fx[i] ); // Q15
684 0 : qr_fx[i] = negate( qr_fx[i] );
685 0 : move16();
686 0 : move16();
687 : }
688 : }
689 :
690 : /* quantize double quaternion */
691 0 : pca_enc_s3_fx( ql_fx, &index[0] );
692 0 : pca_enc_s3_fx( qr_fx, &index[1] );
693 :
694 :
695 : /* write bypass flag to bitstream to indicate active mode */
696 0 : push_next_indice( hMetaData, PCA_MODE_ACTIVE, 1 );
697 :
698 0 : ivas_bitstream_write_int32( hMetaData, index[0], IVAS_PCA_QBITS - 1 );
699 0 : ivas_bitstream_write_int32( hMetaData, index[1], IVAS_PCA_QBITS );
700 :
701 : /* transform */
702 :
703 0 : pca_enc_transform_fx( hPCA, ql_fx, qr_fx, data_fx, input_frame, n_channels );
704 :
705 : /*-----------------------------------------------------------------*
706 : * update state for next frame
707 : *-----------------------------------------------------------------*/
708 :
709 0 : hPCA->prev_bypass_decision = bypass_decision;
710 0 : move16();
711 0 : pca_update_state_fx( hPCA, ql_fx, qr_fx, eigVec_fx16, n_channels );
712 0 : return;
713 : }
|
