;+
; NAME:
;	PLOTORB
; PURPOSE:
;       Plotting satellite trajectories onto Polar Camera images
;       (more generally, generating satellite coordinates for
;	comparison with images)
; CATEGORY:
;	Particle-optical comparison
; CALLING SEQUENCE:
;	PLOTORB, SAT, DATE, START, ENDT, INTERVAL $
;              , OTIME, OLAT, OLONG, OTHER $
;             [, DEVICE, TYPE = TYPE, /NOPLOT, FPALT = FPALT ]
; INPUTS:
;	SAT:	Satellite name in standard 2-line element form
;	DATE:	Date in YYMMDD format
;	START,
;	ENDT:	Start and end times in HHMMSS format
;	INTERVAL: Interval (minutes!) between tabulated points
; OPTIONAL INPUTS:
; KEYWORD PARAMETERS:
;	DEVICE:	Output device (either 'X' (default) or 'PS')
;	TYPE:	Specifies how the satellite position is expressed:
;		0: geocentric equatorial inertial cartesian (km)
;		1: satellite geocentric latitude/longitude/radius
;		2: satellite geodetic   latitude/longitude/alt.
;		3: footprint geocentric latitude/longitude/radius
;		4: footprint geodetic   latitude/longitude/alt.
;		5: footprint eccentric dipole lat/long/MLT
;       /NOPLOT: Just return the results, don't try to plot them
;		onto an image
;       FPALT:  IF TYPE = 3, 4, or 5, this determines the altitude
;		of the footprint
; OUTPUTS:
;	OTIME:	The times of the tabulated positions (HHMMSS)
;	OLAT:	The latitude coordinate (or the GEI X component)
;	OLONG:	The longitude coordinate (or the GEI Y component)
;	OTHER:	The altitude/MLT coordinate (or the GEI Z component)
; OPTIONAL OUTPUTS:
;	None.
; COMMON BLOCKS:
;	None.
; SIDE EFFECTS:
;       A file input.tmp is created to provide input to the orbit
;       calculation program.  A file output.tmp is created by
;	redirecting the output of the orbit calculation program.
; RESTRICTIONS:
;       Elements are available only up to early January 1995 as of
;	this time (950615).
; PROCEDURE:
;	
; EXAMPLE:
;	
; SEE ALSO:
;	
; MODIFICATION HISTORY:
; 	Written by:	DPS, spring 1993
;	Documented by:	DPS, June 1995
;-

PRO plotorb,sat,date,start,endt,interval,otime,olat,olon,other $
	   ,device=device,type=type,noplot=noplot,fpalt=fpalt
;
; This procedure plots the magnetic footprint (at an altitude
; specified by the user; default 110 km) of the specified satellite
; on the given date between the specified start and end times (UT),
; plotting points at the specified intervals.  Output can be directed
; to the specified device.
; 
; Argument formats:
;   sat:	must agree exactly with satellite name in tle.### files
;   date:	YYMMDD
;   start:	HHMMSS (UT) of first footprint point
;   endt:	HHMMSS (UT) of last footprint point; if this time
;               occurs after 24 UT on 'date', specify it as a time
;		> 240000.
;   interval:	minutes between plotted points (can be fractional)
;   device:	device name as accepted by SET_PLOT command
;   noplot:     keyword indicating that instead of plotting the
;               path of the selected satellite on a grid around
;               Rabbit Lake, the path of the satellite is to be
;               returned to the calling routine in the variables
;		ctime, cgdlat, and cgdlon.
;   fpalt:	keyword indicating type of information to return;
;               fpalt >= 0 indicates satellite field line coordinates
;                          at altitude fpalt;
;		fpalt <  0 indicates satellite coordinates;
;
; Algorithm:
; build input file for ~/steele/track/spacetrk/bin/spmod
; use SPAWN command to run spmod with given input file and redirect
;   output to another file (shd the file name be a parameter, or
;   shd it be named to reflect the inputs?)
; select the output device (default to 'X')
; read in the output file and plot the footprint using the MAP calls
; overplot the location and approximate field of view of Rabbit Lake
; if PostScript selected, close the file and send it to the printer
;
; DPS, 28 April 1993

IF N_PARAMS() LT 8 THEN BEGIN
    doc_library,'plotorb'
    RETURN
    ENDIF

; if no plot needed, the following information is not required
IF NOT KEYWORD_SET(noplot) THEN BEGIN
; ascertain user's intentions, if sat = 'UARS'
    resp=''
    IF sat EQ 'UARS' THEN BEGIN
        PRINT,'Do you want the large neighbourhood of Rabbit Lake? (Y/[N])'
        READ,resp
        resp=STRUPCASE(resp)
    ENDIF
ENDIF
; is the footprint altitude specified yet?
IF N_ELEMENTS(fpalt) EQ 0 THEN fpalt=110.
; start by creating spmod input file
OPENW,ip,'~/input.tmp',/GET_LUN
PRINTF,ip,sat
PRINTF,ip,date,start,FORMAT='(I6,".",I6.6)'
PRINTF,ip,date,endt,FORMAT='(I6,".",I6.6)'
PRINTF,ip,interval,FORMAT='(F5.2)'
PRINTF,ip,type,FORMAT='(I1)'
IF type GE 3 THEN BEGIN
    IF N_ELEMENTS(fpalt) EQ 0 THEN fpalt=250.	; reasonable for Eureka
    PRINTF,ip,FLOAT(fpalt),FORMAT='(F4.0)'	; footprint altitude
    ENDIF
