1 files changed, 58 insertions, 0 deletions
diff --git a/log-fourier.hpp b/log-fourier.hpp
new file mode 100644
index 0000000..755f7e9
--- /dev/null
+++ b/log-fourier.hpp
@@ -0,0 +1,58 @@
+#pragma once
+
+#include "types.hpp"
+
+#include <vector>
+#include <tuple>
+
+#include <fftw3.h>
+#include <gsl/gsl_sf_gamma.h>
+
+class LogarithmicFourierTransform {
+private:
+  Complex* a;
+  Complex* â;
+  FFTW_PLAN a_to_â;
+  FFTW_PLAN â_to_a;
+  unsigned N;
+  unsigned pad;
+  Real k;
+  Real Δτ;
+  Real τₛ;
+  Real ωₛ;
+  Real sₛ;
+  std::vector<Real> ts;
+  std::vector<Real> νs;
+  std::vector<Complex> Γs;
+  std::vector<Real> exp1kω;
+  std::vector<Real> exp1kτ;
+  std::vector<Real> expkω;
+  std::vector<Real> expkτ;
+public:
+  Real shift;
+  LogarithmicFourierTransform(unsigned N, Real k, Real Δτ, unsigned pad = 4, Real shift = 0.5);
+  ~LogarithmicFourierTransform();
+  Real τ(unsigned n) const;
+  Real ω(unsigned n) const;
+  Real t(unsigned n) const;
+  Real ν(unsigned n) const;
+  Real s(unsigned n) const;
+  std::vector<Complex> fourier(const std::vector<Real>& c, bool symmetric);
+  std::vector<Real> inverse(const std::vector<Complex>& ĉ);
+};
+
+std::string logFourierFile(std::string prefix, unsigned p, unsigned s, Real λ, Real τ₀, Real β, unsigned log2n, Real Δτ, Real shift);
+
+void logFourierSave(const std::vector<Real>& C, const std::vector<Real>& R, const std::vector<Complex>& Ĉ, const std::vector<Complex>& Ȓ, unsigned p, unsigned s, Real λ, Real τ₀, Real β, unsigned log2n, Real Δτ, Real shift);
+
+bool logFourierLoad(std::vector<Real>& C, std::vector<Real>& R, std::vector<Complex>& Ĉ, std::vector<Complex>& Ȓ, unsigned p, unsigned s, Real λ, Real τ₀, Real β, unsigned log2n, Real Δτ, Real shift);
+
+std::tuple<std::vector<Complex>, std::vector<Complex>> ΣD(LogarithmicFourierTransform& fft, const std::vector<Real>& C, const std::vector<Real>& R, Real β, unsigned p, unsigned s, Real λ);
+
+Real estimateZ(LogarithmicFourierTransform& fft, const std::vector<Real>& C, const std::vector<Complex>& Ĉ, const std::vector<Real>& Ȓ, const std::vector<Complex>& Rt, unsigned p, unsigned s, Real λ, Real τ₀, Real β);
+
+Real energy(const LogarithmicFourierTransform& fft, const std::vector<Real>& C, const std::vector<Real>& R, unsigned p, unsigned s, Real λ, Real β);
+
+Real C0(const LogarithmicFourierTransform& fft, const std::vector<Complex>& Ĉ);
+
+void smooth(std::vector<Real>& C, Real ε);