FIR数字滤波器C语言

你的名字 2022-10-16 05:48 192阅读 0赞

1.单位冲击响应与频响

就如同之前所说的一样，使用下图所示的单位冲击响应，所设计的滤波器，是无法实现的。

![20130529121349327][]

现在，让我们看看其这个滤波器的频响。所谓频响，就是计算其单位冲击响应的离散时间傅里叶变换，

![20130529122006711][]

我们可以看出，这个滤波器的频响的计算结果是实数，并没有虚数部分。也就是，其相位谱一直是0，也就意味着，这个滤波器输入与输出之间没有延迟，这种相位特性称为0延迟相位特性。

但是，这个滤波器无法是无法实现的。我们实际计算一下该滤波器的输入输出，就可以明白。

![20130529123305884][]

这个滤波器在计算的过程中，需要过去的值和未来的值。未来的值是不可预测的，所以，这个滤波器无法实现。为了使得这个滤波器可以实现，我们只好移动其单位冲击响应，使得其不再需要未来的值。比如，就像下面这样移动。

![20130529123551218][]

这样的话，这个滤波器就可以实现了，我们只需要记录其40个过去的输入值和现在的输入值。但是，由于移动了其单位冲击响应，会不会对频响产生什么影响呢，我们来看看。

![20130529123858607][]

为了更好的说明问题，L去代替之前例子中的20。

移动之后频响，我们根据上面式子可以得出两个结论：1，移动不会对幅度谱产生影响。2,，移动会对相位产生一个延迟，这个延迟主要取决于移动的长度，移动的长度越长，延迟越大。但是，这个移动是线性的。

因此，我们把这个移动的相位特性称为，线性相位特性。到这里，我们移动后的，因果的，可实现的滤波器的单位冲击响应，如下所示。

![20130529124432627][]

### 2.窗函数实现的FIR滤波器代码（C语言） ###

**\[cpp\]** [ view plain][view plain] [copy][view plain]

