Index: dynamico_lmdz/simple_physics/phyparam/param/iniorbit.F
===================================================================
--- dynamico_lmdz/simple_physics/phyparam/param/iniorbit.F	(revision 4189)
+++ dynamico_lmdz/simple_physics/phyparam/param/iniorbit.F	(revision 4190)
@@ -1,4 +1,5 @@
       SUBROUTINE iniorbit
      $     (paphelie,pperiheli,pyear_day,pperi_day,pobliq)
+      USE astronomy
       IMPLICIT NONE
 
@@ -32,8 +33,4 @@
 
 c-----------------------------------------------------------------------
-c   Declarations:
-c   -------------
-
-#include "comorbit.h"
 
 c   Arguments:
@@ -56,9 +53,8 @@
       peri_day=pperi_day
 
-      pi=2.*asin(1.)
       PRINT*,'Perihelie en Mkm  ',periheli
       PRINT*,'Aphelise  en Mkm  ',aphelie 
       PRINT*,'obliquite en degres  :',obliquit
-      unitastr=149.597927
+
       e_elips=(aphelie-periheli)/(periheli+aphelie)
       p_elips=0.5*(periheli+aphelie)*(1-e_elips*e_elips)/unitastr
Index: dynamico_lmdz/simple_physics/phyparam/param/iniphyparam.F
===================================================================
--- dynamico_lmdz/simple_physics/phyparam/param/iniphyparam.F	(revision 4189)
+++ dynamico_lmdz/simple_physics/phyparam/param/iniphyparam.F	(revision 4190)
@@ -11,7 +11,7 @@
       use comsaison
       USE geometry_mod, ONLY : longitude,latitude,cell_area
-      USE phys_const, ONLY : pi,rad,g,r,cpp,rcp,dtphys,unjours,mugaz
+      USE phys_const, ONLY : rad,g,r,cpp,rcp,dtphys,unjours,mugaz
       USE planet, ONLY : coefir, coefvis
-      USE planet, ONLY : aphelie, periheli, peri_day, obliquit, year_day
+      USE astronomy
       USE vdif_mod, ONLY : lmixmin, emin_turb 
       IMPLICIT NONE
@@ -197,5 +197,4 @@
       coslon(:)=cos(long(:))
 
-      pi=2.*asin(1.)
       prad=rad
       pg=g
Index: dynamico_lmdz/simple_physics/phyparam/param/orbite.F
===================================================================
--- dynamico_lmdz/simple_physics/phyparam/param/orbite.F	(revision 4189)
+++ dynamico_lmdz/simple_physics/phyparam/param/orbite.F	(revision 4190)
@@ -1,3 +1,4 @@
       SUBROUTINE orbite(pls,pdist_sol,pdecli)
+      USE astronomy
       IMPLICIT NONE
 
@@ -32,6 +33,4 @@
 c   -------------
 
-#include "comorbit.h"
-
 c arguments:
 c ----------
Index: dynamico_lmdz/simple_physics/phyparam/param/paramdef.F
===================================================================
--- dynamico_lmdz/simple_physics/phyparam/param/paramdef.F	(revision 4189)
+++ dynamico_lmdz/simple_physics/phyparam/param/paramdef.F	(revision 4190)
@@ -3,5 +3,5 @@
       USE vdif_mod, ONLY : lmixmin, emin_turb
       USE planet, ONLY : coefir, coefvis
-      USE planet, ONLY : aphelie, periheli, obliquit, peri_day, year_day
+      USE astronomy
       IMPLICIT NONE
 c-----------------------------------------------------------------------
@@ -20,5 +20,5 @@
       character str1*1
       INTEGER ig,ierr
-      real zz,z1,z2,pi
+      real zz,z1,z2
 
       REAL I_mer,I_ter,z0_mer,z0_ter
@@ -29,5 +29,4 @@
 c-----------------------------------------------------------------------
 
-      pi=2.*asin(1.)
       year_day=360.
       periheli=173.
