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.
diaptr.F90 in trunk/NEMO/OPA_SRC/DIA – NEMO

source: trunk/NEMO/OPA_SRC/DIA/diaptr.F90 @ 1397

Last change on this file since 1397 was 1397, checked in by rblod, 15 years ago

Correct diaptr with correct suffix file in mono case, see ticket #361

  • Property svn:eol-style set to native
  • Property svn:keywords set to Id
File size: 46.0 KB
Line 
1MODULE diaptr
2   !!======================================================================
3   !!                       ***  MODULE  diaptr  ***
4   !! Ocean physics:  Computes meridonal transports and zonal means
5   !!=====================================================================
6   !! History :  9.0  !  03-09  (C. Talandier, G. Madec)  Original code
7   !!            9.0  !  06-01  (A. Biastoch)  Allow sub-basins computation
8   !!            9.0  !  03-09  (O. Marti) Add fields
9   !!----------------------------------------------------------------------
10
11   !!----------------------------------------------------------------------
12   !!   dia_ptr      : Poleward Transport Diagnostics module
13   !!   dia_ptr_init : Initialization, namelist read
14   !!   dia_ptr_wri  : Output of poleward fluxes
15   !!   ptr_vjk      : "zonal" sum computation of a "meridional" flux array
16   !!   ptr_tjk      : "zonal" mean computation of a tracer field
17   !!   ptr_vj       : "zonal" and vertical sum computation of a "meridional"
18   !!                : flux array; Generic interface: ptr_vj_3d, ptr_vj_2d
19   !!----------------------------------------------------------------------
20   USE oce           ! ocean dynamics and active tracers
21   USE dom_oce       ! ocean space and time domain
22   USE ldftra_oce    ! ocean active tracers: lateral physics
23   USE lib_mpp
24   USE in_out_manager
25   USE dianam
26   USE phycst
27   USE iom
28   USE ioipsl         
29   USE daymod
30
31   IMPLICIT NONE
32   PRIVATE
33
34   INTERFACE ptr_vj
35      MODULE PROCEDURE ptr_vj_3d, ptr_vj_2d
36   END INTERFACE
37
38   PUBLIC   dia_ptr_init   ! call in opa module
39   PUBLIC   dia_ptr        ! call in step module
40   PUBLIC   ptr_vj         ! call by tra_ldf & tra_adv routines
41   PUBLIC   ptr_vjk        ! call by tra_ldf & tra_adv routines
42
43!!! ** init namelist (namptr)
44   LOGICAL , PUBLIC                 ::   ln_diaptr = .FALSE.   !: Poleward transport flag (T) or not (F)
45   LOGICAL , PUBLIC                 ::   ln_subbas = .FALSE.   !: Atlantic/Pacific/Indian basins calculation
46   LOGICAL , PUBLIC                 ::   ln_diaznl = .FALSE.   !: Add zonal means and meridional stream functions
47   LOGICAL , PUBLIC                 ::   ln_ptrcomp = .FALSE.  !: Add decomposition : overturning (and gyre, soon ...)
48   INTEGER , PUBLIC                 ::   nf_ptr = 15           !: frequency of ptr computation
49   INTEGER , PUBLIC                 ::   nf_ptr_wri = 15       !: frequency of ptr outputs
50
51   REAL(wp), DIMENSION(jpi,jpj), PUBLIC ::   abasin, pbasin, ibasin, dbasin, sbasin !: Sub basin masks
52
53   REAL(wp), PUBLIC, DIMENSION(jpj) ::   pht_adv, pst_adv  !: heat and salt poleward transport: advection
54   REAL(wp), PUBLIC, DIMENSION(jpj) ::   pht_ove_glo, pst_ove_glo, pht_ove_atl, pst_ove_atl, pht_ove_pac, pst_ove_pac, &
55      &        pht_ove_ind, pst_ove_ind, pht_ove_ipc, pst_ove_ipc  !: heat and salt poleward transport: overturning
56   REAL(wp), PUBLIC, DIMENSION(jpj) ::   pht_ldf, pst_ldf      !: heat and salt poleward transport: lateral diffusion
57#if defined key_diaeiv
58   REAL(wp), PUBLIC, DIMENSION(jpj) ::   pht_eiv_glo, pst_eiv_glo, pht_eiv_atl, pst_eiv_atl, pht_eiv_pac, pst_eiv_pac, &
59      &        pht_eiv_ind, pst_eiv_ind, pht_eiv_ipc, pst_eiv_ipc !: heat and salt poleward transport: bolus advection
60#endif
61   REAL(wp), PUBLIC, DIMENSION(jpj) ::   ht_glo,ht_atl,ht_ind,ht_pac,ht_ipc !: heat
62   REAL(wp), PUBLIC, DIMENSION(jpj) ::   st_glo,st_atl,st_ind,st_pac,st_ipc !: salt
63
64   INTEGER :: niter
65   INTEGER :: nidom_ptr
66
67   REAL(wp), DIMENSION(jpj,jpk) ::   tn_jk_glo, sn_jk_glo,  &  !: "zonal" mean temperature and salinity
68      &                              tn_jk_atl, sn_jk_atl,  &
69      &                              tn_jk_pac, sn_jk_pac,  &
70      &                              tn_jk_ind, sn_jk_ind,  &
71      &                              tn_jk_ipc, sn_jk_ipc,  &
72      &                              v_msf_glo       ,  &  !: "meridional" Stream-Function
73      &                              v_msf_atl       ,  & 
74      &                              v_msf_pac       ,  & 
75      &                              v_msf_ind       ,  & 
76      &                              v_msf_ipc       ,  & 
77      &                              surf_jk_glo     ,  &  !: Ocean "zonal" section surface
78      &                              surf_jk_atl     ,  & 
79      &                              surf_jk_pac     ,  & 
80      &                              surf_jk_ind     ,  &         
81      &                              surf_jk_ipc     ,  & 
82      &                              surf_jk_r_glo   ,  &  !: inverse of the ocean "zonal" section surface
83      &                              surf_jk_r_atl   ,  & 
84      &                              surf_jk_r_pac   ,  & 
85      &                              surf_jk_r_ind   ,  &         
86      &                              surf_jk_r_ipc       
87#if defined key_diaeiv
88   REAL(wp), DIMENSION(jpj,jpk) ::   v_msf_eiv_glo, v_msf_eiv_atl, v_msf_eiv_pac,  &
89      &                              v_msf_eiv_ind, v_msf_eiv_ipc !: bolus "meridional" Stream-Function
90#endif
91 
92   !! * Substitutions
93#  include "domzgr_substitute.h90"
94#  include "vectopt_loop_substitute.h90"
95   !!----------------------------------------------------------------------
96   !!   OPA 9.0 , LOCEAN-IPSL (2005)
97   !! $Id$
98   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
99   !!----------------------------------------------------------------------
100
101CONTAINS
102
103   FUNCTION ptr_vj_3d( pva )   RESULT ( p_fval )
104      !!----------------------------------------------------------------------
105      !!                    ***  ROUTINE ptr_vj_3d  ***
106      !!
107      !! ** Purpose :   "zonal" and vertical sum computation of a "meridional"
108      !!      flux array
109      !!
110      !! ** Method  : - i-k sum of pva using the interior 2D vmask (vmask_i).
111      !!      pva is supposed to be a masked flux (i.e. * vmask*e1v*e3v)
112      !!
113      !! ** Action  : - p_fval: i-k-mean poleward flux of pva
114      !!----------------------------------------------------------------------
115      REAL(wp) , INTENT(in), DIMENSION(jpi,jpj,jpk) ::   pva   ! mask flux array at V-point
116      !!
117      INTEGER                  ::   ji, jj, jk   ! dummy loop arguments
118      INTEGER                  ::   ijpj         ! ???
119      REAL(wp), DIMENSION(jpj) ::   p_fval       ! function value
120      !!--------------------------------------------------------------------
121      !
122      ijpj = jpj
123      p_fval(:) = 0.e0
124      DO jk = 1, jpkm1
125         DO jj = 2, jpjm1
126            DO ji = fs_2, fs_jpim1   ! Vector opt.
127               p_fval(jj) = p_fval(jj) + pva(ji,jj,jk) * tmask_i(ji,jj+1) * tmask_i(ji,jj) 
128            END DO
129         END DO
130      END DO
131      !
132#if defined key_mpp_mpi
133      CALL mpp_sum( p_fval, ijpj, ncomm_znl)     !!bug  I presume
134#endif
135      !
136   END FUNCTION ptr_vj_3d
137
138
139   FUNCTION ptr_vj_2d( pva )   RESULT ( p_fval )
140      !!----------------------------------------------------------------------
141      !!                    ***  ROUTINE ptr_vj_2d  ***
142      !!
143      !! ** Purpose :   "zonal" and vertical sum computation of a "meridional"
144      !!      flux array
145      !!
146      !! ** Method  : - i-k sum of pva using the interior 2D vmask (vmask_i).
147      !!      pva is supposed to be a masked flux (i.e. * vmask*e1v*e3v)
148      !!
149      !! ** Action  : - p_fval: i-k-mean poleward flux of pva
150      !!----------------------------------------------------------------------
151      REAL(wp) , INTENT(in), DIMENSION(jpi,jpj) ::   pva   ! mask flux array at V-point
152      !!
153      INTEGER                  ::   ji,jj    ! dummy loop arguments
154      INTEGER                  ::   ijpj     ! ???
155      REAL(wp), DIMENSION(jpj) ::   p_fval   ! function value
156      !!--------------------------------------------------------------------
157      !
158      ijpj = jpj
159      p_fval(:) = 0.e0
160      DO jj = 2, jpjm1
161         DO ji = nldi, nlei   ! No vector optimisation here. Better use a mask ?
162            p_fval(jj) = p_fval(jj) + pva(ji,jj) * tmask_i(ji,jj+1) * tmask_i(ji,jj)
163         END DO
164      END DO
165      !
166#if defined key_mpp_mpi
167      CALL mpp_sum( p_fval, ijpj, ncomm_znl )     !!bug  I presume
168#endif
169      !
170   END FUNCTION ptr_vj_2d
171
172
173   FUNCTION ptr_vjk( pva, bmask )   RESULT ( p_fval )
174      !!----------------------------------------------------------------------
175      !!                    ***  ROUTINE ptr_vjk  ***
176      !!
177      !! ** Purpose :   "zonal" sum computation of a "meridional" flux array
178      !!
179      !! ** Method  : - i-sum of pva using the interior 2D vmask (vmask_i).
180      !!      pva is supposed to be a masked flux (i.e. * vmask*e1v*e3v)
181      !!
182      !! ** Action  : - p_fval: i-k-mean poleward flux of pva
183      !!----------------------------------------------------------------------
184      REAL(wp) , INTENT(in), DIMENSION(jpi,jpj,jpk) ::   pva   ! mask flux array at V-point
185      REAL(wp) , INTENT(in), DIMENSION(jpi,jpj), OPTIONAL :: bmask ! Optional 2D basin mask
186      !!
187      INTEGER                      ::   ji, jj, jk   ! dummy loop arguments
188      INTEGER , DIMENSION (1)      ::   ish
189      INTEGER , DIMENSION (2)      ::   ish2
190      REAL(wp), DIMENSION(jpj*jpk) ::   zwork        ! temporary vector for mpp_sum
191      REAL(wp), DIMENSION(jpj,jpk) ::   p_fval       ! return function value
192      !!--------------------------------------------------------------------
193      !
194      p_fval(:,:) = 0.e0
195      !
196      IF (PRESENT (bmask)) THEN
197         DO jk = 1, jpkm1
198            DO jj = 2, jpjm1
199               DO ji =  nldi, nlei   ! No vector optimisation here. Better use a mask ?
200                  p_fval(jj,jk) = p_fval(jj,jk) + pva(ji,jj,jk) * e1v(ji,jj) * fse3v(ji,jj,jk)   &
201                     &                                          * tmask_i(ji,jj+1) * tmask_i(ji,jj) &
202                     &                                          * bmask(ji,jj)
203               END DO
204            END DO
205         END DO
206      ELSE
207         DO jk = 1, jpkm1
208            DO jj = 2, jpjm1
209               DO ji =  nldi, nlei   ! No vector optimisation here. Better use a mask ?
210                  p_fval(jj,jk) = p_fval(jj,jk) + pva(ji,jj,jk) * e1v(ji,jj) * fse3v(ji,jj,jk)   &
211                     &                                          * tmask_i(ji,jj+1) * tmask_i(ji,jj)
212               END DO
213            END DO
214         END DO
215      END IF
216      !
217#if defined key_mpp_mpi
218      ish(1) = jpj*jpk  ;  ish2(1) = jpj  ;  ish2(2) = jpk
219      zwork(:)= RESHAPE( p_fval, ish )
220      CALL mpp_sum( zwork, jpj*jpk, ncomm_znl )
221      p_fval(:,:)= RESHAPE( zwork, ish2 )
222#endif
223      !
224   END FUNCTION ptr_vjk
225
226   FUNCTION ptr_tjk( pta, bmask )   RESULT ( p_fval )
227      !!----------------------------------------------------------------------
228      !!                    ***  ROUTINE ptr_tjk  ***
229      !!
230      !! ** Purpose :   "zonal" mean computation of a tracer field
231      !!
232      !! ** Method  : - i-sum of mj(pta) using tmask
233      !!      multiplied by the inverse of the surface of the "zonal" ocean
234      !!      section
235      !!
236      !! ** Action  : - p_fval: i-k-mean poleward flux of pta
237      !!----------------------------------------------------------------------
238      REAL(wp) , INTENT(in), DIMENSION(jpi,jpj,jpk) ::   pta    ! tracer flux array at T-point
239      REAL(wp) , INTENT(in), DIMENSION(jpi,jpj), OPTIONAL :: bmask ! Optional 2D basin mask
240      !!
241      INTEGER                     ::   ji, jj, jk   ! dummy loop arguments
242      INTEGER, DIMENSION (1)      ::   ish
243      INTEGER, DIMENSION (2)      ::   ish2
244      REAL(wp),DIMENSION(jpj*jpk) ::   zwork        ! temporary vector for mpp_sum
245      REAL(wp),DIMENSION(jpj,jpk) ::   p_fval       ! return function value
246      !!--------------------------------------------------------------------
247      !
248      p_fval(:,:) = 0.e0
249      IF (PRESENT (bmask)) THEN
250         DO jk = 1, jpkm1
251            DO jj = 2, jpjm1
252               DO ji =  nldi, nlei   ! No vector optimisation here. Better use a mask ?
253                  p_fval(jj,jk) = p_fval(jj,jk) + pta(ji,jj,jk)                  &
254                     &                          * e1t(ji,jj) * fse3t(ji,jj,jk)   &
255                     &                          * tmask_i(ji,jj)                 &
256                     &                          * bmask(ji,jj)
257               END DO
258            END DO
259         END DO
260      ELSE
261         DO jk = 1, jpkm1
262            DO jj = 2, jpjm1
263               DO ji =  nldi, nlei   ! No vector optimisation here. Better use a mask ?
264                  p_fval(jj,jk) = p_fval(jj,jk) + pta(ji,jj,jk)                  &
265                     &                          * e1t(ji,jj) * fse3t(ji,jj,jk)   &
266                     &                          * tmask_i(ji,jj)
267               END DO
268            END DO
269         END DO
270      END IF
271      p_fval(:,:) = p_fval(:,:) * 0.5
272#if defined key_mpp_mpi
273      ish(1) = jpj*jpk   ;   ish2(1) = jpj   ;   ish2(2) = jpk
274      zwork(:)= RESHAPE( p_fval, ish )
275      CALL mpp_sum( zwork, jpj*jpk, ncomm_znl )
276      p_fval(:,:)= RESHAPE(zwork,ish2)
277#endif
278      !
279   END FUNCTION ptr_tjk
280
281
282   SUBROUTINE dia_ptr( kt )
283      !!----------------------------------------------------------------------
284      !!                  ***  ROUTINE dia_ptr  ***
285      !!----------------------------------------------------------------------
286      INTEGER, INTENT(in) ::   kt   ! ocean time step index
287      !!
288      INTEGER  ::   jk, jj, ji               ! dummy loop
289      REAL(wp) ::   zsverdrup,  &              ! conversion from m3/s to Sverdrup
290         &          zpwatt,     &              ! conversion from W    to PW
291         &          zggram                     ! conversion from g    to Pg
292      REAL(wp), DIMENSION(jpi,jpj,jpk) :: vt, vs
293      !!----------------------------------------------------------------------
294
295      IF( kt == nit000 .OR. MOD( kt, nf_ptr ) == 0 )   THEN
296
297         IF ( MOD( kt, nf_ptr ) == 0 ) THEN
298
299            zsverdrup = 1.e-6
300            zpwatt    = 1.e-15
301            zggram    = 1.e-6
302
303            IF ( ln_diaznl ) THEN
304               ! "zonal" mean temperature and salinity at V-points
305               tn_jk_glo(:,:) = ptr_tjk( tn(:,:,:) ) * surf_jk_r_glo(:,:)
306               sn_jk_glo(:,:) = ptr_tjk( sn(:,:,:) ) * surf_jk_r_glo(:,:)
307
308               IF (ln_subbas) THEN
309                  tn_jk_atl(:,:) = ptr_tjk( tn(:,:,:), abasin(:,:) ) * surf_jk_r_atl(:,:)
310                  sn_jk_atl(:,:) = ptr_tjk( sn(:,:,:), abasin(:,:) ) * surf_jk_r_atl(:,:)
311                  tn_jk_pac(:,:) = ptr_tjk( tn(:,:,:), pbasin(:,:) ) * surf_jk_r_pac(:,:)
312                  sn_jk_pac(:,:) = ptr_tjk( sn(:,:,:), pbasin(:,:) ) * surf_jk_r_pac(:,:)
313                  tn_jk_ind(:,:) = ptr_tjk( tn(:,:,:), ibasin(:,:) ) * surf_jk_r_ind(:,:)
314                  sn_jk_ind(:,:) = ptr_tjk( sn(:,:,:), ibasin(:,:) ) * surf_jk_r_ind(:,:)
315                  tn_jk_ipc(:,:) = ptr_tjk( tn(:,:,:), dbasin(:,:) ) * surf_jk_r_ipc(:,:)
316                  sn_jk_ipc(:,:) = ptr_tjk( sn(:,:,:), dbasin(:,:) ) * surf_jk_r_ipc(:,:)
317               ENDIF
318            ENDIF
319
320            !--------------------------------------------------------
321            ! overturning calculation:
322
323            ! horizontal integral and vertical dz
324
325#if defined key_diaeiv
326            v_msf_glo(:,:) = ptr_vjk( vn(:,:,:)+v_eiv(:,:,:) ) 
327            IF( ln_subbas .AND. ln_diaznl ) THEN
328               v_msf_atl(:,:) = ptr_vjk( vn (:,:,:)+v_eiv(:,:,:), abasin(:,:)*sbasin(:,:) ) 
329               v_msf_pac(:,:) = ptr_vjk( vn (:,:,:)+v_eiv(:,:,:), pbasin(:,:)*sbasin(:,:) ) 
330               v_msf_ind(:,:) = ptr_vjk( vn (:,:,:)+v_eiv(:,:,:), ibasin(:,:)*sbasin(:,:) ) 
331               v_msf_ipc(:,:) = ptr_vjk( vn (:,:,:)+v_eiv(:,:,:), dbasin(:,:)*sbasin(:,:) ) 
332            ENDIF
333#else
334            v_msf_glo(:,:) = ptr_vjk( vn(:,:,:) ) 
335            IF( ln_subbas .AND. ln_diaznl ) THEN
336               v_msf_atl(:,:) = ptr_vjk( vn (:,:,:), abasin(:,:)*sbasin(:,:) ) 
337               v_msf_pac(:,:) = ptr_vjk( vn (:,:,:), pbasin(:,:)*sbasin(:,:) ) 
338               v_msf_ind(:,:) = ptr_vjk( vn (:,:,:), ibasin(:,:)*sbasin(:,:) ) 
339               v_msf_ipc(:,:) = ptr_vjk( vn (:,:,:), dbasin(:,:)*sbasin(:,:) ) 
340            ENDIF
341#endif
342
343#if defined key_diaeiv
344            v_msf_eiv_glo(:,:) = ptr_vjk( v_eiv(:,:,:) )
345            IF (ln_subbas ) THEN
346               v_msf_eiv_atl(:,:) = ptr_vjk( v_eiv(:,:,:), abasin(:,:)*sbasin(:,:) )
347               v_msf_eiv_pac(:,:) = ptr_vjk( v_eiv(:,:,:), pbasin(:,:)*sbasin(:,:) )
348               v_msf_eiv_ind(:,:) = ptr_vjk( v_eiv(:,:,:), ibasin(:,:)*sbasin(:,:) )
349               v_msf_eiv_ipc(:,:) = ptr_vjk( v_eiv(:,:,:), dbasin(:,:)*sbasin(:,:) )
350            END IF
351#endif
352
353            ! "Meridional" Stream-Function
354            DO jk = 2,jpk 
355               v_msf_glo(:,jk) = v_msf_glo(:,jk-1) + v_msf_glo(:,jk)
356            END DO
357            v_msf_glo(:,:) = v_msf_glo(:,:) * zsverdrup
358#if defined key_diaeiv
359            ! Bolus "Meridional" Stream-Function
360            DO jk = 2,jpk
361               v_msf_eiv_glo(:,jk) = v_msf_eiv_glo(:,jk-1) + v_msf_eiv_glo(:,jk)
362            END DO
363            v_msf_eiv_glo(:,:) = v_msf_eiv_glo(:,:) * zsverdrup
364            IF ( ln_subbas ) THEN
365               DO jk = 2,jpk
366                  v_msf_eiv_atl(:,jk) = v_msf_eiv_atl(:,jk-1) + v_msf_eiv_atl(:,jk)
367                  v_msf_eiv_pac(:,jk) = v_msf_eiv_pac(:,jk-1) + v_msf_eiv_pac(:,jk)
368                  v_msf_eiv_ind(:,jk) = v_msf_eiv_ind(:,jk-1) + v_msf_eiv_ind(:,jk)
369                  v_msf_eiv_ipc(:,jk) = v_msf_eiv_ipc(:,jk-1) + v_msf_eiv_ipc(:,jk)
370               END DO
371            ENDIF
372#endif
373            !
374            IF( ln_subbas .AND. ln_diaznl ) THEN
375               DO jk = 2,jpk 
376                  v_msf_atl(:,jk) = v_msf_atl(:,jk-1) + v_msf_atl(:,jk)
377                  v_msf_pac(:,jk) = v_msf_pac(:,jk-1) + v_msf_pac(:,jk)
378                  v_msf_ind(:,jk) = v_msf_ind(:,jk-1) + v_msf_ind(:,jk)
379                  v_msf_ipc(:,jk) = v_msf_ipc(:,jk-1) + v_msf_ipc(:,jk)
380               END DO
381               v_msf_atl(:,:) = v_msf_atl(:,:) * zsverdrup
382               v_msf_pac(:,:) = v_msf_pac(:,:) * zsverdrup
383               v_msf_ind(:,:) = v_msf_ind(:,:) * zsverdrup
384               v_msf_ipc(:,:) = v_msf_ipc(:,:) * zsverdrup
385            ENDIF
386
387            ! Transports
388            ! T times V on T points (include bolus velocities)
389#if defined key_diaeiv 
390            DO jj = 1, jpj
391               DO ji = 1, jpi
392                  vt(ji,jj,:) = tn(ji,jj,:) * ( vn(ji,jj,:) + vn(ji,jj-1,:) + u_eiv(ji,jj,:) + u_eiv(ji,jj-1,:) )*0.5
393                  vs(ji,jj,:) = sn(ji,jj,:) * ( vn(ji,jj,:) + vn(ji,jj-1,:) + v_eiv(ji,jj,:) + v_eiv(ji,jj-1,:) )*0.5
394               END DO
395            END DO
396#else
397            DO jj = 1, jpj
398               DO ji = 1, jpi
399                  vt(ji,jj,:) = tn(ji,jj,:) * ( vn(ji,jj,:) + vn(ji,jj-1,:) )*0.5
400                  vs(ji,jj,:) = sn(ji,jj,:) * ( vn(ji,jj,:) + vn(ji,jj-1,:)  )*0.5
401               END DO
402            END DO
403#endif
404
405            ht_glo(:) = SUM( ptr_vjk( vt(:,:,:)), 2 )
406            st_glo(:) = SUM( ptr_vjk( vs(:,:,:)), 2 )
407
408            IF ( ln_subbas ) THEN
409               ht_atl(:) = SUM( ptr_vjk( vt (:,:,:), abasin(:,:)*sbasin(:,:)), 2 )
410               ht_pac(:) = SUM( ptr_vjk( vt (:,:,:), pbasin(:,:)*sbasin(:,:)), 2 )
411               ht_ind(:) = SUM( ptr_vjk( vt (:,:,:), ibasin(:,:)*sbasin(:,:)), 2 )
412               ht_ipc(:) = SUM( ptr_vjk( vt (:,:,:), dbasin(:,:)*sbasin(:,:)), 2 )
413               st_atl(:) = SUM( ptr_vjk( vs (:,:,:), abasin(:,:)*sbasin(:,:)), 2 )
414               st_pac(:) = SUM( ptr_vjk( vs (:,:,:), pbasin(:,:)*sbasin(:,:)), 2 )
415               st_ind(:) = SUM( ptr_vjk( vs (:,:,:), ibasin(:,:)*sbasin(:,:)), 2 )
416               st_ipc(:) = SUM( ptr_vjk( vs (:,:,:), dbasin(:,:)*sbasin(:,:)), 2 )
417            ENDIF
418
419            ! poleward tracer transports:
420            ! overturning components:
421            IF ( ln_ptrcomp ) THEN
422               pht_ove_glo(:) = SUM( v_msf_glo(:,:) * tn_jk_glo(:,:), 2 )   ! SUM over jk
423               pst_ove_glo(:) = SUM( v_msf_glo(:,:) * sn_jk_glo(:,:), 2 ) 
424               IF ( ln_subbas ) THEN
425                  pht_ove_atl(:) = SUM( v_msf_atl(:,:) * tn_jk_atl(:,:), 2 )   ! SUM over jk
426                  pst_ove_atl(:) = SUM( v_msf_atl(:,:) * sn_jk_atl(:,:), 2 ) 
427                  pht_ove_pac(:) = SUM( v_msf_pac(:,:) * tn_jk_pac(:,:), 2 )   ! SUM over jk
428                  pst_ove_pac(:) = SUM( v_msf_pac(:,:) * sn_jk_pac(:,:), 2 ) 
429                  pht_ove_ind(:) = SUM( v_msf_ind(:,:) * tn_jk_ind(:,:), 2 )   ! SUM over jk
430                  pst_ove_ind(:) = SUM( v_msf_ind(:,:) * sn_jk_ind(:,:), 2 ) 
431                  pht_ove_ipc(:) = SUM( v_msf_ipc(:,:) * tn_jk_ipc(:,:), 2 )   ! SUM over jk
432                  pst_ove_ipc(:) = SUM( v_msf_ipc(:,:) * sn_jk_ipc(:,:), 2 ) 
433               END IF
434            END IF
435
436            ! Bolus component
437#if defined key_diaeiv
438            pht_eiv_glo(:) = SUM( v_msf_eiv_glo(:,:) * tn_jk_glo(:,:), 2 )   ! SUM over jk
439            pst_eiv_glo(:) = SUM( v_msf_eiv_glo(:,:) * sn_jk_glo(:,:), 2 )   ! SUM over jk
440            IF ( ln_subbas ) THEN
441               pht_eiv_atl(:) = SUM( v_msf_eiv_glo(:,:) * tn_jk_atl(:,:), 2 )   ! SUM over jk
442               pst_eiv_atl(:) = SUM( v_msf_eiv_glo(:,:) * sn_jk_atl(:,:), 2 )   ! SUM over jk
443               pht_eiv_pac(:) = SUM( v_msf_eiv_pac(:,:) * tn_jk_pac(:,:), 2 )   ! SUM over jk
444               pst_eiv_pac(:) = SUM( v_msf_eiv_pac(:,:) * sn_jk_pac(:,:), 2 )   ! SUM over jk
445               pht_eiv_ind(:) = SUM( v_msf_eiv_ind(:,:) * tn_jk_ind(:,:), 2 )   ! SUM over jk
446               pst_eiv_ind(:) = SUM( v_msf_eiv_ind(:,:) * sn_jk_ind(:,:), 2 )   ! SUM over jk
447               pht_eiv_ipc(:) = SUM( v_msf_eiv_ipc(:,:) * tn_jk_ipc(:,:), 2 )   ! SUM over jk
448               pst_eiv_ipc(:) = SUM( v_msf_eiv_ipc(:,:) * sn_jk_ipc(:,:), 2 )   ! SUM over jk
449            ENDIF
450#endif
451
452            ! conversion in PW and G g
453            zpwatt = zpwatt * rau0 * rcp
454            pht_adv(:) = pht_adv(:) * zpwatt 
455            pht_ldf(:) = pht_ldf(:) * zpwatt
456            pst_adv(:) = pst_adv(:) * zggram
457            pst_ldf(:) = pst_ldf(:) * zggram
458            IF ( ln_ptrcomp ) THEN
459               pht_ove_glo(:) = pht_ove_glo(:) * zpwatt
460               pst_ove_glo(:) = pst_ove_glo(:) * zggram
461            END IF
462#if defined key_diaeiv
463            pht_eiv_glo(:) = pht_eiv_glo(:) * zpwatt
464            pst_eiv_glo(:) = pst_eiv_glo(:) * zggram
465#endif
466            IF( ln_subbas ) THEN
467               ht_atl(:) = ht_atl(:) * zpwatt
468               ht_pac(:) = ht_pac(:) * zpwatt
469               ht_ind(:) = ht_ind(:) * zpwatt
470               ht_ipc(:) = ht_ipc(:) * zpwatt
471               st_atl(:) = st_atl(:) * zggram 
472               st_pac(:) = st_pac(:) * zggram
473               st_ind(:) = st_ind(:) * zggram
474               st_ipc(:) = st_ipc(:) * zggram
475            ENDIF
476         ENDIF
477
478         ! outputs
479         CALL dia_ptr_wri( kt )
480
481      ENDIF
482
483      ! Close the file
484      IF( kt == nitend ) CALL histclo( numptr )
485      !
486   END SUBROUTINE dia_ptr
487
488
489   SUBROUTINE dia_ptr_init
490      !!----------------------------------------------------------------------
491      !!                  ***  ROUTINE dia_ptr_init  ***
492      !!                   
493      !! ** Purpose :   Initialization, namelist read
494      !!----------------------------------------------------------------------
495      NAMELIST/namptr/ ln_diaptr, ln_diaznl, ln_subbas, ln_ptrcomp, nf_ptr, nf_ptr_wri
496      INTEGER :: inum       ! temporary logical unit
497#if defined key_mpp_mpi
498      INTEGER, DIMENSION (1) :: iglo, iloc, iabsf, iabsl, ihals, ihale, idid
499#endif
500      !!----------------------------------------------------------------------
501
502      ! Read Namelist namptr : poleward transport parameters
503      REWIND ( numnam )
504      READ   ( numnam, namptr )
505
506      ! Control print
507      IF(lwp) THEN
508         WRITE(numout,*)
509         WRITE(numout,*) 'dia_ptr_init : poleward transport and msf initialization'
510         WRITE(numout,*) '~~~~~~~~~~~~'
511         WRITE(numout,*) '          Namelist namptr : set ptr parameters'
512         WRITE(numout,*) '             Switch for ptr diagnostic (T) or not (F) ln_diaptr = ', ln_diaptr
513         WRITE(numout,*) '             Atla/Paci/Ind basins computation         ln_subbas = ', ln_subbas
514         WRITE(numout,*) '             Frequency of computation                    nf_ptr = ', nf_ptr
515         WRITE(numout,*) '             Frequency of outputs                    nf_ptr_wri = ', nf_ptr_wri
516      ENDIF
517
518      !
519      ! Define MPI communicator for zonal sum
520      !
521      IF( lk_mpp )  THEN
522         CALL mpp_ini_znl
523      ENDIF
524
525      IF( ln_subbas ) THEN                ! load sub-basin mask
526         CALL iom_open( 'subbasins', inum )
527         CALL iom_get( inum, jpdom_data, 'atlmsk', abasin )      ! Atlantic basin
528         CALL iom_get( inum, jpdom_data, 'pacmsk', pbasin )      ! Pacific basin
529         CALL iom_get( inum, jpdom_data, 'indmsk', ibasin )      ! Indian basin
530         CALL iom_close( inum )
531         dbasin(:,:) = MAX ( pbasin(:,:), ibasin(:,:) )
532         sbasin(:,:) = tmask (:,:,1)
533         WHERE ( gphit (:,:) < -30.e0) sbasin(:,:) = 0.e0
534      ENDIF
535     
536      ! inverse of the ocean "zonal" v-point section
537      surf_jk_glo(:,:) = ptr_tjk( tmask(:,:,:) )
538      surf_jk_r_glo(:,:) = 0.e0
539      WHERE( surf_jk_glo(:,:) /= 0.e0 )   surf_jk_r_glo(:,:) = 1.e0 / surf_jk_glo(:,:)
540     
541      IF (ln_subbas) THEN
542         surf_jk_atl(:,:) = ptr_tjk( tmask (:,:,:), abasin(:,:) )
543         surf_jk_r_atl(:,:) = 0.e0
544         WHERE( surf_jk_atl(:,:) /= 0.e0 )   surf_jk_r_atl(:,:) = 1.e0 / surf_jk_atl(:,:)
545         !
546         surf_jk_pac(:,:) = ptr_tjk( tmask (:,:,:), pbasin(:,:) )
547         surf_jk_r_pac(:,:) = 0.e0
548         WHERE( surf_jk_pac(:,:) /= 0.e0 )   surf_jk_r_pac(:,:) = 1.e0 / surf_jk_pac(:,:)
549         !
550         surf_jk_ind(:,:) = ptr_tjk( tmask (:,:,:), ibasin(:,:) )
551         surf_jk_r_ind(:,:) = 0.e0
552         WHERE( surf_jk_ind(:,:) /= 0.e0 )   surf_jk_r_ind(:,:) = 1.e0 / surf_jk_ind(:,:)
553         !
554         surf_jk_ipc(:,:) = ptr_tjk( tmask (:,:,:), dbasin(:,:) )
555         surf_jk_r_ipc(:,:) = 0.e0
556         WHERE( surf_jk_ipc(:,:) /= 0.e0 )   surf_jk_r_ipc(:,:) = 1.e0 / surf_jk_ipc(:,:)
557      END IF
558
559     
560      !!----------------------------------------------------------------------
561
562#if defined key_mpp_mpi 
563      iglo (1) = jpjglo
564      iloc (1) = nlcj
565      iabsf(1) = njmppt(narea)
566      iabsl(:) = iabsf(:) + iloc(:) - 1
567      ihals(1) = nldj - 1
568      ihale(1) = nlcj - nlej
569      idid (1) = 2
570
571!-$$      IF(lwp) THEN
572!-$$          WRITE(numout,*)
573!-$$          WRITE(numout,*) 'dia_ptr_init :   iloc  = ', iloc
574!-$$          WRITE(numout,*) '~~~~~~~~~~~~     iabsf = ', iabsf
575!-$$          WRITE(numout,*) '                 ihals = ', ihals
576!-$$          WRITE(numout,*) '                 ihale = ', ihale
577!-$$      ENDIF
578
579      CALL flio_dom_set ( jpnj, nproc/jpni, idid, iglo, iloc, iabsf, iabsl, ihals, ihale, 'BOX', nidom_ptr)
580#else
581      nidom_ptr = FLIO_DOM_NONE
582#endif
583     
584   END SUBROUTINE dia_ptr_init
585
586
587   SUBROUTINE dia_ptr_wri( kt )
588      !!---------------------------------------------------------------------
589      !!                ***  ROUTINE dia_ptr_wri  ***
590      !!
591      !! ** Purpose :   output of poleward fluxes
592      !!
593      !! ** Method  :   NetCDF file
594      !!----------------------------------------------------------------------
595      INTEGER, INTENT(in) ::   kt   ! ocean time-step index
596      !!
597      INTEGER, SAVE ::   nhoridz, ndepidzt, ndepidzw
598      INTEGER, SAVE :: ndim  , ndim_atl     , ndim_pac     , ndim_ind     , ndim_ipc
599      INTEGER, SAVE ::         ndim_atl_30  , ndim_pac_30  , ndim_ind_30  , ndim_ipc_30
600      INTEGER, SAVE :: ndim_h, ndim_h_atl_30, ndim_h_pac_30, ndim_h_ind_30, ndim_h_ipc_30
601      INTEGER, SAVE, DIMENSION (jpj*jpk) :: ndex  , ndex_atl     , ndex_pac     , ndex_ind     , ndex_ipc
602      INTEGER, SAVE, DIMENSION (jpj*jpk) ::         ndex_atl_30  , ndex_pac_30  , ndex_ind_30  , ndex_ipc_30
603      INTEGER, SAVE, DIMENSION (jpj)     :: ndex_h, ndex_h_atl_30, ndex_h_pac_30, ndex_h_ind_30, ndex_h_ipc_30
604
605      CHARACTER (len=40)       ::   clhstnam, clop, clop_once, cl_comment   ! temporary names
606      INTEGER                  ::   iline, it, itmod, ji, jj, jk      !
607      REAL(wp)                 ::   zsto, zout, zdt, zjulian   ! temporary scalars
608      REAL(wp), DIMENSION(jpj) ::   zphi, zfoo
609      REAL(wp), DIMENSION(jpj,jpk) :: z_1
610      !!----------------------------------------------------------------------
611
612      ! define time axis
613      it    = kt / nf_ptr
614      itmod = kt - nit000 + 1
615     
616!-$$      IF(lwp) THEN
617!-$$         WRITE(numout,*)
618!-$$         WRITE(numout,*) 'dia_ptr_wri : kt = ', kt, 'it = ', it, ' itmod = ', itmod, ' niter = ', niter
619!-$$         WRITE(numout,*) '~~~~~~~~~~~~'
620!-$$      ENDIF
621
622      ! Initialization
623      ! --------------
624      IF( kt == nit000 ) THEN
625
626         niter = (nit000 - 1) / nf_ptr
627
628!-$$         IF(lwp) THEN
629!-$$            WRITE(numout,*)
630!-$$            WRITE(numout,*) 'dia_ptr_wri : poleward transport and msf writing: initialization , niter = ', niter
631!-$$            WRITE(numout,*) '~~~~~~~~~~~~'
632!-$$         ENDIF
633
634         zdt = rdt
635         IF( nacc == 1 ) zdt = rdtmin
636
637         ! Reference latitude
638         ! ------------------
639         !                                           ! =======================
640         IF( cp_cfg == "orca" ) THEN                 !   ORCA configurations
641            !                                        ! =======================
642
643            IF( jp_cfg == 05  )   iline = 192   ! i-line that passes near the North Pole
644            IF( jp_cfg == 025 )   iline = 384   ! i-line that passes near the North Pole
645            IF( jp_cfg == 1   )   iline =  96   ! i-line that passes near the North Pole
646            IF( jp_cfg == 2   )   iline =  48   ! i-line that passes near the North Pole
647            IF( jp_cfg == 4   )   iline =  24   ! i-line that passes near the North Pole
648            zphi(:) = 0.e0
649            DO ji = mi0(iline), mi1(iline) 
650               zphi(:) = gphiv(ji,:)         ! if iline is in the local domain
651               ! Correct highest latitude for some configurations - will work if domain is parallelized in J ?
652               IF( jp_cfg == 05 ) THEN
653                  DO jj = mj0(jpjdta), mj1(jpjdta) 
654                     zphi( jj ) = zphi(mj0(jpjdta-1)) + (zphi(mj0(jpjdta-1))-zphi(mj0(jpjdta-2)))/2.
655                     zphi( jj ) = MIN( zphi(jj), 90.)
656                  END DO
657               END IF
658               IF( jp_cfg == 1 .OR. jp_cfg == 2 .OR. jp_cfg == 4 ) THEN
659                  DO jj = mj0(jpjdta-1), mj1(jpjdta-1) 
660                     zphi( jj ) = 88.5e0
661                  END DO
662                  DO jj = mj0(jpjdta  ), mj1(jpjdta  ) 
663                     zphi( jj ) = 89.5e0
664                  END DO
665               END IF
666            END DO
667            ! provide the correct zphi to all local domains
668#if defined key_mpp_mpi
669            CALL mpp_sum( zphi, jpj, ncomm_znl )       
670#endif
671
672            !                                        ! =======================
673         ELSE                                        !   OTHER configurations  zjulian = zjulian - adatrj
674            !   set calendar origin to the beginning of the experiment
675            !                                        ! =======================
676            zphi(:) = gphiv(1,:)             ! assume lat/lon coordinate, select the first i-line
677            !
678         ENDIF
679         !
680         ! Work only on westmost processor (will not work if mppini2 is used)
681#if defined key_mpp_mpi
682         IF ( l_znl_root ) THEN 
683#endif
684            !
685            ! OPEN netcdf file
686            ! ----------------
687            ! Define frequency of output and means
688            zsto = nf_ptr * zdt
689            IF( ln_mskland )   THEN    ! put 1.e+20 on land (very expensive!!)
690               clop      = "ave(only(x))"
691               clop_once = "once(only(x))"
692            ELSE                       ! no use of the mask value (require less cpu time)
693               clop      = "ave(x)"       
694               clop_once = "once"
695            ENDIF
696
697            zout = nf_ptr_wri * zdt
698            zfoo(:) = 0.e0
699
700            ! Compute julian date from starting date of the run
701
702            CALL ymds2ju( nyear, nmonth, nday, rdt, zjulian )
703            zjulian = zjulian - adatrj   !   set calendar origin to the beginning of the experiment
704
705            CALL dia_nam( clhstnam, nf_ptr_wri, 'diaptr' )
706            IF(lwp)WRITE( numout,*)" Name of diaptr NETCDF file : ", clhstnam
707
708            ! Horizontal grid : zphi()
709            CALL histbeg(clhstnam, 1, zfoo, jpj, zphi,   &
710               1, 1, 1, jpj, niter, zjulian, zdt*nf_ptr, nhoridz, numptr, domain_id=nidom_ptr)
711            ! Vertical grids : gdept_0, gdepw_0
712            CALL histvert( numptr, "deptht", "Vertical T levels",   &
713               "m", jpk, gdept_0, ndepidzt, "down" )
714            CALL histvert( numptr, "depthw", "Vertical W levels",   &
715               "m", jpk, gdepw_0, ndepidzw, "down" )
716
717            !
718            CALL wheneq ( jpj*jpk, MIN(surf_jk_glo(:,:), 1.e0), 1, 1., ndex  , ndim  )      ! Lat-Depth
719            CALL wheneq ( jpj    , MIN(surf_jk_glo(:,1), 1.e0), 1, 1., ndex_h, ndim_h )     ! Lat
720
721            IF (ln_subbas) THEN
722               z_1 (:,1) = 1.0e0
723               WHERE ( gphit (jpi/2,:) .LT. -30 ) z_1 (:,1) = 0.e0
724               DO jk = 2, jpk
725                  z_1 (:,jk) = z_1 (:,1)
726               END DO
727
728               CALL wheneq ( jpj*jpk, MIN(surf_jk_atl(:,:)         , 1.e0), 1, 1., ndex_atl     , ndim_atl      ) ! Lat-Depth
729               CALL wheneq ( jpj*jpk, MIN(surf_jk_atl(:,:)*z_1(:,:), 1.e0), 1, 1., ndex_atl_30  , ndim_atl_30   ) ! Lat-Depth
730               CALL wheneq ( jpj    , MIN(surf_jk_atl(:,1)*z_1(:,1), 1.e0), 1, 1., ndex_h_atl_30, ndim_h_atl_30 ) ! Lat
731
732               CALL wheneq ( jpj*jpk, MIN(surf_jk_pac(:,:)         , 1.e0), 1, 1., ndex_pac     , ndim_pac      ) ! Lat-Depth
733               CALL wheneq ( jpj*jpk, MIN(surf_jk_pac(:,:)*z_1(:,:), 1.e0), 1, 1., ndex_pac_30  , ndim_pac_30   ) ! Lat-Depth
734               CALL wheneq ( jpj    , MIN(surf_jk_pac(:,1)*z_1(:,1), 1.e0), 1, 1., ndex_h_pac_30, ndim_h_pac_30 ) ! Lat
735
736               CALL wheneq ( jpj*jpk, MIN(surf_jk_ind(:,:)         , 1.e0), 1, 1., ndex_ind     , ndim_ind      ) ! Lat-Depth
737               CALL wheneq ( jpj*jpk, MIN(surf_jk_ind(:,:)*z_1(:,:), 1.e0), 1, 1., ndex_ind_30  , ndim_ind_30   ) ! Lat-Depth
738               CALL wheneq ( jpj    , MIN(surf_jk_ind(:,1)*z_1(:,1), 1.e0), 1, 1., ndex_h_ind_30, ndim_h_ind_30 ) ! Lat
739
740               CALL wheneq ( jpj*jpk, MIN(surf_jk_ipc(:,:)         , 1.e0), 1, 1., ndex_ipc     , ndim_ipc      ) ! Lat-Depth
741               CALL wheneq ( jpj*jpk, MIN(surf_jk_ipc(:,:)*z_1(:,:), 1.e0), 1, 1., ndex_ipc_30  , ndim_ipc_30   ) ! Lat-Depth
742               CALL wheneq ( jpj    , MIN(surf_jk_ipc(:,1)*z_1(:,1), 1.e0), 1, 1., ndex_h_ipc_30, ndim_h_ipc_30 ) ! Lat
743
744            ENDIF
745
746            !
747#if defined key_diaeiv
748            cl_comment = ' (Bolus part included)'
749#else
750            cl_comment = '                      '
751#endif
752            !  Zonal mean T and S
753
754            IF ( ln_diaznl ) THEN
755               CALL histdef( numptr, "zotemglo", "Zonal Mean Temperature","C" ,   &
756                  1, jpj, nhoridz, jpk, 1, jpk, ndepidzt, 32, clop, zsto, zout )
757               CALL histdef( numptr, "zosalglo", "Zonal Mean Salinity","PSU"  ,   &
758                  1, jpj, nhoridz, jpk, 1, jpk, ndepidzt, 32, clop, zsto, zout )
759
760               CALL histdef( numptr, "zosrfglo", "Zonal Mean Surface","m^2"   ,   &
761                  1, jpj, nhoridz, jpk, 1, jpk, ndepidzt, 32, clop_once, zsto, zout )
762
763               IF (ln_subbas) THEN
764                  CALL histdef( numptr, "zotematl", "Zonal Mean Temperature: Atlantic","C" ,   &
765                     1, jpj, nhoridz, jpk, 1, jpk, ndepidzt, 32, clop, zsto, zout )
766                  CALL histdef( numptr, "zosalatl", "Zonal Mean Salinity: Atlantic","PSU"  ,   &
767                     1, jpj, nhoridz, jpk, 1, jpk, ndepidzt, 32, clop, zsto, zout )
768                  CALL histdef( numptr, "zosrfatl", "Zonal Mean Surface: Atlantic","m^2"   ,   &
769                     1, jpj, nhoridz, jpk, 1, jpk, ndepidzt, 32, clop_once, zsto, zout )
770
771                  CALL histdef( numptr, "zotempac", "Zonal Mean Temperature: Pacific","C"  ,   &
772                     1, jpj, nhoridz, jpk, 1, jpk, ndepidzt, 32, clop, zsto, zout )
773                  CALL histdef( numptr, "zosalpac", "Zonal Mean Salinity: Pacific","PSU"   ,   &
774                     1, jpj, nhoridz, jpk, 1, jpk, ndepidzt, 32, clop, zsto, zout )
775                  CALL histdef( numptr, "zosrfpac", "Zonal Mean Surface: Pacific","m^2"    ,   &
776                     1, jpj, nhoridz, jpk, 1, jpk, ndepidzt, 32, clop_once, zsto, zout )
777
778                  CALL histdef( numptr, "zotemind", "Zonal Mean Temperature: Indian","C"   ,   &
779                     1, jpj, nhoridz, jpk, 1, jpk, ndepidzt, 32, clop, zsto, zout )
780                  CALL histdef( numptr, "zosalind", "Zonal Mean Salinity: Indian","PSU"    ,   &
781                     1, jpj, nhoridz, jpk, 1, jpk, ndepidzt, 32, clop, zsto, zout )
782                  CALL histdef( numptr, "zosrfind", "Zonal Mean Surface: Indian","m^2"     ,   &
783                     1, jpj, nhoridz, jpk, 1, jpk, ndepidzt, 32, clop_once, zsto, zout )
784
785                  CALL histdef( numptr, "zotemipc", "Zonal Mean Temperature: Pacific+Indian","C" ,   &
786                     1, jpj, nhoridz, jpk, 1, jpk, ndepidzt, 32, clop, zsto, zout )
787                  CALL histdef( numptr, "zosalipc", "Zonal Mean Salinity: Pacific+Indian","PSU"  ,   &
788                     1, jpj, nhoridz, jpk, 1, jpk, ndepidzt, 32, clop, zsto, zout )
789                  CALL histdef( numptr, "zosrfipc", "Zonal Mean Surface: Pacific+Indian","m^2"   ,   &
790                     1, jpj, nhoridz, jpk, 1, jpk, ndepidzt, 32, clop_once, zsto, zout )
791               ENDIF
792
793            ENDIF
794
795            !  Meridional Stream-Function (Eulerian and Bolus)
796
797            CALL histdef( numptr, "zomsfglo", "Meridional Stream-Function: Global"//TRIM(cl_comment),"Sv" ,   &
798               1, jpj, nhoridz, jpk, 1, jpk, ndepidzw, 32, clop, zsto, zout )
799            IF( ln_subbas .AND. ln_diaznl ) THEN
800               CALL histdef( numptr, "zomsfatl", "Meridional Stream-Function: Atlantic"//TRIM(cl_comment),"Sv" ,   &
801                  1, jpj, nhoridz, jpk, 1, jpk, ndepidzw, 32, clop, zsto, zout )
802               CALL histdef( numptr, "zomsfpac", "Meridional Stream-Function: Pacific"//TRIM(cl_comment),"Sv"  ,   &
803                  1, jpj, nhoridz, jpk, 1, jpk, ndepidzw, 32, clop, zsto, zout )
804               CALL histdef( numptr, "zomsfind", "Meridional Stream-Function: Indian"//TRIM(cl_comment),"Sv"   ,   &
805                  1, jpj, nhoridz, jpk, 1, jpk, ndepidzw, 32, clop, zsto, zout )
806               CALL histdef( numptr, "zomsfipc", "Meridional Stream-Function: Indo-Pacific"//TRIM(cl_comment),"Sv" ,&
807                  1, jpj, nhoridz, jpk, 1, jpk, ndepidzw, 32, clop, zsto, zout )
808            ENDIF
809
810            !  Heat transport
811
812            CALL histdef( numptr, "sophtadv", "Advective Heat Transport"      ,   &
813               "PW", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
814            CALL histdef( numptr, "sophtldf", "Diffusive Heat Transport"      ,   &
815               "PW",1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
816            IF ( ln_ptrcomp ) THEN
817               CALL histdef( numptr, "sophtove", "Overturning Heat Transport"    ,   &
818                  "PW",1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
819            END IF
820            IF( ln_subbas ) THEN
821               CALL histdef( numptr, "sohtatl", "Heat Transport Atlantic"//TRIM(cl_comment),  &
822                  "PW", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
823               CALL histdef( numptr, "sohtpac", "Heat Transport Pacific"//TRIM(cl_comment) ,  &
824                  "PW", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
825               CALL histdef( numptr, "sohtind", "Heat Transport Indian"//TRIM(cl_comment)  ,  &
826                  "PW", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
827               CALL histdef( numptr, "sohtipc", "Heat Transport Pacific+Indian"//TRIM(cl_comment), &
828                  "PW", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
829            ENDIF
830
831
832            !  Salt transport
833
834            CALL histdef( numptr, "sopstadv", "Advective Salt Transport"      ,   &
835               "Giga g/s", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
836            CALL histdef( numptr, "sopstldf", "Diffusive Salt Transport"      ,   &
837               "Giga g/s", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
838            IF ( ln_ptrcomp ) THEN
839               CALL histdef( numptr, "sopstove", "Overturning Salt Transport"    ,   &
840                  "Giga g/s", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
841            END IF
842#if defined key_diaeiv
843            ! Eddy induced velocity
844            CALL histdef( numptr, "zomsfeiv", "Bolus Meridional Stream-Function: global",   &
845               "Sv"      , 1, jpj, nhoridz, jpk, 1, jpk, ndepidzw, 32, clop, zsto, zout )
846            CALL histdef( numptr, "sophteiv", "Bolus Advective Heat Transport",   &
847               "PW"      , 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
848            CALL histdef( numptr, "sopsteiv", "Bolus Advective Salt Transport",   &
849               "Giga g/s", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
850#endif
851            IF( ln_subbas ) THEN
852               CALL histdef( numptr, "sostatl", "Salt Transport Atlantic"//TRIM(cl_comment)      ,  &
853                  "Giga g/s", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
854               CALL histdef( numptr, "sostpac", "Salt Transport Pacific"//TRIM(cl_comment)      ,   &
855                  "Giga g/s", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
856               CALL histdef( numptr, "sostind", "Salt Transport Indian"//TRIM(cl_comment)      ,    &
857                  "Giga g/s", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
858               CALL histdef( numptr, "sostipc", "Salt Transport Pacific+Indian"//TRIM(cl_comment),  &
859                  "Giga g/s", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
860            ENDIF
861
862            CALL histend( numptr )
863
864         END IF
865#if defined key_mpp_mpi
866      END IF
867#endif
868
869#if defined key_mpp_mpi
870      IF( MOD( itmod, nf_ptr ) == 0 .AND. l_znl_root ) THEN
871#else
872      IF( MOD( itmod, nf_ptr ) == 0  ) THEN
873#endif
874         niter = niter + 1
875
876!-$$         IF(lwp) THEN
877!-$$            WRITE(numout,*)
878!-$$            WRITE(numout,*) 'dia_ptr_wri : write Poleward Transports at time-step : kt = ', kt, &
879!-$$               & 'it = ', it, ' itmod = ', itmod, ' niter = ', niter
880!-$$            WRITE(numout,*) '~~~~~~~~~~'
881!-$$            WRITE(numout,*)
882!-$$         ENDIF
883
884         IF (ln_diaznl ) THEN
885            CALL histwrite( numptr, "zosrfglo", niter, surf_jk_glo , ndim, ndex )
886            CALL histwrite( numptr, "zotemglo", niter, tn_jk_glo  , ndim, ndex )
887            CALL histwrite( numptr, "zosalglo", niter, sn_jk_glo  , ndim, ndex )
888
889            IF (ln_subbas) THEN
890               CALL histwrite( numptr, "zosrfatl", niter, surf_jk_atl, ndim_atl, ndex_atl )
891               CALL histwrite( numptr, "zosrfpac", niter, surf_jk_pac, ndim_pac, ndex_pac )
892               CALL histwrite( numptr, "zosrfind", niter, surf_jk_ind, ndim_ind, ndex_ind )
893               CALL histwrite( numptr, "zosrfipc", niter, surf_jk_ipc, ndim_ipc, ndex_ipc )
894
895               CALL histwrite( numptr, "zotematl", niter, tn_jk_atl  , ndim_atl, ndex_atl )
896               CALL histwrite( numptr, "zosalatl", niter, sn_jk_atl  , ndim_atl, ndex_atl )
897               CALL histwrite( numptr, "zotempac", niter, tn_jk_pac  , ndim_pac, ndex_pac )
898               CALL histwrite( numptr, "zosalpac", niter, sn_jk_pac  , ndim_pac, ndex_pac )
899               CALL histwrite( numptr, "zotemind", niter, tn_jk_ind  , ndim_ind, ndex_ind )
900               CALL histwrite( numptr, "zosalind", niter, sn_jk_ind  , ndim_ind, ndex_ind )
901               CALL histwrite( numptr, "zotemipc", niter, tn_jk_ipc  , ndim_ipc, ndex_ipc )
902               CALL histwrite( numptr, "zosalipc", niter, sn_jk_ipc  , ndim_ipc, ndex_ipc )
903            END IF
904         ENDIF
905
906         ! overturning outputs:
907         CALL histwrite( numptr, "zomsfglo", niter, v_msf_glo, ndim, ndex )
908         IF( ln_subbas .AND. ln_diaznl ) THEN
909            CALL histwrite( numptr, "zomsfatl", niter, v_msf_atl , ndim_atl_30, ndex_atl_30 )
910            CALL histwrite( numptr, "zomsfpac", niter, v_msf_pac , ndim_pac_30, ndex_pac_30 )
911            CALL histwrite( numptr, "zomsfind", niter, v_msf_ind , ndim_ind_30, ndex_ind_30 )
912            CALL histwrite( numptr, "zomsfipc", niter, v_msf_ipc , ndim_ipc_30, ndex_ipc_30 )
913         ENDIF
914#if defined key_diaeiv
915         CALL histwrite( numptr, "zomsfeiv", niter, v_msf_eiv_glo, ndim  , ndex   )
916#endif
917
918
919         ! heat transport outputs:
920         IF( ln_subbas ) THEN
921            CALL histwrite( numptr, "sohtatl", niter, ht_atl  , ndim_h_atl_30, ndex_h_atl_30 )
922            CALL histwrite( numptr, "sohtpac", niter, ht_pac  , ndim_h_pac_30, ndex_h_pac_30 )
923            CALL histwrite( numptr, "sohtind", niter, ht_ind  , ndim_h_ind_30, ndex_h_ind_30 )
924            CALL histwrite( numptr, "sohtipc", niter, ht_ipc  , ndim_h_ipc_30, ndex_h_ipc_30 )
925            CALL histwrite( numptr, "sostatl", niter, st_atl  , ndim_h_atl_30, ndex_h_atl_30 )
926            CALL histwrite( numptr, "sostpac", niter, st_pac  , ndim_h_pac_30, ndex_h_pac_30 )
927            CALL histwrite( numptr, "sostind", niter, st_ind  , ndim_h_ind_30, ndex_h_ind_30 )
928            CALL histwrite( numptr, "sostipc", niter, st_ipc  , ndim_h_ipc_30, ndex_h_ipc_30 )
929         ENDIF
930
931         CALL histwrite( numptr, "sophtadv", niter, pht_adv     , ndim_h, ndex_h )
932         CALL histwrite( numptr, "sophtldf", niter, pht_ldf     , ndim_h, ndex_h )
933         CALL histwrite( numptr, "sopstadv", niter, pst_adv     , ndim_h, ndex_h )
934         CALL histwrite( numptr, "sopstldf", niter, pst_ldf     , ndim_h, ndex_h )
935         IF ( ln_ptrcomp ) THEN
936            CALL histwrite( numptr, "sopstove", niter, pst_ove_glo , ndim_h, ndex_h )
937            CALL histwrite( numptr, "sophtove", niter, pht_ove_glo , ndim_h, ndex_h )
938         ENDIF
939#if defined key_diaeiv
940         CALL histwrite( numptr, "sophteiv", niter, pht_eiv_glo  , ndim_h, ndex_h )
941         CALL histwrite( numptr, "sopsteiv", niter, pst_eiv_glo  , ndim_h, ndex_h )
942#endif
943         !
944      ENDIF
945      !
946   END SUBROUTINE dia_ptr_wri
947
948   !!======================================================================
949END MODULE diaptr
Note: See TracBrowser for help on using the repository browser.