;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;; PROCEDURE SOLARS    ;;;
;;; By Darrin Casebeer  ;;;
;;;    Kevin Tubbs      ;;;
;;;    Summer 98        ;;;
;;;;;;;;;;;;;;;;;;;;;;;;;;;

PRO SOLARS

; First delete windows that may be up
; Need to make this more robust in future

WDELETE
WDELETE

WINDOW,0,XSIZE=900,YSIZE=900,TITLE="SOLAR SYSTEM"

G=(6.67d)*10d^(-11d)

Me=(5.976d)*10d^(24d)
Ms=(333480d)*Me

SET_PLOT,'X'
maxi=1000d
minmaxrange=6d*10d^12d
  SURFACE, dist(3,1200),/NODATA,/SAVE,XRANGE=[-MinMaxRange,MinMaxRange], $
  XSTYLE=1,XTITLE='X',YRANGE=[-MinMaxRange,MinMaxRange], YSTYLE=1,$
  YTITLE='Y',CHARSIZE=.5, ZRANGE=[-MinMaxRange,MinMaxRange],$
  ZSTYLE=1,ZTITLE='Z',TITLE='OUTER PLANETS"
CLOSE,10
OPENR,10,"planet.data"
NAME=''
READF,10,NUMPLANETS
FOR J=0,NUMPLANETS-1 DO BEGIN
  READF,10, a,ecc,Eclang,M,NAME
  Eclang=Eclang*2*!Pi/360
  M=M*Me

  a=a*(1.495d)*10d^(11d)
  Rmin=a*((1.0d)-ecc)
  ;PRINT,Rmin,a,M,Ms,G
  Voy=-((G*(M+Ms)*((1.0d)+ecc))/(a*((1.0d)-ecc)))^(0.5d)
  P=((4.0d*!Pi^(2.0d)*a^(3.0d))/(G*(M+Ms)))^(0.5d)
  if STRPOS(NAME,'Halley') GT 0 THEN maxi=5000L ELSE maxi=1000
  dt=P/(maxi)
  PRINT,NAME,"  Period",P," Seconds Intial velocity at Rmin ",Voy," m/s "
  PRINT,"number of steps ",LONG(P/dt)
;  IF J EQ 6 THEN BEGIN
;    WINDOW,2
;    minmaxrange=3d*10d^11d
;    SURFACE, dist(3,1200),/NODATA,/SAVE,XRANGE=[-MinMaxRange,MinMaxRange], $
;    XSTYLE=1,XTITLE='X',YRANGE=[-MinMaxRange,MinMaxRange], YSTYLE=1,$
;    YTITLE='Y',CHARSIZE=.5, ZRANGE=[-MinMaxRange,MinMaxRange],$
;    ZSTYLE=1,ZTITLE='Z',TITLE='INNER PLANETS'
;  ENDIF

 
  Vox=0.0d
  Xo=Rmin
  Yo=0.0d
  X=DBLARR(LONG(maxi))
  Y=DBLARR(LONG(maxi))
  Z=DBLARR(LONG(maxi))
  U=DBLARR(4)
  t=0.0d
  X(0)=Xo
  Y(0)=Yo*COS(Eclang)
  Z(0)=0
  U(0)=Xo
  U(1)=Yo
  U(2)=0.0d
  U(3)=VoY
  FOR I =1L,maxi-1 DO BEGIN
    MVTWOBODY,G,M,Ms,U,t,dt
    X(I)=U(0)*Ms/(M+Ms)
    Y(I)=U(1)*Ms/(M+Ms)
    X(I)=X(I)*COS(Eclang)
    Z(I)=X(I)*(-SIN(Eclang))
    t=t+dt
  PLOTS,X(I),Y(I),Z(I),PSYM=2,COL=92,SYMSIZE=0.1,/T3D
  ENDFOR
ENDFOR
;WRITE_GIF,'solarsystem.gif',TVRD()
PLOTS,0,0,0,PSYM=2,COL=173,SYMSIZE=1.5,/T3D
CLOSE,10
RETURN
END