/* * Copyright(c)1995,97 Mark Olesen <[email protected]> * Queen's Univ at Kingston (Canada) * * Permission to use, copy, modify, and distribute this software for * any purpose without fee is hereby granted, provided that this * entire notice is included in all copies of any software which is * or includes a copy or modification of this software and in all * copies of the supporting documentation for such software. * * THIS SOFTWARE IS BEING PROVIDED "AS IS", WITHOUT ANY EXPRESS OR * IMPLIED WARRANTY. IN PARTICULAR, NEITHER THE AUTHOR NOR QUEEN'S * UNIVERSITY AT KINGSTON MAKES ANY REPRESENTATION OR WARRANTY OF ANY * KIND CONCERNING THE MERCHANTABILITY OF THIS SOFTWARE OR ITS * FITNESS FOR ANY PARTICULAR PURPOSE. * * All of which is to say that you can do what you like with this * source code provided you don't try to sell it as your own and you * include an unaltered copy of this message (including the * copyright). * * It is also implicitly understood that bug fixes and improvements * should make their way back to the general Internet community so * that everyone benefits. * * Changes: * Trivial type modifications by the WebRTC authors. */ /* * File: * WebRtcIsac_Fftn.c * * Public: * WebRtcIsac_Fftn / fftnf (); * * Private: * WebRtcIsac_Fftradix / fftradixf (); * * Descript: * multivariate complex Fourier transform, computed in place * using mixed-radix Fast Fourier Transform algorithm. * * Fortran code by: * RC Singleton, Stanford Research Institute, Sept. 1968 * * translated by f2c (version 19950721). * * int WebRtcIsac_Fftn (int ndim, const int dims[], REAL Re[], REAL Im[], * int iSign, double scaling); * * NDIM = the total number dimensions * DIMS = a vector of array sizes * if NDIM is zero then DIMS must be zero-terminated * * RE and IM hold the real and imaginary components of the data, and return * the resulting real and imaginary Fourier coefficients. Multidimensional * data *must* be allocated contiguously. There is no limit on the number * of dimensions. * * ISIGN = the sign of the complex exponential (ie, forward or inverse FFT) * the magnitude of ISIGN (normally 1) is used to determine the * correct indexing increment (see below). * * SCALING = normalizing constant by which the final result is *divided* * if SCALING == -1, normalize by total dimension of the transform * if SCALING < -1, normalize by the square-root of the total dimension * * example: * tri-variate transform with Re[n1][n2][n3], Im[n1][n2][n3] * * int dims[3] = {n1,n2,n3} * WebRtcIsac_Fftn (3, dims, Re, Im, 1, scaling); * *-----------------------------------------------------------------------* * int WebRtcIsac_Fftradix (REAL Re[], REAL Im[], size_t nTotal, size_t nPass, * size_t nSpan, int iSign, size_t max_factors, * size_t max_perm); * * RE, IM - see above documentation * * Although there is no limit on the number of dimensions, WebRtcIsac_Fftradix() must * be called once for each dimension, but the calls may be in any order. * * NTOTAL = the total number of complex data values * NPASS = the dimension of the current variable * NSPAN/NPASS = the spacing of consecutive data values while indexing the * current variable * ISIGN - see above documentation * * example: * tri-variate transform with Re[n1][n2][n3], Im[n1][n2][n3] * * WebRtcIsac_Fftradix (Re, Im, n1*n2*n3, n1, n1, 1, maxf, maxp); * WebRtcIsac_Fftradix (Re, Im, n1*n2*n3, n2, n1*n2, 1, maxf, maxp); * WebRtcIsac_Fftradix (Re, Im, n1*n2*n3, n3, n1*n2*n3, 1, maxf, maxp); * * single-variate transform, * NTOTAL = N = NSPAN = (number of complex data values), * * WebRtcIsac_Fftradix (Re, Im, n, n, n, 1, maxf, maxp); * * The data can also be stored in a single array with alternating real and * imaginary parts, the magnitude of ISIGN is changed to 2 to give correct * indexing increment, and data [0] and data [1] used to pass the initial * addresses for the sequences of real and imaginary values, * * example: * REAL data [2*NTOTAL]; * WebRtcIsac_Fftradix ( &data[0], &data[1], NTOTAL, nPass, nSpan, 2, maxf, maxp); * * for temporary allocation: * * MAX_FACTORS >= the maximum prime factor of NPASS * MAX_PERM >= the number of prime factors of NPASS. In addition, * if the square-free portion K of NPASS has two or more prime * factors, then MAX_PERM >= (K-1) * * storage in FACTOR for a maximum of 15 prime factors of NPASS. if NPASS * has more than one square-free factor, the product of the square-free * factors must be <= 210 array storage for maximum prime factor of 23 the * following two constants should agree with the array dimensions. * *----------------------------------------------------------------------*/ #include <stdlib.h> #include <math.h> #include "modules/third_party/fft/fft.h" /* double precision routine */ static int WebRtcIsac_Fftradix (double Re[], double Im[], size_t nTotal, size_t nPass, size_t nSpan, int isign, int max_factors, unsigned int max_perm, FFTstr *fftstate); #ifndef M_PI #define M_PI … #endif #ifndef SIN60 #define SIN60 … #define COS72 … #define SIN72 … #endif #define REAL … #define FFTN … #define FFTNS … #define FFTRADIX … #define FFTRADIXS … int WebRtcIsac_Fftns(unsigned int ndim, const int dims[], double Re[], double Im[], int iSign, double scaling, FFTstr *fftstate) { … } /* * singleton's mixed radix routine * * could move allocation out to WebRtcIsac_Fftn(), but leave it here so that it's * possible to make this a standalone function */ static int FFTRADIX (REAL Re[], REAL Im[], size_t nTotal, size_t nPass, size_t nSpan, int iSign, int max_factors, unsigned int max_perm, FFTstr *fftstate) { … } /* ---------------------- end-of-file (c source) ---------------------- */
Codebase Browser
chromium