function [x, y] = findmaxima(in)

npts = length(in);
x = zeros(1,npts/4);
y = x;
wasSmaller = in(2) > in(1);
constStart = 0;
if in(2) == in(1); constStart = 1; end
idx = 1;

for i=2:npts-1
    isSmaller = in(i+1) > in(i);
    isConstant = in(i+1) == in(i);
    
    if constStart==0 & ~isConstant
        if wasSmaller ~= isSmaller
            if wasSmaller
                diff = in(i+1)-in(i-1);
                den = (2*in(i)-in(i-1)-in(i+1));
                x(idx) = i + diff / (2*den);
                y(idx) = in(i) + (diff*diff) / (8*den);
                idx=idx+1;
            end
            wasSmaller = isSmaller;
        end 
    elseif constStart~=0 & isConstant
        continue; % Don't touch wasSmaller!
    elseif constStart~=0 & ~isConstant
        % We haven't touched wasSmaller since the constant zone started!
        if wasSmaller ~= isSmaller
            if wasSmaller
                % put a maximum at the middle point
                x(idx) = 0.5*(constStart+i);
                y(idx) = in(i);
                idx=idx+1;
            end
            wasSmaller = isSmaller;
        end % otherwise it's just a spurious constant value along a sloped region.
        constStart = 0;
    elseif constStart==0 & isConstant
        constStart = i;
    end
end

% delete extra points
x(idx:end) = [];
y(idx:end) = [];

if nargout<2
    x(2,:) = y(1,:);
end