Index: trunk/LMDZ.COMMON/libf/dyn3d_common/dynredem.F90
===================================================================
--- trunk/LMDZ.COMMON/libf/dyn3d_common/dynredem.F	(revision 1422)
+++ trunk/LMDZ.COMMON/libf/dyn3d_common/dynredem.F90	(revision 1508)
@@ -2,646 +2,326 @@
 ! $Id: dynredem.F 1635 2012-07-12 11:37:16Z lguez $
 !
-c
-      SUBROUTINE dynredem0(fichnom,iday_end,phis)
+!
+SUBROUTINE dynredem0(fichnom,iday_end,phis)
 #ifdef CPP_IOIPSL
-      USE IOIPSL
+  USE IOIPSL
 #endif
-      USE infotrac
-      use netcdf95, only: NF95_PUT_VAR
-      use control_mod, only : planet_type
-      USE comvert_mod, ONLY: ap,bp,aps,bps,presnivs,pseudoalt,pa,preff,
-     .			nivsig,nivsigs
-      USE comconst_mod, ONLY: daysec,dtvr,rad,omeg,g,cpp,kappa,pi
-      USE logic_mod, ONLY: fxyhypb,ysinus
-      USE serre_mod, ONLY: clon,clat,grossismx,grossismy,dzoomx,dzoomy,
-     .			taux,tauy
-      USE temps_mod, ONLY: annee_ref,day_ref,itau_dyn,itaufin,
-     .			start_time,hour_ini
-      USE ener_mod, ONLY: etot0,ptot0,ztot0,stot0,ang0
-
-      IMPLICIT NONE
-c=======================================================================
-c Ecriture du fichier de redemarrage sous format NetCDF (initialisation)
-c=======================================================================
-c   Declarations:
-c   -------------
-#include "dimensions.h"
-#include "paramet.h"
-#include "comgeom2.h"
-#include "netcdf.inc"
-#include "iniprint.h"
-
-c   Arguments:
-c   ----------
-      INTEGER iday_end
-      REAL phis(iip1, jjp1)
-      CHARACTER*(*) fichnom
-
-c   Local:
-c   ------
-      INTEGER iq,l
-      INTEGER length
-      PARAMETER (length = 100)
-      REAL tab_cntrl(length) ! tableau des parametres du run
-      INTEGER ierr
-      character*20 modname
-      character*80 abort_message
-
-c   Variables locales pour NetCDF:
-c
-      INTEGER dims2(2), dims3(3), dims4(4)
-      INTEGER idim_index
-      INTEGER idim_rlonu, idim_rlonv, idim_rlatu, idim_rlatv
-      INTEGER idim_s, idim_sig
-      INTEGER idim_tim
-      INTEGER nid,nvarid
-
-      REAL zan0,zjulian,hours
-      INTEGER yyears0,jjour0, mmois0
-      character*30 unites
-
-      character(len=12) :: start_file_type="earth" ! default start file type
-      INTEGER idecal
-
-c-----------------------------------------------------------------------
-      modname='dynredem0'
+  USE infotrac, ONLY: nqtot, tname, ttext
+  USE netcdf, ONLY: NF90_CREATE, NF90_DEF_DIM, NF90_INQ_VARID, NF90_GLOBAL,    &
+                    NF90_CLOSE,  NF90_PUT_ATT, NF90_UNLIMITED, NF90_CLOBBER
+  USE dynredem_mod, ONLY: cre_var, put_var1, put_var2, err, modname, fil
+  use netcdf95, only: NF95_PUT_VAR
+  use control_mod, only : planet_type
+  USE comvert_mod, ONLY: ap,bp,aps,bps,presnivs,pseudoalt,pa,preff, &
+      			nivsig,nivsigs
+  USE comconst_mod, ONLY: daysec,dtvr,rad,omeg,g,cpp,kappa,pi
+  USE logic_mod, ONLY: fxyhypb,ysinus
+  USE serre_mod, ONLY: clon,clat,grossismx,grossismy,dzoomx,dzoomy, &
+      			taux,tauy
+  USE temps_mod, ONLY: annee_ref,day_ref,itau_dyn,itaufin, &
+      			start_time,hour_ini
+  USE ener_mod, ONLY: etot0,ptot0,ztot0,stot0,ang0
+
+  IMPLICIT NONE
+!=======================================================================
+! Writting the NetCDF restart file (initialisation)
+!=======================================================================
+!   Declarations:
+!   -------------
+  include "dimensions.h"
+  include "paramet.h"
+  include "comgeom2.h"
+  include "netcdf.inc"
+  include "iniprint.h"
+
+!===============================================================================
+! Arguments:
+  CHARACTER(LEN=*), INTENT(IN) :: fichnom          !--- FILE NAME
+  INTEGER,          INTENT(IN) :: iday_end         !--- 
+  REAL,             INTENT(IN) :: phis(iip1, jjp1) !--- GROUND GEOPOTENTIAL
+!===============================================================================
+! Local variables:
+  INTEGER :: iq,l
+  INTEGER, PARAMETER :: length=100
+  REAL :: tab_cntrl(length) ! run parameters
+  INTEGER :: ierr
+  CHARACTER(LEN=80) :: abort_message
+
+!   For NetCDF:
+  CHARACTER(LEN=30) :: unites
+  INTEGER :: indexID
+  INTEGER :: rlonuID, rlonvID, rlatuID, rlatvID
+  INTEGER :: sID, sigID, nID, vID, timID
+  INTEGER :: yyears0, jjour0, mmois0
+  REAL    :: zan0, zjulian, hours
+
+  CHARACTER(len=12) :: start_file_type="earth" ! default start file type
+  INTEGER :: idecal
+
+!===============================================================================
+  ! fill dynredem_mod module variables
+  modname='dynredem0'; fil=fichnom
 
 #ifdef CPP_IOIPSL
-      call ymds2ju(annee_ref, 1, iday_end, 0.0, zjulian)
-      call ju2ymds(zjulian, yyears0, mmois0, jjour0, hours)
+  call ymds2ju(annee_ref, 1, iday_end, 0.0, zjulian)
+  call ju2ymds(zjulian, yyears0, mmois0, jjour0, hours)
 #else
 ! set yyears0, mmois0, jjour0 to 0,1,1 (hours is not used)
-      yyears0=0
-      mmois0=1
-      jjour0=1
+  yyears0=0
+  mmois0=1
+  jjour0=1
 #endif        
 
