C+
C NAME:
C	Time2PrecesssionLtd
C PURPOSE:
C	Converts ecliptic coordinates for epoch JEpoch1 to epoch JEpoch2,
C	i.e. corrects for precession between the two epochs.
C	Uses a first order approximation and gives good results only for
C	time differences of several decades or less.
C CATEGORY:
C	Celestial mechanics
C CALLING SEQUENCE:
	subroutine Time2PrecessionLtd(tt1,tt2,lambda,beta)
C CALLS:
C	COSD,SIND,TAND, Time2JEpoch
C INPUTS:
C	tt1		integer			initial epoch
C	tt2		integer			final epoch
C	lambda		double precision	ecliptic longitude for tt1 (deg)
C	beta		double precision	ecliptic latitude  for tt1 (deg)
C OUTPUTS:
C	lambda		double precision	ecliptic longitude for tt2 (deg)
C	beta		double precision	ecliptic latitude  for tt2 (deg)
C PROCEDURE:
C	See R. Green, Spherical Astronomy, Cambridge UP, 1985, Sect. 9.7, p 222
C	A,B,C	precession constants
C MODIFICATION HISTORY:
C-
	    integer		tt1(2)
	    integer		tt2(2)
	    double precision	lambda
	    double precision	beta

	double precision	T
	double precision	A
	double precision	B
	double precision	C
	double precision	D

	double precision	dsind
	double precision	dcosd
	double precision	dtand

	if (abs(beta) .gt. 90.0d0) stop 'PRECESSION_APPROX: invalid BETA'

	call Time2JEpoch(2,tt1,T)
	T = 0.01d0*T						! Time interval in Julian centuries
	A = (1.396971d0 + 0.000309d0*T)*T			! Form. (9.38), p. 224
	B = (0.013056d0 - 0.000010d0*T)*T
	C = 5.1236d0 - 1.1553d0*T 
	if (abs(beta) .eq. 90.0d0) then
	    lambda = 174.8764d0+0.9137d0*T+dsign(1.0d0,beta)*90.0d0 !Form. (9.11a), p. 214
	    D    = (0.47d0-0.0007d0*T)/36.0d0*T			!Form. (9.11b), p. 214
	    beta = beta - D*dsign(1.0d0,beta)
	else
	    D    = lambda - A+B*dtand(beta)*dcosd(lambda+C)
	    beta = beta   - B*dsind(lambda+C)			!JEpoch1 ---> J2000.0
	    lambda = D
	end if

	call Time2JEpoch(0,tt2,T)
	T = 0.01d0*T						!Time interval in Julian centuries
	A = (1.396971d0 + 0.000309d0*T)*T			!Form. (9.38), p. 224
	B = (0.013056d0 - 0.000010d0*T)*T
	C =  5.1236d0   + 0.2416d0*T
	if (abs(beta) .eq. 90.0d0) then
	    lambda = 174.8764d0+0.9137d0*T-dsign(1.0d0,beta)*90.0d0 !Form. (9.11a), p. 214
	    D    = (0.47d0-0.0007d0*T)/36.0d0*T			!Form. (9.11b), p. 214
	    beta = beta+D*dsign(1.0d0,beta)
	else
	    D    = lambda + A-B*dtand(beta)*dcosd(lambda+C)
	    beta = beta   + B*dsind(lambda+C)			!J2000.0 ---> JEpoch2
	    lambda = D
	end if

	lambda = dmod(dmod(lambda,360.0d0)+360.0d0,360.0d0)

	return
	end