gcvgd.f 3.86 KB
      subroutine gcvgd (rkm,gclat,alt,gdlat,ier)
c*
c***********************************************************************
c*
c*          "Copyright [c] CNES 98 - tous droits reserves"
c*          **********************************************
c*
c*PRO MAGLIB
c*
c*VER 99.03.31 - V 1.0
c*VER 01.06.05 - V 2.0
c*VER 03.01.06 - V 2.1
c*
c*AUT adap. CNES - JC KOSIK - janvier 1991
c*AUT port. CISI
c*
c*ROL Theme : Changements de coordonnees
c*ROL         Transformation des composantes geocentriques
c*ROL         en composantes geodesiques.
c*
c*PAR rkm   (I) : distance geocentrique (kilometres)
c*PAR gclat (I) : latitude geocentrique (radians)
c*
c*PAR alt   (O) : altitude (kilometres)
c*PAR gdlat (O) : latitude geodesique (radians)
c*
c*PAR ier   (O) : code de retour
c*
c*NOT ier       : sans objet
c*
c*NOT common    : util
c*
c*INF utilise   : sans objet
c*
c*HST version 1.0 - 99.03.31 - creation de la maglib au CDPP
c*HST version 2.0 - 01.06.05 - correction de commentaires de code
c*HST version 2.1 - 03.01.06 - corrections en compilation avec g77
c
c***********************************************************************
c*
      implicit none
c
c     ---------------------------------
c*FON Declaration identificateur rcs_id
c     ---------------------------------
c
      character rcs_id*100
c
c     --------------------------
c*FON Declaration des parametres
c     --------------------------
c
      double precision rkm
      double precision gclat
      double precision alt, gdlat
      integer ier
c
c     ----------------------------------
c*FON Declaration des variables communes
c     ----------------------------------
c
      double precision pi,dpi,rad,deg,pid,xmu,rayt
c
c*COM pi   : constante pi (obtenue a partir de acos(-1.))
c*COM dpi  : constante 2 * pi
c*COM pid  : constante pi / 2
c*COM rad  : facteur de conversion degres  ----> radians
c*COM deg  : facteur de conversion radians ----> degres
c*COM xmu  : constante de gravitation terrestre (km**3/sec**2)
c*COM rayt : rayon equatorial terrestre (km)
c
      common/util/pi,dpi,rad,deg,pid,xmu,rayt
c
c     ---------------------------------
c*FON Declaration des variables locales
c     ---------------------------------
c
      double precision ep2
      double precision a2,a4,a6,a8,rer
      double precision scl,ccl,s2cl,c2cl,s4cl,c4cl,s8cl,s6cl,dltcl
c*LOC Variables de travail intermediaires
c
      SAVE
c
c     ---------------------------------
c*FON Affectation identificateur rcs_id
c     ---------------------------------
c
      data rcs_id /"
     >$Id$"/
c
c     --------------------------
c*FON Affectation des constantes
c     --------------------------
c
      data ep2 /.0067397d0/
c
c     ******************
c     Debut de programme
c     ******************
c
c     --------------------------------------------------
c*FON Calcul des coefficients de la transformation
c*FON voir Astron. J. vol66 p15, 1961, pour les formules
c     --------------------------------------------------
c
      ier = 0
c
      rer = rkm / rayt
      a2  = ((1.4127348d-8 / rer + .94339131d-8) / rer + .33523288d-2)
     >      / rer
      a4  = (((-1.2545063d-10 / rer + .11760996d-9) /
     >      rer + .11238084d-4) / rer  -.2814244d-5) / rer
      a6  = ((54.939685d-9 / rer - 28.301730d-9) / rer + 3.5435979d-9)
     >      / rer
      a8  = (((320.d0 / rer - 252.d0) / rer + 64.d0) / rer - 5.d0) 
     >      / rer * .98008304d-12
c
      scl   = sin(gclat)
      ccl   = sqrt(1.d0 - scl * scl)
c
      s2cl  = 2.d0 * scl * ccl
      c2cl  = 2.d0 * ccl * ccl - 1.d0
c
      s4cl  = 2.d0 * s2cl * c2cl
      c4cl  = 2.d0 * c2cl * c2cl - 1.d0
c
      s8cl  = 2.d0 * s4cl * c4cl
      s6cl  = s2cl * c4cl + c2cl * s4cl
c
      dltcl = s2cl * a2 + s4cl * a4 + s6cl * a6 + s8cl * a8
c
      gdlat = dltcl + gclat
      alt   = rkm - rayt  /sqrt(1.d0 + ep2 * scl * scl)
c
c     ****************
c     Fin de programme
c     ****************
c
      return
      end