source: readTracFiles_mod.f90 @ 1

Last change on this file since 1 was 1, checked in by dcugnet, 3 years ago

Initial sources.

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