Matlab下各个滤波器传递函数透视图的绘制

来源:互联网 发布:部落冲突矿工升级数据 编辑:程序博客网 时间:2024/05/23 07:25

原文地址:http://blog.csdn.net/ljp1919/article/details/44781215

本文主要绘制常用滤波器的透视图,包括低通,带阻和带通滤波器。而这三类滤波器又各自包含了理想、巴特沃兹和高斯滤波器。如,低通滤波器就可以分为理想低通滤波器、n阶巴特沃兹滤波器和高斯低通滤波器。

第一:低通滤波器

1)理想低通滤波器

传递函数:


透视图结果:


代码:

[cpp] view plain copy
 在CODE上查看代码片派生到我的代码片
  1. a=100;  
  2. b=100;  
  3. U=0:a;  
  4. V=0:b;  
  5. M=length(U);N=length(V);  
  6. D0=10;%W=200;%D0是频带的中心半径;W是频带的宽度  
  7. x1=50;y1=50;  
  8. x0=-50;y0=-50;  
  9. m=fix(M/2); n=fix(N/2);  
  10. H=zeros(M,N);  
  11. n=8;  
  12. for u=1:M  
  13.     for v=1:N  
  14.      a=sqrt((U(u) - 50) .* (U(u)-50) + (V(v) - 50) .* (V(v) - 50));%D(u,v)的值  
  15.             if(a<=D0)%理想陷波器  
  16.                 H(u,v)=1;  
  17.             else  
  18.                 H(u,v)=0;  
  19.             end  
  20.     end  
  21. end  
  22. %在绘制高斯曲面的时候,加上下述代码,显示得美观  
  23. figure;  
  24. surf(U,V,H)  


2)n阶巴特沃兹低通滤波器

传递函数:


透视图结果:


代码:

[cpp] view plain copy
 在CODE上查看代码片派生到我的代码片
  1. a=100;  
  2. b=100;  
  3. U=0:a;  
  4. V=0:b;  
  5. M=length(U);N=length(V);  
  6. D0=10;%W=200;%D0是频带的中心半径;W是频带的宽度  
  7. x1=50;y1=50;  
  8. x0=-50;y0=-50;  
  9. m=fix(M/2); n=fix(N/2);  
  10. H=zeros(M,N);  
  11. n=8;  
  12. for u=1:M  
  13.     for v=1:N  
  14.       a=sqrt((U(u) - 50) .* (U(u)-50) + (V(v) - 50) .* (V(v) - 50));%D(u,v)的值  
  15.       b=1+(a/D0)^2*n;  
  16.       H(u,v)=1/b;%尝试下,非理想滤波器,如巴特沃兹或者高斯  
  17.     end  
  18. end  
  19. figure;  
  20. surf(U,V,H)  
  21. title('n=8')  



3)高斯低通滤波器

传递函数:


透视图结果:


代码:

[cpp] view plain copy
 在CODE上查看代码片派生到我的代码片
  1. a=100;  
  2. b=100;  
  3. U=0:a;  
  4. V=0:b;  
  5. M=length(U);N=length(V);  
  6. D0=10; %D0是频带的中心半径;W是频带的宽度  
  7. x1=50;y1=50;  
  8. x0=-50;y0=-50;  
  9. m=fix(M/2); n=fix(N/2);  
  10. H=zeros(M,N);  
  11.   
  12. for u=1:M  
  13.     for v=1:N  
  14.         D1=((u-m-x0)^2+(v-n-y0).^2)^0.5;  
  15.         D2=((u-m+x0)^2+(v-n+y0).^2)^0.5;  
  16.         D11=((u-m-x1)^2+(v-n-y1).^2)^0.5;  
  17.         D21=((u-m+x1)^2+(v-n+y1).^2)^0.5;  
  18.         %高斯低通曲面  
  19.         H(u,v) = (U(u) - 50) .* (U(u)-50) + (V(v) - 50) .* (V(v) - 50);  
  20.           
  21.     end  
  22. end  
  23. %在绘制高斯曲面的时候,加上下述代码,显示得美观  
  24. fangcha=50;  
  25. H = -H/(2*fangcha);  
  26. H = exp(H) / (sqrt(2*pi) * sqrt(fangcha));  
  27.   
  28. surf(U,V,H)  


第二:带阻滤波器

1)理想带阻滤波器

