lpsig.m

很多matlab的源代码

function [H,P,F,ty] = lpsig(h,N)
%LPSIG	  Amplitude and phase response of Linear Phase Sequences.
%
%         [A,P,F,ty] = LPSIG(h,N) returns type of linear phase sequence h 
%	  A= AMPLITUDE SPECTRUM  vs F (over 0<=F<=1)
%         N=number of frequency intervals in F. Default: N=500 
%         The sequence type is returned in type ty
%	  P=phase, asuming that h is a CAUSAL sequence
%         ty=1 for even sym,odd length   ty=2 for even sym, even length
%         ty=3 for odd sym,odd length    ty=4 for odd sym, even length
% 
%         LPSIG(.) plots results if no output arguments are specified.
%
%         LPSIG (with no input arguments) invokes the following example:
%
%         % Plot the amplitude response of h=[1 2 3 4 4 3 2 1], 
%         >>lpsig([1 2 3 4 4 3 2 1]);


% 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 lpsig,disp('Strike a key to see results of the example')
pause,lpsig([1 2 3 4 4 3 2 1]);return,end

if nargin<2,N=500;end
h=h(:);hf=flipud(h);
if ~any(abs(hf-h)>10*eps), sy='e';      % Even Symmetry
elseif ~any(abs(hf+h)>10*eps), sy='o';  % Odd  Symmetry
else
error('Sequence is not linear phase'),return
end

n=length(h);M=0.5*(n-1);L=fix(0.5*n-1);H=0;

if rem(n,2)==1,  %odd length
if sy=='e', ty=1;H=h(M+1);else, ty=3;end
else
if sy=='e', ty=2;else, ty=4;end
end

%%%NEW
F=(-2*N:2*N)/N;     % Digital Frequency
P=-2*pi*M*F;   % Phase


if ty<3
for k=0:L,H=H+2*h(k+1)*cos(2*(M-k)*pi*F);end
else
for k=0:L,H=H+2*h(k+1)*sin(2*(M-k)*pi*F);end
P=P+0.5*pi;
end


if nargout==0,
subplot(211),plot(F,H),grid,
subplot(212),plot(F,P),grid
end