source: CPL/oasis3-mct/branches/OASIS3-MCT_2.0_branch/lib/mct/mct/m_Merge.F90 @ 4775

Last change on this file since 4775 was 4775, checked in by aclsce, 4 years ago
  • Imported oasis3-mct from Cerfacs svn server (not suppotred anymore).

The version has been extracted from https://oasis3mct.cerfacs.fr/svn/branches/OASIS3-MCT_2.0_branch/oasis3-mct@1818

File size: 102.2 KB
Line 
1!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
2!    Math and Computer Science Division, Argonne National Laboratory   !
3!-----------------------------------------------------------------------
4! CVS m_Merge.F90,v 1.9 2005-11-30 02:58:13 larson Exp
5! CVS MCT_2_8_0
6!BOP -------------------------------------------------------------------
7!
8! !MODULE: m_Merge - Merge flux and state data from multiple sources.
9!
10! !DESCRIPTION:  This module supports {\em merging} of state and flux
11! data from multiple components with overlapping spatial domains for use
12! by another component.  For example, let the vectors ${\bf a}$ and
13! ${\bf b}$ be data from Components $A$ and $B$ that have been
14! interpolated onto the physical grid of another component $C$.  We wish
15! to combine the data from $A$ and $B$ to get a vector ${\bf c}$, which
16! represents the merged data on the grid of component $C$.  This merge
17! process is an element-by-element masked weighted average:
18! $$ c_i = {{{{\prod_{j=1}^J} M_{i}^j} {{\prod_{k=1}^K} F_{i}^k} a_i +
19! {{\prod_{p=1}^P} N_{i}^p} {{\prod_{q=1}^Q} G_{i}^q} b_i} \over
20! {{{\prod_{j=1}^J} M_{i}^j} {{\prod_{k=1}^K} F_{i}^k} +
21! {{\prod_{p=1}^P} N_{i}^p} {{\prod_{q=1}^Q} G_{i}^q}}}, $$
22! Where ${M_{i}^j}$ and ${N_{i}^p}$ are {\em integer masks} (which have
23! value either $0$ or $1$), and ${F_{i}^k}$ and ${G_{i}^q}$ are {\em real
24! masks} (which are in the closed interval $[0,1]$).
25!
26! Currently, we assume that the integer and real masks are stored in
27! the same {\tt GeneralGrid} datatype.  We also assume--and this is of
28! critical importance to the user--that the attributes to be merged are
29! the same for all the inputs and output.  If the user violates this
30! assumption, incorrect merges will occur for any attributes that are
31! present in only some (that is not all) of the inputs.
32!
33! This module supports explicitly the merging data from two, three, and
34! four components.  There is also a routine named {\tt MergeInData} that
35! allows the user to construct other merging schemes.
36!
37! !INTERFACE:
38
39 module m_Merge
40
41!
42! !USES:
43!
44!     No other modules used in the declaration section of this module.
45
46      implicit none
47
48      private   ! except
49
50! !PUBLIC TYPES:
51
52!     None.
53
54! !PUBLIC MEMBER FUNCTIONS:
55
56      public :: MergeTwo      ! Merge Output from two components
57                              ! for use by a third.
58      public :: MergeThree    ! Merge Output from three components
59                              ! for use by a fourth.
60      public :: MergeFour     ! Merge Output from four components
61                              ! for use by a fifth.
62      public :: MergeInData   ! Merge in data from a single component.
63
64    interface MergeTwo ; module procedure &
65         MergeTwoGGSP_, &
66         MergeTwoGGDP_
67    end interface
68    interface MergeThree ; module procedure &
69         MergeThreeGGSP_, &
70         MergeThreeGGDP_
71    end interface
72    interface MergeFour ; module procedure &
73         MergeFourGGSP_, &
74         MergeFourGGDP_
75    end interface
76    interface MergeInData ; module procedure &
77         MergeInDataGGSP_, &
78         MergeInDataGGDP_
79    end interface
80
81! !PUBLIC DATA MEMBERS:
82
83!     None.
84
85! !REVISION HISTORY:
86!       19Jun02 - J.W. Larson <larson@mcs.anl.gov> - Initial version.
87!EOP ___________________________________________________________________
88
89  character(len=*),parameter :: myname='MCT::m_Merge'
90
91 contains
92
93!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
94!    Math and Computer Science Division, Argonne National Laboratory   !
95!BOP -------------------------------------------------------------------
96!
97! !IROUTINE: MergeTwoGGSP_ - Merge Data from Two Sources
98!
99! !DESCRIPTION:  This routine merges {\tt REAL} attribute data from
100! two input {\tt AttrVect} arguments {\tt inAv1} and {\tt inAv2} to
101! a third {\tt AttrVect} {\tt outAv}.  The attributes to be merged are
102! determined entirely by the real attributes of {\tt outAv}.  If
103! {\tt outAv} shares one or more attributes with either of the inputs
104! {\tt inAv1} or {\tt inAv2}, a merge is performed on the individual
105! {\em intersections} of attributes between the pairs $({\tt outAv},
106! {\tt inAv1})$ and $({\tt outAv},{\tt inAv1})$.  Currently, it is assumed
107! that these pairwise intersections are all equal.  This assumption is of
108! critical importance to the user.  If the user violates this
109! assumption, incorrect merges of attributes that are present in some
110! (but not all) of the inputs will result.
111!
112! The merge operatrion is a masked
113! weighted element-by-element sum, as outlined in the following example. 
114! Let the vectors ${\bf a}$ and ${\bf b}$ be data from Components $A$
115! and $B$ that have been interpolated onto the physical grid of another
116! component $C$.  We wish to combine the data from $A$ and $B$ to get
117! a vector ${\bf c}$, which represents the merged data on the grid of
118! component $C$.  The merge relation to obtain the $i$th element of
119! {\bf c} is
120! $$ c_i = {1 \over {W_i}} \bigg\{ {{\prod_{j=1}^J} \kappa_{i}^j}
121! {{\prod_{k=1}^K} \alpha_{i}^k} {a_i} + {{\prod_{l=1}^L} \lambda_{i}^l}
122! {{\prod_{m=1}^M} \beta_{i}^m} {b_i} \bigg\} , $$
123! where
124! $$ {W_i} = {{\prod_{j=1}^J} \kappa_{i}^j} {{\prod_{k=1}^K} \alpha_{i}^k} +
125! {{\prod_{l=1}^L} \lambda_{i}^l} {{\prod_{m=1}^M} \beta_{i}^m}. $$
126! The quantities ${\kappa_{i}^j}$ and ${\lambda_{i}^l}$ are {\em integer
127! masks} (which have value either $0$ or $1$), and ${\alpha_{i}^k}$ and
128! ${\beta_{i}^m}$ are {\em real masks} (which are in the closed interval
129! $[0,1]$).
130!
131! The integer and real masks are stored as attributes to the same input
132! {\tt GeneralGrid} argument {\tt GGrid}.  The mask attribute names are
133! stored as substrings to the colon-separated strings contained in the
134! input {\tt CHARACTER} arguments {\tt iMaskTags1}, {\tt iMaskTags2},
135! {\tt rMaskTags1}, and {\tt rMaskTags2}.  The {\tt LOGICAL} input
136! argument {\tt CheckMasks} governs how the masks are applied.  If
137! ${\tt CheckMasks} = {\tt .TRUE.}$, the entries are checked to ensure
138! they meet the definitions of real and integer masks.  If
139! ${\tt CheckMasks} = {\tt .TRUE.}$ then the masks are multiplied
140! together on an element-by-element basis with no validation of their
141! entries (this option results in slightly higher performance).
142!
143! This routine returns the sume of the masked weights as a diagnostic.
144! This quantity is returned in the output {\tt REAL} array {\tt WeightSum}.
145!
146! The correspondence between the quantities in the above merge relation
147! and the arguments to this routine are summarized in the table.
148! \begin{center}
149! \begin{tabular}{|l|l|l|}\hline
150! {\bf Quantity} & {\bf Stored in} & {\bf Referenced by}  \\
151!  & {\bf Argument}    & {\bf Argument} \\
152! \hline
153! \hline
154! $ {a_i} $ & {\tt inAv1} & \\
155! \hline
156! $ {b_i} $ & {\tt inAv2} & \\
157! \hline
158! $ {c_i} $ & {\tt outAv} & \\
159! \hline
160! $ {\kappa_i^j}, j=1,\ldots,J $ & {\tt GGrid} & {\tt iMaskTags1}\\
161! & & ($J$ items) \\
162! \hline
163! $ {\alpha_i^k}, k=1,\ldots,K $ & {\tt GGrid} & {\tt rMaskTags1}\\
164! & & ($K$ items) \\
165! \hline
166! $ {\lambda_i^l}, l=1,\ldots,L $ & {\tt GGrid} & {\tt iMaskTags2}\\
167! & & ($L$ items) \\
168! \hline
169! $ {\beta_i^m}, m=1,\ldots,M $ & {\tt GGrid} & {\tt rMaskTags2}\\
170! & & ($M$ items) \\
171! \hline
172! $ {W_i} $ & {\tt WeightSum} & \\
173! \hline
174! \end{tabular}
175! \end{center}
176!
177! !INTERFACE:
178
179 subroutine MergeTwoGGSP_(inAv1, iMaskTags1, rMaskTags1, &
180                        inAv2, iMaskTags2, rMaskTags2, &
181                        GGrid, CheckMasks, outAv, WeightSum)
182!
183! !USES:
184!
185      use m_stdio
186      use m_die
187
188      use m_realkinds, only : SP, FP
189
190      use m_List, only : List
191      use m_List, only : List_allocated => allocated
192
193      use m_AttrVect, only : AttrVect
194      use m_AttrVect, only : AttrVect_lsize => lsize
195      use m_AttrVect, only : AttrVect_nRAttr => nRAttr
196
197      use m_GeneralGrid, only : GeneralGrid
198      use m_GeneralGrid, only : GeneralGrid_lsize => lsize
199
200      implicit none
201
202! !INPUT PARAMETERS:
203!
204      type(AttrVect),               intent(IN)    :: inAv1
205      character(len=*),   optional, intent(IN)    :: iMaskTags1
206      character(len=*),   optional, intent(IN)    :: rMaskTags1
207      type(AttrVect),               intent(IN)    :: inAv2
208      character(len=*),   optional, intent(IN)    :: iMaskTags2
209      character(len=*),   optional, intent(IN)    :: rMaskTags2
210      type(GeneralGrid),            intent(IN)    :: GGrid
211      logical,                      intent(IN)    :: CheckMasks
212
213! !INPUT/OUTPUT PARAMETERS:
214!
215      type(AttrVect),               intent(INOUT) :: outAv
216      real(SP),       dimension(:), pointer       :: WeightSum
217
218! !REVISION HISTORY:
219!       19Jun02 - Jay Larson <larson@mcs.anl.gov> - Interface spec.
220!        3Jul02 - Jay Larson <larson@mcs.anl.gov> - Implementation.
221!       10Jul02 - J. Larson <larson@mcs.anl.gov> - Improved argument
222!                 checking.
223!EOP ___________________________________________________________________
224!
225  character(len=*),parameter :: myname_=myname//'::MergeTwoGGSP_'
226
227  integer :: i, j
228  real(FP) :: invWeightSum
229
230       ! Begin argument sanity checks...
231
232       ! Have the input arguments been allocated?
233
234  if(.not.(List_allocated(inAv1%rList))) then
235     write(stderr,'(2a)') myname_, &
236          'ERROR--INPUT argument inAv1 has no real attributes!'
237     call die(myname_)
238  endif
239
240  if(.not.(List_allocated(inAv2%rList))) then
241     write(stderr,'(2a)') myname_, &
242          'ERROR--INPUT argument inAv2 has no real attributes!'
243     call die(myname_)
244  endif
245
246  if(.not.(List_allocated(outaV%rList))) then
247     write(stderr,'(2a)') myname_, &
248          'ERROR--INPUT/OUTPUT argument outAv has no real attributes!'
249     call die(myname_)
250  endif
251
252  if(present(iMaskTags1) .or. present(iMaskTags2)) then
253     if(.not.(List_allocated(GGrid%data%iList))) then
254        write(stderr,'(3a)') myname_, &
255             'ERROR--Integer masking requested, but input argument GGrid ', &
256             'has no integer attributes!'
257        call die(myname_)
258     endif
259  endif
260
261  if(present(rMaskTags1) .or. present(rMaskTags2)) then
262     if(.not.(List_allocated(GGrid%data%rList))) then
263        write(stderr,'(3a)') myname_, &
264             'ERROR--Real masking requested, but input argument GGrid ', &
265             'has no real attributes!'
266        call die(myname_)
267     endif
268  endif
269
270  if(.not.(associated(WeightSum))) then
271     write(stderr,'(2a)') myname_, &
272          'ERROR--INPUT/OUPUT argument WeightSum has not been allocated!'
273     call die(myname_)
274  endif
275
276       ! Do the vector lengths match?
277
278  if(AttrVect_lsize(inAv1) /= AttrVect_lsize(outAv)) then
279     write(stderr,'(2a,2(a,i8))') myname_, &
280          ':: ERROR--Lengths of AttrVect arguments inAv1 and outAv must match.', &
281          'AttrVect_lsize(inAv1) = ',AttrVect_lsize(inAv1), &
282          'AttrVect_lsize(outAv) = ',AttrVect_lsize(outAv)
283     call die(myname_)
284  endif
285
286  if(AttrVect_lsize(inAv2) /= AttrVect_lsize(outAv)) then
287     write(stderr,'(2a,2(a,i8))') myname_, &
288          ':: ERROR--Lengths of AttrVect arguments inAv2 and outAv must match.', &
289          'AttrVect_lsize(inAv2) = ',AttrVect_lsize(inAv2), &
290          'AttrVect_lsize(outAv) = ',AttrVect_lsize(outAv)
291     call die(myname_)
292  endif
293
294  if(AttrVect_lsize(inAv1) /= GeneralGrid_lsize(GGrid)) then
295     write(stderr,'(2a,2(a,i8))') myname_, &
296          ':: ERROR--Lengths of arguments inAv1 and GGrid must match.', &
297          'AttrVect_lsize(inAv1) = ',AttrVect_lsize(inAv1), &
298          'AttrVect_lsize(outAv) = ',GeneralGrid_lsize(GGrid)
299     call die(myname_)
300  endif
301
302  if(AttrVect_lsize(inAv1) /= size(WeightSum)) then
303     write(stderr,'(2a,2(a,i8))') myname_, &
304          ':: ERROR--Lengths of arguments inAv1 and WeightSum must match.', &
305          'AttrVect_lsize(inAv1) = ',AttrVect_lsize(inAv1), &
306          'size(WeightSum) = ',size(WeightSum)
307     call die(myname_)
308  endif
309
310       ! ...end argument sanity checks.
311
312       ! Initialize the elements of WeightSum(:) to zero:
313
314  do i=1,size(WeightSum)
315     WeightSum(i) = 0._FP
316  end do
317
318       ! Process the incoming data one input AttrVect and mask tag
319       ! combination at a time.
320
321       ! First input AttrVect/mask combination...must work through
322       ! all the possible cases for optional arguments iMaskTags1 and
323       ! rMaskTags1.
324
325  if(present(iMaskTags1)) then
326
327     if(present(rMaskTags1)) then ! both real and integer masks
328        call MergeInDataGGSP_(inAv1, iMaskTags1, rMaskTags1, GGrid, &
329                            CheckMasks, outAv, WeightSum)
330     else ! only integer masks
331        call MergeInDataGGSP_(inAv1, iMaskTags=iMaskTags1, GGrid=GGrid, &
332                            CheckMasks=CheckMasks, outAv=outAv, &
333                            WeightSum=WeightSum)
334     endif
335
336  else
337
338     if(present(rMaskTags1)) then ! only real masks
339        call MergeInDataGGSP_(inAv1, rMaskTags=rMaskTags1, GGrid=GGrid, &
340                            CheckMasks=CheckMasks, outAv=outAv, &
341                            WeightSum=WeightSum)
342     else ! no masks at all
343        call MergeInDataGGSP_(inAv1, GGrid=GGrid, &
344                            CheckMasks=CheckMasks, outAv=outAv, &
345                            WeightSum=WeightSum)
346     endif
347
348  endif ! if(present(iMaskTags1))...
349
350       ! Second input AttrVect/mask combination...must work through
351       ! all the possible cases for optional arguments iMaskTags2 and
352       ! rMaskTags2.
353
354  if(present(iMaskTags2)) then
355
356     if(present(rMaskTags2)) then ! both real and integer masks
357        call MergeInDataGGSP_(inAv2, iMaskTags2, rMaskTags2, GGrid, &
358                            CheckMasks, outAv, WeightSum)
359     else ! only integer masks
360        call MergeInDataGGSP_(inAv2, iMaskTags=iMaskTags2, GGrid=GGrid, &
361                            CheckMasks=CheckMasks, outAv=outAv, &
362                            WeightSum=WeightSum)
363     endif
364
365  else
366
367     if(present(rMaskTags2)) then ! only real masks
368        call MergeInDataGGSP_(inAv2, rMaskTags=rMaskTags2, GGrid=GGrid, &
369                            CheckMasks=CheckMasks, outAv=outAv, &
370                            WeightSum=WeightSum)
371     else ! no masks at all
372        call MergeInDataGGSP_(inAv2, GGrid=GGrid, &
373                            CheckMasks=CheckMasks, outAv=outAv, &
374                            WeightSum=WeightSum)
375     endif
376
377  endif ! if(present(iMaskTags2))...
378
379       ! Now we must renormalize the entries in outAv by dividing
380       ! element-by-element by the sums of the merge weights, which
381       ! were accumulated in WeightSum(:)
382
383  do i=1,AttrVect_lsize(outAv)
384
385     if(WeightSum(i) /= 0._FP) then
386        invWeightSum = 1._FP / WeightSum(i)
387     else
388        write(stderr,'(2a,i8,a)') myname_,':: FATAL--WeightSum(', &
389                                  i,') is zero!'
390        call die(myname_)
391     endif
392
393     do j=1,AttrVect_nRAttr(outAv)
394        outAv%rAttr(j,i) = invWeightSum * outAv%rAttr(j,i) 
395     end do
396
397  end do
398
399       ! The merge is now complete.
400
401 end subroutine MergeTwoGGSP_
402
403!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
404!    Math and Computer Science Division, Argonne National Laboratory   !
405!-----------------------------------------------------------------------
406!
407! !IROUTINE: MergeTwoGGDP_ - merge data from two components.
408!
409! !DESCRIPTION:
410! Double precision version of MergeTwoGGSP_
411!
412! !INTERFACE:
413
414 subroutine MergeTwoGGDP_(inAv1, iMaskTags1, rMaskTags1, &
415                        inAv2, iMaskTags2, rMaskTags2, &
416                        GGrid, CheckMasks, outAv, WeightSum)
417!
418! !USES:
419!
420      use m_stdio
421      use m_die
422
423      use m_realkinds, only : DP, FP
424
425      use m_List, only : List
426      use m_List, only : List_allocated => allocated
427
428      use m_AttrVect, only : AttrVect
429      use m_AttrVect, only : AttrVect_lsize => lsize
430      use m_AttrVect, only : AttrVect_nRAttr => nRAttr
431
432      use m_GeneralGrid, only : GeneralGrid
433      use m_GeneralGrid, only : GeneralGrid_lsize => lsize
434
435      implicit none
436
437! !INPUT PARAMETERS:
438!
439      type(AttrVect),               intent(IN)    :: inAv1
440      character(len=*),   optional, intent(IN)    :: iMaskTags1
441      character(len=*),   optional, intent(IN)    :: rMaskTags1
442      type(AttrVect),               intent(IN)    :: inAv2
443      character(len=*),   optional, intent(IN)    :: iMaskTags2
444      character(len=*),   optional, intent(IN)    :: rMaskTags2
445      type(GeneralGrid),            intent(IN)    :: GGrid
446      logical,                      intent(IN)    :: CheckMasks
447
448! !INPUT/OUTPUT PARAMETERS:
449!
450      type(AttrVect),               intent(INOUT) :: outAv
451      real(DP),       dimension(:), pointer       :: WeightSum
452
453! !REVISION HISTORY:
454!       19Jun02 - Jay Larson <larson@mcs.anl.gov> - Interface spec.
455!        3Jul02 - Jay Larson <larson@mcs.anl.gov> - Implementation.
456!       10Jul02 - J. Larson <larson@mcs.anl.gov> - Improved argument
457!                 checking.
458!_______________________________________________________________________
459!
460  character(len=*),parameter :: myname_=myname//'::MergeTwoGGDP_'
461
462  integer :: i, j
463  real(FP) :: invWeightSum
464
465       ! Begin argument sanity checks...
466
467       ! Have the input arguments been allocated?
468
469  if(.not.(List_allocated(inAv1%rList))) then
470     write(stderr,'(2a)') myname_, &
471          'ERROR--INPUT argument inAv1 has no real attributes!'
472     call die(myname_)
473  endif
474
475  if(.not.(List_allocated(inAv2%rList))) then
476     write(stderr,'(2a)') myname_, &
477          'ERROR--INPUT argument inAv2 has no real attributes!'
478     call die(myname_)
479  endif
480
481  if(.not.(List_allocated(outaV%rList))) then
482     write(stderr,'(2a)') myname_, &
483          'ERROR--INPUT/OUTPUT argument outAv has no real attributes!'
484     call die(myname_)
485  endif
486
487  if(present(iMaskTags1) .or. present(iMaskTags2)) then
488     if(.not.(List_allocated(GGrid%data%iList))) then
489        write(stderr,'(3a)') myname_, &
490             'ERROR--Integer masking requested, but input argument GGrid ', &
491             'has no integer attributes!'
492        call die(myname_)
493     endif
494  endif
495
496  if(present(rMaskTags1) .or. present(rMaskTags2)) then
497     if(.not.(List_allocated(GGrid%data%rList))) then
498        write(stderr,'(3a)') myname_, &
499             'ERROR--Real masking requested, but input argument GGrid ', &
500             'has no real attributes!'
501        call die(myname_)
502     endif
503  endif
504
505  if(.not.(associated(WeightSum))) then
506     write(stderr,'(2a)') myname_, &
507          'ERROR--INPUT/OUPUT argument WeightSum has not been allocated!'
508     call die(myname_)
509  endif
510
511       ! Do the vector lengths match?
512
513  if(AttrVect_lsize(inAv1) /= AttrVect_lsize(outAv)) then
514     write(stderr,'(2a,2(a,i8))') myname_, &
515          ':: ERROR--Lengths of AttrVect arguments inAv1 and outAv must match.', &
516          'AttrVect_lsize(inAv1) = ',AttrVect_lsize(inAv1), &
517          'AttrVect_lsize(outAv) = ',AttrVect_lsize(outAv)
518     call die(myname_)
519  endif
520
521  if(AttrVect_lsize(inAv2) /= AttrVect_lsize(outAv)) then
522     write(stderr,'(2a,2(a,i8))') myname_, &
523          ':: ERROR--Lengths of AttrVect arguments inAv2 and outAv must match.', &
524          'AttrVect_lsize(inAv2) = ',AttrVect_lsize(inAv2), &
525          'AttrVect_lsize(outAv) = ',AttrVect_lsize(outAv)
526     call die(myname_)
527  endif
528
529  if(AttrVect_lsize(inAv1) /= GeneralGrid_lsize(GGrid)) then
530     write(stderr,'(2a,2(a,i8))') myname_, &
531          ':: ERROR--Lengths of arguments inAv1 and GGrid must match.', &
532          'AttrVect_lsize(inAv1) = ',AttrVect_lsize(inAv1), &
533          'AttrVect_lsize(outAv) = ',GeneralGrid_lsize(GGrid)
534     call die(myname_)
535  endif
536
537  if(AttrVect_lsize(inAv1) /= size(WeightSum)) then
538     write(stderr,'(2a,2(a,i8))') myname_, &
539          ':: ERROR--Lengths of arguments inAv1 and WeightSum must match.', &
540          'AttrVect_lsize(inAv1) = ',AttrVect_lsize(inAv1), &
541          'size(WeightSum) = ',size(WeightSum)
542     call die(myname_)
543  endif
544
545       ! ...end argument sanity checks.
546
547       ! Initialize the elements of WeightSum(:) to zero:
548
549  do i=1,size(WeightSum)
550     WeightSum(i) = 0._FP
551  end do
552
553       ! Process the incoming data one input AttrVect and mask tag
554       ! combination at a time.
555
556       ! First input AttrVect/mask combination...must work through
557       ! all the possible cases for optional arguments iMaskTags1 and
558       ! rMaskTags1.
559
560  if(present(iMaskTags1)) then
561
562     if(present(rMaskTags1)) then ! both real and integer masks
563        call MergeInDataGGDP_(inAv1, iMaskTags1, rMaskTags1, GGrid, &
564                            CheckMasks, outAv, WeightSum)
565     else ! only integer masks
566        call MergeInDataGGDP_(inAv1, iMaskTags=iMaskTags1, GGrid=GGrid, &
567                            CheckMasks=CheckMasks, outAv=outAv, &
568                            WeightSum=WeightSum)
569     endif
570
571  else
572
573     if(present(rMaskTags1)) then ! only real masks
574        call MergeInDataGGDP_(inAv1, rMaskTags=rMaskTags1, GGrid=GGrid, &
575                            CheckMasks=CheckMasks, outAv=outAv, &
576                            WeightSum=WeightSum)
577     else ! no masks at all
578        call MergeInDataGGDP_(inAv1, GGrid=GGrid, &
579                            CheckMasks=CheckMasks, outAv=outAv, &
580                            WeightSum=WeightSum)
581     endif
582
583  endif ! if(present(iMaskTags1))...
584
585       ! Second input AttrVect/mask combination...must work through
586       ! all the possible cases for optional arguments iMaskTags2 and
587       ! rMaskTags2.
588
589  if(present(iMaskTags2)) then
590
591     if(present(rMaskTags2)) then ! both real and integer masks
592        call MergeInDataGGDP_(inAv2, iMaskTags2, rMaskTags2, GGrid, &
593                            CheckMasks, outAv, WeightSum)
594     else ! only integer masks
595        call MergeInDataGGDP_(inAv2, iMaskTags=iMaskTags2, GGrid=GGrid, &
596                            CheckMasks=CheckMasks, outAv=outAv, &
597                            WeightSum=WeightSum)
598     endif
599
600  else
601
602     if(present(rMaskTags2)) then ! only real masks
603        call MergeInDataGGDP_(inAv2, rMaskTags=rMaskTags2, GGrid=GGrid, &
604                            CheckMasks=CheckMasks, outAv=outAv, &
605                            WeightSum=WeightSum)
606     else ! no masks at all
607        call MergeInDataGGDP_(inAv2, GGrid=GGrid, &
608                            CheckMasks=CheckMasks, outAv=outAv, &
609                            WeightSum=WeightSum)
610     endif
611
612  endif ! if(present(iMaskTags2))...
613
614       ! Now we must renormalize the entries in outAv by dividing
615       ! element-by-element by the sums of the merge weights, which
616       ! were accumulated in WeightSum(:)
617
618  do i=1,AttrVect_lsize(outAv)
619
620     if(WeightSum(i) /= 0._FP) then
621        invWeightSum = 1._FP / WeightSum(i)
622     else
623        write(stderr,'(2a,i8,a)') myname_,':: FATAL--WeightSum(', &
624                                  i,') is zero!'
625        call die(myname_)
626     endif
627
628     do j=1,AttrVect_nRAttr(outAv)
629        outAv%rAttr(j,i) = invWeightSum * outAv%rAttr(j,i) 
630     end do
631
632  end do
633
634       ! The merge is now complete.
635
636 end subroutine MergeTwoGGDP_
637
638!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
639!    Math and Computer Science Division, Argonne National Laboratory   !
640!BOP -------------------------------------------------------------------
641!
642! !IROUTINE: MergeThreeGGSP_ - Merge Data from Three Sources
643!
644! !DESCRIPTION:  This routine merges {\tt REAL} attribute data from
645! three input {\tt AttrVect} arguments {\tt inAv1} , {\tt inAv2}, and
646! {\tt inAv3} to a fourth {\tt AttrVect} {\tt outAv}.  The attributes to
647! be merged are determined entirely by the real attributes of {\tt outAv}.
648!   If {\tt outAv} shares one or more attributes with any of the inputs
649! {\tt inAv1}, {\tt inAv2}, or {\tt inAv3}, a merge is performed on the
650! individual {\em intersections} of attributes between the pairs
651! $({\tt outAv},{\tt inAv1})$,  $({\tt outAv},{\tt inAv2})$,
652! and $({\tt outAv},{\tt inAv3})$.  Currently, it is assumed that these
653! pairwise intersections are all equal.  This assumption is of
654! critical importance to the user.  If the user violates this
655! assumption, incorrect merges of any attributes present only in some
656! (but not all) inputs will result.
657!
658! The merge operatrion is a masked
659! weighted element-by-element sum, as outlined in the following example. 
660! Let the vectors ${\bf a}$,${\bf b}$, and ${\bf c}$ be data from
661! Components $A$, $B$, and $C$ that have been interpolated onto the
662! physical grid of another component $D$.  We wish to combine the data
663! from $A$, $B$ and $C$ to get a vector ${\bf d}$, which represents the
664! merged data on the grid of component $D$.  The merge relation to obtain
665! the $i$th element of ${\bf d}$ is
666! $$ d_i = {1 \over {W_i}} \bigg\{ {{\prod_{j=1}^J} \kappa_{i}^j}
667! {{\prod_{k=1}^K} \alpha_{i}^k} {a_i} + {{\prod_{l=1}^L} \lambda_{i}^l}
668! {{\prod_{m=1}^M} \beta_{i}^m} {b_i} + {{\prod_{p=1}^P} \mu_{i}^p}
669! {{\prod_{q=1}^Q} \gamma_{i}^q} {c_i} \bigg\} , $$
670! where
671! $$ {W_i} = {{\prod_{j=1}^J} \kappa_{i}^j} {{\prod_{k=1}^K} \alpha_{i}^k} +
672! {{\prod_{l=1}^L} \lambda_{i}^l} {{\prod_{m=1}^M} \beta_{i}^m} +
673! {{\prod_{p=1}^P} \mu_{i}^p} {{\prod_{q=1}^Q} \gamma_{i}^q}. $$
674! The quantities ${\kappa_{i}^j}$, ${\lambda_{i}^p}$, and ${\mu_{i}^p}$ are
675! {\em integer masks} (which have value either $0$ or $1$), and
676! ${\alpha_{i}^k}$, ${\beta_{i}^m}$, and ${\gamma_{i}^q}$ are {\em real
677! masks} (which are in the closed interval $[0,1]$).
678!
679! The integer and real masks are stored as attributes to the same input
680! {\tt GeneralGrid} argument {\tt GGrid}.  The mask attribute names are
681! stored as substrings to the colon-separated strings contained in the
682! input {\tt CHARACTER} arguments {\tt iMaskTags1}, {\tt iMaskTags2},
683! {\tt iMaskTags3}, {\tt rMaskTags1}, {\tt rMaskTags2}, and
684! {\tt rMaskTags3}.  The {\tt LOGICAL} input argument {\tt CheckMasks}
685! governs how the masks are applied.  If ${\tt CheckMasks} = {\tt .TRUE.}$,
686! the entries are checked to ensure they meet the definitions of real
687! and integer masks.  If ${\tt CheckMasks} = {\tt .FALSE.}$ then the masks
688! are multiplied together on an element-by-element basis with no validation
689! of their entries (this option results in slightly higher performance). 
690!
691! This routine returns the sum of the masked weights as a diagnostic.
692! This quantity is returned in the output {\tt REAL} array {\tt WeightSum}.
693!
694! The correspondence between the quantities in the above merge relation
695! and the arguments to this routine are summarized in the table.
696! \begin{center}
697! \begin{tabular}{|l|l|l|}\hline
698! {\bf Quantity} & {\bf Stored in} & {\bf Referenced by}  \\
699!  & {\bf Argument}    & {\bf Argument} \\
700! \hline
701! \hline
702! $ {a_i} $ & {\tt inAv1} & \\
703! \hline
704! $ {b_i} $ & {\tt inAv2} & \\
705! \hline
706! $ {c_i} $ & {\tt inAv3} & \\
707! \hline
708! $ {d_i} $ & {\tt outAv} & \\
709! \hline
710! $ {\kappa_i^j}, j=1,\ldots,J $ & {\tt GGrid} & {\tt iMaskTags1}\\
711! & & ($J$ items) \\
712! \hline
713! $ {\alpha_i^k}, k=1,\ldots,K $ & {\tt GGrid} & {\tt rMaskTags1}\\
714! & & ($K$ items) \\
715! \hline
716! $ {\lambda_i^l}, l=1,\ldots,L $ & {\tt GGrid} & {\tt iMaskTags2}\\
717! & & ($L$ items) \\
718! \hline
719! $ {\beta_i^m}, m=1,\ldots,M $ & {\tt GGrid} & {\tt rMaskTags2}\\
720! & & ($M$ items) \\
721! \hline
722! $ {\mu_i^p}, p=1,\ldots,P $ & {\tt GGrid} & {\tt iMaskTags3}\\
723! & & ($L$ items) \\
724! \hline
725! $ {\gamma_i^q}, q=1,\ldots,Q $ & {\tt GGrid} & {\tt rMaskTags3}\\
726! & & ($M$ items) \\
727! \hline
728! $ {W_i} $ & {\tt WeightSum} & \\
729! \hline
730! \end{tabular}
731! \end{center}
732!
733! !INTERFACE:
734
735 subroutine MergeThreeGGSP_(inAv1, iMaskTags1, rMaskTags1, &
736                            inAv2, iMaskTags2, rMaskTags2, &
737                            inAv3, iMaskTags3, rMaskTags3, &
738                            GGrid, CheckMasks, outAv, WeightSum)
739!
740! !USES:
741!
742      use m_stdio
743      use m_die
744
745      use m_realkinds, only : SP, FP
746
747      use m_List, only : List
748      use m_List, only : List_allocated => allocated
749
750      use m_AttrVect, only : AttrVect
751      use m_AttrVect, only : AttrVect_lsize => lsize
752      use m_AttrVect, only : AttrVect_nRAttr => nRAttr
753
754      use m_GeneralGrid, only : GeneralGrid
755      use m_GeneralGrid, only : GeneralGrid_lsize => lsize
756
757      implicit none
758
759! !INPUT PARAMETERS:
760!
761      type(AttrVect),               intent(IN)    :: inAv1
762      character(len=*),   optional, intent(IN)    :: iMaskTags1
763      character(len=*),   optional, intent(IN)    :: rMaskTags1
764      type(AttrVect),               intent(IN)    :: inAv2
765      character(len=*),   optional, intent(IN)    :: iMaskTags2
766      character(len=*),   optional, intent(IN)    :: rMaskTags2
767      type(AttrVect),               intent(IN)    :: inAv3
768      character(len=*),   optional, intent(IN)    :: iMaskTags3
769      character(len=*),   optional, intent(IN)    :: rMaskTags3
770      type(GeneralGrid),            intent(IN)    :: GGrid
771      logical,                      intent(IN)    :: CheckMasks
772
773! !INPUT/OUTPUT PARAMETERS:
774!
775      type(AttrVect),               intent(INOUT) :: outAv
776      real(SP),       dimension(:), pointer       :: WeightSum
777
778! !REVISION HISTORY:
779!       19Jun02 - Jay Larson <larson@mcs.anl.gov> - Interface spec.
780!        3Jul02 - Jay Larson <larson@mcs.anl.gov> - Implementation.
781!       10Jul02 - J. Larson <larson@mcs.anl.gov> - Improved argument
782!                 checking.
783!EOP ___________________________________________________________________
784!
785  character(len=*),parameter :: myname_=myname//'::MergeThreeGGSP_'
786
787  integer :: i, j
788  real(FP) :: invWeightSum
789
790       ! Begin argument sanity checks...
791
792       ! Have the input arguments been allocated?
793
794  if(.not.(List_allocated(inAv1%rList))) then
795     write(stderr,'(2a)') myname_, &
796          'ERROR--INPUT argument inAv1 has no real attributes!'
797     call die(myname_)
798  endif
799
800  if(.not.(List_allocated(inAv2%rList))) then
801     write(stderr,'(2a)') myname_, &
802          'ERROR--INPUT argument inAv2 has no real attributes!'
803     call die(myname_)
804  endif
805
806  if(.not.(List_allocated(inAv3%rList))) then
807     write(stderr,'(2a)') myname_, &
808          'ERROR--INPUT argument inAv3 has no real attributes!'
809     call die(myname_)
810  endif
811
812  if(.not.(List_allocated(outaV%rList))) then
813     write(stderr,'(2a)') myname_, &
814          'ERROR--INPUT/OUTPUT argument outAv has no real attributes!'
815     call die(myname_)
816  endif
817
818  if(present(iMaskTags1) .or. present(iMaskTags2) .or. present(iMaskTags3)) then
819     if(.not.(List_allocated(GGrid%data%iList))) then
820        write(stderr,'(3a)') myname_, &
821             'ERROR--Integer masking requested, but input argument GGrid ', &
822             'has no integer attributes!'
823        call die(myname_)
824     endif
825  endif
826
827  if(present(rMaskTags1) .or. present(rMaskTags2) .or. present(rMaskTags3)) then
828     if(.not.(List_allocated(GGrid%data%rList))) then
829        write(stderr,'(3a)') myname_, &
830             'ERROR--Real masking requested, but input argument GGrid ', &
831             'has no real attributes!'
832        call die(myname_)
833     endif
834  endif
835
836  if(.not.(associated(WeightSum))) then
837     write(stderr,'(2a)') myname_, &
838          'ERROR--INPUT/OUPUT argument WeightSum has not been allocated!'
839     call die(myname_)
840  endif
841
842       ! Do the vector lengths match?
843
844  if(AttrVect_lsize(inAv1) /= AttrVect_lsize(outAv)) then
845     write(stderr,'(2a,2(a,i8))') myname_, &
846          ':: ERROR--Lengths of AttrVect arguments inAv1 and outAv must match.', &
847          'AttrVect_lsize(inAv1) = ',AttrVect_lsize(inAv1), &
848          'AttrVect_lsize(outAv) = ',AttrVect_lsize(outAv)
849     call die(myname_)
850  endif
851
852  if(AttrVect_lsize(inAv2) /= AttrVect_lsize(outAv)) then
853     write(stderr,'(2a,2(a,i8))') myname_, &
854          ':: ERROR--Lengths of AttrVect arguments inAv2 and outAv must match.', &
855          'AttrVect_lsize(inAv2) = ',AttrVect_lsize(inAv2), &
856          'AttrVect_lsize(outAv) = ',AttrVect_lsize(outAv)
857     call die(myname_)
858  endif
859
860  if(AttrVect_lsize(inAv3) /= AttrVect_lsize(outAv)) then
861     write(stderr,'(2a,2(a,i8))') myname_, &
862          ':: ERROR--Lengths of AttrVect arguments inAv3 and outAv must match.', &
863          'AttrVect_lsize(inAv3) = ',AttrVect_lsize(inAv3), &
864          'AttrVect_lsize(outAv) = ',AttrVect_lsize(outAv)
865     call die(myname_)
866  endif
867
868  if(AttrVect_lsize(inAv1) /= GeneralGrid_lsize(GGrid)) then
869     write(stderr,'(2a,2(a,i8))') myname_, &
870          ':: ERROR--Lengths of arguments inAv1 and GGrid must match.', &
871          'AttrVect_lsize(inAv1) = ',AttrVect_lsize(inAv1), &
872          'AttrVect_lsize(outAv) = ',GeneralGrid_lsize(GGrid)
873     call die(myname_)
874  endif
875
876  if(AttrVect_lsize(inAv1) /= size(WeightSum)) then
877     write(stderr,'(2a,2(a,i8))') myname_, &
878          ':: ERROR--Lengths of arguments inAv1 and WeightSum must match.', &
879          'AttrVect_lsize(inAv1) = ',AttrVect_lsize(inAv1), &
880          'size(WeightSum) = ',size(WeightSum)
881     call die(myname_)
882  endif
883
884       ! ...end argument sanity checks.
885
886       ! Initialize the elements of WeightSum(:) to zero:
887
888  do i=1,size(WeightSum)
889     WeightSum(i) = 0._FP
890  end do
891
892       ! Process the incoming data one input AttrVect and mask tag
893       ! combination at a time.
894
895       ! First input AttrVect/mask combination...must work through
896       ! all the possible cases for optional arguments iMaskTags1 and
897       ! rMaskTags1.
898
899  if(present(iMaskTags1)) then
900
901     if(present(rMaskTags1)) then ! both real and integer masks
902        call MergeInDataGGSP_(inAv1, iMaskTags1, rMaskTags1, GGrid, &
903                            CheckMasks, outAv, WeightSum)
904     else ! only integer masks
905        call MergeInDataGGSP_(inAv1, iMaskTags=iMaskTags1, GGrid=GGrid, &
906                            CheckMasks=CheckMasks, outAv=outAv, &
907                            WeightSum=WeightSum)
908     endif
909
910  else
911
912     if(present(rMaskTags1)) then ! only real masks
913        call MergeInDataGGSP_(inAv1, rMaskTags=rMaskTags1, GGrid=GGrid, &
914                            CheckMasks=CheckMasks, outAv=outAv, &
915                            WeightSum=WeightSum)
916     else ! no masks at all
917        call MergeInDataGGSP_(inAv1, GGrid=GGrid, &
918                            CheckMasks=CheckMasks, outAv=outAv, &
919                            WeightSum=WeightSum)
920     endif
921
922  endif ! if(present(iMaskTags1))...
923
924       ! Second input AttrVect/mask combination...must work through
925       ! all the possible cases for optional arguments iMaskTags2 and
926       ! rMaskTags2.
927
928  if(present(iMaskTags2)) then
929
930     if(present(rMaskTags2)) then ! both real and integer masks
931        call MergeInDataGGSP_(inAv2, iMaskTags2, rMaskTags2, GGrid, &
932                            CheckMasks, outAv, WeightSum)
933     else ! only integer masks
934        call MergeInDataGGSP_(inAv2, iMaskTags=iMaskTags2, GGrid=GGrid, &
935                            CheckMasks=CheckMasks, outAv=outAv, &
936                            WeightSum=WeightSum)
937     endif
938
939  else
940
941     if(present(rMaskTags2)) then ! only real masks
942        call MergeInDataGGSP_(inAv2, rMaskTags=rMaskTags2, GGrid=GGrid, &
943                            CheckMasks=CheckMasks, outAv=outAv, &
944                            WeightSum=WeightSum)
945     else ! no masks at all
946        call MergeInDataGGSP_(inAv2, GGrid=GGrid, &
947                            CheckMasks=CheckMasks, outAv=outAv, &
948                            WeightSum=WeightSum)
949     endif
950
951  endif ! if(present(iMaskTags2))...
952
953       ! Third input AttrVect/mask combination...must work through
954       ! all the possible cases for optional arguments iMaskTags3 and
955       ! rMaskTags3.
956
957  if(present(iMaskTags3)) then
958
959     if(present(rMaskTags3)) then ! both real and integer masks
960        call MergeInDataGGSP_(inAv3, iMaskTags3, rMaskTags3, GGrid, &
961                            CheckMasks, outAv, WeightSum)
962     else ! only integer masks
963        call MergeInDataGGSP_(inAv3, iMaskTags=iMaskTags3, GGrid=GGrid, &
964                            CheckMasks=CheckMasks, outAv=outAv, &
965                            WeightSum=WeightSum)
966     endif
967
968  else
969
970     if(present(rMaskTags3)) then ! only real masks
971        call MergeInDataGGSP_(inAv3, rMaskTags=rMaskTags3, GGrid=GGrid, &
972                            CheckMasks=CheckMasks, outAv=outAv, &
973                            WeightSum=WeightSum)
974     else ! no masks at all
975        call MergeInDataGGSP_(inAv3, GGrid=GGrid, &
976                            CheckMasks=CheckMasks, outAv=outAv, &
977                            WeightSum=WeightSum)
978     endif
979
980  endif ! if(present(iMaskTags3))...
981
982       ! Now we must renormalize the entries in outAv by dividing
983       ! element-by-element by the sums of the merge weights, which
984       ! were accumulated in WeightSum(:)
985
986  do i=1,AttrVect_lsize(outAv)
987
988     if(WeightSum(i) /= 0._FP) then
989        invWeightSum = 1._FP / WeightSum(i)
990     else
991        write(stderr,'(2a,i8,a)') myname_,':: FATAL--WeightSum(', &
992                                  i,') is zero!'
993        call die(myname_)
994     endif
995
996     do j=1,AttrVect_nRAttr(outAv)
997        outAv%rAttr(j,i) = invWeightSum * outAv%rAttr(j,i) 
998     end do
999
1000  end do
1001
1002       ! The merge is now complete.
1003
1004 end subroutine MergeThreeGGSP_
1005
1006!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1007!    Math and Computer Science Division, Argonne National Laboratory   !
1008!-----------------------------------------------------------------------
1009!
1010! !IROUTINE: MergeThreeGGDP_ - merge data from three components.
1011!
1012! !DESCRIPTION:
1013! Double precision version of MergeThreeGGSP_
1014!
1015! !INTERFACE:
1016
1017 subroutine MergeThreeGGDP_(inAv1, iMaskTags1, rMaskTags1, &
1018                            inAv2, iMaskTags2, rMaskTags2, &
1019                            inAv3, iMaskTags3, rMaskTags3, &
1020                            GGrid, CheckMasks, outAv, WeightSum)
1021!
1022! !USES:
1023!
1024      use m_stdio
1025      use m_die
1026
1027      use m_realkinds, only : DP, FP
1028
1029      use m_List, only : List
1030      use m_List, only : List_allocated => allocated
1031
1032      use m_AttrVect, only : AttrVect
1033      use m_AttrVect, only : AttrVect_lsize => lsize
1034      use m_AttrVect, only : AttrVect_nRAttr => nRAttr
1035
1036      use m_GeneralGrid, only : GeneralGrid
1037      use m_GeneralGrid, only : GeneralGrid_lsize => lsize
1038
1039      implicit none
1040
1041! !INPUT PARAMETERS:
1042!
1043      type(AttrVect),               intent(IN)    :: inAv1
1044      character(len=*),   optional, intent(IN)    :: iMaskTags1
1045      character(len=*),   optional, intent(IN)    :: rMaskTags1
1046      type(AttrVect),               intent(IN)    :: inAv2
1047      character(len=*),   optional, intent(IN)    :: iMaskTags2
1048      character(len=*),   optional, intent(IN)    :: rMaskTags2
1049      type(AttrVect),               intent(IN)    :: inAv3
1050      character(len=*),   optional, intent(IN)    :: iMaskTags3
1051      character(len=*),   optional, intent(IN)    :: rMaskTags3
1052      type(GeneralGrid),            intent(IN)    :: GGrid
1053      logical,                      intent(IN)    :: CheckMasks
1054
1055! !INPUT/OUTPUT PARAMETERS:
1056!
1057      type(AttrVect),               intent(INOUT) :: outAv
1058      real(DP),       dimension(:), pointer       :: WeightSum
1059
1060! !REVISION HISTORY:
1061!       19Jun02 - Jay Larson <larson@mcs.anl.gov> - Interface spec.
1062!        3Jul02 - Jay Larson <larson@mcs.anl.gov> - Implementation.
1063!       10Jul02 - J. Larson <larson@mcs.anl.gov> - Improved argument
1064!                 checking.
1065!_______________________________________________________________________
1066!
1067  character(len=*),parameter :: myname_=myname//'::MergeThreeGGDP_'
1068
1069  integer :: i, j
1070  real(FP) :: invWeightSum
1071
1072       ! Begin argument sanity checks...
1073
1074       ! Have the input arguments been allocated?
1075
1076  if(.not.(List_allocated(inAv1%rList))) then
1077     write(stderr,'(2a)') myname_, &
1078          'ERROR--INPUT argument inAv1 has no real attributes!'
1079     call die(myname_)
1080  endif
1081
1082  if(.not.(List_allocated(inAv2%rList))) then
1083     write(stderr,'(2a)') myname_, &
1084          'ERROR--INPUT argument inAv2 has no real attributes!'
1085     call die(myname_)
1086  endif
1087
1088  if(.not.(List_allocated(inAv3%rList))) then
1089     write(stderr,'(2a)') myname_, &
1090          'ERROR--INPUT argument inAv3 has no real attributes!'
1091     call die(myname_)
1092  endif
1093
1094  if(.not.(List_allocated(outaV%rList))) then
1095     write(stderr,'(2a)') myname_, &
1096          'ERROR--INPUT/OUTPUT argument outAv has no real attributes!'
1097     call die(myname_)
1098  endif
1099
1100  if(present(iMaskTags1) .or. present(iMaskTags2) .or. present(iMaskTags3)) then
1101     if(.not.(List_allocated(GGrid%data%iList))) then
1102        write(stderr,'(3a)') myname_, &
1103             'ERROR--Integer masking requested, but input argument GGrid ', &
1104             'has no integer attributes!'
1105        call die(myname_)
1106     endif
1107  endif
1108
1109  if(present(rMaskTags1) .or. present(rMaskTags2) .or. present(rMaskTags3)) then
1110     if(.not.(List_allocated(GGrid%data%rList))) then
1111        write(stderr,'(3a)') myname_, &
1112             'ERROR--Real masking requested, but input argument GGrid ', &
1113             'has no real attributes!'
1114        call die(myname_)
1115     endif
1116  endif
1117
1118  if(.not.(associated(WeightSum))) then
1119     write(stderr,'(2a)') myname_, &
1120          'ERROR--INPUT/OUPUT argument WeightSum has not been allocated!'
1121     call die(myname_)
1122  endif
1123
1124       ! Do the vector lengths match?
1125
1126  if(AttrVect_lsize(inAv1) /= AttrVect_lsize(outAv)) then
1127     write(stderr,'(2a,2(a,i8))') myname_, &
1128          ':: ERROR--Lengths of AttrVect arguments inAv1 and outAv must match.', &
1129          'AttrVect_lsize(inAv1) = ',AttrVect_lsize(inAv1), &
1130          'AttrVect_lsize(outAv) = ',AttrVect_lsize(outAv)
1131     call die(myname_)
1132  endif
1133
1134  if(AttrVect_lsize(inAv2) /= AttrVect_lsize(outAv)) then
1135     write(stderr,'(2a,2(a,i8))') myname_, &
1136          ':: ERROR--Lengths of AttrVect arguments inAv2 and outAv must match.', &
1137          'AttrVect_lsize(inAv2) = ',AttrVect_lsize(inAv2), &
1138          'AttrVect_lsize(outAv) = ',AttrVect_lsize(outAv)
1139     call die(myname_)
1140  endif
1141
1142  if(AttrVect_lsize(inAv3) /= AttrVect_lsize(outAv)) then
1143     write(stderr,'(2a,2(a,i8))') myname_, &
1144          ':: ERROR--Lengths of AttrVect arguments inAv3 and outAv must match.', &
1145          'AttrVect_lsize(inAv3) = ',AttrVect_lsize(inAv3), &
1146          'AttrVect_lsize(outAv) = ',AttrVect_lsize(outAv)
1147     call die(myname_)
1148  endif
1149
1150  if(AttrVect_lsize(inAv1) /= GeneralGrid_lsize(GGrid)) then
1151     write(stderr,'(2a,2(a,i8))') myname_, &
1152          ':: ERROR--Lengths of arguments inAv1 and GGrid must match.', &
1153          'AttrVect_lsize(inAv1) = ',AttrVect_lsize(inAv1), &
1154          'AttrVect_lsize(outAv) = ',GeneralGrid_lsize(GGrid)
1155     call die(myname_)
1156  endif
1157
1158  if(AttrVect_lsize(inAv1) /= size(WeightSum)) then
1159     write(stderr,'(2a,2(a,i8))') myname_, &
1160          ':: ERROR--Lengths of arguments inAv1 and WeightSum must match.', &
1161          'AttrVect_lsize(inAv1) = ',AttrVect_lsize(inAv1), &
1162          'size(WeightSum) = ',size(WeightSum)
1163     call die(myname_)
1164  endif
1165
1166       ! ...end argument sanity checks.
1167
1168       ! Initialize the elements of WeightSum(:) to zero:
1169
1170  do i=1,size(WeightSum)
1171     WeightSum(i) = 0._FP
1172  end do
1173
1174       ! Process the incoming data one input AttrVect and mask tag
1175       ! combination at a time.
1176
1177       ! First input AttrVect/mask combination...must work through
1178       ! all the possible cases for optional arguments iMaskTags1 and
1179       ! rMaskTags1.
1180
1181  if(present(iMaskTags1)) then
1182
1183     if(present(rMaskTags1)) then ! both real and integer masks
1184        call MergeInDataGGDP_(inAv1, iMaskTags1, rMaskTags1, GGrid, &
1185                            CheckMasks, outAv, WeightSum)
1186     else ! only integer masks
1187        call MergeInDataGGDP_(inAv1, iMaskTags=iMaskTags1, GGrid=GGrid, &
1188                            CheckMasks=CheckMasks, outAv=outAv, &
1189                            WeightSum=WeightSum)
1190     endif
1191
1192  else
1193
1194     if(present(rMaskTags1)) then ! only real masks
1195        call MergeInDataGGDP_(inAv1, rMaskTags=rMaskTags1, GGrid=GGrid, &
1196                            CheckMasks=CheckMasks, outAv=outAv, &
1197                            WeightSum=WeightSum)
1198     else ! no masks at all
1199        call MergeInDataGGDP_(inAv1, GGrid=GGrid, &
1200                            CheckMasks=CheckMasks, outAv=outAv, &
1201                            WeightSum=WeightSum)
1202     endif
1203
1204  endif ! if(present(iMaskTags1))...
1205
1206       ! Second input AttrVect/mask combination...must work through
1207       ! all the possible cases for optional arguments iMaskTags2 and
1208       ! rMaskTags2.
1209
1210  if(present(iMaskTags2)) then
1211
1212     if(present(rMaskTags2)) then ! both real and integer masks
1213        call MergeInDataGGDP_(inAv2, iMaskTags2, rMaskTags2, GGrid, &
1214                            CheckMasks, outAv, WeightSum)
1215     else ! only integer masks
1216        call MergeInDataGGDP_(inAv2, iMaskTags=iMaskTags2, GGrid=GGrid, &
1217                            CheckMasks=CheckMasks, outAv=outAv, &
1218                            WeightSum=WeightSum)
1219     endif
1220
1221  else
1222
1223     if(present(rMaskTags2)) then ! only real masks
1224        call MergeInDataGGDP_(inAv2, rMaskTags=rMaskTags2, GGrid=GGrid, &
1225                            CheckMasks=CheckMasks, outAv=outAv, &
1226                            WeightSum=WeightSum)
1227     else ! no masks at all
1228        call MergeInDataGGDP_(inAv2, GGrid=GGrid, &
1229                            CheckMasks=CheckMasks, outAv=outAv, &
1230                            WeightSum=WeightSum)
1231     endif
1232
1233  endif ! if(present(iMaskTags2))...
1234
1235       ! Third input AttrVect/mask combination...must work through
1236       ! all the possible cases for optional arguments iMaskTags3 and
1237       ! rMaskTags3.
1238
1239  if(present(iMaskTags3)) then
1240
1241     if(present(rMaskTags3)) then ! both real and integer masks
1242        call MergeInDataGGDP_(inAv3, iMaskTags3, rMaskTags3, GGrid, &
1243                            CheckMasks, outAv, WeightSum)
1244     else ! only integer masks
1245        call MergeInDataGGDP_(inAv3, iMaskTags=iMaskTags3, GGrid=GGrid, &
1246                            CheckMasks=CheckMasks, outAv=outAv, &
1247                            WeightSum=WeightSum)
1248     endif
1249
1250  else
1251
1252     if(present(rMaskTags3)) then ! only real masks
1253        call MergeInDataGGDP_(inAv3, rMaskTags=rMaskTags3, GGrid=GGrid, &
1254                            CheckMasks=CheckMasks, outAv=outAv, &
1255                            WeightSum=WeightSum)
1256     else ! no masks at all
1257        call MergeInDataGGDP_(inAv3, GGrid=GGrid, &
1258                            CheckMasks=CheckMasks, outAv=outAv, &
1259                            WeightSum=WeightSum)
1260     endif
1261
1262  endif ! if(present(iMaskTags3))...
1263
1264       ! Now we must renormalize the entries in outAv by dividing
1265       ! element-by-element by the sums of the merge weights, which
1266       ! were accumulated in WeightSum(:)
1267
1268  do i=1,AttrVect_lsize(outAv)
1269
1270     if(WeightSum(i) /= 0._FP) then
1271        invWeightSum = 1._FP / WeightSum(i)
1272     else
1273        write(stderr,'(2a,i8,a)') myname_,':: FATAL--WeightSum(', &
1274                                  i,') is zero!'
1275        call die(myname_)
1276     endif
1277
1278     do j=1,AttrVect_nRAttr(outAv)
1279        outAv%rAttr(j,i) = invWeightSum * outAv%rAttr(j,i) 
1280     end do
1281
1282  end do
1283
1284       ! The merge is now complete.
1285
1286 end subroutine MergeThreeGGDP_
1287
1288!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1289!    Math and Computer Science Division, Argonne National Laboratory   !
1290!BOP -------------------------------------------------------------------
1291!
1292! !IROUTINE: MergeFourGGSP_ - Merge Data from Four Sources
1293!
1294! !DESCRIPTION:  This routine merges {\tt REAL} attribute data from
1295! four input {\tt AttrVect} arguments {\tt inAv1} , {\tt inAv2},
1296! {\tt inAv3}, and {\tt inAv4} to a fifth {\tt AttrVect} {\tt outAv}.  The
1297! attributes to be merged are determined entirely by the real attributes
1298! of {\tt outAv}.  If {\tt outAv} shares one or more attributes with any of
1299! the inputs {\tt inAv1}, {\tt inAv2}, {\tt inAv3}, or {\tt inAv4}, a merge
1300! is performed on the individual {\em intersections} of attributes between
1301! the pairs $({\tt outAv},{\tt inAv1})$,  $({\tt outAv},{\tt inAv2})$,
1302! $({\tt outAv},{\tt inAv3})$, and $({\tt outAv},{\tt inAv3})$.  Currently,
1303! it is assumed that these pairwise intersections are all equal.  This
1304! assumption is of critical importance to the user.  If the user violates
1305! this assumption, incorrect merges of any attributes present only in some
1306! (but not all) the inputs will result.
1307!
1308! The merge operatrion is a masked
1309! weighted element-by-element sum, as outlined in the following example. 
1310! Let the vectors ${\bf a}$,${\bf b}$, ${\bf c}$ and ${\bf d}$ be data from
1311! Components $A$, $B$, $C$, and $D$ that have been interpolated onto the
1312! physical grid of another component $E$.  We wish to combine the data
1313! from $A$, $B$, $C$, and $D$ to get a vector ${\bf e}$, which represents the
1314! merged data on the grid of component $E$.  The merge relation to obtain
1315! the $i$th element of {\bf e} is
1316! $$ e_i = {1 \over {W_i}} \bigg\{ {{\prod_{j=1}^J} \kappa_{i}^j}
1317! {{\prod_{k=1}^K} \alpha_{i}^k} {a_i} + {{\prod_{l=1}^L} \lambda_{i}^l}
1318! {{\prod_{m=1}^M} \beta_{i}^m} {b_i} + {{\prod_{p=1}^P} \mu_{i}^p}
1319! {{\prod_{q=1}^Q} \gamma_{i}^q} {c_i} + 
1320! {{\prod_{r=1}^R} \nu_{i}^r} {{\prod_{s=1}^S} \delta_{i}^s} {d_i} \bigg\} , $$
1321! where
1322! $$ {W_i} = {{\prod_{j=1}^J} \kappa_{i}^j} {{\prod_{k=1}^K} \alpha_{i}^k} +
1323! {{\prod_{l=1}^L} \lambda_{i}^l} {{\prod_{m=1}^M} \beta_{i}^m} +
1324! {{\prod_{p=1}^P} \mu_{i}^p} {{\prod_{q=1}^Q} \gamma_{i}^q} +
1325! {{\prod_{r=1}^R} \nu_{i}^r} {{\prod_{s=1}^S} \delta_{i}^s}. $$
1326! The quantities ${\kappa_{i}^j}$, ${\lambda_{i}^p}$, ${\mu_{i}^p}$, and
1327! ${\nu_{i}^r}$ are {\em integer masks} (which have value either $0$ or $1$),
1328! and ${\alpha_{i}^k}$, ${\beta_{i}^m}$, ${\gamma_{i}^q}$, and ${\delta_{i}^s}$
1329! are {\em real masks} (which are in the closed interval $[0,1]$).
1330!
1331! The integer and real masks are stored as attributes to the same input
1332! {\tt GeneralGrid} argument {\tt GGrid}.  The mask attribute names are
1333! stored as substrings to the colon-separated strings contained in the
1334! input {\tt CHARACTER} arguments {\tt iMaskTags1}, {\tt iMaskTags2},
1335! {\tt iMaskTags3}, {\tt iMaskTags4}, {\tt rMaskTags1}, and {\tt rMaskTags2},
1336! {\tt rMaskTags3}, and {\tt rMaskTags4}, .  The {\tt LOGICAL} input
1337! argument {\tt CheckMasks} governs how the masks are applied.    If
1338! ${\tt CheckMasks} = {\tt .TRUE.}$, the entries are checked to ensure
1339! they meet the definitions of real and integer masks.  If ${\tt CheckMasks}
1340! = {\tt .FALSE.}$ then the masks are multiplied together on an
1341! element-by-element basis with no validation of their entries (this option
1342! results in slightly higher performance).
1343!
1344! This routine returns the sume of the masked weights as a diagnostic.
1345! This quantity is returned in the output {\tt REAL} array {\tt WeightSum}.
1346!
1347! The correspondence between the quantities in the above merge relation
1348! and the arguments to this routine are summarized in the table.
1349! \begin{center}
1350! \begin{tabular}{|l|l|l|}\hline
1351! {\bf Quantity} & {\bf Stored in} & {\bf Referenced by}  \\
1352!  & {\bf Argument}    & {\bf Argument} \\
1353! \hline
1354! \hline
1355! $ {a_i} $ & {\tt inAv1} & \\
1356! \hline
1357! $ {b_i} $ & {\tt inAv2} & \\
1358! \hline
1359! $ {c_i} $ & {\tt inAv3} & \\
1360! \hline
1361! $ {d_i} $ & {\tt inAv4} & \\
1362! \hline
1363! $ {e_i} $ & {\tt outAv} & \\
1364! \hline
1365! $ {\kappa_i^j}, j=1,\ldots,J $ & {\tt GGrid} & {\tt iMaskTags1}\\
1366! & & ($J$ items) \\
1367! \hline
1368! $ {\alpha_i^k}, k=1,\ldots,K $ & {\tt GGrid} & {\tt rMaskTags1}\\
1369! & & ($K$ items) \\
1370! \hline
1371! $ {\lambda_i^l}, l=1,\ldots,L $ & {\tt GGrid} & {\tt iMaskTags2}\\
1372! & & ($L$ items) \\
1373! \hline
1374! $ {\beta_i^m}, m=1,\ldots,M $ & {\tt GGrid} & {\tt rMaskTags2}\\
1375! & & ($M$ items) \\
1376! \hline
1377! $ {\mu_i^p}, p=1,\ldots,P $ & {\tt GGrid} & {\tt iMaskTags3}\\
1378! & & ($L$ items) \\
1379! \hline
1380! $ {\gamma_i^q}, q=1,\ldots,Q $ & {\tt GGrid} & {\tt rMaskTags3}\\
1381! & & ($M$ items) \\
1382! \hline
1383! $ {\nu_i^r}, r=1,\ldots,R $ & {\tt GGrid} & {\tt iMaskTags4}\\
1384! & & ($L$ items) \\
1385! \hline
1386! $ {\delta_i^s}, s=1,\ldots,S $ & {\tt GGrid} & {\tt rMaskTags4}\\
1387! & & ($M$ items) \\
1388! \hline
1389! $ {W_i} $ & {\tt WeightSum} & \\
1390! \hline
1391! \end{tabular}
1392! \end{center}
1393!
1394! !INTERFACE:
1395
1396 subroutine MergeFourGGSP_(inAv1, iMaskTags1, rMaskTags1, &
1397                           inAv2, iMaskTags2, rMaskTags2, &
1398                           inAv3, iMaskTags3, rMaskTags3, &
1399                           inAv4, iMaskTags4, rMaskTags4, &
1400                           GGrid, CheckMasks, outAv, WeightSum)
1401!
1402! !USES:
1403!
1404      use m_stdio
1405      use m_die
1406
1407      use m_realkinds, only : SP, FP
1408
1409      use m_List, only : List
1410      use m_List, only : List_allocated => allocated
1411
1412      use m_AttrVect, only : AttrVect
1413      use m_AttrVect, only : AttrVect_lsize => lsize
1414      use m_AttrVect, only : AttrVect_nRAttr => nRAttr
1415
1416      use m_GeneralGrid, only : GeneralGrid
1417      use m_GeneralGrid, only : GeneralGrid_lsize => lsize
1418
1419      implicit none
1420
1421! !INPUT PARAMETERS:
1422!
1423      type(AttrVect),               intent(IN)    :: inAv1
1424      character(len=*),   optional, intent(IN)    :: iMaskTags1
1425      character(len=*),   optional, intent(IN)    :: rMaskTags1
1426      type(AttrVect),               intent(IN)    :: inAv2
1427      character(len=*),   optional, intent(IN)    :: iMaskTags2
1428      character(len=*),   optional, intent(IN)    :: rMaskTags2
1429      type(AttrVect),               intent(IN)    :: inAv3
1430      character(len=*),   optional, intent(IN)    :: iMaskTags3
1431      character(len=*),   optional, intent(IN)    :: rMaskTags3
1432      type(AttrVect),               intent(IN)    :: inAv4
1433      character(len=*),   optional, intent(IN)    :: iMaskTags4
1434      character(len=*),   optional, intent(IN)    :: rMaskTags4
1435      type(GeneralGrid),            intent(IN)    :: GGrid
1436      logical,                      intent(IN)    :: CheckMasks
1437
1438! !INPUT/OUTPUT PARAMETERS:
1439!
1440      type(AttrVect),               intent(INOUT) :: outAv
1441      real(SP),       dimension(:), pointer       :: WeightSum
1442
1443! !REVISION HISTORY:
1444!       19Jun02 - Jay Larson <larson@mcs.anl.gov> - Interface spec.
1445!        3Jul02 - Jay Larson <larson@mcs.anl.gov> - Implementation.
1446!       10Jul02 - J. Larson <larson@mcs.anl.gov> - Improved argument
1447!                 checking.
1448!EOP ___________________________________________________________________
1449!
1450  character(len=*),parameter :: myname_=myname//'::MergeFourGGSP_'
1451
1452  integer :: i, j
1453  real(FP) :: invWeightSum
1454
1455       ! Begin argument sanity checks...
1456
1457       ! Have the input arguments been allocated?
1458
1459  if(.not.(List_allocated(inAv1%rList))) then
1460     write(stderr,'(2a)') myname_, &
1461          'ERROR--INPUT argument inAv1 has no real attributes!'
1462     call die(myname_)
1463  endif
1464
1465  if(.not.(List_allocated(inAv2%rList))) then
1466     write(stderr,'(2a)') myname_, &
1467          'ERROR--INPUT argument inAv2 has no real attributes!'
1468     call die(myname_)
1469  endif
1470
1471  if(.not.(List_allocated(inAv3%rList))) then
1472     write(stderr,'(2a)') myname_, &
1473          'ERROR--INPUT argument inAv3 has no real attributes!'
1474     call die(myname_)
1475  endif
1476
1477  if(.not.(List_allocated(inAv4%rList))) then
1478     write(stderr,'(2a)') myname_, &
1479          'ERROR--INPUT argument inAv4 has no real attributes!'
1480     call die(myname_)
1481  endif
1482
1483  if(.not.(List_allocated(outaV%rList))) then
1484     write(stderr,'(2a)') myname_, &
1485          'ERROR--INPUT/OUTPUT argument outAv has no real attributes!'
1486     call die(myname_)
1487  endif
1488
1489  if(present(iMaskTags1) .or. present(iMaskTags2) .or. &
1490                        present(iMaskTags3) .or. present(iMaskTags4)) then
1491     if(.not.(List_allocated(GGrid%data%iList))) then
1492        write(stderr,'(3a)') myname_, &
1493             'ERROR--Integer masking requested, but input argument GGrid ', &
1494             'has no integer attributes!'
1495        call die(myname_)
1496     endif
1497  endif
1498
1499  if(present(rMaskTags1) .or. present(rMaskTags2) .or. &
1500                        present(rMaskTags3) .or. present(rMaskTags4)) then
1501     if(.not.(List_allocated(GGrid%data%rList))) then
1502        write(stderr,'(3a)') myname_, &
1503             'ERROR--Real masking requested, but input argument GGrid ', &
1504             'has no real attributes!'
1505        call die(myname_)
1506     endif
1507  endif
1508
1509  if(.not.(associated(WeightSum))) then
1510     write(stderr,'(2a)') myname_, &
1511          'ERROR--INPUT/OUPUT argument WeightSum has not been allocated!'
1512     call die(myname_)
1513  endif
1514
1515       ! Do the vector lengths match?
1516
1517  if(AttrVect_lsize(inAv1) /= AttrVect_lsize(outAv)) then
1518     write(stderr,'(2a,2(a,i8))') myname_, &
1519          ':: ERROR--Lengths of AttrVect arguments inAv1 and outAv must match.', &
1520          'AttrVect_lsize(inAv1) = ',AttrVect_lsize(inAv1), &
1521          'AttrVect_lsize(outAv) = ',AttrVect_lsize(outAv)
1522     call die(myname_)
1523  endif
1524
1525  if(AttrVect_lsize(inAv2) /= AttrVect_lsize(outAv)) then
1526     write(stderr,'(2a,2(a,i8))') myname_, &
1527          ':: ERROR--Lengths of AttrVect arguments inAv2 and outAv must match.', &
1528          'AttrVect_lsize(inAv2) = ',AttrVect_lsize(inAv2), &
1529          'AttrVect_lsize(outAv) = ',AttrVect_lsize(outAv)
1530     call die(myname_)
1531  endif
1532
1533  if(AttrVect_lsize(inAv3) /= AttrVect_lsize(outAv)) then
1534     write(stderr,'(2a,2(a,i8))') myname_, &
1535          ':: ERROR--Lengths of AttrVect arguments inAv3 and outAv must match.', &
1536          'AttrVect_lsize(inAv3) = ',AttrVect_lsize(inAv3), &
1537          'AttrVect_lsize(outAv) = ',AttrVect_lsize(outAv)
1538     call die(myname_)
1539  endif
1540
1541  if(AttrVect_lsize(inAv4) /= AttrVect_lsize(outAv)) then
1542     write(stderr,'(2a,2(a,i8))') myname_, &
1543          ':: ERROR--Lengths of AttrVect arguments inAv4 and outAv must match.', &
1544          'AttrVect_lsize(inAv4) = ',AttrVect_lsize(inAv4), &
1545          'AttrVect_lsize(outAv) = ',AttrVect_lsize(outAv)
1546     call die(myname_)
1547  endif
1548
1549  if(AttrVect_lsize(inAv1) /= GeneralGrid_lsize(GGrid)) then
1550     write(stderr,'(2a,2(a,i8))') myname_, &
1551          ':: ERROR--Lengths of arguments inAv1 and GGrid must match.', &
1552          'AttrVect_lsize(inAv1) = ',AttrVect_lsize(inAv1), &
1553          'AttrVect_lsize(outAv) = ',GeneralGrid_lsize(GGrid)
1554     call die(myname_)
1555  endif
1556
1557  if(AttrVect_lsize(inAv1) /= size(WeightSum)) then
1558     write(stderr,'(2a,2(a,i8))') myname_, &
1559          ':: ERROR--Lengths of arguments inAv1 and WeightSum must match.', &
1560          'AttrVect_lsize(inAv1) = ',AttrVect_lsize(inAv1), &
1561          'size(WeightSum) = ',size(WeightSum)
1562     call die(myname_)
1563  endif
1564
1565       ! ...end argument sanity checks.
1566
1567       ! Initialize the elements of WeightSum(:) to zero:
1568
1569  do i=1,size(WeightSum)
1570     WeightSum(i) = 0._FP
1571  end do
1572
1573       ! Process the incoming data one input AttrVect and mask tag
1574       ! combination at a time.
1575
1576       ! First input AttrVect/mask combination...must work through
1577       ! all the possible cases for optional arguments iMaskTags1 and
1578       ! rMaskTags1.
1579
1580  if(present(iMaskTags1)) then
1581
1582     if(present(rMaskTags1)) then ! both real and integer masks
1583        call MergeInDataGGSP_(inAv1, iMaskTags1, rMaskTags1, GGrid, &
1584                            CheckMasks, outAv, WeightSum)
1585     else ! only integer masks
1586        call MergeInDataGGSP_(inAv1, iMaskTags=iMaskTags1, GGrid=GGrid, &
1587                            CheckMasks=CheckMasks, outAv=outAv, &
1588                            WeightSum=WeightSum)
1589     endif
1590
1591  else
1592
1593     if(present(rMaskTags1)) then ! only real masks
1594        call MergeInDataGGSP_(inAv1, rMaskTags=rMaskTags1, GGrid=GGrid, &
1595                            CheckMasks=CheckMasks, outAv=outAv, &
1596                            WeightSum=WeightSum)
1597     else ! no masks at all
1598        call MergeInDataGGSP_(inAv1, GGrid=GGrid, &
1599                            CheckMasks=CheckMasks, outAv=outAv, &
1600                            WeightSum=WeightSum)
1601     endif
1602
1603  endif ! if(present(iMaskTags1))...
1604
1605       ! Second input AttrVect/mask combination...must work through
1606       ! all the possible cases for optional arguments iMaskTags2 and
1607       ! rMaskTags2.
1608
1609  if(present(iMaskTags2)) then
1610
1611     if(present(rMaskTags2)) then ! both real and integer masks
1612        call MergeInDataGGSP_(inAv2, iMaskTags2, rMaskTags2, GGrid, &
1613                            CheckMasks, outAv, WeightSum)
1614     else ! only integer masks
1615        call MergeInDataGGSP_(inAv2, iMaskTags=iMaskTags2, GGrid=GGrid, &
1616                            CheckMasks=CheckMasks, outAv=outAv, &
1617                            WeightSum=WeightSum)
1618     endif
1619
1620  else
1621
1622     if(present(rMaskTags2)) then ! only real masks
1623        call MergeInDataGGSP_(inAv2, rMaskTags=rMaskTags2, GGrid=GGrid, &
1624                            CheckMasks=CheckMasks, outAv=outAv, &
1625                            WeightSum=WeightSum)
1626     else ! no masks at all
1627        call MergeInDataGGSP_(inAv2, GGrid=GGrid, &
1628                            CheckMasks=CheckMasks, outAv=outAv, &
1629                            WeightSum=WeightSum)
1630     endif
1631
1632  endif ! if(present(iMaskTags2))...
1633
1634       ! Third input AttrVect/mask combination...must work through
1635       ! all the possible cases for optional arguments iMaskTags3 and
1636       ! rMaskTags3.
1637
1638  if(present(iMaskTags3)) then
1639
1640     if(present(rMaskTags3)) then ! both real and integer masks
1641        call MergeInDataGGSP_(inAv3, iMaskTags3, rMaskTags3, GGrid, &
1642                            CheckMasks, outAv, WeightSum)
1643     else ! only integer masks
1644        call MergeInDataGGSP_(inAv3, iMaskTags=iMaskTags3, GGrid=GGrid, &
1645                            CheckMasks=CheckMasks, outAv=outAv, &
1646                            WeightSum=WeightSum)
1647     endif
1648
1649  else
1650
1651     if(present(rMaskTags3)) then ! only real masks
1652        call MergeInDataGGSP_(inAv3, rMaskTags=rMaskTags3, GGrid=GGrid, &
1653                            CheckMasks=CheckMasks, outAv=outAv, &
1654                            WeightSum=WeightSum)
1655     else ! no masks at all
1656        call MergeInDataGGSP_(inAv3, GGrid=GGrid, &
1657                            CheckMasks=CheckMasks, outAv=outAv, &
1658                            WeightSum=WeightSum)
1659     endif
1660
1661  endif ! if(present(iMaskTags3))...
1662
1663       ! Fourth input AttrVect/mask combination...must work through
1664       ! all the possible cases for optional arguments iMaskTags4 and
1665       ! rMaskTags4.
1666
1667  if(present(iMaskTags4)) then
1668
1669     if(present(rMaskTags4)) then ! both real and integer masks
1670        call MergeInDataGGSP_(inAv4, iMaskTags4, rMaskTags4, GGrid, &
1671                            CheckMasks, outAv, WeightSum)
1672     else ! only integer masks
1673        call MergeInDataGGSP_(inAv4, iMaskTags=iMaskTags4, GGrid=GGrid, &
1674                            CheckMasks=CheckMasks, outAv=outAv, &
1675                            WeightSum=WeightSum)
1676     endif
1677
1678  else
1679
1680     if(present(rMaskTags4)) then ! only real masks
1681        call MergeInDataGGSP_(inAv4, rMaskTags=rMaskTags4, GGrid=GGrid, &
1682                            CheckMasks=CheckMasks, outAv=outAv, &
1683                            WeightSum=WeightSum)
1684     else ! no masks at all
1685        call MergeInDataGGSP_(inAv4, GGrid=GGrid, &
1686                            CheckMasks=CheckMasks, outAv=outAv, &
1687                            WeightSum=WeightSum)
1688     endif
1689
1690  endif ! if(present(iMaskTags4))...
1691
1692       ! Now we must renormalize the entries in outAv by dividing
1693       ! element-by-element by the sums of the merge weights, which
1694       ! were accumulated in WeightSum(:)
1695
1696  do i=1,AttrVect_lsize(outAv)
1697
1698     if(WeightSum(i) /= 0._FP) then
1699        invWeightSum = 1._FP / WeightSum(i)
1700     else
1701        write(stderr,'(2a,i8,a)') myname_,':: FATAL--WeightSum(', &
1702                                  i,') is zero!'
1703        call die(myname_)
1704     endif
1705
1706     do j=1,AttrVect_nRAttr(outAv)
1707        outAv%rAttr(j,i) = invWeightSum * outAv%rAttr(j,i) 
1708     end do
1709
1710  end do
1711
1712       ! The merge is now complete.
1713
1714 end subroutine MergeFourGGSP_
1715
1716!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1717!    Math and Computer Science Division, Argonne National Laboratory   !
1718!-----------------------------------------------------------------------
1719!
1720! !IROUTINE: MergeFourGGDP_ - merge data from four components.
1721!
1722! !DESCRIPTION:
1723! Double precision versions of MergeFourGGSP_
1724!
1725! !INTERFACE:
1726
1727 subroutine MergeFourGGDP_(inAv1, iMaskTags1, rMaskTags1, &
1728                           inAv2, iMaskTags2, rMaskTags2, &
1729                           inAv3, iMaskTags3, rMaskTags3, &
1730                           inAv4, iMaskTags4, rMaskTags4, &
1731                           GGrid, CheckMasks, outAv, WeightSum)
1732!
1733! !USES:
1734!
1735      use m_stdio
1736      use m_die
1737
1738      use m_realkinds, only : DP, FP
1739
1740      use m_List, only : List
1741      use m_List, only : List_allocated => allocated
1742
1743      use m_AttrVect, only : AttrVect
1744      use m_AttrVect, only : AttrVect_lsize => lsize
1745      use m_AttrVect, only : AttrVect_nRAttr => nRAttr
1746
1747      use m_GeneralGrid, only : GeneralGrid
1748      use m_GeneralGrid, only : GeneralGrid_lsize => lsize
1749
1750      implicit none
1751
1752! !INPUT PARAMETERS:
1753!
1754      type(AttrVect),               intent(IN)    :: inAv1
1755      character(len=*),   optional, intent(IN)    :: iMaskTags1
1756      character(len=*),   optional, intent(IN)    :: rMaskTags1
1757      type(AttrVect),               intent(IN)    :: inAv2
1758      character(len=*),   optional, intent(IN)    :: iMaskTags2
1759      character(len=*),   optional, intent(IN)    :: rMaskTags2
1760      type(AttrVect),               intent(IN)    :: inAv3
1761      character(len=*),   optional, intent(IN)    :: iMaskTags3
1762      character(len=*),   optional, intent(IN)    :: rMaskTags3
1763      type(AttrVect),               intent(IN)    :: inAv4
1764      character(len=*),   optional, intent(IN)    :: iMaskTags4
1765      character(len=*),   optional, intent(IN)    :: rMaskTags4
1766      type(GeneralGrid),            intent(IN)    :: GGrid
1767      logical,                      intent(IN)    :: CheckMasks
1768
1769! !INPUT/OUTPUT PARAMETERS:
1770!
1771      type(AttrVect),               intent(INOUT) :: outAv
1772      real(DP),       dimension(:), pointer       :: WeightSum
1773
1774! !REVISION HISTORY:
1775!       19Jun02 - Jay Larson <larson@mcs.anl.gov> - Interface spec.
1776!        3Jul02 - Jay Larson <larson@mcs.anl.gov> - Implementation.
1777!       10Jul02 - J. Larson <larson@mcs.anl.gov> - Improved argument
1778!                 checking.
1779!_______________________________________________________________________
1780!
1781  character(len=*),parameter :: myname_=myname//'::MergeFourGGDP_'
1782
1783  integer :: i, j
1784  real(FP) :: invWeightSum
1785
1786       ! Begin argument sanity checks...
1787
1788       ! Have the input arguments been allocated?
1789
1790  if(.not.(List_allocated(inAv1%rList))) then
1791     write(stderr,'(2a)') myname_, &
1792          'ERROR--INPUT argument inAv1 has no real attributes!'
1793     call die(myname_)
1794  endif
1795
1796  if(.not.(List_allocated(inAv2%rList))) then
1797     write(stderr,'(2a)') myname_, &
1798          'ERROR--INPUT argument inAv2 has no real attributes!'
1799     call die(myname_)
1800  endif
1801
1802  if(.not.(List_allocated(inAv3%rList))) then
1803     write(stderr,'(2a)') myname_, &
1804          'ERROR--INPUT argument inAv3 has no real attributes!'
1805     call die(myname_)
1806  endif
1807
1808  if(.not.(List_allocated(inAv4%rList))) then
1809     write(stderr,'(2a)') myname_, &
1810          'ERROR--INPUT argument inAv4 has no real attributes!'
1811     call die(myname_)
1812  endif
1813
1814  if(.not.(List_allocated(outaV%rList))) then
1815     write(stderr,'(2a)') myname_, &
1816          'ERROR--INPUT/OUTPUT argument outAv has no real attributes!'
1817     call die(myname_)
1818  endif
1819
1820  if(present(iMaskTags1) .or. present(iMaskTags2) .or. &
1821                        present(iMaskTags3) .or. present(iMaskTags4)) then
1822     if(.not.(List_allocated(GGrid%data%iList))) then
1823        write(stderr,'(3a)') myname_, &
1824             'ERROR--Integer masking requested, but input argument GGrid ', &
1825             'has no integer attributes!'
1826        call die(myname_)
1827     endif
1828  endif
1829
1830  if(present(rMaskTags1) .or. present(rMaskTags2) .or. &
1831                        present(rMaskTags3) .or. present(rMaskTags4)) then
1832     if(.not.(List_allocated(GGrid%data%rList))) then
1833        write(stderr,'(3a)') myname_, &
1834             'ERROR--Real masking requested, but input argument GGrid ', &
1835             'has no real attributes!'
1836        call die(myname_)
1837     endif
1838  endif
1839
1840  if(.not.(associated(WeightSum))) then
1841     write(stderr,'(2a)') myname_, &
1842          'ERROR--INPUT/OUPUT argument WeightSum has not been allocated!'
1843     call die(myname_)
1844  endif
1845
1846       ! Do the vector lengths match?
1847
1848  if(AttrVect_lsize(inAv1) /= AttrVect_lsize(outAv)) then
1849     write(stderr,'(2a,2(a,i8))') myname_, &
1850          ':: ERROR--Lengths of AttrVect arguments inAv1 and outAv must match.', &
1851          'AttrVect_lsize(inAv1) = ',AttrVect_lsize(inAv1), &
1852          'AttrVect_lsize(outAv) = ',AttrVect_lsize(outAv)
1853     call die(myname_)
1854  endif
1855
1856  if(AttrVect_lsize(inAv2) /= AttrVect_lsize(outAv)) then
1857     write(stderr,'(2a,2(a,i8))') myname_, &
1858          ':: ERROR--Lengths of AttrVect arguments inAv2 and outAv must match.', &
1859          'AttrVect_lsize(inAv2) = ',AttrVect_lsize(inAv2), &
1860          'AttrVect_lsize(outAv) = ',AttrVect_lsize(outAv)
1861     call die(myname_)
1862  endif
1863
1864  if(AttrVect_lsize(inAv3) /= AttrVect_lsize(outAv)) then
1865     write(stderr,'(2a,2(a,i8))') myname_, &
1866          ':: ERROR--Lengths of AttrVect arguments inAv3 and outAv must match.', &
1867          'AttrVect_lsize(inAv3) = ',AttrVect_lsize(inAv3), &
1868          'AttrVect_lsize(outAv) = ',AttrVect_lsize(outAv)
1869     call die(myname_)
1870  endif
1871
1872  if(AttrVect_lsize(inAv4) /= AttrVect_lsize(outAv)) then
1873     write(stderr,'(2a,2(a,i8))') myname_, &
1874          ':: ERROR--Lengths of AttrVect arguments inAv4 and outAv must match.', &
1875          'AttrVect_lsize(inAv4) = ',AttrVect_lsize(inAv4), &
1876          'AttrVect_lsize(outAv) = ',AttrVect_lsize(outAv)
1877     call die(myname_)
1878  endif
1879
1880  if(AttrVect_lsize(inAv1) /= GeneralGrid_lsize(GGrid)) then
1881     write(stderr,'(2a,2(a,i8))') myname_, &
1882          ':: ERROR--Lengths of arguments inAv1 and GGrid must match.', &
1883          'AttrVect_lsize(inAv1) = ',AttrVect_lsize(inAv1), &
1884          'AttrVect_lsize(outAv) = ',GeneralGrid_lsize(GGrid)
1885     call die(myname_)
1886  endif
1887
1888  if(AttrVect_lsize(inAv1) /= size(WeightSum)) then
1889     write(stderr,'(2a,2(a,i8))') myname_, &
1890          ':: ERROR--Lengths of arguments inAv1 and WeightSum must match.', &
1891          'AttrVect_lsize(inAv1) = ',AttrVect_lsize(inAv1), &
1892          'size(WeightSum) = ',size(WeightSum)
1893     call die(myname_)
1894  endif
1895
1896       ! ...end argument sanity checks.
1897
1898       ! Initialize the elements of WeightSum(:) to zero:
1899
1900  do i=1,size(WeightSum)
1901     WeightSum(i) = 0._FP
1902  end do
1903
1904       ! Process the incoming data one input AttrVect and mask tag
1905       ! combination at a time.
1906
1907       ! First input AttrVect/mask combination...must work through
1908       ! all the possible cases for optional arguments iMaskTags1 and
1909       ! rMaskTags1.
1910
1911  if(present(iMaskTags1)) then
1912
1913     if(present(rMaskTags1)) then ! both real and integer masks
1914        call MergeInDataGGDP_(inAv1, iMaskTags1, rMaskTags1, GGrid, &
1915                            CheckMasks, outAv, WeightSum)
1916     else ! only integer masks
1917        call MergeInDataGGDP_(inAv1, iMaskTags=iMaskTags1, GGrid=GGrid, &
1918                            CheckMasks=CheckMasks, outAv=outAv, &
1919                            WeightSum=WeightSum)
1920     endif
1921
1922  else
1923
1924     if(present(rMaskTags1)) then ! only real masks
1925        call MergeInDataGGDP_(inAv1, rMaskTags=rMaskTags1, GGrid=GGrid, &
1926                            CheckMasks=CheckMasks, outAv=outAv, &
1927                            WeightSum=WeightSum)
1928     else ! no masks at all
1929        call MergeInDataGGDP_(inAv1, GGrid=GGrid, &
1930                            CheckMasks=CheckMasks, outAv=outAv, &
1931                            WeightSum=WeightSum)
1932     endif
1933
1934  endif ! if(present(iMaskTags1))...
1935
1936       ! Second input AttrVect/mask combination...must work through
1937       ! all the possible cases for optional arguments iMaskTags2 and
1938       ! rMaskTags2.
1939
1940  if(present(iMaskTags2)) then
1941
1942     if(present(rMaskTags2)) then ! both real and integer masks
1943        call MergeInDataGGDP_(inAv2, iMaskTags2, rMaskTags2, GGrid, &
1944                            CheckMasks, outAv, WeightSum)
1945     else ! only integer masks
1946        call MergeInDataGGDP_(inAv2, iMaskTags=iMaskTags2, GGrid=GGrid, &
1947                            CheckMasks=CheckMasks, outAv=outAv, &
1948                            WeightSum=WeightSum)
1949     endif
1950
1951  else
1952
1953     if(present(rMaskTags2)) then ! only real masks
1954        call MergeInDataGGDP_(inAv2, rMaskTags=rMaskTags2, GGrid=GGrid, &
1955                            CheckMasks=CheckMasks, outAv=outAv, &
1956                            WeightSum=WeightSum)
1957     else ! no masks at all
1958        call MergeInDataGGDP_(inAv2, GGrid=GGrid, &
1959                            CheckMasks=CheckMasks, outAv=outAv, &
1960                            WeightSum=WeightSum)
1961     endif
1962
1963  endif ! if(present(iMaskTags2))...
1964
1965       ! Third input AttrVect/mask combination...must work through
1966       ! all the possible cases for optional arguments iMaskTags3 and
1967       ! rMaskTags3.
1968
1969  if(present(iMaskTags3)) then
1970
1971     if(present(rMaskTags3)) then ! both real and integer masks
1972        call MergeInDataGGDP_(inAv3, iMaskTags3, rMaskTags3, GGrid, &
1973                            CheckMasks, outAv, WeightSum)
1974     else ! only integer masks
1975        call MergeInDataGGDP_(inAv3, iMaskTags=iMaskTags3, GGrid=GGrid, &
1976                            CheckMasks=CheckMasks, outAv=outAv, &
1977                            WeightSum=WeightSum)
1978     endif
1979
1980  else
1981
1982     if(present(rMaskTags3)) then ! only real masks
1983        call MergeInDataGGDP_(inAv3, rMaskTags=rMaskTags3, GGrid=GGrid, &
1984                            CheckMasks=CheckMasks, outAv=outAv, &
1985                            WeightSum=WeightSum)
1986     else ! no masks at all
1987        call MergeInDataGGDP_(inAv3, GGrid=GGrid, &
1988                            CheckMasks=CheckMasks, outAv=outAv, &
1989                            WeightSum=WeightSum)
1990     endif
1991
1992  endif ! if(present(iMaskTags3))...
1993
1994       ! Fourth input AttrVect/mask combination...must work through
1995       ! all the possible cases for optional arguments iMaskTags4 and
1996       ! rMaskTags4.
1997
1998  if(present(iMaskTags4)) then
1999
2000     if(present(rMaskTags4)) then ! both real and integer masks
2001        call MergeInDataGGDP_(inAv4, iMaskTags4, rMaskTags4, GGrid, &
2002                            CheckMasks, outAv, WeightSum)
2003     else ! only integer masks
2004        call MergeInDataGGDP_(inAv4, iMaskTags=iMaskTags4, GGrid=GGrid, &
2005                            CheckMasks=CheckMasks, outAv=outAv, &
2006                            WeightSum=WeightSum)
2007     endif
2008
2009  else
2010
2011     if(present(rMaskTags4)) then ! only real masks
2012        call MergeInDataGGDP_(inAv4, rMaskTags=rMaskTags4, GGrid=GGrid, &
2013                            CheckMasks=CheckMasks, outAv=outAv, &
2014                            WeightSum=WeightSum)
2015     else ! no masks at all
2016        call MergeInDataGGDP_(inAv4, GGrid=GGrid, &
2017                            CheckMasks=CheckMasks, outAv=outAv, &
2018                            WeightSum=WeightSum)
2019     endif
2020
2021  endif ! if(present(iMaskTags4))...
2022
2023       ! Now we must renormalize the entries in outAv by dividing
2024       ! element-by-element by the sums of the merge weights, which
2025       ! were accumulated in WeightSum(:)
2026
2027  do i=1,AttrVect_lsize(outAv)
2028
2029     if(WeightSum(i) /= 0._FP) then
2030        invWeightSum = 1._FP / WeightSum(i)
2031     else
2032        write(stderr,'(2a,i8,a)') myname_,':: FATAL--WeightSum(', &
2033                                  i,') is zero!'
2034        call die(myname_)
2035     endif
2036
2037     do j=1,AttrVect_nRAttr(outAv)
2038        outAv%rAttr(j,i) = invWeightSum * outAv%rAttr(j,i) 
2039     end do
2040
2041  end do
2042
2043       ! The merge is now complete.
2044
2045 end subroutine MergeFourGGDP_
2046
2047!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
2048!    Math and Computer Science Division, Argonne National Laboratory   !
2049!BOP -------------------------------------------------------------------
2050!
2051! !IROUTINE: MergeInDataGGSP_ - Add Data into a Merge
2052!
2053! !DESCRIPTION:  This routine takes input field data from the input
2054! {\tt AttrVect} argument {\tt inAv}, and merges the real attributes it
2055! shares with the input/output {\tt AttrVect} argument {\tt outAv}. 
2056! The merge is a masked merge of the form
2057! $$ c_i = c_i + {{\prod_{j=1}^J} M_{i}^j} {{\prod_{k=1}^K} F_{i}^k}
2058! a_i , $$
2059! where ${c_i}$ represents one element of one of the real attributes of
2060! {\tt outAv}, and ${a_i}$ represents one element of one of the real
2061! attributes of {\tt inAv}.  The ${M_{i}^j}$ are {\em integer masks} which
2062! have value either $0$ or $1$, and are integer attributes of the input
2063! {\tt GeneralGrid} argument {\tt GGrid}.  The ${F_{i}^k}$ are {\em real
2064! masks} whose values are in the closed interval $[0,1]$, and are real
2065! attributes of the input {\tt GeneralGrid} argument {\tt GGrid}.  The
2066! input {\tt CHARACTER} argument {\tt iMaskTags} is a string of colon-
2067! delimited strings that name the integer attributes in {\tt GGrid}
2068! that are used as the masks ${M_{i}^j}$.  The input {\tt CHARACTER}
2069! argument {\tt rMaskTags} is a string of colon-delimited strings
2070! that name the real attributes in {\tt GGrid} that are used as the
2071! masks ${F_{i}^k}$.  The output {\tt REAL} array {\tt WeightSum} is
2072! used to store a running sum of the product of the masks.  The
2073! {\tt LOGICAL} input argument {\tt CheckMasks} governs how the masks
2074! are applied.  If ${\tt CheckMasks} = {\tt .TRUE.}$, the entries are
2075! checked to ensure they meet the definitions of real and integer masks. 
2076! If ${\tt CheckMasks} = {\tt .FALSE.}$ then the masks are multiplied
2077! together on an element-by-element basis with no validation of their
2078! entries (this option results in slightly higher performance).
2079!
2080! {\tt N.B.:}  The lengths of the {\tt AttrVect} arguments {\tt inAv}
2081! and {\tt outAv} must be equal, and this length must also equal the
2082! lengths of {\tt GGrid} and {\tt WeightSum}.
2083!
2084! {\tt N.B.:}  This algorithm assumes the {\tt AttrVect} argument
2085! {\tt outAv} has been created, and its real attributes have been
2086! initialized.
2087!
2088! {\tt N.B.:}  This algorithm assumes that the array {\tt WeightSum}
2089! has been created and initialized.
2090!
2091! !INTERFACE:
2092
2093 subroutine MergeInDataGGSP_(inAv, iMaskTags, rMaskTags, GGrid, &
2094                             CheckMasks, outAv, WeightSum)
2095!
2096! !USES:
2097!
2098      use m_stdio
2099      use m_die
2100
2101      use m_realkinds, only : SP, FP
2102
2103      use m_String, only : String
2104      use m_String, only : String_clean => clean
2105      use m_String, only : String_ToChar => toChar
2106
2107      use m_List, only : List
2108      use m_List, only : List_init => init
2109      use m_List, only : List_clean => clean
2110      use m_List, only : List_nitem => nitem
2111      use m_List, only : List_get => get
2112      use m_List, only : List_identical => identical
2113      use m_List, only : List_allocated => allocated
2114
2115      use m_AttrVect, only : AttrVect
2116      use m_AttrVect, only : AttrVect_lsize => lsize
2117      use m_AttrVect, only : AttrVect_nRAttr => nRAttr
2118      use m_AttrVect, only : SharedAttrIndexList
2119
2120      use m_GeneralGrid, only : GeneralGrid
2121      use m_GeneralGrid, only : GeneralGrid_lsize => lsize
2122      use m_GeneralGrid, only : GeneralGrid_exportIAttr => exportIAttr
2123      use m_GeneralGrid, only : GeneralGrid_exportRAttr => exportRAttr
2124
2125      implicit none
2126
2127! !INPUT PARAMETERS:
2128!
2129      type(AttrVect),               intent(IN)    :: inAv
2130      character(len=*),   optional, intent(IN)    :: iMaskTags
2131      character(len=*),   optional, intent(IN)    :: rMaskTags
2132      type(GeneralGrid),            intent(IN)    :: GGrid
2133      logical,                      intent(IN)    :: CheckMasks
2134
2135! !INPUT/OUTPUT PARAMETERS:
2136!
2137      type(AttrVect),               intent(INOUT) :: outAv
2138      real(SP),       dimension(:), pointer       :: WeightSum
2139
2140! !REVISION HISTORY:
2141!       19Jun02 - Jay Larson <larson@mcs.anl.gov> - initial verson.
2142!       10Jul02 - J. Larson <larson@mcs.anl.gov> - Improved argument
2143!                 checking.
2144!EOP ___________________________________________________________________
2145!
2146  character(len=*),parameter :: myname_=myname//'::MergeInDataGGSP_'
2147
2148  integer :: i, ierr, j, length
2149  type(String) :: DummStr
2150  type(List) :: iMaskList, rMaskList
2151  integer,  dimension(:), pointer :: iMask,iDummy ! INTEGER mask workspace
2152  real(FP), dimension(:), pointer :: rMask,rDummy ! REAL mask workspace
2153
2154  logical :: RAttrIdentical ! flag to identify identical REAL attribute
2155                            ! lists in inAv and outAv
2156  integer :: NumSharedRAttr ! number of REAL attributes shared by inAv,outAv
2157       ! Cross-index storage for shared REAL attributes of inAv,outAv
2158  integer, dimension(:), pointer :: inAvIndices, outAvIndices
2159
2160       ! Begin argument sanity checks...
2161
2162       ! Have the input arguments been allocated?
2163
2164  if(.not.(List_allocated(inAv%rList))) then
2165     write(stderr,'(2a)') myname_, &
2166          'ERROR--INPUT argument inAv has no real attributes.'
2167     call die(myname_)
2168  endif
2169
2170  if(.not.(List_allocated(outaV%rList))) then
2171     write(stderr,'(2a)') myname_, &
2172          'ERROR--INPUT/OUTPUT argument outAv has no real attributes.'
2173     call die(myname_)
2174  endif
2175
2176  if(present(iMaskTags)) then
2177     if(.not.(List_allocated(GGrid%data%iList))) then
2178        write(stderr,'(3a)') myname_, &
2179             'ERROR--Integer masking requested, but input argument GGrid ', &
2180             'has no integer attributes.'
2181        call die(myname_)
2182     endif
2183  endif
2184
2185  if(present(rMaskTags)) then
2186     if(.not.(List_allocated(GGrid%data%rList))) then
2187        write(stderr,'(3a)') myname_, &
2188             'ERROR--Real masking requested, but input argument GGrid ', &
2189             'has no real attributes.'
2190        call die(myname_)
2191     endif
2192  endif
2193
2194  if(.not.(associated(WeightSum))) then
2195     write(stderr,'(2a)') myname_, &
2196          'ERROR--INPUT/OUPUT argument WeightSum has not been allocated.'
2197     call die(myname_)
2198  endif
2199
2200       ! Do the vector lengths match?
2201
2202  if(AttrVect_lsize(inAv) /= AttrVect_lsize(outAv)) then
2203     write(stderr,'(2a,2(a,i8))') myname_, &
2204          ':: ERROR--Lengths of AttrVect arguments inAv and outAv must match.', &
2205          'AttrVect_lsize(inAv) = ',AttrVect_lsize(inAv), &
2206          'AttrVect_lsize(outAv) = ',AttrVect_lsize(outAv)
2207     call die(myname_)
2208  endif
2209
2210  if(AttrVect_lsize(inAv) /= GeneralGrid_lsize(GGrid)) then
2211     write(stderr,'(2a,2(a,i8))') myname_, &
2212          ':: ERROR--Lengths of arguments inAv and GGrid must match.', &
2213          'AttrVect_lsize(inAv) = ',AttrVect_lsize(inAv), &
2214          'AttrVect_lsize(outAv) = ',GeneralGrid_lsize(GGrid)
2215     call die(myname_)
2216  endif
2217
2218  if(AttrVect_lsize(inAv) /= size(WeightSum)) then
2219     write(stderr,'(2a,2(a,i8))') myname_, &
2220          ':: ERROR--Lengths of arguments inAv and WeightSum must match.', &
2221          'AttrVect_lsize(inAv) = ',AttrVect_lsize(inAv), &
2222          'size(WeightSum) = ',size(WeightSum)
2223     call die(myname_)
2224  endif
2225
2226       ! ...end argument sanity checks.
2227
2228       ! Check for INTEGER masks.  If they are present, retrieve
2229       ! them and combine them into a single integer mask iMask(:)
2230
2231  if(present(iMaskTags)) then
2232
2233       ! allocate two arrays:  iMask (the final product),
2234       ! and iDummy (storage space for each mask as it is retrieved)
2235
2236     allocate(iMask(AttrVect_lsize(inAv)), iDummy(AttrVect_lsize(inAv)), &
2237              stat=ierr)
2238     if(ierr /= 0) then
2239        write(stderr,'(2a,i8)') myname_, &
2240             ':: allocate(iMask(...)...) failed with ierr = ',ierr
2241        call die(myname_)
2242     endif
2243
2244       ! Initialize all the elements of iMask to unity:
2245     iMask = 1
2246
2247       ! turn the colon-delimited string of tags into a List:
2248     call List_init(iMaskList,iMaskTags)
2249     
2250       ! Loop over the items in iMaskList, retrieving each mask
2251       ! into the array iDummy, checking it (if CheckMasks=.TRUE.),
2252       ! and multiplying it element-by-element into the array iMask.
2253
2254     do i=1,List_nitem(iMaskList)
2255       ! grab item as a String
2256        call List_get(DummStr, i, iMaskList) 
2257       ! use this String to identify an INTEGER GeneralGrid attribute
2258       ! for export to iDummy(:)
2259        call GeneralGrid_exportIAttr(GGrid, String_ToChar(DummStr), &
2260                                     iDummy, length)
2261
2262        if(.not.(CheckMasks)) then ! Merely multiply iMask by iDummy:
2263           do j=1,length
2264              iMask(j) = iMask(j) * iDummy(j)
2265           end do
2266        else ! check mask elements and include their effect on iMask
2267           do j=1,length
2268              select case(iDummy(j))
2269              case(0) ! zeroes out iMask(j)
2270                 iMask(j) = 0
2271              case(1) ! leaves iMask(j) untouched
2272              case default ! shut down with an error
2273                 write(stderr,'(5a,i8,a,i8)') myname_, &
2274                      ':: ERROR--illegal mask value (must be 0 or 1).', &
2275                      'Illegal value stored in mask ', &
2276                      String_ToChar(DummStr),'(',j,')=',iDummy(j)
2277                 call die(myname_)
2278              end select
2279           end do
2280        endif ! if(CheckMasks)...
2281       ! clean up dummy String DummStr
2282        call String_clean(DummStr)
2283     end do ! do i=1,List_nitem(iMaskList)...
2284
2285  endif ! if(present(iMaskTags))...
2286
2287       ! Check for REAL masks.  If they are present, retrieve
2288       ! them and combine them into a single real mask rMask(:)
2289
2290  if(present(rMaskTags)) then
2291
2292       ! allocate two arrays:  rMask (the final product),
2293       ! and rDummy (storage space for each mask as it is retrieved)
2294
2295     allocate(rMask(AttrVect_lsize(inAv)), rDummy(AttrVect_lsize(inAv)), &
2296              stat=ierr)
2297     if(ierr /= 0) then
2298        write(stderr,'(2a,i8)') myname_, &
2299             ':: allocate(rMask(...)...) failed with ierr = ',ierr
2300        call die(myname_)
2301     endif
2302
2303       ! Initialize all the elements of rMask to unity:
2304     rMask = 1._FP
2305
2306       ! turn the colon-delimited string of tags into a List:
2307     call List_init(rMaskList,rMaskTags)
2308     
2309       ! Loop over the items in rMaskList, retrieving each mask
2310       ! into the array rDummy, checking it (if CheckMasks=.TRUE.),
2311       ! and multiplying it element-by-element into the array rMask.
2312
2313     do i=1,List_nitem(rMaskList)
2314       ! grab item as a String
2315        call List_get(DummStr, i, rMaskList) 
2316       ! use this String to identify an INTEGER GeneralGrid attribute
2317       ! for export to rDummy(:)
2318        call GeneralGrid_exportRAttr(GGrid, String_ToChar(DummStr), &
2319                                     rDummy, length)
2320
2321        if(.not.(CheckMasks)) then ! Merely multiply rMask by rDummy:
2322           do j=1,length
2323              rMask(j) = rMask(j) * rDummy(j)
2324           end do
2325        else ! check mask elements and include their effect on rMask
2326           do j=1,length
2327              if((iDummy(j) >= 0.) .and. (iDummy(j) <= 1.)) then ! in [0,1]
2328                 rMask(j) = rMask(j) * rDummy(j)
2329              else
2330                 write(stderr,'(5a,i8,a,i8)') myname_, &
2331                      ':: ERROR--illegal mask value (must be in [0.,1.]).', &
2332                      'Illegal value stored in mask ', &
2333                      String_ToChar(DummStr),'(',j,')=',rDummy(j)
2334                 call die(myname_)
2335              endif
2336           end do
2337        endif ! if(CheckMasks)...
2338       ! clean up dummy String DummStr
2339        call String_clean(DummStr)
2340     end do ! do i=1,List_nitem(rMaskList)...
2341
2342  endif ! if(present(rMaskTags))...
2343
2344       ! Now we have (at most) a single INTEGER mask iMask(:) and
2345       ! a single REAL mask rMask(:).  Before we perform the merge,
2346       ! we must tackle one more issue:  are the REAL attributes
2347       ! of inAv and outAv identical and in the same order?  If they
2348       ! are, the merge is a straightforward double loop over the
2349       ! elements and over all the attributes.  If the attribute lists
2350       ! differ, we must cross-reference common attributes, and store
2351       ! their indices.
2352
2353  RAttrIdentical = List_identical(inAv%rList, outAv%rList)
2354  if(.not.(RAttrIdentical)) then 
2355       ! Determine the number of shared REAL attributes NumSharedRAttr,
2356       ! and form cross-index tables inAvIndices, outAvIndices.
2357     call SharedAttrIndexList(inAv, outAv, 'REAL', NumSharedRAttr, &
2358                              inAvIndices, outAvIndices)
2359  endif
2360
2361  if(present(rMaskTags)) then ! REAL masking stored in rMask(:)
2362
2363     if(present(iMaskTags)) then ! also INTEGER mask iMask(:)
2364
2365        if(RAttrIdentical) then ! straight masked multiply
2366           do i=1, AttrVect_lsize(inAv)
2367              do j=1,AttrVect_nRAttr(inAv)
2368                 outAv%rAttr(j,i) = outAv%rAttr(j,i) + &
2369                      rMask(i) * iMask(i) * inAv%rAttr(j,i)
2370              end do ! do j=1,AttrVect_nRAttr(inAv)
2371       ! add in mask contribution to total of merge weights
2372              WeightSum(i) = WeightSum(i) + iMask(i) * rMask(i)
2373           end do ! do i=1,AttrVect_lsize(inAv)...
2374        else ! use previously generated cross-indices
2375           do i=1, AttrVect_lsize(inAv)
2376              do j=1,NumSharedRAttr
2377                 outAv%rAttr(outAVIndices(j),i) = &
2378                                 outAv%rAttr(outAvIndices(j),i) + &
2379                                 rMask(i) * iMask(i) * &
2380                                 inAv%rAttr(inAvIndices(j),i) 
2381              end do ! do j=1,NumSharedRAttr
2382       ! add in mask contribution to total of merge weights
2383              WeightSum(i) = WeightSum(i) + iMask(i) * rMask(i)
2384           end do ! do i=1,AttrVect_lsize(inAv)...
2385        endif ! if(RAttrIdentical)...
2386
2387     else ! rMask(:), but no iMask(:)
2388
2389        if(RAttrIdentical) then ! straight masked multiply
2390           do i=1, AttrVect_lsize(inAv)
2391              do j=1,AttrVect_nRAttr(inAv)
2392                 outAv%rAttr(j,i) = outAv%rAttr(j,i) + &
2393                                    rMask(i) * inAv%rAttr(j,i)
2394              end do ! do j=1,AttrVect_nRAttr(inAv)
2395       ! add in mask contribution to total of merge weights
2396              WeightSum(i) = WeightSum(i) + rMask(i)
2397           end do ! do i=1,AttrVect_lsize(inAv)...
2398        else ! use previously generated cross-indices
2399           do i=1, AttrVect_lsize(inAv)
2400              do j=1,NumSharedRAttr
2401                 outAv%rAttr(outAVIndices(j),i) = &
2402                                 outAv%rAttr(outAvIndices(j),i) + &
2403                                 rMask(i) * inAv%rAttr(inAvIndices(j),i) 
2404              end do ! do j=1,NumSharedRAttr
2405       ! add in mask contribution to total of merge weights
2406              WeightSum(i) = WeightSum(i) + rMask(i)
2407           end do ! do i=1,AttrVect_lsize(inAv)...
2408        endif ! if(RAttrIdentical)
2409
2410     endif ! if(present(iMaskTags))...
2411
2412  else ! No REAL Mask
2413
2414     if(present(iMaskTags)) then ! Have iMask(:), but no rMask(:)
2415
2416        if(RAttrIdentical) then ! straight masked multiply
2417           do i=1, AttrVect_lsize(inAv)
2418              do j=1,AttrVect_nRAttr(inAv)
2419                 outAv%rAttr(j,i) = outAv%rAttr(j,i) + &
2420                                    iMask(i) * inAv%rAttr(j,i)
2421              end do ! do j=1,AttrVect_nRAttr(inAv)
2422       ! add in mask contribution to total of merge weights
2423              WeightSum(i) = WeightSum(i) + iMask(i)
2424           end do ! do i=1,AttrVect_lsize(inAv)...
2425        else ! use previously generated cross-indices
2426           do i=1, AttrVect_lsize(inAv)
2427              do j=1,NumSharedRAttr
2428                 outAv%rAttr(outAVIndices(j),i) = &
2429                                 outAv%rAttr(outAvIndices(j),i) + &
2430                                 iMask(i) * inAv%rAttr(inAvIndices(j),i) 
2431              end do ! do j=1,NumSharedRAttr
2432       ! add in mask contribution to total of merge weights
2433              WeightSum(i) = WeightSum(i) + iMask(i)
2434           end do ! do i=1,AttrVect_lsize(inAv)...
2435        endif ! if(RAttrIdentical)
2436
2437     else ! Neither iMask(:) nor rMask(:)--all elements weighted by unity
2438
2439        if(RAttrIdentical) then ! straight masked multiply
2440           do i=1, AttrVect_lsize(inAv)
2441              do j=1,AttrVect_nRAttr(inAv)
2442                 outAv%rAttr(j,i) = outAv%rAttr(j,i) + inAv%rAttr(j,i)
2443              end do ! do j=1,AttrVect_nRAttr(inAv)
2444       ! add in mask contribution to total of merge weights
2445              WeightSum(i) = WeightSum(i) + 1._FP
2446           end do ! do i=1,AttrVect_lsize(inAv)...
2447        else ! use previously generated cross-indices
2448           do i=1, AttrVect_lsize(inAv)
2449              do j=1,NumSharedRAttr
2450                 outAv%rAttr(outAVIndices(j),i) = &
2451                                 outAv%rAttr(outAvIndices(j),i) + &
2452                                         inAv%rAttr(inAvIndices(j),i) 
2453              end do ! do j=1,NumSharedRAttr
2454       ! add in mask contribution to total of merge weights
2455              WeightSum(i) = WeightSum(i) + 1._FP
2456           end do ! do i=1,AttrVect_lsize(inAv)...
2457        endif ! if(RAttrIdentical)
2458
2459     endif ! if(present(iMaskTags))...
2460
2461  endif ! if(present(rMaskTags))...
2462
2463       ! At this point the merge has been completed.  Now clean
2464       ! up all allocated structures and temporary arrays.
2465
2466  if(present(iMaskTags)) then ! clean up integer mask work space
2467     deallocate(iMask, iDummy, stat=ierr)
2468     if(ierr /= 0) then
2469        write(stderr,'(2a,i8)') myname_, &
2470             ':: deallocate(iMask,...) failed with ierr = ',ierr
2471        call die(myname_)
2472     endif
2473     call List_clean(iMaskList)
2474  endif
2475
2476  if(present(rMaskTags)) then ! clean up real mask work space
2477     deallocate(rMask, rDummy, stat=ierr)
2478     if(ierr /= 0) then
2479        write(stderr,'(2a,i8)') myname_, &
2480             ':: deallocate(rMask,...) failed with ierr = ',ierr
2481        call die(myname_)
2482     endif
2483     call List_clean(rMaskList)
2484  endif
2485
2486  if(.not.(RAttrIdentical)) then ! clean up cross-reference tables
2487     deallocate(inAvIndices, outAvIndices, stat=ierr)
2488     if(ierr /= 0) then
2489        write(stderr,'(2a,i8)') myname_, &
2490             ':: deallocate(inAvIndices,...) failed with ierr = ',ierr
2491        call die(myname_)
2492     endif
2493  endif
2494
2495 end subroutine MergeInDataGGSP_
2496
2497!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
2498!    Math and Computer Science Division, Argonne National Laboratory   !
2499!-----------------------------------------------------------------------
2500!
2501! !IROUTINE: MergeInDataGGDP_ - merge in data from a component.
2502!
2503! !DESCRIPTION:
2504! Double precision version of MergeInDataGGSP_
2505!
2506! !INTERFACE:
2507
2508 subroutine MergeInDataGGDP_(inAv, iMaskTags, rMaskTags, GGrid, &
2509                             CheckMasks, outAv, WeightSum)
2510!
2511! !USES:
2512!
2513      use m_stdio
2514      use m_die
2515
2516      use m_realkinds, only : DP, FP
2517
2518      use m_String, only : String
2519      use m_String, only : String_clean => clean
2520      use m_String, only : String_ToChar => toChar
2521
2522      use m_List, only : List
2523      use m_List, only : List_init => init
2524      use m_List, only : List_clean => clean
2525      use m_List, only : List_nitem => nitem
2526      use m_List, only : List_get => get
2527      use m_List, only : List_identical => identical
2528      use m_List, only : List_allocated => allocated
2529
2530      use m_AttrVect, only : AttrVect
2531      use m_AttrVect, only : AttrVect_lsize => lsize
2532      use m_AttrVect, only : AttrVect_nRAttr => nRAttr
2533      use m_AttrVect, only : SharedAttrIndexList
2534
2535      use m_GeneralGrid, only : GeneralGrid
2536      use m_GeneralGrid, only : GeneralGrid_lsize => lsize
2537      use m_GeneralGrid, only : GeneralGrid_exportIAttr => exportIAttr
2538      use m_GeneralGrid, only : GeneralGrid_exportRAttr => exportRAttr
2539
2540      implicit none
2541
2542! !INPUT PARAMETERS:
2543!
2544      type(AttrVect),               intent(IN)    :: inAv
2545      character(len=*),   optional, intent(IN)    :: iMaskTags
2546      character(len=*),   optional, intent(IN)    :: rMaskTags
2547      type(GeneralGrid),            intent(IN)    :: GGrid
2548      logical,                      intent(IN)    :: CheckMasks
2549
2550! !INPUT/OUTPUT PARAMETERS:
2551!
2552      type(AttrVect),               intent(INOUT) :: outAv
2553      real(DP),       dimension(:), pointer       :: WeightSum
2554
2555! !REVISION HISTORY:
2556!       19Jun02 - Jay Larson <larson@mcs.anl.gov> - initial verson.
2557!       10Jul02 - J. Larson <larson@mcs.anl.gov> - Improved argument
2558!                 checking.
2559!_______________________________________________________________________
2560!
2561  character(len=*),parameter :: myname_=myname//'::MergeInDataGGDP_'
2562
2563  integer :: i, ierr, j, length
2564  type(String) :: DummStr
2565  type(List) :: iMaskList, rMaskList
2566  integer,  dimension(:), pointer :: iMask,iDummy ! INTEGER mask workspace
2567  real(FP), dimension(:), pointer :: rMask,rDummy ! REAL mask workspace
2568
2569  logical :: RAttrIdentical ! flag to identify identical REAL attribute
2570                            ! lists in inAv and outAv
2571  integer :: NumSharedRAttr ! number of REAL attributes shared by inAv,outAv
2572       ! Cross-index storage for shared REAL attributes of inAv,outAv
2573  integer, dimension(:), pointer :: inAvIndices, outAvIndices
2574
2575       ! Begin argument sanity checks...
2576
2577       ! Have the input arguments been allocated?
2578
2579  if(.not.(List_allocated(inAv%rList))) then
2580     write(stderr,'(2a)') myname_, &
2581          'ERROR--INPUT argument inAv has no real attributes.'
2582     call die(myname_)
2583  endif
2584
2585  if(.not.(List_allocated(outaV%rList))) then
2586     write(stderr,'(2a)') myname_, &
2587          'ERROR--INPUT/OUTPUT argument outAv has no real attributes.'
2588     call die(myname_)
2589  endif
2590
2591  if(present(iMaskTags)) then
2592     if(.not.(List_allocated(GGrid%data%iList))) then
2593        write(stderr,'(3a)') myname_, &
2594             'ERROR--Integer masking requested, but input argument GGrid ', &
2595             'has no integer attributes.'
2596        call die(myname_)
2597     endif
2598  endif
2599
2600  if(present(rMaskTags)) then
2601     if(.not.(List_allocated(GGrid%data%rList))) then
2602        write(stderr,'(3a)') myname_, &
2603             'ERROR--Real masking requested, but input argument GGrid ', &
2604             'has no real attributes.'
2605        call die(myname_)
2606     endif
2607  endif
2608
2609  if(.not.(associated(WeightSum))) then
2610     write(stderr,'(2a)') myname_, &
2611          'ERROR--INPUT/OUPUT argument WeightSum has not been allocated.'
2612     call die(myname_)
2613  endif
2614
2615       ! Do the vector lengths match?
2616
2617  if(AttrVect_lsize(inAv) /= AttrVect_lsize(outAv)) then
2618     write(stderr,'(2a,2(a,i8))') myname_, &
2619          ':: ERROR--Lengths of AttrVect arguments inAv and outAv must match.', &
2620          'AttrVect_lsize(inAv) = ',AttrVect_lsize(inAv), &
2621          'AttrVect_lsize(outAv) = ',AttrVect_lsize(outAv)
2622     call die(myname_)
2623  endif
2624
2625  if(AttrVect_lsize(inAv) /= GeneralGrid_lsize(GGrid)) then
2626     write(stderr,'(2a,2(a,i8))') myname_, &
2627          ':: ERROR--Lengths of arguments inAv and GGrid must match.', &
2628          'AttrVect_lsize(inAv) = ',AttrVect_lsize(inAv), &
2629          'AttrVect_lsize(outAv) = ',GeneralGrid_lsize(GGrid)
2630     call die(myname_)
2631  endif
2632
2633  if(AttrVect_lsize(inAv) /= size(WeightSum)) then
2634     write(stderr,'(2a,2(a,i8))') myname_, &
2635          ':: ERROR--Lengths of arguments inAv and WeightSum must match.', &
2636          'AttrVect_lsize(inAv) = ',AttrVect_lsize(inAv), &
2637          'size(WeightSum) = ',size(WeightSum)
2638     call die(myname_)
2639  endif
2640
2641       ! ...end argument sanity checks.
2642
2643       ! Check for INTEGER masks.  If they are present, retrieve
2644       ! them and combine them into a single integer mask iMask(:)
2645
2646  if(present(iMaskTags)) then
2647
2648       ! allocate two arrays:  iMask (the final product),
2649       ! and iDummy (storage space for each mask as it is retrieved)
2650
2651     allocate(iMask(AttrVect_lsize(inAv)), iDummy(AttrVect_lsize(inAv)), &
2652              stat=ierr)
2653     if(ierr /= 0) then
2654        write(stderr,'(2a,i8)') myname_, &
2655             ':: allocate(iMask(...)...) failed with ierr = ',ierr
2656        call die(myname_)
2657     endif
2658
2659       ! Initialize all the elements of iMask to unity:
2660     iMask = 1
2661
2662       ! turn the colon-delimited string of tags into a List:
2663     call List_init(iMaskList,iMaskTags)
2664     
2665       ! Loop over the items in iMaskList, retrieving each mask
2666       ! into the array iDummy, checking it (if CheckMasks=.TRUE.),
2667       ! and multiplying it element-by-element into the array iMask.
2668
2669     do i=1,List_nitem(iMaskList)
2670       ! grab item as a String
2671        call List_get(DummStr, i, iMaskList) 
2672       ! use this String to identify an INTEGER GeneralGrid attribute
2673       ! for export to iDummy(:)
2674        call GeneralGrid_exportIAttr(GGrid, String_ToChar(DummStr), &
2675                                     iDummy, length)
2676
2677        if(.not.(CheckMasks)) then ! Merely multiply iMask by iDummy:
2678           do j=1,length
2679              iMask(j) = iMask(j) * iDummy(j)
2680           end do
2681        else ! check mask elements and include their effect on iMask
2682           do j=1,length
2683              select case(iDummy(j))
2684              case(0) ! zeroes out iMask(j)
2685                 iMask(j) = 0
2686              case(1) ! leaves iMask(j) untouched
2687              case default ! shut down with an error
2688                 write(stderr,'(5a,i8,a,i8)') myname_, &
2689                      ':: ERROR--illegal mask value (must be 0 or 1).', &
2690                      'Illegal value stored in mask ', &
2691                      String_ToChar(DummStr),'(',j,')=',iDummy(j)
2692                 call die(myname_)
2693              end select
2694           end do
2695        endif ! if(CheckMasks)...
2696       ! clean up dummy String DummStr
2697        call String_clean(DummStr)
2698     end do ! do i=1,List_nitem(iMaskList)...
2699
2700  endif ! if(present(iMaskTags))...
2701
2702       ! Check for REAL masks.  If they are present, retrieve
2703       ! them and combine them into a single real mask rMask(:)
2704
2705  if(present(rMaskTags)) then
2706
2707       ! allocate two arrays:  rMask (the final product),
2708       ! and rDummy (storage space for each mask as it is retrieved)
2709
2710     allocate(rMask(AttrVect_lsize(inAv)), rDummy(AttrVect_lsize(inAv)), &
2711              stat=ierr)
2712     if(ierr /= 0) then
2713        write(stderr,'(2a,i8)') myname_, &
2714             ':: allocate(rMask(...)...) failed with ierr = ',ierr
2715        call die(myname_)
2716     endif
2717
2718       ! Initialize all the elements of rMask to unity:
2719     rMask = 1._FP
2720
2721       ! turn the colon-delimited string of tags into a List:
2722     call List_init(rMaskList,rMaskTags)
2723     
2724       ! Loop over the items in rMaskList, retrieving each mask
2725       ! into the array rDummy, checking it (if CheckMasks=.TRUE.),
2726       ! and multiplying it element-by-element into the array rMask.
2727
2728     do i=1,List_nitem(rMaskList)
2729       ! grab item as a String
2730        call List_get(DummStr, i, rMaskList) 
2731       ! use this String to identify an INTEGER GeneralGrid attribute
2732       ! for export to rDummy(:)
2733        call GeneralGrid_exportRAttr(GGrid, String_ToChar(DummStr), &
2734                                     rDummy, length)
2735
2736        if(.not.(CheckMasks)) then ! Merely multiply rMask by rDummy:
2737           do j=1,length
2738              rMask(j) = rMask(j) * rDummy(j)
2739           end do
2740        else ! check mask elements and include their effect on rMask
2741           do j=1,length
2742              if((iDummy(j) >= 0.) .and. (iDummy(j) <= 1.)) then ! in [0,1]
2743                 rMask(j) = rMask(j) * rDummy(j)
2744              else
2745                 write(stderr,'(5a,i8,a,i8)') myname_, &
2746                      ':: ERROR--illegal mask value (must be in [0.,1.]).', &
2747                      'Illegal value stored in mask ', &
2748                      String_ToChar(DummStr),'(',j,')=',rDummy(j)
2749                 call die(myname_)
2750              endif
2751           end do
2752        endif ! if(CheckMasks)...
2753       ! clean up dummy String DummStr
2754        call String_clean(DummStr)
2755     end do ! do i=1,List_nitem(rMaskList)...
2756
2757  endif ! if(present(rMaskTags))...
2758
2759       ! Now we have (at most) a single INTEGER mask iMask(:) and
2760       ! a single REAL mask rMask(:).  Before we perform the merge,
2761       ! we must tackle one more issue:  are the REAL attributes
2762       ! of inAv and outAv identical and in the same order?  If they
2763       ! are, the merge is a straightforward double loop over the
2764       ! elements and over all the attributes.  If the attribute lists
2765       ! differ, we must cross-reference common attributes, and store
2766       ! their indices.
2767
2768  RAttrIdentical = List_identical(inAv%rList, outAv%rList)
2769  if(.not.(RAttrIdentical)) then 
2770       ! Determine the number of shared REAL attributes NumSharedRAttr,
2771       ! and form cross-index tables inAvIndices, outAvIndices.
2772     call SharedAttrIndexList(inAv, outAv, 'REAL', NumSharedRAttr, &
2773                              inAvIndices, outAvIndices)
2774  endif
2775
2776  if(present(rMaskTags)) then ! REAL masking stored in rMask(:)
2777
2778     if(present(iMaskTags)) then ! also INTEGER mask iMask(:)
2779
2780        if(RAttrIdentical) then ! straight masked multiply
2781           do i=1, AttrVect_lsize(inAv)
2782              do j=1,AttrVect_nRAttr(inAv)
2783                 outAv%rAttr(j,i) = outAv%rAttr(j,i) + &
2784                      rMask(i) * iMask(i) * inAv%rAttr(j,i)
2785              end do ! do j=1,AttrVect_nRAttr(inAv)
2786       ! add in mask contribution to total of merge weights
2787              WeightSum(i) = WeightSum(i) + iMask(i) * rMask(i)
2788           end do ! do i=1,AttrVect_lsize(inAv)...
2789        else ! use previously generated cross-indices
2790           do i=1, AttrVect_lsize(inAv)
2791              do j=1,NumSharedRAttr
2792                 outAv%rAttr(outAVIndices(j),i) = &
2793                                 outAv%rAttr(outAvIndices(j),i) + &
2794                                 rMask(i) * iMask(i) * &
2795                                 inAv%rAttr(inAvIndices(j),i) 
2796              end do ! do j=1,NumSharedRAttr
2797       ! add in mask contribution to total of merge weights
2798              WeightSum(i) = WeightSum(i) + iMask(i) * rMask(i)
2799           end do ! do i=1,AttrVect_lsize(inAv)...
2800        endif ! if(RAttrIdentical)...
2801
2802     else ! rMask(:), but no iMask(:)
2803
2804        if(RAttrIdentical) then ! straight masked multiply
2805           do i=1, AttrVect_lsize(inAv)
2806              do j=1,AttrVect_nRAttr(inAv)
2807                 outAv%rAttr(j,i) = outAv%rAttr(j,i) + &
2808                                    rMask(i) * inAv%rAttr(j,i)
2809              end do ! do j=1,AttrVect_nRAttr(inAv)
2810       ! add in mask contribution to total of merge weights
2811              WeightSum(i) = WeightSum(i) + rMask(i)
2812           end do ! do i=1,AttrVect_lsize(inAv)...
2813        else ! use previously generated cross-indices
2814           do i=1, AttrVect_lsize(inAv)
2815              do j=1,NumSharedRAttr
2816                 outAv%rAttr(outAVIndices(j),i) = &
2817                                 outAv%rAttr(outAvIndices(j),i) + &
2818                                 rMask(i) * inAv%rAttr(inAvIndices(j),i) 
2819              end do ! do j=1,NumSharedRAttr
2820       ! add in mask contribution to total of merge weights
2821              WeightSum(i) = WeightSum(i) + rMask(i)
2822           end do ! do i=1,AttrVect_lsize(inAv)...
2823        endif ! if(RAttrIdentical)
2824
2825     endif ! if(present(iMaskTags))...
2826
2827  else ! No REAL Mask
2828
2829     if(present(iMaskTags)) then ! Have iMask(:), but no rMask(:)
2830
2831        if(RAttrIdentical) then ! straight masked multiply
2832           do i=1, AttrVect_lsize(inAv)
2833              do j=1,AttrVect_nRAttr(inAv)
2834                 outAv%rAttr(j,i) = outAv%rAttr(j,i) + &
2835                                    iMask(i) * inAv%rAttr(j,i)
2836              end do ! do j=1,AttrVect_nRAttr(inAv)
2837       ! add in mask contribution to total of merge weights
2838              WeightSum(i) = WeightSum(i) + iMask(i)
2839           end do ! do i=1,AttrVect_lsize(inAv)...
2840        else ! use previously generated cross-indices
2841           do i=1, AttrVect_lsize(inAv)
2842              do j=1,NumSharedRAttr
2843                 outAv%rAttr(outAVIndices(j),i) = &
2844                                 outAv%rAttr(outAvIndices(j),i) + &
2845                                 iMask(i) * inAv%rAttr(inAvIndices(j),i) 
2846              end do ! do j=1,NumSharedRAttr
2847       ! add in mask contribution to total of merge weights
2848              WeightSum(i) = WeightSum(i) + iMask(i)
2849           end do ! do i=1,AttrVect_lsize(inAv)...
2850        endif ! if(RAttrIdentical)
2851
2852     else ! Neither iMask(:) nor rMask(:)--all elements weighted by unity
2853
2854        if(RAttrIdentical) then ! straight masked multiply
2855           do i=1, AttrVect_lsize(inAv)
2856              do j=1,AttrVect_nRAttr(inAv)
2857                 outAv%rAttr(j,i) = outAv%rAttr(j,i) + inAv%rAttr(j,i)
2858              end do ! do j=1,AttrVect_nRAttr(inAv)
2859       ! add in mask contribution to total of merge weights
2860              WeightSum(i) = WeightSum(i) + 1._FP
2861           end do ! do i=1,AttrVect_lsize(inAv)...
2862        else ! use previously generated cross-indices
2863           do i=1, AttrVect_lsize(inAv)
2864              do j=1,NumSharedRAttr
2865                 outAv%rAttr(outAVIndices(j),i) = &
2866                                 outAv%rAttr(outAvIndices(j),i) + &
2867                                         inAv%rAttr(inAvIndices(j),i) 
2868              end do ! do j=1,NumSharedRAttr
2869       ! add in mask contribution to total of merge weights
2870              WeightSum(i) = WeightSum(i) + 1._FP
2871           end do ! do i=1,AttrVect_lsize(inAv)...
2872        endif ! if(RAttrIdentical)
2873
2874     endif ! if(present(iMaskTags))...
2875
2876  endif ! if(present(rMaskTags))...
2877
2878       ! At this point the merge has been completed.  Now clean
2879       ! up all allocated structures and temporary arrays.
2880
2881  if(present(iMaskTags)) then ! clean up integer mask work space
2882     deallocate(iMask, iDummy, stat=ierr)
2883     if(ierr /= 0) then
2884        write(stderr,'(2a,i8)') myname_, &
2885             ':: deallocate(iMask,...) failed with ierr = ',ierr
2886        call die(myname_)
2887     endif
2888     call List_clean(iMaskList)
2889  endif
2890
2891  if(present(rMaskTags)) then ! clean up real mask work space
2892     deallocate(rMask, rDummy, stat=ierr)
2893     if(ierr /= 0) then
2894        write(stderr,'(2a,i8)') myname_, &
2895             ':: deallocate(rMask,...) failed with ierr = ',ierr
2896        call die(myname_)
2897     endif
2898     call List_clean(rMaskList)
2899  endif
2900
2901  if(.not.(RAttrIdentical)) then ! clean up cross-reference tables
2902     deallocate(inAvIndices, outAvIndices, stat=ierr)
2903     if(ierr /= 0) then
2904        write(stderr,'(2a,i8)') myname_, &
2905             ':: deallocate(inAvIndices,...) failed with ierr = ',ierr
2906        call die(myname_)
2907     endif
2908  endif
2909
2910 end subroutine MergeInDataGGDP_
2911
2912 end module m_Merge
Note: See TracBrowser for help on using the repository browser.