diff --git a/Base/Utils/MathFunctions.cpp b/Base/Utils/MathFunctions.cpp index cf92b20107d785af49e8b4cdcd838a0e727fb2a6..c551a5799979f66b6ce7766f0775f977f385ad18 100644 --- a/Base/Utils/MathFunctions.cpp +++ b/Base/Utils/MathFunctions.cpp @@ -285,95 +285,6 @@ complex_t MathFunctions::Bessel_J1_PowSer(const complex_t z) return cj1; } -// ************************************************************************** // -// Fourier transform and convolution -// ************************************************************************** // - -//! @brief simple (and unoptimized) wrapper function -//! for the discrete fast Fourier transformation library (fftw3) - -std::vector<complex_t> MathFunctions::FastFourierTransform(const std::vector<complex_t>& data, - MathFunctions::EFFTDirection ftCase) -{ - double scale(1.); - size_t npx = data.size(); - - fftw_complex* ftData = (fftw_complex*)fftw_malloc(sizeof(fftw_complex) * npx); - fftw_complex* ftResult = (fftw_complex*)fftw_malloc(sizeof(fftw_complex) * npx); - memset(ftData, 0, sizeof(fftw_complex) * npx); - memset(ftResult, 0, sizeof(fftw_complex) * npx); - - for (size_t i = 0; i < npx; i++) { - ftData[i][0] = data[i].real(); - ftData[i][1] = data[i].imag(); - } - - fftw_plan plan; - switch (ftCase) { - case MathFunctions::FORWARD_FFT: - plan = fftw_plan_dft_1d((int)npx, ftData, ftResult, FFTW_FORWARD, FFTW_ESTIMATE); - break; - case MathFunctions::BACKWARD_FFT: - plan = fftw_plan_dft_1d((int)npx, ftData, ftResult, FFTW_BACKWARD, FFTW_ESTIMATE); - scale = 1. / double(npx); - break; - default: - throw std::runtime_error( - "MathFunctions::FastFourierTransform -> Panic! Unknown transform case."); - } - - fftw_execute(plan); - - // saving data for user - std::vector<complex_t> outData; - outData.resize(npx); - for (size_t i = 0; i < npx; i++) - outData[i] = scale * complex_t(ftResult[i][0], ftResult[i][1]); - - fftw_destroy_plan(plan); - fftw_free(ftData); - fftw_free(ftResult); - - return outData; -} - -//! @brief simple (and unoptimized) wrapper function -//! for the discrete fast Fourier transformation library (fftw3); -//! transforms real to complex - -std::vector<complex_t> MathFunctions::FastFourierTransform(const std::vector<double>& data, - MathFunctions::EFFTDirection ftCase) -{ - std::vector<complex_t> cdata; - cdata.resize(data.size()); - for (size_t i = 0; i < data.size(); i++) - cdata[i] = complex_t(data[i], 0); - return MathFunctions::FastFourierTransform(cdata, ftCase); -} - -//! convolution of two real vectors of equal size - -std::vector<complex_t> MathFunctions::ConvolveFFT(const std::vector<double>& signal, - const std::vector<double>& resfunc) -{ - if (signal.size() != resfunc.size()) - throw std::runtime_error("MathFunctions::ConvolveFFT() -> This convolution works only for " - "two vectors of equal size. Use Convolve class instead."); - std::vector<complex_t> fft_signal = - MathFunctions::FastFourierTransform(signal, MathFunctions::FORWARD_FFT); - std::vector<complex_t> fft_resfunc = - MathFunctions::FastFourierTransform(resfunc, MathFunctions::FORWARD_FFT); - - std::vector<complex_t> fft_prod; - fft_prod.resize(fft_signal.size()); - for (size_t i = 0; i < fft_signal.size(); i++) - fft_prod[i] = fft_signal[i] * fft_resfunc[i]; - - std::vector<complex_t> result = - MathFunctions::FastFourierTransform(fft_prod, MathFunctions::BACKWARD_FFT); - return result; -} - // ************************************************************************** // // Random number generators // ************************************************************************** // diff --git a/Base/Utils/MathFunctions.h b/Base/Utils/MathFunctions.h index 291581999943d00b7fe0ca3e789a41d64b848f73..ffd60db19057243b6dc07d15648cdd6efc501e4b 100644 --- a/Base/Utils/MathFunctions.h +++ b/Base/Utils/MathFunctions.h @@ -74,22 +74,6 @@ complex_t Bessel_J1(const complex_t z); //! Complex Bessel function Bessel_J1(x)/x complex_t Bessel_J1c(const complex_t z); -// ************************************************************************** // -// Fourier transform and convolution -// ************************************************************************** // - -// TODO move elsewhere, and rm #include <vector> - -enum EFFTDirection { FORWARD_FFT, BACKWARD_FFT }; - -// TODO: name these two functions differently (SWIG warning 509) -std::vector<complex_t> FastFourierTransform(const std::vector<complex_t>& data, - EFFTDirection tcase); -std::vector<complex_t> FastFourierTransform(const std::vector<double>& data, EFFTDirection tcase); - -std::vector<complex_t> ConvolveFFT(const std::vector<double>& signal, - const std::vector<double>& resfunc); - // ************************************************************************** // // Random number generators // ************************************************************************** // diff --git a/auto/Wrap/doxygenBase.i b/auto/Wrap/doxygenBase.i index 67110f8695246d8a0b098abcfa6cfcf819b24c1f..0070500826cf625f40c73d3c14ff85582fb5fb03 100644 --- a/auto/Wrap/doxygenBase.i +++ b/auto/Wrap/doxygenBase.i @@ -1335,21 +1335,6 @@ Complex Bessel function of the first kind and order 1. Complex Bessel function Bessel_J1(x)/x. "; -%feature("docstring") MathFunctions::FastFourierTransform "std::vector< complex_t > MathFunctions::FastFourierTransform(const std::vector< complex_t > &data, EFFTDirection tcase) - -simple (and unoptimized) wrapper function for the discrete fast Fourier transformation library (fftw3) -"; - -%feature("docstring") MathFunctions::FastFourierTransform "std::vector< complex_t > MathFunctions::FastFourierTransform(const std::vector< double > &data, EFFTDirection tcase) - -simple (and unoptimized) wrapper function for the discrete fast Fourier transformation library (fftw3); transforms real to complex -"; - -%feature("docstring") MathFunctions::ConvolveFFT "std::vector< complex_t > MathFunctions::ConvolveFFT(const std::vector< double > &signal, const std::vector< double > &resfunc) - -convolution of two real vectors of equal size -"; - %feature("docstring") MathFunctions::GenerateUniformRandom "double MathFunctions::GenerateUniformRandom() "; diff --git a/auto/Wrap/libBornAgainBase.py b/auto/Wrap/libBornAgainBase.py index 484e06d51ec520c37150f0d15c8e1bd023099aa8..499a93da852e560ebfdb19fd0a5c6a912ea8af2b 100644 --- a/auto/Wrap/libBornAgainBase.py +++ b/auto/Wrap/libBornAgainBase.py @@ -1872,31 +1872,6 @@ def Bessel_J1c(*args): """ return _libBornAgainBase.Bessel_J1c(*args) -FORWARD_FFT = _libBornAgainBase.FORWARD_FFT - -BACKWARD_FFT = _libBornAgainBase.BACKWARD_FFT - - -def FastFourierTransform(*args): - r""" - FastFourierTransform(vector_complex_t data, MathFunctions::EFFTDirection tcase) -> vector_complex_t - FastFourierTransform(vdouble1d_t data, MathFunctions::EFFTDirection tcase) -> vector_complex_t - std::vector< complex_t > MathFunctions::FastFourierTransform(const std::vector< double > &data, EFFTDirection tcase) - - simple (and unoptimized) wrapper function for the discrete fast Fourier transformation library (fftw3); transforms real to complex - - """ - return _libBornAgainBase.FastFourierTransform(*args) - -def ConvolveFFT(signal, resfunc): - r""" - ConvolveFFT(vdouble1d_t signal, vdouble1d_t resfunc) -> vector_complex_t - std::vector< complex_t > MathFunctions::ConvolveFFT(const std::vector< double > &signal, const std::vector< double > &resfunc) - - convolution of two real vectors of equal size - - """ - return _libBornAgainBase.ConvolveFFT(signal, resfunc) def GenerateUniformRandom(): r""" diff --git a/auto/Wrap/libBornAgainBase_wrap.cpp b/auto/Wrap/libBornAgainBase_wrap.cpp index 052a829cfdf054e6b69a328e4e59257384b05350..b7df226f771af8e7a5f1666d7fb17a83b1789691 100644 --- a/auto/Wrap/libBornAgainBase_wrap.cpp +++ b/auto/Wrap/libBornAgainBase_wrap.cpp @@ -24876,168 +24876,6 @@ fail: } -SWIGINTERN PyObject *_wrap_FastFourierTransform__SWIG_0(PyObject *SWIGUNUSEDPARM(self), Py_ssize_t nobjs, PyObject **swig_obj) { - PyObject *resultobj = 0; - std::vector< complex_t,std::allocator< complex_t > > *arg1 = 0 ; - MathFunctions::EFFTDirection arg2 ; - int res1 = SWIG_OLDOBJ ; - int val2 ; - int ecode2 = 0 ; - std::vector< complex_t,std::allocator< complex_t > > result; - - if ((nobjs < 2) || (nobjs > 2)) SWIG_fail; - { - std::vector< std::complex< double >,std::allocator< std::complex< double > > > *ptr = (std::vector< std::complex< double >,std::allocator< std::complex< double > > > *)0; - res1 = swig::asptr(swig_obj[0], &ptr); - if (!SWIG_IsOK(res1)) { - SWIG_exception_fail(SWIG_ArgError(res1), "in method '" "FastFourierTransform" "', argument " "1"" of type '" "std::vector< complex_t,std::allocator< complex_t > > const &""'"); - } - if (!ptr) { - SWIG_exception_fail(SWIG_ValueError, "invalid null reference " "in method '" "FastFourierTransform" "', argument " "1"" of type '" "std::vector< complex_t,std::allocator< complex_t > > const &""'"); - } - arg1 = ptr; - } - ecode2 = SWIG_AsVal_int(swig_obj[1], &val2); - if (!SWIG_IsOK(ecode2)) { - SWIG_exception_fail(SWIG_ArgError(ecode2), "in method '" "FastFourierTransform" "', argument " "2"" of type '" "MathFunctions::EFFTDirection""'"); - } - arg2 = static_cast< MathFunctions::EFFTDirection >(val2); - result = MathFunctions::FastFourierTransform((std::vector< std::complex< double >,std::allocator< std::complex< double > > > const &)*arg1,arg2); - resultobj = swig::from(static_cast< std::vector< std::complex< double >,std::allocator< std::complex< double > > > >(result)); - if (SWIG_IsNewObj(res1)) delete arg1; - return resultobj; -fail: - if (SWIG_IsNewObj(res1)) delete arg1; - return NULL; -} - - -SWIGINTERN PyObject *_wrap_FastFourierTransform__SWIG_1(PyObject *SWIGUNUSEDPARM(self), Py_ssize_t nobjs, PyObject **swig_obj) { - PyObject *resultobj = 0; - std::vector< double,std::allocator< double > > *arg1 = 0 ; - MathFunctions::EFFTDirection arg2 ; - int res1 = SWIG_OLDOBJ ; - int val2 ; - int ecode2 = 0 ; - std::vector< complex_t,std::allocator< complex_t > > result; - - if ((nobjs < 2) || (nobjs > 2)) SWIG_fail; - { - std::vector< double,std::allocator< double > > *ptr = (std::vector< double,std::allocator< double > > *)0; - res1 = swig::asptr(swig_obj[0], &ptr); - if (!SWIG_IsOK(res1)) { - SWIG_exception_fail(SWIG_ArgError(res1), "in method '" "FastFourierTransform" "', argument " "1"" of type '" "std::vector< double,std::allocator< double > > const &""'"); - } - if (!ptr) { - SWIG_exception_fail(SWIG_ValueError, "invalid null reference " "in method '" "FastFourierTransform" "', argument " "1"" of type '" "std::vector< double,std::allocator< double > > const &""'"); - } - arg1 = ptr; - } - ecode2 = SWIG_AsVal_int(swig_obj[1], &val2); - if (!SWIG_IsOK(ecode2)) { - SWIG_exception_fail(SWIG_ArgError(ecode2), "in method '" "FastFourierTransform" "', argument " "2"" of type '" "MathFunctions::EFFTDirection""'"); - } - arg2 = static_cast< MathFunctions::EFFTDirection >(val2); - result = MathFunctions::FastFourierTransform((std::vector< double,std::allocator< double > > const &)*arg1,arg2); - resultobj = swig::from(static_cast< std::vector< std::complex< double >,std::allocator< std::complex< double > > > >(result)); - if (SWIG_IsNewObj(res1)) delete arg1; - return resultobj; -fail: - if (SWIG_IsNewObj(res1)) delete arg1; - return NULL; -} - - -SWIGINTERN PyObject *_wrap_FastFourierTransform(PyObject *self, PyObject *args) { - Py_ssize_t argc; - PyObject *argv[3] = { - 0 - }; - - if (!(argc = SWIG_Python_UnpackTuple(args, "FastFourierTransform", 0, 2, argv))) SWIG_fail; - --argc; - if (argc == 2) { - int _v; - int res = swig::asptr(argv[0], (std::vector< std::complex< double >,std::allocator< std::complex< double > > >**)(0)); - _v = SWIG_CheckState(res); - if (_v) { - { - int res = SWIG_AsVal_int(argv[1], NULL); - _v = SWIG_CheckState(res); - } - if (_v) { - return _wrap_FastFourierTransform__SWIG_0(self, argc, argv); - } - } - } - if (argc == 2) { - int _v; - int res = swig::asptr(argv[0], (std::vector< double,std::allocator< double > >**)(0)); - _v = SWIG_CheckState(res); - if (_v) { - { - int res = SWIG_AsVal_int(argv[1], NULL); - _v = SWIG_CheckState(res); - } - if (_v) { - return _wrap_FastFourierTransform__SWIG_1(self, argc, argv); - } - } - } - -fail: - SWIG_Python_RaiseOrModifyTypeError("Wrong number or type of arguments for overloaded function 'FastFourierTransform'.\n" - " Possible C/C++ prototypes are:\n" - " MathFunctions::FastFourierTransform(std::vector< complex_t,std::allocator< complex_t > > const &,MathFunctions::EFFTDirection)\n" - " MathFunctions::FastFourierTransform(std::vector< double,std::allocator< double > > const &,MathFunctions::EFFTDirection)\n"); - return 0; -} - - -SWIGINTERN PyObject *_wrap_ConvolveFFT(PyObject *SWIGUNUSEDPARM(self), PyObject *args) { - PyObject *resultobj = 0; - std::vector< double,std::allocator< double > > *arg1 = 0 ; - std::vector< double,std::allocator< double > > *arg2 = 0 ; - int res1 = SWIG_OLDOBJ ; - int res2 = SWIG_OLDOBJ ; - PyObject *swig_obj[2] ; - std::vector< complex_t,std::allocator< complex_t > > result; - - if (!SWIG_Python_UnpackTuple(args, "ConvolveFFT", 2, 2, swig_obj)) SWIG_fail; - { - std::vector< double,std::allocator< double > > *ptr = (std::vector< double,std::allocator< double > > *)0; - res1 = swig::asptr(swig_obj[0], &ptr); - if (!SWIG_IsOK(res1)) { - SWIG_exception_fail(SWIG_ArgError(res1), "in method '" "ConvolveFFT" "', argument " "1"" of type '" "std::vector< double,std::allocator< double > > const &""'"); - } - if (!ptr) { - SWIG_exception_fail(SWIG_ValueError, "invalid null reference " "in method '" "ConvolveFFT" "', argument " "1"" of type '" "std::vector< double,std::allocator< double > > const &""'"); - } - arg1 = ptr; - } - { - std::vector< double,std::allocator< double > > *ptr = (std::vector< double,std::allocator< double > > *)0; - res2 = swig::asptr(swig_obj[1], &ptr); - if (!SWIG_IsOK(res2)) { - SWIG_exception_fail(SWIG_ArgError(res2), "in method '" "ConvolveFFT" "', argument " "2"" of type '" "std::vector< double,std::allocator< double > > const &""'"); - } - if (!ptr) { - SWIG_exception_fail(SWIG_ValueError, "invalid null reference " "in method '" "ConvolveFFT" "', argument " "2"" of type '" "std::vector< double,std::allocator< double > > const &""'"); - } - arg2 = ptr; - } - result = MathFunctions::ConvolveFFT((std::vector< double,std::allocator< double > > const &)*arg1,(std::vector< double,std::allocator< double > > const &)*arg2); - resultobj = swig::from(static_cast< std::vector< std::complex< double >,std::allocator< std::complex< double > > > >(result)); - if (SWIG_IsNewObj(res1)) delete arg1; - if (SWIG_IsNewObj(res2)) delete arg2; - return resultobj; -fail: - if (SWIG_IsNewObj(res1)) delete arg1; - if (SWIG_IsNewObj(res2)) delete arg2; - return NULL; -} - - SWIGINTERN PyObject *_wrap_GenerateUniformRandom(PyObject *SWIGUNUSEDPARM(self), PyObject *args) { PyObject *resultobj = 0; double result; @@ -33838,21 +33676,6 @@ static PyMethodDef SwigMethods[] = { "Complex Bessel function Bessel_J1(x)/x. \n" "\n" ""}, - { "FastFourierTransform", _wrap_FastFourierTransform, METH_VARARGS, "\n" - "FastFourierTransform(vector_complex_t data, MathFunctions::EFFTDirection tcase) -> vector_complex_t\n" - "FastFourierTransform(vdouble1d_t data, MathFunctions::EFFTDirection tcase) -> vector_complex_t\n" - "std::vector< complex_t > MathFunctions::FastFourierTransform(const std::vector< double > &data, EFFTDirection tcase)\n" - "\n" - "simple (and unoptimized) wrapper function for the discrete fast Fourier transformation library (fftw3); transforms real to complex \n" - "\n" - ""}, - { "ConvolveFFT", _wrap_ConvolveFFT, METH_VARARGS, "\n" - "ConvolveFFT(vdouble1d_t signal, vdouble1d_t resfunc) -> vector_complex_t\n" - "std::vector< complex_t > MathFunctions::ConvolveFFT(const std::vector< double > &signal, const std::vector< double > &resfunc)\n" - "\n" - "convolution of two real vectors of equal size \n" - "\n" - ""}, { "GenerateUniformRandom", _wrap_GenerateUniformRandom, METH_NOARGS, "\n" "GenerateUniformRandom() -> double\n" "double MathFunctions::GenerateUniformRandom()\n" @@ -35826,8 +35649,6 @@ SWIG_init(void) { SWIG_addvarlink(globals, "deg", Swig_var_deg_get, Swig_var_deg_set); SWIG_addvarlink(globals, "tesla", Swig_var_tesla_get, Swig_var_tesla_set); SWIG_addvarlink(globals, "gauss", Swig_var_gauss_get, Swig_var_gauss_set); - SWIG_Python_SetConstant(d, "FORWARD_FFT",SWIG_From_int(static_cast< int >(MathFunctions::FORWARD_FFT))); - SWIG_Python_SetConstant(d, "BACKWARD_FFT",SWIG_From_int(static_cast< int >(MathFunctions::BACKWARD_FFT))); #if PY_VERSION_HEX >= 0x03000000 return m; #else