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/UKMO/r5936_CO6_CO5_shelfdiagnostic/NEMOGCM/NEMO/OPA_SRC/DIA – NEMO

source: branches/UKMO/r5936_CO6_CO5_shelfdiagnostic/NEMOGCM/NEMO/OPA_SRC/DIA/dia25h.F90 @ 7104

Last change on this file since 7104 was 7104, checked in by jcastill, 8 years ago

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