calc_en.cpp

实现3GPP的GSM中AMR语音的CODECS。

/* ------------------------------------------------------------------
 * Copyright (C) 2008 PacketVideo
 *
 * Licensed under the Apache License, Version 2.0 (the "License");
 * you may not use this file except in compliance with the License.
 * You may obtain a copy of the License at
 *
 *      http://www.apache.org/licenses/LICENSE-2.0
 *
 * Unless required by applicable law or agreed to in writing, software
 * distributed under the License is distributed on an "AS IS" BASIS,
 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either
 * express or implied.
 * See the License for the specific language governing permissions
 * and limitations under the License.
 * -------------------------------------------------------------------
 */
/****************************************************************************************
Portions of this file are derived from the following 3GPP standard:

    3GPP TS 26.073
    ANSI-C code for the Adaptive Multi-Rate (AMR) speech codec
    Available from http://www.3gpp.org

(C) 2004, 3GPP Organizational Partners (ARIB, ATIS, CCSA, ETSI, TTA, TTC)
Permission to distribute, modify and use this file under the standard license
terms listed above has been obtained from the copyright holder.
****************************************************************************************/
/*
------------------------------------------------------------------------------



 Pathname: ./audio/gsm-amr/c/src/calc_en.c
 Funtions: calc_unfilt_energies
           calc_filt_energies
           calc_target_energy

     Date: 06/16/2000

------------------------------------------------------------------------------
 REVISION HISTORY

 Description:  After review changes

 Description:  Replaced the inlined code for norm_l and saturation checks
    with the original code to fix a bug with 4.75k encode.

 Description: Synchronized file with UMTS version 3.2.0. Updated coding
              template. Removed unnecessary include files.

 Description: Replaced basic_op.h and oper_32b.h with the header files of the
              math functions used by the file.

 Description: Undid optimization in the loop that calculates various dot
              products in calc_filt_energies. Added calls to L_mac function
              instead of doing a straight multiply-accumulate. Updated
              copyright year. Fixed typecasting issue with TI C compiler.

 Description: Add pOverflow as a passed in variable for the EPOC changes.
              Updated code as per review comments regarding pOverflow.

 Description: Changed round function name to pv_round to avoid conflict with
              round function in C standard library.

 Description:  Replaced "int" and/or "char" with OSCL defined types.

 Description: Using inlines from fxp_arithmetic.h .

 Description: Replacing fxp_arithmetic.h with basic_op.h.

 Description:

------------------------------------------------------------------------------
 MODULE DESCRIPTION

 This file contains the functions that calculate the energy coefficients
 for unfiltered and filtered excitation signals, the LTP coding gain, and
 the target energy.

------------------------------------------------------------------------------
*/


/*----------------------------------------------------------------------------
; INCLUDES
----------------------------------------------------------------------------*/

#include "calc_en.h"
#include "typedef.h"
#include "basicop_malloc.h"
#include "l_comp.h"
#include "cnst.h"
#include "log2.h"
#include "basic_op.h"

/*----------------------------------------------------------------------------
; MACROS
; Define module specific macros here
----------------------------------------------------------------------------*/


/*----------------------------------------------------------------------------
; DEFINES
; Include all pre-processor statements here. Include conditional
; compile variables also.
----------------------------------------------------------------------------*/


/*----------------------------------------------------------------------------
; LOCAL FUNCTION DEFINITIONS
; Function Prototype declaration
----------------------------------------------------------------------------*/

/*----------------------------------------------------------------------------
; LOCAL VARIABLE DEFINITIONS
; Variable declaration - defined here and used outside this module
----------------------------------------------------------------------------*/