传递函数:


透视图结果:


代码:

[cpp] view plain copy
 在CODE上查看代码片派生到我的代码片
  1. % %绘制函数剖面线  
  2. % [u,v] = meshgrid(-400:400, -400:400);  
  3. a=100;%图像的尺寸,长  
  4. b=100;%图像的尺寸,宽  
  5. U=0:a;  
  6. V=0:b;  
  7. M=length(U);N=length(V);  
  8. D0=30;%W=200;%D0是频带的中心半径;W是频带的宽度  
  9. x1=50;y1=50;  
  10. x0=-50;y0=-50;  
  11. m=fix(M/2);n=fix(N/2);  
  12. H=zeros(M,N);  
  13. % n=8;  
  14. W=10;%W是频带的宽度  
  15. for u=1:M  
  16.     for v=1:N  
  17. %         D1=((u-m-x0)^2+(v-n-y0).^2)^0.5;  
  18. %         D2=((u-m+x0)^2+(v-n+y0).^2)^0.5;  
  19. %         D11=((u-m-x1)^2+(v-n-y1).^2)^0.5;  
  20. %         D21=((u-m+x1)^2+(v-n+y1).^2)^0.5;  
  21.   
  22. %         % %设计带阻滤波器  
  23. %         H(u,v)=1-exp(-0.5*((D(u,v)^2-D0^2)/(D(u,v)*W))^2);  
  24.             D=sqrt((U(u) - m) .* (U(u)-m) + (V(v) - n) .* (V(v) - n));%D(u,v)的值  
  25.             if((D<D0-W/2) || (D>D0+W/2))%理想带阻滤波器  
  26.                 H(u,v)=1;  
  27.             else  
  28.                 H(u,v)=0;  
  29.             end  
  30.   
  31. %         H(u,v)=1-exp(-0.5*((D1*D2/D0^2)));%尝试下,非理想滤波器,如巴特沃兹或者高斯  
  32. %       a=sqrt((U(u) - 50) .* (U(u)-50) + (V(v) - 50) .* (V(v) - 50));%D(u,v)的值  
  33. %       b=1+(a/D0)^2*n;  
  34. %       H(u,v)=1/b;%尝试下,非理想滤波器,如巴特沃兹或者高斯  
  35.     end  
  36. end  
  37. %在绘制高斯曲面的时候,加上下述代码,显示得美观  
  38. % fangcha=50;  
  39. % H = -H/(2*fangcha);  
  40. % H = exp(H) / (sqrt(2*pi) * sqrt(fangcha));  
  41. figure;  
  42. surf(U,V,H)  
  43. title('理想带阻滤波器:D0=30,W=10')  
  44. % [X,Y] = meshgrid(-2:.2:2, -2:.2:5);                                  
  45. % Z = X .* exp(-X.^2 - Y.^2);                                          
  46. % surf(X,Y,Z)  

2)n阶巴特沃兹带阻滤波器

传递函数:


透视图结果:


代码:

[cpp] view plain copy
 在CODE上查看代码片派生到我的代码片
  1. % %绘制函数剖面线  
  2. % [u,v] = meshgrid(-400:400, -400:400);  
  3. a=100;%图像的尺寸,长  
  4. b=100;%图像的尺寸,宽  
  5. U=0:a;  
  6. V=0:b;  
  7. M=length(U);N=length(V);  
  8. D0=30;%W=200;%D0是频带的中心半径;W是频带的宽度  
  9. x1=50;y1=50;  
  10. x0=-50;y0=-50;  
  11. m=fix(M/2);n=fix(N/2);  
  12. H=zeros(M,N);  
  13. n1=1;%巴特沃兹的阶数  
  14. W=10;%W是频带的宽度  
  15. for u=1:M  
  16.     for v=1:N  
  17. %         D1=((u-m-x0)^2+(v-n-y0).^2)^0.5;  
  18. %         D2=((u-m+x0)^2+(v-n+y0).^2)^0.5;  
  19. %         D11=((u-m-x1)^2+(v-n-y1).^2)^0.5;  
  20. %         D21=((u-m+x1)^2+(v-n+y1).^2)^0.5;  
  21.   
  22. %         % %设计带阻滤波器  
  23. %         H(u,v)=1-exp(-0.5*((D(u,v)^2-D0^2)/(D(u,v)*W))^2);  
  24.             D=sqrt((U(u) - m) .* (U(u)-m) + (V(v) - n) .* (V(v) - n));%D(u,v)的值  
  25. %             if((D<D0-W/2) || (D>D0+W/2))%理想带阻滤波器  
  26. %                 H(u,v)=1;  
  27. %             else  
  28. %                 H(u,v)=0;  
  29. %             end  
  30.             temp=1+(D*W/(D^2-D0^2))^2*n1;  
  31.             H(u,v)=1/temp;  
  32. %         H(u,v)=1-exp(-0.5*((D1*D2/D0^2)));%尝试下,非理想滤波器,如巴特沃兹或者高斯  
  33. %       a=sqrt((U(u) - 50) .* (U(u)-50) + (V(v) - 50) .* (V(v) - 50));%D(u,v)的值  
  34. %       b=1+(a/D0)^2*n;  
  35. %       H(u,v)=1/b;%尝试下,非理想滤波器,如巴特沃兹或者高斯  
  36.     end  
  37. end  
  38. %在绘制高斯曲面的时候,加上下述代码,显示得美观  
  39. % fangcha=50;  
  40. % H = -H/(2*fangcha);  
  41. % H = exp(H) / (sqrt(2*pi) * sqrt(fangcha));  
  42. figure;  
  43. surf(U,V,H)  
  44. title('1阶巴特沃兹带阻滤波器:D0=30,W=10')  
  45. % [X,Y] = meshgrid(-2:.2:2, -2:.2:5);                                  
  46. % Z = X .* exp(-X.^2 - Y.^2);                                          
  47. % surf(X,Y,Z)  
3)高斯带阻滤波器

传递函数:


透视图:


代码:

[cpp] view plain copy
 在CODE上查看代码片派生到我的代码片
  1. % %绘制函数剖面线  
  2. % [u,v] = meshgrid(-400:400, -400:400);  
  3. a=100;%图像的尺寸,长  
  4. b=100;%图像的尺寸,宽  
  5. U=0:a;  
  6. V=0:b;  
  7. M=length(U);N=length(V);  
  8. D0=30;%W=200;%D0是频带的中心半径;W是频带的宽度  
  9. x1=10;y1=10;  
  10. u0=-20;v0=20;%(u0,v0)的值对应移动中心  
  11. m=fix(M/2);n=fix(N/2);  
  12. H=zeros(M,N);  
  13. n1=2;%巴特沃兹的阶数  
  14. W=10;%W是频带的宽度  
  15. for u=1:M  
  16.     for v=1:N  
  17. %         D1=((u-m-u0)^2+(v-n-v0).^2)^0.5;%()  
  18. %         D2=((u-m+u0)^2+(v-n+v0).^2)^0.5;  
  19.           
  20. %          D1=((U(u)-m-u0)^2+(V(v)-n-v0)^2)^0.5;%()  
  21. %         D2=((U(u)-m+u0)^2+(V(v)-n+v0)^2)^0.5;  
  22.   
  23. %         D11=((u-m-x1)^2+(v-n-y1).^2)^0.5;  
  24. %         D21=((u-m+x1)^2+(v-n+y1).^2)^0.5;  
  25.   
  26. %         % %设计带阻滤波器  
  27. %         H(u,v)=1-exp(-0.5*((D(u,v)^2-D0^2)/(D(u,v)*W))^2);  
  28. %             D=sqrt((U(u) - m) .* (U(u)-m) + (V(v) - n) .* (V(v) - n));  
  29.             D=sqrt((U(u) - m) .* (U(u)-m) + (V(v) - n) .* (V(v) - n));%D(u,v)的值  
  30. %             if((D<D0-W/2) || (D>D0+W/2))%理想带阻滤波器  
  31. %                 H(u,v)=1;  
  32. %             else  
  33. %                 H(u,v)=0;  
  34. %             end  
  35.             temp=-0.5*((D^2-D0^2)/(D*W))^2;  
  36.             H(u,v)=1-exp(temp);  
  37. %         H(u,v)=1-exp(-0.5*((D1*D2/D0^2)));%尝试下,非理想滤波器,如巴特沃兹或者高斯  
  38. %       a=sqrt((U(u) - 50) .* (U(u)-50) + (V(v) - 50) .* (V(v) - 50));%D(u,v)的值  
  39. %       b=1+(a/D0)^2*n;  
  40. %       H(u,v)=1/b;%尝试下,非理想滤波器,如巴特沃兹或者高斯  
  41.     end  
  42. end  
  43. %在绘制高斯曲面的时候,加上下述代码,显示得美观  
  44. % fangcha=50;  
  45. % H = -H/(2*fangcha);  
  46. % H = exp(H) / (sqrt(2*pi) * sqrt(fangcha));  
  47. figure;  
  48. surf(U,V,H)  
  49. title('高斯带阻滤波器:D0=30,W=10')  


