       SUBROUTINE dmvtpole (sp, x1, xp, yp, dxqdxp, dxqdyp)

C-----------------------------------------------------------------------C
C                      FA version 1 (fa_v1)                             C
C.......................................................................C
C  ROUTINE TO compute the partial derivative of XQ with respect to the	C
C  pole coordinates: xp and yp						C
C-----------------------------------------------------------------------C
C									C
C  INPUT:  sp: (radian)							C
C	   xp, yp: coordinates of the pole in the terrestrial		C
C		reference system (radian)				C
C	   x1(3,3): transformation matrix from the terrestrial system	C
C		to the celestial system (radian)			C
C									C
C  OUTPUT: dxqdxp(3,3), dxqdyp(3,3): matrices of the partial derivativesC
C		of the transformation matrix XQ with respect to the EOP:C
C		the pole coordinates xp, yp (without dimension)		C
C									C
C   SUBROUTINE INTERFACE:						C
C	caller subroutine: new_transf					C
C	called subroutines: rota, drota, multmat			C
C									C
C-----------------------------------------------------------------------C
C  See notice AM Gontier						C
C-----------------------------------------------------------------------C

        IMPLICIT NONE
        DOUBLE PRECISION x1, xp, yp, sp, dxqdxp, dxqdyp
        DIMENSION x1(3,3), dxqdxp(3,3), dxqdyp(3,3)

C--------------LOCAL VARIABLES------------------------------------------C

        DOUBLE PRECISION xmat1, xmat2, dmat, xmat
        DIMENSION xmat1(3,3), xmat2(3,3), dmat(3,3), xmat(3,3)


         CALL rota (1, yp, xmat1)
         CALL rota (3, -sp, xmat)
         CALL drota (2, xp, dmat)

         CALL multmat (x1, xmat, xmat2, 3, 3, 3)
         CALL multmat (xmat2, xmat1, xmat, 3, 3, 3)
         CALL multmat (xmat, dmat, dxqdxp, 3, 3, 3)

         CALL rota (2, xp, xmat1)
         CALL drota (1, yp, dmat)
         CALL multmat (xmat2, dmat, xmat, 3, 3, 3)
         CALL multmat (xmat, xmat1, dxqdyp, 3, 3, 3)

        END


       SUBROUTINE dangstellair (ak, theta, xpn, xyp, dxqdut1)

C-----------------------------------------------------------------------C
C                      FA version 1 (fa_v1)                             C
C.......................................................................C
C  ROUTINE TO compute the partial derivate of XQ with respect to UT1	C
C-----------------------------------------------------------------------C
C									C
C  INPUT:  ak: constant							C
C	   theta : (radian)						C
C	   xpn(3,3), xyp(3,3): matrices		 			C
C									C
C  OUTPUT: dxqdut1(3,3): matrix of the partial derivatives of the 	C
C		transformation matrix XQ with respect to the EOP UT1	C
C		(without dimension)					C
C									C
C  SUBROUTINE INTERFACE:						C
C	caller subroutine: new_transf					C
C	called subroutines: drota, multmat, scalmat			C
C									C
C-----------------------------------------------------------------------C
C  See notice AM Gontier	pp. 5, 11				C
C-----------------------------------------------------------------------C

        IMPLICIT NONE
        DOUBLE PRECISION xpn, xyp, theta, ak, dxqdut1
        DIMENSION xpn(3,3), xyp(3,3), dxqdut1(3,3)

C---------------LOCAL VARIABLES-----------------------------------------C

        DOUBLE PRECISION dtheta, xmat
        DIMENSION dtheta(3,3), xmat(3,3)


         CALL drota (3, -theta, dtheta)
         CALL multmat (xpn, dtheta, xmat, 3, 3, 3)
         CALL multmat (xmat, xyp, dxqdut1, 3, 3, 3)
         CALL scalmat (-ak, dxqdut1, dxqdut1, 3, 3)

        END


       SUBROUTINE dmatm1 (stheta, x, y, xtheta, xyp, dxqdx, dxqdy)

C-----------------------------------------------------------------------C
C                      FA version 1 (fa_v1)                             C
C.......................................................................C
C  ROUTINE TO compute the partial derivative of the matrix M1, directly	C
C  function of x and y							C
C-----------------------------------------------------------------------C
C									C
C			see notice pp. 6, 12				C
C									C
C  INPUT:  stetha: (arc second)						C
C	   x, y: celestial coordinates of the CEP (arc seconde)		C
C	   xtheta(3,3), xyp(3,3): matrices				C
C									C
C  OUTPUT: dxqdx(3,3), dxqdy(3,3): matrices of the partial derivatives 	C
C		of the transformation matrix XQ with respect to the EOP:C
C		the pole coordinates x, y (without dimension)		C
C									C
C  INCLUDE rad_conv							C
C	   secrad: second of angle in radian				C
C									C
C  SUBROUTINE INTERFACE:						C
C	   caller subroutine: new_transf				C
C	   called subroutines: rota, multmat				C
C									C
C-----------------------------------------------------------------------C
C  See notice AM Gontier						C
C-----------------------------------------------------------------------C

        IMPLICIT NONE
        DOUBLE PRECISION x, y, stheta, xtheta, xyp, dxqdx, dxqdy
        DIMENSION xtheta(3,3), xyp(3,3), dxqdx(3,3), dxqdy(3,3)
        INCLUDE 'rad_conv.i'

