;    unix version.
;    copyright 1999 by Robert J. Nemiroff
;                                      Hard or Soft Copy?
dev = 'n'
read, ' Type of plot (x=window, ps=postscript, gif=gif) ? ', dev
  if (dev eq 'ps') then begin        ; Hardcopy wanted. 
  set_plot, dev
  device, /color 
  endif
;
    if (dev eq 'x') or (dev eq 'gif') then begin
    set_plot, 'x'                    ; X-window display wanted.
    window, 0, retain=2              ; Ability to restore window.
    device, decomposed=0, retain=2
    endif
;                                      Erase previous, save present.
erase & !p.noerase=1
;!FANCY = 1   & !TYPE = 16 + 8 + 4   ; Old IDL commands now obselete.
;                                      Variable ranges.
xmin=-1   &   xmax=1
ymin=-1   &   ymax=1
;                                      Number of tick marks.
!x.ticks = 2 & !y.ticks = 2
!x.tickname = [' ', ' ', ' ']
!y.tickname = [' ', ' ', ' ']
;                                      New IDL commands.
!x.type = 0  & !y.type = 0           ; 0 means linear, 1 means log.
!x.style=1   & !y.style = 1          ; Set plot boundaries in PLOT.
!p.ticklen=-0.02                      ; Small tick marks, thank you.
!p.charsize=2.0                      ; Double sized text.
;
xlow=0.20 & xhigh=0.90               ; Coords on plotter.
ylow=0.20 & yhigh=0.90
!p.clip=[xlow, ylow, xhigh, yhigh]
!p.position = [xlow, ylow, xhigh, yhigh]
!p.region=[xlow, ylow, xhigh, yhigh]
FT = '!6'  
;                                       Redefine arrays.
;
;    Colored version with proper shading.
;    512 x 512 pixel version.
!psym=3
NROW = 11409  &  MCOL = 2
;                                      Define color palate.
;  The zero coordinate gives background color.
;
I = indgen(256)
R = I*0 + 255.
;                                  Star colors:RED
R(0)=0.
R(1:50)=200.
for j=51, 75 do R(j)= 200. + (j-50)*55./25.
R(75:100)=255.
;                                  Land colors: RED
for j=101, 120 do R(j)=150. - (j-100)*7.
for j=121, 140 do R(j)=0.
for j=141, 160 do R(j)=0.
for j=161, 180 do R(j)=(j-160)*10.
for j=181, 200 do R(j)=210. + (j-180)*2.
;                                  Dim star colors.
;for j=201, 254 do R(j)=(j-200)*255./54.
for j=201, 254 do R(j)=(j-200)*200./54. + 55.
R(255)=255.
;                                  Star colors: GREEN
G=I*0+255.
G(0)=0.
for j=1, 100 do G(j)=200.+54.*(SIN(3.14*j/100.))
;                                  Land colors: GREEN
for j=101, 120 do G(j)= 0.
for j=121, 140 do G(j)= (j-120)*9.
for j=141, 160 do G(j)= 180. + (j-140)*3.7
for j=161, 180 do G(j)= 255. - (j-160)*3.7
for j=181, 200 do G(j)= 180. - (j-180)*9.
;                                  Dim star colors.
;for j=201, 254 do G(j)=(j-200)*255./54.
for j=201, 254 do G(j)=(j-200)*200./54. + 55.
G(255)=255.
;                                  Star colors: BLUE
B=I*0+255.
B(0)=0.
for j=1, 100 do B(j)=200. + 54.*(COS(3.14*j/100.))
;                                  Land colors: BLUE
for j=101, 120 do B(j)= 115. + (j-100)*7.
for j=121, 140 do B(j)= 255. 
for j=141, 160 do B(j)= 255. - (j-140)*9.
for j=161, 180 do B(j)= 70. - (j-160)*3.3
for j=181, 200 do B(j)= 0.
;                                  Dim star colors.
;for j=201, 254 do B(j)=(j-200)*255./54.
for j=201, 254 do B(j)=(j-200)*200./54. + 55.
B(255)=255.
;
TVLCT, R, G, B 
;                                      Get and read STAR.DAT.
FILENM = 'star.dat'  ; &  ATTRIB = '/SEQ/VAR/LIST'
FILESPIN = FILENM ;+ ATTRIB
GET_LUN, MLUN                         ;GET AVAILABLE LOGICAL UNIT
OPENR, MLUN, FILESPIN & REC = ''
;                                      Read in blank lines.
STR=''
title=' '
READF, MLUN, '$(a50)', STR
;                                      Read in STAR.DAT parameters.
READF, MLUN, '$(20x, f10.5)', Rs
READF, MLUN, '$(20x, f10.5)', Rad
READF, MLUN, '$(20x, f10.5)', Dobs
READF, MLUN, '$(20x, f10.5)', mcut
READF, MLUN, '$(20x, f10.5)', CthetaRot
READF, MLUN, '$(20x, f10.5)', CphiRot
READF, MLUN, '$(20x, f10.5)', UthetaRot
READF, MLUN, '$(20x, f10.5)', UphiRot
READF, MLUN, '$(20x, a10)', title

