Index: readTracFiles_mod.f90
===================================================================
--- readTracFiles_mod.f90	(revision 12)
+++ readTracFiles_mod.f90	(revision 13)
@@ -1,6 +1,6 @@
 MODULE readTracFiles_mod
 
-  USE strings_mod,    ONLY: msg, testFile,  strFind, strStack, strReduce,  strHead, strCount,   find, maxlen, fmsg, &
-             removeComment, cat, checkList, strIdx,  strParse, strReplace, strTail, reduceExpr, test, get_in, dispTable
+  USE strings_mod,    ONLY: msg, testFile,  strFind, strStack, strReduce,  strHead, strCount, find, fmsg, reduceExpr, &
+             removeComment, cat, checkList, str2int, strParse, strReplace, strTail, strIdx, maxlen, test, dispTable, get_in
   USE trac_types_mod, ONLY: trac_type, isot_type, keys_type
 
@@ -12,9 +12,9 @@
   PUBLIC :: readTracersFiles, indexUpdate, setGeneration             !--- TOOLS ASSOCIATED TO TRACERS  DESCRIPTORS
   PUBLIC :: readIsotopesFile                                         !--- TOOLS ASSOCIATED TO ISOTOPES DESCRIPTORS
-  PUBLIC :: getKey_init, getKey, setDirectKeys                       !--- GET/SET KEYS FROM/TO tracers & isotopes
+  PUBLIC :: getKey_init, getKey, fGetKey, setDirectKeys              !--- GET/SET KEYS FROM/TO tracers & isotopes
 
   PUBLIC :: known_phases, old_phases, nphases, phases_names, &       !--- VARIABLES RELATED TO THE PHASES
-            phases_sep, delPhase, addPhase, &                        !--- + ROUTINES TO ADD/REMOVE PHASE TO/FROM A NAME
-            old2newPhase, new2oldPhase
+            phases_sep, delPhase, addPhase, new2oldPhase,    &       !--- + ROUTINES TO ADD/REMOVE PHASE TO/FROM A NAME
+            old2newName, new2oldName
 
   PUBLIC :: tran0, idxAncestor, ancestor                             !--- GENERATION 0 TRACER + TOOLS FOR GENERATIONS
@@ -30,4 +30,5 @@
   END INTERFACE getKey
 !------------------------------------------------------------------------------------------------------------------------------
+  INTERFACE fGetKey;       MODULE PROCEDURE fgetKeyByIndex_s1, fgetKeyByName_s1; END INTERFACE fGetKey
   INTERFACE tracersSubset; MODULE PROCEDURE trSubset_Indx, trSubset_Name, trSubset_gen0Name; END INTERFACE tracersSubset
   INTERFACE idxAncestor;   MODULE PROCEDURE idxAncestor_1, idxAncestor_m; END INTERFACE idxAncestor
@@ -41,9 +42,9 @@
   !--- SOME PARAMETERS THAT ARE NOT LIKELY TO CHANGE OFTEN
   CHARACTER(LEN=maxlen), SAVE      :: tran0        = 'air'           !--- Default transporting fluid
-  CHARACTER(LEN=maxlen), PARAMETER :: old_phases   = 'vlir'          !--- Old phases for water (no separator)
-  CHARACTER(LEN=maxlen), PARAMETER :: known_phases = 'glsr'          !--- Known phases initials
+  CHARACTER(LEN=maxlen), PARAMETER :: old_phases   = 'vli'           !--- Old phases for water (no separator)
+  CHARACTER(LEN=maxlen), PARAMETER :: known_phases = 'gls'           !--- Known phases initials
   INTEGER,               PARAMETER :: nphases=LEN_TRIM(known_phases) !--- Number of phases
   CHARACTER(LEN=maxlen), SAVE      :: phases_names(nphases) &        !--- Known phases names
-                                = ['gaseous', 'liquid ', 'solid  ', 'cloud  ']
+                                = ['gaseous', 'liquid ', 'solid  ']
   CHARACTER(LEN=1), SAVE :: phases_sep  =  '_'                       !--- Phase separator
   LOGICAL,          SAVE :: tracs_merge = .TRUE.                     !--- Merge/stack tracers lists
@@ -87,6 +88,6 @@
   INTEGER,                      INTENT(OUT) :: fType                  !--- Type of input file found
   TYPE(trac_type), ALLOCATABLE, INTENT(OUT) :: tracs(:)
-  CHARACTER(LEN=maxlen),  ALLOCATABLE ::  s(:), sections(:), trac_files(:)
-  CHARACTER(LEN=maxlen) :: str, fname, mesg, oldH2O, newH2O
+  CHARACTER(LEN=maxlen),  ALLOCATABLE :: s(:), sections(:), trac_files(:)
+  CHARACTER(LEN=maxlen) :: str, fname, mesg
   INTEGER               :: is, nsec, ierr, it, ntrac, ns, ip, ix
   LOGICAL, ALLOCATABLE  :: ll(:), lGen3(:)
@@ -142,8 +143,8 @@
         CALL msg('This file is for air tracers only',           modname, ns == 3 .AND. it == 1)
         CALL msg('This files specifies the transporting fluid', modname, ns == 4 .AND. it == 1)
-        tracs(it)%name = TRIM(s(3))                                  !--- Set %name:   name of the tracer
-        tracs(it)%parent = tran0                                     !--- Set %parent: transporting fluid
-        IF(ns == 4) tracs(it)%parent = s(4)                          !---     default: 'air' or defined in the file
-        tracs(it)%phase = known_phases(1:1)                          !--- Set %phase:  tracer phase (default: "g"azeous)
+        tracs(it)%name   = old2newName(s(3), ip)                     !--- Set %name:   name   of the tracer
+        tracs(it)%parent = tran0                                     !--- Default transporting fluid name
+        IF(ns == 4) tracs(it)%parent = old2newName(s(4))             !--- Set %parent: parent of the tracer
+        tracs(it)%phase = known_phases(ip:ip)                        !--- Set %phase:  tracer phase (default: "g"azeous)
         tracs(it)%component = TRIM(type_trac)                        !--- Set %component: model component name
         tracs(it)%keys%key = ['hadv', 'vadv']                        !--- Set %keys%key
@@ -151,15 +152,6 @@
       END DO
       CLOSE(90)
-      DO ip = 1, nphases                                             !--- Deal with old water names
-        oldH2O = addPhase('H2O', ip, '')
-        newH2O = addPhase('H2O', ip)
-        ix = strIdx(tracs(:)%name, oldH2O)
-        IF(ix == 0) CYCLE
-        tracs(ix)%name  = newH2O                                     !--- Set %name:   name of the tracer
-        WHERE(tracs(:)%parent == oldH2O) tracs(:)%parent = newH2O    !--- Set %parent: transporting fluid
-        tracs(ix)%phase = known_phases(ip:ip)                        !--- Set %phase:  tracer phase
-      END DO
       CALL setGeneration(tracs)                                      !--- Set %iGeneration and %gen0Name
-      WHERE(tracs%iGeneration == 3) tracs%type = 'tag'               !--- Set %type:        'tracer' or 'tag'
+      WHERE(tracs%iGeneration == 2) tracs%type = 'tag'               !--- Set %type:        'tracer' or 'tag'
       IF(test(checkTracers(tracs, fname, fname), lerr)) RETURN       !--- Detect orphans and check phases
       IF(test(checkUnique (tracs, fname, fname), lerr)) RETURN       !--- Detect repeated tracers
@@ -167,22 +159,24 @@
       tracs(:)%keys%name = tracs(:)%name                             !--- Copy tracers names in keys components
     !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
-    CASE(2); IF(test(feedDBase(["tracer.def"],[type_trac]), lerr)) RETURN  !=== SINGLE FILE, COMA-SEPARATED SECTIONS LIST
+    CASE(2); IF(test(feedDBase(["tracer.def"], [type_trac], modname), lerr)) RETURN !=== SINGLE   FILE, MULTIPLE SECTIONS
     !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
