diff --git a/.claude/agents/modernizer.md b/.claude/agents/modernizer.md
new file mode 100644
index 00000000..27e6dead
--- /dev/null
+++ b/.claude/agents/modernizer.md
@@ -0,0 +1,64 @@
+---
+name: modernizer
+description: Refactor a pre-roadmap numerical/ algorithm to reuse shared math/ utilities, add coverage-build infra, and simplify/dedupe tests — preserving existing Q15/Q31 support. Behavior-preserving.
+model: claude-sonnet-4-6
+tools: [Read, Write, Edit, Bash, TodoWrite]
+---
+
+Canonical rules: `AGENTS.md`. You modernize ONE pre-roadmap algorithm at a time so it matches
+current roadmap conventions. Refactor is **behavior-preserving** — no new features, no public-API
+change unless required to remove duplication.
+
+## Workflow
+
+1. Read the target `.hpp`, its `test/Test*.cpp`, and `doc/<domain>/<Name>.md`.
+2. Find duplication: logic that already exists in `numerical/math/`
+   (`CompilerOptimizations.hpp`, `Tolerance.hpp`, `ComplexNumber.hpp`, `QNumber.hpp`,
+   `RecursiveBuffer.hpp`, `Statistics.hpp`, …) or scaffolding repeated across `test/` files
+   (e.g. `CalculateMagnitude`, twiddle-factor mocks — emulate `PowerDensitySpectrumTestSupport.hpp`).
+3. Replace the duplicated logic with the shared utility; delete the local copy.
+4. Add coverage-build infra if missing (per `roadmap/DEPLOYMENT.md`): guarded `extern template`
+   at header bottom + matching `.cpp` (`template class <Name><float, ...>;`) wired via
+   `numerical_add_coverage_sources()`.
+5. Refactor the `test/Test*.cpp` — **mandatory, never skip**. It is a required deliverable, not
+   optional cleanup: apply the **Tests checklist** below on every run, even when the production
+   code needs no change and even when the test already builds green.
+6. Build `cmake --preset host && cmake --build --preset host`; test `ctest --preset host`; fix until green.
+7. Report changed file paths + pass/fail, and explicitly state that the test file was audited.
+
+## Tests checklist — apply every run
+
+The `test/Test*.cpp` is a first-class deliverable of every modernization. Audit and fix:
+- [ ] `StrictMock<...>` only — never a bare mock or `NiceMock`.
+- [ ] Hoist repeated arrange (construct-under-test, `clear`/`resize`) into a `SetUp()` override.
+- [ ] Extract repeated mock-return / `WillOnce(Invoke(...))` tails into fixture helper methods.
+- [ ] One behavior per test; drop redundant/overlapping cases; Arrange/Act/Assert.
+- [ ] `EXPECT_NEAR` + `math::Tolerance<float>()` for float comparisons.
+- [ ] Anonymous-namespace fixture; macros outside; no heap; no comments.
+
+Keep `TYPED_TEST` where the algorithm is multi-type; behavior and assertions stay identical.
+
+## Preserve types — hard rule
+
+Keep existing `Q15`/`Q31` support and its `TYPED_TEST` where the algorithm already has it — do
+**NOT** strip multi-type. The multi-type guard
+(`static_assert(math::is_qnumber<T>::value || std::is_floating_point_v<T>, ...)`) stays for those.
+Float-only migration applies ONLY to algorithms that are already float-only; those follow
+`static_assert(std::is_floating_point_v<T>)` + `TEST_F` on `float`.
+
+## Memory — quick reference
+
+**Forbidden**: `new`/`delete`/`malloc`/`free`, `make_unique`/`make_shared`,
+`std::vector`/`string`/`deque`/`list`/`map`/`set`. Tests too. No recursion.
+
+**Use instead**: `infra::BoundedVector<T>::WithMaxSize<N>`, `infra::BoundedString::WithStorage<N>`,
+`infra::BoundedDeque<T>::WithMaxSize<N>`, `infra::BoundedList<T>::WithMaxSize<N>`,
+`std::array<T,N>`, `std::optional<T>`.
+
+## What NOT to do
+- No behavior/API change beyond removing duplication.
+- No new abstractions except extracting one that is already repeated.
+- Don't touch unrelated algorithms; don't strip `Q15`/`Q31`.
+- No `make_unique` anywhere, including tests.
+
+**Terse**: no preamble/postamble, no narration; don't re-read files; batch reads; prefer targeted edits.
diff --git a/.claude/agents/orchestrator.md b/.claude/agents/orchestrator.md
index 6e081ed2..b313560e 100644
--- a/.claude/agents/orchestrator.md
+++ b/.claude/agents/orchestrator.md
@@ -16,6 +16,7 @@ Triage requests and route to the right specialist via the Agent tool. Do NOT imp
    - **planner** — new algorithm, architectural change, multi-file work
    - **executor** — clear bug fix, small change, existing plan
    - **reviewer** — review existing or recent code
+   - **modernizer** — refactor/dedupe a pre-roadmap algorithm, reuse shared utilities, simplify tests
 
 ## Context to gather
 - Module: `analysis`, `windowing`, `control_analysis`, `controllers`, `dynamics`,
diff --git a/.github/agents/modernizer.agent.md b/.github/agents/modernizer.agent.md
new file mode 100644
index 00000000..2b828869
--- /dev/null
+++ b/.github/agents/modernizer.agent.md
@@ -0,0 +1,67 @@
+---
+description: "Refactor a pre-roadmap numerical/ algorithm to reuse shared math/ utilities, add coverage-build infra, and simplify/dedupe tests — preserving existing Q15/Q31 support. Behavior-preserving."
+tools: [read, edit, search, execute, todo]
+model: "Claude Sonnet 4.6"
+handoffs:
+  - label: "Review Changes"
+    agent: reviewer
+    prompt: "Review the refactoring changes made above against numerical-toolbox project standards."
+---
+
+Canonical rules: `AGENTS.md`. You modernize ONE pre-roadmap algorithm at a time so it matches
+current roadmap conventions. Refactor is **behavior-preserving** — no new features, no public-API
+change unless required to remove duplication.
+
+## Workflow
+
+1. Read the target `.hpp`, its `test/Test*.cpp`, and `doc/<domain>/<Name>.md`.
+2. Find duplication: logic that already exists in `numerical/math/`
+   (`CompilerOptimizations.hpp`, `Tolerance.hpp`, `ComplexNumber.hpp`, `QNumber.hpp`,
+   `RecursiveBuffer.hpp`, `Statistics.hpp`, …) or scaffolding repeated across `test/` files
+   (e.g. `CalculateMagnitude`, twiddle-factor mocks — emulate `PowerDensitySpectrumTestSupport.hpp`).
+3. Replace the duplicated logic with the shared utility; delete the local copy.
+4. Add coverage-build infra if missing (per `roadmap/DEPLOYMENT.md`): guarded `extern template`
+   at header bottom + matching `.cpp` (`template class <Name><float, ...>;`) wired via
+   `numerical_add_coverage_sources()`.
+5. Refactor the `test/Test*.cpp` — **mandatory, never skip**. It is a required deliverable, not
+   optional cleanup: apply the **Tests checklist** below on every run, even when the production
+   code needs no change and even when the test already builds green.
+6. Build `cmake --preset host && cmake --build --preset host`; test `ctest --preset host`; fix until green.
+7. Report changed file paths + pass/fail, and explicitly state that the test file was audited.
+
+## Tests checklist — apply every run
+
+The `test/Test*.cpp` is a first-class deliverable of every modernization. Audit and fix:
+- [ ] `StrictMock<...>` only — never a bare mock or `NiceMock`.
+- [ ] Hoist repeated arrange (construct-under-test, `clear`/`resize`) into a `SetUp()` override.
+- [ ] Extract repeated mock-return / `WillOnce(Invoke(...))` tails into fixture helper methods.
+- [ ] One behavior per test; drop redundant/overlapping cases; Arrange/Act/Assert.
+- [ ] `EXPECT_NEAR` + `math::Tolerance<float>()` for float comparisons.
+- [ ] Anonymous-namespace fixture; macros outside; no heap; no comments.
+
+Keep `TYPED_TEST` where the algorithm is multi-type; behavior and assertions stay identical.
+
+## Preserve types — hard rule
+
+Keep existing `Q15`/`Q31` support and its `TYPED_TEST` where the algorithm already has it — do
+**NOT** strip multi-type. The multi-type guard
+(`static_assert(math::is_qnumber<T>::value || std::is_floating_point_v<T>, ...)`) stays for those.
+Float-only migration applies ONLY to algorithms that are already float-only; those follow
+`static_assert(std::is_floating_point_v<T>)` + `TEST_F` on `float`.
+
+## Memory — quick reference
+
+**Forbidden**: `new`/`delete`/`malloc`/`free`, `make_unique`/`make_shared`,
+`std::vector`/`string`/`deque`/`list`/`map`/`set`. Tests too. No recursion.
+
+**Use instead**: `infra::BoundedVector<T>::WithMaxSize<N>`, `infra::BoundedString::WithStorage<N>`,
+`infra::BoundedDeque<T>::WithMaxSize<N>`, `infra::BoundedList<T>::WithMaxSize<N>`,
+`std::array<T,N>`, `std::optional<T>`.
+
+## What NOT to do
+- No behavior/API change beyond removing duplication.
+- No new abstractions except extracting one that is already repeated.
+- Don't touch unrelated algorithms; don't strip `Q15`/`Q31`.
+- No `make_unique` anywhere, including tests.
+
+**Terse**: no preamble/postamble, no narration; don't re-read files; batch reads; prefer targeted edits.
diff --git a/.github/agents/orchestrator.agent.md b/.github/agents/orchestrator.agent.md
index 63faae21..12faa6ad 100644
--- a/.github/agents/orchestrator.agent.md
+++ b/.github/agents/orchestrator.agent.md
@@ -2,7 +2,7 @@
 description: "Triage development tasks in numerical-toolbox and route to planner, executor, or reviewer. Start here for any new feature, bug fix, or code review."
 tools: [read, search, web, agent]
 model: "Claude Sonnet 4.6"
