gsm.c

mpeg4 video codec mpeg4 video codec

/*
   Written by Mark Podlipec <[e-mail]>.

   Most of this code comes from a GSM 06.10 library by
   Jutta Degener and Carsten Bormann, available via
   <http://www.pobox.com/~jutta/toast.html>.

   That library is distributed with the following copyright:

    Copyright 1992 by Jutta Degener and Carsten Bormann,
    Technische Universitaet Berlin

Any use of this software is permitted provided that this notice is not
removed and that neither the authors nor the Technische Universitaet Berlin
are deemed to have made any representations as to the suitability of this
software for any purpose nor are held responsible for any defects of
this software.  THERE IS ABSOLUTELY NO WARRANTY FOR THIS SOFTWARE.

As a matter of courtesy, the authors request to be informed about uses
this software has found, about bugs in this software, and about any
improvements that may be of general interest.

    Berlin, 15.09.1992
    Jutta Degener
    Carsten Bormann
*/

#include "avcodec.h"

typedef struct XAGSMstate {

	int16_t		dp0[ 280 ];

	int16_t		z1;		/* preprocessing.c, Offset_com. */
	int32_t	L_z2;		/*                  Offset_com. */
	int		mp;		/*                  Preemphasis	*/

	int16_t		u[8];		/* short_term_aly_filter.c	*/
	int16_t		LARpp[2][8]; 	/*                              */
	int16_t		j;		/*                              */

	int16_t            ltp_cut;        /* long_term.c, LTP crosscorr.  */
	int16_t		nrp; /* 40 */	/* long_term.c, synthesis	*/
	int16_t		v[9];		/* short_term.c, synthesis	*/
	int16_t		msr;		/* decoder.c,	Postprocessing	*/

	char		verbose;	/* only used if !NDEBUG		*/
	char		fast;		/* only used if FAST		*/

	char		wav_fmt;	/* only used if WAV49 defined	*/
	unsigned char	frame_index;	/*            odd/even chaining	*/
	unsigned char	frame_chain;	/*   half-byte to carry forward	*/
} XAGSMstate;

static int gsm_decode_init(AVCodecContext * avctx)
{
    XAGSMstate *s = avctx->priv_data;

    memset(s,0,sizeof(XAGSMstate));
    s->nrp = 40;

    return 0;
}

/*   Table 4.3b   Quantization levels of the LTP gain quantizer
 */
/* bc                 0          1        2          3                  */
static const int16_t gsm_QLB[4] = {  3277,    11469,    21299,     32767        };

/*   Table 4.6   Normalized direct mantissa used to compute xM/xmax
 */
/* i                  0      1       2      3      4      5      6      7   */
static const int16_t gsm_FAC[8] = { 18431, 20479, 22527, 24575, 26623, 28671, 30719, 32767 };

#define	MIN_WORD	(-32767 - 1)
#define	MAX_WORD	  32767

#ifdef	SASR		/* flag: >> is a signed arithmetic shift right */
#undef	SASR
#define	SASR(x, by)	((x) >> (by))
#else
#define	SASR(x, by)	((x) >= 0 ? (x) >> (by) : (~(-((x) + 1) >> (by))))
#endif	/* SASR */

/* 
 * #define GSM_MULT_R(a, b) (* word a, word b, !(a == b == MIN_WORD) *)	\
 *	(0x0FFFF & SASR(((int32_t)(a) * (int32_t)(b) + 16384), 15))
 */
#define GSM_MULT_R(a, b) /* word a, word b, !(a == b == MIN_WORD) */	\
	(SASR( ((int32_t)(a) * (int32_t)(b) + 16384), 15 ))

# define GSM_MULT(a,b)	 /* word a, word b, !(a == b == MIN_WORD) */	\
	(SASR( ((int32_t)(a) * (int32_t)(b)), 15 ))

# define GSM_L_MULT(a, b) /* word a, word b */	\
	(((int32_t)(a) * (int32_t)(b)) << 1)