/*
------------------------------------------------------------------------------
 FUNCTION NAME: calc_unfilt_energies
------------------------------------------------------------------------------
 INPUT AND OUTPUT DEFINITIONS

 Inputs:
    res      = LP residual, buffer type Word16
    exc      = LTP excitation (unfiltered), buffer type Word16
    code     = CB innovation (unfiltered), buffer type Word16
    gain_pit = pitch gain,  type Word16
    L_subfr  = Subframe length, type Word16
    frac_en  = energy coefficients (4), fraction part, buffer type Word16
    exp_en   = energy coefficients (4), exponent part, buffer type Word16
    ltpg     = LTP coding gain (log2()), pointer to type Word16
    pOverflow= pointer to value indicating existence of overflow (Flag)

 Outputs:
    frac_en buffer containing new fractional parts of energy coefficients
    exp_en buffer containing new exponential parts of energy coefficients
    ltpg points to new LTP coding gain
    pOverflow = 1 if there is an overflow else it is zero.

 Returns:
    None.

 Global Variables Used:
    None

 Local Variables Needed:
    None

------------------------------------------------------------------------------
 FUNCTION DESCRIPTION

 This function calculates several energy coefficients for unfiltered
 excitation signals and the LTP coding gain

    frac_en[0]*2^exp_en[0] = <res res>    LP residual energy
    frac_en[1]*2^exp_en[1] = <exc exc>    LTP residual energy
    frac_en[2]*2^exp_en[2] = <exc code>   LTP/CB innovation dot product
    frac_en[3]*2^exp_en[3] = <lres lres>  LTP residual energy
    (lres = res - gain_pit*exc)
    ltpg = log2(LP_res_en / LTP_res_en)

------------------------------------------------------------------------------
 REQUIREMENTS

  None.

------------------------------------------------------------------------------
 REFERENCES

 calc_en.c, UMTS GSM AMR speech codec, R99 - Version 3.2.0, March 2, 2001

------------------------------------------------------------------------------
 PSEUDO-CODE

void
calc_unfilt_energies(
    Word16 res[],     // i  : LP residual,                               Q0
    Word16 exc[],     // i  : LTP excitation (unfiltered),               Q0
    Word16 code[],    // i  : CB innovation (unfiltered),                Q13
    Word16 gain_pit,  // i  : pitch gain,                                Q14
    Word16 L_subfr,   // i  : Subframe length

    Word16 frac_en[], // o  : energy coefficients (4), fraction part,    Q15
    Word16 exp_en[],  // o  : energy coefficients (4), exponent part,    Q0
    Word16 *ltpg      // o  : LTP coding gain (log2()),                  Q13
)
{
    Word32 s, L_temp;
    Word16 i, exp, tmp;
    Word16 ltp_res_en, pred_gain;
    Word16 ltpg_exp, ltpg_frac;

    // Compute residual energy
    s = L_mac((Word32) 0, res[0], res[0]);
    for (i = 1; i < L_subfr; i++)
        s = L_mac(s, res[i], res[i]);

    // ResEn := 0 if ResEn < 200.0 (= 400 Q1)
    if (L_sub (s, 400L) < 0)
    {
        frac_en[0] = 0;
        exp_en[0] = -15;
    }
    else
    {
        exp = norm_l(s);
        frac_en[0] = extract_h(L_shl(s, exp));
        exp_en[0] = sub(15, exp);
    }

    // Compute ltp excitation energy
    s = L_mac((Word32) 0, exc[0], exc[0]);
    for (i = 1; i < L_subfr; i++)
        s = L_mac(s, exc[i], exc[i]);

    exp = norm_l(s);
    frac_en[1] = extract_h(L_shl(s, exp));
    exp_en[1] = sub(15, exp);

    // Compute scalar product <exc[],code[]>
    s = L_mac((Word32) 0, exc[0], code[0]);
    for (i = 1; i < L_subfr; i++)
        s = L_mac(s, exc[i], code[i]);

    exp = norm_l(s);
    frac_en[2] = extract_h(L_shl(s, exp));
    exp_en[2] = sub(16-14, exp);

    // Compute energy of LTP residual
    s = 0L;
    for (i = 0; i < L_subfr; i++)
    {
        L_temp = L_mult(exc[i], gain_pit);
        L_temp = L_shl(L_temp, 1);
        tmp = sub(res[i], pv_round(L_temp)); // LTP residual, Q0
        s = L_mac (s, tmp, tmp);
    }

    exp = norm_l(s);
    ltp_res_en = extract_h (L_shl (s, exp));
    exp = sub (15, exp);

    frac_en[3] = ltp_res_en;
    exp_en[3] = exp;

    // calculate LTP coding gain, i.e. energy reduction LP res -> LTP res
    if (ltp_res_en > 0 && frac_en[0] != 0)
    {
        // gain = ResEn / LTPResEn
        pred_gain = div_s (shr (frac_en[0], 1), ltp_res_en);
        exp = sub (exp, exp_en[0]);

        // L_temp = ltpGain * 2^(30 + exp)
        L_temp = L_deposit_h (pred_gain);
        // L_temp = ltpGain * 2^27
        L_temp = L_shr (L_temp, add (exp, 3));

        // Log2 = log2() + 27
        Log2(L_temp, &ltpg_exp, &ltpg_frac);

        // ltpg = log2(LtpGain) * 2^13 --> range: +- 4 = +- 12 dB
        L_temp = L_Comp (sub (ltpg_exp, 27), ltpg_frac);
        *ltpg = pv_round (L_shl (L_temp, 13)); // Q13
    }
    else
    {
        *ltpg = 0;
    }
}


------------------------------------------------------------------------------
 RESOURCES USED [optional]

 When the code is written for a specific target processor the
 the resources used should be documented below.

 HEAP MEMORY USED: x bytes

 STACK MEMORY USED: x bytes

 CLOCK CYCLES: (cycle count equation for this function) + (variable
                used to represent cycle count for each subroutine
                called)
     where: (cycle count variable) = cycle count for [subroutine
                                     name]

------------------------------------------------------------------------------
 CAUTION [optional]
 [State any special notes, constraints or cautions for users of this function]

------------------------------------------------------------------------------
*/

