Index: LMDZ6/branches/Ocean_skin/libf/dyn3d_common/comvert_mod.F90
===================================================================
--- LMDZ6/branches/Ocean_skin/libf/dyn3d_common/comvert_mod.F90	(revision 4013)
+++ LMDZ6/branches/Ocean_skin/libf/dyn3d_common/comvert_mod.F90	(revision 4368)
@@ -10,9 +10,11 @@
 
 PUBLIC :: ap,bp,presnivs,dpres,sig,ds,pa,preff,nivsigs,nivsig, &
-          aps,bps,scaleheight,pseudoalt,disvert_type, pressure_exner
+          aps,bps,scaleheight,pseudoalt,disvert_type, pressure_exner, &
+          presinter
 
 REAL ap(llm+1) ! hybrid pressure contribution at interlayers
 REAL bp (llm+1) ! hybrid sigma contribution at interlayer
 REAL presnivs(llm) ! (reference) pressure at mid-layers
+REAL presinter(llm+1) ! (reference) pressure at interlayers
 REAL dpres(llm)
 REAL sig(llm+1)
Index: LMDZ6/branches/Ocean_skin/libf/dyn3d_common/control_mod.F90
===================================================================
--- LMDZ6/branches/Ocean_skin/libf/dyn3d_common/control_mod.F90	(revision 4013)
+++ LMDZ6/branches/Ocean_skin/libf/dyn3d_common/control_mod.F90	(revision 4368)
@@ -29,5 +29,4 @@
   INTEGER,SAVE :: ip_ebil_dyn
   LOGICAL,SAVE :: offline
-  CHARACTER(len=4),SAVE :: config_inca
   CHARACTER(len=10),SAVE :: planet_type ! planet type ('earth','mars',...)
   LOGICAL,SAVE :: output_grads_dyn ! output dynamics diagnostics in
@@ -38,4 +37,5 @@
   LOGICAL,SAVE :: ok_dyn_ave ! output averaged values of fields in the dynamics
                              ! in NetCDF files dyn_hist*ave.nc
+  LOGICAL,SAVE :: ok_dyn_xios ! xios outputs in dynamics
   LOGICAL,SAVE :: resetvarc  ! allows to reset the variables in sortvarc
 
Index: LMDZ6/branches/Ocean_skin/libf/dyn3d_common/disvert.F90
===================================================================
--- LMDZ6/branches/Ocean_skin/libf/dyn3d_common/disvert.F90	(revision 4013)
+++ LMDZ6/branches/Ocean_skin/libf/dyn3d_common/disvert.F90	(revision 4368)
@@ -11,5 +11,5 @@
   use assert_m, only: assert
   USE comvert_mod, ONLY: ap, bp, aps, bps, nivsigs, nivsig, dpres, presnivs, &
-                         pseudoalt, pa, preff, scaleheight
+                         pseudoalt, pa, preff, scaleheight, presinter
   USE logic_mod, ONLY: ok_strato
 
@@ -35,4 +35,5 @@
 ! dpres(llm)                 !--- PRESSURE DIFFERENCE FOR EACH LAYER
 ! presnivs(llm)              !--- PRESSURE AT EACH MID-LAYER
+! presinter(llm+1)           !--- PRESSURE AT EACH INTERFACE
 ! scaleheight                !--- VERTICAL SCALE HEIGHT            (Earth: 8kms)
 ! nivsig(llm+1)              !--- SIGMA INDEX OF EACH LAYER INTERFACE
@@ -355,4 +356,8 @@
           max(ap(l+1)+bp(l+1)*preff, 1.e-10))
   ENDDO
+  DO l=1, llmp1
+     presinter(l)= ( ap(l)+bp(l)*preff)
+     write(lunout, *)'PRESINTER(', l, ')=', presinter(l)
+  ENDDO
 
   write(lunout, *) trim(modname),': PRESNIVS '
@@ -395,5 +400,5 @@
       ELSE IF(SIGN(f1,f4)/=f1) THEN; x2=xx; f2=f4
       ELSE IF(SIGN(f2,f4)/=f2) THEN; x1=xx; f1=f4
-      ELSE; CALL abort_gcm("ridders",'Algorithm failed (which is odd...')
+      ELSE; CALL abort_gcm("ridders",'Algorithm failed (which is odd...', 1)
       END IF
       IF(ABS(10.*LOG(ABS(x2-x1)))<=1E-5) EXIT       !--- ERROR ON SIG <= 0.01m            
Index: LMDZ6/branches/Ocean_skin/libf/dyn3d_common/grilles_gcm_netcdf_sub.F90
===================================================================
--- LMDZ6/branches/Ocean_skin/libf/dyn3d_common/grilles_gcm_netcdf_sub.F90	(revision 4013)
+++ LMDZ6/branches/Ocean_skin/libf/dyn3d_common/grilles_gcm_netcdf_sub.F90	(revision 4368)
@@ -2,13 +2,17 @@
 ! $Id: $
 !
-! This subroutine creates the file grilles_gcm.nc containg longitudes and
-! latitudes in degrees for grid u and v. This subroutine is called from
-! ce0l. This subroutine corresponds to the first 
-! part in the program create_fausse_var.
+! This subroutine creates the grilles_gcm.nc file, containing:
+! -> longitudes and latitudes in degrees for dynamical grids u, v and scalaire,
+! and the following variables added for INCA (informative anyway)
+! -> vertical levels "presnivs"
+! -> mask (land/sea), area (grid), phis=surface geopotential height = phis/g
 !
+! The subroutine is called in dynphy_lonlat/phylmd/ce0l.F90.
+
 SUBROUTINE grilles_gcm_netcdf_sub(masque,phis)
 
   USE comconst_mod, ONLY: cpp, kappa, g, omeg, daysec, rad, pi
   USE comvert_mod, ONLY: presnivs, preff, pa
+  use netcdf, only: nf90_def_var, nf90_int, nf90_float, nf90_put_var
   
   IMPLICIT NONE
@@ -19,32 +23,35 @@
   INCLUDE "netcdf.inc"
 
-
+!======================== 
   REAL,DIMENSION(iip1,jjp1),INTENT(IN)  :: masque ! masque terre/mer
   REAL,DIMENSION(iip1,jjp1),INTENT(IN)  :: phis   ! geopotentiel au sol
 
-  REAL temp(iim+1,jjm+1)
-  ! Attributs netcdf sortie
+  INTEGER status,i,j
+  
+  ! Attributs netcdf output
   INTEGER ncid_out,rcode_out
-  INTEGER out_lonuid,out_lonvid,out_latuid,out_latvid,out_levid
-  INTEGER out_varid
+  
+  INTEGER out_lonuid,out_lonvid,out_latuid,out_latvid
+  INTEGER out_uid,out_vid,out_tempid
   INTEGER out_lonudim,out_lonvdim
-  INTEGER out_latudim,out_latvdim,out_dim(3)
+  INTEGER out_latudim,out_latvdim,out_dim(2)
   INTEGER out_levdim
-
-  INTEGER start(4),COUNT(4)
-
-  INTEGER status,i,j
-  REAL rlatudeg(jjp1),rlatvdeg(jjm),rlevdeg(llm)
+  !
+  INTEGER :: presnivs_id
+  INTEGER :: mask_id,area_id,phis_id
+  !
+  INTEGER start(2),COUNT(2)
+
+  ! Variables
+  REAL rlatudeg(jjp1),rlatvdeg(jjm),rlev(llm)
   REAL rlonudeg(iip1),rlonvdeg(iip1)
-
-  REAL dlon1(iip1),dlon2(iip1),dlat1(jjp1),dlat2(jjp1)
-  REAL acoslat,dxkm,dykm,resol(iip1,jjp1)
-  REAL,DIMENSION(iip1,jjp1)  :: phis_loc
+  REAL uwnd(iip1,jjp1),vwnd(iip1,jjm),temp(iip1,jjp1) 
+  !
   INTEGER masque_int(iip1,jjp1)
-  INTEGER :: phis_id
-  INTEGER :: area_id
-  INTEGER :: mask_id
-  INTEGER :: presnivs_id
-  
+  REAL :: phis_loc(iip1,jjp1)
+  
+  !======================== 
+  ! CALCULATION of latu, latv, lonu, lonv in deg.
+  ! ---------------------------------------------------
   rad = 6400000
   omeg = 7.272205e-05
@@ -64,4 +71,5 @@
      rlatudeg(j)=rlatu(j)*180./pi
   ENDDO
+  
   DO j=1,jjm
      rlatvdeg(j)=rlatv(j)*180./pi
@@ -72,12 +80,31 @@
      rlonvdeg(i)=rlonv(i)*180./pi + 360.
   ENDDO
-
-
-  !  2 ----- OUVERTURE DE LA SORTIE NETCDF
-  ! ---------------------------------------------------
-  ! CREATION OUTPUT
-  ! ouverture fichier netcdf de sortie out
+ 
+  ! CALCULATION of "false" variables on u, v, s grids
+  ! --------------------------------------------------- 
+   DO i=1,iip1
+     DO j=1,jjp1
+        uwnd(i,j)=MOD(i,2)+MOD(j,2)
+        temp(i,j)=MOD(i,2)+MOD(j,2)
+     ENDDO
+     DO j=1,jjm
+        vwnd(i,j)=MOD(i,2)+MOD(j,2)
+     END DO
+  ENDDO  
+
+  ! CALCULATION of local vars for presnivs, mask, sfc. geopot. height
+  ! --------------------------------------------------- 
+  rlev(:) = presnivs(:)
+  phis_loc(:,:) = phis(:,:)/g
+  masque_int(:,:) = nINT(masque(:,:))
+
+
+  ! OPEN output netcdf file
+  !-------------------------
   status=NF_CREATE('grilles_gcm.nc',IOR(NF_CLOBBER,NF_64BIT_OFFSET),ncid_out)
   CALL handle_err(status)
+  
+  ! DEFINE output dimensions
+  !-------------------------
   status=NF_DEF_DIM(ncid_out,'lonu',iim+1,out_lonudim)
   CALL handle_err(status)
@@ -88,166 +115,118 @@
   status=NF_DEF_DIM(ncid_out,'latv',jjm,out_latvdim)
   CALL handle_err(status)
-
-
-  !   Longitudes en u
-  status=NF_DEF_VAR(ncid_out,'lonu',NF_FLOAT,1,out_lonudim, out_lonuid)
+  !
+  status=NF_DEF_DIM(ncid_out,'lev',llm,out_levdim)
+  CALL handle_err(status)
+  
+  ! DEFINE output variables
+  !-------------------------
+  !   Longitudes on "u" dynamical grid
+  status=NF90_DEF_VAR(ncid_out,'lonu',NF90_FLOAT,out_lonudim, out_lonuid)
   CALL handle_err(status)
   status=NF_PUT_ATT_TEXT(ncid_out,out_lonuid,'units', 12,'degrees_east')
-  status=NF_PUT_ATT_TEXT(ncid_out,out_lonuid,'long_name',9,'Longitude en u')
-
-  !   Longitudes en v
-  status=NF_DEF_VAR(ncid_out,'lonv',NF_FLOAT,1,out_lonvdim, out_lonvid)
+  status=NF_PUT_ATT_TEXT(ncid_out,out_lonuid,'long_name',19,'Longitude on u grid')
+  !   Longitudes on "v" dynamical grid
+  status=NF90_DEF_VAR(ncid_out,'lonv',NF90_FLOAT,out_lonvdim, out_lonvid)
   CALL handle_err(status)
   status=NF_PUT_ATT_TEXT(ncid_out,out_lonvid,'units', 12,'degrees_east')
-  status=NF_PUT_ATT_TEXT(ncid_out,out_lonvid,'long_name', 9,'Longitude en v')
-
-  !   Latitude en u
-  status=NF_DEF_VAR(ncid_out,'latu',NF_FLOAT,1,out_latudim, out_latuid)
+  status=NF_PUT_ATT_TEXT(ncid_out,out_lonvid,'long_name', 19,'Longitude on v grid')
+  !   Latitudes on "u" dynamical grid
+  status=NF90_DEF_VAR(ncid_out,'latu',NF90_FLOAT,out_latudim, out_latuid)
   CALL handle_err(status)
   status=NF_PUT_ATT_TEXT(ncid_out,out_latuid,'units', 13,'degrees_north')
-  status=NF_PUT_ATT_TEXT(ncid_out,out_latuid,'long_name', 8,'Latitude en u')
-
-  !  Latitude en v
-  status=NF_DEF_VAR(ncid_out,'latv',NF_FLOAT,1,out_latvdim, out_latvid)
+  status=NF_PUT_ATT_TEXT(ncid_out,out_latuid,'long_name', 18,'Latitude on u grid')
+  !  Latitudes on "v" dynamical grid
+  status=NF90_DEF_VAR(ncid_out,'latv',NF90_FLOAT,out_latvdim, out_latvid)
   CALL handle_err(status)
   status=NF_PUT_ATT_TEXT(ncid_out,out_latvid,'units', 13,'degrees_north')
-  status=NF_PUT_ATT_TEXT(ncid_out,out_latvid,'long_name', 8,'Latitude en v')
-
-  !   ecriture de la grille u
+  status=NF_PUT_ATT_TEXT(ncid_out,out_latvid,'long_name', 18,'Latitude on v grid')
+  !  "u" lat/lon dynamical grid
   out_dim(1)=out_lonudim
   out_dim(2)=out_latudim
