数字滤波器 C 语言实现库,包含 FIR滤波器、IIR滤波器、自适应滤波器、多速率滤波器 等常用算法,适用于嵌入式系统、信号处理、音频处理等场景。

一、项目结构

DigitalFilters/
├── include/
│   ├── fir_filter.h          ★ FIR滤波器
│   ├── iir_filter.h          ★ IIR滤波器
│   ├── adaptive_filter.h     ★ 自适应滤波器
│   ├── multirate_filter.h   ★ 多速率滤波器
│   ├── filter_design.h      ★ 滤波器设计
│   ├── window_functions.h   ★ 窗函数
│   └── math_utils.h         ★ 数学工具
├── src/
│   ├── fir_filter.c
│   ├── iir_filter.c
│   ├── adaptive_filter.c
│   ├── multirate_filter.c
│   ├── filter_design.c
│   ├── window_functions.c
│   └── math_utils.c
├── examples/
│   ├── fir_lowpass_example.c
│   ├── iir_bandpass_example.c
│   ├── lms_adaptive_example.c
│   └── decimation_example.c
├── tests/
│   ├── test_fir.c
│   ├── test_iir.c
│   └── test_performance.c
└── Makefile

二、核心代码实现

1. 数学工具头文件 (include/math_utils.h)

#ifndef MATH_UTILS_H
#define MATH_UTILS_H

#include <stdint.h>
#include <stddef.h>
#include <math.h>

// 浮点类型定义
typedef float float32_t;
typedef double float64_t;

// 复数结构
typedef struct {
    float32_t real;
    float32_t imag;
} complex_float_t;

// 数学常数
#define PI_F 3.14159265358979323846f
#define TWO_PI_F 6.28318530717958647692f
#define SQRT2_F 1.41421356237309504880f

// 数学函数宏
#define ABS(x) ((x) < 0 ? -(x) : (x))
#define MIN(a, b) ((a) < (b) ? (a) : (b))
#define MAX(a, b) ((a) > (b) ? (a) : (b))
#define CLAMP(x, min, max) (MIN(MAX((x), (min)), (max)))

// 正弦余弦近似(嵌入式优化)
float32_t fast_sin(float32_t x);
float32_t fast_cos(float32_t x);
float32_t fast_sqrt(float32_t x);

// 复数运算
complex_float_t complex_add(complex_float_t a, complex_float_t b);
complex_float_t complex_mul(complex_float_t a, complex_float_t b);
float32_t complex_mag(complex_float_t a);

#endif // MATH_UTILS_H

2. FIR 滤波器头文件 (include/fir_filter.h)

#ifndef FIR_FILTER_H
#define FIR_FILTER_H

#include <stdint.h>
#include <stddef.h>

// FIR滤波器结构
typedef struct {
    float32_t *coeffs;        // 滤波器系数
    float32_t *buffer;        // 延迟线缓冲区
    uint16_t num_taps;        // 抽头数
    uint16_t buffer_index;    // 缓冲区索引
    uint8_t initialized;       // 初始化标志
} fir_filter_t;

// FIR滤波器配置
typedef struct {
    float32_t *coeffs;        // 系数数组
    uint16_t num_taps;        // 抽头数
} fir_config_t;

// 滤波器类型
typedef enum {
    FIR_LOWPASS = 0,
    FIR_HIGHPASS,
    FIR_BANDPASS,
    FIR_BANDSTOP,
    FIR_ALLPASS
} fir_filter_type_t;

// 窗函数类型
typedef enum {
    WINDOW_RECTANGULAR = 0,
    WINDOW_HANN,
    WINDOW_HAMMING,
    WINDOW_BLACKMAN,
    WINDOW_KAISER
} window_type_t;

// 函数声明
int fir_init(fir_filter_t *filter, const fir_config_t *config);
float32_t fir_process(fir_filter_t *filter, float32_t input);
void fir_reset(fir_filter_t *filter);
void fir_set_coeffs(fir_filter_t *filter, const float32_t *coeffs);
void fir_get_coeffs(const fir_filter_t *filter, float32_t *coeffs);
void fir_deinit(fir_filter_t *filter);

// 批量处理函数
void fir_process_block(fir_filter_t *filter, 
                       const float32_t *input, 
                       float32_t *output, 
                       uint16_t block_size);

// 滤波器设计函数
int fir_design_lowpass(fir_filter_t *filter,
                      float32_t cutoff_freq,
                      float32_t sample_rate,
                      uint16_t num_taps,
                      window_type_t window);

int fir_design_highpass(fir_filter_t *filter,
                       float32_t cutoff_freq,
                       float32_t sample_rate,
                       uint16_t num_taps,
                       window_type_t window);

int fir_design_bandpass(fir_filter_t *filter,
                       float32_t low_cutoff,
                       float32_t high_cutoff,
                       float32_t sample_rate,
                       uint16_t num_taps,
                       window_type_t window);

// 对称FIR优化版本(偶数抽头)
float32_t fir_process_symmetrical(fir_filter_t *filter, float32_t input);

// 实时处理宏(高性能)
#define FIR_PROCESS_MACRO(filter, input) \
({ \
    float32_t output = 0.0f; \
    filter.buffer[filter.buffer_index] = input; \
    uint16_t idx = filter.buffer_index; \
    for (uint16_t i = 0; i < filter.num_taps; i++) { \
        output += filter.coeffs[i] * filter.buffer[idx]; \
        if (idx == 0) idx = filter.num_taps - 1; \
        else idx--; \
    } \
    filter.buffer_index = (filter.buffer_index + 1) % filter.num_taps; \
    output; \
})

#endif // FIR_FILTER_H

3. FIR 滤波器实现 (src/fir_filter.c)

#include "fir_filter.h"
#include "math_utils.h"
#include "window_functions.h"
#include <stdlib.h>
#include <string.h>

// 初始化FIR滤波器
int fir_init(fir_filter_t *filter, const fir_config_t *config) {
    if (!filter || !config || config->num_taps == 0) {
        return -1;
    }
    
    // 分配内存
    filter->coeffs = (float32_t *)malloc(config->num_taps * sizeof(float32_t));
    filter->buffer = (float32_t *)malloc(config->num_taps * sizeof(float32_t));
    
    if (!filter->coeffs || !filter->buffer) {
        fir_deinit(filter);
        return -1;
    }
    
    // 复制系数
    memcpy(filter->coeffs, config->coeffs, config->num_taps * sizeof(float32_t));
    
    // 清零缓冲区
    memset(filter->buffer, 0, config->num_taps * sizeof(float32_t));
    
    filter->num_taps = config->num_taps;
    filter->buffer_index = 0;
    filter->initialized = 1;
    
    return 0;
}

