From 7517890d033b4a3f6c92aa58a62a60d53b1f4814 Mon Sep 17 00:00:00 2001 From: Jeff Hill Date: Tue, 10 Feb 2004 21:43:17 +0000 Subject: [PATCH] installed --- src/libCom/math/epicsDSP.h | 42 ++++++++++++++++ src/libCom/math/epicsFFT.cpp | 31 ++++++++++++ src/libCom/math/epicsTemplateFFT.h | 78 ++++++++++++++++++++++++++++++ 3 files changed, 151 insertions(+) create mode 100644 src/libCom/math/epicsDSP.h create mode 100644 src/libCom/math/epicsFFT.cpp create mode 100644 src/libCom/math/epicsTemplateFFT.h diff --git a/src/libCom/math/epicsDSP.h b/src/libCom/math/epicsDSP.h new file mode 100644 index 000000000..f22ba80a6 --- /dev/null +++ b/src/libCom/math/epicsDSP.h @@ -0,0 +1,42 @@ +/*************************************************************************\ +* Copyright (c) 2002 The University of Chicago, as Operator of Argonne +* National Laboratory. +* Copyright (c) 2002 The Regents of the University of California, as +* Operator of Los Alamos National Laboratory. +* Subject to a Software License Agreement found +* in file LICENSE that is included with this distribution. +\*************************************************************************/ + +/* + * Author: Jeff Hill + * Date: 9-29-03 + */ + +#ifndef epicsDSP_h +#define epicsDSP_h + +#ifdef __cplusplus + +#if !defined ( __GNUC__ ) || ( __GNUC__ > 2 || ( __GNUC__ == 2 && __GNUC_MINOR__ >= 95 ) ) + namespace std { + template < class T > class complex; + } +#else + template < class T > class complex; +#endif + +// +// One dimensional FFT computing template function +// +// vec[] Input Sequence replaced by output sequence +// ln Log base two of the number of point in vec +// inverse Inverse fft if true +// +template < class T > +void epicsOneDimFFT ( std::complex < T > vec[], + const unsigned ln, bool inverse ); + +#endif /* __cplusplus */ + +#endif /* epicsDSP_h */ + diff --git a/src/libCom/math/epicsFFT.cpp b/src/libCom/math/epicsFFT.cpp new file mode 100644 index 000000000..ed011cbea --- /dev/null +++ b/src/libCom/math/epicsFFT.cpp @@ -0,0 +1,31 @@ + +/*************************************************************************\ +* Copyright (c) 2002 The University of Chicago, as Operator of Argonne +* National Laboratory. +* Copyright (c) 2002 The Regents of the University of California, as +* Operator of Los Alamos National Laboratory. +* Subject to a Software License Agreement found +* in file LICENSE that is included with this distribution. +\*************************************************************************/ + +// +// Simple fft computing functions +// +// Author: Jeff Hill +// Date: 9-29-03 +// + +#include + +#define epicsExportSharedSymbols +#include +#include + +// template epicsShareFunc void epicsShareAPI epicsOneDimFFT ( +// std::complex < float > vec[], +// const unsigned ln, bool inverse ); + +template epicsShareFunc void epicsShareAPI epicsOneDimFFT ( + std::complex < double > vec[], + const unsigned ln, bool inverse ); + diff --git a/src/libCom/math/epicsTemplateFFT.h b/src/libCom/math/epicsTemplateFFT.h new file mode 100644 index 000000000..f1048d12c --- /dev/null +++ b/src/libCom/math/epicsTemplateFFT.h @@ -0,0 +1,78 @@ +/*************************************************************************\ +* Copyright (c) 2002 The University of Chicago, as Operator of Argonne +* National Laboratory. +* Copyright (c) 2002 The Regents of the University of California, as +* Operator of Los Alamos National Laboratory. +* Subject to a Software License Agreement found +* in file LICENSE that is included with this distribution. +\*************************************************************************/ + +// +// Author: Jeff Hill +// Date: 9-29-03 +// + +#ifndef epicsTemplateFFT_h +#define epicsTemplateFFT_h + +#include + +#include + +// +// Simple one dimensional FFT computing template function +// +// vec[] Input Sequence replaced by output sequence +// ln Log base two of the number of point in vec +// inverse Inverse fft if true +// +template < class T > +void epicsOneDimFFT ( std::complex < T > vec[], + const unsigned ln, bool inverse ) +{ + static const T one = 1.0; + static const T four = 4.0; + static const T PI = four * atan ( one ); + const unsigned n = 1u << ln; + const unsigned nv2 = n >> 1u; + + { + unsigned j = 1u; + + for ( unsigned i = 1; i < n; i++ ) { + if ( i < j ) { + std::complex < T > t = vec[i - 1]; + vec[i - 1] = vec[j - 1]; + vec[j - 1] = t; + } + unsigned k = nv2; + while ( k < j ) { + j = j - k; + k = k >> 1u; + } + j = j + k; + } + } + + const T sign = inverse ? 1.0 : -1.0; + for ( unsigned l = 1; l <= ln; l++ ) { + const unsigned le = 1u << l; + const unsigned le1 = le >> 1u; + const T radians = sign * PI / le1; + const std::complex < T > w = std::polar ( one, radians ); + std::complex < T > u ( 1, 0 ); + + for ( unsigned m = 1u; m <= le1; m++ ) { + for ( unsigned i = m - 1; i < n; i = i + le ) { + unsigned ip = i + le1; + std::complex < T > t0 = vec [ ip ] * u; + vec [ ip ] = vec [ i ] - t0; + vec [ i ] = vec [ i ] + t0; + } + u = u * w; + } + } +} + +#endif // epicsTemplateFFT_h +