;Finds the next node after the given time, and identifies it as ascending
;or descending.
;
;  input
;    t_ - Time to start searching for a node, in GPS microseconds. This is usually
;         a known node time.
;  output
;    which_node - +1 if ascending node, -1 if descending node
;  return
;    Time of next node after t in GPS microseconds
;
;  Note - Will not find a node within 1 minute of the given time, so it is
;         safe to pass the exact time of a node and be assured of finding
;         the next node.
;
;         Won't handle orbits with periods of less than
;         fourteen minutes (no physical orbit will have this problem)
;
;         Function assumes that there is position info available until at least 5 minutes
;         after the next node
function find_next_node,t,which_node=which_node
  ;1. Start by adding a minute to the given time in case the given time was
  ;   exactly on a node
  t_ofs=t
  t_start=60d6
  ;2. Keep track of which hemisphere the start time is in
  z_start=find_node_f(t_start,time_ofs=t_ofs)
  ;3. Set the initial search range to 5 minutes
  range=300d6
  ;4. While the hemisphere of the end of the search range is the same as the hemisphere of the beginning...
  while z_start*find_node_f(t_start+range,time_ofs=t_ofs) gt 0 do begin
  ;   4a. Expand the search range
    range+=300d6 ;expand it by adding 5 minutes
  end
  ;5. Now the hemispheres of the beginning and end of the range are opposite, so there must be a node
  ;   in the range. Use the bisector algorithm to find the node
  t_end=t_start+range
  t_node=rtbis(t_start,t_end,0,"find_node_f",0.01,time_ofs=t_ofs)
  ;6. If the position one second after the node is positive, it's an ascending node, and vice versa
  which_node=find_node_f(t_node+1d6,time_ofs=t_ofs)
  which_node=fix(which_node/abs(which_node))
  return,t_node+t_ofs
end