windo.m

很多matlab的源代码

function w = windo(type,n,a)
% WINDO Generates symmetric windows for digital filter design
%
% 	W = WINDO(TY,N) generates an N-point window of type TY given by
% 	'rect'angle ('boxc'ar), 'bart'lett, 'feje'r, 'blac'kman, 'hamm'ing,
%	'vonh'ann,'par1','par2' (parzen) 'papo'ulis, 'riem'ann, and 'welc'h.
%       W = WINDO(TY,N,A) generates windows requiring additional parameter A 
%       'cosn' or 'lanc'zos   		(A=power, [default A=1])
%       'bell'                		(tapered, taper A<=1, [default A=0.5])
%       'kais'er 'tayl'or or 'mkais'	(beta =A, [default A=2]) 
%       'tuke'y               		(Tapered,taper A<=1, [default A=0.5])
%       'dolp'h ('cheb'yshev) 		(A=ripple dB [default A=60])
%
% TY = 'harris' (with A = 0 to 7) for peak sidelobe level and decay given by
%  0: pk SL: -46.74dB, decay: 30dB/oct	 1: pk SL: -60.95dB, decay: 42dB/oct
%  2: pk SL: -64.19dB, decay: 18dB/oct	 3: pk SL: -70.82dB, decay: 06dB/oct
%  4: pk SL: -71.48dB, decay: 06dB/oct	 5: pk SL: -82.60dB, decay: 30dB/oct 
%  6: pk SL: -93.32dB, decay: 18dB/oct	 7: pk SL: -98.17dB, decay: 06dB/oct 
%
%  NOTE: For windows to be used with the FFT, use WINFFT
%
%        WINDOW (with no input arguments) invokes the following example:
%
%        % Generate and plot a 15 point Dolph window with Ripple=40 dB
%        >>wk=windo('dolph',15,40),dtplot(0:14,wk,'o')


% ADSP Toolbox: Version 2.0 
% For use with "Analog and Digital Signal Processing", 2nd Ed.
% Published by PWS Publishing Co.
%
% Ashok Ambardar, EE Dept. MTU, Houghton, MI 49931, USA
% http://www.ee.mtu/faculty/akambard.html
% e-mail: [e-mail]
% Copyright (c) 1998



if nargin==0,help windo,disp('Strike a key to see results of the example')
pause,nk=windo('dolph',15,40);
 vx=matverch;
 if vx < 4, eval('clg');else,eval('clf');end
dtplot(0:14,nk,'o'),return,end

n1=n-1;nn=fix(n1/2);n0=nn:-1:0;if rem(n,2)==0,n0=n0+.5;end,t=0:nn;m=2*pi*t/n1;
ty=type(1:4);
if ty == 'rect' | ty=='boxc',w=ones(1,n);return,end
if ty == 'hamm',w = .54 - .46*cos(m);
elseif ty == 'hann'|ty=='vonh',w = .5*(1 - cos(m));
elseif ty == 'cosn' ,if nargin<3,a=1;end
b=sin(m/2);w=ones(1,1:nn+1);for j=1:a,w=b.*w;end
elseif ty=='par2' ,w = 1-2*abs((t-0.5*n1)/(n+1));
elseif ty=='blac',w=.42-.5*cos(m)+.08*cos(2*m);
elseif (ty=='bart' | ty=='feje'),w = 2*t/n1;
elseif ty=='papo' ,t1 = 2*n0/n1;w=abs(sin(pi*t1))/pi+(1-abs(t1)).*cos(pi*t1);
elseif ty=='welc' ,w=1-4*((t-0.5*n1)/(n+1)).^2;
elseif ty=='riem' ,t1=2*pi*n0/n1;w =(sin(t1)+eps)./(t1+eps);
elseif ty=='lanc' ,if nargin<3,a=1;end
t1=2*pi*n0/n1;w =(sin(t1)+eps)./(t1+eps);w=(w).^a;
elseif ty=='harr' ,if nargin<3,a=0;end
if a==0,b=[.375 .5 .125 0];end
if a==1,b=[10/32 15/32 6/32 1/32];end
if a==2,b=[.40897 .5 .09103 0];end
if a==3,b=[.42323 .49755 .07922 0];end
if a==4,b=[.4243801 .4973406 .0782793 0];end
if a==5,b=[.338946 .481973 .161054 .018027];end
if a==6,b=[.355768 .487396 .144232 .012604];end
if a==7,b=[.3635819 .4891775 .1365995 .0106411];end
w = b(1)-b(2)*cos(m)+b(3)*cos(2*m)-b(4)*cos(3*m);
elseif ty=='bell' ,if nargin<3,a=0.5;end
w = 0.5*(1-cos(2*pi*t/n1/a)).*(t <= .5*a*n1)+(t > .5*a*n1);
elseif ty=='tuck' ,if nargin<3,a=.5;end
a=1-a;if a==1,w=ones(1:nn+1);else,
w=0.5*(1+cos(pi*(n0-a*nn)/(1-a)/nn)).*(n0 > a*nn)+(n0 <= a*nn);end
elseif ty=='par1' ,z=1-2*t/n1;
w =(2*(1-z).^3).*(t <= 0.25*n1)+(1-6*z.^2+6*z.^3).*(t > 0.25*n1);
elseif (ty=='kais' | ty=='tayl' | ty=='mkai'),if nargin<3,a=2;end
a=a*pi;nn = fix((n-1)/2);n0=nn:-1:0;if rem(n,2)==0,n0=n0+.5;end
beta=a*sqrt(1-4*(n0/(n-1)).^2);b=besin(0,[beta a]);
if ty=='mkai',b=b-1;end,w=b(1:nn+1)/b(nn+2);w=real(w);
elseif ty=='dolp'|ty=='cheb'
if rem(n,2)==0,error('Dolph window requires odd N'),return,end
if nargin<3,a=60;end,nn=fix((n-1)/2);e=10^(-a/20);m=2*nn;b=cosh(acosh(1/e)/m);
f=b*cos(pi*(1:m)/m);d=zeros(1,m);i1=find(abs(f)<=1);i2=find(f>1);i3=find(f<-1);
if ~isempty(i1),d(i1)=cos(m*acos(f(i1)));end
if ~isempty(i2),d(i2)=cosh(m*acosh(f(i2)));end
if ~isempty(i3),d(i3)=((-1)^m)*cosh(m*acosh(-f(i3)));end
w=fft(d);w=abs(w(nn+1:-1:1)/w(1));w(1)=w(1)/2;
else
error('Unrecognized window type !'),return
end
w=[w w(nn+1-rem(n,2):-1:1)];