RESP Make a frequency response of the filter. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% DESCRIPTION: RESP Make a frequency response of the filter. The input filter should be an miir 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 'scale' - spacing of frequencies: 'lin' or 'log' or 'f' - a vector of frequency values OUTPUTS: a - an analysis object The response is calculated as: a0 + a1*exp(s/fs) + ... + an*exp((n-1)s/fs) H(s) = gain * --------------------------------------------- b0 + b1*exp(s/fs) + ... + bm*exp((m-1)s/fs) %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% The following call returns a cell array that contains the default parameter sets names: >> sets = resp(miir, 'Params') The following call returns a parameter list object that contains the default parameters for the jj-th set: >> pl = resp(miir, 'Params', sets{jj}) The following call returns a string that contains the routine CVS version: >> version = resp(miir, 'Version') The following call returns a string that contains the routine category: >> category = resp(miir, 'Category') VERSION: $Id: resp.m,v 1.10 2008/02/24 21:43:59 mauro Exp $ NOTE: Some of the code below is taken from Mathworks's treeplot.m 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 miir object. 0008 % 0009 % The response is returned as a frequency-series in an 0010 % analysis object. 0011 % 0012 % If no outputs are specified, the frequency-series is plotted as 0013 % mag/deg. 0014 % 0015 % CALL: a = resp(filt, pl) 0016 % resp(filt, pl); 0017 % 0018 % PARAMETERS: 'f1' - the start frequency 0019 % 'f2' - the stop frequency 0020 % 'nf' - number of evaluation points 0021 % 'scale' - spacing of frequencies: 'lin' or 'log' 0022 % or 0023 % 'f' - a vector of frequency values 0024 % 0025 % OUTPUTS: 0026 % a - an analysis object 0027 % 0028 % The response is calculated as: 0029 % 0030 % a0 + a1*exp(s/fs) + ... + an*exp((n-1)s/fs) 0031 % H(s) = gain * --------------------------------------------- 0032 % b0 + b1*exp(s/fs) + ... + bm*exp((m-1)s/fs) 0033 % 0034 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0035 % 0036 % The following call returns a cell array that contains the default parameter 0037 % sets names: 0038 % 0039 % >> sets = resp(miir, 'Params') 0040 % 0041 % The following call returns a parameter list object that contains the 0042 % default parameters for the jj-th set: 0043 % 0044 % >> pl = resp(miir, 'Params', sets{jj}) 0045 % 0046 % The following call returns a string that contains the routine CVS version: 0047 % 0048 % >> version = resp(miir, 'Version') 0049 % 0050 % The following call returns a string that contains the routine category: 0051 % 0052 % >> category = resp(miir, 'Category') 0053 % 0054 % 0055 % VERSION: $Id: resp.m,v 1.10 2008/02/24 21:43:59 mauro Exp $ 0056 % 0057 % NOTE: Some of the code below is taken from Mathworks's treeplot.m 0058 % 0059 % HISTORY: 27-08-02 M Hewitson 0060 % Creation 0061 % 0062 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0063 0064 ALGONAME = mfilename; 0065 VERSION = '$Id: resp.m,v 1.10 2008/02/24 21:43:59 mauro Exp $'; 0066 CATEGORY = 'Signal Processing'; 0067 0068 %% 'Params', 'Version' and 'Category' Call 0069 if (nargin == 2 || nargin == 3) && ... 0070 isa(varargin{1}, 'miir') && ... 0071 ischar(varargin{2}) 0072 in = char(varargin{2}); 0073 if strcmp(in, 'Params') 0074 if nargin == 2 0075 varargout{1} = getDefaultPL(); 0076 else 0077 varargout{1} = getDefaultPL(varargin{3}); 0078 end 0079 return 0080 elseif strcmp(in, 'Version') 0081 varargout{1} = VERSION; 0082 return 0083 elseif strcmp(in, 'Category') 0084 varargout{1} = CATEGORY; 0085 return 0086 end 0087 end 0088 0089 0090 %% get input filter 0091 filt = varargin{1}; 0092 if ~isa(filt, 'miir') 0093 error('### first input must be an miir filter object.'); 0094 end 0095 0096 0097 %% Collects parameters 0098 % create empty input if it is not specified 0099 if nargin == 1 0100 pl = plist(); 0101 else 0102 pl = varargin{2}; 0103 end 0104 0105 % parse input parameter list 0106 if ~isempty(find(pl, 'f')) 0107 pl = combine(pl, getDefaultPL('List')); 0108 else 0109 pl = combine(pl, getDefaultPL('Range')); 0110 end 0111 0112 % fill parameter list 0113 if isempty(find(pl, 'f')) 0114 % Compute from frequency range 0115 f1 = find(pl, 'f1'); 0116 f2 = find(pl, 'f2'); 0117 ndata = find(pl, 'nf'); 0118 scale = find(pl, 'scale'); 0119 switch scale 0120 case 'lin' 0121 f = linspace(f1, f2, ndata); 0122 case 'log' 0123 f = logspace(log10(f1), log10(f2), ndata); 0124 end 0125 else 0126 % Compute from frequency list 0127 f = find(pl, 'f'); 0128 % Want to deal with rows 0129 if size(f,1) > 1 0130 f = f.'; 0131 end 0132 ndata = length(f); 0133 end 0134 0135 % create output parameter list, with the frequency array values and the 0136 % input filter 0137 plo = plist('f', f, ... 0138 'filter', filt); 0139 0140 0141 %% compute Laplace vector 0142 s = -1i.*2*pi.*f; 0143 0144 0145 %% Compute filter response 0146 num = zeros(1, ndata); 0147 for n=1:length(filt.a) 0148 num = num + filt.a(n).*exp(s.*(n-1)/filt.fs); 0149 end 0150 denom = zeros(1, ndata); 0151 for n=1:length(filt.b) 0152 denom = denom + filt.b(n).*exp(s.*(n-1)/filt.fs); 0153 end 0154 dresp = num ./ denom; 0155 dresp = dresp .* filt.gain; % apply the gain 0156 0157 % mag = 20*log10(abs(dresp)); 0158 % phase = angle(dresp)*180/pi; 0159 0160 0161 %% Create an analysis object 0162 0163 % create new output fsdata 0164 fs = fsdata(f, dresp, filt.fs); 0165 fs = set(fs, 'name', sprintf('resp(%s)', filt.name)); 0166 0167 % create new output history 0168 h = history(ALGONAME, VERSION, plo); 0169 0170 % make output analysis object 0171 b = ao(fs, h); 0172 0173 % set name 0174 b = setnh(b, 'name', sprintf('resp(%s)', filt.name)); 0175 0176 % Outputs 0177 if nargout == 0 0178 iplot(b) 0179 end 0180 0181 if nargout == 1 0182 varargout{1} = b; 0183 end 0184 if nargout > 1 0185 error('incorrect output arguments'); 0186 end 0187 0188 0189 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0190 function out = getDefaultPL(varargin) 0191 0192 % List of available parameter sets 0193 sets = {'List', 'Range'}; 0194 0195 if nargin == 0 0196 out = sets; 0197 return 0198 end 0199 0200 set = varargin{1}; 0201 0202 switch set 0203 case 'List' 0204 out = plist('f', []); 0205 case 'Range' 0206 out = plist('f1', 0.001,... 0207 'f2', 1,... 0208 'nf', 1000,... 0209 'scale', 'log'); 0210 otherwise 0211 out = plist(); 0212 end 0213 0214 % END 0215