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 @ 508

Last change on this file since 508 was 508, checked in by opalod, 18 years ago

nemo_v1_update_071:RB: add iom for restart and reorganization of restart
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 27.0 KB
Line 
1MODULE diaptr
2   !!======================================================================
3   !!                       ***  MODULE  diaptr  ***
4   !! Ocean physics:  brief description of the purpose of the module
5   !!                 (please no more than 2 lines)
6   !!=====================================================================
7   !! History :  9.0  !  03-09  (C. Talandir, G. Madec)  Original code
8   !!            9.0  !  06-01  (A. Biastoch)  Allow sub-basins computation
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_vtjk     : "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    ! ???
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   INTEGER , PUBLIC                 ::   nf_ptr = 15           !: frequency of ptr computation
47
48   REAL(wp), PUBLIC, DIMENSION(jpj) ::   pht_adv, pst_adv      !: heat and salt poleward transport: advection
49   REAL(wp), PUBLIC, DIMENSION(jpj) ::   pht_ove, pst_ove      !: heat and salt poleward transport: overturning
50   REAL(wp), PUBLIC, DIMENSION(jpj) ::   pht_ldf, pst_ldf      !: heat and salt poleward transport: lateral diffusion
51#if defined key_diaeiv
52   REAL(wp), PUBLIC, DIMENSION(jpj) ::   pht_eiv, pst_eiv      !: heat and salt poleward transport: bolus advection
53#endif
54   REAL(wp), PUBLIC, DIMENSION(jpj) ::   ht_atl,ht_ind,ht_pac  !: heat
55   REAL(wp), PUBLIC, DIMENSION(jpj) ::   st_atl,st_ind,st_pac  !: salt
56
57     
58
59   REAL(wp), DIMENSION(jpj,jpk) ::   tn_jk  , sn_jk  ,  &  !: "zonal" mean temperature and salinity
60      &                              v_msf_atl       ,  &  !: "meridional" Stream-Function
61      &                              v_msf_glo       ,  &  !: "meridional" Stream-Function
62      &                              v_msf_ipc       ,  &  !: "meridional" Stream-Function
63      &                              surf_jk_r             !: inverse of the ocean "zonal" section surface
64#if defined key_diaeiv
65   REAL(wp), DIMENSION(jpj,jpk) ::   v_msf_eiv                  !: bolus "meridional" Stream-Function
66#endif
67   REAL(wp), DIMENSION(jpi,jpj) ::   abasin, pbasin, ibasin     !: return function value
68
69   !! * Substitutions
70#  include "domzgr_substitute.h90"
71#  include "vectopt_loop_substitute.h90"
72   !!----------------------------------------------------------------------
73   !!   OPA 9.0 , LOCEAN-IPSL (2005)
74   !! $Header$
75   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
76   !!----------------------------------------------------------------------
77
78CONTAINS
79
80   FUNCTION ptr_vj_3d( pva )   RESULT ( p_fval )
81      !!----------------------------------------------------------------------
82      !!                    ***  ROUTINE ptr_vj_3d  ***
83      !!
84      !! ** Purpose :   "zonal" and vertical sum computation of a "meridional"
85      !!      flux array
86      !!
87      !! ** Method  : - i-k sum of pva using the interior 2D vmask (vmask_i).
88      !!      pva is supposed to be a masked flux (i.e. * vmask*e1v*e3v)
89      !!
90      !! ** Action  : - p_fval: i-k-mean poleward flux of pva
91      !!----------------------------------------------------------------------
92      REAL(wp) , INTENT(in), DIMENSION(jpi,jpj,jpk) ::   pva   ! mask flux array at V-point
93      !!
94      INTEGER                  ::   ji, jj, jk   ! dummy loop arguments
95      INTEGER                  ::   ijpj         ! ???
96      REAL(wp), DIMENSION(jpj) ::   p_fval       ! function value
97      !!--------------------------------------------------------------------
98      !
99      ijpj = jpj
100      p_fval(:) = 0.e0
101      DO jk = 1, jpkm1
102         DO jj = 2, jpjm1
103            DO ji = fs_2, fs_jpim1   ! Vector opt.
104               p_fval(jj) = p_fval(jj) + pva(ji,jj,jk) * tmask_i(ji,jj+1) * tmask_i(ji,jj) 
105            END DO
106         END DO
107      END DO
108      !
109      IF( lk_mpp )   CALL mpp_sum( p_fval, ijpj )     !!bug  I presume
110      !
111   END FUNCTION ptr_vj_3d
112
113
114   FUNCTION ptr_vj_2d( pva )   RESULT ( p_fval )
115      !!----------------------------------------------------------------------
116      !!                    ***  ROUTINE ptr_vj_2d  ***
117      !!
118      !! ** Purpose :   "zonal" and vertical sum computation of a "meridional"
119      !!      flux array
120      !!
121      !! ** Method  : - i-k sum of pva using the interior 2D vmask (vmask_i).
122      !!      pva is supposed to be a masked flux (i.e. * vmask*e1v*e3v)
123      !!
124      !! ** Action  : - p_fval: i-k-mean poleward flux of pva
125      !!----------------------------------------------------------------------
126      REAL(wp) , INTENT(in), DIMENSION(jpi,jpj) ::   pva   ! mask flux array at V-point
127      !!
128      INTEGER                  ::   ji,jj    ! dummy loop arguments
129      INTEGER                  ::   ijpj     ! ???
130      REAL(wp), DIMENSION(jpj) ::   p_fval   ! function value
131      !!--------------------------------------------------------------------
132      !
133      ijpj = jpj
134      p_fval(:) = 0.e0
135      DO jj = 2, jpjm1
136         DO ji = fs_2, fs_jpim1   ! Vector opt.
137            p_fval(jj) = p_fval(jj) + pva(ji,jj) * tmask_i(ji,jj+1) * tmask_i(ji,jj)
138         END DO
139      END DO
140      !
141      IF( lk_mpp )   CALL mpp_sum( p_fval, ijpj )     !!bug  I presume
142      !
143   END FUNCTION ptr_vj_2d
144
145
146   FUNCTION ptr_vjk( pva )   RESULT ( p_fval )
147      !!----------------------------------------------------------------------
148      !!                    ***  ROUTINE ptr_vjk  ***
149      !!
150      !! ** Purpose :   "zonal" sum computation of a "meridional" flux array
151      !!
152      !! ** Method  : - i-sum of pva using the interior 2D vmask (vmask_i).
153      !!      pva is supposed to be a masked flux (i.e. * vmask*e1v*e3v)
154      !!
155      !! ** Action  : - p_fval: i-k-mean poleward flux of pva
156      !!----------------------------------------------------------------------
157      REAL(wp) , INTENT(in), DIMENSION(jpi,jpj,jpk) ::   pva   ! mask flux array at V-point
158      !!
159      INTEGER                      ::   ji, jj, jk   ! dummy loop arguments
160      INTEGER , DIMENSION (1)      ::   ish
161      INTEGER , DIMENSION (2)      ::   ish2
162      REAL(wp), DIMENSION(jpj*jpk) ::   zwork        ! temporary vector for mpp_sum
163      REAL(wp), DIMENSION(jpj,jpk) ::   p_fval       ! return function value
164      !!--------------------------------------------------------------------
165      !
166      p_fval(:,:) = 0.e0
167      !
168      DO jk = 1, jpkm1
169         DO jj = 2, jpjm1
170           DO ji = fs_2, fs_jpim1
171            p_fval(jj,jk) = p_fval(jj,jk) + pva(ji,jj,jk) * e1v(ji,jj) * fse3v(ji,jj,jk)   &
172               &                                          * tmask_i(ji,jj+1) * tmask_i(ji,jj)
173           END DO
174         END DO
175      END DO
176      !
177      IF(lk_mpp) THEN
178         ish(1) = jpj*jpk  ;  ish2(1) = jpj  ;  ish2(2) = jpk
179         zwork(:)= RESHAPE( p_fval, ish )
180         CALL mpp_sum( zwork, jpj*jpk )
181         p_fval(:,:)= RESHAPE( zwork, ish2 )
182      END IF
183      !
184   END FUNCTION ptr_vjk
185
186
187   FUNCTION ptr_vtjk( pva )   RESULT ( p_fval )
188      !!----------------------------------------------------------------------
189      !!                    ***  ROUTINE ptr_vtjk  ***
190      !!
191      !! ** Purpose :   "zonal" mean computation of a tracer field
192      !!
193      !! ** Method  : - i-sum of mj(pva) using the interior 2D vmask (vmask_i)
194      !!      multiplied by the inverse of the surface of the "zonal" ocean
195      !!      section
196      !!
197      !! ** Action  : - p_fval: i-k-mean poleward flux of pva
198      !!----------------------------------------------------------------------
199      REAL(wp) , INTENT(in), DIMENSION(jpi,jpj,jpk) ::   pva   ! mask flux array at V-point
200      !!
201      INTEGER                     ::   ji, jj, jk   ! dummy loop arguments
202      INTEGER, DIMENSION (1)      ::   ish
203      INTEGER, DIMENSION (2)      ::   ish2
204      REAL(wp),DIMENSION(jpj*jpk) ::   zwork        ! temporary vector for mpp_sum
205      REAL(wp),DIMENSION(jpj,jpk) ::   p_fval       ! return function value
206      !!--------------------------------------------------------------------
207      !
208      p_fval(:,:) = 0.e0
209      DO jk = 1, jpkm1
210         DO jj = 2, jpjm1
211            DO ji = fs_2, fs_jpim1   ! Vector opt.
212               p_fval(jj,jk) = p_fval(jj,jk) + ( pva(ji,jj,jk) + pva(ji,jj+1,jk) )              &
213                  &                          * e1v(ji,jj) * fse3v(ji,jj,jk) * vmask(ji,jj,jk)   &
214                  &                          * tmask_i(ji,jj+1) * tmask_i(ji,jj)
215            END DO
216         END DO
217      END DO
218      p_fval(:,:) = p_fval(:,:) * 0.5
219      IF(lk_mpp) THEN
220         ish(1) = jpj*jpk   ;   ish2(1) = jpj   ;   ish2(2) = jpk
221         zwork(:)= RESHAPE( p_fval, ish )
222         CALL mpp_sum( zwork, jpj*jpk )
223         p_fval(:,:)= RESHAPE(zwork,ish2)
224      END IF
225      !
226   END FUNCTION ptr_vtjk
227
228
229   SUBROUTINE dia_ptr( kt )
230      !!----------------------------------------------------------------------
231      !!                  ***  ROUTINE dia_ptr  ***
232      !!----------------------------------------------------------------------
233      INTEGER, INTENT(in) ::   kt   ! ocean time step index
234      !!
235      INTEGER  ::   jk, jj, ji               ! dummy loop
236      REAL(wp) ::   zsverdrup,  &              ! conversion from m3/s to Sverdrup
237         &          zpwatt,     &              ! conversion from W    to PW
238         &          zggram                     ! conversion from g    to Pg
239
240      REAL(wp), DIMENSION(jpi,jpj,jpk) ::  &
241         v_atl , v_ipc,                    &
242         vt_atl, vt_pac, vt_ind,           &
243         vs_atl, vs_pac, vs_ind,           &
244         zv_eiv
245      INTEGER ::  inum       ! temporary logical unit
246      !!----------------------------------------------------------------------
247
248      IF( kt == nit000 .OR. MOD( kt, nf_ptr ) == 0 )   THEN
249
250         zsverdrup = 1.e-6
251         zpwatt    = 1.e-15
252         zggram    = 1.e-6
253   
254# if defined key_diaeiv
255         zv_eiv(:,:,:) = v_eiv(:,:,:)
256# else
257         zv_eiv(:,:,:) = 0.e0
258# endif
259
260         ! "zonal" mean temperature and salinity at V-points
261         tn_jk(:,:) = ptr_vtjk( tn(:,:,:) ) * surf_jk_r(:,:)
262         sn_jk(:,:) = ptr_vtjk( sn(:,:,:) ) * surf_jk_r(:,:)
263
264         !--------------------------------------------------------
265         ! overturning calculation:
266 
267         IF( ln_subbas ) THEN              ! Basins computation
268
269            IF( kt == nit000 ) THEN                ! load sub-basin mask
270               CALL iom_open( 'subbasins', inum )
271               CALL iom_get( inum, jpdom_data, 'atlmsk', abasin )      ! Atlantic basin
272               CALL iom_get( inum, jpdom_data, 'pacmsk', pbasin )      ! Pacific basin
273               CALL iom_get( inum, jpdom_data, 'indmsk', ibasin )      ! Indian basin
274               CALL iom_close( inum )
275            ENDIF
276
277            ! basin separation:
278            DO jj = 1, jpj
279               DO ji = 1, jpi
280                  ! basin separated velocity
281                  v_atl(ji,jj,:) = (vn(ji,jj,:)+zv_eiv(ji,jj,:))*abasin(ji,jj)
282                  v_ipc(ji,jj,:) = (vn(ji,jj,:)+zv_eiv(ji,jj,:))*(pbasin(ji,jj)+ibasin(ji,jj))
283
284                  ! basin separated T times V on T points
285                  vt_ind(ji,jj,:) = tn(ji,jj,:) *                                 &
286                     &              ( (vn    (ji,jj,:) + vn    (ji,jj-1,:))*0.5   &
287                     &              + (zv_eiv(ji,jj,:) + zv_eiv(ji,jj-1,:))*0.5 ) 
288                  vt_atl(ji,jj,:) = vt_ind(ji,jj,:) * abasin(ji,jj)
289                  vt_pac(ji,jj,:) = vt_ind(ji,jj,:) * pbasin(ji,jj)
290                  vt_ind(ji,jj,:) = vt_ind(ji,jj,:) * ibasin(ji,jj)
291
292                  ! basin separated S times V on T points
293                  vs_ind(ji,jj,:) = sn(ji,jj,:) *                                 &
294                     &              ( (vn    (ji,jj,:) + vn    (ji,jj-1,:))*0.5   &
295                     &              + (zv_eiv(ji,jj,:) + zv_eiv(ji,jj-1,:))*0.5 ) 
296                  vs_atl(ji,jj,:) = vs_ind(ji,jj,:) * abasin(ji,jj)
297                  vs_pac(ji,jj,:) = vs_ind(ji,jj,:) * pbasin(ji,jj)
298                  vs_ind(ji,jj,:) = vs_ind(ji,jj,:) * ibasin(ji,jj)
299               END DO
300            END DO
301
302         ENDIF
303
304         ! horizontal integral and vertical dz
305         v_msf_glo(:,:) = ptr_vjk( vn(:,:,:) ) 
306#if defined key_diaeiv
307         v_msf_eiv(:,:) = ptr_vjk( v_eiv(:,:,:) ) 
308#endif
309         IF( ln_subbas ) THEN
310            v_msf_atl(:,:) = ptr_vjk( v_atl (:,:,:) ) 
311            v_msf_ipc(:,:) = ptr_vjk( v_ipc (:,:,:) ) 
312            ht_atl(:) = SUM( ptr_vjk( vt_atl(:,:,:)), 2 )
313            ht_pac(:) = SUM( ptr_vjk( vt_pac(:,:,:)), 2 )
314            ht_ind(:) = SUM( ptr_vjk( vt_ind(:,:,:)), 2 )
315            st_atl(:) = SUM( ptr_vjk( vs_atl(:,:,:)), 2 )
316            st_pac(:) = SUM( ptr_vjk( vs_pac(:,:,:)), 2 )
317            st_ind(:) = SUM( ptr_vjk( vs_ind(:,:,:)), 2 )
318         ENDIF
319
320         ! poleward tracer transports:
321         ! overturning components:
322         pht_ove(:) = SUM( v_msf_glo(:,:) * tn_jk(:,:), 2 )   ! SUM over jk
323         pst_ove(:) = SUM( v_msf_glo(:,:) * sn_jk(:,:), 2 )   ! SUM over jk
324#if defined key_diaeiv
325         pht_eiv(:) = SUM( v_msf_eiv(:,:) * tn_jk(:,:), 2 )   ! SUM over jk
326         pst_eiv(:) = SUM( v_msf_eiv(:,:) * sn_jk(:,:), 2 )   ! SUM over jk
327#endif
328     
329         ! conversion in PW and G g
330         zpwatt = zpwatt * rau0 * rcp
331         pht_adv(:) = pht_adv(:) * zpwatt 
332         pht_ove(:) = pht_ove(:) * zpwatt
333         pht_ldf(:) = pht_ldf(:) * zpwatt
334         pst_adv(:) = pst_adv(:) * zggram
335         pst_ove(:) = pst_ove(:) * zggram
336         pst_ldf(:) = pst_ldf(:) * zggram
337#if defined key_diaeiv
338         pht_eiv(:) = pht_eiv(:) * zpwatt
339         pst_eiv(:) = pst_eiv(:) * zggram
340#endif
341         IF( ln_subbas ) THEN
342            ht_atl(:) = ht_atl(:) * zpwatt
343            ht_pac(:) = ht_pac(:) * zpwatt
344            ht_ind(:) = ht_ind(:) * zpwatt
345            st_atl(:) = st_atl(:) * zggram 
346            st_pac(:) = st_pac(:) * zggram
347            st_ind(:) = st_ind(:) * zggram
348         ENDIF
349
350         ! "Meridional" Stream-Function
351         DO jk = 2,jpk 
352            v_msf_glo(:,jk) = v_msf_glo(:,jk-1) + v_msf_glo(:,jk)
353         END DO
354         v_msf_glo(:,:) = v_msf_glo(:,:) * zsverdrup
355
356#if defined key_diaeiv
357         ! Bolus "Meridional" Stream-Function
358         DO jk = 2,jpk 
359            v_msf_eiv(:,jk) = v_msf_eiv(:,jk-1) + v_msf_eiv(:,jk)
360         END DO
361         v_msf_eiv(:,:) = v_msf_eiv(:,:) * zsverdrup
362#endif
363
364         IF( ln_subbas ) THEN
365            DO jk = 2,jpk 
366               v_msf_atl(:,jk) = v_msf_atl(:,jk-1) + v_msf_atl(:,jk)
367               v_msf_ipc(:,jk) = v_msf_ipc(:,jk-1) + v_msf_ipc(:,jk)
368            END DO
369            v_msf_atl(:,:) = v_msf_atl(:,:) * zsverdrup
370            v_msf_ipc(:,:) = v_msf_ipc(:,:) * zsverdrup
371         ENDIF
372
373         ! outputs
374         CALL dia_ptr_wri( kt )
375
376      ENDIF
377
378      ! Close the file
379      IF( kt == nitend ) CALL histclo( numptr )
380      !
381   END SUBROUTINE dia_ptr
382
383
384   SUBROUTINE dia_ptr_init
385      !!----------------------------------------------------------------------
386      !!                  ***  ROUTINE dia_ptr_init  ***
387      !!                   
388      !! ** Purpose :   Initialization, namelist read
389      !!----------------------------------------------------------------------
390      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   z_1         ! temporary workspace
391
392      NAMELIST/namptr/ ln_diaptr, ln_subbas, nf_ptr
393      !!----------------------------------------------------------------------
394
395      ! Read Namelist namptr : poleward transport parameters
396      REWIND ( numnam )
397      READ   ( numnam, namptr )
398
399      ! Control print
400      IF(lwp) THEN
401         WRITE(numout,*)
402         WRITE(numout,*) 'dia_ptr_init : poleward transport and msf initialization'
403         WRITE(numout,*) '~~~~~~~~~~~~'
404         WRITE(numout,*) '          Namelist namptr : set ptr parameters'
405         WRITE(numout,*) '             Switch for ptr diagnostic (T) or not (F) ln_diaptr = ', ln_diaptr
406         WRITE(numout,*) '             Atla/Paci/Ind basins computation         ln_subbas = ', ln_subbas
407         WRITE(numout,*) '             Frequency of computation                    nf_ptr = ', nf_ptr
408      ENDIF
409
410      ! inverse of the ocean "zonal" v-point section
411      z_1(:,:,:) = 1.e0
412      surf_jk_r(:,:) = ptr_vtjk( z_1(:,:,:) )
413      WHERE( surf_jk_r(:,:) /= 0.e0 )   surf_jk_r(:,:) = 1.e0 / surf_jk_r(:,:)
414
415   END SUBROUTINE dia_ptr_init
416
417
418   SUBROUTINE dia_ptr_wri( kt )
419      !!---------------------------------------------------------------------
420      !!                ***  ROUTINE dia_ptr_wri  ***
421      !!
422      !! ** Purpose :   output of poleward fluxes
423      !!
424      !! ** Method  :   NetCDF file
425      !!----------------------------------------------------------------------
426      INTEGER, INTENT(in) ::   kt   ! ocean time-step index
427      !!
428      INTEGER, SAVE ::   nhoridz, ndepidzt, ndepidzw, ndex(1)
429
430      CHARACTER (len=40)       ::   clhstnam, clop                   ! temporary names
431      INTEGER                  ::   iline, it, ji                    !
432      REAL(wp)                 ::   zsto, zout, zdt, zmax, zjulian   ! temporary scalars
433      REAL(wp), DIMENSION(jpj) ::   zphi, zfoo
434      !!----------------------------------------------------------------------
435     
436      ! Define frequency of output and means
437      zdt = rdt
438      IF( nacc == 1 ) zdt = rdtmin
439#if defined key_diainstant
440         zsto = nf_ptr * zdt
441         clop = "inst(x)"               ! no use of the mask value (require less cpu time)
442         !!! clop="inst(only(x))"       ! put 1.e+20 on land (very expensive!!)
443#else
444         zsto = zdt
445         clop = "ave(x)"                ! no use of the mask value (require less cpu time)
446         !!! clop="ave(only(x))"        ! put 1.e+20 on land (very expensive!!)
447#endif
448      zout = nf_ptr * zdt
449      zmax = ( nitend - nit000 + 1 ) * zdt
450     
451         
452      ! define time axis
453      it = kt - nit000 + 1
454
455      ! Initialization
456      ! --------------
457      IF( kt == nit000 ) THEN
458     
459      zdt = rdt
460      IF( nacc == 1 ) zdt = rdtmin
461
462         ! Reference latitude
463         ! ------------------
464         !                                           ! =======================
465         IF( cp_cfg == "orca" ) THEN                 !   ORCA configurations
466            !                                        ! =======================
467
468            IF( jp_cfg == 05  )   iline = 192   ! i-line that passes near the North Pole
469            IF( jp_cfg == 025 )   iline = 384   ! i-line that passes near the North Pole
470            IF( jp_cfg == 2   )   iline =  48   ! i-line that passes near the North Pole
471            IF( jp_cfg == 4   )   iline =  24   ! i-line that passes near the North Pole
472            zphi(:) = 0.e0
473            DO ji = mi0(iline), mi1(iline) 
474               zphi(:) = gphiv(ji,:)         ! if iline is in the local domain
475               ! correct highest latitude for ORCA05
476               IF( jp_cfg == 05  ) zphi(jpj) = zphi(jpjm1) + (zphi(jpjm1)-zphi(jpj-2))/2.
477               IF( jp_cfg == 05  ) zphi(jpj) = MIN( zphi(jpj), 90.)
478
479            END DO
480            ! provide the correct zphi to all local domains
481            IF( lk_mpp )   CALL mpp_sum( zphi, jpj )       
482
483            !                                        ! =======================
484         ELSE                                        !   OTHER configurations
485            !                                        ! =======================
486            zphi(:) = gphiv(1,:)             ! assume lat/lon coordinate, select the first i-line
487            !
488         ENDIF
489
490         ! OPEN netcdf file
491         ! ----------------
492         ! Define frequency of output and means
493         zsto = nf_ptr * zdt
494         clop = "ave(x)"
495         zout = nf_ptr * zdt
496         zfoo(:) = 0.e0
497
498         ! Compute julian date from starting date of the run
499
500         CALL ymds2ju( nyear, nmonth, nday, 0.e0, zjulian )
501
502         CALL dia_nam( clhstnam, nf_ptr, 'diaptr' )
503         IF(lwp)WRITE( numout,*)" Name of diaptr NETCDF file ",clhstnam
504
505         ! Horizontal grid : zphi()
506         CALL histbeg(clhstnam, 1, zfoo, jpj, zphi,   &
507            1, 1, 1, jpj, 0, zjulian, zdt, nhoridz, numptr, domain_id=nidom )
508         ! Vertical grids : gdept_0, gdepw_0
509         CALL histvert( numptr, "deptht", "Vertical T levels",   &
510            "m", jpk, gdept_0, ndepidzt )
511         CALL histvert( numptr, "depthw", "Vertical W levels",   &
512            "m", jpk, gdepw_0, ndepidzw )
513         
514         !  Zonal mean T and S
515         
516         CALL histdef( numptr, "zotemglo", "Zonal Mean Temperature","C" ,   &
517            1, jpj, nhoridz, jpk, 1, jpk, ndepidzt, 32, clop, zsto, zout )
518         CALL histdef( numptr, "zosalglo", "Zonal Mean Salinity","PSU"  ,   &
519            1, jpj, nhoridz, jpk, 1, jpk, ndepidzt, 32, clop, zsto, zout )
520
521         !  Meridional Stream-Function (eulerian and bolus)
522         
523         CALL histdef( numptr, "zomsfglo", "Meridional Stream-Function: Global","Sv" ,   &
524            1, jpj, nhoridz, jpk, 1, jpk, ndepidzw, 32, clop, zsto, zout )
525         IF( ln_subbas ) THEN
526            CALL histdef( numptr, "zomsfatl", "Meridional Stream-Function: Atlantic","Sv" ,   &
527               1, jpj, nhoridz, jpk, 1, jpk, ndepidzw, 32, clop, zsto, zout )
528            CALL histdef( numptr, "zomsfipc", "Meridional Stream-Function: Indo-Pacific","Sv" ,&
529               1, jpj, nhoridz, jpk, 1, jpk, ndepidzw, 32, clop, zsto, zout )
530         ENDIF
531
532         !  Heat transport
533
534         CALL histdef( numptr, "sophtadv", "Advective Heat Transport"      ,   &
535            "PW", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
536         CALL histdef( numptr, "sophtldf", "Diffusive Heat Transport"      ,   &
537            "PW",1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
538         CALL histdef( numptr, "sophtove", "Overturning Heat Transport"    ,   &
539            "PW",1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
540         IF( ln_subbas ) THEN
541            CALL histdef( numptr, "sohtatl", "Heat Transport Atlantic"      ,  &
542               "PW", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
543            CALL histdef( numptr, "sohtpac", "Heat Transport Pacific"      ,   &
544               "PW", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
545            CALL histdef( numptr, "sohtind", "Heat Transport Indic"      ,     &
546               "PW", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
547         ENDIF
548
549
550         !  Salt transport
551
552         CALL histdef( numptr, "sopstadv", "Advective Salt Transport"      ,   &
553            "Giga g/s", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
554         CALL histdef( numptr, "sopstldf", "Diffusive Salt Transport"      ,   &
555            "Giga g/s", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
556         CALL histdef( numptr, "sopstove", "Overturning Salt Transport"    ,   &
557            "Giga g/s", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
558
559#if defined key_diaeiv
560         ! Eddy induced velocity
561         CALL histdef( numptr, "zomsfeiv", "Bolus Meridional Stream-Function: global",   &
562            "Sv"      , 1, jpj, nhoridz, jpk, 1, jpk, ndepidzw, 32, clop, zsto, zout )
563         CALL histdef( numptr, "sophteiv", "Bolus Advective Heat Transport",   &
564            "PW"      , 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
565         CALL histdef( numptr, "sopsteiv", "Bolus Advective Salt Transport",   &
566            "Giga g/s", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
567#endif
568         IF( ln_subbas ) THEN
569            CALL histdef( numptr, "sostatl", "Salt Transport Atlantic"      ,    &
570               "Giga g/s", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
571            CALL histdef( numptr, "sostpac", "Salt Transport Pacific"      ,     &
572               "Giga g/s", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
573            CALL histdef( numptr, "sostind", "Salt Transport Indic"      ,       &
574               "Giga g/s", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
575         ENDIF
576         
577
578         CALL histend( numptr )
579
580      ENDIF
581
582      IF( MOD( kt, nf_ptr ) == 0 ) THEN
583
584         IF(lwp) THEN
585            WRITE(numout,*)
586            WRITE(numout,*) 'dia_ptr : write Poleward Transports at time-step : ', kt
587            WRITE(numout,*) '~~~~~~~~'
588            WRITE(numout,*)
589         ENDIF
590
591         ! define time axis
592         it= kt - nit000 + 1
593         ndex(1) = 0
594         CALL histwrite( numptr, "zotemglo", it, tn_jk    , jpj*jpk, ndex )
595         CALL histwrite( numptr, "zosalglo", it, sn_jk    , jpj*jpk, ndex )
596         ! overturning outputs:
597         CALL histwrite( numptr, "zomsfglo", it, v_msf_glo , jpj*jpk, ndex )
598         IF( ln_subbas ) THEN
599            CALL histwrite( numptr, "zomsfatl", it, v_msf_atl , jpj*jpk, ndex )
600            CALL histwrite( numptr, "zomsfipc", it, v_msf_ipc , jpj*jpk, ndex )
601         ENDIF
602         ! heat transport outputs:
603         IF( ln_subbas ) THEN
604            CALL histwrite( numptr, "sohtatl", it, ht_atl  , jpj, ndex )
605            CALL histwrite( numptr, "sohtpac", it, ht_pac  , jpj, ndex )
606            CALL histwrite( numptr, "sohtind", it, ht_ind  , jpj, ndex )
607            CALL histwrite( numptr, "sostatl", it, st_atl  , jpj, ndex )
608            CALL histwrite( numptr, "sostpac", it, st_pac  , jpj, ndex )
609            CALL histwrite( numptr, "sostind", it, st_ind  , jpj, ndex )
610         ENDIF
611
612         CALL histwrite( numptr, "sophtadv", it, pht_adv  , jpj, ndex )
613         CALL histwrite( numptr, "sophtldf", it, pht_ldf  , jpj, ndex )
614         CALL histwrite( numptr, "sophtove", it, pht_ove  , jpj, ndex )
615         CALL histwrite( numptr, "sopstadv", it, pst_adv  , jpj, ndex )
616         CALL histwrite( numptr, "sopstldf", it, pst_ldf  , jpj, ndex )
617         CALL histwrite( numptr, "sopstove", it, pst_ove  , jpj, ndex )
618#if defined key_diaeiv
619         CALL histwrite( numptr, "zomsfeiv", it, v_msf_eiv, jpj*jpk, ndex )
620         CALL histwrite( numptr, "sophteiv", it, pht_eiv  , jpj    , ndex )
621         CALL histwrite( numptr, "sopsteiv", it, pst_eiv  , jpj    , ndex )
622#endif
623 
624      ENDIF
625      !
626   END SUBROUTINE dia_ptr_wri
627
628   !!======================================================================
629END MODULE diaptr
Note: See TracBrowser for help on using the repository browser.