Donate to e Foundation | Murena handsets with /e/OS | Own a part of Murena! Learn more

Commit 3aa85e3b authored by Andy Hung's avatar Andy Hung Committed by Android (Google) Code Review
Browse files

Merge "AudioResampler: Add configurable resampler design"

parents c4f80a51 6bd378fe
Loading
Loading
Loading
Loading
+143 −71
Original line number Diff line number Diff line
@@ -38,6 +38,9 @@

//#define DEBUG_RESAMPLER

// use this for our buffer alignment.  Should be at least 32 bytes.
constexpr size_t CACHE_LINE_SIZE = 64;

namespace android {

/*
@@ -94,7 +97,10 @@ void AudioResamplerDyn<TC, TI, TO>::InBuffer::resize(int CHANNELS, int halfNumCo

    // create new buffer
    TI* state = NULL;
    (void)posix_memalign(reinterpret_cast<void**>(&state), 32, stateCount*sizeof(*state));
    (void)posix_memalign(
            reinterpret_cast<void **>(&state),
            CACHE_LINE_SIZE /* alignment */,
            stateCount * sizeof(*state));
    memset(state, 0, stateCount*sizeof(*state));

    // attempt to preserve state
@@ -185,6 +191,16 @@ AudioResamplerDyn<TC, TI, TO>::AudioResamplerDyn(
    // setSampleRate() for 1:1. (May be removed if precalculated filters are used.)
    mInSampleRate = 0;
    mConstants.set(128, 8, mSampleRate, mSampleRate); // TODO: set better

    // fetch property based resampling parameters
    mPropertyEnableAtSampleRate = property_get_int32(
            "ro.audio.resampler.psd.enable_at_samplerate", mPropertyEnableAtSampleRate);
    mPropertyHalfFilterLength = property_get_int32(
            "ro.audio.resampler.psd.halflength", mPropertyHalfFilterLength);
    mPropertyStopbandAttenuation = property_get_int32(
            "ro.audio.resampler.psd.stopband", mPropertyStopbandAttenuation);
    mPropertyCutoffPercent = property_get_int32(
            "ro.audio.resampler.psd.cutoff_percent", mPropertyCutoffPercent);
}

template<typename TC, typename TI, typename TO>
@@ -215,6 +231,8 @@ void AudioResamplerDyn<TC, TI, TO>::setVolume(float left, float right)
    }
}

// TODO: update to C++11

template<typename T> T max(T a, T b) {return a > b ? a : b;}

