Mercurial > hg > ltpda
view testing/utp_1.1/utps/ao/utp_ao_psd.m @ 52:daf4eab1a51e database-connection-manager tip
Fix. Default password should be [] not an empty string
author | Daniele Nicolodi <nicolodi@science.unitn.it> |
---|---|
date | Wed, 07 Dec 2011 17:29:47 +0100 |
parents | 409a22968d5e |
children |
line wrap: on
line source
% UTP_AO_PSD a set of UTPs for the ao/psd method % % M Hewitson 06-08-08 % % $Id: utp_ao_psd.m,v 1.48 2011/07/21 06:46:48 mauro Exp $ % % <MethodDescription> % % The psd method of the ao class computes the spectral density of time-series AOs. % % </MethodDescription> function results = utp_ao_psd(varargin) % Check the inputs if nargin == 0 addpath(fullfile(fileparts(which(mfilename)), 'reference_files')) % Some keywords class = 'ao'; mthd = 'psd'; results = []; disp('******************************************************'); disp(['**** Running UTPs for ' class '/' mthd]); disp('******************************************************'); % Test AOs [at1,at2,at3,at4,at5,at6,atvec,atmat] = eval(['get_test_objects_' class]); % Exception list for the UTPs: [ple1,ple2,ple3,ple4,ple5,ple6] = get_test_ples(); % Get default window from the preferences prefs = getappdata(0, 'LTPDApreferences'); defaultWinType = char(prefs.getMiscPrefs.getDefaultWindow); % Very-long white noise reference time-series a_n_long = ao(plist('filename', 'ref_whitenoise_20000s_10Hz.xml')); % Split the reference data into different segments nsecs = a_n_long.data.nsecs; Nsegments = 100; ap = a_n_long.split(plist('split_type', 'chunks', 'chunks', Nsegments)); % Run the tests results = [results utp_01]; % getInfo call results = [results utp_02]; % Vector input results = [results utp_03]; % Matrix input results = [results utp_04]; % List input results = [results utp_05]; % Test with mixed input results = [results utp_06]; % Test history is working results = [results utp_07]; % Test against MATLAB's pwelch results = [results utp_08]; % Test output of the data results = [results utp_09]; % Test against MATLAB's pwelch results = [results utp_10]; % Test against MATLAB's pwelch results = [results utp_11(mthd, at1, ple1)]; % Test plotinfo doesn't disappear results = [results utp_12]; % Test "conservation of energy": % - white noise from uniform pdf, random parameters % - no detrending, Rectangular window results = [results utp_13]; % Test "conservation of energy": % - white noise from normal pdf, fixed parameters % - no detrending, Rectangular window results = [results utp_14]; % Test "conservation of energy": % - white noise from uniform pdf, fixed parameters % - no detrending, Rectangular window results = [results utp_15]; % Test "conservation of energy": % - white noise from normal pdf, fixed parameters % - no detrending, Rectangular window results = [results utp_16]; % Test "conservation of energy": % - white noise from uniform pdf, fixed parameters % - no detrending, Rectangular window results = [results utp_17]; % Test units handling: PSD results = [results utp_18]; % Test units handling: ASD results = [results utp_19]; % Test units handling: PS results = [results utp_20]; % Test units handling: AS results = [results utp_21]; % Test "conservation of energy": Welch window + no detrending results = [results utp_22]; % Test "conservation of energy": Welch window + mean detrending results = [results utp_23]; % Test "conservation of energy": Welch window + linear detrending results = [results utp_24]; % Test "conservation of energy": Hanning window + no detrending results = [results utp_25]; % Test "conservation of energy": Hanning window + mean detrending results = [results utp_26]; % Test "conservation of energy": Hanning window + linear detrending results = [results utp_27]; % Test "conservation of energy": Hamming window + no detrending results = [results utp_28]; % Test "conservation of energy": Hamming window + mean detrending results = [results utp_29]; % Test "conservation of energy": Hamming window + linear detrending results = [results utp_30]; % Test "conservation of energy": BH92 window + no detrending results = [results utp_31]; % Test "conservation of energy": BH92 window + mean detrending results = [results utp_32]; % Test "conservation of energy": BH92 window + linear detrending results = [results utp_33]; % Test "conservation of energy": Kaiser200 window + no detrending results = [results utp_34]; % Test "conservation of energy": Kaiser200 window + mean detrending results = [results utp_35]; % Test "conservation of energy": Kaiser200 window + linear detrending results = [results utp_38]; % Test detrending results = [results utp_39]; % Test detrending results = [results utp_40]; % Test detrending results = [results utp_41]; % Test different windows: Rectangular results = [results utp_42]; % Test different windows: BH92 results = [results utp_43]; % Test different windows: Hamming results = [results utp_44]; % Test different windows: Hanning results = [results utp_45]; % Test different windows: Bartlett results = [results utp_46]; % Test different windows: Nuttall3 results = [results utp_47]; % Test different windows: Kaiser psll = [random] results = [results utp_48]; % Test different windows: Kaiser psll = [default] results = [results utp_49]; % Test different windows: Nuttall4 results = [results utp_50]; % Test different windows: SFT3F results = [results utp_51]; % Test number of averages: requested/obtained results = [results utp_52]; % Test number of averages: correct number results = [results utp_53]; % Test number of averages: syntax results = [results utp_54]; % Test "conservation of energy": results = [results utp_60]; % Test 'basic' call: PSD results = [results utp_61]; % Test 'basic' call: ASD results = [results utp_62]; % Test 'basic' call: PS results = [results utp_63]; % Test 'basic' call: AS disp('Done.'); disp('******************************************************'); elseif nargin == 1 % Check for UTP functions if strcmp(varargin{1}, 'isutp') results = 1; else results = 0; end else error('### Incorrect inputs') end %% UTP_01 % <TestDescription> % % Tests that the getInfo call works for this method. % % </TestDescription> function result = utp_01 % <SyntaxDescription> % % Test that the getInfo call works for no sets, all sets, and each set % individually. % % </SyntaxDescription> % <SyntaxCode> try % Call for no sets io(1) = eval([class '.getInfo(''' mthd ''', ''None'')']); % Call for all sets io(2) = eval([class '.getInfo(''' mthd ''')']); % Call for each set for kk=1:numel(io(2).sets) io(kk+2) = eval([class '.getInfo(''' mthd ''', ''' io(2).sets{kk} ''')']); end stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that getInfo call returned an minfo object in all cases. % 2) Check that all plists have the correct parameters. % % </AlgoDescription> % <AlgoCode> atest = true; if stest % check we have minfo objects if isa(io, 'minfo') %%% SET 'None' if ~isempty(io(1).sets), atest = false; end if ~isempty(io(1).plists), atest = false; end %%% Check all Sets if ~any(strcmpi(io(2).sets, 'Default')), atest = false; end if numel(io(2).plists) ~= numel(io(2).sets), atest = false; end %%%%%%%%%% SET 'Default' if io(3).plists.nparams ~= 9, atest = false; end % Check key if ~io(3).plists.isparam('nfft'), atest = false; end if ~io(3).plists.isparam('win'), atest = false; end if ~io(3).plists.isparam('olap'), atest = false; end if ~io(3).plists.isparam('scale'), atest = false; end if ~io(3).plists.isparam('order'), atest = false; end if ~io(3).plists.isparam('navs'), atest = false; end if ~io(3).plists.isparam('times'), atest = false; end if ~io(3).plists.isparam('split'), atest = false; end if ~io(3).plists.isparam('psll'), atest = false; end % Check default value if ~isequal(io(3).plists.find('nfft'), -1), atest = false; end if ~strcmpi(io(3).plists.find('win'), defaultWinType), atest = false; end if ~isequal(io(3).plists.find('olap'), -1), atest = false; end if ~isequal(io(3).plists.find('scale'), 'PSD'), atest = false; end if ~isequal(io(3).plists.find('order'), 0), atest = false; end if ~isequal(io(3).plists.find('navs'), -1), atest = false; end if ~isEmptyDouble(io(3).plists.find('times')), atest = false; end if ~isEmptyDouble(io(3).plists.find('split')), atest = false; end if ~isequal(io(3).plists.find('psll'), 200), atest = false; end % Check options if ~isequal(io(3).plists.getOptionsForParam('nfft'), {-1}), atest = false; end if ~isequal(io(3).plists.getOptionsForParam('win'), specwin.getTypes), atest = false; end if ~isequal(io(3).plists.getOptionsForParam('olap'), {-1}), atest = false; end if ~isequal(io(3).plists.getOptionsForParam('scale'), {'PSD', 'ASD', 'PS', 'AS'}), atest = false; end if ~isequal(io(3).plists.getOptionsForParam('order'), {-1 0 1 2 3 4 5 6 7 8 9}), atest = false; end if ~isequal(io(3).plists.getOptionsForParam('navs'), {-1}), atest = false; end if ~isequal(io(3).plists.getOptionsForParam('times'), {[]}), atest = false; end if ~isequal(io(3).plists.getOptionsForParam('split'), {[]}), atest = false; end if ~isequal(io(3).plists.getOptionsForParam('psll'), {200}), atest = false; end end else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_01 %% UTP_02 % <TestDescription> % % Tests that the psd method works with a vector of AOs as input. % % </TestDescription> function result = utp_02 % <SyntaxDescription> % % Test that the psd method works for a vector of AOs as input. % % </SyntaxDescription> % <SyntaxCode> try avec = [at1 at5 at6]; % Vector output out = psd(avec); % List output [out1, out2, out3] = psd(avec); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that the number of elements in 'out' is the same as in the input. % 2) Check that each output object contains the correct values. % % </AlgoDescription> % <AlgoCode> atest = true; if stest % Check we have the correct number of outputs if numel(out) ~= numel(avec), atest = false; end % Check we have the correct values in the outputs % Vector output for kk = 1:numel(out) if ~eq(out(kk), psd(avec(kk)), ple1), atest = false; end end % List output if ~eq(out1, psd(avec(1)), ple1), atest = false; end if ~eq(out2, psd(avec(2)), ple1), atest = false; end if ~eq(out3, psd(avec(3)), ple1), atest = false; end else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_02 %% UTP_03 % <TestDescription> % % Tests that the psd method works with a matrix of AOs as input. % % </TestDescription> function result = utp_03 % <SyntaxDescription> % % Test that the psd method works for a matrix of AOs as input. % % </SyntaxDescription> % <SyntaxCode> try amat = [at1 at5 at6; at5 at6 at1]; % Vector output out = psd(amat); % List output [out1, out2, out3, out4, out5, out6] = psd(amat); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that the number of elements in 'out' is the same as in the input. % 2) Check that each output object contains the correct values. % % </AlgoDescription> % <AlgoCode> atest = true; if stest % Check we have the correct number of outputs if numel(out) ~= numel(amat), atest = false; end % Check we have the correct values in the outputs % Vector output for kk = 1:numel(out) if ~eq(out(kk), psd(amat(kk)), ple1), atest = false; end end % List output if ~eq(out1, psd(amat(1)), ple1), atest = false; end if ~eq(out2, psd(amat(2)), ple1), atest = false; end if ~eq(out3, psd(amat(3)), ple1), atest = false; end if ~eq(out4, psd(amat(4)), ple1), atest = false; end if ~eq(out5, psd(amat(5)), ple1), atest = false; end if ~eq(out6, psd(amat(6)), ple1), atest = false; end else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_03 %% UTP_04 % <TestDescription> % % Tests that the psd method works with a list of AOs as input. % % </TestDescription> function result = utp_04 % <SyntaxDescription> % % Test that the psd method works for a list of AOs as input. % % </SyntaxDescription> % <SyntaxCode> try % Vector output out = psd(at1,at5,at6); % List output [out1, out2, out3] = psd(at1,at5,at6); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that the number of elements in 'out' is the same as in the input. % 2) Check that each output AO contains the correct data. % % </AlgoDescription> % <AlgoCode> atest = true; if stest % Check we have the correct number of outputs if numel(out) ~= 3, atest = false; end % Check we have the correct values in the outputs % Vector output if ~eq(out(1), psd(at1), ple1), atest = false; end if ~eq(out(2), psd(at5), ple1), atest = false; end if ~eq(out(3), psd(at6), ple1), atest = false; end % List output if ~eq(out1, psd(at1), ple1), atest = false; end if ~eq(out2, psd(at5), ple1), atest = false; end if ~eq(out3, psd(at6), ple1), atest = false; end else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_04 %% UTP_05 % <TestDescription> % % Tests that the psd method works with a mix of different shaped AOs as % input. % % </TestDescription> function result = utp_05 % <SyntaxDescription> % % Test that the psd method works with an input of matrices and vectors % and single AOs. % % </SyntaxDescription> % <SyntaxCode> try % Vector output out = psd(at1,[at5 at6],at5,[at5 at1; at6 at1],at6); % List output [out1, out2, out3, out4, out5, out6, out7, out8, out9] = ... psd(at1,[at5 at6],at5,[at5 at1; at6 at1],at6); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that the number of elements in 'out' is the same as in % input. % 2) Check that each output AO contains the correct data. % % </AlgoDescription> % <AlgoCode> atest = true; if stest % Check we have the correct number of outputs if numel(out) ~= 9, atest = false; end % Check we have the correct values in the outputs % Vector output if ~eq(out(1), psd(at1), ple1), atest = false; end if ~eq(out(2), psd(at5), ple1), atest = false; end if ~eq(out(3), psd(at6), ple1), atest = false; end if ~eq(out(4), psd(at5), ple1), atest = false; end if ~eq(out(5), psd(at5), ple1), atest = false; end if ~eq(out(6), psd(at6), ple1), atest = false; end if ~eq(out(7), psd(at1), ple1), atest = false; end if ~eq(out(8), psd(at1), ple1), atest = false; end if ~eq(out(9), psd(at6), ple1), atest = false; end % List output if ~eq(out1, psd(at1), ple1), atest = false; end if ~eq(out2, psd(at5), ple1), atest = false; end if ~eq(out3, psd(at6), ple1), atest = false; end if ~eq(out4, psd(at5), ple1), atest = false; end if ~eq(out5, psd(at5), ple1), atest = false; end if ~eq(out6, psd(at6), ple1), atest = false; end if ~eq(out7, psd(at1), ple1), atest = false; end if ~eq(out8, psd(at1), ple1), atest = false; end if ~eq(out9, psd(at6), ple1), atest = false; end else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_05 %% UTP_06 % <TestDescription> % % Tests that the psd method properly applies history. % % </TestDescription> function result = utp_06 % <SyntaxDescription> % % Test that the result of applying the psd method can be processed back % to an m-file. % % </SyntaxDescription> % <SyntaxCode> try out = psd(at5); mout = rebuild(out); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that the last entry in the history of 'out' corresponds to % 'psd'. % 2) Check that the re-built object is the same object as 'out'. % % </AlgoDescription> % <AlgoCode> atest = true; if stest % Check the last step in the history of 'out' if ~strcmp(out.hist.methodInfo.mname, 'psd'), atest = false; end % Check the re-built object if ~eq(mout, out, ple2), atest = false; end else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_06 %% UTP_07 % <TestDescription> % % Tests that the psd method agrees with MATLAB's pwelch, within some tolerance, % when configured to use the same parameters. % % </TestDescription> function result = utp_07 % <SyntaxDescription> % % Test that applying psd works on a single AO. % % </SyntaxDescription> % <SyntaxCode> try % Load reference ao a1 = ao('ref_whitenoise_10s_1000Hz.xml'); fs = a1.fs; % Compute PSD Nfft = 2*fs; win = specwin('Hanning', Nfft); pl = plist('Nfft', Nfft, 'Win', win.type, 'order', -1); S1 = a1.psd(pl); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that output agrees with the output of MATLAB's pwelch within tolerance. % % </AlgoDescription> % <AlgoCode> atest = true; TOL = 1.8e-15; if stest % Load reference psd ref = load('ref_psd_10s_1000Hz.txt'); if ~isequal(ref(:,1), S1.x) || any((abs(S1.y - ref(:,2)) ./ S1.y) >= TOL) atest = false; end else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_07 %% UTP_08 % <TestDescription> % % Check that the psd method pass back the output objects to a list of % output variables or to a single variable. % % </TestDescription> function result = utp_08 % <SyntaxDescription> % % Call the method with a list of output variables and with a single output % variable. Additionaly check that the rebuild method works on the output. % % </SyntaxDescription> try % <SyntaxCode> [o1, o2] = psd(at5, at6); o3 = psd(at5, at6); mout1 = rebuild(o1); mout2 = rebuild(o2); mout3 = rebuild(o3); % </SyntaxCode> stest = true; catch err disp(err.message) stest = false; end % <AlgoDescription> % % 1) Check that the output contains the right number of objects % 2) Check that the 'rebuild' method produces the same object as 'out'. % % </AlgoDescription> atest = true; if stest % <AlgoCode> % Check the number of outputs if numel(o1) ~=1, atest = false; end if numel(o2) ~=1, atest = false; end if numel(o3) ~=2, atest = false; end % Check the rebuilding of the object if ~eq(o1, mout1, ple2), atest = false; end if ~eq(o2, mout2, ple2), atest = false; end if ~eq(o3, mout3, ple2), atest = false; end % </AlgoCode> else atest = false; end % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_08 %% UTP_09 % <TestDescription> % % Tests that the psd method agrees with MATLAB's pwelch when % configured to use the same parameters. % % </TestDescription> function result = utp_09 % <SyntaxDescription> % % Test that the applying psd works on a single AO. % % </SyntaxDescription> % <SyntaxCode> try % Load reference ao a1 = ao('ref_whitenoise_3600s_1Hz.xml'); % Compute PSD Nfft = 1000; win = specwin('BH92', Nfft); pl = plist('Nfft', Nfft, 'Win', win.type, 'order', -1); S1 = a1.psd(pl); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that output agrees with the output of MATLAB's pwelch. % % </AlgoDescription> % <AlgoCode> atest = true; TOL = 1e-15; if stest % Load reference psd ref = load('ref_psd_3600s_1Hz.txt'); if ~isequal(ref(:,1), S1.x) || any((abs(S1.y - ref(:,2)) ./ S1.y) >= TOL) atest = false; end else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_09 %% UTP_10 % <TestDescription> % % Tests that the psd method agrees with MATLAB's pwelch when % configured to use the same parameters. % % </TestDescription> function result = utp_10 % <SyntaxDescription> % % Test that the applying psd works on a single AO. % % </SyntaxDescription> % <SyntaxCode> try % Load reference ao a1 = ao('ref_whitenoise_100000s_100mHz.xml'); % Compute PSD Nfft = a1.fs * a1.data.nsecs; psll = 100; win = specwin('Kaiser', Nfft, psll); pl = plist('Nfft', Nfft, 'Win', win.type, 'psll', psll, 'order', -1); S1 = a1.psd(pl); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that output agrees with the output of MATLAB's pwelch. % % </AlgoDescription> % <AlgoCode> atest = true; TOL = 1e-13; if stest % Load reference psd ref = load('ref_psd_100000s_100mHz.txt'); if ~isequal(ref(:,1), S1.x) || any((abs(S1.y - ref(:,2)) ./ S1.y) >= TOL) atest = false; end else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_10 %% UTP_12 % <TestDescription> % % Tests "conservation of energy": % 1) white noise produced from uniform pdf, with a given mean value and % sigma (distribution's 1st and 2nd order) % 2) evaluate the sample mean value m and standard deviation s % 3) psd of the white noise % 4) compares the sqrt(sum(S * df)) with the standard deviation s % % </TestDescription> function result = utp_12 % <SyntaxDescription> % % 1) Prepare the test tsdata: % white noise from uniform distribution + offset % 2) Calculate the statistical parameters % 3) Estimate the psd % % </SyntaxDescription> % <SyntaxCode> try % White noise % Parameters are crucial for the estimate to be correct! noise_type = 'Uniform'; win_type = 'Rectangular'; olap = 0; detrend_order = 0; scale = 'PSD'; [sigma_sample, sigma_psd] = prepare_analyze_noise_energy(win_type, noise_type, scale, ... plist('detrend_order', detrend_order, 'olap', olap)); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that calculated psd energy equals the rms content of the tsdata, estimated by the std of the sample % % </AlgoDescription> % <AlgoCode> atest = true; TOL = 1e-12; if stest if abs(sigma_psd - sigma_sample) / mean([sigma_psd sigma_sample]) >= TOL atest = false; end else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_12 %% UTP_13 % <TestDescription> % % Tests "conservation of energy" with fixed parameters: % 1) white noise produced from normal pdf, with a given mean value and % sigma (distribution's 1st and 2nd order) % 2) evaluate the sample mean value m and standard deviation s % 3) psd of the white noise % 4) compares the sqrt(sum(S * df)) with the standard deviation s % % </TestDescription> function result = utp_13 % <SyntaxDescription> % % 1) Prepare the test tsdata: % white noise from normal distribution + offset % 2) Calculate the statistical parameters % 3) Estimate the psd % % </SyntaxDescription> % <SyntaxCode> try % White noise % Parameters are crucial for the estimate to be correct! noise_type = 'Normal'; win_type = 'Rectangular'; olap = 0; detrend_order = 0; scale = 'PSD'; % Build with fixed parameters fs = 10; nsecs = 3600; sigma_distr = 1; mu_distr = 0; [sigma_sample, sigma_psd] = prepare_analyze_noise_energy(win_type, noise_type, scale, ... plist('detrend_order', detrend_order, 'olap', olap, ... 'fs', fs, 'nsecs', nsecs, 'sigma_distr', sigma_distr, 'mu_distr', mu_distr)); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that calculated psd energy equals the rms content of the tsdata, estimated by the std of the sample % % </AlgoDescription> % <AlgoCode> atest = true; TOL = 1e-12; if stest if abs(sigma_psd - sigma_sample) / mean([sigma_psd sigma_sample]) >= TOL atest = false; end else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_13 %% UTP_14 % <TestDescription> % % Tests "conservation of energy" with fixed parameters: % 1) white noise produced from uniform pdf, with a given mean value and % sigma (distribution's 1st and 2nd order) % 2) evaluate the sample mean value m and standard deviation s % 3) psd of the white noise % 4) compares the sqrt(sum(S * df)) with the standard deviation s % % </TestDescription> function result = utp_14 % <SyntaxDescription> % % 1) Prepare the test tsdata: % white noise from uniform distribution + offset % 2) Calculate the statistical parameters % 3) Estimate the psd % % </SyntaxDescription> % <SyntaxCode> try % White noise % Parameters are crucial for the estimate to be correct! noise_type = 'Uniform'; win_type = 'Rectangular'; olap = 0; detrend_order = 0; scale = 'PSD'; % Build with fixed parameters fs = 1; nsecs = 86400; sigma_distr = 1e-9; mu_distr = 3e-7; [sigma_sample, sigma_psd] = prepare_analyze_noise_energy(win_type, noise_type, scale, ... plist('detrend_order', detrend_order, 'olap', olap, ... 'fs', fs, 'nsecs', nsecs, 'sigma_distr', sigma_distr, 'mu_distr', mu_distr)); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that calculated psd energy equals the rms content of the tsdata, estimated by the std of the sample % % </AlgoDescription> % <AlgoCode> atest = true; TOL = 1e-12; if stest if abs(sigma_psd - sigma_sample) / mean([sigma_psd sigma_sample]) >= TOL atest = false; end else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_14 %% UTP_15 % <TestDescription> % % Tests "conservation of energy" with fixed parameters: % 1) white noise produced from normal pdf, with a given mean value and % sigma (distribution's 1st and 2nd order) % 2) evaluate the sample mean value m and standard deviation s % 3) psd of the white noise % 4) compares the sqrt(sum(S * df)) with the standard deviation s % % </TestDescription> function result = utp_15 % <SyntaxDescription> % % 1) Prepare the test tsdata: % white noise from normal distribution + offset % 2) Calculate the statistical parameters % 3) Estimate the psd % % </SyntaxDescription> % <SyntaxCode> try % White noise % Parameters are crucial for the estimate to be correct! noise_type = 'Normal'; win_type = 'Rectangular'; olap = 0; detrend_order = 0; scale = 'PSD'; % Build with fixed parameters fs = 100; nsecs = 100; sigma_distr = 1e-12; mu_distr = -5e-12; [sigma_sample, sigma_psd] = prepare_analyze_noise_energy(win_type, noise_type, scale, ... plist('detrend_order', detrend_order, 'olap', olap, ... 'fs', fs, 'nsecs', nsecs, 'sigma_distr', sigma_distr, 'mu_distr', mu_distr)); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that calculated psd energy equals the rms content of the tsdata, estimated by the std of the sample % % </AlgoDescription> % <AlgoCode> atest = true; TOL = 1e-12; if stest if abs(sigma_psd - sigma_sample) / mean([sigma_psd sigma_sample]) >= TOL atest = false; end else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_15 %% UTP_16 % <TestDescription> % % Tests "conservation of energy" with fixed parameters: % 1) white noise produced from uniform pdf, with a given mean value and % sigma (distribution's 1st and 2nd order) % 2) evaluate the sample mean value m and standard deviation s % 3) psd of the white noise % 4) compares the sqrt(sum(S * df)) with the standard deviation s % % </TestDescription> function result = utp_16 % <SyntaxDescription> % % 1) Prepare the test tsdata: % white noise from uniform distribution + offset % 2) Calculate the statistical parameters % 3) Estimate the psd % % </SyntaxDescription> % <SyntaxCode> try % White noise % Parameters are crucial for the estimate to be correct! noise_type = 'Uniform'; win_type = 'Rectangular'; olap = 0; detrend_order = 0; scale = 'PSD'; % Build with fixed parameters fs = 0.01; nsecs = 3*86400; sigma_distr = 1e-15; mu_distr = -7.72e-13; [sigma_sample, sigma_psd] = prepare_analyze_noise_energy(win_type, noise_type, scale, ... plist('detrend_order', detrend_order, 'olap', olap, ... 'fs', fs, 'nsecs', nsecs, 'sigma_distr', sigma_distr, 'mu_distr', mu_distr)); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that calculated psd energy equals the rms content of the tsdata, estimated by the std of the sample % % </AlgoDescription> % <AlgoCode> atest = true; TOL = 1e-12; if stest if abs(sigma_psd - sigma_sample) / mean([sigma_psd sigma_sample]) >= TOL atest = false; end else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_16 %% UTP_17 % <TestDescription> % % Tests handling of units: % 1) white noise produced from normal pdf, with a given mean value and % sigma (distribution's 1st and 2nd orders) % 2) PSD of the white noise % 3) compares the units of the input and output % % </TestDescription> function result = utp_17 % <SyntaxDescription> % % 1) Prepare the test tsdata: % white noise from normal distribution + offset % 2) Assign a random unit % 3) PSD of the white noise % % </SyntaxDescription> % <SyntaxCode> try % White noise noise_type = 'Normal'; win_type = 'BH92'; reduce_pts = 10; scale = 'PSD'; olap = specwin(win_type).rov; [a, S, S1, S2] = prepare_analyze_noise_win(win_type, noise_type, scale, ... plist('olap', olap, 'reduce_pts', reduce_pts)); % S1 and S2 are empty in this case stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that (calculated PSD yunits) equals (input units)^2 / Hz % 2) Check that (calculated PSD xunits) equals Hz % % </AlgoDescription> % <AlgoCode> atest = true; if stest if ne(S.yunits, (a.yunits).^2 * unit('Hz^-1')) || ne(S.xunits, unit('Hz')) atest = false; end else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_17 %% UTP_18 % <TestDescription> % % Tests handling of units: % 1) white noise produced from uniform pdf, with a given mean value and % sigma (distribution's 1st and 2nd orders) % 2) ASD of the white noise % 3) compares the units of the input and output % % </TestDescription> function result = utp_18 % <SyntaxDescription> % % 1) Prepare the test tsdata: % white noise from uniform distribution + offset % 2) Assign a random unit % 3) ASD of the white noise % % </SyntaxDescription> % <SyntaxCode> try % White noise noise_type = 'Uniform'; win_type = 'Hamming'; reduce_pts = 4; scale = 'ASD'; olap = specwin(win_type).rov; [a, S, S1, S2] = prepare_analyze_noise_win(win_type, noise_type, scale, ... plist('olap', olap, 'reduce_pts', reduce_pts)); % S1 and S2 are empty in this case stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that (calculated ASD yunits) equals (input units) / Hz^(1/2) % 2) Check that (calculated ASD xunits) equals Hz % % </AlgoDescription> % <AlgoCode> atest = true; if stest if ne(S.yunits, (a.yunits) * unit('Hz^-1/2')) || ne(S.xunits, unit('Hz')) atest = false; end else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_18 %% UTP_19 % <TestDescription> % % Tests handling of units: % 1) white noise produced from normal pdf, with a given mean value and % sigma (distribution's 1st and 2nd orders) % 2) PS of the white noise % 3) compares the units of the input and output % % </TestDescription> function result = utp_19 % <SyntaxDescription> % % 1) Prepare the test tsdata: % white noise from normal distribution + offset % 2) Assign a random unit % 3) PS of the white noise % % </SyntaxDescription> % <SyntaxCode> try % White noise noise_type = 'Normal'; win_type = 'Hanning'; reduce_pts = 5; scale = 'PS'; olap = specwin(win_type).rov; [a, S, S1, S2] = prepare_analyze_noise_win(win_type, noise_type, scale, ... plist('olap', olap, 'reduce_pts', reduce_pts)); % S1 and S2 are empty in this case stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that (calculated PS yunits) equals (input units)^2 % 2) Check that (calculated PS xunits) equals Hz % % </AlgoDescription> % <AlgoCode> atest = true; if stest if ne(S.yunits, (a.yunits).^2) || ne(S.xunits, unit('Hz')) atest = false; end else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_19 %% UTP_20 % <TestDescription> % % Tests handling of units: % 1) white noise produced from uniform distribution, with a given mean value and % sigma (distribution's 1st and 2nd orders) % 2) AS of the white noise % 3) compares the units of the input and output % % </TestDescription> function result = utp_20 % <SyntaxDescription> % % 1) Prepare the test tsdata: % white noise from uniform distribution + offset % 2) Assign a random unit % 3) AS of the white noise % % </SyntaxDescription> % <SyntaxCode> try % White noise noise_type = 'Uniform'; win_type = 'Rectangular'; reduce_pts = 10; scale = 'AS'; olap = specwin(win_type).rov; [a, S, S1, S2] = prepare_analyze_noise_win(win_type, noise_type, scale, ... plist('olap', olap, 'reduce_pts', reduce_pts)); % S1 and S2 are empty in this case stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that (calculated AS yunits) equals (input units) % 2) Check that (calculated AS xunits) equals Hz % % </AlgoDescription> % <AlgoCode> atest = true; if stest if ne(S.yunits, a.yunits) || ne(S.xunits, unit('Hz')) atest = false; end else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_20 %% UTP_21 % <TestDescription> % % Tests "conservation of energy": % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) Calculate the expected level of noise from sample statistics % 3) Calculates the PSD of the individual parts, with: % Welch window and no detrend % 4) Evaluate the mean and standard deviation of the Checks on individual PSDs: % 5) Checks on individual PSDs: mean and standard deviation of the PSD points % 6) Evaluate the expected value, estimated from the std of the % individual segments % 7) Compares with the mean performed on the individual segments % 8) Checks on averaged PSDs: mean and standard deviation of all the PSD points % 9) Compares the grand-mean with the estimated white noise level % </TestDescription> function result = utp_21 % <SyntaxDescription> % % 1) Split the reference data into different segments % 2) Calculates the PSD of the individual parts % % </SyntaxDescription> % <SyntaxCode> try % Calculate the expected level of noise from sample statistics % Evaluate the PSDs on the parts, employing: % - Welch window % - no detrend win_type = 'Welch'; detrend_order = -1; [S_ap, S_sigma_total, S_expected_total] = ... test_psd_energy_prepare_data(win_type, detrend_order, plist); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Evaluate the mean and standard deviation of the Checks on individual PSDs: % 2) Checks on individual PSDs: mean and standard deviation of the PSD points % 3) Evaluate the expected value, estimated from the std of the % individual segments % 4) Compares with the mean performed on the individual segments % 5) Checks on averaged PSDs: mean and standard deviation of all the PSD points % 6) Compares the grand-mean with the estimated white noise level % </AlgoDescription> % <AlgoCode> if stest atest = algo_test_psd_energy(ap, S_ap, S_expected_total, plist); else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_21 %% UTP_22 % <TestDescription> % % Tests "conservation of energy": % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) Calculate the expected level of noise from sample statistics % 3) Calculates the PSD of the individual parts, with: % Welch window and mean detrend % 4) Evaluate the mean and standard deviation of the Checks on individual PSDs: % 5) Checks on individual PSDs: mean and standard deviation of the PSD points % 6) Evaluate the expected value, estimated from the std of the % individual segments % 7) Compares with the mean performed on the individual segments % 8) Checks on averaged PSDs: mean and standard deviation of all the PSD points % 9) Compares the grand-mean with the estimated white noise level % </TestDescription> function result = utp_22 % <SyntaxDescription> % % 1) Split the reference data into different segments % 2) Calculates the PSD of the individual parts % % </SyntaxDescription> % <SyntaxCode> try % Calculate the expected level of noise from sample statistics % Evaluate the PSDs on the parts, employing: % - Welch window % - mean detrend win_type = 'Welch'; detrend_order = 0; [S_ap, S_sigma_total, S_expected_total] = ... test_psd_energy_prepare_data(win_type, detrend_order, plist); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Evaluate the mean and standard deviation of the Checks on individual PSDs: % 2) Checks on individual PSDs: mean and standard deviation of the PSD points % 3) Evaluate the expected value, estimated from the std of the % individual segments % 4) Compares with the mean performed on the individual segments % 5) Checks on averaged PSDs: mean and standard deviation of all the PSD points % 6) Compares the grand-mean with the estimated white noise level % </AlgoDescription> % <AlgoCode> if stest atest = algo_test_psd_energy(ap, S_ap, S_expected_total, plist); else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_22 %% UTP_23 % <TestDescription> % % Tests "conservation of energy": % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) Calculate the expected level of noise from sample statistics % 3) Calculates the PSD of the individual parts, with: % Welch window and linear detrend % 4) Evaluate the mean and standard deviation of the Checks on individual PSDs: % 5) Checks on individual PSDs: mean and standard deviation of the PSD points % 6) Evaluate the expected value, estimated from the std of the % individual segments % 7) Compares with the mean performed on the individual segments % 8) Checks on averaged PSDs: mean and standard deviation of all the PSD points % 9) Compares the grand-mean with the estimated white noise level % </TestDescription> function result = utp_23 % <SyntaxDescription> % % 1) Split the reference data into different segments % 2) Calculates the PSD of the individual parts % % </SyntaxDescription> % <SyntaxCode> try % Calculate the expected level of noise from sample statistics % Evaluate the PSDs on the parts, employing: % - Welch window % - linear detrend win_type = 'Welch'; detrend_order = 1; [S_ap, S_sigma_total, S_expected_total] = ... test_psd_energy_prepare_data(win_type, detrend_order, plist); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Evaluate the mean and standard deviation of the Checks on individual PSDs: % 2) Checks on individual PSDs: mean and standard deviation of the PSD points % 3) Evaluate the expected value, estimated from the std of the % individual segments % 4) Compares with the mean performed on the individual segments % 5) Checks on averaged PSDs: mean and standard deviation of all the PSD points % 6) Compares the grand-mean with the estimated white noise level % </AlgoDescription> % <AlgoCode> if stest atest = algo_test_psd_energy(ap, S_ap, S_expected_total, plist); else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_23 %% UTP_24 % <TestDescription> % % Tests "conservation of energy": % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) Calculate the expected level of noise from sample statistics % 3) Calculates the PSD of the individual parts, with: % Hanning window and no detrend % 4) Evaluate the mean and standard deviation of the Checks on individual PSDs: % 5) Checks on individual PSDs: mean and standard deviation of the PSD points % 6) Evaluate the expected value, estimated from the std of the % individual segments % 7) Compares with the mean performed on the individual segments % 8) Checks on averaged PSDs: mean and standard deviation of all the PSD points % 9) Compares the grand-mean with the estimated white noise level % </TestDescription> function result = utp_24 % <SyntaxDescription> % % 1) Split the reference data into different segments % 2) Calculates the PSD of the individual parts % % </SyntaxDescription> % <SyntaxCode> try % Calculate the expected level of noise from sample statistics % Evaluate the PSDs on the parts, employing: % - Hanning window % - no detrend win_type = 'Hanning'; detrend_order = -1; [S_ap, S_sigma_total, S_expected_total] = ... test_psd_energy_prepare_data(win_type, detrend_order, plist); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Evaluate the mean and standard deviation of the Checks on individual PSDs: % 2) Checks on individual PSDs: mean and standard deviation of the PSD points % 3) Evaluate the expected value, estimated from the std of the % individual segments % 4) Compares with the mean performed on the individual segments % 5) Checks on averaged PSDs: mean and standard deviation of all the PSD points % 6) Compares the grand-mean with the estimated white noise level % </AlgoDescription> % <AlgoCode> if stest atest = algo_test_psd_energy(ap, S_ap, S_expected_total, plist); else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_24 %% UTP_25 % <TestDescription> % % Tests "conservation of energy": % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) Calculate the expected level of noise from sample statistics % 3) Calculates the PSD of the individual parts, with: % Hanning window and mean detrend % 4) Evaluate the mean and standard deviation of the Checks on individual PSDs: % 5) Checks on individual PSDs: mean and standard deviation of the PSD points % 6) Evaluate the expected value, estimated from the std of the % individual segments % 7) Compares with the mean performed on the individual segments % 8) Checks on averaged PSDs: mean and standard deviation of all the PSD points % 9) Compares the grand-mean with the estimated white noise level % </TestDescription> function result = utp_25 % <SyntaxDescription> % % 1) Split the reference data into different segments % 2) Calculates the PSD of the individual parts % % </SyntaxDescription> % <SyntaxCode> try % Calculate the expected level of noise from sample statistics % Evaluate the PSDs on the parts, employing: % - Hanning window % - mean detrend win_type = 'Hanning'; detrend_order = 0; [S_ap, S_sigma_total, S_expected_total] = ... test_psd_energy_prepare_data(win_type, detrend_order, plist); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Evaluate the mean and standard deviation of the Checks on individual PSDs: % 2) Checks on individual PSDs: mean and standard deviation of the PSD points % 3) Evaluate the expected value, estimated from the std of the % individual segments % 4) Compares with the mean performed on the individual segments % 5) Checks on averaged PSDs: mean and standard deviation of all the PSD points % 6) Compares the grand-mean with the estimated white noise level % </AlgoDescription> % <AlgoCode> if stest atest = algo_test_psd_energy(ap, S_ap, S_expected_total, plist); else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_25 %% UTP_26 % <TestDescription> % % Tests "conservation of energy": % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) Calculate the expected level of noise from sample statistics % 3) Calculates the PSD of the individual parts, with: % Hanning window and linear detrend % 4) Evaluate the mean and standard deviation of the Checks on individual PSDs: % 5) Checks on individual PSDs: mean and standard deviation of the PSD points % 6) Evaluate the expected value, estimated from the std of the % individual segments % 7) Compares with the mean performed on the individual segments % 8) Checks on averaged PSDs: mean and standard deviation of all the PSD points % 9) Compares the grand-mean with the estimated white noise level % </TestDescription> function result = utp_26 % <SyntaxDescription> % % 1) Split the reference data into different segments % 2) Calculates the PSD of the individual parts % % </SyntaxDescription> % <SyntaxCode> try % Calculate the expected level of noise from sample statistics % Evaluate the PSDs on the parts, employing: % - Hanning window % - linear detrend win_type = 'Hanning'; detrend_order = 1; [S_ap, S_sigma_total, S_expected_total] = ... test_psd_energy_prepare_data(win_type, detrend_order, plist); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Evaluate the mean and standard deviation of the Checks on individual PSDs: % 2) Checks on individual PSDs: mean and standard deviation of the PSD points % 3) Evaluate the expected value, estimated from the std of the % individual segments % 4) Compares with the mean performed on the individual segments % 5) Checks on averaged PSDs: mean and standard deviation of all the PSD points % 6) Compares the grand-mean with the estimated white noise level % </AlgoDescription> % <AlgoCode> if stest atest = algo_test_psd_energy(ap, S_ap, S_expected_total, plist); else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_26 %% UTP_27 % <TestDescription> % % Tests "conservation of energy": % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) Calculate the expected level of noise from sample statistics % 3) Calculates the PSD of the individual parts, with: % Hamming window and no detrend % 4) Evaluate the mean and standard deviation of the Checks on individual PSDs: % 5) Checks on individual PSDs: mean and standard deviation of the PSD points % 6) Evaluate the expected value, estimated from the std of the % individual segments % 7) Compares with the mean performed on the individual segments % 8) Checks on averaged PSDs: mean and standard deviation of all the PSD points % 9) Compares the grand-mean with the estimated white noise level % </TestDescription> function result = utp_27 % <SyntaxDescription> % % 1) Split the reference data into different segments % 2) Calculates the PSD of the individual parts % % </SyntaxDescription> % <SyntaxCode> try % Calculate the expected level of noise from sample statistics % Evaluate the PSDs on the parts, employing: % - Hamming window % - no detrend win_type = 'Hamming'; detrend_order = -1; [S_ap, S_sigma_total, S_expected_total] = ... test_psd_energy_prepare_data(win_type, detrend_order, plist); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Evaluate the mean and standard deviation of the Checks on individual PSDs: % 2) Checks on individual PSDs: mean and standard deviation of the PSD points % 3) Evaluate the expected value, estimated from the std of the % individual segments % 4) Compares with the mean performed on the individual segments % 5) Checks on averaged PSDs: mean and standard deviation of all the PSD points % 6) Compares the grand-mean with the estimated white noise level % </AlgoDescription> % <AlgoCode> if stest atest = algo_test_psd_energy(ap, S_ap, S_expected_total, plist); else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_27 %% UTP_28 % <TestDescription> % % Tests "conservation of energy": % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) Calculate the expected level of noise from sample statistics % 3) Calculates the PSD of the individual parts, with: % Hamming window and mean detrend % 4) Evaluate the mean and standard deviation of the Checks on individual PSDs: % 5) Checks on individual PSDs: mean and standard deviation of the PSD points % 6) Evaluate the expected value, estimated from the std of the % individual segments % 7) Compares with the mean performed on the individual segments % 8) Checks on averaged PSDs: mean and standard deviation of all the PSD points % 9) Compares the grand-mean with the estimated white noise level % </TestDescription> function result = utp_28 % <SyntaxDescription> % % 1) Split the reference data into different segments % 2) Calculates the PSD of the individual parts % % </SyntaxDescription> % <SyntaxCode> try % Calculate the expected level of noise from sample statistics % Evaluate the PSDs on the parts, employing: % - Hamming window % - mean detrend win_type = 'Hanning'; detrend_order = 0; [S_ap, S_sigma_total, S_expected_total] = ... test_psd_energy_prepare_data(win_type, detrend_order, plist); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Evaluate the mean and standard deviation of the Checks on individual PSDs: % 2) Checks on individual PSDs: mean and standard deviation of the PSD points % 3) Evaluate the expected value, estimated from the std of the % individual segments % 4) Compares with the mean performed on the individual segments % 5) Checks on averaged PSDs: mean and standard deviation of all the PSD points % 6) Compares the grand-mean with the estimated white noise level % </AlgoDescription> % <AlgoCode> if stest atest = algo_test_psd_energy(ap, S_ap, S_expected_total, plist); else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_28 %% UTP_29 % <TestDescription> % % Tests "conservation of energy": % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) Calculate the expected level of noise from sample statistics % 3) Calculates the PSD of the individual parts, with: % Hamming window and linear detrend % 4) Evaluate the mean and standard deviation of the Checks on individual PSDs: % 5) Checks on individual PSDs: mean and standard deviation of the PSD points % 6) Evaluate the expected value, estimated from the std of the % individual segments % 7) Compares with the mean performed on the individual segments % 8) Checks on averaged PSDs: mean and standard deviation of all the PSD points % 9) Compares the grand-mean with the estimated white noise level % </TestDescription> function result = utp_29 % <SyntaxDescription> % % 1) Split the reference data into different segments % 2) Calculates the PSD of the individual parts % % </SyntaxDescription> % <SyntaxCode> try % Calculate the expected level of noise from sample statistics % Evaluate the PSDs on the parts, employing: % - Hamming window % - linear detrend win_type = 'Hamming'; detrend_order = 1; [S_ap, S_sigma_total, S_expected_total] = ... test_psd_energy_prepare_data(win_type, detrend_order, plist); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Evaluate the mean and standard deviation of the Checks on individual PSDs: % 2) Checks on individual PSDs: mean and standard deviation of the PSD points % 3) Evaluate the expected value, estimated from the std of the % individual segments % 4) Compares with the mean performed on the individual segments % 5) Checks on averaged PSDs: mean and standard deviation of all the PSD points % 6) Compares the grand-mean with the estimated white noise level % </AlgoDescription> % <AlgoCode> if stest atest = algo_test_psd_energy(ap, S_ap, S_expected_total, plist); else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_29 %% UTP_30 % <TestDescription> % % Tests "conservation of energy": % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) Calculate the expected level of noise from sample statistics % 3) Calculates the PSD of the individual parts, with: % BH92 window and no detrend % 4) Evaluate the mean and standard deviation of the Checks on individual PSDs: % 5) Checks on individual PSDs: mean and standard deviation of the PSD points % 6) Evaluate the expected value, estimated from the std of the % individual segments % 7) Compares with the mean performed on the individual segments % 8) Checks on averaged PSDs: mean and standard deviation of all the PSD points % 9) Compares the grand-mean with the estimated white noise level % </TestDescription> function result = utp_30 % <SyntaxDescription> % % 1) Split the reference data into different segments % 2) Calculates the PSD of the individual parts % % </SyntaxDescription> % <SyntaxCode> try % Calculate the expected level of noise from sample statistics % Evaluate the PSDs on the parts, employing: % - BH92 window % - no detrend win_type = 'BH92'; detrend_order = -1; [S_ap, S_sigma_total, S_expected_total] = ... test_psd_energy_prepare_data(win_type, detrend_order, plist); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Evaluate the mean and standard deviation of the Checks on individual PSDs: % 2) Checks on individual PSDs: mean and standard deviation of the PSD points % 3) Evaluate the expected value, estimated from the std of the % individual segments % 4) Compares with the mean performed on the individual segments % 5) Checks on averaged PSDs: mean and standard deviation of all the PSD points % 6) Compares the grand-mean with the estimated white noise level % </AlgoDescription> % <AlgoCode> if stest atest = algo_test_psd_energy(ap, S_ap, S_expected_total, plist); else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_30 %% UTP_31 % <TestDescription> % % Tests "conservation of energy": % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) Calculate the expected level of noise from sample statistics % 3) Calculates the PSD of the individual parts, with: % BH92 window and mean detrend % 4) Evaluate the mean and standard deviation of the Checks on individual PSDs: % 5) Checks on individual PSDs: mean and standard deviation of the PSD points % 6) Evaluate the expected value, estimated from the std of the % individual segments % 7) Compares with the mean performed on the individual segments % 8) Checks on averaged PSDs: mean and standard deviation of all the PSD points % 9) Compares the grand-mean with the estimated white noise level % </TestDescription> function result = utp_31 % <SyntaxDescription> % % 1) Split the reference data into different segments % 2) Calculates the PSD of the individual parts % % </SyntaxDescription> % <SyntaxCode> try % Calculate the expected level of noise from sample statistics % Evaluate the PSDs on the parts, employing: % - BH92 window % - mean detrend win_type = 'BH92'; detrend_order = 0; [S_ap, S_sigma_total, S_expected_total] = ... test_psd_energy_prepare_data(win_type, detrend_order, plist); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Evaluate the mean and standard deviation of the Checks on individual PSDs: % 2) Checks on individual PSDs: mean and standard deviation of the PSD points % 3) Evaluate the expected value, estimated from the std of the % individual segments % 4) Compares with the mean performed on the individual segments % 5) Checks on averaged PSDs: mean and standard deviation of all the PSD points % 6) Compares the grand-mean with the estimated white noise level % </AlgoDescription> % <AlgoCode> if stest atest = algo_test_psd_energy(ap, S_ap, S_expected_total, plist); else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_31 %% UTP_32 % <TestDescription> % % Tests "conservation of energy": % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) Calculate the expected level of noise from sample statistics % 3) Calculates the PSD of the individual parts, with: % BH92 window and linear detrend % 4) Evaluate the mean and standard deviation of the Checks on individual PSDs: % 5) Checks on individual PSDs: mean and standard deviation of the PSD points % 6) Evaluate the expected value, estimated from the std of the % individual segments % 7) Compares with the mean performed on the individual segments % 8) Checks on averaged PSDs: mean and standard deviation of all the PSD points % 9) Compares the grand-mean with the estimated white noise level % </TestDescription> function result = utp_32 % <SyntaxDescription> % % 1) Split the reference data into different segments % 2) Calculates the PSD of the individual parts % % </SyntaxDescription> % <SyntaxCode> try % Calculate the expected level of noise from sample statistics % Evaluate the PSDs on the parts, employing: % - BH92 window % - linear detrend win_type = 'BH92'; detrend_order = 1; [S_ap, S_sigma_total, S_expected_total] = ... test_psd_energy_prepare_data(win_type, detrend_order, plist); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Evaluate the mean and standard deviation of the Checks on individual PSDs: % 2) Checks on individual PSDs: mean and standard deviation of the PSD points % 3) Evaluate the expected value, estimated from the std of the % individual segments % 4) Compares with the mean performed on the individual segments % 5) Checks on averaged PSDs: mean and standard deviation of all the PSD points % 6) Compares the grand-mean with the estimated white noise level % </AlgoDescription> % <AlgoCode> if stest atest = algo_test_psd_energy(ap, S_ap, S_expected_total, plist); else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_32 %% UTP_33 % <TestDescription> % % Tests "conservation of energy": % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) Calculate the expected level of noise from sample statistics % 3) Calculates the PSD of the individual parts, with: % Kaiser200 window and no detrend % 4) Evaluate the mean and standard deviation of the Checks on individual PSDs: % 5) Checks on individual PSDs: mean and standard deviation of the PSD points % 6) Evaluate the expected value, estimated from the std of the % individual segments % 7) Compares with the mean performed on the individual segments % 8) Checks on averaged PSDs: mean and standard deviation of all the PSD points % 9) Compares the grand-mean with the estimated white noise level % </TestDescription> function result = utp_33 % <SyntaxDescription> % % 1) Split the reference data into different segments % 2) Calculates the PSD of the individual parts % % </SyntaxDescription> % <SyntaxCode> try % Calculate the expected level of noise from sample statistics % Evaluate the PSDs on the parts, employing: % - Kaiser200 window % - no detrend win_type = 'Kaiser'; psll = 200; detrend_order = -1; [S_ap, S_sigma_total, S_expected_total] = ... test_psd_energy_prepare_data(win_type, detrend_order, plist('psll', psll)); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Evaluate the mean and standard deviation of the Checks on individual PSDs: % 2) Checks on individual PSDs: mean and standard deviation of the PSD points % 3) Evaluate the expected value, estimated from the std of the % individual segments % 4) Compares with the mean performed on the individual segments % 5) Checks on averaged PSDs: mean and standard deviation of all the PSD points % 6) Compares the grand-mean with the estimated white noise level % </AlgoDescription> % <AlgoCode> if stest atest = algo_test_psd_energy(ap, S_ap, S_expected_total, plist); else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_33 %% UTP_34 % <TestDescription> % % Tests "conservation of energy": % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) Calculate the expected level of noise from sample statistics % 3) Calculates the PSD of the individual parts, with: % Kaiser200 window and mean detrend % 4) Evaluate the mean and standard deviation of the Checks on individual PSDs: % 5) Checks on individual PSDs: mean and standard deviation of the PSD points % 6) Evaluate the expected value, estimated from the std of the % individual segments % 7) Compares with the mean performed on the individual segments % 8) Checks on averaged PSDs: mean and standard deviation of all the PSD points % 9) Compares the grand-mean with the estimated white noise level % </TestDescription> function result = utp_34 % <SyntaxDescription> % % 1) Split the reference data into different segments % 2) Calculates the PSD of the individual parts % % </SyntaxDescription> % <SyntaxCode> try % Calculate the expected level of noise from sample statistics % Evaluate the PSDs on the parts, employing: % - Kaiser200 window % - mean detrend win_type = 'Kaiser'; psll = 200; detrend_order = 0; [S_ap, S_sigma_total, S_expected_total] = ... test_psd_energy_prepare_data(win_type, detrend_order, plist('psll', psll)); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Evaluate the mean and standard deviation of the Checks on individual PSDs: % 2) Checks on individual PSDs: mean and standard deviation of the PSD points % 3) Evaluate the expected value, estimated from the std of the % individual segments % 4) Compares with the mean performed on the individual segments % 5) Checks on averaged PSDs: mean and standard deviation of all the PSD points % 6) Compares the grand-mean with the estimated white noise level % </AlgoDescription> % <AlgoCode> if stest atest = algo_test_psd_energy(ap, S_ap, S_expected_total, plist); else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_34 %% UTP_35 % <TestDescription> % % Tests "conservation of energy": % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) Calculate the expected level of noise from sample statistics % 3) Calculates the PSD of the individual parts, with: % Kaiser200 window and linear detrend % 4) Evaluate the mean and standard deviation of the Checks on individual PSDs: % 5) Checks on individual PSDs: mean and standard deviation of the PSD points % 6) Evaluate the expected value, estimated from the std of the % individual segments % 7) Compares with the mean performed on the individual segments % 8) Checks on averaged PSDs: mean and standard deviation of all the PSD points % 9) Compares the grand-mean with the estimated white noise level % </TestDescription> function result = utp_35 % <SyntaxDescription> % % 1) Split the reference data into different segments % 2) Calculates the PSD of the individual parts % % </SyntaxDescription> % <SyntaxCode> try % Calculate the expected level of noise from sample statistics % Evaluate the PSDs on the parts, employing: % - Kaiser200 window % - linear detrend win_type = 'Kaiser'; psll = 200; detrend_order = 1; [S_ap, S_sigma_total, S_expected_total] = ... test_psd_energy_prepare_data(win_type, detrend_order, plist('psll', psll)); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Evaluate the mean and standard deviation of the Checks on individual PSDs: % 2) Checks on individual PSDs: mean and standard deviation of the PSD points % 3) Evaluate the expected value, estimated from the std of the % individual segments % 4) Compares with the mean performed on the individual segments % 5) Checks on averaged PSDs: mean and standard deviation of all the PSD points % 6) Compares the grand-mean with the estimated white noise level % </AlgoDescription> % <AlgoCode> if stest atest = algo_test_psd_energy(ap, S_ap, S_expected_total, plist); else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_35 %% UTP_38 % <TestDescription> % % Tests "conservation of energy": % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) evaluate the sample mean value m and standard deviation s % 3) add a given trend of order n % 4) psd of the noise, with proper detrending % 5) compares the sqrt(sum(S * df)) with the standard deviation s % % </TestDescription> function result = utp_38 % <SyntaxDescription> % % 1) Prepare the test tsdata: % white noise from normal distribution + offset % 2) Calculate the statistical parameters % 3) Add a trend % 4) Estimate the psd % % </SyntaxDescription> % <SyntaxCode> try fs = 10; nsecs = 3600; sigma_distr = 1; trend_0 = 0; trend_1 = 2e-6; trend_2 = 5e-10; % White noise type = 'Normal'; a_n = ao(plist('waveform', 'noise', ... 'type', type, 'fs', fs, 'nsecs', nsecs, 'sigma', sigma_distr)); % Drift signal a_d = ao(plist('tsfcn', [num2str(trend_0) '*t.^0 + ' num2str(trend_1) ' *t + ' num2str(trend_2) ' *t.^2'], ... 'fs', fs, 'nsecs', nsecs)); % Total signal a = a_n + a_d; % Estimate the mean value and the sigma of the white noise time-series data sigma_sample = std(a_n.y, 1); % Evaluate the psd of the white noise time-series data % Parameters are crucial for the estimate to be correct! win = 'Rectangular'; olap = 0; detrend = 0; n_pts = -1; scale = 'PSD'; S_n = a_n.psd(plist('Win', win, 'olap', olap, ... 'Nfft', n_pts, 'order', detrend, 'scale', scale)); % Evaluate the psd of the "drifting" time-series data % Parameters are crucial for the estimate to be correct! win = 'Rectangular'; olap = 0; detrend = 2; n_pts = -1; scale = 'PSD'; S = a.psd(plist('Win', win, 'olap', olap, ... 'Nfft', n_pts, 'order', detrend, 'scale', scale)); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that calculated psd energy equals the rms content of the tsdata, estimated by the std of the sample % % </AlgoDescription> % <AlgoCode> atest = true; TOL = 1e-3; if stest % Integrals of the spectra sigma_psd = sqrt(sum(S.y * S.x(2))); sigma_psd_n = sqrt(sum(S_n.y * S_n.x(2))); if abs(sigma_psd - sigma_sample) / mean([sigma_psd sigma_sample]) >= TOL || ... abs(sigma_psd_n - sigma_sample) / mean([sigma_psd_n sigma_sample]) >= TOL atest = false; end else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_38 %% UTP_39 % <TestDescription> % % Tests "conservation of energy": % 1) white noise produced from uniform pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) evaluate the sample mean value m and standard deviation s % 3) add a given trend of order n % 4) psd of the noise, with proper detrending % 5) compares the sqrt(sum(S * df)) with the standard deviation s % % </TestDescription> function result = utp_39 % <SyntaxDescription> % % 1) Prepare the test tsdata: % white noise from uniform distribution + offset % 2) Calculate the statistical parameters % 3) Add a trend % 4) Estimate the psd % % </SyntaxDescription> % <SyntaxCode> try % Build time-series test data fs = 1; nsecs = 86400; sigma_distr = 1e-9; trend_0 = 3e-7; trend_1 = 1e-7; trend_2 = 2e-14; % White noise type = 'Uniform'; a_n = ao(plist('waveform', 'noise', ... 'type', type, 'fs', fs, 'nsecs', nsecs, 'sigma', sigma_distr)); % Drift signal a_d = ao(plist('tsfcn', [num2str(trend_0) '*t.^0 + ' num2str(trend_1) ' *t + ' num2str(trend_2) ' *t.^2'], ... 'fs', fs, 'nsecs', nsecs)); % Total signal a = a_n + a_d; % Estimate the mean value and the sigma of the white noise time-series data sigma_sample = std(a_n.y, 1); % Evaluate the psd of the white noise time-series data % Parameters are crucial for the estimate to be correct! win = 'Rectangular'; olap = 0; detrend = 0; n_pts = -1; scale = 'PSD'; S_n = a_n.psd(plist('Win', win, 'olap', olap, ... 'Nfft', n_pts, 'order', detrend, 'scale', scale)); % Evaluate the psd of the "drifting" time-series data % Parameters are crucial for the estimate to be correct! win = 'Rectangular'; olap = 0; detrend = 2; n_pts = -1; scale = 'PSD'; S = a.psd(plist('Win', win, 'olap', olap, ... 'Nfft', n_pts, 'order', detrend, 'scale', scale)); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that calculated psd energy equals the rms content of the tsdata, estimated by the std of the sample % % </AlgoDescription> % <AlgoCode> atest = true; TOL = 1e-3; if stest % Integrals of the spectra sigma_psd = sqrt(sum(S.y * S.x(2))); sigma_psd_n = sqrt(sum(S_n.y * S_n.x(2))); if abs(sigma_psd - sigma_sample) / mean([sigma_psd sigma_sample]) >= TOL || ... abs(sigma_psd_n - sigma_sample) / mean([sigma_psd_n sigma_sample]) >= TOL atest = false; end else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_39 %% UTP_40 % <TestDescription> % % Tests "conservation of energy": % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) evaluate the sample mean value m and standard deviation s % 3) add a given trend of order n % 4) psd of the noise, with proper detrending % 5) compares the sqrt(sum(S * df)) with the standard deviation s % % </TestDescription> function result = utp_40 % <SyntaxDescription> % % 1) Prepare the test tsdata: % white noise from normal distribution + offset % 2) Calculate the statistical parameters % 3) Add a trend % 4) Estimate the psd % % </SyntaxDescription> % <SyntaxCode> try % Build time-series test data fs = 100; nsecs = 100; sigma_distr = 1e-12; trend_0 = -5e-12; trend_1 = 4e-13; trend_2 = 1.17e-14; % White noise type = 'Normal'; a_n = ao(plist('waveform', 'noise', ... 'type', type, 'fs', fs, 'nsecs', nsecs, 'sigma', sigma_distr)); % Drift signal a_d = ao(plist('tsfcn', [num2str(trend_0) '*t.^0 + ' num2str(trend_1) ' *t + ' num2str(trend_2) ' *t.^2'], ... 'fs', fs, 'nsecs', nsecs)); % Total signal a = a_n + a_d; % Estimate the mean value and the sigma of the white noise time-series data sigma_sample = std(a_n.y, 1); % Evaluate the psd of the white noise time-series data % Parameters are crucial for the estimate to be correct! win = 'Rectangular'; olap = 0; detrend = 0; n_pts = -1; scale = 'PSD'; S_n = a_n.psd(plist('Win', win, 'olap', olap, ... 'Nfft', n_pts, 'order', detrend, 'scale', scale)); % Evaluate the psd of the "drifting" time-series data % Parameters are crucial for the estimate to be correct! win = 'Rectangular'; olap = 0; detrend = 2; n_pts = -1; scale = 'PSD'; S = a.psd(plist('Win', win, 'olap', olap, ... 'Nfft', n_pts, 'order', detrend, 'scale', scale)); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that calculated psd energy equals the rms content of the tsdata, estimated by the std of the sample % % </AlgoDescription> % <AlgoCode> atest = true; TOL = 1e-3; if stest % Integrals of the spectra sigma_psd = sqrt(sum(S.y * S.x(2))); sigma_psd_n = sqrt(sum(S_n.y * S_n.x(2))); if abs(sigma_psd - sigma_sample) / mean([sigma_psd sigma_sample]) >= TOL || ... abs(sigma_psd_n - sigma_sample) / mean([sigma_psd_n sigma_sample]) >= TOL atest = false; end else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_40 %% UTP_41 % <TestDescription> % % Tests the effect of windowing: % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) psd of the noise, without detrending, Rectangular window % 3) Apply the detrending and the window manually % 4) psd of the noise, without detrending, Rectangular window % 5) compares the to psds % % </TestDescription> function result = utp_41 % <SyntaxDescription> % % 1) Prepare the test tsdata: % white noise from normal distribution + offset % 2) Calculate the statistical parameters % 3) Estimate the psd without detrending, Rectangular window % 4) Manually apply window to the data % 5) Estimate the psd without detrending, Rectangular window % % </SyntaxDescription> % <SyntaxCode> try noise_type = 'Normal'; win_type = 'Rectangular'; scale = 'PSD'; [a, S_name, S_obj, S_man] = prepare_analyze_noise_win(win_type, noise_type, scale, ... plist('win_test', true)); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that calculated psds are identical within a part in 10^12 % % </AlgoDescription> % <AlgoCode> TOL = 1e-12; if stest atest = algo_test_psd_win(S_name, S_obj, S_man, plist('TOL', TOL)); else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_41 %% UTP_42 % <TestDescription> % % Tests the effect of windowing: % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) psd of the noise, without detrending, BH92 window % 3) Apply the detrending and the window manually % 4) psd of the noise, without detrending, Rectangular window % 5) compares the to psds % % </TestDescription> function result = utp_42 % <SyntaxDescription> % % 1) Prepare the test tsdata: % white noise from normal distribution + offset % 2) Calculate the statistical parameters % 3) Estimate the psd without detrending, BH92 window % 4) Manually apply window to the data % 5) Estimate the psd without detrending, Rectangular window % % </SyntaxDescription> % <SyntaxCode> try noise_type = 'Normal'; win_type = 'BH92'; scale = 'PSD'; [a, S_name, S_obj, S_man] = prepare_analyze_noise_win(win_type, noise_type, scale, ... plist('win_test', true)); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that calculated psds are identical within a part in 10^12 % % </AlgoDescription> % <AlgoCode> TOL = 1e-12; if stest atest = algo_test_psd_win(S_name, S_obj, S_man, plist('TOL', TOL)); else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_42 %% UTP_43 % <TestDescription> % % Tests the effect of windowing: % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) psd of the noise, without detrending, Hamming window % 3) Apply the detrending and the window manually % 4) psd of the noise, without detrending, Rectangular window % 5) compares the to psds % % </TestDescription> function result = utp_43 % <SyntaxDescription> % % 1) Prepare the test tsdata: % white noise from normal distribution + offset % 2) Calculate the statistical parameters % 3) Estimate the psd without detrending, Hamming window % 4) Manually apply window to the data % 5) Estimate the psd without detrending, Rectangular window % % </SyntaxDescription> % <SyntaxCode> try noise_type = 'Normal'; win_type = 'Hamming'; scale = 'PSD'; [a, S_name, S_obj, S_man] = prepare_analyze_noise_win(win_type, noise_type, scale, ... plist('win_test', true)); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that calculated psds are identical within a part in 10^12 % % </AlgoDescription> % <AlgoCode> TOL = 1e-12; if stest atest = algo_test_psd_win(S_name, S_obj, S_man, plist('TOL', TOL)); else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_43 %% UTP_44 % <TestDescription> % % Tests the effect of windowing: % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) psd of the noise, without detrending, Hanning window % 3) Apply the detrending and the window manually % 4) psd of the noise, without detrending, Rectangular window % 5) compares the to psds % % </TestDescription> function result = utp_44 % <SyntaxDescription> % % 1) Prepare the test tsdata: % white noise from normal distribution + offset % 2) Calculate the statistical parameters % 3) Estimate the psd without detrending, Hanning window % 4) Manually apply window to the data % 5) Estimate the psd without detrending, Rectangular window % % </SyntaxDescription> % <SyntaxCode> msg = ''; try noise_type = 'Normal'; win_type = 'Hanning'; scale = 'PSD'; [a, S_name, S_obj, S_man] = prepare_analyze_noise_win(win_type, noise_type, scale, ... plist('win_test', true)); stest = true; catch err disp(err.message) msg = err.message; stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that calculated psds are identical within a part in 10^12 % % </AlgoDescription> % <AlgoCode> TOL = 1e-12; if stest [atest, msg] = algo_test_psd_win(S_name, S_obj, S_man, plist('TOL', TOL)); else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename, msg); end % END UTP_44 %% UTP_45 % <TestDescription> % % Tests the effect of windowing: % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) psd of the noise, without detrending, Bartlett window % 3) Apply the detrending and the window manually % 4) psd of the noise, without detrending, Rectangular window % 5) compares the to psds % % </TestDescription> function result = utp_45 % <SyntaxDescription> % % 1) Prepare the test tsdata: % white noise from normal distribution + offset % 2) Calculate the statistical parameters % 3) Estimate the psd without detrending, Bartlett window % 4) Manually apply window to the data % 5) Estimate the psd without detrending, Rectangular window % % </SyntaxDescription> % <SyntaxCode> try noise_type = 'Normal'; win_type = 'Bartlett'; scale = 'PSD'; [a, S_name, S_obj, S_man] = prepare_analyze_noise_win(win_type, noise_type, scale, ... plist('win_test', true)); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that calculated psds are identical within a part in 10^12 % % </AlgoDescription> % <AlgoCode> TOL = 1e-12; if stest atest = algo_test_psd_win(S_name, S_obj, S_man, plist('TOL', TOL)); else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_45 %% UTP_46 % <TestDescription> % % Tests the effect of windowing: % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) psd of the noise, without detrending, Nuttall3 window % 3) Apply the detrending and the window manually % 4) psd of the noise, without detrending, Rectangular window % 5) compares the to psds % % </TestDescription> function result = utp_46 % <SyntaxDescription> % % 1) Prepare the test tsdata: % white noise from normal distribution + offset % 2) Calculate the statistical parameters % 3) Estimate the psd without detrending, Nuttall3 window % 4) Manually apply window to the data % 5) Estimate the psd without detrending, Rectangular window % % </SyntaxDescription> % <SyntaxCode> msg = ''; try noise_type = 'Normal'; win_type = 'Nuttall3'; scale = 'PSD'; [a, S_name, S_obj, S_man] = prepare_analyze_noise_win(win_type, noise_type, scale, ... plist('win_test', true)); stest = true; catch err disp(err.message) msg = err.message; stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that calculated psds are identical within a part in 10^12 % % </AlgoDescription> % <AlgoCode> TOL = 1e-12; if stest [atest, msg] = algo_test_psd_win(S_name, S_obj, S_man, plist('TOL', TOL)); else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename, msg); end % END UTP_46 %% UTP_47 % <TestDescription> % % Tests the effect of windowing: % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) psd of the noise, without detrending, Kaiser psll = random window % 3) Apply the detrending and the window manually % 4) psd of the noise, without detrending, Rectangular window % 5) compares the to psds % % </TestDescription> function result = utp_47 % <SyntaxDescription> % % 1) Prepare the test tsdata: % white noise from normal distribution + offset % 2) Calculate the statistical parameters % 3) Estimate the psd without detrending, Kaiser psll = random window % 4) Manually apply window to the data % 5) Estimate the psd without detrending, Rectangular window % % </SyntaxDescription> % <SyntaxCode> try noise_type = 'Normal'; win_type = 'Kaiser'; psll = utils.math.randelement([10:10:200],1); scale = 'PSD'; [a, S_name, S_obj, S_man] = prepare_analyze_noise_win(win_type, noise_type, scale, ... plist('win_test', true)); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that calculated psds are identical within a part in 10^12 % % </AlgoDescription> % <AlgoCode> TOL = 1e-12; if stest atest = algo_test_psd_win(S_name, S_obj, S_man, plist('TOL', TOL)); else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_47 %% UTP_48 % <TestDescription> % % Tests the effect of windowing: % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) psd of the noise, without detrending, Kaiser psll = default window % 3) Apply the detrending and the window manually % 4) psd of the noise, without detrending, Rectangular window % 5) compares the to psds % % </TestDescription> function result = utp_48 % <SyntaxDescription> % % 1) Prepare the test tsdata: % white noise from normal distribution + offset % 2) Calculate the statistical parameters % 3) Estimate the psd without detrending, Kaiser psll = default window % 4) Manually apply window to the data % 5) Estimate the psd without detrending, Rectangular window % % </SyntaxDescription> % <SyntaxCode> try noise_type = 'Normal'; win_type = 'Kaiser'; scale = 'PSD'; [a, S_name, S_obj, S_man] = prepare_analyze_noise_win(win_type, noise_type, scale, ... plist('win_test', true)); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that calculated psds are identical within a part in 10^12 % % </AlgoDescription> % <AlgoCode> TOL = 1e-12; if stest atest = algo_test_psd_win(S_name, S_obj, S_man, plist('TOL', TOL)); else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_48 %% UTP_49 % <TestDescription> % % Tests the effect of windowing: % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) psd of the noise, without detrending, Nuttall4 window % 3) Apply the detrending and the window manually % 4) psd of the noise, without detrending, Rectangular window % 5) compares the to psds % % </TestDescription> function result = utp_49 % <SyntaxDescription> % % 1) Prepare the test tsdata: % white noise from normal distribution + offset % 2) Calculate the statistical parameters % 3) Estimate the psd without detrending, Nuttall4 window % 4) Manually apply window to the data % 5) Estimate the psd without detrending, Rectangular window % % </SyntaxDescription> % <SyntaxCode> try noise_type = 'Normal'; win_type = 'Nuttall4'; scale = 'PSD'; [a, S_name, S_obj, S_man] = prepare_analyze_noise_win(win_type, noise_type, scale, ... plist('win_test', true)); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that calculated psds are identical within a part in 10^12 % % </AlgoDescription> % <AlgoCode> TOL = 1e-12; if stest atest = algo_test_psd_win(S_name, S_obj, S_man, plist('TOL', TOL)); else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_49 %% UTP_50 % <TestDescription> % % Tests the effect of windowing: % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) psd of the noise, without detrending, SFT3F window % 3) Apply the detrending and the window manually % 4) psd of the noise, without detrending, Rectangular window % 5) compares the to psds % % </TestDescription> function result = utp_50 % <SyntaxDescription> % % 1) Prepare the test tsdata: % white noise from normal distribution + offset % 2) Calculate the statistical parameters % 3) Estimate the psd without detrending, SFT3F window % 4) Manually apply window to the data % 5) Estimate the psd without detrending, Rectangular window % % </SyntaxDescription> % <SyntaxCode> try noise_type = 'Normal'; win_type = 'SFT3F'; scale = 'PSD'; [a, S_name, S_obj, S_man] = prepare_analyze_noise_win(win_type, noise_type, scale, ... plist('win_test', true)); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that calculated psds are identical within a part in 10^12 % % </AlgoDescription> % <AlgoCode> TOL = 1e-12; if stest atest = algo_test_psd_win(S_name, S_obj, S_man, plist('TOL', TOL)); else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_50 %% UTP_51 % <TestDescription> % % Tests the possibility to set the number of averages rather than setting the Nfft: % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) psd of the noise, without detrending, random window, set number of % averages % 3) check the effective number of averages % % </TestDescription> function result = utp_51 % <SyntaxDescription> % % 1) Prepare the test tsdata: % white noise from normal distribution + offset % 2) psd of the noise, without detrending, random window, set number of % averages % % </SyntaxDescription> % <SyntaxCode> try noise_type = 'Normal'; % Evaluate the psd of the white noise time-series data [a, S1, S2, navs] = prepare_analyze_noise_navs(noise_type, plist); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that calculated navs are identical to those requested % % </AlgoDescription> % <AlgoCode> atest = true; if stest % Compare the navs written in the output object with the requested one if ne(navs, S1.data.navs) if ne(find(S1.hist.plistUsed, 'navs'), S1.data.navs) atest = false; end end else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_51 %% UTP_52 % <TestDescription> % % Tests the possibility to set the number of averages rather than setting the Nfft: % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) psd of the noise, without detrending, random window, random navs % 3) get the number of averages % 4) get the nfft used % 5) run psd again, with the nfft used % 6) compare the calculated objects % % </TestDescription> function result = utp_52 % <SyntaxDescription> % % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) psd of the noise, without detrending, random window, random navs % 3) get the number of averages % 4) get the nfft used % 5) run psd again, with the nfft used % % </SyntaxDescription> % <SyntaxCode> try noise_type = 'Normal'; % Evaluate the psd of the white noise time-series data [a, S1, S2, navs] = prepare_analyze_noise_navs(noise_type, plist); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that calculated objects S1 and S2 are identical % % </AlgoDescription> % <AlgoCode> atest = true; if stest % Compare the output objects if ne(S1, S2, ple3) atest = false; end else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_52 %% UTP_53 % <TestDescription> % % Tests the possibility to set the number of averages rather than setting the Nfft: % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) psd of the noise, without detrending, random window, random navs % 3) get the number of averages % 4) get the nfft used % 5) run psd again, with the nfft used % 6) compare navs, nfft, psds % % </TestDescription> function result = utp_53 % <SyntaxDescription> % % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) psd of the noise, without detrending, random window, random navs % 3) get the number of averages % 4) get the nfft used % 5) run psd again, with the nfft used % 6) run psd again, with conflicting parameters, and verify it uses % nfft rather than navs % % </SyntaxDescription> % <SyntaxCode> try noise_type = 'Uniform'; % Evaluate the psd of the white noise time-series data [a, S1, S2, navs] = prepare_analyze_noise_navs(noise_type, plist); npts_3 = fix(find(S1.hist.plistUsed, 'Nfft')/2); % Calculates the psd asking for the number of points AND the window length pl_spec = S1.hist.plistUsed; pl_spec.pset('Nfft', npts_3, 'navs', navs); S3 = a.psd(pl_spec); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that calculated objects S1 and S2 are identical % 2) Check that S3 used different values % % </AlgoDescription> % <AlgoCode> atest = true; if stest % Compare the navs written in the output object with the requested one if ne(S1, S2, ple3) || ... ne(find(S3.hist.plistUsed, 'Nfft'), npts_3) || eq(S3.data.navs, navs) atest = false; end else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_53 %% UTP_54 % <TestDescription> % % Tests "conservation of energy": % 1) white noise produced from normal pdf, with a given mean value and % sigma (distribution's 1st and 2nd order) % 2) evaluate the sample mean value m and standard deviation s % 3) psd of the white noise % 4) compares the sqrt(sum(S * df)) with the standard deviation s % % </TestDescription> function result = utp_54 % <SyntaxDescription> % % 1) Prepare the test tsdata: % white noise from normal distribution + offset % 2) Calculate the statistical parameters % 3) Estimate the psd % % </SyntaxDescription> % <SyntaxCode> try % Array of parameters to pick from fs_list = [0.1;1;10]; nsecs_list = [100:100:10000]'; sigma_distr_list = [1e-6 2e-3 0.25 1:0.1:10]'; mu_distr_list = [1e-6 2e-3 0.25 1:0.1:10]'; % Build time-series test data % Picks the values at random from the list fs = utils.math.randelement(fs_list, 1); nsecs = utils.math.randelement(nsecs_list, 1); sigma_distr = utils.math.randelement(sigma_distr_list, 1); mu_distr = utils.math.randelement(mu_distr_list, 1); % White noise type = 'Normal'; a_n = ao(plist('waveform', 'noise', ... 'type', type, 'fs', fs, 'nsecs', nsecs, 'sigma', sigma_distr)); a_const = ao(mu_distr); a = a_n + a_const; % Estimate the mean value and the sigma of the time-series data sigma_sample = std(a.y, 1); mu_sample = mean(a.y); % Evaluate the psd of the time-series data % Parameters are crucial for the estimate to be correct! win = 'Rectangular'; olap = 0; detrend = 0; n_pts = -1; scale = 'PSD'; S = a.psd(plist('Win', win, 'olap', olap, ... 'Nfft', n_pts, 'order', detrend, 'scale', scale)); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Check that calculated psd energy equals the rms content of the tsdata, estimated by the std of the sample % % </AlgoDescription> % <AlgoCode> atest = true; TOL = 1e-12; if stest % Integral of the spectrum sigma_psd = sqrt(sum(S.y * S.x(2))); if abs(sigma_psd - sigma_sample) / mean([sigma_psd sigma_sample]) >= TOL atest = false; end else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_54 %% UTP_60 % <TestDescription> % % Tests the 'basic' call: % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) psd of the noise, defualt detrending, default window, no averaging % % </TestDescription> function result = utp_60 % <SyntaxDescription> % % 1) Prepare the test tsdata: % white noise from normal distribution + offset % 2) Calculate the statistical parameters % 3) Estimate the psd with default detrending, default window, no averaging % % </SyntaxDescription> % <SyntaxCode> try noise_type = 'Normal'; win_type = char(prefs.getMiscPrefs.getDefaultWindow); scale = 'PSD'; [a, S_name, S_obj, S_man] = prepare_analyze_noise_win(win_type, noise_type, scale, ... plist('win_test', false, 'reduce_pts', false)); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Nothing to check % % </AlgoDescription> % <AlgoCode> if stest atest = true; else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_60 %% UTP_61 % <TestDescription> % % Tests the 'basic' call: % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) psd of the noise, defualt detrending, default window, no averaging % % </TestDescription> function result = utp_61 % <SyntaxDescription> % % 1) Prepare the test tsdata: % white noise from normal distribution + offset % 2) Calculate the statistical parameters % 3) Estimate the psd with default detrending, default window, no averaging % % </SyntaxDescription> % <SyntaxCode> try noise_type = 'Normal'; win_type = char(prefs.getMiscPrefs.getDefaultWindow); scale = 'ASD'; [a, S_name, S_obj, S_man] = prepare_analyze_noise_win(win_type, noise_type, scale, ... plist('win_test', false, 'reduce_pts', false)); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Nothing to check % % </AlgoDescription> % <AlgoCode> if stest atest = true; else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_61 %% UTP_62 % <TestDescription> % % Tests the 'basic' call: % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) psd of the noise, defualt detrending, default window, no averaging % % </TestDescription> function result = utp_62 % <SyntaxDescription> % % 1) Prepare the test tsdata: % white noise from normal distribution + offset % 2) Calculate the statistical parameters % 3) Estimate the psd with default detrending, default window, no averaging % % </SyntaxDescription> % <SyntaxCode> try noise_type = 'Normal'; win_type = char(prefs.getMiscPrefs.getDefaultWindow); scale = 'PS'; [a, S_name, S_obj, S_man] = prepare_analyze_noise_win(win_type, noise_type, scale, ... plist('win_test', false, 'reduce_pts', false)); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Nothing to check % % </AlgoDescription> % <AlgoCode> if stest atest = true; else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_62 %% UTP_63 % <TestDescription> % % Tests the 'basic' call: % 1) white noise produced from normal pdf, with: % a given mean value and sigma (distribution's 1st and 2nd order) % 2) psd of the noise, defualt detrending, default window, no averaging % % </TestDescription> function result = utp_63 % <SyntaxDescription> % % 1) Prepare the test tsdata: % white noise from normal distribution + offset % 2) Calculate the statistical parameters % 3) Estimate the psd with default detrending, default window, no averaging % % </SyntaxDescription> % <SyntaxCode> try noise_type = 'Normal'; win_type = char(prefs.getMiscPrefs.getDefaultWindow); scale = 'AS'; [a, S_name, S_obj, S_man] = prepare_analyze_noise_win(win_type, noise_type, scale, ... plist('win_test', false, 'reduce_pts', false)); stest = true; catch err disp(err.message) stest = false; end % </SyntaxCode> % <AlgoDescription> % % 1) Nothing to check % % </AlgoDescription> % <AlgoCode> if stest atest = true; else atest = false; end % </AlgoCode> % Return a result structure result = utp_prepare_result(atest, stest, dbstack, mfilename); end % END UTP_63 %% Helper function for window call construction function [a, S_name, S_obj, S_w] = prepare_analyze_noise_win(win_type, noise_type, scale, pli) % Array of parameters to pick from fs_list = [0.1;1;2;5;10]; nsecs_list = [20 100 1000:1000:10000]'; sigma_distr_list = [1e-6 2e-3 0.25 1:0.1:10]'; trend_0_list = [1e-6 2e-3 0.25 1:0.1:10]'; % Build time-series test data % Picks the values at random from the list fs = utils.math.randelement(fs_list, 1); nsecs = utils.math.randelement(nsecs_list, 1); sigma_distr = utils.math.randelement(sigma_distr_list, 1); trend_0 = utils.math.randelement(trend_0_list, 1); % Pick units and prefix from those supported unit_list = unit.supportedUnits; % remove the first empty unit '' from the list, because then is it % possible that we add a prefix to an empty unit unit_list = unit_list(2:end); prefix_list = unit.supportedPrefixes; % White noise a_n = ao(plist('waveform', 'noise', ... 'type', noise_type, 'fs', fs, 'nsecs', nsecs, 'sigma', sigma_distr)); % Constant signal a_c = ao(trend_0); % Total signal a = a_n + a_c; % Set units a.setYunits(unit([cell2mat(utils.math.randelement(prefix_list,1)) cell2mat(utils.math.randelement(unit_list,1))])); % Evaluate the psd of the white noise time-series data n_pts_reduction_factor = find(pli, 'reduce_pts'); if isempty(n_pts_reduction_factor) || ~n_pts_reduction_factor n_pts_reduction_factor = 1; % Taking single windows, no reduction in number of points end n_pts = a.nsecs * a.fs/ n_pts_reduction_factor; olap = find(pli, 'olap'); if isempty(olap) olap = 0; %% Should we take the rov instead? end if n_pts_reduction_factor == 1 olap = 0; % This does not matter in the case of single window end detrend_order = 0; switch lower(win_type) case 'kaiser' psll = find(pli, 'psll'); if isempty(psll) psll = find(ao.getInfo('psd').plists, 'psll'); end win = specwin(win_type, fs*nsecs, psll); pl_spec = plist('Win', win_type, 'psll', psll, 'olap', olap, ... 'Nfft', n_pts, 'order', detrend_order, 'scale', scale); otherwise win = specwin(win_type, fs*nsecs); pl_spec = plist('Win', win_type, 'olap', olap, ... 'Nfft', n_pts, 'order', detrend_order, 'scale', scale); end if find(pli, 'win_test') % Makes a copy of the data, and apply the window manually % Apply the detrending a_w = a.detrend(plist('N', detrend_order)); % Apply the window, taking into account the correct normalization a_w.setY(a_w.y .* win.win'./sqrt(win.ws2 / win.len)); % Evaluate the psd of the windowed data S_w = a_w.psd(plist('Win', 'Rectangular', 'olap', olap, ... 'Nfft', n_pts, 'order', -1, 'scale', scale)); % Calls psd applying the detrend and window internally % (passig window object) S_obj = a.psd(pl_spec); else S_obj = ao; S_w = ao; end % Calls psd applying the detrend and window internally % (passig window name) pl_spec.pset('Win', win_type); S_name = a.psd(pl_spec); end %% Helper function for window call construction function [sigma_sample, sigma_psd] = prepare_analyze_noise_energy(win_type, noise_type, scale, pli) % Build time-series test data if find(pli, 'fs') fs = find(pli, 'fs'); nsecs = find(pli, 'nsecs'); sigma_distr = find(pli, 'sigma_distr'); mu_distr = find(pli, 'mu_distr'); else % Array of parameters to pick from fs_list = [0.1;1;10]; nsecs_list = [100:100:10000]; sigma_distr_list = [1e-6 2e-3 0.25 1:0.1:10]'; mu_distr_list = [1e-6 2e-3 0.25 1:0.1:10]'; % Picks the values at random from the list fs = utils.math.randelement(fs_list, 1); nsecs = utils.math.randelement(nsecs_list, 1); sigma_distr = utils.math.randelement(sigma_distr_list, 1); mu_distr = utils.math.randelement(mu_distr_list, 1); end % Pick units and prefix from those supported unit_list = unit.supportedUnits; % remove the first empty unit '' from the list, because then is it % possible that we add a prefix to an empty unit unit_list = unit_list(2:end); prefix_list = unit.supportedPrefixes; % White noise a_n = ao(plist('waveform', 'noise', ... 'type', noise_type, 'fs', fs, 'nsecs', nsecs, 'sigma', sigma_distr)); % Constant signal a_c = ao(mu_distr); % Total signal a = a_n + a_c; % Set units a.setYunits(unit([cell2mat(utils.math.randelement(prefix_list,1)) cell2mat(utils.math.randelement(unit_list,1))])); % Evaluate the psd of the white noise time-series data olap = find(pli, 'olap'); if isempty(olap) olap = 0; %% Should we take the rov instead? end detrend_order = find(pli, 'detrend_order'); if isempty(detrend_order) detrend_order = 0; end n_pts = -1; switch lower(win_type) case 'kaiser' psll = find(pli, 'psll'); if isempty(psll) psll = find(ao.getInfo('psd').plists, 'psll'); end win = specwin(win_type, fs*nsecs, psll); pl_spec = plist('Win', win_type, 'psll', psll, 'olap', olap, ... 'Nfft', n_pts, 'order', detrend_order, 'scale', scale); otherwise win = specwin(win_type, fs*nsecs); pl_spec = plist('Win', win_type, 'olap', olap, ... 'Nfft', n_pts, 'order', detrend_order, 'scale', scale); end % Calls psd applying the detrend and window internally % (passig window name) pl_spec.pset('Win', win_type); S = a.psd(pl_spec); % Evaluate the statistics of the time-series and spectra % Estimate the mean value and the sigma of the time-series data sigma_sample = std(a.y, 1); mu_sample = mean(a.y); % Integral of the spectrum sigma_psd = sqrt(sum(S.y * S.x(2))); end %% Helper function for window call construction, navs option function [a, S, S1, navs] = prepare_analyze_noise_navs(noise_type, pli) % Array of parameters to pick from fs_list = [0.1;1;2;5;10]; nsecs_list = [2000:1000:10000]'; sigma_distr_list = [1e-6 2e-3 0.25 1:0.1:10]'; trend_0_list = [1e-6 2e-3 0.25 1:0.1:10]'; % Build time-series test data % Picks the values at random from the list fs = utils.math.randelement(fs_list, 1); nsecs = utils.math.randelement(nsecs_list, 1); sigma_distr = utils.math.randelement(sigma_distr_list, 1); trend_0 = utils.math.randelement(trend_0_list, 1); % Pick units and prefix from those supported unit_list = unit.supportedUnits; % remove the first empty unit '' from the list, because then is it % possible that we add a prefix to an empty unit unit_list = unit_list(2:end); prefix_list = unit.supportedPrefixes; % White noise a_n = ao(plist('waveform', 'noise', ... 'type', noise_type, 'fs', fs, 'nsecs', nsecs, 'sigma', sigma_distr)); % Constant signal a_c = ao(trend_0); % Total signals a = a_n + a_c; % Set units a.setYunits(unit([cell2mat(utils.math.randelement(prefix_list,1)) cell2mat(utils.math.randelement(unit_list,1))])); % Evaluate the psd of the white noise time-series data olap = 0; detrend_order = 0; n_pts = -1; navs = fix(utils.math.randelement(logspace(0,log10(max(0,a.len/10)),50),1)); % Evaluate the psd of the white noise time-series data % Window win_list = specwin.getTypes; win_type = utils.math.randelement(win_list,1); win_type = win_type{1}; switch lower(win_type) case 'kaiser' psll = utils.math.randelement([0:10:200],1); if psll == 0 psll = find(ao.getInfo('psd').plists, 'psll'); end pl_spec = plist('Win', win_type, 'psll', psll, 'olap', olap, 'order', detrend_order); case 'levelledhanning' pl_spec = plist('Win', 'BH92', 'olap', olap, 'order', detrend_order); otherwise pl_spec = plist('Win', win_type, 'olap', olap, 'order', detrend_order); end % Calls psd asking for the number of averages pl_spec.pset('Nfft', n_pts, 'navs', navs); S = a.psd(pl_spec); % Calls psd asking for the number of points just evaluated pl_spec.pset('Nfft', find(S.hist.plistUsed, 'Nfft')); pl_spec.remove('navs'); S1 = a.psd(pl_spec); end %% Helper function for window call construction, navs option function [S_ap, S_sigma_total, S_expected_total] = test_psd_energy_prepare_data(win_type, detrend_order, pli) % Calculate the expected level of noise from sample statistics fs = a_n_long.fs; S_sigma_total = std(a_n_long.y, 1); S_expected_total = S_sigma_total^2 / fs * 2; olap = 0; n_pts = -1; scale = 'PSD'; switch lower(win_type) case 'kaiser' psll = find(pli, 'psll'); if isempty(psll) psll = find(ao.getInfo('psd').plists, 'psll'); end pl_spec = plist('Win', win_type, 'psll', psll, 'olap', olap, ... 'Nfft', n_pts, 'order', detrend_order, 'scale', scale); case 'levelledhanning' pl_spec = plist('Win', 'BH92', 'olap', olap, ... 'Nfft', n_pts, 'order', detrend_order, 'scale', scale); otherwise pl_spec = plist('Win', win_type, 'olap', olap, ... 'Nfft', n_pts, 'order', detrend_order, 'scale', scale); end % Calls psd applying the detrend and window internally % (passig window name) S_ap = ap.psd(pl_spec); end %% Algorithm test for window comparison function varargout = algo_test_psd_win(S_name, S_obj, S_man, pli) atest = true; msg = ''; pl_def = plist('TOL', 1e-12); TOL = find(parse(pli, pl_def), 'TOL'); % Compare the psd evaluated with the window name / object if ~eq(S_name, S_obj, ple1) atest = false; msg = ['Two PSD objects are not equal: ' lastwarn]; end % Compare the psd evaluated with the two methods res = abs(S_obj - S_man) ./ S_obj; if any(res.y >= TOL) atest = false; msg = 'Two PSDs are not equal'; end if nargout == 1 varargout{1} = atest; elseif nargout == 2 varargout{1} = atest; varargout{2} = msg; else error('Incorrect outputs'); end end %% Algorithm test for energy conservation tests function atest = algo_test_psd_energy(ap, S_ap, S_expected_total, pli) % Checks on individual PSDs: mean and standard deviation of the % PSD points fs = ap.fs; Npoints = length(S_ap(1).y); S_mean = mean(S_ap.y); S_err = std(S_ap.y, 1) / sqrt(Npoints); % Expected value, estimated from the std of the individual segments S_sigma = std(ap.y, 1); S_expected = S_sigma.^2 / fs * 2; % Comparison with the mean performed on the individual segments test_level = abs(S_mean - S_expected) < S_err; sum(test_level) / Nsegments; % TODO: implement hypothesis test here % Checks on averaged PSDs: mean and standard deviation of all the % PSD points S_mean_total = mean(S_mean); S_err_total = std(S_mean) / sqrt(Nsegments); % Compares the grand-mean with the estimated white noise level if sum(abs(S_mean_total - S_expected_total) < S_err_total) atest = true; else atest = false; end end end