X-Git-Url: https://git.creatis.insa-lyon.fr/pubgit/?p=CreaPhase.git;a=blobdiff_plain;f=octave_packages%2Fsignal-1.1.3%2F__power.m;fp=octave_packages%2Fsignal-1.1.3%2F__power.m;h=feb5d5d76cfc88811aa6a0a2a1864e3ab27c7b6b;hp=0000000000000000000000000000000000000000;hb=f5f7a74bd8a4900f0b797da6783be80e11a68d86;hpb=1705066eceaaea976f010f669ce8e972f3734b05 diff --git a/octave_packages/signal-1.1.3/__power.m b/octave_packages/signal-1.1.3/__power.m new file mode 100644 index 0000000..feb5d5d --- /dev/null +++ b/octave_packages/signal-1.1.3/__power.m @@ -0,0 +1,89 @@ +## Copyright (C) 1999 Paul Kienzle +## +## This program is free software; you can redistribute it and/or modify it under +## the terms of the GNU General Public License as published by the Free Software +## Foundation; either version 3 of the License, or (at your option) any later +## version. +## +## This program is distributed in the hope that it will be useful, but WITHOUT +## ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or +## FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more +## details. +## +## You should have received a copy of the GNU General Public License along with +## this program; if not, see . + +## usage: [P, w] = __power (b, a, [, nfft [, Fs]] [, range] [, units]) +## +## Plot the power spectrum of the given ARMA model. +## +## b, a: filter coefficients (b=numerator, a=denominator) +## nfft is number of points at which to sample the power spectrum +## Fs is the sampling frequency of x +## range is 'half' (default) or 'whole' +## units is 'squared' or 'db' (default) +## range and units may be specified any time after the filter, in either +## order +## +## Returns P, the magnitude vector, and w, the frequencies at which it +## is sampled. If there are no return values requested, then plot the power +## spectrum and don't return anything. + +## TODO: consider folding this into freqz --- just one more parameter to +## TODO: distinguish between 'linear', 'log', 'logsquared' and 'squared' + +function [varargout] = __power (b, a, varargin) + if (nargin < 2 || nargin > 6) print_usage; endif + + nfft = []; + Fs = []; + range = []; + range_fact = 1.0; + units = []; + + pos = 0; + for i=1:length(varargin) + arg = varargin{i}; + if strcmp(arg, 'squared') || strcmp(arg, 'db') + units = arg; + elseif strcmp(arg, 'whole') + range = arg; + range_fact = 1.0; + elseif strcmp(arg, 'half') + range = arg; + range_fact = 2.0; + elseif ischar(arg) + usage(usagestr); + elseif pos == 0 + nfft = arg; + pos++; + elseif pos == 1 + Fs = arg; + pos++; + else + usage(usagestr); + endif + endfor + + if isempty(nfft); nfft = 256; endif + if isempty(Fs); Fs = 2; endif + if isempty(range) + range = 'half'; + range_fact = 2.0; + endif + + [P, w] = freqz(b, a, nfft, range, Fs); + + P = (range_fact/Fs)*(P.*conj(P)); + if nargout == 0, + if strcmp(units, 'squared') + plot(w, P, ";;"); + else + plot(w, 10.0*log10(abs(P)), ";;"); + endif + endif + if nargout >= 1, varargout{1} = P; endif + if nargout >= 2, varargout{2} = w; endif + +endfunction +