高频信号处理基础及编程实例-CSDN博客

本文链接：https://blog.csdn.net/chenby186119/article/details/148069031

高频信号处理基础及编程实例

1. 高频信号处理基础理论

1.1 高频信号的特性与挑战

高频信号通常指频率在数百MHz到数十GHz范围内的信号。在这些频率下，信号处理面临独特的挑战：

高频信号的主要特性：

波长较短，易受传输路径影响
存在明显的传播延迟效应
容易受噪声和干扰影响
频率分量丰富，信息容量大
存在寄生效应和阻抗匹配问题

处理挑战：

采样率要求高（奈奎斯特采样定理）
信号衰减与失真严重
计算复杂度高
实时处理要求苛刻
硬件设计复杂（时钟同步、布线等）

1.2 信号采样与量化

采样原理：

奈奎斯特采样定理：采样频率必须至少是信号最高频率的两倍
对于带宽为B的高频信号，理论上最低采样率为2B
实际应用中通常采用更高的采样率（3-4倍最高频率）防止混叠

量化过程：

ADC（模数转换器）将连续幅度信号转换为离散数字值
量化位数决定了动态范围和信噪比
高频系统通常需要高速ADC（GSPS级别采样率）
量化误差会引入噪声和失真

混叠效应：

当采样率不足时，高频分量会被"折叠"到低频区域
解决方法：抗混叠滤波器、带通采样、欠采样技术

1.3 高频数字滤波基础

FIR滤波器：

有限脉冲响应滤波器，具有线性相位特性
结构简单，易于并行实现
计算公式：y[n] = Σ(h[k] × x[n-k])，其中h[k]是滤波器系数

IIR滤波器：

无限脉冲响应滤波器，计算效率更高
可能存在稳定性问题
计算公式：y[n] = Σ(b[k] × x[n-k]) - Σ(a[k] × y[n-k])

高频滤波器设计考量：

过渡带宽度与滤波器阶数的权衡
组延迟与相位响应
计算复杂度与实时处理能力
频率选择性与时域特性的平衡

1.4 频域分析技术

快速傅里叶变换(FFT)：

将时域信号转换到频域的高效算法
复杂度O(N log N)，远优于直接DFT的O(N²)
基本蝶形运算与基2-FFT算法

频谱分析方法：

功率谱密度(PSD)估计
窗函数选择（汉宁窗、黑曼窗、凯撒窗等）
短时傅里叶变换(STFT)处理非平稳信号
平均周期图法降低噪声

高级频域分析：

小波变换提供时频局部化分析
希尔伯特变换计算瞬时频率
双谱分析检测非线性特性
联合时频分析（Wigner-Ville分布等）

1.5 现代高频信号处理算法

自适应滤波：

LMS (Least Mean Square) 算法
RLS (Recursive Least Square) 算法
卡尔曼滤波器跟踪时变信号

谱估计方法：

参数化方法：AR、MA、ARMA模型
高分辨率谱估计：MUSIC、ESPRIT算法
最大熵方法(MEM)

先进调制解调技术：

OFDM（正交频分复用）
扩频技术（直接序列扩频、跳频）
自适应调制与编码

数字下变频(DDC)与数字上变频(DUC)：

I/Q解调与正交处理
CIC (Cascaded Integrator-Comb) 滤波器
多相滤波器与抽取插值技术

2. C/C++实现高频信号处理

2.1 基础信号生成与操作

2.1.1 生成高频正弦信号

cpp

#include <iostream>
#include <vector>
#include <cmath>
#include <complex>

// 生成高频复数正弦信号
std::vector<std::complex<double>> generateComplexSine(
    double frequency,    // 信号频率 (Hz)
    double sampleRate,   // 采样率 (Hz)
    double duration,     // 信号持续时间 (秒)
    double amplitude = 1.0, // 振幅
    double phaseOffset = 0.0) // 相位偏移(弧度)
{
    int numSamples = static_cast<int>(sampleRate * duration);
    std::vector<std::complex<double>> signal(numSamples);
    
    for (int i = 0; i < numSamples; i++) {
        double time = i / sampleRate;
        double angle = 2.0 * M_PI * frequency * time + phaseOffset;
        
        // 使用欧拉公式 e^(jθ) = cos(θ) + j·sin(θ)
        signal[i] = std::complex<double>(
            amplitude * cos(angle), 
            amplitude * sin(angle)
        );
    }
    
    return signal;
}