void calc_unfilt_energies(
    Word16 res[],     /* i  : LP residual,                               Q0  */
    Word16 exc[],     /* i  : LTP excitation (unfiltered),               Q0  */
    Word16 code[],    /* i  : CB innovation (unfiltered),                Q13 */
    Word16 gain_pit,  /* i  : pitch gain,                                Q14 */
    Word16 L_subfr,   /* i  : Subframe length                                */

    Word16 frac_en[], /* o  : energy coefficients (4), fraction part,    Q15 */
    Word16 exp_en[],  /* o  : energy coefficients (4), exponent part,    Q0  */
    Word16 *ltpg,     /* o  : LTP coding gain (log2()),                  Q13 */
    Flag   *pOverflow
)
{
    Word32 s1;      /* Intermediate energy accumulator */
    Word32 s2;      /* Intermediate energy accumulator */
    Word32 s3;      /* Intermediate energy accumulator */
    Word32 s4;      /* Intermediate energy accumulator */
    Word32 L_temp;      /* temporal 32 bits storage */

    Word16 i;       /* index used in all loops */
    Word16 exp;     /* nunmber of '0's or '1's before MSB != 0 */
    Word16 tmp1;        /* temporal storage */
    Word16 tmp2;        /* temporal storage */
    Word16 ltp_res_en;
    Word16 pred_gain;   /* predictor gain */
    Word16 ltpg_exp;    /* LTP gain (exponent) */
    Word16 ltpg_frac;   /* LTP gain (mantissa or fractional part) */

    s1 = 0;
    s2 = 0;
    s3 = 0;
    s4 = 0;

    /*----------------------------------------------------------------------------
    NOTE: Overflow is expected as a result of multiply and accumulated without
        scale down the inputs. This modification is not made at this point
        to have bit exact results with the pre-optimization code. (JT 6/20/00)

    ----------------------------------------------------------------------------*/

    for (i = 0; i < L_subfr; i++)
    {
        tmp1 = res[i];              /* avoid multiple accesses to memory */
        tmp2 = exc[i];

        s1 = amrnb_fxp_mac_16_by_16bb((Word32) tmp1, (Word32) tmp1, s1);   /* Compute residual energy */
        s2 = amrnb_fxp_mac_16_by_16bb((Word32) tmp2, (Word32) tmp2, s2);   /* Compute ltp excitation energy */
        s3 = amrnb_fxp_mac_16_by_16bb((Word32) tmp2, (Word32) code[i], s3);/* Compute scalar product */
        /* <exc[],code[]>         */

        L_temp = L_mult(tmp2, gain_pit, pOverflow);
        L_temp = L_shl(L_temp, 1, pOverflow);
        tmp2   = sub(tmp1, pv_round(L_temp, pOverflow), pOverflow);
        /* LTP residual, Q0 */
        s4     = L_mac(s4, tmp2, tmp2, pOverflow);
        /* Compute energy of LTP residual */
    }
    s1 = s1 << 1;
    s2 = s2 << 1;
    s3 = s3 << 1;

    if (s1 & MIN_32)
    {
        s1 = MAX_32;
        *pOverflow = 1;
    }

    /* ResEn := 0 if ResEn < 200.0 (= 400 Q1) */
    if (s1 < 400L)
    {
        frac_en[0] = 0;
        exp_en[0] = -15;
    }
    else
    {
        exp = norm_l(s1);
        frac_en[0] = (Word16)(L_shl(s1, exp, pOverflow) >> 16);
        exp_en[0] = (15 - exp);
    }

    if (s2 & MIN_32)
    {
        s2 = MAX_32;
        *pOverflow = 1;
    }

    exp = norm_l(s2);
    frac_en[1] = (Word16)(L_shl(s2, exp, pOverflow) >> 16);
    exp_en[1] = sub(15, exp, pOverflow);

    /*  s3 is not always sum of squares */
    exp = norm_l(s3);
    frac_en[2] = (Word16)(L_shl(s3, exp, pOverflow) >> 16);
    exp_en[2]  = 2 - exp;

    exp = norm_l(s4);
    ltp_res_en = (Word16)(L_shl(s4, exp, pOverflow) >> 16);
    exp = sub(15, exp, pOverflow);

    frac_en[3] = ltp_res_en;
    exp_en[3] = exp;

    /* calculate LTP coding gain, i.e. energy reduction LP res -> LTP res */

    if (ltp_res_en > 0 && frac_en[0] != 0)
    {
        /* gain = ResEn / LTPResEn */
        pred_gain = div_s(shr(frac_en[0], 1, pOverflow), ltp_res_en);
        exp = sub(exp, exp_en[0], pOverflow);

        /* L_temp = ltpGain * 2^(30 + exp) */
        L_temp = (Word32) pred_gain << 16;
        /* L_temp = ltpGain * 2^27 */
        L_temp = L_shr(L_temp, (Word16)(exp + 3), pOverflow);

        /* Log2 = log2() + 27 */
        Log2(L_temp, &ltpg_exp, &ltpg_frac, pOverflow);

        /* ltpg = log2(LtpGain) * 2^13 --> range: +- 4 = +- 12 dB */
        L_temp = L_Comp(sub(ltpg_exp, 27, pOverflow), ltpg_frac, pOverflow);
        *ltpg = pv_round(L_shl(L_temp, 13, pOverflow), pOverflow);   /* Q13 */
    }
    else
    {
        *ltpg = 0;
    }

    return;
}

/****************************************************************************/