From f93747c44ab8e1c2d31b458ae1ffc0f85952b67c Mon Sep 17 00:00:00 2001 From: Jean-Marc Valin <Jean-Marc.Valin@csiro.au> Date: Wed, 5 Mar 2008 17:20:30 +1100 Subject: [PATCH] fixed-point: changed find_spectral_pitch() to use single-precision (16-bit) FFT. This involved adding kfft_single.[ch] that redefines kiss_fft a second time with a different prefix. All this is still a bit of a mess now. The mask had to be converted to 16-bit input, but we're still using floats to apply it. --- libcelt/Makefile.am | 8 ++-- libcelt/_kiss_fft_guts.h | 5 +++ libcelt/celt.c | 4 +- libcelt/kfft_single.c | 13 +++++++ libcelt/kfft_single.h | 51 ++++++++++++++++++++++++ libcelt/kiss_fft.c | 7 ++-- libcelt/kiss_fft.h | 4 +- libcelt/kiss_fftr.c | 6 ++- libcelt/pitch.c | 83 ++++++++++++++++++++++++++++++++++------ libcelt/pitch.h | 3 ++ libcelt/psy.c | 8 ++-- libcelt/psy.h | 2 +- 12 files changed, 166 insertions(+), 28 deletions(-) create mode 100644 libcelt/kfft_single.c create mode 100644 libcelt/kfft_single.h diff --git a/libcelt/Makefile.am b/libcelt/Makefile.am index 795c8fb5f..d84b68b67 100644 --- a/libcelt/Makefile.am +++ b/libcelt/Makefile.am @@ -14,8 +14,8 @@ noinst_SCRIPTS = match-test.sh lib_LTLIBRARIES = libcelt.la # Sources for compilation in the library -libcelt_la_SOURCES = bands.c celt.c cwrs.c \ - ecintrin.h entcode.c entdec.c entenc.c header.c kiss_fft.c \ +libcelt_la_SOURCES = bands.c celt.c cwrs.c ecintrin.h \ + entcode.c entdec.c entenc.c header.c kfft_single.c kiss_fft.c \ kiss_fftr.c laplace.c mdct.c modes.c pitch.c psy.c quant_bands.c \ quant_pitch.c rangedec.c rangeenc.c rate.c vq.c @@ -25,8 +25,8 @@ libcelt_la_LDFLAGS = -version-info @CELT_LT_CURRENT@:@CELT_LT_REVISION@:@CELT_LT noinst_HEADERS = _kiss_fft_guts.h arch.h bands.h \ cwrs.h ecintrin.h entcode.h entdec.h entenc.h fixed_generic.h \ - kiss_fft.h kiss_fftr.h laplace.h mdct.h mfrngcod.h mathops.h modes.h \ - os_support.h pgain_table.h pitch.h psy.h \ + kfft_single.h kiss_fft.h kiss_fftr.h laplace.h mdct.h mfrngcod.h \ + mathops.h modes.h os_support.h pgain_table.h pitch.h psy.h \ quant_bands.h quant_pitch.h rate.h stack_alloc.h vq.h noinst_PROGRAMS = testcelt diff --git a/libcelt/_kiss_fft_guts.h b/libcelt/_kiss_fft_guts.h index ddbe16c88..8cc6cba27 100644 --- a/libcelt/_kiss_fft_guts.h +++ b/libcelt/_kiss_fft_guts.h @@ -12,6 +12,9 @@ Redistribution and use in source and binary forms, with or without modification, THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. */ +#ifndef KISS_FFT_GUTS_H +#define KISS_FFT_GUTS_H + #define MIN(a,b) ((a)<(b) ? (a):(b)) #define MAX(a,b) ((a)>(b) ? (a):(b)) @@ -229,3 +232,5 @@ struct kiss_fft_state{ /* a debugging function */ #define pcpx(c)\ fprintf(stderr,"%g + %gi\n",(double)((c)->r),(double)((c)->i) ) + +#endif /* KISS_FFT_GUTS_H */ diff --git a/libcelt/celt.c b/libcelt/celt.c index e1c015f63..7cecb5bb1 100644 --- a/libcelt/celt.c +++ b/libcelt/celt.c @@ -103,7 +103,7 @@ CELTEncoder *celt_encoder_create(const CELTMode *mode) ec_byte_writeinit(&st->buf); ec_enc_init(&st->enc,&st->buf); - st->fft = kiss_fftr_alloc(MAX_PERIOD, 0, 0); + st->fft = pitch_state_alloc(MAX_PERIOD); psydecay_init(&st->psy, MAX_PERIOD/2, st->mode->Fs); st->in_mem = celt_alloc(N*C*sizeof(celt_sig_t)); @@ -130,7 +130,7 @@ void celt_encoder_destroy(CELTEncoder *st) ec_byte_writeclear(&st->buf); - kiss_fft_free(st->fft); + pitch_state_free(st->fft); psydecay_clear(&st->psy); celt_free(st->in_mem); diff --git a/libcelt/kfft_single.c b/libcelt/kfft_single.c new file mode 100644 index 000000000..9416d45d7 --- /dev/null +++ b/libcelt/kfft_single.c @@ -0,0 +1,13 @@ +#ifdef HAVE_CONFIG_H +#include "config.h" +#endif + +#ifdef FIXED_POINT + +#include "kfft_single.h" + +#define SKIP_CONFIG_H +#include "kiss_fft.c" +#include "kiss_fftr.c" + +#endif diff --git a/libcelt/kfft_single.h b/libcelt/kfft_single.h new file mode 100644 index 000000000..7bc081bd8 --- /dev/null +++ b/libcelt/kfft_single.h @@ -0,0 +1,51 @@ +/* (C) 2008 Jean-Marc Valin, CSIRO +*/ +/* + Redistribution and use in source and binary forms, with or without + modification, are permitted provided that the following conditions + are met: + + - Redistributions of source code must retain the above copyright + notice, this list of conditions and the following disclaimer. + + - Redistributions in binary form must reproduce the above copyright + notice, this list of conditions and the following disclaimer in the + documentation and/or other materials provided with the distribution. + + - Neither the name of the Xiph.org Foundation nor the names of its + contributors may be used to endorse or promote products derived from + this software without specific prior written permission. + + THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS + ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT + LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR + A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE FOUNDATION OR + CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, + EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, + PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR + PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF + LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING + NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS + SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. +*/ + +#ifndef KFFT_SINGLE_H +#define KFFT_SINGLE_H + +#ifdef FIXED_POINT + +#ifdef DOUBLE_PRECISION +#undef DOUBLE_PRECISION +#endif + +#ifdef MIXED_PRECISION +#undef MIXED_PRECISION +#endif + +#endif /* FIXED_POINT */ + +#include "kiss_fft.h" +#include "kiss_fftr.h" +#include "_kiss_fft_guts.h" + +#endif /* KFFT_SINGLE_H */ diff --git a/libcelt/kiss_fft.c b/libcelt/kiss_fft.c index d8b769787..b07204455 100644 --- a/libcelt/kiss_fft.c +++ b/libcelt/kiss_fft.c @@ -15,9 +15,10 @@ Redistribution and use in source and binary forms, with or without modification, THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. */ - -#ifdef HAVE_CONFIG_H -#include "config.h" +#ifndef SKIP_CONFIG_H +# ifdef HAVE_CONFIG_H +# include "config.h" +# endif #endif #include "_kiss_fft_guts.h" diff --git a/libcelt/kiss_fft.h b/libcelt/kiss_fft.h index 568ec6f08..f71a83969 100644 --- a/libcelt/kiss_fft.h +++ b/libcelt/kiss_fft.h @@ -36,22 +36,24 @@ extern "C" { #ifdef DOUBLE_PRECISION # define kiss_fft_scalar celt_int32_t # define kiss_twiddle_scalar celt_int32_t +# define KF_SUFFIX _celt_double #else # define kiss_fft_scalar celt_int16_t # define kiss_twiddle_scalar celt_int16_t +# define KF_SUFFIX _celt_single #endif #else # ifndef kiss_fft_scalar /* default is float */ # define kiss_fft_scalar float # define kiss_twiddle_scalar float +# define KF_SUFFIX _celt_single # endif #endif /* This adds a suffix to all the kiss_fft functions so we can easily link with more than one copy of the fft */ -#define KF_SUFFIX _celt_single #define CAT_SUFFIX(a,b) a ## b #define SUF(a,b) CAT_SUFFIX(a, b) diff --git a/libcelt/kiss_fftr.c b/libcelt/kiss_fftr.c index 392f3eb68..1f63365ad 100644 --- a/libcelt/kiss_fftr.c +++ b/libcelt/kiss_fftr.c @@ -16,8 +16,10 @@ Redistribution and use in source and binary forms, with or without modification, THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. */ -#ifdef HAVE_CONFIG_H -#include "config.h" +#ifndef SKIP_CONFIG_H +# ifdef HAVE_CONFIG_H +# include "config.h" +# endif #endif #include "os_support.h" diff --git a/libcelt/pitch.c b/libcelt/pitch.c index c991c949d..c8b3aa092 100644 --- a/libcelt/pitch.c +++ b/libcelt/pitch.c @@ -39,19 +39,73 @@ #include "config.h" #endif +/*#include "_kiss_fft_guts.h" +#include "kiss_fftr.h"*/ +#include "kfft_single.h" + #include <stdio.h> #include <math.h> #include "pitch.h" #include "psy.h" -#include "_kiss_fft_guts.h" -#include "kiss_fftr.h" +#include "os_support.h" + +kiss_fftr_cfg pitch_state_alloc(int max_lag) +{ + return kiss_fftr_alloc_celt_single(max_lag, 0, 0); +} + +void pitch_state_free(kiss_fftr_cfg st) +{ + kiss_fft_free(st); +} + +#ifdef FIXED_POINT +static void normalise16(celt_word16_t *x, int len, celt_word16_t val) +{ + int i; + celt_word16_t maxval = 0; + for (i=0;i<len;i++) + maxval = MAX16(maxval, ABS16(x[i])); + if (maxval > val) + { + int shift = 0; + while (maxval > val) + { + maxval >>= 1; + shift++; + } + if (shift==0) + return; + for (i=0;i<len;i++) + x[i] = SHR16(x[i], shift); + } else { + int shift=0; + if (maxval == 0) + return; + val >>= 1; + while (maxval < val) + { + val >>= 1; + shift++; + } + if (shift==0) + return; + for (i=0;i<len;i++) + x[i] = SHL16(x[i], shift); + } +} +#else +#define normalise16(x,len,val) +#endif + +#define INPUT_SHIFT 15 void find_spectral_pitch(kiss_fftr_cfg fft, struct PsyDecay *decay, const celt_sig_t *x, const celt_sig_t *y, const celt_word16_t *window, int overlap, int lag, int len, int C, int *pitch) { int c, i; celt_word32_t max_corr; - VARDECL(celt_word32_t *X); - VARDECL(celt_word32_t *Y); + VARDECL(celt_word16_t *X); + VARDECL(celt_word16_t *Y); VARDECL(celt_mask_t *curve); int n2; int L2; @@ -59,7 +113,7 @@ void find_spectral_pitch(kiss_fftr_cfg fft, struct PsyDecay *decay, const celt_s SAVE_STACK; n2 = lag/2; L2 = len/2; - ALLOC(X, lag, celt_word32_t); + ALLOC(X, lag, celt_word16_t); ALLOC(curve, n2, celt_mask_t); bitrev = fft->substate->bitrev; @@ -70,8 +124,8 @@ void find_spectral_pitch(kiss_fftr_cfg fft, struct PsyDecay *decay, const celt_s { for (i=0;i<L2;i++) { - X[2*bitrev[i]] += SHR32(x[C*(2*i)+c],1); - X[2*bitrev[i]+1] += SHR32(x[C*(2*i+1)+c],1); + X[2*bitrev[i]] += SHR32(x[C*(2*i)+c],INPUT_SHIFT); + X[2*bitrev[i]+1] += SHR32(x[C*(2*i+1)+c],INPUT_SHIFT); } } /* Applying the window in the bit-reverse domain. It's a bit weird, but it @@ -83,6 +137,8 @@ void find_spectral_pitch(kiss_fftr_cfg fft, struct PsyDecay *decay, const celt_s X[2*bitrev[L2-i-1]] = MULT16_32_Q15(window[2*i+1], X[2*bitrev[L2-i-1]]); X[2*bitrev[L2-i-1]+1] = MULT16_32_Q15(window[2*i], X[2*bitrev[L2-i-1]+1]); } + normalise16(X, lag, 8192); + /*for (i=0;i<lag;i++) printf ("%d ", X[i]);printf ("\n");*/ /* Forward real FFT (in-place) */ kf_work((kiss_fft_cpx*)X, NULL, 1,1, fft->substate->factors,fft->substate, 1, 1, 1); kiss_fftr_twiddles(fft,X); @@ -90,7 +146,7 @@ void find_spectral_pitch(kiss_fftr_cfg fft, struct PsyDecay *decay, const celt_s compute_masking(decay, X, curve, lag); /* Deferred allocation to reduce peak stack usage */ - ALLOC(Y, lag, celt_word32_t); + ALLOC(Y, lag, celt_word16_t); for (i=0;i<lag;i++) Y[i] = 0; /* Sum all channels of the past audio and copy into Y in bit-reverse order */ @@ -98,10 +154,11 @@ void find_spectral_pitch(kiss_fftr_cfg fft, struct PsyDecay *decay, const celt_s { for (i=0;i<n2;i++) { - Y[2*bitrev[i]] += SHR32(y[C*(2*i)+c],1); - Y[2*bitrev[i]+1] += SHR32(y[C*(2*i+1)+c],1); + Y[2*bitrev[i]] += SHR32(y[C*(2*i)+c],INPUT_SHIFT); + Y[2*bitrev[i]+1] += SHR32(y[C*(2*i+1)+c],INPUT_SHIFT); } } + normalise16(Y, lag, 8192); /* Forward real FFT (in-place) */ kf_work((kiss_fft_cpx*)Y, NULL, 1,1, fft->substate->factors,fft->substate, 1, 1, 1); kiss_fftr_twiddles(fft,Y); @@ -113,7 +170,8 @@ void find_spectral_pitch(kiss_fftr_cfg fft, struct PsyDecay *decay, const celt_s celt_word32_t tmp; /*n = 1.f/(1e1+sqrt(sqrt((X[2*i-1]*X[2*i-1] + X[2*i ]*X[2*i ])*(Y[2*i-1]*Y[2*i-1] + Y[2*i ]*Y[2*i ]))));*/ /*n = 1;*/ - n = SIG_SCALING_1/sqrt(1+curve[i]); + /*printf ("%d %d ", X[2*i]*X[2*i]+X[2*i+1]*X[2*i+1], Y[2*i]*Y[2*i]+Y[2*i+1]*Y[2*i+1]);*/ + n = 1.f/sqrt(1+curve[i]); /*printf ("%f ", n);*/ /*n = 1.f/(1+curve[i]);*/ tmp = X[2*i]; @@ -122,8 +180,11 @@ void find_spectral_pitch(kiss_fftr_cfg fft, struct PsyDecay *decay, const celt_s } /*printf ("\n");*/ X[0] = X[1] = 0; + /*for (i=0;i<lag;i++) printf ("%d ", X[i]);printf ("\n");*/ + normalise16(X, lag, 50); /* Inverse half-complex to real FFT gives us the correlation */ kiss_fftri(fft, X, Y); + /*for (i=0;i<lag;i++) printf ("%d ", Y[i]);printf ("\n");*/ /*for (i=0;i<C*lag;i++) printf ("%d %d\n", X[i], xx[i]);*/ diff --git a/libcelt/pitch.h b/libcelt/pitch.h index fe916b90e..062f12621 100644 --- a/libcelt/pitch.h +++ b/libcelt/pitch.h @@ -41,6 +41,9 @@ #include "kiss_fftr.h" #include "psy.h" +kiss_fftr_cfg pitch_state_alloc(int max_lag); +void pitch_state_free(kiss_fftr_cfg st); + /** Find the optimal delay for the pitch prediction. Computation is done in the frequency domain, both to save time and to make it easier to apply psychoacoustic weighting */ diff --git a/libcelt/psy.c b/libcelt/psy.c index 294e2dce7..ed8ec62ec 100644 --- a/libcelt/psy.c +++ b/libcelt/psy.c @@ -126,7 +126,7 @@ static void spreading_func(struct PsyDecay *d, celt_word32_t *psd, celt_mask_t * } /* Compute a marking threshold from the spectrum X. */ -void compute_masking(struct PsyDecay *decay, celt_word32_t *X, celt_mask_t *mask, int len) +void compute_masking(struct PsyDecay *decay, celt_word16_t *X, celt_mask_t *mask, int len) { int i; VARDECL(celt_word32_t *psd); @@ -134,10 +134,10 @@ void compute_masking(struct PsyDecay *decay, celt_word32_t *X, celt_mask_t *mask SAVE_STACK; N=len/2; ALLOC(psd, N, celt_word32_t); - psd[0] = MULT16_16(ROUND(X[0],SIG_SHIFT), ROUND(X[0],SIG_SHIFT)); + psd[0] = MULT16_16(X[0], X[0]); for (i=1;i<N;i++) - psd[i] = ADD32(MULT16_16(ROUND(X[i*2],SIG_SHIFT), ROUND(X[i*2],SIG_SHIFT)), - MULT16_16(ROUND(X[i*2+1],SIG_SHIFT), ROUND(X[i*2+1],SIG_SHIFT))); + psd[i] = ADD32(MULT16_16(X[i*2], X[i*2]), + MULT16_16(X[i*2+1], X[i*2+1])); /* TODO: Do tone masking */ /* Noise masking */ spreading_func(decay, psd, mask, N); diff --git a/libcelt/psy.h b/libcelt/psy.h index 5180ce001..667bebd66 100644 --- a/libcelt/psy.h +++ b/libcelt/psy.h @@ -45,7 +45,7 @@ void psydecay_init(struct PsyDecay *decay, int len, celt_int32_t Fs); void psydecay_clear(struct PsyDecay *decay); /** Compute the masking curve for an input (DFT) spectrum X */ -void compute_masking(struct PsyDecay *decay, celt_word32_t *X, celt_mask_t *mask, int len); +void compute_masking(struct PsyDecay *decay, celt_word16_t *X, celt_mask_t *mask, int len); /** Compute the masking curve for an input (MDCT) spectrum X */ void compute_mdct_masking(struct PsyDecay *decay, celt_word32_t *X, celt_mask_t *mask, int len); -- GitLab