-      !!! AS: idecal is a hack to be able to read planeto starts...
-      !!!     .... while keeping everything OK for LMDZ EARTH
-      if ((planet_type.eq."generic").or.(planet_type.eq."mars")) then
-          write(lunout,*) trim(modname),' : Planeto-like start file'
-          start_file_type="planeto"
-          idecal = 4
-      else
-          write(lunout,*) trim(modname),' : Earth-like start file'
-          idecal = 5
-      endif
-
-      DO l=1,length
-       tab_cntrl(l) = 0.
-      ENDDO
-       tab_cntrl(1)  = REAL(iim)
-       tab_cntrl(2)  = REAL(jjm)
-       tab_cntrl(3)  = REAL(llm)
-       if (start_file_type.eq."earth") then
-         tab_cntrl(4)=REAL(day_ref)
-       else
-         !tab_cntrl(4)=REAL(day_end)
-         tab_cntrl(4)=REAL(iday_end)
-       endif
-       tab_cntrl(5)  = REAL(annee_ref)
-       tab_cntrl(idecal+1)  = rad
-       tab_cntrl(idecal+2)  = omeg
-       tab_cntrl(idecal+3)  = g
-       tab_cntrl(idecal+4)  = cpp
-       tab_cntrl(idecal+5) = kappa
-       tab_cntrl(idecal+6) = daysec
-       tab_cntrl(idecal+7) = dtvr
-       tab_cntrl(idecal+8) = etot0
-       tab_cntrl(idecal+9) = ptot0
-       tab_cntrl(idecal+10) = ztot0
-       tab_cntrl(idecal+11) = stot0
-       tab_cntrl(idecal+12) = ang0
-       tab_cntrl(idecal+13) = pa
-       tab_cntrl(idecal+14) = preff
-c
-c    .....    parametres  pour le zoom      ......   
-
-       tab_cntrl(idecal+15)  = clon
-       tab_cntrl(idecal+16)  = clat
-       tab_cntrl(idecal+17)  = grossismx
-       tab_cntrl(idecal+18)  = grossismy
-c
-      IF ( fxyhypb )   THEN
-       tab_cntrl(idecal+19) = 1.
-       tab_cntrl(idecal+20) = dzoomx
-       tab_cntrl(idecal+21) = dzoomy
-       tab_cntrl(idecal+22) = 0.
-       tab_cntrl(idecal+23) = taux
-       tab_cntrl(idecal+24) = tauy
-      ELSE
-       tab_cntrl(idecal+19) = 0.
-       tab_cntrl(idecal+20) = dzoomx
-       tab_cntrl(idecal+21) = dzoomy
-       tab_cntrl(idecal+22) = 0.
-       tab_cntrl(idecal+23) = 0.
-       tab_cntrl(idecal+24) = 0.
-       IF( ysinus )  tab_cntrl(idecal+22) = 1.
-      ENDIF
-
-      if (start_file_type.eq."earth") then
-       tab_cntrl(idecal+25) = REAL(iday_end)
-       tab_cntrl(idecal+26) = REAL(itau_dyn + itaufin)
-c start_time: start_time of simulation (not necessarily 0.)
-       tab_cntrl(idecal+27) = start_time
-      endif
-
-      if (planet_type=="mars") then ! For Mars only
-        tab_cntrl(29)=hour_ini
-      endif
-c
-c    .........................................................
-c
-c Creation du fichier:
-c
-      ierr = NF_CREATE(fichnom, NF_CLOBBER, nid)
-      IF (ierr.NE.NF_NOERR) THEN
-         write(lunout,*)"dynredem0: Pb d ouverture du fichier "
-     &                  //trim(fichnom)
-         write(lunout,*)' ierr = ', ierr
-         CALL ABORT_GCM("DYNREDEM0", "", 1)
-      ENDIF
-c
-c Preciser quelques attributs globaux:
-c
-      ierr = NF_PUT_ATT_TEXT (nid, NF_GLOBAL, "title", 27,
-     .                       "Fichier demmarage dynamique")
-c
-c Definir les dimensions du fichiers:
-c
-      if (start_file_type.eq."earth") then
-        ierr = NF_DEF_DIM (nid, "index", length, idim_index)
-        ierr = NF_DEF_DIM (nid, "rlonu", iip1, idim_rlonu)
-        ierr = NF_DEF_DIM (nid, "rlatu", jjp1, idim_rlatu)
-        ierr = NF_DEF_DIM (nid, "rlonv", iip1, idim_rlonv)
-        ierr = NF_DEF_DIM (nid, "rlatv", jjm, idim_rlatv)
-        ierr = NF_DEF_DIM (nid, "sigs", llm, idim_s)
-        ierr = NF_DEF_DIM (nid, "sig", llmp1, idim_sig)
-        ierr = NF_DEF_DIM (nid, "temps", NF_UNLIMITED, idim_tim)
-      else
-        ierr = NF_DEF_DIM (nid, "index", length, idim_index)
-        ierr = NF_DEF_DIM (nid, "rlonu", iip1, idim_rlonu)
-        ierr = NF_DEF_DIM (nid, "latitude", jjp1, idim_rlatu)
-        ierr = NF_DEF_DIM (nid, "longitude", iip1, idim_rlonv)
-        ierr = NF_DEF_DIM (nid, "rlatv", jjm, idim_rlatv)
-        ierr = NF_DEF_DIM (nid, "altitude", llm, idim_s)
-        ierr = NF_DEF_DIM (nid, "interlayer", llmp1, idim_sig)
-        ierr = NF_DEF_DIM (nid, "Time", NF_UNLIMITED, idim_tim)
-      endif
-c
-      ierr = NF_ENDDEF(nid) ! sortir du mode de definition
-c
-c Definir et enregistrer certains champs invariants:
-c
-      ierr = NF_REDEF (nid)
-cIM 220306 BEG
-#ifdef NC_DOUBLE
-      ierr = NF_DEF_VAR (nid,"controle",NF_DOUBLE,1,idim_index,nvarid)
-#else
-      ierr = NF_DEF_VAR (nid,"controle",NF_FLOAT,1,idim_index,nvarid)
-#endif
-cIM 220306 END
-      ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 22,
-     .                       "Parametres de controle")
-      ierr = NF_ENDDEF(nid)
-      call NF95_PUT_VAR(nid,nvarid,tab_cntrl)
-c
-      ierr = NF_REDEF (nid)
-cIM 220306 BEG
-#ifdef NC_DOUBLE
-      ierr = NF_DEF_VAR (nid,"rlonu",NF_DOUBLE,1,idim_rlonu,nvarid)
-#else
-      ierr = NF_DEF_VAR (nid,"rlonu",NF_FLOAT,1,idim_rlonu,nvarid)
-#endif
-cIM 220306 END
-      ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 23,
-     .                       "Longitudes des points U")
-      ierr = NF_ENDDEF(nid)
-      call NF95_PUT_VAR(nid,nvarid,rlonu)
-c
-      ierr = NF_REDEF (nid)
-cIM 220306 BEG
-#ifdef NC_DOUBLE
-      ierr = NF_DEF_VAR (nid,"rlatu",NF_DOUBLE,1,idim_rlatu,nvarid)
-#else
-      ierr = NF_DEF_VAR (nid,"rlatu",NF_FLOAT,1,idim_rlatu,nvarid)
-#endif
-cIM 220306 END
-      ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 22,
-     .                       "Latitudes des points U")
-      ierr = NF_ENDDEF(nid)
-      call NF95_PUT_VAR (nid,nvarid,rlatu)
-c
-      ierr = NF_REDEF (nid)
-cIM 220306 BEG
-#ifdef NC_DOUBLE
-      ierr = NF_DEF_VAR (nid,"rlonv",NF_DOUBLE,1,idim_rlonv,nvarid)
-#else
-      ierr = NF_DEF_VAR (nid,"rlonv",NF_FLOAT,1,idim_rlonv,nvarid)
-#endif
-cIM 220306 END
-      ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 23,
-     .                       "Longitudes des points V")
-      ierr = NF_ENDDEF(nid)
-      call NF95_PUT_VAR(nid,nvarid,rlonv)
-c
-      ierr = NF_REDEF (nid)
-cIM 220306 BEG
-#ifdef NC_DOUBLE
-      ierr = NF_DEF_VAR (nid,"rlatv",NF_DOUBLE,1,idim_rlatv,nvarid)
-#else
-      ierr = NF_DEF_VAR (nid,"rlatv",NF_FLOAT,1,idim_rlatv,nvarid)
-#endif
-cIM 220306 END
-      ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 22,
-     .                       "Latitudes des points V")
-      ierr = NF_ENDDEF(nid)
-      call NF95_PUT_VAR(nid,nvarid,rlatv)
-c
-      if (start_file_type.eq."earth") then
-        ierr = NF_REDEF (nid)
-cIM 220306 BEG
-#ifdef NC_DOUBLE
-        ierr = NF_DEF_VAR (nid,"nivsigs",NF_DOUBLE,1,idim_s,nvarid)
-#else
-        ierr = NF_DEF_VAR (nid,"nivsigs",NF_FLOAT,1,idim_s,nvarid)
-#endif
-cIM 220306 END
-        ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 28,
-     .                       "Numero naturel des couches s")
-        ierr = NF_ENDDEF(nid)
-        call NF95_PUT_VAR(nid,nvarid,nivsigs)
-c
-        ierr = NF_REDEF (nid)
-cIM 220306 BEG
-#ifdef NC_DOUBLE
-        ierr = NF_DEF_VAR (nid,"nivsig",NF_DOUBLE,1,idim_sig,nvarid)
-#else
-        ierr = NF_DEF_VAR (nid,"nivsig",NF_FLOAT,1,idim_sig,nvarid)
-#endif
-cIM 220306 END
-        ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 32,
-     .                       "Numero naturel des couches sigma")
-        ierr = NF_ENDDEF(nid)
-        call NF95_PUT_VAR(nid,nvarid,nivsig)
-      endif ! of if (start_file_type.eq."earth")
-c
-      ierr = NF_REDEF (nid)
-cIM 220306 BEG
-#ifdef NC_DOUBLE
-      ierr = NF_DEF_VAR (nid,"ap",NF_DOUBLE,1,idim_sig,nvarid)
-#else
-      ierr = NF_DEF_VAR (nid,"ap",NF_FLOAT,1,idim_sig,nvarid)
-#endif
-cIM 220306 END
-      ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 26,
-     .                       "Coefficient A pour hybride")
-      ierr = NF_ENDDEF(nid)
-      call NF95_PUT_VAR(nid,nvarid,ap)
-c
-      ierr = NF_REDEF (nid)
-cIM 220306 BEG
-#ifdef NC_DOUBLE
-      ierr = NF_DEF_VAR (nid,"bp",NF_DOUBLE,1,idim_sig,nvarid)
-#else
-      ierr = NF_DEF_VAR (nid,"bp",NF_FLOAT,1,idim_sig,nvarid)
-#endif
-cIM 220306 END
-      ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 26,
-     .                       "Coefficient B pour hybride")
-      ierr = NF_ENDDEF(nid)
-      call NF95_PUT_VAR(nid,nvarid,bp)
-c
-      if (start_file_type.ne."earth") then
-        ierr = NF_REDEF (nid)
-cIM 220306 BEG
-#ifdef NC_DOUBLE
-        ierr = NF_DEF_VAR (nid,"aps",NF_DOUBLE,1,idim_s,nvarid)
-#else
-        ierr = NF_DEF_VAR (nid,"aps",NF_FLOAT,1,idim_s,nvarid)
-#endif
-cIM 220306 END
-        ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 37,
-     .                       "Coef AS: hybrid pressure at midlayers")
-        ierr = NF_ENDDEF(nid)
-        call NF95_PUT_VAR(nid,nvarid,aps)
-c
-        ierr = NF_REDEF (nid)
-cIM 220306 BEG
-#ifdef NC_DOUBLE
-        ierr = NF_DEF_VAR (nid,"bps",NF_DOUBLE,1,idim_s,nvarid)
-#else
-        ierr = NF_DEF_VAR (nid,"bps",NF_FLOAT,1,idim_s,nvarid)
-#endif
-cIM 220306 END
-        ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 34,
-     .                       "Coef BS: hybrid sigma at midlayers")
-        ierr = NF_ENDDEF(nid)
-        call NF95_PUT_VAR(nid,nvarid,bps)
-      endif ! of if (start_file_type.ne."earth")
-c
-      ierr = NF_REDEF (nid)
-cIM 220306 BEG
-#ifdef NC_DOUBLE
-      ierr = NF_DEF_VAR (nid,"presnivs",NF_DOUBLE,1,idim_s,nvarid)
-#else
-      ierr = NF_DEF_VAR (nid,"presnivs",NF_FLOAT,1,idim_s,nvarid)
-#endif
-cIM 220306 END
-      ierr = NF_ENDDEF(nid)
-      call NF95_PUT_VAR(nid,nvarid,presnivs)
-c
-      if (start_file_type.ne."earth") then
+  !!! AS: idecal is a hack to be able to read planeto starts...
+  !!!     .... while keeping everything OK for LMDZ EARTH
+  if ((planet_type.eq."generic").or.(planet_type.eq."mars")) then
+    write(lunout,*) trim(modname),' : Planeto-like start file'
+    start_file_type="planeto"
+    idecal = 4
+  else
+    write(lunout,*) trim(modname),' : Earth-like start file'
+    idecal = 5
+  endif
+
+  tab_cntrl(:)  = 0.
+  tab_cntrl(1)  = REAL(iim)
+  tab_cntrl(2)  = REAL(jjm)
+  tab_cntrl(3)  = REAL(llm)
+  if (start_file_type.eq."earth") then
+    tab_cntrl(4)=REAL(day_ref)
+  else
+    !tab_cntrl(4)=REAL(day_end)
+    tab_cntrl(4)=REAL(iday_end)
+  endif
+  tab_cntrl(5)  = REAL(annee_ref)
+  tab_cntrl(idecal+1)  = rad
+  tab_cntrl(idecal+2)  = omeg
+  tab_cntrl(idecal+3)  = g
+  tab_cntrl(idecal+4)  = cpp
+  tab_cntrl(idecal+5) = kappa
+  tab_cntrl(idecal+6) = daysec
+  tab_cntrl(idecal+7) = dtvr
+  tab_cntrl(idecal+8) = etot0
+  tab_cntrl(idecal+9) = ptot0
+  tab_cntrl(idecal+10) = ztot0
+  tab_cntrl(idecal+11) = stot0
+  tab_cntrl(idecal+12) = ang0
+  tab_cntrl(idecal+13) = pa
+  tab_cntrl(idecal+14) = preff
+
+!    .....    parameters for the zoom      ......   
+  tab_cntrl(idecal+15)  = clon
+  tab_cntrl(idecal+16)  = clat
+  tab_cntrl(idecal+17)  = grossismx
+  tab_cntrl(idecal+18)  = grossismy
+!
+  IF ( fxyhypb )   THEN
+    tab_cntrl(idecal+19) = 1.
+    tab_cntrl(idecal+20) = dzoomx
+    tab_cntrl(idecal+21) = dzoomy
+    tab_cntrl(idecal+22) = 0.
+    tab_cntrl(idecal+23) = taux
+    tab_cntrl(idecal+24) = tauy
+  ELSE
+    tab_cntrl(idecal+19) = 0.
+    tab_cntrl(idecal+20) = dzoomx
+    tab_cntrl(idecal+21) = dzoomy
+    tab_cntrl(idecal+22) = 0.
+    tab_cntrl(idecal+23) = 0.
+    tab_cntrl(idecal+24) = 0.
+    IF( ysinus )  tab_cntrl(idecal+22) = 1.
+  ENDIF
+
+  if (start_file_type.eq."earth") then
+    tab_cntrl(idecal+25) = REAL(iday_end)
+    tab_cntrl(idecal+26) = REAL(itau_dyn + itaufin)
+! start_time: start_time of simulation (not necessarily 0.)
+    tab_cntrl(idecal+27) = start_time
+  endif
+
+  if (planet_type=="mars") then ! For Mars only
+    tab_cntrl(29)=hour_ini
+  endif
+
+!--- File creation
+  CALL err(NF90_CREATE(fichnom,NF90_CLOBBER,nid))
+
+!--- Some global attributes
+  CALL err(NF90_PUT_ATT(nid,NF90_GLOBAL,"title","Fichier demarrage dynamique"))
+
+!--- Dimensions
+  if (start_file_type.eq."earth") then
+    CALL err(NF90_DEF_DIM(nid,"index", length, indexID))
+    CALL err(NF90_DEF_DIM(nid,"rlonu", iip1,   rlonuID))
+    CALL err(NF90_DEF_DIM(nid,"rlatu", jjp1,   rlatuID))
+    CALL err(NF90_DEF_DIM(nid,"rlonv", iip1,   rlonvID))
+    CALL err(NF90_DEF_DIM(nid,"rlatv", jjm,    rlatvID))
+    CALL err(NF90_DEF_DIM(nid,"sigs",  llm,        sID))
+    CALL err(NF90_DEF_DIM(nid,"sig",   llmp1,    sigID))
+    CALL err(NF90_DEF_DIM(nid,"temps", NF90_UNLIMITED, timID))
+  else
+    CALL err(NF90_DEF_DIM(nid,"index", length, indexID))
+    CALL err(NF90_DEF_DIM(nid,"rlonu", iip1,   rlonuID))
+    CALL err(NF90_DEF_DIM(nid,"latitude", jjp1,   rlatuID))
+    CALL err(NF90_DEF_DIM(nid,"longitude", iip1,   rlonvID))
+    CALL err(NF90_DEF_DIM(nid,"rlatv", jjm,    rlatvID))
+    CALL err(NF90_DEF_DIM(nid,"altitude",  llm,        sID))
+    CALL err(NF90_DEF_DIM(nid,"interlayer",   llmp1,    sigID))
+    CALL err(NF90_DEF_DIM(nid,"Time", NF90_UNLIMITED, timID))
+  endif
+
+!--- Define and save invariant fields
+  CALL put_var1(nid,"controle","Parametres de controle" ,[indexID],tab_cntrl)
+  CALL put_var1(nid,"rlonu"   ,"Longitudes des points U",[rlonuID],rlonu)
+  CALL put_var1(nid,"rlatu"   ,"Latitudes des points U" ,[rlatuID],rlatu)
+  CALL put_var1(nid,"rlonv"   ,"Longitudes des points V",[rlonvID],rlonv)
+  CALL put_var1(nid,"rlatv"   ,"Latitudes des points V" ,[rlatvID],rlatv)
+  if (start_file_type.eq."earth") then
+    CALL put_var1(nid,"nivsigs" ,"Numero naturel des couches s"    ,[sID]  ,nivsigs)
+    CALL put_var1(nid,"nivsig"  ,"Numero naturel des couches sigma",[sigID],nivsig)
+  endif ! of if (start_file_type.eq."earth")
+  CALL put_var1(nid,"ap"      ,"Coefficient A pour hybride"      ,[sigID],ap)
+  CALL put_var1(nid,"bp"      ,"Coefficient B pour hybride"      ,[sigID],bp)
+  if (start_file_type.ne."earth") then
+    CALL put_var1(nid,"aps","Coef AS: hybrid pressure at midlayers",[sID],aps)
+    CALL put_var1(nid,"bps","Coef BS: hybrid sigma at midlayers",[sID],bps)
+  endif ! of if (start_file_type.eq."earth")
+  CALL put_var1(nid,"presnivs",""                                ,[sID]  ,presnivs)
+  if (start_file_type.ne."earth") then
         ierr = NF_REDEF (nid)
 #ifdef NC_DOUBLE
