Mercurial > hg > ltpda
diff testing/utp_1.1/utps/ao/utp_ao_psd.m @ 44:409a22968d5e default
Add unit tests
| author | Daniele Nicolodi <nicolodi@science.unitn.it> |
|---|---|
| date | Tue, 06 Dec 2011 18:42:11 +0100 |
| parents | |
| children |
line wrap: on
line diff
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/testing/utp_1.1/utps/ao/utp_ao_psd.m Tue Dec 06 18:42:11 2011 +0100 @@ -0,0 +1,4476 @@ +% 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