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

Contents of /MITgcm/pkg/diagnostics/diagstats_global.F

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


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

1 C $Header: $
2 C $Name: $
3
4 #include "DIAG_OPTIONS.h"
5
6 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
7 CBOP 0
8 C !ROUTINE: DIAGSTATS_GLOBAL
9
10 C !INTERFACE:
11 SUBROUTINE DIAGSTATS_GLOBAL(
12 O qtmp1, qtmp2,
13 I undef, nLev, jReg,
14 I ndId, mate, myThid )
15
16 C !DESCRIPTION:
17 C Retrieve averaged model diagnostic
18
19 C !USES:
20 IMPLICIT NONE
21 #include "EEPARAMS.h"
22 #include "SIZE.h"
23 #include "DIAGNOSTICS_SIZE.h"
24 #include "DIAGNOSTICS.h"
25
26 C !INPUT PARAMETERS:
27 C undef ..... Undefined value
28 C nLev .... 2nd Dimension (max Nb of levels) of qtmp1,2 arrays
29 C jReg ..... region Index to be process.
30 C ndId ... diagnostic Id number (in diagnostics long list)
31 C mate ..... counter diagnostic Id number if any ; 0 otherwise
32 C myThid ..... my thread Id number
33 _RL undef
34 INTEGER nLev, jReg, ndId, mate
35 INTEGER myThid
36
37 C !OUTPUT PARAMETERS:
38 C qtmp1 ..... AVERAGED DIAGNOSTIC QUANTITY
39 C qtmp2 ..... working array (used for counter mate statistics)
40 _RL qtmp1(0:nStats,0:nLev)
41 _RL qtmp2(0:nStats,0:nLev)
42 CEOP
43
44 C !LOCAL VARIABLES:
45 INTEGER im, ix, iv
46 PARAMETER ( iv = nStats - 2 , im = nStats - 1 , ix = nStats )
47 INTEGER bi, bj
48 INTEGER i, k, kd, kCnt, klev, kMlev
49
50 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
51
52 C-- Initialize to zero :
53 DO k=0,nLev
54 DO i=0,nStats
55 qtmp1(i,k) = 0.
56 qtmp2(i,k) = 0.
57 ENDDO
58 ENDDO
59
60 klev = kdiag(ndId)
61 IF ( mate.GT.0 ) kMlev = kdiag(mate)
62
63 IF (klev.LE.nLev) THEN
64 C--- Compute global statistics :
65
66 C-- Retrieve tile statistics first
67 DO bj=myByLo(myThid),myByHi(myThid)
68 DO bi=myBxLo(myThid),myBxHi(myThid)
69
70 DO k=1,klev
71 kd = iSdiag(ndId) + k - 1
72 IF ( qSdiag(0,jReg,kd,bi,bj).GT.0. ) THEN
73 IF ( qtmp1(0,k).LE.0. ) THEN
74 DO i=0,nStats
75 qtmp1(i,k) = qSdiag(i,jReg,kd,bi,bj)
76 ENDDO
77 ELSE
78 DO i=0,iv
79 qtmp1(i,k) = qtmp1(i,k) + qSdiag(i,jReg,kd,bi,bj)
80 ENDDO
81 qtmp1(im,k) = MIN( qtmp1(im,k),qSdiag(im,jReg,kd,bi,bj) )
82 qtmp1(ix,k) = MAX( qtmp1(ix,k),qSdiag(ix,jReg,kd,bi,bj) )
83 ENDIF
84 ENDIF
85 ENDDO
86 IF ( mate.GT.0 ) THEN
87 DO k=1,kMlev
88 kd = iSdiag(mate) + k - 1
89 IF ( qSdiag(0,jReg,kd,bi,bj).GT.0. ) THEN
90 IF ( qtmp1(0,k).LE.0. ) THEN
91 DO i=0,nStats
92 qtmp2(i,k) = qSdiag(i,jReg,kd,bi,bj)
93 ENDDO
94 ELSE
95 DO i=0,iv
96 qtmp2(i,k) = qtmp2(i,k) + qSdiag(i,jReg,kd,bi,bj)
97 ENDDO
98 qtmp2(im,k) = MIN( qtmp2(im,k),qSdiag(im,jReg,kd,bi,bj) )
99 qtmp2(ix,k) = MAX( qtmp2(ix,k),qSdiag(ix,jReg,kd,bi,bj) )
100 ENDIF
101 ENDIF
102 ENDDO
103 ENDIF
104
105 C- end tile index loops
106 ENDDO
107 ENDDO
108
109 C-- Global min,max & sum (at each level) over all thread & processors :
110 DO k=1,klev
111 DO i=0,iv
112 _GLOBAL_SUM_R8(qtmp1(i,k),myThid)
113 ENDDO
114 qtmp1(im,k) = -qtmp1(im,k)
115 _GLOBAL_MAX_R8(qtmp1(im,k),myThid)
116 qtmp1(im,k) = -qtmp1(im,k)
117 _GLOBAL_MAX_R8(qtmp1(ix,k),myThid)
118 ENDDO
119 IF ( mate.GT.0 ) THEN
120 DO k=1,kMlev
121 DO i=0,iv
122 _GLOBAL_SUM_R8(qtmp2(i,k),myThid)
123 ENDDO
124 qtmp2(im,k) = -qtmp2(im,k)
125 _GLOBAL_MAX_R8(qtmp2(im,k),myThid)
126 qtmp2(im,k) = -qtmp2(im,k)
127 _GLOBAL_MAX_R8(qtmp2(ix,k),myThid)
128 ENDDO
129 ENDIF
130
131 C-- Vertical integral, min & max :
132 DO k=1,klev
133 IF ( qtmp1(0,0).LE.0. ) THEN
134 DO i=0,nStats
135 qtmp1(i,0) = qtmp1(i,k)
136 ENDDO
137 ELSE
138 DO i=0,iv
139 qtmp1(i,0) = qtmp1(i,0) + qtmp1(i,k)
140 ENDDO
141 qtmp1(im,0) = MIN(qtmp1(im,0),qtmp1(im,k))
142 qtmp1(ix,0) = MAX(qtmp1(ix,0),qtmp1(ix,k))
143 ENDIF
144 ENDDO
145 IF ( mate.GT.0 ) THEN
146 DO k=1,kMlev
147 IF ( qtmp2(0,0).LE.0. ) THEN
148 DO i=0,nStats
149 qtmp2(i,0) = qtmp2(i,k)
150 ENDDO
151 ELSE
152 DO i=0,iv
153 qtmp2(i,0) = qtmp2(i,0) + qtmp2(i,k)
154 ENDDO
155 qtmp2(im,0) = MIN(qtmp2(im,0),qtmp2(im,k))
156 qtmp2(ix,0) = MAX(qtmp2(ix,0),qtmp2(ix,k))
157 ENDIF
158 ENDDO
159 ENDIF
160
161 C-- Average, Standard.Dev.:
162 C- no counter diagnostics => average = Sum / vol :
163 IF ( mate.EQ.0 ) THEN
164 DO k=0,klev
165 IF ( qtmp1(0,k).LE.0. ) THEN
166 DO i=1,nStats
167 qtmp1(i,k) = undef
168 ENDDO
169 ELSE
170 DO i=1,iv
171 qtmp1(i,k) = qtmp1(i,k) / qtmp1(0,k)
172 ENDDO
173 C Variance :
174 qtmp1(iv,k) = qtmp1(iv,k) - qtmp1(1,k)*qtmp1(1,k)
175 C Standard deviation :
176 IF (qtmp1(iv,k).GT.0.) qtmp1(iv,k) = SQRT(qtmp1(iv,k))
177 ENDIF
178 ENDDO
179 ELSE
180 C With counter diagnostics => average = Sum / Sum(counter) :
181 DO k=0,klev
182 kCnt = min(k,kMlev)
183 IF ( qtmp2(0,kCnt).LE.0. ) THEN
184 DO i=1,nStats
185 qtmp1(i,k) = undef
186 ENDDO
187 ELSEIF ( qtmp2(1,kCnt).LE.0. ) THEN
188 DO i=1,iv
189 qtmp1(i,k) = undef
190 ENDDO
191 ELSE
192 DO i=1,iv
193 qtmp1(i,k) = qtmp1(i,k) / qtmp2(1,kCnt)
194 ENDDO
195 C jmc: looks like there is a Pb with how Variance is computed
196 C Variance :
197 qtmp1(iv,k) = qtmp1(iv,k) - qtmp1(1,k)*qtmp1(1,k)
198 C Standard deviation :
199 IF (qtmp1(iv,k).GT.0.) qtmp1(iv,k) = SQRT(qtmp1(iv,k))
200 ENDIF
201 ENDDO
202 ENDIF
203
204 ENDIF
205
206 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
207
208 RETURN
209 END

  ViewVC Help
Powered by ViewVC 1.1.22