line source
+ − % STARTPOLES defines starting poles for fitting procedures ctfit, dtfit.
+ − %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ − % DESCRIPTION
+ − %
+ − % Defines the starting poles for the fitting procedure with ctfit or
+ − % dtfit. Starting poles definition process is different for s-domain
+ − % and z-domain.
+ − % s-domain identification:
+ − % Starting poles can be real or complex in conjugate couples. Real
+ − % poles are chosen on the [-2*pi*f(1),-2*pi*f(end)] intervall. Complex poles can
+ − % be defined with the real and imaginary parts logspaced or linespaced
+ − % on the intervall [2\pi f(1),2\pi f(end)]. Ratio between real and
+ − % imaginary part can be setted by the user.
+ − % z-domain identification:
+ − % Starting poles can be real or come in complex conjugate couples. Real
+ − % poles are chosen on the [-1,1] intervall. Complex poles are
+ − % chosen inside the unit circle as:
+ − % \alfa e^{j\pi\theta}
+ − % where \theta is linespaced inside the intervall [0,2\pi]. In this
+ − % case two different methods are used: the first method define angles
+ − % as \theta = linspace(0,pi,N/2+1), then take out the first element and
+ − % construct the complex conjugate couples. If N is odd the first
+ − % element is added as real pole. With this method the last two
+ − % conjugates poles have a real part very similar to that of the real
+ − % pole. This may generate problems so a second method is implemented in
+ − % which the angle are generated as \theta = linspace(0,pi,N/2+2). Then
+ − % the first and the last elements of the set are taken out and the
+ − % first element is used only for N odd. The last element instead is
+ − % never used. This allow to have well distributed poles on the unit
+ − % circle. The amplitude parameter \alfa can be set by the user.
+ − %
+ − % CALL:
+ − %
+ − % spoles = startpoles(order,f,params)
+ − %
+ − % INPUT:
+ − %
+ − % order: is the function order, ie. the number of desired poles
+ − % f: is the frequency vector in Hz
+ − % params: is a struct with the setting parameters
+ − %
+ − % params.type = 'CONT' --> Output a set of poles for s-domain
+ − % identification
+ − % params.type = 'DISC' --> Output a set of poles for z-domain
+ − % identification
+ − %
+ − % params.spolesopt = 1 --> generate linear spaced real poles on
+ − % the intervall [-2*pi*f(1),-2*pi*f(end)] for s-domain and [-1,1] for z-domain.
+ − % params.spolesopt = 2 --> in case of s-domain generates logspaced
+ − % complex conjugates poles. In case of z-domain generates complex
+ − % conjugates poles with \theta = linspace(0,pi,N/2+1).
+ − % params.spolesopt = 3 --> in case of s-domain generates linespaced
+ − % complex conjugates poles. In case of z-domain generates complex
+ − % conjugates poles with \theta = linspace(0,pi,N/2+2). We advice to
+ − % make use of this option for z-domain identification.
+ − %
+ − % params.pamp = # --> s-domain: set the ratio \alfa/\beta between
+ − % poles real and imaginary parts. Adviced value is 0.01.
+ − % params.pamp = # --> z-domain: set the amplitude of the poles.
+ − % Adviced value is 0.98.
+ − %
+ − % OUTPUT:
+ − %
+ − % spoles: is the set of starting poles
+ − %
+ − %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ − % VERSION: $Id: startpoles.m,v 1.6 2010/04/29 09:00:00 luigi Exp $
+ − % HISTORY: 08-10-2008 L Ferraioli
+ − % Creation
+ − %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ −
+ − function spoles = startpoles(order,f,params)
+ −
+ − % Default input struct
+ − defaultparams = struct('spolesopt',1, 'type','CONT', 'pamp', 0.98);
+ −
+ − names = {'spolesopt', 'type', 'pamp'};
+ −
+ − % collecting input and default params
+ − if ~isempty(params)
+ − for jj=1:length(names)
+ − if isfield(params, names(jj))
+ − defaultparams.(names{1,jj}) = params.(names{1,jj});
+ − end
+ − end
+ − end
+ −
+ − type = defaultparams.type;
+ − spolesopt = defaultparams.spolesopt;
+ − pamp = defaultparams.pamp;
+ −
+ − N = order;
+ −
+ − % switching between continuous and discrete
+ − switch type
+ − case 'CONT'
+ − switch spolesopt
+ −
+ − case 0
+ − disp(' Using external starting poles')
+ −
+ − case 1 % real starting poles
+ − spoles = -1.*2.*pi.*linspace(f(1),f(end),N).';
+ −
+ − case 2 % complex logspaced starting poles
+ − if f(1)==0
+ − bet=2.*pi.*logspace(log10(f(2)),log10(f(end)),N/2);
+ − else
+ − bet=2.*pi.*logspace(log10(f(1)),log10(f(end)),N/2);
+ − end
+ − spoles=[];
+ − for n=1:length(bet)
+ − alf=-bet(n)*pamp;
+ − spoles=[spoles;(alf-1i*bet(n));(alf+1i*bet(n))];
+ − end
+ − if (N-2*floor(N/2)) ~= 0
+ − % rpl = rand(1,1);
+ − % if rpl > 0
+ − % rpl = -1*rpl;
+ − % end
+ − rpl = -1;
+ − spoles = [rpl; spoles];
+ − end
+ −
+ − case 3 % complex linspaced starting poles
+ − bet=linspace(2*pi*f(1),2*pi*f(end),N/2);
+ − spoles=[];
+ − for n=1:length(bet)
+ − alf=-bet(n)*pamp;
+ − spoles=[spoles;(alf-1i*bet(n));(alf+1i*bet(n))];
+ − end
+ − if (N-2*floor(N/2)) ~= 0
+ − % rpl = rand(1,1);
+ − % if rpl > 0
+ − % rpl = -1*rpl;
+ − % end
+ − rpl = -1;
+ − spoles = [rpl; spoles];
+ − end
+ −
+ − end
+ − case 'DISC'
+ − switch spolesopt
+ −
+ − case 0
+ − disp(' Using external starting poles')
+ −
+ − case 1 % real starting poles
+ − spoles = linspace(-0.99,0.99,N).';
+ −
+ − case 2 % complex starting poles
+ − ang = linspace(0,pi,N/2+1);
+ − spoles=[];
+ − for nn=2:length(ang)
+ − spoles=[spoles; exp(1i*ang(nn)); exp(-1i*ang(nn))]; % Taking complex conjugate pairs on the unit circle
+ − end
+ − if (N-2*floor(N/2)) ~= 0
+ − rpl = exp(1i*ang(1));
+ − spoles = [rpl; spoles];
+ − end
+ − spoles = spoles.*pamp; % shifting starting ?poles a little inside the unit circle
+ −
+ − case 3 % complex starting poles
+ − ang = linspace(0,pi,N/2+2);
+ − spoles=[];
+ − for nn=2:length(ang)-1
+ − spoles=[spoles; exp(1i*ang(nn)); exp(-1i*ang(nn))]; % Taking complex conjugate pairs on the unit circle
+ − end
+ − if (N-2*floor(N/2)) ~= 0
+ − rpl = exp(1i*ang(1));
+ − spoles = [rpl; spoles];
+ − end
+ − spoles = spoles.*pamp; % shifting starting ?poles a little inside the unit circle
+ −
+ − end
+ −
+ − end
+ −
+ − end