C-------------LOCAL VARIABLES-------------------------------------------C

        DOUBLE PRECISION xr, yr, x2, x3, x4, y2, xmat, xmat1, xmat2, sthetar
        DIMENSION xmat(3,3), xmat1(3,3), xmat2(3,3)

         xr = x * secrad
         yr = y * secrad
         sthetar = stheta * secrad
         x2 = xr**2
         x3 = x2 * xr
         x4 = x2**2
         y2 = yr**2

         dxqdx(1,1) = -xr - x3 / 2.d0
         dxqdx(1,2) = -yr
         dxqdx(1,3) = 1.d0
         dxqdx(2,1) = 0.d0
         dxqdx(2,2) = 0.d0
         dxqdx(2,3) = 0.d0
         dxqdx(3,1) = -1.d0 - y2 / 2.d0
         dxqdx(3,2) = yr * xr
         dxqdx(3,3) = -xr - x3 / 2.d0

         dxqdy(1,1) = -(x2 * yr) / 4.d0
         dxqdy(1,2) = -xr
         dxqdy(1,3) = 0.d0
         dxqdy(2,1) = 0.d0
         dxqdy(2,2) = -yr - (x2 * yr) / 4.d0
         dxqdy(2,3) =1.d0
         dxqdy(3,1) = -xr * yr
         dxqdy(3,2) = -1.d0 + x2 / 2.d0 + x4 / 8.d0
         dxqdy(3,3) = -yr - (x2 * yr) / 2.d0

         CALL rota (3,sthetar,xmat)
         CALL multmat (dxqdx, xmat, xmat1, 3, 3, 3)
         CALL multmat (dxqdy, xmat, xmat2, 3, 3, 3)
         CALL multmat (xtheta, xyp, xmat, 3, 3, 3)
         CALL multmat (xmat1, xmat, dxqdx, 3, 3, 3)
         CALL multmat (xmat2, xmat, dxqdy, 3, 3, 3)

        END


C  Programm rad_conv.i	C

C-----------------------------------------------------------------------C
C                      FA version 1 (fa_v1)                             C
C.......................................................................C
C  PROGRAMM is used to call the variables for the conversion in radian	C
C  secrad, degrad, hourad						C
C-----------------------------------------------------------------------C
C									C
C  secrad: to convert seconds of angle in radians			C
C  degrad: to convert degree in radians					C
C  hourad: to convert second of time in radians				C
C									C
C  INTERFACE:								C
C	   caller subroutine: arctangente, coord1_XY, dmatm1, dprecnut, C
C			      earth, eop_conv, gravit_delay, mvtpole, 	C
C			      new_transf, prec, precnut, wahr		C
C									C
C-----------------------------------------------------------------------C

        DOUBLE PRECISION secrad, degrad, hourad

        COMMON /rad_conv/ secrad, degrad, hourad


       SUBROUTINE radconv

C-----------------------------------------------------------------------C
C                      FA version 1 (fa_v1)                             C
C.......................................................................C
C  ROUTINE TO initialise variables for the conversion in radian		C
C-----------------------------------------------------------------------C
C									C
C  secrad: to convert seconds of angle in radians			C
C  degrad: to convert degree in radians					C
C  hourad: to convert second of time in radians				C
C									C
C-----------------------------------------------------------------------C

        INCLUDE 'rad_conv.i'
        INCLUDE 'initial.i'

         secrad = tpi / (360.D0 * 60.D0 * 60.D0)
         degrad = tpi / (360.D0)
         hourad = tpi / (24.D0 * 60.D0 * 60.D0)

       END


C  Programm initial.i	C


C-----------------------------------------------------------------------C
C                      FA version 1 (fa_v1)                             C
C.......................................................................C
C  PROGRAMM is used as an include to call variables: vlight, tpi, pi, 	C
C  c2									C
C-----------------------------------------------------------------------C
C									C
C  vlight: light velocity (meters/second)				C
C  tpi: 2 * pi								C
C  pi: pi								C
C  c2: (light velocity)^2						C
C									C
C  INTERFACE:								C
C	  caller subroutines: angstellaire, arctangente, cdgeom, 	C
C			      coord1_XY, dcoordbasel, dcoordsrc, deop, 	C
C			      earth, elev_az, factor, gravit_delay, 	C
C			      ionos_delay, moyer, new_transf, ocean_arg	C
C			      prec, tides, tropos_delay, wahr		C
C									C
C-----------------------------------------------------------------------C
	
        DOUBLE PRECISION vlight, tpi, pi, c2

        COMMON /initial/ vlight, tpi, pi, c2


       SUBROUTINE initiale

C-----------------------------------------------------------------------C
C                      FA version 1 (fa_v1)                             C
C.......................................................................C
C  Routine to initialise variables: vlight, tpi, pi, c2			C
C-----------------------------------------------------------------------C
C									C
C  vlight: light velocity (meters/second)				C
C  tpi: 2 * pi								C
C  pi: pi								C
C  c2: (light velocity)^2						C
C									C
C-----------------------------------------------------------------------C

        INCLUDE 'initial.i'

         vlight = 2.99792458D+08
         tpi = 6.28318530717958648D0
         pi = tpi/2.D0
         c2 = vlight**2

       END