-    CASE(3); IF(test(feedDBase(  trac_files  , sections  ), lerr)) RETURN  !=== MULTIPLE FILES, ONE SECTION EACH FILE
+    CASE(3); IF(test(feedDBase(  trac_files  ,  sections,   modname), lerr)) RETURN !=== MULTIPLE FILES, SINGLE  SECTION
   !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
   END SELECT
   !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 
-  IF(ANY([2,3] == fType) .AND. nsec > 1) THEN
-    IF(tracs_merge) THEN
-      CALL msg('The multiple required sections will be MERGED.',    modname)
-      IF(test(mergeTracers(dBase, tracs), lerr)) RETURN
-    ELSE
-      CALL msg('The multiple required sections will be CUMULATED.', modname)
-      IF(test(cumulTracers(dBase, tracs), lerr)) RETURN
-    END IF
-    WHERE(tracs%gen0Name(1:3) /= 'H2O') tracs%isInPhysics=.TRUE.     !--- Set %isInPhysics: passed to physics
-    CALL setDirectKeys(tracs)                                        !--- Set %iqParent, %iqDescen, %nqDescen, %nqChilds
+  IF(ALL([2,3] /= fType)) RETURN
+
+  IF(nsec  == 1) THEN; 
+    tracs = dBase(1)%trac
+  ELSE IF(tracs_merge) THEN
+    CALL msg('The multiple required sections will be MERGED.',    modname)
+    IF(test(mergeTracers(dBase, tracs), lerr)) RETURN
+  ELSE
+    CALL msg('The multiple required sections will be CUMULATED.', modname)
+    IF(test(cumulTracers(dBase, tracs), lerr)) RETURN
   END IF
+  WHERE(tracs%gen0Name(1:3) /= 'H2O') tracs%isInPhysics=.TRUE.       !--- Set %isInPhysics: passed to physics
+  CALL setDirectKeys(tracs)                                          !--- Set %iqParent, %iqDescen, %nqDescen, %nqChilds
 
 END FUNCTION readTracersFiles
@@ -190,5 +184,5 @@
 
 !==============================================================================================================================
-LOGICAL FUNCTION feedDBase(fnames, snames) RESULT(lerr)
+LOGICAL FUNCTION feedDBase(fnames, snames, modname) RESULT(lerr)
 ! Purpose: Read the sections "snames(is)" (coma-separated list) from each "fnames(is)"
 !   file and create the corresponding tracers set descriptors in the database "dBase":
@@ -200,12 +194,12 @@
   CHARACTER(LEN=*), INTENT(IN)  :: fnames(:)                         !--- Files names
   CHARACTER(LEN=*), INTENT(IN)  :: snames(:)                         !--- Coma-deparated list of sections (one list each file)
-  INTEGER,  ALLOCATABLE :: ndb(:)                                    !--- Nuber of sections for each file
+  CHARACTER(LEN=*), INTENT(IN)  :: modname                           !--- Calling routine name
+  INTEGER,  ALLOCATABLE :: ndb(:)                                    !--- Number of sections for each file
   INTEGER,  ALLOCATABLE :: ixf(:)                                    !--- File index for each section of the expanded list
   LOGICAL,  ALLOCATABLE :: lTg(:)                                    !--- Tagging tracers mask
-  CHARACTER(LEN=maxlen) :: fnm, snm, modname
+  CHARACTER(LEN=maxlen) :: fnm, snm
   INTEGER               :: idb, i
   LOGICAL :: ll
 !------------------------------------------------------------------------------------------------------------------------------
-  modname = 'feedDBase'
   !=== READ THE REQUIRED SECTIONS
   ll = strCount(snames, ',', ndb)                                    !--- Number of sections for each file
@@ -219,4 +213,5 @@
   !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
     fnm = fnames(ixf(idb)); snm = dBase(idb)%name                    !--- FILE AND SECTION NAMES