// 生成多频率正弦信号
std::vector<double> generateMultiFrequencySine(
    const std::vector<double>& frequencies, // 频率数组 (Hz)
    const std::vector<double>& amplitudes,  // 对应振幅数组
    double sampleRate,                     // 采样率 (Hz)
    double duration)                       // 信号持续时间 (秒)
{
    if (frequencies.size() != amplitudes.size()) {
        throw std::invalid_argument("频率和振幅数组长度必须相同");
    }
    
    int numSamples = static_cast<int>(sampleRate * duration);
    std::vector<double> signal(numSamples, 0.0);
    
    for (size_t freqIndex = 0; freqIndex < frequencies.size(); freqIndex++) {
        double frequency = frequencies[freqIndex];
        double amplitude = amplitudes[freqIndex];
        
        for (int i = 0; i < numSamples; i++) {
            double time = i / sampleRate;
            signal[i] += amplitude * sin(2.0 * M_PI * frequency * time);
        }
    }
    
    return signal;
}

// 示例：生成并打印复数正弦信号
void testSignalGeneration() {
    double sampleRate = 1e6;  // 1 MHz采样率
    double frequency = 1e5;   // 100 kHz信号
    double duration = 1e-5;   // 10微秒
    
    auto complexSignal = generateComplexSine(frequency, sampleRate, duration);
    
    std::cout << "复数信号示例 (前10个样本):" << std::endl;
    for (int i = 0; i < std::min(10, static_cast<int>(complexSignal.size())); i++) {
        std::cout << "Sample[" << i << "] = " 
                  << complexSignal[i].real() << " + " 
                  << complexSignal[i].imag() << "i" << std::endl;
    }
}

2.1.2 基本信号操作

cpp

// 对复数信号加入高斯白噪声
std::vector<std::complex<double>> addComplexNoise(
    const std::vector<std::complex<double>>& signal,
    double snr_dB)  // 信噪比(dB)
{
    std::vector<std::complex<double>> noisySignal = signal;
    
    // 计算信号功率
    double signalPower = 0.0;
    for (const auto& sample : signal) {
        signalPower += std::norm(sample);  // |z|^2 = re^2 + im^2
    }
    signalPower /= signal.size();
    
    // 根据SNR计算噪声功率
    double noisePower = signalPower / pow(10.0, snr_dB / 10.0);
    double noiseStdDev = sqrt(noisePower);
    
    // 使用C++11随机数生成器
    std::random_device rd;
    std::mt19937 gen(rd());
    std::normal_distribution<double> dist(0.0, noiseStdDev);
    
    // 为每个样本添加噪声
    for (auto& sample : noisySignal) {
        double noiseReal = dist(gen);
        double noiseImag = dist(gen);
        sample += std::complex<double>(noiseReal, noiseImag);
    }
    
    return noisySignal;
}

// 信号混频 (乘以复指数，实现频率搬移)
std::vector<std::complex<double>> frequencyShift(
    const std::vector<std::complex<double>>& signal,
    double shiftFrequency,  // 频移量 (Hz)
    double sampleRate)      // 采样率 (Hz)
{
    std::vector<std::complex<double>> shiftedSignal(signal.size());
    
    for (size_t i = 0; i < signal.size(); i++) {
        double time = static_cast<double>(i) / sampleRate;
        double angle = 2.0 * M_PI * shiftFrequency * time;
        
        // 复数乘法：信号样本乘以复指数e^(-j·2π·f·t)
        std::complex<double> shifter(cos(angle), -sin(angle));
        shiftedSignal[i] = signal[i] * shifter;
    }
    
    return shiftedSignal;
}

// 信号功率计算
double calculatePower(const std::vector<std::complex<double>>& signal) {
    double power = 0.0;
    for (const auto& sample : signal) {
        power += std::norm(sample);  // |z|^2
    }
    return power / signal.size();
}

// 计算复数信号的包络
std::vector<double> calculateEnvelope(const std::vector<std::complex<double>>& signal) {
    std::vector<double> envelope(signal.size());
    
    for (size_t i = 0; i < signal.size(); i++) {
        envelope[i] = std::abs(signal[i]);  // sqrt(real^2 + imag^2)
    }
    
    return envelope;
}

2.2 快速傅里叶变换(FFT)实现

2.2.1 基于C++的FFT实现

cpp

#include <complex>
#include <vector>
#include <cmath>
#include <stdexcept>

