program smei_pos_tst

	real*8		rx,ry,rz
	integer		tt(2)
	real*8		doy
	character	cstr*80

	print *, 'testing'
	iyr = 2004
	doy = 100.0d0

	call Time2Day8(1,tt,doy)
	call Time2YDoy(1,tt,iyr,tt)

	i = Time2Str('YDOY',tt,cstr)
	print *, cstr(:itrim(cstr))
	call smei_pos(tt, rx,ry,rz)
	print *, rx,ry,rz

	end
	

	subroutine smei_pos(tt,rx,ry,rz)

	implicit double precision (A-H,O-Z)

	    integer	tt(2)
	    real*8	rx
	    real*8	ry
	    real*8	rz

	character	cSay*8		/'smei_pos'/

	include		'openfile.h'
	include		'dirspec.h'
	include		'filparts.h'
	!include	'planet.h'
	!include	'math.h'

	!double precision DE2RA		/0.174532925d-1/
	!double precision E6A		/	   1d-6/
	!double precision PI		/ 3.14159265d0 /
	!double precision PIO2		/ 1.57079633d0 /
	!double precision QO		/	 120d0 /
	!double precision SO		/	  78d0 /
	!double precision TOTHRD		/ 0.66666667d0 /
	!double precision TWOPI		/  6.2831853d0 /
	!double precision X3PIO2		/ 4.71238898d0 /
	!double precision XJ2		/   1.082616d-3/
	!double precision XJ3		/  -0.253881d-5/
	!double precision XJ4		/   -1.65597d-6/
	!double precision XKE		/0.743669161d-1/
	!double precision XKMPER		/   6378.135d0 /
	!double precision XMNPDA		/	1440d0 /
	!double precision AE		/	   1d0 /

	common /E1/ XMO,XNODEO,OMEGAO,EO,XINCL,XNO,XNDT2O,XNDD6O,BSTAR,X,Y,Z,XDOT,YDOT,ZDOT,EPOCH,DS50
	common /C1/ CK2,CK4,E6A,QOMS2T,S,TOTHRD,XJ3,XKE,XKMPER,XMNPDA,AE
	common /C2/ DE2RA,PI,PIO2,TWOPI,X3PIO2

	DATA DE2RA,E6A,PI,PIO2,QO,SO,TOTHRD,TWOPI,X3PIO2,XJ2,XJ3,XJ4,XKE,XKMPER,XMNPDA,AE
     &	    /.174532925d-1,1d-6,3.14159265d0,1.57079633d0,120d0,78d0,0.66666667d0,6.2831853d0,
     &	     4.71238898d0,1.082616d-3,-0.253881d-5,-1.65597d-6,0.743669161d-1,6378.135d0,1440d0,1d0/

	logical		bOpenFile
	integer		Time2Differ
	character	cStr*(FIL__LENGTH)
	integer		uu(2)

	parameter	(NBUF=2000)

	character	cBuf(2,NBUF)*69
	integer		tBuf(2,NBUF)
	integer		iBuf		/0/
	integer		iFlag		/1/
	integer		iBufLast	/1/
	save		iFlag, iBuf, iBufLast, tBuf, cBuf

	if (iBuf .eq. 0) then

	    CK2	  =  0.5d0*XJ2*AE**2d0
	    CK4	  = -0.375d0*XJ4*AE**4d0
	    QOMS2T= ((QO-SO)*AE/XKMPER)**4d0
	    S	  = AE*(1d0+SO/XKMPER)

	    i = iFilePath(cEnvi(:iEnvi)//'EPHEM',1,'sgp','sat27640.txt',cStr)

	    i = 0
	    if (bOpenFile(OPN__TEXT+OPN__REOPEN+OPN__STOP,iU,cStr,i)) then

		read (iU,'(A)') cStr,cStr,cStr,cStr,cStr

		read (iU,'(A)') cStr

		do while (cStr(:13) .ne. '<End of file>')

		    if (iBuf .ge. NBUF) call Say(cSay,'E','NBUF','parameter too small')

		    read (cStr(19:20),'(I2.2)' ) iYr
		    read (cStr(21:32),'(D12.8)') Doy 

		    call Time2Day8(1,uu,Doy)
		    call Time2YDoy(1,uu,2000+iYr,uu)

		    if (iBuf .eq. 1) then

			iBuf = iBuf+1
			tBuf(1,iBuf) = uu(1)
			tBuf(2,iBuf) = uu(2)
			cBuf(1,iBuf) = cStr
			read (iU,'(A)') cBuf(2,iBuf)

		    else if (Time2Differ(uu,tBuf(1,iBuf-1)) .gt. 0) then

			iBuf = iBuf+1
			tBuf(1,iBuf) = uu(1)
			tBuf(2,iBuf) = uu(2)
			cBuf(1,iBuf) = cStr
			read (iU,'(A)') cBuf(2,iBuf)

		    else if (Time2Differ(uu,tBuf(1,iBuf-1)) .lt. 0) then

			write (*,'(1X,A)') cStr(:itrim(cStr))
			call Say(cSay,'E','tle','not chronologically ordered')

		    else

			write (*,'(1X,A)') cStr(:itrim(cStr))
			read (iU,'(A)') cStr
			write (*,'(1X,A)') cStr(:itrim(cStr))
			call Say(cSay,'W','tle','skipping duplicate entry')

		    end if

		    read (iU,'(A)') cStr

		end do

	    end if

	end if

	i = iBufLast
	do while (i .le. iBuf .and. Time2Differ(tBuf(1,i),tt) .lt. 0)
	    i = i+1
	end do

	if (i .eq. 1) then

	    call Time2Delta(tt,tBuf(1,i),uu)
	    call Time2Day8(0,uu,dt1)		! dt1 <= 0 or iBuf=1

	else if (i .gt. iBuf) then

	    i = iBuf
	    call Time2Delta(tt,tBuf(1,i),uu)
	    call Time2Day8(0,uu,dt1)		! dt1 > 0

	else					! tBuf(1,i) >= tt

	    call Time2Delta(tt,tBuf(1,i  ),uu)
	    call Time2Day8(0,uu,dt1)		! dt1 <= 0

	    call Time2Delta(tt,tBuf(1,i-1),uu)
	    call Time2Day8(0,uu,dt2)		! dt2 > 0

	    if (dt2 .lt. -dt1) then
		i = i-1
		dt1 = dt2
	    end if

	end if

	iBufLast = i

	print *, cBuf(1,i)(:itrim(cBuf(1,i)))
	print *, cBuf(2,i)(:itrim(cBuf(2,i)))

	read (cBuf(1,i),'(18X,D14.8,1X,F10.8,2(1X,F6.5,I2))') EPOCH,XNDT2O,XNDD6O,IEXP,BSTAR,IBEXP
	read (cBuf(2,i),'(7X,2(1X,F8.4),1X,F7.7,2(1X,F8.4),1X,F11.8)') XINCL,XNODEO,EO,OMEGAO,XMO,XNO

	print *, EPOCH,XNDT2O,XNDD6O,IEXP,BSTAR,IBEXP
	print *, XINCL,XNODEO,EO,OMEGAO,XMO,XNO

	XNDD6O	= XNDD6O*(10d0**IEXP)
	XNODEO	= XNODEO*DE2RA
	OMEGAO	= OMEGAO*DE2RA
	XMO	= XMO*DE2RA
	XINCL	= XINCL*DE2RA
	TEMP	= TWOPI/XMNPDA/XMNPDA
	XNO	= XNO*TEMP*XMNPDA
	XNDT2O	= XNDT2O*TEMP
	XNDD6O	= XNDD6O*TEMP/XMNPDA

	! INPUT CHECK FOR PERIOD VS EPHEMERIS SELECTED
	! PERIOD GE 225 MINUTES IS DEEP SPACE

	A1	= (XKE/XNO)**TOTHRD
	TEMP	= 1.5d0*CK2*(3d0*COS(XINCL)**2d0-1d0)/(1d0-EO*EO)**1.5d0
	DEL1	= TEMP/(A1*A1)
	AO	= A1*(1d0-DEL1*(0.5d0*TOTHRD+DEL1*(1d0+134d0/81d0*DEL1)))
	DELO	= TEMP/(AO*AO)
	XNODP	= XNO/(1d0+DELO)
	BSTAR	= BSTAR*(10d0**IBEXP)/AE

	print *, iFlag, dt1
	call SGP4D(iFlag,dt1*1440d0)

	print *, X,Y,Z
	rx = X*XKMPER/AE
	ry = Y*XKMPER/AE
	rz = Z*XKMPER/AE

	return
	end