function [points,meshlist]=getpointsmesh(oldpoints,im,cc,fc,kc,omega);
% function [points,meshlist]=getpointsmesh(oldpoints,im,cc,fc,kc,omega);
%	Searches the image row-by-row and looks for one point per row.
%	Point is cacluated to sub-pixel accuracy and triangulated using
%	cc,fc,kc and omega.
%   Appends triangulated point to oldpoints matrix.
%	When finished, returns new points buffer (oldpoints + new points)
%	and a 1xR matrix meshlist.
%	The elements of meshlist are: 0 if no point found
%	along that row, or the index in matrix (points) of the triangulated point
%	found along that row.
%

[r,c]=size(im);

thresh = 50;	% min value for a peak in a row

[xx,yy] = size(oldpoints);

firstndx = yy + 1;
newpoints = zeros(3,r);	% can't generate more points than this
npoints = 0;
meshlist = zeros(1,r);

for i=1:r,
	% process row
	data = im(i,:);
	if (sum(data .* (data > thresh)) > 0),
		% there is at least one pixel above thresh
		data = conv(data,blackman(5));
		data = data(3:length(data));	% shift back over to compensate for filtering
		[y,ndx]=max(data);
		[a,x0,y0]=fitparabola(data(ndx-1:ndx+1));
		pt = [ndx+x0 ; i-1];
		pt = comp_distortion((pt - cc) ./ fc, kc);
		depth = 1 / dot([pt ; 1],omega);
		pt3d = [(pt * depth) ; depth];
		npoints = npoints + 1;
		newpoints(:,npoints) = pt3d;
		meshlist(i) = firstndx + npoints - 1;
	end
end

if npoints > 0,
	newpoints = newpoints(:,1:npoints);
	points = [oldpoints newpoints];
else
	points = oldpoints;
end