template<typename T> T absdiff(T a, T b) {return a > b ? a - b : b - a;}
@@ -223,37 +241,74 @@ template<typename TC, typename TI, typename TO>
void AudioResamplerDyn<TC, TI, TO>::createKaiserFir(Constants &c,
        double stopBandAtten, int inSampleRate, int outSampleRate, double tbwCheat)
{
    TC* buf = NULL;
    static const double atten = 0.9998;   // to avoid ripple overflow
    double fcr;
    double tbw = firKaiserTbw(c.mHalfNumCoefs, stopBandAtten);
    // compute the normalized transition bandwidth
    const double tbw = firKaiserTbw(c.mHalfNumCoefs, stopBandAtten);
    const double halfbw = tbw / 2.;

    (void)posix_memalign(reinterpret_cast<void**>(&buf), 32, (c.mL+1)*c.mHalfNumCoefs*sizeof(TC));
    double fcr; // compute fcr, the 3 dB amplitude cut-off.
    if (inSampleRate < outSampleRate) { // upsample
        fcr = max(0.5*tbwCheat - tbw/2, tbw/2);
        fcr = max(0.5 * tbwCheat - halfbw, halfbw);
    } else { // downsample
        fcr = max(0.5*tbwCheat*outSampleRate/inSampleRate - tbw/2, tbw/2);
        fcr = max(0.5 * tbwCheat * outSampleRate / inSampleRate - halfbw, halfbw);
    }
    // create and set filter
    firKaiserGen(buf, c.mL, c.mHalfNumCoefs, stopBandAtten, fcr, atten);
    c.mFirCoefs = buf;
    if (mCoefBuffer) {
        free(mCoefBuffer);
    createKaiserFir(c, stopBandAtten, fcr);
}
    mCoefBuffer = buf;
#ifdef DEBUG_RESAMPLER

template<typename TC, typename TI, typename TO>
void AudioResamplerDyn<TC, TI, TO>::createKaiserFir(Constants &c,
        double stopBandAtten, double fcr) {
    // compute the normalized transition bandwidth
    const double tbw = firKaiserTbw(c.mHalfNumCoefs, stopBandAtten);
    const int phases = c.mL;
    const int halfLength = c.mHalfNumCoefs;

    // create buffer
    TC *coefs = nullptr;
    int ret = posix_memalign(
            reinterpret_cast<void **>(&coefs),
            CACHE_LINE_SIZE /* alignment */,
            (phases + 1) * halfLength * sizeof(TC));
    LOG_ALWAYS_FATAL_IF(ret != 0, "Cannot allocate buffer memory, ret %d", ret);
    c.mFirCoefs = coefs;
    free(mCoefBuffer);
    mCoefBuffer = coefs;

    // square the computed minimum passband value (extra safety).
    double attenuation =
            computeWindowedSincMinimumPassbandValue(stopBandAtten);
    attenuation *= attenuation;

    // design filter
    firKaiserGen(coefs, phases, halfLength, stopBandAtten, fcr, attenuation);

    // update the design criteria
    mNormalizedCutoffFrequency = fcr;
    mNormalizedTransitionBandwidth = tbw;
    mFilterAttenuation = attenuation;
    mStopbandAttenuationDb = stopBandAtten;
    mPassbandRippleDb = computeWindowedSincPassbandRippleDb(stopBandAtten);

#if 0
    // Keep this debug code in case an app causes resampler design issues.
    const double halfbw = tbw / 2.;
    // print basic filter stats
    printf("L:%d  hnc:%d  stopBandAtten:%lf  fcr:%lf  atten:%lf  tbw:%lf\n",
            c.mL, c.mHalfNumCoefs, stopBandAtten, fcr, atten, tbw);
    // test the filter and report results
    double fp = (fcr - tbw/2)/c.mL;
    double fs = (fcr + tbw/2)/c.mL;
    ALOGD("L:%d  hnc:%d  stopBandAtten:%lf  fcr:%lf  atten:%lf  tbw:%lf\n",
            c.mL, c.mHalfNumCoefs, stopBandAtten, fcr, attenuation, tbw);

    // test the filter and report results.
    // Since this is a polyphase filter, normalized fp and fs must be scaled.
    const double fp = (fcr - halfbw) / phases;
    const double fs = (fcr + halfbw) / phases;

    double passMin, passMax, passRipple;
    double stopMax, stopRipple;
    testFir(buf, c.mL, c.mHalfNumCoefs, fp, fs, /*passSteps*/ 1000, /*stopSteps*/ 100000,

    const int32_t passSteps = 1000;

    testFir(coefs, c.mL, c.mHalfNumCoefs, fp, fs, passSteps, passSteps * c.ML /*stopSteps*/,
            passMin, passMax, passRipple, stopMax, stopRipple);
    printf("passband(%lf, %lf): %.8lf %.8lf %.8lf\n", 0., fp, passMin, passMax, passRipple);
    printf("stopband(%lf, %lf): %.8lf %.3lf\n", fs, 0.5, stopMax, stopRipple);
    ALOGD("passband(%lf, %lf): %.8lf %.8lf %.8lf\n", 0., fp, passMin, passMax, passRipple);
    ALOGD("stopband(%lf, %lf): %.8lf %.3lf\n", fs, 0.5, stopMax, stopRipple);
#endif
}

@@ -304,6 +359,11 @@ void AudioResamplerDyn<TC, TI, TO>::setSampleRate(int32_t inSampleRate)
        mFilterSampleRate = inSampleRate;
        mFilterQuality = getQuality();

        double stopBandAtten;
        double tbwCheat = 1.; // how much we "cheat" into aliasing
        int halfLength;
        double fcr = 0.;

        // Begin Kaiser Filter computation
        //
        // The quantization floor for S16 is about 96db - 10*log_10(#length) + 3dB.
@@ -313,9 +373,16 @@ void AudioResamplerDyn<TC, TI, TO>::setSampleRate(int32_t inSampleRate)
        // 96-98dB
        //

        double stopBandAtten;
        double tbwCheat = 1.; // how much we "cheat" into aliasing
        int halfLength;
        if (mPropertyEnableAtSampleRate >= 0 && mSampleRate >= mPropertyEnableAtSampleRate) {
            // An alternative method which allows allows a greater fcr
            // at the expense of potential aliasing.
            halfLength = mPropertyHalfFilterLength;
            stopBandAtten = mPropertyStopbandAttenuation;
            useS32 = true;
            fcr = mInSampleRate <= mSampleRate
                    ? 0.5 : 0.5 * mSampleRate / mInSampleRate;
            fcr *= mPropertyCutoffPercent / 100.;
        } else {
            if (mFilterQuality == DYN_HIGH_QUALITY) {
                // 32b coefficients, 64 length
                useS32 = true;
@@ -361,6 +428,7 @@ void AudioResamplerDyn<TC, TI, TO>::setSampleRate(int32_t inSampleRate)
                    tbwCheat = 1.01;
                }
            }
        }

        // determine the number of polyphases in the filterbank.
        // for 16b, it is desirable to have 2^(16/2) = 256 phases.
@@ -390,8 +458,12 @@ void AudioResamplerDyn<TC, TI, TO>::setSampleRate(int32_t inSampleRate)

        // create the filter
        mConstants.set(phases, halfLength, inSampleRate, mSampleRate);
        if (fcr > 0.) {
            createKaiserFir(mConstants, stopBandAtten, fcr);
        } else {
            createKaiserFir(mConstants, stopBandAtten,
                    inSampleRate, mSampleRate, tbwCheat);
        }
    } // End Kaiser filter

    // update phase and state based on the new filter.
