New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
dia25h.F90 in branches/2014/dev_r4650_UKMO10_Tidally_Meaned_Diagnostics/NEMOGCM/NEMO/OPA_SRC/DIA – NEMO

source: branches/2014/dev_r4650_UKMO10_Tidally_Meaned_Diagnostics/NEMOGCM/NEMO/OPA_SRC/DIA/dia25h.F90 @ 4756

Last change on this file since 4756 was 4756, checked in by deazer, 10 years ago

Added two new routines, diatmb and dia25h to handle 25hourly and tmb output.
modified diawri to call these routines when logicals are true
logicals are set by new namelist addition set to true in AMM12 cfg and false in reference
default should be false.
additional call in dynspg_ts for barotropic U and V
Created extra fields in field_def.xml and extrae file groups in iodef.xml

File size: 14.2 KB
Line 
1MODULE dia25h 
2   !!======================================================================
3   !!                       ***  MODULE  diaharm  ***
4   !! Harmonic analysis of tidal constituents
5   !!======================================================================
6   !! History :  3.6  !  2014  (E O'Dea)  Original code
7   !!----------------------------------------------------------------------
8   USE oce             ! ocean dynamics and tracers variables
9   USE dom_oce         ! ocean space and time domain
10   USE insitu_tem, ONLY: insitu_t, theta2t
11   USE in_out_manager  ! I/O units
12   USE iom             ! I/0 library
13   USE wrk_nemo        ! working arrays
14#if defined key_zdftke 
15   USE zdftke, ONLY: en
16#endif
17   USE zdf_oce, ONLY: avt, avm
18#if defined key_zdfgls
19   USE zdfgls, ONLY: mxln, en
20#endif
21
22   IMPLICIT NONE
23   PRIVATE
24
25   LOGICAL , PUBLIC ::   ln_dia25h     !:  25h mean output
26   PUBLIC   dia_25h_init               ! routine called by nemogcm.F90
27   PUBLIC   dia_25h                    ! routine called by diawri.F90
28
29  !! * variables for calculating 25-hourly means
30   REAL(wp),SAVE, ALLOCATABLE,   DIMENSION(:,:,:) ::   tn_25h  , sn_25h, insitu_t_25h 
31   REAL(wp),SAVE, ALLOCATABLE,   DIMENSION(:,:)   ::   sshn_25h 
32   REAL(wp),SAVE, ALLOCATABLE,   DIMENSION(:,:,:) ::   un_25h  , vn_25h  , wn_25h
33   REAL(wp),SAVE, ALLOCATABLE,   DIMENSION(:,:,:) ::   avt_25h , avm_25h
34#if defined key_zdfgls || key_zdftke
35   REAL(wp),SAVE, ALLOCATABLE,   DIMENSION(:,:,:) ::   en_25h
36#endif
37#if defined key_zdfgls 
38   REAL(wp),SAVE, ALLOCATABLE,   DIMENSION(:,:,:) ::   mxln_25h
39#endif
40   INTEGER, SAVE :: cnt_25h     ! Counter for 25 hour means
41
42
43
44   !!----------------------------------------------------------------------
45   !! NEMO/OPA 3.6 , NEMO Consortium (2014)
46   !! $Id:$
47   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
48   !!----------------------------------------------------------------------
49CONTAINS
50
51   SUBROUTINE dia_25h_init 
52      !!---------------------------------------------------------------------------
53      !!                  ***  ROUTINE dia_25h_init  ***
54      !!     
55      !! ** Purpose: Initialization of 25h mean namelist
56      !!       
57      !! ** Method : Read namelist
58      !!   History
59      !!   3.6  !  08-14  (E. O'Dea) Routine to initialize dia_25h
60      !!---------------------------------------------------------------------------
61      !!
62      INTEGER ::   ios                 ! Local integer output status for namelist read
63      INTEGER ::   ierror              ! Local integer for memory allocation
64      !
65      NAMELIST/nam_dia25h/ ln_dia25h
66      !!----------------------------------------------------------------------
67      !
68      REWIND ( numnam_ref )              ! Read Namelist nam_dia25h in reference namelist : 25hour mean diagnostics
69      READ   ( numnam_ref, nam_dia25h, IOSTAT=ios, ERR= 901 )
70901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nam_dia25h in reference namelist', lwp )
71
72      REWIND( numnam_cfg )              ! Namelist nam_dia25h in configuration namelist  25hour diagnostics
73      READ  ( numnam_cfg, nam_dia25h, IOSTAT = ios, ERR = 902 )
74902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nam_dia25h in configuration namelist', lwp )
75      IF(lwm) WRITE ( numond, nam_dia25h )
76
77      IF(lwp) THEN                   ! Control print
78         WRITE(numout,*)
79         WRITE(numout,*) 'dia_25h_init : Output 25 hour mean diagnostics'
80         WRITE(numout,*) '~~~~~~~~~~~~'
81         WRITE(numout,*) 'Namelist nam_dia25h : set 25h outputs '
82         WRITE(numout,*) 'Switch for 25h diagnostics (T) or not (F)  ln_dia25h  = ', ln_dia25h
83      ENDIF
84      IF( .NOT. ln_dia25h )   RETURN
85      ! ------------------- !
86      ! 1 - Allocate memory !
87      ! ------------------- !
88      ALLOCATE( tn_25h(jpi,jpj,jpk), STAT=ierror )
89      IF( ierror > 0 ) THEN
90         CALL ctl_stop( 'dia_25h: unable to allocate tn_25h' )   ;   RETURN
91      ENDIF
92      ALLOCATE( sn_25h(jpi,jpj,jpk), STAT=ierror )
93      IF( ierror > 0 ) THEN
94         CALL ctl_stop( 'dia_25h: unable to allocate sn_25h' )   ;   RETURN
95      ENDIF
96      ALLOCATE( insitu_t_25h(jpi,jpj,jpk), STAT=ierror )
97      IF( ierror > 0 ) THEN
98         CALL ctl_stop( 'dia_25h: unable to allocate insitu_t_25h' )   ;   RETURN
99      ENDIF
100      ALLOCATE( un_25h(jpi,jpj,jpk), STAT=ierror )
101      IF( ierror > 0 ) THEN
102         CALL ctl_stop( 'dia_25h: unable to allocate un_25h' )   ;   RETURN
103      ENDIF
104      ALLOCATE( vn_25h(jpi,jpj,jpk), STAT=ierror )
105      IF( ierror > 0 ) THEN
106         CALL ctl_stop( 'dia_25h: unable to allocate vn_25h' )   ;   RETURN
107      ENDIF
108      ALLOCATE( wn_25h(jpi,jpj,jpk), STAT=ierror )
109      IF( ierror > 0 ) THEN
110         CALL ctl_stop( 'dia_25h: unable to allocate wn_25h' )   ;   RETURN
111      ENDIF
112      ALLOCATE( avt_25h(jpi,jpj,jpk), STAT=ierror )
113      IF( ierror > 0 ) THEN
114         CALL ctl_stop( 'dia_25h: unable to allocate avt_25h' )   ;   RETURN
115      ENDIF
116      ALLOCATE( avm_25h(jpi,jpj,jpk), STAT=ierror )
117      IF( ierror > 0 ) THEN
118         CALL ctl_stop( 'dia_25h: unable to allocate avm_25h' )   ;   RETURN
119      ENDIF
120# if defined key_zdfgls || defined key_zdftke
121      ALLOCATE( en_25h(jpi,jpj,jpk), STAT=ierror )
122      IF( ierror > 0 ) THEN
123         CALL ctl_stop( 'dia_25h: unable to allocate en_25h' )   ;   RETURN
124      ENDIF
125#endif
126# if defined key_zdfgls 
127      ALLOCATE( mxln_25h(jpi,jpj,jpk), STAT=ierror )
128      IF( ierror > 0 ) THEN
129         CALL ctl_stop( 'dia_25h: unable to allocate mxln_25h' )   ;   RETURN
130      ENDIF
131#endif
132      ALLOCATE( sshn_25h(jpi,jpj), STAT=ierror )
133      IF( ierror > 0 ) THEN
134         CALL ctl_stop( 'dia_25h: unable to allocate sshn_25h' )   ;   RETURN
135      ENDIF
136      ! ------------------------- !
137      ! 2 - Assign Initial Values !
138      ! ------------------------- !
139      cnt_25h = 1  ! sets the first value of sum at timestep 1 (note - should strictly be at timestep zero so before values used where possible)
140      tn_25h(:,:,:) = tsb(:,:,:,jp_tem)
141      sn_25h(:,:,:) = tsb(:,:,:,jp_sal)
142      CALL theta2t
143      insitu_t_25h(:,:,:) = insitu_t(:,:,:)
144      sshn_25h(:,:) = sshb(:,:)
145      un_25h(:,:,:) = ub(:,:,:)
146      vn_25h(:,:,:) = vb(:,:,:)
147      wn_25h(:,:,:) = wn(:,:,:)
148      avt_25h(:,:,:) = avt(:,:,:)
149      avm_25h(:,:,:) = avm(:,:,:)
150# if defined key_zdfgls || defined key_zdftke
151         en_25h(:,:,:) = en(:,:,:)
152#endif
153# if defined key_zdfgls
154         mxln_25h(:,:,:) = mxln(:,:,:)
155#endif
156#if defined key_lim3 || defined key_lim2
157         CALL ctl_stop('STOP', 'dia_25h not setup yet to do tidemean ice')
158#endif 
159
160      ! -------------------------- !
161      ! 3 - Return to dia_wri      !
162      ! -------------------------- !
163
164
165   END SUBROUTINE dia_25h_init
166
167
168   SUBROUTINE dia_25h( kt ) 
169      !!----------------------------------------------------------------------
170      !!                 ***  ROUTINE dia_25h  ***
171      !!         
172      !!
173      !!--------------------------------------------------------------------
174      !!                   
175      !! ** Purpose :   Write diagnostics with M2/S2 tide removed
176      !!
177      !! ** Method  :   
178      !!      25hr mean outputs for shelf seas
179      !!
180      !! History :
181      !!   ?.0  !  07-04  (A. Hines) New routine, developed from dia_wri_foam
182      !!   3.4  !  02-13  (J. Siddorn) Routine taken from old dia_wri_foam
183      !!   3.6  !  08-14  (E. O'Dea) adapted for VN3.6
184      !!----------------------------------------------------------------------
185      !! * Modules used
186
187      IMPLICIT NONE
188
189      !! * Arguments
190      INTEGER, INTENT( in ) ::   kt      ! ocean time-step index
191
192
193      !! * Local declarations
194      INTEGER ::   ji, jj, jk
195
196      LOGICAL ::   ll_print = .FALSE.    ! =T print and flush numout
197      REAL(wp)                         ::   zsto, zout, zmax, zjulian, zdt, zmdi  ! temporary reals
198      INTEGER                          ::   i_steps                               ! no of timesteps per hour
199      REAL(wp), DIMENSION(jpi,jpj    ) ::   zw2d, un_dm, vn_dm                    ! temporary workspace
200      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zw3d                                  ! temporary workspace
201      REAL(wp), DIMENSION(jpi,jpj,3)   ::   zwtmb                                 ! temporary workspace
202      INTEGER                          ::   nyear0, nmonth0,nday0                 ! start year,month,day
203
204      !!----------------------------------------------------------------------
205
206      ! 0. Initialisation
207      ! -----------------
208      ! Define frequency of summing to create 25 h mean
209      zdt = rdt
210      IF( nacc == 1 ) zdt = rdtmin
211
212      IF( MOD( 3600,INT(zdt) ) == 0 ) THEN
213         i_steps = 3600/INT(zdt)
214      ELSE
215         CALL ctl_stop('STOP', 'dia_wri_tide: timestep must give MOD(3600,rdt) = 0 otherwise no hourly values are possible')
216      ENDIF
217
218#if defined key_lim3 || defined key_lim2
219      CALL ctl_stop('STOP', 'dia_wri_tide not setup yet to do tidemean ice')
220#endif
221
222      ! local variable for debugging
223      ll_print = ll_print .AND. lwp
224
225      ! Sum of 25 hourly instantaneous values to give a 25h mean from 24hours
226      ! every day
227      IF( MOD( kt, i_steps ) == 0  .and. kt .ne. nn_it000 ) THEN
228
229         IF (lwp) THEN
230              WRITE(numout,*) 'dia_wri_tide : Summing instantaneous hourly diagnostics at timestep ',kt
231              WRITE(numout,*) '~~~~~~~~~~~~ '
232         ENDIF
233
234         tn_25h(:,:,:)        = tn_25h(:,:,:) + tsn(:,:,:,jp_tem)
235         sn_25h(:,:,:)        = sn_25h(:,:,:) + tsn(:,:,:,jp_sal)
236         CALL theta2t
237         insitu_t_25h(:,:,:)  = insitu_t_25h(:,:,:) + insitu_t(:,:,:)
238         sshn_25h(:,:)        = sshn_25h(:,:) + sshn (:,:)
239         un_25h(:,:,:)        = un_25h(:,:,:) + un(:,:,:)
240         vn_25h(:,:,:)        = vn_25h(:,:,:) + vn(:,:,:)
241         wn_25h(:,:,:)        = wn_25h(:,:,:) + wn(:,:,:)
242         avt_25h(:,:,:)       = avt_25h(:,:,:) + avt(:,:,:)
243         avm_25h(:,:,:)       = avm_25h(:,:,:) + avm(:,:,:)
244# if defined key_zdfgls || defined key_zdftke
245         en_25h(:,:,:)        = en_25h(:,:,:) + en(:,:,:)
246#endif
247# if defined key_zdfgls
248         mxln_25h(:,:,:)      = mxln_25h(:,:,:) + mxln(:,:,:)
249#endif
250         cnt_25h = cnt_25h + 1
251
252         IF (lwp) THEN
253            WRITE(numout,*) 'dia_tide : Summed the following number of hourly values so far',cnt_25h
254         ENDIF
255
256      ENDIF ! MOD( kt, i_steps ) == 0
257
258         ! Write data for 25 hour mean output streams
259      IF( cnt_25h .EQ. 25 .AND.  MOD( kt, i_steps*24) == 0 .AND. kt .NE. nn_it000 ) THEN
260
261            IF(lwp) THEN
262               WRITE(numout,*) 'dia_wri_tide : Writing 25 hour mean tide diagnostics at timestep', kt
263               WRITE(numout,*) '~~~~~~~~~~~~ '
264            ENDIF
265
266            tn_25h(:,:,:)        = tn_25h(:,:,:) / 25.0_wp
267            sn_25h(:,:,:)        = sn_25h(:,:,:) / 25.0_wp
268            insitu_t_25h(:,:,:)  = insitu_t_25h(:,:,:) / 25.0_wp
269            sshn_25h(:,:)        = sshn_25h(:,:) / 25.0_wp
270            un_25h(:,:,:)        = un_25h(:,:,:) / 25.0_wp
271            vn_25h(:,:,:)        = vn_25h(:,:,:) / 25.0_wp
272            wn_25h(:,:,:)        = wn_25h(:,:,:) / 25.0_wp
273            avt_25h(:,:,:)       = avt_25h(:,:,:) / 25.0_wp
274            avm_25h(:,:,:)       = avm_25h(:,:,:) / 25.0_wp
275# if defined key_zdfgls || defined key_zdftke
276            en_25h(:,:,:)        = en_25h(:,:,:) / 25.0_wp
277#endif
278# if defined key_zdfgls
279            mxln_25h(:,:,:)       = mxln_25h(:,:,:) / 25.0_wp
280#endif
281
282            IF (lwp)  WRITE(numout,*) 'dia_wri_tide : Mean calculated by dividing 25 hour sums and writing output'
283            zmdi=1.e+20 !missing data indicator for masking
284            ! write tracers (instantaneous)
285            zw3d(:,:,:) = tn_25h(:,:,:)*tmask(:,:,:) + zmdi*(1.0-tmask(:,:,:))
286            CALL iom_put("temper25h", zw3d)   ! potential temperature
287            CALL theta2t                                                                    ! calculate insitu temp
288            zw3d(:,:,:) = insitu_t_25h(:,:,:)*tmask(:,:,:) + zmdi*(1.0-tmask(:,:,:))
289            CALL iom_put("tempis25h", zw3d)   ! in-situ temperature
290            zw3d(:,:,:) = sn_25h(:,:,:)*tmask(:,:,:) + zmdi*(1.0-tmask(:,:,:))
291            CALL iom_put( "salin25h", zw3d  )   ! salinity
292            zw2d(:,:) = sshn_25h(:,:)*tmask(:,:,1) + zmdi*(1.0-tmask(:,:,1))
293            CALL iom_put( "ssh25h", zw2d )   ! sea surface
294
295
296            ! Write velocities (instantaneous)
297            zw3d(:,:,:) = un_25h(:,:,:)*umask(:,:,:) + zmdi*(1.0-umask(:,:,:))
298            CALL iom_put("vozocrtx25h", zw3d)    ! i-current
299            zw3d(:,:,:) = vn_25h(:,:,:)*vmask(:,:,:) + zmdi*(1.0-vmask(:,:,:))
300            CALL iom_put("vomecrty25h", zw3d  )   ! j-current
301
302            zw3d(:,:,:) = wn_25h(:,:,:)*tmask(:,:,:) + zmdi*(1.0-tmask(:,:,:))
303            CALL iom_put("vomecrtz25h", zw3d )   ! k-current
304            zw3d(:,:,:) = avt_25h(:,:,:)*tmask(:,:,:) + zmdi*(1.0-tmask(:,:,:))
305            CALL iom_put("avt25h", zw3d )   ! diffusivity
306            zw3d(:,:,:) = avm_25h(:,:,:)*tmask(:,:,:) + zmdi*(1.0-tmask(:,:,:))
307            CALL iom_put("avm25h", zw3d)   ! viscosity
308#if defined key_zdftke || defined key_zdfgls 
309            zw3d(:,:,:) = en_25h(:,:,:)*tmask(:,:,:) + zmdi*(1.0-tmask(:,:,:))
310            CALL iom_put("tke25h", zw3d)   ! tke
311#endif
312#if defined key_zdfgls 
313            zw3d(:,:,:) = mxln_25h(:,:,:)*tmask(:,:,:) + zmdi*(1.0-tmask(:,:,:))
314            CALL iom_put( "mxln25h",zw3d)
315#endif
316
317            ! After the write reset the values to cnt=1 and sum values equal current value
318            tn_25h(:,:,:) = tsn(:,:,:,jp_tem)
319            sn_25h(:,:,:) = tsn(:,:,:,jp_sal)
320            CALL theta2t
321            insitu_t_25h(:,:,:) = insitu_t(:,:,:)
322            sshn_25h(:,:) = sshn (:,:)
323            un_25h(:,:,:) = un(:,:,:)
324            vn_25h(:,:,:) = vn(:,:,:)
325            wn_25h(:,:,:) = wn(:,:,:)
326            avt_25h(:,:,:) = avt(:,:,:)
327            avm_25h(:,:,:) = avm(:,:,:)
328# if defined key_zdfgls || defined key_zdftke
329            en_25h(:,:,:) = en(:,:,:)
330#endif
331# if defined key_zdfgls
332            mxln_25h(:,:,:) = mxln(:,:,:)
333#endif
334            cnt_25h = 1
335            IF (lwp)  WRITE(numout,*) 'dia_wri_tide : After 25hr mean write, reset sum to current value and cnt_25h to one for overlapping average',cnt_25h
336
337      ENDIF !  cnt_25h .EQ. 25 .AND.  MOD( kt, i_steps * 24) == 0 .AND. kt .NE. nn_it000
338
339
340   END SUBROUTINE dia_25h 
341
342   !!======================================================================
343END MODULE dia25h
Note: See TracBrowser for help on using the repository browser.