sbr_hfadj.c

mpeg4 video codec mpeg4 video codec

/*
** FAAD2 - Freeware Advanced Audio (AAC) Decoder including SBR decoding
** Copyright (C) 2003-2004 M. Bakker, Ahead Software AG, http://www.nero.com
**  
** This program is free software; you can redistribute it and/or modify
** it under the terms of the GNU General Public License as published by
** the Free Software Foundation; either version 2 of the License, or
** (at your option) any later version.
** 
** This program is distributed in the hope that it will be useful,
** but WITHOUT ANY WARRANTY; without even the implied warranty of
** MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
** GNU General Public License for more details.
** 
** You should have received a copy of the GNU General Public License
** along with this program; if not, write to the Free Software 
** Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
**
** Any non-GPL usage of this software or parts of this software is strictly
** forbidden.
**
** Commercial non-GPL licensing of this software is possible.
** For more info contact Ahead Software through [e-mail].
**
** $Id: sbr_hfadj.c,v 1.1 2006/02/23 14:38:10 kevin-fu Exp $
**/

/* High Frequency adjustment */

#include "common.h"
#include "structs.h"

#ifdef SBR_DEC

#include "sbr_syntax.h"
#include "sbr_hfadj.h"

#include "sbr_noise.h"


/* static function delcarations */
static void map_noise_data(sbr_info *sbr, sbr_hfadj_info *adj, uint8_t ch);
static void map_sinusoids(sbr_info *sbr, sbr_hfadj_info *adj, uint8_t ch);
static void estimate_current_envelope(sbr_info *sbr, sbr_hfadj_info *adj,
                                      qmf_t Xsbr[MAX_NTSRHFG][64], uint8_t ch);
static void calculate_gain(sbr_info *sbr, sbr_hfadj_info *adj, uint8_t ch);
#ifdef SBR_LOW_POWER
static void calc_gain_groups(sbr_info *sbr, sbr_hfadj_info *adj, real_t *deg, uint8_t ch);
static void aliasing_reduction(sbr_info *sbr, sbr_hfadj_info *adj, real_t *deg, uint8_t ch);
#endif
static void hf_assembly(sbr_info *sbr, sbr_hfadj_info *adj, qmf_t Xsbr[MAX_NTSRHFG][64], uint8_t ch);


void hf_adjustment(sbr_info *sbr, qmf_t Xsbr[MAX_NTSRHFG][64]
#ifdef SBR_LOW_POWER
                   ,real_t *deg /* aliasing degree */
#endif
                   ,uint8_t ch)
{
    ALIGN sbr_hfadj_info adj = {{{0}}};

    map_noise_data(sbr, &adj, ch);
    map_sinusoids(sbr, &adj, ch);

    estimate_current_envelope(sbr, &adj, Xsbr, ch);

    calculate_gain(sbr, &adj, ch);

#ifdef SBR_LOW_POWER
    calc_gain_groups(sbr, &adj, deg, ch);
    aliasing_reduction(sbr, &adj, deg, ch);
#endif

    hf_assembly(sbr, &adj, Xsbr, ch);
}

static void map_noise_data(sbr_info *sbr, sbr_hfadj_info *adj, uint8_t ch)
{
    uint8_t l, i;
    uint32_t m;

    for (l = 0; l < sbr->L_E[ch]; l++)
    {
        for (i = 0; i < sbr->N_Q; i++)
        {
            for (m = sbr->f_table_noise[i]; m < sbr->f_table_noise[i+1]; m++)
            {
                uint8_t k;

                adj->Q_mapped[m - sbr->kx][l] = 0;

                for (k = 0; k < 2; k++)
                {
                    if ((sbr->t_E[ch][l] >= sbr->t_Q[ch][k]) &&
                        (sbr->t_E[ch][l+1] <= sbr->t_Q[ch][k+1]))
                    {
                        adj->Q_mapped[m - sbr->kx][l] =
                            sbr->Q_orig[ch][i][k];
                    }
                }
            }
        }
    }
}