// 处理单个样本
float32_t fir_process(fir_filter_t *filter, float32_t input) {
    if (!filter || !filter->initialized) {
        return 0.0f;
    }
    
    float32_t output = 0.0f;
    uint16_t i;
    
    // 将新样本存入缓冲区
    filter->buffer[filter->buffer_index] = input;
    
    // 卷积计算
    uint16_t buffer_idx = filter->buffer_index;
    for (i = 0; i < filter->num_taps; i++) {
        output += filter->coeffs[i] * filter->buffer[buffer_idx];
        
        // 循环缓冲区索引
        if (buffer_idx == 0) {
            buffer_idx = filter->num_taps - 1;
        } else {
            buffer_idx--;
        }
    }
    
    // 更新缓冲区索引
    filter->buffer_index = (filter->buffer_index + 1) % filter->num_taps;
    
    return output;
}

// 处理数据块
void fir_process_block(fir_filter_t *filter, 
                      const float32_t *input, 
                      float32_t *output, 
                      uint16_t block_size) {
    if (!filter || !input || !output || !filter->initialized) {
        return;
    }
    
    for (uint16_t i = 0; i < block_size; i++) {
        output[i] = fir_process(filter, input[i]);
    }
}

// 对称FIR优化处理(适用于线性相位滤波器)
float32_t fir_process_symmetrical(fir_filter_t *filter, float32_t input) {
    if (!filter || !filter->initialized) {
        return 0.0f;
    }
    
    float32_t output = 0.0f;
    uint16_t half_taps = filter->num_taps / 2;
    
    // 存储新样本
    filter->buffer[filter->buffer_index] = input;
    
    // 对称卷积
    uint16_t idx1 = filter->buffer_index;
    uint16_t idx2 = (filter->buffer_index + half_taps) % filter->num_taps;
    
    for (uint16_t i = 0; i < half_taps; i++) {
        output += filter->coeffs[i] * (filter->buffer[idx1] + filter->buffer[idx2]);
        
        // 更新索引
        idx1 = (idx1 + 1) % filter->num_taps;
        if (idx2 == 0) {
            idx2 = filter->num_taps - 1;
        } else {
            idx2--;
        }
    }
    
    // 处理中心样本(奇数抽头)
    if (filter->num_taps % 2 == 1) {
        output += filter->coeffs[half_taps] * filter->buffer[idx1];
    }
    
    filter->buffer_index = (filter->buffer_index + 1) % filter->num_taps;
    return output;
}

// 重置滤波器状态
void fir_reset(fir_filter_t *filter) {
    if (!filter || !filter->initialized) {
        return;
    }
    
    memset(filter->buffer, 0, filter->num_taps * sizeof(float32_t));
    filter->buffer_index = 0;
}

// 设置滤波器系数
void fir_set_coeffs(fir_filter_t *filter, const float32_t *coeffs) {
    if (!filter || !coeffs || !filter->initialized) {
        return;
    }
    
    memcpy(filter->coeffs, coeffs, filter->num_taps * sizeof(float32_t));
}

// 获取滤波器系数
void fir_get_coeffs(const fir_filter_t *filter, float32_t *coeffs) {
    if (!filter || !coeffs || !filter->initialized) {
        return;
    }
    
    memcpy(coeffs, filter->coeffs, filter->num_taps * sizeof(float32_t));
}

// 设计低通FIR滤波器
int fir_design_lowpass(fir_filter_t *filter,
                     float32_t cutoff_freq,
                     float32_t sample_rate,
                     uint16_t num_taps,
                     window_type_t window) {
    if (!filter || num_taps == 0 || cutoff_freq <= 0 || sample_rate <= 0) {
        return -1;
    }
    
    // 归一化截止频率
    float32_t norm_cutoff = cutoff_freq / sample_rate;
    
    // 生成理想低通滤波器系数
    float32_t *coeffs = (float32_t *)malloc(num_taps * sizeof(float32_t));
    if (!coeffs) {
        return -1;
    }
    
    int center = num_taps / 2;
    float32_t sum = 0.0f;
    
    for (int i = 0; i < num_taps; i++) {
        if (i == center) {
            coeffs[i] = 2.0f * norm_cutoff;
        } else {
            float32_t angle = PI_F * (i - center) * 2.0f * norm_cutoff;
            coeffs[i] = sinf(angle) / (PI_F * (i - center));
        }
        sum += coeffs[i];
    }
    
    // 归一化
    for (int i = 0; i < num_taps; i++) {
        coeffs[i] /= sum;
    }
    
    // 应用窗函数
    apply_window(coeffs, num_taps, window);
    
    // 初始化滤波器
    fir_config_t config = {
        .coeffs = coeffs,
        .num_taps = num_taps
    };
    
    int result = fir_init(filter, &config);
    free(coeffs);
    
    return result;
}

// 设计高通FIR滤波器
int fir_design_highpass(fir_filter_t *filter,
                       float32_t cutoff_freq,
                       float32_t sample_rate,
                       uint16_t num_taps,
                       window_type_t window) {
    if (!filter || num_taps == 0 || cutoff_freq <= 0 || sample_rate <= 0) {
        return -1;
    }
    
    // 先设计低通滤波器
    if (fir_design_lowpass(filter, cutoff_freq, sample_rate, num_taps, window) != 0) {
        return -1;
    }
    
    // 转换为高通:系数取反(除了中心系数)
    int center = num_taps / 2;
    for (int i = 0; i < num_taps; i++) {
        if (i != center) {
            filter->coeffs[i] = -filter->coeffs[i];
        }
    }
    
    return 0;
}

