Pro high_mag_lat, timerange, latm, times, lowlat=lowlat, oval=oval, hc=hc, debug=debug
;+
;NAME:
;      high_mag_lat  
;PURPOSE:
;       To compute magnetic latitude of the orbit and plot the times
;	of high latitude and passage through the ovals
;
;SAMPLE CALLING SEQUENCE:
;       high_mag_lat,['21-may-98 12:00','22-may-98 12:00'],latm,times $
;	  lowlat=45, oval=[50,70]
;INPUTS:
;	timerange - Array of start and stop dates & times
;	lowlat - lower limit for high latitude region (default = 40 degrees mag)
;	oval - array of low and upper magnetic latitude limits for oval 
;	  boundaries (default = [50,70])
;	Assumes North magnetic pole at latitude = 78.5 and longitude = -70
;OUTPUTS
;	latm  - output array of magnetic latitude
;	times - output array of times for latm
;HISTORY:
;        Written 19-May-98 by TDT
;	21-Oct-2002 (TDT)  Fix no-SAA bug, change SAA symbol
;-
version='HIGH_MAG_LAT Vers. 1.1'
clear_utplot
thetam=(90-78.5)/!radeg
phim=-70./!radeg
if (not keyword_set(lowlat)) then lowlat=40.
if (not keyword_set(oval)) then oval=[50.,70.]
if (not keyword_set(debug)) then debug=0

times=timegrid(timerange(0),timerange(1),minutes=1,/string)
saa=gtt_orbit(times,'is_saa')
wsaa=where(saa gt 0,nsaa)
lat=gtt_orbit(times, 'lat')
lon=gtt_orbit(times, 'lon')
phi=lon/!radeg
theta=(90-lat)/!radeg
if (debug) then begin 
  utplot,times,lat & wait,2
  utplot,times,theta & wait,2
endif
dot=cos(theta)*cos(thetam)+sin(theta)*sin(thetam)*(cos(phi)*cos(phim)+sin(phi)*sin(phim))
;utplot,times,dot
latm=90.-acos(dot)*!radeg
if (debug) then begin 
  utplot,times,latm
  outplot,times,lat,line=1  &  wait,2
endif
w=where(abs(latm) gt lowlat) 
highm=0*lat  &  highm(w)=1
;utplot,times(0:120),highm(0:120),yr=[-2,2]
;outplot,times(0:120),latm(0:120)/90.
w=where((abs(latm) gt oval(0)) and (abs(latm) lt oval(1))) 
ovalm=0*lat  &  ovalm(w)=1
!p.title='TRACE High Magnetic Latitude (solid), Oval Passage (dot), SAA (diamonds)'
!p.subtitle='Lower Boundary = '+string(lowlat,'(F6.1)')+' degrees,  Ovals = '+$
  string(oval,'(2F6.1)')+' degrees      Made by '+version
!p.psym=10
utplot,times,highm,yr=[0,1.2]
outplot,times,0.95*ovalm,line=1
if (nsaa gt 0) then outplot,times(wsaa),saa(wsaa)*0.5,psym=4
if  (not keyword_set(hc)) then return
set_plot,'ps'
device,/landscape
utplot,times,highm,yr=[0,1.2]
outplot,times,0.95*ovalm,line=1
if (nsaa gt 0) then outplot,times(wsaa),saa(wsaa)*0.5,psym=4
device,/close
set_plot, 'x
sprint,'idl.ps'
end