Changeset 11993


Ignore:
Timestamp:
2019-11-28T11:20:53+01:00 (6 months ago)
Author:
cetlod
Message:

trunk : undo bad commit. Oups

Location:
NEMO/trunk/src
Files:
25 edited

Legend:

Unmodified
Added
Removed
  • NEMO/trunk/src/OCE/DIA/diaar5.F90

    r11989 r11993  
    7171      INTEGER, INTENT( in ) ::   kt   ! ocean time-step index 
    7272      ! 
    73       INTEGER  ::   ji, jj, jk, iks, ikb                      ! dummy loop arguments 
    74       REAL(wp) ::   zvolssh, zvol, zssh_steric, zztmp, zarho, ztemp, zsal, zmass, zsst 
     73      INTEGER  ::   ji, jj, jk                      ! dummy loop arguments 
     74      REAL(wp) ::   zvolssh, zvol, zssh_steric, zztmp, zarho, ztemp, zsal, zmass 
    7575      REAL(wp) ::   zaw, zbw, zrw 
    7676      ! 
    7777      REAL(wp), ALLOCATABLE, DIMENSION(:,:)     :: zarea_ssh , zbotpres       ! 2D workspace  
    78       REAL(wp), ALLOCATABLE, DIMENSION(:,:)     :: zpe, z2d                   ! 2D workspace  
    79       REAL(wp), ALLOCATABLE, DIMENSION(:,:,:)   :: zrhd , zrhop, ztpot   ! 3D workspace 
     78      REAL(wp), ALLOCATABLE, DIMENSION(:,:)     :: zpe                         ! 2D workspace  
     79      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:)   :: zrhd , zrhop               ! 3D workspace 
    8080      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) :: ztsn                       ! 4D workspace 
    8181 
     
    8686 
    8787      IF( l_ar5 ) THEN  
    88          ALLOCATE( zarea_ssh(jpi,jpj), zbotpres(jpi,jpj), z2d(jpi,jpj) ) 
     88         ALLOCATE( zarea_ssh(jpi,jpj) , zbotpres(jpi,jpj) ) 
    8989         ALLOCATE( zrhd(jpi,jpj,jpk) , zrhop(jpi,jpj,jpk) ) 
    9090         ALLOCATE( ztsn(jpi,jpj,jpk,jpts) ) 
     
    9292      ENDIF 
    9393      ! 
    94       CALL iom_put( 'e2u'      , e2u (:,:) ) 
    95       CALL iom_put( 'e1v'      , e1v (:,:) ) 
    96       CALL iom_put( 'areacello', area(:,:) ) 
    97       ! 
    98       IF( iom_use( 'volcello' ) .OR. iom_use( 'masscello' )  ) THEN   
    99          zrhd(:,:,jpk) = 0._wp        ! ocean volume ; rhd is used as workspace 
    100          DO jk = 1, jpkm1 
    101             zrhd(:,:,jk) = area(:,:) * e3t_n(:,:,jk) * tmask(:,:,jk) 
    102          END DO 
    103          CALL iom_put( 'volcello'  , zrhd(:,:,:)  )  ! WARNING not consistent with CMIP DR where volcello is at ca. 2000 
    104          CALL iom_put( 'masscello' , rau0 * e3t_n(:,:,:) * tmask(:,:,:) )  ! ocean mass 
    105       ENDIF  
    106       ! 
    107       IF( iom_use( 'e3tb' ) )  THEN    ! bottom layer thickness 
    108          DO jj = 1, jpj 
    109             DO ji = 1, jpi 
    110                ikb = mbkt(ji,jj) 
    111                z2d(ji,jj) = e3t_n(ji,jj,ikb) 
    112             END DO 
    113          END DO 
    114          CALL iom_put( 'e3tb', z2d ) 
    115       ENDIF  
    116       ! 
    11794      IF( iom_use( 'voltot' ) .OR. iom_use( 'sshtot' )  .OR. iom_use( 'sshdyn' )  ) THEN     
    11895         !                                         ! total volume of liquid seawater 
    119          zvolssh = glob_sum( 'diaar5', zarea_ssh(:,:) )  
    120          zvol    = vol0 + zvolssh 
     96         zvolssh = SUM( zarea_ssh(:,:) )  
     97         CALL mpp_sum( 'diaar5', zvolssh ) 
     98         zvol = vol0 + zvolssh 
    12199       
    122100         CALL iom_put( 'voltot', zvol               ) 
     
    140118               DO ji = 1, jpi 
    141119                  DO jj = 1, jpj 
    142                      iks = mikt(ji,jj) 
    143                      zbotpres(ji,jj) = zbotpres(ji,jj) + sshn(ji,jj) * zrhd(ji,jj,iks) + riceload(ji,jj) 
     120                     zbotpres(ji,jj) = zbotpres(ji,jj) + sshn(ji,jj) * zrhd(ji,jj,mikt(ji,jj)) + riceload(ji,jj) 
    144121                  END DO 
    145122               END DO 
     
    152129         END IF 
    153130         !                                          
    154          zarho = glob_sum( 'diaar5', area(:,:) * zbotpres(:,:) )  
     131         zarho = SUM( area(:,:) * zbotpres(:,:) )  
     132         CALL mpp_sum( 'diaar5', zarho ) 
    155133         zssh_steric = - zarho / area_tot 
    156134         CALL iom_put( 'sshthster', zssh_steric ) 
     
    169147               DO ji = 1,jpi 
    170148                  DO jj = 1,jpj 
    171                      iks = mikt(ji,jj) 
    172                      zbotpres(ji,jj) = zbotpres(ji,jj) + sshn(ji,jj) * zrhd(ji,jj,iks) + riceload(ji,jj) 
     149                     zbotpres(ji,jj) = zbotpres(ji,jj) + sshn(ji,jj) * zrhd(ji,jj,mikt(ji,jj)) + riceload(ji,jj) 
    173150                  END DO 
    174151               END DO 
     
    178155         END IF 
    179156         !     
    180          zarho = glob_sum( 'diaar5', area(:,:) * zbotpres(:,:) )  
     157         zarho = SUM( area(:,:) * zbotpres(:,:) )  
     158         CALL mpp_sum( 'diaar5', zarho ) 
    181159         zssh_steric = - zarho / area_tot 
    182160         CALL iom_put( 'sshsteric', zssh_steric ) 
     161       
    183162         !                                         ! ocean bottom pressure 
    184163         zztmp = rau0 * grav * 1.e-4_wp               ! recover pressure from pressure anomaly and cover to dbar = 1.e4 Pa 
     
    189168 
    190169      IF( iom_use( 'masstot' ) .OR. iom_use( 'temptot' )  .OR. iom_use( 'saltot' )  ) THEN     
    191           !                                         ! Mean density anomalie, temperature and salinity 
    192           ztsn(:,:,:,:) = 0._wp                    ! ztsn(:,:,1,jp_tem/sal) is used here as 2D Workspace for temperature & salinity 
    193           DO jk = 1, jpkm1 
    194              DO jj = 1, jpj 
    195                 DO ji = 1, jpi 
    196                    zztmp = area(ji,jj) * e3t_n(ji,jj,jk) 
    197                    ztsn(ji,jj,1,jp_tem) = ztsn(ji,jj,1,jp_tem) + zztmp * tsn(ji,jj,jk,jp_tem) 
    198                    ztsn(ji,jj,1,jp_sal) = ztsn(ji,jj,1,jp_sal) + zztmp * tsn(ji,jj,jk,jp_sal) 
    199                 ENDDO 
    200              ENDDO 
    201           ENDDO 
    202  
    203           IF( ln_linssh ) THEN 
     170         !                                         ! Mean density anomalie, temperature and salinity 
     171         ztemp = 0._wp 
     172         zsal  = 0._wp 
     173         DO jk = 1, jpkm1 
     174            DO jj = 1, jpj 
     175               DO ji = 1, jpi 
     176                  zztmp = area(ji,jj) * e3t_n(ji,jj,jk) 
     177                  ztemp = ztemp + zztmp * tsn(ji,jj,jk,jp_tem) 
     178                  zsal  = zsal  + zztmp * tsn(ji,jj,jk,jp_sal) 
     179               END DO 
     180            END DO 
     181         END DO 
     182         IF( ln_linssh ) THEN 
    204183            IF( ln_isfcav ) THEN 
    205184               DO ji = 1, jpi 
    206185                  DO jj = 1, jpj 
    207                      iks = mikt(ji,jj) 
    208                      ztsn(ji,jj,1,jp_tem) = ztsn(ji,jj,1,jp_tem) + zarea_ssh(ji,jj) * tsn(ji,jj,iks,jp_tem)  
    209                      ztsn(ji,jj,1,jp_sal) = ztsn(ji,jj,1,jp_sal) + zarea_ssh(ji,jj) * tsn(ji,jj,iks,jp_sal)  
     186                     ztemp = ztemp + zarea_ssh(ji,jj) * tsn(ji,jj,mikt(ji,jj),jp_tem)  
     187                     zsal  = zsal  + zarea_ssh(ji,jj) * tsn(ji,jj,mikt(ji,jj),jp_sal)  
    210188                  END DO 
    211189               END DO 
    212190            ELSE 
    213                ztsn(:,:,1,jp_tem) = ztsn(:,:,1,jp_tem) + zarea_ssh(:,:) * tsn(:,:,1,jp_tem)  
    214                ztsn(:,:,1,jp_sal) = ztsn(:,:,1,jp_sal) + zarea_ssh(:,:) * tsn(:,:,1,jp_sal)  
     191               ztemp = ztemp + SUM( zarea_ssh(:,:) * tsn(:,:,1,jp_tem) ) 
     192               zsal  = zsal  + SUM( zarea_ssh(:,:) * tsn(:,:,1,jp_sal) ) 
    215193            END IF 
    216194         ENDIF 
    217          ! 
    218          ztemp = glob_sum( 'diaar5', ztsn(:,:,1,jp_tem) ) 
    219          zsal  = glob_sum( 'diaar5', ztsn(:,:,1,jp_sal) ) 
    220          zmass = rau0 * ( zarho + zvol )       
     195         IF( lk_mpp ) THEN   
     196            CALL mpp_sum( 'diaar5', ztemp ) 
     197            CALL mpp_sum( 'diaar5', zsal  ) 
     198         END IF 
     199         ! 
     200         zmass = rau0 * ( zarho + zvol )                 ! total mass of liquid seawater 
     201         ztemp = ztemp / zvol                            ! potential temperature in liquid seawater 
     202         zsal  = zsal  / zvol                            ! Salinity of liquid seawater 
    221203         ! 
    222204         CALL iom_put( 'masstot', zmass ) 
    223          CALL iom_put( 'temptot', ztemp / zvol ) 
    224          CALL iom_put( 'saltot' , zsal  / zvol ) 
    225          ! 
    226       ENDIF      
    227  
    228       IF( ln_teos10 ) THEN        ! ! potential temperature (TEOS-10 case) 
    229          IF( iom_use( 'toce_pot') .OR. iom_use( 'temptot_pot' ) .OR. iom_use( 'sst_pot' )  & 
    230                                   .OR. iom_use( 'ssttot' ) .OR.  iom_use( 'tosmint_pot' ) ) THEN 
    231             ! 
    232             ALLOCATE( ztpot(jpi,jpj,jpk) ) 
    233             ztpot(:,:,jpk) = 0._wp 
    234             ztpot(:,:,:)   = eos_pt_from_ct( tsn(:,:,:,jp_tem), tsn(:,:,:,jp_sal) ) 
    235             ! 
    236             CALL iom_put( 'toce_pot', ztpot(:,:,:) )  ! potential temperature (TEOS-10 case) 
    237             CALL iom_put( 'sst_pot' , ztpot(:,:,1) )  ! surface temperature 
    238             ! 
    239             IF( iom_use( 'temptot_pot' ) ) THEN   ! Output potential temperature in case we use TEOS-10 
    240                z2d(:,:) = 0._wp 
    241                DO jk = 1, jpkm1 
    242                  z2d(:,:) = z2d(:,:) + area(:,:) * e3t_n(:,:,jk) * ztpot(:,:,jk) 
    243                END DO 
    244                ztemp = glob_sum( 'diaar5', z2d(:,:)  )  
    245                CALL iom_put( 'temptot_pot', ztemp / zvol ) 
    246              ENDIF 
    247              ! 
    248              IF( iom_use( 'ssttot' ) ) THEN   ! Output potential temperature in case we use TEOS-10 
    249                zsst = glob_sum( 'diaar5',  area(:,:) * ztpot(:,:,1)  )  
    250                CALL iom_put( 'ssttot', zsst / area_tot ) 
    251              ENDIF 
    252              ! Vertical integral of temperature 
    253              IF( iom_use( 'tosmint_pot') ) THEN 
    254                z2d(:,:) = 0._wp 
    255                DO jk = 1, jpkm1 
    256                   DO jj = 1, jpj 
    257                      DO ji = 1, jpi   ! vector opt. 
    258                         z2d(ji,jj) = z2d(ji,jj) + rau0 * e3t_n(ji,jj,jk) *  ztpot(ji,jj,jk) 
    259                      END DO 
    260                   END DO 
    261                END DO 
    262                CALL iom_put( 'tosmint_pot', z2d )  
    263             ENDIF 
    264             DEALLOCATE( ztpot ) 
    265         ENDIF 
    266       ELSE        
    267          IF( iom_use('ssttot') ) THEN   ! Output sst in case we use EOS-80 
    268             zsst  = glob_sum( 'diaar5', area(:,:) * tsn(:,:,1,jp_tem) ) 
    269             CALL iom_put('ssttot', zsst / area_tot ) 
    270          ENDIF 
     205         CALL iom_put( 'temptot', ztemp ) 
     206         CALL iom_put( 'saltot' , zsal  ) 
     207         ! 
    271208      ENDIF 
    272209 
    273210      IF( iom_use( 'tnpeo' )) THEN     
    274         ! Work done against stratification by vertical mixing 
    275         ! Exclude points where rn2 is negative as convection kicks in here and 
    276         ! work is not being done against stratification 
     211      ! Work done against stratification by vertical mixing 
     212      ! Exclude points where rn2 is negative as convection kicks in here and 
     213      ! work is not being done against stratification 
    277214         ALLOCATE( zpe(jpi,jpj) ) 
    278215         zpe(:,:) = 0._wp 
     
    282219                  DO ji = 1, jpi 
    283220                     IF( rn2(ji,jj,jk) > 0._wp ) THEN 
    284                         zrw = ( gdept_n(ji,jj,jk) - gdepw_n(ji,jj,jk) ) / e3w_n(ji,jj,jk) 
     221                        zrw =   ( gdepw_n(ji,jj,jk  ) - gdept_n(ji,jj,jk) )   & 
     222                           &  / ( gdept_n(ji,jj,jk-1) - gdept_n(ji,jj,jk) ) 
     223!!gm  this can be reduced to :  (depw-dept) / e3w   (NB idem dans bn2 !) 
     224!                        zrw =   ( gdept_n(ji,jj,jk) - gdepw_n(ji,jj,jk) ) / e3w_n(ji,jj,jk) 
     225!!gm end 
    285226                        ! 
    286227                        zaw = rab_n(ji,jj,jk,jp_tem) * (1. - zrw) + rab_n(ji,jj,jk-1,jp_tem)* zrw 
    287228                        zbw = rab_n(ji,jj,jk,jp_sal) * (1. - zrw) + rab_n(ji,jj,jk-1,jp_sal)* zrw 
    288229                        ! 
    289                         zpe(ji, jj) = zpe(ji,jj)   & 
     230                        zpe(ji, jj) = zpe(ji, jj)            & 
    290231                           &        -  grav * (  avt(ji,jj,jk) * zaw * (tsn(ji,jj,jk-1,jp_tem) - tsn(ji,jj,jk,jp_tem) )  & 
    291232                           &                   - avs(ji,jj,jk) * zbw * (tsn(ji,jj,jk-1,jp_sal) - tsn(ji,jj,jk,jp_sal) ) ) 
     
    298239               DO ji = 1, jpi 
    299240                  DO jj = 1, jpj 
    300                      zpe(ji,jj) = zpe(ji,jj) + avt(ji,jj,jk) * MIN(0._wp,rn2(ji,jj,jk)) * rau0 * e3w_n(ji,jj,jk) 
     241                     zpe(ji,jj) = zpe(ji,jj) + avt(ji, jj, jk) * MIN(0._wp,rn2(ji, jj, jk)) * rau0 * e3w_n(ji, jj, jk) 
    301242                  END DO 
    302243               END DO 
    303244            END DO 
    304245         ENDIF 
     246!!gm useless lbc_lnk since the computation above is performed over 1:jpi & 1:jpj 
     247!!gm           CALL lbc_lnk( 'diaar5', zpe, 'T', 1._wp)          
    305248          CALL iom_put( 'tnpeo', zpe ) 
    306249          DEALLOCATE( zpe ) 
     
    308251 
    309252      IF( l_ar5 ) THEN 
    310         DEALLOCATE( zarea_ssh , zbotpres, z2d ) 
     253        DEALLOCATE( zarea_ssh , zbotpres ) 
    311254        DEALLOCATE( zrhd      , zrhop    ) 
    312255        DEALLOCATE( ztsn                 ) 
     
    344287       CALL lbc_lnk( 'diaar5', z2d, 'U', -1. ) 
    345288       IF( cptr == 'adv' ) THEN 
    346           IF( ktra == jp_tem ) CALL iom_put( 'uadv_heattr' , rau0_rcp * z2d )  ! advective heat transport in i-direction 
    347           IF( ktra == jp_sal ) CALL iom_put( 'uadv_salttr' , rau0     * z2d )  ! advective salt transport in i-direction 
     289          IF( ktra == jp_tem ) CALL iom_put( "uadv_heattr" , rau0_rcp * z2d )  ! advective heat transport in i-direction 
     290          IF( ktra == jp_sal ) CALL iom_put( "uadv_salttr" , rau0     * z2d )  ! advective salt transport in i-direction 
    348291       ENDIF 
    349292       IF( cptr == 'ldf' ) THEN 
    350           IF( ktra == jp_tem ) CALL iom_put( 'udiff_heattr' , rau0_rcp * z2d ) ! diffusive heat transport in i-direction 
    351           IF( ktra == jp_sal ) CALL iom_put( 'udiff_salttr' , rau0     * z2d ) ! diffusive salt transport in i-direction 
     293          IF( ktra == jp_tem ) CALL iom_put( "udiff_heattr" , rau0_rcp * z2d ) ! diffusive heat transport in i-direction 
     294          IF( ktra == jp_sal ) CALL iom_put( "udiff_salttr" , rau0     * z2d ) ! diffusive salt transport in i-direction 
    352295       ENDIF 
    353296       ! 
     
    362305       CALL lbc_lnk( 'diaar5', z2d, 'V', -1. ) 
    363306       IF( cptr == 'adv' ) THEN 
    364           IF( ktra == jp_tem ) CALL iom_put( 'vadv_heattr' , rau0_rcp * z2d )  ! advective heat transport in j-direction 
    365           IF( ktra == jp_sal ) CALL iom_put( 'vadv_salttr' , rau0     * z2d )  ! advective salt transport in j-direction 
     307          IF( ktra == jp_tem ) CALL iom_put( "vadv_heattr" , rau0_rcp * z2d )  ! advective heat transport in j-direction 
     308          IF( ktra == jp_sal ) CALL iom_put( "vadv_salttr" , rau0     * z2d )  ! advective salt transport in j-direction 
    366309       ENDIF 
    367310       IF( cptr == 'ldf' ) THEN 
    368           IF( ktra == jp_tem ) CALL iom_put( 'vdiff_heattr' , rau0_rcp * z2d ) ! diffusive heat transport in j-direction 
    369           IF( ktra == jp_sal ) CALL iom_put( 'vdiff_salttr' , rau0     * z2d ) ! diffusive salt transport in j-direction 
     311          IF( ktra == jp_tem ) CALL iom_put( "vdiff_heattr" , rau0_rcp * z2d ) ! diffusive heat transport in j-direction 
     312          IF( ktra == jp_sal ) CALL iom_put( "vdiff_salttr" , rau0     * z2d ) ! diffusive salt transport in j-direction 
    370313       ENDIF 
    371314           
     
    380323      !!---------------------------------------------------------------------- 
    381324      INTEGER  ::   inum 
    382       INTEGER  ::   ik, idep 
     325      INTEGER  ::   ik 
    383326      INTEGER  ::   ji, jj, jk  ! dummy loop indices 
    384327      REAL(wp) ::   zztmp   
    385328      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) ::   zsaldta   ! Jan/Dec levitus salinity 
    386       REAL(wp), ALLOCATABLE, DIMENSION(:,:)     ::   zvol0      
    387329      ! 
    388330      !!---------------------------------------------------------------------- 
     
    398340         IF( dia_ar5_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'dia_ar5_init : unable to allocate arrays' ) 
    399341 
    400          area(:,:) = e1e2t(:,:) 
    401          area_tot  = glob_sum( 'diaar5', area(:,:) ) 
    402  
    403          ALLOCATE( zvol0(jpi,jpj) ) 
    404          zvol0 (:,:) = 0._wp 
     342         area(:,:) = e1e2t(:,:) * tmask_i(:,:) 
     343 
     344         area_tot = SUM( area(:,:) )   ;   CALL mpp_sum( 'diaar5', area_tot ) 
     345 
     346         vol0        = 0._wp 
    405347         thick0(:,:) = 0._wp 
    406348         DO jk = 1, jpkm1 
    407             DO jj = 1, jpj               ! interpolation of salinity at the last ocean level (i.e. the partial step) 
    408                DO ji = 1, jpi 
    409                   idep = tmask(ji,jj,jk) * e3t_0(ji,jj,jk) 
    410                   zvol0 (ji,jj) = zvol0 (ji,jj) +  idep * area(ji,jj) 
    411                   thick0(ji,jj) = thick0(ji,jj) +  idep     
    412                END DO 
    413             END DO 
    414          END DO 
    415          vol0 = glob_sum( 'diaar5', zvol0 ) 
    416          DEALLOCATE( zvol0 ) 
     349            vol0        = vol0        + SUM( area (:,:) * tmask(:,:,jk) * e3t_0(:,:,jk) ) 
     350            thick0(:,:) = thick0(:,:) +    tmask_i(:,:) * tmask(:,:,jk) * e3t_0(:,:,jk) 
     351         END DO 
     352         CALL mpp_sum( 'diaar5', vol0 ) 
    417353 
    418354         IF( iom_use( 'sshthster' ) ) THEN 
    419             ALLOCATE( zsaldta(jpi,jpj,jpk,jpts) ) 
     355            ALLOCATE( zsaldta(jpi,jpj,jpj,jpts) ) 
    420356            CALL iom_open ( 'sali_ref_clim_monthly', inum ) 
    421357            CALL iom_get  ( inum, jpdom_data, 'vosaline' , zsaldta(:,:,:,1), 1  ) 
  • NEMO/trunk/src/OCE/DIA/diahth.F90

    r11989 r11993  
    1111   !!            3.2  !  2009-07  (S. Masson) hc300 bugfix + cleaning + add new diag 
    1212   !!---------------------------------------------------------------------- 
     13#if   defined key_diahth 
     14   !!---------------------------------------------------------------------- 
     15   !!   'key_diahth' :                              thermocline depth diag. 
     16   !!---------------------------------------------------------------------- 
    1317   !!   dia_hth      : Compute varius diagnostics associated with the mixed layer 
    1418   !!---------------------------------------------------------------------- 
     
    2832   PUBLIC   dia_hth_alloc ! routine called by nemogcm.F90 
    2933 
    30    LOGICAL, SAVE  ::   l_hth     !: thermocline-20d depths flag 
     34   LOGICAL , PUBLIC, PARAMETER ::   lk_diahth = .TRUE.    !: thermocline-20d depths flag 
    3135    
    3236   ! note: following variables should move to local variables once iom_put is always used  
    3337   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hth    !: depth of the max vertical temperature gradient [m] 
    3438   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hd20   !: depth of 20 C isotherm                         [m] 
    35    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hd26   !: depth of 26 C isotherm                         [m] 
    3639   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hd28   !: depth of 28 C isotherm                         [m] 
    3740   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   htc3   !: heat content of first 300 m                    [W] 
    38    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   htc7   !: heat content of first 700 m                    [W] 
    39    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   htc20  !: heat content of first 2000 m                   [W] 
    40  
    4141 
    4242   !!---------------------------------------------------------------------- 
     
    5252      !!--------------------------------------------------------------------- 
    5353      ! 
    54       ALLOCATE( hth(jpi,jpj), hd20(jpi,jpj), hd26(jpi,jpj), hd28(jpi,jpj), & 
    55          &      htc3(jpi,jpj), htc7(jpi,jpj), htc20(jpi,jpj), STAT=dia_hth_alloc ) 
     54      ALLOCATE( hth(jpi,jpj), hd20(jpi,jpj), hd28(jpi,jpj), htc3(jpi,jpj), STAT=dia_hth_alloc ) 
    5655      ! 
    5756      CALL mpp_sum ( 'diahth', dia_hth_alloc ) 
     
    8382      INTEGER, INTENT( in ) ::   kt      ! ocean time-step index 
    8483      !! 
    85       INTEGER                      ::   ji, jj, jk            ! dummy loop arguments 
    86       REAL(wp)                     ::   zrho3 = 0.03_wp       ! density     criterion for mixed layer depth 
    87       REAL(wp)                     ::   zrho1 = 0.01_wp       ! density     criterion for mixed layer depth 
    88       REAL(wp)                     ::   ztem2 = 0.2_wp        ! temperature criterion for mixed layer depth 
    89       REAL(wp)                     ::   zztmp, zzdep          ! temporary scalars inside do loop 
    90       REAL(wp)                     ::   zu, zv, zw, zut, zvt  ! temporary workspace 
    91       REAL(wp), DIMENSION(jpi,jpj) ::   zabs2      ! MLD: abs( tn - tn(10m) ) = ztem2  
    92       REAL(wp), DIMENSION(jpi,jpj) ::   ztm2       ! Top of thermocline: tn = tn(10m) - ztem2      
    93       REAL(wp), DIMENSION(jpi,jpj) ::   zrho10_3   ! MLD: rho = rho10m + zrho3       
    94       REAL(wp), DIMENSION(jpi,jpj) ::   zpycn      ! pycnocline: rho = rho10m + (dr/dT)(T,S,10m)*(-0.2 degC) 
    95       REAL(wp), DIMENSION(jpi,jpj) ::   ztinv      ! max of temperature inversion 
    96       REAL(wp), DIMENSION(jpi,jpj) ::   zdepinv    ! depth of temperature inversion 
    97       REAL(wp), DIMENSION(jpi,jpj) ::   zrho0_3    ! MLD rho = rho(surf) = 0.03 
    98       REAL(wp), DIMENSION(jpi,jpj) ::   zrho0_1    ! MLD rho = rho(surf) = 0.01 
    99       REAL(wp), DIMENSION(jpi,jpj) ::   zmaxdzT    ! max of dT/dz 
    100       REAL(wp), DIMENSION(jpi,jpj) ::   zdelr      ! delta rho equivalent to deltaT = 0.2 
     84      INTEGER                          ::   ji, jj, jk            ! dummy loop arguments 
     85      INTEGER                          ::   iid, ilevel           ! temporary integers 
     86      INTEGER, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   ik20, ik28  ! levels 
     87      REAL(wp)                         ::   zavt5 = 5.e-4_wp      ! Kz criterion for the turbocline depth 
     88      REAL(wp)                         ::   zrho3 = 0.03_wp       ! density     criterion for mixed layer depth 
     89      REAL(wp)                         ::   zrho1 = 0.01_wp       ! density     criterion for mixed layer depth 
     90      REAL(wp)                         ::   ztem2 = 0.2_wp        ! temperature criterion for mixed layer depth 
     91      REAL(wp)                         ::   zthick_0, zcoef       ! temporary scalars 
     92      REAL(wp)                         ::   zztmp, zzdep          ! temporary scalars inside do loop 
     93      REAL(wp)                         ::   zu, zv, zw, zut, zvt  ! temporary workspace 
     94      REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   zabs2      ! MLD: abs( tn - tn(10m) ) = ztem2  
     95      REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   ztm2       ! Top of thermocline: tn = tn(10m) - ztem2      
     96      REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   zrho10_3   ! MLD: rho = rho10m + zrho3       
     97      REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   zpycn      ! pycnocline: rho = rho10m + (dr/dT)(T,S,10m)*(-0.2 degC) 
     98      REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   ztinv      ! max of temperature inversion 
     99      REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   zdepinv    ! depth of temperature inversion 
     100      REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   zrho0_3    ! MLD rho = rho(surf) = 0.03 
     101      REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   zrho0_1    ! MLD rho = rho(surf) = 0.01 
     102      REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   zmaxdzT    ! max of dT/dz 
     103      REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   zthick     ! vertical integration thickness  
     104      REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   zdelr      ! delta rho equivalent to deltaT = 0.2 
    101105      !!---------------------------------------------------------------------- 
    102106      IF( ln_timing )   CALL timing_start('dia_hth') 
    103107 
    104108      IF( kt == nit000 ) THEN 
    105          l_hth = .FALSE. 
    106          IF(   iom_use( 'mlddzt'   ) .OR. iom_use( 'mldr0_3'  ) .OR. iom_use( 'mldr0_1'  )    .OR.  &  
    107             &  iom_use( 'mld_dt02' ) .OR. iom_use( 'topthdep' ) .OR. iom_use( 'mldr10_3' )    .OR.  &     
    108             &  iom_use( '20d'      ) .OR. iom_use( '26d'      ) .OR. iom_use( '28d'      )    .OR.  &     
    109             &  iom_use( 'hc300'    ) .OR. iom_use( 'hc700'    ) .OR. iom_use( 'hc2000'   )    .OR.  &     
    110             &  iom_use( 'pycndep'  ) .OR. iom_use( 'tinv'     ) .OR. iom_use( 'depti'    )  ) l_hth = .TRUE. 
    111109         !                                      ! allocate dia_hth array 
    112          IF( l_hth ) THEN  
    113             IF( dia_hth_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'dia_hth : unable to allocate standard arrays' ) 
    114             IF(lwp) WRITE(numout,*) 
    115             IF(lwp) WRITE(numout,*) 'dia_hth : diagnostics of the thermocline depth' 
    116             IF(lwp) WRITE(numout,*) '~~~~~~~ ' 
    117             IF(lwp) WRITE(numout,*) 
    118          ENDIF 
     110         IF( dia_hth_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'dia_hth : unable to allocate standard arrays' ) 
     111 
     112         IF(.NOT. ALLOCATED(ik20) ) THEN 
     113            ALLOCATE(ik20(jpi,jpj), ik28(jpi,jpj), & 
     114               &      zabs2(jpi,jpj),   & 
     115               &      ztm2(jpi,jpj),    & 
     116               &      zrho10_3(jpi,jpj),& 
     117               &      zpycn(jpi,jpj),   & 
     118               &      ztinv(jpi,jpj),   & 
     119               &      zdepinv(jpi,jpj), & 
     120               &      zrho0_3(jpi,jpj), & 
     121               &      zrho0_1(jpi,jpj), & 
     122               &      zmaxdzT(jpi,jpj), & 
     123               &      zthick(jpi,jpj),  & 
     124               &      zdelr(jpi,jpj), STAT=ji) 
     125            CALL mpp_sum('diahth', ji) 
     126            IF( ji /= 0 )   CALL ctl_stop( 'STOP', 'dia_hth : unable to allocate standard ocean arrays' ) 
     127         END IF 
     128 
     129         IF(lwp) WRITE(numout,*) 
     130         IF(lwp) WRITE(numout,*) 'dia_hth : diagnostics of the thermocline depth' 
     131         IF(lwp) WRITE(numout,*) '~~~~~~~ ' 
     132         IF(lwp) WRITE(numout,*) 
    119133      ENDIF 
    120134 
    121       IF( l_hth ) THEN 
    122          ! 
    123          IF( iom_use( 'mlddzt' ) .OR. iom_use( 'mldr0_3' ) .OR. iom_use( 'mldr0_1' ) ) THEN 
    124             ! initialization 
    125             ztinv  (:,:) = 0._wp   
    126             zdepinv(:,:) = 0._wp   
    127             zmaxdzT(:,:) = 0._wp   
    128             DO jj = 1, jpj 
    129                DO ji = 1, jpi 
    130                   zztmp = gdepw_n(ji,jj,mbkt(ji,jj)+1)  
    131                   hth     (ji,jj) = zztmp 
    132                   zabs2   (ji,jj) = zztmp 
    133                   ztm2    (ji,jj) = zztmp 
    134                   zrho10_3(ji,jj) = zztmp 
    135                   zpycn   (ji,jj) = zztmp 
    136                  END DO 
     135      ! initialization 
     136      ztinv  (:,:) = 0._wp   
     137      zdepinv(:,:) = 0._wp   
     138      zmaxdzT(:,:) = 0._wp   
     139      DO jj = 1, jpj 
     140         DO ji = 1, jpi 
     141            zztmp = gdepw_n(ji,jj,mbkt(ji,jj)+1)  
     142            hth     (ji,jj) = zztmp 
     143            zabs2   (ji,jj) = zztmp 
     144            ztm2    (ji,jj) = zztmp 
     145            zrho10_3(ji,jj) = zztmp 
     146            zpycn   (ji,jj) = zztmp 
     147        END DO 
     148      END DO 
     149      IF( nla10 > 1 ) THEN  
     150         DO jj = 1, jpj 
     151            DO ji = 1, jpi 
     152               zztmp = gdepw_n(ji,jj,mbkt(ji,jj)+1)  
     153               zrho0_3(ji,jj) = zztmp 
     154               zrho0_1(ji,jj) = zztmp 
    137155            END DO 
    138             IF( nla10 > 1 ) THEN  
    139                DO jj = 1, jpj 
    140                   DO ji = 1, jpi 
    141                      zztmp = gdepw_n(ji,jj,mbkt(ji,jj)+1)  
    142                      zrho0_3(ji,jj) = zztmp 
    143                      zrho0_1(ji,jj) = zztmp 
    144                   END DO 
    145                END DO 
     156         END DO 
     157      ENDIF 
     158       
     159      ! Preliminary computation 
     160      ! computation of zdelr = (dr/dT)(T,S,10m)*(-0.2 degC) 
     161      DO jj = 1, jpj 
     162         DO ji = 1, jpi 
     163            IF( tmask(ji,jj,nla10) == 1. ) THEN 
     164               zu  =  1779.50 + 11.250 * tsn(ji,jj,nla10,jp_tem) - 3.80   * tsn(ji,jj,nla10,jp_sal)                             & 
     165                  &                                              - 0.0745 * tsn(ji,jj,nla10,jp_tem) * tsn(ji,jj,nla10,jp_tem)   & 
     166                  &                                              - 0.0100 * tsn(ji,jj,nla10,jp_tem) * tsn(ji,jj,nla10,jp_sal) 
     167               zv  =  5891.00 + 38.000 * tsn(ji,jj,nla10,jp_tem) + 3.00   * tsn(ji,jj,nla10,jp_sal)                             & 
     168                  &                                              - 0.3750 * tsn(ji,jj,nla10,jp_tem) * tsn(ji,jj,nla10,jp_tem) 
     169               zut =    11.25 -  0.149 * tsn(ji,jj,nla10,jp_tem) - 0.01   * tsn(ji,jj,nla10,jp_sal) 
     170               zvt =    38.00 -  0.750 * tsn(ji,jj,nla10,jp_tem) 
     171               zw  = (zu + 0.698*zv) * (zu + 0.698*zv) 
     172               zdelr(ji,jj) = ztem2 * (1000.*(zut*zv - zvt*zu)/zw) 
     173            ELSE 
     174               zdelr(ji,jj) = 0._wp 
    146175            ENDIF 
     176         END DO 
     177      END DO 
     178 
     179      ! ------------------------------------------------------------- ! 
     180      ! thermocline depth: strongest vertical gradient of temperature ! 
     181      ! turbocline depth (mixing layer depth): avt = zavt5            ! 
     182      ! MLD: rho = rho(1) + zrho3                                     ! 
     183      ! MLD: rho = rho(1) + zrho1                                     ! 
     184      ! ------------------------------------------------------------- ! 
     185      DO jk = jpkm1, 2, -1   ! loop from bottom to 2 
     186         DO jj = 1, jpj 
     187            DO ji = 1, jpi 
     188               ! 
     189               zzdep = gdepw_n(ji,jj,jk) 
     190               zztmp = ( tsn(ji,jj,jk-1,jp_tem) - tsn(ji,jj,jk,jp_tem) ) / zzdep * tmask(ji,jj,jk)   ! vertical gradient of temperature (dT/dz) 
     191               zzdep = zzdep * tmask(ji,jj,1) 
     192 
     193               IF( zztmp > zmaxdzT(ji,jj) ) THEN                         
     194                  zmaxdzT(ji,jj) = zztmp   ;   hth    (ji,jj) = zzdep                ! max and depth of dT/dz 
     195               ENDIF 
     196                
     197               IF( nla10 > 1 ) THEN  
     198                  zztmp = rhop(ji,jj,jk) - rhop(ji,jj,1)                             ! delta rho(1) 
     199                  IF( zztmp > zrho3 )          zrho0_3(ji,jj) = zzdep                ! > 0.03 
     200                  IF( zztmp > zrho1 )          zrho0_1(ji,jj) = zzdep                ! > 0.01 
     201               ENDIF 
     202 
     203            END DO 
     204         END DO 
     205      END DO 
    147206       
    148             ! Preliminary computation 
    149             ! computation of zdelr = (dr/dT)(T,S,10m)*(-0.2 degC) 
    150             DO jj = 1, jpj 
    151                DO ji = 1, jpi 
    152                   IF( tmask(ji,jj,nla10) == 1. ) THEN 
    153                      zu  =  1779.50 + 11.250 * tsn(ji,jj,nla10,jp_tem) - 3.80   * tsn(ji,jj,nla10,jp_sal)  & 
    154                         &           - 0.0745 * tsn(ji,jj,nla10,jp_tem) * tsn(ji,jj,nla10,jp_tem)   & 
    155                         &           - 0.0100 * tsn(ji,jj,nla10,jp_tem) * tsn(ji,jj,nla10,jp_sal) 
    156                      zv  =  5891.00 + 38.000 * tsn(ji,jj,nla10,jp_tem) + 3.00   * tsn(ji,jj,nla10,jp_sal)  & 
    157                         &           - 0.3750 * tsn(ji,jj,nla10,jp_tem) * tsn(ji,jj,nla10,jp_tem) 
    158                      zut =    11.25 -  0.149 * tsn(ji,jj,nla10,jp_tem) - 0.01   * tsn(ji,jj,nla10,jp_sal) 
    159                      zvt =    38.00 -  0.750 * tsn(ji,jj,nla10,jp_tem) 
    160                      zw  = (zu + 0.698*zv) * (zu + 0.698*zv) 
    161                      zdelr(ji,jj) = ztem2 * (1000.*(zut*zv - zvt*zu)/zw) 
    162                   ELSE 
    163                      zdelr(ji,jj) = 0._wp 
    164                   ENDIF 
    165                END DO 
     207      CALL iom_put( "mlddzt", hth )            ! depth of the thermocline 
     208      IF( nla10 > 1 ) THEN  
     209         CALL iom_put( "mldr0_3", zrho0_3 )   ! MLD delta rho(surf) = 0.03 
     210         CALL iom_put( "mldr0_1", zrho0_1 )   ! MLD delta rho(surf) = 0.01 
     211      ENDIF 
     212 
     213      ! ------------------------------------------------------------- ! 
     214      ! MLD: abs( tn - tn(10m) ) = ztem2                              ! 
     215      ! Top of thermocline: tn = tn(10m) - ztem2                      ! 
     216      ! MLD: rho = rho10m + zrho3                                     ! 
     217      ! pycnocline: rho = rho10m + (dr/dT)(T,S,10m)*(-0.2 degC)       ! 
     218      ! temperature inversion: max( 0, max of tn - tn(10m) )          ! 
     219      ! depth of temperature inversion                                ! 
     220      ! ------------------------------------------------------------- ! 
     221      DO jk = jpkm1, nlb10, -1   ! loop from bottom to nlb10 
     222         DO jj = 1, jpj 
     223            DO ji = 1, jpi 
     224               ! 
     225               zzdep = gdepw_n(ji,jj,jk) * tmask(ji,jj,1) 
     226               ! 
     227               zztmp = tsn(ji,jj,nla10,jp_tem) - tsn(ji,jj,jk,jp_tem)  ! - delta T(10m) 
     228               IF( ABS(zztmp) > ztem2 )      zabs2   (ji,jj) = zzdep   ! abs > 0.2 
     229               IF(     zztmp  > ztem2 )      ztm2    (ji,jj) = zzdep   ! > 0.2 
     230               zztmp = -zztmp                                          ! delta T(10m) 
     231               IF( zztmp >  ztinv(ji,jj) ) THEN                        ! temperature inversion 
     232                  ztinv(ji,jj) = zztmp   ;   zdepinv (ji,jj) = zzdep   ! max value and depth 
     233               ENDIF 
     234 
     235               zztmp = rhop(ji,jj,jk) - rhop(ji,jj,nla10)              ! delta rho(10m) 
     236               IF( zztmp > zrho3        )    zrho10_3(ji,jj) = zzdep   ! > 0.03 
     237               IF( zztmp > zdelr(ji,jj) )    zpycn   (ji,jj) = zzdep   ! > equi. delta T(10m) - 0.2 
     238               ! 
    166239            END DO 
    167  
    168             ! ------------------------------------------------------------- ! 
    169             ! thermocline depth: strongest vertical gradient of temperature ! 
    170             ! turbocline depth (mixing layer depth): avt = zavt5            ! 
    171             ! MLD: rho = rho(1) + zrho3                                     ! 
    172             ! MLD: rho = rho(1) + zrho1                                     ! 
    173             ! ------------------------------------------------------------- ! 
    174             DO jk = jpkm1, 2, -1   ! loop from bottom to 2 
    175                DO jj = 1, jpj 
    176                   DO ji = 1, jpi 
    177                      ! 
    178                      zzdep = gdepw_n(ji,jj,jk) 
    179                      zztmp = ( tsn(ji,jj,jk-1,jp_tem) - tsn(ji,jj,jk,jp_tem) ) & 
    180                             & / zzdep * tmask(ji,jj,jk)   ! vertical gradient of temperature (dT/dz) 
    181                      zzdep = zzdep * tmask(ji,jj,1) 
    182  
    183                      IF( zztmp > zmaxdzT(ji,jj) ) THEN                         
    184                          zmaxdzT(ji,jj) = zztmp    
    185                          hth    (ji,jj) = zzdep                ! max and depth of dT/dz 
    186                      ENDIF 
    187                 
    188                      IF( nla10 > 1 ) THEN  
    189                         zztmp = rhop(ji,jj,jk) - rhop(ji,jj,1)                       ! delta rho(1) 
    190                         IF( zztmp > zrho3 )   zrho0_3(ji,jj) = zzdep                ! > 0.03 
    191                         IF( zztmp > zrho1 )   zrho0_1(ji,jj) = zzdep                ! > 0.01 
    192                      ENDIF 
    193                   END DO 
    194                END DO 
    195             END DO 
    196           
    197             CALL iom_put( 'mlddzt', hth )            ! depth of the thermocline 
    198             IF( nla10 > 1 ) THEN  
    199                CALL iom_put( 'mldr0_3', zrho0_3 )   ! MLD delta rho(surf) = 0.03 
    200                CALL iom_put( 'mldr0_1', zrho0_1 )   ! MLD delta rho(surf) = 0.01 
    201             ENDIF 
    202             ! 
    203          ENDIF 
    204          ! 
    205          IF(  iom_use( 'mld_dt02' ) .OR. iom_use( 'topthdep' ) .OR. iom_use( 'mldr10_3' ) .OR.  &     
    206             &  iom_use( 'pycndep' ) .OR. iom_use( 'tinv'     ) .OR. iom_use( 'depti'    )  ) THEN 
    207             ! ------------------------------------------------------------- ! 
    208             ! MLD: abs( tn - tn(10m) ) = ztem2                              ! 
    209             ! Top of thermocline: tn = tn(10m) - ztem2                      ! 
    210             ! MLD: rho = rho10m + zrho3                                     ! 
    211             ! pycnocline: rho = rho10m + (dr/dT)(T,S,10m)*(-0.2 degC)       ! 
    212             ! temperature inversion: max( 0, max of tn - tn(10m) )          ! 
    213             ! depth of temperature inversion                                ! 
    214             ! ------------------------------------------------------------- ! 
    215             DO jk = jpkm1, nlb10, -1   ! loop from bottom to nlb10 
    216                DO jj = 1, jpj 
    217                   DO ji = 1, jpi 
    218                      ! 
    219                      zzdep = gdepw_n(ji,jj,jk) * tmask(ji,jj,1) 
    220                      ! 
    221                      zztmp = tsn(ji,jj,nla10,jp_tem) - tsn(ji,jj,jk,jp_tem)  ! - delta T(10m) 
    222                      IF( ABS(zztmp) > ztem2 )      zabs2   (ji,jj) = zzdep   ! abs > 0.2 
    223                      IF(     zztmp  > ztem2 )      ztm2    (ji,jj) = zzdep   ! > 0.2 
    224                      zztmp = -zztmp                                          ! delta T(10m) 
    225                      IF( zztmp >  ztinv(ji,jj) ) THEN                        ! temperature inversion 
    226                         ztinv(ji,jj) = zztmp    
    227                         zdepinv (ji,jj) = zzdep   ! max value and depth 
    228                      ENDIF 
    229  
    230                      zztmp = rhop(ji,jj,jk) - rhop(ji,jj,nla10)              ! delta rho(10m) 
    231                      IF( zztmp > zrho3        )    zrho10_3(ji,jj) = zzdep   ! > 0.03 
    232                      IF( zztmp > zdelr(ji,jj) )    zpycn   (ji,jj) = zzdep   ! > equi. delta T(10m) - 0.2 
    233                      ! 
    234                   END DO 
    235                END DO 
    236             END DO 
    237  
    238             CALL iom_put( 'mld_dt02', zabs2    )   ! MLD abs(delta t) - 0.2 
    239             CALL iom_put( 'topthdep', ztm2     )   ! T(10) - 0.2 
    240             CALL iom_put( 'mldr10_3', zrho10_3 )   ! MLD delta rho(10m) = 0.03 
    241             CALL iom_put( 'pycndep' , zpycn    )   ! MLD delta rho equi. delta T(10m) = 0.2 
    242             CALL iom_put( 'tinv'    , ztinv    )   ! max. temp. inv. (t10 ref)  
    243             CALL iom_put( 'depti'   , zdepinv  )   ! depth of max. temp. inv. (t10 ref)  
    244             ! 
    245          ENDIF 
    246   
    247          ! ------------------------------- ! 
    248          !  Depth of 20C/26C/28C isotherm  ! 
    249          ! ------------------------------- ! 
    250          IF( iom_use ('20d') ) THEN  ! depth of the 20 isotherm 
    251             ztem2 = 20. 
    252             CALL dia_hth_dep( ztem2, hd20 )   
    253             CALL iom_put( '20d', hd20 )     
    254          ENDIF 
    255          ! 
    256          IF( iom_use ('26d') ) THEN  ! depth of the 26 isotherm 
    257             ztem2 = 26. 
    258             CALL dia_hth_dep( ztem2, hd26 )   
    259             CALL iom_put( '26d', hd26 )     
    260          ENDIF 
    261          ! 
    262          IF( iom_use ('28d') ) THEN  ! depth of the 28 isotherm 
    263             ztem2 = 28. 
    264             CALL dia_hth_dep( ztem2, hd28 )   
    265             CALL iom_put( '28d', hd28 )     
    266          ENDIF 
    267          
    268          ! ----------------------------- ! 
    269          !  Heat content of first 300 m  ! 
    270          ! ----------------------------- ! 
    271          IF( iom_use ('hc300') ) THEN   
    272             zzdep = 300. 
    273             CALL  dia_hth_htc( zzdep, tsn(:,:,:,jp_tem), htc3 ) 
    274             CALL iom_put( 'hc300', rau0_rcp * htc3 )  ! vertically integrated heat content (J/m2) 
    275          ENDIF 
    276          ! 
    277          ! ----------------------------- ! 
    278          !  Heat content of first 700 m  ! 
    279          ! ----------------------------- ! 
    280          IF( iom_use ('hc700') ) THEN   
    281             zzdep = 700. 
    282             CALL  dia_hth_htc( zzdep, tsn(:,:,:,jp_tem), htc7 ) 
    283             CALL iom_put( 'hc700', rau0_rcp * htc7 )  ! vertically integrated heat content (J/m2) 
    284    
    285          ENDIF 
    286          ! 
    287          ! ----------------------------- ! 
    288          !  Heat content of first 2000 m  ! 
    289          ! ----------------------------- ! 
    290          IF( iom_use ('hc2000') ) THEN   
    291             zzdep = 2000. 
    292             CALL  dia_hth_htc( zzdep, tsn(:,:,:,jp_tem), htc20 ) 
    293             CALL iom_put( 'hc2000', rau0_rcp * htc20 )  ! vertically integrated heat content (J/m2)   
    294          ENDIF 
    295          ! 
    296       ENDIF 
    297  
    298       ! 
    299       IF( ln_timing )   CALL timing_stop('dia_hth') 
    300       ! 
    301    END SUBROUTINE dia_hth 
    302  
    303    SUBROUTINE dia_hth_dep( ptem, pdept ) 
    304       ! 
    305       REAL(wp), INTENT(in) :: ptem 
    306       REAL(wp), DIMENSION(jpi,jpj), INTENT(out) :: pdept      
    307       ! 
    308       INTEGER  :: ji, jj, jk, iid 
    309       REAL(wp) :: zztmp, zzdep 
    310       INTEGER, DIMENSION(jpi,jpj) :: iktem 
    311        
    312       ! --------------------------------------- ! 
    313       ! search deepest level above ptem         ! 
    314       ! --------------------------------------- ! 
    315       iktem(:,:) = 1 
     240         END DO 
     241      END DO 
     242 
     243      CALL iom_put( "mld_dt02", zabs2        )   ! MLD abs(delta t) - 0.2 
     244      CALL iom_put( "topthdep", ztm2         )   ! T(10) - 0.2 
     245      CALL iom_put( "mldr10_3", zrho10_3     )   ! MLD delta rho(10m) = 0.03 
     246      CALL iom_put( "pycndep" , zpycn        )   ! MLD delta rho equi. delta T(10m) = 0.2 
     247      CALL iom_put( "tinv"    , ztinv        )   ! max. temp. inv. (t10 ref)  
     248      CALL iom_put( "depti"   , zdepinv      )   ! depth of max. temp. inv. (t10 ref)  
     249 
     250 
     251      ! ----------------------------------- ! 
     252      ! search deepest level above 20C/28C  ! 
     253      ! ----------------------------------- ! 
     254      ik20(:,:) = 1 
     255      ik28(:,:) = 1 
    316256      DO jk = 1, jpkm1   ! beware temperature is not always decreasing with depth => loop from top to bottom 
    317257         DO jj = 1, jpj 
    318258            DO ji = 1, jpi 
    319259               zztmp = tsn(ji,jj,jk,jp_tem) 
    320                IF( zztmp >= ptem )   iktem(ji,jj) = jk 
     260               IF( zztmp >= 20. )   ik20(ji,jj) = jk 
     261               IF( zztmp >= 28. )   ik28(ji,jj) = jk 
    321262            END DO 
    322263         END DO 
    323264      END DO 
    324265 
    325       ! ------------------------------- ! 
    326       !  Depth of ptem isotherm         ! 
    327       ! ------------------------------- ! 
     266      ! --------------------------- ! 
     267      !  Depth of 20C/28C isotherm  ! 
     268      ! --------------------------- ! 
    328269      DO jj = 1, jpj 
    329270         DO ji = 1, jpi 
    330271            ! 
    331             zzdep = gdepw_n(ji,jj,mbkt(ji,jj)+1)       ! depth of the ocean bottom 
     272            zzdep = gdepw_n(ji,jj,mbkt(ji,jj)+1)       ! depth of the oean bottom 
    332273            ! 
    333             iid = iktem(ji,jj) 
     274            iid = ik20(ji,jj) 
    334275            IF( iid /= 1 ) THEN  
    335                 zztmp =     gdept_n(ji,jj,iid  )   &                     ! linear interpolation 
     276               zztmp =      gdept_n(ji,jj,iid  )   &                     ! linear interpolation 
    336277                  &  + (    gdept_n(ji,jj,iid+1) - gdept_n(ji,jj,iid)                       )   & 
    337278                  &  * ( 20.*tmask(ji,jj,iid+1) - tsn(ji,jj,iid,jp_tem)                       )   & 
    338279                  &  / ( tsn(ji,jj,iid+1,jp_tem) - tsn(ji,jj,iid,jp_tem) + (1.-tmask(ji,jj,1)) ) 
    339                pdept(ji,jj) = MIN( zztmp , zzdep) * tmask(ji,jj,1)       ! bound by the ocean depth 
     280               hd20(ji,jj) = MIN( zztmp , zzdep) * tmask(ji,jj,1)       ! bound by the ocean depth 
    340281            ELSE  
    341                pdept(ji,jj) = 0._wp 
     282               hd20(ji,jj) = 0._wp 
    342283            ENDIF 
    343          END DO 
    344       END DO 
    345       ! 
    346    END SUBROUTINE dia_hth_dep 
    347  
    348  
    349    SUBROUTINE dia_hth_htc( pdep, ptn, phtc ) 
    350       ! 
    351       REAL(wp), INTENT(in) :: pdep     ! depth over the heat content 
    352       REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in) :: ptn    
    353       REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: phtc   
    354       ! 
    355       INTEGER  :: ji, jj, jk, ik 
    356       REAL(wp), DIMENSION(jpi,jpj) :: zthick 
    357       INTEGER , DIMENSION(jpi,jpj) :: ilevel 
    358  
    359  
     284            ! 
     285            iid = ik28(ji,jj) 
     286            IF( iid /= 1 ) THEN  
     287               zztmp =      gdept_n(ji,jj,iid  )   &                     ! linear interpolation 
     288                  &  + (    gdept_n(ji,jj,iid+1) - gdept_n(ji,jj,iid)                       )   & 
     289                  &  * ( 28.*tmask(ji,jj,iid+1) -    tsn(ji,jj,iid,jp_tem)                       )   & 
     290                  &  / (  tsn(ji,jj,iid+1,jp_tem) -    tsn(ji,jj,iid,jp_tem) + (1.-tmask(ji,jj,1)) ) 
     291               hd28(ji,jj) = MIN( zztmp , zzdep ) * tmask(ji,jj,1)      ! bound by the ocean depth 
     292            ELSE  
     293               hd28(ji,jj) = 0._wp 
     294            ENDIF 
     295 
     296         END DO 
     297      END DO 
     298      CALL iom_put( "20d", hd20 )   ! depth of the 20 isotherm 
     299      CALL iom_put( "28d", hd28 )   ! depth of the 28 isotherm 
     300 
     301      ! ----------------------------- ! 
     302      !  Heat content of first 300 m  ! 
     303      ! ----------------------------- ! 
     304 
     305      ! find ilevel with (ilevel+1) the deepest W-level above 300m (we assume we can use e3t_1d to do this search...) 
     306      ilevel   = 0 
     307      zthick_0 = 0._wp 
     308      DO jk = 1, jpkm1                       
     309         zthick_0 = zthick_0 + e3t_1d(jk) 
     310         IF( zthick_0 < 300. )   ilevel = jk 
     311      END DO 
    360312      ! surface boundary condition 
    361        
    362       IF( .NOT. ln_linssh ) THEN   ;   zthick(:,:) = 0._wp       ;   phtc(:,:) = 0._wp                                    
    363       ELSE                         ;   zthick(:,:) = sshn(:,:)   ;   phtc(:,:) = ptn(:,:,1) * sshn(:,:) * tmask(:,:,1)    
     313      IF( ln_linssh ) THEN   ;   zthick(:,:) = sshn(:,:)   ;   htc3(:,:) = tsn(:,:,1,jp_tem) * sshn(:,:) * tmask(:,:,1)   
     314      ELSE                   ;   zthick(:,:) = 0._wp       ;   htc3(:,:) = 0._wp                                    
    364315      ENDIF 
    365       ! 
    366       ilevel(:,:) = 1 
    367       DO jk = 2, jpkm1 
    368          DO jj = 1, jpj 
    369             DO ji = 1, jpi 
    370                IF( ( gdept_n(ji,jj,jk) < pdep ) .AND. ( tmask(ji,jj,jk) == 1 ) ) THEN 
    371                    ilevel(ji,jj) = jk 
    372                    zthick(ji,jj) = zthick(ji,jj) + e3t_n(ji,jj,jk) 
    373                    phtc  (ji,jj) = phtc  (ji,jj) + e3t_n(ji,jj,jk) * ptn(ji,jj,jk) 
    374                ENDIF 
    375             ENDDO 
    376          ENDDO 
    377       ENDDO 
    378       ! 
     316      ! integration down to ilevel 
     317      DO jk = 1, ilevel 
     318         zthick(:,:) = zthick(:,:) + e3t_n(:,:,jk) 
     319         htc3  (:,:) = htc3  (:,:) + e3t_n(:,:,jk) * tsn(:,:,jk,jp_tem) * tmask(:,:,jk) 
     320      END DO 
     321      ! deepest layer 
     322      zthick(:,:) = 300. - zthick(:,:)   !   remaining thickness to reach 300m 
    379323      DO jj = 1, jpj 
    380324         DO ji = 1, jpi 
    381             ik = ilevel(ji,jj) 
    382             zthick(ji,jj) = pdep - zthick(ji,jj)   !   remaining thickness to reach depht pdep 
    383             phtc(ji,jj)   = phtc(ji,jj) + ptn(ji,jj,ik+1) * MIN( e3t_n(ji,jj,ik+1), zthick(ji,jj) ) & 
    384                                                           * tmask(ji,jj,ik+1) 
    385          END DO 
    386       ENDDO 
    387       ! 
    388       ! 
    389    END SUBROUTINE dia_hth_htc 
     325            htc3(ji,jj) = htc3(ji,jj) + tsn(ji,jj,ilevel+1,jp_tem)                  & 
     326               &                      * MIN( e3t_n(ji,jj,ilevel+1), zthick(ji,jj) ) * tmask(ji,jj,ilevel+1) 
     327         END DO 
     328      END DO 
     329      ! from temperature to heat contain 
     330      zcoef = rau0 * rcp 
     331      htc3(:,:) = zcoef * htc3(:,:) 
     332      CALL iom_put( "hc300", htc3 )      ! first 300m heat content 
     333      ! 
     334      IF( ln_timing )   CALL timing_stop('dia_hth') 
     335      ! 
     336   END SUBROUTINE dia_hth 
     337 
     338#else 
     339   !!---------------------------------------------------------------------- 
     340   !!   Default option :                                       Empty module 
     341   !!---------------------------------------------------------------------- 
     342   LOGICAL , PUBLIC, PARAMETER ::   lk_diahth = .FALSE.  !: thermocline-20d depths flag 
     343CONTAINS 
     344   SUBROUTINE dia_hth( kt )         ! Empty routine 
     345      IMPLICIT NONE 
     346      INTEGER, INTENT( in ) :: kt 
     347      WRITE(*,*) 'dia_hth: You should not have seen this print! error?', kt 
     348   END SUBROUTINE dia_hth 
     349#endif 
    390350 
    391351   !!====================================================================== 
  • NEMO/trunk/src/OCE/DIA/diaptr.F90

    r11989 r11993  
    1010   !!            3.6  ! 2014-12  (C. Ethe) use of IOM 
    1111   !!            3.6  ! 2016-06  (T. Graham) Addition of diagnostics for CMIP6 
    12    !!            4.0  ! 2010-08  ( C. Ethe, J. Deshayes ) Improvment 
    1312   !!---------------------------------------------------------------------- 
    1413 
     
    4342 
    4443   !                                  !!** namelist  namptr  ** 
    45    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   hstr_adv, hstr_ldf, hstr_eiv   !: Heat/Salt TRansports(adv, diff, Bolus.) 
    46    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   hstr_ove, hstr_btr, hstr_vtr   !: heat Salt TRansports(overturn, baro, merional) 
    47  
    48    LOGICAL , PUBLIC ::   l_diaptr        !: tracers  trend flag (set from namelist in trdini) 
    49    INTEGER, PARAMETER, PUBLIC ::   nptr = 5  ! (glo, atl, pac, ind, ipc) 
     44   REAL(wp), ALLOCATABLE, SAVE, PUBLIC, DIMENSION(:,:) ::   htr_adv, htr_ldf, htr_eiv   !: Heat TRansports (adv, diff, Bolus.) 
     45   REAL(wp), ALLOCATABLE, SAVE, PUBLIC, DIMENSION(:,:) ::   str_adv, str_ldf, str_eiv   !: Salt TRansports (adv, diff, Bolus.) 
     46   REAL(wp), ALLOCATABLE, SAVE, PUBLIC, DIMENSION(:,:) ::   htr_ove, str_ove   !: heat Salt TRansports ( overturn.) 
     47   REAL(wp), ALLOCATABLE, SAVE, PUBLIC, DIMENSION(:,:) ::   htr_btr, str_btr   !: heat Salt TRansports ( barotropic ) 
     48 
     49   LOGICAL, PUBLIC ::   ln_diaptr   !  Poleward transport flag (T) or not (F) 
     50   LOGICAL, PUBLIC ::   ln_subbas   !  Atlantic/Pacific/Indian basins calculation 
     51   INTEGER, PUBLIC ::   nptr        ! = 1 (l_subbas=F) or = 5 (glo, atl, pac, ind, ipc) (l_subbas=T)  
    5052 
    5153   REAL(wp) ::   rc_sv    = 1.e-6_wp   ! conversion from m3/s to Sverdrup 
    5254   REAL(wp) ::   rc_pwatt = 1.e-15_wp  ! conversion from W    to PW (further x rau0 x Cp) 
    53    REAL(wp) ::   rc_ggram = 1.e-9_wp   ! conversion from g    to Gg  (further x rau0) 
    54  
    55    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: btmsk   ! T-point basin interior masks 
    56    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: btmsk34 ! mask out Southern Ocean (=0 south of 34°S) 
    57  
    58    REAL(wp), TARGET, ALLOCATABLE, SAVE, DIMENSION(:)   :: p_fval1d 
    59    REAL(wp), TARGET, ALLOCATABLE, SAVE, DIMENSION(:,:) :: p_fval2d 
    60  
    61    LOGICAL ::   ll_init = .TRUE.        !: tracers  trend flag (set from namelist in trdini) 
     55   REAL(wp) ::   rc_ggram = 1.e-6_wp   ! conversion from g    to Pg 
     56 
     57   CHARACTER(len=3), ALLOCATABLE, SAVE, DIMENSION(:)     :: clsubb 
     58   REAL(wp),         ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: btmsk   ! T-point basin interior masks 
     59   REAL(wp),         ALLOCATABLE, SAVE, DIMENSION(:,:)   :: btm30   ! mask out Southern Ocean (=0 south of 30°S) 
     60 
     61   REAL(wp), TARGET, ALLOCATABLE, SAVE, DIMENSION(:)     :: p_fval1d 
     62   REAL(wp), TARGET, ALLOCATABLE, SAVE, DIMENSION(:,:)   :: p_fval2d 
     63 
    6264   !! * Substitutions 
    6365#  include "vectopt_loop_substitute.h90" 
     
    6971CONTAINS 
    7072 
    71    SUBROUTINE dia_ptr( kt, pvtr ) 
     73   SUBROUTINE dia_ptr( pvtr ) 
    7274      !!---------------------------------------------------------------------- 
    7375      !!                  ***  ROUTINE dia_ptr  *** 
    7476      !!---------------------------------------------------------------------- 
    75       INTEGER, INTENT( in ) ::   kt   ! ocean time-step index      
    7677      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in), OPTIONAL ::   pvtr   ! j-effective transport 
    7778      ! 
     
    7980      REAL(wp) ::   zsfc,zvfc               ! local scalar 
    8081      REAL(wp), DIMENSION(jpi,jpj)     ::  z2d   ! 2D workspace 
     82      REAL(wp), DIMENSION(jpi,jpj,jpk) ::  z3d   ! 3D workspace 
    8183      REAL(wp), DIMENSION(jpi,jpj,jpk) ::  zmask   ! 3D workspace 
    82       REAL(wp), DIMENSION(jpi,jpj,jpk) ::  z3d    ! 3D workspace 
    8384      REAL(wp), DIMENSION(jpi,jpj,jpk,jpts) ::  zts   ! 3D workspace 
    84       REAL(wp), DIMENSION(jpj)      ::  zvsum, ztsum, zssum   ! 1D workspace 
     85      REAL(wp), DIMENSION(jpj)     ::  vsum   ! 1D workspace 
     86      REAL(wp), DIMENSION(jpj,jpts)     ::  tssum   ! 1D workspace 
     87  
    8588      ! 
    8689      !overturning calculation 
    87       REAL(wp), DIMENSION(jpj,jpk,nptr) :: sjk, r1_sjk, v_msf  ! i-mean i-k-surface and its inverse 
    88       REAL(wp), DIMENSION(jpj,jpk,nptr) :: zt_jk, zs_jk ! i-mean T and S, j-Stream-Function 
    89  
    90       REAL(wp), DIMENSION(jpi,jpj,jpk,nptr)  :: z4d1, z4d2 
    91       REAL(wp), DIMENSION(jpi,jpj,nptr)      :: z3dtr ! i-mean T and S, j-Stream-Function 
     90      REAL(wp), DIMENSION(jpj,jpk,nptr) ::   sjk  , r1_sjk ! i-mean i-k-surface and its inverse 
     91      REAL(wp), DIMENSION(jpj,jpk,nptr) ::   v_msf, sn_jk  , tn_jk ! i-mean T and S, j-Stream-Function 
     92      REAL(wp), DIMENSION(jpi,jpj,jpk) ::  zvn   ! 3D workspace 
     93 
     94 
     95      CHARACTER( len = 12 )  :: cl1 
    9296      !!---------------------------------------------------------------------- 
    9397      ! 
    9498      IF( ln_timing )   CALL timing_start('dia_ptr') 
    9599 
    96       IF( kt == nit000 .AND. ll_init )   CALL dia_ptr_init 
    97       ! 
    98       IF( .NOT. l_diaptr )   RETURN 
    99  
     100      ! 
    100101      IF( PRESENT( pvtr ) ) THEN 
    101          IF( iom_use( 'zomsf' ) ) THEN    ! effective MSF 
    102             DO jn = 1, nptr                                    ! by sub-basins 
    103                z4d1(1,:,:,jn) =  ptr_sjk( pvtr(:,:,:), btmsk34(:,:,jn) )  ! zonal cumulative effective transport excluding closed seas 
    104                DO jk = jpkm1, 1, -1  
    105                   z4d1(1,:,jk,jn) = z4d1(1,:,jk+1,jn) - z4d1(1,:,jk,jn)    ! effective j-Stream-Function (MSF) 
     102         IF( iom_use("zomsfglo") ) THEN    ! effective MSF 
     103            z3d(1,:,:) = ptr_sjk( pvtr(:,:,:) )  ! zonal cumulative effective transport 
     104            DO jk = 2, jpkm1  
     105              z3d(1,:,jk) = z3d(1,:,jk-1) + z3d(1,:,jk)   ! effective j-Stream-Function (MSF) 
     106            END DO 
     107            DO ji = 1, jpi 
     108               z3d(ji,:,:) = z3d(1,:,:) 
     109            ENDDO 
     110            cl1 = TRIM('zomsf'//clsubb(1) ) 
     111            CALL iom_put( cl1, z3d * rc_sv ) 
     112            DO jn = 2, nptr                                    ! by sub-basins 
     113               z3d(1,:,:) =  ptr_sjk( pvtr(:,:,:), btmsk(:,:,jn)*btm30(:,:) )  
     114               DO jk = 2, jpkm1  
     115                  z3d(1,:,jk) = z3d(1,:,jk-1) + z3d(1,:,jk)    ! effective j-Stream-Function (MSF) 
    106116               END DO 
    107117               DO ji = 1, jpi 
    108                   z4d1(ji,:,:,jn) = z4d1(1,:,:,jn) 
    109                ENDDO 
    110             END DO 
    111             CALL iom_put( 'zomsf', z4d1 * rc_sv ) 
    112          ENDIF 
    113          IF(  iom_use( 'sopstove' ) .OR. iom_use( 'sophtove' ) .OR.   & 
    114             & iom_use( 'sopstbtr' ) .OR. iom_use( 'sophtbtr' ) ) THEN 
     118                  z3d(ji,:,:) = z3d(1,:,:) 
     119               ENDDO 
     120               cl1 = TRIM('zomsf'//clsubb(jn) ) 
     121               CALL iom_put( cl1, z3d * rc_sv ) 
     122            END DO 
     123         ENDIF 
     124         IF( iom_use("sopstove") .OR. iom_use("sophtove") .OR. iom_use("sopstbtr") .OR. iom_use("sophtbtr") ) THEN 
    115125            ! define fields multiplied by scalar 
    116126            zmask(:,:,:) = 0._wp 
    117127            zts(:,:,:,:) = 0._wp 
     128            zvn(:,:,:) = 0._wp 
    118129            DO jk = 1, jpkm1 
    119130               DO jj = 1, jpjm1 
     
    123134                     zts(ji,jj,jk,jp_tem) = (tsn(ji,jj,jk,jp_tem)+tsn(ji,jj+1,jk,jp_tem)) * 0.5 * zvfc  !Tracers averaged onto V grid 
    124135                     zts(ji,jj,jk,jp_sal) = (tsn(ji,jj,jk,jp_sal)+tsn(ji,jj+1,jk,jp_sal)) * 0.5 * zvfc 
     136                     zvn(ji,jj,jk)        = vn(ji,jj,jk)         * zvfc 
    125137                  ENDDO 
    126138               ENDDO 
    127139             ENDDO 
    128140         ENDIF 
    129          IF( iom_use( 'sopstove' ) .OR. iom_use( 'sophtove' ) ) THEN 
    130             DO jn = 1, nptr 
    131                sjk(:,:,jn) = ptr_sjk( zmask(:,:,:), btmsk(:,:,jn) ) 
    132                r1_sjk(:,:,jn) = 0._wp 
    133                WHERE( sjk(:,:,jn) /= 0._wp )   r1_sjk(:,:,jn) = 1._wp / sjk(:,:,jn) 
    134                ! i-mean T and S, j-Stream-Function, basin 
    135                zt_jk(:,:,jn) = ptr_sjk( zts(:,:,:,jp_tem), btmsk(:,:,jn) ) * r1_sjk(:,:,jn) 
    136                zs_jk(:,:,jn) = ptr_sjk( zts(:,:,:,jp_sal), btmsk(:,:,jn) ) * r1_sjk(:,:,jn) 
    137                v_msf(:,:,jn) = ptr_sjk( pvtr(:,:,:), btmsk34(:,:,jn) )  
    138                hstr_ove(:,jp_tem,jn) = SUM( v_msf(:,:,jn)*zt_jk(:,:,jn), 2 ) 
    139                hstr_ove(:,jp_sal,jn) = SUM( v_msf(:,:,jn)*zs_jk(:,:,jn), 2 ) 
    140                ! 
    141             ENDDO 
    142             DO jn = 1, nptr 
    143                z3dtr(1,:,jn) = hstr_ove(:,jp_tem,jn) * rc_pwatt  !  (conversion in PW) 
    144                DO ji = 1, jpi 
    145                   z3dtr(ji,:,jn) = z3dtr(1,:,jn) 
    146                ENDDO 
    147             ENDDO 
    148             CALL iom_put( 'sophtove', z3dtr ) 
    149             DO jn = 1, nptr 
    150                z3dtr(1,:,jn) = hstr_ove(:,jp_sal,jn) * rc_ggram !  (conversion in Gg) 
    151                DO ji = 1, jpi 
    152                   z3dtr(ji,:,jn) = z3dtr(1,:,jn) 
    153                ENDDO 
    154             ENDDO 
    155             CALL iom_put( 'sopstove', z3dtr ) 
    156          ENDIF 
    157  
    158          IF( iom_use( 'sopstbtr' ) .OR. iom_use( 'sophtbtr' ) ) THEN 
    159             ! Calculate barotropic heat and salt transport here  
    160             DO jn = 1, nptr 
    161                sjk(:,1,jn) = ptr_sj( zmask(:,:,:), btmsk(:,:,jn) ) 
    162                r1_sjk(:,1,jn) = 0._wp 
    163                WHERE( sjk(:,1,jn) /= 0._wp )   r1_sjk(:,1,jn) = 1._wp / sjk(:,1,jn) 
    164                ! 
    165                zvsum(:) = ptr_sj( pvtr(:,:,:), btmsk34(:,:,jn) ) 
    166                ztsum(:) = ptr_sj( zts(:,:,:,jp_tem), btmsk(:,:,jn) ) 
    167                zssum(:) = ptr_sj( zts(:,:,:,jp_sal), btmsk(:,:,jn) ) 
    168                hstr_btr(:,jp_tem,jn) = zvsum(:) * ztsum(:) * r1_sjk(:,1,jn) 
    169                hstr_btr(:,jp_sal,jn) = zvsum(:) * zssum(:) * r1_sjk(:,1,jn) 
    170                ! 
    171             ENDDO 
    172             DO jn = 1, nptr 
    173                z3dtr(1,:,jn) = hstr_btr(:,jp_tem,jn) * rc_pwatt  !  (conversion in PW) 
    174                DO ji = 1, jpi 
    175                   z3dtr(ji,:,jn) = z3dtr(1,:,jn) 
    176                ENDDO 
    177             ENDDO 
    178             CALL iom_put( 'sophtbtr', z3dtr ) 
    179             DO jn = 1, nptr 
    180                z3dtr(1,:,jn) = hstr_btr(:,jp_sal,jn) * rc_ggram !  (conversion in Gg) 
    181                DO ji = 1, jpi 
    182                   z3dtr(ji,:,jn) = z3dtr(1,:,jn) 
    183                ENDDO 
    184             ENDDO 
    185             CALL iom_put( 'sopstbtr', z3dtr ) 
    186          ENDIF  
     141         IF( iom_use("sopstove") .OR. iom_use("sophtove") ) THEN 
     142             sjk(:,:,1) = ptr_sjk( zmask(:,:,:), btmsk(:,:,1) ) 
     143             r1_sjk(:,:,1) = 0._wp 
     144             WHERE( sjk(:,:,1) /= 0._wp )   r1_sjk(:,:,1) = 1._wp / sjk(:,:,1) 
     145 
     146             ! i-mean T and S, j-Stream-Function, global 
     147             tn_jk(:,:,1) = ptr_sjk( zts(:,:,:,jp_tem) ) * r1_sjk(:,:,1) 
     148             sn_jk(:,:,1) = ptr_sjk( zts(:,:,:,jp_sal) ) * r1_sjk(:,:,1) 
     149             v_msf(:,:,1) = ptr_sjk( zvn(:,:,:) ) 
     150 
     151             htr_ove(:,1) = SUM( v_msf(:,:,1)*tn_jk(:,:,1) ,2 ) 
     152             str_ove(:,1) = SUM( v_msf(:,:,1)*sn_jk(:,:,1) ,2 ) 
     153 
     154             z2d(1,:) = htr_ove(:,1) * rc_pwatt        !  (conversion in PW) 
     155             DO ji = 1, jpi 
     156               z2d(ji,:) = z2d(1,:) 
     157             ENDDO 
     158             cl1 = 'sophtove' 
     159             CALL iom_put( TRIM(cl1), z2d ) 
     160             z2d(1,:) = str_ove(:,1) * rc_ggram        !  (conversion in Gg) 
     161             DO ji = 1, jpi 
     162               z2d(ji,:) = z2d(1,:) 
     163             ENDDO 
     164             cl1 = 'sopstove' 
     165             CALL iom_put( TRIM(cl1), z2d ) 
     166             IF( ln_subbas ) THEN 
     167                DO jn = 2, nptr 
     168                    sjk(:,:,jn) = ptr_sjk( zmask(:,:,:), btmsk(:,:,jn) ) 
     169                    r1_sjk(:,:,jn) = 0._wp 
     170                    WHERE( sjk(:,:,jn) /= 0._wp )   r1_sjk(:,:,jn) = 1._wp / sjk(:,:,jn) 
     171 
     172                    ! i-mean T and S, j-Stream-Function, basin 
     173                    tn_jk(:,:,jn) = ptr_sjk( zts(:,:,:,jp_tem), btmsk(:,:,jn) ) * r1_sjk(:,:,jn) 
     174                    sn_jk(:,:,jn) = ptr_sjk( zts(:,:,:,jp_sal), btmsk(:,:,jn) ) * r1_sjk(:,:,jn) 
     175                    v_msf(:,:,jn) = ptr_sjk( zvn(:,:,:), btmsk(:,:,jn) )  
     176                    htr_ove(:,jn) = SUM( v_msf(:,:,jn)*tn_jk(:,:,jn) ,2 ) 
     177                    str_ove(:,jn) = SUM( v_msf(:,:,jn)*sn_jk(:,:,jn) ,2 ) 
     178 
     179                    z2d(1,:) = htr_ove(:,jn) * rc_pwatt !  (conversion in PW) 
     180                    DO ji = 1, jpi 
     181                        z2d(ji,:) = z2d(1,:) 
     182                    ENDDO 
     183                    cl1 = TRIM('sophtove_'//clsubb(jn)) 
     184                    CALL iom_put( cl1, z2d ) 
     185                    z2d(1,:) = str_ove(:,jn) * rc_ggram        ! (conversion in Gg) 
     186                    DO ji = 1, jpi 
     187                        z2d(ji,:) = z2d(1,:) 
     188                    ENDDO 
     189                    cl1 = TRIM('sopstove_'//clsubb(jn)) 
     190                    CALL iom_put( cl1, z2d ) 
     191                END DO 
     192             ENDIF 
     193         ENDIF 
     194         IF( iom_use("sopstbtr") .OR. iom_use("sophtbtr") ) THEN 
     195         ! Calculate barotropic heat and salt transport here  
     196             sjk(:,1,1) = ptr_sj( zmask(:,:,:), btmsk(:,:,1) ) 
     197             r1_sjk(:,1,1) = 0._wp 
     198             WHERE( sjk(:,1,1) /= 0._wp )   r1_sjk(:,1,1) = 1._wp / sjk(:,1,1) 
     199             
     200            vsum = ptr_sj( zvn(:,:,:), btmsk(:,:,1)) 
     201            tssum(:,jp_tem) = ptr_sj( zts(:,:,:,jp_tem), btmsk(:,:,1) ) 
     202            tssum(:,jp_sal) = ptr_sj( zts(:,:,:,jp_sal), btmsk(:,:,1) ) 
     203            htr_btr(:,1) = vsum * tssum(:,jp_tem) * r1_sjk(:,1,1) 
     204            str_btr(:,1) = vsum * tssum(:,jp_sal) * r1_sjk(:,1,1) 
     205            z2d(1,:) = htr_btr(:,1) * rc_pwatt        !  (conversion in PW) 
     206            DO ji = 2, jpi 
     207               z2d(ji,:) = z2d(1,:) 
     208            ENDDO 
     209            cl1 = 'sophtbtr' 
     210            CALL iom_put( TRIM(cl1), z2d ) 
     211            z2d(1,:) = str_btr(:,1) * rc_ggram        !  (conversion in Gg) 
     212            DO ji = 2, jpi 
     213              z2d(ji,:) = z2d(1,:) 
     214            ENDDO 
     215            cl1 = 'sopstbtr' 
     216            CALL iom_put( TRIM(cl1), z2d ) 
     217            IF( ln_subbas ) THEN 
     218                DO jn = 2, nptr 
     219                    sjk(:,1,jn) = ptr_sj( zmask(:,:,:), btmsk(:,:,jn) ) 
     220                    r1_sjk(:,1,jn) = 0._wp 
     221                    WHERE( sjk(:,1,jn) /= 0._wp )   r1_sjk(:,1,jn) = 1._wp / sjk(:,1,jn) 
     222                    vsum = ptr_sj( zvn(:,:,:), btmsk(:,:,jn)) 
     223                    tssum(:,jp_tem) = ptr_sj( zts(:,:,:,jp_tem), btmsk(:,:,jn) ) 
     224                    tssum(:,jp_sal) = ptr_sj( zts(:,:,:,jp_sal), btmsk(:,:,jn) ) 
     225                    htr_btr(:,jn) = vsum * tssum(:,jp_tem) * r1_sjk(:,1,jn) 
     226                    str_btr(:,jn) = vsum * tssum(:,jp_sal) * r1_sjk(:,1,jn) 
     227                    z2d(1,:) = htr_btr(:,jn) * rc_pwatt !  (conversion in PW) 
     228                    DO ji = 1, jpi 
     229                        z2d(ji,:) = z2d(1,:) 
     230                    ENDDO 
     231                    cl1 = TRIM('sophtbtr_'//clsubb(jn)) 
     232                    CALL iom_put( cl1, z2d ) 
     233                    z2d(1,:) = str_btr(:,jn) * rc_ggram        ! (conversion in Gg) 
     234                    DO ji = 1, jpi 
     235                        z2d(ji,:) = z2d(1,:) 
     236                    ENDDO 
     237                    cl1 = TRIM('sopstbtr_'//clsubb(jn)) 
     238                    CALL iom_put( cl1, z2d ) 
     239               ENDDO 
     240            ENDIF !ln_subbas 
     241         ENDIF !iom_use("sopstbtr....) 
    187242         ! 
    188243      ELSE 
    189244         ! 
    190          zmask(:,:,:) = 0._wp 
    191          zts(:,:,:,:) = 0._wp 
    192          IF( iom_use( 'zotem' ) .OR. iom_use( 'zosal' ) .OR. iom_use( 'zosrf' )  ) THEN    ! i-mean i-k-surface  
     245         IF( iom_use("zotemglo") ) THEN    ! i-mean i-k-surface  
    193246            DO jk = 1, jpkm1 
    194247               DO jj = 1, jpj 
     
    201254               END DO 
    202255            END DO 
    203             ! 
    204256            DO jn = 1, nptr 
    205257               zmask(1,:,:) = ptr_sjk( zmask(:,:,:), btmsk(:,:,jn) ) 
    206                z4d1(:,:,:,jn) = zmask(:,:,:) 
    207             ENDDO 
    208             CALL iom_put( 'zosrf', z4d1 ) 
    209             ! 
    210             DO jn = 1, nptr 
    211                z4d2(1,:,:,jn) = ptr_sjk( zts(:,:,:,jp_tem), btmsk(:,:,jn) ) & 
    212                   &            / MAX( z4d1(1,:,:,jn), 10.e-15 ) 
     258               cl1 = TRIM('zosrf'//clsubb(jn) ) 
     259               CALL iom_put( cl1, zmask ) 
     260               ! 
     261               z3d(1,:,:) = ptr_sjk( zts(:,:,:,jp_tem), btmsk(:,:,jn) ) & 
     262                  &            / MAX( zmask(1,:,:), 10.e-15 ) 
    213263               DO ji = 1, jpi 
    214                   z4d2(ji,:,:,jn) = z4d2(1,:,:,jn) 
    215                ENDDO 
    216             ENDDO 
    217             CALL iom_put( 'zotem', z4d2 ) 
    218             ! 
    219             DO jn = 1, nptr 
    220                z4d2(1,:,:,jn) = ptr_sjk( zts(:,:,:,jp_sal), btmsk(:,:,jn) ) & 
    221                   &            / MAX( z4d1(1,:,:,jn), 10.e-15 ) 
     264                  z3d(ji,:,:) = z3d(1,:,:) 
     265               ENDDO 
     266               cl1 = TRIM('zotem'//clsubb(jn) ) 
     267               CALL iom_put( cl1, z3d ) 
     268               ! 
     269               z3d(1,:,:) = ptr_sjk( zts(:,:,:,jp_sal), btmsk(:,:,jn) ) & 
     270                  &            / MAX( zmask(1,:,:), 10.e-15 ) 
    222271               DO ji = 1, jpi 
    223                   z4d2(ji,:,:,jn) = z4d2(1,:,:,jn) 
    224                ENDDO 
    225             ENDDO 
    226             CALL iom_put( 'zosal', z4d2 ) 
    227             ! 
     272                  z3d(ji,:,:) = z3d(1,:,:) 
     273               ENDDO 
     274               cl1 = TRIM('zosal'//clsubb(jn) ) 
     275               CALL iom_put( cl1, z3d ) 
     276            END DO 
    228277         ENDIF 
    229278         ! 
    230279         !                                ! Advective and diffusive heat and salt transport 
    231          IF( iom_use( 'sophtadv' ) .OR. iom_use( 'sopstadv' ) ) THEN   
    232             !  
    233             DO jn = 1, nptr 
    234                z3dtr(1,:,jn) = hstr_adv(:,jp_tem,jn) * rc_pwatt  !  (conversion in PW) 
     280         IF( iom_use("sophtadv") .OR. iom_use("sopstadv") ) THEN    
     281            z2d(1,:) = htr_adv(:,1) * rc_pwatt        !  (conversion in PW) 
     282            DO ji = 1, jpi 
     283               z2d(ji,:) = z2d(1,:) 
     284            ENDDO 
     285            cl1 = 'sophtadv'                  
     286            CALL iom_put( TRIM(cl1), z2d ) 
     287            z2d(1,:) = str_adv(:,1) * rc_ggram        ! (conversion in Gg) 
     288            DO ji = 1, jpi 
     289               z2d(ji,:) = z2d(1,:) 
     290            ENDDO 
     291            cl1 = 'sopstadv' 
     292            CALL iom_put( TRIM(cl1), z2d ) 
     293            IF( ln_subbas ) THEN 
     294              DO jn=2,nptr 
     295               z2d(1,:) = htr_adv(:,jn) * rc_pwatt        !  (conversion in PW) 
    235296               DO ji = 1, jpi 
    236                   z3dtr(ji,:,jn) = z3dtr(1,:,jn) 
    237                ENDDO 
    238             ENDDO 
    239             CALL iom_put( 'sophtadv', z3dtr ) 
    240             DO jn = 1, nptr 
    241                z3dtr(1,:,jn) = hstr_adv(:,jp_sal,jn) * rc_ggram !  (conversion in Gg) 
     297                 z2d(ji,:) = z2d(1,:) 
     298               ENDDO 
     299               cl1 = TRIM('sophtadv_'//clsubb(jn))                  
     300               CALL iom_put( cl1, z2d ) 
     301               z2d(1,:) = str_adv(:,jn) * rc_ggram        ! (conversion in Gg) 
    242302               DO ji = 1, jpi 
    243                   z3dtr(ji,:,jn) = z3dtr(1,:,jn) 
    244                ENDDO 
    245             ENDDO 
    246             CALL iom_put( 'sopstadv', z3dtr ) 
    247          ENDIF 
    248          ! 
    249          IF( iom_use( 'sophtldf' ) .OR. iom_use( 'sopstldf' ) ) THEN   
    250             !  
    251             DO jn = 1, nptr 
    252                z3dtr(1,:,jn) = hstr_ldf(:,jp_tem,jn) * rc_pwatt  !  (conversion in PW) 
     303                  z2d(ji,:) = z2d(1,:) 
     304               ENDDO 
     305               cl1 = TRIM('sopstadv_'//clsubb(jn))                  
     306               CALL iom_put( cl1, z2d )               
     307              ENDDO 
     308            ENDIF 
     309         ENDIF 
     310         ! 
     311         IF( iom_use("sophtldf") .OR. iom_use("sopstldf") ) THEN    
     312            z2d(1,:) = htr_ldf(:,1) * rc_pwatt        !  (conversion in PW)  
     313            DO ji = 1, jpi 
     314               z2d(ji,:) = z2d(1,:) 
     315            ENDDO 
     316            cl1 = 'sophtldf' 
     317            CALL iom_put( TRIM(cl1), z2d ) 
     318            z2d(1,:) = str_ldf(:,1) * rc_ggram        !  (conversion in Gg) 
     319            DO ji = 1, jpi 
     320               z2d(ji,:) = z2d(1,:) 
     321            ENDDO 
     322            cl1 = 'sopstldf' 
     323            CALL iom_put( TRIM(cl1), z2d ) 
     324            IF( ln_subbas ) THEN 
     325              DO jn=2,nptr 
     326               z2d(1,:) = htr_ldf(:,jn) * rc_pwatt        !  (conversion in PW) 
    253327               DO ji = 1, jpi 
    254                   z3dtr(ji,:,jn) = z3dtr(1,:,jn) 
    255                ENDDO 
    256             ENDDO 
    257             CALL iom_put( 'sophtldf', z3dtr ) 
    258             DO jn = 1, nptr 
    259                z3dtr(1,:,jn) = hstr_ldf(:,jp_sal,jn) * rc_ggram !  (conversion in Gg) 
     328                 z2d(ji,:) = z2d(1,:) 
     329               ENDDO 
     330               cl1 = TRIM('sophtldf_'//clsubb(jn))                  
     331               CALL iom_put( cl1, z2d ) 
     332               z2d(1,:) = str_ldf(:,jn) * rc_ggram        ! (conversion in Gg) 
    260333               DO ji = 1, jpi 
    261                   z3dtr(ji,:,jn) = z3dtr(1,:,jn) 
    262                ENDDO 
    263             ENDDO 
    264             CALL iom_put( 'sopstldf', z3dtr ) 
    265          ENDIF 
    266          ! 
    267          IF( iom_use( 'sophteiv' ) .OR. iom_use( 'sopsteiv' ) ) THEN   
    268             !  
    269             DO jn = 1, nptr 
    270                z3dtr(1,:,jn) = hstr_eiv(:,jp_tem,jn) * rc_pwatt  !  (conversion in PW) 
    271                DO ji = 1, jpi 
    272                   z3dtr(ji,:,jn) = z3dtr(1,:,jn) 
    273                ENDDO 
    274             ENDDO 
    275             CALL iom_put( 'sophteiv', z3dtr ) 
    276             DO jn = 1, nptr 
    277                z3dtr(1,:,jn) = hstr_eiv(:,jp_sal,jn) * rc_ggram !  (conversion in Gg) 
    278                DO ji = 1, jpi 
    279                   z3dtr(ji,:,jn) = z3dtr(1,:,jn) 
    280                ENDDO 
    281             ENDDO 
    282             CALL iom_put( 'sopsteiv', z3dtr ) 
    283          ENDIF 
    284          ! 
    285          IF( iom_use( 'sopstvtr' ) .OR. iom_use( 'sophtvtr' ) ) THEN 
    286             zts(:,:,:,:) = 0._wp 
    287             DO jk = 1, jpkm1 
    288                DO jj = 1, jpjm1 
     334                  z2d(ji,:) = z2d(1,:) 
     335               ENDDO 
     336               cl1 = TRIM('sopstldf_'//clsubb(jn))                  
     337               CALL iom_put( cl1, z2d )               
     338              ENDDO 
     339            ENDIF 
     340         ENDIF 
     341 
     342         IF( iom_use("sophteiv") .OR. iom_use("sopsteiv") ) THEN  
     343            z2d(1,:) = htr_eiv(:,1) * rc_pwatt        !  (conversion in PW)  
     344            DO ji = 1, jpi 
     345               z2d(ji,:) = z2d(1,:) 
     346            ENDDO 
     347            cl1 = 'sophteiv' 
     348            CALL iom_put( TRIM(cl1), z2d ) 
     349            z2d(1,:) = str_eiv(:,1) * rc_ggram        !  (conversion in Gg) 
     350            DO ji = 1, jpi 
     351               z2d(ji,:) = z2d(1,:) 
     352            ENDDO 
     353            cl1 = 'sopsteiv' 
     354            CALL iom_put( TRIM(cl1), z2d ) 
     355            IF( ln_subbas ) THEN 
     356               DO jn=2,nptr 
     357                  z2d(1,:) = htr_eiv(:,jn) * rc_pwatt        !  (conversion in PW) 
    289358                  DO ji = 1, jpi 
    290                      zvfc = e1v(ji,jj) * e3v_n(ji,jj,jk) 
    291                      zts(ji,jj,jk,jp_tem) = (tsn(ji,jj,jk,jp_tem)+tsn(ji,jj+1,jk,jp_tem)) * 0.5 * zvfc  !Tracers averaged onto V grid 
    292                      zts(ji,jj,jk,jp_sal) = (tsn(ji,jj,jk,jp_sal)+tsn(ji,jj+1,jk,jp_sal)) * 0.5 * zvfc 
     359                     z2d(ji,:) = z2d(1,:) 
    293360                  ENDDO 
    294                ENDDO 
    295              ENDDO 
    296              CALL dia_ptr_hst( jp_tem, 'vtr', zts(:,:,:,jp_tem) ) 
    297              CALL dia_ptr_hst( jp_sal, 'vtr', zts(:,:,:,jp_sal) ) 
    298              DO jn = 1, nptr 
    299                 z3dtr(1,:,jn) = hstr_vtr(:,jp_tem,jn) * rc_pwatt  !  (conversion in PW) 
    300                 DO ji = 1, jpi 
    301                    z3dtr(ji,:,jn) = z3dtr(1,:,jn) 
    302                 ENDDO 
    303              ENDDO 
    304              CALL iom_put( 'sophtvtr', z3dtr ) 
    305              DO jn = 1, nptr 
    306                z3dtr(1,:,jn) = hstr_vtr(:,jp_sal,jn) * rc_ggram !  (conversion in Gg) 
    307                DO ji = 1, jpi 
    308                   z3dtr(ji,:,jn) = z3dtr(1,:,jn) 
    309                ENDDO 
    310             ENDDO 
    311             CALL iom_put( 'sopstvtr', z3dtr ) 
    312          ENDIF 
    313          ! 
    314          IF( iom_use( 'uocetr_vsum_cumul' ) ) THEN 
    315             CALL iom_get_var(  'uocetr_vsum_op', z2d ) ! get uocetr_vsum_op from xml 
    316             z2d(:,:) = ptr_ci_2d( z2d(:,:) )   
    317             CALL iom_put( 'uocetr_vsum_cumul', z2d ) 
     361                  cl1 = TRIM('sophteiv_'//clsubb(jn))                  
     362                  CALL iom_put( cl1, z2d ) 
     363                  z2d(1,:) = str_eiv(:,jn) * rc_ggram        ! (conversion in Gg) 
     364                  DO ji = 1, jpi 
     365                     z2d(ji,:) = z2d(1,:) 
     366                  ENDDO 
     367                  cl1 = TRIM('sopsteiv_'//clsubb(jn))  
     368                  CALL iom_put( cl1, z2d )               
     369               ENDDO 
     370            ENDIF 
    318371         ENDIF 
    319372         ! 
     
    331384      !! ** Purpose :   Initialization, namelist read 
    332385      !!---------------------------------------------------------------------- 
    333       INTEGER ::  inum, jn           ! local integers 
    334       !! 
    335       REAL(wp), DIMENSION(jpi,jpj) :: zmsk 
    336       !!---------------------------------------------------------------------- 
    337  
    338       l_diaptr = .FALSE. 
    339       IF(   iom_use( 'zomsf'    ) .OR. iom_use( 'zotem'    ) .OR. iom_use( 'zosal'    ) .OR.  & 
    340          &  iom_use( 'zosrf'    ) .OR. iom_use( 'sopstove' ) .OR. iom_use( 'sophtove' ) .OR.  & 
    341          &  iom_use( 'sopstbtr' ) .OR. iom_use( 'sophtbtr' ) .OR. iom_use( 'sophtadv' ) .OR.  & 
    342          &  iom_use( 'sopstadv' ) .OR. iom_use( 'sophtldf' ) .OR. iom_use( 'sopstldf' ) .OR.  &  
    343          &  iom_use( 'sophteiv' ) .OR. iom_use( 'sopsteiv' ) .OR. iom_use( 'sopstvtr' ) .OR.  & 
    344          &  iom_use( 'sophtvtr' ) .OR. iom_use( 'uocetr_vsum_cumul' ) )  l_diaptr  = .TRUE. 
    345  
    346   
     386      INTEGER ::  jn           ! local integers 
     387      INTEGER ::  inum, ierr   ! local integers 
     388      INTEGER ::  ios          ! Local integer output status for namelist read 
     389      !! 
     390      NAMELIST/namptr/ ln_diaptr, ln_subbas 
     391      !!---------------------------------------------------------------------- 
     392 
     393      REWIND( numnam_ref )              ! Namelist namptr in reference namelist : Poleward transport 
     394      READ  ( numnam_ref, namptr, IOSTAT = ios, ERR = 901) 
     395901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namptr in reference namelist' ) 
     396 
     397      REWIND( numnam_cfg )              ! Namelist namptr in configuration namelist : Poleward transport 
     398      READ  ( numnam_cfg, namptr, IOSTAT = ios, ERR = 902 ) 
     399902   IF( ios >  0 ) CALL ctl_nam ( ios , 'namptr in configuration namelist' ) 
     400      IF(lwm) WRITE ( numond, namptr ) 
     401 
    347402      IF(lwp) THEN                     ! Control print 
    348403         WRITE(numout,*) 
     
    350405         WRITE(numout,*) '~~~~~~~~~~~~' 
    351406         WRITE(numout,*) '   Namelist namptr : set ptr parameters' 
    352          WRITE(numout,*) '      Poleward heat & salt transport (T) or not (F)      l_diaptr  = ', l_diaptr 
     407         WRITE(numout,*) '      Poleward heat & salt transport (T) or not (F)      ln_diaptr  = ', ln_diaptr 
     408         WRITE(numout,*) '      Global (F) or glo/Atl/Pac/Ind/Indo-Pac basins      ln_subbas  = ', ln_subbas 
    353409      ENDIF 
    354410 
    355       IF( l_diaptr ) THEN   
    356          ! 
     411      IF( ln_diaptr ) THEN   
     412         ! 
     413         IF( ln_subbas ) THEN  
     414            nptr = 5            ! Global, Atlantic, Pacific, Indian, Indo-Pacific 
     415            ALLOCATE( clsubb(nptr) ) 
     416            clsubb(1) = 'glo' ;  clsubb(2) = 'atl'  ;  clsubb(3) = 'pac'  ;  clsubb(4) = 'ind'  ;  clsubb(5) = 'ipc' 
     417         ELSE                
     418            nptr = 1       ! Global only 
     419            ALLOCATE( clsubb(nptr) ) 
     420            clsubb(1) = 'glo'  
     421         ENDIF 
     422 
     423         !                                      ! allocate dia_ptr arrays 
    357424         IF( dia_ptr_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'dia_ptr_init : unable to allocate arrays' ) 
    358425 
    359426         rc_pwatt = rc_pwatt * rau0_rcp          ! conversion from K.s-1 to PetaWatt 
    360          rc_ggram = rc_ggram * rau0              ! conversion from m3/s to Gg/s 
    361427 
    362428         IF( lk_mpp )   CALL mpp_ini_znl( numout )     ! Define MPI communicator for zonal sum 
    363429 
    364          btmsk(:,:,1) = tmask_i(:,:)                  
    365          CALL iom_open( 'subbasins', inum,  ldstop = .FALSE.  ) 
    366          CALL iom_get( inum, jpdom_data, 'atlmsk', btmsk(:,:,2) )   ! Atlantic basin 
    367          CALL iom_get( inum, jpdom_data, 'pacmsk', btmsk(:,:,3) )   ! Pacific  basin 
    368          CALL iom_get( inum, jpdom_data, 'indmsk', btmsk(:,:,4) )   ! Indian   basin 
    369          CALL iom_close( inum ) 
    370          btmsk(:,:,5) = MAX ( btmsk(:,:,3), btmsk(:,:,4) )          ! Indo-Pacific basin 
    371          DO jn = 2, nptr 
     430         IF( ln_subbas ) THEN                ! load sub-basin mask 
     431            CALL iom_open( 'subbasins', inum,  ldstop = .FALSE.  ) 
     432            CALL iom_get( inum, jpdom_data, 'atlmsk', btmsk(:,:,2) )   ! Atlantic basin 
     433            CALL iom_get( inum, jpdom_data, 'pacmsk', btmsk(:,:,3) )   ! Pacific  basin 
     434            CALL iom_get( inum, jpdom_data, 'indmsk', btmsk(:,:,4) )   ! Indian   basin 
     435            CALL iom_close( inum ) 
     436            btmsk(:,:,5) = MAX ( btmsk(:,:,3), btmsk(:,:,4) )          ! Indo-Pacific basin 
     437            WHERE( gphit(:,:) < -30._wp)   ;   btm30(:,:) = 0._wp      ! mask out Southern Ocean 
     438            ELSE WHERE                     ;   btm30(:,:) = ssmask(:,:) 
     439            END WHERE 
     440         ENDIF 
     441    
     442         btmsk(:,:,1) = tmask_i(:,:)                                   ! global ocean 
     443       
     444         DO jn = 1, nptr 
    372445            btmsk(:,:,jn) = btmsk(:,:,jn) * tmask_i(:,:)               ! interior domain only 
    373446         END DO 
    374          ! JD : modification so that overturning streamfunction is available in Atlantic at 34S to compare with observations 
    375          WHERE( gphit(:,:)*tmask_i(:,:) < -34._wp) 
    376            zmsk(:,:) = 0._wp      ! mask out Southern Ocean 
    377          ELSE WHERE                   
    378            zmsk(:,:) = ssmask(:,:) 
    379          END WHERE 
    380          btmsk34(:,:,1) = btmsk(:,:,1)                  
    381          DO jn = 2, nptr 
    382             btmsk34(:,:,jn) = btmsk(:,:,jn) * zmsk(:,:)               ! interior domain only 
    383          ENDDO 
    384447 
    385448         ! Initialise arrays to zero because diatpr is called before they are first calculated 
    386449         ! Note that this means diagnostics will not be exactly correct when model run is restarted. 
    387          hstr_adv(:,:,:) = 0._wp            
    388          hstr_ldf(:,:,:) = 0._wp            
    389          hstr_eiv(:,:,:) = 0._wp            
    390          hstr_ove(:,:,:) = 0._wp            
    391          hstr_btr(:,:,:) = 0._wp           ! 
    392          hstr_vtr(:,:,:) = 0._wp           ! 
    393          ! 
    394          ll_init = .FALSE. 
     450         htr_adv(:,:) = 0._wp  ;  str_adv(:,:) =  0._wp  
     451         htr_ldf(:,:) = 0._wp  ;  str_ldf(:,:) =  0._wp  
     452         htr_eiv(:,:) = 0._wp  ;  str_eiv(:,:) =  0._wp  
     453         htr_ove(:,:) = 0._wp  ;   str_ove(:,:) =  0._wp 
     454         htr_btr(:,:) = 0._wp  ;   str_btr(:,:) =  0._wp 
    395455         ! 
    396456      ENDIF  
     
    411471      INTEGER                                        :: jn    ! 
    412472 
    413       ! 
    414473      IF( cptr == 'adv' ) THEN 
    415          IF( ktra == jp_tem )  THEN 
    416              DO jn = 1, nptr 
    417                 hstr_adv(:,jp_tem,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 
    418              ENDDO 
    419          ENDIF 
    420          IF( ktra == jp_sal )  THEN 
    421              DO jn = 1, nptr 
    422                 hstr_adv(:,jp_sal,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 
    423              ENDDO 
    424          ENDIF 
     474         IF( ktra == jp_tem )  htr_adv(:,1) = ptr_sj( pva(:,:,:) ) 
     475         IF( ktra == jp_sal )  str_adv(:,1) = ptr_sj( pva(:,:,:) ) 
    425476      ENDIF 
    426       ! 
    427477      IF( cptr == 'ldf' ) THEN 
    428          IF( ktra == jp_tem )  THEN 
    429              DO jn = 1, nptr 
    430                 hstr_ldf(:,jp_tem,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 
    431              ENDDO 
    432          ENDIF 
    433          IF( ktra == jp_sal )  THEN 
    434              DO jn = 1, nptr 
    435                 hstr_ldf(:,jp_sal,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 
    436              ENDDO 
    437          ENDIF 
     478         IF( ktra == jp_tem )  htr_ldf(:,1) = ptr_sj( pva(:,:,:) ) 
     479         IF( ktra == jp_sal )  str_ldf(:,1) = ptr_sj( pva(:,:,:) ) 
    438480      ENDIF 
    439       ! 
    440481      IF( cptr == 'eiv' ) THEN 
    441          IF( ktra == jp_tem )  THEN 
    442              DO jn = 1, nptr 
    443                 hstr_eiv(:,jp_tem,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 
    444              ENDDO 
    445          ENDIF 
    446          IF( ktra == jp_sal )  THEN 
    447              DO jn = 1, nptr 
    448                 hstr_eiv(:,jp_sal,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 
    449              ENDDO 
    450          ENDIF 
     482         IF( ktra == jp_tem )  htr_eiv(:,1) = ptr_sj( pva(:,:,:) ) 
     483         IF( ktra == jp_sal )  str_eiv(:,1) = ptr_sj( pva(:,:,:) ) 
    451484      ENDIF 
    452485      ! 
    453       IF( cptr == 'vtr' ) THEN 
    454          IF( ktra == jp_tem )  THEN 
    455              DO jn = 1, nptr 
    456                 hstr_vtr(:,jp_tem,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 
    457              ENDDO 
    458          ENDIF 
    459          IF( ktra == jp_sal )  THEN 
    460              DO jn = 1, nptr 
    461                 hstr_vtr(:,jp_sal,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 
    462              ENDDO 
    463          ENDIF 
     486      IF( ln_subbas ) THEN 
     487         ! 
     488         IF( cptr == 'adv' ) THEN 
     489             IF( ktra == jp_tem ) THEN  
     490                DO jn = 2, nptr 
     491                   htr_adv(:,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 
     492                END DO 
     493             ENDIF 
     494             IF( ktra == jp_sal ) THEN  
     495                DO jn = 2, nptr 
     496                   str_adv(:,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 
     497                END DO 
     498             ENDIF 
     499         ENDIF 
     500         IF( cptr == 'ldf' ) THEN 
     501             IF( ktra == jp_tem ) THEN  
     502                DO jn = 2, nptr 
     503                    htr_ldf(:,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 
     504                 END DO 
     505             ENDIF 
     506             IF( ktra == jp_sal ) THEN  
     507                DO jn = 2, nptr 
     508                   str_ldf(:,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 
     509                END DO 
     510             ENDIF 
     511         ENDIF 
     512         IF( cptr == 'eiv' ) THEN 
     513             IF( ktra == jp_tem ) THEN  
     514                DO jn = 2, nptr 
     515                    htr_eiv(:,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 
     516                 END DO 
     517             ENDIF 
     518             IF( ktra == jp_sal ) THEN  
     519                DO jn = 2, nptr 
     520                   str_eiv(:,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 
     521                END DO 
     522             ENDIF 
     523         ENDIF 
     524         ! 
    464525      ENDIF 
    465       ! 
    466526   END SUBROUTINE dia_ptr_hst 
    467527 
     
    476536      ierr(:) = 0 
    477537      ! 
    478       IF( .NOT. ALLOCATED( btmsk ) ) THEN 
    479          ALLOCATE( btmsk(jpi,jpj,nptr)    , btmsk34(jpi,jpj,nptr),   & 
    480             &      hstr_adv(jpj,jpts,nptr), hstr_eiv(jpj,jpts,nptr), & 
    481             &      hstr_ove(jpj,jpts,nptr), hstr_btr(jpj,jpts,nptr), & 
    482             &      hstr_ldf(jpj,jpts,nptr), hstr_vtr(jpj,jpts,nptr), STAT=ierr(1)  ) 
    483             ! 
    484          ALLOCATE( p_fval1d(jpj), p_fval2d(jpj,jpk), Stat=ierr(2)) 
    485          ! 
    486          dia_ptr_alloc = MAXVAL( ierr ) 
    487          CALL mpp_sum( 'diaptr', dia_ptr_alloc ) 
    488       ENDIF 
     538      ALLOCATE( btmsk(jpi,jpj,nptr) ,              & 
     539         &      htr_adv(jpj,nptr) , str_adv(jpj,nptr) ,   & 
     540         &      htr_eiv(jpj,nptr) , str_eiv(jpj,nptr) ,   & 
     541         &      htr_ove(jpj,nptr) , str_ove(jpj,nptr) ,   & 
     542         &      htr_btr(jpj,nptr) , str_btr(jpj,nptr) ,   & 
     543         &      htr_ldf(jpj,nptr) , str_ldf(jpj,nptr) , STAT=ierr(1)  ) 
     544         ! 
     545      ALLOCATE( p_fval1d(jpj), p_fval2d(jpj,jpk), Stat=ierr(2)) 
     546      ! 
     547      ALLOCATE( btm30(jpi,jpj), STAT=ierr(3)  ) 
     548 
     549         ! 
     550      dia_ptr_alloc = MAXVAL( ierr ) 
     551      CALL mpp_sum( 'diaptr', dia_ptr_alloc ) 
    489552      ! 
    490553   END FUNCTION dia_ptr_alloc 
     
    502565      !! ** Action  : - p_fval: i-k-mean poleward flux of pva 
    503566      !!---------------------------------------------------------------------- 
    504       REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk)  ::   pva   ! mask flux array at V-point 
    505       REAL(wp), INTENT(in), DIMENSION(jpi,jpj)      ::   pmsk   ! Optional 2D basin mask 
     567      REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk)       ::   pva   ! mask flux array at V-point 
     568      REAL(wp), INTENT(in), DIMENSION(jpi,jpj), OPTIONAL ::   pmsk   ! Optional 2D basin mask 
    506569      ! 
    507570      INTEGER                  ::   ji, jj, jk   ! dummy loop arguments 
     
    514577      ijpj = jpj 
    515578      p_fval(:) = 0._wp 
    516       DO jk = 1, jpkm1 
    517          DO jj = 2, jpjm1 
    518             DO ji = fs_2, fs_jpim1   ! Vector opt. 
    519                p_fval(jj) = p_fval(jj) + pva(ji,jj,jk) * pmsk(ji,jj) * tmask_i(ji,jj) 
     579      IF( PRESENT( pmsk ) ) THEN  
     580         DO jk = 1, jpkm1 
     581            DO jj = 2, jpjm1 
     582               DO ji = fs_2, fs_jpim1   ! Vector opt. 
     583                  p_fval(jj) = p_fval(jj) + pva(ji,jj,jk) * tmask_i(ji,jj) * pmsk(ji,jj) 
     584               END DO 
    520585            END DO 
    521586         END DO 
    522       END DO 
     587      ELSE 
     588         DO jk = 1, jpkm1 
     589            DO jj = 2, jpjm1 
     590               DO ji = fs_2, fs_jpim1   ! Vector opt. 
     591                  p_fval(jj) = p_fval(jj) + pva(ji,jj,jk) * tmask_i(ji,jj)  
     592               END DO 
     593            END DO 
     594         END DO 
     595      ENDIF 
    523596#if defined key_mpp_mpi 
    524597      CALL mpp_sum( 'diaptr', p_fval, ijpj, ncomm_znl) 
     
    539612      !! ** Action  : - p_fval: i-k-mean poleward flux of pva 
    540613      !!---------------------------------------------------------------------- 
    541       REAL(wp) , INTENT(in), DIMENSION(jpi,jpj) ::   pva   ! mask flux array at V-point 
    542       REAL(wp) , INTENT(in), DIMENSION(jpi,jpj) ::   pmsk   ! Optional 2D basin mask 
     614      REAL(wp) , INTENT(in), DIMENSION(jpi,jpj)           ::   pva   ! mask flux array at V-point 
     615      REAL(wp) , INTENT(in), DIMENSION(jpi,jpj), OPTIONAL ::   pmsk   ! Optional 2D basin mask 
    543616      ! 
    544617      INTEGER                  ::   ji,jj       ! dummy loop arguments 
     
    551624      ijpj = jpj 
    552625      p_fval(:) = 0._wp 
    553       DO jj = 2, jpjm1 
    554          DO ji = fs_2, fs_jpim1   ! Vector opt. 
    555             p_fval(jj) = p_fval(jj) + pva(ji,jj) * pmsk(ji,jj) * tmask_i(ji,jj) 
     626      IF( PRESENT( pmsk ) ) THEN  
     627         DO jj = 2, jpjm1 
     628            DO ji = nldi, nlei   ! No vector optimisation here. Better use a mask ? 
     629               p_fval(jj) = p_fval(jj) + pva(ji,jj) * tmask_i(ji,jj) * pmsk(ji,jj) 
     630            END DO 
    556631         END DO 
    557       END DO 
     632      ELSE 
     633         DO jj = 2, jpjm1 
     634            DO ji = nldi, nlei   ! No vector optimisation here. Better use a mask ? 
     635               p_fval(jj) = p_fval(jj) + pva(ji,jj) * tmask_i(ji,jj) 
     636            END DO 
     637         END DO 
     638      ENDIF 
    558639#if defined key_mpp_mpi 
    559640      CALL mpp_sum( 'diaptr', p_fval, ijpj, ncomm_znl ) 
     
    562643   END FUNCTION ptr_sj_2d 
    563644 
    564    FUNCTION ptr_ci_2d( pva )   RESULT ( p_fval ) 
    565       !!---------------------------------------------------------------------- 
    566       !!                    ***  ROUTINE ptr_ci_2d  *** 
    567       !! 
    568       !! ** Purpose :   "meridional" cumulated sum computation of a j-flux array 
    569       !! 
    570       !! ** Method  : - j cumulated sum of pva using the interior 2D vmask (umask_i). 
    571       !! 
    572       !! ** Action  : - p_fval: j-cumulated sum of pva 
    573       !!---------------------------------------------------------------------- 
    574       REAL(wp) , INTENT(in), DIMENSION(jpi,jpj)  ::   pva   ! mask flux array at V-point 
    575       ! 
    576       INTEGER                  ::   ji,jj,jc       ! dummy loop arguments 
    577       INTEGER                  ::   ijpj        ! ???  
    578       REAL(wp), DIMENSION(jpi,jpj) :: p_fval ! function value 
    579       !!-------------------------------------------------------------------- 
    580       !  
    581       ijpj = jpj  ! ??? 
    582       p_fval(:,:) = 0._wp 
    583       DO jc = 1, jpnj ! looping over all processors in j axis 
    584          DO jj = 2, jpjm1 
    585             DO ji = fs_2, fs_jpim1   ! Vector opt. 
    586                p_fval(ji,jj) = p_fval(ji,jj-1) + pva(ji,jj) * tmask_i(ji,jj) 
    587             END DO 
    588          END DO 
    589          CALL lbc_lnk( 'diaptr', p_fval, 'U', -1. ) 
    590       END DO 
    591       !  
    592    END FUNCTION ptr_ci_2d 
    593  
    594  
    595645 
    596646   FUNCTION ptr_sjk( pta, pmsk )   RESULT ( p_fval ) 
     
    606656      !! 
    607657      IMPLICIT none 
    608       REAL(wp) , INTENT(in), DIMENSION(jpi,jpj,jpk) ::   pta    ! mask flux array at V-point 
    609       REAL(wp) , INTENT(in), DIMENSION(jpi,jpj)    ::   pmsk   ! Optional 2D basin mask 
     658      REAL(wp) , INTENT(in), DIMENSION(jpi,jpj,jpk)           ::   pta    ! mask flux array at V-point 
     659      REAL(wp) , INTENT(in), DIMENSION(jpi,jpj)    , OPTIONAL ::   pmsk   ! Optional 2D basin mask 
    610660      !! 
    611661      INTEGER                           :: ji, jj, jk ! dummy loop arguments 
     
    623673      p_fval(:,:) = 0._wp 
    624674      ! 
    625       DO jk = 1, jpkm1 
    626          DO jj = 2, jpjm1 
    627             DO ji = fs_2, fs_jpim1   ! Vector opt. 
    628                p_fval(jj,jk) = p_fval(jj,jk) + pta(ji,jj,jk) * pmsk(ji,jj) * tmask_i(ji,jj) 
     675      IF( PRESENT( pmsk ) ) THEN  
     676         DO jk = 1, jpkm1 
     677            DO jj = 2, jpjm1 
     678!!gm here, use of tmask_i  ==> no need of loop over nldi, nlei.... 
     679               DO ji =  nldi, nlei   ! No vector optimisation here. Better use a mask ? 
     680                  p_fval(jj,jk) = p_fval(jj,jk) + pta(ji,jj,jk) * pmsk(ji,jj) 
     681               END DO 
    629682            END DO 
    630683         END DO 
    631       END DO 
     684      ELSE  
     685         DO jk = 1, jpkm1 
     686            DO jj = 2, jpjm1 
     687               DO ji =  nldi, nlei   ! No vector optimisation here. Better use a mask ? 
     688                  p_fval(jj,jk) = p_fval(jj,jk) + pta(ji,jj,jk) * tmask_i(ji,jj) 
     689               END DO 
     690            END DO 
     691         END DO 
     692      END IF 
    632693      ! 
    633694#if defined key_mpp_mpi 
  • NEMO/trunk/src/OCE/IOM/iom.F90

    r11989 r11993  
    5656   LOGICAL, PUBLIC, PARAMETER ::   lk_iomput = .FALSE.       !: iom_put flag 
    5757#endif 
    58    PUBLIC iom_init, iom_swap, iom_open, iom_close, iom_setkt, iom_varid, iom_get, iom_get_var 
     58   PUBLIC iom_init, iom_swap, iom_open, iom_close, iom_setkt, iom_varid, iom_get 
    5959   PUBLIC iom_chkatt, iom_getatt, iom_putatt, iom_getszuld, iom_rstput, iom_delay_rst, iom_put 
    6060   PUBLIC iom_use, iom_context_finalize, iom_miss_val 
     
    6262   PRIVATE iom_rp0d, iom_rp1d, iom_rp2d, iom_rp3d 
    6363   PRIVATE iom_g0d, iom_g1d, iom_g2d, iom_g3d, iom_get_123d 
    64    PRIVATE iom_p1d, iom_p2d, iom_p3d, iom_p4d 
     64   PRIVATE iom_p1d, iom_p2d, iom_p3d 
    6565#if defined key_iomput 
    6666   PRIVATE iom_set_domain_attr, iom_set_axis_attr, iom_set_field_attr, iom_set_file_attr, iom_get_file_attr, iom_set_grid_attr 
     
    8383   END INTERFACE 
    8484   INTERFACE iom_put 
    85       MODULE PROCEDURE iom_p0d, iom_p1d, iom_p2d, iom_p3d, iom_p4d 
     85      MODULE PROCEDURE iom_p0d, iom_p1d, iom_p2d, iom_p3d 
    8686   END INTERFACE iom_put 
    8787   
     
    108108      TYPE(xios_date)     :: start_date 
    109109      CHARACTER(len=lc) :: clname 
    110       INTEGER             :: irefyear, irefmonth, irefday 
    111110      INTEGER           :: ji, jkmin 
    112111      LOGICAL :: llrst_context              ! is context related to restart 
     
    140139 
    141140      ! Calendar type is now defined in xml file  
    142       IF (.NOT.(xios_getvar('ref_year' ,irefyear ))) irefyear  = 1900 
    143       IF (.NOT.(xios_getvar('ref_month',irefmonth))) irefmonth = 01 
    144       IF (.NOT.(xios_getvar('ref_day'  ,irefday  ))) irefday   = 01 
    145  
    146141      SELECT CASE ( nleapy )        ! Choose calendar for IOIPSL 
    147       CASE ( 1)   ; CALL xios_define_calendar( TYPE = "Gregorian", time_origin = xios_date(irefyear,irefmonth,irefday,00,00,00), & 
     142      CASE ( 1)   ; CALL xios_define_calendar( TYPE = "Gregorian", time_origin = xios_date(1900,01,01,00,00,00), & 
    148143          &                                    start_date = xios_date(nyear,nmonth,nday,0,0,0) ) 
    149       CASE ( 0)   ; CALL xios_define_calendar( TYPE = "NoLeap"   , time_origin = xios_date(irefyear,irefmonth,irefday,00,00,00), & 
     144      CASE ( 0)   ; CALL xios_define_calendar( TYPE = "NoLeap"   , time_origin = xios_date(1900,01,01,00,00,00), & 
    150145          &                                    start_date = xios_date(nyear,nmonth,nday,0,0,0) ) 
    151       CASE (30)   ; CALL xios_define_calendar( TYPE = "D360"     , time_origin = xios_date(irefyear,irefmonth,irefday,00,00,00), & 
     146      CASE (30)   ; CALL xios_define_calendar( TYPE = "D360"     , time_origin = xios_date(1900,01,01,00,00,00), & 
    152147          &                                    start_date = xios_date(nyear,nmonth,nday,0,0,0) ) 
    153148      END SELECT 
     
    228223          CALL iom_set_axis_attr( "icbcla", class_num ) 
    229224          CALL iom_set_axis_attr( "iax_20C", (/ REAL(20,wp) /) )   ! strange syntaxe and idea... 
    230           CALL iom_set_axis_attr( "iax_26C", (/ REAL(26,wp) /) )   ! strange syntaxe and idea... 
    231225          CALL iom_set_axis_attr( "iax_28C", (/ REAL(28,wp) /) )   ! strange syntaxe and idea... 
    232           CALL iom_set_axis_attr( "basin"  , (/ (REAL(ji,wp), ji=1,5) /) ) 
    233226      ENDIF 
    234227      ! 
     
    13471340   END SUBROUTINE iom_get_123d 
    13481341 
    1349    SUBROUTINE iom_get_var( cdname, z2d) 
    1350       CHARACTER(LEN=*), INTENT(in ) ::   cdname 
    1351       REAL(wp), DIMENSION(jpi,jpj) ::   z2d  
    1352 #if defined key_iomput 
    1353       IF( xios_field_is_active( cdname, at_current_timestep_arg = .TRUE. ) ) THEN 
    1354          z2d(:,:) = 0._wp 
    1355          CALL xios_recv_field( cdname, z2d) 
    1356       ENDIF 
    1357 #else 
    1358       IF( .FALSE. )   WRITE(numout,*) cdname, z2d ! useless test to avoid compilation warnings 
    1359 #endif 
    1360    END SUBROUTINE iom_get_var 
    1361  
    13621342 
    13631343   FUNCTION iom_getszuld ( kiomid )   
     
    17291709   END SUBROUTINE iom_p3d 
    17301710 
    1731    SUBROUTINE iom_p4d( cdname, pfield4d ) 
    1732       CHARACTER(LEN=*)                , INTENT(in) ::   cdname 
    1733       REAL(wp),       DIMENSION(:,:,:,:), INTENT(in) ::   pfield4d 
    1734 #if defined key_iomput 
    1735       CALL xios_send_field(cdname, pfield4d) 
    1736 #else 
    1737       IF( .FALSE. )   WRITE(numout,*) cdname, pfield4d   ! useless test to avoid compilation warnings 
    1738 #endif 
    1739    END SUBROUTINE iom_p4d 
    1740  
    1741  
    17421711#if defined key_iomput 
    17431712   !!---------------------------------------------------------------------- 
     
    20822051      ALLOCATE( zlon(ni*nj) )       ;       zlon(:) = 0._wp 
    20832052      ! 
    2084 !      CALL dom_ngb( -168.53, 65.03, ix, iy, 'T' ) !  i-line that passes through Bering Strait: Reference latitude (used in plots) 
    2085       CALL dom_ngb( 180., 90., ix, iy, 'T' ) !  i-line that passes near the North Pole : Reference latitude (used in plots) 
     2053      CALL dom_ngb( -168.53, 65.03, ix, iy, 'T' ) !  i-line that passes through Bering Strait: Reference latitude (used in plots) 
     2054!      CALL dom_ngb( 180., 90., ix, iy, 'T' ) !  i-line that passes near the North Pole : Reference latitude (used in plots) 
    20862055      CALL iom_set_domain_attr("gznl", ni_glo=jpiglo, nj_glo=jpjglo, ibegin=nimpp+nldi-2, jbegin=njmpp+nldj-2, ni=ni, nj=nj) 
    20872056      CALL iom_set_domain_attr("gznl", data_dim=2, data_ibegin = 1-nldi, data_ni = jpi, data_jbegin = 1-nldj, data_nj = jpj) 
    20882057      CALL iom_set_domain_attr("gznl", lonvalue = zlon,   & 
    20892058         &                             latvalue = RESHAPE(plat(nldi:nlei, nldj:nlej),(/ ni*nj /)))   
    2090       CALL iom_set_zoom_domain_attr("ptr", ibegin=ix-1, jbegin=0, ni=1, nj=jpjglo) 
     2059      CALL iom_set_zoom_domain_attr("znl_T", ibegin=ix-1, jbegin=0, ni=1, nj=jpjglo) 
     2060      CALL iom_set_zoom_domain_attr("znl_W", ibegin=ix-1, jbegin=0, ni=1, nj=jpjglo) 
    20912061      ! 
    20922062      CALL iom_update_file_name('ptr') 
  • NEMO/trunk/src/OCE/LDF/ldftra.F90

    r11989 r11993  
    851851      CALL iom_put( "woce_eiv", zw3d ) 
    852852      ! 
    853       IF( iom_use('weiv_masstr') ) THEN   ! vertical mass transport & its square value 
    854          zw2d(:,:) = rau0 * e1e2t(:,:) 
    855          DO jk = 1, jpk 
    856             zw3d(:,:,jk) = zw3d(:,:,jk) * zw2d(:,:) 
    857          END DO 
    858          CALL iom_put( "weiv_masstr" , zw3d )   
    859       ENDIF 
    860       ! 
    861       IF( iom_use('ueiv_masstr') ) THEN 
    862          zw3d(:,:,:) = 0.e0 
    863          DO jk = 1, jpkm1 
    864             zw3d(:,:,jk) = rau0 * ( psi_uw(:,:,jk+1) - psi_uw(:,:,jk) )  
    865          END DO 
    866          CALL iom_put( "ueiv_masstr", zw3d )                  ! mass transport in i-direction 
    867       ENDIF 
    868853      ! 
    869854      zztmp = 0.5_wp * rau0 * rcp  
     
    885870        CALL iom_put( "ueiv_heattr3d", zztmp * zw3d )                  ! heat transport in i-direction 
    886871      ENDIF 
    887       ! 
    888       IF( iom_use('veiv_masstr') ) THEN 
    889          zw3d(:,:,:) = 0.e0 
    890          DO jk = 1, jpkm1 
    891             zw3d(:,:,jk) = rau0 * ( psi_vw(:,:,jk+1) - psi_vw(:,:,jk) )  
    892          END DO 
    893          CALL iom_put( "veiv_masstr", zw3d )                  ! mass transport in i-direction 
    894       ENDIF 
    895       ! 
    896872      zw2d(:,:)   = 0._wp  
    897873      zw3d(:,:,:) = 0._wp  
     
    909885      CALL iom_put( "veiv_heattr", zztmp * zw3d )                  !  heat transport in j-direction 
    910886      ! 
    911       IF( iom_use( 'sophteiv' ) )   CALL dia_ptr_hst( jp_tem, 'eiv', 0.5 * zw3d ) 
     887      IF( ln_diaptr )  CALL dia_ptr_hst( jp_tem, 'eiv', 0.5 * zw3d ) 
    912888      ! 
    913889      zztmp = 0.5_wp * 0.5 
     
    944920      CALL iom_put( "veiv_salttr", zztmp * zw3d )                  !  salt transport in j-direction 
    945921      ! 
    946       IF( iom_use( 'sopsteiv' ) ) CALL dia_ptr_hst( jp_sal, 'eiv', 0.5 * zw3d ) 
     922      IF( ln_diaptr ) CALL dia_ptr_hst( jp_sal, 'eiv', 0.5 * zw3d ) 
    947923      ! 
    948924      ! 
  • NEMO/trunk/src/OCE/SBC/sbcblk.F90

    r11989 r11993  
    801801      REAL(wp), DIMENSION(jpi,jpj)     ::   zevap, zsnw   ! evaporation and snw distribution after wind blowing (SI3) 
    802802      REAL(wp), DIMENSION(jpi,jpj)     ::   zrhoa 
    803       REAL(wp), DIMENSION(jpi,jpj)     ::   ztmp, ztmp2 
    804803      !!--------------------------------------------------------------------- 
    805804      ! 
     
    914913         qtr_ice_top(:,:,:) = 0._wp  
    915914      END WHERE 
    916       ! 
    917  
    918       IF( iom_use('evap_ao_cea') .OR. iom_use('hflx_evap_cea') ) THEN 
    919          ztmp(:,:) = zevap(:,:) * ( 1._wp - at_i_b(:,:) )  
    920          IF( iom_use('evap_ao_cea'  ) )  CALL iom_put( 'evap_ao_cea'  , ztmp(:,:) * tmask(:,:,1) )   ! ice-free oce evap (cell average) 
    921          IF( iom_use('hflx_evap_cea') )  CALL iom_put( 'hflx_evap_cea', ztmp(:,:) * sst_m(:,:) * rcp * tmask(:,:,1) )   ! heat flux from evap (cell average) 
    922       ENDIF 
    923       IF( iom_use('hflx_rain_cea') ) THEN 
    924          ztmp(:,:) = rcp * ( SUM( (ptsu-rt0) * a_i_b, dim=3 ) + sst_m(:,:) * ( 1._wp - at_i_b(:,:) ) ) 
    925          IF( iom_use('hflx_rain_cea') )  CALL iom_put( 'hflx_rain_cea', ( tprecip(:,:) - sprecip(:,:) ) * ztmp(:,:) )   ! heat flux from rain (cell average) 
    926       ENDIF 
    927       IF( iom_use('hflx_snow_cea') .OR. iom_use('hflx_snow_ao_cea') .OR. iom_use('hflx_snow_ai_cea')  )  THEN 
    928           WHERE( SUM( a_i_b, dim=3 ) > 1.e-10 ) ;   ztmp(:,:) = rcpi * SUM( (ptsu-rt0) * a_i_b, dim=3 ) / SUM( a_i_b, dim=3 ) 
    929           ELSEWHERE                             ;   ztmp(:,:) = rcp * sst_m(:,:)     
    930           ENDWHERE 
    931           ztmp2(:,:) = sprecip(:,:) * ( ztmp(:,:) - rLfus )  
    932           IF( iom_use('hflx_snow_cea')    ) CALL iom_put('hflx_snow_cea'   , ztmp2(:,:) ) ! heat flux from snow (cell average) 
    933           IF( iom_use('hflx_snow_ao_cea') ) CALL iom_put('hflx_snow_ao_cea', ztmp2(:,:) * ( 1._wp - zsnw(:,:) ) ) ! heat flux from snow (over ocean) 
    934           IF( iom_use('hflx_snow_ai_cea') ) CALL iom_put('hflx_snow_ai_cea', ztmp2(:,:) *           zsnw(:,:)   ) ! heat flux from snow (over ice) 
    935       ENDIF 
    936915      ! 
    937916      IF(ln_ctl) THEN 
  • NEMO/trunk/src/OCE/SBC/sbccpl.F90

    r11989 r11993  
    17741774      IF( iom_use('snow_ao_cea') )  CALL iom_put( 'snow_ao_cea' , sprecip(:,:) * ( 1._wp - zsnw(:,:) )                  )  ! Snow over ice-free ocean  (cell average) 
    17751775      IF( iom_use('snow_ai_cea') )  CALL iom_put( 'snow_ai_cea' , sprecip(:,:) *           zsnw(:,:)                    )  ! Snow over sea-ice         (cell average) 
    1776       IF( iom_use('rain_ao_cea') )  CALL iom_put( 'rain_ao_cea' , ( tprecip(:,:) - sprecip(:,:) ) * picefr(:,:)         )  ! liquid precipitation over ocean (cell average) 
    17771776      IF( iom_use('subl_ai_cea') )  CALL iom_put( 'subl_ai_cea' , frcv(jpr_ievp)%z3(:,:,1) * picefr(:,:) * tmask(:,:,1) )  ! Sublimation over sea-ice (cell average) 
    17781777      IF( iom_use('evap_ao_cea') )  CALL iom_put( 'evap_ao_cea' , ( frcv(jpr_tevp)%z3(:,:,1)  & 
     
    18981897#endif 
    18991898      ! outputs 
    1900       IF( srcv(jpr_cal)%laction       ) CALL iom_put('hflx_cal_cea'    , - frcv(jpr_cal)%z3(:,:,1) * rLfus )                      ! latent heat from calving 
    1901       IF( srcv(jpr_icb)%laction       ) CALL iom_put('hflx_icb_cea'    , - frcv(jpr_icb)%z3(:,:,1) * rLfus )                      ! latent heat from icebergs melting 
    1902       IF( iom_use('hflx_rain_cea')    ) CALL iom_put('hflx_rain_cea'   , ( tprecip(:,:) - sprecip(:,:) ) * zcptrain(:,:) )        ! heat flux from rain (cell average) 
    1903       IF( iom_use('hflx_evap_cea')    ) CALL iom_put('hflx_evap_cea'   , ( frcv(jpr_tevp)%z3(:,:,1) - frcv(jpr_ievp)%z3(:,:,1) & 
     1899      IF ( srcv(jpr_cal)%laction       ) CALL iom_put('hflx_cal_cea'    , - frcv(jpr_cal)%z3(:,:,1) * rLfus )                      ! latent heat from calving 
     1900      IF ( srcv(jpr_icb)%laction       ) CALL iom_put('hflx_icb_cea'    , - frcv(jpr_icb)%z3(:,:,1) * rLfus )                      ! latent heat from icebergs melting 
     1901      IF ( iom_use('hflx_rain_cea')    ) CALL iom_put('hflx_rain_cea'   , ( tprecip(:,:) - sprecip(:,:) ) * zcptrain(:,:) )        ! heat flux from rain (cell average) 
     1902      IF ( iom_use('hflx_evap_cea')    ) CALL iom_put('hflx_evap_cea'   , ( frcv(jpr_tevp)%z3(:,:,1) - frcv(jpr_ievp)%z3(:,:,1) & 
    19041903           &                                                              * picefr(:,:) ) * zcptn(:,:) * tmask(:,:,1) )            ! heat flux from evap (cell average) 
    1905       IF( iom_use('hflx_prec_cea')    ) CALL iom_put('hflx_prec_cea'   ,  sprecip(:,:) * ( zcptsnw(:,:) - rLfus ) +  &                    ! heat flux from all precip (cell avg) 
    1906          &                                                               ( tprecip(:,:) - sprecip(:,:) ) * zcptrain(:,:) ) 
    1907       IF( iom_use('hflx_snow_cea')    ) CALL iom_put('hflx_snow_cea'   , sprecip(:,:) * ( zcptsnw(:,:) - rLfus )  )               ! heat flux from snow (cell average) 
    1908       IF( iom_use('hflx_snow_ao_cea') ) CALL iom_put('hflx_snow_ao_cea', sprecip(:,:) * ( zcptsnw(:,:) - rLfus ) & 
     1904      IF ( iom_use('hflx_snow_cea')    ) CALL iom_put('hflx_snow_cea'   , sprecip(:,:) * ( zcptsnw(:,:) - rLfus )  )               ! heat flux from snow (cell average) 
     1905      IF ( iom_use('hflx_snow_ao_cea') ) CALL iom_put('hflx_snow_ao_cea', sprecip(:,:) * ( zcptsnw(:,:) - rLfus ) & 
    19091906           &                                                              * ( 1._wp - zsnw(:,:) )                  )               ! heat flux from snow (over ocean) 
    1910       IF( iom_use('hflx_snow_ai_cea') ) CALL iom_put('hflx_snow_ai_cea', sprecip(:,:) * ( zcptsnw(:,:) - rLfus ) &  
     1907      IF ( iom_use('hflx_snow_ai_cea') ) CALL iom_put('hflx_snow_ai_cea', sprecip(:,:) * ( zcptsnw(:,:) - rLfus ) &  
    19111908           &                                                              *           zsnw(:,:)                    )               ! heat flux from snow (over ice) 
    19121909      ! note: hflx for runoff and iceshelf are done in sbcrnf and sbcisf resp. 
     
    23042301      !                                                      !  CO2 flux from PISCES     !  
    23052302      !                                                      ! ------------------------- ! 
    2306       IF( ssnd(jps_co2)%laction .AND. l_co2cpl )   THEN  
    2307          ztmp1(:,:) = oce_co2(:,:) * 1000.  ! conversion in molC/m2/s 
    2308          CALL cpl_snd( jps_co2, isec, RESHAPE ( ztmp1, (/jpi,jpj,1/) ) , info ) 
    2309       ENDIF 
     2303      IF( ssnd(jps_co2)%laction .AND. l_co2cpl )   CALL cpl_snd( jps_co2, isec, RESHAPE ( oce_co2, (/jpi,jpj,1/) ) , info ) 
    23102304      ! 
    23112305      !                                                      ! ------------------------- ! 
  • NEMO/trunk/src/OCE/SBC/sbcmod.F90

    r11989 r11993  
    244244         fwfisf_b(:,:)   = 0._wp   ;   risf_tsc_b(:,:,:) = 0._wp 
    245245      END IF 
    246       ! 
    247       IF( sbc_ssr_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'sbc_init : unable to allocate sbc_ssr arrays' ) 
    248       IF( .NOT.ln_ssr ) THEN               !* Initialize qrp and erp if no restoring  
    249          qrp(:,:) = 0._wp 
    250          erp(:,:) = 0._wp 
    251       ENDIF 
    252       ! 
    253  
    254246      IF( nn_ice == 0 ) THEN        !* No sea-ice in the domain : ice fraction is always zero 
    255247         IF( nn_components /= jp_iam_opa )   fr_i(:,:) = 0._wp    ! except for OPA in SAS-OPA coupled case 
     
    560552         CALL iom_put( "taum"  , taum       )                   ! wind stress module 
    561553         CALL iom_put( "wspd"  , wndm       )                   ! wind speed  module over free ocean or leads in presence of sea-ice 
    562          CALL iom_put( "qrp", qrp )                             ! heat flux damping 
    563          CALL iom_put( "erp", erp )                             ! freshwater flux damping 
    564554      ENDIF 
    565555      ! 
  • NEMO/trunk/src/OCE/SBC/sbcrnf.F90

    r11989 r11993  
    4343   REAL(wp)                   ::      rn_dep_max        !: depth over which runoffs is spread       (ln_rnf_depth_ini =T) 
    4444   INTEGER                    ::      nn_rnf_depth_file !: create (=1) a runoff depth file or not (=0) 
    45    LOGICAL                    ::   ln_rnf_icb        !: iceberg flux is specified in a file 
    4645   LOGICAL                    ::   ln_rnf_tem        !: temperature river runoffs attribute specified in a file 
    4746   LOGICAL           , PUBLIC ::   ln_rnf_sal        !: salinity    river runoffs attribute specified in a file 
    4847   TYPE(FLD_N)       , PUBLIC ::   sn_rnf            !: information about the runoff file to be read 
    4948   TYPE(FLD_N)                ::   sn_cnf            !: information about the runoff mouth file to be read 
    50    TYPE(FLD_N)                ::   sn_i_rnf        !: information about the iceberg flux file to be read 
    5149   TYPE(FLD_N)                ::   sn_s_rnf          !: information about the salinities of runoff file to be read 
    5250   TYPE(FLD_N)                ::   sn_t_rnf          !: information about the temperatures of runoff file to be read 
     
    6765 
    6866   TYPE(FLD),        ALLOCATABLE, DIMENSION(:) ::   sf_rnf       ! structure: river runoff (file information, fields read) 
    69    TYPE(FLD),        ALLOCATABLE, DIMENSION(:) ::   sf_i_rnf     ! structure: iceberg flux (file information, fields read) 
    7067   TYPE(FLD),        ALLOCATABLE, DIMENSION(:) ::   sf_s_rnf     ! structure: river runoff salinity (file information, fields read)   
    7168   TYPE(FLD),        ALLOCATABLE, DIMENSION(:) ::   sf_t_rnf     ! structure: river runoff temperature (file information, fields read)   
     
    115112      !                                            !-------------------! 
    116113      ! 
    117       ! 
    118       IF( .NOT. l_rnfcpl )  THEN 
    119                             CALL fld_read ( kt, nn_fsbc, sf_rnf   )    ! Read Runoffs data and provide it at kt ( runoffs + iceberg ) 
    120          IF( ln_rnf_icb )   CALL fld_read ( kt, nn_fsbc, sf_i_rnf )    ! idem for iceberg flux if required 
    121       ENDIF 
     114      IF( .NOT. l_rnfcpl )   CALL fld_read ( kt, nn_fsbc, sf_rnf   )    ! Read Runoffs data and provide it at kt 
    122115      IF(   ln_rnf_tem   )   CALL fld_read ( kt, nn_fsbc, sf_t_rnf )    ! idem for runoffs temperature if required 
    123116      IF(   ln_rnf_sal   )   CALL fld_read ( kt, nn_fsbc, sf_s_rnf )    ! idem for runoffs salinity    if required 
     
    125118      IF( MOD( kt - 1, nn_fsbc ) == 0 ) THEN 
    126119         ! 
    127          IF( .NOT. l_rnfcpl ) THEN 
    128              rnf(:,:) = rn_rfact * ( sf_rnf(1)%fnow(:,:,1) ) * tmask(:,:,1)  ! updated runoff value at time step kt 
    129              IF( ln_rnf_icb ) THEN 
    130                 fwficb(:,:) = rn_rfact * ( sf_i_rnf(1)%fnow(:,:,1) ) * tmask(:,:,1)  ! updated runoff value at time step kt 
    131                 CALL iom_put( 'iceberg_cea'  , fwficb(:,:)  )         ! output iceberg flux 
    132                 CALL iom_put( 'hflx_icb_cea' , fwficb(:,:) * rLfus )   ! output Heat Flux into Sea Water due to Iceberg Thermodynamics --> 
    133              ENDIF 
    134          ENDIF 
     120         IF( .NOT. l_rnfcpl )   rnf(:,:) = rn_rfact * ( sf_rnf(1)%fnow(:,:,1) ) * tmask(:,:,1)       ! updated runoff value at time step kt 
    135121         ! 
    136122         !                                                           ! set temperature & salinity content of runoffs 
     
    146132         ELSE                                                        ! use SST as runoffs temperature 
    147133            !CEOD River is fresh water so must at least be 0 unless we consider ice 
    148             rnf_tsc(:,:,jp_tem) = MAX( sst_m(:,:), 0.0_wp ) * rnf(:,:) * r1_rau0 
     134            rnf_tsc(:,:,jp_tem) = MAX(sst_m(:,:),0.0_wp) * rnf(:,:) * r1_rau0 
    149135         ENDIF 
    150136         !                                                           ! use runoffs salinity data 
    151137         IF( ln_rnf_sal )   rnf_tsc(:,:,jp_sal) = ( sf_s_rnf(1)%fnow(:,:,1) ) * rnf(:,:) * r1_rau0 
    152138         !                                                           ! else use S=0 for runoffs (done one for all in the init) 
    153                                          CALL iom_put( 'runoffs'     , rnf(:,:)                         )   ! output runoff mass flux 
     139         IF( iom_use('runoffs') )        CALL iom_put( 'runoffs'     , rnf(:,:)                         )   ! output runoff mass flux 
    154140         IF( iom_use('hflx_rnf_cea') )   CALL iom_put( 'hflx_rnf_cea', rnf_tsc(:,:,jp_tem) * rau0 * rcp )   ! output runoff sensible heat (W/m2) 
    155141      ENDIF 
     
    256242      REAL(wp), DIMENSION(jpi,jpj,2) :: zrnfcl     
    257243      !! 
    258       NAMELIST/namsbc_rnf/ cn_dir            , ln_rnf_depth, ln_rnf_tem, ln_rnf_sal, ln_rnf_icb,   & 
    259          &                 sn_rnf, sn_cnf    , sn_i_rnf, sn_s_rnf    , sn_t_rnf  , sn_dep_rnf,   & 
     244      NAMELIST/namsbc_rnf/ cn_dir            , ln_rnf_depth, ln_rnf_tem, ln_rnf_sal,   & 
     245         &                 sn_rnf, sn_cnf    , sn_s_rnf    , sn_t_rnf  , sn_dep_rnf,   & 
    260246         &                 ln_rnf_mouth      , rn_hrnf     , rn_avt_rnf, rn_rfact,     & 
    261247         &                 ln_rnf_depth_ini  , rn_dep_max  , rn_rnf_max, nn_rnf_depth_file 
     
    313299         IF( sn_rnf%ln_tint ) ALLOCATE( sf_rnf(1)%fdta(jpi,jpj,1,2) ) 
    314300         CALL fld_fill( sf_rnf, (/ sn_rnf /), cn_dir, 'sbc_rnf_init', 'read runoffs data', 'namsbc_rnf', no_print ) 
    315          ! 
    316          IF( ln_rnf_icb ) THEN                      ! Create (if required) sf_i_rnf structure 
    317             IF(lwp) WRITE(numout,*) 
    318             IF(lwp) WRITE(numout,*) '          iceberg flux read in a file' 
    319             ALLOCATE( sf_i_rnf(1), STAT=ierror  ) 
    320             IF( ierror > 0 ) THEN 
    321                CALL ctl_stop( 'sbc_rnf_init: unable to allocate sf_i_rnf structure' )   ;   RETURN 
    322             ENDIF 
    323             ALLOCATE( sf_i_rnf(1)%fnow(jpi,jpj,1)   ) 
    324             IF( sn_i_rnf%ln_tint ) ALLOCATE( sf_i_rnf(1)%fdta(jpi,jpj,1,2) ) 
    325             CALL fld_fill (sf_i_rnf, (/ sn_i_rnf /), cn_dir, 'sbc_rnf_init', 'read iceberg flux data', 'namsbc_rnf' ) 
    326          ELSE 
    327             fwficb(:,:) = 0._wp 
    328          ENDIF 
    329  
    330301      ENDIF 
    331302      ! 
  • NEMO/trunk/src/OCE/SBC/sbcssr.F90

    r11989 r11993  
    3030   PUBLIC   sbc_ssr        ! routine called in sbcmod 
    3131   PUBLIC   sbc_ssr_init   ! routine called in sbcmod 
    32    PUBLIC   sbc_ssr_alloc  ! routine called in sbcmod 
    3332 
    3433   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   erp   !: evaporation damping   [kg/m2/s] 
    3534   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   qrp   !: heat flux damping        [w/m2] 
    36    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   coefice   !: under ice relaxation coefficient 
    3735 
    3836   !                                   !!* Namelist namsbc_ssr * 
     
    4341   LOGICAL         ::   ln_sssr_bnd     ! flag to bound erp term  
    4442   REAL(wp)        ::   rn_sssr_bnd     ! ABS(Max./Min.) value of erp term [mm/day] 
    45    INTEGER         ::   nn_icedmp       ! Control of restoring under ice 
    4643 
    4744   REAL(wp) , ALLOCATABLE, DIMENSION(:) ::   buffer   ! Temporary buffer for exchange 
     
    10097                  END DO 
    10198               END DO 
    102             ENDIF 
    103             ! 
    104             IF( nn_sssr /= 0 .AND. nn_icedmp /= 1 ) THEN 
    105               ! use fraction of ice ( fr_i ) to adjust relaxation under ice if nn_icedmp .ne. 1 
    106               ! n.b. coefice is initialised and fixed to 1._wp if nn_icedmp = 1 
    107                DO jj = 1, jpj 
    108                   DO ji = 1, jpi 
    109                      SELECT CASE ( nn_icedmp ) 
    110                        CASE ( 0 )    ;  coefice(ji,jj) = 1._wp - fr_i(ji,jj)              ! no/reduced damping under ice 
    111                        CASE  DEFAULT ;  coefice(ji,jj) = 1._wp +(nn_icedmp-1)*fr_i(ji,jj) ! reinforced damping (x nn_icedmp) under ice ) 
    112                      END SELECT 
    113                   END DO 
    114                END DO 
     99               CALL iom_put( "qrp", qrp )                             ! heat flux damping 
    115100            ENDIF 
    116101            ! 
     
    120105                  DO ji = 1, jpi 
    121106                     zerp = zsrp * ( 1. - 2.*rnfmsk(ji,jj) )   &      ! No damping in vicinity of river mouths 
    122                         &        *   coefice(ji,jj)            &      ! Optional control of damping under sea-ice 
    123107                        &        * ( sss_m(ji,jj) - sf_sss(1)%fnow(ji,jj,1) ) * tmask(ji,jj,1) 
    124108                     sfx(ji,jj) = sfx(ji,jj) + zerp                 ! salt flux 
     
    126110                  END DO 
    127111               END DO 
     112               CALL iom_put( "erp", erp )                             ! freshwater flux damping 
    128113               ! 
    129114            ELSEIF( nn_sssr == 2 ) THEN                               !* Salinity damping term (volume flux (emp) and associated heat flux (qns) 
     
    133118                  DO ji = 1, jpi                             
    134119                     zerp = zsrp * ( 1. - 2.*rnfmsk(ji,jj) )   &      ! No damping in vicinity of river mouths 
    135                         &        *   coefice(ji,jj)            &      ! Optional control of damping under sea-ice 
    136120                        &        * ( sss_m(ji,jj) - sf_sss(1)%fnow(ji,jj,1) )   & 
    137121                        &        / MAX(  sss_m(ji,jj), 1.e-20   ) * tmask(ji,jj,1) 
     
    142126                  END DO 
    143127               END DO 
     128               CALL iom_put( "erp", erp )                             ! freshwater flux damping 
    144129            ENDIF 
    145130            ! 
     
    169154      CHARACTER(len=100) ::  cn_dir          ! Root directory for location of ssr files 
    170155      TYPE(FLD_N) ::   sn_sst, sn_sss        ! informations about the fields to be read 
    171       NAMELIST/namsbc_ssr/ cn_dir, nn_sstr, nn_sssr, rn_dqdt, rn_deds, sn_sst, & 
    172               & sn_sss, ln_sssr_bnd, rn_sssr_bnd, nn_icedmp 
     156      NAMELIST/namsbc_ssr/ cn_dir, nn_sstr, nn_sssr, rn_dqdt, rn_deds, sn_sst, sn_sss, ln_sssr_bnd, rn_sssr_bnd 
    173157      INTEGER     ::  ios 
    174158      !!---------------------------------------------------------------------- 
     
    198182         WRITE(numout,*) '         flag to bound erp term                 ln_sssr_bnd = ', ln_sssr_bnd 
    199183         WRITE(numout,*) '         ABS(Max./Min.) erp threshold           rn_sssr_bnd = ', rn_sssr_bnd, ' mm/day' 
    200          WRITE(numout,*) '      Cntrl of surface restoration under ice nn_icedmp      = ', nn_icedmp 
    201          WRITE(numout,*) '          ( 0 = no restoration under ice)' 
    202          WRITE(numout,*) '          ( 1 = restoration everywhere  )' 
    203          WRITE(numout,*) '          (>1 = enhanced restoration under ice  )' 
    204       ENDIF 
     184      ENDIF 
     185      ! 
     186      !                            !* Allocate erp and qrp array 
     187      ALLOCATE( qrp(jpi,jpj), erp(jpi,jpj), STAT=ierror ) 
     188      IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_ssr: unable to allocate erp and qrp array' ) 
    205189      ! 
    206190      IF( nn_sstr == 1 ) THEN      !* set sf_sst structure & allocate arrays 
     
    232216      ENDIF 
    233217      ! 
    234       coefice(:,:) = 1._wp         !  Initialise coefice to 1._wp ; will not need to be changed if nn_icedmp=1 
    235218      !                            !* Initialize qrp and erp if no restoring  
    236219      IF( nn_sstr /= 1                   )   qrp(:,:) = 0._wp 
     
    238221      ! 
    239222   END SUBROUTINE sbc_ssr_init 
    240           
    241    INTEGER FUNCTION sbc_ssr_alloc() 
    242       !!---------------------------------------------------------------------- 
    243       !!               ***  FUNCTION sbc_ssr_alloc  *** 
    244       !!---------------------------------------------------------------------- 
    245       sbc_ssr_alloc = 0       ! set to zero if no array to be allocated 
    246       IF( .NOT. ALLOCATED( erp ) ) THEN 
    247          ALLOCATE( qrp(jpi,jpj), erp(jpi,jpj), coefice(jpi,jpj), STAT= sbc_ssr_alloc ) 
    248          ! 
    249          IF( lk_mpp                  )   CALL mpp_sum ( 'sbcssr', sbc_ssr_alloc ) 
    250          IF( sbc_ssr_alloc /= 0 )   CALL ctl_warn('sbc_ssr_alloc: failed to allocate arrays.') 
    251          ! 
    252       ENDIF 
    253    END FUNCTION 
    254223       
    255224   !!====================================================================== 
  • NEMO/trunk/src/OCE/TRA/eosbn2.F90

    r11989 r11993  
    3030   !!   eos_insitu_2d : Compute the in situ density for 2d fields 
    3131   !!   bn2           : Compute the Brunt-Vaisala frequency 
     32   !!   bn2           : compute the Brunt-Vaisala frequency 
    3233   !!   eos_pt_from_ct: compute the potential temperature from the Conservative Temperature 
    3334   !!   eos_rab       : generic interface of in situ thermal/haline expansion ratio  
     
    6667   END INTERFACE 
    6768   ! 
    68    INTERFACE eos_pt_from_ct 
    69       MODULE PROCEDURE eos_pt_from_ct_2d, eos_pt_from_ct_3d 
    70    END INTERFACE 
    71  
    7269   PUBLIC   eos            ! called by step, istate, tranpc and zpsgrd modules 
    7370   PUBLIC   bn2            ! called by step module 
     
    7976 
    8077   !                               !!** Namelist nameos ** 
    81    LOGICAL , PUBLIC ::   ln_TEOS10   
    82    LOGICAL , PUBLIC ::   ln_EOS80   
    83    LOGICAL , PUBLIC ::   ln_SEOS   
     78   LOGICAL , PUBLIC ::   ln_TEOS10 
     79   LOGICAL , PUBLIC ::   ln_EOS80 
     80   LOGICAL , PUBLIC ::   ln_SEOS 
    8481 
    8582   ! Parameters 
     
    939936 
    940937 
    941    FUNCTION eos_pt_from_ct_2d( ctmp, psal ) RESULT( ptmp ) 
     938   FUNCTION eos_pt_from_ct( ctmp, psal ) RESULT( ptmp ) 
    942939      !!---------------------------------------------------------------------- 
    943940      !!                 ***  ROUTINE eos_pt_from_ct  *** 
     
    962959      !!---------------------------------------------------------------------- 
    963960      ! 
    964       IF( ln_timing )   CALL timing_start('eos_pt_from_ct_2d') 
     961      IF( ln_timing )   CALL timing_start('eos_pt_from_ct') 
    965962      ! 
    966963      zdeltaS = 5._wp 
     
    993990      END DO 
    994991      ! 
    995       IF( ln_timing )   CALL timing_stop('eos_pt_from_ct_2d') 
    996       ! 
    997    END FUNCTION eos_pt_from_ct_2d 
    998  
    999     
    1000    FUNCTION eos_pt_from_ct_3d( ctmp, psal ) RESULT( ptmp ) 
    1001       !!---------------------------------------------------------------------- 
    1002       !!                 ***  ROUTINE eos_pt_from_ct  *** 
    1003       !! 
    1004       !! ** Purpose :   Compute pot.temp. from cons. temp. [Celcius] 
    1005       !! 
    1006       !! ** Method  :   rational approximation (5/3th order) of TEOS-10 algorithm 
    1007       !!       checkvalue: pt=20.02391895 Celsius for sa=35.7g/kg, ct=20degC 
    1008       !! 
    1009       !! Reference  :   TEOS-10, UNESCO 
    1010       !!                Rational approximation to TEOS10 algorithm (rms error on WOA13 values: 4.0e-5 degC) 
    1011       !!---------------------------------------------------------------------- 
    1012       REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in) ::   ctmp   ! Cons. Temp [Celcius] 
    1013       REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in) ::   psal   ! salinity   [psu] 
    1014       ! Leave result array automatic rather than making explicitly allocated 
    1015       REAL(wp), DIMENSION(jpi,jpj,jpk) ::   ptmp   ! potential temperature [Celcius] 
    1016       ! 
    1017       INTEGER  ::   ji, jj, jk           ! dummy loop indices 
    1018       REAL(wp) ::   zt , zs , ztm        ! local scalars 
    1019       REAL(wp) ::   zn , zd              ! local scalars 
    1020       REAL(wp) ::   zdeltaS , z1_S0 , z1_T0 
    1021       !!---------------------------------------------------------------------- 
    1022       ! 
    1023       IF ( ln_timing )   CALL timing_start('eos_pt_from_ct_3d') 
    1024       ! 
    1025       zdeltaS = 5._wp 
    1026       z1_S0   = 0.875_wp/35.16504_wp 
    1027       z1_T0   = 1._wp/40._wp 
    1028       ! 
    1029       DO jk = 1, jpkm1 
    1030          DO jj = 1, jpj 
    1031             DO ji = 1, jpi 
    1032                ! 
    1033                zt  = ctmp   (ji,jj,jk) * z1_T0 
    1034                zs  = SQRT( ABS( psal(ji,jj,jk) + zdeltaS ) * r1_S0 ) 
    1035                ztm = tmask(ji,jj,jk) 
    1036                ! 
    1037                zn = ((((-2.1385727895e-01_wp*zt   & 
    1038                   &   - 2.7674419971e-01_wp*zs+1.0728094330_wp)*zt   & 
    1039                   &   + (2.6366564313_wp*zs+3.3546960647_wp)*zs-7.8012209473_wp)*zt   & 
    1040                   &   + ((1.8835586562_wp*zs+7.3949191679_wp)*zs-3.3937395875_wp)*zs-5.6414948432_wp)*zt   & 
    1041                   &   + (((3.5737370589_wp*zs-1.5512427389e+01_wp)*zs+2.4625741105e+01_wp)*zs   & 
    1042                   &      +1.9912291000e+01_wp)*zs-3.2191146312e+01_wp)*zt   & 
    1043                   &   + ((((5.7153204649e-01_wp*zs-3.0943149543_wp)*zs+9.3052495181_wp)*zs   & 
    1044                   &      -9.4528934807_wp)*zs+3.1066408996_wp)*zs-4.3504021262e-01_wp 
    1045                   ! 
    1046                zd = (2.0035003456_wp*zt   & 
    1047                   &   -3.4570358592e-01_wp*zs+5.6471810638_wp)*zt   & 
    1048                   &   + (1.5393993508_wp*zs-6.9394762624_wp)*zs+1.2750522650e+01_wp 
    1049                   ! 
    1050                ptmp(ji,jj,jk) = ( zt / z1_T0 + zn / zd ) * ztm 
    1051                   ! 
    1052             END DO 
    1053          END DO 
    1054       END DO 
    1055       ! 
    1056       IF( ln_timing )   CALL timing_stop('eos_pt_from_ct_3d') 
    1057       ! 
    1058    END FUNCTION eos_pt_from_ct_3d 
     992      IF( ln_timing )   CALL timing_stop('eos_pt_from_ct') 
     993      ! 
     994   END FUNCTION eos_pt_from_ct 
    1059995 
    1060996 
     
    17151651 
    17161652         r1_S0  = 0.875_wp/35.16504_wp   ! Used to convert CT in potential temperature when using bulk formulae (eos_pt_from_ct) 
    1717  
     1653          
    17181654         IF(lwp) THEN 
    17191655            WRITE(numout,*) 
  • NEMO/trunk/src/OCE/TRA/traadv.F90

    r11989 r11993  
    134134      ! 
    135135!!gm ??? 
    136       CALL dia_ptr( kt, zvn )                                    ! diagnose the effective MSF  
     136      IF( ln_diaptr )   CALL dia_ptr( zvn )                                    ! diagnose the effective MSF  
    137137!!gm ??? 
    138138      ! 
    139  
    140139      IF( l_trdtra )   THEN                    !* Save ta and sa trends 
    141140         ALLOCATE( ztrdt(jpi,jpj,jpk), ztrds(jpi,jpj,jpk) ) 
  • NEMO/trunk/src/OCE/TRA/traadv_cen.F90

    r11989 r11993  
    6161      !! ** Action : - update pta  with the now advective tracer trends 
    6262      !!             - send trends to trdtra module for further diagnostcs (l_trdtra=T) 
    63       !!             - poleward advective heat and salt transport (ln_diaptr=T) 
     63      !!             - htr_adv, str_adv : poleward advective heat and salt transport (ln_diaptr=T) 
    6464      !!---------------------------------------------------------------------- 
    6565      INTEGER                              , INTENT(in   ) ::   kt              ! ocean time-step index 
     
    8989      l_hst = .FALSE. 
    9090      l_ptr = .FALSE. 
    91       IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR. ( cdtype == 'TRC' .AND. l_trdtrc ) )       l_trd = .TRUE. 
    92       IF(   cdtype == 'TRA' .AND. ( iom_use( 'sophtadv' ) .OR. iom_use( 'sophtadv' ) )  )    l_ptr = .TRUE.  
     91      IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR. ( cdtype == 'TRC' .AND. l_trdtrc ) )        l_trd = .TRUE. 
     92      IF(   cdtype == 'TRA' .AND. ln_diaptr )                                                 l_ptr = .TRUE.  
    9393      IF(   cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. & 
    94          &                          iom_use("uadv_salttr") .OR. iom_use("vadv_salttr")  ) )  l_hst = .TRUE. 
     94         &                          iom_use("uadv_salttr") .OR. iom_use("vadv_salttr")  ) )   l_hst = .TRUE. 
    9595      ! 
    9696      !                     
  • NEMO/trunk/src/OCE/TRA/traadv_fct.F90

    r11989 r11993  
    6868      !! ** Action : - update pta  with the now advective tracer trends 
    6969      !!             - send trends to trdtra module for further diagnostics (l_trdtra=T) 
    70       !!             - poleward advective heat and salt transport (ln_diaptr=T) 
     70      !!             - htr_adv, str_adv : poleward advective heat and salt transport (ln_diaptr=T) 
    7171      !!---------------------------------------------------------------------- 
    7272      INTEGER                              , INTENT(in   ) ::   kt              ! ocean time-step index 
     
    101101      l_ptr = .FALSE. 
    102102      ll_zAimp = .FALSE. 
    103       IF( ( cdtype == 'TRA' .AND. l_trdtra  ) .OR. ( cdtype =='TRC' .AND. l_trdtrc ) )      l_trd = .TRUE. 
    104       IF(   cdtype == 'TRA' .AND. ( iom_use( 'sophtadv' ) .OR. iom_use( 'sophtadv' ) ) )    l_ptr = .TRUE.  
    105       IF(   cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR.  & 
     103      IF( ( cdtype =='TRA' .AND. l_trdtra  ) .OR. ( cdtype =='TRC' .AND. l_trdtrc ) )       l_trd = .TRUE. 
     104      IF(   cdtype =='TRA' .AND. ln_diaptr )                                                l_ptr = .TRUE.  
     105      IF(   cdtype =='TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR.  & 
    106106         &                         iom_use("uadv_salttr") .OR. iom_use("vadv_salttr")  ) )  l_hst = .TRUE. 
    107107      ! 
  • NEMO/trunk/src/OCE/TRA/traadv_mus.F90

    r11989 r11993  
    6868      !! ** Action : - update pta  with the now advective tracer trends 
    6969      !!             - send trends to trdtra module for further diagnostcs (l_trdtra=T) 
    70       !!             - poleward advective heat and salt transport (ln_diaptr=T) 
     70      !!             - htr_adv, str_adv : poleward advective heat and salt transport (ln_diaptr=T) 
    7171      !! 
    7272      !! References : Estubier, A., and M. Levy, Notes Techn. Pole de Modelisation 
     
    120120      l_ptr = .FALSE. 
    121121      IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR. ( cdtype == 'TRC' .AND. l_trdtrc ) )      l_trd = .TRUE. 
    122       IF(   cdtype == 'TRA' .AND. ( iom_use( 'sophtadv' ) .OR. iom_use( 'sophtadv' ) )  )   l_ptr = .TRUE.  
     122      IF(   cdtype == 'TRA' .AND. ln_diaptr )                                               l_ptr = .TRUE.  
    123123      IF(   cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. & 
    124124         &                          iom_use("uadv_salttr") .OR. iom_use("vadv_salttr")  ) ) l_hst = .TRUE. 
  • NEMO/trunk/src/OCE/TRA/traadv_qck.F90

    r11989 r11993  
    2121   USE trdtra          ! trends manager: tracers  
    2222   USE diaptr          ! poleward transport diagnostics 
    23    USE iom 
    2423   ! 
    2524   USE in_out_manager  ! I/O manager 
     
    8180      !! ** Action : - update pta  with the now advective tracer trends 
    8281      !!             - send trends to trdtra module for further diagnostcs (l_trdtra=T) 
    83       !!             - poleward advective heat and salt transport (ln_diaptr=T) 
     82      !!             - htr_adv, str_adv : poleward advective heat and salt transport (ln_diaptr=T) 
    8483      !! 
    8584      !! ** Reference : Leonard (1979, 1991) 
     
    104103      l_trd = .FALSE. 
    105104      l_ptr = .FALSE. 
    106       IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR. ( cdtype == 'TRC' .AND. l_trdtrc ) )   l_trd = .TRUE. 
    107       IF(   cdtype == 'TRA' .AND. ( iom_use( 'sophtadv' ) .OR. iom_use( 'sophtadv' ) ) ) l_ptr = .TRUE.  
     105      IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR. ( cdtype == 'TRC' .AND. l_trdtrc ) )      l_trd = .TRUE. 
     106      IF(   cdtype == 'TRA' .AND. ln_diaptr )                                              l_ptr = .TRUE.  
    108107      ! 
    109108      ! 
  • NEMO/trunk/src/OCE/TRA/traadv_ubs.F90

    r11989 r11993  
    7979      !! ** Action : - update pta  with the now advective tracer trends 
    8080      !!             - send trends to trdtra module for further diagnostcs (l_trdtra=T) 
    81       !!             - poleward advective heat and salt transport (ln_diaptr=T) 
     81      !!             - htr_adv, str_adv : poleward advective heat and salt transport (ln_diaptr=T) 
    8282      !! 
    8383      !! Reference : Shchepetkin, A. F., J. C. McWilliams, 2005, Ocean Modelling, 9, 347-404.  
     
    111111      l_ptr = .FALSE. 
    112112      IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR. ( cdtype == 'TRC' .AND. l_trdtrc ) )      l_trd = .TRUE. 
    113       IF(   cdtype == 'TRA' .AND. ( iom_use( 'sophtadv' ) .OR. iom_use( 'sophtadv' ) )  )   l_ptr = .TRUE.  
     113      IF(   cdtype == 'TRA' .AND. ln_diaptr )                                               l_ptr = .TRUE.  
    114114      IF(   cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. & 
    115115         &                          iom_use("uadv_salttr") .OR. iom_use("vadv_salttr")  ) ) l_hst = .TRUE. 
  • NEMO/trunk/src/OCE/TRA/trabbc.F90

    r11989 r11993  
    100100      ENDIF 
    101101      ! 
    102       CALL iom_put ( "hfgeou" , rau0_rcp * qgh_trd0(:,:) ) 
    103       ! 
    104102      IF(ln_ctl)   CALL prt_ctl( tab3d_1=tsa(:,:,:,jp_tem), clinfo1=' bbc  - Ta: ', mask1=tmask, clinfo3='tra-ta' ) 
    105103      ! 
  • NEMO/trunk/src/OCE/TRA/traldf_iso.F90

    r11989 r11993  
    124124      l_hst = .FALSE. 
    125125      l_ptr = .FALSE. 
    126       IF( cdtype == 'TRA' .AND. ( iom_use( 'sophtldf' ) .OR. iom_use( 'sopstldf' ) ) )     l_ptr = .TRUE.  
     126      IF( cdtype == 'TRA' .AND. ln_diaptr )                                                 l_ptr = .TRUE.  
    127127      IF( cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. & 
    128128         &                        iom_use("uadv_salttr") .OR. iom_use("vadv_salttr")  ) )   l_hst = .TRUE. 
  • NEMO/trunk/src/OCE/TRA/traldf_lap_blp.F90

    r11989 r11993  
    8989      l_hst = .FALSE. 
    9090      l_ptr = .FALSE. 
    91       IF( cdtype == 'TRA' .AND. ( iom_use( 'sophtldf' ) .OR. iom_use( 'sopstldf' ) ) )     l_ptr = .TRUE. 
     91      IF( cdtype == 'TRA' .AND. ln_diaptr )                                                l_ptr = .TRUE.  
    9292      IF( cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. & 
    9393         &                        iom_use("uadv_salttr") .OR. iom_use("vadv_salttr")  ) )  l_hst = .TRUE. 
  • NEMO/trunk/src/OCE/TRA/traldf_triad.F90

    r11989 r11993  
    110110      l_hst = .FALSE. 
    111111      l_ptr = .FALSE. 
    112       IF( cdtype == 'TRA' .AND. ( iom_use( 'sophtldf' ) .OR. iom_use( 'sopstldf' ) ) )      l_ptr = .TRUE.  
     112      IF( cdtype == 'TRA' .AND. ln_diaptr )                                                 l_ptr = .TRUE.  
    113113      IF( cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. & 
    114114         &                        iom_use("uadv_salttr") .OR. iom_use("vadv_salttr")  ) )   l_hst = .TRUE. 
  • NEMO/trunk/src/OCE/nemogcm.F90

    r11989 r11993  
    479479                           CALL     flo_init    ! drifting Floats 
    480480      IF( ln_diacfl    )   CALL dia_cfl_init    ! Initialise CFL diagnostics 
    481 !                           CALL dia_ptr_init    ! Poleward TRansports initialization 
     481                           CALL dia_ptr_init    ! Poleward TRansports initialization 
    482482                           CALL dia_dct_init    ! Sections tranports 
    483483                           CALL dia_hsb_init    ! heat content, salt content and volume budgets 
  • NEMO/trunk/src/OCE/step.F90

    r11989 r11993  
    154154      IF( l_ldftra_time .OR. l_ldfeiv_time )   CALL ldf_tra( kstp )       !       and/or eiv coeff. 
    155155      IF( l_ldfdyn_time                    )   CALL ldf_dyn( kstp )       ! eddy viscosity coeff.  
     156 
    156157      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
    157158      !  Ocean dynamics : hdiv, ssh, e3, u, v, w 
     
    187188         IF(.NOT.ln_linssh) CALL dom_vvl_sf_nxt( kstp, kcall=2 )  ! after vertical scale factors (update depth average component) 
    188189      ENDIF 
    189                             CALL dyn_zdf    ( kstp )  ! vertical diffusion 
    190       IF( ln_dynspg_ts ) THEN                         ! vertical scale factors and vertical velocity need to be updated 
     190                         CALL dyn_zdf       ( kstp )  ! vertical diffusion 
     191 
     192      IF( ln_dynspg_ts ) THEN                          
    191193                            CALL wzv        ( kstp )              ! now cross-level velocity  
    192194         IF( ln_zad_Aimp )  CALL wAimp      ( kstp )  ! Adaptive-implicit vertical advection partitioning 
    193195      ENDIF 
    194        
    195196 
    196197      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
     
    204205      IF( ln_floats  )   CALL flo_stp ( kstp )        ! drifting Floats 
    205206      IF( ln_diacfl  )   CALL dia_cfl ( kstp )        ! Courant number diagnostics 
    206                          CALL dia_hth ( kstp )        ! Thermocline depth (20 degres isotherm depth) 
     207      IF( lk_diahth  )   CALL dia_hth ( kstp )        ! Thermocline depth (20 degres isotherm depth) 
    207208      IF( ln_diadct  )   CALL dia_dct ( kstp )        ! Transports 
    208209                         CALL dia_ar5 ( kstp )        ! ar5 diag 
    209                          CALL dia_ptr ( kstp )        ! Poleward adv/ldf TRansports diagnostics 
    210210      IF( ln_diaharm )   CALL dia_harm( kstp )        ! Tidal harmonic analysis 
    211211                         CALL dia_wri ( kstp )        ! ocean model: outputs 
     
    243243                         CALL tra_ldf       ( kstp )  ! lateral mixing 
    244244 
     245!!gm : why CALL to dia_ptr has been moved here??? (use trends info?) 
     246      IF( ln_diaptr  )   CALL dia_ptr                 ! Poleward adv/ldf TRansports diagnostics 
     247!!gm 
    245248                         CALL tra_zdf       ( kstp )  ! vertical mixing and after tracer fields 
    246249      IF( ln_zdfnpc  )   CALL tra_npc       ( kstp )  ! update after fields by non-penetrative convection 
  • NEMO/trunk/src/TOP/PISCES/P4Z/p4zflx.F90

    r11989 r11993  
    160160            zfld = zfco2 * chemc(ji,jj,1) * zkgco2(ji,jj)  ! (mol/L) * (m/s) 
    161161            zflu = zh2co3(ji,jj) * zkgco2(ji,jj)                                   ! (mol/L) (m/s) ? 
    162             oce_co2(ji,jj) = ( zfld - zflu ) * tmask(ji,jj,1)  
     162            oce_co2(ji,jj) = ( zfld - zflu ) * rfact2 * e1e2t(ji,jj) * tmask(ji,jj,1) * 1000. 
    163163            ! compute the trend 
    164             tra(ji,jj,1,jpdic) = tra(ji,jj,1,jpdic) + oce_co2(ji,jj) * rfact2 / e3t_n(ji,jj,1) 
     164            tra(ji,jj,1,jpdic) = tra(ji,jj,1,jpdic) + ( zfld - zflu ) * rfact2 / e3t_n(ji,jj,1) * tmask(ji,jj,1) 
    165165 
    166166            ! Compute O2 flux  
     
    174174      IF( iom_use("tcflx") .OR. iom_use("tcflxcum") .OR. kt == nitrst   & 
    175175         &                 .OR. (ln_check_mass .AND. kt == nitend) )    & 
    176          t_oce_co2_flx  = glob_sum( 'p4zflx', oce_co2(:,:) * e1e2t(:,:) * 1000. )                    !  Total Flux of Carbon 
     176         t_oce_co2_flx  = glob_sum( 'p4zflx', oce_co2(:,:) )                    !  Total Flux of Carbon 
    177177      t_oce_co2_flx_cum = t_oce_co2_flx_cum + t_oce_co2_flx       !  Cumulative Total Flux of Carbon 
    178178!      t_atm_co2_flx     = glob_sum( 'p4zflx', satmco2(:,:) * e1e2t(:,:) )       ! Total atmospheric pCO2 
     
    186186 
    187187      IF( lk_iomput .AND. knt == nrdttrc ) THEN 
    188          CALL iom_put( "AtmCo2" , satmco2(:,:) * tmask(:,:,1) )   ! Atmospheric CO2 concentration 
    189          !  
    190188         ALLOCATE( zw2d(jpi,jpj) )   
    191189         IF( iom_use( "Cflx"  ) )  THEN 
    192             zw2d(:,:) = oce_co2(:,:) * 1000.  ! conversion in molC/m2/s 
     190            zw2d(:,:) = oce_co2(:,:) / e1e2t(:,:) * rfact2r 
    193191            CALL iom_put( "Cflx"     , zw2d )  
    194192         ENDIF 
    195193         IF( iom_use( "Oflx"  ) )  THEN 
    196             zw2d(:,:) =  zoflx(:,:) * 1000. 
     194            zw2d(:,:) =  zoflx(:,:) * 1000 * tmask(:,:,1) 
    197195            CALL iom_put( "Oflx" , zw2d ) 
    198196         ENDIF 
     
    205203           CALL iom_put( "Dpco2" ,  zw2d ) 
    206204         ENDIF 
    207          IF( iom_use( "pCO2sea" ) ) THEN 
    208            zw2d(:,:) = ( zh2co3(:,:) / ( chemc(:,:,1) + rtrn ) ) * tmask(:,:,1) 
    209            CALL iom_put( "pCO2sea" ,  zw2d ) 
    210          ENDIF 
    211  
    212205         IF( iom_use( "Dpo2" ) )  THEN 
    213206           zw2d(:,:) = ( atcox * patm(:,:) - atcox * trb(:,:,1,jpoxy) / ( chemo2(:,:,1) + rtrn ) ) * tmask(:,:,1) 
    214207           CALL iom_put( "Dpo2"  , zw2d ) 
    215208         ENDIF 
    216          CALL iom_put( "tcflx"    , t_oce_co2_flx     )   ! molC/s 
    217          CALL iom_put( "tcflxcum" , t_oce_co2_flx_cum )   ! molC 
     209         CALL iom_put( "tcflx"    , t_oce_co2_flx * rfact2r )   ! molC/s 
     210         CALL iom_put( "tcflxcum" , t_oce_co2_flx_cum       )   ! molC 
    218211         ! 
    219212         DEALLOCATE( zw2d ) 
  • NEMO/trunk/src/TOP/PISCES/P4Z/p4zsms.F90

    r11989 r11993  
    3535   INTEGER ::    numco2, numnut, numnit      ! logical unit for co2 budget 
    3636   REAL(wp) ::   alkbudget, no3budget, silbudget, ferbudget, po4budget 
    37    REAL(wp) ::   xfact, xfact1, xfact2, xfact3 
     37   REAL(wp) ::   xfact1, xfact2, xfact3 
    3838 
    3939   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   xnegtr     ! Array used to indicate negative tracer values 
     
    6363      REAL(wp) ::  ztra 
    6464      CHARACTER (len=25) :: charout 
    65       REAL(wp), ALLOCATABLE, DIMENSION(:,:    ) :: zw2d 
    66       REAL(wp), ALLOCATABLE, DIMENSION(:,:,:  ) :: zw3d 
    67       REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) ::   ztrdt   ! 4D workspace 
    68  
    6965      !!--------------------------------------------------------------------- 
    7066      ! 
     
    8985      rfact = r2dttrc 
    9086      ! 
    91       ! trends computation initialisation 
    92       IF( l_trdtrc )  THEN 
    93          ALLOCATE( ztrdt(jpi,jpj,jpk,jp_pisces) )  !* store now fields before applying the Asselin filter 
    94          ztrdt(:,:,:,:)  = trn(:,:,:,:) 
    95       ENDIF 
    96       ! 
    97  
    9887      IF( ( ln_top_euler .AND. kt == nittrc000 )  .OR. ( .NOT.ln_top_euler .AND. kt <= nittrc000 + nn_dttrc ) ) THEN 
    9988         rfactr  = 1. / rfact 
     
    10190         rfact2r = 1. / rfact2 
    10291         xstep = rfact2 / rday         ! Time step duration for biology 
    103          xfact = 1.e+3 * rfact2r 
    10492         IF(lwp) WRITE(numout,*)  
    10593         IF(lwp) WRITE(numout,*) '    Passive Tracer  time step    rfact  = ', rfact, ' rdt = ', rdt 
     
    146134         END DO 
    147135        ! 
    148         IF(  iom_use( 'INTdtAlk' ) .OR. iom_use( 'INTdtDIC' ) .OR. iom_use( 'INTdtFer' ) .OR.  & 
    149           &  iom_use( 'INTdtDIN' ) .OR. iom_use( 'INTdtDIP' ) .OR. iom_use( 'INTdtSil' ) )  THEN 
    150           ! 
    151           ALLOCATE( zw3d(jpi,jpj,jpk), zw2d(jpi,jpj) ) 
    152           zw3d(:,:,jpk) = 0. 
    153           DO jk = 1, jpkm1 
    154               zw3d(:,:,jk) = xnegtr(:,:,jk) * xfact * e3t_n(:,:,jk) * tmask(:,:,jk) 
    155           ENDDO 
    156           ! 
    157           zw2d(:,:) = 0. 
    158           DO jk = 1, jpkm1 
    159              zw2d(:,:) = zw2d(:,:) + zw3d(:,:,jk) * tra(:,:,jk,jptal) 
    160           ENDDO 
    161           CALL iom_put( 'INTdtAlk', zw2d ) 
    162           ! 
    163           zw2d(:,:) = 0. 
    164           DO jk = 1, jpkm1 
    165              zw2d(:,:) = zw2d(:,:) + zw3d(:,:,jk) * tra(:,:,jk,jpdic) 
    166           ENDDO 
    167           CALL iom_put( 'INTdtDIC', zw2d ) 
    168           ! 
    169           zw2d(:,:) = 0. 
    170           DO jk = 1, jpkm1 
    171              zw2d(:,:) = zw2d(:,:) + zw3d(:,:,jk) * rno3 * ( tra(:,:,jk,jpno3) + tra(:,:,jk,jpnh4) ) 
    172           ENDDO 
    173           CALL iom_put( 'INTdtDIN', zw2d ) 
    174           ! 
    175           zw2d(:,:) = 0. 
    176           DO jk = 1, jpkm1 
    177              zw2d(:,:) = zw2d(:,:) + zw3d(:,:,jk) * po4r * tra(:,:,jk,jppo4) 
    178           ENDDO 
    179           CALL iom_put( 'INTdtDIP', zw2d ) 
    180           ! 
    181           zw2d(:,:) = 0. 
    182           DO jk = 1, jpkm1 
    183              zw2d(:,:) = zw2d(:,:) + zw3d(:,:,jk) * tra(:,:,jk,jpfer) 
    184           ENDDO 
    185           CALL iom_put( 'INTdtFer', zw2d ) 
    186           ! 
    187           zw2d(:,:) = 0. 
    188           DO jk = 1, jpkm1 
    189              zw2d(:,:) = zw2d(:,:) + zw3d(:,:,jk) * tra(:,:,jk,jpsil) 
    190           ENDDO 
    191           CALL iom_put( 'INTdtSil', zw2d ) 
    192           ! 
    193           DEALLOCATE( zw3d, zw2d ) 
    194         ENDIF 
    195         ! 
    196136         DO jn = jp_pcs0, jp_pcs1 
    197137            tra(:,:,:,jn) = 0._wp 
     
    204144         ENDIF 
    205145      END DO 
     146 
    206147      ! 
    207148      IF( l_trdtrc ) THEN 
    208149         DO jn = jp_pcs0, jp_pcs1 
    209            ztrdt(:,:,:,jn) = ( trb(:,:,:,jn) - ztrdt(:,:,:,jn) ) * rfact2r  
    210            CALL trd_trc( ztrdt(:,:,:,jn), jn, jptra_sms, kt )   ! save trends 
     150           CALL trd_trc( tra(:,:,:,jn), jn, jptra_sms, kt )   ! save trends 
    211151         END DO 
    212          DEALLOCATE( ztrdt )  
    213152      END IF 
    214153#endif 
Note: See TracChangeset for help on using the changeset viewer.