-        ierr = NF_DEF_VAR(nid,"latitude",NF_DOUBLE,1,idim_rlatu,nvarid)
+        ierr = NF_DEF_VAR(nid,"latitude",NF_DOUBLE,1,rlatuID,vID)
 #else
-        ierr = NF_DEF_VAR(nid,"latitude",NF_FLOAT,1,idim_rlatu,nvarid)
+        ierr = NF_DEF_VAR(nid,"latitude",NF_FLOAT,1,rlatuID,vID)
 #endif
-        ierr =NF_PUT_ATT_TEXT(nid,nvarid,'units',13,"degrees_north")
-        ierr = NF_PUT_ATT_TEXT (nid,nvarid,"long_name", 14,
-     .        "North latitude")
+        ierr =NF_PUT_ATT_TEXT(nid,vID,'units',13,"degrees_north")
+        ierr = NF_PUT_ATT_TEXT (nid,vID,"long_name", 14, &
+              "North latitude")
         ierr = NF_ENDDEF(nid)
-        call NF95_PUT_VAR(nid,nvarid,rlatu*180/pi)
-c
+        call NF95_PUT_VAR(nid,vID,rlatu*180/pi)
+!
         ierr = NF_REDEF (nid)
 #ifdef NC_DOUBLE
-        ierr=NF_DEF_VAR(nid,"longitude",NF_DOUBLE,1,idim_rlonv,nvarid)
+        ierr=NF_DEF_VAR(nid,"longitude",NF_DOUBLE,1,rlonvID,vID)
 #else
