概述
本文分享使用FIR滤波器设计低通滤波器，分享如何使用matlab和arm_dsp库设计滤波器。
正文
matlab设计滤波器系数
-
matlab输入fdatool调出滤波器设计工具
-
输入滤波器各项参数，如下：
-
点击Targets->Generate C header，并选择对应的数据类型即可导出文件。
ARM DSP库FIR滤波器API
这里提供自己整理的驱动文件作为参考，滤波效果还算可以。
#include "fir.h"
/**********************************USER TODO*************************************************/
#define NUM_TAPS 28 //滤波器阶数
#define BLOCK_SIZE 1 //数据处理块大小
//汉明窗 80Hz截止频率 250Hz采样率
static const float32_t firCoeffs32LP[NUM_TAPS + 1] = {
0.0002276871674f, 0.001889977371f, -0.002808483783f, -0.0006533082924f, 0.007655459922f,
-0.008231480606f, -0.006402248517f, 0.02447644062f, -0.01639692485f, -0.02764063515f,
0.0645493865f, -0.02370921522f, -0.116899237f, 0.2843400538f, 0.6392050385f,
0.2843400538f, -0.116899237f, -0.02370921522f, 0.0645493865f, -0.02764063515f,
-0.01639692485f, 0.02447644062f, -0.006402248517f, -0.008231480606f, 0.007655459922f,
-0.0006533082924f, -0.002808483783f, 0.001889977371f, 0.0002276871674f
};
/********************************************************************************************/
static arm_fir_instance_f32 S;
static float32_t firStateF32[BLOCK_SIZE + NUM_TAPS];
/*
@brief Initialization function for the floating-point FIR filter.
@param[in,out] S points to an instance of the floating-point FIR filter structure
@param[in] numTaps number of filter coefficients in the filter
@param[in] pCoeffs points to the filter coefficients buffer
@param[in] pState points to the state buffer
@param[in] blockSize number of samples processed per call
@return none
*/
static void arm_fir_init_f32(arm_fir_instance_f32 * S, uint16_t numTaps, const float32_t * pCoeffs, float32_t * pState, uint32_t blockSize)
{
/* Assign filter taps */
S->numTaps = numTaps;
/* Assign coefficient pointer */
S->pCoeffs = pCoeffs;
/* Clear state buffer. The size is always (blockSize + numTaps - 1) */
memset(pState, 0, (numTaps + (blockSize - 1U)) * sizeof(float32_t));
/* Assign state pointer */
S->pState = pState;
}
/*
@brief Processing function for the floating-point fir filter.
@param[in] S points to an instance of the floating-point fir structure
@param[in] pSrc points to the block of input data
@param[out] pDst points to the block of output data
@param[in] blockSize number of samples to process
@return none
*/
static void arm_fir_f32(const arm_fir_instance_f32 * S, const float32_t * pSrc, float32_t * pDst, uint32_t blockSize)
{
float32_t *pState = S->pState; /* State pointer */
const float32_t *pCoeffs = S->pCoeffs; /* Coefficient pointer */
float32_t *pStateCurnt; /* Points to the current sample of the state */
float32_t *px; /* Temporary pointer for state buffer */
const float32_t *pb; /* Temporary pointer for coefficient buffer */
float32_t acc0; /* Accumulator */
uint32_t numTaps = S->numTaps; /* Number of filter coefficients in the filter */
uint32_t i, tapCnt, blkCnt; /* Loop counters */
float32_t acc1, acc2, acc3, acc4, acc5, acc6, acc7; /* Accumulators */
float32_t x0, x1, x2, x3, x4, x5, x6, x7; /* Temporary variables to hold state values */
float32_t c0; /* Temporary variable to hold coefficient value */
/* S->pState points to state array which contains previous frame (numTaps - 1) samples */
/* pStateCurnt points to the location where the new input data should be written */
pStateCurnt = &(S->pState[(numTaps - 1U)]);
/* Loop unrolling: Compute 8 output values simultaneously.
* The variables acc0 ... acc7 hold output values that are being computed:
*
* acc0 = b[numTaps-1] * x[n-numTaps-1] + b[numTaps-2] * x[n-numTaps-2] + b[numTaps-3] * x[n-numTaps-3] +...+ b[0] * x[0]
* acc1 = b[numTaps-1] * x[n-numTaps] + b[numTaps-2] * x[n-numTaps-1] + b[numTaps-3] * x[n-numTaps-2] +...+ b[0] * x[1]
* acc2 = b[numTaps-1] * x[n-numTaps+1] + b[numTaps-2] * x[n-numTaps] + b[numTaps-3] * x[n-numTaps-1] +...+ b[0] * x[2]
* acc3 = b[numTaps-1] * x[n-numTaps+2] + b[numTaps-2] * x[n-numTaps+1] + b[numTaps-3] * x[n-numTaps] +...+ b[0] * x[3]
*/
blkCnt = blockSize >> 3U;
while (blkCnt > 0U)
{
/* Copy 4 new input samples into the state buffer. */
*pStateCurnt++ = *pSrc++;
*pStateCurnt++ = *pSrc++;
*pStateCurnt++ = *pSrc++;
*pStateCurnt++ = *pSrc++;
/* Set all accumulators to zero */
acc0 = 0.0f;
acc1 = 0.0f;
acc2 = 0.0f;
acc3 = 0.0f;
acc4 = 0.0f;
acc5 = 0.0f;
acc6 = 0.0f;
acc7 = 0.0f;
/* Initialize state pointer */
px = pState;
/* Initialize coefficient pointer */
pb = pCoeffs;
/* This is separated from the others to avoid
* a call to __aeabi_memmove which would be slower
*/
*pStateCurnt++ = *pSrc++;
*pStateCurnt++ = *pSrc++;
*pStateCurnt++ = *pSrc++;
*pStateCurnt++ = *pSrc++;
/* Read the first 7 samples from the state buffer: x[n-numTaps], x[n-numTaps-1], x[n-numTaps-2] */
x0 = *px++;
x1 = *px++;
x2 = *px++;
x3 = *px++;
x4 = *px++;
x5 = *px++;
x6 = *px++;
/* Loop unrolling: process 8 taps at a time. */
tapCnt = numTaps >> 3U;
while (tapCnt > 0U)
{
/* Read the b[numTaps-1] coefficient */
c0 = *(pb++);
/* Read x[n-numTaps-3] sample */
x7 = *(px++);
/* acc0 += b[numTaps-1] * x[n-numTaps] */
acc0 += x0 * c0;
/* acc1 += b[numTaps-1] * x[n-numTaps-1] */
acc1 += x1 * c0;
/* acc2 += b[numTaps-1] * x[n-numTaps-2] */
acc2 += x2 * c0;
/* acc3 += b[numTaps-1] * x[n-numTaps-3] */
acc3 += x3 * c0;
/* acc4 += b[numTaps-1] * x[n-numTaps-4] */
acc4 += x4 * c0;
/* acc1 += b[numTaps-1] * x[n-numTaps-5] */
acc5 += x5 * c0;
/* acc2 += b[numTaps-1] * x[n-numTaps-6] */
acc6 += x6 * c0;
/* acc3 += b[numTaps-1] * x[n-numTaps-7] */
acc7 += x7 * c0;
/* Read the b[numTaps-2] coefficient */
c0 = *(pb++);
/* Read x[n-numTaps-4] sample */
x0 = *(px++);
/* Perform the multiply-accumulate */
acc0 += x1 * c0;
acc1 += x2 * c0;
acc2 += x3 * c0;
acc3 += x4 * c0;
acc4 += x5 * c0;
acc5 += x6 * c0;
acc6 += x7 * c0;
acc7 += x0 * c0;
/* Read the b[numTaps-3] coefficient */
c0 = *(pb++);
/* Read x[n-numTaps-5] sample */
x1 = *(px++);
/* Perform the multiply-accumulates */
acc0 += x2 * c0;
acc1 += x3 * c0;
acc2 += x4 * c0;
acc3 += x5 * c0;
acc4 += x6 * c0;
acc5 += x7 * c0;
acc6 += x0 * c0;
acc7 += x1 * c0;
/* Read the b[numTaps-4] coefficient */
c0 = *(pb++);
/* Read x[n-numTaps-6] sample */
x2 = *(px++);
/* Perform the multiply-accumulates */
acc0 += x3 * c0;
acc1 += x4 * c0;
acc2 += x5 * c0;
acc3 += x6 * c0;
acc4 += x7 * c0;
acc5 += x0 * c0;
acc6 += x1 * c0;
acc7 += x2 * c0;
/* Read the b[numTaps-4] coefficient */
c0 = *(pb++);
/* Read x[n-numTaps-6] sample */
x3 = *(px++);
/* Perform the multiply-accumulates */
acc0 += x4 * c0;
acc1 += x5 * c0;
acc2 += x6 * c0;
acc3 += x7 * c0;
acc4 += x0 * c0;
acc5 += x1 * c0;
acc6 += x2 * c0;
acc7 += x3 * c0;
/* Read the b[numTaps-4] coefficient */
c0 = *(pb++);
/* Read x[n-numTaps-6] sample */
x4 = *(px++);
/* Perform the multiply-accumulates */
acc0 += x5 * c0;
acc1 += x6 * c0;
acc2 += x7 * c0;
acc3 += x0 * c0;
acc4 += x1 * c0;
acc5 += x2 * c0;
acc6 += x3 * c0;
acc7 += x4 * c0;
/* Read the b[numTaps-4] coefficient */
c0 = *(pb++);
/* Read x[n-numTaps-6] sample */
x5 = *(px++);
/* Perform the multiply-accumulates */
acc0 += x6 * c0;
acc1 += x7 * c0;
acc2 += x0 * c0;
acc3 += x1 * c0;
acc4 += x2 * c0;
acc5 += x3 * c0;
acc6 += x4 * c0;
acc7 += x5 * c0;
/* Read the b[numTaps-4] coefficient */
c0 = *(pb++);
/* Read x[n-numTaps-6] sample */
x6 = *(px++);
/* Perform the multiply-accumulates */
acc0 += x7 * c0;
acc1 += x0 * c0;
acc2 += x1 * c0;
acc3 += x2 * c0;
acc4 += x3 * c0;
acc5 += x4 * c0;
acc6 += x5 * c0;
acc7 += x6 * c0;
/* Decrement loop counter */
tapCnt--;
}
/* Loop unrolling: Compute remaining outputs */
tapCnt = numTaps % 0x8U;
while (tapCnt > 0U)
{
/* Read coefficients */
c0 = *(pb++);
/* Fetch 1 state variable */
x7 = *(px++);
/* Perform the multiply-accumulates */
acc0 += x0 * c0;
acc1 += x1 * c0;
acc2 += x2 * c0;
acc3 += x3 * c0;
acc4 += x4 * c0;
acc5 += x5 * c0;
acc6 += x6 * c0;
acc7 += x7 * c0;
/* Reuse the present sample states for next sample */
x0 = x1;
x1 = x2;
x2 = x3;
x3 = x4;
x4 = x5;
x5 = x6;
x6 = x7;
/* Decrement loop counter */
tapCnt--;
}
/* Advance the state pointer by 8 to process the next group of 8 samples */
pState = pState + 8;
/* The results in the 8 accumulators, store in the destination buffer. */
*pDst++ = acc0;
*pDst++ = acc1;
*pDst++ = acc2;
*pDst++ = acc3;
*pDst++ = acc4;
*pDst++ = acc5;
*pDst++ = acc6;
*pDst++ = acc7;
/* Decrement loop counter */
blkCnt--;
}
/* Loop unrolling: Compute remaining output samples */
blkCnt = blockSize % 0x8U;
while (blkCnt > 0U)
{
/* Copy one sample at a time into state buffer */
*pStateCurnt++ = *pSrc++;
/* Set the accumulator to zero */
acc0 = 0.0f;
/* Initialize state pointer */
px = pState;
/* Initialize Coefficient pointer */
pb = pCoeffs;
i = numTaps;
/* Perform the multiply-accumulates */
while (i > 0U)
{
/* acc = b[numTaps-1] * x[n-numTaps-1] + b[numTaps-2] * x[n-numTaps-2] + b[numTaps-3] * x[n-numTaps-3] +...+ b[0] * x[0] */
acc0 += *px++ * *pb++;
i--;
}
/* Store result in destination buffer. */
*pDst++ = acc0;
/* Advance state pointer by 1 for the next sample */
pState = pState + 1U;
/* Decrement loop counter */
blkCnt--;
}
/* Processing is complete.
Now copy the last numTaps - 1 samples to the start of the state buffer.
This prepares the state buffer for the next function call. */
/* Points to the start of the state buffer */
pStateCurnt = S->pState;
/* Loop unrolling: Compute 4 taps at a time */
tapCnt = (numTaps - 1U) >> 2U;
/* Copy data */
while (tapCnt > 0U)
{
*pStateCurnt++ = *pState++;
*pStateCurnt++ = *pState++;
*pStateCurnt++ = *pState++;
*pStateCurnt++ = *pState++;
/* Decrement loop counter */
tapCnt--;
}
/* Calculate remaining number of copies */
tapCnt = (numTaps - 1U) % 0x4U;
/* Copy remaining data */
while (tapCnt > 0U)
{
*pStateCurnt++ = *pState++;
/* Decrement loop counter */
tapCnt--;
}
}
void arm_fir_f32_init(void)
{
arm_fir_init_f32(&S, NUM_TAPS, (float32_t *)&firCoeffs32LP[0], &firStateF32[0], BLOCK_SIZE);
}
void arm_fir_f32_process(float32_t *inputF32, float32_t *outputF32)
{
arm_fir_f32(&S, inputF32, outputF32, BLOCK_SIZE);
}
#ifndef __FIR_H
#define __FIR_H
#include <stdio.h>
#include <stdint.h>
#include <string.h>
typedef float float32_t;
typedef struct
{
uint16_t numTaps; /**< number of filter coefficients in the filter. */
float32_t *pState; /**< points to the state variable array. The array is of length numTaps+blockSize-1. */
const float32_t *pCoeffs; /**< points to the coefficient array. The array is of length numTaps. */
} arm_fir_instance_f32;
void arm_fir_f32_init(void);
void arm_fir_f32_process(float32_t *inputF32, float32_t *outputF32);
#endif
效果
经过matlab将频谱图画出，可以看到50Hz和100Hz的混合信号经过低通滤波后，100Hz噪声基本完全滤除，说明滤波效果还可以。
问：为什么时域信号图上，滤波信号一开始基本没有？
答：因为波形经过 FIR 滤波器后，输出的波形会有一定的延迟。对于线性相位的 FIR，这个群延迟就是一个常数。但是实际应用中这个群延迟是多少呢？ 关于群延迟的数值，滤波器设计工具箱会根据用户的配置计算好。如下：
此实验的滤波器群延迟就为14，因此滤波输出的第15个点才是实际波形的起始点，这个需要在实际项目中特别注意。
总结
高通、带通等其他类型的滤波器设计大体相同，重点是通过matlab输出系数矩阵即可。本文仅提供简单的滤波器设计，实际项目中需要根据需求进行调试优化。
