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

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

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


Revision 1.28 - (hide annotations) (download)
Sun Jan 25 17:00:20 2009 UTC (15 years, 3 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint62, checkpoint61n, checkpoint61o, checkpoint61l, checkpoint61m, checkpoint61j, checkpoint61k, checkpoint61h, checkpoint61i, checkpoint61v, checkpoint61w, checkpoint61t, checkpoint61u, checkpoint61r, checkpoint61s, checkpoint61p, checkpoint61q, checkpoint61z, checkpoint61x, checkpoint61y
Changes since 1.27: +13 -5 lines
fix S/R DIAGNOSTICS_COUNT for "periodic averaging diagnostics".

1 jmc 1.28 C $Header: /u/gcmpack/MITgcm/pkg/diagnostics/diagnostics_utils.F,v 1.27 2008/11/18 21:41:06 jmc Exp $
2 edhill 1.8 C $Name: $
3    
4 edhill 1.9 #include "DIAG_OPTIONS.h"
5    
6 jmc 1.25 C-- File diagnostics_utils.F: General purpose support routines
7     C-- Contents:
8     C-- o GETDIAG
9     C-- o DIAGNOSTICS_COUNT
10 jmc 1.26 C-- o DIAGNOSTICS_GET_POINTERS
11     C-- o DIAGS_GET_PARMS_I (Function)
12 jmc 1.25 C-- o DIAGS_MK_UNITS (Function)
13     C-- o DIAGS_MK_TITLE (Function)
14    
15 edhill 1.12 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
16     CBOP 0
17     C !ROUTINE: GETDIAG
18    
19     C !INTERFACE:
20 jmc 1.18 SUBROUTINE GETDIAG(
21     I levreal, undef,
22     O qtmp,
23 jmc 1.21 I ndId, mate, ip, im, bi, bj, myThid )
24 edhill 1.12
25     C !DESCRIPTION:
26 edhill 1.13 C Retrieve averaged model diagnostic
27 jmc 1.15
28 edhill 1.12 C !USES:
29 jmc 1.18 IMPLICIT NONE
30 molod 1.3 #include "EEPARAMS.h"
31 molod 1.1 #include "SIZE.h"
32 jmc 1.15 #include "DIAGNOSTICS_SIZE.h"
33     #include "DIAGNOSTICS.h"
34 molod 1.3
35 jmc 1.18 C !INPUT PARAMETERS:
36 jmc 1.21 C levreal :: Diagnostic LEVEL
37     C undef :: UNDEFINED VALUE
38     C ndId :: DIAGNOSTIC NUMBER FROM MENU
39     C mate :: counter DIAGNOSTIC NUMBER if any ; 0 otherwise
40     C ip :: pointer to storage array location for diag.
41     C im :: pointer to storage array location for mate
42     C bi :: X-direction tile number
43     C bj :: Y-direction tile number
44     C myThid :: my thread Id number
45 jmc 1.15 _RL levreal
46 molod 1.3 _RL undef
47 jmc 1.21 INTEGER ndId, mate, ip, im
48 jmc 1.18 INTEGER bi,bj, myThid
49 jmc 1.15
50 jmc 1.18 C !OUTPUT PARAMETERS:
51 edhill 1.12 C qtmp ..... AVERAGED DIAGNOSTIC QUANTITY
52 jmc 1.18 _RL qtmp(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
53     CEOP
54 molod 1.11
55 jmc 1.18 C !LOCAL VARIABLES:
56 molod 1.3 _RL factor
57 jmc 1.18 INTEGER i, j, ipnt,ipCt
58     INTEGER lev, levCt, klev
59 molod 1.11
60 jmc 1.21 IF (ndId.GE.1) THEN
61 jmc 1.15 lev = NINT(levreal)
62 jmc 1.21 klev = kdiag(ndId)
63 jmc 1.18 IF (lev.LE.klev) THEN
64 jmc 1.15
65 jmc 1.18 IF ( mate.EQ.0 ) THEN
66     C- No counter diagnostics => average = Sum / ndiag :
67 jmc 1.15
68 jmc 1.21 ipnt = ip + lev - 1
69     factor = FLOAT(ndiag(ip,bi,bj))
70     IF (ndiag(ip,bi,bj).NE.0) factor = 1. _d 0 / factor
71 jmc 1.18
72 jmc 1.27 #ifdef ALLOW_FIZHI
73 jmc 1.18 DO j = 1,sNy+1
74     DO i = 1,sNx+1
75     IF ( qdiag(i,j,ipnt,bi,bj) .LE. undef ) THEN
76     qtmp(i,j) = qdiag(i,j,ipnt,bi,bj)*factor
77     ELSE
78     qtmp(i,j) = undef
79     ENDIF
80     ENDDO
81     ENDDO
82 jmc 1.27 #else /* ALLOW_FIZHI */
83     DO j = 1,sNy+1
84     DO i = 1,sNx+1
85     qtmp(i,j) = qdiag(i,j,ipnt,bi,bj)*factor
86     ENDDO
87     ENDDO
88     #endif /* ALLOW_FIZHI */
89 jmc 1.18
90     ELSE
91     C- With counter diagnostics => average = Sum / counter:
92    
93 jmc 1.21 ipnt = ip + lev - 1
94 jmc 1.18 levCt= MIN(lev,kdiag(mate))
95 jmc 1.21 ipCt = im + levCt - 1
96 jmc 1.18 DO j = 1,sNy+1
97     DO i = 1,sNx+1
98     IF ( qdiag(i,j,ipCt,bi,bj) .NE. 0. ) THEN
99     qtmp(i,j) = qdiag(i,j,ipnt,bi,bj)
100     & / qdiag(i,j,ipCt,bi,bj)
101     ELSE
102     qtmp(i,j) = undef
103     ENDIF
104     ENDDO
105     ENDDO
106 molod 1.1
107 jmc 1.18 ENDIF
108     ENDIF
109     ENDIF
110 molod 1.3
111 jmc 1.15 RETURN
112     END
113 edhill 1.12
114     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
115    
116 jmc 1.15 CBOP 0
117 jmc 1.19 C !ROUTINE: DIAGNOSTICS_COUNT
118     C !INTERFACE:
119 jmc 1.21 SUBROUTINE DIAGNOSTICS_COUNT (chardiag,
120     I biArg, bjArg, myThid)
121 jmc 1.19
122     C !DESCRIPTION:
123     C***********************************************************************
124     C routine to increment the diagnostic counter only
125     C***********************************************************************
126     C !USES:
127     IMPLICIT NONE
128    
129     C == Global variables ===
130     #include "EEPARAMS.h"
131     #include "SIZE.h"
132     #include "DIAGNOSTICS_SIZE.h"
133     #include "DIAGNOSTICS.h"
134    
135     C !INPUT PARAMETERS:
136     C***********************************************************************
137     C Arguments Description
138     C ----------------------
139     C chardiag :: Character expression for diag to increment the counter
140     C biArg :: X-direction tile number, or 0 if called outside bi,bj loops
141     C bjArg :: Y-direction tile number, or 0 if called outside bi,bj loops
142     C myThid :: my thread Id number
143     C***********************************************************************
144     CHARACTER*8 chardiag
145     INTEGER biArg, bjArg
146     INTEGER myThid
147     CEOP
148    
149     C !LOCAL VARIABLES:
150     C ===============
151 jmc 1.21 INTEGER m, n
152     INTEGER bi, bj
153 jmc 1.28 INTEGER ipt, ndId
154 jmc 1.19 c CHARACTER*(MAX_LEN_MBUF) msgBuf
155    
156 jmc 1.28 IF ( biArg.EQ.0 .AND. bjArg.EQ.0 ) THEN
157     bi = myBxLo(myThid)
158     bj = myByLo(myThid)
159     ELSE
160     bi = MIN(biArg,nSx)
161     bj = MIN(bjArg,nSy)
162     ENDIF
163    
164 jmc 1.21 C-- Run through list of active diagnostics to find which counter
165     C to increment (needs to be a valid & active diagnostic-counter)
166 jmc 1.19 DO n=1,nlists
167     DO m=1,nActive(n)
168 jmc 1.21 IF ( chardiag.EQ.flds(m,n) .AND. idiag(m,n).GT.0 ) THEN
169     ipt = idiag(m,n)
170 jmc 1.28 IF (ndiag(ipt,bi,bj).GE.0) THEN
171     ndId = jdiag(m,n)
172     ipt = ipt + pdiag(n,bi,bj)*kdiag(ndId)
173 jmc 1.21 C- Increment the counter for the diagnostic
174     IF ( biArg.EQ.0 .AND. bjArg.EQ.0 ) THEN
175     DO bj=myByLo(myThid), myByHi(myThid)
176     DO bi=myBxLo(myThid), myBxHi(myThid)
177     ndiag(ipt,bi,bj) = ndiag(ipt,bi,bj) + 1
178     ENDDO
179     ENDDO
180     ELSE
181     ndiag(ipt,bi,bj) = ndiag(ipt,bi,bj) + 1
182     ENDIF
183     C- Increment is done
184     ENDIF
185 jmc 1.19 ENDIF
186     ENDDO
187     ENDDO
188    
189 jmc 1.21 RETURN
190 jmc 1.19 END
191    
192     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
193    
194     CBOP 0
195 jmc 1.26 C !ROUTINE: DIAGNOSTICS_GET_POINTERS
196     C !INTERFACE:
197     SUBROUTINE DIAGNOSTICS_GET_POINTERS(
198     I diagName, listId,
199     O ndId, ip,
200     I myThid )
201    
202     C !DESCRIPTION:
203     C *================================================================*
204     C | o Returns the diagnostic Id number and diagnostic
205     C | pointer to storage array for a specified diagnostic.
206     C *================================================================*
207     C | Note: A diagnostics field can be stored multiple times
208     C | (for different output frequency,phase, ...).
209     C | operates in 2 ways:
210     C | o listId =0 => find 1 diagnostics Id & pointer which name matches.
211     C | o listId >0 => find the unique diagnostic Id & pointer with
212     C | the right name and same output time as "listId" output-list
213     C | o return ip=0 if did not find the right diagnostic;
214     C | (ndId <>0 if diagnostic exist but output time does not match)
215     C *================================================================*
216    
217     C !USES:
218     IMPLICIT NONE
219     #include "EEPARAMS.h"
220     #include "SIZE.h"
221     #include "DIAGNOSTICS_SIZE.h"
222     #include "DIAGNOSTICS.h"
223    
224     C !INPUT PARAMETERS:
225     C diagName :: diagnostic identificator name (8 characters long)
226     C listId :: list number that specify the output frequency
227     C myThid :: my Thread Id number
228     C !OUTPUT PARAMETERS:
229     C ndId :: diagnostics Id number (in available diagnostics list)
230     C ip :: diagnostics pointer to storage array
231    
232    
233     CHARACTER*8 diagName
234     INTEGER listId
235     INTEGER ndId, ip
236     INTEGER myThid
237     CEOP
238    
239     C !LOCAL VARIABLES:
240     INTEGER n,m
241    
242     ip = 0
243     ndId = 0
244    
245     IF ( listId.LE.0 ) THEN
246     C-- select the 1rst one which name matches:
247    
248     C- search for this diag. in the active 2D/3D diagnostics list
249     DO n=1,nlists
250     DO m=1,nActive(n)
251     IF ( ip.EQ.0 .AND. diagName.EQ.flds(m,n)
252     & .AND. idiag(m,n).NE.0 ) THEN
253     ip = ABS(idiag(m,n))
254     ndId = jdiag(m,n)
255     ENDIF
256     ENDDO
257     ENDDO
258    
259     ELSEIF ( listId.LE.nlists ) THEN
260     C-- select the unique diagnostic with output-time identical to listId
261    
262     C- search for this diag. in the active 2D/3D diagnostics list
263     DO n=1,nlists
264     IF ( ip.EQ.0
265     & .AND. freq(n) .EQ. freq(listId)
266     & .AND. phase(n).EQ.phase(listId)
267     & .AND. averageFreq(n) .EQ.averageFreq(listId)
268     & .AND. averagePhase(n).EQ.averagePhase(listId)
269     & .AND. averageCycle(n).EQ.averageCycle(listId)
270     & ) THEN
271     DO m=1,nActive(n)
272     IF ( ip.EQ.0 .AND. diagName.EQ.flds(m,n)
273     & .AND. idiag(m,n).NE.0 ) THEN
274     ip = ABS(idiag(m,n))
275     ndId = jdiag(m,n)
276     ENDIF
277     ENDDO
278     ELSEIF ( ip.EQ.0 ) THEN
279     DO m=1,nActive(n)
280     IF ( ip.EQ.0 .AND. diagName.EQ.flds(m,n)
281     & .AND. idiag(m,n).NE.0 ) THEN
282     ndId = jdiag(m,n)
283     ENDIF
284     ENDDO
285     ENDIF
286     ENDDO
287    
288     ELSE
289     STOP 'DIAGNOSTICS_GET_POINTERS: invalid listId number'
290     ENDIF
291    
292     RETURN
293     END
294    
295     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
296    
297     CBOP 0
298     C !ROUTINE: DIAGS_GET_PARMS_I
299    
300     C !INTERFACE:
301     INTEGER FUNCTION DIAGS_GET_PARMS_I(
302     I parName, myThid )
303    
304     C !DESCRIPTION:
305     C *==========================================================*
306     C | FUNCTION DIAGS_GET_PARMS_I
307     C | o Return the value of integer parameter
308     C | from one of the DIAGNOSTICS.h common blocs
309     C *==========================================================*
310    
311     C !USES:
312     IMPLICIT NONE
313     #include "EEPARAMS.h"
314     #include "SIZE.h"
315     #include "DIAGNOSTICS_SIZE.h"
316     #include "DIAGNOSTICS.h"
317    
318     C !INPUT PARAMETERS:
319     C parName :: string used to identify which parameter to get
320     C myThid :: my Thread Id number
321     CHARACTER*(*) parName
322     INTEGER myThid
323     CEOP
324    
325     C !LOCAL VARIABLES:
326     CHARACTER*(MAX_LEN_MBUF) msgBuf
327     INTEGER n
328    
329     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
330    
331     n = LEN(parName)
332     c write(0,'(3A,I4)')
333     c & 'DIAGS_GET_PARMS_I: parName="',parName,'" , length=',n
334    
335     IF ( parName.EQ.'LAST_DIAG_ID' ) THEN
336     DIAGS_GET_PARMS_I = ndiagt
337     ELSE
338     WRITE(msgBuf,'(4A)') 'DIAGS_GET_PARMS_I: ',
339     & ' parName="', parName, '" not known.'
340     CALL PRINT_ERROR( msgBuf, myThid )
341     STOP 'ABNORMAL END: S/R DIAGS_GET_PARMS_I'
342     ENDIF
343    
344     RETURN
345     END
346    
347     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
348    
349     CBOP 0
350 jmc 1.17 C !ROUTINE: DIAGS_MK_UNITS
351    
352     C !INTERFACE:
353 jmc 1.21 CHARACTER*16 FUNCTION DIAGS_MK_UNITS(
354 jmc 1.17 I diagUnitsInPieces, myThid )
355    
356     C !DESCRIPTION:
357     C *==========================================================*
358     C | FUNCTION DIAGS_MK_UNITS
359 jmc 1.21 C | o Return the diagnostic units string (16c) removing
360 jmc 1.17 C | blanks from the input string
361     C *==========================================================*
362    
363     C !USES:
364     IMPLICIT NONE
365     #include "EEPARAMS.h"
366    
367     C !INPUT PARAMETERS:
368 jmc 1.21 C diagUnitsInPieces :: string for diagnostic units: in several
369 jmc 1.17 C pieces, with blanks in between
370     C myThid :: my thread Id number
371     CHARACTER*(*) diagUnitsInPieces
372     INTEGER myThid
373     CEOP
374    
375     C !LOCAL VARIABLES:
376     CHARACTER*(MAX_LEN_MBUF) msgBuf
377     INTEGER i,j,n
378    
379 jmc 1.21 DIAGS_MK_UNITS = ' '
380 jmc 1.17 n = LEN(diagUnitsInPieces)
381 jmc 1.21
382 jmc 1.17 j = 0
383     DO i=1,n
384     IF (diagUnitsInPieces(i:i) .NE. ' ' ) THEN
385     j = j+1
386     IF ( j.LE.16 ) DIAGS_MK_UNITS(j:j) = diagUnitsInPieces(i:i)
387     ENDIF
388     ENDDO
389    
390     IF ( j.GT.16 ) THEN
391     WRITE(msgBuf,'(2A,I4,A)') '**WARNING** ',
392     & 'DIAGS_MK_UNITS: too long (',j,' >16) input string'
393     CALL PRINT_MESSAGE( msgBuf, errorMessageUnit,
394     & SQUEEZE_RIGHT , myThid)
395     WRITE(msgBuf,'(3A)') '**WARNING** ',
396     & 'DIAGS_MK_UNITS: input=', diagUnitsInPieces
397     CALL PRINT_MESSAGE( msgBuf, errorMessageUnit,
398     & SQUEEZE_RIGHT , myThid)
399     ENDIF
400    
401     RETURN
402     END
403 jmc 1.23
404     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
405    
406     CBOP 0
407     C !ROUTINE: DIAGS_MK_TITLE
408    
409     C !INTERFACE:
410     CHARACTER*80 FUNCTION DIAGS_MK_TITLE(
411     I diagTitleInPieces, myThid )
412    
413     C !DESCRIPTION:
414     C *==========================================================*
415     C | FUNCTION DIAGS_MK_TITLE
416     C | o Return the diagnostic title string (80c) removing
417     C | consecutive blanks from the input string
418     C *==========================================================*
419    
420     C !USES:
421     IMPLICIT NONE
422     #include "EEPARAMS.h"
423    
424     C !INPUT PARAMETERS:
425     C diagTitleInPieces :: string for diagnostic units: in several
426     C pieces, with blanks in between
427     C myThid :: my Thread Id number
428     CHARACTER*(*) diagTitleInPieces
429     INTEGER myThid
430     CEOP
431    
432     C !LOCAL VARIABLES:
433     CHARACTER*(MAX_LEN_MBUF) msgBuf
434     LOGICAL flag
435     INTEGER i,j,n
436    
437 molod 1.22 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
438 jmc 1.23
439     DIAGS_MK_TITLE = ' '
440     & //' '
441     n = LEN(diagTitleInPieces)
442    
443     j = 0
444     flag = .FALSE.
445     DO i=1,n
446     IF (diagTitleInPieces(i:i) .NE. ' ' ) THEN
447     IF ( flag ) THEN
448     j = j+1
449     IF (j.LE.80) DIAGS_MK_TITLE(j:j) = ' '
450     ENDIF
451     j = j+1
452     IF ( j.LE.80 ) DIAGS_MK_TITLE(j:j) = diagTitleInPieces(i:i)
453     flag = .FALSE.
454     ELSE
455     flag = j.GE.1
456     ENDIF
457     ENDDO
458    
459     IF ( j.GT.80 ) THEN
460     WRITE(msgBuf,'(2A,I4,A)') '**WARNING** ',
461     & 'DIAGS_MK_TITLE: too long (',j,' >80) input string'
462     CALL PRINT_MESSAGE( msgBuf, errorMessageUnit,
463     & SQUEEZE_RIGHT , myThid)
464     WRITE(msgBuf,'(3A)') '**WARNING** ',
465     & 'DIAGS_MK_TITLE: input=', diagTitleInPieces
466     CALL PRINT_MESSAGE( msgBuf, errorMessageUnit,
467     & SQUEEZE_RIGHT , myThid)
468     ENDIF
469    
470     RETURN
471     END

  ViewVC Help
Powered by ViewVC 1.1.22