MODULE readTracFiles_mod

  USE ioipsl_getin_p_mod, ONLY : getin_p
  USE strings_mod, ONLY: msg, strIdx, dispTable, strHead,  strReduce,  strFind, strStack, removeComment, num2str, str2real, &
       reduceExpr, find, cat, maxlen, checkList, strParse, strReplace, strTail, strCount, maxTableWidth, str2int, str2bool

  IMPLICIT NONE

  PRIVATE

  PUBLIC :: maxlen                                              !--- PARAMETER FOR CASUAL STRING LENGTH
  PUBLIC :: trac_type, tracers, setGeneration, indexUpdate      !--- TRACERS  DESCRIPTION DATABASE + ASSOCIATED TOOLS
  PUBLIC :: testTracersFiles, readTracersFiles                  !--- TRACERS FILES READING ROUTINES
  PUBLIC :: getKeysDBase, setKeysDBase                          !--- TOOLS TO GET/SET THE DATABASE (tracers & isotopes)
  PUBLIC :: addTracer, delTracer                                !--- ADD/REMOVE A TRACER FROM
  PUBLIC :: addKey,    delKey,    getKey,    keys_type          !--- TOOLS TO SET/DEL/GET A KEY   FROM/TO "tracers" / "isotopes"
  PUBLIC :: addPhase,  delPhase,  getPhase,  getiPhase          !--- TOOLS TO SET/DEL/GET A PHASE FROM/TO A TRACER'S NAME
  PUBLIC ::   old_phases, phases_sep, nphases                   !--- VARIABLES RELATED TO THE PHASES
  PUBLIC :: known_phases, phases_names

  PUBLIC :: oldH2OIso, newH2OIso, old2newH2O, new2oldH2O        !--- H2O ISOTOPES BACKWARD COMPATIBILITY (OLD traceur.def)
  PUBLIC :: oldHNO3,   newHNO3                                  !--- HNO3 REPRO   BACKWARD COMPATIBILITY (OLD start.nc)

  PUBLIC :: tran0                                               !--- TRANSPORTING FLUID (USUALLY air)

  !=== FOR ISOTOPES: GENERAL
  PUBLIC :: isot_type, processIsotopes, isoSelect, isoFamilies  !--- ISOTOPES: TYPE, PROCESSING/SELECTION ROUTINES, FAMILIES NAMES
  PUBLIC :: ixIso, nbIso                                        !--- CURRENTLY SELECTED ISOTOPES FAMILY INDEX, NUMBER OF FAMILIES

  !=== FOR ISOTOPES: H2O FAMILY ONLY
  PUBLIC :: iH2O

  !=== FOR ISOTOPES: DEPENDING ON THE SELECTED ISOTOPES FAMILY
  PUBLIC :: isotope, isoKeys                                    !--- SELECTED ISOTOPES DATABASE + ASSOCIATED KEYS
  PUBLIC :: isoName, isoZone, isoPhas                           !--- ISOTOPES AND TAGGING ZONES NAMES AND PHASES
  PUBLIC :: niso,    nzone,   nphas,   ntiso                    !---  " " NUMBERS + ISOTOPES AND TAGGING TRACERS NUMBERS
  PUBLIC :: itZonIso                                            !--- Idx IN isoName(1:niso) = f(tagging idx, isotope idx)
  PUBLIC :: iqIsoPha                                            !--- Idx IN qx(1:nqtot)     = f(isotope idx,   phase idx)
  PUBLIC :: iqWIsoPha                                           !--- Idx IN qx(1:nqtot)     = f(isotope idx,   phase idx) but with normal water first
  PUBLIC :: isoCheck                                            !--- FLAG TO RUN ISOTOPES CHECKING ROUTINES

!------------------------------------------------------------------------------------------------------------------------------
  TYPE :: keys_type                                             !=== TYPE FOR A SET OF KEYS ASSOCIATED TO AN ELEMENT
    CHARACTER(LEN=maxlen)              :: name                  !--- Tracer name
    CHARACTER(LEN=maxlen), ALLOCATABLE :: key(:)                !--- Keys string list
    CHARACTER(LEN=maxlen), ALLOCATABLE :: val(:)                !--- Corresponding values string list
  END TYPE keys_type
!------------------------------------------------------------------------------------------------------------------------------
  TYPE :: trac_type                                             !=== TYPE FOR A SINGLE TRACER NAMED "name"
    CHARACTER(LEN=maxlen) :: name        = ''                   !--- Name of the tracer
    TYPE(keys_type)       :: keys                               !--- <key>=<val> pairs vector
    CHARACTER(LEN=maxlen) :: gen0Name    = ''                   !--- First generation ancestor name
    CHARACTER(LEN=maxlen) :: parent      = ''                   !--- Parent name
    CHARACTER(LEN=maxlen) :: longName    = ''                   !--- Long name (with advection scheme suffix)
    CHARACTER(LEN=maxlen) :: type        = 'tracer'             !--- Type  (so far: 'tracer' / 'tag')
    CHARACTER(LEN=maxlen) :: phase       = 'g'                  !--- Phase ('g'as / 'l'iquid / 's'olid)
    CHARACTER(LEN=maxlen) :: component   = ''                   !--- Coma-separated list of components (Ex: lmdz,inca)
    INTEGER               :: iGeneration = -1                   !--- Generation number (>=0)
    INTEGER               :: iqParent    = 0                    !--- Parent index
    INTEGER,  ALLOCATABLE :: iqDescen(:)                        !--- Descendants index (in growing generation order)
    INTEGER               :: nqDescen    = 0                    !--- Number of descendants (all generations)
    INTEGER               :: nqChildren  = 0                    !--- Number of children  (first generation)
    INTEGER               :: iadv        = 10                   !--- Advection scheme used
    LOGICAL               :: isInPhysics = .TRUE.               !--- "true" tracers: in tr_seri. COUNT(isInPhysics)=nqtottr
    INTEGER               :: iso_iGroup  = 0                    !--- Isotopes group index in isotopes(:)
    INTEGER               :: iso_iName   = 0                    !--- Isotope  name  index in isotopes(iso_iGroup)%trac(:)
    INTEGER               :: iso_iZone   = 0                    !--- Isotope  zone  index in isotopes(iso_iGroup)%zone(:)
    INTEGER               :: iso_iPhase  = 0                    !--- Isotope  phase index in isotopes(iso_iGroup)%phase
  END TYPE trac_type
!------------------------------------------------------------------------------------------------------------------------------
  TYPE :: isot_type                                             !=== TYPE FOR AN ISOTOPES FAMILY DESCENDING ON TRACER "parent"
    CHARACTER(LEN=maxlen)              :: parent                !--- Isotopes family name (parent tracer name ; ex: H2O)
    TYPE(keys_type),       ALLOCATABLE :: keys(:)               !--- Isotopes keys/values pairs list     (length: niso)
    LOGICAL                            :: check=.FALSE.         !--- Flag for checking routines triggering
    CHARACTER(LEN=maxlen), ALLOCATABLE :: trac(:)               !--- Isotopes + tagging tracers list     (length: ntiso)
    CHARACTER(LEN=maxlen), ALLOCATABLE :: zone(:)               !--- Geographic tagging zones names list (length: nzone)
    CHARACTER(LEN=maxlen)              :: phase = 'g'           !--- Phases list: [g][l][s]              (length: nphas)
    INTEGER                            :: niso  = 0             !--- Number of isotopes, excluding tagging tracers
    INTEGER                            :: ntiso = 0             !--- Number of isotopes, including tagging tracers
    INTEGER                            :: nzone = 0             !--- Number of geographic tagging zones
    INTEGER                            :: nphas = 0             !--- Number of phases
    INTEGER,               ALLOCATABLE :: iqIsoPha(:,:)         !--- Idx in "tracers(1:nqtot)" = f(name(1:ntiso),phas)
                                                                !---        (former name: "iqiso"
    INTEGER,               ALLOCATABLE :: iqWIsoPha(:,:)        !--- Idx in "tracers(1:nqtot)" = f([H2O,name(1:ntiso)],phas)
    INTEGER,               ALLOCATABLE :: itZonIso(:,:)         !--- Idx in "trac(1:ntiso)" = f(zone, name(1:niso))
  END TYPE isot_type                                            !---        (former name: "index_trac")
!------------------------------------------------------------------------------------------------------------------------------
  TYPE :: dataBase_type                                         !=== TYPE FOR TRACERS SECTION
    CHARACTER(LEN=maxlen) :: name                               !--- Section name
    TYPE(trac_type), ALLOCATABLE :: trac(:)                     !--- Tracers descriptors
  END TYPE dataBase_type
!------------------------------------------------------------------------------------------------------------------------------
  INTERFACE getKey
    MODULE PROCEDURE &
       getKeyByIndex_s111, getKeyByIndex_sm11, getKeyByIndex_s1m1, getKeyByIndex_smm1, getKeyByIndex_s1mm, getKeyByIndex_smmm, &
       getKeyByIndex_i111, getKeyByIndex_im11, getKeyByIndex_i1m1, getKeyByIndex_imm1, getKeyByIndex_i1mm, getKeyByIndex_immm, &
       getKeyByIndex_r111, getKeyByIndex_rm11, getKeyByIndex_r1m1, getKeyByIndex_rmm1, getKeyByIndex_r1mm, getKeyByIndex_rmmm, &
       getKeyByIndex_l111, getKeyByIndex_lm11, getKeyByIndex_l1m1, getKeyByIndex_lmm1, getKeyByIndex_l1mm, getKeyByIndex_lmmm, &
        getKeyByName_s111,  getKeyByName_sm11,  getKeyByName_s1m1,  getKeyByName_smm1,  getKeyByName_s1mm,  getKeyByName_smmm, &
        getKeyByName_i111,  getKeyByName_im11,  getKeyByName_i1m1,  getKeyByName_imm1,  getKeyByName_i1mm,  getKeyByName_immm, &
        getKeyByName_r111,  getKeyByName_rm11,  getKeyByName_r1m1,  getKeyByName_rmm1,  getKeyByName_r1mm,  getKeyByName_rmmm, &
        getKeyByName_l111,  getKeyByName_lm11,  getKeyByName_l1m1,  getKeyByName_lmm1,  getKeyByName_l1mm,  getKeyByName_lmmm
  END INTERFACE getKey
!------------------------------------------------------------------------------------------------------------------------------
  INTERFACE addKey
    MODULE PROCEDURE addKey_s11, addKey_s1m, addKey_smm, addKey_i11, addKey_i1m, addKey_imm, &
                     addKey_r11, addKey_r1m, addKey_rmm, addKey_l11, addKey_l1m, addKey_lmm
  END INTERFACE addKey
!------------------------------------------------------------------------------------------------------------------------------
  INTERFACE     isoSelect; MODULE PROCEDURE  isoSelectByIndex,  isoSelectByName; END INTERFACE isoSelect
  INTERFACE    old2newH2O; MODULE PROCEDURE  old2newH2O_1,  old2newH2O_m;        END INTERFACE old2newH2O
  INTERFACE    new2oldH2O; MODULE PROCEDURE  new2oldH2O_1,  new2oldH2O_m;        END INTERFACE new2oldH2O
  INTERFACE     addTracer; MODULE PROCEDURE   addTracer_1, addTracer_1def;       END INTERFACE addTracer
  INTERFACE     delTracer; MODULE PROCEDURE   delTracer_1, delTracer_1def;       END INTERFACE delTracer
  INTERFACE      addPhase; MODULE PROCEDURE   addPhase_s1,  addPhase_sm,  addPhase_i1,  addPhase_im; END INTERFACE addPhase
  INTERFACE tracersSubset; MODULE PROCEDURE trSubset_Indx,     trSubset_Name,     trSubset_gen0Name; END INTERFACE tracersSubset
!------------------------------------------------------------------------------------------------------------------------------

  !=== MAIN DATABASE: files sections descriptors
  TYPE(dataBase_type), SAVE, ALLOCATABLE, TARGET :: dBase(:)

  !--- SOME PARAMETERS THAT ARE NOT LIKELY TO CHANGE OFTEN
  CHARACTER(LEN=maxlen), SAVE      :: tran0        = 'air'      !--- Default transporting fluid
  CHARACTER(LEN=maxlen), PARAMETER :: old_phases   = 'vlibfc'    !--- Old phases for water (no separator)
  CHARACTER(LEN=maxlen), PARAMETER :: known_phases = 'glsbfc'    !--- 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    ','blownSnow', 'fracCloud', 'cldVapRat']
  CHARACTER(LEN=1),      SAVE :: phases_sep  =  '_'             !--- Phase separator
  CHARACTER(LEN=maxlen), SAVE :: isoFile = 'isotopes_params.def'!--- Name of the isotopes parameters file

  !--- CORRESPONDANCE BETWEEN OLD AND NEW WATER NAMES
  CHARACTER(LEN=maxlen), SAVE :: oldH2OIso(5) = ['eau',   'HDO',   'O18',   'O17',   'HTO'  ]
  CHARACTER(LEN=maxlen), SAVE :: newH2OIso(5) = ['H216O', 'HDO  ', 'H218O', 'H217O', 'HTO  ']

  !--- CORRESPONDANCE BETWEEN OLD AND NEW HNO3 RELATED SPECIES NAMES (FOR REPROBUS)
  CHARACTER(LEN=maxlen), SAVE ::   oldHNO3(2) = ['HNO3_g ', 'HNO3   ']
  CHARACTER(LEN=maxlen), SAVE ::   newHNO3(2) = ['HNO3   ', 'HNO3tot']

  !=== TRACERS AND ISOTOPES DESCRIPTORS, USED BY getKey
  TYPE(trac_type), ALLOCATABLE, TARGET, SAVE ::  tracers(:)
  TYPE(isot_type), ALLOCATABLE, TARGET, SAVE :: isotopes(:)

  !--- NAMES OF THE ISOTOPES FAMILIES
  CHARACTER(LEN=maxlen),   SAVE, ALLOCATABLE :: isoFamilies(:)

  !=== ALIASES OF VARIABLES FROM SELECTED ISOTOPES FAMILY EMBEDDED IN "isotope" (isotopes(ixIso))
  TYPE(isot_type),         SAVE, POINTER :: isotope             !--- CURRENTLY SELECTED ISOTOPES FAMILY DESCRIPTOR
  INTEGER,                 SAVE          :: ixIso, iH2O=-1      !--- INDEX OF THE SELECTED ISOTOPES FAMILY, H2O FAMILY INDEX
  INTEGER,                 SAVE          :: nbIso               !--- NUMBER OF ISOTOPES FAMILIES
  LOGICAL,                 SAVE          :: isoCheck            !--- FLAG TO TRIGGER THE CHECKING ROUTINES
  TYPE(keys_type),         SAVE, POINTER :: isoKeys(:)          !--- ONE SET OF KEYS FOR EACH ISOTOPE (LISTED IN isoName)
  CHARACTER(LEN=maxlen),   SAVE, POINTER :: isoName(:),   &     !--- ISOTOPES NAMES FOR THE CURRENTLY SELECTED FAMILY
                                            isoZone(:),   &     !--- TAGGING ZONES  FOR THE CURRENTLY SELECTED FAMILY
                                            isoPhas             !--- USED PHASES    FOR THE CURRENTLY SELECTED FAMILY
  INTEGER,                 SAVE          ::  niso, nzone, &     !--- NUMBER OF ISOTOPES, TAGGING ZONES AND PHASES
                                            nphas, ntiso        !--- NUMBER OF PHASES AND ISOTOPES + ISOTOPIC TAGGING TRACERS
  INTEGER,                 SAVE, POINTER ::itZonIso(:,:), &     !--- INDEX IN "isoTrac" AS f(tagging zone idx,  isotope idx)
                                           iqIsoPha(:,:), &     !--- INDEX IN "qx" AS f(isotopic tracer idx, phase idx)
                                           iqWIsoPha(:,:)       !--- INDEX IN "qx" AS f(H2O + isotopic tracer idx, phase idx)

  !=== PARAMETERS FOR DEFAULT BEHAVIOUR
  LOGICAL, PARAMETER :: lTracsMerge = .FALSE.                   !--- Merge/stack tracers lists
  LOGICAL, PARAMETER :: lSortByGen  = .TRUE.                    !--- Sort by growing generation

  CHARACTER(LEN=maxlen) :: modname

CONTAINS

!==============================================================================================================================
!==============================================================================================================================
!=== READ ONE OR SEVERAL TRACER FILES AND FILL A "tr" TRACERS DESCRIPTOR DERIVED TYPE.
!=== THE RETURNED VALUE fType DEPENDS ON WHAT IS FOUND:
!===  0: NO ADEQUATE FILE FOUND ; DEFAULT VALUES MUST BE USED
!===  1: AN "OLD STYLE" TRACERS FILE "traceur.def":
!===    First line: <nb tracers>     Other lines: <hadv> <vadv> <tracer name> [<parent name>]
!===  2: A  "NEW STYLE" TRACERS FILE  "tracer.def" WITH SEVERAL SECTIONS.
!===  3: SEVERAL  "  "  TRACERS FILES "tracer_<component>.def" WITH A SINGLE SECTION IN EACH.
!=== REMARKS:
!===  * EACH SECTION BEGINS WITH A "&<section name> LINE
!===  * DEFAULT VALUES FOR ALL THE SECTIONS OF THE FILE ARE DEFINED IN THE SPECIAL SECTION "&default"
!===  * EACH SECTION LINE HAS THE STRUCTURE:  <name(s)>  <key1>=<value1> <key2>=<value2> ...
!===  * SO FAR, THE DEFINED KEYS ARE: parent, phases, hadv, vadv, type
!===  * <name> AND <parent> CAN BE LISTS OF COMA-SEPARATED TRACERS ; THE ROUTINE EXPAND THESE FACTORIZATIONS.
!=== FUNCTION RETURN VALUE "lerr" IS FALSE IN CASE SOMETHING WENT WRONG.
!=== ABOUT THE KEYS:
!     * The "keys" component (of type keys_type) is in principle enough to store everything we could need.
!     But some variables are stored as direct-access keys to make the code more readable and because they are used often.
!     * Most of the direct-access keys are set in this module, but some are not (longName, iadv and isInPhysicsfor now).
!     * Some of the direct-access keys must be updated (using the routine "setDirectKeys") is a subset of "tracers(:)"
!     is extracted: the indexes are no longer valid for a subset (examples: iqParent, iqDescen).
!     * If you need to convert a %key/%val pair into a direct-access key, add the corresponding line in "setDirectKeys".
!==============================================================================================================================
LOGICAL FUNCTION readTracersFiles(type_trac, tracs, lRepr) RESULT(lerr)
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=*),                               INTENT(IN)  :: type_trac     !--- List of components used
  TYPE(trac_type), ALLOCATABLE, TARGET, OPTIONAL, INTENT(OUT) :: tracs(:)      !--- Tracers descriptor for external storage
  LOGICAL,                              OPTIONAL, INTENT(IN)  :: lRepr         !--- Activate the HNO3 exceptions for REPROBUS
  CHARACTER(LEN=maxlen),  ALLOCATABLE :: s(:), sections(:), trac_files(:)
  CHARACTER(LEN=maxlen) :: str, fname, tname, pname, cname
  INTEGER               :: nsec, ierr, it, ntrac, ns, ip, ix, fType
  INTEGER, ALLOCATABLE  :: iGen(:)
  LOGICAL :: lRep
  TYPE(keys_type), POINTER :: k