// 递归实现的FFT算法
void fft(std::vector<std::complex<double>>& x) {
    const size_t N = x.size();
    if (N <= 1) return;
    
    // 检查N是否为2的幂
    if ((N & (N - 1)) != 0) {
        throw std::invalid_argument("FFT长度必须是2的幂");
    }
    
    // 分割为奇偶序列
    std::vector<std::complex<double>> even(N/2), odd(N/2);
    for (size_t i = 0; i < N/2; i++) {
        even[i] = x[2*i];
        odd[i] = x[2*i + 1];
    }
    
    // 递归计算子问题
    fft(even);
    fft(odd);
    
    // 合并结果
    for (size_t k = 0; k < N/2; k++) {
        std::complex<double> t = std::polar(1.0, -2 * M_PI * k / N) * odd[k];
        x[k] = even[k] + t;
        x[k + N/2] = even[k] - t;
    }
}

// 计算逆FFT
void ifft(std::vector<std::complex<double>>& x) {
    // 对每个元素取共轭
    for (auto& val : x) {
        val = std::conj(val);
    }
    
    // 执行FFT
    fft(x);
    
    // 再次取共轭并除以N
    const size_t N = x.size();
    for (auto& val : x) {
        val = std::conj(val) / static_cast<double>(N);
    }
}

// 迭代实现的FFT (比递归实现更高效)
void fft_iterative(std::vector<std::complex<double>>& x) {
    const size_t N = x.size();
    
    // 检查N是否为2的幂
    if ((N & (N - 1)) != 0) {
        throw std::invalid_argument("FFT长度必须是2的幂");
    }
    
    // 位反转排序
    size_t j = 0;
    for (size_t i = 0; i < N - 1; i++) {
        if (i < j) {
            std::swap(x[i], x[j]);
        }
        size_t k = N / 2;
        while (k <= j) {
            j -= k;
            k /= 2;
        }
        j += k;
    }
    
    // 蝶形运算
    for (size_t len = 2; len <= N; len *= 2) {
        double angle = -2 * M_PI / len;
        std::complex<double> w_m(cos(angle), sin(angle));
        
        for (size_t i = 0; i < N; i += len) {
            std::complex<double> w(1.0, 0.0);
            for (size_t j = 0; j < len/2; j++) {
                std::complex<double> u = x[i + j];
                std::complex<double> t = w * x[i + j + len/2];
                
                x[i + j] = u + t;
                x[i + j + len/2] = u - t;
                
                w *= w_m;
            }
        }
    }
}

2.2.2 功率谱密度计算

cpp

// 使用FFT计算功率谱密度
std::vector<double> calculatePSD(
    const std::vector<std::complex<double>>& signal,
    double sampleRate)
{
    size_t N = signal.size();
    
    // 为了使用FFT，确保N是2的幂
    size_t nextPow2 = 1;
    while (nextPow2 < N) {
        nextPow2 *= 2;
    }
    
    // 复制信号并应用窗函数 (这里使用汉宁窗)
    std::vector<std::complex<double>> windowed(nextPow2, std::complex<double>(0.0, 0.0));
    for (size_t i = 0; i < N; i++) {
        // 汉宁窗: 0.5 * (1 - cos(2π*i/(N-1)))
        double window = 0.5 * (1.0 - cos(2.0 * M_PI * i / (N - 1)));
        windowed[i] = signal[i] * window;
    }
    
    // 执行FFT
    fft_iterative(windowed);
    
    // 计算单边功率谱
    std::vector<double> psd(nextPow2 / 2 + 1);
    
    // 直流分量
    psd[0] = std::norm(windowed[0]) / (sampleRate * N);
    
    // 其他频率分量
    for (size_t i = 1; i < nextPow2 / 2; i++) {
        psd[i] = (std::norm(windowed[i]) + std::norm(windowed[nextPow2 - i])) 
               / (sampleRate * N);
    }
    
    // 奈奎斯特频率
    psd[nextPow2 / 2] = std::norm(windowed[nextPow2 / 2]) / (sampleRate * N);
    
    return psd;
}

// 示例：计算信号的PSD
void testPSDCalculation() {
    double sampleRate = 1e6;   // 1 MHz采样率
    double signalFreq = 1e5;   // 100 kHz信号
    double duration = 1e-3;    // 1毫秒
    
    // 生成信号
    auto signal = generateComplexSine(signalFreq, sampleRate, duration);
    
    // 添加噪声
    auto noisySignal = addComplexNoise(signal, 10.0);  // 10dB SNR
    
    // 计算PSD
    auto psd = calculatePSD(noisySignal, sampleRate);
    
    // 输出部分结果
    std::cout << "PSD结果 (部分数据):" << std::endl;
    for (size_t i = 0; i < std::min(size_t(10), psd.size()); i++) {
        double freq = i * sampleRate / (2.0 * (psd.size() - 1));
        std::cout << "Freq: " << freq << " Hz, PSD: " << 10 * log10(psd[i]) << " dB/Hz" << std::endl;
    }
}

2.3 数字滤波器实现

