Index: LMDZ6/branches/cirrus/libf/dyn3d/check_isotopes.F90
===================================================================
--- LMDZ6/branches/cirrus/libf/dyn3d/check_isotopes.F90	(revision 5202)
+++ LMDZ6/branches/cirrus/libf/dyn3d/check_isotopes.F90	(revision 5203)
@@ -3,4 +3,9 @@
    USE infotrac,    ONLY: nqtot, niso, nphas, isotope, isoCheck, iqIsoPha, isoSelect, &
                           ntiso, iH2O, nzone, tracers, isoName,  itZonIso, getKey
+#ifdef CPP_IOIPSL
+   USE ioipsl,          ONLY: getin
+#else
+   USE ioipsl_getincom, ONLY: getin
+#endif
    IMPLICIT NONE
    include "dimensions.h"
@@ -20,8 +25,6 @@
                       deltaDmin =-999.0, &
                       ridicule  = 1e-12
-   INTEGER, SAVE :: iso_eau, iso_HDO, iso_O18, &
-                             iso_O17, iso_HTO
-   LOGICAL, SAVE :: first=.TRUE.
-   LOGICAL, PARAMETER :: tnat1=.TRUE.
+   INTEGER, SAVE :: iso_eau, iso_O17, iso_O18, iso_HDO, iso_HTO
+   LOGICAL, SAVE :: ltnat1, first=.TRUE.
 
    modname='check_isotopes'
@@ -30,14 +33,16 @@
    IF(niso == 0)        RETURN                   !--- No isotopes => finished
    IF(first) THEN
+      ltnat1 = .TRUE.; CALL getin('tnateq1', ltnat1)
+      ALLOCATE(tnat(niso))
       iso_eau = strIdx(isoName,'H216O')
+      iso_O17 = strIdx(isoName,'H217O')
+      iso_O18 = strIdx(isoName,'H218O')
       iso_HDO = strIdx(isoName,'HDO')
-      iso_O18 = strIdx(isoName,'H218O')
-      iso_O17 = strIdx(isoName,'H217O')
       iso_HTO = strIdx(isoName,'HTO')
-      if (tnat1) then
-              tnat(:)=1.0
-      else
+      IF(ltnat1) THEN
+         tnat(:)=1.0
+      ELSE
          IF(getKey('tnat', tnat)) CALL abort_gcm(modname, 'missing isotopic parameter', 1)
-      endif
+      END IF
       first = .FALSE.
    END IF
@@ -51,5 +56,5 @@
          DO k = 1, llm
             DO i = 1, ip1jmp1
-               IF(ABS(q(i,k,iq)) < borne) CYCLE
+               IF(ABS(q(i,k,iq)) <= borne) CYCLE
                WRITE(msg1,'(s,"(",i0,",",i0,",",i0,") = ",ES12.4)')TRIM(isoName(ixt)),i,k,iq,q(i,k,iq)
                CALL msg(msg1, modname)
Index: LMDZ6/branches/cirrus/libf/dyn3d/dynetat0.F90
===================================================================
--- LMDZ6/branches/cirrus/libf/dyn3d/dynetat0.F90	(revision 5202)
+++ LMDZ6/branches/cirrus/libf/dyn3d/dynetat0.F90	(revision 5203)
@@ -6,9 +6,9 @@
 ! Purpose: Initial state reading.
 !-------------------------------------------------------------------------------
-  USE infotrac,    ONLY: nqtot, tracers, niso, iqIsoPha, iH2O, isoName
+  USE infotrac,    ONLY: nqtot, tracers, niso, iqIsoPha, iH2O, isoName, &
+                         new2oldH2O, newHNO3, oldHNO3, getKey
   USE strings_mod, ONLY: maxlen, msg, strStack, real2str, int2str
   USE netcdf,      ONLY: NF90_OPEN,  NF90_NOWRITE, NF90_INQ_VARID, &
                          NF90_CLOSE, NF90_GET_VAR, NF90_NoErr
-  USE readTracFiles_mod, ONLY: new2oldH2O, newHNO3, oldHNO3, getKey
   USE control_mod, ONLY: planet_type
   USE assert_eq_m, ONLY: assert_eq
@@ -19,4 +19,9 @@
   USE temps_mod, ONLY: annee_ref, day_ini, day_ref, itau_dyn, start_time
   USE ener_mod, ONLY: etot0,ptot0,ztot0,stot0,ang0
+#ifdef CPP_IOIPSL
+  USE IOIPSL,   ONLY: getin
+#else
+  USE ioipsl_getincom, ONLY: getin
+#endif
 
   IMPLICIT NONE
@@ -42,6 +47,5 @@
   INTEGER :: iq, fID, vID, idecal, iqParent, iName, iZone, iPhase
   REAL    :: time, tnat, alpha_ideal, tab_cntrl(length)    !--- RUN PARAMS TABLE
-  LOGICAL :: lSkip, ll
-  LOGICAL,PARAMETER :: tnat1=.TRUE.
+  LOGICAL :: lSkip, ll, ltnat1
 !-------------------------------------------------------------------------------
   modname="dynetat0"
@@ -116,5 +120,5 @@
   var="temps"
   IF(NF90_INQ_VARID(fID,var,vID)/=NF90_NoErr) THEN
-    CALL msg('missing field <temps> ; trying with <Time>', modname)
+    CALL msg('Missing field <temps> ; trying with <Time>', modname)
     var="Time"
     CALL err(NF90_INQ_VARID(fID,var,vID),"inq",var)
@@ -133,4 +137,5 @@
   ll = NF90_INQ_VARID(fID, 'HNO3tot', vID) /= NF90_NoErr                                 !--- DETECT OLD REPRO start.nc FILE
 #endif
+  ltnat1 = .TRUE.; CALL getin('tnateq1', ltnat1)
   DO iq=1,nqtot
     var = tracers(iq)%name
@@ -148,5 +153,5 @@
     !--------------------------------------------------------------------------------------------------------------------------
     ELSE IF(NF90_INQ_VARID(fID, oldVar, vID) == NF90_NoErr) THEN                         !=== TRY WITH ALTERNATE NAME
