From 77e8c94e5b41287e97e36e34a97478400637d102 Mon Sep 17 00:00:00 2001 From: Jaron Kent-Dobias Date: Thu, 8 May 2025 15:59:23 -0300 Subject: Work to debug log-fourier --- log-fourier.cpp | 28 ++++++---------------------- 1 file changed, 6 insertions(+), 22 deletions(-) (limited to 'log-fourier.cpp') diff --git a/log-fourier.cpp b/log-fourier.cpp index 9d1f2cb..bc2dd87 100644 --- a/log-fourier.cpp +++ b/log-fourier.cpp @@ -64,15 +64,13 @@ std::vector LogarithmicFourierTransform::fourier(const std::vector= (pad - 1) * N) { - a[n] = c[pad*N-n-1] * exp((1 - k) * τ(pad*N-n-1)); } else { a[n] = 0; } } FFTW_EXECUTE(a_to_â); for (unsigned n = 0; n < pad*N; n++) { - â[(pad*N / 2 + n) % (pad*N)] *= std::pow(II * σ, II * s(n) - k) * Γ(k - II * s(n)); + â[(pad*N / 2 + n) % (pad*N)] *= std::exp(II * s(n)) * std::pow(II * σ, II * s(n) - k) * Γ(k - II * s(n)); } FFTW_EXECUTE(â_to_a); for (unsigned n = 0; n < N; n++) { @@ -87,30 +85,16 @@ std::vector LogarithmicFourierTransform::inverse(const std::vector c(N); std::vector σs = {1, -1}; for (Real σ : σs) { - for (unsigned n = 0; n < pad * N; n++) { - if (n < N) { - if (σ < 0) { - a[n] = std::conj(ĉ[n]) * std::exp((1 - k) * ω(n)); - } else { - a[n] = ĉ[n] * std::exp((1 - k) * ω(n)); - } - } else if (n >= (pad - 1) * N) { - if (σ < 0) { - a[n] = ĉ[pad*N-n-1] * std::exp((1 - k) * ω(pad*N-n-1)); - } else { - a[n] = std::conj(ĉ[pad*N-n-1]) * std::exp((1 - k) * ω(pad*N-n-1)); - } - } else { - a[n] = 0; - } + for (unsigned n = 0; n < N; n++) { + a[n] = (ĉ[n].real() + II * σ * ĉ[n].imag()) * std::exp((1 - k) * ω(n)); } FFTW_EXECUTE(a_to_â); - for (unsigned n = 0; n < pad*N; n++) { - â[(pad*N / 2 + n) % (pad*N)] *= std::pow(-II * σ, II * s(n) - k) * Γ(k - II * s(n)); + for (unsigned n = 0; n < N; n++) { + â[(N / 2 + n) % N] *= std::pow(-II * σ, II * s(n) - k) * Γ(k - II * s(n)); } FFTW_EXECUTE(â_to_a); for (unsigned n = 0; n < N; n++) { - c[n] += std::exp(-k * τ(n)) * a[(pad - 1)*N+n].real() / (Real)(pad*N) / (2 * M_PI); + c[n] += std::exp(-k * τ(n)) * a[n].real() / (2 * M_PI * N); } } -- cgit v1.2.3-70-g09d2