function [s, z] = seplam_initial(A, B, prtlevel)
% search along lines parallel to the real axis passing through all 2n eigenvalues
% of A and B, looking for ALL real positive s and z on these real lines such that 
% sigma_min(A - zI) = sigma_min(B - zI) = s.  Return the smallest such s and
% the corresponding z.  This ensures that we don't get stuck with a component of
% A's eps-pseudospectrum strictly contained inside a component of B's eps-pseudospectrum,
% for some eps, thus defeating "realimagcheck".

% could make this more efficient
%   can save factor of 2 by looking along lines joining eigenvalues
%   can avoid repeat checks along same line, e.g. if eigenvalues are real
eA = eig(A);
eB = eig(B);
eigtol = 1e-8;
for j=1:length(A)
   [s,z] = linecheck(A, B, eA(j), eigtol);
   [val,indx] = min(s);
   sminA(j) = val;
   zminA(j) = z(indx);
   if prtlevel > 1
       fprintf('\nfound initial upper bound %e on line %d for A', val, j)
   end
end
for j=1:length(B)
   [s,z] = linecheck(A, B, eB(j), eigtol);
   [val,indx] = min(s);
   sminB(j) = val;
   zminB(j) = z(indx);
   if prtlevel > 1
       fprintf('\nfound initial upper bound %e on line %d for B', val, j)
   end
end
[valA,indxA] = min(sminA);
[valB,indxB] = min(sminB);
if valA < valB
   s = valA;
   z = zminA(indxA);
else
   s = valB;
   z = zminB(indxB);
end
% possible that s is inf and z is nan, if no s and z were found on
% any of these lines