PRINTF,ip,FORMAT='()'	; carriage return to autoselect tle file
CLOSE,ip
FREE_LUN,ip
; now spawn a UNIX command to calculate the required footprint
; coordinates
SPAWN,'~steele/track/spacetrk/bin/spmod < ~/input.tmp >! ~/output.tmp',orbres
; read in the coordinates
OPENR,op,'~/output.tmp',/GET_LUN
textline=''
numline=FLTARR(5)
inputs=TRANSPOSE(numline)
REPEAT BEGIN
    READF,op,textline 
    ENDREP UNTIL STRPOS(textline,'Required') GT 0
tlefile=getwrd(textline,5)
edate=0.0D0
READF,op,yr,edate,FORMAT='(18X,I2,F12.8)'
yr=1900+yr
FOR i=0,2 DO READF,op,textline
REPEAT BEGIN
    READF,op,numline
    inputs=[inputs,TRANSPOSE(numline)]
ENDREP UNTIL EOF(op)
CLOSE,op
FREE_LUN,op
; extract file contents from inputs
inputs=inputs(1:*,*)
cdate=inputs(*,0)
otime=inputs(*,1)
olat=inputs(*,2)
olon=inputs(*,3)
olon=olon-(olon GT 180.)*360.
other=inputs(*,4)
; if no plot desired, return satellite path to calling program
IF KEYWORD_SET(noplot) THEN RETURN
; restrict plotted coordinates to points inside lat/lon grid
IF resp NE 'Y' THEN BEGIN
    wedge=-114.
    eedge=-92.
    nedge=62.
    sedge=54.
ENDIF ELSE BEGIN
    wedge=-126.
    eedge=-80.
    nedge=67.
    sedge=49.
ENDELSE
ingrid=WHERE((cgdlon GE wedge) AND (cgdlon LE eedge) AND $
             (cgdlat GE sedge) AND (cgdlat LE nedge),ningrid)
IF ningrid GT 0 THEN BEGIN
    cdate=cdate(ingrid)
    ctime=ctime(ingrid)
    cgdlat=cgdlat(ingrid)
    cgdlon=cgdlon(ingrid)
; set up device
    IF N_ELEMENTS(device) EQ 0 THEN device='X'
    IF STRUPCASE(device) EQ 'PS' THEN psopen,'/tmp/tmp.ps'
; define Rabbit Lake coordinates
    rllat=58.233
    rllon=256.330-360.0
; create map title
    mtitle=STRING(sat,date,start,endt $
          ,FORMAT='(A," near Rabbit Lake, SK, ",I6,", ",I6.6," - ",I6.6," UT")')
; plot map
    MAP_SET,rllat,rllon,/ORTHOGRAPHIC,LIMIT=[sedge,wedge,nedge,eedge] $
           ,/GRID,/LABEL,TITLE=mtitle,LATLAB=-110,LATALIGN=0.0,LONLAB=55 $
           ,LONALIGN=1.0,/NOBORDER,LATDEL=1.0,LONDEL=2.0,XMARGIN=[3,3] $
           ,YMARGIN=[3,3],CHARSIZE=0.75,/CONTINENTS,/USA
    prov
    USERSYM,[1,0,-1,0],[0,1,0,-1],/FILL
    PLOTS,cgdlon,cgdlat,PSYM=-8,NOCLIP=0
    atmin=WHERE(LONG(ctime) MOD 100 EQ 0)
    XYOUTS,cgdlon(atmin),cgdlat(atmin)+0.1,STRTRIM(LONG(ctime(atmin)),2) $
          ,ALIGNMENT=0.5,CHARSIZE=0.75
    ydn2md,yr,FIX(edate),emon,eday
    XYOUTS,rllon,sedge-(nedge-sedge)/22.86,ALIGNMENT=0.5 $
          ,STRING(tlefile,yr,emon,eday $
	  ,FORMAT='(A,": ",I4,"/",I2.2,"/",I2.2)')
; show Rabbit Lake position
    PLOTS,rllon,rllat,PSYM=1
    XYOUTS,rllon-0.4,rllat,'Rabbit Lake',ALIGNMENT=1.
; overplot the horizon at Rabbit Lake
    a=6371.
    nza=4
    hang=!DTOR*15.*(FINDGEN(nza)+1)
    gcrad=hang-ASIN((a*SIN(hang))/(a+fpalt))
    naz=37
    az=FINDGEN(naz)*360./(naz-1)
    hlon=FLTARR(naz,nza)
    hlat=hlon
    FOR i=0,naz-1 DO BEGIN
        rb2ll,rllon,rllat,gcrad,az(i),alon,alat
        hlon(i,*)=alon
        hlat(i,*)=alat
    ENDFOR
    FOR i=0,nza-1 DO BEGIN
        PLOTS,hlon(*,i)-360.,hlat(*,i)
        XYOUTS,hlon(0,i),hlat(0,i)+0.05,ALIGNMENT=0.5,CHARSIZE=0.75 $
              ,STRTRIM(FIX(!RADEG*hang(i)+.001),2)+'!9%!x'
	IF i EQ nza-1 THEN XYOUTS,hlon(0,i),hlat(0,i)+0.25 $
				 ,ALIGNMENT=0.5,CHARSIZE=0.75 $
       			,STRING(fpalt,FORMAT='("(",I3," km altitude)")')
    ENDFOR
; end plot
    IF !D.NAME EQ 'PS' THEN BEGIN
        psclose
        SPAWN,'rmps3 /tmp/tmp.ps b'
    ENDIF
    SPAWN,'rm ~/input.tmp ~/output.tmp'
ENDIF ELSE PRINT,sat+' not visible from Rabbit Lake at these times.'
RETURN
END