!------------------------------------------------------------------------------------------------------------------------------
  lerr = .FALSE.
  modname = 'readTracersFiles'
  IF(.NOT.ALLOCATED(dBase)) ALLOCATE(dBase(0))
  lRep=.FALSE.; IF(PRESENT(lRepr)) lRep = lRepr

  !--- Required sections + corresponding files names (new style single section case) for tests
  lerr = testTracersFiles(modname, type_trac, fType, .FALSE., trac_files, sections); IF(lerr) RETURN
  nsec = SIZE(sections)

  !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  SELECT CASE(fType)                         !--- Set name, component, parent, phase, iGeneration, gen0Name, type
  !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
    CASE(1)                                                          !=== OLD FORMAT "traceur.def"
    !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
      !--- OPEN THE "traceur.def" FILE
      OPEN(90, FILE="traceur.def", FORM='formatted', STATUS='old', POSITION='REWIND', IOSTAT=ierr)

      !--- GET THE TRACERS NUMBER
      READ(90,'(i3)',IOSTAT=ierr)ntrac                               !--- Number of lines/tracers
      lerr = ierr/=0; CALL msg('Invalid format for "'//TRIM(fname)//'"', modname, lerr); IF(lerr) RETURN

      !--- READ THE REMAINING LINES: <hadv> <vadv> <tracer> [<transporting fluid>]
      ALLOCATE(tracers(ntrac))
      DO it = 1, ntrac                                               !=== READ RAW DATA: loop on the line/tracer number
        READ(90,'(a)',IOSTAT=ierr) str
        lerr = ierr>0; CALL msg('Invalid format for "' //TRIM(fname)//'"', modname, lerr); IF(lerr) RETURN
        lerr = ierr<0; CALL msg('Not enough lines in "'//TRIM(fname)//'"', modname, lerr); IF(lerr) RETURN
        lerr = strParse(str, ' ', s, ns)
        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)
        k => tracers(it)%keys

        !=== NAME OF THE TRACER
        tname = old2newH2O(s(3), ip)
        ix = strIdx(oldHNO3, s(3))
        IF(ix /= 0 .AND. lRep) tname = newHNO3(ix)                   !--- Exception for HNO3 (REPROBUS ONLY)
        tracers(it)%name = tname                                     !--- Set the name of the tracer
        CALL addKey('name', tname, k)                                !--- Set the name of the tracer
        tracers(it)%keys%name = tname                                !--- Copy tracers names in keys components

        !=== NAME OF THE COMPONENT
        cname = type_trac                                            !--- Name of the model component
        IF(ANY([(addPhase('H2O', ip), ip = 1, nphases)] == tname)) cname = 'lmdz'
        tracers(it)%component = cname                                !--- Set component
        CALL addKey('component', cname, k)                           !--- Set the name of the model component

        !=== NAME OF THE PARENT
        pname = tran0                                                !--- Default name: default transporting fluid (air)
        IF(ns == 4) THEN
          pname = old2newH2O(s(4))
          ix = strIdx(oldHNO3, s(4))
          IF(ix /= 0 .AND. lRep) pname = newHNO3(ix)                 !--- Exception for HNO3 (REPROBUS ONLY)
        END IF
        tracers(it)%parent = pname                                   !--- Set the parent name
        CALL addKey('parent', pname, k)

        !=== PHASE AND ADVECTION SCHEMES NUMBERS
        tracers(it)%phase = known_phases(ip:ip)                      !--- Set the phase of the tracer (default: "g"azeous)
        CALL addKey('phase', known_phases(ip:ip), k)                 !--- Set the phase of the tracer (default: "g"azeous)
        CALL addKey('hadv', s(1),  k)                                !--- Set the horizontal advection schemes number
        CALL addKey('vadv', s(2),  k)                                !--- Set the vertical   advection schemes number
      END DO
      CLOSE(90)
      lerr = setGeneration(tracers); IF(lerr) RETURN                 !--- Set iGeneration and gen0Name
      lerr = getKey('iGeneration', iGen, tracers(:)%keys)            !--- Generation number
      WHERE(iGen == 2) tracers(:)%type = 'tag'                       !--- Set type:      'tracer' or 'tag'
      DO it = 1, ntrac
        CALL addKey('type', tracers(it)%type, tracers(it)%keys)      !--- Set the type of tracer
      END DO
      lerr = checkTracers(tracers, fname, fname); IF(lerr) RETURN    !--- Detect orphans and check phases
      lerr = checkUnique (tracers, fname, fname); IF(lerr) RETURN    !--- Detect repeated tracers
      CALL sortTracers   (tracers)                                   !--- Sort the tracers
    !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
    CASE(2); lerr=feedDBase(["tracer.def"], [type_trac], modname); IF(lerr) RETURN  !=== SINGLE   FILE, MULTIPLE SECTIONS
    !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
    CASE(3); lerr=feedDBase(  trac_files  ,  sections,   modname); IF(lerr) RETURN  !=== MULTIPLE FILES, SINGLE  SECTION
  !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  END SELECT
  !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  IF(ALL([2,3] /= fType)) RETURN
  IF(nsec == 1) tracers = dBase(1)%trac
  IF(nsec /= 1) THEN
    CALL msg('Multiple sections are MERGED',    modname,      lTracsMerge)
    CALL msg('Multiple sections are CUMULATED', modname, .NOT.lTracsMerge)
    IF(     lTracsMerge) lerr = cumulTracers(dBase, tracers)
    IF(.NOT.lTracsMerge) lerr = cumulTracers(dBase, tracers)
    IF(lerr) RETURN
  END IF
  lerr = indexUpdate(tracers); IF(lerr) RETURN                       !--- Set iqParent, iqDescen, nqDescen, nqChildren
  IF(PRESENT(tracs)) CALL MOVE_ALLOC(FROM=tracers, TO=tracs)
END FUNCTION readTracersFiles
!==============================================================================================================================


!==============================================================================================================================
LOGICAL FUNCTION testTracersFiles(modname, type_trac, fType, lDisp, tracf, sects) RESULT(lerr)
  CHARACTER(LEN=*),                             INTENT(IN)  :: modname, type_trac
  INTEGER,                                      INTENT(OUT) :: fType
  LOGICAL,                            OPTIONAL, INTENT(IN)  :: lDisp
  CHARACTER(LEN=maxlen), ALLOCATABLE, OPTIONAL, INTENT(OUT) :: tracf(:), sects(:)
  CHARACTER(LEN=maxlen), ALLOCATABLE :: trac_files(:), sections(:)
  LOGICAL, ALLOCATABLE :: ll(:)
  LOGICAL :: lD, lFound
  INTEGER :: is, nsec
  lD = .FALSE.; IF(PRESENT(lDisp)) lD = lDisp
  lerr = .FALSE.

  !--- PARSE "type_trac" LIST AND DETERMINE THE TRACERS FILES NAMES (FOR CASE 3: MULTIPLE FILES, SINGLE SECTION PER FILE)
  !--- If type_trac is a scalar (case 1), "sections" and "trac_files" are not usable, but are meaningless for case 1 anyway.
  lerr = strParse(type_trac, '|', sections, n=nsec); IF(lerr) RETURN !--- Parse "type_trac" list
  IF(PRESENT(sects)) sects = sections
  ALLOCATE(trac_files(nsec), ll(nsec))
  DO is=1, nsec
     trac_files(is) = 'tracer_'//TRIM(sections(is))//'.def'
     INQUIRE(FILE=TRIM(trac_files(is)), EXIST=ll(is))
  END DO
  IF(PRESENT(tracf)) tracf = trac_files
  fType = 0
  INQUIRE(FILE='traceur.def', EXIST=lFound); IF(lFound)  fType = 1   !--- OLD STYLE FILE
  INQUIRE(FILE='tracer.def',  EXIST=lFound); IF(lFound)  fType = 2   !--- NEW STYLE ; SINGLE  FILE, SEVERAL SECTIONS
                                             IF(ALL(ll)) fType = 3   !--- NEW STYLE ; SEVERAL FILES, SINGLE SECTION USED
  IF(.NOT.lD) RETURN                                                 !--- NO CHECKING/DISPLAY NEEDED: JUST GET type_trac,fType
  IF(ANY(ll) .AND. fType/=3) THEN                                    !--- MISSING FILES
    lerr = checkList(trac_files, .NOT.ll, 'Failed reading tracers description', 'files', 'missing'); IF(lerr) RETURN
  END IF

  !--- TELLS WHAT WAS IS ABOUT TO BE USED
  CALL msg('Trying to read old-style tracers description file "traceur.def"',                      modname, fType==1)
  CALL msg('Trying to read the new style multi-sections tracers description file "tracer.def"',    modname, fType==2)
  CALL msg('Trying to read the new style single section tracers description files "tracer_*.def"', modname, fType==3)
END FUNCTION testTracersFiles
!==============================================================================================================================

!==============================================================================================================================
LOGICAL FUNCTION feedDBase(fnames, snames, modname) RESULT(lerr)
! Purpose: Read the sections "snames(is)" (pipe-separated list) from each "fnames(is)"
!   file and create the corresponding tracers set descriptors in the database "dBase":
! * dBase(id)%name                : section name
! * dBase(id)%trac(:)%name        : tracers names
! * dBase(id)%trac(it)%keys%key(:): names  of keys associated to tracer dBase(id)%trac(it)%name
! * dBase(id)%trac(it)%keys%val(:): values of keys associated to tracer dBase(id)%trac(it)%name
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=*), INTENT(IN)  :: fnames(:)                         !--- Files names
  CHARACTER(LEN=*), INTENT(IN)  :: snames(:)                         !--- Pipe-deparated list of sections (one list 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
  CHARACTER(LEN=maxlen) :: fnm, snm
  INTEGER               :: idb, i
  LOGICAL :: ll
!------------------------------------------------------------------------------------------------------------------------------
  !=== READ THE REQUIRED SECTIONS
  ll = strCount(snames, '|', ndb)                                    !--- Number of sections for each file
  ALLOCATE(ixf(SUM(ndb)))
  DO i=1, SIZE(fnames)                                               !--- Set name, keys
    lerr = readSections(fnames(i), snames(i), 'default'); IF(lerr) RETURN
    ixf(1+SUM(ndb(1:i-1)):SUM(ndb(1:i))) = i                         !--- File index for each section of the expanded list
  END DO
  !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  DO idb=1,SIZE(dBase)                                               !--- LOOP ON THE LOADED SECTIONS
  !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
    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))])
    lerr = expandSection(dBase(idb)%trac, snm, fnm); IF(lerr) RETURN !--- EXPAND NAMES ;  SET parent, type, component
    lerr = setGeneration(dBase(idb)%trac);           IF(lerr) RETURN !---                 SET iGeneration,  genOName
    lerr = checkTracers (dBase(idb)%trac, snm, fnm); IF(lerr) RETURN !--- CHECK ORPHANS AND PHASES
    lerr = checkUnique  (dBase(idb)%trac, snm, fnm); IF(lerr) RETURN !--- CHECK TRACERS UNIQUENESS
    lerr = expandPhases (dBase(idb)%trac);           IF(lerr) RETURN !--- 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
  !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
END FUNCTION feedDBase
!------------------------------------------------------------------------------------------------------------------------------

!------------------------------------------------------------------------------------------------------------------------------
LOGICAL FUNCTION readSections(fnam,snam,defName) RESULT(lerr)
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=*),           INTENT(IN) :: fnam                     !--- File name
  CHARACTER(LEN=*),           INTENT(IN) :: snam                     !--- Pipe-separated sections list
  CHARACTER(LEN=*), OPTIONAL, INTENT(IN) :: defName                  !--- Special section (default values) name
!------------------------------------------------------------------------------------------------------------------------------
  TYPE(dataBase_type),   ALLOCATABLE :: tdb(:)
  CHARACTER(LEN=maxlen), ALLOCATABLE :: sec(:)
  INTEGER,               ALLOCATABLE ::  ix(:)
  INTEGER :: n0, idb, ndb
  LOGICAL :: ll
!------------------------------------------------------------------------------------------------------------------------------
  n0 = SIZE(dBase) + 1                                               !--- Index for next entry in the database
  CALL readSections_all()                                            !--- Read all the sections of file "fnam"
  ndb= SIZE(dBase)                                                   !--- Current number of sections in the database
  IF(PRESENT(defName)) THEN                                          !--- Add default values to all the tracers
    DO idb=n0,ndb; CALL addDefault(dBase(idb)%trac, defName); END DO !--- and remove the virtual tracer "defName"
  END IF
  ll = strParse(snam, '|', keys = sec)                               !--- Requested sections names
  ix = strIdx(dBase(:)%name, sec(:))                                 !--- Indexes of requested sections in database
  lerr = checkList(sec, ix == 0, 'In file "'//TRIM(fnam)//'"','section(s)', 'missing'); IF(lerr) RETURN
  tdb = dBase(:); dBase = [tdb(1:n0-1),tdb(PACK(ix, MASK=ix/=0))]    !--- Keep requested sections only

CONTAINS

!------------------------------------------------------------------------------------------------------------------------------
SUBROUTINE readSections_all()
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen), ALLOCATABLE ::  s(:), v(:)
  TYPE(trac_type),       ALLOCATABLE :: tt(:)
  TYPE(trac_type)       :: tmp
  CHARACTER(LEN=1024)   :: str, str2
  CHARACTER(LEN=maxlen) :: secn
  INTEGER               :: ierr, n
!------------------------------------------------------------------------------------------------------------------------------
  IF(.NOT.ALLOCATED(dBase)) ALLOCATE(dBase(0))
  OPEN(90, FILE=fnam, FORM='formatted', POSITION='REWIND', STATUS='old')
  DO; str=''
    DO
      READ(90,'(a)', IOSTAT=ierr)str2                                !--- Read a full line
      str=TRIM(str)//' '//TRIM(str2)                                 !--- Append "str" with the current line
      n=LEN_TRIM(str); IF(n == 0) EXIT                               !--- Empty line (probably end of file)
      IF(IACHAR(str(n:n))  /= 92) EXIT                               !--- No "\" continuing line symbol found => end of line
      str = str(1:n-1)                                               !--- Remove the "\" continuing line symbol
    END DO
    str = ADJUSTL(str)                                               !--- Remove the front space
    IF(ierr    /= 0 ) EXIT                                           !--- Finished: error or end of file
    IF(str(1:1)=='#') CYCLE                                          !--- Skip comments lines
    CALL removeComment(str)                                          !--- Skip comments at the end of a line
    IF(LEN_TRIM(str) == 0) CYCLE                                     !--- Empty line (probably end of file)
    IF(str     == '') CYCLE                                          !--- Skip empty line (probably at the end of the file)
    IF(str(1:1)=='&') THEN                                           !=== SECTION HEADER LINE
      ndb  = SIZE(dBase)                                             !--- Number of sections so far
      secn = str(2:LEN_TRIM(str))//' '                               !--- Current section name
      IF(ANY(dBase(:)%name == secn)) CYCLE                           !--- Already known section
      IF(secn(1:7) == 'version')     CYCLE                           !--- Skip the "version" special section
      ndb = ndb + 1                                                  !--- Extend database
      ALLOCATE(tdb(ndb))
      tdb(1:ndb-1)  = dBase
      tdb(ndb)%name = secn
      ALLOCATE(tdb(ndb)%trac(0))
      CALL MOVE_ALLOC(FROM=tdb, TO=dBase)
    ELSE                                                             !=== TRACER LINE
      ll = strParse(str,' ', s, n, v)                                !--- Parse <key>=<val> pairs
      tt = dBase(ndb)%trac(:)
      v(1) = s(1); s(1) = 'name'                                     !--- Convert "name" into a regular key
      tmp%name = v(1); tmp%keys = keys_type(v(1), s(:), v(:))        !--- Set %name and %keys
      dBase(ndb)%trac = [tt(:), tmp]
      DEALLOCATE(tt, tmp%keys%key, tmp%keys%val)
    END IF
  END DO
  CLOSE(90)

END SUBROUTINE readSections_all
!------------------------------------------------------------------------------------------------------------------------------

END FUNCTION readSections
!==============================================================================================================================


!==============================================================================================================================
SUBROUTINE addDefault(t, defName)
!------------------------------------------------------------------------------------------------------------------------------
! Purpose: Add the keys from virtual tracer named "defName" (if any) and remove this virtual tracer.
!------------------------------------------------------------------------------------------------------------------------------
  TYPE(trac_type), ALLOCATABLE, TARGET, INTENT(INOUT) :: t(:)
  CHARACTER(LEN=*),                     INTENT(IN)    :: defName
  INTEGER :: jd, it, k
  TYPE(keys_type), POINTER :: ky
  TYPE(trac_type), ALLOCATABLE :: tt(:)
  jd = strIdx(t(:)%name, defName)
  IF(jd == 0) RETURN
  ky => t(jd)%keys
  DO k = 1, SIZE(ky%key)                                             !--- Loop on the keys of the tracer named "defName"
!   CALL addKey(ky%key(k), ky%val(k), t(:)%keys, .FALSE.)            !--- Add key to all the tracers (no overwriting)
    DO it = 1, SIZE(t); CALL addKey(ky%key(k), ky%val(k), t(it)%keys, .FALSE.); END DO
  END DO
  tt = [t(1:jd-1),t(jd+1:SIZE(t))]; CALL MOVE_ALLOC(FROM=tt, TO=t)   !--- Remove the virtual tracer named "defName"
END SUBROUTINE addDefault
!==============================================================================================================================

!==============================================================================================================================
SUBROUTINE subDefault(t, defName, lSubLocal)
!------------------------------------------------------------------------------------------------------------------------------
! Purpose: Substitute the keys from virtual tracer named "defName" (if any) and remove this virtual tracer.
!          Substitute the keys locally (for the current tracer) if the flag "lSubLocal" is .TRUE.
!------------------------------------------------------------------------------------------------------------------------------
  TYPE(trac_type), ALLOCATABLE, TARGET, INTENT(INOUT) :: t(:)
  CHARACTER(LEN=*),                     INTENT(IN)    :: defName
  LOGICAL,                              INTENT(IN)    :: lSubLocal
  INTEGER :: i0, it, ik
  TYPE(keys_type), POINTER     :: k0, ky
  TYPE(trac_type), ALLOCATABLE :: tt(:)
  i0 = strIdx(t(:)%name, defName)
  IF(i0 == 0) RETURN
  k0 => t(i0)%keys 
  DO it = 1, SIZE(t); IF(it == i0) CYCLE                             !--- Loop on the tracers
    ky => t(it)%keys

    !--- Substitute in the values of <key>=<val> pairs the keys defined in the virtual tracer "defName"
    DO ik = 1, SIZE(k0%key); CALL strReplace(ky%val, k0%key(ik), k0%val(ik), .TRUE.); END DO

    IF(.NOT.lSubLocal) CYCLE
    !--- Substitute in the values of <key>=<val> pairs the keys defined locally (in the current tracer)
    DO ik = 1, SIZE(ky%key); CALL strReplace(ky%val, ky%key(ik), ky%val(ik), .TRUE.); END DO
  END DO
  tt = [t(1:i0-1),t(i0+1:SIZE(t))]; CALL MOVE_ALLOC(FROM=tt, TO=t)   !--- Remove the virtual tracer named "defName"

END SUBROUTINE subDefault
!==============================================================================================================================


!==============================================================================================================================
LOGICAL FUNCTION expandSection(tr, sname, fname) RESULT(lerr)
!------------------------------------------------------------------------------------------------------------------------------
! Purpose: Expand tracers and parents lists in the tracers descriptor "tra".
! Note:  * The following keys are expanded, so are accessible explicitely using "%" operator: "parent" "type".
!        * Default values are provided for these keys because they are necessary.
!------------------------------------------------------------------------------------------------------------------------------
  TYPE(trac_type), ALLOCATABLE, INTENT(INOUT) :: tr(:)                 !--- Tracer derived type vector
  CHARACTER(LEN=*),             INTENT(IN)    :: sname                 !--- Current section name
  CHARACTER(LEN=*), OPTIONAL,   INTENT(IN)    :: fname                 !--- Tracers description file name
  TYPE(trac_type),       ALLOCATABLE :: ttr(:)
  CHARACTER(LEN=maxlen), ALLOCATABLE :: ta(:), pa(:), tname(:), parent(:), dType(:)
  CHARACTER(LEN=maxlen) :: msg1, modname
  INTEGER :: it, nt, iq, nq, itr, ntr, ipr, npr
  LOGICAL :: ll
  modname = 'expandSection'
  lerr = .FALSE.
  nt = SIZE(tr)
  lerr = getKey('name',   tname,  tr(:)%keys);                 IF(lerr) RETURN
  lerr = getKey('parent', parent, tr(:)%keys, def = tran0);    IF(lerr) RETURN
  lerr = getKey('type',   dType,  tr(:)%keys, def = 'tracer'); IF(lerr) RETURN
  nq = 0
  !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  DO it = 1, nt    !=== GET TRACERS NB AFTER EXPANSION + NEEDED KEYS (name, parent, type)
  !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
    !--- Extract useful keys: parent name, type, component name
    tr(it)%component = sname
    CALL addKey('component', sname,  tr(it)%keys)

    !--- Determine the number of tracers and parents ; coherence checking
    ll = strCount( tname(it), ',', ntr)
    ll = strCount(parent(it), ',', npr)

    !--- Tagging tracers only can have multiple parents
    lerr = npr /=1 .AND. TRIM(dType(it)) /= 'tag'
    IF(lerr) THEN
      msg1 = 'Check section "'//TRIM(sname)//'"'
      IF(PRESENT(fname)) msg1 = TRIM(msg1)//' in file "'//TRIM(fname)//'"'
      CALL msg(TRIM(msg1)//': "'//TRIM(tname(it))//'" has several parents but is not a tag', modname); RETURN
    END IF
    nq = nq + ntr*npr                 
  !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  END DO
  !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

  ALLOCATE(ttr(nq))
  iq = 1
  !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  DO it = 1, nt                                                      !=== EXPAND TRACERS AND PARENTS NAMES LISTS
  !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
    ll = strParse( tname(it), ',', ta, ntr)                          !--- Number of tracers
    ll = strParse(parent(it), ',', pa, npr)                          !--- Number of parents
    DO ipr = 1, npr                                                  !--- Loop on parents list elts
      DO itr = 1, ntr                                                !--- Loop on tracers list elts
        ttr(iq)%keys%name = TRIM(ta(itr))
        ttr(iq)%keys%key  = tr(it)%keys%key
        ttr(iq)%keys%val  = tr(it)%keys%val
        ttr(iq)%name      = TRIM(ta(itr))
        ttr(iq)%parent    = TRIM(pa(ipr))
        ttr(iq)%type      = dType(it)
        ttr(iq)%component = sname
        CALL addKey('name',      ta(itr),   ttr(iq)%keys)
        CALL addKey('parent',    pa(ipr),   ttr(iq)%keys)
        CALL addKey('type',      dType(it), ttr(iq)%keys)
        CALL addKey('component', sname,     ttr(iq)%keys)
        iq = iq + 1
      END DO
    END DO
  !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  END DO
  !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  DEALLOCATE(ta,pa)
  CALL MOVE_ALLOC(FROM=ttr, TO=tr)

END FUNCTION expandSection
!==============================================================================================================================


!==============================================================================================================================
LOGICAL FUNCTION setGeneration(tr) RESULT(lerr)
!------------------------------------------------------------------------------------------------------------------------------
! Purpose: Determine, for each tracer of "tr(:)":
!   * iGeneration: the generation number 
!   * gen0Name:    the generation 0 ancestor name
!          Check also for orphan tracers (tracers without parent).
!------------------------------------------------------------------------------------------------------------------------------
  TYPE(trac_type),     INTENT(INOUT) :: tr(:)                        !--- Tracer derived type vector
  INTEGER                            :: iq, jq, ig
  CHARACTER(LEN=maxlen), ALLOCATABLE :: parent(:), tname(:)
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen) :: modname
  modname = 'setGeneration'
  lerr = getKey('name',   tname,  ky=tr(:)%keys); IF(lerr) RETURN
  lerr = getKey('parent', parent, ky=tr(:)%keys); IF(lerr) RETURN
  DO iq = 1, SIZE(tr)
    jq = iq; ig = 0
    DO WHILE(parent(jq) /= tran0)
      jq = strIdx(tname(:), parent(jq))
      lerr = jq == 0
      CALL msg('Orphan tracer "'//TRIM(tname(iq))//'"', modname, lerr); IF(lerr) RETURN
      ig = ig + 1
    END DO
    tr(iq)%gen0Name = tname(jq)
    tr(iq)%iGeneration = ig
    CALL addKey('iGeneration',   ig,  tr(iq)%keys)
    CALL addKey('gen0Name', tname(jq), tr(iq)%keys)
  END DO
END FUNCTION setGeneration
!==============================================================================================================================


!==============================================================================================================================
LOGICAL FUNCTION checkTracers(tr, sname, fname) RESULT(lerr)
!------------------------------------------------------------------------------------------------------------------------------
! Purpose:
!   * check for orphan tracers (without parent)
!   * check wether the phases are known or not (elements of "known_phases")
!------------------------------------------------------------------------------------------------------------------------------
  TYPE(trac_type),            INTENT(IN) :: tr(:)                    !--- Tracer derived type vector
  CHARACTER(LEN=*),           INTENT(IN) :: sname                    !--- Section name
  CHARACTER(LEN=*), OPTIONAL, INTENT(IN) :: fname                    !--- File name
  CHARACTER(LEN=1) :: p
  CHARACTER(LEN=maxlen) :: mesg
  CHARACTER(LEN=maxlen) :: bp(SIZE(tr, DIM=1)), pha                  !--- Bad phases list, phases of current tracer
  CHARACTER(LEN=maxlen), ALLOCATABLE :: tname(:)
  INTEGER,               ALLOCATABLE ::  iGen(:)
  INTEGER :: ip, np, iq, nq
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen) :: modname
  modname = 'checkTracers'
  nq = SIZE(tr,DIM=1)                                                !--- Number of tracers lines
  mesg = 'Check section "'//TRIM(sname)//'"'
  IF(PRESENT(fname)) mesg=TRIM(mesg)//' in file "'//TRIM(fname)//'"'
  lerr = getKey('iGeneration', iGen, tr(:)%keys);               IF(lerr) RETURN
  lerr = getKey('name',       tname, tr(:)%keys);               IF(lerr) RETURN

  !=== CHECK FOR ORPHAN TRACERS
  lerr = checkList(tname, iGen==-1, mesg, 'tracers', 'orphan'); IF(lerr) RETURN

  !=== CHECK PHASES
  DO iq = 1, nq; IF(iGen(iq) /= 0) CYCLE                             !--- Generation O only is checked
    IF(getKey(['phases','phase '], pha, iq, tr(:)%keys, lDisp=.FALSE.)) pha = 'g'   !--- Phase
    np = LEN_TRIM(pha); bp(iq)=' '
    DO ip = 1, np; p = pha(ip:ip); IF(INDEX(known_phases, p) == 0) bp(iq) = TRIM(bp(iq))//p; END DO
    IF(TRIM(bp(iq)) /= '') bp(iq) = TRIM(tname(iq))//': '//TRIM(bp(iq))
  END DO
  lerr = checkList(bp, iGen == 0 .AND. bp /= '', mesg, 'tracers phases', 'unknown')