2.3.1 FIR滤波器

cpp

#include <vector>
#include <complex>
#include <algorithm>

class FIRFilter {
private:
    std::vector<double> coeffs;        // 滤波器系数
    std::vector<std::complex<double>> buffer; // 状态缓冲区
    size_t index;                      // 当前索引位置

public:
    // 构造函数，传入滤波器系数
    FIRFilter(const std::vector<double>& coefficients) : 
        coeffs(coefficients),
        buffer(coefficients.size(), std::complex<double>(0.0, 0.0)),
        index(0) {
    }
    
    // 处理单个样本
    std::complex<double> processSample(std::complex<double> sample) {
        // 保存当前样本到缓冲区
        buffer[index] = sample;
        
        // 计算输出
        std::complex<double> result(0.0, 0.0);
        size_t bufferIndex = index;
        
        for (const auto& coeff : coeffs) {
            result += buffer[bufferIndex] * coeff;
            bufferIndex = (bufferIndex == 0) ? buffer.size() - 1 : bufferIndex - 1;
        }
        
        // 更新索引
        index = (index + 1) % buffer.size();
        
        return result;
    }
    
    // 处理整个信号
    std::vector<std::complex<double>> process(
        const std::vector<std::complex<double>>& input) {
        std::vector<std::complex<double>> output;
        output.reserve(input.size());
        
        // 重置状态
        std::fill(buffer.begin(), buffer.end(), std::complex<double>(0.0, 0.0));
        index = 0;
        
        // 处理每个样本
        for (const auto& sample : input) {
            output.push_back(processSample(sample));
        }
        
        return output;
    }
    
    // 重置滤波器状态
    void reset() {
        std::fill(buffer.begin(), buffer.end(), std::complex<double>(0.0, 0.0));
        index = 0;
    }
    
    // 获取滤波器系数
    const std::vector<double>& getCoefficients() const {
        return coeffs;
    }
};

// 生成低通FIR滤波器系数 (使用窗函数法)
std::vector<double> designLowPassFIR(
    double cutoffFreq,    // 截止频率 (Hz)
    double sampleRate,    // 采样率 (Hz)
    size_t numTaps,       // 滤波器阶数
    const std::string& windowType = "hamming") // 窗类型
{
    std::vector<double> coeffs(numTaps);
    
    // 计算归一化截止频率
    double omega_c = 2.0 * M_PI * cutoffFreq / sampleRate;
    
    // 中心点位置
    int center = numTaps / 2;
    
    // 计算理想低通滤波器的脉冲响应
    for (size_t i = 0; i < numTaps; i++) {
        if (i == center) {
            coeffs[i] = omega_c / M_PI;
        } else {
            double n = i - center;
            coeffs[i] = sin(omega_c * n) / (M_PI * n);
        }
    }
    
    // 应用窗函数
    if (windowType == "hamming") {
        // 汉明窗
        for (size_t i = 0; i < numTaps; i++) {
            double window = 0.54 - 0.46 * cos(2.0 * M_PI * i / (numTaps - 1));
            coeffs[i] *= window;
        }
    } else if (windowType == "blackman") {
        // 布莱克曼窗
        for (size_t i = 0; i < numTaps; i++) {
            double window = 0.42 - 0.5 * cos(2.0 * M_PI * i / (numTaps - 1)) 
                         + 0.08 * cos(4.0 * M_PI * i / (numTaps - 1));
            coeffs[i] *= window;
        }
    } else if (windowType == "hann") {
        // 汉宁窗
        for (size_t i = 0; i < numTaps; i++) {
            double window = 0.5 * (1.0 - cos(2.0 * M_PI * i / (numTaps - 1)));
            coeffs[i] *= window;
        }
    }
    
    return coeffs;
}

2.3.2 IIR滤波器

cpp

class IIRFilter {
private:
    std::vector<double> a_coeffs;  // 反馈系数
    std::vector<double> b_coeffs;  // 前馈系数
    std::vector<std::complex<double>> x_buffer; // 输入缓冲区
    std::vector<std::complex<double>> y_buffer; // 输出缓冲区

public:
    // 构造函数，传入前馈和反馈系数
    IIRFilter(const std::vector<double>& b, const std::vector<double>& a) :
        b_coeffs(b),
        a_coeffs(a),
        x_buffer(b.size(), std::complex<double>(0.0, 0.0)),
        y_buffer(a.size(), std::complex<double>(0.0, 0.0)) {
        
        // 归一化a_0系数
        if (a[0] != 1.0) {
            double a0 = a[0];
            for (auto& coeff : a_coeffs) {
                coeff /= a0;
            }
            for (auto& coeff : b_coeffs) {
                coeff /= a0;
            }
        }
    }
    
