Changeset 01d4d19 for src


Ignore:
Timestamp:
Nov 15, 2018, 3:07:48 AM (13 months ago)
Author:
Paul Brossier <piem@piem.org>
Branches:
feature/autosink, feature/constantq, feature/pitchshift, feature/pydocstrings, feature/timestretch, master
Children:
5999ff2
Parents:
9ef3c6e
Message:

Merge branch 'fix/oddfft' (closes #207)

File:
1 edited

Legend:

Unmodified
Added
Removed
  • src/spectral/fft.c

    r9ef3c6e r01d4d19  
    9090#define aubio_vDSP_vsadd               vDSP_vsadd
    9191#define aubio_vDSP_vsmul               vDSP_vsmul
    92 #define aubio_vDSP_create_fftsetup     vDSP_create_fftsetup
    93 #define aubio_vDSP_destroy_fftsetup    vDSP_destroy_fftsetup
    9492#define aubio_DSPComplex               DSPComplex
    9593#define aubio_DSPSplitComplex          DSPSplitComplex
    96 #define aubio_FFTSetup                 FFTSetup
     94#define aubio_vDSP_DFT_Setup           vDSP_DFT_Setup
     95#define aubio_vDSP_DFT_zrop_CreateSetup vDSP_DFT_zrop_CreateSetup
     96#define aubio_vDSP_DFT_Execute         vDSP_DFT_Execute
     97#define aubio_vDSP_DFT_DestroySetup    vDSP_DFT_DestroySetup
    9798#define aubio_vvsqrt                   vvsqrtf
    9899#else
     
    104105#define aubio_vDSP_vsadd               vDSP_vsaddD
    105106#define aubio_vDSP_vsmul               vDSP_vsmulD
    106 #define aubio_vDSP_create_fftsetup     vDSP_create_fftsetupD
    107 #define aubio_vDSP_destroy_fftsetup    vDSP_destroy_fftsetupD
    108107#define aubio_DSPComplex               DSPDoubleComplex
    109108#define aubio_DSPSplitComplex          DSPDoubleSplitComplex
    110 #define aubio_FFTSetup                 FFTSetupD
     109#define aubio_vDSP_DFT_Setup           vDSP_DFT_SetupD
     110#define aubio_vDSP_DFT_zrop_CreateSetup vDSP_DFT_zrop_CreateSetupD
     111#define aubio_vDSP_DFT_Execute         vDSP_DFT_ExecuteD
     112#define aubio_vDSP_DFT_DestroySetup    vDSP_DFT_DestroySetupD
    111113#define aubio_vvsqrt                   vvsqrt
    112114#endif /* HAVE_AUBIO_DOUBLE */
     
    153155
    154156#elif defined HAVE_ACCELERATE  // using ACCELERATE
    155   int log2fftsize;
    156   aubio_FFTSetup fftSetup;
     157  aubio_vDSP_DFT_Setup fftSetupFwd;
     158  aubio_vDSP_DFT_Setup fftSetupBwd;
    157159  aubio_DSPSplitComplex spec;
    158160  smpl_t *in, *out;
     
    211213
    212214#elif defined HAVE_ACCELERATE  // using ACCELERATE
     215  {
     216    uint_t radix = winsize;
     217    uint_t order = 0;
     218    while ((radix / 2) * 2 == radix) {
     219      radix /= 2;
     220      order++;
     221    }
     222    if (order < 4 || (radix != 1 && radix != 3 && radix != 5 && radix != 15)) {
     223      AUBIO_ERR("fft: vDSP/Accelerate supports FFT with sizes = "
     224          "f * 2 ** n, where n > 4 and f in [1, 3, 5, 15], but requested %d. "
     225          "Use the closest power of two, or try recompiling aubio with "
     226          "--enable-fftw3.\n", winsize);
     227      goto beach;
     228    }
     229  }
    213230  s->winsize = winsize;
    214231  s->fft_size = winsize;
    215232  s->compspec = new_fvec(winsize);
    216   s->log2fftsize = aubio_power_of_two_order(s->fft_size);
    217233  s->in = AUBIO_ARRAY(smpl_t, s->fft_size);
    218234  s->out = AUBIO_ARRAY(smpl_t, s->fft_size);
    219235  s->spec.realp = AUBIO_ARRAY(smpl_t, s->fft_size/2);
    220236  s->spec.imagp = AUBIO_ARRAY(smpl_t, s->fft_size/2);
    221   s->fftSetup = aubio_vDSP_create_fftsetup(s->log2fftsize, FFT_RADIX2);
     237  s->fftSetupFwd = aubio_vDSP_DFT_zrop_CreateSetup(NULL,
     238      s->fft_size, vDSP_DFT_FORWARD);
     239  s->fftSetupBwd = aubio_vDSP_DFT_zrop_CreateSetup(s->fftSetupFwd,
     240      s->fft_size, vDSP_DFT_INVERSE);
    222241
    223242#elif defined HAVE_INTEL_IPP  // using Intel IPP
     
    293312  AUBIO_FREE(s->spec.realp);
    294313  AUBIO_FREE(s->spec.imagp);
    295   aubio_vDSP_destroy_fftsetup(s->fftSetup);
     314  aubio_vDSP_DFT_DestroySetup(s->fftSetupBwd);
     315  aubio_vDSP_DFT_DestroySetup(s->fftSetupFwd);
    296316
    297317#elif defined HAVE_INTEL_IPP  // using Intel IPP
     
    351371  aubio_vDSP_ctoz((aubio_DSPComplex*)s->in, 2, &s->spec, 1, s->fft_size/2);
    352372  // compute the FFT
    353   aubio_vDSP_fft_zrip(s->fftSetup, &s->spec, 1, s->log2fftsize, FFT_FORWARD);
     373  aubio_vDSP_DFT_Execute(s->fftSetupFwd, s->spec.realp, s->spec.imagp,
     374      s->spec.realp, s->spec.imagp);
    354375  // convert from vDSP complex split to [ r0, r1, ..., rN, iN-1, .., i2, i1]
    355376  compspec->data[0] = s->spec.realp[0];
     
    419440  aubio_vDSP_ctoz((aubio_DSPComplex*)s->out, 2, &s->spec, 1, s->fft_size/2);
    420441  // compute the FFT
    421   aubio_vDSP_fft_zrip(s->fftSetup, &s->spec, 1, s->log2fftsize, FFT_INVERSE);
     442  aubio_vDSP_DFT_Execute(s->fftSetupBwd, s->spec.realp, s->spec.imagp,
     443      s->spec.realp, s->spec.imagp);
    422444  // convert result to real output
    423445  aubio_vDSP_ztoc(&s->spec, 1, (aubio_DSPComplex*)output->data, 2, s->fft_size/2);
     
    494516  }
    495517#endif
    496   if (compspec->data[compspec->length/2] < 0) {
    497     spectrum->phas[spectrum->length - 1] = PI;
     518#ifdef HAVE_FFTW3
     519  // for even length only, make sure last element is 0 or PI
     520  if (2 * (compspec->length / 2) == compspec->length) {
     521#endif
     522    if (compspec->data[compspec->length/2] < 0) {
     523      spectrum->phas[spectrum->length - 1] = PI;
     524    } else {
     525      spectrum->phas[spectrum->length - 1] = 0.;
     526    }
     527#ifdef HAVE_FFTW3
    498528  } else {
    499     spectrum->phas[spectrum->length - 1] = 0.;
    500   }
     529    i = spectrum->length - 1;
     530    spectrum->phas[i] = ATAN2(compspec->data[compspec->length-i],
     531        compspec->data[i]);
     532  }
     533#endif
    501534}
    502535
     
    508541        + SQR(compspec->data[compspec->length - i]) );
    509542  }
    510   spectrum->norm[spectrum->length-1] =
    511     ABS(compspec->data[compspec->length/2]);
     543#ifdef HAVE_FFTW3
     544  // for even length, make sure last element is > 0
     545  if (2 * (compspec->length / 2) == compspec->length) {
     546#endif
     547    spectrum->norm[spectrum->length-1] =
     548      ABS(compspec->data[compspec->length/2]);
     549#ifdef HAVE_FFTW3
     550  } else {
     551    i = spectrum->length - 1;
     552    spectrum->norm[i] = SQRT(SQR(compspec->data[i])
     553        + SQR(compspec->data[compspec->length - i]) );
     554  }
     555#endif
    512556}
    513557
Note: See TracChangeset for help on using the changeset viewer.