END FUNCTION checkTracers
!==============================================================================================================================


!==============================================================================================================================
LOGICAL FUNCTION checkUnique(tr, sname, fname) RESULT(lerr)
!------------------------------------------------------------------------------------------------------------------------------
! Purpose: Make sure that tracers are not repeated.
!------------------------------------------------------------------------------------------------------------------------------
  TYPE(trac_type),            INTENT(IN) :: tr(:)                    !--- Tracer derived type vector
  CHARACTER(LEN=*),           INTENT(IN) :: sname                    !--- Section name
  CHARACTER(LEN=*), OPTIONAL, INTENT(IN) :: fname                    !--- File name
!------------------------------------------------------------------------------------------------------------------------------
  INTEGER :: ip, np, iq, nq, k
  LOGICAL, ALLOCATABLE  :: ll(:)
  CHARACTER(LEN=maxlen) :: mesg, phase, tdup(SIZE(tr,DIM=1))
  CHARACTER(LEN=maxlen), ALLOCATABLE :: tname(:), dType(:)
  INTEGER,               ALLOCATABLE :: iGen(:)
  CHARACTER(LEN=1) :: p
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen) :: modname
  modname = 'checkUnique'
  mesg = 'Check section "'//TRIM(sname)//'"'
  IF(PRESENT(fname)) mesg=TRIM(mesg)//' in file "'//TRIM(fname)//'"'
  nq=SIZE(tr,DIM=1); lerr=.FALSE.                                    !--- Number of lines ; error flag
  tdup(:) = ''
  lerr = getKey('name',       tname, tr%keys); IF(lerr) RETURN
  lerr = getKey('type',       dType, tr%keys); IF(lerr) RETURN
  lerr = getKey('iGeneration', iGen, tr%keys); IF(lerr) RETURN
  DO iq = 1, nq
    IF(dType(iq) == 'tag') CYCLE                                     !--- Tags can be repeated
    ll = tname==TRIM(tname(iq))                                      !--- Mask for current tracer name
    IF(COUNT(ll) == 1) CYCLE                                         !--- Tracer is not repeated
    IF(iGen(iq) > 0) THEN
      tdup(iq) = tname(iq)                                           !--- gen>0: MUST be unique
    ELSE
      DO ip = 1, nphases; p = known_phases(ip:ip)                    !--- Loop on known phases
        np = 0
        DO k = 1, nq
          IF(.NOT.ll(k)) CYCLE                                       !--- Skip tracers different from current one
          IF(getKey(['phases','phase '], phase, k, tr%keys, lDisp=.FALSE.)) phase='g'!--- Get current phases
          IF(INDEX(phase, p) /= 0) np = np + 1                       !--- One more appearance of current tracer with phase "p"
        END DO
        IF(np <= 1) CYCLE                                            !--- Regular case: no or a single appearance
        tdup(iq) = TRIM(tdup(iq))//TRIM(phases_names(ip))            !--- Repeated phase
        IF(ANY(tdup(1:iq-1) == tdup(iq))) tdup(iq)=''                !--- Avoid repeating same messages
      END DO
    END IF
    IF(tdup(iq) /= '') tdup(iq)=TRIM(tname(iq))//' in '//TRIM(tdup(iq))//' phase(s)'
  END DO
  lerr = checkList(tdup, tdup/='', mesg, 'tracers', 'duplicated')
END FUNCTION checkUnique
!==============================================================================================================================


!==============================================================================================================================
LOGICAL FUNCTION expandPhases(tr) RESULT(lerr)
!------------------------------------------------------------------------------------------------------------------------------
! 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
!------------------------------------------------------------------------------------------------------------------------------
  TYPE(trac_type),       ALLOCATABLE :: ttr(:)
  INTEGER,               ALLOCATABLE ::  i0(:), iGen(:)
  CHARACTER(LEN=maxlen), ALLOCATABLE :: tname(:), gen0N(:), phase(:), parents(:), dType(:)
  CHARACTER(LEN=maxlen)              ::  nam,     gen0Nm,   pha,      parent
  CHARACTER(LEN=1) :: p
  INTEGER :: ip, np, iq, jq, nq, it, nt, nc, i, n
  LOGICAL :: lTag, lExt
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen) :: modname
  modname = 'expandPhases'
  nq = SIZE(tr, DIM=1)
  nt = 0
  lerr = getKey('name',       tname, tr%keys); IF(lerr) RETURN       !--- Names of the tracers
  lerr = getKey('gen0Name',   gen0N, tr%keys); IF(lerr) RETURN       !--- Names of the tracers of first generation
  lerr = getKey('iGeneration', iGen, tr%keys); IF(lerr) RETURN       !--- Generation number
  lerr = getKey('phases',     phase, tr%keys); IF(lerr) RETURN       !--- Phases names
  lerr = getKey('parent',   parents, tr%keys); IF(lerr) RETURN       !--- Parents names
  lerr = getKey('type',       dType, tr%keys); IF(lerr) RETURN       !--- Tracers types ('tracer' or 'tag')
  DO iq = 1, nq                                                      !--- GET THE NUMBER OF TRACERS
    IF(iGen(iq) /= 0) CYCLE                                          !--- Only deal with generation 0 tracers
    nc = COUNT(gen0N == tname(iq) .AND. iGen /= 0)                   !--- Number of children of tr(iq)
    np = LEN_TRIM(phase(iq))                                         !--- Number of phases   of tr(iq)
    nt = nt + (1+nc) * np                                            !--- Number of tracers after expansion
  END DO
  ALLOCATE(ttr(nt))                                                  !--- Version  of "tr" after phases expansion
  it = 1                                                             !--- Current "ttr(:)" index
  DO iq = 1, nq                                                      !--- Loop on "tr(:)" indexes
    lTag = dType(iq)=='tag'                                          !--- Current tracer is a tag
    i0 = strFind(tname, TRIM(gen0N(iq)), n)                          !--- Indexes of first generation ancestor copies
    np = SUM([( LEN_TRIM(phase(i0(i))), i = 1, n )], 1)              !--- Number of phases for current tracer tr(iq)
    lExt = np > 1                                                    !--- Phase suffix only required if phases number is > 1
    IF(lTag) lExt = lExt .AND. iGen(iq) > 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(iGen(iq) == 0) jq = iq                                      !--- Generation 0: count the current tracer phases only
      pha = phase(jq)                                                !--- Phases list for tr(jq)
      DO ip = 1, LEN_TRIM(pha)                                       !=== LOOP ON PHASES LISTS
        p = pha(ip:ip)
        nam = tname(iq)                                              !--- Tracer name (regular case)
        IF(lTag) nam = TRIM(parents(iq))                             !--- Parent name (tagging case)
        IF(lExt) nam = addPhase(nam, p )                             !--- Phase extension needed
        IF(lTag) nam = TRIM(nam)//'_'//TRIM(tname(iq))               !--- <parent>_<name> for tags
        ttr(it) = tr(iq)                                             !--- Same <key>=<val> pairs
        ttr(it)%name      = TRIM(nam)                                !--- Name with possibly phase suffix
        ttr(it)%keys%name = TRIM(nam)                                !--- Name inside the keys decriptor
        ttr(it)%phase     = p                                        !--- Single phase entry
        CALL addKey('name', nam, ttr(it)%keys)
        CALL addKey('phase', p,  ttr(it)%keys)
        IF(lExt) THEN
          parent = parents(iq); IF(iGen(iq) > 0) parent = addPhase(parent, p)
          gen0Nm =   gen0N(iq); IF(iGen(iq) > 0) gen0Nm = addPhase(gen0Nm, p)
          ttr(it)%parent   = parent
          ttr(it)%gen0Name = gen0Nm
          CALL addKey('parent',   parent, ttr(it)%keys)
          CALL addKey('gen0Name', gen0Nm, ttr(it)%keys)
        END IF
        it = it+1
      END DO
      IF(iGen(iq) == 0) EXIT                                         !--- Break phase loop for gen 0
    END DO
  END DO
  CALL MOVE_ALLOC(FROM=ttr, TO=tr)
  CALL delKey(['phases'],tr)                                         !--- Remove few keys entries

END FUNCTION expandPhases
!==============================================================================================================================


!==============================================================================================================================
SUBROUTINE sortTracers(tr)
!------------------------------------------------------------------------------------------------------------------------------
! Purpose: Sort tracers:
!  * Put water at the beginning of the vector, in the "known_phases" order.
!  * lSortByGen == T: in ascending generations numbers.
!  * lSortByGen == F: tracer + its children sorted by growing generation, one after the other.
!   TO BE ADDED IF NECESSARY: HIGHER MOMENTS AT THE END
!------------------------------------------------------------------------------------------------------------------------------
  TYPE(trac_type), INTENT(INOUT) :: tr(:)                            !--- Tracer derived type vector
!------------------------------------------------------------------------------------------------------------------------------
  TYPE(trac_type),       ALLOCATABLE :: tr2(:)
  INTEGER,               ALLOCATABLE :: iy(:), iz(:)
  INTEGER,               ALLOCATABLE ::  iGen(:)
  CHARACTER(LEN=maxlen), ALLOCATABLE :: tname(:), gen0N(:)
  INTEGER :: ig, ng, iq, jq, ip, nq, n, ix(SIZE(tr)), k
  LOGICAL :: lerr