    // 处理单个样本
    std::complex<double> processSample(std::complex<double> sample) {
        // 更新输入缓冲区
        std::rotate(x_buffer.rbegin(), x_buffer.rbegin() + 1, x_buffer.rend());
        x_buffer[0] = sample;
        
        // 计算输出
        std::complex<double> output = 0.0;
        
        // 前馈部分 (b系数 * 输入)
        for (size_t i = 0; i < b_coeffs.size(); i++) {
            output += b_coeffs[i] * x_buffer[i];
        }
        
        // 反馈部分 (a系数 * 先前输出)
        for (size_t i = 1; i < a_coeffs.size(); i++) {
            output -= a_coeffs[i] * y_buffer[i-1];
        }
        
        // 更新输出缓冲区
        std::rotate(y_buffer.rbegin(), y_buffer.rbegin() + 1, y_buffer.rend());
        y_buffer[0] = output;
        
        return output;
    }
    
    // 处理整个信号
    std::vector<std::complex<double>> process(
        const std::vector<std::complex<double>>& input) {
        std::vector<std::complex<double>> output;
        output.reserve(input.size());
        
        // 重置状态
        std::fill(x_buffer.begin(), x_buffer.end(), std::complex<double>(0.0, 0.0));
        std::fill(y_buffer.begin(), y_buffer.end(), std::complex<double>(0.0, 0.0));
        
        // 处理每个样本
        for (const auto& sample : input) {
            output.push_back(processSample(sample));
        }
        
        return output;
    }
    
    // 重置滤波器状态
    void reset() {
        std::fill(x_buffer.begin(), x_buffer.end(), std::complex<double>(0.0, 0.0));
        std::fill(y_buffer.begin(), y_buffer.end(), std::complex<double>(0.0, 0.0));
    }
};

// 设计二阶带通滤波器系数
std::pair<std::vector<double>, std::vector<double>> designBandpassIIR(
    double centerFreq,   // 中心频率 (Hz)
    double bandwidth,    // 带宽 (Hz)
    double sampleRate,   // 采样率 (Hz)
    double gain = 1.0)   // 增益
{
    // 计算归一化中心频率和带宽
    double omega0 = 2.0 * M_PI * centerFreq / sampleRate;
    double alpha = sin(omega0) * sinh(log(2.0) / 2.0 * bandwidth * M_PI / sampleRate);
    
    // 计算系数
    std::vector<double> b(3);
    std::vector<double> a(3);
    
    b[0] = alpha * gain;
    b[1] = 0.0;
    b[2] = -alpha * gain;
    
    a[0] = 1.0 + alpha;
    a[1] = -2.0 * cos(omega0);
    a[2] = 1.0 - alpha;
    
    return {b, a};
}

2.4 高频信号分析应用

2.4.1 频谱分析器

cpp

#include <iostream>
#include <vector>
#include <complex>
#include <cmath>
#include <algorithm>
#include <string>
#include <fstream>

class SpectrumAnalyzer {
private:
    size_t fftSize;
    double sampleRate;
    std::vector<double> windowCoeffs;
    
    // 计算窗函数
    void createWindow(const std::string& windowType) {
        windowCoeffs.resize(fftSize);
        
        if (windowType == "hamming") {
            for (size_t i = 0; i < fftSize; i++) {
                windowCoeffs[i] = 0.54 - 0.46 * cos(2.0 * M_PI * i / (fftSize - 1));
            }
        } else if (windowType == "hann") {
            for (size_t i = 0; i < fftSize; i++) {
                windowCoeffs[i] = 0.5 * (1.0 - cos(2.0 * M_PI * i / (fftSize - 1)));
            }
        } else if (windowType == "blackman") {
            for (size_t i = 0; i < fftSize; i++) {
                windowCoeffs[i] = 0.42 - 0.5 * cos(2.0 * M_PI * i / (fftSize - 1)) 
                                + 0.08 * cos(4.0 * M_PI * i / (fftSize - 1));
            }
        } else {
            // 默认使用矩形窗 (无窗)
            std::fill(windowCoeffs.begin(), windowCoeffs.end(), 1.0);
        }
    }
    
public:
    SpectrumAnalyzer(size_t fftSize, double sampleRate, 
                    const std::string& windowType = "hamming") :
        fftSize(fftSize),
        sampleRate(sampleRate) {
        
        // 确保FFT大小是2的幂
        if ((fftSize & (fftSize - 1)) != 0) {
            size_t newSize = 1;
            while (newSize < fftSize) {
                newSize *= 2;
            }
            this->fftSize = newSize;
            std::cout << "调整FFT大小为2的幂: " << newSize << std::endl;
        }
        
        // 创建窗函数
        createWindow(windowType);
    }
    