static void map_sinusoids(sbr_info *sbr, sbr_hfadj_info *adj, uint8_t ch)
{
    uint8_t l, i, m, k, k1, k2, delta_S, l_i, u_i;

    if (sbr->bs_frame_class[ch] == FIXFIX)
    {
        sbr->l_A[ch] = -1;
    } else if (sbr->bs_frame_class[ch] == VARFIX) {
        if (sbr->bs_pointer[ch] > 1)
            sbr->l_A[ch] = -1;
        else
            sbr->l_A[ch] = sbr->bs_pointer[ch] - 1;
    } else {
        if (sbr->bs_pointer[ch] == 0)
            sbr->l_A[ch] = -1;
        else
            sbr->l_A[ch] = sbr->L_E[ch] + 1 - sbr->bs_pointer[ch];
    }

    for (l = 0; l < 5; l++)
    {
        for (i = 0; i < 64; i++)
        {
            adj->S_index_mapped[i][l] = 0;
            adj->S_mapped[i][l] = 0;
        }
    }

    for (l = 0; l < sbr->L_E[ch]; l++)
    {
        for (i = 0; i < sbr->N_high; i++)
        {
            for (m = sbr->f_table_res[HI_RES][i]; m < sbr->f_table_res[HI_RES][i+1]; m++)
            {
                uint8_t delta_step = 0;
                if ((l >= sbr->l_A[ch]) || ((sbr->bs_add_harmonic_prev[ch][i]) &&
                    (sbr->bs_add_harmonic_flag_prev[ch])))
                {
                    delta_step = 1;
                }

                if (m == (int32_t)((real_t)(sbr->f_table_res[HI_RES][i+1]+sbr->f_table_res[HI_RES][i])/2.))
                {
                    adj->S_index_mapped[m - sbr->kx][l] =
                        delta_step * sbr->bs_add_harmonic[ch][i];
                } else {
                    adj->S_index_mapped[m - sbr->kx][l] = 0;
                }
            }
        }
    }

    for (l = 0; l < sbr->L_E[ch]; l++)
    {
        for (i = 0; i < sbr->N_high; i++)
        {
            if (sbr->f[ch][l] == 1)
            {
                k1 = i;
                k2 = i + 1;
            } else {
                for (k1 = 0; k1 < sbr->N_low; k1++)
                {
                    if ((sbr->f_table_res[HI_RES][i] >= sbr->f_table_res[LO_RES][k1]) &&
                        (sbr->f_table_res[HI_RES][i+1] <= sbr->f_table_res[LO_RES][k1+1]))
                    {
                        break;
                    }
                }
                for (k2 = 0; k2 < sbr->N_low; k2++)
                {
                    if ((sbr->f_table_res[HI_RES][i+1] >= sbr->f_table_res[LO_RES][k2]) &&
                        (sbr->f_table_res[HI_RES][i+2] <= sbr->f_table_res[LO_RES][k2+1]))
                    {
                        break;
                    }
                }
            }

            l_i = sbr->f_table_res[sbr->f[ch][l]][k1];
            u_i = sbr->f_table_res[sbr->f[ch][l]][k2];

            delta_S = 0;
            for (k = l_i; k < u_i; k++)
            {
                if (adj->S_index_mapped[k - sbr->kx][l] == 1)
                    delta_S = 1;
            }

            for (m = l_i; m < u_i; m++)
            {
                adj->S_mapped[m - sbr->kx][l] = delta_S;
            }
        }
    }
}