// 设计带通FIR滤波器
int fir_design_bandpass(fir_filter_t *filter,
                       float32_t low_cutoff,
                       float32_t high_cutoff,
                       float32_t sample_rate,
                       uint16_t num_taps,
                       window_type_t window) {
    if (!filter || num_taps == 0 || low_cutoff <= 0 || high_cutoff <= 0 || 
        low_cutoff >= high_cutoff || sample_rate <= 0) {
        return -1;
    }
    
    // 设计两个低通滤波器
    fir_filter_t low_filter1, low_filter2;
    if (fir_design_lowpass(&low_filter1, high_cutoff, sample_rate, num_taps, window) != 0 ||
        fir_design_lowpass(&low_filter2, low_cutoff, sample_rate, num_taps, window) != 0) {
        return -1;
    }
    
    // 带通 = 高通(low_cutoff) - 低通(high_cutoff)
    for (int i = 0; i < num_taps; i++) {
        filter->coeffs[i] = low_filter2.coeffs[i] - low_filter1.coeffs[i];
    }
    
    // 释放临时滤波器
    fir_deinit(&low_filter1);
    fir_deinit(&low_filter2);
    
    return 0;
}

// 释放滤波器资源
void fir_deinit(fir_filter_t *filter) {
    if (!filter) {
        return;
    }
    
    if (filter->coeffs) {
        free(filter->coeffs);
        filter->coeffs = NULL;
    }
    
    if (filter->buffer) {
        free(filter->buffer);
        filter->buffer = NULL;
    }
    
    filter->initialized = 0;
    filter->num_taps = 0;
    filter->buffer_index = 0;
}

4. IIR 滤波器头文件 (include/iir_filter.h)

#ifndef IIR_FILTER_H
#define IIR_FILTER_H

#include <stdint.h>
#include <stddef.h>

// IIR滤波器结构
typedef struct {
    float32_t *b_coeffs;        // 分子系数 (feedforward)
    float32_t *a_coeffs;        // 分母系数 (feedback)
    float32_t *x_buffer;        // 输入延迟线
    float32_t *y_buffer;        // 输出延迟线
    uint8_t num_b;              // 分子阶数
    uint8_t num_a;              // 分母阶数
    uint8_t x_index;            // 输入缓冲区索引
    uint8_t y_index;            // 输出缓冲区索引
    uint8_t initialized;        // 初始化标志
} iir_filter_t;

// IIR滤波器配置
typedef struct {
    float32_t *b_coeffs;        // 分子系数
    float32_t *a_coeffs;        // 分母系数
    uint8_t num_b;              // 分子阶数
    uint8_t num_a;              // 分母阶数
} iir_config_t;

// IIR滤波器类型
typedef enum {
    IIR_LOWPASS = 0,
    IIR_HIGHPASS,
    IIR_BANDPASS,
    IIR_BANDSTOP,
    IIR_NOTCH
} iir_filter_type_t;

// 双二阶滤波器结构(二阶节)
typedef struct {
    float32_t b0, b1, b2;       // 分子系数
    float32_t a1, a2;           // 分母系数(a0=1)
    float32_t x1, x2;           // 输入延迟
    float32_t y1, y2;           // 输出延迟
} biquad_t;

// 函数声明
int iir_init(iir_filter_t *filter, const iir_config_t *config);
float32_t iir_process(iir_filter_t *filter, float32_t input);
void iir_reset(iir_filter_t *filter);
void iir_deinit(iir_filter_t *filter);

// 批量处理
void iir_process_block(iir_filter_t *filter,
                      const float32_t *input,
                      float32_t *output,
                      uint16_t block_size);

// 双二阶滤波器(更高效)
int biquad_init(biquad_t *filter, 
                float32_t freq, 
                float32_t q, 
                float32_t gain,
                iir_filter_type_t type,
                float32_t sample_rate);
float32_t biquad_process(biquad_t *filter, float32_t input);
void biquad_reset(biquad_t *filter);

// 级联双二阶滤波器
typedef struct {
    biquad_t *stages;
    uint8_t num_stages;
} cascade_biquad_t;

int cascade_biquad_init(cascade_biquad_t *cascade,
                       const float32_t *freqs,
                       const float32_t *qs,
                       const float32_t *gains,
                       const iir_filter_type_t *types,
                       uint8_t num_stages,
                       float32_t sample_rate);
float32_t cascade_biquad_process(cascade_biquad_t *cascade, float32_t input);
void cascade_biquad_reset(cascade_biquad_t *cascade);
void cascade_biquad_deinit(cascade_biquad_t *cascade);

// 滤波器设计
int iir_design_lowpass(iir_filter_t *filter,
                      float32_t cutoff_freq,
                      float32_t sample_rate,
                      uint8_t order);

int iir_design_highpass(iir_filter_t *filter,
                       float32_t cutoff_freq,
                       float32_t sample_rate,
                       uint8_t order);

// 直接形式II(规范型)实现
float32_t iir_process_df2(iir_filter_t *filter, float32_t input);

#endif // IIR_FILTER_H

5. IIR 滤波器实现 (src/iir_filter.c)

#include "iir_filter.h"
#include "math_utils.h"
#include <stdlib.h>
#include <string.h>

// 初始化IIR滤波器
int iir_init(iir_filter_t *filter, const iir_config_t *config) {
    if (!filter || !config || config->num_b == 0 || config->num_a == 0) {
        return -1;
    }
    
    // 分配内存
    filter->b_coeffs = (float32_t *)malloc(config->num_b * sizeof(float32_t));
    filter->a_coeffs = (float32_t *)malloc(config->num_a * sizeof(float32_t));
    filter->x_buffer = (float32_t *)malloc(config->num_b * sizeof(float32_t));
    filter->y_buffer = (float32_t *)malloc((config->num_a - 1) * sizeof(float32_t));
    
    if (!filter->b_coeffs || !filter->a_coeffs || !filter->x_buffer || !filter->y_buffer) {
        iir_deinit(filter);
        return -1;
    }
    
    // 复制系数
    memcpy(filter->b_coeffs, config->b_coeffs, config->num_b * sizeof(float32_t));
    memcpy(filter->a_coeffs, config->a_coeffs, config->num_a * sizeof(float32_t));
    
    // 清零缓冲区
    memset(filter->x_buffer, 0, config->num_b * sizeof(float32_t));
    memset(filter->y_buffer, 0, (config->num_a - 1) * sizeof(float32_t));
    
    filter->num_b = config->num_b;
    filter->num_a = config->num_a;
    filter->x_index = 0;
    filter->y_index = 0;
    filter->initialized = 1;
    
    return 0;
}

