mirror/Sprout
1
0
Fork 0
mirror of https://github.com/bolero-MURAKAMI/Sprout.git synced 2024-12-23 21:25:49 +00:00
Code Issues Projects Releases Packages Wiki Activity Actions

fix fft implementation

Browse source
This commit is contained in:
bolero-MURAKAMI 2014-04-18 19:33:48 +09:00
parent 07b5f69ebb
commit d1fc657c1f
5 changed files with 72 additions and 17 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

1
sprout/complex.hpp
View file

@ -16,5 +16,6 @@
#include <sprout/complex/hash.hpp>
#include <sprout/complex/nearest.hpp>
#include <sprout/complex/udl.hpp>
#include <sprout/complex/type_traits.hpp>
#endif // #ifndef SPROUT_COMPLEX_HPP

6
sprout/complex/complex.hpp
View file

@ -15,9 +15,9 @@ namespace sprout {
class complex;
namespace detail {
template<typename T>
inline SPROUT_CONSTEXPR T
complex_norm(sprout::complex<T> const& x) {
template<typename Complex>
inline SPROUT_CONSTEXPR typename Complex::value_type
complex_norm(Complex const& x) {
return x.real() * x.real() + x.imag() * x.imag();
}
} // namespace detail

42
sprout/complex/type_traits.hpp Normal file
View file

@ -0,0 +1,42 @@
/*=============================================================================
Copyright (c) 2011-2014 Bolero MURAKAMI
https://github.com/bolero-MURAKAMI/Sprout
Distributed under the Boost Software License, Version 1.0. (See accompanying
file LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
=============================================================================*/
#ifndef SPROUT_COMPLEX_TYPE_TRAITS_HPP
#define SPROUT_COMPLEX_TYPE_TRAITS_HPP
#include <sprout/config.hpp>
#include <sprout/complex/complex.hpp>
#include <sprout/type_traits/integral_constant.hpp>
namespace sprout {
//
// is_complex
//
template<typename T>
struct is_complex
: public sprout::false_type
{};
template<typename T>
struct is_complex<T const>
: public sprout::is_complex<T>
{};
template<typename T>
struct is_complex<T const volatile>
: public sprout::is_complex<T>
{};
template<typename T>
struct is_complex<sprout::complex<T> >
: public sprout::true_type
{};
#if SPROUT_USE_VARIABLE_TEMPLATES
template<typename T>
SPROUT_STATIC_CONSTEXPR bool is_complex_v = sprout::is_complex<T>::value;
#endif // #if SPROUT_USE_VARIABLE_TEMPLATES
} // namespace sprout
#endif // #ifndef SPROUT_COMPLEX_TYPE_TRAITS_HPP

4
sprout/complex/values.hpp
View file

@ -33,12 +33,12 @@ namespace sprout {
return x.imag();
}
template<typename T>
inline SPROUT_CXX14_CONSTEXPR T&
inline SPROUT_CONSTEXPR T&
real(sprout::complex<T>& x) {
return x.real();
}
template<typename T>
inline SPROUT_CXX14_CONSTEXPR T&
inline SPROUT_CONSTEXPR T&
imag(sprout::complex<T>& x) {
return x.imag();
}

36
sprout/numeric/fft/fft.hpp
View file

@ -51,10 +51,14 @@ namespace sprout {
first[j] += first[j1];
value_type x3 = first[j3] - first[j2];
first[j2] += first[j3];
real(first[j1]) = real(x1) - imag(x3);
imag(first[j1]) = imag(x1) + real(x3);
real(first[j3]) = real(x1) + imag(x3);
imag(first[j3]) = imag(x1) - real(x3);
first[j1] = value_type(
real(x1) - imag(x3),
imag(x1) + real(x3)
);
first[j3] = value_type(
real(x1) + imag(x3),
imag(x1) - real(x3)
);
}
}
if (m == size) {
@ -73,12 +77,16 @@ namespace sprout {
first[j2] += first[j3];
elem_type x0r = real(x1) - imag(x3);
elem_type x0i = imag(x1) + real(x3);
real(first[j1]) = w1r * (x0r + x0i);
imag(first[j1]) = w1r * (x0i - x0r);
first[j1] = value_type(
w1r * (x0r + x0i),
w1r * (x0i - x0r)
);
x0r = real(x1) + imag(x3);
x0i = imag(x1) - real(x3);
real(first[j3]) = w1r * (-x0r + x0i);
imag(first[j3]) = w1r * (-x0i - x0r);
first[j3] = value_type(
w1r * (-x0r + x0i),
w1r * (-x0i - x0r)
);
}
}
for (difference_type i = 2 * m; i < size; i += m) {
@ -99,12 +107,16 @@ namespace sprout {
first[j2] += first[j3];
elem_type x0r = real(x1) - imag(x3);
elem_type x0i = imag(x1) + real(x3);
real(first[j1]) = w1r * x0r - w1i * x0i;
imag(first[j1]) = w1r * x0i + w1i * x0r;
first[j1] = value_type(
w1r * x0r - w1i * x0i,
w1r * x0i + w1i * x0r
);
x0r = real(x1) + imag(x3);
x0i = imag(x1) - real(x3);
real(first[j3]) = w3r * x0r - w3i * x0i;
imag(first[j3]) = w3r * x0i + w3i * x0r;
first[j3] = value_type(
w3r * x0r - w3i * x0i,
w3r * x0i + w3i * x0r
);
}
}
}