/[MITgcm]/MITgcm/pkg/diagnostics/diagstats_local.F
ViewVC logotype

Annotation of /MITgcm/pkg/diagnostics/diagstats_local.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph


Revision 1.1 - (hide annotations) (download)
Fri May 20 07:28:52 2005 UTC (18 years, 11 months ago) by jmc
Branch: MAIN
Add new capability: compute & write Global/Regional & per level statistics

1 jmc 1.1 C $Header: $
2     C $Name: $
3    
4     #include "DIAG_OPTIONS.h"
5    
6     CBOP
7     C !ROUTINE: DIAGSTATS_LOCAL
8     C !INTERFACE:
9     SUBROUTINE DIAGSTATS_LOCAL(
10     U statFld,
11     I inpFld,
12     I sizI1,sizI2,sizJ1,sizJ2,sizK,sizTx,sizTy,
13     I iRun,jRun,kIn,biIn,bjIn,
14     I k,bi,bj, region2fill, ndId, parsFld,
15     I myThid)
16    
17     C !DESCRIPTION:
18     C Update array statFld
19     C by adding statistics over the range [1:iRun],[1:jRun]
20     C from input field array inpFld
21     C- note:
22     C a) this S/R should not see DIAGNOSTICS pkg commons blocks (in DIAGNOSTICS.h)
23     C b) for main grid variables, get area & weigting factors (to compute global mean)
24     C from the main common blocks.
25     C c) for other type of grids, call a specifics S/R which include its own
26     C grid common blocks
27    
28     C !USES:
29     IMPLICIT NONE
30    
31     #include "EEPARAMS.h"
32     #include "SIZE.h"
33     #include "DIAGNOSTICS_SIZE.h"
34     #include "PARAMS.h"
35     #include "GRID.h"
36     #include "SURFACE.h"
37     #ifdef ALLOW_FIZHI
38     #include "fizhi_SIZE.h"
39     #include "gridalt_mapping.h"
40     #endif
41    
42     C !INPUT/OUTPUT PARAMETERS:
43     C == Routine Arguments ==
44     C statFld :: cumulative statistics array (updated)
45     C inpFld :: input field array to process (compute stats & add to statFld)
46     C sizI1,sizI2 :: size of inpFld array: 1rst index range (min,max)
47     C sizJ1,sizJ2 :: size of inpFld array: 2nd index range (min,max)
48     C sizK :: size of inpFld array: 3rd dimension
49     C sizTx,sizTy :: size of inpFld array: tile dimensions
50     C iRun,jRun :: range of 1rst & 2nd index
51     C kIn :: level index of inFld array to porcess
52     C biIn,bjIn :: tile indices of inFld array to process
53     C k,bi,bj :: level and tile indices used for weighting (mask,area ...)
54     C region2fill :: indicates whether to compute statistics over this region
55     C ndId :: Diagnostics Id Number (in available diag. list)
56     C parsFld :: parser field with characteristics of the diagnostics
57     C myThid :: my Thread Id number
58     _RL statFld(0:nStats,0:nRegions)
59     INTEGER sizI1,sizI2,sizJ1,sizJ2
60     INTEGER sizK,sizTx,sizTy
61     _RL inpFld(sizI1:sizI2,sizJ1:sizJ2,sizK,sizTx,sizTy)
62     INTEGER iRun, jRun, kIn, biIn, bjIn
63     INTEGER k, bi, bj, ndId
64     INTEGER region2fill(0:nRegions)
65     CHARACTER*16 parsFld
66     INTEGER myThid
67     CEOP
68    
69     C !LOCAL VARIABLES:
70     C i,j :: loop indices
71     INTEGER i, n, km
72     INTEGER im, ix, iv
73     PARAMETER ( iv = nStats - 2 , im = nStats - 1 , ix = nStats )
74     LOGICAL exclSpVal
75     LOGICAL useWeight
76     _RL statLoc(0:nStats)
77     _RL drLoc
78     _RL specialVal
79     _RL getcon
80     EXTERNAL getcon
81    
82     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
83    
84     useWeight = .FALSE.
85     exclSpVal = .FALSE.
86     specialVal = 0.
87     IF ( useFIZHI ) THEN
88     exclSpVal = .TRUE.
89     specialVal = getcon('UNDEF')
90     ENDIF
91    
92     DO n=0,nRegions
93     IF (region2fill(n).NE.0) THEN
94     C--- Compute statistics for this tile, level and region:
95    
96     C- case of a special region (no specific regional mask)
97     IF ( n.EQ.0 ) THEN
98    
99     IF ( parsFld(10:10) .EQ. 'R' ) THEN
100    
101     drLoc = drF(k)
102     IF ( parsFld(9:9).EQ.'L') drLoc = drC(k)
103     IF ( parsFld(9:9).EQ.'U') drLoc = drC(MIN(k+1,Nr))
104     IF ( parsFld(9:9).EQ.'M') useWeight = .TRUE.
105    
106     IF ( parsFld(2:2).EQ.'U' ) THEN
107     CALL DIAGSTATS_R_CALC(
108     O statLoc,
109     I inpFld(sizI1,sizJ1,kIn,biIn,bjIn),
110     I sizI1,sizI2,sizJ1,sizJ2,iRun,jRun,
111     I maskH(1-Olx,1-Oly,bi,bj), maskW(1-Olx,1-Oly,k,bi,bj),
112     I hFacW(1-Olx,1-Oly,k,bi,bj), rAw(1-Olx,1-Oly,bi,bj),
113     I drLoc, specialVal, exclSpVal, useWeight, myThid )
114     c I drLoc, k,bi,bj, parsFld, myThid )
115     ELSEIF ( parsFld(2:2).EQ.'V' ) THEN
116     CALL DIAGSTATS_R_CALC(
117     O statLoc,
118     I inpFld(sizI1,sizJ1,kIn,biIn,bjIn),
119     I sizI1,sizI2,sizJ1,sizJ2,iRun,jRun,
120     I maskH(1-Olx,1-Oly,bi,bj), maskS(1-Olx,1-Oly,k,bi,bj),
121     I hFacS(1-Olx,1-Oly,k,bi,bj), rAs(1-Olx,1-Oly,bi,bj),
122     I drLoc, specialVal, exclSpVal, useWeight, myThid )
123     ELSE
124     CALL DIAGSTATS_R_CALC(
125     O statLoc,
126     I inpFld(sizI1,sizJ1,kIn,biIn,bjIn),
127     I sizI1,sizI2,sizJ1,sizJ2,iRun,jRun,
128     I maskH(1-Olx,1-Oly,bi,bj), maskC(1-Olx,1-Oly,k,bi,bj),
129     I hFacC(1-Olx,1-Oly,k,bi,bj), rA(1-Olx,1-Oly,bi,bj),
130     I drLoc, specialVal, exclSpVal, useWeight, myThid )
131     ENDIF
132    
133     #ifdef ALLOW_FIZHI
134     c ELSEIF ( parsFld(10:10) .EQ. 'L' ) THEN
135     ELSEIF ( parsFld(10:10) .EQ. 'M' ) THEN
136     drLoc = 1. _d 0
137     km = 1 + Nrphys - k
138     CALL DIAGSTATS_R_CALC(
139     O statLoc,
140     I inpFld(sizI1,sizJ1,kIn,biIn,bjIn),
141     I sizI1,sizI2,sizJ1,sizJ2,iRun,jRun,
142     I maskH(1-Olx,1-Oly,bi,bj), maskH(1-Olx,1-Oly,bi,bj),
143     I dpphys(1-Olx,1-Oly,km,bi,bj), rA(1-Olx,1-Oly,bi,bj),
144     I drLoc, specialVal, exclSpVal, useWeight, myThid )
145     #endif
146     #ifdef ALLOW_LAND
147     c ELSEIF ( parsFld(10:10) .EQ. 'G' ) THEN
148     #endif
149     c ELSEIF ( parsFld(10:10) .EQ. 'I' ) THEN
150     c ELSEIF ( parsFld(10:10) .EQ. '1' ) THEN
151     ELSE
152    
153     km = 1
154     IF ( usingPCoords ) km = Nr
155     drLoc = 1. _d 0
156     IF ( parsFld(2:2).EQ.'U' ) THEN
157     CALL DIAGSTATS_R_CALC(
158     O statLoc,
159     I inpFld(sizI1,sizJ1,kIn,biIn,bjIn),
160     I sizI1,sizI2,sizJ1,sizJ2,iRun,jRun,
161     I maskH(1-Olx,1-Oly,bi,bj), maskW(1-Olx,1-Oly,km,bi,bj),
162     I maskW(1-Olx,1-Oly,km,bi,bj),rAw(1-Olx,1-Oly,bi,bj),
163     I drLoc, specialVal, exclSpVal, useWeight, myThid )
164     ELSEIF ( parsFld(2:2).EQ.'V' ) THEN
165     CALL DIAGSTATS_R_CALC(
166     O statLoc,
167     I inpFld(sizI1,sizJ1,kIn,biIn,bjIn),
168     I sizI1,sizI2,sizJ1,sizJ2,iRun,jRun,
169     I maskH(1-Olx,1-Oly,bi,bj), maskS(1-Olx,1-Oly,km,bi,bj),
170     I maskS(1-Olx,1-Oly,km,bi,bj),rAs(1-Olx,1-Oly,bi,bj),
171     I drLoc, specialVal, exclSpVal, useWeight, myThid )
172     ELSE
173     CALL DIAGSTATS_R_CALC(
174     O statLoc,
175     I inpFld(sizI1,sizJ1,kIn,biIn,bjIn),
176     I sizI1,sizI2,sizJ1,sizJ2,iRun,jRun,
177     I maskH(1-Olx,1-Oly,bi,bj), maskH(1-Olx,1-Oly,bi,bj),
178     I maskH(1-Olx,1-Oly,bi,bj), rA(1-Olx,1-Oly,bi,bj),
179     I drLoc, specialVal, exclSpVal, useWeight, myThid )
180     ENDIF
181    
182     ENDIF
183    
184     C Update cumulative statistics array
185     IF ( statLoc(0).GT.0. ) THEN
186     IF ( statFld(0,n).LE.0. ) THEN
187     statFld(im,n) = statLoc(im)
188     statFld(ix,n) = statLoc(ix)
189     ELSE
190     statFld(im,n) = MIN( statFld(im,n), statLoc(im) )
191     statFld(ix,n) = MAX( statFld(ix,n), statLoc(ix) )
192     ENDIF
193     DO i=0,iv
194     statFld(i,n) = statFld(i,n) + statLoc(i)
195     ENDDO
196     ENDIF
197    
198     ENDIF
199    
200     C--- processing region "n" ends here.
201     ENDIF
202     ENDDO
203    
204     RETURN
205     END
206    
207     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
208    
209     CBOP
210     C !ROUTINE: DIAGSTATS_CALC
211     C !INTERFACE:
212     SUBROUTINE DIAGSTATS_R_CALC(
213     U statArr,
214     I inpArr,
215     I sizI1,sizI2,sizJ1,sizJ2, iRun,jRun,
216     I regMask, arrMask, arrhFac, arrArea,
217     I arrDr, specialVal, exclSpVal, useWeight,
218     I myThid )
219     c I arrDr, k,bi,bj, parsFld, myThid )
220    
221     C !DESCRIPTION:
222     C Compute statistics for this tile, level, region
223    
224     C !USES:
225     IMPLICIT NONE
226    
227     #include "EEPARAMS.h"
228     #include "SIZE.h"
229     #include "DIAGNOSTICS_SIZE.h"
230     c #include "PARAMS.h"
231     c #include "GRID.h"
232     c #include "SURFACE.h"
233    
234     C !INPUT/OUTPUT PARAMETERS:
235     C == Routine Arguments ==
236     C statArr :: cumulative statistics array (updated)
237     C inpArr :: input field array to process (compute stats & add to statFld)
238     C sizI1,sizI2 :: size of inpArr array: 1rst index range (min,max)
239     C sizJ1,sizJ2 :: size of inpArr array: 2nd index range (min,max)
240     C iRun,jRun :: range of 1rst & 2nd index to process
241     C regMask :: regional mask
242     C arrMask :: mask for this input array
243     C arrhFac :: weight factor (horizontally varying)
244     C arrArea :: Area weighting factor
245     C arrDr :: uniform weighting factor
246     C specialVal :: special value in input array (to exclude if exclSpVal=T)
247     C exclSpVal :: if T, exclude "specialVal" in input array
248     C useWeight :: use weight factor "arrhFac"
249     Cc k,bi,bj :: level and tile indices used for weighting (mask,area ...)
250     Cc parsFld :: parser field with characteristics of the diagnostics
251     C myThid :: my Thread Id number
252     _RL statArr(0:nStats)
253     INTEGER sizI1,sizI2,sizJ1,sizJ2
254     INTEGER iRun, jRun
255     _RL inpArr (sizI1:sizI2,sizJ1:sizJ2)
256     _RS regMask(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
257     _RS arrMask(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
258     _RS arrhFac(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
259     _RS arrArea(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
260     _RL arrDr
261     _RL specialVal
262     LOGICAL exclSpVal
263     LOGICAL useWeight
264     c INTEGER k, bi, bj
265     c CHARACTER*16 parsFld
266     INTEGER myThid
267     CEOP
268    
269     C !LOCAL VARIABLES:
270     C i,j :: loop indices
271     INTEGER i, j, n
272     INTEGER im, ix, iv
273     PARAMETER ( iv = nStats - 2 , im = nStats - 1 , ix = nStats )
274     _RL tmpVol
275    
276     DO n=0,nStats
277     statArr(n) = 0.
278     ENDDO
279    
280     IF ( exclSpVal ) THEN
281    
282     DO j = 1,jRun
283     DO i = 1,iRun
284     IF (arrMask(i,j).NE.0. .AND. inpArr(i,j).NE.specialVal) THEN
285     IF ( statArr(0).EQ.0. ) THEN
286     statArr(im) = inpArr(i,j)
287     statArr(ix) = inpArr(i,j)
288     ELSE
289     statArr(im) = MIN(inpArr(i,j),statArr(im))
290     statArr(ix) = MAX(inpArr(i,j),statArr(ix))
291     ENDIF
292     IF ( useWeight ) THEN
293     tmpVol = arrDr*arrhFac(i,j)*arrArea(i,j)
294     ELSE
295     tmpVol = arrDr*arrArea(i,j)
296     ENDIF
297     statArr(0) = statArr(0) + tmpVol
298     statArr(1) = statArr(1) + tmpVol*inpArr(i,j)
299     statArr(2) = statArr(2) + tmpVol*inpArr(i,j)*inpArr(i,j)
300     ENDIF
301     ENDDO
302     ENDDO
303    
304     ELSE
305    
306     DO j = 1,jRun
307     DO i = 1,iRun
308     c IF ( regMask(i,j).NE.0. .AND. arrMask(i,j).NE.0. ) THEN
309     IF ( arrMask(i,j).NE.0. ) THEN
310     IF ( statArr(0).EQ.0. ) THEN
311     statArr(im) = inpArr(i,j)
312     statArr(ix) = inpArr(i,j)
313     ELSE
314     statArr(im) = MIN(inpArr(i,j),statArr(im))
315     statArr(ix) = MAX(inpArr(i,j),statArr(ix))
316     ENDIF
317     IF ( useWeight ) THEN
318     tmpVol = arrDr*arrhFac(i,j)*arrArea(i,j)
319     ELSE
320     tmpVol = arrDr*arrArea(i,j)
321     ENDIF
322     statArr(0) = statArr(0) + tmpVol
323     statArr(1) = statArr(1) + tmpVol*inpArr(i,j)
324     statArr(2) = statArr(2) + tmpVol*inpArr(i,j)*inpArr(i,j)
325     ENDIF
326     ENDDO
327     ENDDO
328    
329     ENDIF
330    
331     RETURN
332     END

  ViewVC Help
Powered by ViewVC 1.1.22