// 处理单个样本(直接形式I)
float32_t iir_process(iir_filter_t *filter, float32_t input) {
    if (!filter || !filter->initialized) {
        return 0.0f;
    }
    
    float32_t output = 0.0f;
    uint8_t i;
    
    // 存储输入样本
    filter->x_buffer[filter->x_index] = input;
    
    // 计算输出:y[n] = sum(b_k * x[n-k]) - sum(a_k * y[n-k])
    uint8_t x_idx = filter->x_index;
    uint8_t y_idx = filter->y_index;
    
    // 分子部分(前馈)
    for (i = 0; i < filter->num_b; i++) {
        output += filter->b_coeffs[i] * filter->x_buffer[x_idx];
        if (x_idx == 0) {
            x_idx = filter->num_b - 1;
        } else {
            x_idx--;
        }
    }
    
    // 分母部分(反馈)
    for (i = 0; i < filter->num_a - 1; i++) {
        output -= filter->a_coeffs[i + 1] * filter->y_buffer[y_idx];
        if (y_idx == 0) {
            y_idx = filter->num_a - 2;
        } else {
            y_idx--;
        }
    }
    
    // 除以a0(通常a0=1,这里省略除法以提高效率)
    output /= filter->a_coeffs[0];
    
    // 存储输出样本
    filter->y_buffer[filter->y_index] = output;
    
    // 更新索引
    filter->x_index = (filter->x_index + 1) % filter->num_b;
    filter->y_index = (filter->y_index + 1) % (filter->num_a - 1);
    
    return output;
}

// 直接形式II(规范型)实现
float32_t iir_process_df2(iir_filter_t *filter, float32_t input) {
    if (!filter || !filter->initialized) {
        return 0.0f;
    }
    
    float32_t w = input;
    float32_t output = 0.0f;
    uint8_t i;
    
    // 计算反馈部分
    for (i = 1; i < filter->num_a; i++) {
        w -= filter->a_coeffs[i] * filter->y_buffer[(filter->y_index + i - 1) % (filter->num_a - 1)];
    }
    
    // 计算前馈部分
    for (i = 0; i < filter->num_b; i++) {
        output += filter->b_coeffs[i] * filter->y_buffer[(filter->y_index + i) % (filter->num_b - 1)];
    }
    
    // 更新延迟线
    filter->y_buffer[filter->y_index] = w;
    filter->y_index = (filter->y_index + 1) % (filter->num_a - 1);
    
    return output;
}

// 双二阶滤波器处理
float32_t biquad_process(biquad_t *filter, float32_t input) {
    if (!filter) {
        return 0.0f;
    }
    
    // 直接形式II转置
    float32_t output = filter->b0 * input + filter->x1;
    filter->x1 = filter->b1 * input - filter->a1 * output + filter->x2;
    filter->x2 = filter->b2 * input - filter->a2 * output;
    
    return output;
}

// 设计双二阶低通滤波器
int biquad_init(biquad_t *filter, 
                float32_t freq, 
                float32_t q, 
                float32_t gain,
                iir_filter_type_t type,
                float32_t sample_rate) {
    if (!filter || freq <= 0 || q <= 0 || sample_rate <= 0) {
        return -1;
    }
    
    float32_t omega = TWO_PI_F * freq / sample_rate;
    float32_t sn = sinf(omega);
    float32_t cs = cosf(omega);
    float32_t alpha = sn / (2.0f * q);
    float32_t A = powf(10.0f, gain / 40.0f);
    float32_t beta = sqrtf(A + A);
    
    switch (type) {
        case IIR_LOWPASS:
            filter->b0 = (1.0f - cs) / 2.0f;
            filter->b1 = 1.0f - cs;
            filter->b2 = (1.0f - cs) / 2.0f;
            filter->a1 = -2.0f * cs;
            filter->a2 = 1.0f - alpha;
            break;
            
        case IIR_HIGHPASS:
            filter->b0 = (1.0f + cs) / 2.0f;
            filter->b1 = -(1.0f + cs);
            filter->b2 = (1.0f + cs) / 2.0f;
            filter->a1 = -2.0f * cs;
            filter->a2 = 1.0f - alpha;
            break;
            
        case IIR_BANDPASS:
            filter->b0 = alpha;
            filter->b1 = 0.0f;
            filter->b2 = -alpha;
            filter->a1 = -2.0f * cs;
            filter->a2 = 1.0f - alpha;
            break;
            
        case IIR_NOTCH:
            filter->b0 = 1.0f;
            filter->b1 = -2.0f * cs;
            filter->b2 = 1.0f;
            filter->a1 = -2.0f * cs;
            filter->a2 = 1.0f - alpha;
            break;
            
        default:
            return -1;
    }
    
    // 归一化
    float32_t a0 = 1.0f + alpha;
    filter->b0 /= a0;
    filter->b1 /= a0;
    filter->b2 /= a0;
    filter->a1 /= a0;
    filter->a2 /= a0;
    
    // 清零延迟
    filter->x1 = filter->x2 = filter->y1 = filter->y2 = 0.0f;
    
    return 0;
}

// 级联双二阶滤波器处理
float32_t cascade_biquad_process(cascade_biquad_t *cascade, float32_t input) {
    if (!cascade || !cascade->stages) {
        return 0.0f;
    }
    
    float32_t output = input;
    for (uint8_t i = 0; i < cascade->num_stages; i++) {
        output = biquad_process(&cascade->stages[i], output);
    }
    
    return output;
}

// 重置IIR滤波器
void iir_reset(iir_filter_t *filter) {
    if (!filter || !filter->initialized) {
        return;
    }
    
    memset(filter->x_buffer, 0, filter->num_b * sizeof(float32_t));
    memset(filter->y_buffer, 0, (filter->num_a - 1) * sizeof(float32_t));
    filter->x_index = 0;
    filter->y_index = 0;
}

// 释放IIR滤波器资源
void iir_deinit(iir_filter_t *filter) {
    if (!filter) {
        return;
    }
    
    if (filter->b_coeffs) {
        free(filter->b_coeffs);
        filter->b_coeffs = NULL;
    }
    
    if (filter->a_coeffs) {
        free(filter->a_coeffs);
        filter->a_coeffs = NULL;
    }
    
    if (filter->x_buffer) {
        free(filter->x_buffer);
        filter->x_buffer = NULL;
    }
    
    if (filter->y_buffer) {
        free(filter->y_buffer);
        filter->y_buffer = NULL;
    }
    
    filter->initialized = 0;
}

6. 自适应滤波器 (include/adaptive_filter.h)