-  status=NF_DEF_VAR(ncid_out,'grille_u',NF_FLOAT,2,out_dim, out_varid)
-  CALL handle_err(status)
-  status=NF_PUT_ATT_TEXT(ncid_out,out_varid,'units', 6,'Kelvin')
-  status=NF_PUT_ATT_TEXT(ncid_out,out_varid,'long_name', 16,'Grille aux point u')
-
-  !   ecriture de la grille v
+  status=NF90_DEF_VAR(ncid_out,'grille_u',NF90_FLOAT,out_dim, out_uid)
+  CALL handle_err(status)
+  status=NF_PUT_ATT_TEXT(ncid_out,out_uid,'units', 3,'m/s')
+  status=NF_PUT_ATT_TEXT(ncid_out,out_uid,'long_name', 21,'u-wind dynamical grid')
+  !  "v" lat/lon dynamical grid
   out_dim(1)=out_lonvdim
   out_dim(2)=out_latvdim
-  status=NF_DEF_VAR(ncid_out,'grille_v',NF_FLOAT,2,out_dim, out_varid)
-  CALL handle_err(status)
-  status=NF_PUT_ATT_TEXT(ncid_out,out_varid,'units', 6,'Kelvin')
-  status=NF_PUT_ATT_TEXT(ncid_out,out_varid,'long_name', 16,'Grille aux point v')
-
-  !   ecriture de la grille u
+  status=NF90_DEF_VAR(ncid_out,'grille_v',NF90_FLOAT,out_dim, out_vid)
+  CALL handle_err(status)
+  status=NF_PUT_ATT_TEXT(ncid_out,out_vid,'units', 3,'m/s')
+  status=NF_PUT_ATT_TEXT(ncid_out,out_vid,'long_name', 21,'v-wind dynamical grid')
+  !  "s" (scalar) lat/lon dynamical grid
   out_dim(1)=out_lonvdim
   out_dim(2)=out_latudim
-  status=NF_DEF_VAR(ncid_out,'grille_s',NF_FLOAT,2,out_dim, out_varid)
-  CALL handle_err(status)
-  status=NF_PUT_ATT_TEXT(ncid_out,out_varid,'units', 6,'Kelvin')
-  status=NF_PUT_ATT_TEXT(ncid_out,out_varid,'long_name',16,'Grille aux point u')
-
-  status=NF_ENDDEF(ncid_out)
-  write(*,*) "COUCOU 6"
-  CALL handle_err(status)
-  ! 5) ----- FERMETURE DES FICHIERS NETCDF------------------
-  ! --------------------------------------------------------
-  ! 3-b- Ecriture de la grille pour la sortie
-  ! rajoute l'ecriture de la grille
-
-#ifdef NC_DOUBLE
-  status=NF_PUT_VARA_DOUBLE(ncid_out,out_lonuid,1,iim+1,rlonudeg)
-  CALL handle_err(status)
-  status=NF_PUT_VARA_DOUBLE(ncid_out,out_lonvid,1,iim+1,rlonvdeg)
-  CALL handle_err(status)
-  status=NF_PUT_VARA_DOUBLE(ncid_out,out_latuid,1,jjm+1,rlatudeg)
-  CALL handle_err(status)
-  status=NF_PUT_VARA_DOUBLE(ncid_out,out_latvid,1,jjm,rlatvdeg)
-  CALL handle_err(status)
-#else
-  status=NF_PUT_VARA_REAL(ncid_out,out_lonuid,1,iim+1,rlonudeg)
-  CALL handle_err(status)
-  status=NF_PUT_VARA_REAL(ncid_out,out_lonvid,1,iim+1,rlonvdeg)
-  CALL handle_err(status)
-  status=NF_PUT_VARA_REAL(ncid_out,out_latuid,1,jjm+1,rlatudeg)
-  CALL handle_err(status)
-  status=NF_PUT_VARA_REAL(ncid_out,out_latvid,1,jjm,rlatvdeg)
-  CALL handle_err(status)
-#endif
-
-  start(1)=1
-  start(2)=1
-  start(3)=1
-  start(4)=1
-
-  COUNT(1)=iim+1
-  COUNT(2)=jjm+1
-  COUNT(3)=1
-  COUNT(4)=1
-
-  DO j=1,jjm+1
-     DO i=1,iim+1
-        temp(i,j)=MOD(i,2)+MOD(j,2)
-     ENDDO
-  ENDDO
-
-#ifdef NC_DOUBLE
-  status=NF_PUT_VARA_DOUBLE(ncid_out,out_varid,start, count,temp)
-  CALL handle_err(status)
-#else
-  status=NF_PUT_VARA_REAL(ncid_out,out_varid,start, count,temp)
-  CALL handle_err(status)
-#endif
-
-  ! On re-ouvre le fichier pour rajouter 4 nouvelles variables necessaire pour INCA
-! lev - phis - aire - mask
-! rlevdeg(:) = presnivs
-  rlevdeg(:) = presnivs(:)
-  phis_loc(:,:) = phis(:,:)/g
-
-! niveaux de pression verticaux
-  status = NF_REDEF (ncid_out)
-  CALL handle_err(status)
-  status=NF_DEF_DIM(ncid_out,'lev',llm,out_levdim)
-  CALL handle_err(status)
-  status=NF_DEF_VAR(ncid_out,'presnivs',NF_FLOAT,1,out_levdim,&
-                    presnivs_id)
-  CALL handle_err(status)
-  
-! fields
+  status=NF90_DEF_VAR(ncid_out,'grille_s',NF90_FLOAT,out_dim, out_tempid)
+  CALL handle_err(status)
+  status=NF_PUT_ATT_TEXT(ncid_out,out_tempid,'units', 6,'Kelvin')
+  status=NF_PUT_ATT_TEXT(ncid_out,out_tempid,'long_name',21,'scalar dynamical grid')
+  !
+  ! for INCA :
+  ! vertical levels "presnivs"
+  status=NF90_DEF_VAR(ncid_out,'presnivs',NF90_FLOAT,out_levdim, presnivs_id)
+  CALL handle_err(status)
+  status=NF_PUT_ATT_TEXT(ncid_out,presnivs_id,'units', 2,'Pa')
+  status=NF_PUT_ATT_TEXT(ncid_out,presnivs_id,'long_name',15,'Vertical levels')
+  ! surface geopotential height: named "phis" as the sfc geopotential, is actually phis/g
   out_dim(1)=out_lonvdim
   out_dim(2)=out_latudim
-
-  status = nf_def_var(ncid_out,'phis',NF_FLOAT,2,out_dim,phis_id)
-  CALL handle_err(status)
-  status = nf_def_var(ncid_out,'aire',NF_FLOAT,2,out_dim,area_id)
-  CALL handle_err(status)
-  status = nf_def_var(ncid_out,'mask',NF_INT  ,2,out_dim,mask_id)
-  CALL handle_err(status)
-
-  status=NF_ENDDEF(ncid_out)
-  CALL handle_err(status)
-
-  ! ecriture des variables
-#ifdef NC_DOUBLE
-  status=NF_PUT_VARA_DOUBLE(ncid_out,presnivs_id,1,llm,rlevdeg)
-  CALL handle_err(status)
-#else
-  status=NF_PUT_VARA_REAL(ncid_out,out_levid,1,llm,rlevdeg)
-  CALL handle_err(status)
-#endif
-
-  start(1)=1
-  start(2)=1
-  start(3)=1
-  start(4)=0
+  status = nf90_def_var(ncid_out,'phis',NF90_FLOAT,out_dim,phis_id)
+  CALL handle_err(status)
+  status=NF_PUT_ATT_TEXT(ncid_out,phis_id,'units', 1,'m')
+  status=NF_PUT_ATT_TEXT(ncid_out,phis_id,'long_name',27,'surface geopotential height')
+  ! gridcell area
+  status = nf90_def_var(ncid_out,'aire',NF90_FLOAT,out_dim,area_id)
+  CALL handle_err(status)
+  status=NF_PUT_ATT_TEXT(ncid_out,area_id,'units', 2,'m2')
+  status=NF_PUT_ATT_TEXT(ncid_out,area_id,'long_name',13,'gridcell area')
+  ! land-sea mask (nearest integer approx)
+  status = nf90_def_var(ncid_out,'mask',NF90_INT,out_dim,mask_id)
+  CALL handle_err(status)
+  status=NF_PUT_ATT_TEXT(ncid_out,mask_id,'long_name',27,'land-sea mask (nINT approx)')
+  
+  ! END the 'define' mode in netCDF file
+  status=NF_ENDDEF(ncid_out)  
+  CALL handle_err(status)
+  
+  ! WRITE the variables
+  !-------------------------
+  ! 1D : lonu, lonv,latu,latv ; INCA : presnivs
+  status=NF90_PUT_VAR(ncid_out,out_lonuid,rlonudeg,[1],[iip1])
+  CALL handle_err(status)
+  status=NF90_PUT_VAR(ncid_out,out_lonvid,rlonvdeg,[1],[iip1])
+  CALL handle_err(status)
+  status=NF90_PUT_VAR(ncid_out,out_latuid,rlatudeg,[1],[jjp1])
+  CALL handle_err(status)
+  status=NF90_PUT_VAR(ncid_out,out_latvid,rlatvdeg,[1],[jjm])
+  CALL handle_err(status)
+  status=NF90_PUT_VAR(ncid_out,presnivs_id,rlev,[1],[llm])
+  CALL handle_err(status)
+
+  ! 2D : grille_u,grille_v,grille_s ; INCA: phis,aire,mask
+  start(:)=1 
   COUNT(1)=iip1
-  COUNT(2)=jjp1
-  COUNT(3)=1
-  COUNT(4)=0
-
-  status = nf_put_vara_double(ncid_out, phis_id,start,count, phis_loc)
-  CALL handle_err(status)
-  status = nf_put_vara_double(ncid_out, area_id,start,count, aire)
-  CALL handle_err(status)
-  masque_int(:,:) = nINT(masque(:,:))
-  status = nf_put_vara_int(ncid_out, mask_id,start,count,masque_int)
-  CALL handle_err(status)
-  
-  ! fermeture du fichier netcdf
+  
+  COUNT(2)=jjp1  ! for "u" and "s" grids
+  status=NF90_PUT_VAR(ncid_out,out_uid,uwnd,start, count)
+  CALL handle_err(status)
+  COUNT(2)=jjm  ! for "v" grid
+  status=NF90_PUT_VAR(ncid_out,out_vid,vwnd,start, count)
+  CALL handle_err(status) 
+  COUNT(2)=jjp1  ! as "s" grid, for all the following vars
+  status=NF90_PUT_VAR(ncid_out,out_tempid,temp,start, count)
+  CALL handle_err(status)
+  status = nf90_put_var(ncid_out, phis_id, phis_loc,start,count)
+  CALL handle_err(status) 
+  status = nf90_put_var(ncid_out, area_id, aire,start,count)
+  CALL handle_err(status) 
+  status = nf90_put_var(ncid_out, mask_id,masque_int,start,count)
+  CALL handle_err(status)
+  
+  ! CLOSE netcdf file
   CALL ncclos(ncid_out,rcode_out)
+  write(*,*) "END grilles_gcm_netcdf_sub OK"
 
 END SUBROUTINE grilles_gcm_netcdf_sub
-
 
 
Index: LMDZ6/branches/Ocean_skin/libf/dyn3d_common/infotrac.F90
===================================================================
--- LMDZ6/branches/Ocean_skin/libf/dyn3d_common/infotrac.F90	(revision 4013)
+++ LMDZ6/branches/Ocean_skin/libf/dyn3d_common/infotrac.F90	(revision 4368)
@@ -1,98 +1,132 @@
-! $Id$
+!$Id$
 !
 MODULE infotrac
 