-        ierr=NF_DEF_VAR(nid,"longitude",NF_FLOAT,1,idim_rlonv,nvarid)
+        ierr=NF_DEF_VAR(nid,"longitude",NF_FLOAT,1,rlonvID,vID)
 #endif
-        ierr = NF_PUT_ATT_TEXT (nid,nvarid,"long_name", 14,
-     .        "East longitude")
-        ierr = NF_PUT_ATT_TEXT(nid,nvarid,'units',12,"degrees_east")
+        ierr = NF_PUT_ATT_TEXT (nid,vID,"long_name", 14, &
+              "East longitude")
+        ierr = NF_PUT_ATT_TEXT(nid,vID,'units',12,"degrees_east")
         ierr = NF_ENDDEF(nid)
-        call NF95_PUT_VAR(nid,nvarid,rlonv*180/pi)
-c
+        call NF95_PUT_VAR(nid,vID,rlonv*180/pi)
+!
         ierr = NF_REDEF (nid)
 #ifdef NC_DOUBLE
-        ierr = NF_DEF_VAR (nid, "altitude", NF_DOUBLE, 1,
-     .       idim_s,nvarid)
+        ierr = NF_DEF_VAR (nid, "altitude", NF_DOUBLE, 1, &
+             sID,vID)
 #else
