0001 function varargout = ltpda_linedetect(varargin)
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021
0022
0023
0024
0025
0026 if nargin == 1
0027 in = char(varargin{1});
0028 if strcmp(in, 'Params')
0029 varargout{1} = getDefaultPL();
0030 return
0031 end
0032 end
0033
0034
0035 invars = {};
0036 for j=1:nargin
0037 if isa(varargin{j}, 'ao')
0038 invars = [invars cellstr(inputname(j))];
0039 end
0040 end
0041
0042 ALGONAME = mfilename;
0043 VERSION = '$Id: ltpda_linedetect.m,v 1.5 2007/06/12 12:32:14 hewitson Exp $';
0044
0045
0046 as = [];
0047 ps = [];
0048 queries = [];
0049 for j=1:nargin
0050 if isa(varargin{j}, 'ao')
0051 as = [as varargin{j}];
0052 end
0053 if isa(varargin{j}, 'plist')
0054 ps = [ps varargin{j}];
0055 end
0056 end
0057
0058 Na = length(as);
0059 if isempty(as)
0060 error('### Please input at least one AO.');
0061 end
0062
0063
0064 if ~isempty(ps)
0065 pl = combine(ps, getDefaultPL());
0066 else
0067 pl = getDefaultPL();
0068 end
0069
0070
0071 N = find(pl, 'N');
0072 fsearch = find(pl, 'fsearch');
0073 thresh = find(pl, 'thresh');
0074 bw = find(pl, 'bw');
0075 hc = find(pl, 'hc');
0076
0077
0078 bo = [];
0079 for j=1:Na
0080
0081 a = as(j);
0082 d = a.data;
0083
0084 if isa(d, 'fsdata')
0085
0086
0087 nf = ltpda_nfest(a, pl);
0088
0089
0090 r = a./nf;
0091
0092
0093 lines = findLines(d.xx, get(r.data, 'f'), get(r.data, 'xx'), thresh, N, fsearch);
0094
0095
0096 fs = fsdata([lines(:).f], [lines(:).a], d.fs);
0097
0098 else
0099 error('### I can only find lines in frequency-series AOs.');
0100 end
0101
0102
0103
0104
0105 h = history(ALGONAME, VERSION, pl, a.hist);
0106 h = set(h, 'invars', invars);
0107
0108
0109 b = ao(fs, h);
0110
0111
0112
0113 if isempty(invars{j})
0114 n1 = a.name;
0115 else
0116 n1 = invars{j};
0117 end
0118 b = set(b, 'name', sprintf('lines(%s)', n1));
0119
0120
0121 bo = [bo b];
0122
0123 end
0124
0125
0126 varargout{1} = bo;
0127
0128
0129
0130 function lines = findLines(axx, f, naxx, thresh, N, fsearch)
0131
0132
0133 l = 0;
0134 pmax = 0;
0135 pmaxi = 0;
0136 line = [];
0137 idx = find( f>=fsearch(1) & f<=fsearch(2) );
0138 for j=1:length(idx)
0139 v = naxx(idx(j));
0140 if v > thresh
0141 if v > pmax
0142 pmax = v;
0143 pmaxi = idx(j);
0144 end
0145 else
0146 if pmax > 0
0147
0148 fidx = pmaxi;
0149 l = l+1;
0150 line(l).idx = fidx;
0151 line(l).f = f(fidx);
0152 line(l).a = axx(fidx);
0153 end
0154 pmax = 0;
0155 end
0156 end
0157
0158
0159 lines = [];
0160 if ~isempty(line)
0161 [bl, lidx] = sort([line.a], 'descend');
0162 lidxs = lidx(1:min([N length(lidx)]));
0163 lines = line(lidxs);
0164 disp(sprintf(' + found %d lines.', length([lines.f])));
0165 else
0166 disp(' + found 0 lines.');
0167 end
0168
0169
0170
0171
0172
0173 function pl = getDefaultPL()
0174
0175 pl = plist();
0176 pl = append(pl, param('N', 10));
0177 pl = append(pl, param('fsearch', [0 1e10]));
0178 pl = append(pl, param('thresh', 2));
0179 pl = append(pl, param('bw', 20));
0180 pl = append(pl, param('hc', 0.8));
0181
0182
0183
0184
0185