#ifndef ADAPTIVE_FILTER_H
#define ADAPTIVE_FILTER_H

#include <stdint.h>
#include <stddef.h>

// 自适应滤波器类型
typedef enum {
    ADAPTIVE_LMS = 0,      // 最小均方
    ADAPTIVE_NLMS,        // 归一化LMS
    ADAPTIVE_RLS,         // 递归最小二乘
    ADAPTIVE_APA          // 仿射投影算法
} adaptive_algorithm_t;

// LMS自适应滤波器结构
typedef struct {
    float32_t *weights;      // 权系数
    float32_t *buffer;      // 输入缓冲区
    float32_t step_size;    // 步长
    float32_t leakage;      // 泄漏因子
    uint16_t order;         // 滤波器阶数
    uint16_t buffer_index;  // 缓冲区索引
    uint8_t initialized;    // 初始化标志
} lms_filter_t;

// NLMS自适应滤波器结构
typedef struct {
    float32_t *weights;      // 权系数
    float32_t *buffer;      // 输入缓冲区
    float32_t step_size;    // 步长
    float32_t eps;          // 正则化参数
    uint16_t order;         // 滤波器阶数
    uint16_t buffer_index;  // 缓冲区索引
    uint8_t initialized;    // 初始化标志
} nlms_filter_t;

// RLS自适应滤波器结构
typedef struct {
    float32_t *weights;      // 权系数
    float32_t *p_matrix;    // 逆相关矩阵
    float32_t *buffer;      // 输入缓冲区
    float32_t lambda;       // 遗忘因子
    float32_t delta;        // 初始化参数
    uint16_t order;         // 滤波器阶数
    uint16_t buffer_index;  // 缓冲区索引
    uint8_t initialized;    // 初始化标志
} rls_filter_t;

// 函数声明
int lms_init(lms_filter_t *filter, uint16_t order, float32_t step_size);
float32_t lms_process(lms_filter_t *filter, float32_t input, float32_t desired);
void lms_reset(lms_filter_t *filter);
void lms_deinit(lms_filter_t *filter);

int nlms_init(nlms_filter_t *filter, uint16_t order, float32_t step_size, float32_t eps);
float32_t nlms_process(nlms_filter_t *filter, float32_t input, float32_t desired);
void nlms_reset(nlms_filter_t *filter);
void nlms_deinit(nlms_filter_t *filter);

int rls_init(rls_filter_t *filter, uint16_t order, float32_t lambda, float32_t delta);
float32_t rls_process(rls_filter_t *filter, float32_t input, float32_t desired);
void rls_reset(rls_filter_t *filter);
void rls_deinit(rls_filter_t *filter);

// 噪声消除器
typedef struct {
    lms_filter_t adaptive_filter;
    float32_t *noise_reference;
    uint16_t buffer_size;
} noise_canceller_t;

int noise_canceller_init(noise_canceller_t *nc, uint16_t order, float32_t step_size, uint16_t buffer_size);
float32_t noise_canceller_process(noise_canceller_t *nc, float32_t signal_plus_noise, float32_t noise_reference);
void noise_canceller_deinit(noise_canceller_t *nc);

#endif // ADAPTIVE_FILTER_H

7. 自适应滤波器实现 (src/adaptive_filter.c)

#include "adaptive_filter.h"
#include "math_utils.h"
#include <stdlib.h>
#include <string.h>

// 初始化LMS滤波器
int lms_init(lms_filter_t *filter, uint16_t order, float32_t step_size) {
    if (!filter || order == 0 || step_size <= 0) {
        return -1;
    }
    
    filter->weights = (float32_t *)malloc(order * sizeof(float32_t));
    filter->buffer = (float32_t *)malloc(order * sizeof(float32_t));
    
    if (!filter->weights || !filter->buffer) {
        lms_deinit(filter);
        return -1;
    }
    
    memset(filter->weights, 0, order * sizeof(float32_t));
    memset(filter->buffer, 0, order * sizeof(float32_t));
    
    filter->step_size = step_size;
    filter->leakage = 1.0f;
    filter->order = order;
    filter->buffer_index = 0;
    filter->initialized = 1;
    
    return 0;
}

// LMS处理
float32_t lms_process(lms_filter_t *filter, float32_t input, float32_t desired) {
    if (!filter || !filter->initialized) {
        return 0.0f;
    }
    
    // 存储输入样本
    filter->buffer[filter->buffer_index] = input;
    
    // 计算输出
    float32_t output = 0.0f;
    uint16_t idx = filter->buffer_index;
    for (uint16_t i = 0; i < filter->order; i++) {
        output += filter->weights[i] * filter->buffer[idx];
        if (idx == 0) {
            idx = filter->order - 1;
        } else {
            idx--;
        }
    }
    
    // 计算误差
    float32_t error = desired - output;
    
    // 更新权系数(带泄漏)
    idx = filter->buffer_index;
    for (uint16_t i = 0; i < filter->order; i++) {
        filter->weights[i] = filter->leakage * filter->weights[i] + 
                           filter->step_size * error * filter->buffer[idx];
        if (idx == 0) {
            idx = filter->order - 1;
        } else {
            idx--;
        }
    }
    
    // 更新缓冲区索引
    filter->buffer_index = (filter->buffer_index + 1) % filter->order;
    
    return output;
}

// 初始化NLMS滤波器
int nlms_init(nlms_filter_t *filter, uint16_t order, float32_t step_size, float32_t eps) {
    if (!filter || order == 0 || step_size <= 0) {
        return -1;
    }
    
    filter->weights = (float32_t *)malloc(order * sizeof(float32_t));
    filter->buffer = (float32_t *)malloc(order * sizeof(float32_t));
    
    if (!filter->weights || !filter->buffer) {
        nlms_deinit(filter);
        return -1;
    }
    
    memset(filter->weights, 0, order * sizeof(float32_t));
    memset(filter->buffer, 0, order * sizeof(float32_t));
    
    filter->step_size = step_size;
    filter->eps = eps;
    filter->order = order;
    filter->buffer_index = 0;
    filter->initialized = 1;
    
    return 0;
}