    // 分析信号并返回功率谱
    std::vector<double> analyze(const std::vector<std::complex<double>>& signal, 
                               bool dBScale = true) {
        // 准备输入数据
        std::vector<std::complex<double>> windowed(fftSize, std::complex<double>(0, 0));
        
        // 应用窗函数 (处理信号长度可能小于FFT大小的情况)
        size_t dataLength = std::min(signal.size(), fftSize);
        for (size_t i = 0; i < dataLength; i++) {
            windowed[i] = signal[i] * windowCoeffs[i];
        }
        
        // 执行FFT
        std::vector<std::complex<double>> fftResult = windowed;
        fft_iterative(fftResult);
        
        // 计算功率谱
        std::vector<double> powerSpectrum(fftSize / 2 + 1);
        
        // 直流分量
        powerSpectrum[0] = std::norm(fftResult[0]) / fftSize;
        
        // 其他频率分量
        for (size_t i = 1; i < fftSize / 2; i++) {
            powerSpectrum[i] = (std::norm(fftResult[i]) + std::norm(fftResult[fftSize - i])) 
                             / fftSize;
        }
        
        // 奈奎斯特频率
        powerSpectrum[fftSize / 2] = std::norm(fftResult[fftSize / 2]) / fftSize;
        
        // 转换为dB刻度
        if (dBScale) {
            for (auto& power : powerSpectrum) {
                // 防止对0或负数取对数
                if (power <= 0) {
                    power = -100; // 设置一个较低的dB值
                } else {
                    power = 10.0 * log10(power);
                }
            }
        }
        
        return powerSpectrum;
    }
    
    // 获取频率轴
    std::vector<double> getFrequencyAxis() {
        std::vector<double> freqAxis(fftSize / 2 + 1);
        
        for (size_t i = 0; i <= fftSize / 2; i++) {
            freqAxis[i] = i * sampleRate / fftSize;
        }
        
        return freqAxis;
    }
    
    // 保存谱图到CSV文件
    void saveSpectrumToCSV(const std::vector<double>& spectrum, 
                          const std::string& filename) {
        std::ofstream file(filename);
        if (!file.is_open()) {
            std::cerr << "无法打开文件: " << filename << std::endl;
            return;
        }
        
        std::vector<double> freqAxis = getFrequencyAxis();
        
        file << "Frequency(Hz),Power(dB)\n";
        for (size_t i = 0; i < spectrum.size(); i++) {
            file << freqAxis[i] << "," << spectrum[i] << "\n";
        }
        
        file.close();
    }
};

// 示例：使用频谱分析器分析信号
void testSpectrumAnalyzer() {
    // 参数设置
    double sampleRate = 10e6;    // 10 MHz采样率
    size_t fftSize = 4096;       // 4K FFT点数
    
    // 创建分析器
    SpectrumAnalyzer analyzer(fftSize, sampleRate, "hamming");
    
    // 生成测试信号 (包含多个频率分量)
    std::vector<double> freqs = {1e6, 2.5e6, 3.3e6};  // 1MHz, 2.5MHz, 3.3MHz
    std::vector<double> amps = {1.0, 0.5, 0.25};      // 不同幅度
    
    // 创建复数信号
    std::vector<std::complex<double>> signal;
    signal.reserve(fftSize);
    
    for (size_t i = 0; i < fftSize; i++) {
        double time = i / sampleRate;
        std::complex<double> sample(0.0, 0.0);
        
        for (size_t j = 0; j < freqs.size(); j++) {
            double angle = 2.0 * M_PI * freqs[j] * time;
            sample += std::complex<double>(
                amps[j] * cos(angle),
                amps[j] * sin(angle)
            );
        }
        
        signal.push_back(sample);
    }
    
    // 添加噪声
    auto noisySignal = addComplexNoise(signal, 20.0);  // 20dB SNR
    
    // 分析信号
    auto spectrum = analyzer.analyze(noisySignal, true);
    
    // 保存结果
    analyzer.saveSpectrumToCSV(spectrum, "spectrum_analysis.csv");
    
    std::cout << "频谱分析完成，结果已保存到 spectrum_analysis.csv" << std::endl;
}

2.4.2 数字下变频(DDC)实现

cpp

class DigitalDownConverter {
private:
    double inputSampleRate;
    double outputSampleRate;
    double centerFrequency;
    size_t decimationFactor;
    
    std::vector<double> firCoeffs;
    FIRFilter *lowPassFilter;
    