+67 −0
Original line number Diff line number Diff line
@@ -55,6 +55,39 @@ public:
    virtual size_t resample(int32_t* out, size_t outFrameCount,
            AudioBufferProvider* provider);

    // Make available key design criteria for testing
    int getHalfLength() const {
        return mConstants.mHalfNumCoefs;
    }

    const TC *getFilterCoefs() const {
        return mConstants.mFirCoefs;
    }

    int getPhases() const {
        return mConstants.mL;
    }

    double getStopbandAttenuationDb() const {
        return mStopbandAttenuationDb;
    }

    double getPassbandRippleDb() const {
        return mPassbandRippleDb;
    }

    double getNormalizedTransitionBandwidth() const {
        return mNormalizedTransitionBandwidth;
    }

    double getFilterAttenuation() const {
        return mFilterAttenuation;
    }

    double getNormalizedCutoffFrequency() const {
        return mNormalizedCutoffFrequency;
    }

private:

    class Constants { // stores the filter constants.
@@ -112,6 +145,8 @@ private:
    void createKaiserFir(Constants &c, double stopBandAtten,
            int inSampleRate, int outSampleRate, double tbwCheat);

    void createKaiserFir(Constants &c, double stopBandAtten, double fcr);

    template<int CHANNELS, bool LOCKED, int STRIDE>
    size_t resample(TO* out, size_t outFrameCount, AudioBufferProvider* provider);

@@ -127,6 +162,38 @@ private:
            int32_t mFilterSampleRate; // designed filter sample rate.
        src_quality mFilterQuality;    // designed filter quality.
              void* mCoefBuffer;       // if a filter is created, this is not null

    // Property selected design parameters.
              // This will enable fixed high quality resampling.

              // 32 char PROP_NAME_MAX limit enforced before Android O

              // Use for sample rates greater than or equal to this value.
              // Set to non-negative to enable, negative to disable.
              int32_t mPropertyEnableAtSampleRate = 48000;
                      // "ro.audio.resampler.psd.enable_at_samplerate"

              // Specify HALF the resampling filter length.
              // Set to a value which is a multiple of 4.
              int32_t mPropertyHalfFilterLength = 32;
                      // "ro.audio.resampler.psd.halflength"

              // Specify the stopband attenuation in positive dB.
              // Set to a value greater or equal to 20.
              int32_t mPropertyStopbandAttenuation = 90;
                      // "ro.audio.resampler.psd.stopband"

              // Specify the cutoff frequency as a percentage of Nyquist.
              // Set to a value between 50 and 100.
              int32_t mPropertyCutoffPercent = 100;
                      // "ro.audio.resampler.psd.cutoff_percent"

    // Filter creation design parameters, see setSampleRate()
             double mStopbandAttenuationDb = 0.;
             double mPassbandRippleDb = 0.;
             double mNormalizedTransitionBandwidth = 0.;
             double mFilterAttenuation = 0.;
             double mNormalizedCutoffFrequency = 0.;
};

} // namespace android
+78 −29
Original line number Diff line number Diff line
@@ -546,8 +546,9 @@ static void testFir(const T* coef, int L, int halfNumCoef,
        }
        wstart += wstep;
    }
    // renormalize - this is only needed for integer filter types
    double norm = 1./((1ULL<<(sizeof(T)*8-1))*L);
    // renormalize - this is needed for integer filter types, use 1 for float or double.
    constexpr int64_t integralShift = std::is_integral<T>::value ? (sizeof(T) * 8 - 1) : 0;
    const double norm = 1. / (L << integralShift);

    firMin = fmin * norm;
    firMax = fmax * norm;
