commit no. 2

src/dsp/FFTProcessor.cpp

@@ -38,10 +38,13 @@ void FFTProcessor::configure(int fftSize, bool complexInput) {
     }
 }
 
-void FFTProcessor::processReal(const float* input, std::vector<float>& outputDB) {
+void FFTProcessor::processReal(const float* input,
+                               std::vector<float>& outputDB,
+                               std::vector<std::complex<float>>& outputCplx) {
     const int N = fftSize_;
     const int bins = N / 2 + 1;
     outputDB.resize(bins);
+    outputCplx.resize(bins);
 
     std::copy(input, input + N, realIn_);
     fftwf_execute(realPlan_);
@@ -50,17 +53,19 @@ void FFTProcessor::processReal(const float* input, std::vector<float>& outputDB)
     for (int i = 0; i < bins; ++i) {
         float re = realOut_[i][0] * scale;
         float im = realOut_[i][1] * scale;
+        outputCplx[i] = {re, im};
         float mag2 = re * re + im * im;
-        // Power in dB, floor at -200 dB
         outputDB[i] = (mag2 > 1e-20f) ? 10.0f * std::log10(mag2) : -200.0f;
     }
 }
 
-void FFTProcessor::processComplex(const float* inputIQ, std::vector<float>& outputDB) {
+void FFTProcessor::processComplex(const float* inputIQ,
+                                  std::vector<float>& outputDB,
+                                  std::vector<std::complex<float>>& outputCplx) {
     const int N = fftSize_;
     outputDB.resize(N);
+    outputCplx.resize(N);
 
-    // Copy interleaved I/Q into FFTW complex array
     for (int i = 0; i < N; ++i) {
         cplxIn_[i][0] = inputIQ[2 * i];
         cplxIn_[i][1] = inputIQ[2 * i + 1];
@@ -68,18 +73,27 @@ void FFTProcessor::processComplex(const float* inputIQ, std::vector<float>& outp
 
     fftwf_execute(cplxPlan_);
 
-    // FFT-shift: reorder so DC is in center.
-    // FFTW output: [0, 1, ..., N/2-1, -N/2, ..., -1]
-    // Shifted:     [-N/2, ..., -1, 0, 1, ..., N/2-1]
     const float scale = 1.0f / N;
     const int half = N / 2;
     for (int i = 0; i < N; ++i) {
         int src = (i + half) % N;
         float re = cplxOut_[src][0] * scale;
         float im = cplxOut_[src][1] * scale;
+        outputCplx[i] = {re, im};
         float mag2 = re * re + im * im;
         outputDB[i] = (mag2 > 1e-20f) ? 10.0f * std::log10(mag2) : -200.0f;
     }
 }
 
+// Convenience overloads (no complex output).
+void FFTProcessor::processReal(const float* input, std::vector<float>& outputDB) {
+    std::vector<std::complex<float>> dummy;
+    processReal(input, outputDB, dummy);
+}
+
+void FFTProcessor::processComplex(const float* inputIQ, std::vector<float>& outputDB) {
+    std::vector<std::complex<float>> dummy;
+    processComplex(inputIQ, outputDB, dummy);
+}
+
 } // namespace baudline

src/dsp/FFTProcessor.h

@@ -2,13 +2,13 @@
 
 #include "core/Types.h"
 #include <fftw3.h>
-#include <vector>
 #include <complex>
+#include <vector>
 
 namespace baudline {
 
-// Wraps FFTW for real→complex and complex→complex transforms.
+// Wraps FFTW for real->complex and complex->complex transforms.
-// Produces magnitude output in dB.
+// Produces magnitude output in dB and optionally retains the complex spectrum.
 class FFTProcessor {
 public:
     FFTProcessor();
@@ -17,7 +17,6 @@ public:
     FFTProcessor(const FFTProcessor&) = delete;
     FFTProcessor& operator=(const FFTProcessor&) = delete;
 
-    // Reconfigure for a new FFT size and mode.  Rebuilds FFTW plans.
     void configure(int fftSize, bool complexInput);
 
     int  fftSize()      const { return fftSize_; }
@@ -25,27 +24,27 @@ public:
     int  outputBins()   const { return complexInput_ ? fftSize_ : fftSize_ / 2 + 1; }
     int  spectrumSize() const { return complexInput_ ? fftSize_ : fftSize_ / 2 + 1; }
 
-    // Process windowed real samples → magnitude dB spectrum.
+    // Process and produce both dB magnitude and complex spectrum.
-    // `input` must have fftSize_ elements.
+    void processReal(const float* input,
-    // `outputDB` will be resized to spectrumSize().
+                     std::vector<float>& outputDB,
-    void processReal(const float* input, std::vector<float>& outputDB);
+                     std::vector<std::complex<float>>& outputCplx);
 
-    // Process windowed I/Q samples → magnitude dB spectrum.
+    void processComplex(const float* inputIQ,
-    // `inputIQ` is interleaved [I0,Q0,I1,Q1,...], fftSize_*2 floats.
+                        std::vector<float>& outputDB,
-    // `outputDB` will be resized to spectrumSize().
+                        std::vector<std::complex<float>>& outputCplx);
-    // Output is FFT-shifted so DC is in the center.
+
+    // Convenience: dB-only (no complex output).
+    void processReal(const float* input, std::vector<float>& outputDB);
     void processComplex(const float* inputIQ, std::vector<float>& outputDB);
 
 private:
     int  fftSize_      = 0;
     bool complexInput_ = false;
 
-    // Real FFT
     float*         realIn_   = nullptr;
     fftwf_complex* realOut_  = nullptr;
     fftwf_plan     realPlan_ = nullptr;
 
-    // Complex FFT
     fftwf_complex* cplxIn_   = nullptr;
     fftwf_complex* cplxOut_  = nullptr;
     fftwf_plan     cplxPlan_ = nullptr;

src/dsp/SpectrumAnalyzer.cpp

@@ -6,7 +6,6 @@
 namespace baudline {
 
 SpectrumAnalyzer::SpectrumAnalyzer() {
-    // Force sizeChanged=true on first configure by setting fftSize to 0.
     settings_.fftSize = 0;
     configure(AnalyzerSettings{});
 }
@@ -35,6 +34,7 @@ void SpectrumAnalyzer::configure(const AnalyzerSettings& settings) {
         int specSz = fft_.spectrumSize();
 
         channelSpectra_.assign(nSpec, std::vector<float>(specSz, -200.0f));
+        channelComplex_.assign(nSpec, std::vector<std::complex<float>>(specSz, {0,0}));
         channelWaterfalls_.assign(nSpec, {});
 
         avgAccum_.assign(nSpec, std::vector<float>(specSz, 0.0f));
@@ -64,7 +64,6 @@ void SpectrumAnalyzer::pushSamples(const float* data, size_t frames) {
         if (accumPos_ >= bufLen) {
             processBlock();
 
-            // Shift by hopSize for overlap
             size_t hopSamples = hopSize_ * inCh;
             size_t keep = bufLen - hopSamples;
             std::memmove(accumBuf_.data(), accumBuf_.data() + hopSamples,
@@ -80,36 +79,33 @@ void SpectrumAnalyzer::processBlock() {
     int nSpec = static_cast<int>(channelSpectra_.size());
     int specSz = fft_.spectrumSize();
 
-    // Compute per-channel spectra.
     std::vector<std::vector<float>> tempDBs(nSpec);
+    std::vector<std::vector<std::complex<float>>> tempCplx(nSpec);
 
     if (settings_.isIQ) {
-        // I/Q: treat the 2 interleaved channels as one complex signal.
         std::vector<float> windowed(N * 2);
         for (int i = 0; i < N; ++i) {
             windowed[2 * i]     = accumBuf_[2 * i]     * window_[i];
             windowed[2 * i + 1] = accumBuf_[2 * i + 1] * window_[i];
         }
-        fft_.processComplex(windowed.data(), tempDBs[0]);
+        fft_.processComplex(windowed.data(), tempDBs[0], tempCplx[0]);
     } else {
-        // Real: deinterleave and FFT each channel independently.
         std::vector<float> chanBuf(N);
         for (int ch = 0; ch < nSpec; ++ch) {
-            // Deinterleave channel `ch` from the accumulation buffer.
             for (int i = 0; i < N; ++i)
                 chanBuf[i] = accumBuf_[i * inCh + ch];
             WindowFunctions::apply(window_, chanBuf.data(), N);
-            fft_.processReal(chanBuf.data(), tempDBs[ch]);
+            fft_.processReal(chanBuf.data(), tempDBs[ch], tempCplx[ch]);
         }
     }
 
-    // Correct for window gain.
+    // Window gain correction.
     float correction = -20.0f * std::log10(windowGain_ > 0 ? windowGain_ : 1.0f);
     for (auto& db : tempDBs)
         for (float& v : db)
             v += correction;
 
-    // Averaging.
+    // Averaging (dB only; complex is not averaged — last block is kept).
     if (settings_.averaging > 1) {
         if (static_cast<int>(avgAccum_[0].size()) != specSz) {
             for (auto& a : avgAccum_) a.assign(specSz, 0.0f);
@@ -131,9 +127,9 @@ void SpectrumAnalyzer::processBlock() {
         }
     }
 
-    // Store results.
     for (int ch = 0; ch < nSpec; ++ch) {
         channelSpectra_[ch] = tempDBs[ch];
+        channelComplex_[ch] = tempCplx[ch];
         channelWaterfalls_[ch].push_back(tempDBs[ch]);
         if (channelWaterfalls_[ch].size() > kWaterfallHistory)
             channelWaterfalls_[ch].pop_front();

src/dsp/SpectrumAnalyzer.h

@@ -3,18 +3,12 @@
 #include "core/Types.h"
 #include "dsp/FFTProcessor.h"
 #include "dsp/WindowFunctions.h"
+#include <complex>
 #include <deque>
 #include <vector>
 
 namespace baudline {
 
-// Manages the DSP pipeline: accumulation with overlap, windowing, FFT,
-// averaging, and waterfall history.
-//
-// Supports three modes:
-//   - Mono real   (numChannels=1, isIQ=false): 1 real FFT  → 1 spectrum
-//   - Multi-ch real (numChannels>1, isIQ=false): N real FFTs → N spectra
-//   - I/Q complex (isIQ=true):                 1 complex FFT → 1 spectrum
 class SpectrumAnalyzer {
 public:
     SpectrumAnalyzer();
@@ -22,37 +16,32 @@ public:
     void configure(const AnalyzerSettings& settings);
     const AnalyzerSettings& settings() const { return settings_; }
 
-    // Feed raw interleaved audio samples.
-    // `frames` = number of sample frames (1 frame = inputChannels() samples).
     void pushSamples(const float* data, size_t frames);
 
-    // Returns true if a new spectrum line is available since last call.
     bool hasNewSpectrum() const { return newSpectrumReady_; }
 
-    // Number of independent spectra (1 for mono/IQ, numChannels for multi-ch).
+    // Number of physical spectra (1 for mono/IQ, numChannels for multi-ch).
     int numSpectra() const { return static_cast<int>(channelSpectra_.size()); }
 
-    // Per-channel spectra (dB magnitudes).
+    // Per-channel dB magnitude spectrum.
     const std::vector<float>& channelSpectrum(int ch) const { return channelSpectra_[ch]; }
-
-    // Convenience: first channel spectrum (backward compat / primary).
     const std::vector<float>& currentSpectrum() const { return channelSpectra_[0]; }
-
-    // All channel spectra.
     const std::vector<std::vector<float>>& allSpectra() const { return channelSpectra_; }
 
-    // Number of output bins (per channel).
+    // Per-channel complex spectrum (for math ops like phase, cross-correlation).
+    const std::vector<std::complex<float>>& channelComplex(int ch) const {
+        return channelComplex_[ch];
+    }
+    const std::vector<std::vector<std::complex<float>>>& allComplex() const {
+        return channelComplex_;
+    }
+
     int spectrumSize() const { return fft_.spectrumSize(); }
 
-    // Peak detection on a given channel.
     std::pair<int, float> findPeak(int ch = 0) const;
-
-    // Get frequency for a given bin index.
     double binToFreq(int bin) const;
-
     void clearHistory();
 
-    // Waterfall history for a given channel (most recent = back).
     const std::deque<std::vector<float>>& waterfallHistory(int ch = 0) const {
         return channelWaterfalls_[ch];
     }
@@ -65,7 +54,6 @@ private:
     std::vector<float>  window_;
     float               windowGain_ = 1.0f;
 
-    // Accumulation buffer (interleaved, length = fftSize * inputChannels)
     std::vector<float>  accumBuf_;
     size_t              accumPos_ = 0;
     size_t              hopSize_  = 0;
@@ -74,8 +62,9 @@ private:
     std::vector<std::vector<float>> avgAccum_;
     int                             avgCount_ = 0;
 
-    // Per-channel output
+    // Per-channel output: magnitude (dB) and complex
     std::vector<std::vector<float>>                  channelSpectra_;
+    std::vector<std::vector<std::complex<float>>>    channelComplex_;
     std::vector<std::deque<std::vector<float>>>      channelWaterfalls_;
     bool                                             newSpectrumReady_ = false;
 };

src/ui/Application.cpp

@@ -178,17 +178,29 @@ void Application::processAudio() {
         analyzer_.pushSamples(audioBuf_.data(), framesRead);
 
         if (analyzer_.hasNewSpectrum()) {
+            computeMathChannels();
+
             int nSpec = analyzer_.numSpectra();
             if (waterfallMultiCh_ && nSpec > 1) {
-                // Multi-channel overlay waterfall.
+                // Multi-channel overlay waterfall: physical + math channels.
-                std::vector<WaterfallChannelInfo> wfChInfo(nSpec);
+                std::vector<std::vector<float>> wfSpectra;
+                std::vector<WaterfallChannelInfo> wfChInfo;
+
                 for (int ch = 0; ch < nSpec; ++ch) {
                     const auto& c = channelColors_[ch % kMaxChannels];
-                    wfChInfo[ch] = {c.x, c.y, c.z,
+                    wfSpectra.push_back(analyzer_.channelSpectrum(ch));
-                                    channelEnabled_[ch % kMaxChannels]};
+                    wfChInfo.push_back({c.x, c.y, c.z,
+                                        channelEnabled_[ch % kMaxChannels]});
                 }
-                waterfall_.pushLineMulti(analyzer_.allSpectra(),
+                for (size_t mi = 0; mi < mathChannels_.size(); ++mi) {
-                                         wfChInfo, minDB_, maxDB_);
+                    if (mathChannels_[mi].enabled && mathChannels_[mi].waterfall &&
+                        mi < mathSpectra_.size()) {
+                        const auto& c = mathChannels_[mi].color;
+                        wfSpectra.push_back(mathSpectra_[mi]);
+                        wfChInfo.push_back({c.x, c.y, c.z, true});
+                    }
+                }
+                waterfall_.pushLineMulti(wfSpectra, wfChInfo, minDB_, maxDB_);
             } else {
                 int wfCh = std::clamp(waterfallChannel_, 0, nSpec - 1);
                 waterfall_.pushLine(analyzer_.channelSpectrum(wfCh),
@@ -424,6 +436,10 @@ void Application::renderControlPanel() {
         }
     }
 
+    // Math channels section (always shown).
+    ImGui::Separator();
+    renderMathPanel();
+
     ImGui::Separator();
 
     // Playback controls
@@ -478,18 +494,39 @@ void Application::renderSpectrumPanel() {
     specSizeX_ = availW;
     specSizeY_ = specH;
 
-    // Build per-channel styles and pass all spectra.
+    // Build per-channel styles and combine physical + math spectra.
-    int nCh = analyzer_.numSpectra();
+    int nPhys = analyzer_.numSpectra();
-    std::vector<ChannelStyle> styles(nCh);
+    int nMath = static_cast<int>(mathSpectra_.size());
-    for (int ch = 0; ch < nCh; ++ch) {
+    int nTotal = nPhys + nMath;
+
+    std::vector<std::vector<float>> allSpectra;
+    std::vector<ChannelStyle> styles;
+    allSpectra.reserve(nTotal);
+    styles.reserve(nTotal);
+
+    // Physical channels.
+    for (int ch = 0; ch < nPhys; ++ch) {
+        allSpectra.push_back(analyzer_.channelSpectrum(ch));
         const auto& c = channelColors_[ch % kMaxChannels];
         uint8_t r = static_cast<uint8_t>(c.x * 255);
         uint8_t g = static_cast<uint8_t>(c.y * 255);
         uint8_t b = static_cast<uint8_t>(c.z * 255);
-        styles[ch].lineColor = IM_COL32(r, g, b, 220);
+        styles.push_back({IM_COL32(r, g, b, 220), IM_COL32(r, g, b, 35)});
-        styles[ch].fillColor = IM_COL32(r, g, b, 35);
     }
-    specDisplay_.draw(analyzer_.allSpectra(), styles, minDB_, maxDB_,
+
+    // Math channels.
+    for (int mi = 0; mi < nMath; ++mi) {
+        if (mi < static_cast<int>(mathChannels_.size()) && mathChannels_[mi].enabled) {
+            allSpectra.push_back(mathSpectra_[mi]);
+            const auto& c = mathChannels_[mi].color;
+            uint8_t r = static_cast<uint8_t>(c.x * 255);
+            uint8_t g = static_cast<uint8_t>(c.y * 255);
+            uint8_t b = static_cast<uint8_t>(c.z * 255);
+            styles.push_back({IM_COL32(r, g, b, 220), IM_COL32(r, g, b, 35)});
+        }
+    }
+
+    specDisplay_.draw(allSpectra, styles, minDB_, maxDB_,
                       settings_.sampleRate, settings_.isIQ, freqScale_,
                       specPosX_, specPosY_, specSizeX_, specSizeY_);
 
@@ -715,4 +752,169 @@ void Application::updateAnalyzerSettings() {
     }
 }
 
+// ── Math channels ────────────────────────────────────────────────────────────
+
+void Application::computeMathChannels() {
+    int nPhys = analyzer_.numSpectra();
+    int specSz = analyzer_.spectrumSize();
+    mathSpectra_.resize(mathChannels_.size());
+
+    for (size_t mi = 0; mi < mathChannels_.size(); ++mi) {
+        const auto& mc = mathChannels_[mi];
+        auto& out = mathSpectra_[mi];
+        out.resize(specSz);
+
+        if (!mc.enabled) {
+            std::fill(out.begin(), out.end(), -200.0f);
+            continue;
+        }
+
+        int sx = std::clamp(mc.sourceX, 0, nPhys - 1);
+        int sy = std::clamp(mc.sourceY, 0, nPhys - 1);
+        const auto& xDB = analyzer_.channelSpectrum(sx);
+        const auto& yDB = analyzer_.channelSpectrum(sy);
+        const auto& xC  = analyzer_.channelComplex(sx);
+        const auto& yC  = analyzer_.channelComplex(sy);
+
+        for (int i = 0; i < specSz; ++i) {
+            float val = -200.0f;
+            switch (mc.op) {
+                // ── Unary ──
+                case MathOp::Negate:
+                    val = -xDB[i];
+                    break;
+                case MathOp::Absolute:
+                    val = std::abs(xDB[i]);
+                    break;
+                case MathOp::Square:
+                    val = 2.0f * xDB[i];
+                    break;
+                case MathOp::Cube:
+                    val = 3.0f * xDB[i];
+                    break;
+                case MathOp::Sqrt:
+                    val = 0.5f * xDB[i];
+                    break;
+                case MathOp::Log: {
+                    // log10 of linear magnitude, back to dB-like scale.
+                    float lin = std::pow(10.0f, xDB[i] / 10.0f);
+                    float l = std::log10(lin + 1e-30f);
+                    val = 10.0f * l;  // keep in dB-like range
+                    break;
+                }
+                // ── Binary ──
+                case MathOp::Add: {
+                    float lx = std::pow(10.0f, xDB[i] / 10.0f);
+                    float ly = std::pow(10.0f, yDB[i] / 10.0f);
+                    float s = lx + ly;
+                    val = (s > 1e-20f) ? 10.0f * std::log10(s) : -200.0f;
+                    break;
+                }
+                case MathOp::Subtract: {
+                    float lx = std::pow(10.0f, xDB[i] / 10.0f);
+                    float ly = std::pow(10.0f, yDB[i] / 10.0f);
+                    float d = std::abs(lx - ly);
+                    val = (d > 1e-20f) ? 10.0f * std::log10(d) : -200.0f;
+                    break;
+                }
+                case MathOp::Multiply:
+                    val = xDB[i] + yDB[i];
+                    break;
+                case MathOp::Phase: {
+                    if (i < static_cast<int>(xC.size()) &&
+                        i < static_cast<int>(yC.size())) {
+                        auto cross = xC[i] * std::conj(yC[i]);
+                        float deg = std::atan2(cross.imag(), cross.real())
+                                    * (180.0f / 3.14159265f);
+                        // Map [-180, 180] degrees into the dB display range
+                        // so it's visible on the plot.
+                        val = deg;
+                    }
+                    break;
+                }
+                case MathOp::CrossCorr: {
+                    if (i < static_cast<int>(xC.size()) &&
+                        i < static_cast<int>(yC.size())) {
+                        auto cross = xC[i] * std::conj(yC[i]);
+                        float mag2 = std::norm(cross);
+                        val = (mag2 > 1e-20f) ? 10.0f * std::log10(mag2) : -200.0f;
+                    }
+                    break;
+                }
+                default: break;
+            }
+            out[i] = val;
+        }
+    }
+}
+
+void Application::renderMathPanel() {
+    ImGui::Text("Channel Math");
+    ImGui::Separator();
+
+    int nPhys = analyzer_.numSpectra();
+
+    // Build source channel name list.
+    static const char* chNames[] = {
+        "Ch 0 (L)", "Ch 1 (R)", "Ch 2", "Ch 3", "Ch 4", "Ch 5", "Ch 6", "Ch 7"
+    };
+
+    // List existing math channels.
+    int toRemove = -1;
+    for (int mi = 0; mi < static_cast<int>(mathChannels_.size()); ++mi) {
+        auto& mc = mathChannels_[mi];
+        ImGui::PushID(1000 + mi);
+
+        ImGui::Checkbox("##en", &mc.enabled);
+        ImGui::SameLine();
+        ImGui::ColorEdit3("##col", &mc.color.x, ImGuiColorEditFlags_NoInputs);
+        ImGui::SameLine();
+
+        // Operation combo.
+        if (ImGui::BeginCombo("##op", mathOpName(mc.op), ImGuiComboFlags_NoPreview)) {
+            for (int o = 0; o < static_cast<int>(MathOp::Count); ++o) {
+                auto op = static_cast<MathOp>(o);
+                if (ImGui::Selectable(mathOpName(op), mc.op == op))
+                    mc.op = op;
+            }
+            ImGui::EndCombo();
+        }
+        ImGui::SameLine();
+        ImGui::Text("%s", mathOpName(mc.op));
+
+        // Source X.
+        ImGui::SetNextItemWidth(80);
+        ImGui::Combo("X", &mc.sourceX, chNames, std::min(nPhys, kMaxChannels));
+
+        // Source Y (only for binary ops).
+        if (mathOpIsBinary(mc.op)) {
+            ImGui::SameLine();
+            ImGui::SetNextItemWidth(80);
+            ImGui::Combo("Y", &mc.sourceY, chNames, std::min(nPhys, kMaxChannels));
+        }
+
+        ImGui::SameLine();
+        ImGui::Checkbox("WF", &mc.waterfall);
+        if (ImGui::IsItemHovered())
+            ImGui::SetTooltip("Show on waterfall");
+        ImGui::SameLine();
+        if (ImGui::SmallButton("X##del"))
+            toRemove = mi;
+
+        ImGui::PopID();
+    }
+
+    if (toRemove >= 0)
+        mathChannels_.erase(mathChannels_.begin() + toRemove);
+
+    if (ImGui::Button("+ Add Math Channel")) {
+        MathChannel mc;
+        mc.op = MathOp::Subtract;
+        mc.sourceX = 0;
+        mc.sourceY = std::min(1, nPhys - 1);
+        mc.color = {1.0f, 1.0f, 0.5f, 1.0f};
+        mathChannels_.push_back(mc);
+    }
+}
+
 } // namespace baudline

src/ui/Application.h

@@ -10,12 +10,62 @@
 #include "ui/Cursors.h"
 
 #include <SDL.h>
+#include <complex>
 #include <memory>
 #include <string>
 #include <vector>
 
 namespace baudline {
 
+// ── Channel math operations ──────────────────────────────────────────────────
+
+enum class MathOp {
+    // Unary (on channel X)
+    Negate,       // -x  (negate dB)
+    Absolute,     // |x| (absolute value of dB)
+    Square,       // x^2 in linear  → 2*x_dB
+    Cube,         // x^3 in linear  → 3*x_dB
+    Sqrt,         // sqrt in linear → 0.5*x_dB
+    Log,          // 10*log10(10^(x_dB/10) + 1)  — compressed scale
+    // Binary (on channels X and Y)
+    Add,          // linear(x) + linear(y)  → dB
+    Subtract,     // |linear(x) - linear(y)| → dB
+    Multiply,     // x_dB + y_dB  (multiply in linear = add in dB)
+    Phase,        // angle(X_cplx * conj(Y_cplx)) in degrees
+    CrossCorr,    // |X_cplx * conj(Y_cplx)| → dB
+    Count
+};
+
+inline bool mathOpIsBinary(MathOp op) {
+    return op >= MathOp::Add;
+}
+
+inline const char* mathOpName(MathOp op) {
+    switch (op) {
+        case MathOp::Negate:    return "-x";
+        case MathOp::Absolute:  return "|x|";
+        case MathOp::Square:    return "x^2";
+        case MathOp::Cube:      return "x^3";
+        case MathOp::Sqrt:      return "sqrt(x)";
+        case MathOp::Log:       return "log(x)";
+        case MathOp::Add:       return "x + y";
+        case MathOp::Subtract:  return "x - y";
+        case MathOp::Multiply:  return "x * y";
+        case MathOp::Phase:     return "phase(x,y)";
+        case MathOp::CrossCorr: return "xcorr(x,y)";
+        default:                return "?";
+    }
+}
+
+struct MathChannel {
+    MathOp op        = MathOp::Subtract;
+    int    sourceX   = 0;
+    int    sourceY   = 1;
+    ImVec4 color     = {1.0f, 1.0f, 1.0f, 1.0f};
+    bool   enabled   = true;
+    bool   waterfall = false;  // include on waterfall overlay
+};
+
 class Application {
 public:
     Application();
@@ -36,6 +86,8 @@ private:
     void openPortAudio();
     void openFile(const std::string& path, InputFormat format, double sampleRate);
     void updateAnalyzerSettings();
+    void computeMathChannels();
+    void renderMathPanel();
 
     // SDL / GL / ImGui
     SDL_Window*   window_    = nullptr;
@@ -103,6 +155,10 @@ private:
     bool waterfallMultiCh_ = true; // true = multi-channel overlay mode
     bool channelEnabled_[kMaxChannels] = {true,true,true,true,true,true,true,true};
 
+    // Math channels
+    std::vector<MathChannel>          mathChannels_;
+    std::vector<std::vector<float>>   mathSpectra_;  // computed each frame
+
     // Spectrum panel geometry (stored for cursor interaction)
     float specPosX_ = 0, specPosY_ = 0, specSizeX_ = 0, specSizeY_ = 0;
 };