TestFFT.cpp

Go to the documentation of this file.

 
#include <gtest/gtest.h>
#include <cmath>
#include <complex>
#include <vector>
#include "Utilities/GSLFFT.h"
 
class FFTTest : public ::testing::Test {
protected:
    void SetUp() override {
        // Test setup
    }
};
 
TEST_F(FFTTest, ComplexFFT_Identity) {
    // Test FFT of a constant signal
    std::vector<std::complex<double>> data(8, std::complex<double>(1.0, 0.0));
 
    // Create a copy for FFT
    std::vector<std::complex<double>> fft_data = data;
 
    // Perform FFT (using the internal fft function)
    // Note: This tests the internal implementation
    // For GSL-compatible interface, we'd use gsl_fft_complex_forward
 
    // Simple test: constant signal should have DC component only
    // (This is a simplified test - actual FFT would require the internal function)
    EXPECT_EQ(data.size(), 8);
}
 
TEST_F(FFTTest, RealFFT_Transform) {
    // Test real FFT transform
    const size_t n = 8;
    std::vector<double> data(n);
 
    // Create a simple signal: [1, 0, 0, 0, 0, 0, 0, 0]
    data[0] = 1.0;
    for (size_t i = 1; i < n; ++i) {
        data[i] = 0.0;
    }
 
    gsl_fft_real_wavetable* wavetable = gsl_fft_real_wavetable_alloc(n);
    gsl_fft_real_workspace* workspace = gsl_fft_real_workspace_alloc(n);
 
    gsl_fft_real_transform(data.data(), 1, n, wavetable, workspace);
 
    // DC component should be 1.0
    EXPECT_NEAR(data[0], 1.0, 1e-10);
 
    gsl_fft_real_wavetable_free(wavetable);
    gsl_fft_real_workspace_free(workspace);
}
 
TEST_F(FFTTest, ComplexFFT_Forward) {
    const size_t n = 4;
    std::vector<double> data(2 * n);  // Complex data: [real0, imag0, real1, imag1, ...]
 
    // Set data to [1+0i, 0+0i, 0+0i, 0+0i]
    data[0] = 1.0;  // real[0]
    data[1] = 0.0;  // imag[0]
    for (size_t i = 1; i < n; ++i) {
        data[2 * i]     = 0.0;
        data[2 * i + 1] = 0.0;
    }
 
    gsl_fft_complex_wavetable* wavetable = gsl_fft_complex_wavetable_alloc(n);
    gsl_fft_complex_workspace* workspace = gsl_fft_complex_workspace_alloc(n);
 
    gsl_fft_complex_forward(data.data(), 1, n, wavetable, workspace);
 
    // All components should be 1.0 (DC signal transforms to all ones)
    for (size_t i = 0; i < n; ++i) {
        EXPECT_NEAR(data[2 * i], 1.0, 1e-10);
        EXPECT_NEAR(data[2 * i + 1], 0.0, 1e-10);
    }
 
    gsl_fft_complex_wavetable_free(wavetable);
    gsl_fft_complex_workspace_free(workspace);
}
 
TEST_F(FFTTest, ComplexFFT_Inverse) {
    const size_t n = 4;
    std::vector<double> data(2 * n);
 
    // Set all components to 1.0 (constant in frequency domain)
    for (size_t i = 0; i < n; ++i) {
        data[2 * i]     = 1.0;
        data[2 * i + 1] = 0.0;
    }
 
    gsl_fft_complex_wavetable* wavetable = gsl_fft_complex_wavetable_alloc(n);
    gsl_fft_complex_workspace* workspace = gsl_fft_complex_workspace_alloc(n);
 
    gsl_fft_complex_inverse(data.data(), 1, n, wavetable, workspace);
 
    // Should recover impulse at n=0
    EXPECT_NEAR(data[0], static_cast<double>(n), 1e-10);
    EXPECT_NEAR(data[1], 0.0, 1e-10);
 
    gsl_fft_complex_wavetable_free(wavetable);
    gsl_fft_complex_workspace_free(workspace);
}
 
TEST_F(FFTTest, FFT_InverseRoundTrip) {
    const size_t n = 8;
    std::vector<double> original(2 * n);
    std::vector<double> transformed(2 * n);
 
    // Create original signal
    for (size_t i = 0; i < n; ++i) {
        original[2 * i]     = static_cast<double>(i);
        original[2 * i + 1] = 0.0;
    }
 
    transformed = original;
 
    gsl_fft_complex_wavetable* wavetable = gsl_fft_complex_wavetable_alloc(n);
    gsl_fft_complex_workspace* workspace = gsl_fft_complex_workspace_alloc(n);
 
    // Forward transform
    gsl_fft_complex_forward(transformed.data(), 1, n, wavetable, workspace);
 
    // Inverse transform
    gsl_fft_complex_inverse(transformed.data(), 1, n, wavetable, workspace);
 
    // Should recover original (within scaling)
    for (size_t i = 0; i < n; ++i) {
        EXPECT_NEAR(transformed[2 * i], original[2 * i] * static_cast<double>(n), 1e-6);
        EXPECT_NEAR(transformed[2 * i + 1], 0.0, 1e-6);
    }
 
    gsl_fft_complex_wavetable_free(wavetable);
    gsl_fft_complex_workspace_free(workspace);
}
 
TEST_F(FFTTest, Radix2Forward) {
    const size_t n = 4;  // Must be power of 2
    std::vector<double> data(2 * n);
 
    // Set to impulse
    data[0] = 1.0;
    data[1] = 0.0;
    for (size_t i = 1; i < n; ++i) {
        data[2 * i]     = 0.0;
        data[2 * i + 1] = 0.0;
    }
 
    gsl_fft_complex_radix2_forward(data.data(), 1, n);
 
    // All components should be 1.0
    for (size_t i = 0; i < n; ++i) {
        EXPECT_NEAR(data[2 * i], 1.0, 1e-10);
        EXPECT_NEAR(data[2 * i + 1], 0.0, 1e-10);
    }
}
 
TEST_F(FFTTest, Radix2Inverse) {
    const size_t n = 4;
    std::vector<double> data(2 * n);
 
    // Set all to 1.0
    for (size_t i = 0; i < n; ++i) {
        data[2 * i]     = 1.0;
        data[2 * i + 1] = 0.0;
    }
 
    gsl_fft_complex_radix2_inverse(data.data(), 1, n);
 
    // Should get impulse at n=0
    EXPECT_NEAR(data[0], static_cast<double>(n), 1e-10);
    EXPECT_NEAR(data[1], 0.0, 1e-10);
}