@@ -557,9 +558,12 @@ static void testFir(const T* coef, int L, int halfNumCoef,
 * evaluates the |H(f)| lowpass band characteristics.
 *
 * This function tests the lowpass characteristics for the overall polyphase filter,
 * and is used to verify the design.  For this case, fp should be set to the
 * and is used to verify the design.
 *
 * For a polyphase filter (L > 1), typically fp should be set to the
 * passband normalized frequency from 0 to 0.5 for the overall filter (thus it
 * is the designed polyphase bank value / L).  Likewise for fs.
 * Similarly the stopSteps should be L * passSteps for equivalent accuracy.
 *
 * @param coef is the designed polyphase filter banks
 *
@@ -609,6 +613,74 @@ static void testFir(const T* coef, int L, int halfNumCoef,
    stopRipple = -20.*log10(fmax); // stopband ripple/attenuation
}

/*
 * Estimate the windowed sinc minimum passband value.
 *
 * This is the minimum value for a windowed sinc filter in its passband,
 * which is identical to the scaling required not to cause overflow of a 0dBFS signal.
 * The actual value used to attenuate the filter amplitude should be slightly
 * smaller than this (suggest squaring) as this is just an estimate.
 *
 * As a windowed sinc has a passband ripple commensurate to the stopband attenuation
 * due to Gibb's phenomenon from truncating the sinc, we derive this value from
 * the design stopbandAttenuationDb (a positive value).
 */
static inline double computeWindowedSincMinimumPassbandValue(
        double stopBandAttenuationDb) {
    return 1. - pow(10. /* base */, stopBandAttenuationDb * (-1. / 20.));
}

/*
 * Compute the windowed sinc passband ripple from stopband attenuation.
 *
 * As a windowed sinc has an passband ripple commensurate to the stopband attenuation
 * due to Gibb's phenomenon from truncating the sinc, we derive this value from
 * the design stopbandAttenuationDb (a positive value).
 */
static inline double computeWindowedSincPassbandRippleDb(
        double stopBandAttenuationDb) {
    return -20. * log10(computeWindowedSincMinimumPassbandValue(stopBandAttenuationDb));
}

/*
 * Kaiser window Beta value
 *
 * Formula 3.2.5, 3.2.7, Vaidyanathan, _Multirate Systems and Filter Banks_, p. 48
 * Formula 7.75, Oppenheim and Schafer, _Discrete-time Signal Processing, 3e_, p. 542
 *
 * See also: http://melodi.ee.washington.edu/courses/ee518/notes/lec17.pdf
 *
 * Kaiser window and beta parameter
 *
 *         | 0.1102*(A - 8.7)                         A > 50
 *  Beta = | 0.5842*(A - 21)^0.4 + 0.07886*(A - 21)   21 < A <= 50
 *         | 0.                                       A <= 21
 *
 * with A is the desired stop-band attenuation in positive dBFS
 *
 *    30 dB    2.210
 *    40 dB    3.384
 *    50 dB    4.538
 *    60 dB    5.658
 *    70 dB    6.764
 *    80 dB    7.865
 *    90 dB    8.960
 *   100 dB   10.056
 *
 * For some values of stopBandAttenuationDb the function may be computed
 * at compile time.
 */
static inline constexpr double computeBeta(double stopBandAttenuationDb) {
    if (stopBandAttenuationDb > 50.) {
        return 0.1102 * (stopBandAttenuationDb - 8.7);
    }
    const double offset = stopBandAttenuationDb - 21.;
    if (offset > 0.) {
        return 0.5842 * pow(offset, 0.4) + 0.07886 * offset;
    }
    return 0.;
}

/*
 * Calculates the overall polyphase filter based on a windowed sinc function.
 *
@@ -642,31 +714,8 @@ static void testFir(const T* coef, int L, int halfNumCoef,
template <typename T>
static inline void firKaiserGen(T* coef, int L, int halfNumCoef,
        double stopBandAtten, double fcr, double atten) {
    //
    // Formula 3.2.5, 3.2.7, Vaidyanathan, _Multirate Systems and Filter Banks_, p. 48
    // Formula 7.75, Oppenheim and Schafer, _Discrete-time Signal Processing, 3e_, p. 542
    //
    // See also: http://melodi.ee.washington.edu/courses/ee518/notes/lec17.pdf
    //
    // Kaiser window and beta parameter
    //
    //         | 0.1102*(A - 8.7)                         A > 50
    //  beta = | 0.5842*(A - 21)^0.4 + 0.07886*(A - 21)   21 <= A <= 50
    //         | 0.                                       A < 21
    //
    // with A is the desired stop-band attenuation in dBFS
    //
    //    30 dB    2.210
    //    40 dB    3.384
    //    50 dB    4.538
    //    60 dB    5.658
    //    70 dB    6.764
    //    80 dB    7.865
    //    90 dB    8.960
    //   100 dB   10.056

    const int N = L * halfNumCoef; // non-negative half
    const double beta = 0.1102 * (stopBandAtten - 8.7); // >= 50dB always
    const double beta = computeBeta(stopBandAtten);
    const double xstep = (2. * M_PI) * fcr / L;
    const double xfrac = 1. / N;
    const double yscale = atten * L / (I0(beta) * M_PI);
@@ -696,9 +745,9 @@ static inline void firKaiserGen(T* coef, int L, int halfNumCoef,
                sg.advance();
            }

            if (is_same<T, int16_t>::value) { // int16_t needs noise shaping
            if (std::is_same<T, int16_t>::value) { // int16_t needs noise shaping
                *coef++ = static_cast<T>(toint(y, 1ULL<<(sizeof(T)*8-1), err));
            } else if (is_same<T, int32_t>::value) {
            } else if (std::is_same<T, int32_t>::value) {
                *coef++ = static_cast<T>(toint(y, 1ULL<<(sizeof(T)*8-1)));
            } else { // assumed float or double
                *coef++ = static_cast<T>(y);
+2 −2
Original line number Diff line number Diff line
@@ -14,8 +14,8 @@ mm

echo "waiting for device"
adb root && adb wait-for-device remount
adb push $OUT/system/lib/libaudioresampler.so /system/lib
adb push $OUT/system/lib64/libaudioresampler.so /system/lib64
adb push $OUT/system/lib/libaudioprocessing.so /system/lib
adb push $OUT/system/lib64/libaudioprocessing.so /system/lib64
adb push $OUT/data/nativetest/resampler_tests/resampler_tests /data/nativetest/resampler_tests/resampler_tests
adb push $OUT/data/nativetest64/resampler_tests/resampler_tests /data/nativetest64/resampler_tests/resampler_tests

+84 −0
Original line number Diff line number Diff line
@@ -29,6 +29,7 @@
#include <unistd.h>

#include <iostream>
#include <memory>
#include <utility>
#include <vector>

@@ -37,6 +38,8 @@
#include <media/AudioBufferProvider.h>

#include <media/AudioResampler.h>
#include "../AudioResamplerDyn.h"
#include "../AudioResamplerFirGen.h"
#include "test_utils.h"

template <typename T>
@@ -242,6 +245,60 @@ void testStopbandDownconversion(size_t channels,
    delete resampler;
}

void testFilterResponse(
        size_t channels, unsigned inputFreq, unsigned outputFreq)
{
    // create resampler
    using ResamplerType = android::AudioResamplerDyn<float, float, float>;
    std::unique_ptr<ResamplerType> rdyn(
            static_cast<ResamplerType *>(
                    android::AudioResampler::create(
                            AUDIO_FORMAT_PCM_FLOAT,
                            channels,
                            outputFreq,
                            android::AudioResampler::DYN_HIGH_QUALITY)));
    rdyn->setSampleRate(inputFreq);

    // get design parameters
    const int phases = rdyn->getPhases();
    const int halfLength = rdyn->getHalfLength();
    const float *coefs = rdyn->getFilterCoefs();
    const double fcr = rdyn->getNormalizedCutoffFrequency();
    const double tbw = rdyn->getNormalizedTransitionBandwidth();
    const double attenuation = rdyn->getFilterAttenuation();
    const double stopbandDb = rdyn->getStopbandAttenuationDb();
    const double passbandDb = rdyn->getPassbandRippleDb();
    const double fp = fcr - tbw / 2;
    const double fs = fcr + tbw / 2;

    printf("inputFreq:%d outputFreq:%d design"
            " phases:%d halfLength:%d"
            " fcr:%lf fp:%lf fs:%lf tbw:%lf"
            " attenuation:%lf stopRipple:%.lf passRipple:%lf"
            "\n",
            inputFreq, outputFreq,
            phases, halfLength,
            fcr, fp, fs, tbw,
            attenuation, stopbandDb, passbandDb);

    // verify design parameters
    constexpr int32_t passSteps = 1000;
    double passMin, passMax, passRipple, stopMax, stopRipple;
    android::testFir(coefs, phases, halfLength, fp / phases, fs / phases,
            passSteps, phases * passSteps /* stopSteps */,
            passMin, passMax, passRipple,
            stopMax, stopRipple);
    printf("inputFreq:%d outputFreq:%d verify"
            " passMin:%lf passMax:%lf passRipple:%lf stopMax:%lf stopRipple:%lf"
            "\n",
            inputFreq, outputFreq,
            passMin, passMax, passRipple, stopMax, stopRipple);

    ASSERT_GT(stopRipple, 60.);  // enough stopband attenuation
    ASSERT_LT(passRipple, 0.2);  // small passband ripple
    ASSERT_GT(passMin, 0.99);    // we do not attenuate the signal (ideally 1.)
}

/* Buffer increment test
 *
 * We compare a reference output, where we consume and process the entire
@@ -484,3 +541,30 @@ TEST(audioflinger_resampler, stopbandresponse_float_multichannel) {
    }
}

TEST(audioflinger_resampler, filterresponse) {
    std::vector<int> inSampleRates{
        8000,
        11025,
        12000,
        16000,
        22050,
        24000,
        32000,
        44100,
        48000,
        88200,
        96000,
        176400,
        192000,
    };
    std::vector<int> outSampleRates{
        48000,
        96000,
    };

    for (int outSampleRate : outSampleRates) {
        for (int inSampleRate : inSampleRates) {
            testFilterResponse(2 /* channels */, inSampleRate, outSampleRate);
        }
    }
}