!  tr2 is introduced in order to cope with a bug in gfortran 4.8.5 compiler
!------------------------------------------------------------------------------------------------------------------------------
  lerr = getKey('iGeneration', iGen, tr%keys); IF(lerr) RETURN       !--- Generation number
  nq = SIZE(tr)
  DO ip = nphases, 1, -1
    lerr = getKey('name',     tname, tr%keys); IF(lerr) RETURN       !--- Names of the tracers of first generation
    iq = strIdx(tname, addPhase('H2O', ip))
    IF(iq == 0) CYCLE
    tr2 = tr(:)
    tr = [tr2(iq), tr2(1:iq-1), tr2(iq+1:nq)]
  END DO
  IF(lSortByGen) THEN
    iq = 1
    ng = MAXVAL(iGen, MASK=.TRUE., DIM=1)                            !--- Number of generations
    DO ig = 0, ng                                                    !--- Loop on generations
      iy = PACK([(k, k=1, nq)], MASK=iGen(:) == ig)                  !--- Generation ig tracers indexes
      n = SIZE(iy)
      ix(iq:iq+n-1) = iy                                             !--- Stack growing generations idxs
      iq = iq + n
    END DO
  ELSE
    lerr = getKey('gen0Name',   gen0N, tr%keys); IF(lerr) RETURN     !--- Names of the tracers    iq = 1
    DO jq = 1, nq                                                    !--- Loop on generation 0 tracers
      IF(iGen(jq) /= 0) CYCLE                                        !--- Skip generations /= 0
      ix(iq) = jq                                                    !--- Generation 0 ancestor index first
      iq = iq + 1                                                    !--- Next "iq" for next generations tracers
      iy = strFind(gen0N(:), TRIM(tname(jq)))                        !--- Indices of "tr(jq)" children in "tr(:)"
      ng = MAXVAL(iGen(iy), 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(iGen(iy), ig, n)                                   !--- Indices of the tracers "tr(iy(:))" of generation "ig"
        ix(iq:iq+n-1) = iy(iz)                                       !--- Same indexes in "tr(:)"
        iq = iq + n
      END DO
    END DO
  END IF
  tr = tr(ix)                                                        !--- Reorder the tracers
END SUBROUTINE sortTracers
!==============================================================================================================================


!==============================================================================================================================
LOGICAL FUNCTION mergeTracers(sections, tr) RESULT(lerr)
  TYPE(dataBase_type),  TARGET, INTENT(IN)  :: sections(:)
  TYPE(trac_type), ALLOCATABLE, INTENT(OUT) ::       tr(:)
  TYPE(trac_type), POINTER ::   t1(:),   t2(:)
  TYPE(keys_type), POINTER ::   k1(:),   k2(:)
  INTEGER,     ALLOCATABLE :: ixct(:), ixck(:)
  INTEGER :: is, ik, ik1, ik2, nk2, i1, i2, nt2
  CHARACTER(LEN=maxlen) :: s1, v1, v2, tnam, knam, modname
  CHARACTER(LEN=maxlen), ALLOCATABLE :: keys(:), n1(:), n2(:)
  modname = 'mergeTracers'
  lerr = .FALSE.
  keys = ['parent     ', 'type       ', 'iGeneration']               !--- Mandatory keys
  t1 => sections(1)%trac(:); k1 => t1(:)%keys                        !--- Alias: first tracers section, corresponding keys
  lerr = getKey('name', n1, k1); IF(lerr) RETURN                     !--- Names of the tracers
  tr = t1
  !----------------------------------------------------------------------------------------------------------------------------
  DO is=2,SIZE(sections)                                             !=== SEVERAL SECTIONS: MERGE THEM
  !----------------------------------------------------------------------------------------------------------------------------
    t2  => sections(is)%trac(:)                                      !--- Alias: current tracers section
    k2  => t2(:)%keys
    lerr = getKey('name', n2, k2); IF(lerr) RETURN                   !--- Names of the tracers
    nt2  = SIZE(t2(:), DIM=1)                                        !--- Number of tracers in section
    ixct = strIdx(n1(:), n2(:))                                      !--- Indexes of common tracers
    tr = [tr, PACK(t2, MASK= ixct==0)]                               !--- Append with new tracers
    IF( ALL(ixct == 0) ) CYCLE                                       !--- No common tracers => done
    CALL msg('Tracers defined in previous sections and duplicated in "'//TRIM(sections(is)%name)//'":', modname)
    CALL msg(n1(PACK(ixct, MASK = ixct/=0)), modname, nmax=128)      !--- Display duplicates (the 128 first at most)
    !--------------------------------------------------------------------------------------------------------------------------
    DO i2=1,nt2; tnam = TRIM(t2(i2)%name)                            !=== LOOP ON COMMON TRACERS
    !--------------------------------------------------------------------------------------------------------------------------
      i1 = ixct(i2); IF(i1 == 0) CYCLE                               !--- Idx in t1(:) ; skip new tracers

      !=== CHECK WETHER ESSENTIAL KEYS ARE IDENTICAL OR NOT
      s1=' of "'//TRIM(tnam)//'" in "'//TRIM(sections(is)%name)//'" not matching previous value'
      DO ik = 1, SIZE(keys)
        lerr = getKey(keys(ik), v1, i1, k1)
        lerr = getKey(keys(ik), v2, i2, k2)
        lerr = v1 /= v2; CALL msg(TRIM(keys(ik))//TRIM(s1), modname, lerr); IF(lerr) RETURN
      END DO

      !=== GET THE INDICES IN tr(i2)%keys%key(:) OF THE KEYS ALSO PRESENT IN tr(i1)%keys%key(:)
      nk2  =   SIZE(k2(i2)%key(:))                                   !--- Keys number in current section
      ixck = strIdx(k1(i1)%key(:), k2(i2)%key(:))                    !--- Common keys indexes
      !--- APPEND THE NEW KEYS PAIRS IN tr(i1)%keys%key(:)
      tr(i1)%keys%key = [ tr(i1)%keys%key, PACK(tr(i2)%keys%key, MASK = ixck==0)]
      tr(i1)%keys%val = [ tr(i1)%keys%val, PACK(tr(i2)%keys%val, MASK = ixck==0)]

      !=== KEEP TRACK OF THE COMPONENTS NAMES: COMA-SEPARATED LIST
      lerr = getKey('component', v1, i1, k1)
      lerr = getKey('component', v2, i2, k2)
      tr(i1)%component = TRIM(v1)//','//TRIM(v2)
      CALL addKey('component', TRIM(v1)//','//TRIM(v2), tr(i1)%keys)

      !=== FOR TRACERS COMMON TO PREVIOUS AND CURRENT SECTIONS: CHECK WETHER SOME KEYS HAVE DIFFERENT VALUES ; KEEP OLD ONE
      DO ik2 = 1, nk2                                                !--- Collect the corresponding indices
        ik1 = ixck(ik2); IF(ik1 == 0) CYCLE
        IF(k1(i1)%val(ik1) == k2(i2)%val(ik2)) ixck(ik2)=0
      END DO
      IF(ALL(ixck==0)) CYCLE                                         !--- No identical keys with /=values => nothing to display
      CALL msg('Key(s)'//TRIM(s1), modname)                          !--- Display the  keys with /=values (names list)
      DO ik2 = 1, nk2                                                !--- Loop on keys found in both t1(:) and t2(:)
        knam = k2(i2)%key(ik2)                                       !--- Name of the current key
        ik1 = ixck(ik2)                                              !--- Corresponding index in t1(:)
        IF(ik1 == 0) CYCLE                                           !--- New keys are skipped
        v1 = k1(i1)%val(ik1); v2 = k2(i2)%val(ik2)                   !--- Key values in t1(:) and t2(:)
        CALL msg(' * '//TRIM(knam)//'='//TRIM(v2)//' ; previous value kept:'//TRIM(v1), modname)
      END DO
      !------------------------------------------------------------------------------------------------------------------------
    END DO
    !--------------------------------------------------------------------------------------------------------------------------
  END DO
  CALL sortTracers(tr)

END FUNCTION mergeTracers
!==============================================================================================================================

!==============================================================================================================================
LOGICAL FUNCTION cumulTracers(sections, tr, lRename) RESULT(lerr)
  TYPE(dataBase_type),  TARGET, INTENT(IN)  :: sections(:)
  TYPE(trac_type), ALLOCATABLE, INTENT(OUT) ::       tr(:)
  LOGICAL,            OPTIONAL, INTENT(IN)  :: lRename               !--- .TRUE.: add a section suffix to identical names
  CHARACTER(LEN=maxlen)  :: tnam_new, modname
  CHARACTER(LEN=maxlen), ALLOCATABLE :: tname(:), parent(:), comp(:)
  INTEGER :: iq, jq, is
  modname = 'cumulTracers'
  lerr = .FALSE.
  tr = [( sections(is)%trac(:), is = 1, SIZE(sections) )]            !--- Concatenated tracers vector
  IF(PRESENT(lRename)) THEN; IF(lRename) RETURN; END IF              !--- No renaming: finished
  lerr = getKey('name',     tname, tr%keys); IF(lerr) RETURN         !--- Names
  lerr = getKey('parent',  parent, tr%keys); IF(lerr) RETURN         !--- Parents
  lerr = getKey('component', comp, tr%keys); IF(lerr) RETURN         !--- Component name
  !----------------------------------------------------------------------------------------------------------------------------
  DO iq = 1, SIZE(tr); IF(COUNT(tname == tname(iq)) == 1) CYCLE      !=== LOOP ON TRACERS
  !----------------------------------------------------------------------------------------------------------------------------
    tnam_new = TRIM(tname(iq))//'_'//TRIM(comp(iq))                  !--- Same with section extension
    CALL addKey('name', tnam_new, tr(iq)%keys)                       !--- Modify tracer name
    tr(iq)%name = TRIM(tnam_new)                                     !--- Modify tracer name
    !--------------------------------------------------------------------------------------------------------------------------
    DO jq = 1, SIZE(tr); IF(parent(jq) /= tname(iq)) CYCLE           !=== LOOP ON TRACERS PARENTS
    !--------------------------------------------------------------------------------------------------------------------------
      CALL addKey('parent', tnam_new, tr(jq)%keys)                   !--- Modify tracer name
      tr(jq)%parent = TRIM(tnam_new)                                 !--- Modify tracer name
    !--------------------------------------------------------------------------------------------------------------------------
    END DO
  !----------------------------------------------------------------------------------------------------------------------------
  END DO
  !----------------------------------------------------------------------------------------------------------------------------
  CALL sortTracers(tr)
END FUNCTION cumulTracers
!==============================================================================================================================


!==============================================================================================================================
LOGICAL  FUNCTION  dispTraSection(message, sname, modname) RESULT(lerr)
  CHARACTER(LEN=*), INTENT(IN) :: message, sname, modname
  CHARACTER(LEN=maxlen), ALLOCATABLE :: s(:,:), n(:), tmp(:)
  CHARACTER(LEN=maxlen) :: p
  INTEGER :: idb, iq, nq
  idb = strIdx(dBase(:)%name, sname); IF(idb == 0) RETURN
  nq = SIZE(dBase(idb)%trac)
  p = ''
  CALL append(['iq'],     .TRUE. ); IF(lerr) RETURN
  CALL append(['name'],   .TRUE. ); IF(lerr) RETURN
  CALL append(['phases','phase '], .FALSE., 'pha'); IF(lerr) RETURN
  CALL append(['hadv'],   .TRUE. ); IF(lerr) RETURN
  CALL append(['vadv'],   .TRUE. ); IF(lerr) RETURN
  CALL append(['parent'], .FALSE.); IF(lerr) RETURN
  CALL append(['iGen'],   .FALSE.); IF(lerr) RETURN
  CALL msg(TRIM(message)//':', modname)
  lerr = dispTable(p, n, s, nColMax=maxTableWidth, nHead=2, sub=modname); IF(lerr) RETURN

CONTAINS

SUBROUTINE append(nam, lMandatory, snam)
! Test whether key named "nam(:)" is available.
!  * yes: - get its value for all species in "tmp(:)" and append table "s(:,:)" with it
!         - append titles list with "nam(1)" (or, if specified, "snam", usually a short name).
!  * no:  return to calling routine with an error flag if the required key is mandatory
  CHARACTER(LEN=*),           INTENT(IN) :: nam(:)
  LOGICAL,                    INTENT(IN) :: lMandatory
  CHARACTER(LEN=*), OPTIONAL, INTENT(IN) :: snam
  INTEGER :: m
  CHARACTER(LEN=maxlen), ALLOCATABLE :: n0(:)
  CHARACTER(LEN=maxlen) :: nm
  CHARACTER(LEN=maxlen) :: tmp2(nq)

  lerr = .FALSE.
  IF(nam(1) == 'iq') THEN
    tmp2 = num2str([(iq, iq=1, nq)])
    tmp = tmp2
  ELSE
    lerr = getKey(nam, tmp, dBase(idb)%trac(:)%keys, lDisp=lMandatory)
  END IF
  IF(lerr) THEN; lerr = lMandatory; RETURN; END IF
  nm = nam(1); IF(PRESENT(snam)) nm = snam
  p = TRIM(p)//'s'
  IF(ALLOCATED(s)) THEN; s = cat(s, tmp); ELSE; ALLOCATE(s(nq,1)); s(:,1) = tmp; END IF
  IF(ALLOCATED(n)) THEN; m = SIZE(n); ALLOCATE(n0(m+1)); n0(1:m)=n; n0(m+1)=nm; CALL MOVE_ALLOC(FROM=n0, TO=n)
  ELSE; n=nam(1:1); END IF
END SUBROUTINE append

END FUNCTION dispTraSection
!==============================================================================================================================


!==============================================================================================================================
!=== CREATE TRACER(S) ALIAS: SCALAR/VECTOR FROM NAME(S) OR INDICE(S) ==========================================================
!==============================================================================================================================
LOGICAL FUNCTION aliasTracer(tname, trac, alias) RESULT(lerr)                  !=== TRACER NAMED "tname" - SCALAR
  CHARACTER(LEN=*),         INTENT(IN)  :: tname
  TYPE(trac_type), TARGET,  INTENT(IN)  :: trac(:)
  TYPE(trac_type), POINTER, INTENT(OUT) :: alias
  INTEGER :: it
  CHARACTER(LEN=maxlen), ALLOCATABLE :: tnames(:)
  alias => NULL()
  lerr = getKey('name', tnames, trac(:)%keys)
  it = strIdx(tnames, tname)
  lerr = it /= 0; IF(.NOT.lerr) alias => trac(it)
END FUNCTION aliasTracer
!==============================================================================================================================
LOGICAL FUNCTION trSubset_Indx(trac, idx, alias) RESULT(lerr)                  !=== TRACERS WITH INDICES "idx(:)" - VECTOR
  TYPE(trac_type), ALLOCATABLE, INTENT(IN)  ::  trac(:)
  INTEGER,                      INTENT(IN)  ::   idx(:)
  TYPE(trac_type), ALLOCATABLE, INTENT(OUT) :: alias(:)
  alias = trac(idx)
  lerr = indexUpdate(alias)
END FUNCTION trSubset_Indx
!------------------------------------------------------------------------------------------------------------------------------
LOGICAL FUNCTION trSubset_Name(trac, tname, alias) RESULT(lerr)                !=== TRACERS NAMED "tname(:)" - VECTOR
  TYPE(trac_type), ALLOCATABLE, INTENT(IN)  ::  trac(:)
  CHARACTER(LEN=*),             INTENT(IN)  :: tname(:)
  TYPE(trac_type), ALLOCATABLE, INTENT(OUT) :: alias(:)
  CHARACTER(LEN=maxlen), ALLOCATABLE :: tnames(:)
  lerr = getKey('name', tnames, trac(:)%keys)
  alias = trac(strIdx(tnames, tname))
  lerr = indexUpdate(alias)
END FUNCTION trSubset_Name
!==============================================================================================================================
LOGICAL FUNCTION trSubset_gen0Name(trac, gen0Nm, alias) RESULT(lerr)           !=== TRACERS OF COMMON 1st GENERATION ANCESTOR
  TYPE(trac_type), ALLOCATABLE, INTENT(IN)  :: trac(:)
  CHARACTER(LEN=*),             INTENT(IN)  :: gen0Nm
  TYPE(trac_type), ALLOCATABLE, INTENT(OUT) :: alias(:)
  CHARACTER(LEN=maxlen), ALLOCATABLE :: gen0N(:)
  lerr = getKey('gen0Name', gen0N, trac(:)%keys)
  alias = trac(strFind(delPhase(gen0N), gen0Nm))
  lerr = indexUpdate(alias)
END FUNCTION trSubset_gen0Name
!==============================================================================================================================


!==============================================================================================================================
!=== UPDATE THE INDEXES iqParent, iqDescend AND iGeneration IN THE TRACERS DESCRIPTOR LIST "tr" (USEFULL FOR SUBSETS) =========
!==============================================================================================================================
LOGICAL FUNCTION indexUpdate(tr) RESULT(lerr)
  TYPE(trac_type), INTENT(INOUT) :: tr(:)
  INTEGER :: iq, jq, nq, ig, nGen
  INTEGER,               ALLOCATABLE :: iqDescen(:), ix(:), iy(:)
  CHARACTER(LEN=maxlen), ALLOCATABLE :: tnames(:), parent(:)
  INTEGER,       DIMENSION(SIZE(tr)) :: iqParent, iGen
  lerr = getKey('name',   tnames, tr%keys); IF(lerr) RETURN          !--- Names
  lerr = getKey('parent', parent, tr%keys); IF(lerr) RETURN          !--- Parents
  nq = SIZE(tr)

  !=== iqParent, iGeneration
  DO iq = 1, nq; iGen(iq) = 0; jq = iq
    iqParent(iq) = strIdx(tnames, parent(iq))
    DO; jq = strIdx(tnames, parent(jq)); IF(jq == 0) EXIT; iGen(iq) = iGen(iq) + 1; END DO
    CALL addKey('iqParent', iqParent(iq), tr(iq)%keys)
    CALL addKey('iqGeneration', iGen(iq), tr(iq)%keys)
    tr(iq)%iqParent = iqParent(iq)
    tr(iq)%iGeneration =  iGen(iq)
  END DO

  !=== nqChildren, iqDescen, nqDescen
  nGen = MAXVAL(iGen, MASK=.TRUE.)
  DO iq = 1, nq
    ix = [iq]; ALLOCATE(iqDescen(0))
    DO ig = iGen(iq)+1, nGen
      iy = find(iqParent, ix); iqDescen = [iqDescen, iy]; ix = iy
      IF(ig /= iGen(iq)+1) CYCLE
      CALL addKey('nqChildren', SIZE(iqDescen), tr(iq)%keys)
      tr(iq)%nqChildren = SIZE(iqDescen)
    END DO
    CALL addKey('iqDescen', strStack(num2str(iqDescen)), tr(iq)%keys)
    CALL addKey('nqDescen',             SIZE(iqDescen),  tr(iq)%keys)
    tr(iq)%iqDescen =      iqDescen
    tr(iq)%nqDescen = SIZE(iqDescen)
    DEALLOCATE(iqDescen)
  END DO
END FUNCTION indexUpdate
!==============================================================================================================================
  
  
!==============================================================================================================================
!=== READ FILE "fnam" TO APPEND THE "dBase" TRACERS DATABASE WITH AS MUCH SECTIONS AS PARENTS NAMES IN "isot(:)%parent":   ====
!===  * Each section dBase(i)%name contains the isotopes "dBase(i)%trac(:)" descending on "dBase(i)%name"="iso(i)%parent"  ====
!===  * For each isotopes family, the <key>=<val> vector of each tracer is moved into the isotopes descriptor "isot"       ====
!=== NOTES:                                                                                                                ====
!===  * Most of the "isot" components have been defined in the calling routine (processIsotopes):                          ====
!===      parent,  nzone, zone(:),  niso, keys(:)%name,  ntiso, trac(:),  nphas, phas,  iqIsoPha(:,:),  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             ====
!===  * In case keys are found both in the "params" section and the "*.def" file, the later value is retained              ====
!===  * On each isotope line, defined keys can be used for other keys defintions (single level depth substitution)         ====
!===  * The routine gives an error if a required isotope is not available in the database stored in "fnam"                 ====
!==============================================================================================================================
LOGICAL FUNCTION readIsotopesFile(fnam, isot) RESULT(lerr)
  CHARACTER(LEN=*),        INTENT(IN)    :: fnam                     !--- Input file name
  TYPE(isot_type), TARGET, INTENT(INOUT) :: isot(:)                  !--- Isotopes descriptors (field %parent must be defined!)
  LOGICAL :: lFound
  INTEGER :: is, iis, it, idb, ndb, nb0
  CHARACTER(LEN=maxlen), ALLOCATABLE :: vals(:)
  CHARACTER(LEN=maxlen)              :: modname
  TYPE(trac_type),           POINTER ::   tt(:), t
  TYPE(dataBase_type),   ALLOCATABLE ::  tdb(:)
  modname = 'readIsotopesFile'

  !--- THE INPUT FILE MUST BE PRESENT
  INQUIRE(FILE=TRIM(fnam), EXIST=lFound); lerr = .NOT.lFound
  CALL msg('Missing isotopes parameters file "'//TRIM(fnam)//'"', modname, lerr); IF(lerr) RETURN

  !--- READ THE FILE SECTIONS, ONE EACH ISOTOPES FAMILY
  nb0 = SIZE(dBase, DIM=1)+1                                         !--- Next database element index
  lerr = readSections(fnam,strStack(isot(:)%parent,'|')); IF(lerr) RETURN !--- Read sections, one each isotopes family
  ndb = SIZE(dBase, DIM=1)                                           !--- Current database size
  DO idb = nb0, ndb
    iis = idb-nb0+1

    !--- GET FEW GLOBAL KEYS FROM "def" FILES AND ADD THEM TO THE 'params' SECTION
!    CALL addKeysFromDef(dBase(idb)%trac, 'params')

    !--- SUBSTITUTE THE KEYS DEFINED IN THE 'params' VIRTUAL TRACER ; SUBSTITUTE LOCAL KEYS ; REMOVE 'params' VIRTUAL TRACER
    CALL subDefault(dBase(idb)%trac, 'params', .TRUE.)

    tt => dBase(idb)%trac

    !--- REDUCE THE EXPRESSIONS TO OBTAIN SCALARS AND TRANSFER THEM TO THE "isot" ISOTOPES DESCRIPTORS VECTOR
    DO it = 1, SIZE(dBase(idb)%trac)
      t => dBase(idb)%trac(it)
      is = strIdx(isot(iis)%keys(:)%name, t%name)                    !--- Index in "isot(iis)%keys(:)%name" of isotope "t%name"
      IF(is == 0) CYCLE
      lerr = ANY(reduceExpr(t%keys%val, vals)); IF(lerr) RETURN      !--- Reduce expressions ; detect non-numerical elements
      isot(iis)%keys(is)%key = t%keys%key
      isot(iis)%keys(is)%val = vals
    END DO

    !--- CHECK FOR MISSING ISOTOPES (NO KEYS ALLOCATED)
    lerr = checkList(isot(iis)%keys(:)%name, .NOT.[( ALLOCATED(isot(iis)%keys(is)%key), is=1, SIZE(isot(iis)%keys) )], &
                     'Check file "'//TRIM(fnam)//'" in section "'//TRIM(dBase(idb)%name)//'"', 'isotopes', 'missing')
    IF(lerr) RETURN
  END DO

  !--- CLEAN THE DATABASE ENTRIES
  IF(nb0 == 1) THEN
    DEALLOCATE(dBase); ALLOCATE(dBase(0))
  ELSE
    ALLOCATE(tdb(nb0-1)); tdb(1:nb0-1)=dBase(1:nb0-1); CALL MOVE_ALLOC(FROM=tdb, TO=dBase)
  END IF

  lerr = dispIsotopes()

CONTAINS

!------------------------------------------------------------------------------------------------------------------------------
LOGICAL FUNCTION dispIsotopes() RESULT(lerr)
  INTEGER :: ik, nk, ip, it, nt
  CHARACTER(LEN=maxlen) :: prf
  CHARACTER(LEN=maxlen), ALLOCATABLE :: ttl(:), val(:,:)
  CALL msg('Isotopes parameters read from file "'//TRIM(fnam)//'":', modname)
  DO ip = 1, SIZE(isot)                                              !--- Loop on isotopes families
    nk = SIZE(isot(ip)%keys(1)%key)                                  !--- Same keys for each isotope
    nt = SIZE(isot(ip)%keys)                                         !--- Number of isotopes
    prf = 'i'//REPEAT('s',nk+1)                                      !--- Profile for table printing
    ALLOCATE(ttl(nk+2), val(nt,nk+1))
    ttl(1:2) = ['it  ','name']; ttl(3:nk+2) = isot(ip)%keys(1)%key(:)!--- Titles line with keys names
    val(:,1) = isot(ip)%keys(:)%name                                 !--- Values table 1st column: isotopes names  
    DO ik = 1, nk
      DO it = 1, nt
        val(it,ik+1) = isot(ip)%keys(it)%val(ik)                     !--- Other columns: keys values
      END DO
    END DO
    lerr = dispTable(prf, ttl, val, cat([(it,it=1,nt)]), rFmt='(EN8.4)', nColMax=maxTableWidth, nHead=2, sub=modname)
    CALL msg('Problem with the table content', modname, lerr); IF(lerr) RETURN
    DEALLOCATE(ttl, val)
  END DO        
END FUNCTION dispIsotopes
!------------------------------------------------------------------------------------------------------------------------------

END FUNCTION readIsotopesFile
!==============================================================================================================================


!==============================================================================================================================
!=== IF ISOTOPES (2ND GENERATION TRACERS) ARE DETECTED:                                                                     ===
!===    * COMPUTE MOST OF THE RELATED QUANTITIES ("isot" COMPONENTS).                                                       ===
!===    * COMPUTE FEW ISOTOPES-DEDICATED "trac" COMPONENTS                                                                  ===
!===    * CALL readIsotopesFile TO GET PHYSICAL QUANTITIES (<key>=<val> PAIRS)                                              ===
!===      NOTE: THIS IS DONE HERE (IN A ROUTINE CALLED BY THE DYNAMIC), BECAUSE THE DYNAMIC NEEDS FEW PHYSICAL PARAMETERS.  ===
!==============================================================================================================================
LOGICAL FUNCTION processIsotopes(iNames) RESULT(lerr)
  CHARACTER(LEN=maxlen), TARGET, OPTIONAL, INTENT(IN)  :: iNames(:)
  CHARACTER(LEN=maxlen), ALLOCATABLE :: p(:), str(:)                 !--- Temporary storage
  CHARACTER(LEN=maxlen), ALLOCATABLE :: tname(:), parent(:), dType(:), phase(:), gen0N(:)
  CHARACTER(LEN=maxlen) :: iName, modname
  CHARACTER(LEN=1)   :: ph                                           !--- Phase
  INTEGER, ALLOCATABLE ::  iGen(:)
  INTEGER :: ic, ip, iq, it, iz
  LOGICAL, ALLOCATABLE :: ll(:)                                      !--- Mask
  TYPE(trac_type), POINTER   ::  t(:), t1
  TYPE(isot_type), POINTER   ::  i

  lerr = .FALSE.
  modname = 'readIsotopesFile'

  t => tracers

  lerr = getKey('name',       tname, t%keys); IF(lerr) RETURN       !--- Names
  lerr = getKey('parent',    parent, t%keys); IF(lerr) RETURN       !--- Parents
  lerr = getKey('type',       dType, t%keys); IF(lerr) RETURN       !--- Tracer type
  lerr = getKey('phase',      phase, t%keys); IF(lerr) RETURN       !--- Phase
  lerr = getKey('gen0Name',   gen0N, t%keys); IF(lerr) RETURN       !--- 1st generation ancestor name
  lerr = getKey('iGeneration', iGen, t%keys); IF(lerr) RETURN       !--- Generation number

  !--- GET FROM "tracers" THE FULL LIST OF AVAILABLE ISOTOPES FAMILIES
  p = PACK(delPhase(parent), MASK = dType=='tracer' .AND. iGen==1)
  CALL strReduce(p, nbIso)

  !--- CHECK WHETHER NEEDED ISOTOPES FAMILIES "iNames" ARE AVAILABLE OR NOT
  IF(PRESENT(iNames)) THEN
    DO it = 1, SIZE(iNames)
      lerr = ALL(p /= iNames(it))
      CALL msg('No isotopes family "'//TRIM(iNames(it))//'" found among tracers', modname, lerr); IF(lerr) RETURN
    END DO
    p = iNames; nbIso = SIZE(p)
  END IF
  IF(ALLOCATED(isotopes)) DEALLOCATE(isotopes)
  ALLOCATE(isotopes(nbIso), isoFamilies(nbIso))

  IF(nbIso==0) RETURN                                                !=== NO ISOTOPES: FINISHED

  !--- ISOTOPES RELATED VARIABLES ; NULL OR EMPTY IF NO ISOTOPES
  isotopes(:)%parent = p
  DO ic = 1, SIZE(p)                                                 !--- Loop on isotopes families
    i => isotopes(ic)
    iname = i%parent                                                 !--- Current isotopes family name (parent tracer name)

    !=== Isotopes children of tracer "iname": mask, names, number (same for each phase of "iname")
    ll = dType=='tracer' .AND. delPhase(parent) == iname .AND. phase == 'g'
    str = PACK(delPhase(tname), MASK = ll)                           !--- Effectively found isotopes of "iname"
    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 = dType=='tag'    .AND. delPhase(gen0N) == iname .AND. iGen == 2
    i%zone = PACK(strTail(tname,'_',.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 children of tracer "iname" (same for each phase of "iname")
    !    NOTE: One might like to create a similar variable for 2nd generation tagging tracers (tagging the gen1 tracers)
    str = PACK(delPhase(tname), MASK=ll)
    CALL strReduce(str)
    i%ntiso = i%niso + SIZE(str)                                     !--- Number of isotopes + their geographic tracers [ntiso]
    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"
    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 => tracers(iq)
      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(i%trac, strHead(delPhase(t1%name),'_',.TRUE.)) !--- Current isotope       idx in effective isotopes list
      t1%iso_iZone  = strIdx(i%zone,          strTail(t1%name, '_',.TRUE.)) !--- 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
    i%iqIsoPha = 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 iq (index in dyn array, size nqtot) from the water and 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
    i%iqWIsoPha = RESHAPE( [( [strIdx(t%name,   addPhase('H2O',i%phase(ip:ip))), i%iqIsoPha(:,ip)], ip=1,i%nphas)], &
                         [1+i%ntiso, i%nphas] )
    !=== Table used to get ix (index in tagging tracers isotopes list, size ntiso) from the zone and isotope indexes
    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

  !=== READ PHYSICAL PARAMETERS FROM isoFile FILE
!  lerr = readIsotopesFile(isoFile, isotopes); IF(lerr) RETURN! on commente pour ne pas chercher isotopes_params.def

  !--- GET THE isoCheck ENTRY FROM THE *.DEF FILES (MIGHT BE CHANGED TO A CLASS-DEPENDANT KEYWORD)
  CALL getin_p('ok_iso_verif', isoCheck, .TRUE.)

  !=== CHECK CONSISTENCY
  lerr = testIsotopes(); IF(lerr) RETURN

  !=== SELECT WATER ISOTOPES CLASS OR, IF UNFOUND, THE FIRST ISOTOPES CLASS
  IF(.NOT.isoSelect('H2O', .TRUE.)) THEN; iH2O = ixIso; ELSE; lerr = isoSelect(1, .TRUE.); END IF

  !=== COLLECT THE NAMES OF THE ISOTOPES FAMILIES
  isoFamilies = isotopes(:)%parent

CONTAINS

!------------------------------------------------------------------------------------------------------------------------------
LOGICAL FUNCTION testIsotopes() RESULT(lerr)     !--- MAKE SURE MEMBERS OF AN ISOTOPES FAMILY ARE PRESENT IN THE SAME PHASES
!------------------------------------------------------------------------------------------------------------------------------
  INTEGER :: ix, it, ip, np, iz, nz, npha, nzon
  TYPE(isot_type), POINTER :: i
  DO ix = 1, nbIso
    i => isotopes(ix)
    !--- Check whether each isotope and tagging isotopic tracer is present in the same number of phases
    DO it = 1, i%ntiso; npha = i%nphas
      np = SUM([(COUNT(tracers(:)%name == addPhase(i%trac(it), i%phase(ip:ip))), ip=1, npha)])
      lerr = np /= npha
      CALL msg(TRIM(num2str(np))//' phases instead of '//TRIM(num2str(npha))//' for '//TRIM(i%trac(it)), modname, lerr)
      IF(lerr) RETURN
    END DO
    DO it = 1, i%niso; nzon = i%nzone
      nz = SUM([(COUNT(i%trac == TRIM(i%trac(it))//'_'//i%zone(iz)), iz=1, nzon)])
      lerr = nz /= nzon
      CALL msg(TRIM(num2str(nz))//' tagging zones instead of '//TRIM(num2str(nzon))//' for '//TRIM(i%trac(it)), modname, lerr)
      IF(lerr) RETURN
    END DO
  END DO
END FUNCTION testIsotopes
!------------------------------------------------------------------------------------------------------------------------------

END FUNCTION processIsotopes
!==============================================================================================================================


!==============================================================================================================================
!=== THE ROUTINE isoSelect IS USED TO SWITCH FROM AN ISOTOPE FAMILY TO ANOTHER: ISOTOPES DEPENDENT PARAMETERS ARE UPDATED
!     Single generic "isoSelect" routine, using the predefined index of the family (fast version) or its name (first call).
!==============================================================================================================================
LOGICAL FUNCTION isoSelectByName(iName, lVerbose) RESULT(lerr)
   IMPLICIT NONE
   CHARACTER(LEN=*),  INTENT(IN) :: iName
   LOGICAL, OPTIONAL, INTENT(IN) :: lVerbose
   INTEGER :: iIso
   LOGICAL :: lV
   lV = .FALSE.; IF(PRESENT(lVerbose)) lV = lVerbose
   iIso = strIdx(isotopes(:)%parent, iName)
   lerr = iIso == 0
   IF(lerr) THEN
      niso = 0; ntiso = 0; nzone = 0; nphas = 0; isoCheck=.FALSE.
      CALL msg('no isotope family named "'//TRIM(iName)//'"', ll=lV)
      RETURN
   END IF
   lerr = isoSelectByIndex(iIso, lV)
END FUNCTION isoSelectByName
!==============================================================================================================================
LOGICAL FUNCTION isoSelectByIndex(iIso, lVerbose) RESULT(lerr)
   IMPLICIT NONE
   INTEGER,           INTENT(IN) :: iIso
   LOGICAL, OPTIONAL, INTENT(IN) :: lVerbose
   LOGICAL :: lV
   lv = .FALSE.; IF(PRESENT(lVerbose)) lv = lVerbose
   lerr = .FALSE.
   IF(iIso == ixIso) RETURN                                          !--- Nothing to do if the index is already OK
   lerr = iIso<=0 .OR. iIso>SIZE(isotopes)
   CALL msg('Inconsistent isotopes family index '//TRIM(num2str(iIso))//': should be > 0 and <= '&
          //TRIM(num2str(SIZE(isotopes)))//'"', ll = lerr .AND. lV)
   IF(lerr) RETURN
   ixIso = iIso                                                      !--- Update currently selected family index
   isotope  => isotopes(ixIso)                                       !--- Select corresponding component
   isoKeys  => isotope%keys;     niso     = isotope%niso
   isoName  => isotope%trac;     ntiso    = isotope%ntiso
   isoZone  => isotope%zone;     nzone    = isotope%nzone
   isoPhas  => isotope%phase;    nphas    = isotope%nphas
   itZonIso => isotope%itZonIso; isoCheck = isotope%check
   iqIsoPha => isotope%iqIsoPha
   iqWIsoPha => isotope%iqWIsoPha
END FUNCTION isoSelectByIndex
!==============================================================================================================================


!==============================================================================================================================
!=== ADD THE <key>=<val> PAIR TO THE "ky[(:)]" KEY[S] DESCRIPTOR[S] OR THE <key>=<val(:)> PAIRS TO THE "ky(:)" KEYS DESCRIPTORS
!==============================================================================================================================
SUBROUTINE addKey_s11(key, sval, ky, lOverWrite)
  CHARACTER(LEN=*),  INTENT(IN)    :: key, sval
  TYPE(keys_type),   INTENT(INOUT) :: ky
  LOGICAL, OPTIONAL, INTENT(IN)    :: lOverWrite
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen), ALLOCATABLE :: k(:), v(:)
  INTEGER :: iky, nky
  LOGICAL :: lo
  lo=.TRUE.; IF(PRESENT(lOverWrite)) lo=lOverWrite
  IF(.NOT.ALLOCATED(ky%key)) THEN
    ALLOCATE(ky%key(1)); ky%key(1)=key
    ALLOCATE(ky%val(1)); ky%val(1)=sval
    RETURN
  END IF
  iky = strIdx(ky%key,key)
  IF(iky == 0) THEN
    nky = SIZE(ky%key)
    ALLOCATE(k(nky+1)); k(1:nky) = ky%key; k(nky+1) = key;  ky%key = k
    ALLOCATE(v(nky+1)); v(1:nky) = ky%val; v(nky+1) = sval; ky%val = v
  ELSE IF(lo) THEN
    ky%key(iky) = key; ky%val(iky) = sval
  END IF
END SUBROUTINE addKey_s11
!==============================================================================================================================
SUBROUTINE addKey_i11(key, ival, ky, lOverWrite)
  CHARACTER(LEN=*),  INTENT(IN)    :: key
  INTEGER,           INTENT(IN)    :: ival
  TYPE(keys_type),   INTENT(INOUT) :: ky
  LOGICAL, OPTIONAL, INTENT(IN)    :: lOverWrite
!------------------------------------------------------------------------------------------------------------------------------
  CALL addKey_s11(key, num2str(ival), ky, lOverWrite)
END SUBROUTINE addKey_i11
!==============================================================================================================================
SUBROUTINE addKey_r11(key, rval, ky, lOverWrite)
  CHARACTER(LEN=*),  INTENT(IN)    :: key
  REAL,              INTENT(IN)    :: rval
  TYPE(keys_type),   INTENT(INOUT) :: ky
  LOGICAL, OPTIONAL, INTENT(IN)    :: lOverWrite
!------------------------------------------------------------------------------------------------------------------------------
  CALL addKey_s11(key, num2str(rval), ky, lOverWrite)
END SUBROUTINE addKey_r11
!==============================================================================================================================
SUBROUTINE addKey_l11(key, lval, ky, lOverWrite)
  CHARACTER(LEN=*),  INTENT(IN)    :: key
  LOGICAL,           INTENT(IN)    :: lval
  TYPE(keys_type),   INTENT(INOUT) :: ky
  LOGICAL, OPTIONAL, INTENT(IN)    :: lOverWrite
!------------------------------------------------------------------------------------------------------------------------------
  CALL addKey_s11(key, num2str(lval), ky, lOverWrite)
END SUBROUTINE addKey_l11
!==============================================================================================================================
!==============================================================================================================================
SUBROUTINE addKey_s1m(key, sval, ky, lOverWrite)
  CHARACTER(LEN=*),  INTENT(IN)    :: key, sval
  TYPE(keys_type),   INTENT(INOUT) :: ky(:)
  LOGICAL, OPTIONAL, INTENT(IN)    :: lOverWrite
!------------------------------------------------------------------------------------------------------------------------------
  INTEGER :: itr
  DO itr = 1, SIZE(ky); CALL addKey_s11(key, sval, ky(itr), lOverWrite); END DO
END SUBROUTINE addKey_s1m
!==============================================================================================================================
SUBROUTINE addKey_i1m(key, ival, ky, lOverWrite)
  CHARACTER(LEN=*),  INTENT(IN)    :: key
  INTEGER,           INTENT(IN)    :: ival
  TYPE(keys_type),   INTENT(INOUT) :: ky(:)
  LOGICAL, OPTIONAL, INTENT(IN)    :: lOverWrite
!------------------------------------------------------------------------------------------------------------------------------
  INTEGER :: itr
  DO itr = 1, SIZE(ky); CALL addKey_s11(key, num2str(ival), ky(itr), lOverWrite); END DO
END SUBROUTINE addKey_i1m
!==============================================================================================================================
SUBROUTINE addKey_r1m(key, rval, ky, lOverWrite)
  CHARACTER(LEN=*),  INTENT(IN)    :: key
  REAL,              INTENT(IN)    :: rval
  TYPE(keys_type),   INTENT(INOUT) :: ky(:)
  LOGICAL, OPTIONAL, INTENT(IN)    :: lOverWrite
!------------------------------------------------------------------------------------------------------------------------------
  INTEGER :: itr
  DO itr = 1, SIZE(ky); CALL addKey_s11(key, num2str(rval), ky(itr), lOverWrite); END DO
END SUBROUTINE addKey_r1m
!==============================================================================================================================
SUBROUTINE addKey_l1m(key, lval, ky, lOverWrite)
  CHARACTER(LEN=*),  INTENT(IN)    :: key
  LOGICAL,           INTENT(IN)    :: lval
  TYPE(keys_type),   INTENT(INOUT) :: ky(:)
  LOGICAL, OPTIONAL, INTENT(IN)    :: lOverWrite
!------------------------------------------------------------------------------------------------------------------------------
  INTEGER :: itr
  DO itr = 1, SIZE(ky); CALL addKey_s11(key, num2str(lval), ky(itr), lOverWrite); END DO
END SUBROUTINE addKey_l1m
!==============================================================================================================================
!==============================================================================================================================
SUBROUTINE addKey_smm(key, sval, ky, lOverWrite)
  CHARACTER(LEN=*),  INTENT(IN)    :: key, sval(:)
  TYPE(keys_type),   INTENT(INOUT) :: ky(:)
  LOGICAL, OPTIONAL, INTENT(IN)    :: lOverWrite
!------------------------------------------------------------------------------------------------------------------------------
  INTEGER :: itr
  DO itr = 1, SIZE(ky); CALL addKey_s11(key, sval(itr), ky(itr), lOverWrite); END DO
END SUBROUTINE addKey_smm
!==============================================================================================================================
SUBROUTINE addKey_imm(key, ival, ky, lOverWrite)
  CHARACTER(LEN=*),  INTENT(IN)    :: key
  INTEGER,           INTENT(IN)    :: ival(:)
  TYPE(keys_type),   INTENT(INOUT) :: ky(:)
  LOGICAL, OPTIONAL, INTENT(IN)    :: lOverWrite
!------------------------------------------------------------------------------------------------------------------------------
  INTEGER :: itr
  DO itr = 1, SIZE(ky); CALL addKey_s11(key, num2str(ival(itr)), ky(itr), lOverWrite); END DO
END SUBROUTINE addKey_imm
!==============================================================================================================================
SUBROUTINE addKey_rmm(key, rval, ky, lOverWrite)
  CHARACTER(LEN=*),  INTENT(IN)    :: key
  REAL,              INTENT(IN)    :: rval(:)
  TYPE(keys_type),   INTENT(INOUT) :: ky(:)
  LOGICAL, OPTIONAL, INTENT(IN)    :: lOverWrite
!------------------------------------------------------------------------------------------------------------------------------
  INTEGER :: itr
  DO itr = 1, SIZE(ky); CALL addKey_s11(key, num2str(rval(itr)), ky(itr), lOverWrite); END DO
END SUBROUTINE addKey_rmm
!==============================================================================================================================
SUBROUTINE addKey_lmm(key, lval, ky, lOverWrite)
  CHARACTER(LEN=*),  INTENT(IN)    :: key
  LOGICAL,           INTENT(IN)    :: lval(:)
  TYPE(keys_type),   INTENT(INOUT) :: ky(:)
  LOGICAL, OPTIONAL, INTENT(IN)    :: lOverWrite
!------------------------------------------------------------------------------------------------------------------------------
  INTEGER :: itr
  DO itr = 1, SIZE(ky); CALL addKey_s11(key, num2str(lval(itr)), ky(itr), lOverWrite); END DO
END SUBROUTINE addKey_lmm
!==============================================================================================================================


!==============================================================================================================================
!=== OVERWRITE THE KEYS OF THE TRACER NAMED "tr0" WITH THE VALUES FOUND IN THE *.def FILES, IF ANY. ===========================
!==============================================================================================================================
!SUBROUTINE addKeysFromDef(t, tr0)
!  TYPE(trac_type), ALLOCATABLE, INTENT(INOUT) :: t(:)
!  CHARACTER(LEN=*),             INTENT(IN)    :: tr0
!------------------------------------------------------------------------------------------------------------------------------
!  CHARACTER(LEN=maxlen) :: val
!  INTEGER               :: ik, jd
!  jd = strIdx(t%name, tr0)
!  IF(jd == 0) RETURN
!  DO ik = 1, SIZE(t(jd)%keys%key)
!    CALL getin_p(t(jd)%keys%key(ik), val, '*none*')
!    IF(val /= '*none*') CALL addKey(t(jd)%keys%key(ik), val, t(jd)%keys, .TRUE.)
!  END DO
!END SUBROUTINE addKeysFromDef
!==============================================================================================================================


!==============================================================================================================================
!=== REMOVE THE KEYS NAMED "keyn(:)" FROM EITHER THE "itr"th OR ALL THE KEYS DESCRIPTORS OF "ky(:)" ===========================
!==============================================================================================================================
SUBROUTINE delKey_1(itr, keyn, ky)
  INTEGER,          INTENT(IN)    :: itr
  CHARACTER(LEN=*), INTENT(IN)    :: keyn(:)
  TYPE(trac_type),  INTENT(INOUT) :: ky(:)
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen), ALLOCATABLE :: k(:), v(:)
  LOGICAL,               ALLOCATABLE :: ll(:)
  INTEGER :: iky
  IF(itr<=0 .OR. itr>SIZE(ky, DIM=1)) RETURN                          !--- Index is out of range
  ll = [( ALL(keyn/=ky(itr)%keys%key(iky)), iky=1, SIZE(ky(itr)%keys%key) )]
  k = PACK(ky(itr)%keys%key, MASK=ll); CALL MOVE_ALLOC(FROM=k, TO=ky(itr)%keys%key)
  v = PACK(ky(itr)%keys%val, MASK=ll); CALL MOVE_ALLOC(FROM=v, TO=ky(itr)%keys%val)
END SUBROUTINE delKey_1
!==============================================================================================================================
SUBROUTINE delKey(keyn, ky)
  CHARACTER(LEN=*), INTENT(IN)    :: keyn(:)
  TYPE(trac_type),  INTENT(INOUT) :: ky(:)
!------------------------------------------------------------------------------------------------------------------------------
  INTEGER :: iky
  DO iky = 1, SIZE(ky); CALL delKey_1(iky, keyn, ky); END DO
END SUBROUTINE delKey
!==============================================================================================================================


!==============================================================================================================================
!===   INTERNAL FUNCTION: GET THE VALUE OF A KEY FOR "itr"th TRACER FROM A "keys_type" DERIVED TYPE AND RETURN THE RESULT   ===
!===   IF keyn CONTAINS SEVERAL ELEMENTS, TRY WITH EACH ELEMENT ONE AFTER THE OTHER                                         ===
!==============================================================================================================================
CHARACTER(LEN=maxlen) FUNCTION fgetKeyIdx(itr, keyn, ky, lerr) RESULT(val)
  INTEGER,                    INTENT(IN)  :: itr
  CHARACTER(LEN=*),           INTENT(IN)  :: keyn(:)
  TYPE(keys_type),            INTENT(IN)  :: ky(:)
  LOGICAL,          OPTIONAL, INTENT(OUT) :: lerr
!------------------------------------------------------------------------------------------------------------------------------
  INTEGER :: ik
  LOGICAL :: ler
  ler = .TRUE.
  DO ik = 1, SIZE(keyn)
    CALL getKeyIdx(keyn(ik)); IF(.NOT.ler) EXIT
  END DO
  IF(PRESENT(lerr)) lerr = ler

CONTAINS

SUBROUTINE getKeyIdx(keyn)
  CHARACTER(LEN=*), INTENT(IN) :: keyn
!------------------------------------------------------------------------------------------------------------------------------
  INTEGER :: iky
  iky = 0; val = ''
  ler = itr <= 0 .OR. itr > SIZE(ky); IF(ler) RETURN
  iky = strIdx(ky(itr)%key(:), keyn)
  ler = iky == 0;                     IF(ler) RETURN
  val = ky(itr)%val(iky)
END SUBROUTINE getKeyIdx

END FUNCTION fgetKeyIdx
!==============================================================================================================================


!==============================================================================================================================
!===                                          GET KEYS VALUES FROM TRACERS INDICES                                          ===
!==============================================================================================================================
!=== TRY TO GET THE KEY NAMED "key" FOR THE "itr"th TRACER IN:                                                              ===
!===  * ARGUMENT "ky" DATABASE IF SPECIFIED ; OTHERWISE:                                                                    ===
!===  * IN INTERNAL TRACERS DATABASE "tracers(:)%keys" (IF UNFOUND, IN INTERNAL ISOTOPES DATABASE "isotope%keys(:)")        ===
!=== THE RETURNED VALUE (STRING, AN INTEGER, A REAL OR A LOGICAL) CAN BE EITHER:                                            ===
!===  * A SCALAR                                                                                                            ===
!===  * A VECTOR, WHICH IS THE FOUND VALUE PARSED WITH THE SEPARATOR ","                                                    ===
!===                                                                                                                        ===
!=== SYNTAX:       lerr = getKeyByIndex_{sirl}{1m}{1m}1(keyn[(:)], val[(:)], itr[, ky(:)]          [, def][, lDisp])        ===
!==============================================================================================================================
!=== IF "itr" IS NOT PRESENT, VALUES FOR ALL THE TRACERS OF THE SELECTED DATABASE ARE STORED IN THE VECTOR "val(:)"         ===
!=== THE NAME OF THE TRACERS FOUND IN THE EFFECTIVELY USED DATABASE CAN BE RETURNED OPTIONALLY IN "nam(:)"                  ===
!=== SYNTAX        lerr = getKeyByIndex_{sirl}{1m}mm   (keyn[(:)], val (:)      [, ky(:)][, nam(:)][, def][, lDisp])        ===
!==============================================================================================================================
LOGICAL FUNCTION getKeyByIndex_s111(keyn, val, itr, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),          INTENT(IN)  :: keyn
  CHARACTER(LEN=maxlen),     INTENT(OUT) :: val
  INTEGER,                   INTENT(IN)  :: itr
  TYPE(keys_type), OPTIONAL, INTENT(IN)  :: ky(:)
  CHARACTER(LEN=*),OPTIONAL, INTENT(IN)  :: def
  LOGICAL,         OPTIONAL, INTENT(IN)  :: lDisp
  lerr = getKeyByIndex_sm11([keyn], val, itr, ky, def, lDisp)
END FUNCTION getKeyByIndex_s111
!==============================================================================================================================
LOGICAL FUNCTION getKeyByIndex_i111(keyn, val, itr, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),          INTENT(IN)  :: keyn
  INTEGER,                   INTENT(OUT) :: val
  INTEGER,                   INTENT(IN)  :: itr
  TYPE(keys_type), OPTIONAL, INTENT(IN)  :: ky(:)
  INTEGER,         OPTIONAL, INTENT(IN)  :: def
  LOGICAL,         OPTIONAL, INTENT(IN)  :: lDisp
  lerr = getKeyByIndex_im11([keyn], val, itr, ky, def, lDisp)
END FUNCTION getKeyByIndex_i111
!==============================================================================================================================
LOGICAL FUNCTION getKeyByIndex_r111(keyn, val, itr, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),          INTENT(IN)  :: keyn
  REAL   ,                   INTENT(OUT) :: val
  INTEGER,                   INTENT(IN)  :: itr
  TYPE(keys_type), OPTIONAL, INTENT(IN)  :: ky(:)
  REAL,            OPTIONAL, INTENT(IN)  :: def
  LOGICAL,         OPTIONAL, INTENT(IN)  :: lDisp
  lerr = getKeyByIndex_rm11([keyn], val, itr, ky, def, lDisp)
END FUNCTION getKeyByIndex_r111
!==============================================================================================================================
LOGICAL FUNCTION getKeyByIndex_l111(keyn, val, itr, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),          INTENT(IN)  :: keyn
  LOGICAL,                   INTENT(OUT) :: val
  INTEGER,                   INTENT(IN)  :: itr
  TYPE(keys_type), OPTIONAL, INTENT(IN)  :: ky(:)
  LOGICAL,         OPTIONAL, INTENT(IN)  :: def
  LOGICAL,         OPTIONAL, INTENT(IN)  :: lDisp
  lerr = getKeyByIndex_lm11([keyn], val, itr, ky, def, lDisp)
END FUNCTION getKeyByIndex_l111
!==============================================================================================================================
!==============================================================================================================================
LOGICAL FUNCTION getKeyByIndex_sm11(keyn, val, itr, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),          INTENT(IN)  :: keyn(:)
  CHARACTER(LEN=maxlen),     INTENT(OUT) :: val
  INTEGER,                   INTENT(IN)  :: itr
  TYPE(keys_type), OPTIONAL, INTENT(IN)  :: ky(:)
  CHARACTER(LEN=*),OPTIONAL, INTENT(IN)  :: def
  LOGICAL,         OPTIONAL, INTENT(IN)  :: lDisp
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen) :: s
  LOGICAL :: lD
  lD = .TRUE.; IF(PRESENT(lDisp)) lD = lDisp
  s = 'key "'//TRIM(strStack(keyn, '/'))//'" for tracer nr. '//TRIM(num2str(itr))
  lerr = .TRUE.
  IF(lerr .AND. PRESENT(ky))         val = fgetKey(ky)                                   !--- "ky"
  IF(lerr .AND. ALLOCATED (tracers)) val = fgetKey(tracers(:)%keys)                      !--- "tracers"
  IF(lerr .AND. ASSOCIATED(isotope)) val = fgetKey(isotope%keys(:))                      !--- "isotope"
  IF(lerr .AND. PRESENT(def)) THEN
     val = def; lerr = .NOT.PRESENT(def)
     CALL msg('Using defaut value of '//TRIM(s)//': '//TRIM(def), modname, lD)
  END IF
  CALL msg('No '//TRIM(s)//' or out of bounds index', modname, lD .AND. lerr)

CONTAINS

CHARACTER(LEN=maxlen) FUNCTION fgetKey(ky) RESULT(val)
  TYPE(keys_type),  INTENT(IN)  :: ky(:)
  lerr = SIZE(ky) == 0; IF(lerr) RETURN
  val = fgetKeyIdx(itr, keyn(:), ky, lerr)
END FUNCTION fgetKey

END FUNCTION getKeyByIndex_sm11
!==============================================================================================================================
LOGICAL FUNCTION getKeyByIndex_im11(keyn, val, itr, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),          INTENT(IN)  :: keyn(:)
  INTEGER,                   INTENT(OUT) :: val
  INTEGER,                   INTENT(IN)  :: itr
  TYPE(keys_type), OPTIONAL, INTENT(IN)  :: ky(:)
  INTEGER,         OPTIONAL, INTENT(IN)  :: def
  LOGICAL,         OPTIONAL, INTENT(IN)  :: lDisp
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen) :: sval, s
  IF(     PRESENT(def)) lerr = getKeyByIndex_sm11(keyn, sval, itr, ky, num2str(def), lDisp)
  IF(.NOT.PRESENT(def)) lerr = getKeyByIndex_sm11(keyn, sval, itr, ky, lDisp=lDisp)
  IF(lerr) RETURN
  val = str2int(sval)
  lerr = val == -HUGE(1)
  s = 'key "'//TRIM(strStack(keyn, '/'))//'" of tracer nr. '//TRIM(num2str(itr))//' is not'
  CALL msg(TRIM(s)//' an integer: '//TRIM(sval), modname, lerr)
END FUNCTION getKeyByIndex_im11
!==============================================================================================================================
LOGICAL FUNCTION getKeyByIndex_rm11(keyn, val, itr, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),          INTENT(IN)  :: keyn(:)
  REAL   ,                   INTENT(OUT) :: val
  INTEGER,                   INTENT(IN)  :: itr
  TYPE(keys_type), OPTIONAL, INTENT(IN)  :: ky(:)
  REAL,            OPTIONAL, INTENT(IN)  :: def
  LOGICAL,         OPTIONAL, INTENT(IN)  :: lDisp
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen) :: sval, s
  IF(     PRESENT(def)) lerr = getKeyByIndex_sm11(keyn, sval, itr, ky, num2str(def), lDisp)
  IF(.NOT.PRESENT(def)) lerr = getKeyByIndex_sm11(keyn, sval, itr, ky, lDisp=lDisp)
  IF(lerr) RETURN
  val = str2real(sval)
  lerr = val == -HUGE(1.)
  s = 'key "'//TRIM(strStack(keyn, '/'))//'" of tracer nr. '//TRIM(num2str(itr))//' is not'
  CALL msg(TRIM(s)//' a real: '//TRIM(sval), modname, lerr)
END FUNCTION getKeyByIndex_rm11
!==============================================================================================================================
LOGICAL FUNCTION getKeyByIndex_lm11(keyn, val, itr, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),          INTENT(IN)  :: keyn(:)
  LOGICAL,                   INTENT(OUT) :: val
  INTEGER,                   INTENT(IN)  :: itr
  TYPE(keys_type), OPTIONAL, INTENT(IN)  :: ky(:)
  LOGICAL,         OPTIONAL, INTENT(IN)  :: def
  LOGICAL,         OPTIONAL, INTENT(IN)  :: lDisp
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen) :: sval, s
  INTEGER               :: ival
  IF(     PRESENT(def)) lerr = getKeyByIndex_sm11(keyn, sval, itr, ky, num2str(def), lDisp)
  IF(.NOT.PRESENT(def)) lerr = getKeyByIndex_sm11(keyn, sval, itr, ky, lDisp=lDisp)
  IF(lerr) RETURN
  ival = str2bool(sval)
  lerr = ival == -1
  s = 'key "'//TRIM(strStack(keyn, '/'))//'" of tracer nr. '//TRIM(num2str(itr))//' is not'
  CALL msg(TRIM(s)//' a boolean: '//TRIM(sval), modname, lerr)
  IF(.NOT.lerr) val = ival == 1
END FUNCTION getKeyByIndex_lm11
!==============================================================================================================================
!==============================================================================================================================
LOGICAL FUNCTION getKeyByIndex_s1m1(keyn, val, itr, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),                   INTENT(IN)  :: keyn
  CHARACTER(LEN=maxlen), ALLOCATABLE, INTENT(OUT) :: val(:)
  INTEGER,                            INTENT(IN)  :: itr
  TYPE(keys_type),          OPTIONAL, INTENT(IN)  :: ky(:)
  CHARACTER(LEN=*),         OPTIONAL, INTENT(IN)  :: def
  LOGICAL,                  OPTIONAL, INTENT(IN)  :: lDisp
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen)              :: sval
  lerr = getKeyByIndex_sm11([keyn], sval, itr, ky, def, lDisp); IF(lerr) RETURN
  lerr = strParse(sval, ',', val)
  CALL msg('can''t parse '//TRIM(sval), modname, lerr); IF(lerr) RETURN
END FUNCTION getKeyByIndex_s1m1
!==============================================================================================================================
LOGICAL FUNCTION getKeyByIndex_i1m1(keyn, val, itr, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),          INTENT(IN)  :: keyn
  INTEGER,      ALLOCATABLE, INTENT(OUT) :: val(:)
  INTEGER,                   INTENT(IN)  :: itr
  TYPE(keys_type), OPTIONAL, INTENT(IN)  :: ky(:)
  INTEGER,         OPTIONAL, INTENT(IN)  :: def
  LOGICAL,         OPTIONAL, INTENT(IN)  :: lDisp
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen)              :: sval, s
  CHARACTER(LEN=maxlen), ALLOCATABLE :: svals(:)
  IF(     PRESENT(def)) lerr = getKeyByIndex_sm11([keyn], sval, itr, ky, num2str(def), lDisp)
  IF(.NOT.PRESENT(def)) lerr = getKeyByIndex_sm11([keyn], sval, itr, ky, lDisp=lDisp)
  IF(lerr) RETURN
  lerr = strParse(sval, ',', svals)
  CALL msg('can''t parse '//TRIM(sval), modname, lerr); IF(lerr) RETURN
  val = str2int(svals)
  lerr = ANY(val == -HUGE(1))
  s = 'key "'//TRIM(keyn)//'" of tracer nr. '//TRIM(num2str(itr))//' is not'
  CALL msg(TRIM(s)//' an integer: '//TRIM(sval), modname, lerr)
END FUNCTION getKeyByIndex_i1m1
!==============================================================================================================================
LOGICAL FUNCTION getKeyByIndex_r1m1(keyn, val, itr, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),          INTENT(IN)  :: keyn
  REAL,         ALLOCATABLE, INTENT(OUT) :: val(:)
  INTEGER,                   INTENT(IN)  :: itr
  TYPE(keys_type), OPTIONAL, INTENT(IN)  :: ky(:)
  REAL,            OPTIONAL, INTENT(IN)  :: def
  LOGICAL,         OPTIONAL, INTENT(IN)  :: lDisp
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen)              :: sval, s
  CHARACTER(LEN=maxlen), ALLOCATABLE :: svals(:)
  IF(     PRESENT(def)) lerr = getKeyByIndex_sm11([keyn], sval, itr, ky, num2str(def), lDisp)
  IF(.NOT.PRESENT(def)) lerr = getKeyByIndex_sm11([keyn], sval, itr, ky, lDisp=lDisp)
  lerr = strParse(sval, ',', svals)
  CALL msg('can''t parse '//TRIM(sval), modname, lerr); IF(lerr) RETURN
  val = str2real(svals)
  lerr = ANY(val == -HUGE(1.))
  s = 'key "'//TRIM(keyn)//'" of tracer nr. '//TRIM(num2str(itr))//' is not'
  CALL msg(TRIM(s)//' a real: '//TRIM(sval), modname, lerr)
END FUNCTION getKeyByIndex_r1m1
!==============================================================================================================================
LOGICAL FUNCTION getKeyByIndex_l1m1(keyn, val, itr, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),          INTENT(IN)  :: keyn
  LOGICAL,      ALLOCATABLE, INTENT(OUT) :: val(:)
  INTEGER,                   INTENT(IN)  :: itr
  TYPE(keys_type), OPTIONAL, INTENT(IN)  :: ky(:)
  LOGICAL,         OPTIONAL, INTENT(IN)  :: def
  LOGICAL,         OPTIONAL, INTENT(IN)  :: lDisp
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen)              :: sval, s
  CHARACTER(LEN=maxlen), ALLOCATABLE :: svals(:)
  INTEGER,               ALLOCATABLE :: ivals(:)
  IF(     PRESENT(def)) lerr = getKeyByIndex_sm11([keyn], sval, itr, ky, num2str(def), lDisp)
  IF(.NOT.PRESENT(def)) lerr = getKeyByIndex_sm11([keyn], sval, itr, ky, lDisp=lDisp)
  lerr = strParse(sval, ',', svals)
  CALL msg('can''t parse '//TRIM(sval), modname, lerr); IF(lerr) RETURN
  ivals = str2bool(svals)
  lerr = ANY(ivals == -1)
  s = 'key "'//TRIM(keyn)//'" of tracer nr. '//TRIM(num2str(itr))//' is not'
  CALL msg(TRIM(s)//' a boolean: '//TRIM(sval), modname, lerr)
  IF(.NOT.lerr) val = ivals == 1
END FUNCTION getKeyByIndex_l1m1
!==============================================================================================================================
!==============================================================================================================================
LOGICAL FUNCTION getKeyByIndex_smm1(keyn, val, itr, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),                   INTENT(IN)  :: keyn(:)
  CHARACTER(LEN=maxlen), ALLOCATABLE, INTENT(OUT) :: val(:)
  INTEGER,                            INTENT(IN)  :: itr
  TYPE(keys_type),          OPTIONAL, INTENT(IN)  :: ky(:)
  CHARACTER(LEN=*),         OPTIONAL, INTENT(IN)  :: def
  LOGICAL,                  OPTIONAL, INTENT(IN)  :: lDisp
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen) :: sval
  lerr = getKeyByIndex_sm11(keyn, sval, itr, ky, def, lDisp); IF(lerr) RETURN
  lerr = strParse(sval, ',', val)
  CALL msg('can''t parse '//TRIM(sval), modname, lerr); IF(lerr) RETURN
END FUNCTION getKeyByIndex_smm1
!==============================================================================================================================
LOGICAL FUNCTION getKeyByIndex_imm1(keyn, val, itr, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),          INTENT(IN)  :: keyn(:)
  INTEGER,      ALLOCATABLE, INTENT(OUT) :: val(:)
  INTEGER,                   INTENT(IN)  :: itr
  TYPE(keys_type), OPTIONAL, INTENT(IN)  :: ky(:)
  INTEGER,         OPTIONAL, INTENT(IN)  :: def
  LOGICAL,         OPTIONAL, INTENT(IN)  :: lDisp
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen)              :: sval, s
  CHARACTER(LEN=maxlen), ALLOCATABLE :: svals(:)
  IF(     PRESENT(def)) lerr = getKeyByIndex_sm11(keyn, sval, itr, ky, num2str(def), lDisp)
  IF(.NOT.PRESENT(def)) lerr = getKeyByIndex_sm11(keyn, sval, itr, ky, lDisp=lDisp)
  IF(lerr) RETURN
  lerr = strParse(sval, ',', svals)
  CALL msg('can''t parse '//TRIM(sval), modname, lerr); IF(lerr) RETURN
  val = str2int(svals)
  lerr = ANY(val == -HUGE(1))
  s = 'key "'//TRIM(strStack(keyn, '/'))//'" of tracer nr. '//TRIM(num2str(itr))//' is not'
  CALL msg(TRIM(s)//' an integer: '//TRIM(sval), modname, lerr)
END FUNCTION getKeyByIndex_imm1
!==============================================================================================================================
LOGICAL FUNCTION getKeyByIndex_rmm1(keyn, val, itr, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),          INTENT(IN)  :: keyn(:)
  REAL,         ALLOCATABLE, INTENT(OUT) :: val(:)
  INTEGER,                   INTENT(IN)  :: itr
  TYPE(keys_type), OPTIONAL, INTENT(IN)  :: ky(:)
  REAL,            OPTIONAL, INTENT(IN)  :: def
  LOGICAL,         OPTIONAL, INTENT(IN)  :: lDisp
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen)              :: sval, s
  CHARACTER(LEN=maxlen), ALLOCATABLE :: svals(:)
  IF(     PRESENT(def)) lerr = getKeyByIndex_sm11(keyn, sval, itr, ky, num2str(def), lDisp)
  IF(.NOT.PRESENT(def)) lerr = getKeyByIndex_sm11(keyn, sval, itr, ky, lDisp=lDisp)
  IF(lerr) RETURN
  lerr = strParse(sval, ',', svals)
  CALL msg('can''t parse '//TRIM(sval), modname, lerr); IF(lerr) RETURN
  val = str2real(svals)
  lerr = ANY(val == -HUGE(1.))
  s = 'key "'//TRIM(strStack(keyn, '/'))//'" of tracer nr. '//TRIM(num2str(itr))//' is not'
  CALL msg(TRIM(s)//' a real: '//TRIM(sval), modname, lerr)
END FUNCTION getKeyByIndex_rmm1
!==============================================================================================================================
LOGICAL FUNCTION getKeyByIndex_lmm1(keyn, val, itr, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),          INTENT(IN)  :: keyn(:)
  LOGICAL,      ALLOCATABLE, INTENT(OUT) :: val(:)
  INTEGER,                   INTENT(IN)  :: itr
  TYPE(keys_type), OPTIONAL, INTENT(IN)  :: ky(:)
  LOGICAL,         OPTIONAL, INTENT(IN)  :: def
  LOGICAL,         OPTIONAL, INTENT(IN)  :: lDisp
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen)              :: sval, s
  CHARACTER(LEN=maxlen), ALLOCATABLE :: svals(:)
  INTEGER,               ALLOCATABLE :: ivals(:)
  IF(     PRESENT(def)) lerr = getKeyByIndex_sm11(keyn, sval, itr, ky, num2str(def), lDisp)
  IF(.NOT.PRESENT(def)) lerr = getKeyByIndex_sm11(keyn, sval, itr, ky, lDisp=lDisp)
  IF(lerr) RETURN
  lerr = strParse(sval, ',', svals)
  CALL msg('can''t parse '//TRIM(sval), modname, lerr); IF(lerr) RETURN
  ivals = str2bool(svals)
  lerr = ANY(ivals == -1)
  s = 'key "'//TRIM(strStack(keyn, '/'))//'" of tracer nr. '//TRIM(num2str(itr))//' is not'
  CALL msg(TRIM(s)//' a boolean: '//TRIM(sval), modname, lerr)
  IF(.NOT.lerr) val = ivals == 1
END FUNCTION getKeyByIndex_lmm1
!==============================================================================================================================
!==============================================================================================================================
LOGICAL FUNCTION getKeyByIndex_s1mm(keyn, val, ky, nam, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),                             INTENT(IN)  :: keyn
  CHARACTER(LEN=maxlen),           ALLOCATABLE, INTENT(OUT) :: val(:)
  TYPE(keys_type),       OPTIONAL,              INTENT(IN)  ::  ky(:)
  CHARACTER(LEN=maxlen), OPTIONAL, ALLOCATABLE, INTENT(OUT) :: nam(:)
  CHARACTER(LEN=*),      OPTIONAL,              INTENT(IN)  :: def
  LOGICAL,               OPTIONAL,              INTENT(IN)  :: lDisp
  lerr = getKeyByIndex_smmm([keyn], val, ky, nam, def, lDisp)
END FUNCTION getKeyByIndex_s1mm
!==============================================================================================================================
LOGICAL FUNCTION getKeyByIndex_i1mm(keyn, val, ky, nam, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),                             INTENT(IN)  :: keyn
  INTEGER,                         ALLOCATABLE, INTENT(OUT) :: val(:)
  TYPE(keys_type),       OPTIONAL,              INTENT(IN)  ::  ky(:)
  CHARACTER(LEN=maxlen), OPTIONAL, ALLOCATABLE, INTENT(OUT) :: nam(:)
  INTEGER,               OPTIONAL,              INTENT(IN)  :: def
  LOGICAL,               OPTIONAL,              INTENT(IN)  :: lDisp
  lerr = getKeyByIndex_immm([keyn], val, ky, nam, def, lDisp)
END FUNCTION getKeyByIndex_i1mm
!==============================================================================================================================
LOGICAL FUNCTION getKeyByIndex_r1mm(keyn, val, ky, nam, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),                             INTENT(IN)  :: keyn
  REAL,                            ALLOCATABLE, INTENT(OUT) :: val(:)
  TYPE(keys_type),       OPTIONAL,              INTENT(IN)  ::  ky(:)
  CHARACTER(LEN=maxlen), OPTIONAL, ALLOCATABLE, INTENT(OUT) :: nam(:)
  REAL,                  OPTIONAL,              INTENT(IN)  :: def
  LOGICAL,               OPTIONAL,              INTENT(IN)  :: lDisp
  lerr = getKeyByIndex_rmmm([keyn], val, ky, nam, def, lDisp)
END FUNCTION getKeyByIndex_r1mm
!==============================================================================================================================
LOGICAL FUNCTION getKeyByIndex_l1mm(keyn, val, ky, nam, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),                             INTENT(IN)  :: keyn
  LOGICAL,                         ALLOCATABLE, INTENT(OUT) :: val(:)
  TYPE(keys_type),       OPTIONAL,              INTENT(IN)  ::  ky(:)
  CHARACTER(LEN=maxlen), OPTIONAL, ALLOCATABLE, INTENT(OUT) :: nam(:)
  LOGICAL,               OPTIONAL,              INTENT(IN)  :: def
  LOGICAL,               OPTIONAL,              INTENT(IN)  :: lDisp
  lerr = getKeyByIndex_lmmm([keyn], val, ky, nam, def, lDisp)
END FUNCTION getKeyByIndex_l1mm
!==============================================================================================================================
!==============================================================================================================================
LOGICAL FUNCTION getKeyByIndex_smmm(keyn, val, ky, nam, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),                             INTENT(IN)  :: keyn(:)
  CHARACTER(LEN=maxlen),           ALLOCATABLE, INTENT(OUT) ::  val(:)
  TYPE(keys_type),       OPTIONAL,              INTENT(IN)  ::   ky(:)
  CHARACTER(LEN=maxlen), OPTIONAL, ALLOCATABLE, INTENT(OUT) ::  nam(:)
  CHARACTER(LEN=*),      OPTIONAL,              INTENT(IN)  ::  def
  LOGICAL,               OPTIONAL,              INTENT(IN)  :: lDisp
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen) :: s
  CHARACTER(LEN=maxlen), ALLOCATABLE :: tname(:)
  INTEGER :: iq, nq(3), k
  LOGICAL :: lD, l(3)
  lD = .TRUE.; IF(PRESENT(lDisp)) lD = lDisp
  s = 'key "'//TRIM(strStack(keyn, '/'))//'"'
  lerr = .TRUE.
  IF(PRESENT(ky)) THEN;                 val = fgetKey(ky)                                !--- "ky"
  ELSE IF(ALLOCATED(tracers)) THEN;     val = fgetKey(tracers(:)%keys)                   !--- "tracers"
     IF(lerr .AND. ASSOCIATED(isotope)) val = fgetKey(isotope%keys(:))                   !--- "isotope"
  END IF
  IF(.NOT.lerr) THEN; IF(PRESENT(nam)) nam = tname; RETURN; END IF
  IF(.NOT.PRESENT(def)) THEN; CALL msg('No '//TRIM(s)//' found', modname, lD); RETURN; END IF

  !--- DEFAULT VALUE
  l = [PRESENT(ky), ALLOCATED(tracers), ASSOCIATED(isotope)]; nq(:) = 0
  IF(l(1)) nq(1) = SIZE(ky)
  IF(l(2)) nq(2) = SIZE(tracers)
  IF(l(3)) nq(3) = SIZE(isotope%keys)
  DO k = 1, 3; IF(l(k) .AND. nq(k) /= 0) THEN; val = [(def, iq = 1, nq(k))]; EXIT; END IF; END DO
  lerr = k == 4
  CALL msg('Using defaut value for '//TRIM(s)//': '//TRIM(def), modname, lD .AND. .NOT.lerr)
  CALL msg('No '//TRIM(s), modname, lD .AND. lerr)

CONTAINS

FUNCTION fgetKey(ky) RESULT(val)
  CHARACTER(LEN=maxlen), ALLOCATABLE :: val(:)
  TYPE(keys_type),       INTENT(IN)  :: ky(:)
  LOGICAL :: ler(SIZE(ky))
  INTEGER :: iq
  lerr = SIZE(ky) == 0; IF(lerr) RETURN
  tname = ky%name
  val = [(fgetKeyIdx(iq, keyn(:), ky, ler(iq)), iq = 1, SIZE(ky))]
  lerr = ANY(ler)
END FUNCTION fgetKey

END FUNCTION getKeyByIndex_smmm
!==============================================================================================================================
LOGICAL FUNCTION getKeyByIndex_immm(keyn, val, ky, nam, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),                             INTENT(IN)  :: keyn(:)
  INTEGER,                         ALLOCATABLE, INTENT(OUT) ::  val(:)
  TYPE(keys_type),       OPTIONAL,              INTENT(IN)  ::   ky(:)
  CHARACTER(LEN=maxlen), OPTIONAL, ALLOCATABLE, INTENT(OUT) ::  nam(:)
  INTEGER,               OPTIONAL,              INTENT(IN)  ::  def
  LOGICAL,               OPTIONAL,              INTENT(IN)  :: lDisp
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen) :: s
  CHARACTER(LEN=maxlen), ALLOCATABLE :: svals(:), tname(:)
  LOGICAL,               ALLOCATABLE ::    ll(:)
  IF(     PRESENT(def)) lerr = getKeyByIndex_smmm(keyn, svals, ky, tname, num2str(def), lDisp)
  IF(.NOT.PRESENT(def)) lerr = getKeyByIndex_smmm(keyn, svals, ky, tname, lDisp=lDisp)
  IF(lerr) RETURN
  val = str2int(svals)
  ll = val == -HUGE(1)
  lerr = ANY(ll); IF(.NOT.lerr) THEN; IF(PRESENT(nam)) nam = tname; RETURN; END IF
  s = 'key "'//TRIM(strStack(keyn, '/'))//'" of '//TRIM(strStack(PACK(tname, MASK=ll)))//' is not'
  CALL msg(TRIM(s)//' an integer: '//TRIM(strStack(svals, MASK=ll)), modname, lerr)
  IF(.NOT.lerr .AND. PRESENT(nam)) nam = tname
END FUNCTION getKeyByIndex_immm
!==============================================================================================================================
LOGICAL FUNCTION getKeyByIndex_rmmm(keyn, val, ky, nam, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),                             INTENT(IN)  :: keyn(:)
  REAL,                            ALLOCATABLE, INTENT(OUT) ::  val(:)
  TYPE(keys_type),       OPTIONAL,              INTENT(IN)  ::   ky(:)
  CHARACTER(LEN=maxlen), OPTIONAL, ALLOCATABLE, INTENT(OUT) ::  nam(:)
  REAL,                  OPTIONAL,              INTENT(IN)  ::  def
  LOGICAL,               OPTIONAL,              INTENT(IN)  :: lDisp
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen) :: s
  CHARACTER(LEN=maxlen), ALLOCATABLE :: svals(:), tname(:)
  LOGICAL,               ALLOCATABLE ::    ll(:)
  IF(     PRESENT(def)) lerr = getKeyByIndex_smmm(keyn, svals, ky, tname, num2str(def), lDisp)
  IF(.NOT.PRESENT(def)) lerr = getKeyByIndex_smmm(keyn, svals, ky, tname, lDisp=lDisp)
  IF(lerr) RETURN
  val = str2real(svals)
  ll = val == -HUGE(1.)
  lerr = ANY(ll); IF(.NOT.lerr) THEN; IF(PRESENT(nam)) nam = tname; RETURN; END IF
  s = 'key "'//TRIM(strStack(keyn, '/'))//'" of '//TRIM(strStack(PACK(tname, MASK=ll)))//' is not a'
  CALL msg(TRIM(s)//' a real: '//TRIM(strStack(svals, MASK=ll)), modname)
END FUNCTION getKeyByIndex_rmmm
!==============================================================================================================================
LOGICAL FUNCTION getKeyByIndex_lmmm(keyn, val, ky, nam, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),                             INTENT(IN)  :: keyn(:)
  LOGICAL,                         ALLOCATABLE, INTENT(OUT) ::  val(:)
  TYPE(keys_type),       OPTIONAL,              INTENT(IN)  ::   ky(:)
  CHARACTER(LEN=maxlen), OPTIONAL, ALLOCATABLE, INTENT(OUT) ::  nam(:)
  LOGICAL,               OPTIONAL,              INTENT(IN)  ::  def
  LOGICAL,               OPTIONAL,              INTENT(IN)  :: lDisp
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen) :: s
  CHARACTER(LEN=maxlen), ALLOCATABLE :: svals(:), tname(:)
  LOGICAL,               ALLOCATABLE ::    ll(:)
  INTEGER,               ALLOCATABLE :: ivals(:)
  IF(     PRESENT(def)) lerr = getKeyByIndex_smmm(keyn, svals, ky, tname, num2str(def), lDisp)
  IF(.NOT.PRESENT(def)) lerr = getKeyByIndex_smmm(keyn, svals, ky, tname, lDisp=lDisp)
  IF(lerr) RETURN
  ivals = str2bool(svals)
  ll = ivals == -1
  lerr = ANY(ll); IF(.NOT.lerr) THEN; val = ivals == 1; IF(PRESENT(nam)) nam = tname; RETURN; END IF
  s = 'key "'//TRIM(strStack(keyn, '/'))//'" of '//TRIM(strStack(tname, MASK=ll))//' is not'
  CALL msg(TRIM(s)//' a boolean: '//TRIM(strStack(svals, MASK=ll)), modname)
END FUNCTION getKeyByIndex_lmmm
!==============================================================================================================================



!==============================================================================================================================
!===                                           GET KEYS VALUES FROM TRACERS NAMES                                           ===
!==============================================================================================================================
!=== TRY TO GET THE KEY NAMED "key" FOR THE TRACER NAMED "tname" IN:                                                        ===
!===  * ARGUMENT "ky" DATABASE IF SPECIFIED ; OTHERWISE:                                                                    ===
!===  * IN INTERNAL TRACERS DATABASE "tracers(:)%keys" (IF UNFOUND, IN INTERNAL ISOTOPES DATABASE "isotope%keys(:)")        ===
!=== THE RETURNED VALUE (STRING, AN INTEGER, A REAL OR A LOGICAL) CAN BE EITHER:                                            ===
!===  * A SCALAR                                                                                                            ===
!===  * A VECTOR, WHICH IS THE FOUND VALUE PARSED WITH THE SEPARATOR ","                                                    ===
!===                                                                                                                        ===
!=== SYNTAX:       lerr = getKeyByName_{sirl}{1m}{1m}1(keyn[(:)], val[(:)], tname  [, ky(:)][, def][, lDisp])               ===
!==============================================================================================================================
!=== IF "tname(:)" IS A VECTOR, THE RETURNED VALUES (ONE EACH "tname(:)" ELEMENT) ARE STORED IN THE VECTOR "val(:)"         ===
!===                                                                                                                        ===
!=== SYNTAX        lerr = getKeyByName_{sirl}{1m}mm   (keyn[(:)], val (:), tname(:)[, ky(:)][, def][, lDisp])               ===
!==============================================================================================================================
LOGICAL FUNCTION getKeyByName_s111(keyn, val, tname, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),          INTENT(IN)  :: keyn, tname
  CHARACTER(LEN=maxlen),     INTENT(OUT) :: val
  TYPE(keys_type), OPTIONAL, INTENT(IN)  :: ky(:)
  CHARACTER(LEN=*),OPTIONAL, INTENT(IN)  :: def
  LOGICAL,         OPTIONAL, INTENT(IN)  :: lDisp
  lerr = getKeyByName_sm11([keyn], val, tname, ky, def, lDisp)
END FUNCTION getKeyByName_s111
!==============================================================================================================================
LOGICAL FUNCTION getKeyByName_i111(keyn, val, tname, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),          INTENT(IN)  :: keyn, tname
  INTEGER,                   INTENT(OUT) :: val
  TYPE(keys_type), OPTIONAL, INTENT(IN)  :: ky(:)
  INTEGER,         OPTIONAL, INTENT(IN)  :: def
  LOGICAL,         OPTIONAL, INTENT(IN)  :: lDisp
  lerr = getKeyByName_im11([keyn], val, tname, ky, def, lDisp)
END FUNCTION getKeyByName_i111
!==============================================================================================================================
LOGICAL FUNCTION getKeyByName_r111(keyn, val, tname, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),          INTENT(IN)  :: keyn, tname
  REAL   ,                   INTENT(OUT) :: val
  TYPE(keys_type), OPTIONAL, INTENT(IN)  :: ky(:)
  REAL,            OPTIONAL, INTENT(IN)  :: def
  LOGICAL,         OPTIONAL, INTENT(IN)  :: lDisp
  lerr = getKeyByName_rm11([keyn], val, tname, ky, def, lDisp)
END FUNCTION getKeyByName_r111
!==============================================================================================================================
LOGICAL FUNCTION getKeyByName_l111(keyn, val, tname, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),          INTENT(IN)  :: keyn, tname
  LOGICAL,                   INTENT(OUT) :: val
  TYPE(keys_type), OPTIONAL, INTENT(IN)  :: ky(:)
  LOGICAL,         OPTIONAL, INTENT(IN)  :: def
  LOGICAL,         OPTIONAL, INTENT(IN)  :: lDisp
  lerr = getKeyByName_lm11([keyn], val, tname, ky, def, lDisp)
END FUNCTION getKeyByName_l111
!==============================================================================================================================
!==============================================================================================================================
LOGICAL FUNCTION getKeyByName_sm11(keyn, val, tname, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),          INTENT(IN)  :: keyn(:), tname
  CHARACTER(LEN=maxlen),     INTENT(OUT) :: val
  TYPE(keys_type), OPTIONAL, INTENT(IN)  :: ky(:)
  CHARACTER(LEN=*),OPTIONAL, INTENT(IN)  :: def
  LOGICAL,         OPTIONAL, INTENT(IN)  :: lDisp
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen) :: s, tnam
  LOGICAL :: lD
  lD = .TRUE.; IF(PRESENT(lDisp)) lD = lDisp
  s = 'key "'//TRIM(strStack(keyn, '/'))//'" for "'//TRIM(tname)//'"'
  lerr = .TRUE.
  tnam = strHead(delPhase(tname),'_',.TRUE.)                                             !--- Remove phase and tag
  IF(lerr .AND. PRESENT(ky))         val = fgetKey(ky)                                   !--- "ky"
  IF(lerr .AND. ALLOCATED (tracers)) val = fgetKey(tracers(:)%keys)                      !--- "tracers"
  IF(lerr .AND. ASSOCIATED(isotope)) val = fgetKey(isotope%keys(:))                      !--- "isotope"
  IF(lerr .AND. PRESENT(def)) THEN
     val = def; lerr = .NOT.PRESENT(def)
     CALL msg('Using defaut value of '//TRIM(s)//': '//TRIM(def), modname, lD)
  END IF
  CALL msg('No '//TRIM(s)//' or out of bounds index', modname, lD .AND. lerr)

CONTAINS

 CHARACTER(LEN=maxlen) FUNCTION fgetKey(ky) RESULT(val)
  TYPE(keys_type),  INTENT(IN)  :: ky(:)
  lerr = SIZE(ky) == 0
  IF(lerr) RETURN
           val = fgetKeyIdx(strIdx(ky%name, tname), [keyn], ky, lerr)
  IF(lerr) val = fgetKeyIdx(strIdx(ky%name, tnam ), [keyn], ky, lerr)

END FUNCTION fgetKey

END FUNCTION getKeyByName_sm11
!==============================================================================================================================
LOGICAL FUNCTION getKeyByName_im11(keyn, val, tname, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),          INTENT(IN)  :: keyn(:), tname
  INTEGER,                   INTENT(OUT) :: val
  TYPE(keys_type), OPTIONAL, INTENT(IN)  :: ky(:)
  INTEGER,         OPTIONAL, INTENT(IN)  :: def
  LOGICAL,         OPTIONAL, INTENT(IN)  :: lDisp
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen) :: sval, s
  IF(     PRESENT(def)) lerr = getKeyByName_sm11(keyn, sval, tname, ky, num2str(def), lDisp)
  IF(.NOT.PRESENT(def)) lerr = getKeyByName_sm11(keyn, sval, tname, ky, lDisp=lDisp)
  IF(lerr) RETURN
  val = str2int(sval)
  lerr = val == -HUGE(1)
  s = 'key "'//TRIM(strStack(keyn, '/'))//'" of '//TRIM(tname)//' is not'
  CALL msg(TRIM(s)//' an integer: '//TRIM(sval), modname, lerr)
END FUNCTION getKeyByName_im11
!==============================================================================================================================
LOGICAL FUNCTION getKeyByName_rm11(keyn, val, tname, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),          INTENT(IN)  :: keyn(:), tname
  REAL   ,                   INTENT(OUT) :: val
  TYPE(keys_type), OPTIONAL, INTENT(IN)  :: ky(:)
  REAL,            OPTIONAL, INTENT(IN)  :: def
  LOGICAL,         OPTIONAL, INTENT(IN)  :: lDisp
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen) :: sval, s
  IF(     PRESENT(def)) lerr = getKeyByName_sm11(keyn, sval, tname, ky, num2str(def), lDisp)
  IF(.NOT.PRESENT(def)) lerr = getKeyByName_sm11(keyn, sval, tname, ky, lDisp=lDisp)
  IF(lerr) RETURN
  val = str2real(sval)
  lerr = val == -HUGE(1.)
  s = 'key "'//TRIM(strStack(keyn, '/'))//'" of '//TRIM(tname)//' is not'
  CALL msg(TRIM(s)//' a real: '//TRIM(sval), modname, lerr)
END FUNCTION getKeyByName_rm11
!==============================================================================================================================
LOGICAL FUNCTION getKeyByName_lm11(keyn, val, tname, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),          INTENT(IN)  :: keyn(:), tname
  LOGICAL,                   INTENT(OUT) :: val
  TYPE(keys_type), OPTIONAL, INTENT(IN)  :: ky(:)
  LOGICAL,         OPTIONAL, INTENT(IN)  :: def
  LOGICAL,         OPTIONAL, INTENT(IN)  :: lDisp
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen) :: sval, s
  INTEGER               :: ival
  IF(     PRESENT(def)) lerr = getKeyByName_sm11(keyn, sval, tname, ky, num2str(def), lDisp)
  IF(.NOT.PRESENT(def)) lerr = getKeyByName_sm11(keyn, sval, tname, ky, lDisp=lDisp)
  IF(lerr) RETURN
  ival = str2bool(sval)
  lerr = ival == -1
  s = 'key "'//TRIM(strStack(keyn, '/'))//'" of '//TRIM(tname)//' is not'
  CALL msg(TRIM(s)//' a boolean: '//TRIM(sval), modname, lerr)
  IF(.NOT.lerr) val = ival == 1
END FUNCTION getKeyByName_lm11
!==============================================================================================================================
!==============================================================================================================================
LOGICAL FUNCTION getKeyByName_s1m1(keyn, val, tname, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),                   INTENT(IN)  :: keyn, tname
  CHARACTER(LEN=maxlen), ALLOCATABLE, INTENT(OUT) :: val(:)
  TYPE(keys_type),          OPTIONAL, INTENT(IN)  :: ky(:)
  CHARACTER(LEN=*),         OPTIONAL, INTENT(IN)  :: def
  LOGICAL,                  OPTIONAL, INTENT(IN)  :: lDisp
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen)              :: sval
  lerr = getKeyByName_sm11([keyn], sval, tname, ky, def, lDisp); IF(lerr) RETURN
  lerr = strParse(sval, ',', val)
  CALL msg('can''t parse '//TRIM(sval), modname, lerr); IF(lerr) RETURN
END FUNCTION getKeyByName_s1m1
!==============================================================================================================================
LOGICAL FUNCTION getKeyByName_i1m1(keyn, val, tname, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),          INTENT(IN)  :: keyn, tname
  INTEGER,      ALLOCATABLE, INTENT(OUT) :: val(:)
  TYPE(keys_type), OPTIONAL, INTENT(IN)  :: ky(:)
  INTEGER,         OPTIONAL, INTENT(IN)  :: def
  LOGICAL,         OPTIONAL, INTENT(IN)  :: lDisp
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen)              :: sval, s
  CHARACTER(LEN=maxlen), ALLOCATABLE :: svals(:)
  IF(     PRESENT(def)) lerr = getKeyByName_sm11([keyn], sval, tname, ky, num2str(def), lDisp)
  IF(.NOT.PRESENT(def)) lerr = getKeyByName_sm11([keyn], sval, tname, ky, lDisp=lDisp)
  IF(lerr) RETURN
  lerr = strParse(sval, ',', svals)
  CALL msg('can''t parse '//TRIM(sval), modname, lerr); IF(lerr) RETURN
  val = str2int(svals)
  lerr = ANY(val == -HUGE(1))
  s = 'key "'//TRIM(keyn)//'" of '//TRIM(tname)//' is not'
  CALL msg(TRIM(s)//' an integer: '//TRIM(sval), modname, lerr)
END FUNCTION getKeyByName_i1m1
!==============================================================================================================================
LOGICAL FUNCTION getKeyByName_r1m1(keyn, val, tname, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),          INTENT(IN)  :: keyn, tname
  REAL,         ALLOCATABLE, INTENT(OUT) :: val(:)
  TYPE(keys_type), OPTIONAL, INTENT(IN)  :: ky(:)
  REAL,            OPTIONAL, INTENT(IN)  :: def
  LOGICAL,         OPTIONAL, INTENT(IN)  :: lDisp
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen)              :: sval, s
  CHARACTER(LEN=maxlen), ALLOCATABLE :: svals(:)
  IF(     PRESENT(def)) lerr = getKeyByName_sm11([keyn], sval, tname, ky, num2str(def), lDisp)
  IF(.NOT.PRESENT(def)) lerr = getKeyByName_sm11([keyn], sval, tname, ky, lDisp=lDisp)
  IF(lerr) RETURN
  lerr = strParse(sval, ',', svals)
  CALL msg('can''t parse '//TRIM(sval), modname, lerr); IF(lerr) RETURN
  val = str2real(svals)
  lerr = ANY(val == -HUGE(1.))
  s = 'key "'//TRIM(keyn)//'" of '//TRIM(tname)//' is not'
  CALL msg(TRIM(s)//' a real: '//TRIM(sval), modname, lerr)
END FUNCTION getKeyByName_r1m1
!==============================================================================================================================
LOGICAL FUNCTION getKeyByName_l1m1(keyn, val, tname, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),          INTENT(IN)  :: keyn, tname
  LOGICAL,      ALLOCATABLE, INTENT(OUT) :: val(:)
  TYPE(keys_type), OPTIONAL, INTENT(IN)  :: ky(:)
  LOGICAL,         OPTIONAL, INTENT(IN)  :: def
  LOGICAL,         OPTIONAL, INTENT(IN)  :: lDisp
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen)              :: sval, s
  CHARACTER(LEN=maxlen), ALLOCATABLE :: svals(:)
  INTEGER,               ALLOCATABLE :: ivals(:)
  IF(     PRESENT(def)) lerr = getKeyByName_sm11([keyn], sval, tname, ky, num2str(def), lDisp)
  IF(.NOT.PRESENT(def)) lerr = getKeyByName_sm11([keyn], sval, tname, ky, lDisp=lDisp)
  IF(lerr) RETURN
  lerr = strParse(sval, ',', svals)
  CALL msg('can''t parse '//TRIM(sval), modname, lerr); IF(lerr) RETURN
  ivals = str2bool(svals)
  lerr = ANY(ivals == -1)
  s = 'key "'//TRIM(keyn)//'" of '//TRIM(tname)//' is not'
  CALL msg(TRIM(s)//' a boolean: '//TRIM(sval), modname, lerr)
  IF(.NOT.lerr) val = ivals == 1
END FUNCTION getKeyByName_l1m1
!==============================================================================================================================
!==============================================================================================================================
LOGICAL FUNCTION getKeyByName_smm1(keyn, val, tname, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),                   INTENT(IN)  :: keyn(:), tname
  CHARACTER(LEN=maxlen), ALLOCATABLE, INTENT(OUT) :: val(:)
  TYPE(keys_type),          OPTIONAL, INTENT(IN)  :: ky(:)
  CHARACTER(LEN=*),         OPTIONAL, INTENT(IN)  :: def
  LOGICAL,                  OPTIONAL, INTENT(IN)  :: lDisp
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen) :: sval
  lerr = getKeyByName_sm11(keyn, sval, tname, ky, def, lDisp); IF(lerr) RETURN
  lerr = strParse(sval, ',', val)
  CALL msg('can''t parse '//TRIM(sval), modname, lerr); IF(lerr) RETURN
END FUNCTION getKeyByName_smm1
!==============================================================================================================================
LOGICAL FUNCTION getKeyByName_imm1(keyn, val, tname, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),          INTENT(IN)  :: keyn(:), tname
  INTEGER,      ALLOCATABLE, INTENT(OUT) :: val(:)
  TYPE(keys_type), OPTIONAL, INTENT(IN)  :: ky(:)
  INTEGER,         OPTIONAL, INTENT(IN)  :: def
  LOGICAL,         OPTIONAL, INTENT(IN)  :: lDisp
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen)              :: sval, s
  CHARACTER(LEN=maxlen), ALLOCATABLE :: svals(:)
  IF(     PRESENT(def)) lerr = getKeyByName_sm11(keyn, sval, tname, ky, num2str(def), lDisp)
  IF(.NOT.PRESENT(def)) lerr = getKeyByName_sm11(keyn, sval, tname, ky, lDisp=lDisp)
  IF(lerr) RETURN
  lerr = strParse(sval, ',', svals)
  CALL msg('can''t parse '//TRIM(sval), modname, lerr); IF(lerr) RETURN
  val = str2int(svals)
  lerr = ANY(val == -HUGE(1))
  s = 'key "'//TRIM(strStack(keyn, '/'))//'" of '//TRIM(tname)//' is not'
  CALL msg(TRIM(s)//' an integer: '//TRIM(sval), modname, lerr)
END FUNCTION getKeyByName_imm1
!==============================================================================================================================
LOGICAL FUNCTION getKeyByName_rmm1(keyn, val, tname, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),          INTENT(IN)  :: keyn(:), tname
  REAL,         ALLOCATABLE, INTENT(OUT) :: val(:)
  TYPE(keys_type), OPTIONAL, INTENT(IN)  :: ky(:)
  REAL,            OPTIONAL, INTENT(IN)  :: def
  LOGICAL,         OPTIONAL, INTENT(IN)  :: lDisp
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen)              :: sval, s
  CHARACTER(LEN=maxlen), ALLOCATABLE :: svals(:)
  IF(     PRESENT(def)) lerr = getKeyByName_sm11(keyn, sval, tname, ky, num2str(def), lDisp)
  IF(.NOT.PRESENT(def)) lerr = getKeyByName_sm11(keyn, sval, tname, ky, lDisp=lDisp)
  IF(lerr) RETURN
  lerr = strParse(sval, ',', svals)
  CALL msg('can''t parse '//TRIM(sval), modname, lerr); IF(lerr) RETURN
  val = str2real(svals)
  lerr = ANY(val == -HUGE(1.))
  s = 'key "'//TRIM(strStack(keyn, '/'))//'" of '//TRIM(tname)//' is not'
  CALL msg(TRIM(s)//' a real: '//TRIM(sval), modname, lerr)
END FUNCTION getKeyByName_rmm1
!==============================================================================================================================
LOGICAL FUNCTION getKeyByName_lmm1(keyn, val, tname, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),          INTENT(IN)  :: keyn(:), tname
  LOGICAL,      ALLOCATABLE, INTENT(OUT) :: val(:)
  TYPE(keys_type), OPTIONAL, INTENT(IN)  :: ky(:)
  LOGICAL,         OPTIONAL, INTENT(IN)  :: def
  LOGICAL,         OPTIONAL, INTENT(IN)  :: lDisp
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen)              :: sval, s
  CHARACTER(LEN=maxlen), ALLOCATABLE :: svals(:)
  INTEGER,               ALLOCATABLE :: ivals(:)
  IF(     PRESENT(def)) lerr = getKeyByName_sm11(keyn, sval, tname, ky, num2str(def), lDisp)
  IF(.NOT.PRESENT(def)) lerr = getKeyByName_sm11(keyn, sval, tname, ky, lDisp=lDisp)
  IF(lerr) RETURN
  lerr = strParse(sval, ',', svals)
  CALL msg('can''t parse '//TRIM(sval), modname, lerr); IF(lerr) RETURN
  ivals = str2bool(svals)
  lerr = ANY(ivals == -1)
  s = 'key "'//TRIM(strStack(keyn, '/'))//'" of '//TRIM(tname)//' is not'
  CALL msg(TRIM(s)//' a boolean: '//TRIM(sval), modname, lerr)
  IF(.NOT.lerr) val = ivals == 1
END FUNCTION getKeyByName_lmm1
!==============================================================================================================================
!==============================================================================================================================
LOGICAL FUNCTION getKeyByName_s1mm(keyn, val, tname, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),                   INTENT(IN)  :: keyn, tname(:)
  CHARACTER(LEN=maxlen), ALLOCATABLE, INTENT(OUT) :: val(:)
  TYPE(keys_type),       OPTIONAL,    INTENT(IN)  ::  ky(:)
  CHARACTER(LEN=*),      OPTIONAL,    INTENT(IN)  :: def
  LOGICAL,               OPTIONAL,    INTENT(IN)  :: lDisp
  lerr = getKeyByName_smmm([keyn], val, tname, ky, def, lDisp)
END FUNCTION getKeyByName_s1mm
!==============================================================================================================================
LOGICAL FUNCTION getKeyByName_i1mm(keyn, val, tname, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),          INTENT(IN)  :: keyn, tname(:)
  INTEGER,      ALLOCATABLE, INTENT(OUT) :: val(:)
  TYPE(keys_type), OPTIONAL, INTENT(IN)  ::  ky(:)
  INTEGER,         OPTIONAL, INTENT(IN)  :: def
  LOGICAL,         OPTIONAL, INTENT(IN)  :: lDisp
  lerr = getKeyByName_immm([keyn], val, tname, ky, def, lDisp)
END FUNCTION getKeyByName_i1mm
!==============================================================================================================================
LOGICAL FUNCTION getKeyByName_r1mm(keyn, val, tname, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),          INTENT(IN)  :: keyn, tname(:)
  REAL,         ALLOCATABLE, INTENT(OUT) :: val(:)
  TYPE(keys_type), OPTIONAL, INTENT(IN)  ::  ky(:)
  REAL,            OPTIONAL, INTENT(IN)  :: def
  LOGICAL,         OPTIONAL, INTENT(IN)  :: lDisp
  lerr = getKeyByName_rmmm([keyn], val, tname, ky, def, lDisp)
END FUNCTION getKeyByName_r1mm
!==============================================================================================================================
LOGICAL FUNCTION getKeyByName_l1mm(keyn, val, tname, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),          INTENT(IN)  :: keyn, tname(:)
  LOGICAL,      ALLOCATABLE, INTENT(OUT) :: val(:)
  TYPE(keys_type), OPTIONAL, INTENT(IN)  ::  ky(:)
  LOGICAL,         OPTIONAL, INTENT(IN)  :: def
  LOGICAL,         OPTIONAL, INTENT(IN)  :: lDisp
  lerr = getKeyByName_lmmm([keyn], val, tname, ky, def, lDisp)
END FUNCTION getKeyByName_l1mm
!==============================================================================================================================
!==============================================================================================================================
LOGICAL FUNCTION getKeyByName_smmm(keyn, val, tname, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),                   INTENT(IN)  :: keyn(:), tname(:)
  CHARACTER(LEN=maxlen), ALLOCATABLE, INTENT(OUT) ::  val(:)
  TYPE(keys_type),       OPTIONAL,    INTENT(IN)  ::   ky(:)
  CHARACTER(LEN=*),      OPTIONAL,    INTENT(IN)  ::   def
  LOGICAL,               OPTIONAL,    INTENT(IN)  :: lDisp
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen) :: s
  INTEGER :: iq, nq
  LOGICAL :: lD
  nq = SIZE(tname); ALLOCATE(val(nq))
  lD = .TRUE.; IF(PRESENT(lDisp)) lD = lDisp
  s = 'key "'//TRIM(strStack(keyn, '/'))//'"'
  lerr = .TRUE.
  IF(PRESENT(ky)) THEN;                 val = fgetKey(ky)                                !--- "ky"
  ELSE IF(ALLOCATED(tracers)) THEN;     val = fgetKey(tracers(:)%keys)                   !--- "tracers"
     IF(lerr .AND. ASSOCIATED(isotope)) val = fgetKey(isotope%keys(:))                   !--- "isotope"
  END IF
  IF(.NOT.PRESENT(def)) THEN; CALL msg('No '//TRIM(s)//' found', modname, lD); RETURN; END IF

  !--- DEFAULT VALUE
  val = [(def, iq = 1, SIZE(tname))]
  CALL msg('Using defaut value for '//TRIM(s)//': '//TRIM(def), modname, lD)

CONTAINS

FUNCTION fgetKey(ky) RESULT(val)
  CHARACTER(LEN=maxlen), ALLOCATABLE :: val(:)
  TYPE(keys_type),       INTENT(IN)  :: ky(:)
  LOGICAL,               ALLOCATABLE :: ler(:)
  lerr = SIZE(ky) == 0; IF(lerr) RETURN
  ALLOCATE(ler(SIZE(tname)))
  val = [(fgetKeyIdx(strIdx(ky(:)%name, tname(iq)), keyn, ky, ler(iq)), iq = 1, SIZE(tname))]
  lerr = ANY(ler)
END FUNCTION fgetKey

END FUNCTION getKeyByName_smmm
!==============================================================================================================================
LOGICAL FUNCTION getKeyByName_immm(keyn, val, tname, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),          INTENT(IN)  :: keyn(:), tname(:)
  INTEGER,      ALLOCATABLE, INTENT(OUT) ::  val(:)
  TYPE(keys_type), OPTIONAL, INTENT(IN)  ::   ky(:)
  INTEGER,         OPTIONAL, INTENT(IN)  ::  def
  LOGICAL,         OPTIONAL, INTENT(IN)  :: lDisp
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen) :: s
  CHARACTER(LEN=maxlen), ALLOCATABLE :: svals(:)
  LOGICAL,               ALLOCATABLE ::    ll(:)
  IF(     PRESENT(def)) lerr = getKeyByName_smmm(keyn, svals, tname, ky, num2str(def), lDisp)
  IF(.NOT.PRESENT(def)) lerr = getKeyByName_smmm(keyn, svals, tname, ky, lDisp=lDisp)
  IF(lerr) RETURN
  val = str2int(svals)
  ll = val == -HUGE(1)
  lerr = ANY(ll); IF(.NOT.lerr) RETURN
  s = 'key "'//TRIM(strStack(keyn, '/'))//'" of '//TRIM(strStack(tname, MASK=ll))//' is not'
  CALL msg(TRIM(s)//' an integer: '//TRIM(strStack(svals, MASK=ll)), modname)
END FUNCTION getKeyByName_immm
!==============================================================================================================================
LOGICAL FUNCTION getKeyByName_rmmm(keyn, val, tname, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),          INTENT(IN)  :: keyn(:), tname(:)
  REAL,         ALLOCATABLE, INTENT(OUT) ::  val(:)
  TYPE(keys_type), OPTIONAL, INTENT(IN)  ::   ky(:)
  REAL,            OPTIONAL, INTENT(IN)  ::  def
  LOGICAL,         OPTIONAL, INTENT(IN)  :: lDisp
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen) :: s
  CHARACTER(LEN=maxlen), ALLOCATABLE :: svals(:)
  LOGICAL,               ALLOCATABLE ::    ll(:)
  IF(     PRESENT(def)) lerr = getKeyByName_smmm(keyn, svals, tname, ky, num2str(def), lDisp)
  IF(.NOT.PRESENT(def)) lerr = getKeyByName_smmm(keyn, svals, tname, ky, lDisp=lDisp)
  IF(lerr) RETURN
  val = str2real(svals)
  ll = val == -HUGE(1.)
  lerr = ANY(ll); IF(.NOT.lerr) RETURN
  s = 'key "'//TRIM(strStack(keyn, '/'))//'" of '//TRIM(strStack(tname, MASK=ll))//' is not'
  CALL msg(TRIM(s)//' a real: '//TRIM(strStack(svals, MASK=ll)), modname)
END FUNCTION getKeyByName_rmmm
!==============================================================================================================================
LOGICAL FUNCTION getKeyByName_lmmm(keyn, val, tname, ky, def, lDisp) RESULT(lerr)
  CHARACTER(LEN=*),          INTENT(IN)  :: keyn(:), tname(:)
  LOGICAL,      ALLOCATABLE, INTENT(OUT) ::  val(:)
  TYPE(keys_type), OPTIONAL, INTENT(IN)  ::   ky(:)
  LOGICAL,         OPTIONAL, INTENT(IN)  ::  def
  LOGICAL,         OPTIONAL, INTENT(IN)  :: lDisp
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen) :: s
  CHARACTER(LEN=maxlen), ALLOCATABLE :: svals(:)
  LOGICAL,               ALLOCATABLE ::    ll(:)
  INTEGER,               ALLOCATABLE :: ivals(:)
  IF(     PRESENT(def)) lerr = getKeyByName_smmm(keyn, svals, tname, ky, num2str(def), lDisp)
  IF(.NOT.PRESENT(def)) lerr = getKeyByName_smmm(keyn, svals, tname, ky, lDisp=lDisp)
  IF(lerr) RETURN
  ivals = str2bool(svals)
  ll = ivals == -1
  lerr = ANY(ll); IF(.NOT.lerr) THEN; val = ivals == 1; RETURN; END IF
  s = 'key "'//TRIM(strStack(keyn, '/'))//'" of '//TRIM(strStack(tname, MASK=ll))//' is not'
  CALL msg(TRIM(s)//' a boolean: '//TRIM(strStack(svals, MASK=ll)), modname)
END FUNCTION getKeyByName_lmmm
!==============================================================================================================================


!==============================================================================================================================
!=== ROUTINES TO GET OR PUT TRACERS AND ISOTOPES DATABASES, SINCE tracers AND isotopes ARE PRIVATE VARIABLES ==================
!==============================================================================================================================
SUBROUTINE setKeysDBase(tracers_, isotopes_, isotope_)
  TYPE(trac_type), OPTIONAL, INTENT(IN) ::  tracers_(:)
  TYPE(isot_type), OPTIONAL, INTENT(IN) :: isotopes_(:)
  TYPE(isot_type), OPTIONAL, INTENT(IN) :: isotope_
!------------------------------------------------------------------------------------------------------------------------------
  TYPE(isot_type), ALLOCATABLE :: iso(:)
  INTEGER :: ix, nbIso
  IF(PRESENT( tracers_)) THEN;  tracers =  tracers_; ELSE; ALLOCATE( tracers(0)); END IF
  IF(PRESENT(isotopes_)) THEN; isotopes = isotopes_; ELSE; ALLOCATE(isotopes(0)); END IF
  IF(PRESENT(isotope_ )) THEN
    ix = strIdx(isotopes(:)%parent, isotope_%parent)
    IF(ix /= 0) THEN
      isotopes(ix) = isotope_
    ELSE
      nbIso = SIZE(isotopes); ALLOCATE(iso(nbIso+1)); iso(1:nbIso) = isotopes; iso(nbIso+1) = isotope_
      CALL MOVE_ALLOC(FROM=iso, TO=isotopes)
    END IF
  END IF
END SUBROUTINE setKeysDBase
!==============================================================================================================================
SUBROUTINE getKeysDBase(tracers_, isotopes_, isotope_)
  TYPE(trac_type), OPTIONAL, ALLOCATABLE, INTENT(OUT) ::  tracers_(:)
  TYPE(isot_type), OPTIONAL, ALLOCATABLE, INTENT(OUT) :: isotopes_(:)
  TYPE(isot_type), OPTIONAL,              INTENT(OUT) :: isotope_
!------------------------------------------------------------------------------------------------------------------------------
  INTEGER :: ix
  IF(PRESENT( tracers_)) THEN;  tracers_ =  tracers; ELSE; ALLOCATE( tracers_(0)); END IF
  IF(PRESENT(isotopes_)) THEN; isotopes_ = isotopes; ELSE; ALLOCATE(isotopes_(0)); END IF
  IF(PRESENT(isotope_ )) THEN; ix = strIdx(isotopes(:)%parent, isotope%parent); IF(ix /= 0) isotope_=isotopes(ix); END IF
END SUBROUTINE getKeysDBase
!==============================================================================================================================


!==============================================================================================================================
!=== REMOVE, IF ANY, OR ADD THE PHASE SUFFIX OF THE TRACER NAMED "s" ==========================================================
!==============================================================================================================================
ELEMENTAL CHARACTER(LEN=maxlen) FUNCTION delPhase(s) RESULT(out)
  CHARACTER(LEN=*), INTENT(IN) :: s
!------------------------------------------------------------------------------------------------------------------------------
  INTEGER :: ix, ip, ns
  out = s; ns = LEN_TRIM(s)
  IF(s == '' .OR. ns<=2) RETURN                                                !--- Empty string or LEN(name)<=2: nothing to do
  IF(s(1:3)=='H2O' .AND. INDEX(old_phases,s(4:4))/=0 .AND. (ns==4 .OR. s(5:5)=='_')) THEN
    out='H2O'//s(5:ns)                                                         !--- H2O<phase>[_<iso>][_<tag>]
  ELSE IF(s(ns-1:ns-1)==phases_sep .AND. INDEX(known_phases,s(ns:ns))/=0) THEN
    out = s(1:ns-2); RETURN                                                    !--- <var><phase_sep><phase>
  ELSE; DO ip = 1, nphases; ix = INDEX(s, phases_sep//known_phases(ip:ip)//'_'); IF(ix /= 0) EXIT; END DO
    IF(ip /= nphases+1) out = s(1:ix-1)//s(ix+2:ns)                            !--- <var><phase_sep><phase>_<tag>
  END IF
END FUNCTION delPhase
!==============================================================================================================================
CHARACTER(LEN=maxlen) FUNCTION addPhase_s1(s,pha) RESULT(out)
  CHARACTER(LEN=*),           INTENT(IN) :: s
  CHARACTER(LEN=1),           INTENT(IN) :: pha
!------------------------------------------------------------------------------------------------------------------------------
  INTEGER :: l, i
  out = s
  IF(s == '') RETURN                                                           !--- Empty string: nothing to do
  i = INDEX(s, '_')                                                            !--- /=0 for <var>_<tag> tracers names
  l = LEN_TRIM(s)
  IF(i == 0) out =  TRIM(s)//phases_sep//pha                                   !--- <var>       => return <var><sep><pha>
  IF(i /= 0) out = s(1:i-1)//phases_sep//pha//'_'//s(i+1:l)                    !--- <var>_<tag> => return <var><sep><pha>_<tag>
END FUNCTION addPhase_s1
!==============================================================================================================================
FUNCTION addPhase_sm(s,pha) RESULT(out)
  CHARACTER(LEN=*),           INTENT(IN) :: s(:)
  CHARACTER(LEN=1),           INTENT(IN) :: pha
  CHARACTER(LEN=maxlen),     ALLOCATABLE :: out(:)
!------------------------------------------------------------------------------------------------------------------------------
  INTEGER :: k
  out = [( addPhase_s1(s(k), pha), k=1, SIZE(s) )]
END FUNCTION addPhase_sm
!==============================================================================================================================
CHARACTER(LEN=maxlen) FUNCTION addPhase_i1(s,ipha,phases) RESULT(out)
  CHARACTER(LEN=*),           INTENT(IN) :: s
  INTEGER,                    INTENT(IN) :: ipha
  CHARACTER(LEN=*), OPTIONAL, INTENT(IN) :: phases
!------------------------------------------------------------------------------------------------------------------------------
  out = s
  IF(s == '') RETURN                                                           !--- Empty string: nothing to do
  IF(ipha == 0 .OR. ipha > nphases) RETURN                                     !--- Absurd index: no phase to add
  IF(     PRESENT(phases)) out = addPhase_s1(s,       phases(ipha:ipha))
  IF(.NOT.PRESENT(phases)) out = addPhase_s1(s, known_phases(ipha:ipha))
END FUNCTION addPhase_i1
!==============================================================================================================================
FUNCTION addPhase_im(s,ipha,phases) RESULT(out)
  CHARACTER(LEN=*),           INTENT(IN) :: s(:)
  INTEGER,                    INTENT(IN) :: ipha
  CHARACTER(LEN=*), OPTIONAL, INTENT(IN) :: phases
  CHARACTER(LEN=maxlen),     ALLOCATABLE :: out(:)
!------------------------------------------------------------------------------------------------------------------------------
  INTEGER :: k
  IF(     PRESENT(phases)) out = [( addPhase_i1(s(k), ipha,       phases), k=1, SIZE(s) )]
  IF(.NOT.PRESENT(phases)) out = [( addPhase_i1(s(k), ipha, known_phases), k=1, SIZE(s) )]
END FUNCTION addPhase_im
!==============================================================================================================================


!==============================================================================================================================
!=== APPEND TRACERS DATABASE "tracs" WITH TRACERS/KEYS "tname"/"keys" ; SAME FOR INTERNAL DATABASE "tracers" ==================
!==============================================================================================================================
LOGICAL FUNCTION addTracer_1(tname, keys, tracs) RESULT(lerr)
  CHARACTER(LEN=*),             INTENT(IN)    :: tname
  TYPE(keys_type),              INTENT(IN)    ::  keys
  TYPE(trac_type), ALLOCATABLE, INTENT(INOUT) :: tracs(:)
  TYPE(trac_type), ALLOCATABLE :: tr(:)
  CHARACTER(LEN=maxlen), ALLOCATABLE :: tnames(:)
  INTEGER :: nt, ix
  IF(ALLOCATED(tracs)) THEN
     lerr = getKey('name', tnames, ky=tracs(:)%keys); IF(lerr) RETURN
     nt = SIZE(tracs)
     ix = strIdx(tnames, tname)
     CALL msg('Modifying existing tracer "'//TRIM(tname)//'"', modname, ix /= 0)
     CALL msg('Appending with tracer "'    //TRIM(tname)//'"', modname, ix == 0)
     IF(ix == 0) THEN
        ix = nt+1; ALLOCATE(tr(nt+1)); tr(1:nt) = tracs(1:nt); CALL MOVE_ALLOC(FROM=tr, TO=tracs)
     END IF
  ELSE
     CALL msg('Creating a tracer descriptor with tracer "'//TRIM(tname)//'"', modname)
     ix = 1; ALLOCATE(tracs(1))
  END IF
  CALL addKey('name', tname, tracs(ix)%keys)
  tracs(ix)%name = tname
  tracs(ix)%keys = keys

END FUNCTION addTracer_1
!==============================================================================================================================
LOGICAL FUNCTION addTracer_1def(tname, keys) RESULT(lerr)
  CHARACTER(LEN=*),             INTENT(IN)    :: tname
  TYPE(keys_type),              INTENT(IN)    ::  keys
  lerr = addTracer_1(tname, keys, tracers)
END FUNCTION addTracer_1def
!==============================================================================================================================


!==============================================================================================================================
LOGICAL FUNCTION delTracer_1(tname, tracs) RESULT(lerr)
  CHARACTER(LEN=*),                     INTENT(IN)    :: tname
  TYPE(trac_type), ALLOCATABLE, TARGET, INTENT(INOUT) :: tracs(:)
  TYPE(trac_type), ALLOCATABLE :: tr(:)
  CHARACTER(LEN=maxlen), ALLOCATABLE :: tnames(:)
  INTEGER :: nt, ix
  lerr = .NOT.ALLOCATED(tracs)
  CALL msg('Can''t remove tracer "'//TRIM(tname)//'" from an empty tracers descriptor', modname, lerr); IF(lerr) RETURN
  nt = SIZE(tracs)
  lerr = getKey('name', tnames, ky=tracs(:)%keys); IF(lerr) RETURN
  ix = strIdx(tnames, tname)
  CALL msg('Removing tracer "'             //TRIM(tname)//'"', modname, ix /= 0)
  CALL msg('Can''t remove unknown tracer "'//TRIM(tname)//'"', modname, ix == 0)
  IF(ix /= 0) THEN
     ALLOCATE(tr(nt-1)); tr(1:ix-1) = tracs(1:ix-1); tr(ix:nt-1) = tracs(ix+1:nt); CALL MOVE_ALLOC(FROM=tr, TO=tracs)
  END IF
END FUNCTION delTracer_1
!==============================================================================================================================
LOGICAL FUNCTION delTracer_1def(tname) RESULT(lerr)
  CHARACTER(LEN=*), INTENT(IN) :: tname
  lerr = delTracer(tname, tracers)
END FUNCTION delTracer_1def
!==============================================================================================================================


!==============================================================================================================================
!=== GET PHASE INDEX IN THE POSSIBLE PHASES LIST OR IN A SPECIFIED LIST ("phases") ============================================
!==============================================================================================================================
INTEGER FUNCTION getiPhase(tname, phases) RESULT(iPhase)
  CHARACTER(LEN=*),           INTENT(IN)  :: tname
  CHARACTER(LEN=*), OPTIONAL, INTENT(IN)  :: phases
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen) :: phase
  IF(     PRESENT(phases)) phase = getPhase(tname,       phases, iPhase)
  IF(.NOT.PRESENT(phases)) phase = getPhase(tname, known_phases, iPhase)
END FUNCTION getiPhase
!==============================================================================================================================
CHARACTER(LEN=1) FUNCTION getPhase(tname, phases, iPhase) RESULT(phase)
  CHARACTER(LEN=*),           INTENT(IN)  :: tname
  CHARACTER(LEN=*), OPTIONAL, INTENT(IN)  :: phases
  INTEGER,          OPTIONAL, INTENT(OUT) :: iPhase
!------------------------------------------------------------------------------------------------------------------------------
  INTEGER :: ip
  phase = TRIM(strHead(strTail(tname, phases_sep), '_', .TRUE.))     !--- <nam><sep><pha>[_<tag>] -> <pha>[_<tag>] -> <pha>
  IF(     PRESENT(phases)) ip = INDEX(      phases, phase)
  IF(.NOT.PRESENT(phases)) ip = INDEX(known_phases, phase)
  IF(ip == 0) phase = 'g'
  IF(PRESENT(iPhase)) iPhase = ip
END FUNCTION getPhase
!==============================================================================================================================


!==============================================================================================================================
!============ CONVERT WATER-DERIVED NAMES FROM FORMER TO CURRENT CONVENTION ; OTHER NAMES ARE LEFT UNTOUCHED ==================
!======= NAMES STRUCTURE: H2O[<phase>][_<isotope>][_<tag>] (<phase> from "old_phases", <isotope> from "oldH2OIso") ============
!==============================================================================================================================
CHARACTER(LEN=maxlen) FUNCTION old2newH2O_1(oldName, iPhase) RESULT(newName)
  CHARACTER(LEN=*),  INTENT(IN)  :: oldName
  INTEGER, OPTIONAL, INTENT(OUT) :: iPhase
!------------------------------------------------------------------------------------------------------------------------------
  CHARACTER(LEN=maxlen), ALLOCATABLE :: tmp(:)
  INTEGER :: ix, ip, nt
  LOGICAL :: lerr
  newName = oldName
  IF(PRESENT(iPhase)) iPhase = 1                                               !--- Default: gaseous phase
  lerr = strParse(oldName, '_', tmp, nt)                                       !--- Parsing: 1 up to 3 elements.
  ip = strIdx( [('H2O'//old_phases(ix:ix), ix=1, nphases)], tmp(1) )           !--- Phase index
  IF(ip /= 0 .AND. PRESENT(iPhase)) iPhase = ip                                !--- Returning phase index
  IF(ip == 0 .AND. tmp(1) /= 'H2O')   RETURN                                   !--- Not an old-style water-related species
  IF(nt == 1) THEN
    newName = addPhase('H2O',ip)                                               !=== WATER WITH OR WITHOUT PHASE
  ELSE
    ix = strIdx(oldH2OIso, tmp(2))                                             !--- Index in the known isotopes list
    IF(ix /= 0) newName = newH2OIso(ix)                                        !--- Move to new isotope name
    IF(ip /= 0) newName = addPhase(newName, ip)                                !--- Add phase to isotope name
    IF(nt == 3) newName = TRIM(newName)//'_'//TRIM(tmp(3))                     !=== WATER ISOTOPE OR TAGGING TRACER
  END IF
END FUNCTION old2newH2O_1
!==============================================================================================================================
FUNCTION old2newH2O_m(oldName) RESULT(newName)
  CHARACTER(LEN=*), INTENT(IN) :: oldName(:)
  CHARACTER(LEN=maxlen)        :: newName(SIZE(oldName))
!------------------------------------------------------------------------------------------------------------------------------
  INTEGER :: i
  newName = [(old2newH2O_1(oldName(i)), i=1, SIZE(oldName))]
END FUNCTION old2newH2O_m
!==============================================================================================================================


!==============================================================================================================================
!============ CONVERT WATER-DERIVED NAMES FROM CURRENT TO FORMER CONVENTION ; OTHER NAMES ARE LEFT UNTOUCHED ==================
!==== NAMES STRUCTURE: <var>[<phase_sep><phase>][_<tag>] (<phase> from "known_phases", <var> = 'H2O' or from "newH2OIso") =====
!==============================================================================================================================
CHARACTER(LEN=maxlen) FUNCTION new2oldH2O_1(newName, iPhase) RESULT(oldName)
  CHARACTER(LEN=*),  INTENT(IN)  :: newName
  INTEGER, OPTIONAL, INTENT(OUT) :: iPhase
!------------------------------------------------------------------------------------------------------------------------------
  INTEGER :: ix, ip
  CHARACTER(LEN=maxlen) :: var
  oldName = newName
  ip = getiPhase(newName)                                                      !--- Phase index
  IF(PRESENT(iPhase)) iPhase = MAX(ip, 1)                                      !--- Return phase index ; default: 1 (gazeous)
  var = TRIM(strHead(newName, phases_sep, .TRUE.))                             !--- Variable without phase and tag
  ix = strIdx(newH2OIso, var)                                                  !--- Index in the known H2O isotopes list
  IF(ix == 0 .AND. var /= 'H2O') RETURN                                        !--- Neither H2O nor an H2O isotope => finished
  oldName = 'H2O'
  IF(ip /= 0) oldName = TRIM(oldName)//old_phases(ip:ip)                       !--- Add phase if needed
  IF(ix /= 0) oldName = TRIM(oldName)//'_'//oldH2OIso(ix)                      !--- H2O isotope name
  IF(newName /= addPhase(var, ip)) &
    oldName = TRIM(oldName)//strTail(newName, '_', .TRUE.)                     !--- Add the tag suffix
  IF(ip == 0 .AND. ix /= 0) oldName = strTail(oldName, '_')                    !--- Isotope with no phase: remove 'H2O_' prefix
END FUNCTION new2oldH2O_1
!==============================================================================================================================
FUNCTION new2oldH2O_m(newName) RESULT(oldName)
  CHARACTER(LEN=*), INTENT(IN) :: newName(:)
  CHARACTER(LEN=maxlen)        :: oldName(SIZE(newName))
!------------------------------------------------------------------------------------------------------------------------------
  INTEGER :: i
  oldName = [(new2oldH2O_1(newName(i)), i=1, SIZE(newName))]
END FUNCTION new2oldH2O_m
!==============================================================================================================================

END MODULE readTracFiles_mod