1.  \#include <stdio.h>
2.  \#include <math.h>
3.  \#include <malloc.h>
4.  \#include <string.h>
5.  
6.  
7.  \#define pi (3.1415926)
8.  
9.  /\*-------------Win Type----------------\*/
10. \#define Hamming (1)
11. 
12. 
13. 
14. double Input\_Data\[\] =
15. \{
16. 0.000000 , 0.896802 , 1.538842 , 1.760074 , 1.538842 , 1.000000 , 0.363271 , -0.142040 , -0.363271 , -0.278768,
17. 0.000000 , 0.278768 , 0.363271 , 0.142020 , -0.363271 , -1.000000 , -1.538842 , -1.760074 , -1.538842 , -0.896802,
18. 0.000000 , 0.896802 , 1.538842 , 1.760074 , 1.538842 , 1.000000 , 0.363271 , -0.142040 , -0.363271 , -0.278768,
19. 0.000000 , 0.278768 , 0.363271 , 0.142020 , -0.363271 , -1.000000 , -1.538842 , -1.760074 , -1.538842 , -0.896802,
20. 0.000000 , 0.896802 , 1.538842 , 1.760074 , 1.538842 , 1.000000 , 0.363271 , -0.142040 , -0.363271 , -0.278768,
21. 0.000000 , 0.278768 , 0.363271 , 0.142020 , -0.363271 , -1.000000 , -1.538842 , -1.760074 , -1.538842 , -0.896802,
22. 0.000000 , 0.896802 , 1.538842 , 1.760074 , 1.538842 , 1.000000 , 0.363271 , -0.142040 , -0.363271 , -0.278768,
23. 0.000000 , 0.278768 , 0.363271 , 0.142020 , -0.363271 , -1.000000 , -1.538842 , -1.760074 , -1.538842 , -0.896802,
24. 0.000000 , 55
25. \};
26. 
27. 
28. 
29. 
30. double sinc(double n)
31. \{
32. if(0==n) return (double)1;
33. else return (double)(sin(pi\*n)/(pi\*n));
34. \}
35. 
36. int Unit\_Impulse\_Response(int N,double w\_c,
37. int Win\_Type,
38. double \*Output\_Data)
39. \{
40. signed int Count = 0;
41. 
42. for(Count = -(N-1)/2;Count <= (N-1)/2;Count++)
43. \{
44. \*(Output\_Data+Count+((N-1)/2)) = (w\_c/pi)\*sinc((w\_c/pi)\*(double)(Count));
45. //printf("%d %lf ",Count+((N-1)/2)+1,\*(Output\_Data+Count+((N-1)/2)));
46. //if(Count%4 == 0) printf("\\n");
47. \}
48. 
49. 
50. switch (Win\_Type)
51. \{
52. 
53. case Hamming: printf("Hamming \\n");
54. for(Count = -(N-1)/2;Count <= (N-1)/2;Count++)
55. \{
56. \*(Output\_Data+Count+((N-1)/2)) \*= (0.54 + 0.46 \* cos((2\*pi\*Count)/(N-1)));
57. //printf("%d %lf ",Count+((N-1)/2)+1,\*(Output\_Data+Count+((N-1)/2)));
58. //if(((Count+1)%5 == 0)&&(Count != -20)) printf("\\n");
59. \}
60. break;
61. 
62. 
63. default: printf("default Hamming \\n");
64. for(Count = -(N-1)/2;Count <= (N-1)/2;Count++)
65. \{
66. \*(Output\_Data+Count+((N-1)/2)) \*= (0.54 + 0.46 \* cos((2\*pi\*Count)/(N-1)));
67. //printf("%d %lf ",Count+((N-1)/2)+1,\*(Output\_Data+Count+((N-1)/2)));
68. //if(((Count+1)%5 == 0)&&(Count != -20)) printf("\\n");
69. \}
70. 
71. break;
72. \}
73. 
74. return (int)1;
75. \}
76. 
77. 
78. void Save\_Input\_Date (double Scand,
79. int Depth,
80. double \*Input\_Data)
81. \{
82. int Count;
83. 
84. for(Count = 0 ; Count < Depth-1 ; Count++)
85. \{
86. \*(Input\_Data + Count) = \*(Input\_Data + Count + 1);
87. \}
88. 
89. \*(Input\_Data + Depth-1) = Scand;
90. \}
91. 
92. 
93. 
94. double Real\_Time\_FIR\_Filter(double \*b,
95. int b\_Lenth,
96. double \*Input\_Data)
97. \{
98. int Count;
99. double Output\_Data = 0;
100. 
101. Input\_Data += b\_Lenth - 1;
102. 
103. for(Count = 0; Count < b\_Lenth ;Count++)
104. \{
105. Output\_Data += (\*(b + Count)) \*
106. (\*(Input\_Data - Count));
107. \}
108. 
109. return (double)Output\_Data;
110. \}
111. 
112. 
113. 
114. 
115. 
116. int main(void)
117. \{
118. double w\_p = pi/10;
119. double w\_s = pi/5;
120. double w\_c = (w\_s + w\_p)/2;
121. printf("w\_c = %f \\n" , w\_c);
122. 
123. int N = 0;
124. N = (int) ((6.6\*pi)/(w\_s - w\_p) + 0.5);
125. if(0 == N%2) N++;
126. printf("N = %d \\n" , N);
127. 
128. double \*Impulse\_Response;
129. Impulse\_Response = (double \*) malloc(sizeof(double)\*N);
130. memset(Impulse\_Response,
131. 0,
132. sizeof(double)\*N);
133. 
134. Unit\_Impulse\_Response(N,w\_c,
135. Hamming,
136. Impulse\_Response);
137. 
138. double \*Input\_Buffer;
139. double Output\_Data = 0;
140. Input\_Buffer = (double \*) malloc(sizeof(double)\*N);
141. memset(Input\_Buffer,
142. 0,
143. sizeof(double)\*N);
144. int Count = 0;
145. 
146. FILE \*fs;
147. fs=fopen("FIR\_Data.txt","w");
148. 
149. while(1)
150. \{
151. if(Input\_Data\[Count\] == 55) break;
152. 
153. Save\_Input\_Date (Input\_Data\[Count\],
154. N,
155. Input\_Buffer);
156. 
157. Output\_Data = Real\_Time\_FIR\_Filter(Impulse\_Response,
158. N,
159. Input\_Buffer);
160. 
161. 
162. fprintf(fs,"%lf,",Output\_Data);
163. //if(((Count+1)%5 == 0)&&(Count != 0)) fprintf(fs,"\\r\\n"); 
164. 
165. Count++;
166. \}
167. 
168. /\*---------------------display-------------------------------- 
169.  for(Count = 0; Count < N;Count++)
170.  \{ 
171.  printf("%d %lf ",Count,\*(Input\_Buffer+Count));
172.  if(((Count+1)%5 == 0)&&(Count != 0)) printf("\\n"); 
173.  \} 
174.  \*/
175. 
176. fclose(fs);
177. printf("Finish \\n");
178. return (int)0;
179. \}