// NLMS处理
float32_t nlms_process(nlms_filter_t *filter, float32_t input, float32_t desired) {
    if (!filter || !filter->initialized) {
        return 0.0f;
    }
    
    // 存储输入样本
    filter->buffer[filter->buffer_index] = input;
    
    // 计算输出
    float32_t output = 0.0f;
    uint16_t idx = filter->buffer_index;
    for (uint16_t i = 0; i < filter->order; i++) {
        output += filter->weights[i] * filter->buffer[idx];
        if (idx == 0) {
            idx = filter->order - 1;
        } else {
            idx--;
        }
    }
    
    // 计算误差
    float32_t error = desired - output;
    
    // 计算输入功率
    float32_t power = filter->eps;
    idx = filter->buffer_index;
    for (uint16_t i = 0; i < filter->order; i++) {
        power += filter->buffer[idx] * filter->buffer[idx];
        if (idx == 0) {
            idx = filter->order - 1;
        } else {
            idx--;
        }
    }
    
    // 更新权系数
    float32_t mu_norm = filter->step_size / power;
    idx = filter->buffer_index;
    for (uint16_t i = 0; i < filter->order; i++) {
        filter->weights[i] += mu_norm * error * filter->buffer[idx];
        if (idx == 0) {
            idx = filter->order - 1;
        } else {
            idx--;
        }
    }
    
    // 更新缓冲区索引
    filter->buffer_index = (filter->buffer_index + 1) % filter->order;
    
    return output;
}

// 噪声消除器处理
float32_t noise_canceller_process(noise_canceller_t *nc, float32_t signal_plus_noise, float32_t noise_reference) {
    if (!nc || !nc->adaptive_filter.initialized) {
        return 0.0f;
    }
    
    // 使用自适应滤波器估计噪声
    float32_t estimated_noise = lms_process(&nc->adaptive_filter, noise_reference, signal_plus_noise);
    
    // 从混合信号中减去估计的噪声
    float32_t clean_signal = signal_plus_noise - estimated_noise;
    
    return clean_signal;
}

// 释放LMS滤波器资源
void lms_deinit(lms_filter_t *filter) {
    if (!filter) {
        return;
    }
    
    if (filter->weights) {
        free(filter->weights);
        filter->weights = NULL;
    }
    
    if (filter->buffer) {
        free(filter->buffer);
        filter->buffer = NULL;
    }
    
    filter->initialized = 0;
}

8. 多速率滤波器 (include/multirate_filter.h)

#ifndef MULTIRATE_FILTER_H
#define MULTIRATE_FILTER_H

#include <stdint.h>
#include <stddef.h>
#include "fir_filter.h"

// 抽取器结构
typedef struct {
    fir_filter_t antialiasing_filter;  // 抗混叠滤波器
    uint16_t decimation_factor;        // 抽取因子
    uint16_t phase;                   // 相位
    uint8_t initialized;              // 初始化标志
} decimator_t;

// 插值器结构
typedef struct {
    fir_filter_t interpolation_filter; // 插值滤波器
    uint16_t interpolation_factor;     // 插值因子
    uint16_t phase;                   // 相位
    uint8_t initialized;              // 初始化标志
} interpolator_t;

// 多相滤波器结构
typedef struct {
    fir_filter_t **polyphase_filters;  // 多相滤波器组
    uint16_t num_phases;              // 相位数
    uint16_t current_phase;           // 当前相位
    uint8_t initialized;              // 初始化标志
} polyphase_filter_t;

// 分数倍采样率转换器
typedef struct {
    interpolator_t interpolator;      // 插值器
    decimator_t decimator;            // 抽取器
    float32_t rate_ratio;             // 采样率比
    uint8_t initialized;              // 初始化标志
} resampler_t;

// 函数声明
int decimator_init(decimator_t *decimator, 
                   uint16_t decimation_factor,
                   uint16_t filter_order,
                   float32_t cutoff_freq,
                   float32_t sample_rate);
float32_t decimator_process(decimator_t *decimator, float32_t input);
void decimator_process_block(decimator_t *decimator,
                            const float32_t *input,
                            float32_t *output,
                            uint16_t input_size,
                            uint16_t *output_size);
void decimator_deinit(decimator_t *decimator);

int interpolator_init(interpolator_t *interpolator,
                      uint16_t interpolation_factor,
                      uint16_t filter_order,
                      float32_t cutoff_freq,
                      float32_t sample_rate);
float32_t interpolator_process(interpolator_t *interpolator, float32_t input);
void interpolator_process_block(interpolator_t *interpolator,
                               const float32_t *input,
                               float32_t *output,
                               uint16_t input_size,
                               uint16_t *output_size);
void interpolator_deinit(interpolator_t *interpolator);

int resampler_init(resampler_t *resampler,
                   float32_t input_rate,
                   float32_t output_rate,
                   uint16_t filter_order);
float32_t resampler_process(resampler_t *resampler, float32_t input);
void resampler_deinit(resampler_t *resampler);

// 高效的多相实现
int polyphase_filter_init(polyphase_filter_t *filter,
                         const float32_t *coeffs,
                         uint16_t num_coeffs,
                         uint16_t num_phases);
float32_t polyphase_filter_process(polyphase_filter_t *filter, float32_t input);
void polyphase_filter_deinit(polyphase_filter_t *filter);

#endif // MULTIRATE_FILTER_H

9. 窗函数 (include/window_functions.h)

#ifndef WINDOW_FUNCTIONS_H
#define WINDOW_FUNCTIONS_H

#include <stdint.h>
#include <stddef.h>

// 窗函数类型
typedef enum {
    WINDOW_RECTANGULAR = 0,
    WINDOW_HANN,
    WINDOW_HAMMING,
    WINDOW_BLACKMAN,
    WINDOW_KAISER,
    WINDOW_GAUSSIAN,
    WINDOW_TRIANGULAR
} window_type_t;

// 函数声明
void apply_window(float32_t *data, uint16_t length, window_type_t window);
void rectangular_window(float32_t *data, uint16_t length);
void hann_window(float32_t *data, uint16_t length);
void hamming_window(float32_t *data, uint16_t length);
void blackman_window(float32_t *data, uint16_t length);
void kaiser_window(float32_t *data, uint16_t length, float32_t beta);
void gaussian_window(float32_t *data, uint16_t length, float32_t sigma);

// 计算窗函数的参数
float32_t calculate_kaiser_beta(float32_t ripple_db);
uint16_t calculate_filter_order(float32_t transition_width, float32_t sampling_freq, window_type_t window);

#endif // WINDOW_FUNCTIONS_H

