source: readTracFiles_mod.f90 @ 13

Last change on this file since 13 was 13, checked in by dcugnet, 3 years ago
  • Rename "iTraPha" and "iZonIso": now named "iqTraPha" and "itZonIso".
  • Add the possibility of a cloud H2O tracer (nqo=4) for ok_ice_sursat=y and 4 H2O tracers in traceur.def (commit from O. Boucher on the LMDZ trunk branch, r4075).
  • Simplify the naming (old and new conventions) operations using old2newName and new2oldName routines.
  • Remove the useless "old2newPhase".
  • Add a function form "fGetKey" for string keys, using either a string tracer name or its index.
  • Function addPhase accepts now a phase index instead of the string phase itself.
  • The fType==1 case (old traceur.def files) now really complies with the old conventions, especially for isotopes with name structure <parent[phase]>_<isotope>.
  • Add the calling routine name "modname" to some routines to improve displayed messages accuracy.
  • Propagates the new generations counting convention (from 0 and not from 1).
  • Few fixes about the %isInPhysics and %component frames.
  • Few fixes about the code imposing the water to be at the beginning of the tracers table.
File size: 97.2 KB
Line 
1MODULE readTracFiles_mod
2
3  USE strings_mod,    ONLY: msg, testFile,  strFind, strStack, strReduce,  strHead, strCount, find, fmsg, reduceExpr, &
4             removeComment, cat, checkList, str2int, strParse, strReplace, strTail, strIdx, maxlen, test, dispTable, get_in
5  USE trac_types_mod, ONLY: trac_type, isot_type, keys_type
6
7  IMPLICIT NONE
8
9  PRIVATE
10
11  PUBLIC :: initIsotopes, maxlen, trac_type, isot_type, keys_type
12  PUBLIC :: readTracersFiles, indexUpdate, setGeneration             !--- TOOLS ASSOCIATED TO TRACERS  DESCRIPTORS
13  PUBLIC :: readIsotopesFile                                         !--- TOOLS ASSOCIATED TO ISOTOPES DESCRIPTORS
14  PUBLIC :: getKey_init, getKey, fGetKey, setDirectKeys              !--- GET/SET KEYS FROM/TO tracers & isotopes
15
16  PUBLIC :: known_phases, old_phases, nphases, phases_names, &       !--- VARIABLES RELATED TO THE PHASES
17            phases_sep, delPhase, addPhase, new2oldPhase,    &       !--- + ROUTINES TO ADD/REMOVE PHASE TO/FROM A NAME
18            old2newName, new2oldName
19
20  PUBLIC :: tran0, idxAncestor, ancestor                             !--- GENERATION 0 TRACER + TOOLS FOR GENERATIONS
21
22!------------------------------------------------------------------------------------------------------------------------------
23  TYPE :: dataBase_type                                              !=== TYPE FOR TRACERS SECTION
24    CHARACTER(LEN=maxlen)  :: name                                   !--- Section name
25    TYPE(trac_type), ALLOCATABLE :: trac(:)                          !--- Tracers descriptors
26  END TYPE dataBase_type
27!------------------------------------------------------------------------------------------------------------------------------
28  INTERFACE getKey
29    MODULE PROCEDURE getKeyByName_s1, getKeyByName_i1, getKeyByName_r1, getKeyByName_sm, getKeyByName_im, getKeyByName_rm
30  END INTERFACE getKey
31!------------------------------------------------------------------------------------------------------------------------------
32  INTERFACE fGetKey;       MODULE PROCEDURE fgetKeyByIndex_s1, fgetKeyByName_s1; END INTERFACE fGetKey
33  INTERFACE tracersSubset; MODULE PROCEDURE trSubset_Indx, trSubset_Name, trSubset_gen0Name; END INTERFACE tracersSubset
34  INTERFACE idxAncestor;   MODULE PROCEDURE idxAncestor_1, idxAncestor_m; END INTERFACE idxAncestor
35  INTERFACE    ancestor;   MODULE PROCEDURE    ancestor_1,    ancestor_m; END INTERFACE    ancestor
36  INTERFACE    addPhase;   MODULE PROCEDURE   addPhase_s1,   addPhase_sm,   addPhase_i1,   addPhase_im; END INTERFACE addPhase
37!------------------------------------------------------------------------------------------------------------------------------
38
39  !=== MAIN DATABASE: files sections descriptors
40  TYPE(dataBase_type), SAVE, ALLOCATABLE, TARGET :: dBase(:)
41
42  !--- SOME PARAMETERS THAT ARE NOT LIKELY TO CHANGE OFTEN
43  CHARACTER(LEN=maxlen), SAVE      :: tran0        = 'air'           !--- Default transporting fluid
44  CHARACTER(LEN=maxlen), PARAMETER :: old_phases   = 'vli'           !--- Old phases for water (no separator)
45  CHARACTER(LEN=maxlen), PARAMETER :: known_phases = 'gls'           !--- Known phases initials
46  INTEGER,               PARAMETER :: nphases=LEN_TRIM(known_phases) !--- Number of phases
47  CHARACTER(LEN=maxlen), SAVE      :: phases_names(nphases) &        !--- Known phases names
48                                = ['gaseous', 'liquid ', 'solid  ']
49  CHARACTER(LEN=1), SAVE :: phases_sep  =  '_'                       !--- Phase separator
50  LOGICAL,          SAVE :: tracs_merge = .TRUE.                     !--- Merge/stack tracers lists
51  LOGICAL,          SAVE :: lSortByGen  = .TRUE.                     !--- Sort by growing generation
52
53  !=== LOCAL COPIES OF THE TRACERS AND ISOTOPES DESCRIPTORS, USED BY getKey (INITIALIZED IN getKey_init)
54  TYPE(trac_type), ALLOCATABLE, TARGET, SAVE ::  tracers(:)
55  TYPE(isot_type), ALLOCATABLE, TARGET, SAVE :: isotopes(:)
56
57  CHARACTER(LEN=maxlen) :: modname
58
59CONTAINS
60
61!==============================================================================================================================
62!==============================================================================================================================
63!=== READ ONE OR SEVERAL TRACER FILES AND FILL A "tr" TRACERS DESCRIPTOR DERIVED TYPE.
64!=== THE RETURN VALUE fType DEPENDS ON WHAT IS FOUND:
65!===  0: NO ADEQUATE FILE FOUND ; DEFAULT VALUES MUST BE USED
66!===  1: AN "OLD STYLE" TRACERS FILE "traceur.def":
67!===    First line: <nb tracers>     Other lines: <hadv> <vadv> <tracer name> [<parent name>]
68!===  2: A  "NEW STYLE" TRACERS FILE  "tracer.def" WITH SEVERAL SECTIONS.
69!===  3: SEVERAL  "  "  TRACERS FILES "tracer_<component>.def" WITH A SINGLE SECTION IN EACH.
70!=== REMARKS:
71!===  * EACH SECTION BEGINS WITH A "&<section name> LINE
72!===  * DEFAULT VALUES FOR ALL THE SECTIONS OF THE FILE ARE DEFINED IN THE SPECIAL SECTION "&default"
73!===  * EACH SECTION LINE HAS THE STRUCTURE:  <name(s)>  <key1>=<value1> <key2>=<value2> ...
74!===  * SO FAR, THE DEFINED KEYS ARE: parent, phases, hadv, vadv, type
75!===  * <name> AND <parent> CAN BE LISTS OF COMA-SEPARATED TRACERS ; THE ROUTINE EXPAND THESE FACTORIZATIONS.
76!=== FUNCTION RETURN VALUE "lerr" IS FALSE IN CASE SOMETHING WENT WRONG.
77!=== ABOUT THE KEYS:
78!     * The "keys" component (of type keys_type) is in principle enough to store everything we could need.
79!     But some variables are stored as direct-access keys to make the code more readable and because they are used often.
80!     * Most of the direct-access keys are set in this module, but some are not (longName, iadv, isAdvected for now).
81!     * Some of the direct-access keys must be updated (using the routine "setDirectKeys") is a subset of "tracers(:)"
82!     is extracted: the indexes are no longer valid for a subset (examples: iqParent, iqDescen).
83!     * If you need to convert a %key/%val pair into a direct-access key, add the corresponding line in "setDirectKeys".
84!==============================================================================================================================
85LOGICAL FUNCTION readTracersFiles(type_trac, fType, tracs) RESULT(lerr)
86!------------------------------------------------------------------------------------------------------------------------------
87  CHARACTER(LEN=*),             INTENT(IN)  :: type_trac              !--- List of components used
88  INTEGER,                      INTENT(OUT) :: fType                  !--- Type of input file found
89  TYPE(trac_type), ALLOCATABLE, INTENT(OUT) :: tracs(:)
90  CHARACTER(LEN=maxlen),  ALLOCATABLE :: s(:), sections(:), trac_files(:)
91  CHARACTER(LEN=maxlen) :: str, fname, mesg
92  INTEGER               :: is, nsec, ierr, it, ntrac, ns, ip, ix
93  LOGICAL, ALLOCATABLE  :: ll(:), lGen3(:)
94!------------------------------------------------------------------------------------------------------------------------------
95  lerr = .FALSE.
96  modname = 'readTracersFiles'
97  IF(.NOT.ALLOCATED(dBase)) ALLOCATE(dBase(0))
98
99  !--- Required sections + corresponding files names (new style single section case)
100  IF(test(strParse(type_trac, ',', sections), lerr)) RETURN          !--- Parse "type_trac" list
101
102  nsec = SIZE(sections, DIM=1)
103  ALLOCATE(trac_files(nsec)); DO is=1, nsec; trac_files(is) = 'tracer_'//TRIM(sections(is))//'.def'; END DO
104
105  !--- LOOK AT AVAILABLE FILES
106  ll = .NOT.testFile(trac_files)
107  fType = 0
108  IF(.NOT.testFile('traceur.def') .AND. nsec==1) fType = 1           !--- OLD STYLE FILE
109  IF(.NOT.testFile('tracer.def'))                fType = 2           !--- NEW STYLE ; SINGLE  FILE, SEVERAL SECTIONS
110  IF(ALL(ll))                                    fType = 3           !--- NEW STYLE ; SEVERAL FILES, SINGLE SECTION USED
111  IF(ANY(ll) .AND. fType/=3) THEN                                    !--- MISSING FILES
112    IF(test(checkList(trac_files, .NOT.ll, 'Failed reading tracers description', 'files', 'missing'), lerr)) RETURN
113  END IF
114
115  !--- CHECK WHETHER type_trac AND FILE TYPE ARE COMPATIBLE
116  IF(test(fmsg('No multiple sections for the old format "traceur.def"', ll = SIZE(sections)>1 .AND. fType==1), lerr)) RETURN
117
118  !--- TELLS WHAT WAS IS ABOUT TO BE USED
119  IF (fmsg('No adequate tracers description file(s) found ; default values will be used',          modname, fType==0)) RETURN
120  CALL msg('Trying to read old-style tracers description file "traceur.def"',                      modname, fType==1)
121  CALL msg('Trying to read the new style multi-sections tracers description file "tracer.def"',    modname, fType==2)
122  CALL msg('Trying to read the new style single section tracers description files "tracer_*.def"', modname, fType==3)
123
124  !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
125  SELECT CASE(fType)                         !--- Set %name, %genOName, %parent, %type, %phase, %component, %iGeneration, %keys
126  !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
127    CASE(1)                                                               !=== OLD FORMAT "traceur.def"
128    !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
129      !--- OPEN THE "traceur.def" FILE
130      OPEN(90, FILE="traceur.def", FORM='formatted', STATUS='old', IOSTAT=ierr)
131
132      !--- GET THE TRACERS NUMBER
133      READ(90,'(i3)',IOSTAT=ierr)ntrac                               !--- Number of lines/tracers
134      IF(test(fmsg('Invalid format for "'//TRIM(fname)//'"', modname, ierr /= 0), lerr)) RETURN
135
136      !--- READ THE REMAINING LINES: <hadv> <vadv> <tracer> [<transporting fluid>]
137      ALLOCATE(tracs(ntrac))
138      DO it=1,ntrac                                                  !=== READ RAW DATA: loop on the line/tracer number
139        READ(90,'(a)',IOSTAT=ierr) str
140        IF(test(fmsg('Invalid format for "' //TRIM(fname)//'"', modname, ierr>0), lerr)) RETURN
141        IF(test(fmsg('Not enough lines in "'//TRIM(fname)//'"', modname, ierr<0), lerr)) RETURN
142        ll = strParse(str, ' ', s, n=ns)
143        CALL msg('This file is for air tracers only',           modname, ns == 3 .AND. it == 1)
144        CALL msg('This files specifies the transporting fluid', modname, ns == 4 .AND. it == 1)
145        tracs(it)%name   = old2newName(s(3), ip)                     !--- Set %name:   name   of the tracer
146        tracs(it)%parent = tran0                                     !--- Default transporting fluid name
147        IF(ns == 4) tracs(it)%parent = old2newName(s(4))             !--- Set %parent: parent of the tracer
148        tracs(it)%phase = known_phases(ip:ip)                        !--- Set %phase:  tracer phase (default: "g"azeous)
149        tracs(it)%component = TRIM(type_trac)                        !--- Set %component: model component name
150        tracs(it)%keys%key = ['hadv', 'vadv']                        !--- Set %keys%key
151        tracs(it)%keys%val = s(1:2)                                  !--- Set %keys%val
152      END DO
153      CLOSE(90)
154      CALL setGeneration(tracs)                                      !--- Set %iGeneration and %gen0Name
155      WHERE(tracs%iGeneration == 2) tracs%type = 'tag'               !--- Set %type:        'tracer' or 'tag'
156      IF(test(checkTracers(tracs, fname, fname), lerr)) RETURN       !--- Detect orphans and check phases
157      IF(test(checkUnique (tracs, fname, fname), lerr)) RETURN       !--- Detect repeated tracers
158      CALL sortTracers  (tracs)                                      !--- Sort the tracers
159      tracs(:)%keys%name = tracs(:)%name                             !--- Copy tracers names in keys components
160    !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
161    CASE(2); IF(test(feedDBase(["tracer.def"], [type_trac], modname), lerr)) RETURN !=== SINGLE   FILE, MULTIPLE SECTIONS
162    !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
163    CASE(3); IF(test(feedDBase(  trac_files  ,  sections,   modname), lerr)) RETURN !=== MULTIPLE FILES, SINGLE  SECTION
164  !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
165  END SELECT
166  !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
167
168  IF(ALL([2,3] /= fType)) RETURN
169
170  IF(nsec  == 1) THEN;
171    tracs = dBase(1)%trac
172  ELSE IF(tracs_merge) THEN
173    CALL msg('The multiple required sections will be MERGED.',    modname)
174    IF(test(mergeTracers(dBase, tracs), lerr)) RETURN
175  ELSE
176    CALL msg('The multiple required sections will be CUMULATED.', modname)
177    IF(test(cumulTracers(dBase, tracs), lerr)) RETURN
178  END IF
179  WHERE(tracs%gen0Name(1:3) /= 'H2O') tracs%isInPhysics=.TRUE.       !--- Set %isInPhysics: passed to physics
180  CALL setDirectKeys(tracs)                                          !--- Set %iqParent, %iqDescen, %nqDescen, %nqChilds
181
182END FUNCTION readTracersFiles
183!==============================================================================================================================
184
185!==============================================================================================================================
186LOGICAL FUNCTION feedDBase(fnames, snames, modname) RESULT(lerr)
187! Purpose: Read the sections "snames(is)" (coma-separated list) from each "fnames(is)"
188!   file and create the corresponding tracers set descriptors in the database "dBase":
189! * dBase(id)%name                : section name
190! * dBase(id)%trac(:)%name        : tracers names
191! * dBase(id)%trac(it)%keys%key(:): names  of keys associated to tracer dBase(id)%trac(it)%name
192! * dBase(id)%trac(it)%keys%val(:): values of keys associated to tracer dBase(id)%trac(it)%name
193!------------------------------------------------------------------------------------------------------------------------------
194  CHARACTER(LEN=*), INTENT(IN)  :: fnames(:)                         !--- Files names
195  CHARACTER(LEN=*), INTENT(IN)  :: snames(:)                         !--- Coma-deparated list of sections (one list each file)
196  CHARACTER(LEN=*), INTENT(IN)  :: modname                           !--- Calling routine name
197  INTEGER,  ALLOCATABLE :: ndb(:)                                    !--- Number of sections for each file
198  INTEGER,  ALLOCATABLE :: ixf(:)                                    !--- File index for each section of the expanded list
199  LOGICAL,  ALLOCATABLE :: lTg(:)                                    !--- Tagging tracers mask
200  CHARACTER(LEN=maxlen) :: fnm, snm
201  INTEGER               :: idb, i
202  LOGICAL :: ll
203!------------------------------------------------------------------------------------------------------------------------------
204  !=== READ THE REQUIRED SECTIONS
205  ll = strCount(snames, ',', ndb)                                    !--- Number of sections for each file
206  ALLOCATE(ixf(SUM(ndb)))
207  DO i=1, SIZE(fnames)                                               !--- Set %name, %keys
208    IF(test(readSections(fnames(i), snames(i), 'default'), lerr)) RETURN
209    ixf(1+SUM(ndb(1:i-1)):SUM(ndb(1:i))) = i                         !--- File index for each section of the expanded list
210  END DO
211  !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
212  DO idb=1,SIZE(dBase)                                               !--- LOOP ON THE LOADED SECTIONS
213  !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
214    fnm = fnames(ixf(idb)); snm = dBase(idb)%name                    !--- FILE AND SECTION NAMES
215    lerr = ANY([(dispTraSection('RAW CONTENT OF SECTION "'//TRIM(snm)//'"', snm, modname), idb=1, SIZE(dBase))])
216    IF(test(expandSection(dBase(idb)%trac, snm, fnm), lerr)) RETURN  !--- EXPAND NAMES ;  set %parent, %type, %component
217    CALL setGeneration   (dBase(idb)%trac)                           !---                 set %iGeneration,   %genOName
218    IF(test(checkTracers (dBase(idb)%trac, snm, fnm), lerr)) RETURN  !--- CHECK ORPHANS AND PHASES
219    IF(test(checkUnique  (dBase(idb)%trac, snm, fnm), lerr)) RETURN  !--- CHECK TRACERS UNIQUENESS
220    CALL expandPhases    (dBase(idb)%trac)                           !--- EXPAND PHASES ; set %phase
221    CALL sortTracers     (dBase(idb)%trac)                           !--- SORT TRACERS
222    lerr = ANY([(dispTraSection('EXPANDED CONTENT OF SECTION "'//TRIM(snm)//'"', snm, modname), idb=1, SIZE(dBase))])
223  !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
224  END DO
225  !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
226END FUNCTION feedDBase
227!------------------------------------------------------------------------------------------------------------------------------
228
229!------------------------------------------------------------------------------------------------------------------------------
230LOGICAL FUNCTION readSections(fnam,snam,defName) RESULT(lerr)
231!------------------------------------------------------------------------------------------------------------------------------
232  CHARACTER(LEN=*),           INTENT(IN) :: fnam                     !--- File name
233  CHARACTER(LEN=*),           INTENT(IN) :: snam                     !--- Coma-separated sections list
234  CHARACTER(LEN=*), OPTIONAL, INTENT(IN) :: defName                  !--- Special section (default values) name
235!------------------------------------------------------------------------------------------------------------------------------
236  TYPE(dataBase_type),   ALLOCATABLE :: tdb(:)
237  CHARACTER(LEN=maxlen), ALLOCATABLE :: sec(:)
238  INTEGER,               ALLOCATABLE ::  ix(:)
239  INTEGER :: n0, idb, ndb, i, j
240  LOGICAL :: ll
241!------------------------------------------------------------------------------------------------------------------------------
242  n0 = SIZE(dBase) + 1                                               !--- Index for next entry in the database
243  CALL readSections_all()                                            !--- Read all the sections of file "fnam"
244  ndb= SIZE(dBase)                                                   !--- Current number of sections in the database
245  IF(PRESENT(defName)) THEN                                          !--- Add default values to all the tracers
246    DO idb=n0,ndb; CALL addDefault(dBase(idb)%trac, defName); END DO !--- and remove the virtual tracer "defName"
247  END IF
248  ll = strParse(snam, ',', keys = sec)                               !--- Requested sections names
249  ix = strIdx(dBase(:)%name, sec(:))                                 !--- Indexes of requested sections in database
250  IF(test(checkList(sec, ix == 0, 'In file "'//TRIM(fnam)//'"','section(s)', 'missing'), lerr)) RETURN
251  tdb = dBase(:); dBase = [tdb(1:n0-1),tdb(PACK(ix, MASK=ix/=0))]    !--- Keep requested sections only
252
253CONTAINS
254
255!------------------------------------------------------------------------------------------------------------------------------
256SUBROUTINE readSections_all()
257!------------------------------------------------------------------------------------------------------------------------------
258  CHARACTER(LEN=maxlen), ALLOCATABLE ::  s(:), v(:)
259  TYPE(trac_type),       ALLOCATABLE :: tt(:)
260  TYPE(trac_type)       :: tmp
261  CHARACTER(LEN=1024)   :: str
262  CHARACTER(LEN=maxlen) :: secn
263  INTEGER               :: ierr, n
264!------------------------------------------------------------------------------------------------------------------------------
265  IF(.NOT.ALLOCATED(dBase)) ALLOCATE(dBase(0))
266  OPEN(90, FILE=fnam, FORM='formatted', STATUS='old')
267  DO; READ(90,'(a)', IOSTAT=ierr)str
268    IF(ierr    /= 0 ) EXIT                                           !--- Finished: error or end of file
269    IF(str(1:1)=='#') CYCLE                                          !--- Skip comments lines
270    CALL removeComment(str)                                          !--- Skip comments at the end of a line
271    IF(str     == '') CYCLE                                          !--- Skip empty line (probably at the end of the file)
272    IF(str(1:1)=='&') THEN                                           !=== SECTION HEADER LINE
273      ndb  = SIZE(dBase)                                             !--- Number of sections so far
274      secn = str(2:LEN_TRIM(str))//' '                               !--- Current section name
275      IF(ANY(dBase(:)%name == secn)) CYCLE                           !--- Already known section
276      IF(secn(1:7) == 'version')     CYCLE                           !--- Skip the "version" special section
277      ndb = ndb + 1                                                  !--- Extend database
278      ALLOCATE(tdb(ndb))
279      tdb(1:ndb-1)  = dBase
280      tdb(ndb)%name = secn
281      ALLOCATE(tdb(ndb)%trac(0))
282      CALL MOVE_ALLOC(FROM=tdb, TO=dBase)
283    ELSE                                                             !=== TRACER LINE
284      ll = strParse(str,' ', keys = s, vals = v, n = n)              !--- Parse <key>=<val> pairs
285      tt = dBase(ndb)%trac(:)
286      tmp%name = s(1); tmp%keys = keys_type(s(1), s(2:n), v(2:n))    !--- Set %name and %keys
287      dBase(ndb)%trac = [tt(:), tmp]
288      DEALLOCATE(tt)
289!      dBase(ndb)%trac = [dBase(ndb)%trac(:), tra(name=s(1), keys=keys_type(s(1), s(2:n), v(2:n)))]
290    END IF
291  END DO
292  CLOSE(90)
293
294END SUBROUTINE readSections_all
295!------------------------------------------------------------------------------------------------------------------------------
296
297END FUNCTION readSections
298!==============================================================================================================================
299
300
301!==============================================================================================================================
302SUBROUTINE addDefault(t, defName)
303!------------------------------------------------------------------------------------------------------------------------------
304! Purpose: Add the keys from virtual tracer named "defName" (if any) and remove this virtual tracer.
305!------------------------------------------------------------------------------------------------------------------------------
306  TYPE(trac_type), ALLOCATABLE, TARGET, INTENT(INOUT) :: t(:)
307  CHARACTER(LEN=*),                     INTENT(IN)    :: defName
308  INTEGER :: jd, it, k
309  TYPE(keys_type), POINTER :: ky
310  TYPE(trac_type), ALLOCATABLE :: tt(:)
311  jd = strIdx(t(:)%name, defName)
312  IF(jd == 0) RETURN
313  ky => t(jd)%keys
314  DO k = 1, SIZE(ky%key)                                             !--- Loop on the keys of the tracer named "defName"
315    CALL addKey_m(ky%key(k), ky%val(k), t(:)%keys)                   !--- Add key to all the tracers (no overwriting)
316  END DO
317  tt = [t(1:jd-1),t(jd+1:SIZE(t))]; CALL MOVE_ALLOC(FROM=tt, TO=t)   !--- Remove the virtual tracer named "defName"
318END SUBROUTINE addDefault
319!==============================================================================================================================
320
321!==============================================================================================================================
322SUBROUTINE subDefault(t, defName, lSubLocal)
323!------------------------------------------------------------------------------------------------------------------------------
324! Purpose: Substitute the keys from virtual tracer named "defName" (if any) and remove this virtual tracer.
325!          Substitute the keys locally (for the current tracer) if the flag "lSubLocal" is .TRUE.
326!------------------------------------------------------------------------------------------------------------------------------
327  TYPE(trac_type), ALLOCATABLE, TARGET, INTENT(INOUT) :: t(:)
328  CHARACTER(LEN=*),                     INTENT(IN)    :: defName
329  LOGICAL,                              INTENT(IN)    :: lSubLocal
330  INTEGER :: i0, it, ik
331  TYPE(keys_type), POINTER     :: k0, ky
332  TYPE(trac_type), ALLOCATABLE :: tt(:)
333  i0 = strIdx(t(:)%name, defName)
334  IF(i0 == 0) RETURN
335  k0 => t(i0)%keys
336  DO it = 1, SIZE(t); IF(it == i0) CYCLE                             !--- Loop on the tracers
337    ky => t(it)%keys
338
339    !--- Substitute in the values of <key>=<val> pairs the keys defined in the virtual tracer "defName"
340    DO ik = 1, SIZE(k0%key); CALL strReplace(ky%val, k0%key(ik), k0%val(ik), .TRUE.); END DO
341
342    IF(.NOT.lSubLocal) CYCLE
343    !--- Substitute in the values of <key>=<val> pairs the keys defined locally (in the current tracer)
344    DO ik = 1, SIZE(ky%key); CALL strReplace(ky%val, ky%key(ik), ky%val(ik), .TRUE.); END DO
345  END DO
346  tt = [t(1:i0-1),t(i0+1:SIZE(t))]; CALL MOVE_ALLOC(FROM=tt, TO=t)   !--- Remove the virtual tracer named "defName"
347
348END SUBROUTINE subDefault
349!==============================================================================================================================
350
351
352!==============================================================================================================================
353LOGICAL FUNCTION expandSection(tr, sname, fname) RESULT(lerr)
354!------------------------------------------------------------------------------------------------------------------------------
355! Purpose: Expand tracers and parents lists in the tracers descriptor "tra".
356! Note:  * The following keys are expanded, so are accessible explicitely using "%" operator: "parent" "type".
357!        * Default values are provided for these keys because they are necessary.
358!------------------------------------------------------------------------------------------------------------------------------
359  TYPE(trac_type), ALLOCATABLE, INTENT(INOUT) :: tr(:)                 !--- Tracer derived type vector
360  CHARACTER(LEN=*),             INTENT(IN)    :: sname
361  CHARACTER(LEN=*), OPTIONAL,   INTENT(IN)    :: fname
362  TYPE(trac_type),       ALLOCATABLE :: ttr(:)
363  CHARACTER(LEN=maxlen), ALLOCATABLE :: ta(:), pa(:)
364  CHARACTER(LEN=maxlen) :: msg1, modname
365  INTEGER :: it, nt, iq, nq, itr, ntr, ipr, npr, i
366  LOGICAL :: ll
367  modname = 'expandSection'
368  lerr = .FALSE.
369  nt = SIZE(tr)
370  nq = 0
371  !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
372  DO it = 1, nt    !=== GET TRACERS NB AFTER EXPANSION + NEEDED KEYS (name, parent, type)
373  !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
374    !--- Extract useful keys: parent name, type, component name
375    tr(it)%parent    = fgetKey(it, 'parent', tr(:)%keys,  tran0  )
376    tr(it)%type      = fgetKey(it, 'type'  , tr(:)%keys, 'tracer')
377    tr(it)%component = sname
378
379    !--- Determine the number of tracers and parents ; coherence checking
380    ll = strCount(tr(it)%name,   ',', ntr)
381    ll = strCount(tr(it)%parent, ',', npr)
382
383    !--- Tagging tracers only can have multiple parents
384    IF(test(npr/=1 .AND. TRIM(tr(it)%type)/='tag', lerr)) THEN
385      msg1 = 'Check section "'//TRIM(sname)//'"'
386      IF(PRESENT(fname)) msg1=TRIM(msg1)//' in file "'//TRIM(fname)//'"'
387      CALL msg(TRIM(msg1)//': "'//TRIM(tr(it)%name)//'" has several parents but is not a tag', modname); RETURN
388    END IF
389    nq = nq + ntr*npr                 
390  !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
391  END DO
392  !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
393  CALL delKey(['parent','type  '], tr)
394
395  ALLOCATE(ttr(nq))
396  iq = 1
397  !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
398  DO it = 1, nt                                                      !=== EXPAND TRACERS AND PARENTS NAMES LISTS
399  !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
400    ll = strParse(tr(it)%name,   ',', ta, n=ntr)                     !--- Number of tracers
401    ll = strParse(tr(it)%parent, ',', pa, n=npr)                     !--- Number of parents
402    DO ipr=1,npr                                                     !--- Loop on parents list elts
403      DO itr=1,ntr                                                   !--- Loop on tracers list elts
404        i = iq+itr-1+(ipr-1)*ntr
405        ttr(i)%name   = TRIM(ta(itr))
406        ttr(i)%parent = TRIM(pa(ipr))
407        ttr(i)%keys   = keys_type(ta(itr), tr(it)%keys%key, tr(it)%keys%val)
408      END DO
409    END DO
410    ttr(iq:iq+ntr*npr-1)%type      = tr(it)%type                     !--- Duplicating type
411    ttr(iq:iq+ntr*npr-1)%component = tr(it)%component                !--- Duplicating type
412    iq = iq + ntr*npr
413  !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
414  END DO
415  !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
416  DEALLOCATE(ta,pa)
417  CALL MOVE_ALLOC(FROM=ttr, TO=tr)
418
419END FUNCTION expandSection
420!==============================================================================================================================
421
422!==============================================================================================================================
423SUBROUTINE setGeneration(tr)
424!------------------------------------------------------------------------------------------------------------------------------
425! Purpose: Determine, for each tracer of "tr(:)":
426!   * %iGeneration: the generation number
427!   * %gen0Name:    the generation 0 ancestor name
428!------------------------------------------------------------------------------------------------------------------------------
429  TYPE(trac_type),     INTENT(INOUT) :: tr(:)                        !--- Tracer derived type vector
430  INTEGER                            :: iq, nq, ig
431  LOGICAL,               ALLOCATABLE :: lg(:)
432  CHARACTER(LEN=maxlen), ALLOCATABLE :: prn(:)
433!------------------------------------------------------------------------------------------------------------------------------
434  tr(:)%iGeneration = -1                                             !--- error if -1
435  nq = SIZE(tr, DIM=1)                                               !--- Number of tracers lines
436  lg = tr(:)%parent == tran0                                         !--- Flag for generation 0 tracers
437  WHERE(lg) tr(:)%iGeneration = 0                                    !--- Generation 0 tracers
438
439  !=== Determine generation for each tracer
440  ig=-1; prn = [tran0]
441  DO                                                                 !--- Update current generation flag
442    IF(ig/=-1) prn = PACK( tr(:)%name, MASK=tr(:)%iGeneration == ig)
443    lg(:) = [(ANY(prn(:) == tr(iq)%parent), iq=1, nq)]               !--- Current generation tracers flag
444    IF( ALL( .NOT. lg ) ) EXIT                                       !--- Empty current generation
445    ig = ig+1; WHERE(lg) tr(:)%iGeneration = ig
446  END DO
447  tr(:)%gen0Name = ancestor(tr)                                      !--- First generation ancestor name
448
449END SUBROUTINE setGeneration
450!==============================================================================================================================
451
452!==============================================================================================================================
453LOGICAL FUNCTION checkTracers(tr, sname, fname) RESULT(lerr)
454!------------------------------------------------------------------------------------------------------------------------------
455! Purpose:
456!   * check for orphan tracers (without known parent)
457!   * check wether the phases are known or not ("g"aseous, "l"iquid or "s"olid so far)
458!------------------------------------------------------------------------------------------------------------------------------
459  TYPE(trac_type),            INTENT(IN) :: tr(:)                    !--- Tracer derived type vector
460  CHARACTER(LEN=*),           INTENT(IN) :: sname                    !--- Section name
461  CHARACTER(LEN=*), OPTIONAL, INTENT(IN) :: fname                    !--- File name
462  CHARACTER(LEN=maxlen) :: mesg
463  CHARACTER(LEN=maxlen) :: bp(SIZE(tr, DIM=1)), pha                  !--- Bad phases list, phases of current tracer
464  CHARACTER(LEN=1) :: p
465  INTEGER :: ip, np, iq, nq
466!------------------------------------------------------------------------------------------------------------------------------
467  nq = SIZE(tr,DIM=1)                                                !--- Number of tracers lines
468  mesg = 'Check section "'//TRIM(sname)//'"'
469  IF(PRESENT(fname)) mesg=TRIM(mesg)//' in file "'//TRIM(fname)//'"'
470
471  !=== CHECK FOR ORPHAN TRACERS
472  IF(test(checkList(tr%name, tr%iGeneration==-1, mesg, 'tracers', 'orphan'), lerr)) RETURN
473
474  !=== CHECK PHASES
475  DO iq=1,nq; IF(tr(iq)%iGeneration/=0) CYCLE                        !--- Generation O only is checked
476    pha = fgetKey(iq, 'phases', tr(:)%keys, 'g')                     !--- Phases
477    np = LEN_TRIM(pha); bp(iq)=' '
478    DO ip=1,np; p = pha(ip:ip); IF(INDEX(known_phases,p)==0) bp(iq) = TRIM(bp(iq))//p; END DO
479    IF(TRIM(bp(iq)) /= '') bp(iq) = TRIM(tr(iq)%name)//': '//TRIM(bp(iq))
480  END DO
481  lerr = checkList(bp, tr%iGeneration==0 .AND. bp/='', mesg, 'tracers phases', 'unknown')
482END FUNCTION checkTracers
483!==============================================================================================================================
484
485!==============================================================================================================================
486LOGICAL FUNCTION checkUnique(tr, sname, fname) RESULT(lerr)
487!------------------------------------------------------------------------------------------------------------------------------
488! Purpose: Make sure that tracers are not repeated.
489!------------------------------------------------------------------------------------------------------------------------------
490  TYPE(trac_type),            INTENT(IN) :: tr(:)                    !--- Tracer derived type vector
491  CHARACTER(LEN=*),           INTENT(IN) :: sname                    !--- Section name
492  CHARACTER(LEN=*), OPTIONAL, INTENT(IN) :: fname                    !--- File name
493!------------------------------------------------------------------------------------------------------------------------------
494  INTEGER :: ip, np, iq, nq, k
495  LOGICAL, ALLOCATABLE  :: ll(:)
496  CHARACTER(LEN=maxlen) :: mesg, tnam, tdup(SIZE(tr,DIM=1))
497  CHARACTER(LEN=1)      :: p
498!------------------------------------------------------------------------------------------------------------------------------
499  mesg = 'Check section "'//TRIM(sname)//'"'
500  IF(PRESENT(fname)) mesg=TRIM(mesg)//' in file "'//TRIM(fname)//'"'
501  nq=SIZE(tr,DIM=1); lerr=.FALSE.                                    !--- Number of lines ; error flag
502  tdup(:) = ''
503  DO iq=1,nq; IF(tr(iq)%type == 'tag') CYCLE                         !--- Tags can be repeated
504    tnam = TRIM(tr(iq)%name)
505    ll = tr(:)%name==TRIM(tnam)                                      !--- Mask for current tracer name
506    IF(COUNT(ll)==1 ) CYCLE                                          !--- Tracer is not repeated
507    IF(tr(iq)%iGeneration>0) THEN
508      tdup(iq) = tnam                                                !--- gen>0: MUST be unique
509    ELSE
510      DO ip=1,nphases; p=known_phases(ip:ip)                         !--- Loop on known phases
511        !--- Number of appearances of the current tracer with known phase "p"
512        np = COUNT( PACK( [(INDEX(fgetKey(k, 'phases', tr(:)%keys, 'g'),p), k=1, nq)] /=0 , MASK=ll ) )
513        IF(np <=1) CYCLE
514        tdup(iq) = TRIM(tdup(iq))//TRIM(phases_names(ip))
515        IF(ANY(tdup(1:iq-1) == tdup(iq))) tdup(iq)=''                !--- Avoid repeating same messages
516      END DO
517    END IF
518    IF(tdup(iq) /= '') tdup(iq)=TRIM(tnam)//' in '//TRIM(tdup(iq))//' phase(s)'
519  END DO
520  lerr = checkList(tdup, tdup/='', mesg, 'tracers', 'duplicated')
521END FUNCTION checkUnique
522!==============================================================================================================================
523
524!==============================================================================================================================
525SUBROUTINE expandPhases(tr)
526!------------------------------------------------------------------------------------------------------------------------------
527! Purpose: Expand the phases in the tracers descriptor "tr". Phases are not repeated for a tracer, thanks to "checkUnique".
528!------------------------------------------------------------------------------------------------------------------------------
529  TYPE(trac_type), ALLOCATABLE, INTENT(INOUT) :: tr(:)               !--- Tracer derived type vector
530!------------------------------------------------------------------------------------------------------------------------------
531  TYPE(trac_type), ALLOCATABLE :: ttr(:)
532  INTEGER,   ALLOCATABLE ::  i0(:)
533  CHARACTER(LEN=maxlen)  :: nam, pha, trn
534  CHARACTER(LEN=1) :: p
535  INTEGER :: ip, np, iq, jq, nq, it, nt, nc, i, n
536  LOGICAL :: lTg, lEx
537!------------------------------------------------------------------------------------------------------------------------------
538  nq = SIZE(tr, DIM=1)
539  nt = 0
540  DO iq = 1, nq                                                      !--- GET THE NUMBER OF TRACERS
541    IF(tr(iq)%iGeneration /= 0) CYCLE                                !--- Only deal with generation 0 tracers
542    nc = COUNT(tr(:)%gen0Name==tr(iq)%name .AND. tr%iGeneration/=0)  !--- Number of childs of tr(iq)
543    tr(iq)%phase = fgetKey(iq, 'phases', tr(:)%keys)                 !--- Phases list      of tr(iq)
544    np = LEN_TRIM(tr(iq)%phase)                                      !--- Number of phases of tr(iq)
545    nt = nt + (1+nc) * np                                            !--- Number of tracers after expansion
546  END DO
547  ALLOCATE(ttr(nt))                                                  !--- Version  of "tr" after phases expansion
548  it = 1                                                             !--- Current "ttr(:)" index
549  DO iq = 1, nq                                                      !--- Loop on "tr(:)" indexes
550    lTg = tr(iq)%type=='tag'                                         !--- Current tracer is a tag
551    i0 = strFind(tr(:)%name, TRIM(tr(iq)%gen0Name), n)               !--- Indexes of first generation ancestor copies
552    np = SUM([( LEN_TRIM(tr(i0(i))%phase),i=1,n )], 1)               !--- Number of phases for current tracer tr(iq)
553    lEx = np>1                                                       !--- Phase suffix only required if phases number is > 1
554    IF(lTg) lEx = lEx .AND. tr(iq)%iGeneration>0                     !--- No phase suffix for generation 0 tags
555    DO i=1,n                                                         !=== LOOP ON GENERATION 0 ANCESTORS
556      jq = i0(i)                                                     !--- tr(jq): ith tracer with same gen 0 ancestor as tr(iq)
557      IF(tr(iq)%iGeneration==0) jq=iq                                !--- Generation 0: count the current tracer phases only
558      pha = tr(jq)%phase                                             !--- Phases list for tr(jq)
559      DO ip = 1, LEN_TRIM(pha)                                       !=== LOOP ON PHASES LISTS
560        p = pha(ip:ip)
561        trn = TRIM(tr(iq)%name); nam = trn                           !--- Tracer name (regular case)
562        IF(lTg) nam = TRIM(tr(iq)%parent)                            !--- Parent name (tagging case)
563        IF(lEx) nam = addPhase(nam, p )                              !--- Phase extension needed
564        IF(lTg) nam = TRIM(nam)//'_'//TRIM(trn)                      !--- <parent>_<name> for tags
565        ttr(it) = tr(iq)                                             !--- Same <key>=<val> pairs
566        ttr(it)%name      = TRIM(nam)                                !--- Name with possibly phase suffix
567        ttr(it)%keys%name = TRIM(nam)                                !--- Name inside the keys decriptor
568        ttr(it)%phase     = p                                        !--- Single phase entry
569        IF(lEx .AND. tr(iq)%iGeneration>0) THEN
570          ttr(it)%parent   = addPhase(ttr(it)%parent,   p)
571          ttr(it)%gen0Name = addPhase(ttr(it)%gen0Name, p)
572        END IF
573        it = it+1
574      END DO
575      IF(tr(iq)%iGeneration==0) EXIT                                 !--- Break phase loop for gen 0
576    END DO
577  END DO
578  CALL MOVE_ALLOC(FROM=ttr, TO=tr)
579  CALL delKey(['phases'],tr)                                         !--- Remove few keys entries
580
581END SUBROUTINE expandPhases
582!==============================================================================================================================
583
584!==============================================================================================================================
585SUBROUTINE sortTracers(tr)
586!------------------------------------------------------------------------------------------------------------------------------
587! Purpose: Sort tracers:
588!  * Put water at the beginning of the vector, in the "known_phases" order.
589!  * lGrowGen == T: in ascending generations numbers.
590!  * lGrowGen == F: tracer + its childs sorted by growing generation, one after the other.
591!   TO BE ADDED IF NECESSARY: HIGHER MOMENTS AT THE END
592!------------------------------------------------------------------------------------------------------------------------------
593  TYPE(trac_type), ALLOCATABLE, INTENT(INOUT) :: tr(:)         !--- Tracer derived type vector
594!  TYPE(trac_type), ALLOCATABLE :: ttr(:)
595  INTEGER,         ALLOCATABLE :: iy(:), iz(:)
596  INTEGER :: ig, ng, iq, jq, ip, nq, n, ix(SIZE(tr)), k
597!------------------------------------------------------------------------------------------------------------------------------
598  nq = SIZE(tr)
599  DO ip = nphases, 1, -1
600    iq = strIdx(tr(:)%name, addPhase('H2O', ip))
601    IF(iq == 0) CYCLE
602    tr = tr([iq, 1:iq-1, iq+1:nq])
603!    tr(:)%name = nam
604  END DO
605  IF(lSortByGen) THEN
606    iq = 1
607    ng = MAXVAL(tr(:)%iGeneration, MASK=.TRUE., DIM=1)               !--- Number of generations
608    DO ig = 0, ng                                                    !--- Loop on generations
609      iy = PACK([(k, k=1, nq)], MASK=tr(:)%iGeneration==ig)          !--- Generation ig tracers indexes
610      n = SIZE(iy)
611      ix(iq:iq+n-1) = iy                                             !--- Stack growing generations idxs
612      iq = iq + n
613    END DO
614  ELSE
615    iq = 1
616    DO jq = 1, nq                                                    !--- Loop on generation 0 tracers
617      IF(tr(jq)%iGeneration /= 0) CYCLE                              !--- Skip generations /= 0
618      ix(iq) = jq                                                    !--- Generation 0 ancestor index first
619      iq = iq + 1                                                    !--- Next "iq" for next generations tracers
620      iy = strFind(tr(:)%gen0Name, TRIM(tr(jq)%name))                !--- Indexes of "tr(jq)" childs in "tr(:)"
621      ng = MAXVAL(tr(iy)%iGeneration, MASK=.TRUE., DIM=1)            !--- Number of generations of the "tr(jq)" family
622      DO ig = 1, ng                                                  !--- Loop   on generations of the "tr(jq)" family
623        iz = find(tr(iy)%iGeneration, ig, n)                         !--- Indexes of the tracers "tr(iy(:))" of generation "ig"
624        ix(iq:iq+n-1) = iy(iz)                                       !--- Same indexes in "tr(:)"
625        iq = iq + n
626      END DO
627    END DO
628  END IF
629  tr = tr(ix)                                                        !--- Reorder the tracers
630END SUBROUTINE sortTracers
631!==============================================================================================================================
632
633!==============================================================================================================================
634LOGICAL FUNCTION mergeTracers(sections, tr) RESULT(lerr)
635  TYPE(dataBase_type),  TARGET, INTENT(IN)  :: sections(:)
636  TYPE(trac_type), ALLOCATABLE, INTENT(OUT) ::       tr(:)
637  TYPE(trac_type), POINTER ::   t1(:),   t2(:)
638  INTEGER,     ALLOCATABLE :: ixct(:), ixck(:)
639  INTEGER :: is, k1, k2, nk2, i1, i2, nt2
640  CHARACTER(LEN=maxlen) :: s1, v1, v2, tnam, knam, modname
641  modname = 'mergeTracers'
642  lerr = .FALSE.
643  t1 => sections(1)%trac(:)                                          !--- Alias: first tracers section
644  tr = t1
645  !----------------------------------------------------------------------------------------------------------------------------
646  DO is=2,SIZE(sections)                                             !=== SEVERAL SECTIONS: MERGE THEM
647  !----------------------------------------------------------------------------------------------------------------------------
648    t2  => sections(is)%trac(:)                                      !--- Alias: current tracers section
649    nt2  = SIZE(t2(:), DIM=1)                                        !--- Number of tracers in section
650    ixct = strIdx(t1(:)%name, t2(:)%name)                            !--- Indexes of common tracers
651    tr = [tr, PACK(t2, MASK= ixct==0)]                               !--- Append with new tracers
652    IF( ALL(ixct == 0) ) CYCLE                                       !--- No common tracers => done
653    CALL msg('Tracers defined in previous sections and duplicated in "'//TRIM(sections(is)%name)//'":', modname)
654    CALL msg(t1(PACK(ixct, MASK = ixct/=0))%name, modname, nmax=128) !--- Display duplicates (the 128 first at most)
655    !--------------------------------------------------------------------------------------------------------------------------
656    DO i2=1,nt2; tnam = TRIM(t2(i2)%name)                            !=== LOOP ON COMMON TRACERS
657    !--------------------------------------------------------------------------------------------------------------------------
658      i1 = ixct(i2); IF(i1 == 0) CYCLE                               !--- Idx in t1(:) ; skip new tracers
659
660      !=== CHECK WETHER ESSENTIAL KEYS ARE IDENTICAL OR NOT
661      s1=' of "'//TRIM(tnam)//'" in "'//TRIM(sections(is)%name)//'" not matching previous value'
662     
663      IF(test(fmsg('Parent name'//TRIM(s1), modname, t1(i1)%parent      /= t2(i2)%parent),      lerr)) RETURN
664      IF(test(fmsg('Type'       //TRIM(s1), modname, t1(i1)%type        /= t2(i2)%type),        lerr)) RETURN
665      IF(test(fmsg('Generation' //TRIM(s1), modname, t1(i1)%iGeneration /= t2(i2)%iGeneration), lerr)) RETURN
666
667      !=== APPEND <key>=<val> PAIRS NOT PREVIOULSLY DEFINED
668      nk2  = SIZE(t2(i2)%keys%key(:))                                !--- Keys number in current section
669      ixck = strIdx(t1(i1)%keys%key(:), t2(i2)%keys%key(:))          !--- Common keys indexes
670
671      !=== APPEND NEW KEYS
672      tr(i1)%keys%key = [ tr(i1)%keys%key, PACK(tr(i2)%keys%key, MASK = ixck==0)]
673      tr(i1)%keys%val = [ tr(i1)%keys%val, PACK(tr(i2)%keys%val, MASK = ixck==0)]
674
675      !--- KEEP TRACK OF THE COMPONENTS NAMES
676      tr(i1)%component = TRIM(tr(i1)%component)//','//TRIM(tr(i2)%component)
677
678      !--- SELECT COMMON TRACERS WITH DIFFERING KEYS VALUES (PREVIOUS VALUE IS KEPT)
679      DO k2=1,nk2
680        k1 = ixck(k2); IF(k1 == 0) CYCLE
681        IF(t1(i1)%keys%val(k1) == t2(i2)%keys%val(k2)) ixck(k2)=0
682      END DO
683      IF(ALL(ixck==0)) CYCLE                                         !--- No identical keys with /=values
684
685      !--- DISPLAY INFORMATION: OLD VALUES ARE KEPT FOR THE KEYS FOUND IN PREVIOUS AND CURRENT SECTIONS
686      CALL msg('Key(s)'//TRIM(s1), modname)
687      DO k2 = 1, nk2                                                 !--- Loop on keys found in both t1(:) and t2(:)
688        knam = t2(i2)%keys%key(k2)                                   !--- Name of the current key
689        k1 = ixck(k2)                                                !--- Corresponding index in t1(:)
690        IF(k1 == 0) CYCLE                                            !--- New keys are skipped
691        v1 = t1(i1)%keys%val(k1); v2 = t2(i2)%keys%val(k2)           !--- Key values in t1(:) and t2(:)
692        CALL msg(' * '//TRIM(knam)//'='//TRIM(v2)//' ; previous value kept:'//TRIM(v1), modname)
693      END DO
694      !------------------------------------------------------------------------------------------------------------------------
695    END DO
696    !--------------------------------------------------------------------------------------------------------------------------
697  END DO
698  CALL sortTracers(tr)
699
700END FUNCTION mergeTracers
701!==============================================================================================================================
702
703!==============================================================================================================================
704LOGICAL FUNCTION cumulTracers(sections, tr) RESULT(lerr)
705  TYPE(dataBase_type),  TARGET, INTENT(IN)  :: sections(:)
706  TYPE(trac_type), ALLOCATABLE, INTENT(OUT) ::       tr(:)
707  TYPE(trac_type), POINTER     :: t1(:), t2(:)
708  INTEGER,   ALLOCATABLE :: nt(:)
709  CHARACTER(LEN=maxlen)  :: tnam, tnam_new
710  INTEGER :: iq, nq, is, ns, nsec
711  lerr = .FALSE.                                                     !--- Can't fail ; kept to match "mergeTracer" interface.
712  nsec =  SIZE(sections)
713  tr = [(      sections(is)%trac(:) , is=1, nsec )]                  !--- Concatenated tracers vector
714  nt = [( SIZE(sections(is)%trac(:)), is=1, nsec )]                  !--- Number of tracers in each section
715  !----------------------------------------------------------------------------------------------------------------------------
716  DO is=1, nsec                                                      !=== LOOP ON SECTIONS
717  !----------------------------------------------------------------------------------------------------------------------------
718    t1 => sections(is)%trac(:)
719    !--------------------------------------------------------------------------------------------------------------------------
720    DO iq=1, nt(is)                                                  !=== LOOP ON TRACERS
721    !--------------------------------------------------------------------------------------------------------------------------
722      tnam = TRIM(t1(iq)%name)                                       !--- Original name
723      IF(COUNT(t1%name == tnam) == 1) CYCLE                          !--- Current tracer is not duplicated: finished
724      tnam_new = TRIM(tnam)//'_'//TRIM(sections(is)%name)            !--- Same with section extension
725      nq = SUM(nt(1:is-1))                                           !--- Number of tracers in previous sections
726      ns = nt(is)                                                    !--- Number of tracers in the current section
727      tr(iq + nq)%name = TRIM(tnam_new)                              !--- Modify tracer name
728      WHERE(tr(1+nq:ns+nq)%parent==tnam) tr(1+nq:ns+nq)%parent=tnam_new  !--- Modify parent name
729    !--------------------------------------------------------------------------------------------------------------------------
730    END DO
731  !----------------------------------------------------------------------------------------------------------------------------
732  END DO
733  !----------------------------------------------------------------------------------------------------------------------------
734  CALL sortTracers(tr)
735END FUNCTION cumulTracers
736!==============================================================================================================================
737
738!==============================================================================================================================
739SUBROUTINE setDirectKeys(tr)
740  TYPE(trac_type), INTENT(INOUT) :: tr(:)
741
742  !--- Update %iqParent, %iqDescen, %nqDescen, %nqChilds
743  CALL indexUpdate(tr)
744
745  !--- Extract some direct-access keys
746!  DO iq = 1, SIZE(tr)
747!    tr(iq)%keys%<key> = getKey_prv(it, "<key>", tr%keys, <default_value> )
748!  END DO
749END SUBROUTINE setDirectKeys
750!==============================================================================================================================
751
752!==============================================================================================================================
753LOGICAL FUNCTION dispTraSection(message, sname, modname) RESULT(lerr)
754  CHARACTER(LEN=*), INTENT(IN) :: message, sname, modname
755  INTEGER :: idb, iq, nq
756  INTEGER, ALLOCATABLE :: hadv(:), vadv(:)
757  CHARACTER(LEN=maxlen), ALLOCATABLE :: phas(:)
758  TYPE(trac_type), POINTER :: tm(:)
759  lerr = .FALSE.
760  idb = strIdx(dBase(:)%name, sname); IF(idb == 0) RETURN
761  tm => dBase(idb)%trac
762  nq = SIZE(tm)
763  !--- BEWARE ! Can't always already use the "getKeyByName" functions.
764  !             Names must first include the phases for tracers defined on multiple lines.
765  hadv = str2int([(fgetKey(iq, 'hadv',  tm(:)%keys, '10'), iq=1, nq)])
766  vadv = str2int([(fgetKey(iq, 'vadv',  tm(:)%keys, '10'), iq=1, nq)])
767  phas =         [(fgetKey(iq, 'phases',tm(:)%keys, 'g' ), iq=1, nq)]
768  CALL msg(TRIM(message)//':', modname)
769  IF(tm(1)%parent == '') THEN
770    IF(test(dispTable('iiiss',   ['iq    ','hadv  ','vadv  ','name  ','phase '], cat(tm%name, phas), cat([(iq, iq=1, nq)], &
771                                            hadv,    vadv),                 sub=modname), lerr)) RETURN
772  ELSE
773    IF(test(dispTable('iiissis', ['iq    ','hadv  ','vadv  ','name  ','parent','igen  ','phase '], cat(tm%name, tm%parent, &
774           tm%phase), cat([(iq, iq=1, nq)], hadv,    vadv, tm%iGeneration), sub=modname), lerr)) RETURN
775  END IF
776END FUNCTION dispTraSection
777!==============================================================================================================================
778!==============================================================================================================================
779
780
781!==============================================================================================================================
782!== CREATE A SCALAR ALIAS OF THE COMPONENT OF THE TRACERS DESCRIPTOR "t" NAMED "tname" ========================================
783!==============================================================================================================================
784FUNCTION aliasTracer(tname, t) RESULT(out)
785  TYPE(trac_type),         POINTER    :: out
786  CHARACTER(LEN=*),        INTENT(IN) :: tname
787  TYPE(trac_type), TARGET, INTENT(IN) :: t(:)
788  INTEGER :: it
789  it = strIdx(t(:)%name, tname)
790  out => NULL(); IF(it /= 0) out => t(it)
791END FUNCTION aliasTracer
792!------------------------------------------------------------------------------------------------------------------------------
793
794
795!==============================================================================================================================
796!=== FROM A LIST OF INDEXES OR NAMES, CREATE A SUBSET OF THE TRACERS DESCRIPTORS LIST "trac" ==================================
797!==============================================================================================================================
798FUNCTION trSubset_Indx(trac,idx) RESULT(out)
799  TYPE(trac_type), ALLOCATABLE             ::  out(:)
800  TYPE(trac_type), ALLOCATABLE, INTENT(IN) :: trac(:)
801  INTEGER,                      INTENT(IN) ::  idx(:)
802  out = trac(idx)
803  CALL indexUpdate(out)
804END FUNCTION trSubset_Indx
805!------------------------------------------------------------------------------------------------------------------------------
806FUNCTION trSubset_Name(trac,nam) RESULT(out)
807  TYPE(trac_type), ALLOCATABLE             ::  out(:)
808  TYPE(trac_type), ALLOCATABLE, INTENT(IN) :: trac(:)
809  CHARACTER(LEN=*),             INTENT(IN) ::  nam(:)
810  out = trac(strIdx(trac(:)%name, nam))
811  CALL indexUpdate(out)
812END FUNCTION trSubset_Name
813!------------------------------------------------------------------------------------------------------------------------------
814
815
816!==============================================================================================================================
817!=== CREATE THE SUBSET OF THE TRACERS DESCRIPTORS LIST "trac" HAVING THE FIRST GENERATION ANCESTOR NAMED "nam" ================
818!==============================================================================================================================
819FUNCTION trSubset_gen0Name(trac,nam) RESULT(out)
820  TYPE(trac_type), ALLOCATABLE             ::  out(:)
821  TYPE(trac_type), ALLOCATABLE, INTENT(IN) :: trac(:)
822  CHARACTER(LEN=*),             INTENT(IN) ::  nam
823  out = trac(strFind(delPhase(trac(:)%gen0Name), nam))
824  CALL indexUpdate(out)
825END FUNCTION trSubset_gen0Name
826!------------------------------------------------------------------------------------------------------------------------------
827
828
829!==============================================================================================================================
830!=== UPDATE THE INDEXES iqParent, iqDescend AND iGeneration IN THE TRACERS DESCRIPTOR LIST "tr" (USEFULL FOR SUBSETS) =========
831!==============================================================================================================================
832SUBROUTINE indexUpdate(tr)
833  TYPE(trac_type), INTENT(INOUT) :: tr(:)
834  INTEGER :: iq, ig, ng, ngen
835  INTEGER, ALLOCATABLE :: ix(:)
836  tr(:)%iqParent = strIdx( tr(:)%name, tr(:)%parent )                !--- Parent index
837  ngen = MAXVAL(tr(:)%iGeneration, MASK=.TRUE.)
838  DO iq = 1, SIZE(tr)
839    ng = tr(iq)%iGeneration                                          !--- Generation of the current tracer
840    ix = idxAncestor(tr, igen = ng); ix = PACK(ix, ix/=0)            !--- Indexes of the tracers with ancestor tr(iq)
841    !--- Childs indexes in growing generation order
842    tr(iq)%iqDescen = [( PACK(ix, MASK = tr(ix)%iGeneration == ig), ig = ng+1, ngen)]
843    tr(iq)%nqDescen =     SUM(  [( COUNT(tr(ix)%iGeneration == ig), ig = ng+1, ngen)] )
844    tr(iq)%nqChilds =              COUNT(tr(ix)%iGeneration == ng+1)
845  END DO
846END SUBROUTINE indexUpdate
847!------------------------------------------------------------------------------------------------------------------------------
848 
849 
850!==============================================================================================================================
851!=== READ FILE "fnam" TO APPEND THE "dBase" TRACERS DATABASE WITH AS MUCH SECTIONS AS PARENTS NAMES IN "isot(:)%parent":   ====
852!===  * Each section dBase(i)%name contains the isotopes "dBase(i)%trac(:)" descending on "dBase(i)%name"="iso(i)%parent"  ====
853!===  * For each isotopes class, the <key>=<val> vector of each tracer is moved into the isotopes descriptor "isot"        ====
854!=== NOTES:                                                                                                                ====
855!===  * Most of the "isot" components have been defined in the calling routine (initIsotopes):                             ====
856!===      parent,  nzone, zone(:),  niso, keys(:)%name,  ntiso, trac(:),  nphas, phas,  iqTraPha(:,:),  itZonPhi(:,:)      ====
857!===  * Same syntax for isotopes file and "tracer.def": a tracers section contains one line for each of its isotopes       ====
858!===  * Each tracers section can contain a "params" virtual isotope line of isotopes parameters default values             ====
859!===  * In case keys are found both in the "params" section and the "*.def" file, the later value is retained              ====
860!===  * On each isotope line, defined keys can be used for other keys defintions (single level depth substitution)         ====
861!===  * The routine gives an error if a required isotope is not available in the database stored in "fnam"                 ====
862!==============================================================================================================================
863LOGICAL FUNCTION readIsotopesFile(fnam, isot) RESULT(lerr)
864  CHARACTER(LEN=*),        INTENT(IN)    :: fnam                     !--- Input file name
865  TYPE(isot_type), TARGET, INTENT(INOUT) :: isot(:)                  !--- Isotopes descriptors (field "prnt" must be defined !)
866  INTEGER :: ik, is, it, idb, nk0, i, iis
867  INTEGER :: nk, ns, nt, ndb, nb0, i0
868  CHARACTER(LEN=maxlen), POINTER     :: k(:), v(:), k0(:), v0(:)
869  CHARACTER(LEN=maxlen), ALLOCATABLE :: vals(:)
870  CHARACTER(LEN=maxlen)              :: val, modname
871  TYPE(keys_type),           POINTER ::   ky(:)
872  TYPE(trac_type),           POINTER ::   tt(:), t
873  TYPE(dataBase_type),   ALLOCATABLE ::  tdb(:)
874  LOGICAL,               ALLOCATABLE :: liso(:)
875  modname = 'readIsotopesFile'
876
877  !--- THE INPUT FILE MUST BE PRESENT
878  IF(test(fmsg('Missing isotopes parameters file "'//TRIM(fnam)//'"', modname, testFile(fnam)),lerr)) RETURN
879
880  !--- READ THE FILE SECTIONS, ONE EACH PARENT TRACER
881  nb0 = SIZE(dBase, DIM=1)+1                                         !--- Next database element index
882  IF(test(readSections(fnam,strStack(isot(:)%parent,',')),lerr)) RETURN !--- Read sections, one each parent tracer
883  ndb = SIZE(dBase, DIM=1)                                           !--- Current database size
884  DO idb = nb0, ndb
885   iis = idb-nb0+1
886
887    !--- GET FEW GLOBAL KEYS FROM "def" FILES AND ADD THEM TO THE 'params' SECTION
888    CALL addKeysFromDef(dBase(idb)%trac, 'params')
889
890    !--- SUBSTITUTE THE KEYS DEFINED IN THE 'params' VIRTUAL TRACER ; SUBSTITUTE LOCAL KEYS ; REMOVE 'params' VIRTUAL TRACER
891    CALL subDefault(dBase(idb)%trac, 'params', .TRUE.)
892
893    tt => dBase(idb)%trac
894
895    !--- REDUCE THE EXPRESSIONS TO OBTAIN SCALARS AND TRANSFER THEM TO THE "isot" ISOTOPES DESCRIPTORS VECTOR
896    DO it = 1, SIZE(dBase(idb)%trac)
897      t => dBase(idb)%trac(it)
898      is = strIdx(isot(iis)%keys(:)%name, t%name)                    !--- Index in "isot(iis)%keys(:)%name" of isotope "t%name"
899      IF(is == 0) CYCLE
900      liso = reduceExpr(t%keys%val, vals)                            !--- Reduce expressions (for substituted variables)
901      IF(test(ANY(liso), lerr)) RETURN                               !--- Some non-numerical elements were found
902      isot(iis)%keys(is)%key = PACK(t%keys%key, MASK=.NOT.liso)
903      isot(iis)%keys(is)%val = PACK(  vals,     MASK=.NOT.liso)
904    END DO
905
906    !--- CHECK FOR MISSING ISOTOPES (NO KEYS ALLOCATED)
907    liso = [( ALLOCATED(isot(iis)%keys(is)%key), is=1, SIZE(isot(iis)%keys) )]
908    IF(test(checkList(isot(iis)%keys(:)%name, .NOT.liso, &
909      'Check file "'//TRIM(fnam)//'" in section "'//TRIM(dBase(idb)%name)//'"', 'isotopes', 'missing'),lerr)) RETURN
910  END DO
911
912  !--- CLEAN THE DATABASE ENTRIES
913  IF(nb0 == 1) THEN
914    DEALLOCATE(dBase); ALLOCATE(dBase(0))
915  ELSE
916    ALLOCATE(tdb(nb0-1)); tdb(1:nb0-1)=dBase(1:nb0-1); CALL MOVE_ALLOC(FROM=tdb, TO=dBase)
917  END IF
918
919  !--- GET THE isoCheck ENTRY FROM THE *.DEF FILES (MIGHT BE CHANGED TO A CLASS-DEPENDANT KEYWORD)
920  CALL get_in('ok_iso_verif', isot(strIdx(isot%parent, 'H2O'))%check, .FALSE.)
921
922  lerr = dispIsotopes(isot, 'Isotopes parameters read from file', modname)
923
924END FUNCTION readIsotopesFile
925!==============================================================================================================================
926
927!==============================================================================================================================
928!=== IF ISOTOPES (2ND GENERATION TRACERS) ARE DETECTED:                                                                     ===
929!===    * COMPUTE MOST OF THE RELATED QUANTITIES ("isot" COMPONENTS).                                                       ===
930!===    * COMPUTE FEW ISOTOPES-DEDICATED "trac" COMPONENTS                                                                  ===
931!===    * CALL readIsotopesFile TO GET PHYSICAL QUANTITIES (<key>=<val> PAIRS)                                              ===
932!===      NOTE: THIS IS DONE HERE (IN A ROUTINE CALLED BY THE DYNAMIC), BECAUSE THE DYNAMIC NEEDS FEW PHYSICAL PARAMETERS.  ===
933!==============================================================================================================================
934SUBROUTINE initIsotopes(trac, isot)
935  TYPE(trac_type), ALLOCATABLE, TARGET, INTENT(INOUT) :: trac(:)
936  TYPE(isot_type), ALLOCATABLE, TARGET, INTENT(INOUT) :: isot(:)
937  CHARACTER(LEN=maxlen), ALLOCATABLE :: p(:), str(:)                 !--- Temporary storage
938  CHARACTER(LEN=maxlen) :: iname
939  CHARACTER(LEN=1)   :: ph                                           !--- Phase
940  INTEGER :: nbIso, ic, ip, iq, it, iz
941  LOGICAL, ALLOCATABLE :: ll(:)                                      !--- Mask
942  TYPE(trac_type), POINTER   ::  t(:), t1
943  TYPE(isot_type), POINTER   ::  i
944
945  t => trac
946
947  p = PACK(delPhase(t%parent), MASK = t%type=='tracer' .AND. t%iGeneration==1) !--- Parents of generation 1 isotopes
948  CALL strReduce(p, nbIso)
949  ALLOCATE(isot(nbIso))
950
951  IF(nbIso==0) RETURN                                                !=== NO ISOTOPES: FINISHED
952
953  !--- ISOTOPES RELATED VARIABLES ; NULL OR EMPTY IF NO ISOTOPES
954  isot(:)%parent = p
955  DO ic = 1, SIZE(p)                                                 !--- Loop on isotopes classes
956    i => isot(ic)
957    iname = i%parent                                                 !--- Current isotopes class name (parent tracer name)
958
959    !=== Isotopes childs of tracer "iname": mask, names, number (same for each phase of "iname")
960    ll = t(:)%type=='tracer' .AND. delPhase(t(:)%parent) == iname .AND. t(:)%phase == 'g'
961    str = PACK(delPhase(t(:)%name), MASK = ll)                       !--- Effectively found isotopes of "iname"
962    i%niso = SIZE(str)                                               !--- Number of "effectively found isotopes of "iname"
963    ALLOCATE(i%keys(i%niso))
964    FORALL(it = 1:i%niso) i%keys(it)%name = str(it)
965
966    !=== Geographic tagging tracers descending on tracer "iname": mask, names, number
967    ll = t(:)%type=='tag'    .AND. delPhase(t(:)%gen0Name) == iname .AND. t(:)%iGeneration == 2
968    i%zone = PACK(strTail(t(:)%name,'_',lFirst=.TRUE.), MASK = ll)   !--- Tagging zones names  for isotopes category "iname"
969    CALL strReduce(i%zone)
970    i%nzone = SIZE(i%zone)                                           !--- Tagging zones number for isotopes category "iname"
971
972    !=== Geographic tracers of the isotopes childs of tracer "iname" (same for each phase of "iname")
973    !    NOTE: One might like to create a similar variable for 2nd generation tagging tracers (tagging the gen1 tracers)
974    str = PACK(delPhase(t(:)%name), MASK=ll)
975    CALL strReduce(str)
976    i%ntiso = i%niso + SIZE(str)                                     !--- Number of isotopes + their geographic tracers [ntraciso]
977    ALLOCATE(i%trac(i%ntiso))
978    FORALL(it = 1:i%niso) i%trac(it) = i%keys(it)%name
979    FORALL(it = i%niso+1:i%ntiso) i%trac(it) = str(it-i%niso)
980
981    !=== Phases for tracer "iname"
982    i%phase = ''
983    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
984    i%nphas = LEN_TRIM(i%phase)                                       !--- Equal to "nqo" for water
985
986    !=== Tables giving the index in a table of effectively found items for each dynamical tracer (1<=iq<=nqtot)
987    DO iq = 1, SIZE(t)
988      t1 => trac(iq)
989      IF(delPhase(t1%gen0Name)/=iname .OR. t1%iGeneration==0) CYCLE  !--- Only deal with tracers descending on "iname"
990      t1%iso_iGroup = ic                                             !--- Isotopes family       idx in list "isotopes(:)%parent"
991      t1%iso_iName  = strIdx(i%trac, delPhase(strHead(t1%name,'_'))) !--- Current isotope       idx in effective isotopes list
992      t1%iso_iZone  = strIdx(i%zone,          strTail(t1%name,'_') ) !--- Current isotope zone  idx in effective zones    list
993      t1%iso_iPhase =  INDEX(i%phase,TRIM(t1%phase))                 !--- Current isotope phase idx in effective phases   list
994      IF(t1%iGeneration /= 2) t1%iso_iZone = 0                       !--- Skip possible generation 1 tagging tracers
995    END DO
996
997    !=== Table used to get iq (index in dyn array, size nqtot) from the isotope and phase indexes ; the full isotopes list
998    !    (including tagging tracers) is sorted this way:  iso1, iso2, ..., iso1_zone1, iso2_zone1, ..., iso1_zoneN, iso2_zoneN
999    i%iqTraPha = RESHAPE( [( (strIdx(t%name,  addPhase(i%trac(it),i%phase(ip:ip))),    it=1, i%ntiso), ip=1, i%nphas)], &
1000                         [i%ntiso, i%nphas] )
1001    !=== Table used to get ix (index in tagging tracers isotopes list, size ntiso) from the zone and isotope indexes
1002    i%itZonIso = RESHAPE( [( (strIdx(i%trac(:), TRIM(i%trac(it))//'_'//TRIM(i%zone(iz))), iz=1, i%nzone), it=1, i%niso )], &
1003                         [i%nzone, i%niso] )
1004  END DO
1005
1006  !=== READ PHYSICAL PARAMETERS FROM "isotopes_params.def" FILE
1007  !    DONE HERE, AND NOT ONLY IN "infotrac_phy", BECAUSE SOME PHYSICAL PARAMS ARE NEEDED FOR RESTARTS (tnat AND alpha_ideal)
1008  IF(readIsotopesFile('isotopes_params.def',isot)) CALL abort_gcm(modname,'Problem when reading isotopes parameters',1)
1009
1010END SUBROUTINE initIsotopes
1011!==============================================================================================================================
1012
1013
1014!==============================================================================================================================
1015LOGICAL FUNCTION dispIsotopes(ides, message, modname) RESULT(lerr)
1016  TYPE(isot_type),  INTENT(IN) :: ides(:)                            !--- Isotopes descriptor vector
1017  CHARACTER(LEN=*), INTENT(IN) :: message                            !--- Message to display
1018  CHARACTER(LEN=*), INTENT(IN) :: modname                            !--- Calling subroutine name
1019  INTEGER :: ik, nk, ip, it, nt
1020  CHARACTER(LEN=maxlen) :: prf
1021  CHARACTER(LEN=maxlen), ALLOCATABLE :: ttl(:), val(:,:)
1022  CALL msg(TRIM(message)//':', modname)
1023  DO ip = 1, SIZE(ides)                                              !--- Loop on parents tracers
1024    nk = SIZE(ides(ip)%keys(1)%key)                                  !--- Same keys for each isotope
1025    nt = SIZE(ides(ip)%keys)                                         !--- Number of isotopes
1026    prf = 'i'//REPEAT('s',nk+1)                                      !--- Profile for table printing
1027    ALLOCATE(ttl(nk+2), val(nt,nk+1))
1028    ttl(1:2) = ['iq  ','name']; ttl(3:nk+2) = ides(ip)%keys(1)%key(:)!--- Titles line with keys names
1029    val(:,1) = ides(ip)%keys(:)%name                                 !--- Values table 1st column: isotopes names 
1030    DO ik = 1, nk
1031      DO it = 1, nt
1032        val(it,ik+1) = ides(ip)%keys(it)%val(ik)                     !--- Other columns: keys values
1033      END DO
1034    END DO
1035    IF(test(fmsg('Problem with the table content', modname, dispTable(prf, ttl, val, cat([(it,it=1,nt)]), rFmt='(EN8.4)',     &
1036       sub=modname)), lerr)) RETURN
1037    DEALLOCATE(ttl, val)
1038  END DO       
1039END FUNCTION dispIsotopes
1040!==============================================================================================================================
1041
1042
1043!==============================================================================================================================
1044SUBROUTINE addKey_1(key, val, ky, lOverWrite)
1045!------------------------------------------------------------------------------------------------------------------------------
1046! Purpose: Add the <key>=<val> pair in the "ky" keys descriptor.
1047!------------------------------------------------------------------------------------------------------------------------------
1048  CHARACTER(LEN=*),  INTENT(IN)    :: key, val
1049  TYPE(keys_type),   INTENT(INOUT) :: ky
1050  LOGICAL, OPTIONAL, INTENT(IN)    :: lOverWrite
1051  CHARACTER(LEN=maxlen), ALLOCATABLE :: k(:), v(:)
1052  INTEGER :: iky, nky
1053  LOGICAL :: lo
1054!------------------------------------------------------------------------------------------------------------------------------
1055  lo=.FALSE.; IF(PRESENT(lOverWrite)) lo=lOverWrite
1056  iky = strIdx(ky%key,key)
1057  IF(iky == 0) THEN
1058    nky = SIZE(ky%key)
1059    IF(nky == 0) THEN; ky%key = [key]; ky%val = [val]; ELSE; ky%key = [ky%key, key]; ky%val = [ky%val, val]; END IF
1060  ELSE IF(lo) THEN                                                   !--- Overwriting
1061    ky%key(iky) = key; ky%val(iky) = val
1062  END IF
1063END SUBROUTINE addKey_1
1064!==============================================================================================================================
1065SUBROUTINE addKey_m(key, val, ky, lOverWrite)
1066!------------------------------------------------------------------------------------------------------------------------------
1067! Purpose: Add the <key>=<val> pair in all the components of the "ky" keys descriptor.
1068!------------------------------------------------------------------------------------------------------------------------------
1069  CHARACTER(LEN=*),  INTENT(IN)    :: key, val
1070  TYPE(keys_type),   INTENT(INOUT) :: ky(:)
1071  LOGICAL, OPTIONAL, INTENT(IN)    :: lOverWrite
1072  INTEGER :: itr
1073  LOGICAL :: lo
1074!------------------------------------------------------------------------------------------------------------------------------
1075  lo=.FALSE.; IF(PRESENT(lOverWrite)) lo=lOverWrite
1076  DO itr = 1, SIZE(ky); CALL addKey_1(key, val, ky(itr), lo); END DO
1077END SUBROUTINE addKey_m
1078!==============================================================================================================================
1079SUBROUTINE addKeysFromDef(t, tr0)
1080!------------------------------------------------------------------------------------------------------------------------------
1081! Purpose: The values of the keys of the tracer named "tr0" are overwritten by the values found in the *.def files, if any.
1082!------------------------------------------------------------------------------------------------------------------------------
1083  TYPE(trac_type), ALLOCATABLE, INTENT(INOUT) :: t(:)
1084  CHARACTER(LEN=*),             INTENT(IN)    :: tr0
1085  CHARACTER(LEN=maxlen) :: val
1086  INTEGER               :: ik, jd
1087  jd = strIdx(t%name, tr0)
1088  IF(jd == 0) RETURN
1089  DO ik = 1, SIZE(t(jd)%keys%key)
1090    CALL get_in(t(jd)%keys%key(ik), val, '*none*')
1091    IF(val /= '*none*') CALL addKey_1(t(jd)%keys%key(ik), val, t(jd)%keys, .TRUE.)
1092  END DO
1093END SUBROUTINE addKeysFromDef
1094!==============================================================================================================================
1095SUBROUTINE delKey_1(itr, keyn, ky)
1096!------------------------------------------------------------------------------------------------------------------------------
1097! Purpose: Internal routine.
1098!   Remove <key>=<val> pairs in the "itr"th component of the "ky" keys descriptor.
1099!------------------------------------------------------------------------------------------------------------------------------
1100  INTEGER,          INTENT(IN)    :: itr
1101  CHARACTER(LEN=*), INTENT(IN)    :: keyn(:)
1102  TYPE(trac_type),  INTENT(INOUT) :: ky(:)
1103  CHARACTER(LEN=maxlen), ALLOCATABLE :: k(:), v(:)
1104  LOGICAL,               ALLOCATABLE :: ll(:)
1105  INTEGER :: iky
1106!------------------------------------------------------------------------------------------------------------------------------
1107  IF(itr<=0 .OR. itr>SIZE(ky, DIM=1)) RETURN                          !--- Index is out of range
1108  ll = [( ALL(keyn/=ky(itr)%keys%key(iky)), iky=1, SIZE(ky(itr)%keys%key) )]
1109  k = PACK(ky(itr)%keys%key, MASK=ll); CALL MOVE_ALLOC(FROM=k, TO=ky(itr)%keys%key)
1110  v = PACK(ky(itr)%keys%val, MASK=ll); CALL MOVE_ALLOC(FROM=v, TO=ky(itr)%keys%val)
1111END SUBROUTINE delKey_1
1112!==============================================================================================================================
1113SUBROUTINE delKey(keyn, ky)
1114!------------------------------------------------------------------------------------------------------------------------------
1115! Purpose: Internal routine.
1116!   Remove <key>=<val> pairs in all the components of the "t" tracers descriptor.
1117!------------------------------------------------------------------------------------------------------------------------------
1118  CHARACTER(LEN=*), INTENT(IN)    :: keyn(:)
1119  TYPE(trac_type),  INTENT(INOUT) :: ky(:)
1120  INTEGER :: iky
1121!------------------------------------------------------------------------------------------------------------------------------
1122  DO iky = 1, SIZE(ky); CALL delKey_1(iky, keyn, ky); END DO
1123END SUBROUTINE delKey
1124!==============================================================================================================================
1125
1126
1127!==============================================================================================================================
1128!=== PUBLIC ROUTINES: GET A KEY FROM A <key>=<val> LIST ; VECTORS, TRACER AND DATABASE VERSIONS ===============================
1129!=== BEWARE !!! IF THE "ky" ARGUMENT IS NOT PRESENT, THEN THE VARIABLES "tracers" AND "isotopes" ARE USED. ====================
1130!===     THEY ARE LOCAL TO THIS MODULE, SO MUST MUST BE INITIALIZED FIRST USING the "getKey_init" ROUTINE  ====================
1131!==============================================================================================================================
1132SUBROUTINE getKey_init(tracers_, isotopes_)
1133  TYPE(trac_type), OPTIONAL, INTENT(IN) ::  tracers_(:)
1134  TYPE(isot_type), OPTIONAL, INTENT(IN) :: isotopes_(:)
1135  IF(PRESENT( tracers_))  tracers =  tracers_
1136  IF(PRESENT(isotopes_)) isotopes = isotopes_
1137END SUBROUTINE getKey_init
1138!==============================================================================================================================
1139CHARACTER(LEN=maxlen) FUNCTION fgetKeyByIndex_s1(itr, keyn, ky, def_val) RESULT(val)
1140!------------------------------------------------------------------------------------------------------------------------------
1141! Purpose: Internal function ; get a string formatted key value (returned argument) from its name and the tracer index.
1142!------------------------------------------------------------------------------------------------------------------------------
1143  INTEGER,                    INTENT(IN) :: itr
1144  CHARACTER(LEN=*),           INTENT(IN) :: keyn
1145  TYPE(keys_type),            INTENT(IN) :: ky(:)
1146  CHARACTER(LEN=*), OPTIONAL, INTENT(IN) :: def_val
1147!------------------------------------------------------------------------------------------------------------------------------
1148  INTEGER :: iky
1149  iky = 0;  IF(itr >  0 .AND. itr <= SIZE(ky)) iky = strIdx(ky(itr)%key(:), keyn)
1150  val = ''; IF(iky /= 0) val = ky(itr)%val(iky)                      !--- Key was found
1151  IF(PRESENT(def_val) .AND. iky == 0) val = def_val                  !--- Default value from arguments
1152END FUNCTION fgetKeyByIndex_s1
1153!==============================================================================================================================
1154CHARACTER(LEN=maxlen) FUNCTION fgetKeyByName_s1(tname, keyn, ky, def_val, lerr) RESULT(val)
1155!------------------------------------------------------------------------------------------------------------------------------
1156! Purpose: Internal function ; get a string formatted key value (returned argument) from its name and the tracer name.
1157!------------------------------------------------------------------------------------------------------------------------------
1158  CHARACTER(LEN=*),           INTENT(IN)  :: tname, keyn
1159  TYPE(keys_type),            INTENT(IN)  :: ky(:)
1160  CHARACTER(LEN=*), OPTIONAL, INTENT(IN)  :: def_val
1161  LOGICAL,          OPTIONAL, INTENT(OUT) :: lerr
1162!------------------------------------------------------------------------------------------------------------------------------
1163  INTEGER :: iky, itr
1164  val = ''; iky = 0
1165  itr = strIdx(ky(:)%name, tname)                                    !--- Get the index of the wanted tracer
1166  IF(PRESENT(lerr)) lerr = itr==0; IF(itr == 0) RETURN
1167  IF(itr >  0 .AND. itr <= SIZE(ky)) iky = strIdx(ky(itr)%key(:), keyn)
1168  IF(iky /= 0) val = ky(itr)%val(iky)                                !--- Key was found
1169  IF(PRESENT(def_val) .AND. iky == 0) val = def_val                  !--- Default value from arguments
1170END FUNCTION fgetKeyByName_s1
1171!==============================================================================================================================
1172LOGICAL FUNCTION getKeyByName_s1(keyn, val, tname, ky) RESULT(lerr)
1173  !--- Purpose: Get the value of the key named "keyn" for the tracer named "tnam".
1174  !     * "ky" unspecified: try in "tracers" for "tnam" with phase and tagging suffixes, then in "isotopes" without.
1175  !     * "ky"   specified: try in "ky"      for "tnam" with phase and tagging suffixes, then without.
1176  !    The returned error code is always .FALSE.: an empty string is returned when the key hasn't been found.
1177  CHARACTER(LEN=*),          INTENT(IN)  :: keyn
1178  CHARACTER(LEN=maxlen),     INTENT(OUT) :: val
1179  CHARACTER(LEN=*),          INTENT(IN)  :: tname
1180  TYPE(keys_type), OPTIONAL, INTENT(IN)  :: ky(:)
1181  CHARACTER(LEN=maxlen) :: tnam
1182  INTEGER, ALLOCATABLE  :: is(:)
1183  INTEGER :: i, itr
1184  tnam = delPhase(strHead(tname,'_',.FALSE.))                        !--- Remove tag and phase
1185  IF(PRESENT(ky)) THEN
1186    val = fgetKeyByName_s1(tname, keyn, ky, lerr=lerr)               !--- "ky" and "tname"
1187    IF(val /= '' .OR. lerr)      RETURN
1188    val = fgetKeyByName_s1(tnam,  keyn, ky, lerr=lerr)               !--- "ky" and "tnam"
1189  ELSE
1190    IF(.NOT.ALLOCATED(tracers))  RETURN
1191    val = fgetKeyByName_s1(tname, keyn, tracers(:)%keys, lerr=lerr)  !--- "tracers" and "tname"
1192    IF(val /= ''.AND..NOT.lerr)  RETURN
1193    IF(.NOT.ALLOCATED(isotopes)) RETURN
1194    IF(SIZE(isotopes) == 0)      RETURN
1195    !--- Search the "is" isotopes class index of the isotope named "tnam"
1196    is = find([(ANY(isotopes(i)%keys(:)%name == tnam), i=1, SIZE(isotopes))])
1197    IF(test(SIZE(is) == 0,lerr)) RETURN
1198    val = fgetKeyByName_s1(tname, keyn, isotopes(is(1))%keys(:),lerr=lerr)!--- "isotopes" and "tnam"
1199  END IF
1200END FUNCTION getKeyByName_s1
1201!==============================================================================================================================
1202LOGICAL FUNCTION getKeyByName_sm(keyn, val, tname, ky) RESULT(lerr)
1203  CHARACTER(LEN=*),                   INTENT(IN)  :: keyn
1204  CHARACTER(LEN=maxlen), ALLOCATABLE, INTENT(OUT) ::   val(:)
1205  CHARACTER(LEN=*), TARGET, OPTIONAL, INTENT(IN)  :: tname(:)
1206  TYPE(keys_type),  TARGET, OPTIONAL, INTENT(IN)  ::    ky(:)
1207  TYPE(keys_type),           POINTER :: k(:)
1208  CHARACTER(LEN=maxlen), ALLOCATABLE :: n(:)
1209  INTEGER :: iq, nq
1210  IF(test(.NOT.(PRESENT(tname).OR.PRESENT(ky)), lerr)) RETURN
1211  IF(PRESENT(ky   )) nq = SIZE(ky%name)
1212  IF(PRESENT(tname)) nq = SIZE(  tname)
1213  ALLOCATE(val(nq))
1214  IF(PRESENT(tname)) THEN
1215    IF(     PRESENT(ky)) lerr = ANY([(getKeyByName_s1(keyn, val(iq), tname(iq),   ky), iq=1, nq)])
1216    IF(.NOT.PRESENT(ky)) lerr = ANY([(getKeyByName_s1(keyn, val(iq), tname(iq)      ), iq=1, nq)])
1217  ELSE;                  lerr = ANY([(getKeyByName_s1(keyn, val(iq), ky(iq)%name, ky), iq=1, nq)])
1218  END IF
1219END FUNCTION getKeyByName_sm
1220!==============================================================================================================================
1221LOGICAL FUNCTION getKeyByName_i1(keyn, val, tname, ky) RESULT(lerr)
1222  CHARACTER(LEN=*),          INTENT(IN)  :: keyn
1223  INTEGER,                   INTENT(OUT) :: val
1224  CHARACTER(LEN=*),          INTENT(IN)  :: tname
1225  TYPE(keys_type), OPTIONAL, INTENT(IN)  :: ky(:)
1226  CHARACTER(LEN=maxlen) :: sval
1227  INTEGER :: ierr
1228  IF(     PRESENT(ky)) lerr = getKeyByName_s1(keyn, sval, tname, ky)
1229  IF(.NOT.PRESENT(ky)) lerr = getKeyByName_s1(keyn, sval, tname)
1230  IF(test(fmsg('key "'//TRIM(keyn)//'" or tracer "'//TRIM(tname)//'" is missing',        modname, lerr), lerr)) RETURN
1231  READ(sval, *, IOSTAT=ierr) val
1232  IF(test(fmsg('key "'//TRIM(keyn)//'" of tracer "'//TRIM(tname)//'" is not an integer', modname, lerr), lerr)) RETURN
1233END FUNCTION getKeyByName_i1
1234!==============================================================================================================================
1235LOGICAL FUNCTION getKeyByName_im(keyn, val, tname, ky) RESULT(lerr)
1236  CHARACTER(LEN=*),                   INTENT(IN)  :: keyn
1237  INTEGER,               ALLOCATABLE, INTENT(OUT) ::   val(:)
1238  CHARACTER(LEN=*), TARGET, OPTIONAL, INTENT(IN)  :: tname(:)
1239  TYPE(keys_type),  TARGET, OPTIONAL, INTENT(IN)  ::    ky(:)
1240  TYPE(keys_type),           POINTER :: k(:)
1241  CHARACTER(LEN=maxlen), ALLOCATABLE :: n(:)
1242  INTEGER :: iq, nq
1243  IF(test(.NOT.(PRESENT(tname).OR.PRESENT(ky)), lerr)) RETURN
1244  IF(PRESENT(ky   )) nq = SIZE(ky%name)
1245  IF(PRESENT(tname)) nq = SIZE(  tname)
1246  ALLOCATE(val(nq))
1247  IF(PRESENT(tname)) THEN
1248    IF(     PRESENT(ky)) lerr = ANY([(getKeyByName_i1(keyn, val(iq), tname(iq),   ky), iq=1, nq)])
1249    IF(.NOT.PRESENT(ky)) lerr = ANY([(getKeyByName_i1(keyn, val(iq), tname(iq)      ), iq=1, nq)])
1250  ELSE;                  lerr = ANY([(getKeyByName_i1(keyn, val(iq), ky(iq)%name, ky), iq=1, nq)])
1251  END IF
1252END FUNCTION getKeyByName_im
1253!==============================================================================================================================
1254LOGICAL FUNCTION getKeyByName_r1(keyn, val, tname, ky) RESULT(lerr)
1255  CHARACTER(LEN=*),          INTENT(IN)  :: keyn
1256  REAL,                      INTENT(OUT) :: val
1257  CHARACTER(LEN=*),          INTENT(IN)  :: tname
1258  TYPE(keys_type), OPTIONAL, INTENT(IN)  :: ky(:)
1259  CHARACTER(LEN=maxlen) :: sval
1260  INTEGER :: ierr
1261  IF(     PRESENT(ky)) lerr = getKeyByName_s1(keyn, sval, tname, ky)
1262  IF(.NOT.PRESENT(ky)) lerr = getKeyByName_s1(keyn, sval, tname)
1263  IF(test(fmsg('key "'//TRIM(keyn)//'" or tracer "'//TRIM(tname)//'" is missing',    modname, lerr), lerr)) RETURN
1264  READ(sval, *, IOSTAT=ierr) val
1265  IF(test(fmsg('key "'//TRIM(keyn)//'" of tracer "'//TRIM(tname)//'" is not a real', modname, lerr), lerr)) RETURN
1266END FUNCTION getKeyByName_r1
1267!==============================================================================================================================
1268LOGICAL FUNCTION getKeyByName_rm(keyn, val, tname, ky) RESULT(lerr)
1269  CHARACTER(LEN=*),                   INTENT(IN)  :: keyn
1270  REAL,                  ALLOCATABLE, INTENT(OUT) ::   val(:)
1271  CHARACTER(LEN=*), TARGET, OPTIONAL, INTENT(IN)  :: tname(:)
1272  TYPE(keys_type),  TARGET, OPTIONAL, INTENT(IN)  ::    ky(:)
1273  TYPE(keys_type),           POINTER :: k(:)
1274  CHARACTER(LEN=maxlen), ALLOCATABLE :: n(:)
1275  INTEGER :: iq, nq
1276  IF(test(.NOT.(PRESENT(tname).OR.PRESENT(ky)), lerr)) RETURN
1277  IF(PRESENT(ky   )) nq = SIZE(ky%name)
1278  IF(PRESENT(tname)) nq = SIZE(  tname)
1279  ALLOCATE(val(nq))
1280  IF(PRESENT(tname)) THEN
1281    IF(     PRESENT(ky)) lerr = ANY([(getKeyByName_r1(keyn, val(iq), tname(iq),   ky), iq=1, nq)])
1282    IF(.NOT.PRESENT(ky)) lerr = ANY([(getKeyByName_r1(keyn, val(iq), tname(iq)      ), iq=1, nq)])
1283  ELSE;                  lerr = ANY([(getKeyByName_r1(keyn, val(iq), ky(iq)%name, ky), iq=1, nq)])
1284  END IF
1285END FUNCTION getKeyByName_rm
1286!==============================================================================================================================
1287
1288
1289!==============================================================================================================================
1290!=== REMOVE, IF ANY, OR ADD THE PHASE SUFFIX OF THE TRACER NAMED "s" ==========================================================
1291!==============================================================================================================================
1292ELEMENTAL CHARACTER(LEN=maxlen) FUNCTION delPhase(s) RESULT(out)
1293  CHARACTER(LEN=*), INTENT(IN) :: s
1294  INTEGER :: l, i, ix
1295  CHARACTER(LEN=maxlen) :: sh, st
1296  out = s
1297  IF(s == '') RETURN                                                           !--- Empty string: nothing to do
1298
1299  !--- Special case: old phases for water, no phases separator
1300  i = INDEX(s,'_'); sh = s; IF(i/=0) sh=s(1:i-1); st='H2O'; IF(i/=0) st='H2O_'//s(i+1:LEN_TRIM(s))
1301  IF(ANY([('H2O'//old_phases(ix:ix), ix=1, nphases)] == sh)) THEN; out=st; RETURN; END IF
1302
1303  !--- Index of found phase in "known_phases"
1304  ix = MAXLOC( [(i, i=1,nphases)], MASK=[( INDEX(s, phases_sep//known_phases(i:i))/=0, i=1, nphases)], DIM=1 )
1305  IF(ix == 0) RETURN                                                           !--- No phase pattern found
1306  i = INDEX(s, phases_sep//known_phases(ix:ix))                                !--- Index of <sep><pha> pattern in "str"
1307  l = LEN_TRIM(s)
1308  IF(i == l-1) THEN                                                            !--- <var><sep><pha>       => return <var>
1309    out = s(1:l-2)
1310  ELSE IF(s(i+2:i+2) == '_') THEN                                              !--- <var><sep><pha>_<tag> => return <var>_<tag>
1311    out = s(1:i-1)//s(i+2:l)
1312  END IF
1313END FUNCTION delPhase
1314!------------------------------------------------------------------------------------------------------------------------------
1315CHARACTER(LEN=maxlen) FUNCTION addPhase_s1(s,pha,ph_sep) RESULT(out)
1316  CHARACTER(LEN=*),           INTENT(IN) :: s
1317  CHARACTER(LEN=1),           INTENT(IN) :: pha
1318  CHARACTER(LEN=*), OPTIONAL, INTENT(IN) :: ph_sep
1319  CHARACTER(LEN=1) :: psep
1320  INTEGER :: l, i
1321  out = s
1322  IF(s == '') RETURN                                                           !--- Empty string: nothing to do
1323  psep = phases_sep; IF(PRESENT(ph_sep)) psep = ph_sep
1324  i = INDEX(s, '_')                                                            !--- /=0 for <var>_<tag> tracers names
1325  l = LEN_TRIM(s)
1326  IF(i == 0) out =  TRIM(s)//TRIM(psep)//pha                                   !--- <var>       => return <var><sep><pha>
1327  IF(i /= 0) out = s(1:i-1)//TRIM(psep)//pha//'_'//s(i+1:l)                    !--- <var>_<tag> => return <var><sep><pha>_<tag>
1328END FUNCTION addPhase_s1
1329!------------------------------------------------------------------------------------------------------------------------------
1330FUNCTION addPhase_sm(s,pha,ph_sep) RESULT(out)
1331  CHARACTER(LEN=*),           INTENT(IN) :: s(:)
1332  CHARACTER(LEN=1),           INTENT(IN) :: pha
1333  CHARACTER(LEN=*), OPTIONAL, INTENT(IN) :: ph_sep
1334  CHARACTER(LEN=maxlen),     ALLOCATABLE :: out(:)
1335  CHARACTER(LEN=1) :: psep
1336  INTEGER :: k
1337  psep = phases_sep; IF(PRESENT(ph_sep)) psep = ph_sep
1338  out = [( addPhase_s1(s(k), pha, psep), k=1, SIZE(s) )]
1339END FUNCTION addPhase_sm
1340!------------------------------------------------------------------------------------------------------------------------------
1341CHARACTER(LEN=maxlen) FUNCTION addPhase_i1(s,ipha,ph_sep) RESULT(out)
1342  CHARACTER(LEN=*),           INTENT(IN) :: s
1343  INTEGER,                    INTENT(IN) :: ipha
1344  CHARACTER(LEN=*), OPTIONAL, INTENT(IN) :: ph_sep
1345  CHARACTER(LEN=1) :: psep
1346  out = s
1347  IF(s == '') RETURN                                                           !--- Empty string: nothing to do
1348  psep = phases_sep; IF(PRESENT(ph_sep)) psep = ph_sep
1349  IF(psep == '') out = addPhase_s1(s,   old_phases(ipha:ipha), psep)
1350  IF(psep /= '') out = addPhase_s1(s, known_phases(ipha:ipha), psep)
1351END FUNCTION addPhase_i1
1352!------------------------------------------------------------------------------------------------------------------------------
1353FUNCTION addPhase_im(s,ipha,ph_sep) RESULT(out)
1354  CHARACTER(LEN=*),           INTENT(IN) :: s(:)
1355  INTEGER,                    INTENT(IN) :: ipha
1356  CHARACTER(LEN=*), OPTIONAL, INTENT(IN) :: ph_sep
1357  CHARACTER(LEN=maxlen),     ALLOCATABLE :: out(:)
1358  CHARACTER(LEN=1) :: psep
1359  INTEGER :: k
1360  psep = phases_sep; IF(PRESENT(ph_sep)) psep = ph_sep
1361  out = [( addPhase_i1(s(k), ipha, psep), k=1, SIZE(s) )]
1362END FUNCTION addPhase_im
1363!------------------------------------------------------------------------------------------------------------------------------
1364
1365
1366INTEGER FUNCTION getiPhase(tname, lPhase) RESULT(iphase)
1367  CHARACTER(LEN=*),  INTENT(IN)  :: tname
1368  LOGICAL, OPTIONAL, INTENT(OUT) :: lPhase
1369  CHARACTER(LEN=maxlen) :: s1
1370  INTEGER :: ip
1371  IF(PRESENT(lPhase)) lPhase = .TRUE.
1372
1373  !--- Old tracer name descending on water: H2O[v][l][i][_<isotope>][_<tag>]
1374  iphase = strIdx([('H2O'//old_phases(ip:ip), ip=1, nphases)], tname(1:MIN(4,LEN_TRIM(tname))))
1375  IF(iphase /= 0) RETURN
1376
1377  !--- New tracer name: <name>[_<phase>][_<tag>]
1378  iphase = INDEX(known_phases, TRIM(strHead(strTail(tname, phases_sep, .TRUE.), '_', .TRUE.)))
1379  IF(iphase /= 0) RETURN
1380
1381  !---Default case: 1 (gaseous phase)
1382  iphase = 1
1383  IF(PRESENT(lPhase)) lPhase = .FALSE.
1384END FUNCTION getiPhase 
1385
1386!------------------------------------------------------------------------------------------------------------------------------
1387CHARACTER(LEN=1) FUNCTION new2oldPhase(np) RESULT(op)
1388  CHARACTER(LEN=1), INTENT(IN) :: np
1389  op = old_phases(INDEX(known_phases,np):INDEX(known_phases,np))
1390END FUNCTION new2oldPhase
1391!------------------------------------------------------------------------------------------------------------------------------
1392
1393!------------------------------------------------------------------------------------------------------------------------------
1394CHARACTER(LEN=maxlen) FUNCTION old2newName(oldName, iPhase) RESULT(newName)
1395  !--- Convert an old style name into a new one.
1396  !    Only usable with old style "traceur.def" files, in which only water isotopes are allowed.
1397  !    In these files, H2O descendants names are: H2O<phase>[_<isotope>][_<tag>], with:
1398  !    phase = v, l or i ; isotope = eau, HDO, O18, O17 or HTO.
1399  CHARACTER(LEN=*),  INTENT(IN)  :: oldName
1400  INTEGER, OPTIONAL, INTENT(OUT) :: iPhase
1401  CHARACTER(LEN=maxlen) :: oldIso(5) = ['eau',     'HDO',     'O18',     'O17',     'HTO'    ]
1402  CHARACTER(LEN=maxlen) :: newIso(5) = ['H2[16]O', 'H[2]HO ', 'H2[18]O', 'H2[17]O', 'H[3]HO ']
1403  CHARACTER(LEN=maxlen), ALLOCATABLE :: tmp(:)
1404  INTEGER :: ix, ip, it, nt
1405  LOGICAL :: lPhase, lerr
1406  ip = getiPhase(oldName, lPhase)                                              !--- Get the phase ; lPhase==T: phase is needed
1407  IF(PRESENT(iPhase)) iPhase = ip
1408  IF(.NOT.lPhase) THEN; newName = oldName ; RETURN; END IF                     !--- Not a water descendant
1409  newName = addPhase('H2O', ip)
1410  lerr = strParse(oldName, '_', tmp, n=nt)
1411  IF(nt == 1) RETURN                                                           !--- H2O with phase
1412  ix = strIdx(oldIso, tmp(2))
1413  newName = tmp(2); IF(ix /= 0) newName = newIso(ix)                           !--- Isotope name
1414  IF(lPhase)  newName = addPhase(newName, ip)                                  !--- Phase is needed
1415  IF(nt == 3) newName = TRIM(newName)//'_'//TRIM(tmp(3))                       !--- Tagging tracer
1416END FUNCTION old2newName
1417!------------------------------------------------------------------------------------------------------------------------------
1418
1419!------------------------------------------------------------------------------------------------------------------------------
1420CHARACTER(LEN=maxlen) FUNCTION new2oldName(newName, iPhase) RESULT(oldName)
1421  !--- Convert a new style name into an old one.
1422  !    Only convertable names are water descendants names H2O_<phase>, <isotope>_<phase>, <isotope>_<phase>_<tag>, with:
1423  !    phase = g, l or s ; isotope = H2[16]O, H[2]O, H2<[18]O, H2[17]O or H[3]O.
1424  CHARACTER(LEN=*),  INTENT(IN)    :: newName
1425  INTEGER, OPTIONAL, INTENT(OUT)   :: iPhase
1426  CHARACTER(LEN=maxlen) :: oldIso(5) = ['eau',     'HDO',     'O18',     'O17',     'HTO'    ]
1427  CHARACTER(LEN=maxlen) :: newIso(5) = ['H2[16]O', 'H[2]HO ', 'H2[18]O', 'H2[17]O', 'H[3]HO '], tag
1428  INTEGER :: ix, ip, it, nt
1429  LOGICAL :: lPhase, lH2O
1430  lH2O = newName(1:MIN(3,LEN_TRIM(newName)))=='H2O'
1431  ix = strIdx(newIso, strHead(strHead(newName,'_',.TRUE.),phases_sep,.TRUE.))  !--- Isotope index
1432  IF(ix == 0 .AND. .NOT.lH2O) THEN; oldName=newName; RETURN; END IF            !--- Not a water descendant
1433  ip = getiPhase(newName, lPhase)                                              !--- Get the phase ; lPhase==T: phase is needed
1434  oldName = 'H2O'; IF(lPhase) oldName = addPhase('H2O', ip, '')                !--- H2O with phase
1435  IF(ix == 0) RETURN
1436  oldName = TRIM(oldName)//'_'//oldIso(ix)                                     !--- Isotope
1437  tag = strTail(delPhase(newName), TRIM(newIso(ix)))
1438  IF(tag /= delPhase(newName) .AND. tag /= '') oldName = TRIM(oldName)//tag    !--- Tagging tracer
1439END FUNCTION new2oldName
1440!------------------------------------------------------------------------------------------------------------------------------
1441
1442
1443!==============================================================================================================================
1444!=== GET THE NAME(S) OF THE ANCESTOR(S) OF TRACER(S) "tname" AT GENERATION "igen"  IN THE TRACERS DESCRIPTORS LIST "tr" =======
1445!==============================================================================================================================
1446CHARACTER(LEN=maxlen) FUNCTION ancestor_1(t, tname, igen) RESULT(out)
1447  TYPE(trac_type),   INTENT(IN) :: t(:)
1448  CHARACTER(LEN=*),  INTENT(IN) :: tname
1449  INTEGER, OPTIONAL, INTENT(IN) :: igen
1450  INTEGER :: ig, ix
1451  ig = 0; IF(PRESENT(igen)) ig = igen
1452  ix = idxAncestor_1(t, tname, ig)
1453  out = ''; IF(ix /= 0) out = t(ix)%name
1454END FUNCTION ancestor_1
1455!------------------------------------------------------------------------------------------------------------------------------
1456FUNCTION ancestor_m(t, tname, igen) RESULT(out)
1457  CHARACTER(LEN=maxlen), ALLOCATABLE     ::   out(:)
1458  TYPE(trac_type),            INTENT(IN) ::     t(:)
1459  CHARACTER(LEN=*), OPTIONAL, INTENT(IN) :: tname(:)
1460  INTEGER,          OPTIONAL, INTENT(IN) :: igen
1461  INTEGER, ALLOCATABLE :: ix(:)
1462  INTEGER :: ig
1463  ig = 0; IF(PRESENT(igen)) ig = igen
1464  IF(     PRESENT(tname)) ix = idxAncestor_m(t, tname,     ig)
1465  IF(.NOT.PRESENT(tname)) ix = idxAncestor_m(t, t(:)%name, ig)
1466  ALLOCATE(out(SIZE(ix))); out(:) = ''
1467  WHERE(ix /= 0) out = t(ix)%name
1468END FUNCTION ancestor_m
1469!==============================================================================================================================
1470
1471
1472!==============================================================================================================================
1473!=== GET THE INDEX(ES) OF THE ANCESTOR(S) OF TRACER(S) "tname" AT GENERATION "igen"  IN THE TRACERS DESCRIPTORS LIST "tr" =====
1474!==============================================================================================================================
1475INTEGER FUNCTION idxAncestor_1(t, tname, igen) RESULT(out)
1476! Return the name of the generation "igen" (>=0) ancestor of "tname"
1477  TYPE(trac_type),   INTENT(IN) :: t(:)
1478  CHARACTER(LEN=*),  INTENT(IN) :: tname
1479  INTEGER, OPTIONAL, INTENT(IN) :: igen
1480  INTEGER :: ig
1481  ig = 0; IF(PRESENT(igen)) ig = igen
1482  out = strIdx(t(:)%name, tname)
1483  IF(out == 0)                 RETURN            !--- Tracer not found
1484  IF(t(out)%iGeneration <= ig) RETURN            !--- Tracer has a lower generation number than asked generation 'igen"
1485  DO WHILE(t(out)%iGeneration > ig); out = strIdx(t(:)%name, t(out)%parent); END DO
1486END FUNCTION idxAncestor_1
1487!------------------------------------------------------------------------------------------------------------------------------
1488FUNCTION idxAncestor_m(t, tname, igen) RESULT(out)
1489  INTEGER,          ALLOCATABLE          ::   out(:)
1490  TYPE(trac_type),            INTENT(IN) ::     t(:)
1491  CHARACTER(LEN=*), OPTIONAL, INTENT(IN) :: tname(:)
1492  INTEGER,          OPTIONAL, INTENT(IN) :: igen
1493  INTEGER :: ig, ix
1494  ig = 0; IF(PRESENT(igen)) ig = igen
1495  IF(     PRESENT(tname)) out = [(idxAncestor_1(t, tname(ix),  ig), ix=1, SIZE(tname))]
1496  IF(.NOT.PRESENT(tname)) out = [(idxAncestor_1(t, t(ix)%name, ig), ix=1, SIZE(t))]
1497END FUNCTION idxAncestor_m
1498!==============================================================================================================================
1499
1500
1501END MODULE readTracFiles_mod
Note: See TracBrowser for help on using the repository browser.