-      CALL msg('Tracer <'//TRIM(var)//'> is missing => initialized to <'//TRIM(oldVar)//'>', modname)
+      CALL msg('Missing tracer <'//TRIM(var)//'> => initialized to <'//TRIM(oldVar)//'>', modname)
       CALL err(NF90_GET_VAR(fID,vID,q(:,:,:,iq)),"get",oldVar)
     !--------------------------------------------------------------------------------------------------------------------------
@@ -156,16 +161,16 @@
       iqParent = tracers(iq)%iqParent
       IF(tracers(iq)%iso_iZone == 0) THEN
-         if (tnat1) then
-                 tnat=1.0
-                 alpha_ideal=1.0
-                 write(*,*) 'attention dans dynetat0: les alpha_ideal sont a 1'
-         else
+         IF(ltnat1) THEN
+            tnat = 1.0
+            alpha_ideal = 1.0
+            CALL msg(' !!!  Beware: alpha_ideal put to 1  !!!', modname)
+         ELSE
           IF(getKey('tnat', tnat, isoName(iName)) .OR. getKey('alpha', alpha_ideal, isoName(iName))) &
             CALL abort_gcm(TRIM(modname), 'missing isotopic parameters', 1)
-         endif
-         CALL msg('Tracer <'//TRIM(var)//'> is missing => initialized with a simplified Rayleigh distillation law.', modname)
+         END IF
+         CALL msg('Missing tracer <'//TRIM(var)//'> => initialized with a simplified Rayleigh distillation law.', modname)
          q(:,:,:,iq) = q(:,:,:,iqParent)*tnat*(q(:,:,:,iqParent)/30.e-3)**(alpha_ideal-1.)
       ELSE
-         CALL msg('Tracer <'//TRIM(var)//'> is missing => initialized to its parent isotope concentration.', modname)
+         CALL msg('Missing tracer <'//TRIM(var)//'> => initialized to its parent isotope concentration.', modname)
          ! Camille 9 mars 2023: attention!! seuls les tags qui correspondent à
          ! izone=izone_init (définie dans isotrac_mod) sont initialisés comme
@@ -181,5 +186,5 @@
     !--------------------------------------------------------------------------------------------------------------------------
     ELSE                                                                                 !=== MISSING: SET TO 0
-      CALL msg('Tracer <'//TRIM(var)//'> is missing => initialized to zero', modname)
+      CALL msg('Missing tracer <'//TRIM(var)//'> => initialized to zero', modname)
       q(:,:,:,iq)=0.
     !--------------------------------------------------------------------------------------------------------------------------
Index: LMDZ6/branches/cirrus/libf/dyn3d/iniacademic.F90
===================================================================
--- LMDZ6/branches/cirrus/libf/dyn3d/iniacademic.F90	(revision 5202)
+++ LMDZ6/branches/cirrus/libf/dyn3d/iniacademic.F90	(revision 5203)
@@ -5,5 +5,5 @@
 
   USE filtreg_mod, ONLY: inifilr
-  USE infotrac,    ONLY: nqtot, niso, iqIsoPha, tracers, getKey, isoName
+  USE infotrac,    ONLY: nqtot, niso, iqIsoPha, tracers, getKey, isoName, addPhase
   USE control_mod, ONLY: day_step,planet_type
   use exner_hyb_m, only: exner_hyb
@@ -21,5 +21,4 @@
   USE temps_mod, ONLY: annee_ref, day_ini, day_ref
   USE ener_mod, ONLY: etot0,ptot0,ztot0,stot0,ang0
-  USE readTracFiles_mod, ONLY: addPhase
   use netcdf, only : NF90_NOWRITE,NF90_OPEN,NF90_NOERR,NF90_INQ_VARID
   use netcdf, only : NF90_CLOSE, NF90_GET_VAR
@@ -80,5 +79,5 @@
 
   REAL zdtvr, tnat, alpha_ideal
-  LOGICAL,PARAMETER :: tnat1=.true.
+  LOGICAL :: ltnat1
   
   character(len=*),parameter :: modname="iniacademic"
@@ -309,4 +308,5 @@
         ! bulk initialization of tracers
         if (planet_type=="earth") then
+           ltnat1 = .TRUE.; CALL getin('tnateq1', ltnat1)
            ! Earth: first two tracers will be water
            do iq=1,nqtot
@@ -322,12 +322,12 @@
               iqParent = tracers(iq)%iqParent
               IF(tracers(iq)%iso_iZone == 0) THEN
-                 if (tnat1) then
-                         tnat=1.0
-                         alpha_ideal=1.0
-                         write(*,*) 'Attention dans iniacademic: alpha_ideal=1'
-                 else
+                 IF(ltnat1) THEN
+                    tnat = 1.0
+                    alpha_ideal = 1.0
+                    WRITE(lunout, *)'In '//TRIM(modname)//': !!!  Beware: alpha_ideal put to 1  !!!'
+                 ELSE
                     IF(getKey('tnat', tnat, isoName(iName)) .OR. getKey('alpha', alpha_ideal, isoName(iName))) &
                     CALL abort_gcm(TRIM(modname), 'missing isotopic parameters', 1)
-                 endif
+                 END IF
                  q(:,:,iq) = q(:,:,iqParent)*tnat*(q(:,:,iqParent)/30.e-3)**(alpha_ideal-1.)
               ELSE !IF(tracers(iq)%iso_iZone == 0) THEN
Index: LMDZ6/branches/cirrus/libf/dyn3d/qminimum.F
===================================================================
--- LMDZ6/branches/cirrus/libf/dyn3d/qminimum.F	(revision 5202)
+++ LMDZ6/branches/cirrus/libf/dyn3d/qminimum.F	(revision 5203)
@@ -4,7 +4,6 @@
       SUBROUTINE qminimum( q,nqtot,deltap )
 
-      USE infotrac, ONLY: niso, ntiso,iqIsoPha, tracers
+      USE infotrac, ONLY: niso, ntiso, iqIsoPha, tracers, addPhase
       USE strings_mod, ONLY: strIdx
-      USE readTracFiles_mod, ONLY: addPhase
       IMPLICIT none
 c
Index: LMDZ6/branches/cirrus/libf/dyn3d_common/infotrac.F90
===================================================================
--- LMDZ6/branches/cirrus/libf/dyn3d_common/infotrac.F90	(revision 5202)
+++ LMDZ6/branches/cirrus/libf/dyn3d_common/infotrac.F90	(revision 5203)
@@ -7,4 +7,5 @@
         delPhase, niso, getKey, isot_type, processIsotopes,  isotope, maxTableWidth, iqIsoPha, nphas, ixIso, isoPhas, &
         addPhase, iH2O, addKey, isoSelect, testTracersFiles, isoKeys, indexUpdate,  iqWIsoPha, nbIso, ntiso, isoName, isoCheck
+   USE readTracFiles_mod, ONLY: new2oldH2O, newHNO3, oldHNO3
    IMPLICIT NONE
 
@@ -16,4 +17,6 @@
    PUBLIC :: nqtot,   nbtr,   nqo,   nqCO2,   nqtottr      !--- Main dimensions
    PUBLIC :: conv_flg, pbl_flg                             !--- Convection & boundary layer activation keys
+   PUBLIC :: new2oldH2O, newHNO3, oldHNO3                  !--- For backwards compatibility in dynetat0
+   PUBLIC :: addPhase, delPhase                            !--- Add/remove the phase from the name of a tracer
 
    !=== FOR ISOTOPES: General
@@ -21,12 +24,12 @@
    PUBLIC :: isoSelect, ixIso                              !--- Isotopes family selection tool + selected family index
    !=== FOR ISOTOPES: Specific to water
-   PUBLIC :: iH2O                                          !--- H2O isotopes class index
+   PUBLIC :: iH2O                                          !--- Value of "ixIso" for "H2O" isotopes class
    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 :: isotope                                       !--- Selected isotopes database (argument of getKey)
+   PUBLIC :: isoKeys, isoName, isoZone, isoPhas            !--- Isotopes keys & names, tagging zones names, phases
+   PUBLIC ::    niso,   ntiso,   nzone,   nphas            !--- Number of   "   "
+   PUBLIC :: itZonIso                                      !--- index "it" in "isoName(1:niso)" = f(tagging idx, isotope idx)
+   PUBLIC :: iqIsoPha                                      !--- index "iq" in "qx"              = f(isotope idx,   phase idx)
    PUBLIC :: isoCheck                                      !--- Run isotopes checking routines
    !=== FOR BOTH TRACERS AND ISOTOPES
@@ -36,5 +39,5 @@
 !  |--------------------+-----------------------+-----------------+---------------+----------------------------|
 !  | water in different |    water tagging      |  water isotopes | other tracers | additional tracers moments |
-!  | phases: H2O_[gls]  |      isotopes         |                 |               |  for higher order schemes  |
+!  | phases: H2O_[glsrb]|      isotopes         |                 |               |  for higher order schemes  |
 !  |--------------------+-----------------------+-----------------+---------------+----------------------------|
 !  |                    |                       |                 |               |                            |
@@ -61,9 +64,11 @@
 !  |-------------+------------------------------------------------------+-------------+------------------------+
 !  | name        | Name (short)                                         | tname       |                        |
+!  | keys        | key/val pairs accessible with "getKey" routine       | /           |                        |
 !  | 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]              |
+!  | phase       | Phases list ("g"as / "l"iquid / "s"olid              |             | [g|l|s|r|b]            |
+!  |             |              "r"(cloud) / "b"lowing)                 | /           |                        |
 !  | component   | Name(s) of the merged/cumulated section(s)           | /           | coma-separated names   |
 !  | iGeneration | Generation (>=1)                                     | /           |                        |
@@ -72,8 +77,7 @@
 !  | 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  |
+!  | 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                 |
@@ -91,5 +95,5 @@
 !  | 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  |
+!  | phase  | nphas  | Phases                     list + number         |                    | [g|l|s|r|b] 1:5 |
 !  | 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         |
@@ -99,14 +103,14 @@
 
    !=== 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
+   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
+   CHARACTER(LEN=maxlen), SAVE :: type_trac                     !--- Keyword for tracers type(s)
 
    !=== VARIABLES FOR INCA
-   INTEGER,               SAVE, ALLOCATABLE :: conv_flg(:), &   !--- Convection     activation ; needed for INCA        (nbtr)
-                                                pbl_flg(:)      !--- Boundary layer activation ; needed for INCA        (nbtr)
+   INTEGER, DIMENSION(:), SAVE, ALLOCATABLE :: &
+                    conv_flg, pbl_flg                           !--- Convection / boundary layer activation (nbtr)
 
 CONTAINS
@@ -115,5 +119,5 @@
    USE control_mod, ONLY: planet_type
 #ifdef REPROBUS
-   USE CHEM_REP,    ONLY: Init_chem_rep_trac
+   USE CHEM_REP, ONLY: Init_chem_rep_trac
 #endif
    IMPLICIT NONE
@@ -160,6 +164,5 @@
    TYPE(trac_type), ALLOCATABLE, TARGET :: ttr(:)
    TYPE(trac_type), POINTER             :: t1, t(:)
-   CHARACTER(LEN=maxlen),   ALLOCATABLE :: types_trac(:)  !--- Keyword for tracers type(s), parsed version
-
+   CHARACTER(LEN=maxlen),   ALLOCATABLE :: types_trac(:)             !--- Keywords for tracers type(s), parsed version
    CHARACTER(LEN=*), PARAMETER :: modname="init_infotrac"
 !------------------------------------------------------------------------------------------------------------------------------
@@ -171,13 +174,12 @@
    descrq(10:20) = ['VL1','VLP','FH1','FH2','VLH','   ','PPM','PPS','PPP','   ','SLP']
    descrq(30)    =  'PRA'
-    
-   CALL msg('type_trac = "'//TRIM(type_trac)//'"', modname)
 
    lerr=strParse(type_trac, '|', types_trac, n=nt)
    IF (nt .GT. 1) THEN
       IF (nt .GT. 2) CALL abort_gcm(modname, 'you need to modify type_trac, this version is not supported by lmdz', 1)
-      if (nt .EQ. 2) type_trac=types_trac(2)
+      IF (nt .EQ. 2) type_trac=types_trac(2)
    ENDIF
 
+   CALL msg('type_trac = "'//TRIM(type_trac)//'"', modname)
 
    
@@ -213,19 +215,20 @@
 
 !==============================================================================================================================
-! 1) Get the numbers of: true (first order only) tracers "nqtrue", water tracers "nqo" (vapor/liquid/solid)
-!==============================================================================================================================
-   texp = type_trac                                                  !=== EXPANDED VERSION OF "type_trac", WITH "|" SEPARATOR
+! 0) DETERMINE THE TYPE OF THE INPUT TRACERS DESCRIPTION FILE AND READ IT
+!==============================================================================================================================
+   texp = type_trac                                                            !=== EXPANDED (WITH "|" SEPARATOR) "type_trac"
    IF(texp == 'inco') texp = 'co2i|inca'
    IF(texp /= 'lmdz') texp = 'lmdz|'//TRIM(texp)
-
-   !=== DETERMINE THE TYPE OF THE INPUT TRACERS DESCRIPTION FILE
    IF(testTracersFiles(modname, texp, fType, .TRUE.)) CALL abort_gcm(modname, 'problem with tracers file(s)',1)
    ttp = type_trac; IF(fType /= 1) ttp = texp
-
    IF(readTracersFiles(ttp, lRepr=type_trac=='repr')) 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"
+
+!==============================================================================================================================
+! 1) Get various numbers: "nqtrue" (first order only tracers), "nqo" (water phases), 'nbtr' (tracers passed to physics), etc.
+!==============================================================================================================================
+   !---------------------------------------------------------------------------------------------------------------------------
+   IF(fType == 0) CALL abort_gcm(modname, 'Missing "traceur.def", "tracer.def" or "tracer_<keyword>.def tracers file.',1)
+   !---------------------------------------------------------------------------------------------------------------------------
+   IF(fType == 1 .AND. ANY(['inca', 'inco'] == type_trac)) THEN      !=== FOUND OLD STYLE INCA "traceur.def"
    !---------------------------------------------------------------------------------------------------------------------------
 #ifdef INCA
@@ -264,14 +267,13 @@
    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
+   nqtrue = SIZE(tracers)                                                                               !--- "true" tracers
+   nqo    =      COUNT(tracers(:)%component == 'lmdz' .AND. delPhase(tracers(:)%name)     == 'H2O')     !--- Water phases
+   nbtr = nqtrue-COUNT(tracers(:)%component == 'lmdz' .AND. delPhase(tracers(:)%gen0Name) == 'H2O')     !--- Passed to phytrac
+   nqCO2  =      COUNT( [type_trac == 'inco', type_trac == 'co2i'] )
 #ifdef INCA
-      nqINCA = COUNT(tracers(:)%component == 'inca')
-#endif
-      lerr = getKey('hadv', hadv, ky=tracers(:)%keys)
-      lerr = getKey('vadv', vadv, ky=tracers(:)%keys)
+   nqINCA =      COUNT(tracers(:)%component == 'inca')
+#endif
+   IF(getKey('hadv', hadv, ky=tracers(:)%keys)) CALL abort_gcm(modname, 'missing key "hadv"', 1)
+   IF(getKey('vadv', vadv, ky=tracers(:)%keys)) CALL abort_gcm(modname, 'missing key "vadv"', 1)
    !---------------------------------------------------------------------------------------------------------------------------
    END IF
@@ -279,8 +281,7 @@
 
 #ifdef REPROBUS
-   !--- Transfert the number of tracers to Reprobus
    CALL Init_chem_rep_trac(nbtr, nqo, tracers(:)%name)
-
-#endif
+#endif
+
 !==============================================================================================================================
 ! 2) Calculate nqtot, number of tracers needed (greater if advection schemes 20 or 30 have been chosen).
@@ -332,10 +333,9 @@
       IF(iad == -1) CALL abort_gcm(modname, msg1, 1)
 
-      !--- SET FIELDS %longName, %iadv, %isAdvected, %isInPhysics
+      !--- 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...
+      t1%isInPhysics= delPhase(t1%gen0Name) /= 'H2O' .OR. t1%component /= 'lmdz' !=== MORE EXCEPTIONS ? CO2i, SURSAT CLOUD H2O
       ttr(iq)       = t1
 
@@ -347,4 +347,5 @@
       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)%gen0Name    = [ (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) ]
@@ -356,9 +357,9 @@
    CALL MOVE_ALLOC(FROM=ttr, TO=tracers)
 
-   !--- SET FIELDS %iqParent, %nqChildren, %iGeneration, %iqDescen, %nqDescen
+   !--- SET FIELDS iqParent, iqDescen, nqDescen, nqChildren
    IF(indexUpdate(tracers)) CALL abort_gcm(modname, 'problem with tracers indices update', 1)
 
    !=== TEST ADVECTION SCHEME
-   DO iq=1,nqtot ; t1 => tracers(iq); iad = t1%iadv
+   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)
@@ -405,10 +406,8 @@
 #endif
    t => tracers
-   CALL msg('Information stored in infotrac :', modname)
-
-   IF(dispTable('isssssssssiiiiiiiii', ['iq  ', 'name', 'lNam', 'g0Nm', 'prnt', 'type', 'phas', 'comp',    &
-                'isPh', 'isAd', 'iadv', 'iGen', 'iqPr', 'nqDe', 'nqCh', 'iGrp', 'iNam', 'iZon', 'iPha'],   &
-      cat(t%name, t%longName, t%gen0Name, t%parent, t%type, t%phase, t%component, bool2str(t%isInPhysics), &
-                                                                                  bool2str(t%isAdvected)), &
+   CALL msg('Information stored in '//TRIM(modname)//': ', modname)
+   IF(dispTable('isssssssiiiiiiiii', ['iq  ', 'name', 'lNam', 'g0Nm', 'prnt', 'type', 'phas', 'comp',      &
+                              'iAdv', 'iGen', 'iqPr', 'nqDe', 'nqCh', 'iGrp', 'iNam', 'iZon', 'iPha'],     &
+      cat(t%name, t%longName, t%gen0Name, t%parent, t%type, t%phase, t%component),                         &
       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))   &
Index: LMDZ6/branches/cirrus/libf/dyn3d_common/iso_verif_dyn.F
===================================================================
--- LMDZ6/branches/cirrus/libf/dyn3d_common/iso_verif_dyn.F	(revision 5202)
+++ LMDZ6/branches/cirrus/libf/dyn3d_common/iso_verif_dyn.F	(revision 5203)
@@ -64,4 +64,9 @@
         function iso_verif_aberrant_nostop
      :           (x,iso,q,err_msg)
+#ifdef CPP_IOIPSL
+        USE IOIPSL, ONLY: getin
+#else
+        USE ioipsl_getincom, ONLY: getin
+#endif
         USE infotrac, ONLY: isoName, getKey
         implicit none
@@ -77,13 +82,21 @@
         parameter (qmin=1e-11)
         parameter (deltaDmax=200.0,deltaDmin=-999.9)
+        LOGICAL, SAVE :: ltnat1
+        LOGICAL, SAVE :: lFirst=.TRUE.
 
         ! output
         integer iso_verif_aberrant_nostop
 
+        IF(lFirst) THEN
+           ltnat1 = .TRUE.; CALL getin('tnateq1', ltnat1)
+           lFirst = .FALSE.
+        END IF
         iso_verif_aberrant_nostop=0
 
         ! verifier que HDO est raisonable
          if (q.gt.qmin) then
-             IF(getKey('tnat', tnat, isoName(iso))) THEN
+             IF(ltnat1) THEN
+                tnat = 1.0
+             ELSE IF(getKey('tnat', tnat, isoName(iso))) THEN
                   err_msg = 'Missing isotopic parameter "tnat"'
                   iso_verif_aberrant_nostop=1
Index: LMDZ6/branches/cirrus/libf/dyn3dmem/check_isotopes_loc.F90
===================================================================
--- LMDZ6/branches/cirrus/libf/dyn3dmem/check_isotopes_loc.F90	(revision 5202)
+++ LMDZ6/branches/cirrus/libf/dyn3dmem/check_isotopes_loc.F90	(revision 5203)
@@ -4,4 +4,9 @@
    USE infotrac,    ONLY: nqtot, niso, nphas, isotope, isoCheck, iqIsoPha, isoSelect, &
                           ntiso, iH2O, nzone, tracers, isoName,  itZonIso, getKey
+#ifdef CPP_IOIPSL
+   USE ioipsl,          ONLY: getin
+#else
+   USE ioipsl_getincom, ONLY: getin
+#endif
    IMPLICIT NONE
    include "dimensions.h"
@@ -21,8 +26,8 @@
                       deltaDmin =-999.0, &
                       ridicule  = 1e-12
-   INTEGER, SAVE :: iso_eau, iso_HDO, iso_O18, & !--- OpenMP shared variables
-                             iso_O17, iso_HTO
+   INTEGER, SAVE :: iso_eau, iso_HDO, iso_O18, iso_O17, iso_HTO
+!$OMP THREADPRIVATE(iso_eau, iso_HDO, iso_O18, iso_O17, iso_HTO)
+   LOGICAL       :: ltnat1
    LOGICAL, SAVE :: first=.TRUE.
-   LOGICAL, PARAMETER :: tnat1=.TRUE.
 !$OMP THREADPRIVATE(first)
 
@@ -33,14 +38,16 @@
    IF(first) THEN
 !$OMP MASTER
+      ltnat1 = .TRUE.; CALL getin('tnateq1', ltnat1)
+      ALLOCATE(tnat(niso))
       iso_eau = strIdx(isoName,'H216O')
+      iso_O17 = strIdx(isoName,'H217O')
+      iso_O18 = strIdx(isoName,'H218O')
       iso_HDO = strIdx(isoName,'HDO')
-      iso_O18 = strIdx(isoName,'H218O')
-      iso_O17 = strIdx(isoName,'H217O')
       iso_HTO = strIdx(isoName,'HTO')
-      if (tnat1) then
-              tnat(:)=1.0
-      else
+      IF(ltnat1) THEN
+         tnat(:)=1.0
+      ELSE
          IF(getKey('tnat', tnat)) CALL abort_gcm(modname, 'missing isotopic parameter', 1)
-      endif
+      END IF
 !$OMP END MASTER
 !$OMP BARRIER
@@ -57,5 +64,5 @@
          DO k = 1, llm
             DO i = ijb, ije
-               IF(ABS(q(i,k,iq))<=borne) CYCLE
+               IF(ABS(q(i,k,iq)) <= borne) CYCLE
                WRITE(msg1,'(s,"(",i0,",",i0,",",i0,") = ",ES12.4)')TRIM(isoName(ixt)),i,k,iq,q(i,k,iq)
                CALL msg(msg1, modname)
Index: LMDZ6/branches/cirrus/libf/dyn3dmem/dynetat0_loc.F90
===================================================================
--- LMDZ6/branches/cirrus/libf/dyn3dmem/dynetat0_loc.F90	(revision 5202)
+++ LMDZ6/branches/cirrus/libf/dyn3dmem/dynetat0_loc.F90	(revision 5203)
@@ -7,9 +7,9 @@
 !-------------------------------------------------------------------------------
   USE parallel_lmdz
-  USE infotrac,    ONLY: nqtot, tracers, niso, iqIsoPha, iH2O, isoName
+  USE infotrac,    ONLY: nqtot, tracers, niso, iqIsoPha, iH2O, isoName, &
+                         new2oldH2O, newHNO3, oldHNO3, getKey
   USE strings_mod, ONLY: maxlen, msg, strStack, real2str, int2str, strIdx
   USE netcdf,      ONLY: NF90_OPEN,  NF90_NOWRITE, NF90_INQUIRE_DIMENSION, NF90_INQ_VARID, &
                          NF90_CLOSE, NF90_GET_VAR, NF90_INQUIRE_VARIABLE,  NF90_NoErr
-  USE readTracFiles_mod, ONLY: new2oldH2O, newHNO3, oldHNO3, getKey
   USE control_mod, ONLY: planet_type
   USE assert_eq_m, ONLY: assert_eq
@@ -20,4 +20,9 @@
   USE temps_mod, ONLY: annee_ref, day_ini, day_ref, itau_dyn, start_time
   USE ener_mod, ONLY: etot0,ptot0,ztot0,stot0,ang0
+#ifdef CPP_IOIPSL
+  USE IOIPSL,   ONLY: getin
+#else
+  USE ioipsl_getincom, ONLY: getin
+#endif
 
   IMPLICIT NONE
@@ -47,6 +52,5 @@
   REAL,             ALLOCATABLE :: ucov_glo(:,:),    q_glo(:,:), phis_glo(:)
   REAL,             ALLOCATABLE :: teta_glo(:,:)
-  LOGICAL :: lSkip, ll
-  LOGICAL,PARAMETER :: tnat1=.TRUE.
+  LOGICAL :: lSkip, ll, ltnat1
 !-------------------------------------------------------------------------------
   modname="dynetat0_loc"
@@ -158,4 +162,5 @@
   ll = NF90_INQ_VARID(fID, 'HNO3tot', vID) /= NF90_NoErr                                 !--- DETECT OLD REPRO start.nc FILE
 #endif
+  ltnat1 = .TRUE.; CALL getin('tnateq1', ltnat1)
   DO iq=1,nqtot
     var = tracers(iq)%name
@@ -173,5 +178,5 @@
     !--------------------------------------------------------------------------------------------------------------------------
     ELSE IF(NF90_INQ_VARID(fID, oldVar, vID) == NF90_NoErr) THEN                         !=== TRY WITH ALTERNATE NAME
-      CALL msg('Tracer <'//TRIM(var)//'> is missing => initialized to <'//TRIM(oldVar)//'>', modname)
+      CALL msg('Missing tracer <'//TRIM(var)//'> => initialized to <'//TRIM(oldVar)//'>', modname)
       CALL get_var2(oldVar, q_glo); q(ijb_u:ije_u,:,iq)=q_glo(ijb_u:ije_u,:)
     !--------------------------------------------------------------------------------------------------------------------------
@@ -181,18 +186,18 @@
       iqParent = tracers(iq)%iqParent
       IF(tracers(iq)%iso_iZone == 0) THEN
-         if (tnat1) then
-                 tnat=1.0
-                 alpha_ideal=1.0
-                 write(*,*) 'attention dans dynetat0: les alpha_ideal sont a 1'
-         else
+         IF(ltnat1) THEN
+            tnat = 1.0
+            alpha_ideal = 1.0
+            CALL msg(' !!!  Beware: alpha_ideal put to 1  !!!', modname)
+         ELSE
           IF(getKey('tnat', tnat, isoName(iName)) .OR. getKey('alpha', alpha_ideal, isoName(iName))) &
             CALL abort_gcm(TRIM(modname), 'missing isotopic parameters', 1)
-         endif
-         CALL msg('Tracer <'//TRIM(var)//'> is missing => initialized with a simplified Rayleigh distillation law.', modname)
+         END IF
+         CALL msg('Missing tracer <'//TRIM(var)//'> => initialized with a simplified Rayleigh distillation law.', modname)
          q(ijb_u:ije_u,:,iq) = q(ijb_u:ije_u,:,iqParent)*tnat*(q(ijb_u:ije_u,:,iqParent)/30.e-3)**(alpha_ideal-1.)
          ! Camille 9 mars 2023: point de vigilence: initialisation incohérente
          ! avec celle de xt_ancien dans la physiq.
       ELSE
-         CALL msg('Tracer <'//TRIM(var)//'> is missing => initialized to its parent isotope concentration.', modname)
+         CALL msg('Missing tracer <'//TRIM(var)//'> => initialized to its parent isotope concentration.', modname)
          ! Camille 9 mars 2023: attention!! seuls les tags qui correspondent à
          ! izone=izone_init (définie dans isotrac_mod) sont initialisés comme
@@ -208,5 +213,5 @@
     !--------------------------------------------------------------------------------------------------------------------------
     ELSE                                                                                 !=== MISSING: SET TO 0
-      CALL msg('Tracer <'//TRIM(var)//'> is missing => initialized to zero', modname)
+      CALL msg('missing tracer <'//TRIM(var)//'> => initialized to zero', modname)
       q(ijb_u:ije_u,:,iq)=0.
     !--------------------------------------------------------------------------------------------------------------------------
Index: LMDZ6/branches/cirrus/libf/dyn3dmem/iniacademic_loc.F90
===================================================================
--- LMDZ6/branches/cirrus/libf/dyn3dmem/iniacademic_loc.F90	(revision 5202)
+++ LMDZ6/branches/cirrus/libf/dyn3dmem/iniacademic_loc.F90	(revision 5203)
@@ -5,5 +5,5 @@
 
   USE filtreg_mod, ONLY: inifilr
-  USE infotrac,    ONLY: nqtot, niso, iqIsoPha, tracers, getKey, isoName
+  USE infotrac,    ONLY: nqtot, niso, iqIsoPha, tracers, getKey, isoName, addPhase
   USE control_mod, ONLY: day_step,planet_type
   use exner_hyb_m, only: exner_hyb
@@ -22,5 +22,4 @@
   USE temps_mod, ONLY: annee_ref, day_ini, day_ref
   USE ener_mod, ONLY: etot0,ptot0,ztot0,stot0,ang0
-  USE readTracFiles_mod, ONLY: addPhase
   use netcdf, only : NF90_NOWRITE,NF90_OPEN,NF90_NOERR,NF90_INQ_VARID
   use netcdf, only : NF90_CLOSE, NF90_GET_VAR
@@ -85,5 +84,5 @@
 
   REAL zdtvr, tnat, alpha_ideal
-  LOGICAL,PARAMETER :: tnat1=.true.
+  LOGICAL :: ltnat1
   
   character(len=*),parameter :: modname="iniacademic"
@@ -311,4 +310,5 @@
         ! bulk initialization of tracers
         if (planet_type=="earth") then
+           ltnat1 = .TRUE.; CALL getin('tnateq1', ltnat1)
            ! Earth: first two tracers will be water
            do iq=1,nqtot
@@ -324,12 +324,12 @@
               iqParent = tracers(iq)%iqParent
               IF(tracers(iq)%iso_iZone == 0) THEN
-                 if (tnat1) then
-                         tnat=1.0
-                         alpha_ideal=1.0
-                         write(*,*) 'Attention dans iniacademic: alpha_ideal=1'
-                 else
+                 IF(ltnat1) THEN
+                    tnat = 1.0
+                    alpha_ideal = 1.0
+                    WRITE(lunout, *) 'In '//TRIM(modname)//': !!!  Beware: alpha_ideal put to 1  !!!'
+                 ELSE
                     IF(getKey('tnat', tnat, isoName(iName)) .OR. getKey('alpha', alpha_ideal, isoName(iName))) &
                     CALL abort_gcm(TRIM(modname), 'missing isotopic parameters', 1)
-                 endif
+                 END IF
                  q(ijb_u:ije_u,:,iq) = q(ijb_u:ije_u,:,iqParent)*tnat*(q(ijb_u:ije_u,:,iqParent)/30.e-3)**(alpha_ideal-1.)
               ELSE !IF(tracers(iq)%iso_iZone == 0) THEN
Index: LMDZ6/branches/cirrus/libf/dyn3dmem/qminimum_loc.F
===================================================================
--- LMDZ6/branches/cirrus/libf/dyn3dmem/qminimum_loc.F	(revision 5202)
+++ LMDZ6/branches/cirrus/libf/dyn3dmem/qminimum_loc.F	(revision 5203)
@@ -4,8 +4,7 @@
       SUBROUTINE qminimum_loc( q,nqtot,deltap )
       USE parallel_lmdz
-      USE infotrac, ONLY: niso, ntiso, iqIsoPha, tracers, 
+      USE infotrac, ONLY: niso, ntiso, iqIsoPha, tracers, addPhase,
      &                    isoCheck, min_qParent
       USE strings_mod, ONLY: strIdx
-      USE readTracFiles_mod, ONLY: addPhase
       IMPLICIT none
 c
Index: LMDZ6/branches/cirrus/libf/phylmd/infotrac_phy.F90
===================================================================
--- LMDZ6/branches/cirrus/libf/phylmd/infotrac_phy.F90	(revision 5202)
+++ LMDZ6/branches/cirrus/libf/phylmd/infotrac_phy.F90	(revision 5203)
@@ -7,4 +7,5 @@
         delPhase, niso, getKey, isot_type, processIsotopes,  isotope, maxTableWidth, iqIsoPha, nphas, ixIso, isoPhas, &
         addPhase, iH2O, addKey, isoSelect, testTracersFiles, isoKeys, indexUpdate,  iqWIsoPha, nbIso, ntiso, isoName, isoCheck
+   USE readTracFiles_mod, ONLY: new2oldH2O
    IMPLICIT NONE
 
@@ -16,24 +17,24 @@
    PUBLIC :: nqtot,   nbtr,   nqo,   nqCO2,   nqtottr      !--- Main dimensions
    PUBLIC :: conv_flg, pbl_flg                             !--- Convection & boundary layer activation keys
-#ifdef CPP_StratAer
+   PUBLIC :: new2oldH2O                                    !--- For backwards compatibility in phyetat0
+   PUBLIC :: addPhase, delPhase                            !--- Add/remove the phase from the name of a tracer
+#if defined CPP_StratAer || defined REPROBUS
    PUBLIC :: nbtr_bin, nbtr_sulgas                         !--- Number of aerosols bins and sulfur gases for StratAer model
    PUBLIC :: id_OCS_strat, id_SO2_strat, id_H2SO4_strat, id_BIN01_strat, id_TEST_strat
 #endif
 
-   !=== FOR WATER
-   PUBLIC :: ivap, iliq, isol
    !=== 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 :: iH2O                                          !--- Value of "ixIso" for "H2O" isotopes class
+   PUBLIC :: ivap, iliq, isol
    !=== 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 :: iqWIsoPha                                      !--- Same as iqIsoPha but with normal water phases
-
+   PUBLIC :: isotope                                       !--- Selected isotopes database (argument of getKey)
+   PUBLIC :: isoKeys, isoName, isoZone, isoPhas            !--- Isotopes keys & names, tagging zones names, phases
+   PUBLIC ::    niso,   ntiso,   nzone,   nphas            !--- Number of   "   "
+   PUBLIC :: itZonIso                                      !--- Index "it" in "isoName(1:niso)" = f(tagging idx, isotope idx)
+   PUBLIC :: iqIsoPha                                      !--- Index "iq" in "qx"              = f(isotope idx,   phase idx)
+   PUBLIC :: iqWIsoPha                                     !--- Same as iqIsoPha but with normal water phases
    PUBLIC :: isoCheck                                      !--- Run isotopes checking routines
    !=== FOR BOTH TRACERS AND ISOTOPES
@@ -43,5 +44,5 @@
 !  |--------------------+-----------------------+-----------------+---------------+----------------------------|
 !  | water in different |    water tagging      |  water isotopes | other tracers | additional tracers moments |
-!  | phases: H2O_[gls]  |      isotopes         |                 |               |  for higher order schemes  |
+!  | phases: H2O_[glsrb]|      isotopes         |                 |               |  for higher order schemes  |
 !  |--------------------+-----------------------+-----------------+---------------+----------------------------|
 !  |                    |                       |                 |               |                            |
@@ -68,9 +69,11 @@
 !  |-------------+------------------------------------------------------+-------------+------------------------+
 !  | name        | Name (short)                                         | tname       |                        |
+!  | keys        | key/val pairs accessible with "getKey" routine       | /           |                        |
 !  | 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]              |
+!  | phase       | Phases list ("g"as / "l"iquid / "s"olid              |             | [g|l|s|r|b]            |
+!  |             |              "r"(cloud) / "b"lowing)                 | /           |                        |
 !  | component   | Name(s) of the merged/cumulated section(s)           | /           | coma-separated names   |
 !  | iGeneration | Generation (>=1)                                     | /           |                        |
@@ -79,7 +82,6 @@
 !  | 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       | /           |                        |
-!  | isAdvected  | advected tracers flag (.TRUE. if iadv >= 0)          | /           | nqtrue  .TRUE. values  |
-!  | isInPhysics | tracers not extracted from the main table in physics | /           | nqtottr .TRUE. values  |
+!  | 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                 |
@@ -97,14 +99,18 @@
 !  | 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  |
+!  | phase  | nphas  | Phases                     list + number         |                    | [g|l|s|r|b] 1:5 |
 !  | iqIsoPha        | Index in "qx"           = f(name(1:ntiso)),phas) | iqiso              | 1:nqtot         |
-!  | iqWIsoPha       | Index in "qx"           = f(name(1:ntiso)),phas) | iqiso              | 1:nqtot         |
+!  | iqWIsoPha       | Index in "qx"       = f(name(1:ntiso+nqo)),phas) |   /                | 1:nqtot         |
 !  | itZonIso        | Index in "trac(1:ntiso)"= f(zone, name(1:niso))  | index_trac         | 1:ntiso         |
 !  +-----------------+--------------------------------------------------+--------------------+-----------------+
 
+   !=== INDICES FOR WATER
+   INTEGER, SAVE :: ivap, iliq, isol
+!$OMP THREADPRIVATE(ivap, iliq, isol)
+
    !=== 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
+   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)
@@ -112,14 +118,10 @@
 !$OMP THREADPRIVATE(nqtot, nbtr, nqo, nqtottr, nqCO2, type_trac)
 
-   !=== INDICES OF WATER
-   INTEGER,               SAVE :: ivap,iliq,isol ! Indices for vap, liq and ice
-!$OMP THREADPRIVATE(ivap,iliq,isol)
-
    !=== VARIABLES FOR INCA
-   INTEGER,               SAVE, ALLOCATABLE :: conv_flg(:), &   !--- Convection     activation ; needed for INCA        (nbtr)
-                                                pbl_flg(:)      !--- Boundary layer activation ; needed for INCA        (nbtr)
+   INTEGER, DIMENSION(:), SAVE, ALLOCATABLE :: &
+                    conv_flg, pbl_flg                           !--- Convection / boundary layer activation (nbtr)
 !$OMP THREADPRIVATE(conv_flg, pbl_flg)
 
-#ifdef CPP_StratAer
+#if defined CPP_StratAer || defined REPROBUS
   !=== SPECIFIC TO STRATOSPHERIC AEROSOLS (CK/OB)
   INTEGER, SAVE ::  nbtr_bin, nbtr_sulgas         !--- number of aerosols bins and sulfur gases for StratAer model
@@ -173,15 +175,14 @@
    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, texp, ttp                          !--- String for messages and expanded tracers type
+   CHARACTER(LEN=maxlen) :: msg1, texp, ttp                          !--- Strings for messages and expanded tracers type
    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 :: iq, jq, nt, im, nm, k                                 !--- Indexes and temporary variables
+   INTEGER :: iq, jq, nt, im, nm                                     !--- Indexes and temporary variables
    LOGICAL :: lerr, lInit
    TYPE(trac_type), ALLOCATABLE, TARGET :: ttr(:)
    TYPE(trac_type), POINTER             :: t1, t(:)
-   CHARACTER(LEN=maxlen),   ALLOCATABLE :: types_trac(:)  !--- Keyword for tracers type(s), parsed version
-   
+   CHARACTER(LEN=maxlen),   ALLOCATABLE :: types_trac(:)             !--- Keywords for tracers type(s), parsed version
    CHARACTER(LEN=*), PARAMETER :: modname="init_infotrac_phy"
 !------------------------------------------------------------------------------------------------------------------------------
@@ -202,5 +203,4 @@
    ENDIF
 
-   
    CALL msg('type_trac = "'//TRIM(type_trac)//'"', modname)
    lInit = .NOT.ALLOCATED(tracers)
@@ -243,13 +243,10 @@
 
 !==============================================================================================================================
-! 1) Get the numbers of: true (first order only) tracers "nqtrue", water tracers "nqo" (vapor/liquid/solid)
-!==============================================================================================================================
-   texp = type_trac                                                  !=== EXPANDED VERSION OF "type_trac", WITH "|" SEPARATOR
+! 0) DETERMINE THE TYPE OF THE INPUT TRACERS DESCRIPTION FILE AND READ IT
+!==============================================================================================================================
+   texp = type_trac                                                            !=== EXPANDED (WITH "|" SEPARATOR) "type_trac"
    IF(texp == 'inco') texp = 'co2i|inca'
    IF(texp /= 'lmdz') texp = 'lmdz|'//TRIM(texp)
-
-   !=== DETERMINE THE TYPE OF THE INPUT TRACERS DESCRIPTION FILE
    IF(testTracersFiles(modname, texp, fType, lInit)) CALL abort_physic(modname, 'problem with tracers file(s)',1)
-
    ttp = type_trac; IF(fType /= 1) ttp = texp
 
@@ -262,8 +259,11 @@
 !##############################################################################################################################
 
-   !---------------------------------------------------------------------------------------------------------------------------
-   IF(fType == 0) CALL abort_physic(modname, 'Missing tracers file: "traceur.def", "tracer.def" or "tracer_<keyword>.def file.',1)
-   !---------------------------------------------------------------------------------------------------------------------------
-   IF(fType == 1 .AND. ANY(['inca','inco']==type_trac) .AND. lInit) THEN  !=== FOUND OLD STYLE INCA "traceur.def"
+!==============================================================================================================================
+! 1) Get various numbers: "nqtrue" (first order only tracers), "nqo" (water phases), 'nbtr' (tracers passed to physics), etc.
+!==============================================================================================================================
+   !---------------------------------------------------------------------------------------------------------------------------
+   IF(fType == 0) CALL abort_physic(modname, 'Missing "traceur.def", "tracer.def" or "tracer_<keyword>.def tracers file.',1)
+   !---------------------------------------------------------------------------------------------------------------------------
+   IF(fType == 1 .AND. ANY(['inca', 'inco'] == type_trac) .AND. lInit) THEN  !=== FOUND OLD STYLE INCA "traceur.def"
    !---------------------------------------------------------------------------------------------------------------------------
 #ifdef INCA
@@ -281,5 +281,5 @@
       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      :nqtrue)%name      = [('CO2     ', iq=1, nqCO2), solsym_inca]
       ttr(1+nqo+nqCO2:nqtrue)%parent    = tran0
       ttr(1+nqo+nqCO2:nqtrue)%phase     = 'g'
@@ -302,21 +302,19 @@
    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
+   nqtrue = SIZE(tracers)                                                                               !--- "true" tracers
+   nqo    =      COUNT(tracers(:)%component == 'lmdz' .AND. delPhase(tracers(:)%name)     == 'H2O')     !--- Water phases
+   nbtr = nqtrue-COUNT(tracers(:)%component == 'lmdz' .AND. delPhase(tracers(:)%gen0Name) == 'H2O')     !--- Passed to phytrac
+   nqCO2  =      COUNT( [type_trac == 'inco', type_trac == 'co2i'] )
 #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
+   nqINCA =      COUNT(tracers(:)%component == 'inca')
+#endif
+   IF(getKey('hadv', hadv, ky=tracers(:)%keys)) CALL abort_physic(modname, 'missing key "hadv"', 1)
+   IF(getKey('vadv', vadv, ky=tracers(:)%keys)) CALL abort_physic(modname, 'missing key "vadv"', 1)
+   !---------------------------------------------------------------------------------------------------------------------------
+   END IF
+   !---------------------------------------------------------------------------------------------------------------------------
+
 #ifdef REPROBUS
-   CALL Init_chem_rep_trac(nbtr, nqo, tracers(:)%name)
+   CALL Init_chem_rep_trac(nbtr, nqo, tracers(:)%name)                         !--- Transfert the number of tracers to Reprobus
 #endif
 
@@ -359,9 +357,8 @@
       IF(iad == -1) CALL abort_physic(modname, msg1, 1)
 
-      !--- SET FIELDS %longName, %isAdvected, %isInPhysics
+      !--- SET FIELDS longName, isAdvected, isInPhysics
       t1%longName   = t1%name; IF(iad > 0) t1%longName=TRIM(t1%name)//descrq(iad)
       t1%isAdvected = iad >= 0
-      t1%isInPhysics= delPhase(t1%gen0Name) /= 'H2O' &
-                          .OR. t1%component /= 'lmdz' !=== OTHER EXCEPTIONS TO BE ADDED: CO2i, SURSATURATED WATER CLOUD...
+      t1%isInPhysics= delPhase(t1%gen0Name) /= 'H2O' .OR. t1%component /= 'lmdz' !=== MORE EXCEPTIONS ? CO2i, SURSAT CLOUD H2O
       ttr(iq)       = t1
 
@@ -372,8 +369,10 @@
       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)%isAdvected  = [(.FALSE., im=1, nm) ]
+      ttr(jq+1:jq+nm)%name        = [ (TRIM(t1%name)    //'-'//TRIM(suff(im)), im=1, nm) ]
+      ttr(jq+1:jq+nm)%gen0Name    = [ (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)%isAdvected  = [ (.FALSE., im=1, nm) ]
+      ttr(jq+1:jq+nm)%isInPhysics = [ (.FALSE., im=1, nm) ]
       jq = jq + nm
    END DO
@@ -381,6 +380,6 @@
    CALL MOVE_ALLOC(FROM=ttr, TO=tracers)
 
-   !--- SET FIELDS %iqParent, %nqChildren, %iGeneration, %iqDescen, %nqDescen
-   IF(indexUpdate(tracers)) CALL abort_gcm(modname, 'problem when processing isotopes parameters', 1)
+   !--- SET FIELDS iqParent, iqDescen, nqDescen, nqChildren, iGeneration
+   IF(indexUpdate(tracers)) CALL abort_physic(modname, 'problem with tracers indices update', 1)
 
 !##############################################################################################################################
@@ -404,6 +403,6 @@
 !##############################################################################################################################
    !--- Convection / boundary layer activation for all tracers
-   IF (.NOT.ALLOCATED(conv_flg)) ALLOCATE(conv_flg(nbtr)); conv_flg(1:nbtr) = 1
-   IF (.NOT.ALLOCATED(pbl_flg)) ALLOCATE( pbl_flg(nbtr));  pbl_flg(1:nbtr) = 1
+   IF(.NOT.ALLOCATED(conv_flg)) ALLOCATE(conv_flg(nbtr)); conv_flg(1:nbtr) = 1
+   IF(.NOT.ALLOCATED( pbl_flg)) ALLOCATE( pbl_flg(nbtr));  pbl_flg(1:nbtr) = 1
 
    !--- Note: nqtottr can differ from nbtr when nmom/=0
@@ -424,8 +423,9 @@
 #endif
    t => tracers
-   CALL msg('Information stored in infotrac_phy :', modname)
-   IF(dispTable('issssssssiiiiiiii', ['iq  ', 'name', 'lNam', 'g0Nm', 'prnt', 'type', 'phas', 'comp',      &
-                              'isPh', 'iGen', 'iqPr', 'nqDe', 'nqCh', 'iGrp', 'iNam', 'iZon', 'iPha'],     &
-      cat(t%name, t%longName, t%gen0Name, t%parent, t%type, t%phase, t%component, bool2str(t%isInPhysics)),&
+   CALL msg('Information stored in '//TRIM(modname)//': ', modname)
+   IF(dispTable('isssssssssiiiiiiii', ['iq  ', 'name', 'lNam', 'g0Nm', 'prnt', 'type', 'phas', 'comp',     &
+                       'isPh', 'isAd', 'iGen', 'iqPr', 'nqDe', 'nqCh', 'iGrp', 'iNam', 'iZon', 'iPha'],    &
+      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%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))   &
Index: LMDZ6/branches/cirrus/libf/phylmd/lmdz_lscp.F90
===================================================================
--- LMDZ6/branches/cirrus/libf/phylmd/lmdz_lscp.F90	(revision 5202)
+++ LMDZ6/branches/cirrus/libf/phylmd/lmdz_lscp.F90	(revision 5203)
@@ -1526,8 +1526,4 @@
                 znebprecipcld(i)=0.0
             ENDIF
-        !IF ( ((1-zfice(i))*zoliq(i) .GT. 0.) .AND. (zt(i) .LE. 233.15) ) THEN
-        !print*,'WARNING LEA OLIQ A <-40°C '
-        !print*,'zt,Tbef,oliq,oice,cldfraliq,icefrac,rneb',zt(i),Tbef(i),(1-zfice(i))*zoliq(i),zfice(i)*zoliq(i),cldfraliq(i,k),zfice(i),rneb(i,k)
-        !ENDIF 
         ENDDO
 
Index: LMDZ6/branches/cirrus/libf/phylmd/lmdz_wake.F90
===================================================================
--- LMDZ6/branches/cirrus/libf/phylmd/lmdz_wake.F90	(revision 5202)
+++ LMDZ6/branches/cirrus/libf/phylmd/lmdz_wake.F90	(revision 5203)
@@ -1,5 +1,8 @@
 MODULE lmdz_wake
 
-! $Id$
+  USE lmdz_wake_ini, ONLY: CPPKEY_IOPHYS_WK
+
+  IMPLICIT NONE; PRIVATE
+  PUBLIC wake
 
 CONTAINS
@@ -353,7 +356,5 @@
 !
 IF (first_call) THEN
-!!#define IOPHYS_WK
-#undef IOPHYS_WK
-#ifdef IOPHYS_WK
+IF (CPPKEY_IOPHYS_WK) THEN
   IF (phys_sub) THEN
     call iophys_ini(dtimesub)
@@ -361,5 +362,5 @@
     call iophys_ini(dtime)
   ENDIF
-#endif
+END IF
   first_call = .false.
 ENDIF   !(first_call)
@@ -771,7 +772,7 @@
 
   END DO
-#ifdef IOPHYS_WK
+IF (CPPKEY_IOPHYS_WK) THEN
   IF (.not.phys_sub) CALL iophys_ecrit('wape_a',1,'wape_a','J/kg',wape)
-#endif
+END IF
 
   ! 2.2 Prognostic variable update
@@ -1007,5 +1008,5 @@
 !!--------------------------------------------------------
     ELSEIF (iflag_wk_pop_dyn == 3) THEN
-#ifdef IOPHYS_WK
+IF (CPPKEY_IOPHYS_WK) THEN
     IF (phys_sub) THEN
      CALL iophys_ecrit('ptop',1,'ptop','Pa',ptop)
@@ -1018,5 +1019,5 @@
      CALL iophys_ecrit('irad_wk',1,'irad_wk','m',irad_wk)
     ENDIF
-#endif
+END IF
   !
      CALL wake_popdyn_3 ( klon, klev, phys_sub, wk_adv, dtimesub, wgen, &
@@ -1072,5 +1073,5 @@
 !!--------------------------------------------------------
 
-#ifdef IOPHYS_WK
+IF (CPPKEY_IOPHYS_WK) THEN
     IF (phys_sub) THEN
      CALL iophys_ecrit('wdensa',1,'wdensa','m',wdens)
@@ -1079,5 +1080,5 @@
      CALL iophys_ecrit('asigmawa',1,'asigmawa','m',asigmaw)
     ENDIF
-#endif
+END IF
     ! calcul de la difference de vitesse verticale poche - zone non perturbee
     ! IM 060208 differences par rapport au code initial; init. a 0 dp_deltomg
@@ -1669,5 +1670,5 @@
           END DO
 
-#ifdef IOPHYS_WK
+IF (CPPKEY_IOPHYS_WK) THEN
     IF (phys_sub) THEN
      CALL iophys_ecrit('wdensb',1,'wdensb','m',wdens)
@@ -1702,5 +1703,5 @@
      CALL iophys_ecrit('d_asig_bnd2',1,'d_asig_bnd2','m',d_asig_bnd2)
     ENDIF
-#endif
+END IF
         ENDIF ! (iflag_wk_pop_dyn == 3)
       ENDIF ! (iflag_wk_pop_dyn >= 2)
@@ -1866,7 +1867,7 @@
   !
 
-#ifdef IOPHYS_WK
+IF (CPPKEY_IOPHYS_WK) THEN
     IF (.not.phys_sub) CALL iophys_ecrit('wape_b',1,'wape_b','J/kg',wape)
-#endif
+END IF
   IF (prt_level>=10) THEN
     PRINT *, 'wake-5, sigmaw(igout), cstar(igout), wape(igout), ptop(igout) ', &
@@ -2009,7 +2010,7 @@
     END IF
   END DO
-#ifdef IOPHYS_WK
+IF (CPPKEY_IOPHYS_WK) THEN
   IF (.not.phys_sub) CALL iophys_ecrit('wape2_a',1,'wape2_a','J/kg',wape2)
-#endif
+END IF
 
 
@@ -2041,7 +2042,7 @@
     END DO
   END IF
-#ifdef IOPHYS_WK
+IF (CPPKEY_IOPHYS_WK) THEN
   IF (.not.phys_sub) CALL iophys_ecrit('wape2_b',1,'wape2_b','J/kg',wape2)
-#endif
+END IF
 
 
@@ -2087,7 +2088,7 @@
         gwake(i) = .TRUE.
       END IF
-#ifdef IOPHYS_WK
+IF (CPPKEY_IOPHYS_WK) THEN
   IF (.not.phys_sub) CALL iophys_ecrit('cstar2',1,'cstar2','J/kg',cstar2)
-#endif
+END IF
     END IF  ! (ok_qx_qw(i))
   END DO
@@ -2123,5 +2124,5 @@
   END DO
     IF (iflag_wk_pop_dyn >= 3) THEN
-#ifdef IOPHYS_WK
+IF (CPPKEY_IOPHYS_WK) THEN
       IF (.not.phys_sub) THEN
        CALL iophys_ecrit('fip',1,'fip','J/kg',fip)
@@ -2163,5 +2164,5 @@
        CALL iophys_ecrit('d_asig_bnd2',1,'d_asig_bnd2','m',d_asig_bnd2)
       ENDIF  ! (.not.phys_sub)
-#endif
+END IF
     ENDIF  ! (iflag_wk_pop_dyn >= 3)
   ! Limitation de sigmaw
@@ -2268,7 +2269,7 @@
                       wape(igout),wape2(igout),ktopw(igout),OK_qx_qw(igout)
   ENDIF
-#ifdef IOPHYS_WK
+IF (CPPKEY_IOPHYS_WK) THEN
   IF (.not.phys_sub) CALL iophys_ecrit('wape_c',1,'wape_c','J/kg',wape)
-#endif
+END IF
 
 
@@ -3255,7 +3256,7 @@
           d_sig_spread(i) =  d_sig_spread(i)*dtimesub 
           d_sigmaw(i) =  d_sig_gen(i) + d_sig_death(i) + d_sig_col(i) + d_sig_spread(i)
-#ifdef IOPHYS_WK
+IF (CPPKEY_IOPHYS_WK) THEN
           IF (phys_sub) call iophys_ecrit('d_sigmaw0',1,'d_sigmaw0','',d_sigmaw)
-#endif
+END IF
 
          
@@ -3266,5 +3267,5 @@
           d_sigmaw(i) = d_sigmaw_targ
 !!          d_sigmaw(i) = max(d_sigmaw(i), sigmad-sigmaw(i))
-#ifdef IOPHYS_WK
+IF (CPPKEY_IOPHYS_WK) THEN
           IF (phys_sub) THEN
              call iophys_ecrit('tauwk_inv',1,'tau_wk_inv_min','',tau_wk_inv_min)
@@ -3276,5 +3277,5 @@
              call iophys_ecrit('d_sig_bnd',1,'d_sig_bnd','',d_sig_bnd)
           ENDIF
-#endif
+END IF
           d_asig_death(i) = - asigmaw(i)/tau_prime(i)
           d_asig_aicol(i) = (agfl(i)*iwdens(i) + igfl(i)*awdens(i))*cstar(i)*is_wk(i)
@@ -3287,7 +3288,7 @@
           d_asig_spread(i) =  d_asig_spread(i)*dtimesub 
           d_asigmaw(i) =  d_sig_gen(i) + d_asig_death(i) + d_asig_aicol(i) + d_asig_iicol(i) + d_asig_spread(i)
-#ifdef IOPHYS_WK
+IF (CPPKEY_IOPHYS_WK) THEN
           IF (phys_sub) call iophys_ecrit('d_asigmaw0',1,'d_asigmaw0','',d_asigmaw)
-#endif
+END IF
 
           d_sigmaw_targ = min(max(d_asigmaw(i),-asigmaw(i)), sigmaw(i)-asigmaw(i))
@@ -3295,5 +3296,5 @@
           d_asig_bnd(i) = d_sigmaw_targ - d_asigmaw(i)
           d_asigmaw(i) = d_sigmaw_targ
-#ifdef IOPHYS_WK
+IF (CPPKEY_IOPHYS_WK) THEN
           IF (phys_sub) THEN
              call iophys_ecrit('d_asigmaw',1,'d_asigmaw','',d_asigmaw)
@@ -3304,5 +3305,5 @@
              call iophys_ecrit('d_asig_bnd',1,'d_asig_bnd','',d_asig_bnd)
           ENDIF
-#endif
+END IF
           d_dens_gen(i) = wgen(i)
           d_dens_death(i) = - iwdens(i)*tau_wk_inv_min 
@@ -3318,5 +3319,5 @@
           d_dens_bnd(i) = d_wdens_targ - d_wdens(i)
           d_wdens(i) = d_wdens_targ
-#ifdef IOPHYS_WK
+IF (CPPKEY_IOPHYS_WK) THEN
     IF (phys_sub) THEN 
         call iophys_ecrit('d_wdens',1,'d_wdens','',d_wdens)
@@ -3325,5 +3326,5 @@
         call iophys_ecrit('d_dens_col',1,'d_dens_col','',d_dens_col)
     ENDIF
-#endif
+END IF
 
           d_adens_death(i) = -awdens(i)/tau_prime(i)
@@ -3335,5 +3336,5 @@
           d_adens_acol(i)  =   d_adens_acol(i)*dtimesub
           d_awdens(i) =   d_dens_gen(i) + d_adens_death(i) + d_adens_icol(i) + d_adens_acol(i)     
-#ifdef IOPHYS_WK
+IF (CPPKEY_IOPHYS_WK) THEN
     IF (phys_sub) THEN 
         call iophys_ecrit('d_awdens',1,'d_awdens','',d_awdens)
@@ -3342,5 +3343,5 @@
         call iophys_ecrit('d_adens_acol',1,'d_adens_acol','',d_adens_acol)
     ENDIF
-#endif
+END IF
           d_wdens_targ = min(max(d_awdens(i),-awdens(i)), wdens(i)-awdens(i))
 !!          d_dens_bnd(i) = d_dens_bnd(i) + d_wdens_targ - d_wdens(i)
Index: LMDZ6/branches/cirrus/libf/phylmd/lmdz_wake_ini.F90
===================================================================
--- LMDZ6/branches/cirrus/libf/phylmd/lmdz_wake_ini.F90	(revision 5202)
+++ LMDZ6/branches/cirrus/libf/phylmd/lmdz_wake_ini.F90	(revision 5203)
@@ -81,4 +81,12 @@
   !$OMP THREADPRIVATE(wk_int_delta_t_min)
 
+! CPP key used only in this module for debugging purposes. jyg 09/24
+!!#define IOPHYS_WK
+#ifdef IOPHYS_WK
+  LOGICAL, PARAMETER :: CPPKEY_IOPHYS_WK = .TRUE.
+#else
+  LOGICAL, PARAMETER :: CPPKEY_IOPHYS_WK = .FALSE.
+#endif
+
 
 
Index: LMDZ6/branches/cirrus/libf/phylmd/phyetat0_mod.F90
===================================================================
--- LMDZ6/branches/cirrus/libf/phylmd/phyetat0_mod.F90	(revision 5202)
+++ LMDZ6/branches/cirrus/libf/phylmd/phyetat0_mod.F90	(revision 5203)
@@ -33,6 +33,6 @@
   USE geometry_mod,     ONLY: longitude_deg, latitude_deg
   USE iostart,          ONLY: close_startphy, get_field, get_var, open_startphy
-  USE infotrac_phy,     ONLY: nqtot, nbtr, type_trac, tracers
-  USE readTracFiles_mod,ONLY: maxlen, new2oldH2O
+  USE infotrac_phy,     ONLY: nqtot, nbtr, type_trac, tracers, new2oldH2O
+  USE strings_mod,      ONLY: maxlen
   USE traclmdz_mod,     ONLY: traclmdz_from_restart
   USE carbon_cycle_mod, ONLY: carbon_cycle_init, carbon_cycle_cpl, carbon_cycle_tr, carbon_cycle_rad, co2_send, RCO2_glo
Index: LMDZ6/branches/cirrus/libf/phylmd/phys_local_var_mod.F90
===================================================================
--- LMDZ6/branches/cirrus/libf/phylmd/phys_local_var_mod.F90	(revision 5202)
+++ LMDZ6/branches/cirrus/libf/phylmd/phys_local_var_mod.F90	(revision 5203)
@@ -400,6 +400,6 @@
       REAL,ALLOCATABLE,SAVE,DIMENSION(:) :: tpot, tpote, ue, uq, uwat, ve, vq, vwat, zxffonte
 !$OMP THREADPRIVATE(tpot, tpote, ue, uq, uwat, ve, vq, vwat, zxffonte)
-      REAL,ALLOCATABLE,SAVE,DIMENSION(:) :: zxustartlic, zxrhoslic, zxqsaltlic
-!$OMP THREADPRIVATE(zxustartlic, zxrhoslic, zxqsaltlic)
+      REAL,ALLOCATABLE,SAVE,DIMENSION(:) :: zxustartlic, zxrhoslic, zxqsaltlic, tempsmoothlic
+!$OMP THREADPRIVATE(zxustartlic, zxrhoslic, zxqsaltlic, tempsmoothlic)
       REAL,ALLOCATABLE,SAVE,DIMENSION(:) :: zxfqcalving
 !$OMP THREADPRIVATE(zxfqcalving)
@@ -1051,6 +1051,6 @@
       ALLOCATE(zxrunofflic(klon), runoff_diag(klon))
       runoff_diag(:)=0.
-      ALLOCATE(zxustartlic(klon), zxrhoslic(klon), zxqsaltlic(klon))
-      zxustartlic(:)=0. ; zxrhoslic(:)=0. ; zxqsaltlic(:)=0.
+      ALLOCATE(zxustartlic(klon), zxrhoslic(klon), zxqsaltlic(klon), tempsmoothlic(klon))
+      zxustartlic(:)=0. ; zxrhoslic(:)=0. ; zxqsaltlic(:)=0. ; tempsmoothlic(:)=0.
       ALLOCATE(rain_lsc(klon))
       ALLOCATE(rain_num(klon))
@@ -1480,5 +1480,5 @@
 ! SN runoff_diag
       DEALLOCATE(zxrunofflic, runoff_diag)
-      DEALLOCATE(zxustartlic, zxrhoslic, zxqsaltlic)
+      DEALLOCATE(zxustartlic, zxrhoslic, zxqsaltlic, tempsmoothlic)
       DEALLOCATE(zxtsol, snow_lsc, zxfqfonte, zxqsurf)
       DEALLOCATE(rain_lsc)
Index: LMDZ6/branches/cirrus/libf/phylmd/physiq_mod.F90
===================================================================
--- LMDZ6/branches/cirrus/libf/phylmd/physiq_mod.F90	(revision 5202)
+++ LMDZ6/branches/cirrus/libf/phylmd/physiq_mod.F90	(revision 5203)
@@ -39,6 +39,5 @@
     USE ioipsl_getin_p_mod, ONLY : getin_p
     USE indice_sol_mod
-    USE infotrac_phy, ONLY: nqtot, nbtr, nqo, tracers, type_trac
-    USE readTracFiles_mod, ONLY: addPhase
+    USE infotrac_phy, ONLY: nqtot, nbtr, nqo, tracers, type_trac, addPhase
     USE strings_mod,  ONLY: strIdx
     USE iophy
@@ -1270,4 +1269,7 @@
     !albedo SB <<<
 
+    !--Lea Raillard qs_ini 
+    REAL, dimension(klon,klev) :: qs_ini
+
     !--OB variables for mass fixer (hard coded for now)
     REAL qql1(klon),qql2(klon),corrqql
@@ -2461,4 +2463,6 @@
        ENDDO
     ENDDO
+    ! Lea Raillard qs_ini for cloud phase param.
+    qs_ini(:,:)=qs_seri(:,:)
     !
     !--OB water mass fixer 
@@ -3858,5 +3862,5 @@
 
     CALL lscp(klon,klev,phys_tstep,missing_val,paprs,pplay, &
-         t_seri, q_seri,qs_ancien,ptconv,ratqs, &
+         t_seri, q_seri,qs_ini,ptconv,ratqs, &
          d_t_lsc, d_q_lsc, d_ql_lsc, d_qi_lsc, rneb, rneblsvol, & 
          pfraclr, pfracld, cldfraliq, sigma2_icefracturb, mean_icefracturb,  &
Index: LMDZ6/branches/cirrus/libf/phylmd/surf_landice_mod.F90
===================================================================
--- LMDZ6/branches/cirrus/libf/phylmd/surf_landice_mod.F90	(revision 5202)
+++ LMDZ6/branches/cirrus/libf/phylmd/surf_landice_mod.F90	(revision 5203)
@@ -36,5 +36,5 @@
     USE cpl_mod,          ONLY : cpl_send_landice_fields
     USE calcul_fluxs_mod
-    USE phys_local_var_mod, ONLY : zxrhoslic, zxustartlic, zxqsaltlic
+    USE phys_local_var_mod, ONLY : zxrhoslic, zxustartlic, zxqsaltlic, tempsmoothlic
     USE phys_output_var_mod, ONLY : snow_o,zfra_o
 #ifdef ISO   
@@ -185,4 +185,5 @@
 
     REAL,DIMENSION(klon) :: alb1,alb2
+    REAL                 :: time_tempsmooth,coef_tempsmooth
     REAL,DIMENSION(klon) :: precip_totsnow, evap_totsnow
     REAL, DIMENSION (klon,6) :: alb6
@@ -230,4 +231,8 @@
            PRINT*, 'alb_nir_sno_lic',alb_nir_sno_lic
   
+  DO j=1,knon
+       i = knindex(j)
+       tempsmoothlic(i) = temp_air(j)
+  ENDDO
   firstcall=.false.
   ENDIF
@@ -445,5 +450,5 @@
     z0h(1:knon) = z0h_landice
 else
-    ! parameterization of z0=f(T) following measurements in Adelie Land by Amory et al 2018
+    ! parameterization of z0=f(T) following measurements in Adelie Land by Amory et al 2017
     coefa = 0.1658 !0.1862 !Ant
     coefb = -50.3869 !-55.7718 !Ant
@@ -456,12 +461,17 @@
     coefc = log(z03/z02)/(ta3-ta2)
     coefd = log(z03)-coefc*ta3
+    time_tempsmooth=2.*86400.
+    coef_tempsmooth=min(1.,dtime/time_tempsmooth)
+    !coef_tempsmooth=0.
     do j=1,knon 
-      if (temp_air(j) .lt. ta1) then
+      i=knindex(j)
+      tempsmoothlic(i)=temp_air(j)*coef_tempsmooth+tempsmoothlic(i)*(1.-coef_tempsmooth)
+      if (tempsmoothlic(i) .lt. ta1) then
         z0m(j) = z01
-      else if (temp_air(j).ge.ta1 .and. temp_air(j).lt.ta2) then
-        z0m(j) = exp(coefa*temp_air(j) + coefb)
-      else if (temp_air(j).ge.ta2 .and. temp_air(j).lt.ta3) then
+      else if (tempsmoothlic(i).ge.ta1 .and. tempsmoothlic(i).lt.ta2) then
+        z0m(j) = exp(coefa*tempsmoothlic(i) + coefb)
+      else if (tempsmoothlic(i).ge.ta2 .and. tempsmoothlic(i).lt.ta3) then
         ! if st > 0, melting induce smooth surface
-        z0m(j) = exp(coefc*temp_air(j) + coefd)
+        z0m(j) = exp(coefc*tempsmoothlic(i) + coefd)
       else
         z0m(j) = z03
Index: LMDZ6/branches/cirrus/libf/phylmdiso/isotopes_mod.F90
===================================================================
--- LMDZ6/branches/cirrus/libf/phylmdiso/isotopes_mod.F90	(revision 5202)
+++ LMDZ6/branches/cirrus/libf/phylmdiso/isotopes_mod.F90	(revision 5203)
@@ -161,5 +161,6 @@
    !=== Local variables:
    INTEGER :: ixt
-
+   LOGICAL :: ltnat1
+   CHARACTER(LEN=maxlen) :: modname, sxt
  
    !--- For H2[17]O
@@ -170,10 +171,7 @@
    LOGICAL, PARAMETER ::   ok_nocinsat = .FALSE. ! if T: no sursaturation effect for ice
    LOGICAL, PARAMETER :: Rdefault_smow = .FALSE. ! if T: Rdefault=smow; if F: nul
-   LOGICAL, PARAMETER :: tnat1 = .TRUE. ! If T: all tnats are 1.
 
    !--- For [3]H
    INTEGER :: iessai
-
-   CHARACTER(LEN=maxlen) :: modname, sxt
 
    modname = 'iso_init'
@@ -265,4 +263,5 @@
    IF(ANY(isoName == 'HTO')) &
    CALL get_in('ok_prod_nucl_tritium', ok_prod_nucl_tritium, .FALSE., .FALSE.)
+   CALL get_in('tnateq1', ltnat1, .TRUE.)
 
    ! Ocean composition
@@ -295,9 +294,5 @@
        tkcin1(ixt) = 0.0005016
        tkcin2(ixt) = 0.0014432
-       if (tnat1) then 
-               tnat(ixt)=1 
-       else 
-               tnat(ixt)=0.
-       endif
+       tnat(ixt) = 0.0; IF(ltnat1) tnat(ixt)=1 
        toce(ixt)=4.0E-19 ! rapport T/H = 0.2 TU Dreisigacker and Roether 1978
        tcorr(ixt)=1.
@@ -318,9 +313,6 @@
        tkcin1(ixt) = tkcin1_O18*fac_kcin
        tkcin2(ixt) = tkcin2_O18*fac_kcin
-       if (tnat1) then 
-               tnat(ixt)=1 
-       else 
-               tnat(ixt)=0.004/100. ! O17 représente 0.004% de l'oxygène
-       endif
+       tnat(ixt)=0.004/100. ! O17 représente 0.004% de l'oxygène
+       IF(ltnat1) tnat(ixt)=1 
        toce(ixt)=tnat(ixt)*(1.0+deltaO18_oce/1000.0)**pente_MWL
        tcorr(ixt)=1.0+fac_enrichoce18*pente_MWL ! donné par Amaelle            
@@ -338,9 +330,5 @@
        tkcin1(ixt) = tkcin1_O18
        tkcin2(ixt) = tkcin2_O18
-       if (tnat1) then
-               tnat(ixt)=1 
-       else 
-               tnat(ixt)=2005.2E-6
-       endif
+       tnat(ixt)=2005.2E-6; IF(ltnat1) tnat(ixt)=1 
        toce(ixt)=tnat(ixt)*(1.0+deltaO18_oce/1000.0)
        tcorr(ixt)=1.0+fac_enrichoce18
@@ -362,9 +350,5 @@
        tkcin1(ixt) = tkcin1_O18*fac_kcin
        tkcin2(ixt) = tkcin2_O18*fac_kcin
-       if (tnat1) then 
-               tnat(ixt)=1
-       else
-               tnat(ixt)=155.76E-6
-       endif
+       tnat(ixt)=155.76E-6; IF(ltnat1) tnat(ixt)=1
        toce(ixt)=tnat(ixt)*(1.0+pente_MWL*deltaO18_oce/1000.0)
        tcorr(ixt)=1.0+fac_enrichoce18*pente_MWL          
Index: LMDZ6/branches/cirrus/libf/phylmdiso/isotopes_routines_mod.F90
===================================================================
--- LMDZ6/branches/cirrus/libf/phylmdiso/isotopes_routines_mod.F90	(revision 5202)
+++ LMDZ6/branches/cirrus/libf/phylmdiso/isotopes_routines_mod.F90	(revision 5203)
@@ -16419,6 +16419,6 @@
    USE isotopes_mod,      ONLY: isoName,iso_HDO,iso_eau
    USE phyetat0_get_mod,  ONLY: phyetat0_get, phyetat0_srf
-   USE readTracFiles_mod, ONLY: new2oldH2O
-   USE strings_mod,       ONLY: strIdx, strTail, maxlen, msg, int2str
+   USE infotrac_phy,      ONLY: new2oldH2O
+   USE strings_mod,       ONLY: strIdx, strHead, strTail, maxlen, msg, int2str
 #ifdef ISOVERIF
    USE isotopes_verif_mod
@@ -16459,6 +16459,5 @@
       outiso = isoName(ixt)
       oldIso = strTail(new2oldH2O(outiso), '_')            !--- Remove "H2O_" from "H2O_<iso>[_<tag>]"
-      i = INDEX(outiso, '_', .TRUE.)
-      oldIso2 = outiso(1:i-1)//outiso(i+1:LEN_TRIM(outiso)) ! CR 2023: on ajoute cette possibilité aussi, elle correspond au cas le plus récent.
+      oldIso2= TRIM(strHead(outiso,'_'))//strTail(outiso,'_') ! CR 2023: most recent possibility
 !      write(*,*) 'tmp 16541:'
 !      write(*,*) 'outiso=',outiso
Index: LMDZ6/branches/cirrus/libf/phylmdiso/isotrac_mod.F90
===================================================================
--- LMDZ6/branches/cirrus/libf/phylmdiso/isotrac_mod.F90	(revision 5202)
+++ LMDZ6/branches/cirrus/libf/phylmdiso/isotrac_mod.F90	(revision 5203)
@@ -3,7 +3,6 @@
 
 MODULE isotrac_mod
-  USE infotrac_phy,      ONLY: niso, ntiso, nzone
-  USE readTracFiles_mod, ONLY: delPhase
-  USE isotopes_mod,      ONLY: ridicule, get_in
+  USE infotrac_phy, ONLY: niso, ntiso, nzone, delPhase
+  USE isotopes_mod, ONLY: ridicule, get_in
 
   IMPLICIT NONE
Index: LMDZ6/branches/cirrus/libf/phylmdiso/phyetat0_mod.F90
===================================================================
--- LMDZ6/branches/cirrus/libf/phylmdiso/phyetat0_mod.F90	(revision 5202)
+++ LMDZ6/branches/cirrus/libf/phylmdiso/phyetat0_mod.F90	(revision 5203)
@@ -41,6 +41,6 @@
   USE geometry_mod,     ONLY: longitude_deg, latitude_deg
   USE iostart,          ONLY: close_startphy, get_field, get_var, open_startphy
-  USE infotrac_phy,     ONLY: nqtot, nbtr, type_trac, tracers
-  USE readTracFiles_mod,ONLY: maxlen, new2oldH2O
+  USE infotrac_phy,     ONLY: nqtot, nbtr, type_trac, tracers, new2oldH2O
+  USE strings_mod,      ONLY: maxlen
   USE traclmdz_mod,     ONLY: traclmdz_from_restart
   USE carbon_cycle_mod, ONLY: carbon_cycle_init, carbon_cycle_cpl, carbon_cycle_tr, carbon_cycle_rad, co2_send, RCO2_glo
Index: LMDZ6/branches/cirrus/libf/phylmdiso/physiq_mod.F90
===================================================================
--- LMDZ6/branches/cirrus/libf/phylmdiso/physiq_mod.F90	(revision 5202)
+++ LMDZ6/branches/cirrus/libf/phylmdiso/physiq_mod.F90	(revision 5203)
@@ -39,6 +39,5 @@
     USE ioipsl_getin_p_mod, ONLY : getin_p
     USE indice_sol_mod
-    USE infotrac_phy, ONLY: nqtot, nbtr, nqo, tracers, type_trac,ivap,iliq,isol
-    USE readTracFiles_mod, ONLY: addPhase
+    USE infotrac_phy, ONLY: nqtot, nbtr, nqo, tracers, type_trac,addPhase, ivap, iliq, isol
     USE strings_mod,  ONLY: strIdx
     USE iophy
@@ -1375,4 +1374,7 @@
 !$OMP THREADPRIVATE(SFRWL)
     !albedo SB <<<
+
+    !--Lea Raillard qs_ini 
+    REAL, dimension(klon,klev) :: qs_ini
 
     !--OB variables for mass fixer (hard coded for now)
@@ -5090,5 +5092,5 @@
 
     CALL lscp(klon,klev,phys_tstep,missing_val,paprs,pplay, &
-         t_seri, q_seri,qs_ancien,ptconv,ratqs, &
+         t_seri, q_seri,qs_ini,ptconv,ratqs, &
          d_t_lsc, d_q_lsc, d_ql_lsc, d_qi_lsc, rneb, rneblsvol, & 
          pfraclr, pfracld, cldfraliq, sigma2_icefracturb, mean_icefracturb,  &
@@ -7172,4 +7174,7 @@
        ENDDO
     ENDDO
+    
+    ! Lea Raillard qs_ini for cloud phase param.
+    qs_ini(:,:)=qs_seri(:,:)
 
     ! C Risi: dispatcher les isotopes dans les xt_seri