10. 窗函数实现 (src/window_functions.c)

#include "window_functions.h"
#include "math_utils.h"
#include <math.h>

// 应用窗函数
void apply_window(float32_t *data, uint16_t length, window_type_t window) {
    if (!data || length == 0) {
        return;
    }
    
    switch (window) {
        case WINDOW_RECTANGULAR:
            rectangular_window(data, length);
            break;
        case WINDOW_HANN:
            hann_window(data, length);
            break;
        case WINDOW_HAMMING:
            hamming_window(data, length);
            break;
        case WINDOW_BLACKMAN:
            blackman_window(data, length);
            break;
        case WINDOW_KAISER:
            kaiser_window(data, length, 8.6f);  // 默认beta值
            break;
        case WINDOW_GAUSSIAN:
            gaussian_window(data, length, 0.5f);  // 默认sigma
            break;
        default:
            rectangular_window(data, length);
            break;
    }
}

// 矩形窗
void rectangular_window(float32_t *data, uint16_t length) {
    for (uint16_t i = 0; i < length; i++) {
        data[i] = 1.0f;
    }
}

// 汉宁窗
void hann_window(float32_t *data, uint16_t length) {
    for (uint16_t i = 0; i < length; i++) {
        float32_t factor = 2.0f * PI_F * i / (length - 1);
        data[i] = 0.5f * (1.0f - cosf(factor));
    }
}

// 汉明窗
void hamming_window(float32_t *data, uint16_t length) {
    for (uint16_t i = 0; i < length; i++) {
        float32_t factor = 2.0f * PI_F * i / (length - 1);
        data[i] = 0.54f - 0.46f * cosf(factor);
    }
}

// 布莱克曼窗
void blackman_window(float32_t *data, uint16_t length) {
    for (uint16_t i = 0; i < length; i++) {
        float32_t factor = 2.0f * PI_F * i / (length - 1);
        data[i] = 0.42f - 0.5f * cosf(factor) + 0.08f * cosf(2.0f * factor);
    }
}

// 凯泽窗
void kaiser_window(float32_t *data, uint16_t length, float32_t beta) {
    float32_t i0_beta = bessel_i0(beta);
    
    for (uint16_t i = 0; i < length; i++) {
        float32_t arg = beta * sqrtf(1.0f - powf(2.0f * i / (length - 1) - 1.0f, 2.0f));
        data[i] = bessel_i0(arg) / i0_beta;
    }
}

// 高斯窗
void gaussian_window(float32_t *data, uint16_t length, float32_t sigma) {
    float32_t scale = -0.5f / (sigma * sigma);
    
    for (uint16_t i = 0; i < length; i++) {
        float32_t x = (i - length / 2.0f) / (length / 2.0f);
        data[i] = expf(scale * x * x);
    }
}

// 贝塞尔函数I0(用于凯泽窗)
float32_t bessel_i0(float32_t x) {
    float32_t sum = 1.0f;
    float32_t term = 1.0f;
    float32_t half_x = x / 2.0f;
    
    for (uint16_t k = 1; k <= 20; k++) {
        term *= (half_x * half_x) / (k * k);
        sum += term;
    }
    
    return sum;
}

// 计算凯泽窗的beta参数
float32_t calculate_kaiser_beta(float32_t ripple_db) {
    if (ripple_db > 50) {
        return 0.1102f * (ripple_db - 8.7f);
    } else if (ripple_db >= 21) {
        return 0.5842f * powf(ripple_db - 21.0f, 0.4f) + 0.07886f * (ripple_db - 21.0f);
    } else {
        return 0.0f;
    }
}

// 计算滤波器阶数
uint16_t calculate_filter_order(float32_t transition_width, float32_t sampling_freq, window_type_t window) {
    float32_t normalized_transition = transition_width / sampling_freq;
    
    switch (window) {
        case WINDOW_RECTANGULAR:
            return (uint16_t)(1.8f / normalized_transition);
        case WINDOW_HANN:
            return (uint16_t)(6.2f / normalized_transition);
        case WINDOW_HAMMING:
            return (uint16_t)(6.6f / normalized_transition);
        case WINDOW_BLACKMAN:
            return (uint16_t)(11.0f / normalized_transition);
        case WINDOW_KAISER:
            return (uint16_t)(10.0f / normalized_transition);
        default:
            return (uint16_t)(6.2f / normalized_transition);
    }
}

11. 示例程序 (examples/fir_lowpass_example.c)

#include <stdio.h>
#include <stdlib.h>
#include "fir_filter.h"
#include "window_functions.h"

#define SAMPLE_RATE 48000.0f
#define CUTOFF_FREQ 5000.0f
#define NUM_TAPS 101

int main() {
    printf("FIR低通滤波器示例\n");
    printf("采样率: %.0f Hz\n", SAMPLE_RATE);
    printf("截止频率: %.0f Hz\n", CUTOFF_FREQ);
    printf("抽头数: %d\n\n", NUM_TAPS);
    
    // 创建滤波器
    fir_filter_t filter;
    fir_config_t config;
    
    // 设计滤波器
    if (fir_design_lowpass(&filter, CUTOFF_FREQ, SAMPLE_RATE, NUM_TAPS, WINDOW_HAMMING) != 0) {
        printf("滤波器设计失败!\n");
        return -1;
    }
    
    printf("滤波器设计成功!\n");
    printf("滤波器系数:\n");
    
    // 打印前10个系数
    for (int i = 0; i < 10; i++) {
        printf("h[%d] = %.6f\n", i, filter.coeffs[i]);
    }
    printf("...\n");
    
    // 测试信号:包含1kHz和10kHz的正弦波
    printf("\n测试信号处理:\n");
    printf("输入: 1kHz + 10kHz 正弦波\n");
    printf("期望: 只保留1kHz成分\n\n");
    
    float32_t input_signal[1000];
    float32_t output_signal[1000];
    
    // 生成测试信号
    for (int i = 0; i < 1000; i++) {
        float32_t t = (float32_t)i / SAMPLE_RATE;
        input_signal[i] = sinf(2.0f * PI_F * 1000.0f * t) +  // 1kHz
                        sinf(2.0f * PI_F * 10000.0f * t);   // 10kHz
    }
    
    // 处理信号
    fir_process_block(&filter, input_signal, output_signal, 1000);
    
    // 显示前10个输出样本
    printf("前10个输出样本:\n");
    for (int i = 0; i < 10; i++) {
        printf("out[%d] = %.6f\n", i, output_signal[i]);
    }
    
    // 计算频率响应
    printf("\n频率响应(幅度):\n");
    for (float32_t freq = 0; freq <= SAMPLE_RATE/2; freq += 1000) {
        float32_t response = 0.0f;
        for (int i = 0; i < filter.num_taps; i++) {
            float32_t phase = 2.0f * PI_F * freq * i / SAMPLE_RATE;
            response += filter.coeffs[i] * cosf(phase);
        }
        printf("%.0f Hz: %.4f\n", freq, response);
    }
    
    // 清理资源
    fir_deinit(&filter);
    
    return 0;
}