第三:陷波器

1)理想陷波器

传递函数:


透视图:


代码:

[cpp] view plain copy
 在CODE上查看代码片派生到我的代码片
  1. % %绘制函数剖面线  
  2. % [u,v] = meshgrid(-400:400, -400:400);  
  3. a=100;%图像的尺寸,长  
  4. b=100;%图像的尺寸,宽  
  5. U=0:a;  
  6. V=0:b;  
  7. M=length(U);N=length(V);  
  8. D0=10;%W=200;%D0是频带的中心半径;W是频带的宽度  
  9. x1=10;y1=10;  
  10. u0=-30;v0=30;%(u0,v0)的值对应移动中心  
  11. m=fix(M/2);n=fix(N/2);  
  12. H=zeros(M,N);  
  13. n1=2;%巴特沃兹的阶数  
  14. W=10;%W是频带的宽度  
  15. for u=1:M  
  16.     for v=1:N  
  17. %         D1=((u-m-u0)^2+(v-n-v0).^2)^0.5;%()  
  18. %         D2=((u-m+u0)^2+(v-n+v0).^2)^0.5;  
  19.           
  20.          D1=((U(u)-m-u0)^2+(V(v)-n-v0)^2)^0.5;%()  
  21.         D2=((U(u)-m+u0)^2+(V(v)-n+v0)^2)^0.5;  
  22.   
  23. %         D11=((u-m-x1)^2+(v-n-y1).^2)^0.5;  
  24. %         D21=((u-m+x1)^2+(v-n+y1).^2)^0.5;  
  25.   
  26. %         % %设计带阻滤波器  
  27. %         H(u,v)=1-exp(-0.5*((D(u,v)^2-D0^2)/(D(u,v)*W))^2);  
  28. %             D=sqrt((U(u) - m) .* (U(u)-m) + (V(v) - n) .* (V(v) - n));  
  29. %             D=sqrt((U(u) - m) .* (U(u)-m) + (V(v) - n) .* (V(v) - n));%D(u,v)的值  
  30.             if((D1<=D0) || (D2<=D0))%理想带阻滤波器  
  31.                 H(u,v)=0;  
  32.             else  
  33.                 H(u,v)=1;  
  34.             end  
  35. %             temp=-0.5*((D^2-D0^2)/(D*W))^2;  
  36. %             H(u,v)=1-exp(temp);  
  37. %         H(u,v)=1-exp(-0.5*((D1*D2/D0^2)));%尝试下,非理想滤波器,如巴特沃兹或者高斯  
  38. %       a=sqrt((U(u) - 50) .* (U(u)-50) + (V(v) - 50) .* (V(v) - 50));%D(u,v)的值  
  39. %       b=1+(a/D0)^2*n;  
  40. %       H(u,v)=1/b;%尝试下,非理想滤波器,如巴特沃兹或者高斯  
  41.     end  
  42. end  
  43. %在绘制高斯曲面的时候,加上下述代码,显示得美观  
  44. % fangcha=50;  
  45. % H = -H/(2*fangcha);  
  46. % H = exp(H) / (sqrt(2*pi) * sqrt(fangcha));  
  47. figure;  
  48. surf(U,V,H)  
  49. title('理想带阻滤波器:D0=10,u0=-30;v0=30')  
  50. % [X,Y] = meshgrid(-2:.2:2, -2:.2:5);                                  
  51. % Z = X .* exp(-X.^2 - Y.^2);                                          
  52. % surf(X,Y,Z)  


2)n阶巴特沃兹陷波器:

传递函数:


透视图:

代码:

