function orbits = nscat_25_orbits(aE, aW, dE, dW)

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Determine which orbits crossed region (using inventory of NSCAT files  with 
% [orbit time 1st equatorial crossing] )
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% load NSCAT_inv_new
load INV
T1 = INV(:,2);				% 1st data Time 

% First, subset the inventory based on time limits.

% Second, subset the inventory based on 1st equatorial crossing longitude.
% This requires different logic for two cases where eastern boundary
% determined by eqlonlim may be greater or less than western boundary, i.e.
% whether or not zero longitude is crossed by the equatorial limits
% determined for ascending and descending passes.

eqX1 = INV(:,3);	% Longitude of 1st rev equatorial crossing

% Case 1 : neither ascending nor descending limits cross zero longitude
if ((aW < aE) & (dW < dE))
  idx =  find((eqX1 >= aW & eqX1 <= aE ) | (eqX1 >= dW & eqX1 <= dE));

% Case 2: ascending limits do not cross zero longitude; descending limits do 

elseif ((aW < aE) & (dW > dE))
  idx =  find((eqX1 >= aW & eqX1 <= aE ) | (eqX1 >= dW) | (eqX1 <= dE));

% Case 3 : ascending limits cross zero longitude; descending do not

elseif ((aW > aE) & (dW < dE))
  idx =  find((eqX1 >= aW) | (eqX1 <= aE) | (eqX1 >= dW & eqX1 <= dE));

% Case 4 : both ascending and descending limits cross zero longitude

else
  idx =  find((eqX1 >= aW) | (eqX1 <= aE) | (eqX1 >= dW) | (eqX1 <= dE));
end

% YO JAMES : USE FULL INVENTORY FOR TESTING IF IT WORKS
%% SUBSET
% subset_inv = INV(idx,:);
%% DON'T SUBSET
%subset_inv = INV;	% USE FULL INVENTORY FOR TEST ONLY

% subset inv : filenames 69 to 84 for JD 117 123
% Full inventory : orbits 446 to 496 cover JD 119 to 123
% time_to_GUI = subset_inv(:,2);

orbits = INV(idx,:);

return;