-! nqtot : total number of tracers and higher order of moment, water vapor and liquid included
-  INTEGER, SAVE :: nqtot
-!CR: on ajoute le nombre de traceurs de l eau
-  INTEGER, SAVE :: nqo
-
-! nbtr : number of tracers not including higher order of moment or water vapor or liquid
-!        number of tracers used in the physics
-  INTEGER, SAVE :: nbtr
-
-! CRisi: on retranche les isotopes des traceurs habituels
-! On fait un tableaux d'indices des traceurs qui passeront dans phytrac
-  INTEGER, SAVE :: nqtottr
-  INTEGER, ALLOCATABLE, DIMENSION(:), SAVE :: itr_indice
-
-! CRisi: nb traceurs peres= directement advectes par l'air
-  INTEGER, SAVE :: nqperes
-
-! ThL: nb traceurs INCA
-  INTEGER, SAVE :: nqINCA
-
-! ThL: nb traceurs CO2
-  INTEGER, SAVE :: nqCO2
-
-! Name variables
-  INTEGER,PARAMETER :: tname_lenmax=128 
-  CHARACTER(len=tname_lenmax), ALLOCATABLE, DIMENSION(:), SAVE :: tname ! tracer short name for restart and diagnostics
-  CHARACTER(len=tname_lenmax+3), ALLOCATABLE, DIMENSION(:), SAVE :: ttext ! tracer long name for diagnostics
-
-! iadv  : index of trasport schema for each tracer
-  INTEGER, ALLOCATABLE, DIMENSION(:), SAVE    :: iadv
-
-! niadv : vector keeping the coorspondance between all tracers(nqtot) treated in the 
-!         dynamic part of the code and the tracers (nbtr+2) used in the physics part of the code. 
-  INTEGER, ALLOCATABLE, DIMENSION(:), SAVE    :: niadv ! equivalent dyn / physique
-
-! CRisi: tableaux de fils
-  INTEGER, ALLOCATABLE, DIMENSION(:), SAVE    :: nqfils
-  INTEGER, ALLOCATABLE, DIMENSION(:), SAVE    :: nqdesc ! nombres de fils + nombre de tous les petits fils sur toutes les generations
-  INTEGER, SAVE :: nqdesc_tot
-  INTEGER, ALLOCATABLE, DIMENSION(:,:), SAVE    :: iqfils
-  INTEGER, ALLOCATABLE, DIMENSION(:), SAVE    :: iqpere
-  REAL :: qperemin,masseqmin,ratiomin ! MVals et CRisi
-  PARAMETER (qperemin=1e-16,masseqmin=1e-16,ratiomin=1e-16) ! MVals
-
-! conv_flg(it)=0 : convection desactivated for tracer number it 
-  INTEGER, ALLOCATABLE, DIMENSION(:), SAVE  :: conv_flg
-! pbl_flg(it)=0  : boundary layer diffusion desactivaded for tracer number it 
-  INTEGER, ALLOCATABLE, DIMENSION(:), SAVE  :: pbl_flg
-
-  CHARACTER(len=4),SAVE :: type_trac
-  CHARACTER(len=8),DIMENSION(:),ALLOCATABLE, SAVE :: solsym
-
-! CRisi: cas particulier des isotopes
-  LOGICAL,SAVE :: ok_isotopes,ok_iso_verif,ok_isotrac,ok_init_iso
-  INTEGER :: niso_possibles   
-  PARAMETER ( niso_possibles=5)
-  REAL, DIMENSION (niso_possibles),SAVE :: tnat,alpha_ideal
-  LOGICAL, DIMENSION(niso_possibles),SAVE ::  use_iso
-  INTEGER, ALLOCATABLE, DIMENSION(:,:), SAVE ::  iqiso ! donne indice iq en fn de (ixt,phase) 
-  INTEGER, ALLOCATABLE, DIMENSION(:), SAVE ::  iso_num ! donne numero iso entre 1 et niso_possibles en fn de nqtot
-  INTEGER, ALLOCATABLE, DIMENSION(:), SAVE ::  iso_indnum ! donne numero iso entre 1 et niso effectif en fn de nqtot
-  INTEGER, ALLOCATABLE, DIMENSION(:), SAVE ::  zone_num ! donne numero de la zone de tracage en fn de nqtot
-  INTEGER, ALLOCATABLE, DIMENSION(:), SAVE ::  phase_num ! donne numero de la zone de tracage en fn de nqtot
-  INTEGER, DIMENSION(niso_possibles), SAVE :: indnum_fn_num ! donne indice entre entre 1 et niso en fonction du numero d isotope entre 1 et niso_possibles
-  INTEGER, ALLOCATABLE, DIMENSION(:,:), SAVE ::  index_trac ! numero ixt en fn izone, indnum entre 1 et niso
-  INTEGER,SAVE :: niso,ntraceurs_zone,ntraciso
-
-#ifdef CPP_StratAer
-!--CK/OB for stratospheric aerosols
-  INTEGER, SAVE :: nbtr_bin
-  INTEGER, SAVE :: nbtr_sulgas
-  INTEGER, SAVE :: id_OCS_strat
-  INTEGER, SAVE :: id_SO2_strat
-  INTEGER, SAVE :: id_H2SO4_strat
-  INTEGER, SAVE :: id_BIN01_strat
-  INTEGER, SAVE :: id_TEST_strat
-#endif
- 
+   USE       strings_mod, ONLY: msg, fmsg, maxlen, cat, dispTable, int2str, bool2str, strStack, strParse
+   USE readTracFiles_mod, ONLY: trac_type, readTracersFiles, tracers, setGeneration, itZonIso, nbIso, tran0, delPhase, &
+                        getKey, isot_type, readIsotopesFile, isotope, maxTableWidth, iqIsoPha, ntiso, ixIso, addPhase, &
+                   indexUpdate, isoSelect, isoPhas, isoZone, isoName, isoKeys, iH2O, isoCheck, nphas, nzone, niso
+   IMPLICIT NONE
+
+   PRIVATE
+
+   !=== FOR TRACERS:
+   PUBLIC :: init_infotrac                                 !--- Initialization of the tracers
+   PUBLIC :: tracers, type_trac, types_trac                !--- Full tracers database, tracers type keyword
+   PUBLIC :: nqtot,   nbtr,   nqo,   nqCO2,   nqtottr      !--- Main dimensions
+   PUBLIC :: conv_flg, pbl_flg                             !--- Convection & boundary layer activation keys
+
+   !=== FOR ISOTOPES: General
+   PUBLIC :: isot_type, nbIso                              !--- Derived type, full isotopes families database + nb of families
+   PUBLIC :: isoSelect, ixIso                              !--- Isotopes family selection tool + selected family index
+   !=== FOR ISOTOPES: Specific to water
+   PUBLIC :: iH2O                                          !--- H2O isotopes class index
+   PUBLIC :: min_qParent, min_qMass, min_ratio             !--- Min. values for various isotopic quantities
+   !=== FOR ISOTOPES: Depending on the selected isotopes family
+   PUBLIC :: isotope, isoKeys                              !--- Selected isotopes database + associated keys (cf. getKey)
+   PUBLIC :: isoName, isoZone, isoPhas                     !--- Isotopes and tagging zones names, phases
+   PUBLIC :: niso,    nzone,   nphas,   ntiso              !---  " " numbers + isotopes & tagging tracers number
+   PUBLIC :: itZonIso                                      !--- idx "it" (in "isoName(1:niso)") = function(tagging idx, isotope idx)
+   PUBLIC :: iqIsoPha                                      !--- idx "iq" (in "qx") = function(isotope idx, phase idx) + aliases
+   PUBLIC :: isoCheck                                      !--- Run isotopes checking routines
+   !=== FOR BOTH TRACERS AND ISOTOPES
+   PUBLIC :: getKey                                        !--- Get a key from "tracers" or "isotope"
+
+!=== CONVENTIONS FOR TRACERS NUMBERS:
+!  |--------------------+-----------------------+-----------------+---------------+----------------------------|
+!  | water in different |    water tagging      |  water isotopes | other tracers | additional tracers moments |
+!  | phases: H2O_[gls]  |      isotopes         |                 |               |  for higher order schemes  |
+!  |--------------------+-----------------------+-----------------+---------------+----------------------------|
+!  |                    |                       |                 |               |                            |
+!  |<--     nqo      -->|<-- nqo*niso* nzone -->|<-- nqo*niso  -->|<--  nbtr   -->|<--        (nmom)        -->|         
+!  |                    |                                         |                                            |
+!  |                    |<-- nqo*niso*(nzone+1)  =   nqo*ntiso -->|<--    nqtottr = nbtr + nmom             -->|
+!  |                                                                              = nqtot - nqo*(ntiso+1)      |
+!  |                                                                                                           |
+!  |<--                        nqtrue  =  nbtr + nqo*(ntiso+1)                 -->|                            |
+!  |                                                                                                           |
+!  |<--                        nqtot   =  nqtrue + nmom                                                     -->|
+!  |                                                                                                           |
+!  |-----------------------------------------------------------------------------------------------------------|
+!  NOTES FOR THIS TABLE:
+!  * Used "niso", "nzone" and "ntiso" are components of "isotopes(ip)" for water (isotopes(ip)%parent == 'H2O'),
+!    since water is so far the sole tracers family, except passive CO2, removed from the main tracers table.
+!  * For water, "nqo" is equal to the more general field "isotopes(ip)%nphas".
+!  * "niso", "nzone", "ntiso", "nphas" are defined for other isotopic tracers families, if any.
+!
+!=== DERIVED TYPE EMBEDDING MOST OF THE TRACERS-RELATED QUANTITIES (LENGTH: nqtot)
+!    Each entry is accessible using "%" sign.
+!  |-------------+------------------------------------------------------+-------------+------------------------+
+!  |  entry      | Meaning                                              | Former name | Possible values        |
+!  |-------------+------------------------------------------------------+-------------+------------------------+
+!  | name        | Name (short)                                         | tname       |                        |
+!  | gen0Name    | Name of the 1st generation ancestor                  | /           |                        |
+!  | parent      | Name of the parent                                   | /           |                        |
+!  | longName    | Long name (with adv. scheme suffix) for outputs      | ttext       |                        |
+!  | type        | Type (so far: tracer or tag)                         | /           | tracer,tag             |
+!  | phase       | Phases list ("g"as / "l"iquid / "s"olid)             | /           | [g][l][s]              |
+!  | component   | Name(s) of the merged/cumulated section(s)           | /           | coma-separated names   |
+!  | iGeneration | Generation (>=1)                                     | /           |                        |
+!  | iqParent    | Index of the parent tracer                           | iqpere      | 1:nqtot                |
+!  | iqDescen    | Indexes of the childs       (all generations)        | iqfils      | 1:nqtot                |
+!  | nqDescen    | Number of the descendants   (all generations)        | nqdesc      | 1:nqtot                |
+!  | nqChildren  | Number of childs            (1st generation only)    | nqfils      | 1:nqtot                |
+!  | keys        | key/val pairs accessible with "getKey" routine       | /           |                        |
+!  | iadv        | Advection scheme number                              | iadv        | 1,2,10-20(exc.15,19),30|
+!  | isAdvected  | advected tracers flag (.TRUE. if iadv >= 0)          | /           | nqtrue  .TRUE. values  |
+!  | isInPhysics | tracers not extracted from the main table in physics | /           | nqtottr .TRUE. values  |
+!  | iso_iGroup  | Isotopes group index in isotopes(:)                  | /           | 1:nbIso                |
+!  | iso_iName   | Isotope  name  index in isotopes(iso_iGroup)%trac(:) | iso_indnum  | 1:niso                 |
+!  | iso_iZone   | Isotope  zone  index in isotopes(iso_iGroup)%zone(:) | zone_num    | 1:nzone                |
+!  | iso_iPhas   | Isotope  phase index in isotopes(iso_iGroup)%phas(:) | phase_num   | 1:nphas                |
+!  +-------------+------------------------------------------------------+-------------+------------------------+
+!
+!=== DERIVED TYPE EMBEDDING MOST OF THE ISOTOPES-RELATED QUANTITIES (LENGTH: nbIso, NUMBER OF ISOTOPES FAMILIES)
+!    Each entry is accessible using "%" sign.
+!  |-----------------+--------------------------------------------------+--------------------+-----------------+
+!  |  entry | length | Meaning                                          |    Former name     | Possible values |
+!  |-----------------+--------------------------------------------------+--------------------+-----------------+
+!  | parent          | Parent tracer (isotopes family name)             |                    |                 |
+!  | keys   | niso   | Isotopes keys/values pairs list + number         |                    |                 |
+!  | trac   | ntiso  | Isotopes + tagging tracers list + number         | / | ntraciso       |                 |
+!  | zone   | nzone  | Geographic tagging zones   list + number         | / | ntraceurs_zone |                 |
+!  | phase  | nphas  | Phases                     list + number         |                    | [g][l][s], 1:3  |
+!  | iqIsoPha        | Index in "qx"           = f(name(1:ntiso)),phas) | iqiso              | 1:nqtot         |
+!  | itZonIso        | Index in "trac(1:ntiso)"= f(zone, name(1:niso))  | index_trac         | 1:ntiso         |
+!  +-----------------+--------------------------------------------------+--------------------+-----------------+
+
+   REAL, PARAMETER :: min_qParent = 1.e-30, min_qMass = 1.e-18, min_ratio = 1.e-16 ! MVals et CRisi
+
+   !=== DIMENSIONS OF THE TRACERS TABLES AND OTHER SCALAR VARIABLES
+   INTEGER,               SAVE :: nqtot,  &                     !--- Tracers nb in dynamics (incl. higher moments + H2O)
+                                  nbtr,   &                     !--- Tracers nb in physics  (excl. higher moments + H2O)
+                                  nqo,    &                     !--- Number of water phases
+                                  nqtottr, &                    !--- Number of tracers passed to phytrac (TO BE DELETED ?)
+                                  nqCO2                         !--- Number of tracers of CO2  (ThL)
+   CHARACTER(LEN=maxlen), SAVE :: type_trac                     !--- Keyword for tracers type(s)
+   CHARACTER(LEN=maxlen), SAVE, ALLOCATABLE :: types_trac(:)    !--- Keyword for tracers type(s), parsed version
+
+   !=== VARIABLES FOR INCA
+   INTEGER,               SAVE, ALLOCATABLE :: conv_flg(:), &   !--- Convection     activation ; needed for INCA        (nbtr)
+                                                pbl_flg(:)      !--- Boundary layer activation ; needed for INCA        (nbtr)
+
 CONTAINS
 
-  SUBROUTINE infotrac_init
-    USE control_mod, ONLY: planet_type, config_inca
+SUBROUTINE init_infotrac
+   USE control_mod, ONLY: planet_type
 #ifdef REPROBUS
-    USE CHEM_REP, ONLY : Init_chem_rep_trac
-#endif
-    IMPLICIT NONE
-!=======================================================================
+   USE CHEM_REP,    ONLY: Init_chem_rep_trac
+#endif
+   IMPLICIT NONE
+!==============================================================================================================================
 !
 !   Auteur:  P. Le Van /L. Fairhead/F.Hourdin
 !   -------
-!   Modif special traceur F.Forget 05/94
-!   Modif M-A Filiberti 02/02 lecture de traceur.def
+!
+!   Modifications:
+!   --------------
+!   05/94: F.Forget      Modif special traceur
+!   02/02: M-A Filiberti Lecture de traceur.def
+!   01/22: D. Cugnet     Nouveaux tracer.def et tracer_*.def + encapsulation (types trac_type et isot_type)
 !
 !   Objet:
@@ -100,882 +134,285 @@
 !   GCM LMD nouvelle grille
 !
-!=======================================================================
+!==============================================================================================================================
 !   ... modification de l'integration de q ( 26/04/94 ) ....
-!-----------------------------------------------------------------------
-! Declarations
-
-    INCLUDE "dimensions.h"
-    INCLUDE "iniprint.h"
-
+!------------------------------------------------------------------------------------------------------------------------------
+! Declarations:
+   INCLUDE "dimensions.h"
+   INCLUDE "iniprint.h"
+
+!------------------------------------------------------------------------------------------------------------------------------
 ! Local variables
