source: trunk/NEMO/OPA_SRC/DIA/diawri_dimg.h90 @ 1057

Last change on this file since 1057 was 1057, checked in by rblod, 13 years ago

ctlopn with dimg diagnostics and additonnal stuff related to dimg

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 15.2 KB
Line 
1  !!----------------------------------------------------------------------
2  !!                        ***  diawri_dimg.h90  ***
3  !!----------------------------------------------------------------------
4  !!   OPA 9.0 , LOCEAN-IPSL (2005)
5  !! $Id$
6  !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt
7  !!----------------------------------------------------------------------
8
9  SUBROUTINE dia_wri (kt, kindic)
10    !!----------------------------------------------------------------------
11    !!           *** routine dia_wri ***
12    !!
13    !! ** Purpose : output dynamics and tracer fields on direct access file
14    !!              suitable for MPP computing
15    !!
16    !! ** define key : 'key_dimgout'
17    !!
18    !! **  Method : Default is to cumulate the values over the interval between
19    !!      2 output, and each nwrite time-steps  the mean value is  computed
20    !!      and written to the direct access file.
21    !!     If 'key_diainstant' is defined, no mean values are computed and the
22    !!     instantaneous fields are dump.
23    !!     kindic is 0 or >0 in normal condition. When < 0 it indicates an error
24    !!     condition and instantaneous file output is forced.
25    !!       Each processor creates its own file with its local data
26    !!     Merging all the files is performed off line by a dedicated program
27    !!
28    !! ** Arguments :
29    !!     kt      : time-step number
30    !!     kindinc :  error condition indicator : >=0 :  OK, < 0 : error.
31    !!
32    !! ** Naming convention for files
33    !!
34    !! {cexper}_{var}_y----m--d--.dimg
35    !!   cexper is the name of the experience, given in the namelist
36    !!   var can be either U, V, T, S, KZ, SSH, ...
37    !!   var can also be 2D, which means that each level of the file is a 2D field as described below
38    !!    y----m--d--  is the date at the time of the dump
39    !!    For mpp output, each processor dumps its own memory, on appropriate record range
40    !!    (direct access : for level jk of a klev field on proc narea irec = 1+ klev*(narea -1) + jk )
41    !!    To be tested with a lot of procs !!!!
42    !!
43    !!  level 1:  utau(:,:) * umask(:,:,1) zonal stress in N.m-2
44    !!  level 2:  vtau(:,:) * vmask(:,:,1) meridional stress in N. m-2
45    !!  level 3:  qsr + qns                total heat flux (W/m2)
46    !!  level 4:  emp (:,:)               E-P flux (mm/day)
47    !!  level 5:  tb  (:,:,1)-sst          model SST -forcing sst (degree C)
48    !!  level 6:  bsfb(:,:)         streamfunction (m**3/s)
49    !!  level 7:  qsr (:,:)         solar flux (W/m2)
50    !!  level 8:  qrp (:,:)                relax component of T flux.
51    !!  level 9:  erp (:,:)                relax component of S flux
52    !!  level 10: hmld(:,:)                turbocline depth
53    !!  level 11: hmlp(:,:)                mixed layer depth
54    !!  level 12: fr_i(:,:)                ice fraction (between 0 and 1)
55    !!  level 13: sst(:,:)                 the observed SST we relax to.
56    !!  level 14: qct(:,:)                 equivalent flux due to treshold SST
57    !!  level 15: fbt(:,:)                 feedback term .
58    !!  level 16: emps(:,:)                concentration/dilution water flux
59    !!  level 17: fsalt(:,:)               Ice=>ocean net freshwater
60    !!  level 18: gps(:,:)                 the surface pressure (m).
61    !!  level 19: spgu(:,:)                the surface pressure gradient in X direction.
62    !!  level 20: spgv(:,:)                the surface pressure gradient in Y direction.
63    !!
64    !! History
65    !!      original  : 91-03 ()
66    !!      additions : 91-11 (G. Madec)
67    !!      additions : 92-06 (M. Imbard) correction restart file
68    !!      additions : 92-07 (M. Imbard) split into diawri and rstwri
69    !!      additions : 93-03 (M. Imbard) suppress writibm
70    !!      additions : 94-12 (M. Imbard) acces direct files
71    !!      additions : 97-2002 ( Clipper Group ) dimg files
72    !!                  dec 2003 ( J.M. Molines) f90, mpp output for OPA9.0
73    !!   9.0  !  05-05  (S. Theetten) add emps fsalt move gps spgu spgv 2 lines below
74    !!   9.0  !  05-11  (V. Garnier) Surface pressure gradient organization
75    !!----------------------------------------------------------------------
76    !! * modules used
77    USE lib_mpp
78
79    !! * Arguments
80    INTEGER ,INTENT(in) :: kt, kindic
81
82    !! * local declarations
83    INTEGER :: inbsel, jk
84!!  INTEGER :: iwrite
85    INTEGER :: iyear,imon,iday
86    INTEGER, SAVE :: nmoyct
87
88#if defined key_diainstant
89    LOGICAL, PARAMETER :: ll_dia_inst=.TRUE.  !: for instantaneous output
90#else
91    LOGICAL, PARAMETER :: ll_dia_inst=.FALSE. !: for average output
92#endif
93
94    REAL(wp), SAVE, DIMENSION (jpi,jpj,jpk) ::  um , vm   ! used to compute mean u, v fields
95    REAL(wp), SAVE, DIMENSION (jpi,jpj,jpk) ::  wm        ! used to compute mean w fields
96    REAL(wp), SAVE, DIMENSION (jpi,jpj,jpk) ::  avtm      ! used to compute mean kz fields
97    REAL(wp), SAVE, DIMENSION (jpi,jpj,jpk) ::  tm , sm   ! used to compute mean t, s fields
98    REAL(wp), SAVE, DIMENSION (jpi,jpj,jpk) ::  fsel      ! used to compute mean 2d fields
99    REAL(wp) :: zdtj
100    !
101    CHARACTER(LEN=80) :: clname
102    CHARACTER(LEN=80) :: cltext
103    CHARACTER(LEN=80) :: clmode
104    CHARACTER(LEN= 4) :: clver
105    !
106    !  Initialization
107    !  ---------------
108    !
109#ifdef key_diaspr
110    inbsel = 20
111#else
112    inbsel = 17
113#endif
114#if defined key_flx_core
115    inbsel = 23
116#endif
117
118    IF( inbsel >  jpk) THEN
119       IF( lwp) WRITE(numout,*)  &
120            ' STOP inbsel =',inbsel,' is larger than jpk=',jpk
121       STOP
122    ENDIF
123
124
125    iyear = ndastp/10000
126    imon = (ndastp-iyear*10000)/100
127    iday = ndastp - imon*100 - iyear*10000
128
129    !     
130    !! dimg format V1.0 should start with the 4 char. string '@!01'
131    !!
132    clver='@!01'
133    !
134    IF( .NOT. ll_dia_inst ) THEN
135       !
136       !! * Mean output section
137       !! ----------------------
138       !
139       IF( kt == nit000 .AND. lwp ) WRITE(numout,*) &
140            'THE OUTPUT FILES CONTAINS THE AVERAGE OF EACH FIELD'
141       !
142       IF( kt == nit000  ) THEN
143          ! reset arrays for average computation
144          nmoyct = 0
145          !
146          um(:,:,:) = 0._wp
147          vm(:,:,:) = 0._wp
148          wm(:,:,:) = 0._wp
149          avtm(:,:,:) = 0._wp
150          tm(:,:,:) = 0._wp
151          sm(:,:,:) = 0._wp
152          fsel(:,:,:) = 0._wp
153          !
154       ENDIF
155
156       !  cumulate values
157       !  ---------------
158
159       nmoyct = nmoyct+1
160       !
161       um(:,:,:)=um(:,:,:) + un (:,:,:)
162       vm(:,:,:)=vm(:,:,:) + vn (:,:,:)
163       wm(:,:,:)=wm(:,:,:) + wn (:,:,:)
164       avtm(:,:,:)=avtm(:,:,:) + avt (:,:,:)
165       tm(:,:,:)=tm(:,:,:) + tn (:,:,:)
166       sm(:,:,:)=sm(:,:,:) + sn (:,:,:)
167       !
168       fsel(:,:,1 ) = fsel(:,:,1 ) + utau(:,:) * umask(:,:,1)
169       fsel(:,:,2 ) = fsel(:,:,2 ) + vtau(:,:) * vmask(:,:,1)
170       fsel(:,:,3 ) = fsel(:,:,3 ) + qsr (:,:) + qns  (:,:)
171       fsel(:,:,4 ) = fsel(:,:,4 ) + emp (:,:)
172       fsel(:,:,5 ) = fsel(:,:,5 ) + tb  (:,:,1) - sf_sst(1)%fnow(:,:)
173#if ! defined key_dynspg_rl
174       fsel(:,:,6 ) = fsel(:,:,6 ) + sshn(:,:)    ! SSH
175#else
176       fsel(:,:,6 ) = fsel(:,:,6 ) + bsfn(:,:)    ! BSF
177#endif
178       fsel(:,:,7 ) = fsel(:,:,7 ) + qsr(:,:)
179       fsel(:,:,8 ) = fsel(:,:,8 ) + qrp (:,:)
180       fsel(:,:,9 ) = fsel(:,:,9 ) + erp (:,:)
181       fsel(:,:,10) = fsel(:,:,10) + hmld(:,:)
182       fsel(:,:,11) = fsel(:,:,11) + hmlp(:,:)
183       fsel(:,:,12) = fsel(:,:,12) + fr_i(:,:)
184       fsel(:,:,13) = fsel(:,:,13) + sf_sst(1)%fnow(:,:)
185       !        fsel(:,:,14) = fsel(:,:,14) + qct(:,:)
186       !        fsel(:,:,15) = fsel(:,:,15) + fbt(:,:)
187       fsel(:,:,16) = fsel(:,:,16) + emps(:,:)
188#ifdef key_diaspr   
189       fsel(:,:,18) = fsel(:,:,18) + gps(:,:)/g
190#endif
191       !
192       ! Output of dynamics and tracer fields and selected fields (numwri)
193       ! -----------------------------------------------------------------
194       !
195       !
196       zdtj=rdt/86400.   ! time step in days
197       WRITE(clmode,'(f5.1,a)' ) nwrite*zdtj,' days average'
198
199       !       iwrite=NINT(adatrj/rwrite)
200       !      IF (abs(adatrj-iwrite*rwrite) < zdtj/2.      &
201
202       IF(  ( MOD (kt-nit000+1,nwrite) ==  0 )          &
203            &   .OR.       kindic <   0            &
204            &   .OR. ( kt == 1 .AND. kindic > 0)  ) THEN
205          ! it is time to make a dump on file
206          ! compute average
207          um(:,:,:) = um(:,:,:) / nmoyct
208          vm(:,:,:) = vm(:,:,:) / nmoyct
209          wm(:,:,:) = wm(:,:,:) / nmoyct
210          avtm(:,:,:) = avtm(:,:,:) / nmoyct
211          tm(:,:,:) = tm(:,:,:) / nmoyct
212          sm(:,:,:) = sm(:,:,:) / nmoyct
213          !
214          fsel(:,:,:) = fsel(:,:,:) / nmoyct
215          !
216          ! note : the surface pressure is not averaged, but rather
217          ! computed from the averaged gradients.
218          !
219#ifdef key_diaspr
220          fsel(:,:,18)= gps(:,:)/g
221          fsel(:,:,19)= spgu(:,:)
222          fsel(:,:,20)= spgv(:,:)
223#endif
224          ! mask mean field with tmask except utau vtau (1,2)
225          DO jk=3,inbsel
226            fsel(:,:,jk)=fsel(:,:,jk)*tmask(:,:,1)
227          END DO
228       ENDIF
229       !
230    ELSE   ! ll_dia_inst true
231       !
232       !! * Instantaneous output section
233       !! ------------------------------
234       !
235       IF( kt == nit000 .AND. lwp ) WRITE(numout,*) &
236            'THE OUTPUT FILES CONTAINS INSTANTANEOUS VALUES OF EACH FIELD'
237       !
238       zdtj=rdt/86400.   ! time step in days
239       !  iwrite=NINT(adatrj/rwrite)
240       clmode='instantaneous'
241       !     IF (abs(adatrj-iwrite*rwrite) < zdtj/2.  &
242       IF (  ( MOD (kt-nit000+1,nwrite) ==  0 )          &
243            &   .OR.       kindic <   0            &
244            &   .OR. ( kt == 1 .AND. kindic > 0)  ) THEN
245          !
246          ! transfer wp arrays to sp arrays for dimg files
247          fsel(:,:,:) = 0._wp
248          !
249          fsel(:,:,1 ) = utau(:,:) * umask(:,:,1)
250          fsel(:,:,2 ) = vtau(:,:) * vmask(:,:,1)
251          fsel(:,:,3 ) = (qsr (:,:) + qns (:,:)) * tmask(:,:,1)
252          fsel(:,:,4 ) = emp (:,:) * tmask(:,:,1)
253          fsel(:,:,5 ) = (tb  (:,:,1) - sf_sst(1)%fnow(:,:) ) *tmask(:,:,1)
254
255#if ! defined key_dynspg_rl
256          fsel(:,:,6 ) = sshn(:,:)
257#else
258          fsel(:,:,6 ) = bsfn(:,:)
259#endif
260          fsel(:,:,7 ) = qsr (:,:) * tmask(:,:,1)
261          fsel(:,:,8 ) = qrp (:,:) * tmask(:,:,1)
262          fsel(:,:,9 ) = erp (:,:) * tmask(:,:,1)
263          fsel(:,:,10) = hmld(:,:) * tmask(:,:,1)
264          fsel(:,:,11) = hmlp(:,:) * tmask(:,:,1)
265          fsel(:,:,12) = fr_i(:,:) * tmask(:,:,1)
266          fsel(:,:,13) = sf_sst(1)%fnow(:,:)
267          !         fsel(:,:,14) =  qct(:,:)
268          !         fsel(:,:,15) =  fbt(:,:)
269          fsel(:,:,16) =  emps(:,:) * tmask(:,:,1)
270#ifdef key_diaspr           
271          fsel(:,:,18) =      gps(:,:) /g
272          fsel(:,:,19) =      spgu(:,:)
273          fsel(:,:,20) =      spgv(:,:)
274#endif
275          !
276          !         qct(:,:) = 0._wp
277       ENDIF
278    ENDIF
279    !
280    ! Opening of the datrj.out file with the absolute time in day of each dump
281    ! this file gives a record of the dump date for post processing ( ASCII file )
282    !
283    IF(  ( MOD (kt-nit000+1,nwrite) ==  0 )          &
284         &   .OR.       kindic <   0            &
285         &   .OR. ( kt == 1 .AND. kindic > 0)  ) THEN
286
287       IF( lwp) WRITE(numout,*)'Days since the begining of the run :',adatrj
288
289       !! * U section
290
291       WRITE(clname,9000) TRIM(cexper),'U',iyear,imon,iday
292       cltext=TRIM(cexper)//' U(m/s) '//TRIM(clmode)
293       IF ( kindic < 0 )   cltext=TRIM(cexper)//' U(m/s)  instantaneous (explosion)'
294       !
295       IF( ll_dia_inst) THEN
296          CALL dia_wri_dimg(clname, cltext, un, jpk, 'T')
297
298       ELSE
299          IF( kindic ==  -3 ) THEN
300             ! ... in case of explosion on umax, dump instantateous u field instead of mean.
301             CALL dia_wri_dimg(clname, cltext, un, jpk, 'T')
302          ELSE
303             CALL dia_wri_dimg(clname, cltext, um, jpk, 'T')
304          ENDIF
305       ENDIF
306
307       !! * V section
308
309       WRITE(clname,9000) TRIM(cexper),'V',iyear,imon,iday
310       cltext=TRIM(cexper)//' V(m/s) '//TRIM(clmode)
311       !
312       IF( ll_dia_inst) THEN
313          CALL dia_wri_dimg(clname, cltext, vn, jpk, 'T')
314       ELSE
315          CALL dia_wri_dimg(clname, cltext, vm, jpk, 'T')
316       ENDIF
317       !
318
319       !! * KZ section
320
321       WRITE(clname,9000) TRIM(cexper),'KZ',iyear,imon,iday
322       cltext=TRIM(cexper)//' KZ(m2/s) '//TRIM(clmode)
323
324       IF( ll_dia_inst) THEN
325          CALL dia_wri_dimg(clname, cltext, avt, jpk, 'W')
326       ELSE
327          CALL dia_wri_dimg(clname, cltext, avtm, jpk, 'W')
328       ENDIF
329       !
330
331       !! * W section
332
333       WRITE(clname,9000) TRIM(cexper),'W',iyear,imon,iday
334       cltext=TRIM(cexper)//' W(m/s) '//TRIM(clmode)
335
336       IF( ll_dia_inst) THEN
337          CALL dia_wri_dimg(clname, cltext, wn, jpk, 'W')
338       ELSE
339          CALL dia_wri_dimg(clname, cltext, wm, jpk, 'W')
340       ENDIF
341
342       !! * T section
343
344       WRITE(clname,9000) TRIM(cexper),'T',iyear,imon,iday
345       cltext=TRIM(cexper)//' T (DegC) '//TRIM(clmode)
346
347       IF( ll_dia_inst) THEN
348          CALL dia_wri_dimg(clname, cltext, tn, jpk, 'T')
349       ELSE
350          CALL dia_wri_dimg(clname, cltext, tm, jpk, 'T')
351       ENDIF
352       !
353
354       !! * S section
355
356       WRITE(clname,9000) TRIM(cexper),'S',iyear,imon,iday
357       cltext=TRIM(cexper)//' S (PSU) '//TRIM(clmode)
358
359       IF( ll_dia_inst) THEN
360          CALL dia_wri_dimg(clname, cltext, sn, jpk, 'T')
361       ELSE
362          CALL dia_wri_dimg(clname, cltext, sm, jpk, 'T')
363       ENDIF
364       !
365
366       !! * 2D section
367
368       WRITE(clname,9000) TRIM(cexper),'2D',iyear,imon,iday
369       cltext='2D fields '//TRIM(clmode)
370
371       IF( ll_dia_inst) THEN
372          CALL dia_wri_dimg(clname, cltext, fsel, inbsel, '2')
373       ELSE
374          CALL dia_wri_dimg(clname, cltext, fsel, inbsel, '2')
375       ENDIF
376
377       IF( lk_mpp )   CALL mppsync   ! synchronization in mpp
378
379       !! * Log message in numout
380
381       IF( lwp)WRITE(numout,*) ' '
382       IF( lwp)WRITE(numout,*) ' **** WRITE in numwri ',kt
383
384       IF( lwp .AND.        ll_dia_inst) WRITE(numout,*) '    instantaneous fields'
385       IF( lwp .AND. .NOT.  ll_dia_inst) WRITE(numout,*) '    average fields with ',nmoyct,'pdt'
386       !
387       !
388       !! * Reset cumulating arrays  and counter to 0 after writing
389       !
390       IF( .NOT. ll_dia_inst ) THEN
391          nmoyct = 0
392          !
393          um(:,:,:) = 0._wp
394          vm(:,:,:) = 0._wp
395          wm(:,:,:) = 0._wp
396          tm(:,:,:) = 0._wp
397          sm(:,:,:) = 0._wp
398          fsel(:,:,:) = 0._wp
399          avtm(:,:,:) = 0._wp
400       ENDIF
401    ENDIF
402    !
4039000 FORMAT(a,"_",a,"_y",i4.4,"m",i2.2,"d",i2.2,".dimgproc")
404
405  END SUBROUTINE dia_wri
406
407  SUBROUTINE dia_wri_state ( cdfile_name)
408    !!-------------------------------------------------------------------
409    !!        ***     ROUTINE dia_wri_state  ***
410    !!
411    !! ** Purpose :   Dummy routine for compatibility with IOIPSL output
412    !!
413    !! ** History :
414    !!      9.O  ! 03-06  (J.M. Molines ) dimgout
415    !!--------------------------------------------------------------------
416    !! * Arguments
417    CHARACTER (len=*), INTENT(in) ::   cdfile_name   ! name of the file created
418    !!--------------------------------------------------------------------
419
420    IF( lwp) WRITE(numout,*) 'dia_wri_state: Dummy call', cdfile_name
421    IF( lwp) WRITE(numout,*) '-------------'
422    IF( lwp) WRITE(numout,*)
423
424  END SUBROUTINE dia_wri_state
Note: See TracBrowser for help on using the repository browser.