123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367 |
- /*********************************************************************************
- ** ITU-T G.722.1 (2005-05) - Fixed point implementation for main body and Annex C
- ** > Software Release 2.1 (2008-06)
- ** (Simple repackaging; no change from 2005-05 Release 2.0 code)
- **
- ** © 2004 Polycom, Inc.
- **
- ** All rights reserved.
- **
- *********************************************************************************/
- /*********************************************************************************
- * Filename: dct_type_iv_a.c
- *
- * Purpose: Discrete Cosine Transform, Type IV used for MLT
- *
- * The basis functions are
- *
- * cos(PI*(t+0.5)*(k+0.5)/block_length)
- *
- * for time t and basis function number k. Due to the symmetry of the expression
- * in t and k, it is clear that the forward and inverse transforms are the same.
- *
- *********************************************************************************/
- /*********************************************************************************
- Include files
- *********************************************************************************/
- #include "defs.h"
- #include "count.h"
- #include "dct4_a.h"
- /*********************************************************************************
- External variable declarations
- *********************************************************************************/
- extern Word16 anal_bias[DCT_LENGTH];
- extern Word16 dct_core_a[DCT_LENGTH_DIV_32][DCT_LENGTH_DIV_32];
- extern cos_msin_t a_cos_msin_2 [DCT_LENGTH_DIV_32];
- extern cos_msin_t a_cos_msin_4 [DCT_LENGTH_DIV_16];
- extern cos_msin_t a_cos_msin_8 [DCT_LENGTH_DIV_8];
- extern cos_msin_t a_cos_msin_16[DCT_LENGTH_DIV_4];
- extern cos_msin_t a_cos_msin_32[DCT_LENGTH_DIV_2];
- extern cos_msin_t a_cos_msin_64[DCT_LENGTH];
- extern cos_msin_t *a_cos_msin_table[];
- /*********************************************************************************
- Function: dct_type_iv_a
- Syntax: void dct_type_iv_a (input, output, dct_length)
- Word16 input[], output[], dct_length;
- Description: Discrete Cosine Transform, Type IV used for MLT
- Design Notes:
-
- WMOPS: | 24kbit | 32kbit
- -------|--------------|----------------
- AVG | 1.14 | 1.14
- -------|--------------|----------------
- MAX | 1.14 | 1.14
- -------|--------------|----------------
-
- 14kHz | 24kbit | 32kbit | 48kbit
- -------|--------------|----------------|----------------
- AVG | 2.57 | 2.57 | 2.57
- -------|--------------|----------------|----------------
- MAX | 2.57 | 2.57 | 2.57
- -------|--------------|----------------|----------------
- *********************************************************************************/
- void dct_type_iv_a (Word16 *input,Word16 *output,Word16 dct_length)
- {
- Word16 buffer_a[MAX_DCT_LENGTH], buffer_b[MAX_DCT_LENGTH], buffer_c[MAX_DCT_LENGTH];
- Word16 *in_ptr, *in_ptr_low, *in_ptr_high, *next_in_base;
- Word16 *out_ptr_low, *out_ptr_high, *next_out_base;
- Word16 *out_buffer, *in_buffer, *buffer_swap;
- Word16 in_val_low, in_val_high;
- Word16 out_val_low, out_val_high;
- Word16 in_low_even, in_low_odd;
- Word16 in_high_even, in_high_odd;
- Word16 out_low_even, out_low_odd;
- Word16 out_high_even, out_high_odd;
- Word16 *pair_ptr;
- Word16 cos_even, cos_odd, msin_even, msin_odd;
- Word16 neg_cos_odd;
- Word16 neg_msin_even;
- Word32 sum;
- Word16 set_span, set_count, set_count_log, pairs_left, sets_left;
- Word16 i,k;
- Word16 index;
- cos_msin_t **table_ptr_ptr, *cos_msin_ptr;
-
- Word16 temp;
- Word32 acca;
- Word16 dct_length_log;
- /*++++++++++++++++++++++++++++++++++++++++++++++++++++++*/
- /* Do the sum/difference butterflies, the first part of */
- /* converting one N-point transform into N/2 two-point */
- /* transforms, where N = 1 << DCT_LENGTH_LOG. = 64/128 */
- /*++++++++++++++++++++++++++++++++++++++++++++++++++++++*/
- test();
- if (dct_length==DCT_LENGTH)
- {
- dct_length_log = DCT_LENGTH_LOG;
- /* Add bias offsets */
- for (i=0;i<dct_length;i++)
- {
- input[i] = add(input[i],anal_bias[i]);
- move16();
- }
- }
- else
- dct_length_log = MAX_DCT_LENGTH_LOG;
- index = 0L;
- move16();
- in_buffer = input;
- move16();
- out_buffer = buffer_a;
- move16();
- temp = sub(dct_length_log,2);
- for (set_count_log=0;set_count_log<=temp;set_count_log++)
- {
- /*===========================================================*/
- /* Initialization for the loop over sets at the current size */
- /*===========================================================*/
- /* set_span = 1 << (DCT_LENGTH_LOG - set_count_log); */
- set_span = shr_nocheck(dct_length,set_count_log);
- set_count = shl_nocheck(1,set_count_log);
- in_ptr = in_buffer;
- move16();
- next_out_base = out_buffer;
- move16();
- /*=====================================*/
- /* Loop over all the sets of this size */
- /*=====================================*/
- for (sets_left=set_count;sets_left>0;sets_left--)
- {
- /*||||||||||||||||||||||||||||||||||||||||||||*/
- /* Set up output pointers for the current set */
- /*||||||||||||||||||||||||||||||||||||||||||||*/
- out_ptr_low = next_out_base;
- next_out_base = next_out_base + set_span;
- out_ptr_high = next_out_base;
- /*||||||||||||||||||||||||||||||||||||||||||||||||||*/
- /* Loop over all the butterflies in the current set */
- /*||||||||||||||||||||||||||||||||||||||||||||||||||*/
- do
- {
- in_val_low = *in_ptr++;
- in_val_high = *in_ptr++;
- // blp: addition of two 16bits vars, there's no way
- // they'll overflow a 32bit var
- //acca = L_add(in_val_low,in_val_high);
- acca = (in_val_low + in_val_high);
- acca = L_shr_nocheck(acca,1);
- out_val_low = extract_l(acca);
- acca = L_sub(in_val_low,in_val_high);
- acca = L_shr_nocheck(acca,1);
- out_val_high = extract_l(acca);
- *out_ptr_low++ = out_val_low;
- *--out_ptr_high = out_val_high;
- test();
- } while (out_ptr_low < out_ptr_high);
- } /* End of loop over sets of the current size */
- /*============================================================*/
- /* Decide which buffers to use as input and output next time. */
- /* Except for the first time (when the input buffer is the */
- /* subroutine input) we just alternate the local buffers. */
- /*============================================================*/
- in_buffer = out_buffer;
- move16();
- if (out_buffer == buffer_a)
- out_buffer = buffer_b;
- else
- out_buffer = buffer_a;
- index = add(index,1);
- } /* End of loop over set sizes */
- /*++++++++++++++++++++++++++++++++*/
- /* Do N/2 two-point transforms, */
- /* where N = 1 << DCT_LENGTH_LOG */
- /*++++++++++++++++++++++++++++++++*/
- pair_ptr = in_buffer;
- move16();
- buffer_swap = buffer_c;
- move16();
- temp = sub(dct_length_log,1);
- temp = shl_nocheck(1,temp);
- for (pairs_left=temp; pairs_left > 0; pairs_left--)
- {
- for ( k=0; k<CORE_SIZE; k++ )
- {
- #if PJ_HAS_INT64
- /* blp: danger danger! not really compatible but faster */
- pj_int64_t sum64=0;
- move32();
-
- for ( i=0; i<CORE_SIZE; i++ )
- {
- sum64 += L_mult(pair_ptr[i], dct_core_a[i][k]);
- }
- sum = L_saturate(sum64);
- #else
- sum=0L;
- move32();
- for ( i=0; i<CORE_SIZE; i++ )
- {
- sum = L_mac(sum, pair_ptr[i],dct_core_a[i][k]);
- }
- #endif
- buffer_swap[k] = itu_round(sum);
- }
- /* address arithmetic */
- pair_ptr += CORE_SIZE;
- buffer_swap += CORE_SIZE;
- }
- for (i=0;i<dct_length;i++)
- {
- in_buffer[i] = buffer_c[i];
- move16();
- }
-
- table_ptr_ptr = a_cos_msin_table;
- /*++++++++++++++++++++++++++++++*/
- /* Perform rotation butterflies */
- /*++++++++++++++++++++++++++++++*/
- temp = sub(dct_length_log,2);
- for (set_count_log = temp; set_count_log >= 0; set_count_log--)
- {
- /*===========================================================*/
- /* Initialization for the loop over sets at the current size */
- /*===========================================================*/
- /* set_span = 1 << (DCT_LENGTH_LOG - set_count_log); */
- set_span = shr_nocheck(dct_length,set_count_log);
- set_count = shl_nocheck(1,set_count_log);
- next_in_base = in_buffer;
- move16();
- test();
- if (set_count_log == 0)
- {
- next_out_base = output;
- }
- else
- {
- next_out_base = out_buffer;
- }
- /*=====================================*/
- /* Loop over all the sets of this size */
- /*=====================================*/
- for (sets_left = set_count; sets_left > 0;sets_left--)
- {
- /*|||||||||||||||||||||||||||||||||||||||||*/
- /* Set up the pointers for the current set */
- /*|||||||||||||||||||||||||||||||||||||||||*/
- in_ptr_low = next_in_base;
- move16();
- temp = shr_nocheck(set_span,1);
- /* address arithmetic */
- in_ptr_high = in_ptr_low + temp;
- next_in_base += set_span;
- out_ptr_low = next_out_base;
- next_out_base += set_span;
- out_ptr_high = next_out_base;
- cos_msin_ptr = *table_ptr_ptr;
- /*||||||||||||||||||||||||||||||||||||||||||||||||||||||*/
- /* Loop over all the butterfly pairs in the current set */
- /*||||||||||||||||||||||||||||||||||||||||||||||||||||||*/
- do
- {
- /* address arithmetic */
- in_low_even = *in_ptr_low++;
- in_low_odd = *in_ptr_low++;
- in_high_even = *in_ptr_high++;
- in_high_odd = *in_ptr_high++;
- cos_even = cos_msin_ptr[0].cosine;
- move16();
- msin_even = cos_msin_ptr[0].minus_sine;
- move16();
- cos_odd = cos_msin_ptr[1].cosine;
- move16();
- msin_odd = cos_msin_ptr[1].minus_sine;
- move16();
- cos_msin_ptr += 2;
- sum = 0L;
- sum=L_mac(sum,cos_even,in_low_even);
- neg_msin_even = negate(msin_even);
- sum=L_mac(sum,neg_msin_even,in_high_even);
- out_low_even = itu_round(sum);
- sum = 0L;
- sum=L_mac(sum,msin_even,in_low_even);
- sum=L_mac(sum,cos_even,in_high_even);
- out_high_even= itu_round(sum);
- sum = 0L;
- sum=L_mac(sum,cos_odd,in_low_odd);
- sum=L_mac(sum,msin_odd,in_high_odd);
- out_low_odd= itu_round(sum);
- sum = 0L;
- sum=L_mac(sum,msin_odd,in_low_odd);
- neg_cos_odd = negate(cos_odd);
- sum=L_mac(sum,neg_cos_odd,in_high_odd);
- out_high_odd= itu_round(sum);
- *out_ptr_low++ = out_low_even;
- *--out_ptr_high = out_high_even;
- *out_ptr_low++ = out_low_odd;
- *--out_ptr_high = out_high_odd;
- test();
- } while (out_ptr_low < out_ptr_high);
- } /* End of loop over sets of the current size */
- /*=============================================*/
- /* Swap input and output buffers for next time */
- /*=============================================*/
- buffer_swap = in_buffer;
- in_buffer = out_buffer;
- out_buffer = buffer_swap;
- table_ptr_ptr++;
- }
- }
|