-        ierr = NF_DEF_VAR (nid, "altitude", NF_FLOAT, 1,
-     .       idim_s,nvarid)
+        ierr = NF_DEF_VAR (nid, "altitude", NF_FLOAT, 1, &
+             sID,vID)
 #endif
-        ierr = NF_PUT_ATT_TEXT(nid,nvarid,"long_name",10,"pseudo-alt")
-        ierr = NF_PUT_ATT_TEXT (nid,nvarid,'units',2,"km")
-        ierr = NF_PUT_ATT_TEXT (nid,nvarid,'positive',2,"up")
+        ierr = NF_PUT_ATT_TEXT(nid,vID,"long_name",10,"pseudo-alt")
+        ierr = NF_PUT_ATT_TEXT (nid,vID,'units',2,"km")
+        ierr = NF_PUT_ATT_TEXT (nid,vID,'positive',2,"up")
         ierr = NF_ENDDEF(nid)
-        call NF95_PUT_VAR(nid,nvarid,pseudoalt)
-      endif ! of if (start_file_type.ne."earth")
-c
-c Coefficients de passage cov. <-> contra. <--> naturel
-c
-      ierr = NF_REDEF (nid)
-      dims2(1) = idim_rlonu
-      dims2(2) = idim_rlatu
-cIM 220306 BEG
-#ifdef NC_DOUBLE
-      ierr = NF_DEF_VAR (nid,"cu",NF_DOUBLE,2,dims2,nvarid)
-#else
-      ierr = NF_DEF_VAR (nid,"cu",NF_FLOAT,2,dims2,nvarid)
-#endif
-cIM 220306 END
-      ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 29,
-     .                       "Coefficient de passage pour U")
-      ierr = NF_ENDDEF(nid)
-      call NF95_PUT_VAR(nid,nvarid,cu)
-c
-      ierr = NF_REDEF (nid)
-      dims2(1) = idim_rlonv
-      dims2(2) = idim_rlatv
-cIM 220306 BEG
-#ifdef NC_DOUBLE
-      ierr = NF_DEF_VAR (nid,"cv",NF_DOUBLE,2,dims2,nvarid)
-#else
-      ierr = NF_DEF_VAR (nid,"cv",NF_FLOAT,2,dims2,nvarid)
-#endif
-cIM 220306 END
-      ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 29,
-     .                       "Coefficient de passage pour V")
-      ierr = NF_ENDDEF(nid)
-      call NF95_PUT_VAR(nid,nvarid,cv)
-c
-c Aire de chaque maille:
-c
-      ierr = NF_REDEF (nid)
-      dims2(1) = idim_rlonv
-      dims2(2) = idim_rlatu
-cIM 220306 BEG
-#ifdef NC_DOUBLE
-      ierr = NF_DEF_VAR (nid,"aire",NF_DOUBLE,2,dims2,nvarid)
-#else
-      ierr = NF_DEF_VAR (nid,"aire",NF_FLOAT,2,dims2,nvarid)
-#endif
-cIM 220306 END
-      ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 22,
-     .                       "Aires de chaque maille")
-      ierr = NF_ENDDEF(nid)
-      call NF95_PUT_VAR(nid,nvarid,aire)
-c
-c Geopentiel au sol:
-c
-      ierr = NF_REDEF (nid)
-      dims2(1) = idim_rlonv
-      dims2(2) = idim_rlatu
-cIM 220306 BEG
-#ifdef NC_DOUBLE
-      ierr = NF_DEF_VAR (nid,"phisinit",NF_DOUBLE,2,dims2,nvarid)
-#else
-      ierr = NF_DEF_VAR (nid,"phisinit",NF_FLOAT,2,dims2,nvarid)
-#endif
-cIM 220306 END
-      ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 19,
-     .                       "Geopotentiel au sol")
-      ierr = NF_ENDDEF(nid)
-      call NF95_PUT_VAR(nid,nvarid,phis)
-c
-c Definir les variables pour pouvoir les enregistrer plus tard:
-c
-      ierr = NF_REDEF (nid) ! entrer dans le mode de definition
-c
-      if (start_file_type.eq."earth") then
-cIM 220306 BEG
-#ifdef NC_DOUBLE
-        ierr = NF_DEF_VAR (nid,"temps",NF_DOUBLE,1,idim_tim,nvarid)
-#else
-        ierr = NF_DEF_VAR (nid,"temps",NF_FLOAT,1,idim_tim,nvarid)
-#endif
-cIM 220306 END
-      else ! start_file_type=="planeto"
-#ifdef NC_DOUBLE
-        ierr = NF_DEF_VAR (nid,"Time",NF_DOUBLE,1,idim_tim,nvarid)
-#else
-        ierr = NF_DEF_VAR (nid,"Time",NF_FLOAT,1,idim_tim,nvarid)
-#endif
-      endif ! of if (start_file_type.eq."earth")
-      ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 19,
-     .                       "Temps de simulation")
-      write(unites,200)yyears0,mmois0,jjour0
-200   format('days since ',i4,'-',i2.2,'-',i2.2,' 00:00:00')
-      ierr = NF_PUT_ATT_TEXT (nid, nvarid, "units", 30,
-     .                         unites)
-
-c
-      dims4(1) = idim_rlonu
-      dims4(2) = idim_rlatu
-      dims4(3) = idim_s
-      dims4(4) = idim_tim
-cIM 220306 BEG
-#ifdef NC_DOUBLE
-      ierr = NF_DEF_VAR (nid,"ucov",NF_DOUBLE,4,dims4,nvarid)
-#else
-      ierr = NF_DEF_VAR (nid,"ucov",NF_FLOAT,4,dims4,nvarid)
-#endif
-cIM 220306 END
-      ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 9,
-     .                       "Vitesse U")
-c
-      dims4(1) = idim_rlonv
-      dims4(2) = idim_rlatv
-      dims4(3) = idim_s
-      dims4(4) = idim_tim
-cIM 220306 BEG
-#ifdef NC_DOUBLE
-      ierr = NF_DEF_VAR (nid,"vcov",NF_DOUBLE,4,dims4,nvarid)
-#else
-      ierr = NF_DEF_VAR (nid,"vcov",NF_FLOAT,4,dims4,nvarid)
-#endif
-cIM 220306 END
-      ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 9,
-     .                       "Vitesse V")
-c
-      dims4(1) = idim_rlonv
-      dims4(2) = idim_rlatu
-      dims4(3) = idim_s
-      dims4(4) = idim_tim
-cIM 220306 BEG
-#ifdef NC_DOUBLE
-      ierr = NF_DEF_VAR (nid,"teta",NF_DOUBLE,4,dims4,nvarid)
-#else
-      ierr = NF_DEF_VAR (nid,"teta",NF_FLOAT,4,dims4,nvarid)
-#endif
-cIM 220306 END
-      ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 11,
-     .                       "Temperature")
-c
-      dims4(1) = idim_rlonv
-      dims4(2) = idim_rlatu
-      dims4(3) = idim_s
-      dims4(4) = idim_tim
-      IF(nqtot.GE.1) THEN
-      DO iq=1,nqtot
-cIM 220306 BEG
-#ifdef NC_DOUBLE
-      ierr = NF_DEF_VAR (nid,tname(iq),NF_DOUBLE,4,dims4,nvarid)
-#else
-      ierr = NF_DEF_VAR (nid,tname(iq),NF_FLOAT,4,dims4,nvarid)
-#endif
-cIM 220306 END
-      ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 12,ttext(iq))
-      ENDDO
-      ENDIF
-c
-      dims4(1) = idim_rlonv
-      dims4(2) = idim_rlatu
-      dims4(3) = idim_s
-      dims4(4) = idim_tim
-cIM 220306 BEG
-#ifdef NC_DOUBLE
-      ierr = NF_DEF_VAR (nid,"masse",NF_DOUBLE,4,dims4,nvarid)
-#else
-      ierr = NF_DEF_VAR (nid,"masse",NF_FLOAT,4,dims4,nvarid)
-#endif
-cIM 220306 END
-      ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 12,
-     .                       "C est quoi ?")
-c
-      dims3(1) = idim_rlonv
-      dims3(2) = idim_rlatu
-      dims3(3) = idim_tim
-cIM 220306 BEG
-#ifdef NC_DOUBLE
-      ierr = NF_DEF_VAR (nid,"ps",NF_DOUBLE,3,dims3,nvarid)
-#else
-      ierr = NF_DEF_VAR (nid,"ps",NF_FLOAT,3,dims3,nvarid)
-#endif
-cIM 220306 END
-      ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 15,
-     .                       "Pression au sol")
-c
-      ierr = NF_ENDDEF(nid) ! sortir du mode de definition
-      ierr = NF_CLOSE(nid) ! fermer le fichier
-
-      write(lunout,*)'dynredem0: iim,jjm,llm,iday_end',
-     &               iim,jjm,llm,iday_end
-      write(lunout,*)'dynredem0: rad,omeg,g,cpp,kappa',
-     &        rad,omeg,g,cpp,kappa
-
-      RETURN
-      END
+        call NF95_PUT_VAR(nid,vID,pseudoalt)
+        CALL err(NF_REDEF(nid))
+  endif ! of if (start_file_type.ne."earth")
+
+! covariant <-> contravariant <-> natural conversion coefficients
+  CALL put_var2(nid,"cu","Coefficient de passage pour U",[rlonuID,rlatuID],cu)
+  CALL put_var2(nid,"cv","Coefficient de passage pour V",[rlonvID,rlatvID],cv)
+  CALL put_var2(nid,"aire","Aires de chaque maille"     ,[rlonvID,rlatuID],aire)
+  CALL put_var2(nid,"phisinit","Geopotentiel au sol"    ,[rlonvID,rlatuID],phis)
+
+
+! Define variables that will be stored later:
+  WRITE(unites,"('days since ',i4,'-',i2.2,'-',i2.2,' 00:00:00')"),&
+               yyears0,mmois0,jjour0
+  IF (planet_type.eq."earth") THEN
+    CALL cre_var(nid,"temps","Temps de simulation",[timID],unites)
+  ELSE
+    CALL cre_var(nid,"Time","Temps de simulation",[timID],unites)
+  ENDIF
+
+  CALL cre_var(nid,"ucov" ,"Vitesse U"  ,[rlonuID,rlatuID,sID,timID])
+  CALL cre_var(nid,"vcov" ,"Vitesse V"  ,[rlonvID,rlatvID,sID,timID])
+  CALL cre_var(nid,"teta" ,"Temperature",[rlonvID,rlatuID,sID,timID])
+
+  IF(nqtot.GE.1) THEN
+    DO iq=1,nqtot
+      CALL cre_var(nid,tname(iq),ttext(iq),[rlonvID,rlatuID,sID,timID])
+    END DO
+  ENDIF
+
+  CALL cre_var(nid,"masse","Masse d air"    ,[rlonvID,rlatuID,sID,timID])
+  CALL cre_var(nid,"ps"   ,"Pression au sol",[rlonvID,rlatuID    ,timID])
+
+  CALL err(NF90_CLOSE (nid)) ! close file
+
+  WRITE(lunout,*)TRIM(modname)//': iim,jjm,llm,iday_end',iim,jjm,llm,iday_end
+  WRITE(lunout,*)TRIM(modname)//': rad,omeg,g,cpp,kappa',rad,omeg,g,cpp,kappa
+
+END SUBROUTINE dynredem0
 
 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
 
