;+
;PROCEDURE: plot_fa_crossing
;
;PURPOSE:
;
; Plots magnetic footprint of spacecraft across earth. (FLAT,FLNG)
; Shows the night/day terminator, auroral ovals, apogee and perigee
; footprints, etc. (See examples below.)
;
; Keywords are grouped into categories: TIMESPAN, VIEW, DISPLAY,
; AURORAL ZONE, MISCELLANEOUS.
;
; TIMESPAN KEYWORDS:
;
; ORBIT The orbit to plot. If unset, will plot over interval
; containing present time and show craft position right now
; (unless TMIN and TMAX are set). ORBIT cannot be greater than
; the last orbit listed in the predicted orbit almanac file.
; TMIN, TMAX Time points to be included on the chart. Use these if
; you want to display timespan in distant future.
; Should not be more than a couple hours apart. Format
; of TMIN and TMAX must be the type of string accepted
; by str_to_time() or a double float in seconds since
; 1970. ORBIT must not be set if these keywords are to
; be used. These times will be labeled on the map as t1
; and t2. (Good for showing AOS and LOS taken from
; contact schedules.)
; XMARK Set this to a time (string or double float) to have that point
; labeled on the map as a big X. If none of ORBIT,
; TMIN, TMAX are set then this time will be used as the
; reference time for which to create the plot. (Good
; for showing conjunctions.)
;
; VIEW KEYWORDS:
;
; KIRUNA, POKER, WALLOPS, MCMURDO, CANBERRA, SANTIAGO, BERKELEY:
; Set any of these for an overhead view of the station.
; VIEWPOINT This allows arbitrary views of the Earth. Should be a
; 3-element array designating lattitude, longitude, and
; rotation in degrees.
; WHOLE If set, will not confine plot to polar regions.
; SOUTH Set to view earth from directly over South Pole.
; MAGPOLE View from above the magnetic pole with magnetic local
; noon at the top of the plot. Gridlines are still
; geographic. Poles are not symmetric because of the
; eccentric dipole. (When this keyword is not set
; the display defaults to above the geographic poles,
; with geographic local noon at the top.)
; ZOOM Enlarges the map. 1 is normal, 2 is twice as big.
; See WINSIZE keyword for zooming in postscript mode.
;
; DISPLAY KEYWORDS:
;
; WINSIZE The width of the plot window in pixels. Default width is
; 640. Height is scaled automatically so that lattitude
; lines will be circular when output on a Tek printer:
; H = W * 1.031. In PostScript mode (POST keyword set),
; WINSIZE acts like a magnification factor; WINSIZE=10
; is normal.
; PC Set this to nonzero if you are using a windowing system
; that does not provide backing store to retain hidden
; windows. This keyword also sets the VECTOR_FONTS keyword.
; GREY Set this keyword to output the plot in greyscale.
; (Run @startup to restore colortable.)
; FILL Fills oceans and continents with solid color to make a
; pretty plot. Setting this keyword when printing to
; a color printer is discouraged.
; VECTOR_FONTS Disable switch to device (hardware) font from
; (default) Hershey vector-drawn fonts. Device fonts
; may be prettier, but vector fonts are device
; independent.
; POST Set this to a filename to direct the ouput to an 8-bit color
; postscript file instead of the graphics window. When
; viewing the postscript file, be sure your viewer is
; set to 8 bits, not 24, and that it is displaying
; "perfect colors". (Do not add ".ps" to the name.)
; GIF Captures the image on the output graphics window to a
; GIF file. (Add ".gif" to the name yourself.)
;
; AURORAL ZONE KEYWORDS:
;
; ACTIVITY Set the Activity Index of aurora. ( 0-6, default 3)
; SSCZONE Show the auroral zone used by the Satellite Situation
; Center instead of the one from Holzworth and Meng.
; (THIS IS NOT PERFECTED YET.) See auroral_zone_ssc.pro
; and notes below.
;
; MISCELLANEOUS KEYWORDS:
;
; POLAR Plots the magnetic footprint of the POLAR spacecraft
; and information about its closest approach to that of
; FAST. Good for finding conjunctions or displaying
; known ones. (For a description of how POLAR orbit data
; is obtained, see the documentation for
; get_po_orbit.pro.)
; DRAG_PROP If set, the orbit propagator will include the effects
; of atmospheric drag. The orbit track may be
; inaccurate by a few degrees for propagations as short
; as a few weeks. This keyword minimizes this error.
;ALMANAC_INFO Prints a line under the plot title telling the orbit
; file used, the last epoch in that orbit file, and
; whether or not drag was included in the orbit propagation.
;
; EXAMPLES: Show what FAST is doing now, and make it pretty:
;
; IDL> plot_fa_crossing, /fill
;
; Show what FAST is doing for Christmas, make output
; greyscale, and copy output to a GIF file:
;
; IDL> plot_fa_crossing, xmark='98-12-25/00:00:00', $
; gif='~/images/northpole.gif', /grey
;
; Show the polar crossing of orbit 3000 from a
; bird's-eye view over Canberra, and send the output to
; an 8-bit color postscript file:
;
; IDL> plot_fa_crossing, /can, post='~/images/orb3000', $
; orbit=3000
;
;NOTES ON THE
;AURORAL ZONE: The auroral zone is drawn using a procedure by Jim
; Clemens (suggested by Holzworth & Meng, and corrected
; by Mike Temerin). The source for this procedure is
; "Mathematical Representation of the Auroral Oval",
; R.H. Holzworth and C.-I. Meng, Geophysical Research
; Letters, Vol.2, No.9, Sept 1975.
;
; The mathematical representation of the auroral ovals
; is found by fitting Feldstein statistical ovals in
; corrected geomagnetic coordinates to a 7-parameter
; Fourier series:
;
; theta = A1 + A2 cos(phi + A3)
; + A4 cos(2phi + 2A5)
; + A6 cos(3phi + 3A7)
;
; theta = corrected geomagnetic co-lattitude
; phi = 2(pi)(MLT)/(24hrs)
;
; The best fit constants are found for each value of
; Q=0..6, where Q is the activity index which describes
; how quiet (0) or active (6) the Feldstein auroral
; ovals are. The characteristic radius A1 of the ovals
; increases monotonically with Q. (Source: Holzworth and
; Meng.)
;
; The function azonloc.pro generates the southern
; auroral oval by reflection of the northern oval
; through the magnetic equator.
;
; The formula in Holzworth and Meng gives the ovals in
; corrected geomagnetic coordinates. The procedure
; transform_mag_geo.pro converts them to geograhic
; coordinates using the eccentric dipole model.
;
;CREATED BY: J.Rauchleiba 96-12-20
;-
pro plot_fa_crossing, $
ORBIT=orbit, $
TMIN=time1, $
TMAX=time2, $
XMARK=xmarkin, $
ZOOM=zoom, $
DRAG_PROP=drag_prop, $
ALMANAC_INFO=almanac_info, $
WHOLE=whole, $
VIEWPOINT=pov, $
KIRUNA=kir, $
POKER=pok, $
WALLOPS=wal, $
MCMURDO=mcm, $
CANBERRA=can, $
SANTIAGO=san, $
BERKELEY=ber, $
SOUTH=south, $
MAGPOLE=magpole, $
SSCZONE=ssc, $
ACTIVITY=act, $
WINSIZE=win, $
PC=pc, $
GREY=grey, $
FILL=fill, $
VECTOR_FONTS=vdf, $
POST=psfile, $
GIF=giffile, $
POLAR=polar, $
USER_OVERPLOT = uop
thm_init
loadct2,43
cols=get_colors()
window,10,xsize=500,ysize=500
; FIDDLE WITH THE COLORS
;n_colors = !d.n_colors ; Number of colors already loaded
n_colors = !d.table_size ; Number of colors already loaded
black = 0 ; Assuming loadct2, 39
white = n_colors - 1
blue = fix(.33*float(n_colors))
if keyword_set(grey) then begin ; Make Greyscale colortable
g1 = indgen(256)
g2 = indgen(256)
g3 = indgen(256)
g1 = congrid(g1, n_colors) ; Sample n_colors of new table
g2 = congrid(g2, n_colors) ; so size of new table equal to
g3 = congrid(g3, n_colors) ; size of old table.
tvlct, g1, g2, g3
col_seas = white ; Seas
col_cont = fix(n_colors*4/5) ; Continents
col_fast = fix(n_colors/4) ; Spacecraft
col_term = fix(n_colors/2) ; Terminator line
col_oval = fix(n_colors/3) ; Auroral Ovals
col_tags = black ; Various Labels
erase, white ; PS device non-supportive
endif else begin
col_seas = blue
col_cont = fix(.81*float(n_colors))
col_fast = 6
if keyword_set(fill) then col_term=5 else col_term=blue
col_oval = 4
if keyword_set(fill) then col_tags=white else col_tags=black
;; PS device does not support erase command.
;; Must use polyfill if /FILL in this case.
if !d.name EQ 'PS' AND keyword_set(fill) then begin
polyfill, [0,1,1,0,0], [0,0,1,1,0], /normal, col=black
endif else begin
if keyword_set(fill) then erase, black else erase, white
endelse
endelse
; FONT SETUP
oldfont = !p.font
if NOT keyword_set(vdf) and NOT keyword_set(pc) then !p.font = 0
; VIEWPOINT AND PLOT LIMITS
case 1 of
keyword_set(kir): pov = [67.883,21.067,0 ]
keyword_set(pok): pov = [64.8,-147.85,0 ]
keyword_set(wal): pov = [37.9,-75.4,0]
keyword_set(mcm): pov = [-77.85,166.62,0]
keyword_set(can): pov = [-35.35,149.167,0]
keyword_set(san): pov = [-33.433,-70.667,0]
keyword_set(ber): pov = [37.9,-122.1,0]
else:
endcase
if keyword_set(whole) then latlim=[-90, 90] else latlim = [40, -40]
if keyword_set(south) then begin
lat_ref = -90.
endif else begin
if keyword_set(pov) then lat_ref = pov(0) else lat_ref = 90.
endelse
;;ns_sign = (lat_ref/abs(lat_ref) + 1)/2 ; 1:0 :: N:S
; GET THE MAGNETIC PRINT DATA OF FAST
;get_fa_orbit, tmin, tmax, /all, orbit_file=orbit_file, drag_prop=drag_prop,/no_sto, struc=fa_orbit, status=status
fa_orbitrange,orbit
trange=timerange()
orbits=fa_time_to_orbit(trange)
if orbits[0] eq orbits[1] then orbit_string = string(orbits[0]) else orbit_string=string(orbit[0])+'-'+string(orbit[1])
fa_k0_load,'orb'
;if status NE 0 then message, 'Error returned by get_fa_orbit.'
get_data, 'flat', data=tmp1
get_data, 'flng', data=tmp2
clock = tmp1.x
flat = tmp1.y
flng = tmp2.y
tmin = trange(0)
tmax = trange(1)
;clock = fa_orbit.TIME
;flat = fa_orbit.FLAT
;flng = fa_orbit.FLNG
Mtheta = flat
Mphi = flng
; Set time_N to when craft nearest N/S pole if orbit keyword set
; or median of tmin and tmax if they are specified
; or right now if no time specification.
; middle_time will be set to either of the latter two if orbit not set.
; This variable used by both the terminator and auroral oval generators.
if keyword_set(orbit) then begin
dummy = min(abs(lat_ref - Flat), sub_N)
time_N = clock(sub_N)
endif else time_N=middle_time
; GET LAT AND LNG OF SUN
center_time = time_N
if (!VERSION.RELEASE LE '5.4') then begin
hr_min = (str_sep((str_sep(time_to_str(center_time),'/'))(1),':'))(0:1)
endif else begin
hr_min = (strsplit((strsplit(time_to_str(center_time),'/',/EXTRACT))(1),':',/EXTRACT))(0:1)
endelse
gmt_hrs = float(hr_min(0)) + float(hr_min(1))/60.0
noon_long = (!pi - gmt_hrs*(!pi/12.)) MOD (2*!pi)
t0 = str_to_time('96-12-21/0:00') ;
ang = (center_time - t0)/(365.25d*24d*3600d)*2d*!pi
tilt = .410152
noon_lat = -(tilt * cos(ang))
; HANDLE SPECIAL POLE VIEWS
if keyword_set(magpole) then begin ; MAG pole view
if keyword_set(south) then begin
pov = [-75.3, 118.6, 0.0] ; override
rotcoeff = -1
rot_add = 0.0
endif else begin
pov = [82.7, -92.0, 0.0] ; override
rotcoeff = +1
rot_add = !pi
endelse
mnoon_long = noon_long -pov(1)*!dtor ; Sm pole is 92 deg W of 0.
rotate = (mnoon_long)*rotcoeff + rot_add
rerange, rotate, /deg
pov(2) = rotate
endif else if NOT keyword_set(pov) then begin ; GEO pole view
if keyword_set(south) then begin
pov = [-90., 0., 0.]
rotcoeff = -1
rot_add = 0.0
endif else begin
pov = [90., 0., 0.]
rotcoeff = +1
rot_add = !pi
endelse
rotate = noon_long*rotcoeff + rot_add
rerange, rotate, /deg
pov(2) = rotate
endif
; MAKE THE MAP (This sets up the plotting axes.)
if keyword_set(orbit) $
then maptitle = 'ORBIT ' + strtrim(string(orbit), 2)+' '+time_to_str(tmin) $
else maptitle = time_to_str(tmin) + ' ' + time_to_str(tmax)
maptitle = '!3' + maptitle
if keyword_set(zoom) then scale = 81.25e6/float(zoom)/(win/640.)
map_set, /ortho, pov(0), pov(1), pov(2), /noerase, /noborder, scale=scale, $
horizon=keyword_set(fill), e_horizon={fill:keyword_set(fill), color:blue}
xyouts, .5, .975, /norm, align=.5, color=col_tags, charsize=1.3, maptitle
map_continents, fill=fill, color=col_cont
map_grid, glinestyle=0, color=0, label=2, latlab=-20, latalign=0.5, $
lats = [-90,-75,-60,-45,-30,-15,0,15,30,45,60,75,90], $
lonlab=.75*lat_ref, lonalign=0.5, $
lons = [-180,-135,-90,-45,0,45,90,135,135,180]
; If ORBIT PROPAGATION INFO IS DESIRED
if keyword_set(almanac_info) then begin
if (!VERSION.RELEASE LE '5.4') then begin
ofpieces = str_sep(orbit_file, '/')
endif else begin
ofpieces = strsplit(orbit_file, '/', /EXTRACT)
endelse
ofname = ofpieces(n_elements(ofpieces) - 1)
if keyword_set(drag_prop) then drag_inclu='with drag.' $
else drag_inclu='without drag.'
last_epoch = time_to_str(find_last_epoch(orbit_file))
propinfo = 'Propagation: '+ofname+' ('+last_epoch+') '+drag_inclu
xyouts, .29, .955, propinfo, /norm, color=col_tags, charsize=.90
endif
; SHOW SOME GROUND STATIONS AND OTHER POINTS OF INTEREST
citylng = [-147.85,21.067,-70.667,-75.4,166.62,149.167,-122.2,-92.0,118.6]
citylat = [64.8,67.883,-33.433,37.9,-77.85,-35.35,37.9,82.7,-75.3]
citytag = ['!3 Poker',' Kiruna',' Santiago',' Wallops',' McMurdo',$
' Canberra',' Berkeley',' !6N!3',' !6S!3']
oplot, citylng, citylat, psym=2, color=black
xyouts, citylng, citylat, citytag, color=black
;;xyouts, -147.85, 64.8, '!3*Poker', color=black
;;xyouts, 21.067, 67.883, '*Kiruna', color=black
;;xyouts, -70.667, -33.433, '*Santiago', color=black
;;xyouts, -75.4, 37.9, '*Wallops', color=black
;;xyouts, 166.62, -77.85, '*McMurdo', color=black
;;xyouts, 149.167, -35.35, '*Canberra', color=black
;;xyouts, -122.2, 37.9, '*Berkeley', color=black
;;xyouts, -92.0, 82.7, '*!6N!3', color=black ; Magnetic North
;;xyouts, 118.6, -75.3, '*!6S!3', color=black ; Magnetic South
; GET POSITION OF FAST IN GSE TO FIND WHEN IN SHADOW
;;get_data, 'fa_pos', data=tmp
get_data,'r',data=tmp
;pos_arr = fa_orbit.fa_pos/6372.1 ; Position in GEI, km --> Re
store_data, 'fa_pos_re', data={x:clock, y:tmp.y/6372.1}
coord_trans, 'fa_pos_re', 'fa_pos_gse', 'GEIGSE'
get_data, 'fa_pos_gse', data=gse_stc
store_data, 'fa_pos_re', /delete
store_data, 'fa_pos_gse', /delete
; Array shade holds indices of data that are in shadow
; Craft in shade if x < 0 and y^2 + z^2 < Re^2 (1)
shade = where( gse_stc.y(*,0) LT 0 AND $
(gse_stc.y(*,1)*gse_stc.y(*,1) + gse_stc.y(*,2)*gse_stc.y(*,2)) LE 1 )
n_datapts = n_elements(gse_stc.x)
n_shadepts = n_elements(shade)
if shade(0) NE -1 then begin ; If there are points in shade
; The index to the shade array of the last shade point before sunrise:
shade_break = (where(shade NE indgen(n_shadepts)+shade(0))) (0) - 1
; data_break is index of last data point in shade before sunrise.
; If the data pts ref'd by shade not all adjacent,
; set data_break to the index of the last point in the continuous segment,
; otherwise just set it to last point ref'd by shade.
; (If there was no break in the shade then shade_break equals -2)
if shade_break NE -2 then data_break = shade(shade_break) $
else data_break = shade(n_shadepts - 1)
; Assume, if shaded (shade(0) ne -1), data_break is either the last point
; in continuous shade before sunlight (shade_break EQ -2), or is the last
; point in the first interval of shade where the shade is broken up into
; two intervals (shade_break NE -2). data_break is never the last point
; in the data array unless daybreak happens just before timespan. In this
; case, data_break = n_datapts - 1 AND shade_break EQ -2.
endif
; PLOT THE PATH OF THE CRAFT IN N AND S POLAR REGIONS (UNLESS /WHOLE).
oplot, Mphi, Mtheta, min_value=latlim(0), color=col_fast, linestyle=0, thick=3
; If we didn't already plot all the necessary points with the last plot command
if not keyword_set(whole) then $
oplot, Mphi, Mtheta, max_value=latlim(1), color=col_fast, linestyle=0, thick=3
;if defined(uop) then begin ; user-defined overplot routine
if keyword_set(uop) then begin ; user-defined overplot routine
call_procedure, uop
endif
; OVERPLOT A BLACK LINE WHERE CRAFT IS IN SHADOW
if shade(0) NE -1 then begin ; if there are some points in shade
oplot, Mphi(shade(0):data_break), Mtheta(shade(0):data_break), $
min_value=latlim(0), color=black, linestyle=0, thick=2
if not keyword_set(whole) then begin ; if we didn't aready plot all
oplot, Mphi(shade(0):data_break), Mtheta(shade(0):data_break), $
max_value=latlim(1), color=black, linestyle=0, thick=2
endif
; We must plot the remaining shady points if there was a break in the shade
if data_break LT (n_datapts - 1) AND shade_break NE -2 then begin
last_shade = shade(n_elements(shade) - 1)
oplot, Mphi(shade(shade_break+1):last_shade), $
Mtheta(shade(shade_break+1):last_shade), $
min_value=latlim(0), color=black, linestyle=0, thick=2
if not keyword_set(whole) then begin
oplot, Mphi(shade(shade_break+1):last_shade), $
Mtheta(shade(shade_break+1):last_shade), $
max_value=latlim(1),color=black,linestyle=0,thick=2
endif
endif
endif
if 0 then begin ; this is no longer needed
; Show print of craft right now if "now" is between tmin, tmax.
; Write text to plot window telling UT, ORBIT, FLAT, FLNG right now.
if (plot_craft) then begin
dt = min( abs(clock - now_time), pos_ind) ; pos_ind=nowdataindex
oplot, [Flng(pos_ind)], [Flat(pos_ind)], psym=6, syms=3, col=col_fast
;;get_data, 'ORBIT', data=orbit_nums
get_data
; orbit_nums = fa_orbit.ORBIT
orbit_nums = fa_orbit.ORBIT
current_orb = orbit_nums(pos_ind)
if (!VERSION.RELEASE LE '5.4') then begin
xyouts, .8, .95, /norm, color=col_tags, $
'CURRENT FAST INFO' + $
'!C------------' + $
'!CUT: ' + (str_sep(time_to_str(now_time), '/'))(1) + $
'!CORBIT: ' + strtrim(string(current_orb), 2) + $
'!CFLAT: ' + string(format='(F6.1)', Flat(pos_ind)) + $
'!CFLNG: ' + string(format='(F6.1)', Flng(pos_ind))
endif else begin
xyouts, .8, .95, /norm, color=col_tags, $
'CURRENT FAST INFO' + $
'!C------------' + $
'!CUT: ' + (strsplit(time_to_str(now_time), '/', /EXTRACT))(1) + $
'!CORBIT: ' + strtrim(string(current_orb), 2) + $
'!CFLAT: ' + string(format='(F6.1)', Flat(pos_ind)) + $
'!CFLNG: ' + string(format='(F6.1)', Flng(pos_ind))
endelse
endif
endif
;; ADD TIME TICKS, LABELS TO PATH
label_foot_ticks, time=clock, latitude=Flat, longitude=Flng, interval=300, $
latlim=45, color=col_tags
; LABEL TIME1 AND TIME2 IF SUPPLIED EXPLICITLY BY USER
if keyword_set(time1s) AND keyword_set(time2s) AND NOT keyword_set(orbit) $
then begin
dt1 = min( abs(clock - time1s), time1_ind )
dt2 = min( abs(clock - time2s), time2_ind )
xyouts, Flng(time1_ind), Flat(time1_ind), '.t1', $
align=0.0, color=col_tags
xyouts, Flng(time2_ind), Flat(time2_ind), '.t2', $
align=0.0, color=col_tags
endif
; PUT AN X MARK AT XMARK
if keyword_set(xmark) then begin
dtx = min( abs(clock - xmark), x_ind )
oplot, [Flng(x_ind)], [Flat(x_ind)], psym=7, symsize=3.0, color=col_tags
endif
; PLOT THE TERMINATOR
terminator, time=time_N, tlat, tlng ; tlat, tlng as yet undefined
; The use of max_value and min_value below is a patch to remove ugly line
if pov(0) EQ 90 then min_term = 10 else min_term = -90
if pov(0) EQ -90 then max_term = -10 else max_term = 90
oplot, tlng, tlat, thick=3, color=col_term, min_val=min_term, max_val=max_term
; APOGEE, PERIGEE
;;get_data, 'ALT', data=alt
;;get_data, 'LAT', data=Glat
get_data,'alt',data=tmp & alt=tmp.y
get_data,'ilat',data=tmp & ilat=tmp.y
;alt = fa_orbit.ALT
;Glat = fa_orbit.LAT
max_alt = max(alt, app_sub)
min_lat = min(alt, per_sub)
xyouts, Flng(app_sub), Flat(app_sub), 'A', /data, color=col_tags
xyouts, Flng(per_sub), Flat(per_sub), 'P', /data, color=col_tags
;apogee_lat = Glat(app_sub) ; Geographic lattitude of apogee
apogee_lat = ilat(app_sub) ; Invariant lattitude of apogee
perigee_lat = ilat(per_sub) ; Invariant lattitude of apogee
; AURORAL OVALS
IF n_elements(act) NE 1 then act = 3 ; Default activity index
; Get Corrected Geomagnetic Lattitudes for all MLT
mlt = findgen(241)/10. ; [.1, .2, ..., 23.9, 24.0]
if keyword_set(ssc) then begin
nelat = auroral_zone_ssc(mlt, poleward=nplat)
selat = auroral_zone_ssc(mlt, poleward=splat, /south)
endif else begin
nelat = auroral_zone(mlt,act,/lat) ; Corr. Geomag. Lat for all MLT
nplat = auroral_zone(mlt,act,/lat,/pole) ; Northern poleward
selat = auroral_zone(mlt,act,/lat,/south) ; Southern equatorward
splat = auroral_zone(mlt,act,/lat,/south,/pole) ; Southern poleward
endelse
if (!VERSION.RELEASE LE '5.4') then begin
date_time = str_sep(time_to_str(time_N), '/') ; time_N is reference time
year = fix(strmid(date_time(0), 0, 4)) ; 4-digit integer
N_hms = fix(str_sep(date_time(1),':')) ; [hh,mm,ss]
UT_hrs = N_hms(0) + N_hms(1)/60. ; Hours into UT day
endif else begin
date_time = strsplit(time_to_str(time_N), '/', /EXTRACT) ; time_N is reference time
year = fix(strmid(date_time(0), 0, 4)) ; 4-digit integer
N_hms = fix(strsplit(date_time(1),':', /EXTRACT)) ; [hh,mm,ss]
UT_hrs = N_hms(0) + N_hms(1)/60. ; Hours into UT day
endelse
; Convert MLT to MLNG using UT and get MLNG of Sun
sun_lng = [!pi - UT_hrs*(!pi/12.)] ; LNG of Sun
sun_lat = [0.]
mag_to_geo, sun_lat, sun_lng, /mag ; MLNG of Sun
azon_mlng = mlt*(!pi/12.) + sun_lng(0) - !pi ; MLNG of oval
nelon = azon_mlng
nplon = nelon
selon = nelon
splon = nelon
; Convert auroral ovals from MAG to GEO
transform_mag_geo, nelat, nelon, tnelat, tnelon, year=year
transform_mag_geo, nplat, nplon, tnplat, tnplon, year=year
transform_mag_geo, selat, selon, tselat, tselon, year=year
transform_mag_geo, splat, splon, tsplat, tsplon, year=year
rerange, tnelon, tnelat, /deg
rerange, tnplon, tnplat, /deg
rerange, tselon, tselat, /deg
rerange, tsplon, tsplat, /deg
; Plot the ovals
oplot, tnelon, tnelat, color=col_oval, thick=2
oplot, tnplon, tnplat, color=col_oval, thick=2
oplot, tselon, tselat, color=col_oval, thick=2
oplot, tsplon, tsplat, color=col_oval, thick=2
; Put a dot on the ovals at MLT=12, 24
dot_lng = [tnelon(240), tnplon(240), tselon(240), tsplon(240), $
tnelon(120), tnplon(120), tselon(120), tsplon(120)]
dot_lat = [tnelat(240), tnplat(240), tselat(240), tsplat(240), $
tnelat(120), tnplat(120), tselat(120), tsplat(120)]
oplot, psym=4, symsize=.5, dot_lng, dot_lat, color=black
; COLOR-CODED LEGEND, LABELS, ETC.
leg_loc = convert_coord(.02, .95, /norm, /to_device)
line_xpos = make_array(8, /int, value=leg_loc(0))
line_ypos = intarr(8)
if !d.name EQ 'PS' then linewidth=500 else linewidth=12
for p=0,7 do line_ypos(p) = leg_loc(1) - p*linewidth
colors = [col_tags, col_tags, $
col_fast, col_fast, $
col_tags, col_tags, $
col_term, col_oval]
ds = string("260B) ;" degree symbol
legend = ['LEGEND', $
'----------', $
; '[] FAST Position Now', $
'+ FAST Path', $
'A Apogee (' + string(format='(F5.1)', apogee_lat) + ds + ' ilat)', $
'P Perigee(' + string(format='(F5.1)', perigee_lat) + ds + ' ilat)', $
'- Terminator', $
'- Auroral Oval']
xyouts, line_xpos, line_ypos, /device, legend, color=colors
; Tell when craft enters and exits shade. Must have shade(0) NE -1
; True transition from shade to sunlight if data_break NE n_datapts-1
; True transition from sunlight to shade if shade(0) NE 0.
if shade(0) NE -1 then begin
if (!VERSION.RELEASE LE '5.4') then begin
if data_break NE (n_datapts - 1) $
then ecl_ex = 'Eclipse Ext: '+(str_sep(time_to_str(gse_stc.x(data_break)), '/'))(1) $
else ecl_ex = ''
if shade(0) NE 0 $
then ecl_en = 'Eclipse Ent: '+(str_sep(time_to_str(gse_stc.x(shade(0))), '/'))(1) $
else ecl_en = ''
endif else begin
if data_break NE (n_datapts - 1) $
then ecl_ex = 'Eclipse Ext: '+(strsplit(time_to_str(gse_stc.x(data_break)), '/', /EXTRACT))(1) $
else ecl_ex = ''
if shade(0) NE 0 $
then ecl_en = 'Eclipse Ent: '+(strsplit(time_to_str(gse_stc.x(shade(0))), '/', /EXTRACT))(1) $
else ecl_en = ''
endelse
ecl_label = [ecl_ex, ecl_en]
leg_loc = convert_coord(.98, .05, /norm, /to_device)
line_xpos = make_array(2, /int, value=leg_loc(0))
line_ypos = intarr(2)
if !d.name EQ 'PS' then linewidth=500 else linewidth=12
for p=0,1 do line_ypos(p) = leg_loc(1) - p*linewidth
xyouts, line_xpos, line_ypos, /device, ecl_label, color=col_tags, align=1.0
endif
; Write TMIN and TMAX, or XMARK
if keyword_set(time1) AND keyword_set(time2) then begin
if (!VERSION.RELEASE LE '5.4') then begin
t12_label = ['t1 = '+(str_sep(time_to_str(time1s),'/'))(1), $
't2 = '+(str_sep(time_to_str(time2s),'/'))(1) ]
endif else begin
t12_label = ['t1 = '+(strsplit(time_to_str(time1s),'/',/EXTRACT))(1), $
't2 = '+(strsplit(time_to_str(time2s),'/',/EXTRACT))(1) ]
endelse
leg_loc = convert_coord(.02, .05, /norm, /to_device)
line_xpos = make_array(2, /int, value=leg_loc(0))
line_ypos = intarr(2)
if !d.name EQ 'PS' then linewidth=500 else linewidth=12
for p=0,1 do line_ypos(p) = leg_loc(1) - p*linewidth
xyouts, line_xpos, line_ypos, /device, t12_label, color=col_tags
endif
if keyword_set(xmark) then xyouts, .02, .05, /normal, col=col_tags, $
'X = '+time_to_str(xmark)
; Add POLAR spacecraft footprint if requested
if keyword_set(polar) then begin
fast_flat = double(flat*!dtor)
fast_flng = double(flng*!dtor)
fast_times = clock
plot_po_overlay, tag_color=col_tags, $
time=clock, FLAT=polar_flat, FLNG=polar_flng
;;get_data, 'FLAT', data=po_flat_stc
;;get_data, 'FLNG', data=po_flng_stc
polar_flat = polar_flat*!dtor
polar_flng = polar_flng*!dtor
;; Find closest approach of footprints
fx = cos(fast_flat)*cos(fast_flng)
px = cos(polar_flat)*cos(polar_flng)
fy = cos(fast_flat)*sin(fast_flng)
py = cos(polar_flat)*sin(polar_flng)
fz = sin(fast_flat)
pz = sin(polar_flat)
linesep = min(sqrt((fx-px)^2 + (fy-py)^2 + (fz-pz)^2), cji)
;;angsep = asin(linesep/2d)*!radeg
;;print, 'Geocent Ang Sep: ', time_to_str(fast_times(cji)), angsep, ds
;; Label Conjunction
if (!VERSION.RELEASE LE '5.4') then begin
cnj_time = (str_sep(time_to_str(fast_times(cji)), '/'))(1)
endif else begin
cnj_time = (strsplit(time_to_str(fast_times(cji)), '/', /EXTRACT))(1)
endelse
d_lat = string((fast_flat(cji) - polar_flat(cji))*!radeg, format='(F5.1)')
d_lng = string((fast_flng(cji) - polar_flng(cji))*!radeg, format='(F5.1)')
cnj_loc = convert_coord(.02, .17, /normal, /to_device)
cnj_xpos = make_array(4, /int, value=cnj_loc(0))
cnj_ypos = intarr(4)
for p=0,3 do cnj_ypos(p) = cnj_loc(1) - p*linewidth
cnj_stg = ['CONJUNCTION:', cnj_time, '!4D!3Lat '+d_lat, '!4D!3Lng '+d_lng]
xyouts, cnj_xpos, cnj_ypos, /device, cnj_stg, color=col_tags
;;oplot, [fast_flng(cji),polar_flng(cji)]*!radeg, $
;; [fast_flat(cji),polar_flat(cji)]*!radeg,color=col_oval
endif
; Close postscript file if opened
if keyword_set(psfile) then pclose
; Write to GIF file if desired
if keyword_set(giffile) AND (!d.name EQ 'X' OR !d.name EQ 'Z') then begin
tvlct, /get, roy, gee, biv
write_gif, giffile, tvrd(), roy, gee, biv
endif
; Reset to original font system and plot window
!p.font = oldfont
if keyword_set(orig_plot_pos) then !p.position=orig_plot_pos
; Message
print, 'Plot_fa_crossing graphical interface: PFCGUI.PRO'
end