12. Makefile

CC = gcc
CFLAGS = -Wall -O2 -std=c99 -fPIC
AR = ar
ARFLAGS = rcs

# 目录
INCLUDE_DIR = include
SRC_DIR = src
EXAMPLES_DIR = examples
TESTS_DIR = tests
BUILD_DIR = build
LIB_DIR = lib

# 目标
LIBRARY = $(LIB_DIR)/libdigitalfilters.a
SHARED_LIB = $(LIB_DIR)/libdigitalfilters.so

# 源文件
SRCS = $(wildcard $(SRC_DIR)/*.c)
OBJS = $(SRCS:$(SRC_DIR)/%.c=$(BUILD_DIR)/%.o)

# 示例程序
EXAMPLES = fir_lowpass_example iir_bandpass_example lms_adaptive_example decimation_example

.PHONY: all clean directories library examples test

all: directories library examples

directories:
	mkdir -p $(BUILD_DIR) $(LIB_DIR)

library: $(LIBRARY) $(SHARED_LIB)

$(LIBRARY): $(OBJS)
	$(AR) $(ARFLAGS) $@ $^

$(SHARED_LIB): $(OBJS)
	$(CC) -shared -o $@ $^

$(BUILD_DIR)/%.o: $(SRC_DIR)/%.c
	$(CC) $(CFLAGS) -I$(INCLUDE_DIR) -c $< -o $@

examples: $(EXAMPLES)

fir_lowpass_example: $(BUILD_DIR)/fir_lowpass_example.o $(LIBRARY)
	$(CC) -o $@ $^ -lm

iir_bandpass_example: $(BUILD_DIR)/iir_bandpass_example.o $(LIBRARY)
	$(CC) -o $@ $^ -lm

lms_adaptive_example: $(BUILD_DIR)/lms_adaptive_example.o $(LIBRARY)
	$(CC) -o $@ $^ -lm

decimation_example: $(BUILD_DIR)/decimation_example.o $(LIBRARY)
	$(CC) -o $@ $^ -lm

$(BUILD_DIR)/%.o: $(EXAMPLES_DIR)/%.c
	$(CC) $(CFLAGS) -I$(INCLUDE_DIR) -c $< -o $@

test: $(TESTS_DIR)/test_fir.c $(LIBRARY)
	$(CC) $(CFLAGS) -I$(INCLUDE_DIR) -o test_fir $(TESTS_DIR)/test_fir.c $(LIBRARY) -lm
	./test_fir

clean:
	rm -rf $(BUILD_DIR) $(LIB_DIR) $(EXAMPLES) test_fir

install: library
	cp $(LIBRARY) /usr/local/lib/
	cp $(SHARED_LIB) /usr/local/lib/
	cp -r $(INCLUDE_DIR)/*.h /usr/local/include/

uninstall:
	rm -f /usr/local/lib/libdigitalfilters.*
	rm -f /usr/local/include/fir_filter.h
	rm -f /usr/local/include/iir_filter.h
	rm -f /usr/local/include/adaptive_filter.h
	rm -f /usr/local/include/multirate_filter.h
	rm -f /usr/local/include/window_functions.h
	rm -f /usr/local/include/math_utils.h

help:
	@echo "可用目标:"
	@echo "  all        - 编译所有内容"
	@echo "  library    - 编译静态和动态库"
	@echo "  examples   - 编译示例程序"
	@echo "  test       - 运行测试"
	@echo "  clean      - 清理编译文件"
	@echo "  install    - 安装库和头文件"
	@echo "  uninstall  - 卸载库和头文件"

参考代码 数字滤波器的C语言实现 www.youwenfan.com/contentcsv/112860.html

三、使用说明

1. 基本用法

#include "fir_filter.h"

int main() {
    // 1. 创建滤波器
    fir_filter_t filter;
    fir_config_t config;
    
    // 2. 设计低通滤波器
    fir_design_lowpass(&filter, 1000.0f, 8000.0f, 51, WINDOW_HAMMING);
    
    // 3. 处理信号
    float32_t input = 1.0f;
    float32_t output = fir_process(&filter, input);
    
    // 4. 清理
    fir_deinit(&filter);
    
    return 0;
}

2. 滤波器类型对比

滤波器类型 优点 缺点 适用场景
FIR 线性相位、稳定、设计简单 阶数高、延迟大 音频处理、图像处理
IIR 阶数低、效率高 非线性相位、可能不稳定 实时处理、嵌入式
自适应 自动跟踪变化 收敛慢、计算复杂 噪声消除、回声消除
多速率 改变采样率 设计复杂 采样率转换、多速率系统

3. 性能优化建议

  1. 使用对称FIR:减少50%乘法运算
  2. 直接形式II:减少内存使用
  3. 定点运算:使用Q格式提高嵌入式性能
  4. SIMD指令:使用NEON/SSE加速
  5. 查表法:预先计算三角函数等

4. 嵌入式移植

// 嵌入式版本(无动态内存分配)
static float32_t coeffs[51];
static float32_t buffer[51];
static fir_filter_t filter = {
    .coeffs = coeffs,
    .buffer = buffer,
    .num_taps = 51,
    .buffer_index = 0,
    .initialized = 1
};

// 初始化
void init_filter(void) {
    // 手动设置系数
    coeffs[0] = 0.001f;
    // ... 设置其他系数
}

// 处理
float32_t process_sample(float32_t input) {
    return fir_process(&filter, input);
}
Logo

智能硬件社区聚焦AI智能硬件技术生态,汇聚嵌入式AI、物联网硬件开发者,打造交流分享平台,同步全国赛事资讯、开展 OPC 核心人才招募,助力技术落地与开发者成长。

更多推荐