# define GSM_L_ADD(a, b)	\
	( (a) <  0 ? ( (b) >= 0 ? (a) + (b)	\
		 : (utmp = (uint32_t)-((a) + 1) + (uint32_t)-((b) + 1)) \
		   >= MAX_LONGWORD ? MIN_LONGWORD : -(int32_t)utmp-2 )   \
	: ((b) <= 0 ? (a) + (b)   \
	          : (utmp = (uint32_t)(a) + (uint32_t)(b)) >= MAX_LONGWORD \
		    ? MAX_LONGWORD : utmp))

/*
 * # define GSM_ADD(a, b)	\
 * 	((ltmp = (int32_t)(a) + (int32_t)(b)) >= MAX_WORD \
 * 	? MAX_WORD : ltmp <= MIN_WORD ? MIN_WORD : ltmp)
 */
/* Nonportable, but faster: */

#define	GSM_ADD(a, b)	\
	((uint32_t)((ltmp = (int32_t)(a) + (int32_t)(b)) - MIN_WORD) > \
		MAX_WORD - MIN_WORD ? (ltmp > 0 ? MAX_WORD : MIN_WORD) : ltmp)

# define GSM_SUB(a, b)	\
	((ltmp = (int32_t)(a) - (int32_t)(b)) >= MAX_WORD \
	? MAX_WORD : ltmp <= MIN_WORD ? MIN_WORD : ltmp)

# define GSM_ABS(a)	((a) < 0 ? ((a) == MIN_WORD ? MAX_WORD : -(a)) : (a))

/****************/
#define saturate(x)     \
        ((x) < MIN_WORD ? MIN_WORD : (x) > MAX_WORD ? MAX_WORD: (x))

/****************/
static int16_t gsm_sub (
int16_t a,
int16_t b)
{
        int32_t diff = (int32_t)a - (int32_t)b;
        return saturate(diff);
}

/****************/
static int16_t gsm_asr (
int16_t a,
int n)
{
        if (n >= 16) return -(a < 0);
        if (n <= -16) return 0;
        if (n < 0) return a << -n;

#       ifdef   SASR
                return a >> n;
#       else
                if (a >= 0) return a >> n;
                else return -(int16_t)( -(uint16_t)a >> n );
#       endif
}

/****************/
static int16_t gsm_asl (
int16_t a,
int n)
{
        if (n >= 16) return 0;
        if (n <= -16) return -(a < 0);
        if (n < 0) return gsm_asr(a, -n);
        return a << n;
}


/*
 * Copyright 1992 by Jutta Degener and Carsten Bormann, Technische
 * Universitaet Berlin.  See the accompanying file "COPYRIGHT" for
 * details.  THERE IS ABSOLUTELY NO WARRANTY FOR THIS SOFTWARE.
 */

/**** 4.2.17 */
static void RPE_grid_positioning(
int16_t            Mc,             /* grid position        IN      */
register int16_t   * xMp,          /* [0..12]              IN      */
register int16_t   * ep)           /* [0..39]              OUT     */
/*
 *  This procedure computes the reconstructed long term residual signal
 *  ep[0..39] for the LTP analysis filter.  The inputs are the Mc
 *  which is the grid position selection and the xMp[0..12] decoded
 *  RPE samples which are upsampled by a factor of 3 by inserting zero
 *  values.
 */
{
        int     i = 13;

        assert(0 <= Mc && Mc <= 3);

        switch (Mc) {
                case 3: *ep++ = 0;
                case 2:  do {
                                *ep++ = 0;
                case 1:         *ep++ = 0;
                case 0:         *ep++ = *xMp++;
                         } while (--i);
        }
        while (++Mc < 4) *ep++ = 0;

        /*

        int i, k;
        for (k = 0; k <= 39; k++) ep[k] = 0;
        for (i = 0; i <= 12; i++) {
                ep[ Mc + (3*i) ] = xMp[i];
        }
        */
}