;print, Rs, Rad, Dobs, mcut, CthetaRot, CphiRot, UthetaRot, UphiRot
;
;    Generates lensed star-field image figure: outer image to Einstein ring.
;                                     Redefine symbols.
!noeras=1
ANG=(9.24*FINDGEN(40)) / !RADEG
USERSYM, COS(ANG)*1., SIN(ANG)*1., FILL=1
!P.COLOR=255
polyfill, [-2.13, -2.13, -2.1, -2.1], [-2.13, -2.1, -2.1, -2.13]
!p.psym=8
;                                      Define constants.
;                                      Get outer star file 1.
FILENM = 'bhsky1.out'  
FILESPIN = FILENM 
GET_LUN, MLUN                         ;GET AVAILABLE LOGICAL UNIT
OPENR, MLUN, FILESPIN & REC = ''
READF, MLUN, NROW
;
X=FLTARR(2)
Y=FLTARR(2)
X(0)=-2
Y(0)=-2
;                                      Set axis color by plotting dummy point.
X(1)=-2
Y(1)=-2
!P.COLOR=255
if (dev eq 'ps') then !p.color=0
plot, x, y, xrange=[xmin, xmax], yrange=[ymin, ymax]
;
IF NROW EQ 0 THEN GOTO, STAR2
READF, MLUN, AS, BS, PS, CS
X(1)=AS
Y(1)=BS
!P.COLOR=CS
oplot, x, y, symsiz=PS/4.0
;
IF NROW EQ 1 THEN GOTO, STAR2
FOR K=0, NROW-2 DO BEGIN
READF, MLUN, AS, BS, PS, CS
X(1)=AS
Y(1)=BS
!P.COLOR=CS
oplot, x, y, symsiz=PS/4.0
ENDFOR
;
STAR2:
CLOSE, MLUN & FREE_LUN, MLUN
;
;    Generates lensed star-field image figure: inner image to Einstein ring.
;                                      Get inner star file 2.
FILENM = 'bhsky2.out'  
FILESPIN = FILENM 
GET_LUN, MLUN                         ;GET AVAILABLE LOGICAL UNIT
OPENR, MLUN, FILESPIN & REC = ''
READF, MLUN, NROW
;
X=FLTARR(2)
Y=FLTARR(2)
X(0)=-2
Y(0)=-2
;
IF NROW EQ 0 THEN GOTO, THENDER
READF, MLUN, AS, BS, PS, CS
X(1)=AS
Y(1)=BS
!P.COLOR=CS
oplot, x, y, symsiz=PS/4.0
;
IF NROW EQ 1 THEN GOTO, THENDER
FOR K=0, NROW-2 DO BEGIN
READF, MLUN, AS, BS, PS, CS
X(1)=AS
Y(1)=BS
!P.COLOR=CS
oplot, x, y, symsiz=PS/4.0
ENDFOR
;
CLOSE, MLUN & FREE_LUN, MLUN

if (dev eq 'gif') then begin
title=strcompress(title, /remove_all) + '.gif'
print, ' Generating output file: ' + title
image24=TVRD(True=1)
image2D=Color_Quan(image24, 1, R, G, B, cube=6, colors=256)
write_gif, title, image2D, R, G, B 
endif

;
THENDER:
STOP
END