-    INTEGER, ALLOCATABLE, DIMENSION(:) :: hadv  ! index of horizontal trasport schema
-    INTEGER, ALLOCATABLE, DIMENSION(:) :: vadv  ! index of vertical trasport schema
-
-    INTEGER, ALLOCATABLE, DIMENSION(:) :: hadv_inca  ! index of horizontal trasport schema
-    INTEGER, ALLOCATABLE, DIMENSION(:) :: vadv_inca  ! index of vertical trasport schema
-
-    INTEGER, ALLOCATABLE, DIMENSION(:) :: conv_flg_inca
-    INTEGER, ALLOCATABLE, DIMENSION(:) :: pbl_flg_inca
-    CHARACTER(len=8), ALLOCATABLE, DIMENSION(:) :: solsym_inca
-
-    CHARACTER(len=tname_lenmax), ALLOCATABLE, DIMENSION(:) :: tnom_0  ! tracer short name
-    CHARACTER(len=tname_lenmax), ALLOCATABLE, DIMENSION(:) :: tnom_transp ! transporting fluid short name: CRisi
-    CHARACTER(len=3), DIMENSION(30) :: descrq
-    CHARACTER(len=1), DIMENSION(3)  :: txts
-    CHARACTER(len=2), DIMENSION(9)  :: txtp
-    CHARACTER(len=tname_lenmax)               :: str1,str2
-  
-    INTEGER :: nqtrue  ! number of tracers read from tracer.def, without higer order of moment
-    INTEGER :: iq, new_iq, iiq, jq, ierr,itr
-    INTEGER :: ifils,ipere,generation ! CRisi
-    LOGICAL :: continu,nouveau_traceurdef
-    INTEGER :: IOstatus ! gestion de la retrocompatibilite de traceur.def
-    CHARACTER(len=2*tname_lenmax+1) :: tchaine    
-
-    character(len=*),parameter :: modname="infotrac_init"
-
-!-----------------------------------------------------------------------
+   INTEGER, ALLOCATABLE :: hadv(:), vadv(:)                          !--- Horizontal/vertical transport scheme number
+#ifdef INCA
+   INTEGER, ALLOCATABLE :: had (:), hadv_inca(:), conv_flg_inca(:), &!--- Variables specific to INCA
+                           vad (:), vadv_inca(:),  pbl_flg_inca(:)
+   CHARACTER(LEN=8), ALLOCATABLE :: solsym_inca(:)                   !--- Tracers names for INCA
+   INTEGER :: nqINCA
+#endif
+   CHARACTER(LEN=2)      ::   suff(9)                                !--- Suffixes for schemes of order 3 or 4 (Prather)
+   CHARACTER(LEN=3)      :: descrq(30)                               !--- Advection scheme description tags
+   CHARACTER(LEN=maxlen) :: msg1                                     !--- String for messages
+   INTEGER :: fType                                                  !--- Tracers description file type ; 0: none
+                                                                     !--- 1/2/3: "traceur.def"/"tracer.def"/"tracer_*.def"
+   INTEGER :: nqtrue                                                 !--- Tracers nb from tracer.def (no higher order moments)
+   INTEGER :: iad                                                    !--- Advection scheme number
+   INTEGER :: ic, iq, jq, it, nt, im, nm, iz, k                      !--- Indexes and temporary variables
+   LOGICAL :: lerr, ll, lRepr
+   CHARACTER(LEN=1) :: p
+   TYPE(trac_type), ALLOCATABLE, TARGET :: ttr(:)
+   TYPE(trac_type), POINTER             :: t1, t(:)
+   INTEGER :: ierr
+
+   CHARACTER(LEN=*), PARAMETER :: modname="init_infotrac"
+!------------------------------------------------------------------------------------------------------------------------------
 ! Initialization :
-!
-    txts=(/'x','y','z'/)
-    txtp=(/'x ','y ','z ','xx','xy','xz','yy','yz','zz'/)
-
-    descrq(14)='VLH'
-    descrq(10)='VL1'
-    descrq(11)='VLP'
-    descrq(12)='FH1'
-    descrq(13)='FH2'
-    descrq(16)='PPM'
-    descrq(17)='PPS'
-    descrq(18)='PPP'
-    descrq(20)='SLP'
-    descrq(30)='PRA'
+!------------------------------------------------------------------------------------------------------------------------------
+   suff          = ['x ','y ','z ','xx','xy','xz','yy','yz','zz']
+   descrq( 1: 2) = ['LMV','BAK']
+   descrq(10:20) = ['VL1','VLP','FH1','FH2','VLH','   ','PPM','PPS','PPP','   ','SLP']
+   descrq(30)    =  'PRA'
     