static void estimate_current_envelope(sbr_info *sbr, sbr_hfadj_info *adj,
                                      qmf_t Xsbr[MAX_NTSRHFG][64], uint8_t ch)
{
    uint8_t m, l, j, k, k_l, k_h, p;
    real_t nrg, div;

    if (sbr->bs_interpol_freq == 1)
    {
        for (l = 0; l < sbr->L_E[ch]; l++)
        {
            uint8_t i, l_i, u_i;

            l_i = sbr->t_E[ch][l];
            u_i = sbr->t_E[ch][l+1];

            div = (real_t)(u_i - l_i);

            for (m = 0; m < sbr->M; m++)
            {
                nrg = 0;

                for (i = l_i + sbr->tHFAdj; i < u_i + sbr->tHFAdj; i++)
                {
                    nrg += MUL_R(QMF_RE(Xsbr[i][m + sbr->kx]), QMF_RE(Xsbr[i][m + sbr->kx]))
#ifndef SBR_LOW_POWER
                        + MUL_R(QMF_IM(Xsbr[i][m + sbr->kx]), QMF_IM(Xsbr[i][m + sbr->kx]))
#endif
                        ;
                }

                sbr->E_curr[ch][m][l] = nrg / div;
#ifdef SBR_LOW_POWER
                sbr->E_curr[ch][m][l] *= 2;
#endif
            }
        }
    } else {
        for (l = 0; l < sbr->L_E[ch]; l++)
        {
            for (p = 0; p < sbr->n[sbr->f[ch][l]]; p++)
            {
                k_l = sbr->f_table_res[sbr->f[ch][l]][p];
                k_h = sbr->f_table_res[sbr->f[ch][l]][p+1];

                for (k = k_l; k < k_h; k++)
                {
                    uint8_t i, l_i, u_i;
                    nrg = 0.0;

                    l_i = sbr->t_E[ch][l];
                    u_i = sbr->t_E[ch][l+1];

                    div = (real_t)((u_i - l_i)*(k_h - k_l));

                    for (i = l_i + sbr->tHFAdj; i < u_i + sbr->tHFAdj; i++)
                    {
                        for (j = k_l; j < k_h; j++)
                        {
                            nrg += MUL_R(QMF_RE(Xsbr[i][j]), QMF_RE(Xsbr[i][j]))
#ifndef SBR_LOW_POWER
                                + MUL_R(QMF_IM(Xsbr[i][j]), QMF_IM(Xsbr[i][j]))
#endif
                                ;
                        }
                    }

                    sbr->E_curr[ch][k - sbr->kx][l] = nrg / div;
#ifdef SBR_LOW_POWER
                    sbr->E_curr[ch][k - sbr->kx][l] *= 2;
#endif
                }
            }
        }
    }
}


#define EPS (1e-12)

#define ONE (1)


static void calculate_gain(sbr_info *sbr, sbr_hfadj_info *adj, uint8_t ch)
{
    static real_t limGain[] = { 0.5, 1.0, 2.0, 1e10 };
    uint8_t m, l, k, i;

    ALIGN real_t Q_M_lim[64];
    ALIGN real_t G_lim[64];
    ALIGN real_t G_boost;
    ALIGN real_t S_M[64];
    ALIGN uint8_t table_map_res_to_m[64];


    for (l = 0; l < sbr->L_E[ch]; l++)
    {
        real_t delta = (l == sbr->l_A[ch] || l == sbr->prevEnvIsShort[ch]) ? 0 : 1;

        for (i = 0; i < sbr->n[sbr->f[ch][l]]; i++)
        {
            for (m = sbr->f_table_res[sbr->f[ch][l]][i]; m < sbr->f_table_res[sbr->f[ch][l]][i+1]; m++)
            {
                table_map_res_to_m[m - sbr->kx] = i;
            }
        }

        for (k = 0; k < sbr->N_L[sbr->bs_limiter_bands]; k++)
        {
            real_t G_max;
            real_t den = 0;
            real_t acc1 = 0;
            real_t acc2 = 0;

            for (m = sbr->f_table_lim[sbr->bs_limiter_bands][k];
                 m < sbr->f_table_lim[sbr->bs_limiter_bands][k+1]; m++)
            {
                acc1 += sbr->E_orig[ch][table_map_res_to_m[m]][l];
                acc2 += sbr->E_curr[ch][m][l];
            }

            G_max = ((EPS + acc1)/(EPS + acc2)) * limGain[sbr->bs_limiter_gains];
            G_max = min(G_max, 1e10);

            for (m = sbr->f_table_lim[sbr->bs_limiter_bands][k];
                 m < sbr->f_table_lim[sbr->bs_limiter_bands][k+1]; m++)
            {
                real_t d, Q_M, G;
                real_t div2;

                div2 = adj->Q_mapped[m][l] / (1 + adj->Q_mapped[m][l]);
                Q_M = sbr->E_orig[ch][table_map_res_to_m[m]][l] * div2;

                /* 12-Nov: Changed S_mapped to S_index_mapped */
                if (adj->S_index_mapped[m][l] == 0)
                {
                    S_M[m] = 0;
                } else {
                    real_t div;

                    div = adj->S_index_mapped[m][l] / (1. + adj->Q_mapped[m][l]);
                    S_M[m] = sbr->E_orig[ch][table_map_res_to_m[m]][l] * div;
                }