From 6103b583bf284987d6171d7079e73b7f54f07111 Mon Sep 17 00:00:00 2001 From: Christoph Ruegg Date: Fri, 28 Oct 2016 09:06:07 +0200 Subject: [PATCH] FFT-MKL: use complex complex mode for conjugate even storage, simplify --- src/NativeProviders/MKL/fft.cpp | 23 +++++++--------------- src/Numerics/IntegralTransforms/Fourier.cs | 10 ++++++++++ 2 files changed, 17 insertions(+), 16 deletions(-) diff --git a/src/NativeProviders/MKL/fft.cpp b/src/NativeProviders/MKL/fft.cpp index 072c4ca3..fccf5ce2 100644 --- a/src/NativeProviders/MKL/fft.cpp +++ b/src/NativeProviders/MKL/fft.cpp @@ -6,11 +6,6 @@ #include #include "mkl_dfti.h" -inline MKL_LONG fft_free(DFTI_DESCRIPTOR_HANDLE* handle) -{ - return DftiFreeDescriptor(handle); -} - template inline MKL_LONG fft_create_1d(DFTI_DESCRIPTOR_HANDLE* handle, const MKL_LONG n, const Precision forward_scale, const Precision backward_scale, const DFTI_CONFIG_VALUE precision, const DFTI_CONFIG_VALUE domain) { @@ -18,6 +13,7 @@ inline MKL_LONG fft_create_1d(DFTI_DESCRIPTOR_HANDLE* handle, const MKL_LONG n, DFTI_DESCRIPTOR_HANDLE descriptor = *handle; if (0 == status) status = DftiSetValue(descriptor, DFTI_FORWARD_SCALE, forward_scale); if (0 == status) status = DftiSetValue(descriptor, DFTI_BACKWARD_SCALE, backward_scale); + if (0 == status) status = DftiSetValue(descriptor, DFTI_CONJUGATE_EVEN_STORAGE, DFTI_COMPLEX_COMPLEX); if (0 == status) status = DftiCommitDescriptor(descriptor); return status; } @@ -29,21 +25,16 @@ inline MKL_LONG fft_create_md(DFTI_DESCRIPTOR_HANDLE* handle, MKL_LONG dimension DFTI_DESCRIPTOR_HANDLE descriptor = *handle; if (0 == status) status = DftiSetValue(descriptor, DFTI_FORWARD_SCALE, forward_scale); if (0 == status) status = DftiSetValue(descriptor, DFTI_BACKWARD_SCALE, backward_scale); + if (0 == status) status = DftiSetValue(descriptor, DFTI_CONJUGATE_EVEN_STORAGE, DFTI_COMPLEX_COMPLEX); if (0 == status) status = DftiCommitDescriptor(descriptor); return status; } -template -inline MKL_LONG fft_compute(const DFTI_DESCRIPTOR_HANDLE handle, Data x[], FFT fft) -{ - return fft(handle, x); -} - extern "C" { DLLEXPORT MKL_LONG x_fft_free(DFTI_DESCRIPTOR_HANDLE* handle) { - return fft_free(handle); + return DftiFreeDescriptor(handle); } DLLEXPORT MKL_LONG z_fft_create(DFTI_DESCRIPTOR_HANDLE* handle, const MKL_LONG n, const double forward_scale, const double backward_scale) @@ -68,21 +59,21 @@ extern "C" { DLLEXPORT MKL_LONG z_fft_forward(const DFTI_DESCRIPTOR_HANDLE handle, MKL_Complex16 x[]) { - return fft_compute(handle, x, DftiComputeForward); + return DftiComputeForward(handle, x); } DLLEXPORT MKL_LONG c_fft_forward(const DFTI_DESCRIPTOR_HANDLE handle, MKL_Complex8 x[]) { - return fft_compute(handle, x, DftiComputeForward); + return DftiComputeForward(handle, x); } DLLEXPORT MKL_LONG z_fft_backward(const DFTI_DESCRIPTOR_HANDLE handle, MKL_Complex16 x[]) { - return fft_compute(handle, x, DftiComputeBackward); + return DftiComputeBackward(handle, x); } DLLEXPORT MKL_LONG c_fft_backward(const DFTI_DESCRIPTOR_HANDLE handle, MKL_Complex8 x[]) { - return fft_compute(handle, x, DftiComputeBackward); + return DftiComputeBackward(handle, x); } } diff --git a/src/Numerics/IntegralTransforms/Fourier.cs b/src/Numerics/IntegralTransforms/Fourier.cs index b853ca30..5e7f383f 100644 --- a/src/Numerics/IntegralTransforms/Fourier.cs +++ b/src/Numerics/IntegralTransforms/Fourier.cs @@ -110,6 +110,16 @@ namespace MathNet.Numerics.IntegralTransforms } } + public static void ForwardReal(double[] samples, FourierOptions options = FourierOptions.Default) + { + if (real.Length != imaginary.Length) + { + throw new ArgumentException(Resources.ArgumentArraysSameLength); + } + + int length = real.Length/2; + } + /// /// Applies the forward Fast Fourier Transform (FFT) to multiple dimensional sample data. ///