RESP Make a frequency response of the filter. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% DESCRIPTION: RESP Make a frequency response of the filter. The input filter should be an mfir object. The response is returned as a frequency-series in an analysis object. If no outputs are specified, the frequency-series is plotted as mag/deg. CALL: a = resp(filt, pl); resp(filt, pl); PARAMETERS: f1 - the start frequency f2 - the stop frequency nf - number of evaluation points or f - a vector of frequency values OUTPUTS: a - an analysis object The response is calculated as: H(s) = gain * (a0 + a1*exp(s/fs) + ... + an*exp((n-1)s/fs)) VERSION: $Id: resp.m,v 1.4 2007/06/07 08:51:26 ingo Exp $ HISTORY: 27-08-02 M Hewitson Creation %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0001 function [varargout] = resp(varargin) 0002 % RESP Make a frequency response of the filter. 0003 % 0004 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0005 % 0006 % DESCRIPTION: RESP Make a frequency response of the filter. 0007 % The input filter should be an mfir object. 0008 % 0009 % The response is returned as a frequency-series in an analysis object. 0010 % 0011 % If no outputs are specified, the frequency-series is plotted as 0012 % mag/deg. 0013 % 0014 % CALL: a = resp(filt, pl); 0015 % resp(filt, pl); 0016 % 0017 % PARAMETERS: f1 - the start frequency 0018 % f2 - the stop frequency 0019 % nf - number of evaluation points 0020 % or 0021 % f - a vector of frequency values 0022 % 0023 % OUTPUTS: a - an analysis object 0024 % 0025 % The response is calculated as: 0026 % 0027 % H(s) = gain * (a0 + a1*exp(s/fs) + ... + an*exp((n-1)s/fs)) 0028 % 0029 % VERSION: $Id: resp.m,v 1.4 2007/06/07 08:51:26 ingo Exp $ 0030 % 0031 % 0032 % HISTORY: 27-08-02 M Hewitson 0033 % Creation 0034 % 0035 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0036 0037 0038 ALGONAME = mfilename; 0039 VERSION = '$Id: resp.m,v 1.4 2007/06/07 08:51:26 ingo Exp $'; 0040 0041 % get input filter 0042 filt = varargin{1}; 0043 if ~isa(filt, 'mfir') 0044 error('### first input must be an mfir filter object.'); 0045 end 0046 0047 0048 % create empty input if it is not specified 0049 if nargin == 1 0050 pl = plist(); 0051 else 0052 pl = varargin{2}; 0053 end 0054 0055 % create empty output parameter list 0056 plo = plist(); 0057 0058 % parse input parameter list 0059 f1 = find(pl, 'f1'); 0060 f2 = find(pl, 'f2'); 0061 ndata = find(pl, 'nf'); 0062 f = find(pl, 'f'); 0063 0064 % fill parameter list 0065 if isempty(f) 0066 0067 % we look for the other parameter set 0068 if isempty(f1) 0069 f1 = 0; 0070 end 0071 plo = append(plo, param('f1', f1)); 0072 % we look for the other parameter set 0073 if isempty(f2) 0074 f2 = filt.fs/2; 0075 end 0076 plo = append(plo, param('f2', f2)); 0077 % we look for the other parameter set 0078 if isempty(ndata) 0079 ndata = 1000; 0080 end 0081 plo = append(plo, param('nf', ndata)); 0082 0083 % compute f from this 0084 f = linspace(f1, f2, ndata); 0085 else 0086 plo = append(plo, param('f', f)); 0087 0088 ndata = length(f); 0089 0090 end 0091 0092 % append input filter to output parameter list 0093 plo = append(plo, param('filter', filt)); 0094 0095 % compute Laplace vector 0096 s = -1i.*2*pi.*f; 0097 0098 % Compute filter response 0099 num = zeros(1, ndata); 0100 for n=1:length(filt.a) 0101 num = num + filt.a(n).*exp(s.*(n-1)/filt.fs); 0102 end 0103 dresp = num .* filt.g; % apply the gain 0104 0105 % num = zeros(1, ndata); 0106 % for n=1:length(filt.a) 0107 % num = num + filt.a(n).*exp(s.*(n-1)/filt.fs); 0108 % end 0109 % denom = zeros(1, ndata); 0110 % for n=1:length(filt.b) 0111 % denom = denom + filt.b(n).*exp(s.*(n-1)/filt.fs); 0112 % end 0113 % dresp = num ./ denom; 0114 % dresp = dresp .* filt.g; % apply the gain 0115 0116 % mag = 20*log10(abs(dresp)); 0117 % phase = angle(dresp)*180/pi; 0118 0119 %----------------------------------------------------------- 0120 % Create an analysis object 0121 0122 % create new output fsdata 0123 fs = fsdata(f, dresp, filt.fs); 0124 fs = set(fs, 'name', sprintf('resp(%s)', filt.name)); 0125 0126 % create new output history 0127 h = history(ALGONAME, VERSION, plo); 0128 0129 % make output analysis object 0130 b = ao(fs, h); 0131 0132 % set name 0133 b = set(b, 'name', sprintf('resp(%s)', filt.name)); 0134 0135 % Outputs 0136 if nargout == 0 0137 plot(b) 0138 end 0139 0140 if nargout == 1 0141 varargout{1} = b; 0142 end 0143 if nargout > 1 0144 error('incorrect output arguments'); 0145 end 0146 0147