+    lerr = ANY([(dispTraSection('RAW CONTENT OF SECTION "'//TRIM(snm)//'"', snm, modname), idb=1, SIZE(dBase))])
     IF(test(expandSection(dBase(idb)%trac, snm, fnm), lerr)) RETURN  !--- EXPAND NAMES ;  set %parent, %type, %component
     CALL setGeneration   (dBase(idb)%trac)                           !---                 set %iGeneration,   %genOName
@@ -225,11 +220,8 @@
     CALL expandPhases    (dBase(idb)%trac)                           !--- EXPAND PHASES ; set %phase
     CALL sortTracers     (dBase(idb)%trac)                           !--- SORT TRACERS
+    lerr = ANY([(dispTraSection('EXPANDED CONTENT OF SECTION "'//TRIM(snm)//'"', snm, modname), idb=1, SIZE(dBase))])
   !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
   END DO
   !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
-
-  !=== DISPLAY BASIC INFORMATION
-  lerr = ANY([( dispTraSection('Expanded list for section "'//TRIM(dBase(idb)%name)//'"', dBase(idb)%name, modname), &
-                idb=1, SIZE(dBase) )])
 END FUNCTION feedDBase
 !------------------------------------------------------------------------------------------------------------------------------
@@ -406,14 +398,16 @@
   DO it = 1, nt                                                      !=== EXPAND TRACERS AND PARENTS NAMES LISTS
   !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
-    ll = strParse(tr(it)%name, ',', ta, n=ntr)                       !--- Number of tracers
+    ll = strParse(tr(it)%name,   ',', ta, n=ntr)                     !--- Number of tracers
     ll = strParse(tr(it)%parent, ',', pa, n=npr)                     !--- Number of parents
     DO ipr=1,npr                                                     !--- Loop on parents list elts
       DO itr=1,ntr                                                   !--- Loop on tracers list elts
         i = iq+itr-1+(ipr-1)*ntr
-        ttr(i)%name = TRIM(ta(itr)); ttr(i)%parent = pa(ipr)
-        ttr(i)%keys = keys_type(ta(itr), tr(it)%keys%key, tr(it)%keys%val)
+        ttr(i)%name   = TRIM(ta(itr))
+        ttr(i)%parent = TRIM(pa(ipr))
+        ttr(i)%keys   = keys_type(ta(itr), tr(it)%keys%key, tr(it)%keys%val)
       END DO
     END DO
-    ttr(iq:iq+ntr*npr-1)%type = tr(it)%type                          !--- Duplicating type
+    ttr(iq:iq+ntr*npr-1)%type      = tr(it)%type                     !--- Duplicating type
+    ttr(iq:iq+ntr*npr-1)%component = tr(it)%component                !--- Duplicating type
     iq = iq + ntr*npr
   !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
@@ -440,6 +434,6 @@
   tr(:)%iGeneration = -1                                             !--- error if -1
   nq = SIZE(tr, DIM=1)                                               !--- Number of tracers lines
-  lg = tr(:)%parent == tran0                                         !--- First generation tracers flag
-  WHERE(lg) tr(:)%iGeneration = 0                                    !--- First generation tracers
+  lg = tr(:)%parent == tran0                                         !--- Flag for generation 0 tracers
+  WHERE(lg) tr(:)%iGeneration = 0                                    !--- Generation 0 tracers
 
   !=== Determine generation for each tracer
@@ -511,6 +505,6 @@
     ll = tr(:)%name==TRIM(tnam)                                      !--- Mask for current tracer name
     IF(COUNT(ll)==1 ) CYCLE                                          !--- Tracer is not repeated
-    IF(tr(iq)%iGeneration>1) THEN
-      tdup(iq) = tnam                                                !--- gen>1: MUST be unique
+    IF(tr(iq)%iGeneration>0) THEN
+      tdup(iq) = tnam                                                !--- gen>0: MUST be unique
     ELSE
       DO ip=1,nphases; p=known_phases(ip:ip)                         !--- Loop on known phases
@@ -531,5 +525,5 @@
 SUBROUTINE expandPhases(tr)
 !------------------------------------------------------------------------------------------------------------------------------
-! Purpose: Expand the phases in the tracers descriptor "tr".
+! Purpose: Expand the phases in the tracers descriptor "tr". Phases are not repeated for a tracer, thanks to "checkUnique".
 !------------------------------------------------------------------------------------------------------------------------------
   TYPE(trac_type), ALLOCATABLE, INTENT(INOUT) :: tr(:)               !--- Tracer derived type vector
@@ -538,4 +532,5 @@
   INTEGER,   ALLOCATABLE ::  i0(:)
   CHARACTER(LEN=maxlen)  :: nam, pha, trn
+  CHARACTER(LEN=1) :: p
   INTEGER :: ip, np, iq, jq, nq, it, nt, nc, i, n
   LOGICAL :: lTg, lEx
@@ -544,38 +539,39 @@
   nt = 0
   DO iq = 1, nq                                                      !--- GET THE NUMBER OF TRACERS
-    IF(tr(iq)%iGeneration /= 1) CYCLE
-    nc = COUNT(tr(:)%gen0Name==tr(iq)%name .AND. tr%iGeneration/=1)  !--- Number of childs of tr(iq)
+    IF(tr(iq)%iGeneration /= 0) CYCLE                                !--- Only deal with generation 0 tracers
+    nc = COUNT(tr(:)%gen0Name==tr(iq)%name .AND. tr%iGeneration/=0)  !--- Number of childs of tr(iq)
     tr(iq)%phase = fgetKey(iq, 'phases', tr(:)%keys)                 !--- Phases list      of tr(iq)
     np = LEN_TRIM(tr(iq)%phase)                                      !--- Number of phases of tr(iq)
     nt = nt + (1+nc) * np                                            !--- Number of tracers after expansion
   END DO
-  ALLOCATE(ttr(nt))
+  ALLOCATE(ttr(nt))                                                  !--- Version  of "tr" after phases expansion
   it = 1                                                             !--- Current "ttr(:)" index
   DO iq = 1, nq                                                      !--- Loop on "tr(:)" indexes
     lTg = tr(iq)%type=='tag'                                         !--- Current tracer is a tag
     i0 = strFind(tr(:)%name, TRIM(tr(iq)%gen0Name), n)               !--- Indexes of first generation ancestor copies
-    np = SUM( [( LEN_TRIM(tr(i0(i))%phase),i=1,n )],1)               !--- Number of phases for current tracer tr(iq)
-    lEx = np>1                                                       !--- Need of a phase suffix
-    IF(lTg) lEx=lEx.AND.tr(iq)%iGeneration>1                         !--- No phase suffix for first generation tags
-    DO i=1,n                                                         !=== LOOP ON FIRST GENERATION ANCESTORS
-      jq=i0(i)                                                       !--- tr(jq): ith copy of 1st gen. ancestor of tr(iq)
-      IF(tr(iq)%iGeneration==1) jq=iq                                !--- Generation 1: current tracer phases only
+    np = SUM([( LEN_TRIM(tr(i0(i))%phase),i=1,n )], 1)               !--- Number of phases for current tracer tr(iq)
+    lEx = np>1                                                       !--- Phase suffix only required if phases number is > 1
+    IF(lTg) lEx = lEx .AND. tr(iq)%iGeneration>0                     !--- No phase suffix for generation 0 tags
+    DO i=1,n                                                         !=== LOOP ON GENERATION 0 ANCESTORS
+      jq = i0(i)                                                     !--- tr(jq): ith tracer with same gen 0 ancestor as tr(iq)
+      IF(tr(iq)%iGeneration==0) jq=iq                                !--- Generation 0: count the current tracer phases only
       pha = tr(jq)%phase                                             !--- Phases list for tr(jq)
-      DO ip=1,LEN_TRIM(pha)                                          !=== LOOP ON PHASES LISTS
-        trn=TRIM(tr(iq)%name); nam=trn                               !--- Tracer name (regular case)
+      DO ip = 1, LEN_TRIM(pha)                                       !=== LOOP ON PHASES LISTS
+        p = pha(ip:ip)
+        trn = TRIM(tr(iq)%name); nam = trn                           !--- Tracer name (regular case)
         IF(lTg) nam = TRIM(tr(iq)%parent)                            !--- Parent name (tagging case)
-        IF(lEx) nam = addPhase(nam, ip)                              !--- Phase extension needed
+        IF(lEx) nam = addPhase(nam, p )                              !--- Phase extension needed
         IF(lTg) nam = TRIM(nam)//'_'//TRIM(trn)                      !--- <parent>_<name> for tags
         ttr(it) = tr(iq)                                             !--- Same <key>=<val> pairs
-        ttr(it)%name = TRIM(nam)                                     !--- Name with possibly phase suffix
+        ttr(it)%name      = TRIM(nam)                                !--- Name with possibly phase suffix
         ttr(it)%keys%name = TRIM(nam)                                !--- Name inside the keys decriptor
-        ttr(it)%phase = pha(ip:ip)                                   !--- Single phase entry
-        IF(lEx.AND.tr(iq)%iGeneration>1) THEN
-          ttr(it)%parent   = addPhase(ttr(it)%parent,   ip)          !--- Modify parent name
-          ttr(it)%gen0Name = addPhase(ttr(it)%gen0Name, ip)          !--- Modify generation 0 ancestor name
+        ttr(it)%phase     = p                                        !--- Single phase entry
+        IF(lEx .AND. tr(iq)%iGeneration>0) THEN
+          ttr(it)%parent   = addPhase(ttr(it)%parent,   p)
+          ttr(it)%gen0Name = addPhase(ttr(it)%gen0Name, p)
         END IF
-        it=it+1
+        it = it+1
       END DO
-      IF(tr(iq)%iGeneration==1) EXIT                                 !--- Break phase loop for gen 1
+      IF(tr(iq)%iGeneration==0) EXIT                                 !--- Break phase loop for gen 0
     END DO
   END DO
@@ -590,22 +586,23 @@
 !------------------------------------------------------------------------------------------------------------------------------
 ! Purpose: Sort tracers:
-!  * Put water at first places, in the "known_phases" order.
+!  * Put water at the beginning of the vector, in the "known_phases" order.
 !  * lGrowGen == T: in ascending generations numbers.
 !  * lGrowGen == F: tracer + its childs sorted by growing generation, one after the other.
 !   TO BE ADDED IF NECESSARY: HIGHER MOMENTS AT THE END
 !------------------------------------------------------------------------------------------------------------------------------
-  TYPE(trac_type), ALLOCATABLE, INTENT(INOUT) :: tr(:)               !--- Tracer derived type vector
+  TYPE(trac_type), ALLOCATABLE, INTENT(INOUT) :: tr(:)         !--- Tracer derived type vector
+!  TYPE(trac_type), ALLOCATABLE :: ttr(:)
+  INTEGER,         ALLOCATABLE :: iy(:), iz(:)
   INTEGER :: ig, ng, iq, jq, ip, nq, n, ix(SIZE(tr)), k
-  INTEGER, ALLOCATABLE :: iy(:), iz(:)
 !------------------------------------------------------------------------------------------------------------------------------
   nq = SIZE(tr)
-  iy = [(k, k=1, nq)]
   DO ip = nphases, 1, -1
     iq = strIdx(tr(:)%name, addPhase('H2O', ip))
-    IF(iq/=0) iy = [iq, iy(1:iq-1), iy(iq+1:nq)]
-  END DO
-  tr = tr(iy)                                                        !--- Water displaces at first positions
-  iq = 1
+    IF(iq == 0) CYCLE
+    tr = tr([iq, 1:iq-1, iq+1:nq])
+!    tr(:)%name = nam
+  END DO
   IF(lSortByGen) THEN
+    iq = 1
     ng = MAXVAL(tr(:)%iGeneration, MASK=.TRUE., DIM=1)               !--- Number of generations
     DO ig = 0, ng                                                    !--- Loop on generations
@@ -616,11 +613,12 @@
     END DO
   ELSE
-    DO jq = 1, nq                                                    !--- Loop on first generation tracers
-      IF(tr(jq)%iGeneration /= 1) CYCLE                              !--- Skip generations >= 1
-      ix(iq) = jq                                                    !--- First generation ancestor index first
-      iq = iq + 1
+    iq = 1
+    DO jq = 1, nq                                                    !--- Loop on generation 0 tracers
+      IF(tr(jq)%iGeneration /= 0) CYCLE                              !--- Skip generations /= 0
+      ix(iq) = jq                                                    !--- Generation 0 ancestor index first
+      iq = iq + 1                                                    !--- Next "iq" for next generations tracers
       iy = strFind(tr(:)%gen0Name, TRIM(tr(jq)%name))                !--- Indexes of "tr(jq)" childs in "tr(:)"
-      ng = MAXVAL(tr(iy)%iGeneration, MASK=.TRUE., DIM=1)            !--- Generations number of the "tr(jq)" family
-      DO ig = 2, ng                                                  !--- Loop   on generations for the tr(jq) family
+      ng = MAXVAL(tr(iy)%iGeneration, MASK=.TRUE., DIM=1)            !--- Number of generations of the "tr(jq)" family
+      DO ig = 1, ng                                                  !--- Loop   on generations of the "tr(jq)" family
         iz = find(tr(iy)%iGeneration, ig, n)                         !--- Indexes of the tracers "tr(iy(:))" of generation "ig"
         ix(iq:iq+n-1) = iy(iz)                                       !--- Same indexes in "tr(:)"
@@ -757,4 +755,5 @@
   INTEGER :: idb, iq, nq
   INTEGER, ALLOCATABLE :: hadv(:), vadv(:)
+  CHARACTER(LEN=maxlen), ALLOCATABLE :: phas(:)
   TYPE(trac_type), POINTER :: tm(:)
   lerr = .FALSE.
@@ -762,9 +761,17 @@
   tm => dBase(idb)%trac
   nq = SIZE(tm)
-  IF(test(getKeyByName_im('hadv', hadv, tm(:)%name, tm(:)%keys),lerr)) RETURN
-  IF(test(getKeyByName_im('vadv', vadv, tm(:)%name, tm(:)%keys),lerr)) RETURN
+  !--- BEWARE ! Can't always already use the "getKeyByName" functions.
+  !             Names must first include the phases for tracers defined on multiple lines.
+  hadv = str2int([(fgetKey(iq, 'hadv',  tm(:)%keys, '10'), iq=1, nq)])
+  vadv = str2int([(fgetKey(iq, 'vadv',  tm(:)%keys, '10'), iq=1, nq)])
+  phas =         [(fgetKey(iq, 'phases',tm(:)%keys, 'g' ), iq=1, nq)]
   CALL msg(TRIM(message)//':', modname)
-  IF(test(dispTable('iiissis', ['iq        ','hadv      ','vadv      ','short name','parent    ','igen      ','phase     '], &
-    cat(tm(:)%name,  tm(:)%parent, tm(:)%phase), cat([(iq, iq=1, nq)],  hadv,  vadv, tm(:)%iGeneration)), lerr)) RETURN
+  IF(tm(1)%parent == '') THEN
+    IF(test(dispTable('iiiss',   ['iq    ','hadv  ','vadv  ','name  ','phase '], cat(tm%name, phas), cat([(iq, iq=1, nq)], &
+                                            hadv,    vadv),                 sub=modname), lerr)) RETURN
+  ELSE
+    IF(test(dispTable('iiissis', ['iq    ','hadv  ','vadv  ','name  ','parent','igen  ','phase '], cat(tm%name, tm%parent, &
+           tm%phase), cat([(iq, iq=1, nq)], hadv,    vadv, tm%iGeneration), sub=modname), lerr)) RETURN
+  END IF
 END FUNCTION dispTraSection
 !==============================================================================================================================
@@ -847,5 +854,5 @@
 !=== NOTES:                                                                                                                ====
 !===  * Most of the "isot" components have been defined in the calling routine (initIsotopes):                             ====
-!===      parent,  nzone, zone(:),  niso, keys(:)%name,  ntiso, trac(:),  nphas, phas,  iTraPha(:,:),  iZonPhi(:,:)        ====
+!===      parent,  nzone, zone(:),  niso, keys(:)%name,  ntiso, trac(:),  nphas, phas,  iqTraPha(:,:),  itZonPhi(:,:)      ====
 !===  * Same syntax for isotopes file and "tracer.def": a tracers section contains one line for each of its isotopes       ====
 !===  * Each tracers section can contain a "params" virtual isotope line of isotopes parameters default values             ====
@@ -909,4 +916,8 @@
     ALLOCATE(tdb(nb0-1)); tdb(1:nb0-1)=dBase(1:nb0-1); CALL MOVE_ALLOC(FROM=tdb, TO=dBase)
   END IF
+
+  !--- GET THE isoCheck ENTRY FROM THE *.DEF FILES (MIGHT BE CHANGED TO A CLASS-DEPENDANT KEYWORD)
+  CALL get_in('ok_iso_verif', isot(strIdx(isot%parent, 'H2O'))%check, .FALSE.)
+
   lerr = dispIsotopes(isot, 'Isotopes parameters read from file', modname)
 
@@ -930,9 +941,9 @@
   LOGICAL, ALLOCATABLE :: ll(:)                                      !--- Mask
   TYPE(trac_type), POINTER   ::  t(:), t1
-  TYPE(isot_type), POINTER   ::  s
+  TYPE(isot_type), POINTER   ::  i
 
   t => trac
 
-  p = PACK(delPhase(t%parent), MASK = t%type=='tracer' .AND. t%iGeneration==2) !--- Parents of 2nd generation isotopes
+  p = PACK(delPhase(t%parent), MASK = t%type=='tracer' .AND. t%iGeneration==1) !--- Parents of generation 1 isotopes
   CALL strReduce(p, nbIso)
   ALLOCATE(isot(nbIso))
@@ -943,19 +954,19 @@
   isot(:)%parent = p
   DO ic = 1, SIZE(p)                                                 !--- Loop on isotopes classes
-    s => isot(ic)
-    iname = s%parent                                                 !--- Current isotopes class name (parent tracer name)
+    i => isot(ic)
+    iname = i%parent                                                 !--- Current isotopes class name (parent tracer name)
 
     !=== Isotopes childs of tracer "iname": mask, names, number (same for each phase of "iname")
     ll = t(:)%type=='tracer' .AND. delPhase(t(:)%parent) == iname .AND. t(:)%phase == 'g'
     str = PACK(delPhase(t(:)%name), MASK = ll)                       !--- Effectively found isotopes of "iname"
-    s%niso = SIZE(str)                                               !--- Number of "effectively found isotopes of "iname"
-    ALLOCATE(s%keys(s%niso))
-    FORALL(it = 1:s%niso) s%keys(it)%name = str(it)
+    i%niso = SIZE(str)                                               !--- Number of "effectively found isotopes of "iname"
+    ALLOCATE(i%keys(i%niso))
+    FORALL(it = 1:i%niso) i%keys(it)%name = str(it)
 
     !=== Geographic tagging tracers descending on tracer "iname": mask, names, number
-    ll = t(:)%type=='tag'    .AND. delPhase(t(:)%gen0Name) == iname .AND. t(:)%iGeneration == 3
-    s%zone = PACK(strTail(t(:)%name,'_',lFirst=.TRUE.), MASK = ll)   !--- Tagging zones names  for isotopes category "iname"
-    CALL strReduce(s%zone)
-    s%nzone = SIZE(s%zone)                                           !--- Tagging zones number for isotopes category "iname"
+    ll = t(:)%type=='tag'    .AND. delPhase(t(:)%gen0Name) == iname .AND. t(:)%iGeneration == 2
+    i%zone = PACK(strTail(t(:)%name,'_',lFirst=.TRUE.), MASK = ll)   !--- Tagging zones names  for isotopes category "iname"
+    CALL strReduce(i%zone)
+    i%nzone = SIZE(i%zone)                                           !--- Tagging zones number for isotopes category "iname"
 
     !=== Geographic tracers of the isotopes childs of tracer "iname" (same for each phase of "iname")
@@ -963,33 +974,32 @@
     str = PACK(delPhase(t(:)%name), MASK=ll)
     CALL strReduce(str)
-    s%ntiso = s%niso + SIZE(str)                                     !--- Number of isotopes + their geographic tracers [ntraciso]
-    ALLOCATE(s%trac(s%ntiso))
-    FORALL(it = 1:s%niso) s%trac(it) = s%keys(it)%name
-    FORALL(it = s%niso+1:s%ntiso) s%trac(it) = str(it-s%niso)
+    i%ntiso = i%niso + SIZE(str)                                     !--- Number of isotopes + their geographic tracers [ntraciso]
+    ALLOCATE(i%trac(i%ntiso))
+    FORALL(it = 1:i%niso) i%trac(it) = i%keys(it)%name
+    FORALL(it = i%niso+1:i%ntiso) i%trac(it) = str(it-i%niso)
 
     !=== Phases for tracer "iname"
-    s%phase = ''
-    DO ip = 1, nphases; ph = known_phases(ip:ip); IF(strIdx(t%name,addPhase(iname, ph)) /= 0) s%phase = TRIM(s%phase)//ph; END DO
-    s%nphas = LEN_TRIM(s%phase)                                       !--- Equal to "nqo" for water
+    i%phase = ''
+    DO ip = 1, nphases; ph = known_phases(ip:ip); IF(strIdx(t%name,addPhase(iname, ph)) /= 0) i%phase = TRIM(i%phase)//ph; END DO
+    i%nphas = LEN_TRIM(i%phase)                                       !--- Equal to "nqo" for water
 
     !=== Tables giving the index in a table of effectively found items for each dynamical tracer (1<=iq<=nqtot)
     DO iq = 1, SIZE(t)
       t1 => trac(iq)
-      IF(delPhase(t1%gen0Name) /= iname) CYCLE                       !--- Only deal with tracers descending on "iname"
+      IF(delPhase(t1%gen0Name)/=iname .OR. t1%iGeneration==0) CYCLE  !--- Only deal with tracers descending on "iname"
       t1%iso_iGroup = ic                                             !--- Isotopes family       idx in list "isotopes(:)%parent"
-      t1%iso_iName  = strIdx(s%trac, delPhase(strHead(t1%name,'_'))) !--- Current isotope       idx in effective isotopes list
-      t1%iso_iZone  = strIdx(s%zone,          strTail(t1%name,'_') ) !--- Current isotope zone  idx in effective zones    list
-      t1%iso_iPhase =  INDEX(s%phase,TRIM(t1%phase))                 !--- Current isotope phase idx in effective phases   list
-      IF(t1%iGeneration /= 3) t1%iso_iZone = 0                       !--- Skip possible generation 2 tagging tracers
+      t1%iso_iName  = strIdx(i%trac, delPhase(strHead(t1%name,'_'))) !--- Current isotope       idx in effective isotopes list
+      t1%iso_iZone  = strIdx(i%zone,          strTail(t1%name,'_') ) !--- Current isotope zone  idx in effective zones    list
+      t1%iso_iPhase =  INDEX(i%phase,TRIM(t1%phase))                 !--- Current isotope phase idx in effective phases   list
+      IF(t1%iGeneration /= 2) t1%iso_iZone = 0                       !--- Skip possible generation 1 tagging tracers
     END DO
 
     !=== Table used to get iq (index in dyn array, size nqtot) from the isotope and phase indexes ; the full isotopes list
     !    (including tagging tracers) is sorted this way:  iso1, iso2, ..., iso1_zone1, iso2_zone1, ..., iso1_zoneN, iso2_zoneN
-    s%iTraPha = RESHAPE( [( (strIdx(t(:)%name,  addPhase(s%trac(it),s%phase(ip:ip))),    it=1, s%ntiso), ip=1, s%nphas)], &
-                         [s%ntiso, s%nphas] )
-
+    i%iqTraPha = RESHAPE( [( (strIdx(t%name,  addPhase(i%trac(it),i%phase(ip:ip))),    it=1, i%ntiso), ip=1, i%nphas)], &
+                         [i%ntiso, i%nphas] )
     !=== Table used to get ix (index in tagging tracers isotopes list, size ntiso) from the zone and isotope indexes
-    s%iZonIso = RESHAPE( [( (strIdx(s%trac(:), TRIM(s%trac(it))//'_'//TRIM(s%zone(iz))), iz=1, s%nzone), it=1, s%niso )], &
-                         [s%nzone, s%niso] )
+    i%itZonIso = RESHAPE( [( (strIdx(i%trac(:), TRIM(i%trac(it))//'_'//TRIM(i%zone(iz))), iz=1, i%nzone), it=1, i%niso )], &
+                         [i%nzone, i%niso] )
   END DO
 
@@ -1023,6 +1033,6 @@
       END DO
     END DO
-    IF(test(fmsg('Problem with the table content', modname, dispTable(prf, ttl, val, cat([(it,it=1,nt)]), rFmt='(EN8.4)')), &
-       lerr)) RETURN
+    IF(test(fmsg('Problem with the table content', modname, dispTable(prf, ttl, val, cat([(it,it=1,nt)]), rFmt='(EN8.4)',     &
+       sub=modname)), lerr)) RETURN
     DEALLOCATE(ttl, val)
   END DO        
@@ -1078,6 +1088,6 @@
   IF(jd == 0) RETURN
   DO ik = 1, SIZE(t(jd)%keys%key)
-    CALL get_in(t(jd)%keys%key(ik), val, 'zzzz')
-    IF(val /= 'zzzz') CALL addKey_1(t(jd)%keys%key(ik), val, t(jd)%keys, .TRUE.)
+    CALL get_in(t(jd)%keys%key(ik), val, '*none*')
+    IF(val /= '*none*') CALL addKey_1(t(jd)%keys%key(ik), val, t(jd)%keys, .TRUE.)
   END DO
 END SUBROUTINE addKeysFromDef
@@ -1127,7 +1137,7 @@
 END SUBROUTINE getKey_init
 !==============================================================================================================================
-CHARACTER(LEN=maxlen) FUNCTION fgetKey(itr, keyn, ky, def_val) RESULT(out)
-!------------------------------------------------------------------------------------------------------------------------------
-! Purpose: Internal function ; get a key value in string format (this is the returned argument).
+CHARACTER(LEN=maxlen) FUNCTION fgetKeyByIndex_s1(itr, keyn, ky, def_val) RESULT(val)
+!------------------------------------------------------------------------------------------------------------------------------
+! Purpose: Internal function ; get a string formatted key value (returned argument) from its name and the tracer index.
 !------------------------------------------------------------------------------------------------------------------------------
   INTEGER,                    INTENT(IN) :: itr
@@ -1136,9 +1146,27 @@
   CHARACTER(LEN=*), OPTIONAL, INTENT(IN) :: def_val
 !------------------------------------------------------------------------------------------------------------------------------
-  INTEGER :: ik
-  ik = 0; IF(itr>0 .AND. itr<=SIZE(ky)) ik = strIdx(ky(itr)%key(:), keyn)
-  out = '';              IF(ik /= 0) out = ky(itr)%val(ik)           !--- Key was found
-  IF(PRESENT(def_val) .AND. ik == 0) out = def_val                   !--- Default value from arguments
-END FUNCTION fgetKey
+  INTEGER :: iky
+  iky = 0;  IF(itr >  0 .AND. itr <= SIZE(ky)) iky = strIdx(ky(itr)%key(:), keyn)
+  val = ''; IF(iky /= 0) val = ky(itr)%val(iky)                      !--- Key was found
+  IF(PRESENT(def_val) .AND. iky == 0) val = def_val                  !--- Default value from arguments
+END FUNCTION fgetKeyByIndex_s1
+!==============================================================================================================================
+CHARACTER(LEN=maxlen) FUNCTION fgetKeyByName_s1(tname, keyn, ky, def_val, lerr) RESULT(val)
+!------------------------------------------------------------------------------------------------------------------------------
+! Purpose: Internal function ; get a string formatted key value (returned argument) from its name and the tracer name.
+!------------------------------------------------------------------------------------------------------------------------------
+  CHARACTER(LEN=*),           INTENT(IN)  :: tname, keyn
+  TYPE(keys_type),            INTENT(IN)  :: ky(:)
+  CHARACTER(LEN=*), OPTIONAL, INTENT(IN)  :: def_val
+  LOGICAL,          OPTIONAL, INTENT(OUT) :: lerr
+!------------------------------------------------------------------------------------------------------------------------------
+  INTEGER :: iky, itr
+  val = ''; iky = 0
+  itr = strIdx(ky(:)%name, tname)                                    !--- Get the index of the wanted tracer
+  IF(PRESENT(lerr)) lerr = itr==0; IF(itr == 0) RETURN
+  IF(itr >  0 .AND. itr <= SIZE(ky)) iky = strIdx(ky(itr)%key(:), keyn)
+  IF(iky /= 0) val = ky(itr)%val(iky)                                !--- Key was found
+  IF(PRESENT(def_val) .AND. iky == 0) val = def_val                  !--- Default value from arguments
+END FUNCTION fgetKeyByName_s1
 !==============================================================================================================================
 LOGICAL FUNCTION getKeyByName_s1(keyn, val, tname, ky) RESULT(lerr)
@@ -1151,102 +1179,108 @@
   CHARACTER(LEN=*),          INTENT(IN)  :: tname
   TYPE(keys_type), OPTIONAL, INTENT(IN)  :: ky(:)
-  INTEGER :: is
-  lerr = .FALSE.
+  CHARACTER(LEN=maxlen) :: tnam
+  INTEGER, ALLOCATABLE  :: is(:)
+  INTEGER :: i, itr
+  tnam = delPhase(strHead(tname,'_',.FALSE.))                        !--- Remove tag and phase
   IF(PRESENT(ky)) THEN
-    val = getKeyByName_prv(keyn, tname, ky);              IF(val /= '') RETURN !--- "ky" and "tnam"
-    val = getKeyByName_prv(keyn, delPhase(strHead(tname,'_')), ky)             !--- "ky" and "tnam" without phase
+    val = fgetKeyByName_s1(tname, keyn, ky, lerr=lerr)               !--- "ky" and "tname"
+    IF(val /= '' .OR. lerr)      RETURN
+    val = fgetKeyByName_s1(tnam,  keyn, ky, lerr=lerr)               !--- "ky" and "tnam"
   ELSE
     IF(.NOT.ALLOCATED(tracers))  RETURN
-    val = getKeyByName_prv(keyn, tname, tracers(:)%keys); IF(val /= '') RETURN !--- "tracers" and "tnam"
+    val = fgetKeyByName_s1(tname, keyn, tracers(:)%keys, lerr=lerr)  !--- "tracers" and "tname"
+    IF(val /= ''.AND..NOT.lerr)  RETURN
     IF(.NOT.ALLOCATED(isotopes)) RETURN
     IF(SIZE(isotopes) == 0)      RETURN
-    DO is = 1, SIZE(isotopes); IF(strIdx(isotopes(is)%keys(:)%name, delPhase(strHead(tname,'_'))) /= 0) EXIT; END DO
-    IF(is /= 0) val = getKeyByName_prv(keyn, tname, isotopes(is)%keys(:))      !--- "isotopes" and "tnam" without phase
+    !--- Search the "is" isotopes class index of the isotope named "tnam"
+    is = find([(ANY(isotopes(i)%keys(:)%name == tnam), i=1, SIZE(isotopes))])
+    IF(test(SIZE(is) == 0,lerr)) RETURN
+    val = fgetKeyByName_s1(tname, keyn, isotopes(is(1))%keys(:),lerr=lerr)!--- "isotopes" and "tnam"
   END IF
-
-CONTAINS
-
-FUNCTION getKeyByName_prv(keyn, tname, ky) RESULT(val)
-  CHARACTER(LEN=maxlen)         :: val
-  CHARACTER(LEN=*), INTENT(IN)  :: keyn
-  CHARACTER(LEN=*), INTENT(IN)  :: tname
-  TYPE(keys_type),  INTENT(IN)  :: ky(:)
-  INTEGER :: itr, iky
-  val = ''; iky = 0
-  itr = strIdx(ky(:)%name, tname);                 IF(itr==0) RETURN           !--- Get the index of the wanted tracer
-  IF(itr /= 0) iky = strIdx(ky(itr)%key(:), keyn); IF(iky==0) RETURN           !--- Wanted key    index
-  val = ky(itr)%val(iky)
-END FUNCTION getKeyByName_prv
-
 END FUNCTION getKeyByName_s1
 !==============================================================================================================================
-LOGICAL FUNCTION getKeyByName_sm(keyn, val, tnam, ky) RESULT(lerr)
+LOGICAL FUNCTION getKeyByName_sm(keyn, val, tname, ky) RESULT(lerr)
   CHARACTER(LEN=*),                   INTENT(IN)  :: keyn
-  CHARACTER(LEN=maxlen), ALLOCATABLE, INTENT(OUT) ::  val(:)
-  CHARACTER(LEN=*), TARGET, OPTIONAL, INTENT(IN)  :: tnam(:)
-  TYPE(keys_type),  TARGET, OPTIONAL, INTENT(IN)  ::   ky(:)
-  TYPE(keys_type),       POINTER :: k(:)
-  CHARACTER(LEN=maxlen), POINTER :: n(:)
-  INTEGER :: iq
-  k => tracers(:)%keys; IF(PRESENT(ky  )) k => ky
-  n => k(:)%name;       IF(PRESENT(tnam)) n => tnam
-  ALLOCATE(val(SIZE(n)))
-  lerr = ANY([(getKeyByName_s1(keyn, val(iq), n(iq), k), iq=1, SIZE(n))])
+  CHARACTER(LEN=maxlen), ALLOCATABLE, INTENT(OUT) ::   val(:)
+  CHARACTER(LEN=*), TARGET, OPTIONAL, INTENT(IN)  :: tname(:)
+  TYPE(keys_type),  TARGET, OPTIONAL, INTENT(IN)  ::    ky(:)
+  TYPE(keys_type),           POINTER :: k(:)
+  CHARACTER(LEN=maxlen), ALLOCATABLE :: n(:)
+  INTEGER :: iq, nq
+  IF(test(.NOT.(PRESENT(tname).OR.PRESENT(ky)), lerr)) RETURN
+  IF(PRESENT(ky   )) nq = SIZE(ky%name)
+  IF(PRESENT(tname)) nq = SIZE(  tname)
+  ALLOCATE(val(nq))
+  IF(PRESENT(tname)) THEN
+    IF(     PRESENT(ky)) lerr = ANY([(getKeyByName_s1(keyn, val(iq), tname(iq),   ky), iq=1, nq)])
+    IF(.NOT.PRESENT(ky)) lerr = ANY([(getKeyByName_s1(keyn, val(iq), tname(iq)      ), iq=1, nq)])
+  ELSE;                  lerr = ANY([(getKeyByName_s1(keyn, val(iq), ky(iq)%name, ky), iq=1, nq)])
+  END IF
 END FUNCTION getKeyByName_sm
 !==============================================================================================================================
-LOGICAL FUNCTION getKeyByName_i1(keyn, val, tnam, ky) RESULT(lerr)
+LOGICAL FUNCTION getKeyByName_i1(keyn, val, tname, ky) RESULT(lerr)
   CHARACTER(LEN=*),          INTENT(IN)  :: keyn
   INTEGER,                   INTENT(OUT) :: val
-  CHARACTER(LEN=*),          INTENT(IN)  :: tnam
+  CHARACTER(LEN=*),          INTENT(IN)  :: tname
   TYPE(keys_type), OPTIONAL, INTENT(IN)  :: ky(:)
   CHARACTER(LEN=maxlen) :: sval
   INTEGER :: ierr
-  IF(     PRESENT(ky)) lerr = getKeyByName_s1(keyn, sval, tnam, ky)
-  IF(.NOT.PRESENT(ky)) lerr = getKeyByName_s1(keyn, sval, tnam)
-  IF(test(fmsg('key "'//TRIM(keyn)//'" or tracer "'//TRIM(tnam)//'" is missing',        modname, lerr), lerr)) RETURN
+  IF(     PRESENT(ky)) lerr = getKeyByName_s1(keyn, sval, tname, ky)
+  IF(.NOT.PRESENT(ky)) lerr = getKeyByName_s1(keyn, sval, tname)
+  IF(test(fmsg('key "'//TRIM(keyn)//'" or tracer "'//TRIM(tname)//'" is missing',        modname, lerr), lerr)) RETURN
   READ(sval, *, IOSTAT=ierr) val
-  IF(test(fmsg('key "'//TRIM(keyn)//'" of tracer "'//TRIM(tnam)//'" is not an integer', modname, lerr), lerr)) RETURN
+  IF(test(fmsg('key "'//TRIM(keyn)//'" of tracer "'//TRIM(tname)//'" is not an integer', modname, lerr), lerr)) RETURN
 END FUNCTION getKeyByName_i1
 !==============================================================================================================================
-LOGICAL FUNCTION getKeyByName_im(keyn, val, tnam, ky) RESULT(lerr)
+LOGICAL FUNCTION getKeyByName_im(keyn, val, tname, ky) RESULT(lerr)
   CHARACTER(LEN=*),                   INTENT(IN)  :: keyn
-  INTEGER,               ALLOCATABLE, INTENT(OUT) ::  val(:)
-  CHARACTER(LEN=*), TARGET, OPTIONAL, INTENT(IN)  :: tnam(:)
-  TYPE(keys_type),  TARGET, OPTIONAL, INTENT(IN)  ::   ky(:)
-  TYPE(keys_type),       POINTER :: k(:)
-  CHARACTER(LEN=maxlen), POINTER :: n(:)
-  INTEGER :: iq
-  k => tracers(:)%keys; IF(PRESENT(ky  )) k => ky
-  n => k(:)%name;       IF(PRESENT(tnam)) n => tnam
-  ALLOCATE(val(SIZE(n)))
-  lerr = ANY([(getKeyByName_i1(keyn, val(iq), n(iq), k), iq=1, SIZE(n))])
+  INTEGER,               ALLOCATABLE, INTENT(OUT) ::   val(:)
+  CHARACTER(LEN=*), TARGET, OPTIONAL, INTENT(IN)  :: tname(:)
+  TYPE(keys_type),  TARGET, OPTIONAL, INTENT(IN)  ::    ky(:)
+  TYPE(keys_type),           POINTER :: k(:)
+  CHARACTER(LEN=maxlen), ALLOCATABLE :: n(:)
+  INTEGER :: iq, nq
+  IF(test(.NOT.(PRESENT(tname).OR.PRESENT(ky)), lerr)) RETURN
+  IF(PRESENT(ky   )) nq = SIZE(ky%name)
+  IF(PRESENT(tname)) nq = SIZE(  tname)
+  ALLOCATE(val(nq))
+  IF(PRESENT(tname)) THEN
+    IF(     PRESENT(ky)) lerr = ANY([(getKeyByName_i1(keyn, val(iq), tname(iq),   ky), iq=1, nq)])
+    IF(.NOT.PRESENT(ky)) lerr = ANY([(getKeyByName_i1(keyn, val(iq), tname(iq)      ), iq=1, nq)])
+  ELSE;                  lerr = ANY([(getKeyByName_i1(keyn, val(iq), ky(iq)%name, ky), iq=1, nq)])
+  END IF
 END FUNCTION getKeyByName_im
 !==============================================================================================================================
-LOGICAL FUNCTION getKeyByName_r1(keyn, val, tnam, ky) RESULT(lerr)
+LOGICAL FUNCTION getKeyByName_r1(keyn, val, tname, ky) RESULT(lerr)
   CHARACTER(LEN=*),          INTENT(IN)  :: keyn
   REAL,                      INTENT(OUT) :: val
-  CHARACTER(LEN=*),          INTENT(IN)  :: tnam
+  CHARACTER(LEN=*),          INTENT(IN)  :: tname
   TYPE(keys_type), OPTIONAL, INTENT(IN)  :: ky(:)
   CHARACTER(LEN=maxlen) :: sval
   INTEGER :: ierr
-  IF(     PRESENT(ky)) lerr = getKeyByName_s1(keyn, sval, tnam, ky)
-  IF(.NOT.PRESENT(ky)) lerr = getKeyByName_s1(keyn, sval, tnam)
-  IF(test(fmsg('key "'//TRIM(keyn)//'" or tracer "'//TRIM(tnam)//'" is missing',    modname, lerr), lerr)) RETURN
+  IF(     PRESENT(ky)) lerr = getKeyByName_s1(keyn, sval, tname, ky)
+  IF(.NOT.PRESENT(ky)) lerr = getKeyByName_s1(keyn, sval, tname)
+  IF(test(fmsg('key "'//TRIM(keyn)//'" or tracer "'//TRIM(tname)//'" is missing',    modname, lerr), lerr)) RETURN
   READ(sval, *, IOSTAT=ierr) val
-  IF(test(fmsg('key "'//TRIM(keyn)//'" of tracer "'//TRIM(tnam)//'" is not a real', modname, lerr), lerr)) RETURN
+  IF(test(fmsg('key "'//TRIM(keyn)//'" of tracer "'//TRIM(tname)//'" is not a real', modname, lerr), lerr)) RETURN
 END FUNCTION getKeyByName_r1
 !==============================================================================================================================
-LOGICAL FUNCTION getKeyByName_rm(keyn, val, tnam, ky) RESULT(lerr)
+LOGICAL FUNCTION getKeyByName_rm(keyn, val, tname, ky) RESULT(lerr)
   CHARACTER(LEN=*),                   INTENT(IN)  :: keyn
-  REAL,                  ALLOCATABLE, INTENT(OUT) ::  val(:)
-  CHARACTER(LEN=*), TARGET, OPTIONAL, INTENT(IN)  :: tnam(:)
-  TYPE(keys_type),  TARGET, OPTIONAL, INTENT(IN)  ::   ky(:)
-  TYPE(keys_type),       POINTER :: k(:)
-  CHARACTER(LEN=maxlen), POINTER :: n(:)
-  INTEGER :: iq
-  k => tracers(:)%keys; IF(PRESENT(ky  )) k => ky
-  n => k(:)%name;       IF(PRESENT(tnam)) n => tnam
-  ALLOCATE(val(SIZE(n)))
-  lerr = ANY([(getKeyByName_r1(keyn, val(iq), n(iq), k), iq=1, SIZE(n))])
+  REAL,                  ALLOCATABLE, INTENT(OUT) ::   val(:)
+  CHARACTER(LEN=*), TARGET, OPTIONAL, INTENT(IN)  :: tname(:)
+  TYPE(keys_type),  TARGET, OPTIONAL, INTENT(IN)  ::    ky(:)
+  TYPE(keys_type),           POINTER :: k(:)
+  CHARACTER(LEN=maxlen), ALLOCATABLE :: n(:)
+  INTEGER :: iq, nq
+  IF(test(.NOT.(PRESENT(tname).OR.PRESENT(ky)), lerr)) RETURN
+  IF(PRESENT(ky   )) nq = SIZE(ky%name)
+  IF(PRESENT(tname)) nq = SIZE(  tname)
+  ALLOCATE(val(nq))
+  IF(PRESENT(tname)) THEN
+    IF(     PRESENT(ky)) lerr = ANY([(getKeyByName_r1(keyn, val(iq), tname(iq),   ky), iq=1, nq)])
+    IF(.NOT.PRESENT(ky)) lerr = ANY([(getKeyByName_r1(keyn, val(iq), tname(iq)      ), iq=1, nq)])
+  ELSE;                  lerr = ANY([(getKeyByName_r1(keyn, val(iq), ky(iq)%name, ky), iq=1, nq)])
+  END IF
 END FUNCTION getKeyByName_rm
 !==============================================================================================================================
@@ -1329,13 +1363,81 @@
 !------------------------------------------------------------------------------------------------------------------------------
 
-CHARACTER(LEN=1) FUNCTION old2newPhase(op) RESULT(np)
-  CHARACTER(LEN=1), INTENT(IN) :: op
-  np = known_phases(INDEX(old_phases,op):INDEX(old_phases,op))
-END FUNCTION old2newPhase
-
+
+INTEGER FUNCTION getiPhase(tname, lPhase) RESULT(iphase)
+  CHARACTER(LEN=*),  INTENT(IN)  :: tname
+  LOGICAL, OPTIONAL, INTENT(OUT) :: lPhase
+  CHARACTER(LEN=maxlen) :: s1
+  INTEGER :: ip
+  IF(PRESENT(lPhase)) lPhase = .TRUE.
+
+  !--- Old tracer name descending on water: H2O[v][l][i][_<isotope>][_<tag>]
+  iphase = strIdx([('H2O'//old_phases(ip:ip), ip=1, nphases)], tname(1:MIN(4,LEN_TRIM(tname))))
+  IF(iphase /= 0) RETURN
+
+  !--- New tracer name: <name>[_<phase>][_<tag>]
+  iphase = INDEX(known_phases, TRIM(strHead(strTail(tname, phases_sep, .TRUE.), '_', .TRUE.)))
+  IF(iphase /= 0) RETURN
+
+  !---Default case: 1 (gaseous phase)
+  iphase = 1
+  IF(PRESENT(lPhase)) lPhase = .FALSE.
+END FUNCTION getiPhase  
+
+!------------------------------------------------------------------------------------------------------------------------------
 CHARACTER(LEN=1) FUNCTION new2oldPhase(np) RESULT(op)
   CHARACTER(LEN=1), INTENT(IN) :: np
   op = old_phases(INDEX(known_phases,np):INDEX(known_phases,np))
 END FUNCTION new2oldPhase
+!------------------------------------------------------------------------------------------------------------------------------
+
+!------------------------------------------------------------------------------------------------------------------------------
+CHARACTER(LEN=maxlen) FUNCTION old2newName(oldName, iPhase) RESULT(newName)
+  !--- Convert an old style name into a new one.
+  !    Only usable with old style "traceur.def" files, in which only water isotopes are allowed.
+  !    In these files, H2O descendants names are: H2O<phase>[_<isotope>][_<tag>], with:
+  !    phase = v, l or i ; isotope = eau, HDO, O18, O17 or HTO.
+  CHARACTER(LEN=*),  INTENT(IN)  :: oldName
+  INTEGER, OPTIONAL, INTENT(OUT) :: iPhase
+  CHARACTER(LEN=maxlen) :: oldIso(5) = ['eau',     'HDO',     'O18',     'O17',     'HTO'    ]
+  CHARACTER(LEN=maxlen) :: newIso(5) = ['H2[16]O', 'H[2]HO ', 'H2[18]O', 'H2[17]O', 'H[3]HO ']
+  CHARACTER(LEN=maxlen), ALLOCATABLE :: tmp(:)
+  INTEGER :: ix, ip, it, nt
+  LOGICAL :: lPhase, lerr
+  ip = getiPhase(oldName, lPhase)                                              !--- Get the phase ; lPhase==T: phase is needed
+  IF(PRESENT(iPhase)) iPhase = ip
+  IF(.NOT.lPhase) THEN; newName = oldName ; RETURN; END IF                     !--- Not a water descendant
+  newName = addPhase('H2O', ip)
+  lerr = strParse(oldName, '_', tmp, n=nt)
+  IF(nt == 1) RETURN                                                           !--- H2O with phase
+  ix = strIdx(oldIso, tmp(2))
+  newName = tmp(2); IF(ix /= 0) newName = newIso(ix)                           !--- Isotope name
+  IF(lPhase)  newName = addPhase(newName, ip)                                  !--- Phase is needed
+  IF(nt == 3) newName = TRIM(newName)//'_'//TRIM(tmp(3))                       !--- Tagging tracer
+END FUNCTION old2newName
+!------------------------------------------------------------------------------------------------------------------------------
+
+!------------------------------------------------------------------------------------------------------------------------------
+CHARACTER(LEN=maxlen) FUNCTION new2oldName(newName, iPhase) RESULT(oldName)
+  !--- Convert a new style name into an old one.
+  !    Only convertable names are water descendants names H2O_<phase>, <isotope>_<phase>, <isotope>_<phase>_<tag>, with:
+  !    phase = g, l or s ; isotope = H2[16]O, H[2]O, H2<[18]O, H2[17]O or H[3]O.
+  CHARACTER(LEN=*),  INTENT(IN)    :: newName
+  INTEGER, OPTIONAL, INTENT(OUT)   :: iPhase
+  CHARACTER(LEN=maxlen) :: oldIso(5) = ['eau',     'HDO',     'O18',     'O17',     'HTO'    ]
+  CHARACTER(LEN=maxlen) :: newIso(5) = ['H2[16]O', 'H[2]HO ', 'H2[18]O', 'H2[17]O', 'H[3]HO '], tag
+  INTEGER :: ix, ip, it, nt
+  LOGICAL :: lPhase, lH2O
+  lH2O = newName(1:MIN(3,LEN_TRIM(newName)))=='H2O'
+  ix = strIdx(newIso, strHead(strHead(newName,'_',.TRUE.),phases_sep,.TRUE.))  !--- Isotope index
+  IF(ix == 0 .AND. .NOT.lH2O) THEN; oldName=newName; RETURN; END IF            !--- Not a water descendant
+  ip = getiPhase(newName, lPhase)                                              !--- Get the phase ; lPhase==T: phase is needed
+  oldName = 'H2O'; IF(lPhase) oldName = addPhase('H2O', ip, '')                !--- H2O with phase
+  IF(ix == 0) RETURN
+  oldName = TRIM(oldName)//'_'//oldIso(ix)                                     !--- Isotope
+  tag = strTail(delPhase(newName), TRIM(newIso(ix)))
+  IF(tag /= delPhase(newName) .AND. tag /= '') oldName = TRIM(oldName)//tag    !--- Tagging tracer
+END FUNCTION new2oldName
+!------------------------------------------------------------------------------------------------------------------------------
+
 
 !==============================================================================================================================