/**** 4.2.16 */
static void APCM_inverse_quantization (
register int16_t   * xMc,  /* [0..12]                      IN      */
int16_t            mant,
int16_t            exp,
register int16_t   * xMp)  /* [0..12]                      OUT     */
/* 
 *  This part is for decoding the RPE sequence of coded xMc[0..12]
 *  samples to obtain the xMp[0..12] array.  Table 4.6 is used to get
 *  the mantissa of xmaxc (FAC[0..7]).
 */
{
        int     i;
        int16_t    temp, temp1, temp2, temp3;
        int32_t        ltmp;

        assert( mant >= 0 && mant <= 7 ); 

        temp1 = gsm_FAC[ mant ];        /* see 4.2-15 for mant */
        temp2 = gsm_sub( 6, exp );      /* see 4.2-15 for exp  */
        temp3 = gsm_asl( 1, gsm_sub( temp2, 1 ));

        for (i = 13; i--;) {

                assert( *xMc <= 7 && *xMc >= 0 );       /* 3 bit unsigned */

                /* temp = gsm_sub( *xMc++ << 1, 7 ); */
                temp = (*xMc++ << 1) - 7;               /* restore sign   */
                assert( temp <= 7 && temp >= -7 );      /* 4 bit signed   */

                temp <<= 12;                            /* 16 bit signed  */
                temp = GSM_MULT_R( temp1, temp );
                temp = GSM_ADD( temp, temp3 );
                *xMp++ = gsm_asr( temp, temp2 );
        }
}


/**** 4.12.15 */
static void APCM_quantization_xmaxc_to_exp_mant (
int16_t            xmaxc,          /* IN   */
int16_t            * exp_out,      /* OUT  */
int16_t            * mant_out)    /* OUT  */
{
  int16_t    exp, mant;

  /* Compute exponent and mantissa of the decoded version of xmaxc
   */

        exp = 0;
        if (xmaxc > 15) exp = SASR(xmaxc, 3) - 1;
        mant = xmaxc - (exp << 3);

        if (mant == 0) {
                exp  = -4;
                mant = 7;
        }
        else {
                while (mant <= 7) {
                        mant = mant << 1 | 1;
                        exp--;
                }
                mant -= 8;
        }

        assert( exp  >= -4 && exp <= 6 );
        assert( mant >= 0 && mant <= 7 );

        *exp_out  = exp;
        *mant_out = mant;
}

static void Gsm_RPE_Decoding (
XAGSMstate        * S,
int16_t            xmaxcr,
int16_t            Mcr,
int16_t            * xMcr,  /* [0..12], 3 bits             IN      */
int16_t            * erp)   /* [0..39]                     OUT     */

{
        int16_t    exp, mant;
        int16_t    xMp[ 13 ];

        APCM_quantization_xmaxc_to_exp_mant( xmaxcr, &exp, &mant );
        APCM_inverse_quantization( xMcr, mant, exp, xMp );
        RPE_grid_positioning( Mcr, xMp, erp );

}


/*
 *  4.3 FIXED POINT IMPLEMENTATION OF THE RPE-LTP DECODER
 */

static void Postprocessing(
XAGSMstate	* S,
register int16_t 	* s)
{
  register int		k;
  register int16_t		msr = S->msr;
  register int32_t	ltmp;	/* for GSM_ADD */
  register int16_t		tmp;

  for (k = 160; k--; s++) 
  {
    tmp = GSM_MULT_R( msr, 28180 );
    msr = GSM_ADD(*s, tmp);  	   /* Deemphasis 	     */
    *s  = GSM_ADD(msr, msr) & 0xFFF8;  /* Truncation & Upscaling */
  }
  S->msr = msr;
}

/**** 4.3.2 */
static void Gsm_Long_Term_Synthesis_Filtering (
XAGSMstate        * S,
int16_t                    Ncr,
int16_t                    bcr,
register int16_t           * erp,     /* [0..39]                    IN */
register int16_t           * drp)     /* [-120..-1] IN, [-120..40] OUT */

