Mercurial > hg > ltpda
diff m-toolbox/test/test_ao_zDomainFit_1.m @ 0:f0afece42f48
Import.
| author | Daniele Nicolodi <nicolodi@science.unitn.it> |
|---|---|
| date | Wed, 23 Nov 2011 19:22:13 +0100 |
| parents | |
| children |
line wrap: on
line diff
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/m-toolbox/test/test_ao_zDomainFit_1.m Wed Nov 23 19:22:13 2011 +0100 @@ -0,0 +1,199 @@ +% Test script for zDomainFit +% +% L. Ferraioli 02-12-08 +% +% $Id: test_ao_zDomainFit_1.m,v 1.7 2009/12/02 16:50:49 luigi Exp $ +% +%% Building a frequency response + +% Create a frequency-series AO +% pl_data = plist('fsfcn', '0.01./(0.01+f)', 'f1', 1e-6, 'f2', 5, 'nf', 1000); +pl_data = plist('fsfcn', '(1e-3./(f).^2 + 1e3./(0.001+f) + 1e5.*f.^2).*1e-10', 'f1', 1e-6, 'f2', 5, 'nf', 100); +a = ao(pl_data); +a.setName; + +% iplot(a) + +%% Fitting 1 +% Check if Mean Square Error is lower than FITTOL and its relative +% variation is lower than MSEVARTOL + +% Fitting parameter list +pl_fit = plist('FS',10,... + 'AutoSearch','on',... + 'StartPoles',[],... + 'StartPolesOpt','clog',... + 'maxiter',40,... + 'minorder',10,... + 'maxorder',25,... + 'weights',[],... + 'weightparam','abs',... + 'CONDTYPE','MSE',... + 'FITTOL',1e-4,... % check if MSE is lower than 1e-4 + 'MSEVARTOL',1e-2,... + 'Plot','off',... + 'ForceStability','off',... + 'CheckProgress','off'); + +% Do fit +mod = zDomainFit(a, pl_fit); + +%% Fitting 2 +% Check if Rresiduals log difference is larger than FITTOL and MSE relative +% variation is lower than MSEVARTOL + +% Fitting parameter list +pl_fit = plist('FS',10,... + 'AutoSearch','on',... + 'StartPoles',[],... + 'StartPolesOpt','clog',... + 'maxiter',40,... + 'minorder',10,... + 'maxorder',25,... + 'weights',[],... + 'weightparam','abs',... + 'CONDTYPE','RLD',... + 'FITTOL',1,... % check if RLD is larger than 1 + 'MSEVARTOL',1e-2,... + 'Plot','off',... + 'ForceStability','off',... + 'CheckProgress','off'); + +% Do fit +mod = zDomainFit(a, pl_fit); + +%% Fitting 3 +% Check if Rresiduals spectral flatness is larger than FITTOL and MSE relative +% variation is lower than MSEVARTOL + +% Fitting parameter list +pl_fit = plist('FS',10,... + 'AutoSearch','on',... + 'StartPoles',[],... + 'StartPolesOpt','clog',... + 'maxiter',40,... + 'minorder',10,... + 'maxorder',25,... + 'weights',[],... + 'weightparam','abs',... + 'CONDTYPE','RSF',... + 'FITTOL',0.03,... % check if RSF is larger than 0.03 + 'MSEVARTOL',1e-2,... + 'Plot','off',... + 'ForceStability','off',... + 'CheckProgress','off'); + +% Do fit +mod = zDomainFit(a, pl_fit); + +%% Comparison + +resp = mod.procinfo.find('FIT_RESP'); +resids = mod.procinfo.find('FIT_RESIDUALS'); +mse = mod.procinfo.find('FIT_MSE'); +iplot(a,resp,abs(resids)) + +%% Testing the ability of reproducing models + +clear all +% Define a pratial fraction miir object +f = logspace(-6,log10(5),30); +fs = 10; +res = [0.7 0.2+0.01i 0.2-0.01i]; +poles = [0.5 0.1+0.07i 0.1-0.07i]; +tmod(3,1) = miir; + +for ii = 1:3 + tmod(ii,1) = miir(res(ii),[1 -poles(ii)],fs); +end + +btmod = filterbank(plist('filters',tmod,'type','parallel')); + +rtmod = resp(btmod.filters,plist('bank','parallel','f',f.')); + +% do the fit on the response +fs = 10; +plstd2 = plist('FS',fs,... + 'AutoSearch','off',... + 'StartPoles',[],... + 'StartPolesOpt','clog',... + 'maxiter',60,... + 'minorder',3,... + 'maxorder',3,... + 'weights',[],... + 'weightparam','abs',... + 'CONDTYPE','MSE',... + 'FITTOL',1e-3,... % check if MSE is lower than 1e-3 + 'MSEVARTOL',1e-2,... + 'Plot','off',... + 'ForceStability','off',... + 'CheckProgress','off'); +out = zDomainFit(rtmod,plstd2); + +% Check that coefficients are equal within a tolerance +for kk=1:numel(out.filters) + if any(abs(tmod(kk).a - out.filters(kk).a)>1e-6) + disp([num2str(kk) 'numerator coefficients are not consistent']); + else + disp([num2str(kk) 'numerator coefficients are consistent']) + end + if any(abs(tmod(kk).b - out.filters(kk).b)>1e-6) + disp([num2str(kk) 'denominator coefficients are not consistent']); + else + disp([num2str(kk) 'denominator coefficients are consistent']) + end + +end + + +%% help example + +pl = plist('fsfcn', '(1e-3./(2.*pi.*1i.*f).^2 + 1e3./(0.001+2.*pi.*1i.*f) + 1e5.*(2.*pi.*1i.*f).^2).*1e-10', 'f1', 1e-6, 'f2', 5, 'nf', 100); +a = ao(pl); +a.setName; + +% Fitting parameter list +pl_fit = plist('FS',10,... + 'AutoSearch','on',... + 'StartPoles',[],... + 'StartPolesOpt','clog',... + 'maxiter',50,... + 'minorder',15,... + 'maxorder',30,... + 'weights',[],... + 'weightparam','abs',... + 'CONDTYPE','MSE',... + 'FITTOL',1e-2,... + 'MSEVARTOL',1e-1,... + 'Plot','on',... + 'ForceStability','on',... + 'CheckProgress','off'); + +% Do fit +mod = zDomainFit(a, pl_fit); + +%% help test + +pl = plist('fsfcn', '(1e-3./(2.*pi.*1i.*f).^2 + 1e3./(0.001+2.*pi.*1i.*f) + 1e5.*(2.*pi.*1i.*f).^2).*1e-10', 'f1', 1e-6, 'f2', 5, 'nf', 100); +a = ao(pl); +a.setName; + +% Fitting parameter list +pl_fit = plist('FS',10,... + 'AutoSearch','on',... + 'StartPoles',[],... + 'StartPolesOpt','clog',... + 'maxiter',50,... + 'minorder',15,... + 'maxorder',30,... + 'weights',[],... + 'weightparam','abs',... + 'CONDTYPE','MSE',... + 'FITTOL',1e-2,... + 'MSEVARTOL',1e-1,... + 'Plot','on',... + 'ForceStability','on',... + 'CheckProgress','off'); + +% Do fit +mod = zDomainFit(a, pl_fit);