双二阶IIR滤波器——直接I型、直接II型C语言实现

双二阶滤波器

双二阶滤波器的生成

双二阶滤波器可以利用MATLAB的滤波器设计工具箱生成，根据目的，在Filter Designer界面中设置滤波器类型、设计方法、阶数、采样频率、截止频率等参数，如下所示：

在“当前滤波器信息”一栏中可以看到当前滤波器的结构信息、阶数、二阶节数、稳定性，当前界面显示的就是直接I型IIR滤波器，阶数为4，所以共有两个二阶节，且是稳定的。

滤波器结构转换

直接I型和直接II型滤波器可相互转换，具体步骤如下：
Edit->Convert Structure->选择期望的滤波器结构类型

可选项中的前两项分别为直接I型和直接II型，点击OK即可完成结构转换。

导出

得到期望的滤波器之后就可以导出了，为了得到双二阶滤波器中每个二阶节的参数，导出步骤如下所示：
File->Export->Export to中选择Coefficient File(ASCII)->Options中选择Decimal，点击Export后选择存放位置，保存后可得包含双二阶滤波器所有信息的fcf文件，如下所示：

% Generated by MATLAB(R) 9.6 and Signal Processing Toolbox 8.2.
% Generated on: 18-Aug-2021 16:19:34

% Coefficient Format: Decimal

% Discrete-Time IIR Filter (real)                                                                
% -------------------------------                                                                
% Filter Structure    : Direct-Form I, Second-Order Sections                                     
% Number of Sections  : 2                                                                        
% Stable              : Yes                                                                      
% Linear Phase        : No                                                                       

                                                                                                
SOS Matrix:                                                                                      
1  2  1  1  -1.637122299309369921616053034085780382156  0.756463440348990046047106261539738625288
1  2  1  1  -1.396616129407874051793214675853960216045  0.498425116157480918577249440204468555748
                                                                                                 
Scale Values:                                                                                    
0.029835285259905086618914538121316581964                                                        
0.02545224668740174792103125867015478434   

差分方程为：

与SOS Matrix的参数对应关系依次为：
b0 b1 b2 a0 a1 a2

C语言实现

直接I型

float filterDirectI(float DataIn, float *pFilterCoefs, float *pScaleVals, float *xDelay,float *yDelay, int numStages)
{
	float a1, a2;
	float b0, b1, b2;
	float y, input, output;
	float scaleVal;
	int indexStage;
	int numCoefsStage;
	int indexTemp;
	float temp[2];

	numCoefsStage = 6;
	input = DataIn;
	for (indexStage = 0; indexStage < numStages; indexStage++)
	{
		indexTemp = numCoefsStage * indexStage;
		// y[n] = b0*x[n] + b1*x[n-1] + b2*x[n-2] - a1*y[n-1] - a2*y[n-2]
		b0 = pFilterCoefs[indexTemp];
		b1 = pFilterCoefs[indexTemp + 1];
		b2 = pFilterCoefs[indexTemp + 2];
		a1 = pFilterCoefs[indexTemp + 4];
		a2 = pFilterCoefs[indexTemp + 5];
		// 两个系数
		scaleVal = pScaleVals[indexStage];
		//直接I型差分方程实现形式：
		//每个二阶节都需要四个延时因子，分为输入延时因子和输出延时因子
		//可以直接合为一个式子
		temp[indexStage] = scaleVal * (b0 * input + b1 * xDelay[indexStage * 2] + b2 * xDelay[indexStage * 2 + 1]);
		y = temp[indexStage] - a1 * yDelay[indexStage * 2] - a2 * yDelay[indexStage * 2 + 1];
 
		//延时因子的移位
		xDelay[indexStage * 2 + 1] = xDelay[indexStage * 2];
		xDelay[indexStage * 2] = input;
		yDelay[indexStage * 2 + 1] = yDelay[indexStage * 2];
		yDelay[indexStage * 2] = y;

		//进入下一个二阶节需要重新赋值，即下一个二阶节的输入是当前二阶节的输出
		input = y;
	}
	output = y;
	return output;
}

注：只考虑了两个二阶节的情况，更高阶的未验证。

直接II型

float  filter_IIR_BiquadCascade(float DataIn, float *pFilterCoefs, float *pScaleVals, float *pDelay, uint16_t numStages)
{
	/* Input:
	 * DataIn：                  输入数据（x[n]）
	 * pFilterCoefs:滤波器系数矩阵，一共有6*numStage 个数据
	 * pScaleVals：      幅度因子，一共有numStages个数据，每个二阶节滤波完成后乘以相应位置上的Scale
	 * pDelay：时延因子，用于存放滤波过程中的中间变量
	 * numStages:二阶节的个数 */
	float a1, a2;
	float b0, b1, b2;
	float y, input, output;
	float scaleVal;
	uint16_t indexStage;
	uint16_t numCoefsStage;
	uint16_t indexTemp;

	numCoefsStage = 6;
	input = DataIn;
	for (indexStage = 0; indexStage < numStages; indexStage++)
	{
		indexTemp = numCoefsStage * indexStage;
		// y[n] = b0*x[n] + b1*x[n-1] + b2*x[n-2] + a1*y[n-1] + a2*y[n-2]
		b0 = pFilterCoefs[indexTemp];
		b1 = pFilterCoefs[indexTemp + 1];
		b2 = pFilterCoefs[indexTemp + 2];
		a1 = pFilterCoefs[indexTemp + 4];
		a2 = pFilterCoefs[indexTemp + 5];
		// 两个系数
		scaleVal = pScaleVals[indexStage];

		pDelay[indexTemp] = scaleVal * input - a1 * pDelay[indexTemp + 1] - a2 * pDelay[indexTemp + 2];
		y = b0 * pDelay[indexTemp] + b1 * pDelay[indexTemp + 1] + b2 * pDelay[indexTemp + 2];

		pDelay[indexTemp + 2] = pDelay[indexTemp + 1];
		pDelay[indexTemp + 1] = pDelay[indexTemp];

		input = y;
	}
	output = y;
	return output;

}

经过与MATLAB filter函数对比验证，函数能够实现两种IIR滤波器的滤波功能。