function [Nsealed, xl, yl,  xr, yr] = sealedVolumes(tup, tdown, xL, yL, xR, yR, Z)
%UNTITLED2 Summary of this function goes here
%   Detailed explanation goes here

tup = pi - deg2rad(tup); tdown = pi + deg2rad(tdown); %angular limits
thL = zeros(1, length(xL)); thR = zeros(1, length(xR));
thLmin = 10000; thRmax = -10000;

for ii = 1 : length(xL)
   % if yL(ii) > 0
        thL(ii) =  mod(atan2(yL(ii), xL(ii)), 2*pi);
        if thL(ii) >= tup && thL(ii) < thLmin
           thLmin = thL(ii);
           IDl = ii;
        end
   % end
end
for ii = 1 : length(xR)
    if yR(ii) < 0
        thR(ii) = mod(atan2(yR(ii), xR(ii)), 2*pi);
        if thR(ii) <= tdown && thR(ii) > thRmax
           thRmax = thR(ii);
           IDr = ii;
        end
    end
end

%noti IDl, IDr
 if IDr < IDl, Nsealed = IDl - IDr;
 else Nsealed = Z - (Idr - IDl);
 end

 xl=xL(IDl); yl=yL(IDl);
 xr=xR(IDr); yr=yR(IDr);

 end