Index: dynamico_lmdz/simple_physics/phyparam/param/phyparam.F
===================================================================
--- dynamico_lmdz/simple_physics/phyparam/param/phyparam.F	(revision 4189)
+++ dynamico_lmdz/simple_physics/phyparam/param/phyparam.F	(revision 4190)
@@ -12,6 +12,7 @@
       USE phys_const
       USE planet
+      USE astronomy
       USE vdif_mod, ONLY : vdif
-c
+c     
       IMPLICIT NONE
 c=======================================================================
Index: dynamico_lmdz/simple_physics/phyparam/param/solarlong.F
===================================================================
--- dynamico_lmdz/simple_physics/phyparam/param/solarlong.F	(revision 4189)
+++ 	(revision )
@@ -1,100 +1,0 @@
-      SUBROUTINE solarlong(pday,psollong)
-      IMPLICIT NONE
-
-c=======================================================================
-c
-c   Objet:
-c   ------
-c
-c      Calcul de la distance soleil-planete et de la declinaison
-c   en fonction du jour de l'annee.
-c
-c
-c   Methode:
-c   --------
-c
-c      Calcul complet de l'elipse
-c
-c   Interface:
-c   ----------
-c
-c      Uncommon comprenant les parametres orbitaux.
-c
-c   Arguments:
-c   ----------
-c
-c   Input:
-c   ------
-c   pday          jour de l'annee (le jour 0 correspondant a l'equinoxe)
-c   lwrite        clef logique pour sorties de controle
-c
-c   Output:
-c   -------
-c   pdist_sol     distance entre le soleil et la planete
-c                 ( en unite astronomique pour utiliser la constante 
-c                  solaire terrestre 1370 Wm-2 )
-c   pdecli        declinaison ( en radians )
-c
-c=======================================================================
-c-----------------------------------------------------------------------
-c   Declarations:
-c   -------------
-
-#include "comorbit.h"
-
-c arguments:
-c ----------
-
-      REAL pday,pdist_sol,pdecli,psollong
-      LOGICAL lwrite
-
-c Local:
-c ------
-
-      REAL zanom,xref,zx0,zdx,zteta,zz
-      INTEGER iter
-
-
-c-----------------------------------------------------------------------
-c calcul de l'angle polaire et de la distance au soleil :
-c -------------------------------------------------------
-
-c  calcul de l'zanomalie moyenne
-
-      zz=(pday-peri_day)/year_day
-      zanom=2.*pi*(zz-nint(zz))
-      xref=abs(zanom)
-
-c  resolution de l'equation horaire  zx0 - e * sin (zx0) = xref
-c  methode de Newton
-
-      zx0=xref+e_elips*sin(xref)
-      DO 110 iter=1,10
-         zdx=-(zx0-e_elips*sin(zx0)-xref)/(1.-e_elips*cos(zx0))
-         if(abs(zdx).le.(1.e-7)) goto 120
-         zx0=zx0+zdx
-110   continue
-120   continue
-      zx0=zx0+zdx
-      if(zanom.lt.0.) zx0=-zx0
-
-c zteta est la longitude solaire
-
-      zteta=2.*atan(sqrt((1.+e_elips)/(1.-e_elips))*tan(zx0/2.))
-
-      psollong=zteta-timeperi
-
-      IF(psollong.LT.0.) psollong=psollong+2.*pi
-      IF(psollong.GT.2.*pi) psollong=psollong-2.*pi
-c-----------------------------------------------------------------------
-c   sorties eventuelles:
-c   ---------------------
-
-      IF (lwrite) THEN
-         PRINT*,'jour de l"annee   :',pday
-         PRINT*,'distance au soleil (en unite astronomique) :',pdist_sol
-         PRINT*,'declinaison (en degres) :',pdecli*180./pi
-      ENDIF
-
-      RETURN
-      END
Index: dynamico_lmdz/simple_physics/phyparam/physics/astronomy.F90
===================================================================
--- dynamico_lmdz/simple_physics/phyparam/physics/astronomy.F90	(revision 4190)
+++ dynamico_lmdz/simple_physics/phyparam/physics/astronomy.F90	(revision 4190)
@@ -0,0 +1,87 @@
+MODULE astronomy
+  IMPLICIT NONE
+  SAVE
+
+  REAL :: aphelie, periheli, year_day, peri_day, obliquit, &
+       timeperi, e_elips,p_elips
+  REAL, PARAMETER :: unitastr=149.597927, & ! millions of km
+       pi=2.*ASIN(1.)
+CONTAINS
+ 
+  SUBROUTINE solarlong(pday,psollong)
+
+!=======================================================================
+!      Calcul de la distance soleil-planete et de la declinaison
+!   en fonction du jour de l'annee.
+!
+!   Methode:
+!   --------
+!      Calcul complet de l'ellipse
+!
+!   Input:
+!   ------
+!   pday          jour de l'annee (le jour 0 correspondant a l'equinoxe)
+!   lwrite        clef logique pour sorties de controle
+!
+!   Output:
+!   -------
+!   pdist_sol     distance entre le soleil et la planete
+!                 ( en unite astronomique pour utiliser la constante 
+!                  solaire terrestre 1370 Wm-2 )
+!   pdecli        declinaison ( en radians )
+!
+!=======================================================================
+
+    USE planet
+    REAL, INTENT(IN) :: pday
+    REAL, INTENT(OUT) :: psollong
+    LOGICAL, PARAMETER ::  lwrite=.TRUE.
+    REAL pdist_sol, pdecli
+
+! Local:
+! ------
+    REAL zanom,xref,zx0,zdx,zteta,zz
+    INTEGER iter
+
+!--------------------------------------------------------
+! calcul de l'angle polaire et de la distance au soleil :
+! -------------------------------------------------------
+
+!  calcul de l'zanomalie moyenne
+
+    zz=(pday-peri_day)/year_day
+    zanom=2.*pi*(zz-nint(zz))
+    xref=abs(zanom)
+
+!  resolution de l'equation horaire  zx0 - e * sin (zx0) = xref
+!  methode de Newton
+    
+    zx0=xref+e_elips*sin(xref)
+    DO iter=1,10
+       zdx=-(zx0-e_elips*sin(zx0)-xref)/(1.-e_elips*cos(zx0))
+       zx0=zx0+zdx
+    END DO
+    zx0=zx0+zdx
+    if(zanom.lt.0.) zx0=-zx0
+    
+    ! zteta est la longitude solaire
+    
+    zteta=2.*atan(sqrt((1.+e_elips)/(1.-e_elips))*tan(zx0/2.))
+    
+    psollong=zteta-timeperi
+    
+    IF(psollong.LT.0.) psollong=psollong+2.*pi
+    IF(psollong.GT.2.*pi) psollong=psollong-2.*pi
+    !-----------------------------------------------------------------------
+    !   sorties eventuelles:
+    !   ---------------------
+    
+    IF (lwrite) THEN
+       PRINT*,'jour de l"annee   :',pday
+       PRINT*,'distance au soleil (en unite astronomique) :',pdist_sol
+       PRINT*,'declinaison (en degres) :',pdecli*180./pi
+    ENDIF
+    
+  END SUBROUTINE solarlong
+  
+END MODULE astronomy
Index: dynamico_lmdz/simple_physics/phyparam/physics/phys_const.F90
===================================================================
--- dynamico_lmdz/simple_physics/phyparam/physics/phys_const.F90	(revision 4189)
+++ dynamico_lmdz/simple_physics/phyparam/physics/phys_const.F90	(revision 4190)
@@ -1,5 +1,5 @@
 MODULE phys_const
   
-  REAL pi,rad,g,r,cpp,rcp,dtphys,unjours,mugaz,omeg
+  REAL rad,g,r,cpp,rcp,dtphys,unjours,mugaz,omeg
 
 END MODULE phys_const
Index: dynamico_lmdz/simple_physics/phyparam/physics/planet.F90
===================================================================
--- dynamico_lmdz/simple_physics/phyparam/physics/planet.F90	(revision 4189)
+++ dynamico_lmdz/simple_physics/phyparam/physics/planet.F90	(revision 4190)
@@ -3,6 +3,5 @@
   SAVE
 
-  REAL :: aphelie, periheli, year_day, peri_day, obliquit, &
-       coefvis, coefir
+  REAL :: coefvis, coefir
   
 END MODULE planet
