Commit 99f7595a authored by JÜRGEN Gerstacker's avatar JÜRGEN Gerstacker
Browse files

added new 32bit calculations

parent 5c74348a
Loading
Loading
Loading
Loading
+223 −13
Original line number Diff line number Diff line
@@ -12,8 +12,12 @@
#include "prot_fx.h"
#include "ivas_prot_fx.h"

//#define FIX_ISSUE_1966
//#define DEBUG_ISSUE_1966
#define FIX_ISSUE_1966
#define DEBUG_ISSUE_1966
#ifdef FIX_ISSUE_1966
#include <math.h>
#include <limits.h>
#endif

#ifdef DEBUG_ISSUE_1966
extern int32_t frame;
@@ -29,7 +33,15 @@ static void getEnvelope( const Word16 nSamples, const Word32 *powerSpec, Word16
#else
static void getEnvelope( const Word16 nSamples, const Word32 *powerSpec, Word16 F0, Word32 *envelope, Word32 *smoothedSpectrum );
#endif
static void GetF0( Word16 const nSamples, Word16 const nSamplesCore, Word32 const *const powerSpectrum, Word32 const pitchLag, Word16 *const pOrigF0, Word16 *const pF0 );
static void GetF0( Word16 const nSamples, Word16 const nSamplesCore, Word32 const *const powerSpectrum, Word32 const pitchLag, Word16 *const pOrigF0,
#ifdef FIX_ISSUE_1966
    Word32 *const pOrigF0_32,
#endif
    Word16 *const pF0
#ifdef FIX_ISSUE_1966
    , Word32 *const pF0_32
#endif
);
static void findStrongestHarmonics( const Word16 nSamples, const Word32 *powerSpectrum, const Word16 F0, const Word16 nTotalHarmonics, Word16 *pHarmonicIndexes, Word16 *pnHarmonics );
static void CorrectF0( const Word16 *pHarmonicIndexes, const Word16 nHarmonics, Word16 *pF0 );
#ifdef FIX_ISSUE_1966
@@ -38,8 +50,34 @@ static void findCandidates( const Word16 nSamples, const Word32 *MDCTSpectrum, c
#else
static void findCandidates( const Word16 nSamples, const Word32 *MDCTSpectrum, const Word16 MDCTSpectrum_exp, Word16 *thresholdModificationNew, Word16 floorPowerSpectrum );
#endif
static void modifyThreshold( Word16 i, Word16 F0, Word16 threshold, Word16 *thresholdModification );
static void modifyThresholds( Word16 F0, Word16 origF0, Word16 *thresholdModification );
static void modifyThreshold( Word16 i, Word16 F0,
#ifdef FIX_ISSUE_1966
                             Word32 F0_32,
#endif
                             Word16 threshold,
#ifdef FIX_ISSUE_1966
                             Word32 threshold_32,
#endif
                             Word16 *thresholdModification
#ifdef FIX_ISSUE_1966
                             ,
                             Word32 *thresholdModification_32
#endif
);
static void modifyThresholds( Word16 F0,
#ifdef FIX_ISSUE_1966
                              Word32 F0_32,
#endif
                              Word16 origF0,
#ifdef FIX_ISSUE_1966
                              Word32 origF0_32,
#endif
                              Word16 *thresholdModification
#ifdef FIX_ISSUE_1966
                              ,
                              Word32 *thresholdModification_32
#endif
);
#ifdef FIX_ISSUE_1966
static void RefineThresholdsUsingPitch_fx( const Word16 nSamples, const Word16 nSamplesCore, const Word32 powerSpectrum[], const Word32 lastPitchLag, const Word32 currentPitchLag, 
                                           Word16 *pF0, Word32 *pF0_32, Word16 *thresholdModification, Word32 *thresholdModification32 );
@@ -54,6 +92,44 @@ static void findTonalComponents_fx( Word16 *indexOfTonalPeak, Word16 *lowerIndex
static void findTonalComponents_fx( Word16 *indexOfTonalPeak, Word16 *lowerIndex, Word16 *upperIndex, Word16 *numIndexes, Word16 nSamples, const Word32 *powerSpectrum, const Word16 powerSpectrum_e, Word16 F0, Word16 *thresholdModification, Word16 element_mode );
#endif

#ifdef FIX_ISSUE_1966
static inline int32_t double_to_q31_with_exp(double x, int *out_exp) {
    // Zerlege x in Mantisse und Exponent: x = m * 2^e, m in [0.5, 1) oder (-1, -0.5]
    int e;
    double m = frexp(x, &e); // x = m * 2^e, mit |m| in (0.5, 1]

    // Normalisiere so, dass |m| in [0.5, 1). Wir möchten später m in [-1,1) skalieren.
    // Wir wählen hier eine Darstellung mit Mantisse in [-1, 1).
    // Skaliere Mantisse weiter nach [-1, 1) falls nötig:
    // Keep m wie von frexp (|m| in (0.5,1]). Wir skalieren bei Bedarf:
    // Wir wollen final m_q31 = round(m * 2^31) und e beibehalten, aber sicherstellen
    // dass der resultierende Wert im Q31-Bereich liegt.

    // In Q31 kann Mantisse m_q31 = round(m * 2^31).
    // Grenzen beachten: |m| < 1, also m_q31 liegt in (-2^31, 2^31-1). Wir saturieren ggf.
    double scaled = m * (double)(1ULL << 31); // 2^31
    // Rundung
    if (scaled >= 0.0)
        scaled += 0.5;
    else
        scaled -= 0.5;

    // Saturation in int32_t
    if (scaled > (double)INT32_MAX) scaled = (double)INT32_MAX;
    if (scaled < (double)INT32_MIN) scaled = (double)INT32_MIN;

    int32_t m_q31 = (int32_t)scaled;

    if (out_exp) *out_exp = e;

    return m_q31;
}
// Optionaler Helper: Mantisse aus Q31 zurück zu double (nur zur Überprüfung)
static inline double q31_to_double(int32_t q31) {
    // q31 → double in [-1, 1)
    return ((double)q31) / (double)(1ULL << 31);
}
#endif

/*-------------------------------------------------------------------*
 * DetectTonalComponents()
@@ -444,7 +520,14 @@ static void GetF0(
    /*i   -  Qx */                   /*is justed handed over and given back*/
    Word32 /*int*/ const pitchLag,   /*i   -  Q16*/
    Word16 /*short*/ *const pOrigF0, /*o   -  Q10*/
    Word16 /*short*/ *const pF0 )    /*o   -  Q10*/
#ifdef FIX_ISSUE_1966
    Word32 *const pOrigF0_32,
#endif
    Word16 /*short*/ *const pF0    /*o   -  Q10*/
#ifdef FIX_ISSUE_1966
    , Word32 *const pF0_32
#endif
)
{
    Word16 /*short*/ tmpPitchLag;
    Word16 /*short*/ rgiStrongHarmonics[MAX_PEAKS_FROM_PITCH]; /*Q0*/
@@ -467,9 +550,31 @@ static void GetF0(
        /**pF0 = nSamplesCore/tmpPitchLag;*/
        BASOP_Util_Divide_MantExp( nSamplesCore, 0, tmpPitchLag, -( 1 /*division by 2*/ + 4 /*accommodate accuracy-prevention-leftshift*/ ), pF0, &tmp ); /*pF0 is Q15*/
        move16();
#ifdef FIX_ISSUE_1966
        int tmp_32 = 0;
        double k = nSamplesCore / (tmpPitchLag*0.5f*0.0625f);
        *pF0_32 = double_to_q31_with_exp(k, &tmp_32);
#ifdef DEBUG_ISSUE_1966
        static float err_max=0.f;
        double f1 = *pF0 * pow(2.f, -15+tmp);
        float tmp2= err_max;
        err_max= fmax( err_max, (f1-k)/k );
        if (err_max != tmp2) 
            printf("# frame=%d fx32=0x%08x fx16=0x%04x err=%0.5f,%f,   exp: tmp=%d, tmp_32=%d \n",
                            frame,    *pF0_32,    *pF0,    (f1-k)/k,err_max,     tmp,      tmp_32 );
#endif
        //assert(tmp_32 == tmp);
        if (tmp_32 != tmp) printf("\n\n ERROR L:%d, tmp,tmp_32 = %d,%d\n\n\n", __LINE__, tmp,tmp_32);
#endif
        *pF0 = shr_sat( *pF0, sub( 5, tmp ) ); /*Q10 without scalingfactor*/
#ifdef FIX_ISSUE_1966
        *pF0_32 = L_shr_sat( *pF0_32, sub( 5, tmp_32 ) ); /*Q10 without scalingfactor*/
#endif
        move16();
        *pOrigF0 = *pF0; /*Q10*/
#ifdef FIX_ISSUE_1966
        *pOrigF0_32 = *pF0_32; /*Q26*/
#endif
        move16();
        tmp = 2 * LAST_HARMONIC_POS_TO_CHECK;
        if ( LT_16( nSamples, 2 * LAST_HARMONIC_POS_TO_CHECK ) )
@@ -477,8 +582,28 @@ static void GetF0(
            move16();
            tmp = nSamples;
        }
#ifdef FIX_ISSUE_1966
        double k2 = ((Word32)tmp<<15) / (*pF0_32*pow(2.,-16+5));
        Word32 nTotalHarmonics_32 = double_to_q31_with_exp(k2, &tmp_32);
#endif
        BASOP_Util_Divide_MantExp( tmp, 15, *pF0, 5, &nTotalHarmonics, &tmp );
#ifdef DEBUG_ISSUE_1966
        if (tmp != tmp_32) printf("tmp=%d  tmp_32=%d\n", tmp,tmp_32 );
#endif
        assert( tmp == tmp_32 );
#ifdef DEBUG_ISSUE_1966__
        printf("# frame=%d k2=%f fx16=%f,  \n",
                        frame, k2,     nTotalHarmonics*pow(2.,-15+tmp ) );
#endif
        nTotalHarmonics = shl( nTotalHarmonics, sub( tmp, 15 ) );
#ifdef FIX_ISSUE_1966
        nTotalHarmonics_32 = L_shl( nTotalHarmonics_32, L_sub( tmp_32, 31 ) );
#ifdef DEBUG_ISSUE_1966___
        if (nTotalHarmonics != nTotalHarmonics_32)
            printf("# nTotalHarmonics=%d nTotalHarmonics_32=%d\n",
                                      nTotalHarmonics,      nTotalHarmonics_32 );
#endif
#endif

        /* Get in rgiStrongHarmonics all i for which i*F0 are the strongest harmonics */
        findStrongestHarmonics( nSamples, powerSpectrum, *pF0, nTotalHarmonics, rgiStrongHarmonics, &nStrongHarmonics );
@@ -741,8 +866,19 @@ static void CorrectF0(
static void modifyThreshold(
    Word16 /*short*/ i,                        /*I   - Q0 */
    Word16 /*short*/ F0,                       /*I   - Q10*/
#ifdef FIX_ISSUE_1966
    Word32           F0_32,                    /*I   - Qi26*/
#endif
    Word16 /*short*/ threshold,                /*I   - Q10*/
    Word16 * /*short*/ thresholdModification ) /*I/O - Q10*/
#ifdef FIX_ISSUE_1966
    Word32           threshold_32,             /*I   - Q26*/
#endif
    Word16 * /*short*/ thresholdModification   /*I   - Q10*/
#ifdef FIX_ISSUE_1966
    ,
    Word32 * /*short*/ thresholdModification_32
#endif
) /*I/O - Q10*/
{
    Word32 harmonic;
    Word16 fractional /*Q15*/;
@@ -775,8 +911,19 @@ static void modifyThreshold(

static void modifyThresholds(
    Word16 /*short*/ F0,                       /*I   - Q10*/
#ifdef FIX_ISSUE_1966
    Word32           F0_32,                    /*I   - Qi26*/
#endif
    Word16 /*short*/ origF0,                   /*I   - Q10*/
    Word16 * /*short*/ thresholdModification ) /*I/O - Q10*/
#ifdef FIX_ISSUE_1966
    Word32           origF0_32,                /*I   - Qi26*/
#endif
    Word16 * /*short*/ thresholdModification
#ifdef FIX_ISSUE_1966
    ,
    Word32 * /*short*/ thresholdModification_32
#endif
) /*I/O - Q10*/
{
    Word16 /*int*/ i, /*int*/ nHarmonics;
    Word16 tmp, tmpM, tmpE;
@@ -789,7 +936,20 @@ static void modifyThresholds(

            FOR( i = 1; i <= nHarmonics; i++ )
            {
                modifyThreshold( i, origF0, 717 /*0.7f Q10*/ /*0.7f in Q10*/, thresholdModification );
                modifyThreshold( i, origF0,
#ifdef DEBUG_ISSUE_1966
                                 origF0_32,
#endif
                                 717 /*0.7f Q10*/ /*0.7f in Q10*/,
#ifdef DEBUG_ISSUE_1966
                                 717 * 0x10000,
#endif
                                 thresholdModification
#ifdef DEBUG_ISSUE_1966
                                 ,
                                 thresholdModification_32
#endif
               );
            }
        }
        IF( F0 > 0 )
@@ -802,11 +962,37 @@ static void modifyThresholds(

            FOR( i = tmp; i > 0; i-- )
            {
                modifyThreshold( i, origF0, 358 /*0.35f Q10*/, thresholdModification );
                modifyThreshold( i, origF0,
#ifdef DEBUG_ISSUE_1966
                                 origF0_32,
#endif
                                 358 /*0.35f Q10*/,
#ifdef DEBUG_ISSUE_1966
                                 358 * 0x10000,
#endif
                                 thresholdModification
#ifdef DEBUG_ISSUE_1966
                                 ,
                                 thresholdModification_32
#endif
                );
            }
            FOR( i = 1; i <= nHarmonics; i++ )
            {
                modifyThreshold( i, F0, 358 /*0.35f Q10*/, thresholdModification );
                modifyThreshold( i, F0,
#ifdef DEBUG_ISSUE_1966
                                 origF0_32,
#endif
                                 358 /*0.35f Q10*/,
#ifdef DEBUG_ISSUE_1966
                                 358 * 0x10000,
#endif
                                 thresholdModification
#ifdef DEBUG_ISSUE_1966
                                 ,
                                 thresholdModification_32
#endif
                );
            }
        }
    }
@@ -1184,6 +1370,9 @@ static void RefineThresholdsUsingPitch_fx(
{
    Word16 pitchIsStable;
    Word16 origF0;
#ifdef FIX_ISSUE_1966
    Word32 origF0_32;
#endif
    Word32 L_tmp;

    /*pitchIsStable = (fabs(lastPitchLag-currentPitchLag) < 0.25f);*/
@@ -1198,9 +1387,30 @@ static void RefineThresholdsUsingPitch_fx(

    IF( pitchIsStable )
    {
        GetF0( nSamples, nSamplesCore, powerSpectrum, lastPitchLag, &origF0, pF0 );
        GetF0( nSamples, nSamplesCore, powerSpectrum, lastPitchLag, &origF0,
#ifdef FIX_ISSUE_1966
                                                                    &origF0_32,
#endif
                                                                    pF0
#ifdef FIX_ISSUE_1966
                                                                    ,pF0_32
#endif
        );

        modifyThresholds( *pF0, origF0, thresholdModification );
        modifyThresholds( *pF0,
#ifdef FIX_ISSUE_1966
                          *pF0_32,
#endif
                          origF0,
#ifdef FIX_ISSUE_1966
                          origF0_32,
#endif
                          thresholdModification
#ifdef FIX_ISSUE_1966
                          ,
                          thresholdModification32
#endif
        );
    }
    ELSE
    {