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

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

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


Revision 1.20 - (hide annotations) (download)
Thu Aug 25 21:54:54 2005 UTC (18 years, 8 months ago) by jmc
Branch: MAIN
Changes since 1.19: +15 -138 lines
do vertical interpolation in dedicated S/R DIAGNOSTICS_INTERP_VERT

1 jmc 1.20 C $Header: /u/gcmpack/MITgcm/pkg/diagnostics/diagnostics_out.F,v 1.19 2005/07/18 18:35:19 molod Exp $
2 jmc 1.1 C $Name: $
3    
4     #include "DIAG_OPTIONS.h"
5    
6     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
7     CBOP 0
8     C !ROUTINE: DIAGNOSTICS_OUT
9    
10     C !INTERFACE:
11     SUBROUTINE DIAGNOSTICS_OUT(
12 jmc 1.15 I listId,
13 jmc 1.1 I myIter,
14 edhill 1.14 I myTime,
15 jmc 1.1 I myThid )
16    
17     C !DESCRIPTION:
18     C Write output for diagnostics fields.
19 jmc 1.15
20 jmc 1.1 C !USES:
21 jmc 1.3 IMPLICIT NONE
22 jmc 1.1 #include "SIZE.h"
23     #include "EEPARAMS.h"
24     #include "PARAMS.h"
25 edhill 1.7 #include "GRID.h"
26 jmc 1.3 #include "DIAGNOSTICS_SIZE.h"
27     #include "DIAGNOSTICS.h"
28 jmc 1.1
29     #ifdef ALLOW_FIZHI
30     #include "fizhi_SIZE.h"
31     #else
32 jmc 1.3 INTEGER Nrphys
33     PARAMETER (Nrphys=0)
34 jmc 1.1 #endif
35    
36    
37     C !INPUT PARAMETERS:
38 jmc 1.15 C listId :: Diagnostics list number being written
39 jmc 1.3 C myIter :: current iteration number
40 jmc 1.15 C myTime :: current time of simulation (s)
41 jmc 1.3 C myThid :: my Thread Id number
42 edhill 1.14 _RL myTime
43 jmc 1.15 INTEGER listId, myIter, myThid
44 jmc 1.1 CEOP
45    
46 jmc 1.3 C !LOCAL VARIABLES:
47 jmc 1.15 C i,j,k :: loop indices
48     C md :: field number in the list "listId".
49     C ndId :: diagnostics Id number (in available diagnostics list)
50     C mate :: counter mate Id number (in available diagnostics list)
51     C ip :: diagnostics pointer to storage array
52     C im :: counter-mate pointer to storage array
53     INTEGER i, j, k
54     INTEGER bi, bj
55     INTEGER md, ndId, ip, im
56     INTEGER mate, mVec
57 jmc 1.3 CHARACTER*8 parms1
58     CHARACTER*3 mate_index
59 jmc 1.1 _RL qtmp1(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr+Nrphys,nSx,nSy)
60     _RL undef, getcon
61 jmc 1.3 EXTERNAL getcon
62     INTEGER ILNBLNK
63     EXTERNAL ILNBLNK
64     INTEGER ilen
65 jmc 1.20 INTEGER nlevsout
66 jmc 1.1
67 jmc 1.6 INTEGER ioUnit
68 jmc 1.11 CHARACTER*(MAX_LEN_FNAM) fn
69 jmc 1.1 CHARACTER*(MAX_LEN_MBUF) suff
70 jmc 1.3 CHARACTER*(MAX_LEN_MBUF) msgBuf
71     LOGICAL glf
72 jmc 1.1 #ifdef ALLOW_MNC
73 jmc 1.3 INTEGER ii
74 jmc 1.1 CHARACTER*(MAX_LEN_FNAM) diag_mnc_bn
75 jmc 1.3 INTEGER CW_DIMS, NLEN
76     PARAMETER ( CW_DIMS = 10 )
77     PARAMETER ( NLEN = 80 )
78     INTEGER dim(CW_DIMS), ib(CW_DIMS), ie(CW_DIMS)
79     CHARACTER*(NLEN) dn(CW_DIMS)
80 edhill 1.7 CHARACTER*(NLEN) d_cw_name
81 jmc 1.3 CHARACTER*(NLEN) dn_blnk
82 jmc 1.20 #ifdef DIAG_MNC_COORD_NEEDSWORK
83     CHARACTER*(5) ctmp
84 edhill 1.7 _RS ztmp(Nr+Nrphys)
85 jmc 1.20 #endif
86 jmc 1.1 #endif /* ALLOW_MNC */
87    
88 jmc 1.3 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
89    
90 jmc 1.6 ioUnit= standardMessageUnit
91 jmc 1.1 undef = getcon('UNDEF')
92     glf = globalFiles
93     WRITE(suff,'(I10.10)') myIter
94 jmc 1.15 ilen = ILNBLNK(fnames(listId))
95     WRITE( fn, '(A,A,A)' ) fnames(listId)(1:ilen),'.',suff(1:10)
96 jmc 1.1
97     #ifdef ALLOW_MNC
98     IF (useMNC .AND. diag_mnc) THEN
99     DO i = 1,MAX_LEN_FNAM
100     diag_mnc_bn(i:i) = ' '
101     ENDDO
102     DO i = 1,NLEN
103     dn_blnk(i:i) = ' '
104     ENDDO
105 jmc 1.15 WRITE( diag_mnc_bn, '(A)' ) fnames(listId)(1:ilen)
106 jmc 1.1
107     C Update the record dimension by writing the iteration number
108     CALL MNC_CW_SET_UDIM(diag_mnc_bn, -1, myThid)
109 edhill 1.14 CALL MNC_CW_RL_W_S('D',diag_mnc_bn,0,0,'T',myTime,myThid)
110 jmc 1.1 CALL MNC_CW_SET_UDIM(diag_mnc_bn, 0, myThid)
111    
112     dn(1)(1:NLEN) = dn_blnk(1:NLEN)
113 jmc 1.15 WRITE(dn(1),'(a,i6.6)') 'Zmd', nlevels(listId)
114     dim(1) = nlevels(listId)
115 jmc 1.1 ib(1) = 1
116 jmc 1.15 ie(1) = nlevels(listId)
117 jmc 1.1
118     CALL MNC_CW_ADD_GNAME('diag_levels', 1,
119     & dim, dn, ib, ie, myThid)
120 edhill 1.7 CALL MNC_CW_ADD_VNAME('diag_levels', 'diag_levels',
121 jmc 1.1 & 0,0, myThid)
122 edhill 1.7 CALL MNC_CW_ADD_VATTR_TEXT('diag_levels','description',
123     & 'Idicies of vertical levels within the source arrays',
124 jmc 1.1 & myThid)
125    
126 edhill 1.9 CALL MNC_CW_RL_W('D',diag_mnc_bn,0,0,
127 jmc 1.15 & 'diag_levels', levs(1,listId), myThid)
128 jmc 1.1
129 edhill 1.7 CALL MNC_CW_DEL_VNAME('diag_levels', myThid)
130 jmc 1.1 CALL MNC_CW_DEL_GNAME('diag_levels', myThid)
131 edhill 1.7
132 edhill 1.16 #ifdef DIAG_MNC_COORD_NEEDSWORK
133     C This part has been placed in an #ifdef because, as its currently
134     C written, it will only work with variables defined on a dynamics
135     C grid. As we start using diagnostics for physics grids, ice
136     C levels, land levels, etc. the different vertical coordinate
137     C dimensions will have to be taken into account.
138    
139 edhill 1.7 C Now define: Zmdxxxxxx, Zudxxxxxx, Zldxxxxxx
140     ctmp(1:5) = 'mul '
141     DO i = 1,3
142     dn(1)(1:NLEN) = dn_blnk(1:NLEN)
143 jmc 1.15 WRITE(dn(1),'(3a,i6.6)') 'Z',ctmp(i:i),'d',nlevels(listId)
144 edhill 1.7 CALL MNC_CW_ADD_GNAME(dn(1), 1, dim, dn, ib, ie, myThid)
145     CALL MNC_CW_ADD_VNAME(dn(1), dn(1), 0,0, myThid)
146 edhill 1.10
147     C The following three ztmp() loops should eventually be modified
148     C to reflect the fractional nature of levs(j,l) -- they should
149     C do something like:
150     C ztmp(j) = rC(INT(FLOOR(levs(j,l))))
151     C + ( rC(INT(FLOOR(levs(j,l))))
152     C + rC(INT(CEIL(levs(j,l)))) )
153     C / ( levs(j,l) - FLOOR(levs(j,l)) )
154     C for averaged levels.
155     IF (i .EQ. 1) THEN
156 jmc 1.15 DO j = 1,nlevels(listId)
157     ztmp(j) = rC(NINT(levs(j,listId)))
158 edhill 1.10 ENDDO
159     CALL MNC_CW_ADD_VATTR_TEXT(dn(1),'description',
160     & 'Dimensional coordinate value at the mid point',
161     & myThid)
162     ELSEIF (i .EQ. 2) THEN
163 jmc 1.15 DO j = 1,nlevels(listId)
164     ztmp(j) = rF(NINT(levs(j,listId)) + 1)
165 edhill 1.10 ENDDO
166     CALL MNC_CW_ADD_VATTR_TEXT(dn(1),'description',
167     & 'Dimensional coordinate value at the upper point',
168     & myThid)
169     ELSEIF (i .EQ. 3) THEN
170 jmc 1.15 DO j = 1,nlevels(listId)
171     ztmp(j) = rF(NINT(levs(j,listId)))
172 edhill 1.10 ENDDO
173     CALL MNC_CW_ADD_VATTR_TEXT(dn(1),'description',
174     & 'Dimensional coordinate value at the lower point',
175     & myThid)
176     ENDIF
177 edhill 1.7 CALL MNC_CW_RS_W('D',diag_mnc_bn,0,0, dn(1), ztmp, myThid)
178     CALL MNC_CW_DEL_VNAME(dn(1), myThid)
179     CALL MNC_CW_DEL_GNAME(dn(1), myThid)
180     ENDDO
181 edhill 1.16 #endif /* DIAG_MNC_COORD_NEEDSWORK */
182 edhill 1.7
183 jmc 1.1 ENDIF
184     #endif /* ALLOW_MNC */
185    
186 jmc 1.15 DO md = 1,nfields(listId)
187     ndId = jdiag(md,listId)
188     parms1 = gdiag(ndId)(1:8)
189     IF ( idiag(md,listId).NE.0 .AND. parms1(5:5).NE.'D' ) THEN
190 jmc 1.3 C-- Start processing 1 Fld :
191    
192 jmc 1.15 ip = ABS(idiag(md,listId))
193     im = mdiag(md,listId)
194     IF ( ndiag(ip,1,1).EQ.0 ) THEN
195 jmc 1.3 C- Empty diagnostics case :
196    
197     _BEGIN_MASTER( myThid )
198     WRITE(msgBuf,'(A,I10)')
199     & '- WARNING - from DIAGNOSTICS_OUT at iter=', myIter
200 jmc 1.15 CALL PRINT_MESSAGE( msgBuf, errorMessageUnit,
201 jmc 1.3 & SQUEEZE_RIGHT, myThid)
202     WRITE(msgBuf,'(A,I4,3A,I3,2A)')
203 jmc 1.15 & '- WARNING - diag.#',ndId, ' : ',flds(md,listId),
204     & ' (#',md,' ) in outp.Stream: ',fnames(listId)
205     CALL PRINT_MESSAGE( msgBuf, errorMessageUnit,
206 jmc 1.3 & SQUEEZE_RIGHT, myThid)
207     WRITE(msgBuf,'(A,I2,A)')
208 jmc 1.15 & '- WARNING - has not been filled (ndiag=',
209     & ndiag(ip,1,1), ' )'
210     CALL PRINT_MESSAGE( msgBuf, errorMessageUnit,
211 jmc 1.3 & SQUEEZE_RIGHT, myThid)
212     WRITE(msgBuf,'(A)')
213     & 'WARNING DIAGNOSTICS_OUT => write ZEROS instead'
214 jmc 1.15 CALL PRINT_MESSAGE( msgBuf, errorMessageUnit,
215 jmc 1.3 & SQUEEZE_RIGHT, myThid)
216     _END_MASTER( myThid )
217     DO bj = myByLo(myThid), myByHi(myThid)
218     DO bi = myBxLo(myThid), myBxHi(myThid)
219 jmc 1.15 DO k = 1,nlevels(listId)
220 jmc 1.3 DO j = 1-OLy,sNy+OLy
221     DO i = 1-OLx,sNx+OLx
222     qtmp1(i,j,k,bi,bj) = 0. _d 0
223     ENDDO
224     ENDDO
225     ENDDO
226     ENDDO
227     ENDDO
228    
229     ELSE
230     C- diagnostics is not empty :
231    
232 jmc 1.15 IF ( myThid.EQ.1 ) WRITE(ioUnit,'(A,I3,3A,I8,2A)')
233     & ' Computing Diagnostic # ', ndId, ' ', cdiag(ndId),
234     & ' Counter:',ndiag(ip,1,1),' Parms: ',gdiag(ndId)
235 jmc 1.3
236 jmc 1.6 IF ( parms1(5:5).EQ.'C' ) THEN
237     C Check for Mate of a Counter Diagnostic
238     C --------------------------------------
239     mate_index = parms1(6:8)
240     READ (mate_index,'(I3)') mate
241 jmc 1.15 IF ( myThid.EQ.1 ) WRITE(ioUnit,'(3A,I3,2A)')
242     & ' use Counter Mate for ', cdiag(ndId),
243     & ' Diagnostic # ',mate, ' ', cdiag(mate)
244    
245 jmc 1.6 ELSE
246     mate = 0
247 jmc 1.3
248     C Check for Mate of a Vector Diagnostic
249     C -------------------------------------
250     IF ( parms1(1:1).EQ.'U' .OR. parms1(1:1).EQ.'V' ) THEN
251     mate_index = parms1(6:8)
252 jmc 1.6 READ (mate_index,'(I3)') mVec
253 jmc 1.15 IF ( im.GT.0 .AND. ndiag(MAX(1,im),1,1).GT.0 ) THEN
254     IF ( myThid.EQ.1 ) WRITE(ioUnit,'(3A,I3,3A)')
255     & ' Vector Mate for ', cdiag(ndId),
256     & ' Diagnostic # ',mVec, ' ', cdiag(mVec),
257     & ' exists '
258 jmc 1.3 ELSE
259 jmc 1.15 IF ( myThid.EQ.1 ) WRITE(ioUnit,'(3A,I3,3A)')
260     & ' Vector Mate for ', cdiag(ndId),
261     & ' Diagnostic # ',mVec, ' ', cdiag(mVec),
262     & ' not enabled'
263 jmc 1.3 ENDIF
264     ENDIF
265 jmc 1.6 ENDIF
266 jmc 1.3
267 jmc 1.6 DO bj = myByLo(myThid), myByHi(myThid)
268     DO bi = myBxLo(myThid), myBxHi(myThid)
269 jmc 1.15 DO k = 1,nlevels(listId)
270 jmc 1.6 CALL GETDIAG(
271 jmc 1.15 I levs(k,listId),undef,
272 jmc 1.6 O qtmp1(1-OLx,1-OLy,k,bi,bj),
273 jmc 1.15 I ndId,mate,ip,im,bi,bj,myThid)
274 jmc 1.3 ENDDO
275 jmc 1.6 ENDDO
276     ENDDO
277 jmc 1.1
278 jmc 1.3 C- end of empty diag / not empty block
279     ENDIF
280 jmc 1.1
281 molod 1.17 nlevsout = nlevels(listId)
282    
283     C-----------------------------------------------------------------------
284 jmc 1.20 C Check to see if we need to interpolate before output
285 molod 1.17 C-----------------------------------------------------------------------
286 jmc 1.20 IF ( fflags(listId)(2:2).EQ.'P' ) THEN
287     C- Do vertical interpolation:
288     CALL DIAGNOSTICS_INTERP_VERT(
289     I listId, md, ndId, ip, im,
290     U nlevsout,
291     U qtmp1,
292     I undef,
293     I myTime, myIter, myThid )
294     ENDIF
295 molod 1.17
296 jmc 1.1 #ifdef ALLOW_MDSIO
297 jmc 1.3 C Prepare for mdsio optionality
298     IF (diag_mdsio) THEN
299 jmc 1.15 IF (fflags(listId)(1:1) .EQ. ' ') THEN
300 edhill 1.13 C This is the old default behavior
301 jmc 1.15 CALL MDSWRITEFIELD_NEW(fn,writeBinaryPrec,glf,'RL',
302 molod 1.17 & Nr+Nrphys,nlevsout,qtmp1,md,myIter,myThid)
303 jmc 1.15 ELSEIF (fflags(listId)(1:1) .EQ. 'R') THEN
304 edhill 1.13 C Force it to be 32-bit precision
305 jmc 1.15 CALL MDSWRITEFIELD_NEW(fn,precFloat32,glf,'RL',
306 molod 1.17 & Nr+Nrphys,nlevsout,qtmp1,md,myIter,myThid)
307 jmc 1.15 ELSEIF (fflags(listId)(1:1) .EQ. 'D') THEN
308 edhill 1.13 C Force it to be 64-bit precision
309 jmc 1.15 CALL MDSWRITEFIELD_NEW(fn,precFloat64,glf,'RL',
310 molod 1.17 & Nr+Nrphys,nlevsout,qtmp1,md,myIter,myThid)
311 edhill 1.13 ENDIF
312 jmc 1.3 ENDIF
313 jmc 1.1 #endif /* ALLOW_MDSIO */
314    
315     #ifdef ALLOW_MNC
316 jmc 1.3 IF (useMNC .AND. diag_mnc) THEN
317 jmc 1.1
318 jmc 1.3 _BEGIN_MASTER( myThid )
319 jmc 1.1
320 jmc 1.3 DO ii = 1,CW_DIMS
321 edhill 1.7 d_cw_name(1:NLEN) = dn_blnk(1:NLEN)
322 jmc 1.3 dn(ii)(1:NLEN) = dn_blnk(1:NLEN)
323     ENDDO
324    
325 edhill 1.7 C Note that the "d_cw_name" variable is a hack that hides a
326     C subtlety within MNC. Basically, each MNC-wrapped file is
327     C caching its own concept of what each "grid name" (that is, a
328     C dimension group name) means. So one cannot re-use the same
329     C "grid" name for different collections of dimensions within a
330 jmc 1.15 C given file. By appending the "ndId" values to each name, we
331 edhill 1.7 C guarantee uniqueness within each MNC-produced file.
332 jmc 1.15 WRITE(d_cw_name,'(a,i6.6)') 'd_cw_',ndId
333 edhill 1.7
334 edhill 1.5 C XY dimensions
335     dim(1) = sNx + 2*OLx
336     dim(2) = sNy + 2*OLy
337     ib(1) = OLx + 1
338     ib(2) = OLy + 1
339 jmc 1.15 IF (gdiag(ndId)(2:2) .EQ. 'M') THEN
340 edhill 1.5 dn(1)(1:2) = 'X'
341     ie(1) = OLx + sNx
342     dn(2)(1:2) = 'Y'
343     ie(2) = OLy + sNy
344 jmc 1.15 ELSEIF (gdiag(ndId)(2:2) .EQ. 'U') THEN
345 edhill 1.5 dn(1)(1:3) = 'Xp1'
346     ie(1) = OLx + sNx + 1
347     dn(2)(1:2) = 'Y'
348     ie(2) = OLy + sNy
349 jmc 1.15 ELSEIF (gdiag(ndId)(2:2) .EQ. 'V') THEN
350 edhill 1.5 dn(1)(1:2) = 'X'
351     ie(1) = OLx + sNx
352     dn(2)(1:3) = 'Yp1'
353     ie(2) = OLy + sNy + 1
354 jmc 1.15 ELSEIF (gdiag(ndId)(2:2) .EQ. 'Z') THEN
355 edhill 1.5 dn(1)(1:3) = 'Xp1'
356     ie(1) = OLx + sNx + 1
357     dn(2)(1:3) = 'Yp1'
358     ie(2) = OLy + sNy + 1
359     ENDIF
360    
361 jmc 1.3 C Z is special since it varies
362 jmc 1.15 WRITE(dn(3),'(a,i6.6)') 'Zd', nlevels(listId)
363     IF ( (gdiag(ndId)(10:10) .EQ. 'R')
364     & .AND. (gdiag(ndId)(9:9) .EQ. 'M') ) THEN
365     WRITE(dn(3),'(a,i6.6)') 'Zmd', nlevels(listId)
366 edhill 1.7 ENDIF
367 jmc 1.15 IF ( (gdiag(ndId)(10:10) .EQ. 'R')
368     & .AND. (gdiag(ndId)(9:9) .EQ. 'L') ) THEN
369     WRITE(dn(3),'(a,i6.6)') 'Zld', nlevels(listId)
370 edhill 1.7 ENDIF
371 jmc 1.15 IF ( (gdiag(ndId)(10:10) .EQ. 'R')
372     & .AND. (gdiag(ndId)(9:9) .EQ. 'U') ) THEN
373     WRITE(dn(3),'(a,i6.6)') 'Zud', nlevels(listId)
374 edhill 1.7 ENDIF
375 jmc 1.3 dim(3) = Nr+Nrphys
376     ib(3) = 1
377 jmc 1.15 ie(3) = nlevels(listId)
378 jmc 1.1
379 edhill 1.5 C Time dimension
380     dn(4)(1:1) = 'T'
381     dim(4) = -1
382     ib(4) = 1
383     ie(4) = 1
384    
385 edhill 1.7 CALL MNC_CW_ADD_GNAME(d_cw_name, 4,
386 jmc 1.1 & dim, dn, ib, ie, myThid)
387 jmc 1.15 CALL MNC_CW_ADD_VNAME(cdiag(ndId), d_cw_name,
388 jmc 1.1 & 4,5, myThid)
389 jmc 1.15 CALL MNC_CW_ADD_VATTR_TEXT(cdiag(ndId),'description',
390     & tdiag(ndId),myThid)
391     CALL MNC_CW_ADD_VATTR_TEXT(cdiag(ndId),'units',
392     & udiag(ndId),myThid)
393 jmc 1.1
394 jmc 1.15 IF ((fflags(listId)(1:1) .EQ. ' ')
395     & .OR. (fflags(listId)(1:1) .EQ. 'R')) THEN
396 edhill 1.13 CALL MNC_CW_RL_W('R',diag_mnc_bn,0,0,
397 jmc 1.15 & cdiag(ndId), qtmp1, myThid)
398     ELSEIF (fflags(listId)(1:1) .EQ. 'D') THEN
399 edhill 1.13 CALL MNC_CW_RL_W('D',diag_mnc_bn,0,0,
400 jmc 1.15 & cdiag(ndId), qtmp1, myThid)
401 edhill 1.13 ENDIF
402    
403 jmc 1.15 CALL MNC_CW_DEL_VNAME(cdiag(ndId), myThid)
404 edhill 1.7 CALL MNC_CW_DEL_GNAME(d_cw_name, myThid)
405 jmc 1.1
406 jmc 1.3 _END_MASTER( myThid )
407 jmc 1.1
408 jmc 1.3 ENDIF
409 jmc 1.1 #endif /* ALLOW_MNC */
410    
411 jmc 1.15 C-- end of Processing Fld # md
412 jmc 1.3 ENDIF
413     ENDDO
414 jmc 1.1
415 jmc 1.15 RETURN
416 jmc 1.3 END
417 jmc 1.15
418 jmc 1.1 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|

  ViewVC Help
Powered by ViewVC 1.1.22