/*
 *  This procedure uses the bcr and Ncr parameter to realize the
 *  long term synthesis filtering.  The decoding of bcr needs
 *  table 4.3b.
 */
{
        register int32_t       ltmp;   /* for ADD */
        register int            k;
        int16_t                    brp, drpp, Nr;

        /*  Check the limits of Nr.
         */
        Nr = Ncr < 40 || Ncr > 120 ? S->nrp : Ncr;
        S->nrp = Nr;
        assert(Nr >= 40 && Nr <= 120);

        /*  Decoding of the LTP gain bcr
         */
        brp = gsm_QLB[ bcr ];

        /*  Computation of the reconstructed short term residual 
         *  signal drp[0..39]
         */
        assert(brp != MIN_WORD);

        for (k = 0; k <= 39; k++) {
                drpp   = GSM_MULT_R( brp, drp[ k - Nr ] );
                drp[k] = GSM_ADD( erp[k], drpp );
        }

        /*
         *  Update of the reconstructed short term residual signal
         *  drp[ -1..-120 ]
         */

        for (k = 0; k <= 119; k++) drp[ -120 + k ] = drp[ -80 + k ];
}

static void Short_term_synthesis_filtering (
XAGSMstate *S,
register int16_t   *rrp,  /* [0..7]       IN      */
register int    k,      /* k_end - k_start      */
register int16_t   *wt,   /* [0..k-1]     IN      */
register int16_t   *sr)   /* [0..k-1]     OUT     */
{
        register int16_t           * v = S->v;
        register int            i;
        register int16_t           sri, tmp1, tmp2;
        register int32_t       ltmp;   /* for GSM_ADD  & GSM_SUB */

        while (k--) {
                sri = *wt++;
                for (i = 8; i--;) {

                        /* sri = GSM_SUB( sri, gsm_mult_r( rrp[i], v[i] ) );
                         */
                        tmp1 = rrp[i];
                        tmp2 = v[i];
                        tmp2 =  ( tmp1 == MIN_WORD && tmp2 == MIN_WORD
                                ? MAX_WORD
                                : 0x0FFFF & (( (int32_t)tmp1 * (int32_t)tmp2
                                             + 16384) >> 15)) ;

                        sri  = GSM_SUB( sri, tmp2 );

                        /* v[i+1] = GSM_ADD( v[i], gsm_mult_r( rrp[i], sri ) );
                         */
                        tmp1  = ( tmp1 == MIN_WORD && sri == MIN_WORD
                                ? MAX_WORD
                                : 0x0FFFF & (( (int32_t)tmp1 * (int32_t)sri
                                             + 16384) >> 15)) ;

                        v[i+1] = GSM_ADD( v[i], tmp1);
                }
                *sr++ = v[0] = sri;
        }
}

/* 4.2.8 */

static void Decoding_of_the_coded_Log_Area_Ratios (
int16_t    * LARc,         /* coded log area ratio [0..7]  IN      */
int16_t    * LARpp)        /* out: decoded ..                      */
{
        register int16_t   temp1 /* , temp2 */;
        register long   ltmp;   /* for GSM_ADD */

        /*  This procedure requires for efficient implementation
         *  two tables.
         *
         *  INVA[1..8] = integer( (32768 * 8) / real_A[1..8])
         *  MIC[1..8]  = minimum value of the LARc[1..8]
         */

        /*  Compute the LARpp[1..8]
         */

        /*      for (i = 1; i <= 8; i++, B++, MIC++, INVA++, LARc++, LARpp++) {
         *
         *              temp1  = GSM_ADD( *LARc, *MIC ) << 10;
         *              temp2  = *B << 1;
         *              temp1  = GSM_SUB( temp1, temp2 );
         *
         *              assert(*INVA != MIN_WORD);
         *
         *              temp1  = GSM_MULT_R( *INVA, temp1 );
         *              *LARpp = GSM_ADD( temp1, temp1 );
         *      }
         */

#undef  STEP
#define STEP( B, MIC, INVA )    \
                temp1    = GSM_ADD( *LARc++, MIC ) << 10;       \
                temp1    = GSM_SUB( temp1, B << 1 );            \
                temp1    = GSM_MULT_R( INVA, temp1 );           \
                *LARpp++ = GSM_ADD( temp1, temp1 );

        STEP(      0,  -32,  13107 );
        STEP(      0,  -32,  13107 );
        STEP(   2048,  -16,  13107 );
        STEP(  -2560,  -16,  13107 );