% 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 endend