-      SUBROUTINE dynredem1(fichnom,time,
-     .                     vcov,ucov,teta,q,masse,ps)
-      USE infotrac
-      USE control_mod, only : planet_type
-      use netcdf, only: NF90_get_VAR
-      use netcdf95, only: NF95_PUT_VAR
-      USE temps_mod, ONLY: itaufin,itau_dyn
+SUBROUTINE dynredem1(fichnom,time,vcov,ucov,teta,q,masse,ps)
+!
+!-------------------------------------------------------------------------------
+! Purpose: Write the NetCDF restart file (append).
+!-------------------------------------------------------------------------------
+  USE infotrac, ONLY: nqtot, tname, type_trac
+  USE control_mod, only : planet_type
+  USE netcdf,   ONLY: NF90_OPEN,  NF90_NOWRITE, NF90_GET_VAR, NF90_INQ_VARID,  &
+                      NF90_CLOSE, NF90_WRITE,   NF90_PUT_VAR, NF90_NoErr
+  use netcdf95, only: NF95_PUT_VAR
+  USE temps_mod, ONLY: itaufin,itau_dyn
+  USE dynredem_mod, ONLY: dynredem_write_u, dynredem_write_v, dynredem_read_u, &
+                          err, modname, fil, msg
  
-      IMPLICIT NONE
-c=================================================================
-c  Ecriture du fichier de redemarrage sous format NetCDF
-c=================================================================
-#include "dimensions.h"
-#include "paramet.h"
-#include "netcdf.inc"
-#include "comgeom.h"
-#include "iniprint.h"
-
-
-      INTEGER l
-      REAL vcov(iip1,jjm,llm),ucov(iip1, jjp1,llm) 
-      REAL teta(iip1, jjp1,llm)                   
-      REAL ps(iip1, jjp1),masse(iip1, jjp1,llm)                   
-      REAL q(iip1, jjp1, llm, nqtot)
-      CHARACTER*(*) fichnom
-     
-      REAL time
-      INTEGER nid, nvarid, nid_trac, nvarid_trac
-      REAL trac_tmp(ip1jmp1,llm)      
-      INTEGER ierr, ierr_file 
-      INTEGER iq
-      INTEGER length
-      PARAMETER (length = 100)
-      REAL tab_cntrl(length) ! tableau des parametres du run
-      character(len=*),parameter :: modname='dynredem1'
-      character*80 abort_message
-c
-      INTEGER nb
-      SAVE nb
-      DATA nb / 0 /
-      character(len=12) :: start_file_type="earth" ! default start file type
-
-      if ((planet_type.eq."generic").or.(planet_type.eq."mars")) then
-          write(lunout,*) trim(modname),' : Planeto-like start file'
-          start_file_type="planeto"
-      else
-          write(lunout,*) trim(modname),' : Earth-like start file'
-      endif 
-
-      ierr = NF_OPEN(fichnom, NF_WRITE, nid)
-      IF (ierr .NE. NF_NOERR) THEN
-         PRINT*, "dynredem1: Pb. d ouverture "//trim(fichnom)
-         CALL abort_gcm("dynredem1", "", 1)
-      ENDIF
-
-c  Ecriture/extension de la coordonnee temps
-
-      nb = nb + 1
-      if (start_file_type.eq."earth") then
-        ierr = NF_INQ_VARID(nid, "temps", nvarid)
+  IMPLICIT NONE
+  include "dimensions.h"
+  include "paramet.h"
+  include "netcdf.inc"
+  include "comgeom.h"
+  include "iniprint.h"
+!===============================================================================
+! Arguments:
+  CHARACTER(LEN=*), INTENT(IN) :: fichnom              !-- FILE NAME
+  REAL, INTENT(IN)    ::  time                         !-- TIME
+  REAL, INTENT(IN)    ::  vcov(iip1,jjm, llm)          !-- V COVARIANT WIND
+  REAL, INTENT(IN)    ::  ucov(iip1,jjp1,llm)          !-- U COVARIANT WIND
+  REAL, INTENT(IN)    ::  teta(iip1,jjp1,llm)          !-- POTENTIAL TEMPERATURE
+  REAL, INTENT(INOUT) ::     q(iip1,jjp1,llm,nqtot)    !-- TRACERS
+  REAL, INTENT(IN)    :: masse(iip1,jjp1,llm)          !-- MASS PER CELL
+  REAL, INTENT(IN)    ::    ps(iip1,jjp1)              !-- GROUND PRESSURE
+!===============================================================================
+! Local variables:
+  INTEGER :: l, iq, nid, vID, ierr, nid_trac, vID_trac
+  INTEGER,SAVE :: nb=0
+  INTEGER, PARAMETER :: length=100
+  REAL               :: tab_cntrl(length) ! tableau des parametres du run
+  CHARACTER(LEN=256) :: var, dum
+  LOGICAL            :: lread_inca
+  CHARACTER(LEN=80) :: abort_message
+  CHARACTER(len=12) :: start_file_type="earth" ! default start file type
+
+  ! fill dynredem_mod module variables
+  modname='dynredem1'; fil=fichnom
+
+  if ((planet_type.eq."generic").or.(planet_type.eq."mars")) then
+      write(lunout,*) trim(modname),' : Planeto-like start file'
+      start_file_type="planeto"
+  else
+      write(lunout,*) trim(modname),' : Earth-like start file'
+  endif 
+
+  CALL err(NF90_OPEN(fil,NF90_WRITE,nid),"open",fil)
+
+!--- Write/extend time coordinate
+  nb = nb + 1
+  if (start_file_type.eq."earth") then
+        ierr = NF_INQ_VARID(nid, "temps", vID)
         IF (ierr .NE. NF_NOERR) THEN
           write(lunout,*) NF_STRERROR(ierr)