[cpp] view plain copy
 在CODE上查看代码片派生到我的代码片
  1. a=100;%图像的尺寸,长  
  2. b=100;%图像的尺寸,宽  
  3. U=0:a;  
  4. V=0:b;  
  5. M=length(U);N=length(V);  
  6. D0=20;%W=200;%D0是频带的中心半径;W是频带的宽度  
  7. x1=10;y1=10;  
  8. u0=-30;v0=30;%(u0,v0)的值对应移动中心  
  9. m=fix(M/2);n=fix(N/2);  
  10. H=zeros(M,N);  
  11. n1=2;%巴特沃兹的阶数  
  12. W=10;%W是频带的宽度  
  13. for u=1:M  
  14.     for v=1:N  
  15.         D1=((U(u)-m-u0)^2+(V(v)-n-v0)^2)^0.5;  
  16.         D2=((U(u)-m+u0)^2+(V(v)-n+v0)^2)^0.5;  
  17.             temp=1+(D0^2/(D1*D2))^n1;  
  18.             H(u,v)=1/temp;  
  19.     end  
  20. end  
  21. figure;  
  22. surf(U,V,H)  
  23. title('2阶巴特沃兹陷波器:D0=20,u0=-30;v0=30')  


3)高斯陷波器

传递函数:


透视图:


代码:

[cpp] view plain copy
 在CODE上查看代码片派生到我的代码片
  1. % %绘制函数剖面线  
  2. % [u,v] = meshgrid(-400:400, -400:400);  
  3. a=100;%图像的尺寸,长  
  4. b=100;%图像的尺寸,宽  
  5. U=0:a;  
  6. V=0:b;  
  7. M=length(U);N=length(V);  
  8. D0=10;%W=200;%D0是频带的中心半径;W是频带的宽度  
  9. x1=10;y1=10;  
  10. u0=-20;v0=20;%(u0,v0)的值对应移动中心  
  11. m=fix(M/2);n=fix(N/2);  
  12. H=zeros(M,N);  
  13. n1=2;%巴特沃兹的阶数  
  14. W=10;%W是频带的宽度  
  15. for u=1:M  
  16.     for v=1:N  
  17. %         D1=((u-m-u0)^2+(v-n-v0).^2)^0.5;%()  
  18. %         D2=((u-m+u0)^2+(v-n+v0).^2)^0.5;  
  19.           
  20.          D1=((U(u)-m-u0)^2+(V(v)-n-v0)^2)^0.5;%()  
  21.         D2=((U(u)-m+u0)^2+(V(v)-n+v0)^2)^0.5;  
  22. %         D11=((u-m-x1)^2+(v-n-y1).^2)^0.5;  
  23. %         D21=((u-m+x1)^2+(v-n+y1).^2)^0.5;  
  24.   
  25. %         % %设计带阻滤波器  
  26. %         H(u,v)=1-exp(-0.5*((D(u,v)^2-D0^2)/(D(u,v)*W))^2);  
  27. %             D=sqrt((U(u) - m) .* (U(u)-m) + (V(v) - n) .* (V(v) - n));%D(u,v)的值  
  28. %             if((D<D0-W/2) || (D>D0+W/2))%理想带阻滤波器  
  29. %                 H(u,v)=1;  
  30. %             else  
  31. %                 H(u,v)=0;  
  32. %             end  
  33.             temp=-0.5*(D1*D2/D0^2);  
  34.             H(u,v)=1-exp(temp);  
  35. %         H(u,v)=1-exp(-0.5*((D1*D2/D0^2)));%尝试下,非理想滤波器,如巴特沃兹或者高斯  
  36. %       a=sqrt((U(u) - 50) .* (U(u)-50) + (V(v) - 50) .* (V(v) - 50));%D(u,v)的值  
  37. %       b=1+(a/D0)^2*n;  
  38. %       H(u,v)=1/b;%尝试下,非理想滤波器,如巴特沃兹或者高斯  
  39.     end  
  40. end  
  41. %在绘制高斯曲面的时候,加上下述代码,显示得美观  
  42. % fangcha=50;  
  43. % H = -H/(2*fangcha);  
  44. % H = exp(H) / (sqrt(2*pi) * sqrt(fangcha));  
  45. figure;  
  46. surf(U,V,H)  
  47. title('高斯陷波器:D0=10,u0=-20;v0=20')  
  48. % [X,Y] = meshgrid(-2:.2:2, -2:.2:5);                                  
  49. % Z = X .* exp(-X.^2 - Y.^2);                                          
  50. % surf(X,Y,Z)  
0 0
原创粉丝点击