单位冲击响应与频响

news/2024/11/22 22:28:36/

[数字信号处理]单位冲击响应与频响以及FIR实现代码(C语言)

分类: 数字信号处理 993人阅读 评论(6) 收藏 举报
解卷绕 FIR 滤波器 数字信号处理 C语言

目录(?)[+]

1.单位冲击响应与频响

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


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


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

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


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


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


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

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

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


2.窗函数实现的FIR滤波器代码(C语言)

[cpp] view plain copy
  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.频响的问题

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

   

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

        

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

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

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


4.实际的滤波效果

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


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


      其输出波形是如下。

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

       博客地址: http://blog.csdn.net/thnh169/

http://www.ppmy.cn/news/475649.html

相关文章

带你读懂——频率响应与采样频率之间的关系

频响范围 频率响应&#xff1a;不同频率下的输入信号经过系统后响应之后的输出信号增益。大白话就是&#xff0c;输入信号频率是 x x x Hz&#xff0c;幅值为 y y y mg&#xff0c;观察此时的输出信号幅值为 A y Ay Ay mg&#xff0c;此时升高或降低了 A A A倍。 电压增益计算…

EQ频响曲线绘制和DRC特性曲线绘制

目录 1 EQ 系数计算和频响曲线绘制 1.1 基本流程 1.2 EQ参数输入 1.3 滤波器系数计算 1.4 频率响应计算 1.5 曲线绘制 2.DRC特性曲线绘制 2.1 基本流程 2.2 参数输入 2.3 增益计算 2.4 静态特性曲线绘制 3 Android 工程实现 1 EQ 系数计算和频响曲线绘制 1.1 基本…

机械阻抗法与频响分析

目录 1. 什么是频响分析2. 机械阻抗和导纳的概念3. 集中参数元件的阻抗和导纳3.1 阻尼元件的阻抗和导纳3.2 机械阻抗网络图的建立3.3 机械阻抗的串并联计算 4. 单自由度紫铜的频响特性和导纳曲线4.1 无阻尼系统4.1.1 自由系统4.2 约束系统 4.2 有阻尼系统4.2.1 自由系统4.2.2 约…

LOTO示波器 实测 开环增益频响曲线/电源环路响应稳定性

一般我们用的电源系统/控制系统或者信号处理系统都可以简单理解成负反馈控制系统。最典型的&#xff0c;运放组成的信号放大电路就是这样的系统。本文以最简单的运放信号放大电路为例&#xff0c;演示如何使用LOTO示波器测量控制系统的开环增益频响曲线&#xff0c;以及演示电源…

ANSYS APDL谐响应分析——悬臂梁的频响函数计算以及幅值、角度(相位)、分贝计算

问题描述 研究一根悬臂梁&#xff0c;材质为钢材。长度 L 2 L2 L2 米&#xff1b;截面为矩形&#xff0c;矩形的长度为 H 5 c m H 5cm H5cm&#xff0c;宽度为 B 2 c m B 2cm B2cm 。 建模思路&#xff1a; 先建立节点&#xff0c;然后用节点生成单元。使用n命令&…

什么是冲激函数、时域卷积、冲激响应以及频响曲线

1&#xff0c; 线性时不变系统和非线性时变系统 1.1 线性时不变系统 线性时不变系统就是同时满足线性系统和时不变系统。 1. 线性系统 系统的输入输出之间满足线性叠加原理的系统称为线性系统。即 y 1 ( n ) T [ a 1 x 1 ( n ) ] , y 2 ( n ) T [ a 2 x 2 ( n ) ] y_1(n…

LOTO示波器如何测试阻抗的频响曲线

LOTO示波器如何测试阻抗的频响曲线 模块的输入输出端口&#xff0c;在电路分析上&#xff0c;一般简单表征为电阻来进行计算和分析。但多数情况下&#xff0c;这些端口并不是纯电阻的特性&#xff0c;更精确一些&#xff0c;它可能是电阻电容以及电感的组合&#xff0c;表现为非…

汽车振动响应分析-频响函数法(附程序)

文章目录 汽车的四自由度振动模型频响函数法MATLAB程序 汽车的四自由度振动模型 m&#xff1a;车身的等效质量&#xff1b; m_1&#xff1a;前轴的等效质量(包含电机的质量在内)&#xff1b; m_2&#xff1a;后轴的等效质量&#xff1b; m_3&#xff1a;人体和座椅的等效质量&a…