### 3.频响的问题 ###

按照上面程序，参数如下设定。

![20130529125948369][]

运行程序，我们就实现了一个FIR滤波器。我们使用Matlab做出其频响。

![20130529125322923][]

好了，这里可以看出，从其幅度特性看来，我们确实实现了一个低通滤波器。但是，相位特性就比较奇怪（为了方便看出问题，我已经进行了解卷绕，至于什么是解卷绕，为什么要解卷绕，之后会说）。

那么，问题来了！按照道理来说，这个FIR滤波器应该是拥有线性相位特性的，但是为什么这里的线性相位特性确不是一条直线！在接近于阻带起始频率的地方，有什么会震荡?这个问题，之后再解决吧。

[ \[数字信号处理\]相位特性解卷绕][Link 1] <-----------戳我

### 4.实际的滤波效果 ###

为了验证效果，我们输入实际的信号看看。这里，我们选择采样周期为10ms，那么，其通带频率和阻带起始频率为

![20130529131705108][]

为了验证其性质，我选择了1KHz和3KHz的频率混合，最终输出。输入的波形如下。

![20130529131708174][]

其输出波形是如下。

![20130529131711159][]

红色的“+”是Matlab计算的结果，黑的o是我用C语言实现的滤波器的计算结果，蓝的\*号是1KHz的信号，也就是希望的输出。可以看出，这个滤波器有一定的延迟，但是滤波效果还是不错。

from： [http://blog.csdn.net/thnh169/][http_blog.csdn.net_thnh169]

[20130529121349327]: /images/20221014/e451e82da5894731a67988cd3f0c4d7d.png
[20130529122006711]: /images/20221014/3fa0bf19bbf54015b93c7522eb6e7e1a.png
[20130529123305884]: https://img-blog.csdn.net/20130529123305884
[20130529123551218]: https://img-blog.csdn.net/20130529123551218
[20130529123858607]: https://img-blog.csdn.net/20130529123858607
[20130529124432627]: https://img-blog.csdn.net/20130529124432627
[view plain]: http://blog.csdn.net/thnh169/article/details/8988740#
[20130529125948369]: https://img-blog.csdn.net/20130529125948369
[20130529125322923]: https://img-blog.csdn.net/20130529125322923
[Link 1]: http://blog.csdn.net/thnh169/article/details/9007125
[20130529131705108]: https://img-blog.csdn.net/20130529131705108
[20130529131708174]: https://img-blog.csdn.net/20130529131708174
[20130529131711159]: https://img-blog.csdn.net/20130529131711159
[http_blog.csdn.net_thnh169]: http://blog.csdn.net/thnh169/