Index: /LMDZ4/branches/LMDZ4-dev/libf/dyn3d/calfis.F
===================================================================
--- /LMDZ4/branches/LMDZ4-dev/libf/dyn3d/calfis.F	(revision 1200)
+++ /LMDZ4/branches/LMDZ4-dev/libf/dyn3d/calfis.F	(revision 1201)
@@ -1,10 +1,9 @@
 !
-! $Header$
+! $Id$
 !
 C
 C
       SUBROUTINE calfis(lafin,
-     $                  rdayvrai,
-     $                  heure,
+     $                  jD_cur, jH_cur,
      $                  pucov,
      $                  pvcov,
@@ -102,5 +101,5 @@
 c    -----------
       LOGICAL  lafin
-      REAL heure
+
 
       REAL pvcov(iip1,jjm,llm)
@@ -170,5 +169,6 @@
       DATA firstcal/.true./
       SAVE firstcal,debut
-      REAL rdayvrai
+!      REAL rdayvrai
+      REAL :: jD_cur, jH_cur
 c
 c-----------------------------------------------------------------------
@@ -445,6 +445,6 @@
      .             debut,
      .             lafin,
-     .             rdayvrai,
-     .             heure,
+     .             jD_cur,
+     .             jH_cur,
      .             dtphys,
      .             zplev,
Index: /LMDZ4/branches/LMDZ4-dev/libf/dyn3d/conf_gcm.F
===================================================================
--- /LMDZ4/branches/LMDZ4-dev/libf/dyn3d/conf_gcm.F	(revision 1200)
+++ /LMDZ4/branches/LMDZ4-dev/libf/dyn3d/conf_gcm.F	(revision 1201)
@@ -38,4 +38,5 @@
 #include "serre.h"
 #include "comdissnew.h"
+#include "temps.h"
 #include "comconst.h"
 
@@ -111,4 +112,12 @@
       planet_type="earth"
       CALL getin('planet_type',planet_type)
+
+!Config  Key  = calend
+!Config  Desc = type de calendrier utilise
+!Config  Def  = earth_360d
+!Config  Help = valeur possible: earth_360d, earth_365d, earth_366d
+!Config         
+      calend = 'earth_360d'
+      CALL getin('calend', calend)
 
 !Config  Key  = dayref
@@ -576,4 +585,5 @@
       write(lunout,*)' Configuration des parametres du gcm: '
       write(lunout,*)' planet_type = ', planet_type
+      write(lunout,*)' calend = ', calend
       write(lunout,*)' dayref = ', dayref
       write(lunout,*)' anneeref = ', anneeref
@@ -762,4 +772,5 @@
       write(lunout,*)' Configuration des parametres du gcm: '
       write(lunout,*)' planet_type = ', planet_type
+      write(lunout,*)' calend = ', calend
       write(lunout,*)' dayref = ', dayref
       write(lunout,*)' anneeref = ', anneeref
Index: /LMDZ4/branches/LMDZ4-dev/libf/dyn3d/gcm.F
===================================================================
--- /LMDZ4/branches/LMDZ4-dev/libf/dyn3d/gcm.F	(revision 1200)
+++ /LMDZ4/branches/LMDZ4-dev/libf/dyn3d/gcm.F	(revision 1201)
@@ -113,5 +113,4 @@
       real time_step, t_wrt, t_ops
 
-      REAL rdayvrai,rdaym_ini,rday_ecri
       LOGICAL first
 
@@ -135,8 +134,8 @@
       character (len=20) :: modname
       character (len=80) :: abort_message
-
-C Calendrier
-      LOGICAL true_calendar
-      PARAMETER (true_calendar = .false.)
+! locales pour gestion du temps
+      INTEGER :: an, mois, jour
+      REAL :: heure
+
 
 c-----------------------------------------------------------------------
@@ -163,15 +162,4 @@
 
 
-c-----------------------------------------------------------------------
-c   Choix du calendrier
-c   -------------------
-
-#ifdef CPP_IOIPSL
-      if (true_calendar) then
-        call ioconf_calendar('gregorian')
-      else
-        call ioconf_calendar('360d')
-      endif
-#endif
 c----------------------------------------------------------------------
 c  lecture des fichiers gcm.def ou run.def
@@ -197,4 +185,26 @@
       endif
 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+c-----------------------------------------------------------------------
+c   Choix du calendrier
+c   -------------------
+
+c      calend = 'earth_365d'
+
+#ifdef CPP_IOIPSL
+      if (calend == 'earth_360d') then
+        call ioconf_calendar('360d')
+        write(lunout,*)'CALENDRIER CHOISI: Terrestre a 360 jours/an'
+      else if (calend == 'earth_365d') then
+        call ioconf_calendar('noleap')
+        write(lunout,*)'CALENDRIER CHOISI: Terrestre a 365 jours/an'
+      else if (calend == 'earth_366d') then
+        call ioconf_calendar('gregorian')
+        write(lunout,*)'CALENDRIER CHOISI: Terrestre bissextile'
+      else
+        abort_message = 'Mauvais choix de calendrier'
+        call abort_gcm(modname,abort_message,1)
+      endif
+#endif
+c-----------------------------------------------------------------------
 
       IF (config_inca /= 'none') THEN
@@ -303,5 +313,5 @@
           annee_ref = anneeref
           day_ref = dayref
-          day_ini = dayref
+          day_ini = 1
           itau_dyn = 0
           itau_phy = 0
@@ -314,8 +324,20 @@
       endif
 
+      mois = 1
+      heure = 0.
+      call ymds2ju(annee_ref, mois, day_ref, heure, jD_ref)
+      jH_ref = jD_ref - int(jD_ref)
+      jD_ref = int(jD_ref)
+
 #ifdef CPP_IOIPSL
       call ioconf_startdate(annee_ref,0,day_ref, 0.)
 #endif
 
+      write(lunout,*)'DEBUG'
+      write(lunout,*)'annee_ref, mois, day_ref, heure, jD_ref'
+      write(lunout,*)annee_ref, mois, day_ref, heure, jD_ref
+      call ju2ymds(jD_ref+jH_ref,an, mois, jour, heure)
+      write(lunout,*)'jD_ref+jH_ref,an, mois, jour, heure'
+      write(lunout,*)jD_ref+jH_ref,an, mois, jour, heure
 
 c  nombre d'etats dans les fichiers demarrage et histoire
@@ -391,4 +413,10 @@
       WRITE(lunout,300)day_ini,day_end
  300  FORMAT('1'/,15x,'run du jour',i7,2x,'au jour',i7//)
+      call ju2ymds(jD_ref+day_ini-1,an, mois, jour, heure)
+      write (lunout,301)jour, mois, an
+      call ju2ymds(jD_ref+day_end-1,an, mois, jour, heure)
+      write (lunout,302)jour, mois, an
+ 301  FORMAT('1'/,15x,'run du ', i2,'/',i2,'/',i4)
+ 302  FORMAT('1'/,15x,'    au ', i2,'/',i2,'/',i4)
 
       if (planet_type.eq."earth") then
Index: /LMDZ4/branches/LMDZ4-dev/libf/dyn3d/ini_paramLMDZ_dyn.h
===================================================================
--- /LMDZ4/branches/LMDZ4-dev/libf/dyn3d/ini_paramLMDZ_dyn.h	(revision 1200)
+++ /LMDZ4/branches/LMDZ4-dev/libf/dyn3d/ini_paramLMDZ_dyn.h	(revision 1201)
@@ -2,7 +2,7 @@
       dt_cum = dtvr*day_step
 
-      zan = annee_ref
-      dayref = day_ref
-      CALL ymds2ju(zan, 1, dayref, 0.0, zjulian)
+!      zan = annee_ref
+!      dayref = day_ref
+!      CALL ymds2ju(zan, 1, dayref, 0.0, zjulian)
       tau0 = itau_dyn
 c
@@ -15,5 +15,5 @@
      .                 iip1,rlong, jjp1,rlatg,
      .                 1,1,1,1,
-     .                 tau0, zjulian, dt_cum,
+     .                 tau0, jD_ref+jH_ref , dt_cum,
      .                 thoriid, nid_ctesGCM)
 c
@@ -134,5 +134,5 @@
 c
          CALL histdef(nid_ctesGCM, "true_calendar", 
-     ."Choix du calendrier: 1=gregorien ,0=calen. a 360 j",
+     ."Choix du calendrier",
      .                "-",iip1,jjp1,thoriid, 1,1,1, -99, 32, 
      .                "once", dt_cum,dt_cum)
Index: /LMDZ4/branches/LMDZ4-dev/libf/dyn3d/leapfrog.F
===================================================================
--- /LMDZ4/branches/LMDZ4-dev/libf/dyn3d/leapfrog.F	(revision 1200)
+++ /LMDZ4/branches/LMDZ4-dev/libf/dyn3d/leapfrog.F	(revision 1201)
@@ -1,2 +1,4 @@
+!
+! $Id$
 !
 c
@@ -112,6 +114,6 @@
 c
       INTEGER itau,itaufinp1,iav
-      INTEGER*4  iday ! jour julien
-      REAL       time ! Heure de la journee en fraction d'1 jour
+!      INTEGER  iday ! jour julien
+      REAL       time 
 
       REAL  SSUM
@@ -125,10 +127,16 @@
       real time_step, t_wrt, t_ops
 
-      REAL rdayvrai,rdaym_ini
+!      REAL rdayvrai,rdaym_ini
+! jD_cur: jour julien courant
+! jH_cur: heure julienne courante
+      REAL :: jD_cur, jH_cur
+      INTEGER :: an, mois, jour
+      REAL :: secondes
+
       LOGICAL first,callinigrads
 cIM : pour sortir les param. du modele dans un fis. netcdf 110106
       save first
       data first/.true./
-      real dt_cum, zjulian
+      real dt_cum
       character*10 infile
       integer zan, tau0, thoriid
@@ -167,8 +175,4 @@
       character*80 abort_message
 
-C Calendrier
-      LOGICAL true_calendar
-      PARAMETER (true_calendar = .false.)
-
       logical dissip_conservative
       save dissip_conservative
@@ -193,10 +197,10 @@
 
       itau = 0
-      iday = day_ini+itau/day_step
-      time = FLOAT(itau-(iday-day_ini)*day_step)/day_step+time_0
-         IF(time.GT.1.) THEN
-          time = time-1.
-          iday = iday+1
-         ENDIF
+c$$$      iday = day_ini+itau/day_step
+c$$$      time = FLOAT(itau-(iday-day_ini)*day_step)/day_step+time_0
+c$$$         IF(time.GT.1.) THEN
+c$$$          time = time-1.
+c$$$          iday = iday+1
+c$$$         ENDIF
 
 
@@ -214,4 +218,8 @@
 
    1  CONTINUE
+
+      jD_cur = jD_ref + (day_ini - 1) + int (itau * dtvr / daysec) 
+      jH_cur = jH_ref +                                                 &
+     &          (itau * dtvr / daysec - int(itau * dtvr / daysec)) 
 
 
@@ -284,7 +292,8 @@
       CALL geopot  ( ip1jmp1, teta  , pk , pks,  phis  , phi   )
 
+      time = jD_cur + jH_cur
       CALL caldyn 
      $  ( itau,ucov,vcov,teta,ps,masse,pk,pkf,phis ,
-     $    phi,conser,du,dv,dteta,dp,w, pbaru,pbarv, time+iday-day_ini )
+     $    phi,conser,du,dv,dteta,dp,w, pbaru,pbarv, time )
 
 
@@ -345,7 +354,12 @@
          CALL exner_hyb(  ip1jmp1, ps, p,alpha,beta,pks, pk, pkf )
 
