Changeset aadd27a for src


Ignore:
Timestamp:
Nov 7, 2007, 4:30:29 PM (17 years ago)
Author:
Paul Brossier <piem@piem.org>
Branches:
feature/autosink, feature/cnn, feature/cnn_org, feature/constantq, feature/crepe, feature/crepe_org, feature/pitchshift, feature/pydocstrings, feature/timestretch, fix/ffmpeg5, master, pitchshift, sampler, timestretch, yinfft+
Children:
8b2dc90
Parents:
f1eda56
Message:

fft.c,fft.h: refactorised, no more mfft, only one fft

Location:
src
Files:
2 edited

Legend:

Unmodified
Added
Removed

  • src/fft.c

    rf1eda56 raadd27a  
    1919
    2020#include "aubio_priv.h"
    21 #include "sample.h"
     21#include "fvec.h"
     22#include "cvec.h"
    2223#include "mathutils.h"
    2324#include "fft.h"
    … …  
    4142
    4243struct _aubio_fft_t {
     44  uint_t winsize;
     45  uint_t channels;
    4346  uint_t fft_size;
    44   uint_t channels;
    45   real_t    *in, *out;
    46   fft_data_t   *specdata;
     47  real_t *in, *out;
    4748  fftw_plan   pfw, pbw;
     49  fft_data_t * specdata;     /* complex spectral data */
     50  fvec_t * compspec;
    4851};
    4952
    50 static void aubio_fft_getspectrum(fft_data_t * spectrum, smpl_t *norm, smpl_t * phas, uint_t size);
    51 
    52 aubio_fft_t * new_aubio_fft(uint_t size) {
     53aubio_fft_t * new_aubio_fft(uint_t winsize, uint_t channels) {
    5354  aubio_fft_t * s = AUBIO_NEW(aubio_fft_t);
     55  s->winsize  = winsize;
     56  s->channels = channels;
    5457  /* allocate memory */
    55   s->in       = AUBIO_ARRAY(real_t,size);
    56   s->out      = AUBIO_ARRAY(real_t,size);
     58  s->in       = AUBIO_ARRAY(real_t,winsize);
     59  s->out      = AUBIO_ARRAY(real_t,winsize);
     60  s->compspec = new_fvec(winsize,channels);
    5761  /* create plans */
    5862#ifdef HAVE_COMPLEX_H
    59   s->fft_size = size/2+1;
     63  s->fft_size = winsize/2+1;
    6064  s->specdata = (fft_data_t*)fftw_malloc(sizeof(fft_data_t)*s->fft_size);
    61   s->pfw = fftw_plan_dft_r2c_1d(size, s->in,  s->specdata, FFTW_ESTIMATE);
    62   s->pbw = fftw_plan_dft_c2r_1d(size, s->specdata, s->out, FFTW_ESTIMATE);
     65  s->pfw = fftw_plan_dft_r2c_1d(winsize, s->in,  s->specdata, FFTW_ESTIMATE);
     66  s->pbw = fftw_plan_dft_c2r_1d(winsize, s->specdata, s->out, FFTW_ESTIMATE);
    6367#else
    64   s->fft_size = size;
     68  s->fft_size = winsize;
    6569  s->specdata = (fft_data_t*)fftw_malloc(sizeof(fft_data_t)*s->fft_size);
    66   s->pfw = fftw_plan_r2r_1d(size, s->in,  s->specdata, FFTW_R2HC, FFTW_ESTIMATE);
    67   s->pbw = fftw_plan_r2r_1d(size, s->specdata, s->out, FFTW_HC2R, FFTW_ESTIMATE);
     70  s->pfw = fftw_plan_r2r_1d(winsize, s->in,  s->specdata, FFTW_R2HC, FFTW_ESTIMATE);
     71  s->pbw = fftw_plan_r2r_1d(winsize, s->specdata, s->out, FFTW_HC2R, FFTW_ESTIMATE);
    6872#endif
    6973  return s;
    … …  
    7276void del_aubio_fft(aubio_fft_t * s) {
    7377  /* destroy data */
     78  del_fvec(s->compspec);
    7479  fftw_destroy_plan(s->pfw);
    7580  fftw_destroy_plan(s->pbw);
    … …  
    8085}
    8186
    82 void aubio_fft_do(const aubio_fft_t * s,
    83     const smpl_t * data, fft_data_t * spectrum, const uint_t size) {
    84   uint_t i;
    85   for (i=0;i<size;i++) s->in[i] = data[i];
    86   fftw_execute(s->pfw);
    87   for (i=0; i < s->fft_size; i++) spectrum[i] = s->specdata[i];
     87void aubio_fft_do(aubio_fft_t * s, fvec_t * input, cvec_t * spectrum) {
     88  aubio_fft_do_complex(s, input, s->compspec);
     89  aubio_fft_get_spectrum(s->compspec, spectrum);
    8890}
    8991
    90 void aubio_fft_rdo(const aubio_fft_t * s,
    91     const fft_data_t * spectrum, smpl_t * data, const uint_t size) {
    92   uint_t i;
    93   const smpl_t renorm = 1./(smpl_t)size;
    94   for (i=0; i < s->fft_size; i++) s->specdata[i] = spectrum[i];
    95   fftw_execute(s->pbw);
    96   for (i=0;i<size;i++) data[i] = s->out[i]*renorm;
     92void aubio_fft_rdo(aubio_fft_t * s, cvec_t * spectrum, fvec_t * output) {
     93  aubio_fft_get_realimag(spectrum, s->compspec);
     94  aubio_fft_rdo_complex(s, s->compspec, output);
    9795}
    9896
    99 #ifdef HAVE_COMPLEX_H
    100 
    101 void aubio_fft_getnorm(smpl_t * norm, fft_data_t * spectrum, uint_t size) {
    102   uint_t i;
    103   for (i=0;i<size/2+1;i++) norm[i] = ABSC(spectrum[i]);
    104 }
    105 
    106 void aubio_fft_getphas(smpl_t * phas, fft_data_t * spectrum, uint_t size) {
    107   uint_t i;
    108   for (i=0;i<size/2+1;i++) phas[i] = ARGC(spectrum[i]);
    109 }
    110 
    111 void aubio_fft_getspectrum(fft_data_t * spectrum, smpl_t *norm, smpl_t * phas, uint_t size) {
    112   uint_t j;
    113   for (j=0; j<size/2+1; j++) {
    114     spectrum[j]  = CEXPC(I*phas[j]);
    115     spectrum[j] *= norm[j];
     97void aubio_fft_do_complex(aubio_fft_t * s, fvec_t * input, fvec_t * compspec) {
     98  uint_t i, j;
     99  for (i = 0; i < s->channels; i++) {
     100    for (j=0; j < s->winsize; j++) {
     101      s->in[j] = input->data[i][j];
     102    }
     103    fftw_execute(s->pfw);
     104#if HAVE_COMPLEX_H
     105    compspec->data[i][0] = REAL(s->specdata[0]);
     106    for (j = 1; j < s->fft_size -1 ; j++) {
     107      compspec->data[i][j] = REAL(s->specdata[j]);
     108      compspec->data[i][compspec->length - j] = IMAG(s->specdata[j]);
     109    }
     110    compspec->data[i][s->fft_size-1] = REAL(s->specdata[s->fft_size-1]);
     111#else
     112    for (j = 0; j < s->fft_size; j++) {
     113      compspec->data[i][j] = s->specdata[j];
     114    }
     115#endif
    116116  }
    117117}
    118118
     119void aubio_fft_rdo_complex(aubio_fft_t * s, fvec_t * compspec, fvec_t * output) {
     120  uint_t i, j;
     121  const smpl_t renorm = 1./(smpl_t)s->winsize;
     122  for (i = 0; i < compspec->channels; i++) {
     123#if HAVE_COMPLEX_H
     124    s->specdata[0] = compspec->data[i][0];
     125    for (j=1; j < s->fft_size - 1; j++) {
     126      s->specdata[j] = compspec->data[i][j] +
     127        I * compspec->data[i][compspec->length - j];
     128    }
     129    s->specdata[s->fft_size - 1] = compspec->data[i][s->fft_size - 1];
    119130#else
    120 
    121 void aubio_fft_getnorm(smpl_t * norm, fft_data_t * spectrum, uint_t size) {
    122   uint_t i;
    123   norm[0] = spectrum[0];
    124   for (i=1;i<size/2;i++) norm[i] = SQRT((SQR(spectrum[i]) + SQR(spectrum[size-i])));
    125   norm[size/2] = spectrum[size/2];
    126 }
    127 
    128 void aubio_fft_getphas(smpl_t * phas, fft_data_t * spectrum, uint_t size) {
    129   uint_t i;
    130   phas[0] = 0;
    131   for (i=1;i<size/2+1;i++) phas[i] = atan2f(spectrum[size-i] , spectrum[i]);
    132   phas[size/2] = 0;
    133 }
    134 
    135 void aubio_fft_getspectrum(fft_data_t * spectrum, smpl_t *norm, smpl_t * phas, uint_t size) {
    136   uint_t j;
    137   for (j=0; j<size/2+1; j++) {
    138     spectrum[j]       = norm[j]*COS(phas[j]);
    139   }
    140   for (j=1; j<size/2+1; j++) {
    141     spectrum[size-j]  = norm[j]*SIN(phas[j]);
     131    for (j=0; j < s->fft_size; j++) {
     132      s->specdata[j] = compspec->data[i][j];
     133    }
     134#endif
     135    fftw_execute(s->pbw);
     136    for (j = 0; j < output->length; j++) {
     137      output->data[i][j] = s->out[j]*renorm;
     138    }
    142139  }
    143140}
    144141
    145 #endif
    146 
    147 /* new interface aubio_mfft */
    148 struct _aubio_mfft_t {
    149         aubio_fft_t * fft;      /* fftw interface */
    150         fft_data_t ** spec;     /* complex spectral data */
    151         uint_t winsize;
    152         uint_t channels;
    153 };
    154 
    155 aubio_mfft_t * new_aubio_mfft(uint_t winsize, uint_t channels){
    156   uint_t i;
    157   aubio_mfft_t * fft = AUBIO_NEW(aubio_mfft_t);
    158   fft->winsize       = winsize;
    159   fft->channels      = channels;
    160   fft->fft           = new_aubio_fft(winsize);
    161   fft->spec          = AUBIO_ARRAY(fft_data_t*,channels);
    162   for (i=0; i < channels; i++)
    163     fft->spec[i] = AUBIO_ARRAY(fft_data_t,winsize);
    164   return fft;
     142void aubio_fft_get_spectrum(fvec_t * compspec, cvec_t * spectrum) {
     143  aubio_fft_get_phas(compspec, spectrum);
     144  aubio_fft_get_norm(compspec, spectrum);
    165145}
    166146
    167 /* execute stft */
    168 void aubio_mfft_do (aubio_mfft_t * fft,fvec_t * in,cvec_t * fftgrain){
    169   uint_t i=0;
    170   /* execute stft */
    171   for (i=0; i < fft->channels; i++) {
    172     aubio_fft_do (fft->fft,in->data[i],fft->spec[i],fft->winsize);
    173     /* put norm and phase into fftgrain */
    174     aubio_fft_getnorm(fftgrain->norm[i], fft->spec[i], fft->winsize);
    175     aubio_fft_getphas(fftgrain->phas[i], fft->spec[i], fft->winsize);
     147void aubio_fft_get_realimag(cvec_t * spectrum, fvec_t * compspec) {
     148  aubio_fft_get_imag(spectrum, compspec);
     149  aubio_fft_get_real(spectrum, compspec);
     150}
     151
     152void aubio_fft_get_phas(fvec_t * compspec, cvec_t * spectrum) {
     153  uint_t i, j;
     154  for (i = 0; i < spectrum->channels; i++) {
     155    spectrum->phas[i][0] = 0.;
     156    for (j=1; j < spectrum->length - 1; j++) {
     157      spectrum->phas[i][j] = atan2f(compspec->data[i][compspec->length-j],
     158          compspec->data[i][j]);
     159    }
     160    spectrum->phas[i][spectrum->length-1] = 0.;
    176161  }
    177162}
    178163
    179 /* execute inverse fourier transform */
    180 void aubio_mfft_rdo(aubio_mfft_t * fft,cvec_t * fftgrain, fvec_t * out){
    181   uint_t i=0;
    182   for (i=0; i < fft->channels; i++) {
    183     aubio_fft_getspectrum(fft->spec[i],fftgrain->norm[i],fftgrain->phas[i],fft->winsize);
    184     aubio_fft_rdo(fft->fft,fft->spec[i],out->data[i],fft->winsize);
     164void aubio_fft_get_norm(fvec_t * compspec, cvec_t * spectrum) {
     165  uint_t i, j = 0;
     166  for (i = 0; i < spectrum->channels; i++) {
     167    spectrum->norm[i][0] = compspec->data[i][0];
     168    for (j=1; j < spectrum->length - 1; j++) {
     169      spectrum->norm[i][j] = SQRT(SQR(compspec->data[i][j])
     170          + SQR(compspec->data[i][compspec->length - j]) );
     171    }
     172    spectrum->norm[i][spectrum->length-1] = compspec->data[i][compspec->length/2];
    185173  }
    186174}
    187175
    188 void del_aubio_mfft(aubio_mfft_t * fft) {
    189   uint_t i;
    190   for (i=0; i < fft->channels; i++)
    191     AUBIO_FREE(fft->spec[i]);
    192   AUBIO_FREE(fft->spec);
    193   del_aubio_fft(fft->fft);
    194   AUBIO_FREE(fft);       
     176void aubio_fft_get_imag(cvec_t * spectrum, fvec_t * compspec) {
     177  uint_t i, j;
     178  for (i = 0; i < compspec->channels; i++) {
     179    for (j = 1; j < compspec->length / 2 + 1; j++) {
     180      compspec->data[i][compspec->length - j] =
     181        spectrum->norm[i][j]*SIN(spectrum->phas[i][j]);
     182    }
     183  }
    195184}
     185
     186void aubio_fft_get_real(cvec_t * spectrum, fvec_t * compspec) {
     187  uint_t i, j;
     188  for (i = 0; i < compspec->channels; i++) {
     189    for (j = 0; j< compspec->length / 2 + 1; j++) {
     190      compspec->data[i][j] =
     191        spectrum->norm[i][j]*COS(spectrum->phas[i][j]);
     192    }
     193  }
     194}

  • src/fft.h

    rf1eda56 raadd27a  
    6767
    6868  \param size length of the FFT
     69  \param channels number of channels
    6970
    7071*/
    71 aubio_fft_t * new_aubio_fft(uint_t size);
     72aubio_fft_t * new_aubio_fft(uint_t size, uint_t channels);
    7273/** delete FFT object
    7374
    … …  
    7677*/
    7778void del_aubio_fft(aubio_fft_t * s);
     79
    7880/** compute forward FFT
    7981
    8082  \param s fft object as returned by new_aubio_fft
    81   \param data input signal
     83  \param input input signal
    8284  \param spectrum output spectrum
    83   \param size length of the input vector
    8485
    8586*/
    86 void aubio_fft_do (const aubio_fft_t *s, const smpl_t * data,
    87     fft_data_t * spectrum, const uint_t size);
     87void aubio_fft_do (aubio_fft_t *s, fvec_t * input, cvec_t * spectrum);
    8888/** compute backward (inverse) FFT
    8989
    9090  \param s fft object as returned by new_aubio_fft
    9191  \param spectrum input spectrum
    92   \param data output signal
    93   \param size length of the input vector
     92  \param output output signal
    9493
    9594*/
    96 void aubio_fft_rdo(const aubio_fft_t *s, const fft_data_t * spectrum,
    97     smpl_t * data, const uint_t size);
    98 /** compute norm vector from input spectrum
     95void aubio_fft_rdo (aubio_fft_t *s, cvec_t * spectrum, fvec_t * output);
    9996
    100   \param norm magnitude vector output
    101   \param spectrum spectral data input
    102   \param size size of the vectors
     97/** compute forward FFT
     98
     99  \param s fft object as returned by new_aubio_fft
     100  \param input real input signal
     101  \param compspec complex output fft real/imag
    103102
    104103*/
    105 void aubio_fft_getnorm(smpl_t * norm, fft_data_t * spectrum, uint_t size);
    106 /** compute phase vector from input spectrum
    107  
    108   \param phase phase vector output
    109   \param spectrum spectral data input
    110   \param size size of the vectors
     104void aubio_fft_do_complex (aubio_fft_t *s, fvec_t * input, fvec_t * compspec);
     105/** compute backward (inverse) FFT from real/imag
     106
     107  \param s fft object as returned by new_aubio_fft
     108  \param compspec real/imag input fft array
     109  \param output real output array
    111110
    112111*/
    113 void aubio_fft_getphas(smpl_t * phase, fft_data_t * spectrum, uint_t size);
     112void aubio_fft_rdo_complex (aubio_fft_t *s, fvec_t * compspec, fvec_t * output);
    114113
    115 /** FFT object (using cvec)
     114/** convert real/imag spectrum to norm/phas spectrum
    116115
    117   This object works similarly as aubio_fft_t, except the spectral data is
    118   stored in a cvec_t as two vectors, magnitude and phase.
     116  \param compspec real/imag input fft array
     117  \param spectrum cvec norm/phas output array
    119118
    120119*/
    121 typedef struct _aubio_mfft_t aubio_mfft_t;
     120void aubio_fft_get_spectrum(fvec_t * compspec, cvec_t * spectrum);
     121/** convert real/imag spectrum to norm/phas spectrum
    122122
    123 /** create new FFT computation object
    124 
    125   \param winsize length of the FFT
    126   \param channels number of channels
     123  \param compspec real/imag input fft array
     124  \param spectrum cvec norm/phas output array
    127125
    128126*/
    129 aubio_mfft_t * new_aubio_mfft(uint_t winsize, uint_t channels);
    130 /** compute forward FFT
     127void aubio_fft_get_realimag(cvec_t * spectrum, fvec_t * compspec);
    131128
    132   \param fft fft object as returned by new_aubio_mfft
    133   \param in input signal
    134   \param fftgrain output spectrum
     129/** compute phas spectrum from real/imag parts
     130
     131  \param compspec real/imag input fft array
     132  \param spectrum cvec norm/phas output array
    135133
    136134*/
    137 void aubio_mfft_do (aubio_mfft_t * fft,fvec_t * in,cvec_t * fftgrain);
    138 /** compute backward (inverse) FFT
     135void aubio_fft_get_phas(fvec_t * compspec, cvec_t * spectrum);
     136/** compute imaginary part from the norm/phas cvec
    139137
    140   \param fft fft object as returned by new_aubio_mfft
    141   \param fftgrain input spectrum (cvec)
    142   \param out output signal
     138  \param spectrum norm/phas input array
     139  \param compspec real/imag output fft array
    143140
    144141*/
    145 void aubio_mfft_rdo(aubio_mfft_t * fft,cvec_t * fftgrain, fvec_t * out);
    146 /** delete FFT object
     142void aubio_fft_get_imag(cvec_t * spectrum, fvec_t * compspec);
    147143
    148   \param fft fft object as returned by new_aubio_mfft
     144/** compute norm component from real/imag parts
     145
     146  \param compspec real/imag input fft array
     147  \param spectrum cvec norm/phas output array
    149148
    150149*/
    151 void del_aubio_mfft(aubio_mfft_t * fft);
     150void aubio_fft_get_norm(fvec_t * compspec, cvec_t * spectrum);
     151/** compute real part from norm/phas components
    152152
     153  \param spectrum norm/phas input array
     154  \param compspec real/imag output fft array
     155
     156*/
     157void aubio_fft_get_real(cvec_t * spectrum, fvec_t * compspec);
    153158
    154159#ifdef __cplusplus
    … …  
    156161#endif
    157162
    158 #endif
     163#endif // FFT_H_
Note: See TracChangeset for help on using the changeset viewer.