subroutine xyobs1980(mjd,dpsi,deps,X,Y)

C-------------------------------------------------------------------------C
C									  C
C       Conversion of observed dpsi, deps to X, Y                         C
C									  C
C-------------------------------------------------------------------------C
C	 								  C
C	References :                                                      C
C		Lieske et al. 1977, A&A 58, pp. 1-16                      C
C	 	Capitaine et al. 2000, A&A 355, pp. 398-405               C
C	 	Capitaine et al. 2003, A&A 400, pp. 1145-1154             C
C               IERS Conventions 2000, chap. 5                            C
C	 								  C
C-------------------------------------------------------------------------C
C	 								  C
C	Input :                                                           C
C		mjd : time in mjd                                         C
C		dpsi, deps : celestial pole offsets in arcsec             C
C                            refered to the IAU1980 nutation              C
C	Output :                                                          C
C               X, Y : CIP coordinates in the ICRF in arcsec	          C
C	 								  C
C-------------------------------------------------------------------------C
C	 								  C
C	Subroutine interface : iau_NUT80 (SOFA routine)		 	  C
C	 								  C
C-------------------------------------------------------------------------C
C	 								  C
C	Author : Sebastien Lambert (Sebastien.Lambert@obspm.fr)		  C
C	Last modified : 22/11/02	 				  C
C	 								  C
C-------------------------------------------------------------------------C

	implicit none

C	Input variables

	double precision               :: mjd,dpsi,deps

C 	Output variables

	double precision               :: X,Y
	
C	Other variables
	
	double precision               :: depsobs,dpsiobs
	double precision               :: deps80,dpsi80,deps1,dpsi1
	double precision               :: psiA,omeA,epsA,chiA
	double precision               :: psi,omega
	double precision               :: Xb,Yb
	double precision               :: jc,eps0
        double precision, parameter    :: a2r=4.84813681109535993D-6
        double precision, parameter    :: r2a=2.06264806247096355D+5

C	Julian centuries from J2000.0

	jc=(mjd-51544.5D0)/36525.D0

C	Lieske et al. (1977) expressions in arcsec and converted in radians

	eps0=84381.448D0

	psiA=5038.7784D0*jc
     .      -1.07259D0*jc*jc
     .      -0.001147D0*jc*jc*jc

        omeA=eps0
     .      +0.05127D0*jc*jc
     .      -0.007726D0*jc*jc*jc

        epsA=eps0
     .      -46.8150D0*jc
     .      -0.00059D0*jc*jc
     .      +0.001813D0*jc*jc*jc

        chiA=10.5526D0*jc
     .      -2.38064D0*jc*jc
     .      -0.001125D0*jc*jc*jc

     	eps0=eps0*a2r
	psiA=psiA*a2r
	omeA=omeA*a2r
	epsA=epsA*a2r
	chiA=chiA*a2r

C	Input values of celestial pole offsets in arcsec

	dpsiobs=dpsi
	depsobs=deps

C	IAU1980 model for the date (output in radians)

	call iau_NUT80(2400000.5D0,mjd,dpsi80,deps80)

C	Observed dpsi, deps in radians

	dpsiobs=dpsiobs*a2r
	depsobs=depsobs*a2r

	dpsiobs=dpsi80+dpsiobs
	depsobs=deps80+depsobs

C	Chapter 5, equation (22)

	dpsi1=(dpsiobs*dsin(epsA)*dcos(chiA)-depsobs*dsin(chiA))/dsin(omeA)
	deps1=dpsiobs*dsin(epsA)*dsin(chiA)+depsobs*dcos(chiA)

C	Chapter 5, equation (21)

	omega=omeA+deps1
	psi=psiA+dpsi1

C	Chapter 5, equation (20)

	Xb=dsin(omega)*dsin(psi)
	Yb=-dsin(eps0)*dcos(omega)+dcos(eps0)*dsin(omega)*dcos(psi)

C	Chapter 5, equation (18)

	X=Xb
	Y=Yb

C	Conversion radians in arcsec

	X=X*r2a
	Y=Y*r2a

	end