-
-    ! Coherence test between parameter type_trac, config_inca and preprocessing keys
-    IF (type_trac=='inca') THEN
-       WRITE(lunout,*) 'You have chosen to couple with INCA chemistry model : type_trac=', &
-            type_trac,' config_inca=',config_inca
-       IF (config_inca/='aero' .AND. config_inca/='aeNP' .AND. config_inca/='chem') THEN
-          WRITE(lunout,*) 'Incoherence between type_trac and config_inca. Model stops. Modify run.def'
-          CALL abort_gcm('infotrac_init','Incoherence between type_trac and config_inca',1)
-       ENDIF
+   CALL msg('type_trac = "'//TRIM(type_trac)//'"', modname)
+   IF(strParse(type_trac, '|', types_trac, n=nt)) CALL abort_gcm(modname,'can''t parse "type_trac = '//TRIM(type_trac)//'"',1)
+
+   !---------------------------------------------------------------------------------------------------------------------------
+   DO it = 1, nt                                                          !--- nt>1=> "type_trac": coma-separated keywords list
+   !---------------------------------------------------------------------------------------------------------------------------
+      !--- MESSAGE ABOUT THE CHOSEN CONFIGURATION
+      msg1 = 'For type_trac = "'//TRIM(types_trac(it))//'":'
+      SELECT CASE(types_trac(it))
+         CASE('inca'); CALL msg(TRIM(msg1)//' coupling with INCA chemistry model',        modname)
+         CASE('inco'); CALL msg(TRIM(msg1)//' coupling jointly with INCA and CO2 cycle',  modname)
+         CASE('repr'); CALL msg(TRIM(msg1)//' coupling with REPROBUS chemistry model',    modname)
+         CASE('co2i'); CALL msg(TRIM(msg1)//' you have chosen to run with CO2 cycle',     modname)
+         CASE('coag'); CALL msg(TRIM(msg1)//' tracers are treated for COAGULATION tests', modname)
+         CASE('lmdz'); CALL msg(TRIM(msg1)//' tracers are treated in LMDZ only',          modname)
+         CASE DEFAULT; CALL abort_gcm(modname,'type_trac='//TRIM(types_trac(it))//' not possible yet.',1)
+      END SELECT
+
+      !--- COHERENCE TEST BETWEEN "type_trac" AND PREPROCESSING KEYS
+      SELECT CASE(types_trac(it))
+         CASE('inca', 'inco')
 #ifndef INCA
-       WRITE(lunout,*) 'To run this option you must add cpp key INCA and compile with INCA code'
-       CALL abort_gcm('infotrac_init','You must compile with cpp key INCA',1)
-#endif
-    ELSE IF (type_trac=='repr') THEN
-       WRITE(lunout,*) 'You have chosen to couple with REPROBUS chemestry model : type_trac=', type_trac
+            CALL abort_gcm(modname, 'You must add cpp key INCA and compile with INCA code', 1)
+#endif
+         CASE('repr')
 #ifndef REPROBUS
-       WRITE(lunout,*) 'To run this option you must add cpp key REPROBUS and compile with REPRPBUS code'
-       CALL abort_gcm('infotrac_init','You must compile with cpp key REPROBUS',1)
-#endif
-    ELSE IF (type_trac == 'co2i') THEN
-       WRITE(lunout,*) 'You have chosen to run with CO2 cycle: type_trac=', type_trac
-    ELSE IF (type_trac == 'coag') THEN
-       WRITE(lunout,*) 'Tracers are treated for COAGULATION tests : type_trac=', type_trac
+            CALL abort_gcm(modname, 'You must add cpp key REPROBUS and compile with REPROBUS code', 1)
+#endif
+         CASE('coag')
 #ifndef CPP_StratAer
-       WRITE(lunout,*) 'To run this option you must add cpp key StratAer and compile with StratAer code'
-       CALL abort_gcm('infotrac_init','You must compile with cpp key StratAer',1)
-#endif
-    ELSE IF (type_trac == 'lmdz') THEN
-       WRITE(lunout,*) 'Tracers are treated in LMDZ only : type_trac=', type_trac
-    ELSE IF (type_trac == 'inco') THEN ! ThL
-       WRITE(lunout,*) 'Using jointly INCA and CO2 cycle: type_trac =', type_trac
-       IF (config_inca/='aero' .AND. config_inca/='aeNP' .AND. config_inca/='chem') THEN
-          WRITE(lunout,*) 'Incoherence between type_trac and config_inca. Model stops. Modify run.def'
-          CALL abort_gcm('infotrac_init','Incoherence between type_trac and config_inca',1)
-       ENDIF
-#ifndef INCA
-       WRITE(lunout,*) 'To run this option you must add cpp key INCA and compilewith INCA code'
-       CALL abort_gcm('infotrac_init','You must compile with cpp key INCA',1)
-#endif    
-    ELSE
-       WRITE(lunout,*) 'type_trac=',type_trac,' not possible. Model stops'
-       CALL abort_gcm('infotrac_init','bad parameter',1)
-    ENDIF
-
-    ! Test if config_inca is other then none for run without INCA
-    IF (type_trac/='inca' .AND. type_trac/='inco' .AND. config_inca/='none') THEN
-       WRITE(lunout,*) 'config_inca will now be changed to none as you do not couple with INCA model'
-       config_inca='none'
-    ENDIF
-
-!-----------------------------------------------------------------------
-!
-! 1) Get the true number of tracers + water vapor/liquid
-!    Here true tracers (nqtrue) means declared tracers (only first order)
-!
-!-----------------------------------------------------------------------
-    IF (type_trac == 'lmdz' .OR. type_trac == 'repr' .OR. type_trac == 'coag' .OR. type_trac == 'co2i') THEN
-       IF (type_trac=='co2i') THEN
-          nqCO2 = 1
-       ELSE
-          nqCO2 = 0
-       ENDIF
-       OPEN(90,file='traceur.def',form='formatted',status='old', iostat=ierr)
-       IF(ierr.EQ.0) THEN
-          WRITE(lunout,*) trim(modname),': Open traceur.def : ok'
-          READ(90,*) nqtrue
-          write(lunout,*) 'nqtrue=',nqtrue
-       ELSE 
-          WRITE(lunout,*) trim(modname),': Failed opening traceur.def'
-          CALL abort_gcm(modname,"file traceur.def not found!",1)
-       ENDIF
-!jyg<
-!!       if ( planet_type=='earth') then
-!!         ! For Earth, water vapour & liquid tracers are not in the physics
-!!         nbtr=nqtrue-2
-!!       else
-!!         ! Other planets (for now); we have the same number of tracers
-!!         ! in the dynamics than in the physics
-!!         nbtr=nqtrue
-!!       endif
-!>jyg
-    ELSE ! type_trac=inca or inco
-       IF (type_trac=='inco') THEN
-          nqCO2 = 1
-       ELSE
-          nqCO2 = 0
-       ENDIF
-!jyg<
-       ! The traceur.def file is used to define the number "nqo" of water phases 
-       ! present in the simulation. Default : nqo = 2.
-       OPEN(90,file='traceur.def',form='formatted',status='old', iostat=ierr)
-       IF(ierr.EQ.0) THEN
-          WRITE(lunout,*) trim(modname),': Open traceur.def : ok'
-          READ(90,*) nqo
-       ELSE 
-          WRITE(lunout,*) trim(modname),': Failed opening traceur.def'
-          CALL abort_gcm(modname,"file traceur.def not found!",1)
-       ENDIF
-       IF (nqo /= 2 .AND. nqo /= 3 ) THEN
-          IF (nqo == 4 .AND. type_trac=='inco') THEN ! ThL
-             WRITE(lunout,*) trim(modname),': you are coupling with INCA, and also using CO2i.'
-             nqo = 3    ! A ameliorier... je force 3 traceurs eau...  ThL
-             WRITE(lunout,*) trim(modname),': nqo = ',nqo
-          ELSE
-          WRITE(lunout,*) trim(modname),': nqo=',nqo, ' is not allowed. Only 2 or 3 water phases allowed'
-          CALL abort_gcm('infotrac_init','Bad number of water phases',1)
-          ENDIF
-       ENDIF
-       ! nbtr has been read from INCA by init_const_lmdz() in gcm.F 
+            CALL abort_gcm(modname, 'You must add cpp key StratAer and compile with StratAer code', 1)
+#endif
+      END SELECT
+
+   !---------------------------------------------------------------------------------------------------------------------------
+   END DO
+   !---------------------------------------------------------------------------------------------------------------------------
+
+   nqCO2 = COUNT( [ANY(types_trac == 'inco') .OR. (ANY(types_trac == 'co2i') .AND. ANY(types_trac == 'inca'))] )
+
+!==============================================================================================================================
+! 1) Get the numbers of: true (first order only) tracers "nqtrue", water tracers "nqo" (vapor/liquid/solid)
+!==============================================================================================================================
+   lRepr = ANY(types_trac(:) == 'repr')
+   IF(readTracersFiles(type_trac, fType, lRepr)) CALL abort_gcm(modname, 'problem with tracers file(s)',1)
+   !---------------------------------------------------------------------------------------------------------------------------
+   IF(fType == 0) CALL abort_gcm(modname, 'Missing tracers file: "traceur.def", "tracer.def" or "tracer_<keyword>.def file.',1)
+   !---------------------------------------------------------------------------------------------------------------------------
+   IF(fType == 1 .AND. ANY(['inca','inco'] == type_trac)) THEN       !=== FOUND OLD STYLE INCA "traceur.def" (single type_trac)
+   !---------------------------------------------------------------------------------------------------------------------------
 #ifdef INCA
-       CALL Init_chem_inca_trac(nqINCA)
-#else
-       nqINCA=0
-#endif
-       nbtr=nqINCA+nqCO2
-       nqtrue=nbtr+nqo
-       WRITE(lunout,*) trim(modname),': nqo = ',nqo
-       WRITE(lunout,*) trim(modname),': nbtr = ',nbtr
-       WRITE(lunout,*) trim(modname),': nqtrue = ',nqtrue
-       WRITE(lunout,*) trim(modname),': nqCO2 = ',nqCO2
-       WRITE(lunout,*) trim(modname),': nqINCA = ',nqINCA
-       ALLOCATE(hadv_inca(nqINCA), vadv_inca(nqINCA), conv_flg_inca(nqINCA), pbl_flg_inca(nqINCA), solsym_inca(nqINCA))
-    ENDIF   ! type_trac 'inca' ou 'inco'
-!>jyg
-
-    IF ((planet_type=="earth").and.(nqtrue < 2)) THEN
-       WRITE(lunout,*) trim(modname),': nqtrue=',nqtrue, ' is not allowed. 2 tracers is the minimum'
-       CALL abort_gcm('infotrac_init','Not enough tracers',1)
-    ENDIF
-    
-!jyg<
-       
-!
-! Allocate variables depending on nqtrue
-!
-    ALLOCATE(tnom_0(nqtrue), hadv(nqtrue), vadv(nqtrue),tnom_transp(nqtrue))
-
-
-!-----------------------------------------------------------------------
-! 2)     Choix  des schemas d'advection pour l'eau et les traceurs
-!
-!     iadv = 1    schema  transport type "humidite specifique LMD"
-!     iadv = 2    schema   amont
-!     iadv = 14   schema  Van-leer + humidite specifique 
-!                            Modif F.Codron
-!     iadv = 10   schema  Van-leer (retenu pour l'eau vapeur et liquide)
-!     iadv = 11   schema  Van-Leer pour hadv et version PPM (Monotone) pour vadv
-!     iadv = 12   schema  Frederic Hourdin I
-!     iadv = 13   schema  Frederic Hourdin II
-!     iadv = 16   schema  PPM Monotone(Collela & Woodward 1984)
-!     iadv = 17   schema  PPM Semi Monotone (overshoots autorises)
-!     iadv = 18   schema  PPM Positif Defini (overshoots undershoots autorises)
-!     iadv = 20   schema  Slopes
-!     iadv = 30   schema  Prather
-!
-!        Dans le tableau q(ij,l,iq) : iq = 1  pour l'eau vapeur
-!                                     iq = 2  pour l'eau liquide
-!       Et eventuellement             iq = 3,nqtot pour les autres traceurs
-!
-!        iadv(1): choix pour l'eau vap. et  iadv(2) : choix pour l'eau liq.
-!------------------------------------------------------------------------
-!
-!    Get choice of advection schema from file tracer.def or from INCA
-!---------------------------------------------------------------------
-    IF (type_trac == 'lmdz' .OR. type_trac == 'repr' .OR. type_trac == 'coag' .OR. type_trac == 'co2i') THEN
-
-          ! Continue to read tracer.def
-          DO iq=1,nqtrue
-
-             write(*,*) 'infotrac 237: iq=',iq
-             ! CRisi: ajout du nom du fluide transporteur
-             ! mais rester retro compatible
-             READ(90,'(I2,X,I2,X,A)',IOSTAT=IOstatus) hadv(iq),vadv(iq),tchaine
-             write(lunout,*) 'iq,hadv(iq),vadv(iq)=',iq,hadv(iq),vadv(iq)
-             write(lunout,*) 'tchaine=',trim(tchaine)
-             write(*,*) 'infotrac 238: IOstatus=',IOstatus
-             if (IOstatus.ne.0) then
-                CALL abort_gcm('infotrac_init','Pb dans la lecture de traceur.def',1)
-             endif
-             ! Y-a-t-il 1 ou 2 noms de traceurs? -> On regarde s'il y a un
-             ! espace ou pas au milieu de la chaine.
-             continu=.true.
-             nouveau_traceurdef=.false.
-             iiq=1
-             do while (continu)
-                if (tchaine(iiq:iiq).eq.' ') then
-                  nouveau_traceurdef=.true.
-                  continu=.false.
-                else if (iiq.lt.LEN_TRIM(tchaine)) then
-                  iiq=iiq+1
-                else
-                  continu=.false.
-                endif
-             enddo
-             write(*,*) 'iiq,nouveau_traceurdef=',iiq,nouveau_traceurdef
-             if (nouveau_traceurdef) then
-                write(lunout,*) 'C''est la nouvelle version de traceur.def'
-                tnom_0(iq)=tchaine(1:iiq-1)
-                tnom_transp(iq)=tchaine(iiq+1:)
-             else
-                write(lunout,*) 'C''est l''ancienne version de traceur.def'
-                write(lunout,*) 'On suppose que les traceurs sont tous d''air'
-                tnom_0(iq)=tchaine
-                tnom_transp(iq) = 'air'
-             endif
-             write(lunout,*) 'tnom_0(iq)=<',trim(tnom_0(iq)),'>'
-             write(lunout,*) 'tnom_transp(iq)=<',trim(tnom_transp(iq)),'>'
-
-          ENDDO!DO iq=1,nqtrue
-
-          CLOSE(90)  
-
-       WRITE(lunout,*) trim(modname),': Valeur de traceur.def :'
-       WRITE(lunout,*) trim(modname),': nombre total de traceurs ',nqtrue
-       DO iq=1,nqtrue
-          WRITE(lunout,*) hadv(iq),vadv(iq),' ',trim(tnom_0(iq)),' ',trim(tnom_transp(iq))
-       END DO
-
-       IF ( planet_type=='earth') THEN
-         !CR: nombre de traceurs de l eau
-         IF (tnom_0(3) == 'H2Oi') THEN
-            nqo=3
-         ELSE
-            nqo=2
-         ENDIF
-         ! For Earth, water vapour & liquid tracers are not in the physics
-         nbtr=nqtrue-nqo
-       ELSE
-         ! Other planets (for now); we have the same number of tracers
-         ! in the dynamics than in the physics
-         nbtr=nqtrue
-       ENDIF
-
-#ifdef CPP_StratAer
-       IF (type_trac == 'coag') THEN
-         nbtr_bin=0
-         nbtr_sulgas=0
-         DO iq=1,nqtrue
-           IF (tnom_0(iq)(1:3)=='BIN') THEN !check if tracer name contains 'BIN'
-             nbtr_bin=nbtr_bin+1
-           ENDIF
-           IF (tnom_0(iq)(1:3)=='GAS') THEN !check if tracer name contains 'GAS'
-             nbtr_sulgas=nbtr_sulgas+1
-           ENDIF
-         ENDDO
-         print*,'nbtr_bin=',nbtr_bin
-         print*,'nbtr_sulgas=',nbtr_sulgas
-         DO iq=1,nqtrue
-           IF (tnom_0(iq)=='GASOCS') THEN
-             id_OCS_strat=iq-nqo
-           ENDIF
-           IF (tnom_0(iq)=='GASSO2') THEN
-             id_SO2_strat=iq-nqo
-           ENDIF
-           IF (tnom_0(iq)=='GASH2SO4') THEN
-             id_H2SO4_strat=iq-nqo
-           ENDIF
-           IF (tnom_0(iq)=='BIN01') THEN
-             id_BIN01_strat=iq-nqo
-           ENDIF
-           IF (tnom_0(iq)=='GASTEST') THEN
-             id_TEST_strat=iq-nqo
-           ENDIF
-         ENDDO
-         print*,'id_OCS_strat  =',id_OCS_strat
-         print*,'id_SO2_strat  =',id_SO2_strat
-         print*,'id_H2SO4_strat=',id_H2SO4_strat
-         print*,'id_BIN01_strat=',id_BIN01_strat
-       ENDIF
-#endif
-
-    ENDIF  ! (type_trac == 'lmdz' .OR. type_trac == 'repr' .OR. type_trac = 'coag' .OR. type_trac = 'co2i')
-!jyg<
-!
-
-! Transfert number of tracers to Reprobus
-    IF (type_trac == 'repr') THEN
+      nqo = SIZE(tracers) - nqCO2
+      CALL Init_chem_inca_trac(nqINCA)                               !--- Get nqINCA from INCA
+      nbtr = nqINCA + nqCO2                                          !--- Number of tracers passed to phytrac
+      nqtrue = nbtr + nqo                                            !--- Total number of "true" tracers
+      IF(ALL([2,3] /= nqo)) CALL abort_gcm(modname, 'Only 2 or 3 water phases allowed ; found nqo='//TRIM(int2str(nqo)), 1)
+      ALLOCATE(hadv(nqtrue), hadv_inca(nqINCA), conv_flg_inca(nqINCA), solsym_inca(nqINCA))
+      ALLOCATE(vadv(nqtrue), vadv_inca(nqINCA),  pbl_flg_inca(nqINCA))
+      CALL init_transport(solsym_inca, conv_flg_inca, pbl_flg_inca, hadv_inca, vadv_inca)
+      ALLOCATE(ttr(nqtrue))
+      ttr(1:nqo+nqCO2)                    = tracers
+      ttr(1    :      nqo   )%component   = 'lmdz'
+      ttr(1+nqo:nqCO2+nqo   )%component   = 'co2i'
+      ttr(1+nqo+nqCO2:nqtrue)%component   = 'inca'
+      ttr(1+nqo      :nqtrue)%name        = [('CO2     ', k=1, nqCO2), solsym_inca]
+      ttr(1+nqo+nqCO2:nqtrue)%parent      = tran0
+      ttr(1+nqo+nqCO2:nqtrue)%phase       = 'g'
+      lerr = getKey('hadv', had, ky=tracers(:)%keys)
+      lerr = getKey('vadv', vad, ky=tracers(:)%keys)
+      hadv(1:nqo) = had(:); hadv(nqo+1:nqtrue) = hadv_inca
+      vadv(1:nqo) = vad(:); vadv(nqo+1:nqtrue) = vadv_inca
+      CALL MOVE_ALLOC(FROM=ttr, TO=tracers)
+      CALL setGeneration(tracers)                                    !--- SET FIELDS %iGeneration, %gen0Name
+      DEALLOCATE(had, hadv_inca, vad, vadv_inca, conv_flg_inca, pbl_flg_inca, solsym_inca)
+#endif
+   !---------------------------------------------------------------------------------------------------------------------------
+   ELSE                                                              !=== OTHER CASES (OLD OR NEW FORMAT, NO INCA MODULE)
+   !---------------------------------------------------------------------------------------------------------------------------
+      nqo    =        COUNT(delPhase(tracers(:)%name)     == 'H2O' &
+                               .AND. tracers(:)%component == 'lmdz') !--- Number of water phases
+      nqtrue = SIZE(tracers)                                         !--- Total number of "true" tracers
+      nbtr   = nqtrue-COUNT(delPhase(tracers(:)%gen0Name) == 'H2O' &
+                               .AND. tracers(:)%component == 'lmdz') !--- Number of tracers passed to phytrac
+#ifdef INCA
+      nqINCA = COUNT(tracers(:)%component == 'inca')
+#endif
+      lerr = getKey('hadv', hadv, ky=tracers(:)%keys)
+      lerr = getKey('vadv', vadv, ky=tracers(:)%keys)
+   !---------------------------------------------------------------------------------------------------------------------------
+   END IF
+   !---------------------------------------------------------------------------------------------------------------------------
+
+   !--- Transfert the number of tracers to Reprobus
 #ifdef REPROBUS
-       CALL Init_chem_rep_trac(nbtr,nqo,tnom_0)
-#endif
-    ENDIF
-!
-! Allocate variables depending on nbtr
-!
-    ALLOCATE(conv_flg(nbtr), pbl_flg(nbtr), solsym(nbtr))
-    conv_flg(:) = 1 ! convection activated for all tracers
-    pbl_flg(:)  = 1 ! boundary layer activated for all tracers
-
-    IF (type_trac == 'inca' .OR. type_trac == 'inco') THEN   ! config_inca='aero' ou 'chem'
-!>jyg
-! le module de chimie fournit les noms des traceurs
-! et les schemas d'advection associes. excepte pour ceux lus 
-! dans traceur.def 
-
-          DO iq=1,nqo+nqCO2
-
-             write(*,*) 'infotrac 237: iq=',iq
-             ! CRisi: ajout du nom du fluide transporteur
-             ! mais rester retro compatible
-             READ(90,'(I2,X,I2,X,A)',IOSTAT=IOstatus) hadv(iq),vadv(iq),tchaine
-             write(lunout,*) 'iq,hadv(iq),vadv(iq)=',iq,hadv(iq),vadv(iq)
-             write(lunout,*) 'tchaine=',trim(tchaine)
-             write(*,*) 'infotrac 238: IOstatus=',IOstatus
-             if (IOstatus.ne.0) then
-                CALL abort_gcm('infotrac_init','Pb dans la lecture de traceur.def',1)
-             endif
-             ! Y-a-t-il 1 ou 2 noms de traceurs? -> On regarde s'il y a un
-             ! espace ou pas au milieu de la chaine.
-             continu=.true.
-             nouveau_traceurdef=.false.
-             iiq=1
-
-             do while (continu)
-                if (tchaine(iiq:iiq).eq.' ') then
-                  nouveau_traceurdef=.true.
-                  continu=.false.
-                else if (iiq.lt.LEN_TRIM(tchaine)) then
-                  iiq=iiq+1
-                else
-                  continu=.false.
-                endif
-             enddo
-
-             write(*,*) 'iiq,nouveau_traceurdef=',iiq,nouveau_traceurdef
-
-             if (nouveau_traceurdef) then
-                write(lunout,*) 'C''est la nouvelle version de traceur.def'
-                tnom_0(iq)=tchaine(1:iiq-1)
-                tnom_transp(iq)=tchaine(iiq+1:)
-             else
-                write(lunout,*) 'C''est l''ancienne version de traceur.def'
-                write(lunout,*) 'On suppose que les traceurs sont tous d''air'
-                tnom_0(iq)=tchaine
-                tnom_transp(iq) = 'air'
-             endif
-
-             write(lunout,*) 'tnom_0(iq)=<',trim(tnom_0(iq)),'>'
-             write(lunout,*) 'tnom_transp(iq)=<',trim(tnom_transp(iq)),'>'
-
-          ENDDO  !DO iq=1,nqo
-          CLOSE(90)  
-
- 
+   CALL Init_chem_rep_trac(nbtr, nqo, tracers(:)%name)
+#endif
+
+!==============================================================================================================================
+! 2) Calculate nqtot, number of tracers needed (greater if advection schemes 20 or 30 have been chosen).
+!==============================================================================================================================
+   DO iq = 1, nqtrue
+      IF( hadv(iq)<20 .OR. (ANY(hadv(iq)==[20,30]) .AND. hadv(iq)==vadv(iq)) ) CYCLE
+      WRITE(msg1,'("The choice hadv=",i0,", vadv=",i0,a)')hadv(iq),vadv(iq),' for "'//TRIM(tracers(iq)%name)//'" is not available'
+      CALL abort_gcm(modname, TRIM(msg1), 1)
+   END DO
+   nqtot =    COUNT( hadv< 20 .AND. vadv< 20 ) &                     !--- No additional tracer
+         +  4*COUNT( hadv==20 .AND. vadv==20 ) &                     !--- 3  additional tracers
+         + 10*COUNT( hadv==30 .AND. vadv==30 )                       !--- 9  additional tracers
+
+   !--- More tracers due to the choice of advection scheme => assign total number of tracers
+   IF( nqtot /= nqtrue ) THEN
+      CALL msg('The choice of advection scheme for one or more tracers makes it necessary to add tracers')
+      CALL msg('The number of true tracers is '//TRIM(int2str(nqtrue)))
+      CALL msg('The total number of tracers needed is '//TRIM(int2str(nqtot)))
+   END IF
+
+!==============================================================================================================================
+! 3) Determine the advection scheme choice for water and tracers "iadv" and the fields long name, isAdvected.
+!     iadv = 1    "LMDZ-specific humidity transport" (for H2O vapour)          LMV
+!     iadv = 2    backward                           (for H2O liquid)          BAK
+!     iadv = 14   Van-Leer + specific humidity, modified by Francis Codron     VLH
+!     iadv = 10   Van-Leer (chosen for vapour and liquid water)                VL1
+!     iadv = 11   Van-Leer for hadv and PPM version (Monotonic) for vadv       VLP
+!     iadv = 12   Frederic Hourdin I                                           FH1
+!     iadv = 13   Frederic Hourdin II                                          FH2
+!     iadv = 16   Monotonic         PPM (Collela & Woodward 1984)              PPM
+!     iadv = 17   Semi-monotonic    PPM (overshoots allowed)                   PPS
+!     iadv = 18   Definite positive PPM (overshoots and undershoots allowed)   PPP
+!     iadv = 20   Slopes                                                       SLP
+!     iadv = 30   Prather                                                      PRA
+!
+!        In array q(ij,l,iq) : iq = 1/2[/3]    for vapour/liquid[/ice] water
+!        And optionaly:        iq = 3[4],nqtot for other tracers
+!==============================================================================================================================
+   ALLOCATE(ttr(nqtot))
+   jq = nqtrue+1; tracers(:)%iadv = -1
+   DO iq = 1, nqtrue
+      t1 => tracers(iq)
+
+      !--- VERIFY THE CHOICE OF ADVECTION SCHEME
+      iad = -1
+      IF(hadv(iq)     ==    vadv(iq)    ) iad = hadv(iq)
+      IF(hadv(iq)==10 .AND. vadv(iq)==16) iad = 11
+      WRITE(msg1,'("Bad choice of advection scheme for ",a,": hadv = ",i0,", vadv = ",i0)')TRIM(t1%name), hadv(iq), vadv(iq)
+      IF(iad == -1) CALL abort_gcm(modname, msg1, 1)
+
+      !--- SET FIELDS %longName, %iadv, %isAdvected, %isInPhysics
+      t1%longName   = t1%name; IF(iad > 0) t1%longName=TRIM(t1%name)//descrq(iad)
+      t1%iadv       = iad
+      t1%isAdvected = iad >= 0
+      t1%isInPhysics= delPhase(t1%gen0Name) /= 'H2O' &
+                          .OR. t1%component /= 'lmdz' !=== OTHER EXCEPTIONS TO BE ADDED: CO2i, SURSATURATED WATER CLOUD...
+      ttr(iq)       = t1
+
+      !--- DEFINE THE HIGHER ORDER TRACERS, IF ANY
+      nm = 0
+      IF(iad == 20) nm = 3                                           !--- 2nd order scheme
+      IF(iad == 30) nm = 9                                           !--- 3rd order scheme
+      IF(nm == 0) CYCLE                                              !--- No higher moments
+      ttr(jq+1:jq+nm)             = t1
+      ttr(jq+1:jq+nm)%name        = [(TRIM(t1%name)    //'-'//TRIM(suff(im)), im=1, nm) ]
+      ttr(jq+1:jq+nm)%parent      = [(TRIM(t1%parent)  //'-'//TRIM(suff(im)), im=1, nm) ]
+      ttr(jq+1:jq+nm)%longName    = [(TRIM(t1%longName)//'-'//TRIM(suff(im)), im=1, nm) ]
+      ttr(jq+1:jq+nm)%iadv        = [(-iad,    im=1, nm) ]
+      ttr(jq+1:jq+nm)%isAdvected  = [(.FALSE., im=1, nm) ]
+      jq = jq + nm
+   END DO
+   DEALLOCATE(hadv, vadv)
+   CALL MOVE_ALLOC(FROM=ttr, TO=tracers)
+
+   !--- SET FIELDS %iqParent, %nqChildren, %iGeneration, %iqDescen, %nqDescen
+   CALL indexUpdate(tracers)
+
+   !=== TEST ADVECTION SCHEME
+   DO iq=1,nqtot ; t1 => tracers(iq); iad = t1%iadv
+
+      !--- ONLY TESTED VALUES FOR TRACERS FOR NOW:               iadv = 14, 10 (and 0 for non-transported tracers)
+      IF(ALL([10,14,0] /= iad)) &
+         CALL abort_gcm(modname, 'Not tested for iadv='//TRIM(int2str(iad))//' ; 10 or 14 only are allowed !', 1)
+
+      !--- ONLY TESTED VALUES FOR PARENTS HAVING CHILDS FOR NOW: iadv = 14, 10 (PARENTS: TRACERS OF GENERATION 1)
+      IF(ALL([10,14] /= iad) .AND. t1%iGeneration == 1 .AND. ANY(tracers(:)%iGeneration > 1)) &
+         CALL abort_gcm(modname, 'iadv='//TRIM(int2str(iad))//' not implemented for parents ; 10 or 14 only are allowed !', 1)
+
+      !--- ONLY TESTED VALUES FOR CHILDS FOR NOW:                iadv = 10     (CHILDS:  TRACERS OF GENERATION GREATER THAN 1)
+      IF(fmsg('WARNING ! iadv='//TRIM(int2str(iad))//' not implemented for childs. Setting iadv=10 for "'//TRIM(t1%name)//'"',&
+         modname, iad /= 10 .AND. t1%iGeneration > 1)) t1%iadv = 10
+
+      !--- ONLY VALID SCHEME NUMBER FOR WATER VAPOUR:            iadv = 14
+      ll = t1%name /= addPhase('H2O','g')
+      IF(fmsg('WARNING ! iadv=14 is valid for water vapour only. Setting iadv=10 for "'//TRIM(t1%name)//'".', &
+         modname, iad == 14 .AND. ll))                 t1%iadv = 10
+   END DO
+
+   !=== READ PHYSICAL PARAMETERS FOR ISOTOPES ; DONE HERE BECAUSE dynetat0 AND iniacademic NEED "tnat" AND "alpha_ideal"
+   niso = 0; nzone = 0; nphas = nqo; ntiso = 0; isoCheck = .FALSE.
+   IF(readIsotopesFile()) CALL abort_gcm(modname, 'Problem when reading isotopes parameters', 1)
+
+   !--- Convection / boundary layer activation for all tracers
+   ALLOCATE(conv_flg(nbtr)); conv_flg(1:nbtr) = 1
+   ALLOCATE( pbl_flg(nbtr));  pbl_flg(1:nbtr) = 1
+
+   !--- Note: nqtottr can differ from nbtr when nmom/=0
+   nqtottr = nqtot - COUNT(delPhase(tracers%gen0Name) == 'H2O' .AND. tracers%component == 'lmdz')
+   IF(COUNT(tracers%iso_iName == 0) - COUNT(delPhase(tracers%name) == 'H2O' .AND. tracers%component == 'lmdz') /= nqtottr) &
+      CALL abort_gcm(modname, 'pb dans le calcul de nqtottr', 1)
+
+   !=== DISPLAY THE RESULTS
+   IF(prt_level > 1) THEN
+      CALL msg('nqo    = '//TRIM(int2str(nqo)),    modname)
+      CALL msg('nbtr   = '//TRIM(int2str(nbtr)),   modname)
+      CALL msg('nqtrue = '//TRIM(int2str(nqtrue)), modname)
+      CALL msg('nqtot  = '//TRIM(int2str(nqtot)),  modname)
+      CALL msg('niso   = '//TRIM(int2str(niso)),   modname)
+      CALL msg('ntiso  = '//TRIM(int2str(ntiso)),  modname)
 #ifdef INCA
-       CALL init_transport( &
-            hadv_inca, &
-            vadv_inca, &
-            conv_flg_inca, &
-            pbl_flg_inca,  &
-            solsym_inca)
-
-       conv_flg(1+nqCO2:nbtr) = conv_flg_inca
-       pbl_flg(1+nqCO2:nbtr) = pbl_flg_inca
-       solsym(1+nqCO2:nbtr) = solsym_inca
-
-       IF (type_trac == 'inco') THEN
-          conv_flg(1:nqCO2) = 1
-          pbl_flg(1:nqCO2) = 1
-          solsym(1:nqCO2) = 'CO2'
-       ENDIF
-#endif
-
-!jyg<
-       DO iq = nqo+nqCO2+1, nqtrue
-          hadv(iq) = hadv_inca(iq-nqo-nqCO2)
-          vadv(iq) = vadv_inca(iq-nqo-nqCO2) 
-          tnom_0(iq)=solsym_inca(iq-nqo-nqCO2)
-          tnom_transp(iq) = 'air'
-       END DO
-
-    ENDIF ! (type_trac == 'inca' or 'inco')
-
-!-----------------------------------------------------------------------
-!
-! 3) Verify if advection schema 20 or 30 choosen
-!    Calculate total number of tracers needed: nqtot
-!    Allocate variables depending on total number of tracers
-!-----------------------------------------------------------------------
-    new_iq=0
-    DO iq=1,nqtrue
-       ! Add tracers for certain advection schema
-       IF (hadv(iq)<20 .AND. vadv(iq)<20 ) THEN
-          new_iq=new_iq+1  ! no tracers added
-       ELSE IF (hadv(iq)==20 .AND. vadv(iq)==20 ) THEN
-          new_iq=new_iq+4  ! 3 tracers added
-       ELSE IF (hadv(iq)==30 .AND. vadv(iq)==30 ) THEN
-          new_iq=new_iq+10 ! 9 tracers added
-       ELSE
-          WRITE(lunout,*) trim(modname),': This choice of advection schema is not available',iq,hadv(iq),vadv(iq)
-          CALL abort_gcm('infotrac_init','Bad choice of advection schema - 1',1)
-       ENDIF
-    END DO
-    
-    IF (new_iq /= nqtrue) THEN
-       ! The choice of advection schema imposes more tracers
-       ! Assigne total number of tracers
-       nqtot = new_iq
-
-       WRITE(lunout,*) trim(modname),': The choice of advection schema for one or more tracers'
-       WRITE(lunout,*) 'makes it necessary to add tracers'
-       WRITE(lunout,*) trim(modname)//': ',nqtrue,' is the number of true tracers'
-       WRITE(lunout,*) trim(modname)//': ',nqtot, ' is the total number of tracers needed'
-
-    ELSE
-       ! The true number of tracers is also the total number
-       nqtot = nqtrue
-    ENDIF
-
-!
-! Allocate variables with total number of tracers, nqtot
-!
-    ALLOCATE(tname(nqtot), ttext(nqtot))
-    ALLOCATE(iadv(nqtot), niadv(nqtot))
-
-!-----------------------------------------------------------------------
-!
-! 4) Determine iadv, long and short name
-!
-!-----------------------------------------------------------------------
-    new_iq=0
-    DO iq=1,nqtrue
-       new_iq=new_iq+1
-
-       ! Verify choice of advection schema
-       IF (hadv(iq)==vadv(iq)) THEN
-          iadv(new_iq)=hadv(iq)
-       ELSE IF (hadv(iq)==10 .AND. vadv(iq)==16) THEN
-          iadv(new_iq)=11
-       ELSE
-          WRITE(lunout,*)trim(modname),': This choice of advection schema is not available',iq,hadv(iq),vadv(iq)
-
-          CALL abort_gcm('infotrac_init','Bad choice of advection schema - 2',1)
-       ENDIF
-      
-       str1=tnom_0(iq)
-       tname(new_iq)= tnom_0(iq)
-       IF (iadv(new_iq)==0) THEN
-          ttext(new_iq)=trim(str1)
-       ELSE
-          ttext(new_iq)=trim(tnom_0(iq))//descrq(iadv(new_iq))
-       ENDIF
-
-       ! schemas tenant compte des moments d'ordre superieur
-       str2=ttext(new_iq)
-       IF (iadv(new_iq)==20) THEN
-          DO jq=1,3
-             new_iq=new_iq+1
-             iadv(new_iq)=-20
-             ttext(new_iq)=trim(str2)//txts(jq)
-             tname(new_iq)=trim(str1)//txts(jq)
-          END DO
-       ELSE IF (iadv(new_iq)==30) THEN
-          DO jq=1,9
-             new_iq=new_iq+1
-             iadv(new_iq)=-30
-             ttext(new_iq)=trim(str2)//txtp(jq)
-             tname(new_iq)=trim(str1)//txtp(jq)
-          END DO
-       ENDIF
-    END DO
-
-!
-! Find vector keeping the correspodence between true and total tracers
-!
-    niadv(:)=0
-    iiq=0
-    DO iq=1,nqtot
-       IF(iadv(iq).GE.0) THEN
-          ! True tracer
-          iiq=iiq+1
-          niadv(iiq)=iq
-       ENDIF
-    END DO
-
-
-    WRITE(lunout,*) trim(modname),': Information stored in infotrac :'
-    WRITE(lunout,*) trim(modname),': iadv  niadv tname  ttext :'
-
-    DO iq=1,nqtot
-       WRITE(lunout,*) iadv(iq),niadv(iq), ' ',trim(tname(iq)),' ',trim(ttext(iq))
-    END DO
-
-!
-! Test for advection schema. 
-! This version of LMDZ only garantees iadv=10 and iadv=14 (14 only for water vapour) .
-!
-    DO iq=1,nqtot
-       IF (iadv(iq)/=10 .AND. iadv(iq)/=14 .AND. iadv(iq)/=0) THEN
-          WRITE(lunout,*)trim(modname),' STOP : The option iadv=',iadv(iq),' is not tested in this version of LMDZ'
-          CALL abort_gcm('infotrac_init','In this version only iadv=10 and iadv=14 is tested!',1)
-       ELSE IF (iadv(iq)==14 .AND. iq/=1) THEN
-          WRITE(lunout,*)trim(modname),'STOP : The option iadv=',iadv(iq),' is not tested in this version of LMDZ'
-          CALL abort_gcm('infotrac_init','In this version iadv=14 is only permitted for water vapour!',1)
-       ENDIF
-    END DO
-
-
-! CRisi: quels sont les traceurs fils et les traceurs peres.
-! initialiser tous les tableaux d'indices lies aux traceurs familiaux
-! + verifier que tous les peres sont ecrits en premieres positions
-    ALLOCATE(nqfils(nqtot),nqdesc(nqtot))    
-    ALLOCATE(iqfils(nqtot,nqtot))   
-    ALLOCATE(iqpere(nqtot))
-    nqperes=0
-    nqfils(:)=0
-    nqdesc(:)=0
-    iqfils(:,:)=0
-    iqpere(:)=0
-    nqdesc_tot=0    
-    DO iq=1,nqtot
-      if (tnom_transp(iq) == 'air') then
-        ! ceci est un traceur pere
-        WRITE(lunout,*) 'Le traceur',iq,', appele ',trim(tnom_0(iq)),', est un pere'
-        nqperes=nqperes+1
-        iqpere(iq)=0
-      else !if (tnom_transp(iq) == 'air') then
-        ! ceci est un fils. Qui est son pere?
-        WRITE(lunout,*) 'Le traceur',iq,', appele ',trim(tnom_0(iq)),', est un fils'
-        continu=.true.
-        ipere=1
-        do while (continu)           
-          if (tnom_transp(iq) == tnom_0(ipere)) then
-            ! Son pere est ipere
-            WRITE(lunout,*) 'Le traceur',iq,'appele ', &
-      &          trim(tnom_0(iq)),' est le fils de ',ipere,'appele ',trim(tnom_0(ipere))
-            if (iq.eq.ipere) then
-                CALL abort_gcm('infotrac_init','Un fils est son propre pere',1)
-            endif
-            nqfils(ipere)=nqfils(ipere)+1  
-            iqfils(nqfils(ipere),ipere)=iq
-            iqpere(iq)=ipere         
-            continu=.false.
-          else !if (tnom_transp(iq) == tnom_0(ipere)) then
-            ipere=ipere+1 
-            if (ipere.gt.nqtot) then
-                WRITE(lunout,*) 'Le traceur',iq,'appele ', &
-      &          trim(tnom_0(iq)),', est orphelin.'
-                CALL abort_gcm('infotrac_init','Un traceur est orphelin',1)
-            endif !if (ipere.gt.nqtot) then
-          endif !if (tnom_transp(iq) == tnom_0(ipere)) then
-        enddo !do while (continu)
-      endif !if (tnom_transp(iq) == 'air') then
-    enddo !DO iq=1,nqtot
-    WRITE(lunout,*) 'infotrac: nqperes=',nqperes    
-    WRITE(lunout,*) 'nqfils=',nqfils
-    WRITE(lunout,*) 'iqpere=',iqpere
-    WRITE(lunout,*) 'iqfils=',iqfils
-
-! Calculer le nombre de descendants a partir de iqfils et de nbfils
-    DO iq=1,nqtot    
-      generation=0
-      continu=.true.
-      ifils=iq
-      do while (continu)
-        ipere=iqpere(ifils)
-        if (ipere.gt.0) then
-         nqdesc(ipere)=nqdesc(ipere)+1   
-         nqdesc_tot=nqdesc_tot+1      
-         iqfils(nqdesc(ipere),ipere)=iq
-         ifils=ipere
-         generation=generation+1
-        else !if (ipere.gt.0) then
-         continu=.false.
-        endif !if (ipere.gt.0) then
-      enddo !do while (continu)    
-      WRITE(lunout,*) 'Le traceur ',iq,', appele ',trim(tnom_0(iq)),' est un traceur de generation: ',generation
-    enddo !DO iq=1,nqtot
-    WRITE(lunout,*) 'infotrac: nqdesc=',nqdesc
-    WRITE(lunout,*) 'iqfils=',iqfils
-    WRITE(lunout,*) 'nqdesc_tot=',nqdesc_tot
-
-! Interdire autres schemas que 10 pour les traceurs fils, et autres schemas
-! que 10 et 14 si des peres ont des fils
-    do iq=1,nqtot
-      if (iqpere(iq).gt.0) then
-        ! ce traceur a un pere qui n'est pas l'air
-        ! Seul le schema 10 est autorise
-        if (iadv(iq)/=10) then
-           WRITE(lunout,*)trim(modname),' STOP : The option iadv=',iadv(iq),' is not implemented for sons'
-          CALL abort_gcm('infotrac_init','Sons should be advected by scheme 10',1)
-        endif
-        ! Le traceur pere ne peut etre advecte que par schema 10 ou 14:
-        IF (iadv(iqpere(iq))/=10 .AND. iadv(iqpere(iq))/=14) THEN
-          WRITE(lunout,*)trim(modname),' STOP : The option iadv=',iadv(iq),' is not implemented for fathers'
-          CALL abort_gcm('infotrac_init','Fathers should be advected by scheme 10 ou 14',1)
-        endif !IF (iadv(iqpere(iq))/=10 .AND. iadv(iqpere(iq))/=14) THEN
-     endif !if (iqpere(iq).gt.0) the
-    enddo !do iq=1,nqtot
-
-
-
-! detecter quels sont les traceurs isotopiques parmi des traceurs
-    call infotrac_isoinit(tnom_0,nqtrue)
-
-!    if (ntraciso.gt.0) then
-! le 18 sep 2020: on enleve la condition ntraciso.gt.0 car nqtottr doit etre
-! connu meme si il n'y a pas d'isotopes!
-        write(lunout,*) 'infotrac 702: nbtr,ntraciso=',nbtr,ntraciso
-! retrancher les traceurs isotopiques de la liste des traceurs qui passent dans
-! phytrac
-        nbtr=nbtr-nqo*ntraciso
-
-! faire un tableau d'indice des traceurs qui passeront dans phytrac
-        nqtottr=nqtot-nqo*(1+ntraciso)
-        write(lunout,*) 'infotrac 704: nqtottr,nqtot,nqo=',nqtottr,nqtot,nqo
-        ! Rq: nqtottr n'est pas forcement egal a nbtr dans le cas ou new_iq /= nqtrue
-        ALLOCATE (itr_indice(nqtottr))  
-        itr_indice(:)=0  
-        itr=0
-        do iq=nqo+1, nqtot
-          if (iso_num(iq).eq.0) then
-            itr=itr+1
-            write(*,*) 'itr=',itr
-            itr_indice(itr)=iq
-          endif !if (iso_num(iq).eq.0) then
-        enddo
-        if (itr.ne.nqtottr) then
-            CALL abort_gcm('infotrac_init','pb dans le calcul de nqtottr',1)
-        endif
-        write(lunout,*) 'itr_indice=',itr_indice
-!    endif !if (ntraciso.gt.0) then
-
-!-----------------------------------------------------------------------
-! Finalize :
-!
-    DEALLOCATE(tnom_0, hadv, vadv,tnom_transp)
-
-    WRITE(lunout,*) 'infotrac init fin'
-
-  END SUBROUTINE infotrac_init
-
-  SUBROUTINE infotrac_isoinit(tnom_0,nqtrue)
-
-#ifdef CPP_IOIPSL
-  use IOIPSL
-#else
-  ! if not using IOIPSL, we still need to use (a local version of) getin
-  use ioipsl_getincom
-#endif
-  implicit none
- 
-    ! inputs
-    INTEGER,INTENT(IN) :: nqtrue
-    CHARACTER(len=*),INTENT(IN) :: tnom_0(nqtrue)
-    
-    ! locals    
-    CHARACTER(len=3), DIMENSION(niso_possibles) :: tnom_iso
-    INTEGER, ALLOCATABLE,DIMENSION(:,:) :: nb_iso,nb_traciso
-    INTEGER, ALLOCATABLE,DIMENSION(:) :: nb_isoind
-    INTEGER :: ntraceurs_zone_prec,iq,phase,ixt,iiso,izone
-    CHARACTER(len=tname_lenmax) :: tnom_trac
-    INCLUDE "iniprint.h"
-
-    tnom_iso=(/'eau','HDO','O18','O17','HTO'/)
-
-    ALLOCATE(nb_iso(niso_possibles,nqo))
-    ALLOCATE(nb_isoind(nqo))
-    ALLOCATE(nb_traciso(niso_possibles,nqo))
-    ALLOCATE(iso_num(nqtot))
-    ALLOCATE(iso_indnum(nqtot))
-    ALLOCATE(zone_num(nqtot))
-    ALLOCATE(phase_num(nqtot))
-     
-    iso_num(:)=0
-    iso_indnum(:)=0
-    zone_num(:)=0
-    phase_num(:)=0
-    indnum_fn_num(:)=0
-    use_iso(:)=.false.  
-    nb_iso(:,:)=0  
-    nb_isoind(:)=0     
-    nb_traciso(:,:)=0
-    niso=0
-    ntraceurs_zone=0  
-    ntraceurs_zone_prec=0
-    ntraciso=0
-
-    do iq=nqo+1,nqtot
-!       write(lunout,*) 'infotrac 569: iq,tnom_0(iq)=',iq,tnom_0(iq)
-       do phase=1,nqo    
-        do ixt= 1,niso_possibles   
-         tnom_trac=trim(tnom_0(phase))//'_'
-         tnom_trac=trim(tnom_trac)//trim(tnom_iso(ixt))
-!         write(*,*) 'phase,ixt,tnom_trac=',phase,ixt,tnom_trac      
-         IF (tnom_0(iq) == tnom_trac) then
-!          write(lunout,*) 'Ce traceur est un isotope'
-          nb_iso(ixt,phase)=nb_iso(ixt,phase)+1   
-          nb_isoind(phase)=nb_isoind(phase)+1   
-          iso_num(iq)=ixt
-          iso_indnum(iq)=nb_isoind(phase)
-          indnum_fn_num(ixt)=iso_indnum(iq)
-          phase_num(iq)=phase
-!          write(lunout,*) 'iso_num(iq)=',iso_num(iq)
-!          write(lunout,*) 'iso_indnum(iq)=',iso_indnum(iq)
-!          write(lunout,*) 'indnum_fn_num(ixt)=',indnum_fn_num(ixt)
-!          write(lunout,*) 'phase_num(iq)=',phase_num(iq)
-          goto 20
-         else if (iqpere(iq).gt.0) then          
-          if (tnom_0(iqpere(iq)) == tnom_trac) then
-!           write(lunout,*) 'Ce traceur est le fils d''un isotope'
-           ! c'est un traceur d'isotope
-           nb_traciso(ixt,phase)=nb_traciso(ixt,phase)+1
-           iso_num(iq)=ixt
-           iso_indnum(iq)=indnum_fn_num(ixt)
-           zone_num(iq)=nb_traciso(ixt,phase)
-           phase_num(iq)=phase
-!           write(lunout,*) 'iso_num(iq)=',iso_num(iq)
-!           write(lunout,*) 'phase_num(iq)=',phase_num(iq)
-!           write(lunout,*) 'zone_num(iq)=',zone_num(iq)
-           goto 20
-          endif !if (tnom_0(iqpere(iq)) == trim(tnom_0(phase))//trim(tnom_iso(ixt))) then
-         endif !IF (tnom_0(iq) == trim(tnom_0(phase))//trim(tnom_iso(ixt))) then
-        enddo !do ixt= niso_possibles
-       enddo !do phase=1,nqo
-  20   continue
-      enddo !do iq=1,nqtot
-
-!      write(lunout,*) 'iso_num=',iso_num
-!      write(lunout,*) 'iso_indnum=',iso_indnum
-!      write(lunout,*) 'zone_num=',zone_num  
-!      write(lunout,*) 'phase_num=',phase_num
-!      write(lunout,*) 'indnum_fn_num=',indnum_fn_num
-
-      do ixt= 1,niso_possibles  
-
-        if (nb_iso(ixt,1).eq.1) then
-          ! on verifie que toutes les phases ont le meme nombre de
-          ! traceurs
-          do phase=2,nqo
-            if (nb_iso(ixt,phase).ne.nb_iso(ixt,1)) then
-!              write(lunout,*) 'ixt,phase,nb_iso=',ixt,phase,nb_iso(ixt,phase)
-              CALL abort_gcm('infotrac_init','Phases must have same number of isotopes',1)
-            endif
-          enddo !do phase=2,nqo
-
-          niso=niso+1
-          use_iso(ixt)=.true.
-          ntraceurs_zone=nb_traciso(ixt,1)
-
-          ! on verifie que toutes les phases ont le meme nombre de
-          ! traceurs
-          do phase=2,nqo
-            if (nb_traciso(ixt,phase).ne.ntraceurs_zone) then
-              write(lunout,*) 'ixt,phase,nb_traciso=',ixt,phase,nb_traciso(ixt,phase)
-              write(lunout,*) 'ntraceurs_zone=',ntraceurs_zone
-              CALL abort_gcm('infotrac_init','Phases must have same number of tracers',1)
-            endif  
-          enddo  !do phase=2,nqo
-          ! on verifie que tous les isotopes ont le meme nombre de
-          ! traceurs
-          if (ntraceurs_zone_prec.gt.0) then               
-            if (ntraceurs_zone.eq.ntraceurs_zone_prec) then
-              ntraceurs_zone_prec=ntraceurs_zone
-            else !if (ntraceurs_zone.eq.ntraceurs_zone_prec) then
-              write(*,*) 'ntraceurs_zone_prec,ntraceurs_zone=',ntraceurs_zone_prec,ntraceurs_zone   
-              CALL abort_gcm('infotrac_init', &
-               &'Isotope tracers are not well defined in traceur.def',1)           
-            endif !if (ntraceurs_zone.eq.ntraceurs_zone_prec) then
-           endif !if (ntraceurs_zone_prec.gt.0) then
-
-        else if (nb_iso(ixt,1).ne.0) then
-           WRITE(lunout,*) 'nqo,ixt=',nqo,ixt
-           WRITE(lunout,*) 'nb_iso(ixt,1)=',nb_iso(ixt,1)    
-           CALL abort_gcm('infotrac_init','Isotopes are not well defined in traceur.def',1)     
-        endif   !if (nb_iso(ixt,1).eq.1) then       
-    enddo ! do ixt= niso_possibles
-
-    ! dimensions isotopique:
-    ntraciso=niso*(ntraceurs_zone+1)
-!    WRITE(lunout,*) 'niso=',niso
-!    WRITE(lunout,*) 'ntraceurs_zone,ntraciso=',ntraceurs_zone,ntraciso    
- 
-    ! flags isotopiques:
-    if (niso.gt.0) then
-        ok_isotopes=.true.
-    else
-        ok_isotopes=.false.
-    endif 
-!    WRITE(lunout,*) 'ok_isotopes=',ok_isotopes
- 
-    if (ok_isotopes) then
-        ok_iso_verif=.false.
-        call getin('ok_iso_verif',ok_iso_verif)
-        ok_init_iso=.false.
-        call getin('ok_init_iso',ok_init_iso)
-        tnat=(/1.0,155.76e-6,2005.2e-6,0.004/100.,0.0/)
-        alpha_ideal=(/1.0,1.01,1.006,1.003,1.0/)
-    endif !if (ok_isotopes) then  
-!    WRITE(lunout,*) 'ok_iso_verif=',ok_iso_verif
-!    WRITE(lunout,*) 'ok_init_iso=',ok_init_iso
-
-    if (ntraceurs_zone.gt.0) then
-        ok_isotrac=.true.
-    else
-        ok_isotrac=.false.
-    endif   
-!    WRITE(lunout,*) 'ok_isotrac=',ok_isotrac
-
-    ! remplissage du tableau iqiso(ntraciso,phase)
-    ALLOCATE(iqiso(ntraciso,nqo))   
-    iqiso(:,:)=0     
-    do iq=1,nqtot
-        if (iso_num(iq).gt.0) then
-          ixt=iso_indnum(iq)+zone_num(iq)*niso
-          iqiso(ixt,phase_num(iq))=iq
-        endif
-    enddo
-!    WRITE(lunout,*) 'iqiso=',iqiso
-
-    ! replissage du tableau index_trac(ntraceurs_zone,niso)
-    ALLOCATE(index_trac(ntraceurs_zone,niso))  
-    if (ok_isotrac) then
-        do iiso=1,niso
-          do izone=1,ntraceurs_zone 
-             index_trac(izone,iiso)=iiso+izone*niso
-          enddo
-        enddo
-    else !if (ok_isotrac) then      
-        index_trac(:,:)=0.0
-    endif !if (ok_isotrac) then
-!    write(lunout,*) 'index_trac=',index_trac    
-
-! Finalize :
-    DEALLOCATE(nb_iso)
-
-  END SUBROUTINE infotrac_isoinit
+      CALL msg('nqCO2  = '//TRIM(int2str(nqCO2)),  modname)
+      CALL msg('nqINCA = '//TRIM(int2str(nqINCA)), modname)
+#endif
+   END IF
+   t => tracers
+   CALL msg('Information stored in infotrac :', modname)
+   IF(dispTable('isssssssssiiiiiiiii', &
+      ['iq    ', 'name  ', 'lName ', 'gen0N ', 'parent', 'type  ', 'phase ', 'compon', 'isPhy ', 'isAdv ', &
+       'iadv  ', 'iGen  ', 'iqPar ', 'nqDes ', 'nqChld', 'iGroup', 'iName ', 'iZone ', 'iPhase'],          &
+      cat(t%name, t%longName, t%gen0Name, t%parent, t%type, t%phase, t%component, bool2str(t%isInPhysics), &
+                                                                                  bool2str(t%isAdvected)), &
+      cat([(iq, iq=1, nqtot)], t%iadv, t%iGeneration, t%iqParent, t%nqDescen, t%nqChildren, t%iso_iGroup,  &
+                  t%iso_iName, t%iso_iZone, t%iso_iPhase), nColMax=maxTableWidth, nHead=2, sub=modname))   &
+      CALL abort_gcm(modname, "problem with the tracers table content", 1)
+   IF(niso > 0) THEN
+      CALL msg('Where, for isotopes family "'//TRIM(isotope%parent)//'":', modname)
+      CALL msg('  isoKeys%name = '//strStack(isoKeys%name), modname)
+      CALL msg('  isoName = '//strStack(isoName),      modname)
+      CALL msg('  isoZone = '//strStack(isoZone),      modname)
+      CALL msg('  isoPhas = '//TRIM(isoPhas),          modname)
+   ELSE
+      CALL msg('No isotopes identified.', modname)
+   END IF
+   CALL msg('end', modname)
+
+END SUBROUTINE init_infotrac
 
 END MODULE infotrac
Index: LMDZ6/branches/Ocean_skin/libf/dyn3d_common/initdynav.F90
===================================================================
--- LMDZ6/branches/Ocean_skin/libf/dyn3d_common/initdynav.F90	(revision 4013)
+++ LMDZ6/branches/Ocean_skin/libf/dyn3d_common/initdynav.F90	(revision 4368)
@@ -6,5 +6,5 @@
   USE IOIPSL
 #endif
-  USE infotrac, ONLY : nqtot, ttext
+  USE infotrac, ONLY : nqtot
   use com_io_dyn_mod, only : histaveid,histvaveid,histuaveid, &
        dynhistave_file,dynhistvave_file,dynhistuave_file
@@ -158,5 +158,6 @@
 
   !        DO iq=1,nqtot
-  !          call histdef(histaveid, ttext(iq), ttext(iq), '-', &
+  !          call histdef(histaveid, tracers(iq)%name, &
+  !                                  tracers(iq)%longName, '-',  &
   !                  iip1, jjp1, thoriid, llm, 1, llm, zvertiid, &
   !                  32, 'ave(X)', t_ops, t_wrt)
Index: LMDZ6/branches/Ocean_skin/libf/dyn3d_common/inithist.F
===================================================================
--- LMDZ6/branches/Ocean_skin/libf/dyn3d_common/inithist.F	(revision 4013)
+++ LMDZ6/branches/Ocean_skin/libf/dyn3d_common/inithist.F	(revision 4368)
@@ -7,5 +7,5 @@
        USE IOIPSL
 #endif
-       USE infotrac, ONLY : nqtot, ttext
+       USE infotrac, ONLY : nqtot
        use com_io_dyn_mod, only : histid,histvid,histuid,               &
      &                        dynhist_file,dynhistv_file,dynhistu_file
@@ -157,5 +157,6 @@
 !
 !        DO iq=1,nqtot
-!          call histdef(histid, ttext(iq),  ttext(iq), '-',
+!          call histdef(histid, tracers(iq)%name, 
+!                               tracers(iq)%longName, '-',
 !     .             iip1, jjp1, thoriid, llm, 1, llm, zvertiid,
 !     .             32, 'inst(X)', t_ops, t_wrt)
Index: LMDZ6/branches/Ocean_skin/libf/dyn3d_common/iso_verif_dyn.F
===================================================================
--- LMDZ6/branches/Ocean_skin/libf/dyn3d_common/iso_verif_dyn.F	(revision 4013)
+++ LMDZ6/branches/Ocean_skin/libf/dyn3d_common/iso_verif_dyn.F	(revision 4368)
@@ -64,5 +64,5 @@
         function iso_verif_aberrant_nostop
      :           (x,iso,q,err_msg)
-        USE infotrac
+        USE infotrac, ONLY: isoName, getKey
         implicit none
         
@@ -74,5 +74,5 @@
         ! locals
         real qmin,deltaD
-        real deltaDmax,deltaDmin
+        real deltaDmax,deltaDmin,tnat
         parameter (qmin=1e-11)
         parameter (deltaDmax=200.0,deltaDmin=-999.9)
@@ -85,5 +85,10 @@
         ! verifier que HDO est raisonable
          if (q.gt.qmin) then
-             deltaD=(x/q/tnat(iso)-1)*1000
+             IF(getKey('tnat', tnat, isoName(iso))) THEN
+                  err_msg = 'Missing isotopic parameter "tnat"'
+                  iso_verif_aberrant_nostop=1
+                  RETURN
+             END IF
+             deltaD=(x/q/tnat-1)*1000
              if ((deltaD.gt.deltaDmax).or.(deltaD.lt.deltaDmin)) then
                   write(*,*) 'erreur detectee par iso_verif_aberrant:'
Index: LMDZ6/branches/Ocean_skin/libf/dyn3d_common/writedynav.F90
===================================================================
--- LMDZ6/branches/Ocean_skin/libf/dyn3d_common/writedynav.F90	(revision 4013)
+++ LMDZ6/branches/Ocean_skin/libf/dyn3d_common/writedynav.F90	(revision 4368)
@@ -6,5 +6,5 @@
   USE ioipsl
 #endif
-  USE infotrac, ONLY : nqtot, ttext
+  USE infotrac, ONLY : nqtot
   use com_io_dyn_mod, only : histaveid, histvaveid, histuaveid
   USE comconst_mod, ONLY: cpp
@@ -106,6 +106,6 @@
 
   !  DO iq=1, nqtot
-  !       call histwrite(histaveid, ttext(iq), itau_w, q(:, :, iq), &
-  !                   iip1*jjp1*llm, ndexu)
+  !       call histwrite(histaveid, tracers(iq)%longName, itau_w, &
+  !                   q(:, :, iq), iip1*jjp1*llm, ndexu)
   ! enddo
 
Index: LMDZ6/branches/Ocean_skin/libf/dyn3d_common/writehist.F
===================================================================
--- LMDZ6/branches/Ocean_skin/libf/dyn3d_common/writehist.F	(revision 4013)
+++ LMDZ6/branches/Ocean_skin/libf/dyn3d_common/writehist.F	(revision 4368)
@@ -7,5 +7,5 @@
       USE ioipsl
 #endif
-      USE infotrac, ONLY : nqtot, ttext
+      USE infotrac, ONLY : nqtot
       use com_io_dyn_mod, only : histid,histvid,histuid
       USE temps_mod, ONLY: itau_dyn
@@ -100,6 +100,6 @@
 C
 !        DO iq=1,nqtot
-!          call histwrite(histid, ttext(iq), itau_w, q(:,:,iq), 
-!     .                   iip1*jjp1*llm, ndexu)
+!          call histwrite(histid, tracers(iq)%longName, itau_w, 
+!     .                   q(:,:,iq), iip1*jjp1*llm, ndexu)
 !        enddo
 !C
