Mercurial > hg > ltpda
view m-toolbox/test/utils/test_ndeigcsd.m @ 7:1e91f84a4be8 database-connection-manager
Make ltpda_up.retrieve work with java.sql.Connection objects
| author | Daniele Nicolodi <nicolodi@science.unitn.it> |
|---|---|
| date | Mon, 05 Dec 2011 16:20:06 +0100 |
| parents | f0afece42f48 |
| children |
line wrap: on
line source
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % TEST ndeigcsd % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % HISTORY: 23-04-2009 L Ferraioli % Creation %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% 2 Dim Test %%% Define startinf TFs fs = 10; f = [logspace(-6,log10(fs/2),2000)]'; % Model Stefano TF11 dRes11 = [2.44554138162509e-011 - 1.79482547894083e-011i; 2.44554138162509e-011 + 1.79482547894083e-011i; 2.66402334803101e-009 + 1.1025122049153e-009i; 2.66402334803101e-009 - 1.1025122049153e-009i; -7.3560293387644e-009; -1.82811618589835e-009 - 1.21803627800855e-009i; -1.82811618589835e-009 + 1.21803627800855e-009i; 1.16258677367555e-009; 1.65216557639319e-016; -1.78092396888606e-016; -2.80420398962379e-017; 9.21305973049041e-013 - 8.24686706827269e-014i; 9.21305973049041e-013 + 8.24686706827269e-014i; 5.10730060739905e-010 - 3.76571756625722e-011i; 5.10730060739905e-010 + 3.76571756625722e-011i; 3.45893698149735e-009; 3.98139182134446e-014 - 8.25503935419059e-014i; 3.98139182134446e-014 + 8.25503935419059e-014i; -1.40595719147164e-011]; dPoles11 = [0.843464045655194 - 0.0959986292915475i; 0.843464045655194 + 0.0959986292915475i; 0.953187595424927 - 0.0190043625473383i; 0.953187595424927 + 0.0190043625473383i; 0.967176277937188; 0.995012027005247 - 0.00268322602801729i; 0.995012027005247 + 0.00268322602801729i; 0.996564761885673; 0.999999366165445; 0.999981722418555; 0.999921882627659; 0.999624431675213 - 0.000813407848742761i; 0.999624431675213 + 0.000813407848742761i; 0.997312006278751 - 0.00265611346834941i; 0.997312006278751 + 0.00265611346834941i; 0.990516544257531; 0.477796923118318 - 0.311064085401834i; 0.477796923118318 + 0.311064085401834i; 0]; dDTerms11 = 0; % Model Stefano TF12 dRes12 = [1.44258422208796e-017 + 7.07359428613009e-019i; 1.44258422208796e-017 - 7.07359428613009e-019i; -3.4918408053655e-021 - 1.05662874569329e-021i; -3.4918408053655e-021 + 1.05662874569329e-021i; -7.61773292876976e-021; 4.84357724603939e-020 + 2.38824204294595e-019i; 4.84357724603939e-020 - 2.38824204294595e-019i; -4.07088520945753e-020 - 2.31474543846105e-019i; -4.07088520945753e-020 + 2.31474543846105e-019i; 8.73316588658882e-023; -5.21840635377469e-020; 1.8461911504859e-023; 5.20105247464461e-020; -4.68960092394415e-022; -1.44261407664171e-017 + 6.8922564526833e-019i; -1.44261407664171e-017 - 6.8922564526833e-019i; 3.13688133935426e-022]; dPoles12 = [0.477546340377332 - 0.310830571032376i; 0.477546340377332 + 0.310830571032376i; 0.99790715414307 - 0.0028490561287024i; 0.99790715414307 + 0.0028490561287024i; 0.998014205354671 ; 0.999585354543332 - 0.000780408757425194i; 0.999585354543332 + 0.000780408757425194i; 0.99966003029931 - 0.000830944038363768i; 0.99966003029931 + 0.000830944038363768i; 0.999962770401331 ; 0.999981881865521 ; 0.999999365763457 ; 0.999981706320212 ; 0.99992421574188 ; 0.477898460791003 + 0.311001926610074i; 0.477898460791003 - 0.311001926610074i; 0]; dDTerms12 = 0; % Model Stefano Tf21 dRes21 = [-1.80035241582968e-016 + 1.99543917791863e-015i; -1.80035241582968e-016 - 1.99543917791863e-015i; -1.85590889333759e-013 - 1.23844418827409e-014i; -1.85590889333759e-013 + 1.23844418827409e-014i; 5.03656596876842e-013 ; -2.62470963499904e-013 + 2.30024232938878e-012i; -2.62470963499904e-013 - 2.30024232938878e-012i; -9.83780507870955e-018 ; 3.40426735130194e-021 ; 9.78322351492755e-018 ; -1.65010934542937e-020 ; 2.60918565203438e-015 + 1.0546609464659e-015i; 2.60918565203438e-015 - 1.0546609464659e-015i; 3.18105585405455e-014 + 2.48839990780042e-013i; 3.18105585405455e-014 - 2.48839990780042e-013i; 3.23021641947666e-013 ; 4.81265000078114e-016 - 3.18269170053848e-017i; 4.81265000078114e-016 + 3.18269170053848e-017i; 5.16260024128201e-018]; dPoles21 = [0.872004077421604 - 0.110344282822693i; 0.872004077421604 + 0.110344282822693i; 0.956884129232757 - 0.0225532091775074i; 0.956884129232757 + 0.0225532091775074i; 0.966514825697177 ; 0.995140550419744 - 0.000755127639524413i; 0.995140550419744 + 0.000755127639524413i; 0.999981802194393 ; 0.99999936576546 ; 0.999981722418555 ; 0.999921882627659 ; 0.999624431675213 + 0.000813407848742761i; 0.999624431675213 - 0.000813407848742761i; 0.997312006278751 + 0.00265611346834941i; 0.997312006278751 - 0.00265611346834941i; 0.990516544257531 ; 0.477796923118318 + 0.311064085401834i; 0.477796923118318 - 0.311064085401834i; 0]; dDTerms21 = 0; % Model Stefano Tf22 dRes22 = [1.1284521501259e-014; -3.72133611555879e-014 - 2.08232683444075e-014i; -3.72133611555879e-014 + 2.08232683444075e-014i; 9.84930639106637e-014 - 1.46640810672565e-013i; 9.84930639106637e-014 + 1.46640810672565e-013i; 2.51323684013671e-014 ; -5.64078525288305e-014 ; -1.62476406586366e-014 ; -1.08424815979566e-011 + 8.32328079357669e-012i; -1.08424815979566e-011 - 8.32328079357669e-012i; 2.41831559776112e-011]; dPoles22 = [0.988511243978897; 0.997305870640646 + 0.00211760900132725i; 0.997305870640646 - 0.00211760900132725i; 0.999626453270255 + 0.0008125673525946i; 0.999626453270255 - 0.0008125673525946i; 0.999999366366222 ; 0.999981706320212 ; 0.99992421574188 ; 0.477898460791003 - 0.311001926610074i; 0.477898460791003 + 0.311001926610074i; 0]; dDTerms22 = 0; % response calculation pfparams.type = 'disc'; pfparams.freq = f; pfparams.fs = 10; % TF11 pfparams.res = dRes11; pfparams.pol = dPoles11; pfparams.dterm = dDTerms11; pfr = utils.math.pfresp(pfparams); mtf11 = pfr.resp; % TF12 pfparams.res = dRes12; pfparams.pol = dPoles12; pfparams.dterm = dDTerms12; pfr = utils.math.pfresp(pfparams); mtf12 = pfr.resp; % TF21 pfparams.res = dRes21; pfparams.pol = dPoles21; pfparams.dterm = dDTerms21; pfr = utils.math.pfresp(pfparams); mtf21 = pfr.resp; % TF22 pfparams.res = dRes22; pfparams.pol = dPoles22; pfparams.dterm = dDTerms22; pfr = utils.math.pfresp(pfparams); mtf22 = pfr.resp; % CSD calculation with Papoulis Method Style % csd11 = mtf11.*conj(mtf11)+mtf12.*conj(mtf12); % csd12 = mtf11.*conj(mtf21)+mtf12.*conj(mtf22); % csd22 = mtf22.*conj(mtf22)+mtf21.*conj(mtf21); % csd21 = conj(csd12); % CSD = zeros(2,2,length(f)); % for hh = 1:length(f) % CSD(:,:,hh) = [csd11(hh) csd12(hh);csd21(hh) csd22(hh)]; % end % CSD = zeros(2,2,length(f)); % for hh = 1:length(f) % CSD(:,:,hh) = [csd11(hh) 0;0 csd22(hh)]; % end %% Get CSD Model fs = 10; f = [logspace(-6,log10(fs/2),2000)]'; % currPath = cd; % cd MDC1 % assuming to start from ..\ltpda\software\m-toolbox\test % [TF,CSD] = mdc1_tf_models(plist('f',f,'fs',fs)); % cd(currPath) % get noise psd shapes % parameters names parnames = {'DH','S11','S1D','SD1','SDD','Dh1','Dh2','w1','w2',... 'TH','Th1','Th2'}; % nominal params values parvalues = {1,1,0,-1e-4,1,1,1,-1.3e-6,-2.0e-6,0.35,0.25,0.28}; % Noise model N = matrix(plist('built-in','mdc3_ifo2ifo_noise')); % evalueat model [rw,cl]=size(N.objs); for ii=1:rw for jj=1:cl N.objs(ii,jj).setParams(parnames,parvalues); N.objs(ii,jj).setXvals(f); sp(ii,jj)=eval(N.objs(ii,jj)); if ii==jj sp(ii,jj)= abs(sp(ii,jj)); end end end csd11 = sp(1,1).y; csd12 = sp(1,2).y; csd21 = sp(2,1).y; csd22 = sp(2,2).y; CSD = zeros(2,2,length(f)); for hh = 1:length(f) CSD(:,:,hh) = [csd11(hh) csd12(hh);csd21(hh) csd22(hh)]; end %% Eigcsd 2-dim h = utils.math.ndeigcsd(CSD,'MTD','PAP'); %% Test eigenshuffle [l,m,npts] = size(CSD); % Finding suppression k = min(sqrt(csd11./csd22)); if k>=1 suppr = floor(k); else n=0; while k<1 k=k*10; n=n+1; end k = floor(k); suppr = k*10^(-n); end supmat = [1 0;0 suppr]; % isupmat = [1 0;0 1/suppr]; isupmat = inv(supmat); PP = CSD; for ii=1:npts PP(:,:,ii)=supmat*CSD(:,:,npts)*supmat; end [V,D] = eigenshuffle(CSD); h = ones(l,m,npts); for ii=1:npts HH = isupmat*V(:,:,ii)*sqrt(diag(D(:,ii))); h(:,:,ii) = HH; end %% Plot CSD % CSD calculation ncsd11 = squeeze(h(1,1,:).*conj(h(1,1,:))+h(1,2,:).*conj(h(1,2,:))); ncsd12 = squeeze(h(1,1,:).*conj(h(2,1,:))+h(1,2,:).*conj(h(2,2,:))); ncsd22 = squeeze(h(2,2,:).*conj(h(2,2,:))+h(2,1,:).*conj(h(2,1,:))); ncsd21 = conj(csd12); figure() loglog(f,abs(csd11),'k') hold on grid on loglog(f,abs(ncsd11),'r') figure() loglog(f,abs(csd12),'k') hold on grid on loglog(f,abs(ncsd12),'r') figure() semilogx(f,angle(csd12),'k') hold on grid on semilogx(f,angle(ncsd12),'r') figure() loglog(f,abs(csd22),'k') hold on grid on loglog(f,abs(ncsd22),'r') %% plot TFs figure() loglog(f,abs(mtf11),'k') hold on grid on loglog(f,abs(squeeze(h(1,1,:))),'r') figure() loglog(f,abs(mtf12),'k') hold on grid on loglog(f,abs(squeeze(h(1,2,:))),'r') figure() loglog(f,abs(mtf21),'k') hold on grid on loglog(f,abs(squeeze(h(2,1,:))),'r') figure() loglog(f,abs(mtf22),'k') hold on grid on loglog(f,abs(squeeze(h(2,2,:))),'r') %% plot TFs phase figure() semilogx(f,angle(mtf11),'k') hold on grid on semilogx(f,angle(squeeze(h(1,1,:))),'r') figure() semilogx(f,angle(mtf12),'k') hold on grid on semilogx(f,angle(squeeze(h(1,2,:))),'r') figure() semilogx(f,angle(mtf21),'k') hold on grid on semilogx(f,angle(squeeze(h(2,1,:))),'r') figure() semilogx(f,angle(mtf22),'k') hold on grid on semilogx(f,angle(squeeze(h(2,2,:))),'r') %% 3 Dim Test % sampling frequency fs = 1; % Hz % frequency vector f = logspace(-6,log10(0.5),300); f = f.'; %%% Set the 3 channel model % d = [1 -1.1 0.5 0.02]; % % h11 = miir(1e-1.*[1 -0.5],d,fs); % h12 = miir(1e-4.*[0 -0.5],d,fs); % h13 = miir(1e-2.*[1 0.2],d,fs); % % h21 = miir(1e-3.*[0 0.4],d,fs); % h22 = miir(1e-5.*[1 -0.6],d,fs); % h23 = miir(1e-6.*[1 0.3],d,fs); % % h31 = miir(1e-4.*[0 0.1 -0.2],d,fs); % h32 = miir(1e-7.*[0 -0.6],d,fs); % h33 = miir(1e-5.*[1 0.05 0.3],d,fs); h11 = pzmodel(plist('GAIN', [0.01], 'POLES', [pz(9.9999999999999995e-007,NaN) pz(0.01,NaN) pz(0.02,NaN)], 'ZEROS', [pz(0.0001,NaN) pz(0.002,NaN)])); h12 = pzmodel(plist('GAIN', [0.0001], 'POLES', [pz(0.10000000000000001,NaN) pz(0.0030000000000000001,NaN) pz(0.0050000000000000001,NaN) pz(1.0000000000000001e-005,NaN)], 'ZEROS', [pz(0.00050000000000000001,NaN) pz(0.001,NaN)])); h13 = pzmodel(plist('GAIN', [0.00001], 'POLES', [pz(5.0000000000000004e-006,NaN) pz(0.02,NaN)], 'ZEROS', pz(0.00040000000000000002,NaN))); h21 = pzmodel(plist('GAIN', [0.00001], 'POLES', [pz(5.0000000000000004e-006,NaN) pz(0.10000000000000001,NaN) pz(0.050000000000000003,NaN)], 'ZEROS', [pz(0.002,NaN) pz(0.00020000000000000001,NaN)])); h22 = pzmodel(plist('GAIN', [0.001], 'POLES', [pz(0.01,NaN) pz(9.9999999999999995e-008,NaN) pz(0.002,NaN) pz(0.001,NaN)], 'ZEROS', [pz(1.0000000000000001e-005,NaN) pz(2.0000000000000002e-005,NaN) pz(5.0000000000000002e-005,NaN) pz(0.20000000000000001,NaN)])); h23 = pzmodel(plist('GAIN', [0.0001], 'POLES', [pz(0.01,NaN) pz(0.029999999999999999,NaN) pz(0.10000000000000001,NaN) pz(5.0000000000000002e-005,NaN)], 'ZEROS', [pz(0.0001,NaN) pz(0.0050000000000000001,NaN)])); h31 = pzmodel(plist('GAIN', [0.001], 'POLES', [pz(0.01,NaN) pz(0.029999999999999999,NaN) pz(0.10000000000000001,NaN) pz(0.00050000000000000001,NaN)], 'ZEROS', [pz(0.0001,NaN) pz(0.0050000000000000001,NaN)])); h32 = pzmodel(plist('GAIN', [0.00001], 'POLES', [pz(0.01,NaN) pz(0.029999999999999999,NaN) pz(0.10000000000000001,NaN) pz(0.00050000000000000001,NaN) pz(0.00029999999999999997,NaN)], 'ZEROS', [pz(0.002,NaN) pz(0.059999999999999998,NaN)])); h33 = pzmodel(plist('GAIN', [0.05], 'POLES', [pz(0.0001,NaN) pz(9.9999999999999995e-008,NaN) pz(0.002,NaN)], 'ZEROS', [pz(3.0000000000000001e-005,NaN) pz(5.0000000000000002e-005,NaN) pz(0.10000000000000001,NaN)])); %%% TFs resps % filters response plresp = plist('f',f); rh11 = resp(h11,plresp); rh12 = resp(h12,plresp); rh13 = resp(h13,plresp); rh21 = resp(h21,plresp); rh22 = resp(h22,plresp); rh23 = resp(h23,plresp); rh31 = resp(h31,plresp); rh32 = resp(h32,plresp); rh33 = resp(h33,plresp); % iplot(rh11,rh12,rh13) % iplot(rh21,rh22,rh23) % iplot(rh31,rh32,rh33) % iplot(rh11,rh22,rh33) %%% Spectra calculation with Kay style method % Note that the definition of cross-spectral matrix follows that of s M % Kay, Modern Spectral Estimation fs = 1; % Hz % csd11 G11 = rh11.y.*conj(rh11.y) + rh12.y.*conj(rh12.y) + rh13.y.*conj(rh13.y); G11 = ao(plist('xvals', f, 'yvals', G11, 'fs', fs, 'dtype', 'fsdata', 'yunits', unit('m^2').*unit('Hz^-1'))); G11.setName; % csd12 G12 = conj(rh11.y).*rh21.y + conj(rh12.y).*rh22.y + conj(rh13.y).*rh23.y; G12 = ao(plist('xvals', f, 'yvals', G12, 'fs', fs, 'dtype', 'fsdata', 'yunits', unit('m^2').*unit('Hz^-1'))); G12.setName; % csd13 G13 = conj(rh11.y).*rh31.y + conj(rh12.y).*rh32.y + conj(rh13.y).*rh33.y; G13 = ao(plist('xvals', f, 'yvals', G13, 'fs', fs, 'dtype', 'fsdata', 'yunits', unit('m^2').*unit('Hz^-1'))); G13.setName; % csd21 G21 = conj(G12); G21.setName; % csd22 G22 = rh21.y.*conj(rh21.y) + rh22.y.*conj(rh22.y) + rh23.y.*conj(rh23.y); G22 = ao(plist('xvals', f, 'yvals', G22, 'fs', fs, 'dtype', 'fsdata', 'yunits', unit('m^2').*unit('Hz^-1'))); G22.setName; % csd23 G23 = conj(rh21.y).*rh31.y + conj(rh22.y).*rh32.y + conj(rh23.y).*rh33.y; G23 = ao(plist('xvals', f, 'yvals', G23, 'fs', fs, 'dtype', 'fsdata', 'yunits', unit('m^2').*unit('Hz^-1'))); G23.setName; % csd31 G31 = conj(G13); G31.setName; % csd32 G32 = conj(G23); G32.setName; % csd33 G33 = conj(rh31.y).*rh31.y + conj(rh32.y).*rh32.y + conj(rh33.y).*rh33.y; G33 = ao(plist('xvals', f, 'yvals', G33, 'fs', fs, 'dtype', 'fsdata', 'yunits', unit('m^2').*unit('Hz^-1'))); G33.setName; % iplot(G11,G12,G13,G22,G23,G33) % iplot(G11,G22,G33) %%% Build CSD matrix CSD = zeros(3,3,length(f)); for hh = 1:length(f) CSD(:,:,hh) = [G11.y(hh) G12.y(hh) G13.y(hh);G21.y(hh) G22.y(hh) G23.y(hh);G31.y(hh) G32.y(hh) G33.y(hh)]; end % for hh = 1:length(f) % CSD(:,:,hh) = [G11.y(hh) 0 0;0 G22.y(hh) 0;0 0 G33.y(hh)]; % end %%% 3-dim Eigcsd with kay method stule h = utils.math.ndeigcsd(CSD,'MTD','KAY'); % %%% Eigenshuffle % % [l,m,kk] = size(CSD); % suppr = ones(l,l); % for ii = 2:l % k = ones(l,1); % % for jj = ii-1:-1:1 % k(jj) = min(sqrt(CSD(jj,jj,:)./CSD(ii,ii,:))); % if k(jj)>=1 % suppr(jj,ii) = floor(k(jj)); % else % n=0; % while k(jj)<1 % k(jj)=k(jj)*10; % n=n+1; % end % k(jj) = floor(k(jj)); % suppr(jj,ii) = k(jj)*10^(-n); % % suppr(ii) = suppr(ii)*suppr(ii-1); % end % end % % csuppr(ii) = prod(suppr(:,ii)); % end % csuppr = prod(suppr,2); % supmat = diag(csuppr); % % ssup = rot90(rot90(supmat)); % ssup = supmat; % % ssup = eye(3); % % ssup = supmat*supmat.'; % issup = inv(ssup); % % for jj = 1:kk % nCSD(:,:,jj) = ssup*CSD(:,:,jj)*ssup; % end % % [Vseq,Dseq] = eigenshuffle(nCSD); % % % get h = V*sqrt(D); % [nn,mm,kk] = size(Vseq); % h = zeros(nn,mm,kk); % for dd = 1:kk % % h(:,:,dd) = rot90(issup*Vseq(:,:,dd)*sqrt(diag(Dseq(:,dd)))); % h(:,:,dd) = issup*Vseq(:,:,dd)*sqrt(diag(Dseq(:,dd))); % end %% Build TF AOs nh11 = ao(plist('xvals', f, 'yvals', squeeze(h(1,1,:)), 'fs', fs, 'dtype', 'fsdata')); nh11.setName; nh12 = ao(plist('xvals', f, 'yvals', squeeze(h(1,2,:)), 'fs', fs, 'dtype', 'fsdata')); nh12.setName; nh13 = ao(plist('xvals', f, 'yvals', squeeze(h(1,3,:)), 'fs', fs, 'dtype', 'fsdata')); nh13.setName; nh21 = ao(plist('xvals', f, 'yvals', squeeze(h(2,1,:)), 'fs', fs, 'dtype', 'fsdata')); nh21.setName; nh22 = ao(plist('xvals', f, 'yvals', squeeze(h(2,2,:)), 'fs', fs, 'dtype', 'fsdata')); nh22.setName; nh23 = ao(plist('xvals', f, 'yvals', squeeze(h(2,3,:)), 'fs', fs, 'dtype', 'fsdata')); nh23.setName; nh31 = ao(plist('xvals', f, 'yvals', squeeze(h(3,1,:)), 'fs', fs, 'dtype', 'fsdata')); nh31.setName; nh32 = ao(plist('xvals', f, 'yvals', squeeze(h(3,2,:)), 'fs', fs, 'dtype', 'fsdata')); nh32.setName; nh33 = ao(plist('xvals', f, 'yvals', squeeze(h(3,3,:)), 'fs', fs, 'dtype', 'fsdata')); nh33.setName; iplot(rh11,nh11) iplot(rh12,nh12) iplot(rh13,nh13) iplot(rh21,nh21) iplot(rh22,nh22) iplot(rh23,nh23) iplot(rh31,nh31) iplot(rh32,nh32) iplot(rh33,nh33) %% Build CSD AOs % Note that the definition of cross-spectral matrix follows that of s M % Kay, Modern Spectral Estimation % csd11 H11 = nh11.y.*conj(nh11.y) + nh12.y.*conj(nh12.y) + nh13.y.*conj(nh13.y); H11 = ao(plist('xvals', f, 'yvals', H11, 'fs', fs, 'dtype', 'fsdata', 'yunits', unit('m^2').*unit('Hz^-1'))); H11.setName; % csd12 H12 = conj(nh11.y).*nh21.y + conj(nh12.y).*nh22.y + conj(nh13.y).*nh23.y; H12 = ao(plist('xvals', f, 'yvals', H12, 'fs', fs, 'dtype', 'fsdata', 'yunits', unit('m^2').*unit('Hz^-1'))); H12.setName; % csd13 H13 = conj(nh11.y).*nh31.y + conj(nh12.y).*nh32.y + conj(nh13.y).*nh33.y; H13 = ao(plist('xvals', f, 'yvals', H13, 'fs', fs, 'dtype', 'fsdata', 'yunits', unit('m^2').*unit('Hz^-1'))); H13.setName; % csd21 H21 = conj(H12); H21.setName; % csd22 H22 = nh21.y.*conj(nh21.y) + nh22.y.*conj(nh22.y) + nh23.y.*conj(nh23.y); H22 = ao(plist('xvals', f, 'yvals', H22, 'fs', fs, 'dtype', 'fsdata', 'yunits', unit('m^2').*unit('Hz^-1'))); H22.setName; % csd23 H23 = conj(nh21.y).*nh31.y + conj(nh22.y).*nh32.y + conj(nh23.y).*nh33.y; H23 = ao(plist('xvals', f, 'yvals', H23, 'fs', fs, 'dtype', 'fsdata', 'yunits', unit('m^2').*unit('Hz^-1'))); H23.setName; % csd31 H31 = conj(H13); H31.setName; % csd32 H32 = conj(H23); H32.setName; % csd33 H33 = conj(nh31.y).*nh31.y + conj(nh32.y).*nh32.y + conj(nh33.y).*nh33.y; H33 = ao(plist('xvals', f, 'yvals', H33, 'fs', fs, 'dtype', 'fsdata', 'yunits', unit('m^2').*unit('Hz^-1'))); H33.setName; iplot(G11,H11) iplot(G12,H12) iplot(G13,H13) iplot(G22,H22) iplot(G23,H23) iplot(G33,H33)