@@ -649,6 +329,6 @@
           CALL abort_gcm(modname,abort_message,ierr)
         ENDIF
-      else
-        ierr = NF_INQ_VARID(nid,"Time", nvarid)
+ else
+        ierr = NF_INQ_VARID(nid,"Time", vID)
         IF (ierr .NE. NF_NOERR) THEN
           write(lunout,*) NF_STRERROR(ierr)
@@ -656,145 +336,54 @@
           CALL abort_gcm(modname,abort_message,ierr)
         ENDIF
-      endif ! of if (start_file_type.eq."earth") 
-      call NF95_PUT_VAR(nid,nvarid,time,start=(/nb/))
-      write(lunout,*) "dynredem1: Enregistrement pour ", nb, time
-
-c
-c  Re-ecriture du tableau de controle, itaufin n'est plus defini quand
-c  on passe dans dynredem0
-      ierr = NF_INQ_VARID (nid, "controle", nvarid)
-      IF (ierr .NE. NF_NOERR) THEN
-         abort_message="dynredem1: Le champ <controle> est absent"
-         ierr = 1
-         CALL abort_gcm(modname,abort_message,ierr)
-      ENDIF
-      ierr = NF90_GET_VAR(nid, nvarid, tab_cntrl)
-      if (start_file_type=="earth") then
-        tab_cntrl(31) = REAL(itau_dyn + itaufin)
-      else
-        tab_cntrl(31) = 0
-      endif
-      call NF95_PUT_VAR(nid,nvarid,tab_cntrl)
-
-c  Ecriture des champs
-c
-      ierr = NF_INQ_VARID(nid, "ucov", nvarid)
-      IF (ierr .NE. NF_NOERR) THEN
-         abort_message="Variable ucov n est pas definie"
-         ierr=1
-         CALL abort_gcm(modname,abort_message,ierr)
-      ENDIF
-      call NF95_PUT_VAR(nid,nvarid,ucov,start=(/1,1,1,nb/))
-
-      ierr = NF_INQ_VARID(nid, "vcov", nvarid)
-      IF (ierr .NE. NF_NOERR) THEN
-         abort_message="Variable vcov n est pas definie"
-         ierr=1
-         CALL abort_gcm(modname,abort_message,ierr)
-      ENDIF
-      call NF95_PUT_VAR(nid,nvarid,vcov,start=(/1,1,1,nb/))
-
-      ierr = NF_INQ_VARID(nid, "teta", nvarid)
-      IF (ierr .NE. NF_NOERR) THEN
-         abort_message="Variable teta n est pas definie"
-         ierr=1
-         CALL abort_gcm(modname,abort_message,ierr)
-      ENDIF
-      call NF95_PUT_VAR(nid,nvarid,teta,start=(/1,1,1,nb/))
-
-      IF (type_trac == 'inca') THEN
-! Ajout Anne pour lecture valeurs traceurs dans un fichier start_trac.nc
-         ierr_file = NF_OPEN ("start_trac.nc", NF_NOWRITE,nid_trac)
-         IF (ierr_file .NE.NF_NOERR) THEN
-            write(lunout,*)'dynredem1: Pb d''ouverture du fichier',
-     &                     ' start_trac.nc'
-            write(lunout,*)' ierr = ', ierr_file 
-         ENDIF
-      END IF
-
-      IF(nqtot.GE.1) THEN
-      do iq=1,nqtot 
-
-         IF (type_trac /= 'inca') THEN
-            ierr = NF_INQ_VARID(nid, tname(iq), nvarid)
-            IF (ierr .NE. NF_NOERR) THEN
-               abort_message="Variable  tname(iq) n est pas definie"
-               ierr=1
-               CALL abort_gcm(modname,abort_message,ierr)
-            ENDIF
-            call NF95_PUT_VAR(nid,nvarid,q(:,:,:,iq),start=(/1,1,1,nb/))
-        ELSE ! type_trac = inca
-! lecture de la valeur du traceur dans start_trac.nc
-           IF (ierr_file .ne. 2) THEN
-             ierr = NF_INQ_VARID (nid_trac, tname(iq), nvarid_trac)
-             IF (ierr .NE. NF_NOERR) THEN
-                write(lunout,*) "dynredem1: ",trim(tname(iq)),
-     &                          " est absent de start_trac.nc"
-                ierr = NF_INQ_VARID(nid, tname(iq), nvarid)
-                IF (ierr .NE. NF_NOERR) THEN
-                   abort_message="dynredem1: Variable "//
-     &                     trim(tname(iq))//" n est pas definie"
-                   ierr=1
-                   CALL abort_gcm(modname,abort_message,ierr)
-                ENDIF
-                call NF95_PUT_VAR(nid,nvarid,q(:,:,:,iq))
-                
-             ELSE
-                write(lunout,*) "dynredem1: ",trim(tname(iq)),
-     &              " est present dans start_trac.nc"
-               ierr = NF90_GET_VAR(nid_trac, nvarid_trac, trac_tmp)
-                IF (ierr .NE. NF_NOERR) THEN
-                   abort_message="dynredem1: Lecture echouee pour"//
-     &                    trim(tname(iq))
-                   ierr=1
-                   CALL abort_gcm(modname,abort_message,ierr)
-                ENDIF
-                ierr = NF_INQ_VARID(nid, tname(iq), nvarid)
-                IF (ierr .NE. NF_NOERR) THEN
-                   abort_message="dynredem1: Variable "//
-     &                trim(tname(iq))//" n est pas definie"
-                   ierr=1
-                   CALL abort_gcm(modname,abort_message,ierr)
-                ENDIF
-                call NF95_PUT_VAR(nid, nvarid, trac_tmp)
-               
-             ENDIF ! IF (ierr .NE. NF_NOERR)
-! fin lecture du traceur
-          ELSE                  ! si il n'y a pas de fichier start_trac.nc
-!             print *, 'il n y a pas de fichier start_trac'
-             ierr = NF_INQ_VARID(nid, tname(iq), nvarid)
-             IF (ierr .NE. NF_NOERR) THEN
-                abort_message="dynredem1: Variable "//
-     &                trim(tname(iq))//" n est pas definie"
-                   ierr=1
-                   CALL abort_gcm(modname,abort_message,ierr)
-             ENDIF
-             call NF95_PUT_VAR(nid,nvarid,q(:,:,:,iq),
-     &                  start=(/1,1,1,nb/))
-          ENDIF ! (ierr_file .ne. 2)
-       END IF   !type_trac
-      
-      ENDDO
-      ENDIF
-c
-      ierr = NF_INQ_VARID(nid, "masse", nvarid)
-      IF (ierr .NE. NF_NOERR) THEN
-         abort_message="dynredem1: Variable masse n est pas definie"
-         ierr=1
-         CALL abort_gcm(modname,abort_message,ierr)
-      ENDIF
-      call NF95_PUT_VAR(nid,nvarid,masse,start=(/1,1,1,nb/))
-c
-      ierr = NF_INQ_VARID(nid, "ps", nvarid)
-      IF (ierr .NE. NF_NOERR) THEN
-         abort_message="dynredem1: Variable ps n est pas definie"
-         ierr=1
-         CALL abort_gcm(modname,abort_message,ierr)
-      ENDIF
-      call NF95_PUT_VAR(nid,nvarid,ps,start=(/1,1,nb/))
-
-      ierr = NF_CLOSE(nid)
-c
-      RETURN
-      END
-
+  endif ! of if (start_file_type.eq."earth") 
+  call NF95_PUT_VAR(nid,vID,time,start=(/nb/))
+  WRITE(lunout,*)TRIM(modname)//": Saving for ", nb, time
+
+!--- Rewrite control table (itaufin undefined in dynredem0)
+  var="controle"
+  CALL err(NF90_INQ_VARID(nid,var,vID),"inq",var)
+  CALL err(NF90_GET_VAR(nid,vID,tab_cntrl),"get",var)
+  if (start_file_type=="earth") then
+    tab_cntrl(31) = REAL(itau_dyn + itaufin)
+  else
+    tab_cntrl(31) = 0
+  endif
+  CALL err(NF90_INQ_VARID(nid,var,vID),"inq",var)
+  CALL err(NF90_PUT_VAR(nid,vID,tab_cntrl),"put",var)
+
+!--- Save fields
+  CALL dynredem_write_u(nid,"ucov" ,ucov ,llm, nb)
+  CALL dynredem_write_v(nid,"vcov" ,vcov ,llm, nb)
+  CALL dynredem_write_u(nid,"teta" ,teta ,llm, nb)
+  CALL dynredem_write_u(nid,"masse",masse,llm, nb)
+  CALL dynredem_write_u(nid,"ps"   ,ps   ,1, nb)
+
+!--- Tracers in file "start_trac.nc" (added by Anne)
+  lread_inca=.FALSE.; fil="start_trac.nc"
+  IF(type_trac=='inca') INQUIRE(FILE=fil,EXIST=lread_inca)
+  IF(lread_inca) CALL err(NF90_OPEN(fil,NF90_NOWRITE,nid_trac),"open")
+
+!--- Save tracers
+  IF(nqtot.GE.1) THEN
+    DO iq=1,nqtot
+      var=tname(iq); ierr=-1
+      IF(lread_inca) THEN                  !--- Possibly read from "start_trac.nc"
+        fil="start_trac.nc"
+        ierr=NF90_INQ_VARID(nid_trac,var,vID_trac)
+        dum='inq'; IF(ierr==NF90_NoErr) dum='fnd'
+        WRITE(lunout,*)msg(dum,var)
+
+
+        IF(ierr==NF90_NoErr) CALL dynredem_read_u(nid_trac,var,q(:,:,:,iq),llm)
+      END IF ! of IF(lread_inca)
+      fil=fichnom
+      CALL dynredem_write_u(nid,var,q(:,:,:,iq),llm,nb)
+    END DO ! of DO iq=1,nqtot
+  ENDIF ! of IF(nqtot.GE.1)
+
+  CALL err(NF90_CLOSE(nid),"close")
+  fil="start_trac.nc"
+  IF(lread_inca) CALL err(NF90_CLOSE(nid_trac),"close")
+
+END SUBROUTINE dynredem1
+