    double phase;
    double phaseIncrement;
    
public:
    DigitalDownConverter(double inputRate, double outputRate, 
                        double centerFreq, size_t numTaps = 51) :
        inputSampleRate(inputRate),
        outputSampleRate(outputRate),
        centerFrequency(centerFreq),
        phase(0.0)
    {
        // 计算抽取因子
        decimationFactor = static_cast<size_t>(std::round(inputRate / outputRate));
        
        // 确保输出采样率与抽取因子匹配
        outputSampleRate = inputSampleRate / decimationFactor;
        
        // 设置相位增量
        phaseIncrement = 2.0 * M_PI * centerFrequency / inputSampleRate;
        
        // 设计低通滤波器，截止频率为输出采样率的0.4倍(防止混叠)
        double cutoffFreq = 0.4 * outputSampleRate;
        firCoeffs = designLowPassFIR(cutoffFreq, inputSampleRate, numTaps, "blackman");
        
        // 创建滤波器
        lowPassFilter = new FIRFilter(firCoeffs);
    }
    
    ~DigitalDownConverter() {
        delete lowPassFilter;
    }
    
    // 处理输入信号
    std::vector<std::complex<double>> process(
        const std::vector<std::complex<double>>& input) {
        
        std::vector<std::complex<double>> mixedSignal;
        mixedSignal.reserve(input.size());
        
        // 1. 混频器 - 移动中心频率到DC
        for (const auto& sample : input) {
            // 使用复数乘法实现频移
            std::complex<double> localOsc(cos(phase), -sin(phase));
            mixedSignal.push_back(sample * localOsc);
            
            // 更新相位
            phase += phaseIncrement;
            if (phase >= 2.0 * M_PI) {
                phase -= 2.0 * M_PI;  // 保持相位在[0, 2π)范围内
            }
        }
        
        // 2. 低通滤波
        std::vector<std::complex<double>> filteredSignal = 
            lowPassFilter->process(mixedSignal);
        
        // 3. 抽取 (取每隔decimationFactor个样本)
        std::vector<std::complex<double>> outputSignal;
        outputSignal.reserve(filteredSignal.size() / decimationFactor + 1);
        
        for (size_t i = 0; i < filteredSignal.size(); i += decimationFactor) {
            outputSignal.push_back(filteredSignal[i]);
        }
        
        return outputSignal;
    }
    
    // 重置DDC状态
    void reset() {
        phase = 0.0;
        if (lowPassFilter) {
            lowPassFilter->reset();
        }
    }
    
    // 获取当前输出采样率
    double getOutputSampleRate() const {
        return outputSampleRate;
    }
};

// 示例：使用DDC处理高频信号
void testDDC() {
    // 参数设置
    double inputSampleRate = 50e6;     // 50 MHz输入采样率
    double outputSampleRate = 5e6;     // 5 MHz输出采样率
    double signalFreq = 12.5e6;        // 12.5 MHz信号频率
    double centerFreq = 12.5e6;        // 下变频中心频率
    double duration = 1e-4;            // 0.1毫秒
    
    // 生成测试信号
    size_t numSamples = static_cast<size_t>(inputSampleRate * duration);
    std::vector<std::complex<double>> inputSignal = 
        generateComplexSine(signalFreq, inputSampleRate, duration);
    
    // 添加噪声
    auto noisySignal = addComplexNoise(inputSignal, 15.0);  // 15dB SNR
    
    // 创建DDC
    DigitalDownConverter ddc(inputSampleRate, outputSampleRate, centerFreq);
    
    // 处理信号
    auto outputSignal = ddc.process(noisySignal);
    
    std::cout << "DDC处理完成。" << std::endl;
    std::cout << "输入样本数: " << noisySignal.size() << std::endl;
    std::cout << "输出样本数: " << outputSignal.size() << std::endl;
    std::cout << "抽取比: " << static_cast<double>(noisySignal.size()) / outputSignal.size() << std::endl;
    std::cout << "输出采样率: " << ddc.getOutputSampleRate() / 1e6 << " MHz" << std::endl;
    
    // 分析输出信号频谱
    size_t fftSize = 2048;
    SpectrumAnalyzer analyzer(fftSize, ddc.getOutputSampleRate(), "hamming");
    auto spectrum = analyzer.analyze(
        std::vector<std::complex<double>>(
            outputSignal.begin(), 
            outputSignal.begin() + std::min(outputSignal.size(), fftSize)
        )
    );
    
    // 保存结果
    analyzer.saveSpectrumToCSV(spectrum, "ddc_output_spectrum.csv");
    
    std::cout << "DDC输出信号频谱已保存到 ddc_output_spectrum.csv" << std::endl;
}