-agents: [planner, executor, reviewer]
+agents: [planner, executor, reviewer, modernizer]
 handoffs:
   - label: "Plan Implementation"
     agent: planner
@@ -13,6 +13,9 @@ handoffs:
   - label: "Review Code"
     agent: reviewer
     prompt: "Review the code changes described above against numerical-toolbox project standards."
+  - label: "Modernize Legacy Algorithm"
+    agent: modernizer
+    prompt: "Refactor the pre-roadmap algorithm described above to reuse shared utilities and simplify tests, following numerical-toolbox conventions."
 ---
 
 Triage requests and route to the right specialist. Do NOT implement or plan yourself.
@@ -26,6 +29,7 @@ Triage requests and route to the right specialist. Do NOT implement or plan your
    - **planner** — new algorithm, architectural change, multi-file work
    - **executor** — clear bug fix, small change, existing plan
    - **reviewer** — review existing or recent code
+   - **modernizer** — refactor/dedupe a pre-roadmap algorithm, reuse shared utilities, simplify tests
 
 ## Context to gather
 - Module: `analysis`, `windowing`, `control_analysis`, `controllers`, `dynamics`,
diff --git a/TESTING.md b/TESTING.md
index a26e2846..e39d152b 100644
--- a/TESTING.md
+++ b/TESTING.md
@@ -72,7 +72,7 @@ Canonical rules still apply ([AGENTS.md](AGENTS.md), [testing.instructions.md](.
 ## Per-family detail
 
 ### 1. Signal transforms — `analysis/`
-`FastFourierTransformRadix2Impl`, `RealFft`, `DiscreteCosineTransform`, `GoertzelAlgorithm`,
+`FastFourierTransformRadix2Impl`, `RealFastFourierTransform`, `DiscreteCosineTransform`, `GoertzelAlgorithm`,
 `ConvolutionCorrelation`, `PowerDensitySpectrum`, `SignalDetectors`, `windowing/`.
 
 - **M1 accuracy** — known transform pairs: δ[n] → flat spectrum; single sinusoid → single bin at its
@@ -80,7 +80,7 @@ Canonical rules still apply ([AGENTS.md](AGENTS.md), [testing.instructions.md](.
 - **M7 Parseval / energy** — `Σ|x|² ≈ (1/N)·Σ|X|²`; assert residual near 0.
 - **M1 linearity** — `F(a·x + b·y) = a·F(x) + b·F(y)`.
 - **M1 round-trip** — `Inverse(Forward(x)) ≈ x`; assert reconstruction RMSE.
-- **M7 symmetry** — real input ⇒ conjugate-symmetric spectrum (RealFft): `X[N-k] = conj(X[k])`.
+- **M7 symmetry** — real input ⇒ conjugate-symmetric spectrum (RealFastFourierTransform): `X[N-k] = conj(X[k])`.
 - **Convolution** — matches the direct sum; `x * δ = x`; commutativity; output length.
 - **PSD** — non-negative; total power = signal variance; spectral peak at the tone's frequency.
 - **Goertzel** — single-bin magnitude equals the full-FFT bin.
diff --git a/doc/analysis/README.md b/doc/analysis/README.md
index 9e7606c7..b324d768 100644
--- a/doc/analysis/README.md
+++ b/doc/analysis/README.md
@@ -7,7 +7,7 @@ Signal analysis algorithms for frequency-domain decomposition and spectral estim
 | Algorithm                                               | Description                                                                                      |
 |---------------------------------------------------------|--------------------------------------------------------------------------------------------------|
 | [Fast Fourier Transform](FastFourierTransform.md)       | Efficient computation of the Discrete Fourier Transform using the Cooley-Tukey radix-2 algorithm |
-| [Real-Input FFT](RealFft.md)                            | Length-N real FFT via even/odd split into two N/2-point complex DFTs, halving compute and memory |
+| [Real-Input FFT](RealFastFourierTransform.md)                            | Length-N real FFT via even/odd split into two N/2-point complex DFTs, halving compute and memory |
 | [Power Spectral Density](PowerDensitySpectrum.md)       | Estimation of signal power distribution across frequencies using Welch's method                  |
 | [Discrete Cosine Transform](DiscreteCosineTransform.md) | Real-valued frequency decomposition via cosine basis functions, computed through FFT             |
 | [Signal Detectors](SignalDetectors.md)                  | Peak hold, zero-crossing counter, and RMS envelope detectors for real-time signal monitoring     |
diff --git a/doc/analysis/RealFft.md b/doc/analysis/RealFastFourierTransform.md
similarity index 100%
rename from doc/analysis/RealFft.md
rename to doc/analysis/RealFastFourierTransform.md
diff --git a/numerical/analysis/CMakeLists.txt b/numerical/analysis/CMakeLists.txt
index 540ff3fb..39645b5b 100644
--- a/numerical/analysis/CMakeLists.txt
+++ b/numerical/analysis/CMakeLists.txt
@@ -19,16 +19,17 @@ target_sources(numerical.analysis PRIVATE
     FastFourierTransformRadix2Impl.hpp
     GoertzelAlgorithm.hpp
     PowerDensitySpectrum.hpp
-    RealFft.hpp
+    RealFastFourierTransform.hpp
     SignalDetectors.hpp
 )
 
 numerical_add_coverage_sources(numerical.analysis
     ConvolutionCorrelation.cpp
+    DiscreteCosineTransform.cpp
     FastFourierTransformRadix2Impl.cpp
     GoertzelAlgorithm.cpp
     PowerDensitySpectrum.cpp
-    RealFft.cpp
+    RealFastFourierTransform.cpp
     SignalDetectors.cpp
 )
 
diff --git a/numerical/analysis/DiscreteCosineTransform.cpp b/numerical/analysis/DiscreteCosineTransform.cpp
new file mode 100644
index 00000000..a0e2998f
--- /dev/null
+++ b/numerical/analysis/DiscreteCosineTransform.cpp
@@ -0,0 +1,9 @@
+#include "numerical/analysis/DiscreteCosineTransform.hpp"
+#include "numerical/math/QNumber.hpp"
+
+namespace analysis
+{
+    template class DiscreteConsineTransform<float, 8>;
+    template class DiscreteConsineTransform<math::Q15, 8>;
+    template class DiscreteConsineTransform<math::Q31, 8>;
+}
diff --git a/numerical/analysis/DiscreteCosineTransform.hpp b/numerical/analysis/DiscreteCosineTransform.hpp
index ee8d96f1..72a5d022 100644
--- a/numerical/analysis/DiscreteCosineTransform.hpp
+++ b/numerical/analysis/DiscreteCosineTransform.hpp
@@ -6,8 +6,10 @@
 
 #include "infra/util/BoundedVector.hpp"
 #include "numerical/analysis/FastFourierTransform.hpp"
+#include "numerical/math/CompilerOptimizations.hpp"
 #include "numerical/math/ComplexNumber.hpp"
 #include <cmath>
+#include <numbers>
 
 namespace analysis
 {
@@ -15,7 +17,7 @@ namespace analysis
     class DiscreteConsineTransform
     {
         static_assert((Length & (Length - 1)) == 0, "DiscreteConsineTransform size must be a power of 2");
-        static_assert(math::is_qnumber<QNumberType>::value || std::is_floating_point<QNumberType>::value,
+        static_assert(math::is_qnumber<QNumberType>::value || std::is_floating_point_v<QNumberType>,
             "DiscreteConsineTransform can only be instantiated with math::QNumber types or floating point.");
 
     public:
@@ -30,7 +32,6 @@ namespace analysis
         FastFourierTransform<QNumberType>& fft;
         typename infra::BoundedVector<QNumberType>::template WithMaxSize<Length> output;
         typename VectorComplex::template WithMaxSize<Length> complexBuffer;
-        static constexpr float PI = 3.14159265358979323846f;
     };
 
     // Implementation //
@@ -54,7 +55,7 @@ namespace analysis
 
         for (std::size_t k = 1; k < Length; ++k)
         {
-            float angle = -k * PI / (2.0f * Length);
+            float angle = -static_cast<float>(k) * std::numbers::pi_v<float> / (2.0f * Length);
             float scale = 2.0f / std::sqrt(static_cast<float>(Length));
 
             float real = math::ToFloat(fftResult[k].Real());
@@ -73,7 +74,7 @@ namespace analysis
 
         for (std::size_t k = 1; k < Length; ++k)
         {
-            float angle = k * PI / (2.0f * Length);
+            float angle = static_cast<float>(k) * std::numbers::pi_v<float> / (2.0f * static_cast<float>(Length));
             float scale = std::sqrt(static_cast<float>(Length)) / 2.0f;
 
             float value = math::ToFloat(input[k]) * scale;
@@ -82,4 +83,10 @@ namespace analysis
 
         return fft.Inverse(complexBuffer);
     }
+
+#ifdef NUMERICAL_TOOLBOX_COVERAGE_BUILD
+    extern template class DiscreteConsineTransform<float, 8>;
+    extern template class DiscreteConsineTransform<math::Q15, 8>;
+    extern template class DiscreteConsineTransform<math::Q31, 8>;
+#endif
 }
diff --git a/numerical/analysis/FastFourierTransformRadix2Impl.hpp b/numerical/analysis/FastFourierTransformRadix2Impl.hpp
index 5388ec15..a6f615ba 100644
--- a/numerical/analysis/FastFourierTransformRadix2Impl.hpp
+++ b/numerical/analysis/FastFourierTransformRadix2Impl.hpp
@@ -4,10 +4,9 @@
 #pragma GCC optimize("O3", "fast-math")
 #endif
 
-#include "numerical/math/CompilerOptimizations.hpp"
-
 #include "infra/util/BoundedVector.hpp"
 #include "numerical/analysis/FastFourierTransform.hpp"
+#include "numerical/math/CompilerOptimizations.hpp"
 #include "numerical/math/ComplexNumber.hpp"
 
 namespace analysis
@@ -17,12 +16,13 @@ namespace analysis
         : public FastFourierTransform<QNumberType>
     {
         static_assert((Length & (Length - 1)) == 0, "FastFourierTransformRadix2Impl size must be a power of 2");
+        static_assert(math::is_qnumber<QNumberType>::value || std::is_floating_point_v<QNumberType>, "QNumberType must be a floating-point or Q-number type");
 
     public:
         using VectorComplex = typename FastFourierTransform<QNumberType>::VectorComplex;
         using VectorReal = typename FastFourierTransform<QNumberType>::VectorReal;
 
-        explicit FastFourierTransformRadix2Impl(TwiddleFactors<QNumberType, Length / 2>& twinddleFactors);
+        explicit FastFourierTransformRadix2Impl(TwiddleFactors<QNumberType, Length / 2>& twiddleFactors);
 
         VectorComplex& Forward(VectorReal& input) override;
         VectorReal& Inverse(VectorComplex& input) override;
@@ -37,14 +37,14 @@ namespace analysis
         const std::size_t log2_n = FastFourierTransform<QNumberType>::Log2(Length);
         const std::size_t radix = 2;
         const std::size_t radixBits = FastFourierTransform<QNumberType>::Log2(radix);
-        TwiddleFactors<QNumberType, Length / 2>& twinddleFactors;
+        TwiddleFactors<QNumberType, Length / 2>& twiddleFactors;
         typename infra::BoundedVector<math::Complex<QNumberType>>::template WithMaxSize<Length> frequencyDomain;
         typename infra::BoundedVector<QNumberType>::template WithMaxSize<Length> timeDomain;
     };
 
     template<typename QNumberType, std::size_t Length>
-    FastFourierTransformRadix2Impl<QNumberType, Length>::FastFourierTransformRadix2Impl(TwiddleFactors<QNumberType, Length / 2>& twinddleFactors)
-        : twinddleFactors(twinddleFactors)
+    FastFourierTransformRadix2Impl<QNumberType, Length>::FastFourierTransformRadix2Impl(TwiddleFactors<QNumberType, Length / 2>& twiddleFactors)
+        : twiddleFactors(twiddleFactors)
     {}
 
     template<typename QNumberType, std::size_t Length>
@@ -73,7 +73,7 @@ namespace analysis
                 {
                     math::Complex<QNumberType>& a = frequencyDomain[j];
                     math::Complex<QNumberType>& b = frequencyDomain[j + halfStep];
-                    math::Complex twiddle = twinddleFactors[k * stepFactor];
+                    math::Complex twiddle = twiddleFactors[k * stepFactor];
 
                     math::Complex temp = b * twiddle;
                     b = (a - temp);
diff --git a/numerical/analysis/GoertzelAlgorithm.hpp b/numerical/analysis/GoertzelAlgorithm.hpp
index 0c43fbd5..911c02bf 100644
--- a/numerical/analysis/GoertzelAlgorithm.hpp
+++ b/numerical/analysis/GoertzelAlgorithm.hpp
@@ -35,9 +35,9 @@ namespace analysis
         T cosine;
         T sine;
         std::size_t blockSize;
-        T s1;
-        T s2;
-        std::size_t sampleCount;
+        T s1{ T{ 0 } };
+        T s2{ T{ 0 } };
+        std::size_t sampleCount{ 0 };
     };
 
     template<typename T>
@@ -46,9 +46,6 @@ namespace analysis
         , cosine{ std::cos(T{ 2 } * std::numbers::pi_v<T> * static_cast<T>(k) / static_cast<T>(blockLength)) }
         , sine{ std::sin(T{ 2 } * std::numbers::pi_v<T> * static_cast<T>(k) / static_cast<T>(blockLength)) }
         , blockSize{ blockLength }
-        , s1{ T{ 0 } }
-        , s2{ T{ 0 } }
-        , sampleCount{ 0 }
     {}
 
     template<typename T>
diff --git a/numerical/analysis/PowerDensitySpectrum.hpp b/numerical/analysis/PowerDensitySpectrum.hpp
index 5e732627..0c8b50ad 100644
--- a/numerical/analysis/PowerDensitySpectrum.hpp
+++ b/numerical/analysis/PowerDensitySpectrum.hpp
@@ -16,7 +16,7 @@
 
 namespace analysis
 {
-    template<typename QNumberType, std::size_t SegmentSize, typename Fft, typename TwindleFactor, std::size_t Overlap>
+    template<typename QNumberType, std::size_t SegmentSize, typename Fft, typename TwiddleFactor, std::size_t Overlap>
     class PowerSpectralDensity
     {
         static_assert(math::is_qnumber_v<QNumberType> ||
@@ -33,72 +33,89 @@ namespace analysis
             "Only values supported by overlap are: 0%, 25%, 50% and 75%.");
 
     public:
-        explicit PowerSpectralDensity(windowing::Window<QNumberType>& window, QNumberType samplingTimeInSeconds)
-            : window(window)
-            , samplingTimeInSeconds(samplingTimeInSeconds)
-            , frequencyResolution(QNumberType(samplingTimeInSeconds * segmentSizeInverted))
-        {}
-
         using VectorReal = typename FastFourierTransform<QNumberType>::VectorReal;
         using VectorComplex = typename FastFourierTransform<QNumberType>::VectorComplex;
 
-        OPTIMIZE_FOR_SPEED VectorReal& Calculate(const VectorReal& input)
-        {
-            really_assert(input.size() >= SegmentSize);
-
-            ResetOutput();
-
-            for (std::size_t i = 0, segmentCount = 0; i + SegmentSize <= input.size(); i += step, ++segmentCount)
-            {
-                ResetSegment();
-
-                for (std::size_t j = 0; j < SegmentSize; ++j)
-                    segment[j] = (QNumberType(input[i + j] * window(j, SegmentSize)));
+        explicit PowerSpectralDensity(windowing::Window<QNumberType>& window, QNumberType samplingTimeInSeconds);
 
-                auto& spectrum = fft.Forward(segment);
-
-                for (std::size_t k = 0; k <= SegmentSize / 2; ++k)
-                    y[k] += QNumberType(math::ToFloat(MagnitudeSquared(spectrum[k])) / static_cast<float>(SegmentSize));
-            }
-
-            for (std::size_t i = 0; i < y.size(); ++i)
-                y[i] *= frequencyResolution;
-
-            return y;
-        }
+        OPTIMIZE_FOR_SPEED VectorReal& Calculate(const VectorReal& input);
 
     private:
-        QNumberType MagnitudeSquared(const math::Complex<QNumberType>& data) const
-        {
-            return data.Real() * data.Real() + data.Imaginary() * data.Imaginary();
-        }
-
-        void ResetOutput()
-        {
-            y.clear();
-            y.resize(SegmentSize / 2 + 1);
-        }
-
-        void ResetSegment()
-        {
-            segment.clear();
-            segment.resize(SegmentSize);
-        }
+        QNumberType MagnitudeSquared(const math::Complex<QNumberType>& data) const;
+        void ResetOutput();
+        void ResetSegment();
 
     private:
         static constexpr QNumberType segmentSizeInverted = QNumberType(1.0f / static_cast<float>(SegmentSize));
         static constexpr std::size_t overlapSize = (SegmentSize * Overlap) / 100;
         static constexpr std::size_t step = SegmentSize - overlapSize;
         QNumberType frequencyResolution;
-        QNumberType scale;
         windowing::Window<QNumberType>& window;
         QNumberType samplingTimeInSeconds;
-        TwindleFactor twindleFacors;
-        Fft fft{ twindleFacors };
+        TwiddleFactor twiddleFactors;
+        Fft fft{ twiddleFactors };
         typename infra::BoundedVector<QNumberType>::template WithMaxSize<SegmentSize> segment;
         typename FastFourierTransform<QNumberType>::VectorReal::template WithMaxSize<SegmentSize / 2 + 1> y;
     };
 
+    /// Implementation ///
+
+    template<typename QNumberType, std::size_t SegmentSize, typename Fft, typename TwiddleFactor, std::size_t Overlap>
+    PowerSpectralDensity<QNumberType, SegmentSize, Fft, TwiddleFactor, Overlap>::PowerSpectralDensity(
+        windowing::Window<QNumberType>& window, QNumberType samplingTimeInSeconds)
+        : window(window)
+        , samplingTimeInSeconds(samplingTimeInSeconds)
+        , frequencyResolution(QNumberType(samplingTimeInSeconds * segmentSizeInverted))
+    {}
+
+    template<typename QNumberType, std::size_t SegmentSize, typename Fft, typename TwiddleFactor, std::size_t Overlap>
+    OPTIMIZE_FOR_SPEED typename PowerSpectralDensity<QNumberType, SegmentSize, Fft, TwiddleFactor, Overlap>::VectorReal&
+    PowerSpectralDensity<QNumberType, SegmentSize, Fft, TwiddleFactor, Overlap>::Calculate(const VectorReal& input)
+    {
+        really_assert(input.size() >= SegmentSize);
+
+        ResetOutput();
+
+        for (std::size_t i = 0; i + SegmentSize <= input.size(); i += step)
+        {
+            ResetSegment();
+
+            for (std::size_t j = 0; j < SegmentSize; ++j)
+                segment[j] = QNumberType(input[i + j] * window(j, SegmentSize));
+
+            auto& spectrum = fft.Forward(segment);
+
+            for (std::size_t k = 0; k <= SegmentSize / 2; ++k)
+                y[k] += QNumberType(math::ToFloat(MagnitudeSquared(spectrum[k])) / static_cast<float>(SegmentSize));
+        }
+
+        for (std::size_t i = 0; i < y.size(); ++i)
+            y[i] *= frequencyResolution;
+
+        return y;
+    }
+
+    template<typename QNumberType, std::size_t SegmentSize, typename Fft, typename TwiddleFactor, std::size_t Overlap>
+    QNumberType PowerSpectralDensity<QNumberType, SegmentSize, Fft, TwiddleFactor, Overlap>::MagnitudeSquared(
+        const math::Complex<QNumberType>& data) const
+    {
+        return data.Real() * data.Real() + data.Imaginary() * data.Imaginary();
+    }
+
+    template<typename QNumberType, std::size_t SegmentSize, typename Fft, typename TwiddleFactor, std::size_t Overlap>
+    void PowerSpectralDensity<QNumberType, SegmentSize, Fft, TwiddleFactor, Overlap>::ResetOutput()
+    {
+        y.clear();
+        y.resize(SegmentSize / 2 + 1);
+    }
+
+    template<typename QNumberType, std::size_t SegmentSize, typename Fft, typename TwiddleFactor, std::size_t Overlap>
+    void PowerSpectralDensity<QNumberType, SegmentSize, Fft, TwiddleFactor, Overlap>::ResetSegment()
+    {
+        segment.clear();
+        segment.resize(SegmentSize);
+    }
+
 #ifdef NUMERICAL_TOOLBOX_COVERAGE_BUILD
     extern template class PowerSpectralDensity<float, 512, test::FftStub<float, 512>, test::TwiddleFactorsStub<float, 256>, 50>;
 #endif
diff --git a/numerical/analysis/RealFastFourierTransform.cpp b/numerical/analysis/RealFastFourierTransform.cpp
new file mode 100644
index 00000000..058aacb7
--- /dev/null
+++ b/numerical/analysis/RealFastFourierTransform.cpp
@@ -0,0 +1,6 @@
+#include "numerical/analysis/RealFastFourierTransform.hpp"
+
+namespace analysis
+{
+    template class RealFastFourierTransform<float, 16>;
+}
diff --git a/numerical/analysis/RealFft.hpp b/numerical/analysis/RealFastFourierTransform.hpp
similarity index 75%
rename from numerical/analysis/RealFft.hpp
rename to numerical/analysis/RealFastFourierTransform.hpp
index e493972b..5b11f3c6 100644
--- a/numerical/analysis/RealFft.hpp
+++ b/numerical/analysis/RealFastFourierTransform.hpp
@@ -8,23 +8,21 @@
 #include "numerical/analysis/FastFourierTransform.hpp"
 #include "numerical/math/CompilerOptimizations.hpp"
 #include "numerical/math/ComplexNumber.hpp"
-#include <cmath>
-#include <numbers>
 
 namespace analysis
 {
     template<typename T, std::size_t N>
-    class RealFft
+    class RealFastFourierTransform
     {
-        static_assert(std::is_floating_point_v<T>, "RealFft supports floating-point types");
-        static_assert(N >= 4 && (N & (N - 1)) == 0, "RealFft length N must be an even power of two");
+        static_assert(std::is_floating_point_v<T>, "RealFastFourierTransform supports floating-point types");
+        static_assert(N >= 4 && (N & (N - 1)) == 0, "RealFastFourierTransform length N must be an even power of two");
 
     public:
         using Complex = math::Complex<T>;
         using VectorComplex = infra::BoundedVector<Complex>;
         using VectorReal = infra::BoundedVector<T>;
 
-        explicit RealFft(FastFourierTransform<T>& engine);
+        explicit RealFastFourierTransform(FastFourierTransform<T>& engine, TwiddleFactors<T, N / 2>& twiddleFactors);
 
         OPTIMIZE_FOR_SPEED VectorComplex& Forward(VectorReal& realInput);
         OPTIMIZE_FOR_SPEED VectorReal& Inverse(VectorComplex& halfSpectrum);
@@ -33,35 +31,30 @@ namespace analysis
         static constexpr std::size_t half{ N / 2 };
 
         FastFourierTransform<T>& engine;
+        TwiddleFactors<T, half>& twiddleFactors;
         typename VectorReal::template WithMaxSize<half> evenBuf;
         typename VectorReal::template WithMaxSize<half> oddBuf;
-        typename VectorComplex::template WithMaxSize<half> twiddles;
         typename VectorComplex::template WithMaxSize<half> eCopy;
         typename VectorComplex::template WithMaxSize<half + 1> spectrum;
         typename VectorReal::template WithMaxSize<N> timeDomain;
     };
 
+    /// Implementation ///
+
     template<typename T, std::size_t N>
-    RealFft<T, N>::RealFft(FastFourierTransform<T>& engine)
+    RealFastFourierTransform<T, N>::RealFastFourierTransform(FastFourierTransform<T>& engine, TwiddleFactors<T, N / 2>& twiddleFactors)
         : engine{ engine }
+        , twiddleFactors{ twiddleFactors }
     {
         evenBuf.resize(half);
         oddBuf.resize(half);
-        twiddles.resize(half);
         eCopy.resize(half);
         spectrum.resize(half + 1);
         timeDomain.resize(N);
-
-        twiddles[0] = Complex{ T(1), T(0) };
-        for (std::size_t k{ 1 }; k < half; ++k)
-        {
-            T angle{ -T(2) * std::numbers::pi_v<T> * static_cast<T>(k) / static_cast<T>(N) };
-            twiddles[k] = Complex{ std::cos(angle), std::sin(angle) };
-        }
     }
 
     template<typename T, std::size_t N>
-    OPTIMIZE_FOR_SPEED typename RealFft<T, N>::VectorComplex& RealFft<T, N>::Forward(VectorReal& realInput)
+    OPTIMIZE_FOR_SPEED typename RealFastFourierTransform<T, N>::VectorComplex& RealFastFourierTransform<T, N>::Forward(VectorReal& realInput)
     {
         for (std::size_t n{ 0 }; n < half; ++n)
         {
@@ -80,7 +73,7 @@ namespace analysis
 
         for (std::size_t k{ 1 }; k < half; ++k)
         {
-            Complex twidO{ twiddles[k] * O[k] };
+            Complex twidO{ twiddleFactors[k] * O[k] };
             spectrum[k] = Complex{ eCopy[k].Real() + twidO.Real(), eCopy[k].Imaginary() + twidO.Imaginary() };
         }
 
@@ -88,7 +81,7 @@ namespace analysis
     }
 
     template<typename T, std::size_t N>
-    OPTIMIZE_FOR_SPEED typename RealFft<T, N>::VectorReal& RealFft<T, N>::Inverse(VectorComplex& halfSpectrum)
+    OPTIMIZE_FOR_SPEED typename RealFastFourierTransform<T, N>::VectorReal& RealFastFourierTransform<T, N>::Inverse(VectorComplex& halfSpectrum)
     {
         typename VectorComplex::template WithMaxSize<half> E;
         typename VectorComplex::template WithMaxSize<half> O;
@@ -104,7 +97,7 @@ namespace analysis
             Complex xNk{ halfSpectrum[half - k].Real(), -halfSpectrum[half - k].Imaginary() };
             Complex sum{ T(0.5) * (xk.Real() + xNk.Real()), T(0.5) * (xk.Imaginary() + xNk.Imaginary()) };
             Complex diff{ T(0.5) * (xk.Real() - xNk.Real()), T(0.5) * (xk.Imaginary() - xNk.Imaginary()) };
-            Complex twConj{ twiddles[k].Real(), -twiddles[k].Imaginary() };
+            Complex twConj{ twiddleFactors[k].Real(), -twiddleFactors[k].Imaginary() };
             E[k] = sum;
             O[k] = twConj * diff;
         }
@@ -127,6 +120,6 @@ namespace analysis
     }
 
 #ifdef NUMERICAL_TOOLBOX_COVERAGE_BUILD
-    extern template class RealFft<float, 16>;
+    extern template class RealFastFourierTransform<float, 16>;
 #endif
 }
diff --git a/numerical/analysis/RealFft.cpp b/numerical/analysis/RealFft.cpp
deleted file mode 100644
index 9e0b5750..00000000
--- a/numerical/analysis/RealFft.cpp
+++ /dev/null
@@ -1,6 +0,0 @@
-#include "numerical/analysis/RealFft.hpp"
-
-namespace analysis
-{
-    template class RealFft<float, 16>;
-}
diff --git a/numerical/analysis/SignalDetectors.hpp b/numerical/analysis/SignalDetectors.hpp
index 0f82ddae..1b68b2cc 100644
--- a/numerical/analysis/SignalDetectors.hpp
+++ b/numerical/analysis/SignalDetectors.hpp
@@ -17,21 +17,10 @@ namespace analysis
         static_assert(std::is_floating_point_v<T>, "PeakHold supports floating-point types only");
 
     public:
-        explicit PeakHold(T decay = T{ 1 })
-            : decay{ decay }
-        {}
-
-        OPTIMIZE_FOR_SPEED T Update(T x)
-        {
-            T m{ std::abs(x) };
-            peak = (m > peak * decay) ? m : peak * decay;
-            return peak;
-        }
-
-        void Reset()
-        {
-            peak = T{ 0 };
-        }
+        explicit PeakHold(T decay = T{ 1 });
+
+        OPTIMIZE_FOR_SPEED T Update(T x);
+        void Reset();
 
     private:
         T peak{ T{ 0 } };
@@ -44,29 +33,11 @@ namespace analysis
         static_assert(std::is_floating_point_v<T>, "ZeroCrossingCounter supports floating-point types only");
 
     public:
-        explicit ZeroCrossingCounter(T hysteresis = T{ 0 })
-            : hysteresis{ hysteresis }
-        {}
-
-        OPTIMIZE_FOR_SPEED bool Update(T x)
-        {
-            bool crossed{ (previous < T{ 0 } ? x > T{ 0 } : x < T{ 0 }) && (std::abs(x) > hysteresis) };
-            if (crossed)
-                ++count;
-            previous = x;
-            return crossed;
-        }
-
-        uint32_t Count() const
-        {
-            return count;
-        }
-
-        void Reset()
-        {
-            previous = T{ 0 };
-            count = 0;
-        }
+        explicit ZeroCrossingCounter(T hysteresis = T{ 0 });
+
+        OPTIMIZE_FOR_SPEED bool Update(T x);
+        uint32_t Count() const;
+        void Reset();
 
     private:
         T previous{ T{ 0 } };
@@ -80,26 +51,83 @@ namespace analysis
         static_assert(std::is_floating_point_v<T>, "RmsEnvelope supports floating-point types only");
 
     public:
-        explicit RmsEnvelope(T alpha)
-            : alpha{ alpha }
-        {}
+        explicit RmsEnvelope(T alpha);
 
-        OPTIMIZE_FOR_SPEED T Update(T x)
-        {
-            meanSquare += alpha * (x * x - meanSquare);
-            return std::sqrt(meanSquare);
-        }
-
-        void Reset(T value = T{ 0 })
-        {
-            meanSquare = value;
-        }
+        OPTIMIZE_FOR_SPEED T Update(T x);
+        void Reset(T value = T{ 0 });
 
     private:
         T alpha;
         T meanSquare{ T{ 0 } };
     };
 
+    /// Implementation ///
+
+    template<typename T>
+    PeakHold<T>::PeakHold(T decay)
+        : decay{ decay }
+    {}
+
+    template<typename T>
+    OPTIMIZE_FOR_SPEED T PeakHold<T>::Update(T x)
+    {
+        T m{ std::abs(x) };
+        peak = (m > peak * decay) ? m : peak * decay;
+        return peak;
+    }
+
+    template<typename T>
+    void PeakHold<T>::Reset()
+    {
+        peak = T{ 0 };
+    }
+
+    template<typename T>
+    ZeroCrossingCounter<T>::ZeroCrossingCounter(T hysteresis)
+        : hysteresis{ hysteresis }
+    {}
+
+    template<typename T>
+    OPTIMIZE_FOR_SPEED bool ZeroCrossingCounter<T>::Update(T x)
+    {
+        bool crossed{ (previous < T{ 0 } ? x > T{ 0 } : x < T{ 0 }) && (std::abs(x) > hysteresis) };
+        if (crossed)
+            ++count;
+        previous = x;
+        return crossed;
+    }
+
+    template<typename T>
+    uint32_t ZeroCrossingCounter<T>::Count() const
+    {
+        return count;
+    }
+
+    template<typename T>
+    void ZeroCrossingCounter<T>::Reset()
+    {
+        previous = T{ 0 };
+        count = 0;
+    }
+
+    template<typename T>
+    RmsEnvelope<T>::RmsEnvelope(T alpha)
+        : alpha{ alpha }
+    {}
+
+    template<typename T>
+    OPTIMIZE_FOR_SPEED T RmsEnvelope<T>::Update(T x)
+    {
+        meanSquare += alpha * (x * x - meanSquare);
+        return std::sqrt(meanSquare);
+    }
+
+    template<typename T>
+    void RmsEnvelope<T>::Reset(T value)
+    {
+        meanSquare = value;
+    }
+
 #ifdef NUMERICAL_TOOLBOX_COVERAGE_BUILD
     extern template class PeakHold<float>;
     extern template class ZeroCrossingCounter<float>;
diff --git a/numerical/analysis/test/CMakeLists.txt b/numerical/analysis/test/CMakeLists.txt
index 8458ad3b..85403209 100644
--- a/numerical/analysis/test/CMakeLists.txt
+++ b/numerical/analysis/test/CMakeLists.txt
@@ -15,6 +15,6 @@ target_sources(numerical.analysis_test PRIVATE
     TestFastFourierTransformRadix2Impl.cpp
     TestGoertzelAlgorithm.cpp
     TestPowerDensitySpectrum.cpp
-    TestRealFft.cpp
+    TestRealFastFourierTransform.cpp
     TestSignalDetectors.cpp
 )
diff --git a/numerical/analysis/test/TestDiscreteCosineTransform.cpp b/numerical/analysis/test/TestDiscreteCosineTransform.cpp
index b9acd300..2adc34b6 100644
--- a/numerical/analysis/test/TestDiscreteCosineTransform.cpp
+++ b/numerical/analysis/test/TestDiscreteCosineTransform.cpp
@@ -22,11 +22,32 @@ namespace
         using VectorReal = typename analysis::FastFourierTransform<T>::VectorReal;
         using VectorComplex = typename analysis::FastFourierTransform<T>::VectorComplex;
 
-        MockFft<T> mockFft;
+        ::testing::StrictMock<MockFft<T>> mockFft;
         std::optional<analysis::DiscreteConsineTransform<T, Length>> dct;
         typename VectorReal::template WithMaxSize<Length> input;
         typename VectorReal::template WithMaxSize<Length> output;
         typename VectorComplex::template WithMaxSize<Length> fftOutput;
+
+        void SetUp() override
+        {
+            dct.emplace(mockFft);
+            input.clear();
+            input.resize(Length);
+        }
+
+        typename VectorComplex::template WithMaxSize<Length>& EmptyForwardResult()
+        {
+            fftOutput.clear();
+            fftOutput.resize(fftOutput.max_size());
+            return fftOutput;
+        }
+
+        typename VectorReal::template WithMaxSize<Length>& EmptyInverseResult()
+        {
+            output.clear();
+            output.resize(output.max_size());
+            return output;
+        }
     };
 
     using TestedTypes = ::testing::Types<float, math::Q15, math::Q31>;
@@ -35,10 +56,6 @@ namespace
 
 TYPED_TEST(TestDiscreteConsineTransform, forward_transform_calls_fft_with_preprocessed_data)
 {
-    this->dct.emplace(this->mockFft);
-    this->input.clear();
-    this->input.resize(TestDiscreteConsineTransform<TypeParam>::Length);
-
     for (size_t i = 0; i < this->input.size(); ++i)
         this->input[i] = TypeParam(0.1f);
 
@@ -48,9 +65,7 @@ TYPED_TEST(TestDiscreteConsineTransform, forward_transform_calls_fft_with_prepro
                 for (size_t i = 0; i < input.size() / 2; ++i)
                     EXPECT_EQ(math::ToFloat(input[i]), math::ToFloat(this->input[2 * i]));
 
-                this->fftOutput.clear();
-                this->fftOutput.resize(this->fftOutput.max_size());
-                return this->fftOutput;
+                return this->EmptyForwardResult();
             }));
 
     this->dct->Forward(this->input);
@@ -58,10 +73,6 @@ TYPED_TEST(TestDiscreteConsineTransform, forward_transform_calls_fft_with_prepro
 
 TYPED_TEST(TestDiscreteConsineTransform, inverse_transform_calls_fft_with_preprocessed_data)
 {
-    this->dct.emplace(this->mockFft);
-    this->input.clear();
-    this->input.resize(TestDiscreteConsineTransform<TypeParam>::Length);
-
     for (size_t i = 0; i < this->input.size(); ++i)
         this->input[i] = TypeParam(0.1f);
 
@@ -70,9 +81,7 @@ TYPED_TEST(TestDiscreteConsineTransform, inverse_transform_calls_fft_with_prepro
             {
                 EXPECT_LE(std::abs(math::ToFloat(input[0].Real())), 1.0f);
 
-                this->output.clear();
-                this->output.resize(this->output.max_size());
-                return this->output;
+                return this->EmptyInverseResult();
             }));
 
     this->dct->Inverse(this->input);
@@ -80,17 +89,11 @@ TYPED_TEST(TestDiscreteConsineTransform, inverse_transform_calls_fft_with_prepro
 
 TYPED_TEST(TestDiscreteConsineTransform, input_size_matches_fft_size)
 {
-    this->dct.emplace(this->mockFft);
-    this->input.clear();
-    this->input.resize(TestDiscreteConsineTransform<TypeParam>::Length);
-
     EXPECT_CALL(this->mockFft, Forward(::testing::_))
         .WillOnce(::testing::Invoke([this](typename analysis::FastFourierTransform<TypeParam>::VectorReal& input) -> typename analysis::FastFourierTransform<TypeParam>::VectorComplex&
             {
                 EXPECT_EQ(input.size(), TestDiscreteConsineTransform<TypeParam>::Length);
-                this->fftOutput.clear();
-                this->fftOutput.resize(this->fftOutput.max_size());
-                return this->fftOutput;
+                return this->EmptyForwardResult();
             }));
 
     this->dct->Forward(this->input);
@@ -98,10 +101,6 @@ TYPED_TEST(TestDiscreteConsineTransform, input_size_matches_fft_size)
 
 TYPED_TEST(TestDiscreteConsineTransform, scaling_stays_within_fixed_point_range)
 {
-    this->dct.emplace(this->mockFft);
-    this->input.clear();
-    this->input.resize(TestDiscreteConsineTransform<TypeParam>::Length);
-
     std::fill(this->input.begin(), this->input.end(), TypeParam(0.9999f));
 
     EXPECT_CALL(this->mockFft, Forward(::testing::_))
@@ -110,9 +109,7 @@ TYPED_TEST(TestDiscreteConsineTransform, scaling_stays_within_fixed_point_range)
                 for (const auto& value : input)
                     EXPECT_LE(std::abs(math::ToFloat(value)), 1.0f);
 
-                this->fftOutput.clear();
-                this->fftOutput.resize(this->fftOutput.max_size());
-                return this->fftOutput;
+                return this->EmptyForwardResult();
             }));
 
     this->dct->Forward(this->input);
diff --git a/numerical/analysis/test/TestFastFourierTransformRadix2Impl.cpp b/numerical/analysis/test/TestFastFourierTransformRadix2Impl.cpp
index 7e61f5bc..0c57cfa5 100644
--- a/numerical/analysis/test/TestFastFourierTransformRadix2Impl.cpp
+++ b/numerical/analysis/test/TestFastFourierTransformRadix2Impl.cpp
@@ -46,6 +46,11 @@ namespace
         using VectorComplex = typename analysis::FastFourierTransform<T>::VectorComplex;
         using VectorReal = typename analysis::FastFourierTransform<T>::VectorReal;
 
+        void SetUp() override
+        {
+            fft.emplace(twiddleFactors);
+        }
+
         MockTwiddleFactors<T> twiddleFactors{};
         std::optional<analysis::FastFourierTransformRadix2Impl<T, Length>> fft;
         typename VectorReal::template WithMaxSize<Length> timeDomain;
@@ -58,7 +63,6 @@ namespace
 
 TYPED_TEST(TestFastFourierTransform, zero_input_produces_zero_output)
 {
-    this->fft.emplace(this->twiddleFactors);
     this->timeDomain.clear();
     this->timeDomain.resize(TestFastFourierTransform<TypeParam>::Length);
 
@@ -70,7 +74,6 @@ TYPED_TEST(TestFastFourierTransform, zero_input_produces_zero_output)
 
 TYPED_TEST(TestFastFourierTransform, dc_signal_appears_in_zero_frequency_bin)
 {
-    this->fft.emplace(this->twiddleFactors);
     this->timeDomain.clear();
     this->timeDomain.resize(TestFastFourierTransform<TypeParam>::Length);
     std::fill(this->timeDomain.begin(), this->timeDomain.end(), TypeParam(0.1f));
@@ -87,7 +90,6 @@ TYPED_TEST(TestFastFourierTransform, dc_signal_appears_in_zero_frequency_bin)
 
 TYPED_TEST(TestFastFourierTransform, forward_and_inverse_transform_recovers_original_signal)
 {
-    this->fft.emplace(this->twiddleFactors);
     constexpr std::array<float, TestFastFourierTransform<TypeParam>::Length> signal = {
         0.1f, 0.07f, 0.0f, -0.07f, -0.1f, -0.07f, 0.0f, 0.07f
     };
@@ -105,7 +107,6 @@ TYPED_TEST(TestFastFourierTransform, forward_and_inverse_transform_recovers_orig
 
 TYPED_TEST(TestFastFourierTransform, nyquist_frequency_detection)
 {
-    this->fft.emplace(this->twiddleFactors);
     this->timeDomain.clear();
 
     for (size_t i = 0; i < TestFastFourierTransform<TypeParam>::Length; ++i)
diff --git a/numerical/analysis/test/TestRealFft.cpp b/numerical/analysis/test/TestRealFastFourierTransform.cpp
similarity index 81%
rename from numerical/analysis/test/TestRealFft.cpp
rename to numerical/analysis/test/TestRealFastFourierTransform.cpp
index 4dc9b0e1..4f338307 100644
--- a/numerical/analysis/test/TestRealFft.cpp
+++ b/numerical/analysis/test/TestRealFastFourierTransform.cpp
@@ -1,5 +1,5 @@
 #include "numerical/analysis/FastFourierTransformRadix2Impl.hpp"
-#include "numerical/analysis/RealFft.hpp"
+#include "numerical/analysis/RealFastFourierTransform.hpp"
 #include "numerical/math/Tolerance.hpp"
 #include "gmock/gmock.h"
 #include <array>
@@ -29,19 +29,41 @@ namespace
         std::array<math::Complex<float>, 4> factors{};
     };
 
-    class TestRealFft : public ::testing::Test
+    class TwiddleFactors16 : public analysis::TwiddleFactors<float, 8>
+    {
+    public:
+        TwiddleFactors16()
+        {
+            for (std::size_t k{ 0 }; k < 8; ++k)
+            {
+                float angle{ -2.0f * std::numbers::pi_v<float> * static_cast<float>(k) / 16.0f };
+                factors[k] = math::Complex<float>{ std::cos(angle), std::sin(angle) };
+            }
+        }
+
+        math::Complex<float>& operator[](std::size_t n) override
+        {
+            return factors[n];
+        }
+
+    private:
+        std::array<math::Complex<float>, 8> factors{};
+    };
+
+    class TestRealFastFourierTransform : public ::testing::Test
     {
     public:
         static constexpr std::size_t N{ 16 };
 
         TwiddleFactors8 twiddleFactors{};
+        TwiddleFactors16 realTwiddleFactors{};
         analysis::FastFourierTransformRadix2Impl<float, N / 2> engine{ twiddleFactors };
-        analysis::RealFft<float, N> rfft{ engine };
+        analysis::RealFastFourierTransform<float, N> rfft{ engine, realTwiddleFactors };
         typename infra::BoundedVector<float>::template WithMaxSize<N> input{};
     };
 }
 
-TEST_F(TestRealFft, dc_input_has_only_dc_bin)
+TEST_F(TestRealFastFourierTransform, dc_input_has_only_dc_bin)
 {
     input.resize(N, 1.0f);
 
@@ -55,7 +77,7 @@ TEST_F(TestRealFft, dc_input_has_only_dc_bin)
     }
 }
 
-TEST_F(TestRealFft, single_real_sinusoid_hits_one_bin)
+TEST_F(TestRealFastFourierTransform, single_real_sinusoid_hits_one_bin)
 {
     input.resize(N);
     for (std::size_t n{ 0 }; n < N; ++n)
@@ -73,7 +95,7 @@ TEST_F(TestRealFft, single_real_sinusoid_hits_one_bin)
     }
 }
 
-TEST_F(TestRealFft, matches_reference_complex_fft)
+TEST_F(TestRealFastFourierTransform, matches_reference_complex_fft)
 {
     static constexpr std::array<float, N> signal{ 0.1f, 0.3f, -0.2f, 0.5f, 0.4f, -0.1f, 0.0f, 0.2f,
         -0.3f, 0.6f, 0.1f, -0.4f, 0.2f, 0.0f, -0.1f, 0.3f };
@@ -98,7 +120,7 @@ TEST_F(TestRealFft, matches_reference_complex_fft)
     }
 }
 
-TEST_F(TestRealFft, nyquist_and_dc_are_real)
+TEST_F(TestRealFastFourierTransform, nyquist_and_dc_are_real)
 {
     input.resize(N);
     for (std::size_t n{ 0 }; n < N; ++n)
@@ -110,7 +132,7 @@ TEST_F(TestRealFft, nyquist_and_dc_are_real)
     EXPECT_NEAR(spectrum[N / 2].Imaginary(), 0.0f, math::Tolerance<float>());
 }
 
-TEST_F(TestRealFft, inverse_is_left_inverse)
+TEST_F(TestRealFastFourierTransform, inverse_is_left_inverse)
 {
     static constexpr std::array<float, N> signal{ 0.5f, -0.3f, 0.1f, 0.8f, -0.6f, 0.2f, 0.0f, -0.4f,
         0.7f, -0.1f, 0.3f, -0.5f, 0.4f, 0.1f, -0.2f, 0.6f };
@@ -125,7 +147,7 @@ TEST_F(TestRealFft, inverse_is_left_inverse)
         EXPECT_NEAR(recovered[n], signal[n], 1e-2f);
 }
 
-TEST_F(TestRealFft, linearity_holds)
+TEST_F(TestRealFastFourierTransform, linearity_holds)
 {
     static constexpr float a{ 2.0f };
     static constexpr float b{ 0.5f };
@@ -171,7 +193,7 @@ TEST_F(TestRealFft, linearity_holds)
     }
 }
 
-TEST_F(TestRealFft, impulse_has_flat_magnitude)
+TEST_F(TestRealFastFourierTransform, impulse_has_flat_magnitude)
 {
     input.resize(N, 0.0f);
     input[0] = 1.0f;
diff --git a/roadmap/README.md b/roadmap/README.md
index 15ac5a60..dc847c42 100644
--- a/roadmap/README.md
+++ b/roadmap/README.md
@@ -96,7 +96,7 @@ with `EXPECT_NEAR` and `math::Tolerance<float>()` (or an explicit tolerance).
 `FeedbackLinearization` (40) · `BacksteppingControl` (41) · `ModelReferenceAdaptiveControl` (47)
 
 ### `analysis`
-`SignalDetectors` (5) · `ConvolutionCorrelation` (11) · `GoertzelAlgorithm` (13) · `RealFft` (25) · `HilbertTransform` (37) · `DiscreteWaveletTransform` (38)
+`SignalDetectors` (5) · `ConvolutionCorrelation` (11) · `GoertzelAlgorithm` (13) · `RealFastFourierTransform` (25) · `HilbertTransform` (37) · `DiscreteWaveletTransform` (38)
 
 ### `estimators/offline`
 `PolynomialFitting` (12) · `TotalLeastSquares` (44) · `DynamicParameterIdentification` (M22)
diff --git a/roadmap/analysis/DiscreteWaveletTransform/explanation.md b/roadmap/analysis/DiscreteWaveletTransform/explanation.md
index 3890eef9..ab29e17d 100644
--- a/roadmap/analysis/DiscreteWaveletTransform/explanation.md
+++ b/roadmap/analysis/DiscreteWaveletTransform/explanation.md
@@ -31,5 +31,5 @@ Trans. Pattern Analysis and Machine Intelligence*, 11(7), 674–693, 1989; I. Da
 on Wavelets*, SIAM, 1992.
 
 ## See also
-`ConvolutionCorrelation` (the FIR filtering underneath), `FastFourierTransform` / `RealFft`
+`ConvolutionCorrelation` (the FIR filtering underneath), `FastFourierTransform` / `RealFastFourierTransform`
 (alternative time-frequency view), `math::RecursiveBuffer`.
diff --git a/roadmap/analysis/HilbertTransform/explanation.md b/roadmap/analysis/HilbertTransform/explanation.md
index 4cc913d5..695e1ec8 100644
--- a/roadmap/analysis/HilbertTransform/explanation.md
+++ b/roadmap/analysis/HilbertTransform/explanation.md
@@ -33,5 +33,5 @@ S. L. Marple, "Computing the Discrete-Time Analytic Signal via FFT," *IEEE Trans
 47(9), 2600–2603, 1999.
 
 ## See also
-`RealFft` / `FastFourierTransform` (FFT back-end), `ConvolutionCorrelation` (FIR realisation),
+`RealFastFourierTransform` / `FastFourierTransform` (FFT back-end), `ConvolutionCorrelation` (FIR realisation),
 `SignalDetectors` (cheap non-coherent envelope).
diff --git a/roadmap/analysis/HilbertTransform/implementation.md b/roadmap/analysis/HilbertTransform/implementation.md
index 0b74c6a1..d5be3545 100644
--- a/roadmap/analysis/HilbertTransform/implementation.md
+++ b/roadmap/analysis/HilbertTransform/implementation.md
@@ -71,7 +71,7 @@ function InstantaneousFrequency(phaseNow, phasePrev, Ts):
   streams continuously but only approximates the ideal `−90°` shift across its passband.
 - Instantaneous **phase must be unwrapped** before differencing, else frequency spikes by `2π`.
 - FIR Hilbert taps are antisymmetric with every even tap zero — exploit that to halve the MACs.
-- Reuses `FastFourierTransform` / `RealFft` for the block back-end.
+- Reuses `FastFourierTransform` / `RealFastFourierTransform` for the block back-end.
 - Float-only: `static_assert(std::is_floating_point_v<T>)`; the generic `T` signature keeps a
   `Q15`/`Q31` specialisation cheap to add later.