-           rdaym_ini  = itau * dtvr / daysec
-           rdayvrai   = rdaym_ini  + day_ini
-
+!           rdaym_ini  = itau * dtvr / daysec
+!           rdayvrai   = rdaym_ini  + day_ini
+           jD_cur = jD_ref + (day_ini - 1) + int (itau * dtvr / daysec) 
+           jH_cur = jH_ref +                                            &
+     &              (itau * dtvr / daysec - int(itau * dtvr / daysec)) 
+!         write(lunout,*)'itau, jD_cur = ', itau, jD_cur, jH_cur
+!         call ju2ymds(jD_cur+jH_cur, an, mois, jour, secondes)
+!         write(lunout,*)'current date = ',an, mois, jour, secondes 
 
 c rajout debug
@@ -379,5 +393,5 @@
 #endif
 ! #endif of #ifdef CPP_IOIPSL
-         CALL calfis( lafin ,rdayvrai,time  ,
+         CALL calfis( lafin , jD_cur, jH_cur,
      $               ucov,vcov,teta,q,masse,ps,p,pk,phis,phi ,
      $               du,dv,dteta,dq,
@@ -507,10 +521,10 @@
             IF(forward. OR. leapf) THEN
               itau= itau + 1
-              iday= day_ini+itau/day_step
-              time= FLOAT(itau-(iday-day_ini)*day_step)/day_step+time_0
-                IF(time.GT.1.) THEN
-                  time = time-1.
-                  iday = iday+1
-                ENDIF
+c$$$              iday= day_ini+itau/day_step
+c$$$              time= FLOAT(itau-(iday-day_ini)*day_step)/day_step+time_0
+c$$$                IF(time.GT.1.) THEN
+c$$$                  time = time-1.
+c$$$                  iday = iday+1
+c$$$                ENDIF
             ENDIF
 
@@ -632,11 +646,11 @@
 
              itau =  itau + 1
-             iday = day_ini+itau/day_step
-             time = FLOAT(itau-(iday-day_ini)*day_step)/day_step+time_0
-
-                  IF(time.GT.1.) THEN
-                   time = time-1.
-                   iday = iday+1
-                  ENDIF
+c$$$             iday = day_ini+itau/day_step
+c$$$             time = FLOAT(itau-(iday-day_ini)*day_step)/day_step+time_0
+c$$$
+c$$$                  IF(time.GT.1.) THEN
+c$$$                   time = time-1.
+c$$$                   iday = iday+1
+c$$$                  ENDIF
 
                forward =  .FALSE.
Index: /LMDZ4/branches/LMDZ4-dev/libf/dyn3d/sortvarc.F
===================================================================
--- /LMDZ4/branches/LMDZ4-dev/libf/dyn3d/sortvarc.F	(revision 1200)
+++ /LMDZ4/branches/LMDZ4-dev/libf/dyn3d/sortvarc.F	(revision 1201)
@@ -129,6 +129,7 @@
       ang   = SSUM(     llm,  angl, 1 )
 
-      rday = FLOAT(INT ( day_ini + time ))
+c      rday = FLOAT(INT ( day_ini + time ))
 c
+       rday = FLOAT(INT(time-jD_ref-jH_ref))
       IF(ptot0.eq.0.)  THEN
          PRINT 3500, itau, rday, heure,time
@@ -156,6 +157,6 @@
       RETURN
 
-3500   FORMAT('0'10(1h*),4x,'pas'i7,5x,'jour'f5.0,'heure'f5.1,4x 
-     *   ,'date',f10.5,4x,10(1h*))
+3500   FORMAT('0'10(1h*),4x,'pas'i7,5x,'jour'f9.0,'heure'f5.1,4x 
+     *   ,'date',f14.4,4x,10(1h*))
 4000   FORMAT(10x,'masse',4x,'rmsdpdt',7x,'energie',2x,'enstrophie'
      * ,2x,'entropie',3x,'rmsv',4x,'mt.ang',/,'GLOB  '
Index: /LMDZ4/branches/LMDZ4-dev/libf/dyn3d/temps.h
===================================================================
--- /LMDZ4/branches/LMDZ4-dev/libf/dyn3d/temps.h	(revision 1200)
+++ /LMDZ4/branches/LMDZ4-dev/libf/dyn3d/temps.h	(revision 1201)
@@ -1,4 +1,4 @@
 !
-! $Header$
+! $Id$
 !
 !  ATTENTION!!!!: ce fichier include est compatible format fixe/format libre
@@ -8,14 +8,17 @@
 !
 !
+! jD_ref = jour julien de la date de reference (lancement de l'experience)
+! hD_ref = "heure" julienne de la date de reference
 !-----------------------------------------------------------------------
 ! INCLUDE 'temps.h'
 
       COMMON/temps/itaufin, dt, day_ini, day_end, annee_ref, day_ref,   &
-     &             itau_dyn, itau_phy
+     &             itau_dyn, itau_phy, jD_ref, jH_ref, calend
 
       INTEGER   itaufin
-      INTEGER(kind=4) itau_dyn, itau_phy
-      INTEGER(kind=4) day_ini, day_end, annee_ref, day_ref
-      REAL      dt
+      INTEGER itau_dyn, itau_phy
+      INTEGER day_ini, day_end, annee_ref, day_ref
+      REAL      dt, jD_ref, jH_ref
+      CHARACTER (len=10) :: calend
 
 !-----------------------------------------------------------------------
Index: /LMDZ4/branches/LMDZ4-dev/libf/dyn3d/write_paramLMDZ_dyn.h
===================================================================
--- /LMDZ4/branches/LMDZ4-dev/libf/dyn3d/write_paramLMDZ_dyn.h	(revision 1200)
+++ /LMDZ4/branches/LMDZ4-dev/libf/dyn3d/write_paramLMDZ_dyn.h	(revision 1201)
@@ -107,9 +107,12 @@
      .               zx_tmp_2d,iip1*jjp1,ndex2d)
 c
-      IF(true_calendar) THEN
-       zx_tmp_2d(1:iip1,1:jjp1)=1.
-      ELSE
-       zx_tmp_2d(1:iip1,1:jjp1)=0.
-      ENDIF
+      if (calend == 'earth_360d') then
+        zx_tmp_2d(1:iip1,1:jjp1)=1.
+      else if (calend == 'earth_365d') then
+        zx_tmp_2d(1:iip1,1:jjp1)=2.
+      else if (calend == 'earth_366d') then
+        zx_tmp_2d(1:iip1,1:jjp1)=3.
+      endif
+
       CALL histwrite(nid_ctesGCM, "true_calendar", itau_w,
      .               zx_tmp_2d,iip1*jjp1,ndex2d)
Index: /LMDZ4/branches/LMDZ4-dev/libf/dyn3dpar/calfis_p.F
===================================================================
--- /LMDZ4/branches/LMDZ4-dev/libf/dyn3dpar/calfis_p.F	(revision 1200)
+++ /LMDZ4/branches/LMDZ4-dev/libf/dyn3dpar/calfis_p.F	(revision 1201)
@@ -1,9 +1,9 @@
 !
-! $Header$
+! $Id$
 !
 C
 C
       SUBROUTINE calfis_p(lafin,
-     $                  rdayvrai,
+     $                  jD_cur, jH_cur,
      $                  heure,
      $                  pucov,
@@ -209,5 +209,5 @@
       SAVE firstcal,debut
 c$OMP THREADPRIVATE(firstcal,debut)
-      REAL rdayvrai
+      REAL :: jD_cur, jH_cur
       
       REAL,SAVE,dimension(1:iim,1:llm):: du_send,du_recv,dv_send,dv_recv
@@ -334,5 +334,5 @@
 c   ---------------
 c
-
+n
       DO iq=1,nqtot
          iiq=niadv(iq)
@@ -665,6 +665,6 @@
      .             debut,
      .             lafin,
-     .             rdayvrai,
-     .             heure,
+     .             jD_cur,
+     .             jH_cur,
      .             dtphys,
      .             zplev_omp,
Index: /LMDZ4/branches/LMDZ4-dev/libf/dyn3dpar/conf_gcm.F
===================================================================
--- /LMDZ4/branches/LMDZ4-dev/libf/dyn3dpar/conf_gcm.F	(revision 1200)
+++ /LMDZ4/branches/LMDZ4-dev/libf/dyn3dpar/conf_gcm.F	(revision 1201)
@@ -44,4 +44,5 @@
 !#include "clesphys.h"
 #include "iniprint.h"
+#include "temps.h"
 #include "comconst.h"
 
@@ -122,4 +123,12 @@
       CALL getin('planet_type',planet_type)
 
+!Config  Key  = calend
+!Config  Desc = type de calendrier utilise
+!Config  Def  = earth_360d
+!Config  Help = valeur possible: earth_360d, earth_365d, earth_366d
+!Config         
+      calend = 'earth_360d'
+      CALL getin('calend', calend)
+
 !Config  Key  = dayref
 !Config  Desc = Jour de l'etat initial
@@ -588,4 +597,5 @@
       write(lunout,*)' Configuration des parametres du gcm: '
       write(lunout,*)' planet_type = ', planet_type
+      write(lunout,*)' calend = ', calend
       write(lunout,*)' dayref = ', dayref
       write(lunout,*)' anneeref = ', anneeref
@@ -812,4 +822,5 @@
       write(lunout,*)' Configuration des parametres du gcm: '
       write(lunout,*)' planet_type = ', planet_type
+      write(lunout,*)' calend = ', calend
       write(lunout,*)' dayref = ', dayref
       write(lunout,*)' anneeref = ', anneeref
Index: /LMDZ4/branches/LMDZ4-dev/libf/dyn3dpar/gcm.F
===================================================================
--- /LMDZ4/branches/LMDZ4-dev/libf/dyn3dpar/gcm.F	(revision 1200)
+++ /LMDZ4/branches/LMDZ4-dev/libf/dyn3dpar/gcm.F	(revision 1201)
@@ -1,4 +1,4 @@
 !
-! $Header$
+! $Id$
 !
 c
@@ -113,6 +113,4 @@
       real time_step, t_wrt, t_ops
 
-c      REAL rdayvrai,rdaym_ini,rday_ecri
-c      LOGICAL first
 
       LOGICAL call_iniphys
@@ -133,10 +131,10 @@
 
       character (len=80) :: dynhist_file, dynhistave_file
-      character (len=20) ::modname
-      character (len=80) ::abort_message
-
-C Calendrier
-      LOGICAL true_calendar
-      PARAMETER (true_calendar = .false.)
+      character (len=20) :: modname
+      character (len=80) :: abort_message
+! locales pour gestion du temps
+      INTEGER :: an, mois, jour
+      REAL :: heure
+
 
 c-----------------------------------------------------------------------
@@ -165,15 +163,5 @@
 
 
-c-----------------------------------------------------------------------
-c   Choix du calendrier
-c   -------------------
-
-#ifdef CPP_IOIPSL
-      if (true_calendar) then
-        call ioconf_calendar('gregorian')
-      else
-        call ioconf_calendar('360d')
-      endif
-#endif
+
 c----------------------------------------------------------------------
 c  lecture des fichiers gcm.def ou run.def
@@ -221,4 +209,26 @@
 #endif
       endif ! of if (planet_type.eq."earth")
+
+c-----------------------------------------------------------------------
+c   Choix du calendrier
+c   -------------------
+
+c      calend = 'earth_365d'
+
+#ifdef CPP_IOIPSL
+      if (calend == 'earth_360d') then
+        call ioconf_calendar('360d')
+        write(lunout,*)'CALENDRIER CHOISI: Terrestre a 360 jours/an'
+      else if (calend == 'earth_365d') then
+        call ioconf_calendar('noleap')
+        write(lunout,*)'CALENDRIER CHOISI: Terrestre a 365 jours/an'
+      else if (calend == 'earth_366d') then
+        call ioconf_calendar('gregorian')
+        write(lunout,*)'CALENDRIER CHOISI: Terrestre bissextile'
+      else
+        abort_message = 'Mauvais choix de calendrier'
+        call abort_gcm(modname,abort_message,1)
+      endif
+#endif
 
       IF (config_inca /= 'none') THEN
@@ -306,5 +316,5 @@
       if (annee_ref .ne. anneeref .or. day_ref .ne. dayref) then
         write(lunout,*)
-     .  ' Attention les dates initiales lues dans le fichier'
+     .  'GCM: Attention les dates initiales lues dans le fichier'
         write(lunout,*)
      .  ' restart ne correspondent pas a celles lues dans '
@@ -312,14 +322,14 @@
         if (raz_date .ne. 1) then
           write(lunout,*)
-     .    ' On garde les dates du fichier restart'
+     .    'GCM: On garde les dates du fichier restart'
         else
           annee_ref = anneeref
           day_ref = dayref
-          day_ini = dayref
+          day_ini = 1
           itau_dyn = 0
           itau_phy = 0
           time_0 = 0.
           write(lunout,*)
-     .   ' On reinitialise a la date lue dans gcm.def'
+     .   'GCM: On reinitialise a la date lue dans gcm.def'
         endif
       ELSE
@@ -327,7 +337,20 @@
       endif
 
-#ifdef CPP_IOIPSL 
-      call ioconf_startdate(annee_ref,0,day_ref,0.) 
-#endif 
+      mois = 1
+      heure = 0.
+      call ymds2ju(annee_ref, mois, day_ref, heure, jD_ref)
+      jH_ref = jD_ref - int(jD_ref)
+      jD_ref = int(jD_ref)
+
+#ifdef CPP_IOIPSL
+      call ioconf_startdate(annee_ref,0,day_ref, 0.)
+#endif
+
+      write(lunout,*)'DEBUG'
+      write(lunout,*)'annee_ref, mois, day_ref, heure, jD_ref'
+      write(lunout,*)annee_ref, mois, day_ref, heure, jD_ref
+      call ju2ymds(jD_ref+jH_ref,an, mois, jour, heure)
+      write(lunout,*)'jD_ref+jH_ref,an, mois, jour, heure'
+      write(lunout,*)jD_ref+jH_ref,an, mois, jour, heure
 
 c  nombre d'etats dans les fichiers demarrage et histoire
@@ -413,4 +436,10 @@
       WRITE(lunout,300)day_ini,day_end
  300  FORMAT('1'/,15x,'run du jour',i7,2x,'au jour',i7//)
+      call ju2ymds(jD_ref+day_ini-1,an, mois, jour, heure)
+      write (lunout,301)jour, mois, an
+      call ju2ymds(jD_ref+day_end-1,an, mois, jour, heure)
+      write (lunout,302)jour, mois, an
+ 301  FORMAT('1'/,15x,'run du ', i2,'/',i2,'/',i4)
+ 302  FORMAT('1'/,15x,'    au ', i2,'/',i2,'/',i4)
 
 !#ifdef CPP_IOIPSL
Index: /LMDZ4/branches/LMDZ4-dev/libf/dyn3dpar/leapfrog_p.F
===================================================================
--- /LMDZ4/branches/LMDZ4-dev/libf/dyn3dpar/leapfrog_p.F	(revision 1200)
+++ /LMDZ4/branches/LMDZ4-dev/libf/dyn3dpar/leapfrog_p.F	(revision 1201)
@@ -1,4 +1,4 @@
 ! 
-! $Header$
+! $Id$
 !
 c
@@ -120,6 +120,6 @@
 c
       INTEGER itau,itaufinp1,iav
-      INTEGER*4  iday ! jour julien
-      REAL       time ! Heure de la journee en fraction d'1 jour
+!      INTEGER  iday ! jour julien
+      REAL       time 
 
       REAL  SSUM
@@ -134,5 +134,10 @@
       real time_step, t_wrt, t_ops
 
-      REAL rdayvrai,rdaym_ini
+! jD_cur: jour julien courant
+! jH_cur: heure julienne courante
+      REAL :: jD_cur, jH_cur
+      INTEGER :: an, mois, jour
+      REAL :: secondes
+
       LOGICAL first,callinigrads
 
@@ -162,7 +167,4 @@
       character*80 abort_message
 
-C Calendrier
-      LOGICAL true_calendar
-      PARAMETER (true_calendar = .false.)
 
       logical,PARAMETER :: dissip_conservative=.TRUE.
@@ -207,10 +209,10 @@
 
       itau = 0
-      iday = day_ini+itau/day_step
-      time = FLOAT(itau-(iday-day_ini)*day_step)/day_step+time_0
-         IF(time.GT.1.) THEN
-          time = time-1.
-          iday = iday+1
-         ENDIF
+c$$$      iday = day_ini+itau/day_step
+c$$$      time = FLOAT(itau-(iday-day_ini)*day_step)/day_step+time_0
+c$$$         IF(time.GT.1.) THEN
+c$$$          time = time-1.
+c$$$          iday = iday+1
+c$$$         ENDIF
 
 c Allocate variables depending on dynamic variable nqtot
@@ -241,10 +243,7 @@
    1  CONTINUE
 
-c$OMP MASTER
-
-      CALL barrier
-      
-c$OMP END MASTER
-c$OMP BARRIER
+      jD_cur = jD_ref + (day_ini - 1) + int (itau * dtvr / daysec) 
+      jH_cur = jH_ref +                                                 &
+     &          (itau * dtvr / daysec - int(itau * dtvr / daysec)) 
 
 
@@ -551,7 +550,8 @@
 c$OMP BARRIER
 !      CALL FTRACE_REGION_BEGIN("caldyn")
+      time = jD_cur + jH_cur 
       CALL caldyn_p 
      $  ( itau,ucov,vcov,teta,ps,masse,pk,pkf,phis ,
-     $    phi,conser,du,dv,dteta,dp,w, pbaru,pbarv, time+iday-day_ini )
+     $    phi,conser,du,dv,dteta,dp,w, pbaru,pbarv, time )
 
 !      CALL FTRACE_REGION_END("caldyn")
@@ -684,7 +684,8 @@
          CALL exner_hyb_p(  ip1jmp1, ps, p,alpha,beta,pks, pk, pkf )
 c$OMP BARRIER
-           rdaym_ini  = itau * dtvr / daysec
-           rdayvrai   = rdaym_ini  + day_ini
-
+           jD_cur = jD_ref + (day_ini -1) + int (itau * dtvr / daysec) 
+           jH_cur = jH_ref +                                            &
+     &              (itau * dtvr / daysec - int(itau * dtvr / daysec)) 
+         call ju2ymds(jD_cur+jH_cur, an, mois, jour, secondes)
 
 c rajout debug
@@ -773,5 +774,5 @@
 cc$OMP BARRIER
 !        CALL FTRACE_REGION_BEGIN("calfis")
-        CALL calfis_p(lafin ,rdayvrai,time  ,
+        CALL calfis_p(lafin ,jD_cur, jH_cur,
      $               ucov,vcov,teta,q,masse,ps,p,pk,phis,phi ,
      $               du,dv,dteta,dq,
@@ -1274,10 +1275,10 @@
             IF(forward. OR. leapf) THEN
               itau= itau + 1
-              iday= day_ini+itau/day_step
-              time= FLOAT(itau-(iday-day_ini)*day_step)/day_step+time_0
-                IF(time.GT.1.) THEN
-                  time = time-1.
-                  iday = iday+1
-                ENDIF
+c$$$              iday= day_ini+itau/day_step
+c$$$              time= FLOAT(itau-(iday-day_ini)*day_step)/day_step+time_0
+c$$$                IF(time.GT.1.) THEN
+c$$$                  time = time-1.
+c$$$                  iday = iday+1
+c$$$                ENDIF
             ENDIF
 
@@ -1444,11 +1445,11 @@
 
              itau =  itau + 1
-             iday = day_ini+itau/day_step
-             time = FLOAT(itau-(iday-day_ini)*day_step)/day_step+time_0
-
-                  IF(time.GT.1.) THEN
-                   time = time-1.
-                   iday = iday+1
-                  ENDIF
+c$$$             iday = day_ini+itau/day_step
+c$$$             time = FLOAT(itau-(iday-day_ini)*day_step)/day_step+time_0
+c$$$
+c$$$                  IF(time.GT.1.) THEN
+c$$$                   time = time-1.
+c$$$                   iday = iday+1
+c$$$                  ENDIF
 
                forward =  .FALSE.
Index: /LMDZ4/branches/LMDZ4-dev/libf/dyn3dpar/sortvarc.F
===================================================================
--- /LMDZ4/branches/LMDZ4-dev/libf/dyn3dpar/sortvarc.F	(revision 1200)
+++ /LMDZ4/branches/LMDZ4-dev/libf/dyn3dpar/sortvarc.F	(revision 1201)
@@ -129,6 +129,7 @@
       ang   = SSUM(     llm,  angl, 1 )
 
-      rday = FLOAT(INT ( day_ini + time ))
+c      rday = FLOAT(INT ( day_ini + time ))
 c
+       rday = FLOAT(INT(time-jD_ref-jH_ref))
       IF(ptot0.eq.0.)  THEN
          PRINT 3500, itau, rday, heure,time
@@ -156,6 +157,6 @@
       RETURN
 
-3500   FORMAT('0'10(1h*),4x,'pas'i7,5x,'jour'f5.0,'heure'f5.1,4x 
-     *   ,'date',f10.5,4x,10(1h*))
+3500   FORMAT('0'10(1h*),4x,'pas'i7,5x,'jour'f9.0,'heure'f5.1,4x 
+     *   ,'date',f14.4,4x,10(1h*))
 4000   FORMAT(10x,'masse',4x,'rmsdpdt',7x,'energie',2x,'enstrophie'
      * ,2x,'entropie',3x,'rmsv',4x,'mt.ang',/,'GLOB  '
Index: /LMDZ4/branches/LMDZ4-dev/libf/dyn3dpar/temps.h
===================================================================
--- /LMDZ4/branches/LMDZ4-dev/libf/dyn3dpar/temps.h	(revision 1200)
+++ /LMDZ4/branches/LMDZ4-dev/libf/dyn3dpar/temps.h	(revision 1201)
@@ -1,4 +1,4 @@
 !
-! $Header$
+! $Id$
 !
 !  ATTENTION!!!!: ce fichier include est compatible format fixe/format libre
@@ -8,13 +8,17 @@
 !
 !
+! jD_ref = jour julien de la date de reference (lancement de l'experience)
+! hD_ref = "heure" julienne de la date de reference
 !-----------------------------------------------------------------------
 ! INCLUDE 'temps.h'
 
       COMMON/temps/itaufin, dt, day_ini, day_end, annee_ref, day_ref,   &
-     &             itau_dyn, itau_phy
+     &             itau_dyn, itau_phy, jD_ref, jH_ref, calend
 
       INTEGER   itaufin
-      INTEGER(kind=4) itau_dyn, itau_phy
-      INTEGER(kind=4) day_ini, day_end, annee_ref, day_ref
-      REAL      dt
-!$OMP THREADPRIVATE(/temps/)
+      INTEGER itau_dyn, itau_phy
+      INTEGER day_ini, day_end, annee_ref, day_ref
+      REAL      dt, jD_ref, jH_ref
+      CHARACTER (len=10) :: calend
+
+!-----------------------------------------------------------------------
Index: /LMDZ4/branches/LMDZ4-dev/libf/phylmd/hgardfou.F
===================================================================
--- /LMDZ4/branches/LMDZ4-dev/libf/phylmd/hgardfou.F	(revision 1200)
+++ /LMDZ4/branches/LMDZ4-dev/libf/phylmd/hgardfou.F	(revision 1201)
@@ -1,3 +1,4 @@
 !
+! $Id$
       SUBROUTINE hgardfou (t,tsol,text)
       use dimphy
@@ -12,4 +13,6 @@
       REAL t(klon,klev), tsol(klon,nbsrf)
       CHARACTER*(*) text
+      character (len=20) :: modname = 'hgardfou'
+      character (len=80) :: abort_message
 C
       INTEGER i, k, nsrf
@@ -124,6 +127,6 @@
 c
       IF (.NOT. ok) THEN
-         PRINT*, 'hgardfou s arrete ', text
-         CALL abort
+         abort_message= 'hgardfou s arrete '//text
+         CALL abort_gcm (modname,abort_message,1)
       ENDIF
 
Index: /LMDZ4/branches/LMDZ4-dev/libf/phylmd/iniorbit.F
===================================================================
--- /LMDZ4/branches/LMDZ4-dev/libf/phylmd/iniorbit.F	(revision 1201)
+++ /LMDZ4/branches/LMDZ4-dev/libf/phylmd/iniorbit.F	(revision 1201)
@@ -0,0 +1,104 @@
+      SUBROUTINE iniorbit
+     $     (paphelie,pperiheli,pyear_day,pperi_day,pobliq)
+      IMPLICIT NONE
+
+c=======================================================================
+c
+c   Auteur:
+c   -------
+c     Frederic Hourdin      22 Fevrier 1991
+c
+c   Objet:
+c   ------
+c    Initialisation du sous programme orbite qui calcule
+c    a une date donnee de l'annee de duree year_day commencant
+c    a l'equinoxe de printemps et dont le perihelie se situe
+c    a la date peri_day, la distance au soleil et la declinaison.
+c
+c   Interface:
+c   ----------
+c   - Doit etre appele avant d'utiliser orbite.
+c   - initialise une partie du common planete.h
+c
+c   Arguments:
+c   ----------
+c
+c   Input:
+c   ------
+c   aphelie       \   aphelie et perihelie de l'orbite
+c   periheli      /   en millions de kilometres.
+c
+c=======================================================================
+
+c-----------------------------------------------------------------------
+c   Declarations:
+c   -------------
+
+#include "planete.h"
+#include "YOMCST.h"
+
+c   Arguments:
+c   ----------
+
+      REAL paphelie,pperiheli,pyear_day,pperi_day,pobliq
+
+c   Local:
+c   ------
+
+      REAL zxref,zanom,zz,zx0,zdx, pi
+      INTEGER iter
+
+c-----------------------------------------------------------------------
+
+      pi=2.*asin(1.)
+
+      aphelie =paphelie
+      periheli=pperiheli
+      year_day=pyear_day
+      obliquit=pobliq
+      peri_day=pperi_day
+
+      PRINT*,'Perihelie en Mkm  ',periheli
+      PRINT*,'Aphelie  en Mkm   ',aphelie 
+      PRINT*,'obliquite en degres  :',obliquit
+      PRINT*,'Jours dans l annee : ',year_day
+      PRINT*,'Date perihelie : ',peri_day
+      unitastr=149.597870
+      e_elips=(aphelie-periheli)/(periheli+aphelie)
+      p_elips=0.5*(periheli+aphelie)*(1-e_elips*e_elips)/unitastr
+
+      print*,'e_elips',e_elips
+      print*,'p_elips',p_elips
+
+c-----------------------------------------------------------------------
+c calcul de l'angle polaire et de la distance au soleil :
+c -------------------------------------------------------
+
+c  calcul de l'zanomalie moyenne
+
+      zz=(year_day-pperi_day)/year_day
+      zanom=2.*pi*(zz-nint(zz))
+      zxref=abs(zanom)
+      PRINT*,'zanom  ',zanom
+
+c  resolution de l'equation horaire  zx0 - e * sin (zx0) = zxref
+c  methode de Newton
+
+      zx0=zxref+R_ecc*sin(zxref)
+      DO 110 iter=1,100
+         zdx=-(zx0-R_ecc*sin(zx0)-zxref)/(1.-R_ecc*cos(zx0))
+         if(abs(zdx).le.(1.e-12)) goto 120
+         zx0=zx0+zdx
+110   continue
+120   continue
+      zx0=zx0+zdx
+      if(zanom.lt.0.) zx0=-zx0
+      PRINT*,'zx0   ',zx0
+
+c zteta est la longitude solaire
+
+      timeperi=2.*atan(sqrt((1.+R_ecc)/(1.-R_ecc))*tan(zx0/2.))
+      PRINT*,'longitude solaire du perihelie timeperi = ',timeperi
+
+      RETURN
+      END
Index: /LMDZ4/branches/LMDZ4-dev/libf/phylmd/physiq.F
===================================================================
--- /LMDZ4/branches/LMDZ4-dev/libf/phylmd/physiq.F	(revision 1200)
+++ /LMDZ4/branches/LMDZ4-dev/libf/phylmd/physiq.F	(revision 1201)
@@ -4,5 +4,5 @@
 
       SUBROUTINE physiq (nlon,nlev,
-     .            debut,lafin,rjourvrai,gmtime,pdtphys,
+     .            debut,lafin,jD_cur, jH_cur,pdtphys,
      .            paprs,pplay,pphi,pphis,presnivs,clesphy0,
      .            u,v,t,qx,
@@ -60,6 +60,6 @@
 c debut---input-L-variable logique indiquant le premier passage
 c lafin---input-L-variable logique indiquant le dernier passage
-c rjour---input-R-numero du jour de l'experience
-c gmtime--input-R-temps universel dans la journee (0 a 86400 s)
+c jD_cur--input-R-jour courant a l'appel de la physique (jour julien)
+c jH_cur--input-R-heure courante a l'appel de la physique (jour julien)
 c pdtphys-input-R-pas d'integration pour la physique (seconde)
 c paprs---input-R-pression pour chaque inter-couche (en Pa)
@@ -186,6 +186,6 @@
       INTEGER nlon
       INTEGER nlev
-      REAL rjourvrai
-      REAL gmtime
+      REAL :: jD_cur, jH_cur
+
       REAL pdtphys
       LOGICAL debut, lafin
@@ -703,5 +703,16 @@
 c Conditions aux limites
 c
-      INTEGER julien
+!
+! Gestion calendrier
+!
+      REAL :: jD_1jan, jH_1jan
+      INTEGER :: year_cur, mth_cur, day_cur, days_elapsed
+      REAL :: hour, day_since_equinox
+! Date de l'equinoxe de printemps
+      INTEGER, parameter :: mth_eq=3, day_eq=21
+      REAL :: jD_eq
+
+      LOGICAL, parameter :: new_orbit = .true.
+
 c
       INTEGER lmt_pas
@@ -1087,4 +1098,6 @@
 c$OMP THREADPRIVATE(first)
 
+      integer iunit
+
       logical, save::  read_climoz ! read ozone climatology
       integer, save:: ncid_climoz ! NetCDF file containing ozone climatology
@@ -1114,7 +1127,7 @@
          write(lunout,*) 'DEBUT DE PHYSIQ !!!!!!!!!!!!!!!!!!!!'
          write(lunout,*)
-     s 'nlon,klev,nqtot,debut,lafin,rjourvrai,gmtime,pdtphys'
+     s 'nlon,klev,nqtot,debut,lafin, jD_cur, jH_cur,pdtphys'
          write(lunout,*)
-     s  nlon,klev,nqtot,debut,lafin,rjourvrai,gmtime,pdtphys 
+     s  nlon,klev,nqtot,debut,lafin, jD_cur, jH_cur,pdtphys 
 
          write(lunout,*) 'papers, play, phi, u, v, t, omega'
@@ -1184,5 +1197,13 @@
 
 c======================================================================
-      xjour = rjourvrai
+! Gestion calendrier
+!
+      call ju2ymds(jD_cur+jH_cur, year_cur, mth_cur, day_cur, hour)
+      call ymds2ju(year_cur, 1, 1, 0., jD_1jan)
+      jH_1jan = jD_1jan - int (jD_1jan)
+      jD_1jan = int (jD_1jan) 
+      xjour = jD_cur - jD_1jan
+      days_elapsed = jD_cur - jD_1jan
+
 c
 c Si c'est le debut, il faut initialiser plusieurs choses
@@ -1472,6 +1493,5 @@
 
 cXXXPB Positionner date0 pour initialisation de ORCHIDEE
-      date0 = zjulian
-C      date0 = day_ini
+      date0 = jD_ref 
       WRITE(*,*) 'physiq date0 : ',date0
 c
@@ -1491,7 +1511,8 @@
          CALL VTe(VTphysiq)
          CALL VTb(VTinca)
-         iii = MOD(NINT(xjour),360)
-         calday = FLOAT(iii) + gmtime
-         WRITE(lunout,*) 'initial time ', xjour, calday
+!         iii = MOD(NINT(xjour),360)
+!         calday = FLOAT(iii) + jH_cur
+         calday = FLOAT(days_elapsed) + jH_cur
+         WRITE(lunout,*) 'initial time chemini', days_elapsed, calday
 
          CALL chemini( 
@@ -1534,7 +1555,4 @@
 !
       itap   = itap + 1
-      julien = MOD(NINT(xjour),360)
-      if (julien .eq. 0) julien = 360
-
 !
 ! Update fraction of the sub-surfaces (pctsrf) and 
@@ -1542,7 +1560,6 @@
 ! on the surface fraction.
 !
-      CALL change_srf_frac(itap, dtime, julien, 
+      CALL change_srf_frac(itap, dtime, days_elapsed+1, 
      *     pctsrf, falb1, falb2, ftsol, u10m, v10m, pbl_tke)
-
 
 ! Tendances bidons pour les processus qui n'affectent pas certaines
@@ -1697,5 +1714,6 @@
          if (read_climoz) then
 C           Ozone climatology from a NetCDF file
-            call regr_pr_av(ncid_climoz, "tro3", julien, press_climoz,
+            call regr_pr_av(ncid_climoz, "tro3", days_elapsed+1, 
+     &           press_climoz,
      $           paprs, wo)
 !           Convert from mole fraction of ozone to column density of ozone in a
@@ -1707,5 +1725,6 @@
 C           "zmasse" changes a little.)
          else
-            CALL ozonecm(real(julien), rlat, paprs, wo)
+            CALL ozonecm(real(days_elapsed+1), rlat, paprs, wo)
+
          end if
       ENDIF
@@ -1749,20 +1768,28 @@
 ! doit donc etre placé avant radlwsw et pbl_surface
 
+! calcul selon la routine utilisee pour les planetes
+      if (new_orbit) then
+        call ymds2ju(year_cur, mth_eq, day_eq,0., jD_eq)
+        day_since_equinox = (jD_cur + jH_cur) - jD_eq 
+!        day_since_equinox = (jD_cur) - jD_eq 
+        call solarlong(day_since_equinox, zlongi, dist)
+      else      
+! calcul selon la routine utilisee pour l'AR4
 !   choix entre calcul de la longitude solaire vraie ou valeur fixee a
 !   solarlong0
-
-      if (solarlong0<-999.) then
-         CALL orbite(FLOAT(julien),zlongi,dist)
-      else
-         zlongi=solarlong0  ! longitude solaire vraie
-         dist=1.            ! distance au soleil / moyenne
+        if (solarlong0<-999.) then
+           CALL orbite(FLOAT(days_elapsed+1),zlongi,dist)
+        else
+           zlongi=solarlong0  ! longitude solaire vraie
+           dist=1.            ! distance au soleil / moyenne
+        endif
       endif
-
-      if(prt_level.ge.1) print*,'Longitude solaire ',zlongi,solarlong0
+      if(prt_level.ge.1)                                                &
+     &    write(lunout,*)'Longitude solaire ',zlongi,solarlong0,dist
 
 !  Avec ou sans cycle diurne
       IF (cycle_diurne) THEN
         zdtime=dtime*FLOAT(radpas) ! pas de temps du rayonnement (s)
-        CALL zenang(zlongi,gmtime,zdtime,rlat,rlon,rmu0,fract)
+        CALL zenang(zlongi,jH_cur,zdtime,rlat,rlon,rmu0,fract)
       ELSE
         CALL angle(zlongi, rlat, fract, rmu0)
@@ -1797,5 +1824,5 @@
 
       CALL pbl_surface( 
-     e     dtime,     date0,     itap,    julien,
+     e     dtime,     date0,     itap,    days_elapsed+1,
      e     debut,     lafin,
      e     rlon,      rlat,      rugoro,  rmu0,     
@@ -2094,6 +2121,8 @@
 
             if (itop_con(i).gt.klev-3) then
-               print*,'La convection monte trop haut '
-               print*,'itop_con(,',i,',)=',itop_con(i)
+              if(prt_level >= 9) then
+                write(lunout,*)'La convection monte trop haut '
+                write(lunout,*)'itop_con(,',i,',)=',itop_con(i)
+              endif
             endif
           ENDDO
@@ -2635,5 +2664,5 @@
          IF (.NOT. aerosol_couple)
      &        CALL readaerosol_optic(
-     &        debut, new_aod, flag_aerosol, rjourvrai, pdtphys,
+     &        debut, new_aod, flag_aerosol, jD_cur-jD_ref, pdtphys,
      &        pplay, paprs, t_seri, rhcl, 
      &        mass_solu_aero, mass_solu_aero_pi, 
@@ -2764,5 +2793,5 @@
          CALL VTe(VTphysiq)
          CALL VTb(VTinca)
-         calday = FLOAT(julien) + gmtime
+         calday = FLOAT(days_elapsed + 1) + jH_cur
 
          IF (config_inca == 'aero') THEN
@@ -2775,6 +2804,6 @@
 
          CALL chemhook_begin (calday,
-     $                          julien,
-     $                          gmtime,
+     $                          days_elapsed,
+     $                          jH_cur,
      $                          pctsrf(1,1),
      $                          rlat,
@@ -3106,5 +3135,5 @@
       IF (is_sequential) THEN
       
-        CALL aaam_bud (27,klon,klev,rjourvrai,gmtime,
+        CALL aaam_bud (27,klon,klev,jD_cur-jD_ref,jH_cur,
      C                 ra,rg,romega,
      C                 rlat,rlon,pphis,
@@ -3133,5 +3162,5 @@
 
       call phytrac (
-     I     itap,     julien,    gmtime,   debut,
+     I     itap,     days_elapsed+1,    jH_cur,   debut,
      I     lafin,    dtime,     u, v,     t,
      I     paprs,    pplay,     pmfu,     pmfd, 
@@ -3351,7 +3380,7 @@
       write(lunout,*) 'FIN DE PHYSIQ !!!!!!!!!!!!!!!!!!!!'
       write(lunout,*)
-     s 'nlon,klev,nqtot,debut,lafin,rjourvrai,gmtime,pdtphys pct tlos'
+     s 'nlon,klev,nqtot,debut,lafin,jD_cur, jH_cur, pdtphys pct tlos'
       write(lunout,*)
-     s  nlon,klev,nqtot,debut,lafin,rjourvrai,gmtime,pdtphys,
+     s  nlon,klev,nqtot,debut,lafin, jD_cur, jH_cur ,pdtphys,
      s  pctsrf(igout,is_ter), pctsrf(igout,is_lic),pctsrf(igout,is_oce),
      s  pctsrf(igout,is_sic)
@@ -3447,4 +3476,5 @@
       ENDIF
       
+!      first=.false.
 
       RETURN
Index: /LMDZ4/branches/LMDZ4-dev/libf/phylmd/solarlong.F
===================================================================
--- /LMDZ4/branches/LMDZ4-dev/libf/phylmd/solarlong.F	(revision 1201)
+++ /LMDZ4/branches/LMDZ4-dev/libf/phylmd/solarlong.F	(revision 1201)
@@ -0,0 +1,131 @@
+      SUBROUTINE solarlong(pday,psollong,pdist_sol)
+
+      USE ioipsl
+
+      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 "planete.h"
+#include "YOMCST.h"
+      include 'iniprint.h'
+
+c arguments:
+c ----------
+
+      REAL pday,pdist_sol,pdecli,psollong
+      LOGICAL lwrite
+
+c Local:
+c ------
+
+      REAL zanom,xref,zx0,zdx,zteta,zz,pi
+      INTEGER iter
+      REAL :: pyear_day,pperi_day
+      REAL :: jD_eq, jD_peri
+
+c-----------------------------------------------------------------------
+c calcul de l'angle polaire et de la distance au soleil :
+c -------------------------------------------------------
+
+c   Initialisation eventuelle:
+      if(.not.unitastr.gt.1.e-4) then
+        call ioget_calendar(pyear_day)
+        call ymds2ju(2000, 3, 21, 0., jD_eq)
+        call ymds2ju(2001, 1, 4, 0., jD_peri)
+        pperi_day = jD_peri - jD_eq
+        pperi_day = R_peri + 180.
+        write(lunout,*)' Number of days in a year = ',pyear_day
+c         call iniorbit(249.22,206.66,669.,485.,25.2)
+         call iniorbit(152.59,146.61,pyear_day,pperi_day,R_incl)
+      endif
+
+c  calcul de l'zanomalie moyenne
+
+      zz=(pday-peri_day)/year_day
+      pi=2.*asin(1.)
+      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)
+      zx0=xref+R_ecc*sin(xref)
+      DO 110 iter=1,10
+!         zdx=-(zx0-e_elips*sin(zx0)-xref)/(1.-e_elips*cos(zx0))
+         zdx=-(zx0-R_ecc*sin(zx0)-xref)/(1.-R_ecc*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.))
+      zteta=2.*atan(sqrt((1.+R_ecc)/(1.-R_ecc))*tan(zx0/2.))
+
+      psollong=zteta-timeperi
+
+      IF(psollong.LT.0.) psollong=psollong+2.*pi
+      IF(psollong.GT.2.*pi) psollong=psollong-2.*pi
+
+      psollong = psollong * 180. / pi
+
+c distance soleil
+
+      pdist_sol = (1-R_ecc*R_ecc)
+     &      /(1+R_ecc*COS(pi/180.*(psollong-(R_peri+180.0))))
+!      pdist_sol = (1-e_elips*e_elips)
+!     &      /(1+e_elips*COS(pi/180.*(psollong-(R_peri+180.0))))
+c-----------------------------------------------------------------------
+c   sorties eventuelles:
+c   ---------------------
+
+c     IF (lwrite) THEN
+c        PRINT*,'jour de l"annee   :',pday
+c        PRINT*,'distance au soleil (en unite astronomique) :',pdist_sol
+c        PRINT*,'declinaison (en degres) :',pdecli*180./pi
+c     ENDIF
+
+      RETURN
+      END