2.5 实时高频信号处理系统

2.5.1 缓冲区管理

cpp

#include <vector>
#include <complex>
#include <mutex>
#include <condition_variable>
#include <thread>
#include <atomic>

// 线程安全的循环缓冲区
template<typename T>
class CircularBuffer {
private:
    std::vector<T> buffer;
    size_t head;  // 写入位置
    size_t tail;  // 读取位置
    size_t count; // 当前元素数量
    size_t capacity;
    
    std::mutex mutex;
    std::condition_variable notEmpty;
    std::condition_variable notFull;
    
public:
    CircularBuffer(size_t size) :
        buffer(size),
        head(0),
        tail(0),
        count(0),
        capacity(size) {
    }
    
    // 写入单个元素 (阻塞直到有空间)
    void write(const T& item) {
        std::unique_lock<std::mutex> lock(mutex);
        
        // 等待直到缓冲区不满
        notFull.wait(lock, [this]() { return count < capacity; });
        
        buffer[head] = item;
        head = (head + 1) % capacity;
        ++count;
        
        // 通知可能在等待的读取者
        notEmpty.notify_one();
    }
    
    // 读取单个元素 (阻塞直到有数据)
    T read() {
        std::unique_lock<std::mutex> lock(mutex);
        
        // 等待直到缓冲区非空
        notEmpty.wait(lock, [this]() { return count > 0; });
        
        T item = buffer[tail];
        tail = (tail + 1) % capacity;
        --count;
        
        // 通知可能在等待的写入者
        notFull.notify_one();
        
        return item;
    }
    
    // 尝试写入单个元素 (非阻塞)
    bool tryWrite(const T& item) {
        std::unique_lock<std::mutex> lock(mutex, std::try_to_lock);
        if (!lock.owns_lock() || count >= capacity) {
            return false;
        }
        
        buffer[head] = item;
        head = (head + 1) % capacity;
        ++count;
        
        notEmpty.notify_one();
        return true;
    }
    
    // 尝试读取单个元素 (非阻塞)
    bool tryRead(T& item) {
        std::unique_lock<std::mutex> lock(mutex, std::try_to_lock);
        if (!lock.owns_lock() || count == 0) {
            return false;
        }
        
        item = buffer[tail];
        tail = (tail + 1) % capacity;
        --count;
        
        notFull.notify_one();
        return true;
    }
    
    // 批量写入 (最多写入size个元素，返回实际写入数量)
    size_t writeMultiple(const std::vector<T>& items, size_t maxItems = 0) {
        std::unique_lock<std::mutex> lock(mutex);
        
        if (items.empty()) {
            return 0;
        }
        
        size_t itemsToWrite = maxItems > 0 ? std::min(items.size(), maxItems) : items.size();
        
        // 等待直到有足够空间
        notFull.wait(lock, [this, itemsToWrite]() { 
            return capacity - count >= itemsToWrite; 
        });
        
        // 写入数据
        for (size_t i = 0; i < itemsToWrite; ++i) {
            buffer[head] = items[i];
            head = (head + 1) % capacity;
            ++count;
        }
        
        // 通知读取者
        if (itemsToWrite == 1) {
            notEmpty.notify_one();
        } else {
            notEmpty.notify_all();
        }
        
        return itemsToWrite;
    }
    
    // 批量读取 (尝试读取maxItems个元素，返回实际读取数量)
    size_t readMultiple(std::vector<T>& items, size_t maxItems) {
        std::unique_lock<std::mutex> lock(mutex);
        
        if (count == 0 || maxItems == 0) {
            return 0;
        }
        
        size_t itemsToRead = std::min(count, maxItems);
        items.resize(itemsToRead);
        
        // 读取数据
        for (size_t i = 0; i < itemsToRead; ++i) {
            items[i] = buffer[tail];
            tail = (tail + 1) % capacity;
            --count;
        }
        
        // 通知写入者
        if (itemsToRead == 1) {
            notFull.notify_one();
        } else {
            notFull.notify_all();
        }
        
        return itemsToRead;
    }
    
    // 获取当前元素数量
    size_t size() const {
        std::lock_guard<std::mutex> lock(mutex);
        return count;
    }
    
    // 检查缓冲区是否为空
    bool empty() const {
        std::lock_guard<std::mutex> lock(mutex);
        return count == 0;
    }
    
    // 检查缓冲区是否已满
    bool full() const {
        std::lock_guard<std::mutex> lock(mutex);
        return count == capacity;
    }
    
    // 重置缓冲区
    void clear() {
        std: