Index: /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/CONFIG/ORCA2_LIM/EXP00/iodef.xml
===================================================================
--- /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/CONFIG/ORCA2_LIM/EXP00/iodef.xml	(revision 3315)
+++ /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/CONFIG/ORCA2_LIM/EXP00/iodef.xml	(revision 3316)
@@ -19,5 +19,5 @@
       <group id="grid_T" axis_ref="none" grid_ref="grid_T">
 	<field id="toce"         description="temperature"                               unit="degC" axis_ref="deptht"   />
-       	<field id="soce"         description="salinity"                                  unit="psu"  axis_ref="deptht"   />
+   <field id="soce"         description="salinity"                                  unit="psu"  axis_ref="deptht"   />
 	<field id="sst"          description="sea surface temperature"                   unit="degC"                     />
 	<field id="sst2"         description="square of sea surface temperature"         unit="degC2"                    />
@@ -117,5 +117,29 @@
     <field id="Mmx_v"          description="Mm current barotrope along j-axis harmonic real part "      unit="m/s"        />
     <field id="Mmy_v"          description="Mm current barotrope along j-axis harmonic imaginary part " unit="m/s"        />
-     </group>
+    <!-- variables available with key_trdtra -->
+    <field id="ttrd_xad"     description="temperature-trend: i-advection"           unit="degC/s" axis_ref="deptht" />
+    <field id="strd_xad"     description="salinity   -trend: i-advection"           unit="psu/s"  axis_ref="deptht" />
+    <field id="ttrd_yad"     description="temperature-trend: j-advection"           unit="degC/s" axis_ref="deptht" />
+    <field id="strd_yad"     description="salinity   -trend: j-advection"           unit="psu/s"  axis_ref="deptht" />
+    <field id="ttrd_zad"     description="temperature-trend: k-advection"           unit="degC/s" axis_ref="deptht" />
+    <field id="strd_zad"     description="salinity   -trend: k-advection"           unit="psu/s"  axis_ref="deptht" />
+    <field id="ttrd_sad"     description="temperature-trend: surface adv. (no-vvl)" unit="degC/s"                   />
+    <field id="strd_sad"     description="salinity   -trend: surface adv. (no-vvl)" unit="psu/s"                    />
+    <field id="ttrd_ldf"     description="temperature-trend: lateral  diffusion"    unit="degC/s" axis_ref="deptht" />
+    <field id="strd_ldf"     description="salinity   -trend: lateral  diffusion"    unit="psu/s"  axis_ref="deptht" />
+    <field id="ttrd_zdf"     description="temperature-trend: vertical diffusion"    unit="degC/s" axis_ref="deptht" />
+    <field id="strd_zdf"     description="salinity   -trend: vertical diffusion"    unit="psu/s"  axis_ref="deptht" />
+    <field id="ttrd_zdfp"    description="temperature-trend: pure vert. diffusion"  unit="degC/s" axis_ref="deptht" />
+    <field id="strd_zdfp"    description="salinity   -trend: pure vert. diffusion"  unit="psu/s"  axis_ref="deptht" />
+    <field id="ttrd_dmp"     description="temperature-trend: interior restoring"    unit="degC/s" axis_ref="deptht" />
+    <field id="strd_dmp"     description="salinity   -trend: interior restoring"    unit="psu/s"  axis_ref="deptht" />
+    <field id="ttrd_bbl"     description="temperature-trend: bottom boundary layer" unit="degC/s" axis_ref="deptht" />
+    <field id="strd_bbl"     description="salinity   -trend: bottom boundary layer" unit="psu/s"  axis_ref="deptht" />
+    <field id="ttrd_npc"     description="temperature-trend: non-penetrative conv." unit="degC/s" axis_ref="deptht" />
+    <field id="strd_npc"     description="salinity   -trend: non-penetrative conv." unit="psu/s"  axis_ref="deptht" />
+    <field id="ttrd_qns"     description="temperature-trend: non-solar heating"     unit="degC/s"                   />
+    <field id="ttrd_qsr"     description="temperature-trend: solar penetr. heating" unit="degC/s"  axis_ref="deptht" />
+    <field id="ttrd_bbc"     description="temperature-trend: geothermal heating"    unit="degC/s"  axis_ref="deptht" />
+      </group>
 
       <!-- SBC -->
@@ -215,5 +239,17 @@
 	<field id="ueiv_heattr"  description="ocean bolus heat transport along i-axis"     unit="W"    axis_ref="none" />
 	<field id="udiff_heattr" description="ocean diffusion heat transport along i-axis" unit="W"    axis_ref="none" />
-     </group>
+   <!-- variables available with key_trddyn -->
+   <field id="utrd_hpg"     description="i-trend: hydrostatic pressure gradient"      unit="m/s^2"                />
+   <field id="utrd_spg"     description="i-trend: surface     pressure gradient"      unit="m/s^2"                />
+   <field id="utrd_keg"     description="i-trend: KE gradient         or hor. adv."   unit="m/s^2"                />
+   <field id="utrd_rvo"     description="i-trend: relative  vorticity or metric term" unit="m/s^2"                />
+   <field id="utrd_pvo"     description="i-trend: planetary vorticity"                unit="m/s^2"                />
+   <field id="utrd_zad"     description="i-trend: vertical  advection"                unit="m/s^2"                />
+   <field id="utrd_udx"     description="i-trend: U.dx[U]"                            unit="m/s^2"                />
+   <field id="utrd_ldf"     description="i-trend: lateral   diffusion"                unit="m/s^2"                />
+   <field id="utrd_zdf"     description="i-trend: vertical  diffusion"                unit="m/s^2"                />
+   <field id="utrd_tau"     description="i-trend: wind stress "                       unit="m/s^2" axis_ref="none"/>
+   <field id="utrd_bfr"     description="i-trend: bottom friction"                    unit="m/s^2"                />   
+      </group>
       
       <!-- V grid -->
@@ -233,4 +269,16 @@
 	<field id="veiv_heattr"  description="ocean bolus heat transport along j-axis"     unit="W"    axis_ref="none" />
 	<field id="vdiff_heattr" description="ocean diffusion heat transport along j-axis" unit="W"    axis_ref="none" />
+   <!-- variables available with key_trddyn -->
+   <field id="vtrd_hpg"     description="j-trend: hydrostatic pressure gradient"      unit="m/s^2"                />
+   <field id="vtrd_spg"     description="j-trend: surface     pressure gradient"      unit="m/s^2"                />
+   <field id="vtrd_keg"     description="j-trend: KE gradient         or hor. adv."   unit="m/s^2"                />
+   <field id="vtrd_rvo"     description="j-trend: relative  vorticity or metric term" unit="m/s^2"                />
+   <field id="vtrd_pvo"     description="j-trend: planetary vorticity"                unit="m/s^2"                />
+   <field id="vtrd_zad"     description="j-trend: vertical  advection"                unit="m/s^2"                />
+   <field id="vtrd_vdy"     description="i-trend: V.dx[V]"                            unit="m/s^2"                />
+   <field id="vtrd_ldf"     description="j-trend: lateral   diffusion"                unit="m/s^2"                />
+   <field id="vtrd_zdf"     description="j-trend: vertical  diffusion"                unit="m/s^2"                />
+   <field id="vtrd_tau"     description="j-trend: wind stress "                       unit="m/s^2" axis_ref="none"/>
+   <field id="vtrd_bfr"     description="j-trend: bottom friction"                    unit="m/s^2"                /> 
       </group>
       
Index: /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/CONFIG/ORCA2_LIM/EXP00/namelist
===================================================================
--- /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/CONFIG/ORCA2_LIM/EXP00/namelist	(revision 3315)
+++ /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/CONFIG/ORCA2_LIM/EXP00/namelist	(revision 3316)
@@ -814,5 +814,13 @@
 !              !       or mixed-layer trends or barotropic vorticity    ("key_trdmld" or     "key_trdvor")
 !-----------------------------------------------------------------------
-   nn_trd      = 365       !  time step frequency dynamics and tracers trends
+   ln_3D_trd_d = .FALSE.   ! (T) 3D momentum             trends or (F) not
+   ln_3D_trd_t = .FALSE.   ! (T) 3D tracer               trends or (F) not
+   ln_ML_trd_d = .FALSE.   ! (T) 2D tracer   trends averaged over the mixed layer 
+   ln_ML_trd_t = .FALSE.   ! (T) 2D momentum trends averaged over the mixed layer 
+   ln_PE_trd   = .FALSE.   ! (T) 3D Potential Energy     trends or (F) not
+   ln_KE_trd   = .FALSE.   ! (T) 3D Kinetic   Energy     trends or (F) not
+   ln_vor_trd  = .FALSE.   ! (T) 3D barotropic vorticity trends or (F) not
+   ln_glo_trd  = .FALSE.   ! (T) global domain averaged diag for T, T^2, KE, and PE
+   nn_trd      = 365       !  print frequency (ln_glo_trd=T) (unit=time step)
    nn_ctls     =   0       !  control surface type in mixed-layer trends (0,1 or n<jpk)
    rn_ucf      =   1.      !  unit conversion factor (=1 -> /seconds ; =86400. -> /day)
Index: /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/DYN/dynadv_cen2.F90
===================================================================
--- /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/DYN/dynadv_cen2.F90	(revision 3315)
+++ /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/DYN/dynadv_cen2.F90	(revision 3316)
@@ -103,5 +103,5 @@
          zfu_uw(:,:,:) = ua(:,:,:) - zfu_uw(:,:,:)
          zfv_vw(:,:,:) = va(:,:,:) - zfv_vw(:,:,:)
-         CALL trd_mod( zfu_uw, zfv_vw, jpdyn_trd_had, 'DYN', kt )
+         CALL trd_mod( zfu_uw, zfv_vw, jpdyn_trd_keg, 'DYN', kt )
          zfu_t(:,:,:) = ua(:,:,:)
          zfv_t(:,:,:) = va(:,:,:)
Index: /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/DYN/dynadv_ubs.F90
===================================================================
--- /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/DYN/dynadv_ubs.F90	(revision 3315)
+++ /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/DYN/dynadv_ubs.F90	(revision 3316)
@@ -196,5 +196,5 @@
          zfu_uw(:,:,:) = ua(:,:,:) - zfu_uw(:,:,:)
          zfv_vw(:,:,:) = va(:,:,:) - zfv_vw(:,:,:)
-         CALL trd_mod( zfu_uw, zfv_vw, jpdyn_trd_had, 'DYN', kt )
+         CALL trd_mod( zfu_uw, zfv_vw, jpdyn_trd_keg, 'DYN', kt )
          zfu_t(:,:,:) = ua(:,:,:)
          zfv_t(:,:,:) = va(:,:,:)
Index: /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/DYN/dynnxt.F90
===================================================================
--- /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/DYN/dynnxt.F90	(revision 3315)
+++ /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/DYN/dynnxt.F90	(revision 3316)
@@ -17,4 +17,5 @@
    !!            3.3  !  2010-09  (D. Storkey, E.O'Dea) Bug fix for BDY module
    !!            3.3  !  2011-03  (P. Oddo) Bug fix for time-splitting+(BDY-OBC) and not VVL
+	!!            3.4  !  2012-02  (G. Madec) add the diagnostic of the time filter trends
    !!-------------------------------------------------------------------------
   
@@ -25,4 +26,5 @@
    USE dom_oce         ! ocean space and time domain
    USE sbc_oce         ! Surface boundary condition: ocean fields
+   USE trdmod_oce      ! ocean trends
    USE phycst          ! physical constants
    USE dynspg_oce      ! type of surface pressure gradient
@@ -38,12 +40,13 @@
    USE bdyvol          ! ocean open boundary condition (bdy_vol routines)
    USE in_out_manager  ! I/O manager
+   USE iom             ! I/O manager library
    USE lbclnk          ! lateral boundary condition (or mpp link)
    USE lib_mpp         ! MPP library
    USE wrk_nemo        ! Memory Allocation
    USE prtctl          ! Print control
+   USE timing          ! Timing
 #if defined key_agrif
    USE agrif_opa_interp
 #endif
-   USE timing          ! Timing
 
    IMPLICIT NONE
@@ -81,5 +84,5 @@
       !!             at the local domain boundaries through lbc_lnk call,
       !!             at the one-way open boundaries (lk_obc=T),
-      !!             at the AGRIF zoom     boundaries (lk_agrif=T)
+      !!             at the AGRIF zoom   boundaries (lk_agrif=T)
       !!
       !!              * Apply the time filter applied and swap of the dynamics
@@ -101,12 +104,12 @@
       REAL(wp) ::   z2dt         ! temporary scalar
 #endif
-      REAL(wp) ::   zue3a, zue3n, zue3b, zuf, zec   ! local scalars
-      REAL(wp) ::   zve3a, zve3n, zve3b, zvf        !   -      -
-      REAL(wp), POINTER, DIMENSION(:,:,:) ::  ze3u_f, ze3v_f 
+      REAL(wp) ::   zue3a, zue3n, zue3b, zuf, zec      ! local scalars
+      REAL(wp) ::   zve3a, zve3n, zve3b, zvf, z1_2dt   !   -      -
+      REAL(wp), POINTER, DIMENSION(:,:,:) ::  ze3u_f, ze3v_f, zua, zva 
       !!----------------------------------------------------------------------
       !
       IF( nn_timing == 1 )  CALL timing_start('dyn_nxt')
       !
-      CALL wrk_alloc( jpi,jpj,jpk, ze3u_f, ze3v_f )
+      CALL wrk_alloc( jpi,jpj,jpk, ze3u_f, ze3v_f, zua, zva )
       !
       IF( kt == nit000 ) THEN
@@ -186,4 +189,16 @@
 # endif
 #endif
+
+      IF( ln_3D_trd_d ) THEN             ! 3D output: total momentum trends a prepare the atf trend computation
+         z1_2dt = 1._wp / (2. * rdt)            ! Euler or leap-frog time step 
+         IF( neuler == 0 .AND. kt == nit000 )  z1_2dt = 1._wp / rdt
+         zua(:,:,:) = ( ua(:,:,:) - ub(:,:,:) ) * z1_2dt
+         zva(:,:,:) = ( va(:,:,:) - vb(:,:,:) ) * z1_2dt
+         CALL iom_put( "utrd_tot", zua )        ! total momentum trends (but the asselin time filter)
+         CALL iom_put( "vtrd_tot", zva )
+         zua(:,:,:) = un(:,:,:)                 ! save the before velocity before the asselin filter
+         zva(:,:,:) = vn(:,:,:)                 ! (caution: there is a shift by 1 timestep in the
+         !                                      !  computation of the asselin filter trends)
+      ENDIF
 
       ! Time filter and swap of dynamics arrays
@@ -276,8 +291,15 @@
       ENDIF
 
+      IF( ln_3D_trd_d ) THEN             ! 3D output: asselin filter trends on momentum
+        zua(:,:,:) = ( ub(:,:,:) - zua(:,:,:) ) * z1_2dt
+        zva(:,:,:) = ( vb(:,:,:) - zva(:,:,:) ) * z1_2dt
+        CALL iom_put( "utrd_atf", zua )        ! asselin filter trends on momentum 
+        CALL iom_put( "vtrd_atf", zva )
+      ENDIF
+
       IF(ln_ctl)   CALL prt_ctl( tab3d_1=un, clinfo1=' nxt  - Un: ', mask1=umask,   &
          &                       tab3d_2=vn, clinfo2=' Vn: '       , mask2=vmask )
       ! 
-      CALL wrk_dealloc( jpi,jpj,jpk, ze3u_f, ze3v_f )
+      CALL wrk_dealloc( jpi,jpj,jpk, ze3u_f, ze3v_f, zua, zva )
       !
       IF( nn_timing == 1 )  CALL timing_stop('dyn_nxt')
Index: /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/DYN/dynvor.F90
===================================================================
--- /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/DYN/dynvor.F90	(revision 3315)
+++ /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/DYN/dynvor.F90	(revision 3316)
@@ -15,4 +15,5 @@
    !!            3.2  ! 2009-04  (R. Benshila)  vvl: correction of een scheme
    !!            3.3  ! 2010-10  (C. Ethe, G. Madec) reorganisation of initialisation phase
+	!!            3.5  ! 2012-02  (G. Madec) suppress jpdyn_trd_dat vorticity trend simplification
    !!----------------------------------------------------------------------
 
@@ -115,5 +116,4 @@
             ztrdv(:,:,:) = va(:,:,:) - ztrdv(:,:,:)
             CALL trd_mod( ztrdu, ztrdv, jpdyn_trd_pvo, 'DYN', kt )
-            CALL trd_mod( ztrdu, ztrdv, jpdyn_trd_dat, 'DYN', kt )
          ELSE
             CALL vor_ene( kt, ntot, ua, va )                ! total vorticity
@@ -134,5 +134,4 @@
             ztrdv(:,:,:) = va(:,:,:) - ztrdv(:,:,:)
             CALL trd_mod( ztrdu, ztrdv, jpdyn_trd_pvo, 'DYN', kt )
-            CALL trd_mod( ztrdu, ztrdv, jpdyn_trd_dat, 'DYN', kt )
          ELSE
             CALL vor_ens( kt, ntot, ua, va )                ! total vorticity
@@ -153,5 +152,4 @@
             ztrdv(:,:,:) = va(:,:,:) - ztrdv(:,:,:)
             CALL trd_mod( ztrdu, ztrdv, jpdyn_trd_pvo, 'DYN', kt )
-            CALL trd_mod( ztrdu, ztrdv, jpdyn_trd_dat, 'DYN', kt )
          ELSE
             CALL vor_mix( kt )                               ! total vorticity (mix=ens-ene)
@@ -172,5 +170,4 @@
             ztrdv(:,:,:) = va(:,:,:) - ztrdv(:,:,:)
             CALL trd_mod( ztrdu, ztrdv, jpdyn_trd_pvo, 'DYN', kt )
-            CALL trd_mod( ztrdu, ztrdv, jpdyn_trd_dat, 'DYN', kt )
          ELSE
             CALL vor_een( kt, ntot, ua, va )                ! total vorticity
Index: /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/TRA/traldf.F90
===================================================================
--- /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/TRA/traldf.F90	(revision 3315)
+++ /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/TRA/traldf.F90	(revision 3316)
@@ -35,6 +35,6 @@
    PRIVATE
 
-   PUBLIC   tra_ldf         ! called by step.F90 
-   PUBLIC   tra_ldf_init    ! called by opa.F90 
+   PUBLIC   tra_ldf        ! called by step.F90 
+   PUBLIC   tra_ldf_init   ! called by opa.F90 
    !
    INTEGER ::   nldf = 0   ! type of lateral diffusion used defined from ln_traldf_... namlist logicals)
Index: /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/TRD/trdicp.F90
===================================================================
--- /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/TRD/trdicp.F90	(revision 3315)
+++ /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/TRD/trdicp.F90	(revision 3316)
@@ -5,4 +5,5 @@
    !!=====================================================================
    !! History :  1.0  !  2004-08 (C. Talandier) New trends organization
+   !!            3.5  !  2012-02 (G. Madec)  add 3D tracer zdf trend output using iom
    !!----------------------------------------------------------------------
 #if  defined key_trdtra   ||   defined key_trddyn   ||   defined key_esopa
@@ -22,10 +23,11 @@
    USE ldfdyn_oce      ! ocean dynamics: lateral physics
    USE zdf_oce         ! ocean vertical physics
-   USE in_out_manager  ! I/O manager
-   USE lib_mpp         ! distibuted memory computing library
+   USE zdfddm          ! ocean vertical physics: double diffusion
    USE eosbn2          ! equation of state
    USE phycst          ! physical constants
+   USE lib_mpp         ! distibuted memory computing library
+   USE in_out_manager  ! I/O manager
+   USE iom             ! I/O manager library
    USE wrk_nemo        ! Memory allocation
-
 
    IMPLICIT NONE
@@ -44,4 +46,5 @@
 #  include "domzgr_substitute.h90"
 #  include "vectopt_loop_substitute.h90"
+#  include "zdfddm_substitute.h90"
    !!----------------------------------------------------------------------
    !! NEMO/OPA 3.3 , NEMO Consortium (2010)
@@ -55,11 +58,11 @@
       !!                  ***  ROUTINE trd_2d  ***
       !! 
-      !! ** Purpose : verify the basin averaged properties of tracers and/or
-      !!              momentum equations at every time step frequency nn_trd.
-      !!----------------------------------------------------------------------
-      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) ::   ptrd2dx   ! Temperature or U trend 
-      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) ::   ptrd2dy   ! Salinity    or V trend
-      INTEGER                     , INTENT(in   ) ::   ktrd      ! tracer trend index
-      CHARACTER(len=3)            , INTENT(in   ) ::   ctype     ! momentum ('DYN') or tracers ('TRA') trends
+      !! ** Purpose :   compute and print the domain averaged properties of tracers 
+      !!              and/or momentum equations at each nn_trd time step.
+      !!----------------------------------------------------------------------
+      REAL(wp), DIMENSION(:,:), INTENT(inout) ::   ptrd2dx   ! Temperature or U trend 
+      REAL(wp), DIMENSION(:,:), INTENT(inout) ::   ptrd2dy   ! Salinity    or V trend
+      INTEGER                 , INTENT(in   ) ::   ktrd      ! tracer trend index
+      CHARACTER(len=3)        , INTENT(in   ) ::   ctype     ! momentum ('DYN') or tracers ('TRA') trends
       !!
       INTEGER ::   ji, jj   ! loop indices
@@ -123,8 +126,8 @@
       !!              momentum equations at every time step frequency nn_trd.
       !!----------------------------------------------------------------------
-      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   ptrd3dx   ! Temperature or U trend 
-      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   ptrd3dy   ! Salinity    or V trend
-      INTEGER,                          INTENT(in   ) ::   ktrd      ! momentum or tracer trend index
-      CHARACTER(len=3),                 INTENT(in   ) ::   ctype     ! momentum ('DYN') or tracers ('TRA') trends
+      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   ptrd3dx   ! Temperature or U trend 
+      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   ptrd3dy   ! Salinity    or V trend
+      INTEGER,                    INTENT(in   ) ::   ktrd      ! momentum or tracer trend index
+      CHARACTER(len=3),           INTENT(in   ) ::   ctype     ! momentum ('DYN') or tracers ('TRA') trends
       !!
       INTEGER  ::   ji, jj, jk   ! dummy loop indices
@@ -253,5 +256,4 @@
       !! ** Purpose :  write dynamic trends in ocean.output 
       !!----------------------------------------------------------------------
-      !
       INTEGER, INTENT(in) ::   kt   ! ocean time-step index
       !
@@ -261,4 +263,6 @@
       !!----------------------------------------------------------------------
 
+      IF( ln_3D_trd_t .AND. ln_3D_trd_d )   RETURN            ! do nothing if 3D output with IOM
+
       CALL wrk_alloc( jpi, jpj, jpk, zkx, zky, zkz, zkepe )
 
@@ -266,5 +270,5 @@
       ! -------------------
 
-      IF( MOD(kt,nn_trd) == 0 .OR. kt == nit000 .OR. kt == nitend ) THEN
+      IF( MOD( kt, nn_trd ) == 0 .OR. kt == nit000 .OR. kt == nitend ) THEN
 
          ! I.1 Conversion potential energy - kinetic energy
@@ -449,5 +453,38 @@
       !!----------------------------------------------------------------------
       INTEGER, INTENT(in) ::   kt   ! ocean time-step index
-      !!----------------------------------------------------------------------
+      !
+      INTEGER  ::   jk   ! loop indices
+      REAL(wp), POINTER, DIMENSION(:,:,:)  ::   zwt, zws, ztrdt, ztrds   ! 3D workspace
+      !!----------------------------------------------------------------------
+
+
+      IF( ln_3D_trd_t ) THEN      ! 3D output: treat the vertical diffusion trends (if iso)
+         !
+         CALL wrk_alloc( jpi, jpj, jpk, zwt, zws, ztrdt, ztrds )
+         !
+         IF( ln_traldf_iso ) THEN      ! iso-neutral diffusion : re-compute the PURE vertical diffusive trend
+            !                                 !  zdf trends using now field (called after the swap)
+            zwt(:,:, 1 ) = 0._wp   ;   zws(:,:, 1 ) = 0._wp            ! vertical diffusive fluxes
+            zwt(:,:,jpk) = 0._wp   ;   zws(:,:,jpk) = 0._wp
+            DO jk = 2, jpk
+               zwt(:,:,jk) =   avt(:,:,jk) * ( tsn(:,:,jk-1,jp_tem) - tsn(:,:,jk,jp_tem) ) / fse3w(:,:,jk) * tmask(:,:,jk)
+               zws(:,:,jk) = fsavs(:,:,jk) * ( tsn(:,:,jk-1,jp_sal) - tsn(:,:,jk,jp_sal) ) / fse3w(:,:,jk) * tmask(:,:,jk)
+            END DO
+            !
+            ztrdt(:,:,jpk) = 0._wp   ;   ztrds(:,:,jpk) = 0._wp
+            DO jk = 1, jpkm1
+               ztrdt(:,:,jk) = ( zwt(:,:,jk) - zwt(:,:,jk+1) ) / fse3t(:,:,jk)
+               ztrds(:,:,jk) = ( zws(:,:,jk) - zws(:,:,jk+1) ) / fse3t(:,:,jk)
+            END DO
+            CALL iom_put( "ttrd_zdfp", ztrdt )        ! PURE vertical diffusion (no isoneutral contribution)
+            CALL iom_put( "strd_zdfp", ztrds )
+         ENDIF
+         !
+         CALL wrk_dealloc( jpi, jpj, jpk, zwt, zws, ztrdt, ztrds )
+         !
+         RETURN                     ! do nothing else if 3D output with IOM
+         !
+      ENDIF
+
 
       ! I. Tracers trends
Index: /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/TRD/trdicp_oce.F90
===================================================================
--- /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/TRD/trdicp_oce.F90	(revision 3315)
+++ /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/TRD/trdicp_oce.F90	(revision 3316)
@@ -4,8 +4,6 @@
    !! Ocean trends :   set tracer and momentum trend variables
    !!======================================================================
-
-   !!----------------------------------------------------------------------
-   !!   'key_trdtra'   or                         tracer trends diagnostics
-   !!   'key_trddyn'                            momentum trends diagnostics
+   !! History :  1.0  !  2004-08 (C. Talandier) New trends organization
+   !!            3.5  !  2012-02 (G. Madec)  add 3D tracer zdf trend output using iom
    !!----------------------------------------------------------------------
    USE par_oce                 ! ocean parameters
Index: /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/TRD/trdmod.F90
===================================================================
--- /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/TRD/trdmod.F90	(revision 3315)
+++ /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/TRD/trdmod.F90	(revision 3316)
@@ -6,9 +6,12 @@
    !! History :  1.0  !  2004-08  (C. Talandier) Original code
    !!             -   !  2005-04  (C. Deltel)    Add Asselin trend in the ML budget
-   !!            3.3  ! 2010-10  (C. Ethe, G. Madec) reorganisation of initialisation phase
+   !!            3.3  !  2010-10  (C. Ethe, G. Madec) reorganisation of initialisation phase
+   !!            3.5  !  2012-02  (G. Madec) add 3D trends output for T, S, U, V, PE and KE
    !!----------------------------------------------------------------------
 #if  defined key_trdtra || defined key_trddyn || defined key_trdmld || defined key_trdvor || defined key_esopa
    !!----------------------------------------------------------------------
-   !!   trd_mod          : Call the trend to be computed
+   !!   trd_mod          : manage the type of trend diagnostics
+   !!   trd_3Diom        : output 3D momentum and/or tracer trends using IOM
+   !!   trd_budget       : domain averaged budget of trends (including kinetic energy and tracer variance trends)
    !!   trd_mod_init     : Initialization step
    !!----------------------------------------------------------------------
@@ -24,7 +27,7 @@
    USE trdmld                  ! ocean active mixed layer tracers trends 
    USE in_out_manager          ! I/O manager
+   USE iom					   ! I/O manager library
    USE lib_mpp                 ! MPP library
    USE wrk_nemo                ! Memory allocation
-
 
    IMPLICIT NONE
@@ -44,5 +47,4 @@
    !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
    !!----------------------------------------------------------------------
-
 CONTAINS
 
@@ -51,99 +53,40 @@
       !!                  ***  ROUTINE trd_mod  ***
       !! 
-      !! ** Purpose : Dispatch all trends computation, e.g. vorticity, mld or 
-      !!              integral constraints
-      !!----------------------------------------------------------------------
-      !
+      !! ** Purpose :   Dispatch all trends computation, e.g. 3D output, integral
+      !!                constraints, barotropic vorticity, kinetic enrgy, 
+      !!                potential energy, and/or mixed layer budget.
+      !!----------------------------------------------------------------------
       REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   ptrdx   ! Temperature or U trend 
       REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   ptrdy   ! Salinity    or V trend
+      INTEGER                   , INTENT(in   ) ::   ktrd    ! tracer trend index
       CHARACTER(len=3)          , INTENT(in   ) ::   ctype   ! momentum or tracers trends type 'DYN'/'TRA'
       INTEGER                   , INTENT(in   ) ::   kt      ! time step
-      INTEGER                   , INTENT(in   ) ::   ktrd    ! tracer trend index
       !!
       INTEGER ::   ji, jj   ! dummy loop indices
-      REAL(wp), POINTER, DIMENSION(:,:)  :: ztswu, ztswv, ztbfu, ztbfv, z2dx, z2dy 
-      !!----------------------------------------------------------------------
-
-      CALL wrk_alloc( jpi, jpj, ztswu, ztswv, ztbfu, ztbfv, z2dx, z2dy )
-
-      z2dx(:,:) = 0._wp   ;   z2dy(:,:) = 0._wp                            ! initialization of workspace arrays
-
-      IF( neuler == 0 .AND. kt == nit000    ) THEN   ;   r2dt =      rdt   ! = rdtra (restart with Euler time stepping)
-      ELSEIF(               kt <= nit000 + 1) THEN   ;   r2dt = 2. * rdt   ! = 2 rdttra (leapfrog)
-      ENDIF
-
-      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
-      ! I. Integral Constraints Properties for momentum and/or tracers trends
-      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
-
-      IF( ( mod(kt,nn_trd) == 0 .OR. kt == nit000 .OR. kt == nitend) )   THEN
-         !
-         IF( lk_trdtra .AND. ctype == 'TRA' )   THEN       ! active tracer trends
-            SELECT CASE ( ktrd )
-            CASE ( jptra_trd_ldf )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpt_ldf, ctype )   ! lateral diff
-            CASE ( jptra_trd_zdf )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpt_zdf, ctype )   ! vertical diff (Kz)
-            CASE ( jptra_trd_bbc )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpt_bbc, ctype )   ! bottom boundary cond
-            CASE ( jptra_trd_bbl )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpt_bbl, ctype )   ! bottom boundary layer
-            CASE ( jptra_trd_npc )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpt_npc, ctype )   ! static instability mixing
-            CASE ( jptra_trd_dmp )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpt_dmp, ctype )   ! damping
-            CASE ( jptra_trd_qsr )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpt_qsr, ctype )   ! penetrative solar radiat.
-            CASE ( jptra_trd_nsr )   ;   z2dx(:,:) = ptrdx(:,:,1)   
-                                         z2dy(:,:) = ptrdy(:,:,1)
-                                         CALL trd_icp( z2dx , z2dy , jpicpt_nsr, ctype )   ! non solar radiation
-            CASE ( jptra_trd_xad )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpt_xad, ctype )   ! x- horiz adv
-            CASE ( jptra_trd_yad )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpt_yad, ctype )   ! y- horiz adv
-            CASE ( jptra_trd_zad )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpt_zad, ctype )   ! z- vertical adv 
-                                         CALL trd_icp( ptrdx, ptrdy, jpicpt_zad, ctype )   
-                                         ! compute the surface flux condition wn(:,:,1)*tsn(:,:,1,jp_tem)
-                                         z2dx(:,:) = wn(:,:,1)*tsn(:,:,1,jp_tem)/fse3t(:,:,1)
-                                         z2dy(:,:) = wn(:,:,1)*tsn(:,:,1,jp_sal)/fse3t(:,:,1)
-                                         CALL trd_icp( z2dx , z2dy , jpicpt_zl1, ctype )   ! 1st z- vertical adv 
-            END SELECT
-         END IF
-
-         IF( lk_trddyn .AND. ctype == 'DYN' )   THEN       ! momentum trends 
-            !
-            SELECT CASE ( ktrd )
-            CASE ( jpdyn_trd_hpg )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpd_hpg, ctype )   ! hydrost. pressure grad
-            CASE ( jpdyn_trd_keg )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpd_keg, ctype )   ! KE gradient 
-            CASE ( jpdyn_trd_rvo )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpd_rvo, ctype )   ! relative vorticity 
-            CASE ( jpdyn_trd_pvo )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpd_pvo, ctype )   ! planetary vorticity
-            CASE ( jpdyn_trd_ldf )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpd_ldf, ctype )   ! lateral diffusion 
-            CASE ( jpdyn_trd_had )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpd_had, ctype )   ! horizontal advection 
-            CASE ( jpdyn_trd_zad )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpd_zad, ctype )   ! vertical advection 
-            CASE ( jpdyn_trd_spg )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpd_spg, ctype )   ! surface pressure grad.
-            CASE ( jpdyn_trd_dat )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpd_dat, ctype )   ! damping term
-            CASE ( jpdyn_trd_zdf )                                                         ! vertical diffusion 
-               ! subtract surface forcing/bottom friction trends 
-               ! from vertical diffusive momentum trends
-               ztswu(:,:) = 0._wp   ;   ztswv(:,:) = 0._wp
-               ztbfu(:,:) = 0._wp   ;   ztbfv(:,:) = 0._wp 
-               DO jj = 2, jpjm1   
-                  DO ji = fs_2, fs_jpim1   ! vector opt.
-                     ! save the surface forcing momentum fluxes
-                     ztswu(ji,jj) = utau(ji,jj) / ( fse3u(ji,jj,1)*rau0 )
-                     ztswv(ji,jj) = vtau(ji,jj) / ( fse3v(ji,jj,1)*rau0 )
-                     ! bottom friction contribution now handled explicitly
-                     ptrdx(ji,jj,1) = ptrdx(ji,jj,1) - ztswu(ji,jj)
-                     ptrdy(ji,jj,1) = ptrdy(ji,jj,1) - ztswv(ji,jj)
-                  END DO
-               END DO
-               !
-               CALL trd_icp( ptrdx, ptrdy, jpicpd_zdf, ctype )   
-               CALL trd_icp( ztswu, ztswv, jpicpd_swf, ctype )                               ! wind stress forcing term
-               ! bottom friction contribution now handled explicitly
-            CASE ( jpdyn_trd_bfr )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpd_bfr, ctype )     ! bottom friction term
-            END SELECT
-            !
-         END IF
-         !
-      END IF
-
-      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
-      ! II. Vorticity trends
-      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
-
-      IF( lk_trdvor .AND. ctype == 'DYN' )   THEN
-         !
+      REAL(wp), POINTER, DIMENSION(:,:) ::   ztswu, ztswv    ! 2D workspace 
+      !!----------------------------------------------------------------------
+
+      CALL wrk_alloc( jpi, jpj, ztswu, ztswv )
+
+      IF( neuler == 0 .AND. kt == nit000    ) THEN   ;   r2dt =      rdt      ! = rdtra (restart with Euler time stepping)
+      ELSEIF(               kt <= nit000 + 1) THEN   ;   r2dt = 2. * rdt      ! = 2 rdttra (leapfrog)
+      ENDIF
+
+      !                                            !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
+      IF( ln_3D_trd_d .OR. ln_3D_trd_t ) THEN      !   3D output of momentum and/or tracers trends using IOM interface
+         !                                         !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
+         CALL trd_3Diom ( ptrdx, ptrdy, ktrd, ctype, kt )
+         !
+      ENDIF
+         !                                         !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
+      IF( ln_glo_trd ) THEN                        ! I. Integral Constraints Properties for momentum and/or tracers trends
+         !                                         !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
+         CALL trd_budget( ptrdx, ptrdy, ktrd, ctype, kt )
+         !
+      ENDIF
+
+      !                                            !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
+      IF( lk_trdvor .AND. ctype == 'DYN' ) THEN    ! II. Vorticity trends
+         !                                         !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
          SELECT CASE ( ktrd ) 
          CASE ( jpdyn_trd_hpg )   ;   CALL trd_vor_zint( ptrdx, ptrdy, jpvor_prg )   ! Hydrostatique Pressure Gradient 
@@ -152,28 +95,19 @@
          CASE ( jpdyn_trd_pvo )   ;   CALL trd_vor_zint( ptrdx, ptrdy, jpvor_pvo )   ! Planetary Vorticity Term 
          CASE ( jpdyn_trd_ldf )   ;   CALL trd_vor_zint( ptrdx, ptrdy, jpvor_ldf )   ! Horizontal Diffusion 
-         CASE ( jpdyn_trd_had )   ;   CALL ctl_warn('Vorticity for horizontal advection trend never checked')   
          CASE ( jpdyn_trd_zad )   ;   CALL trd_vor_zint( ptrdx, ptrdy, jpvor_zad )   ! Vertical Advection 
          CASE ( jpdyn_trd_spg )   ;   CALL trd_vor_zint( ptrdx, ptrdy, jpvor_spg )   ! Surface Pressure Grad. 
-         CASE ( jpdyn_trd_dat )   ;   CALL trd_vor_zint( ptrdx, ptrdy, jpvor_bev )   ! Beta V  
          CASE ( jpdyn_trd_zdf )                                                      ! Vertical Diffusion 
-            ! subtract surface forcing/bottom friction trends 
-            ! from vertical diffusive momentum trends
             ztswu(:,:) = 0.e0   ;   ztswv(:,:) = 0.e0
-            ztbfu(:,:) = 0.e0   ;   ztbfv(:,:) = 0.e0  
-            DO jj = 2, jpjm1   
+            DO jj = 2, jpjm1                                                             ! wind stress trends
                DO ji = fs_2, fs_jpim1   ! vector opt.
-                  ! save the surface forcing momentum fluxes
-                  ztswu(ji,jj) = utau(ji,jj) / ( fse3u(ji,jj,1)*rau0 )
-                  ztswv(ji,jj) = vtau(ji,jj) / ( fse3v(ji,jj,1)*rau0 )
-                  !
-                  ptrdx(ji,jj,1     ) = ptrdx(ji,jj,1     ) - ztswu(ji,jj)
-                  ptrdy(ji,jj,1     ) = ptrdy(ji,jj,1     ) - ztswv(ji,jj)
+                  ztswu(ji,jj) = ( utau_b(ji,jj) + utau(ji,jj) ) / ( fse3u(ji,jj,1) * rau0 )
+                  ztswv(ji,jj) = ( vtau_b(ji,jj) + vtau(ji,jj) ) / ( fse3v(ji,jj,1) * rau0 )
                END DO
             END DO
             !
-            CALL trd_vor_zint( ptrdx, ptrdy, jpvor_zdf )   
-            CALL trd_vor_zint( ztswu, ztswv, jpvor_swf )                               ! Wind stress forcing term
+            CALL trd_vor_zint( ptrdx, ptrdy, jpvor_zdf )                             ! zdf trend including surf./bot. stresses 
+            CALL trd_vor_zint( ztswu, ztswv, jpvor_swf )                             ! surface wind stress 
          CASE ( jpdyn_trd_bfr )
-            CALL trd_vor_zint( ptrdx, ptrdy, jpvor_bfr )                               ! Bottom friction term
+            CALL trd_vor_zint( ptrdx, ptrdy, jpvor_bfr )                             ! Bottom stress
          END SELECT
          !
@@ -184,6 +118,5 @@
       !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
 
-      IF( lk_trdmld .AND. ctype == 'TRA' )   THEN
-         
+      IF( lk_trdmld .AND. ctype == 'TRA' )   THEN   
          !-----------------------------------------------------------------------------------------------
          ! W.A.R.N.I.N.G :
@@ -198,47 +131,217 @@
 
          SELECT CASE ( ktrd )
-         CASE ( jptra_trd_xad )   ;   CALL trd_mld_zint( ptrdx, ptrdy, jpmld_xad, '3D' )   ! merid. advection
-         CASE ( jptra_trd_yad )   ;   CALL trd_mld_zint( ptrdx, ptrdy, jpmld_yad, '3D' )   ! zonal  advection
-         CASE ( jptra_trd_zad )   ;   CALL trd_mld_zint( ptrdx, ptrdy, jpmld_zad, '3D' )   ! vertical advection
-         CASE ( jptra_trd_ldf )   ;   CALL trd_mld_zint( ptrdx, ptrdy, jpmld_ldf, '3D' )   ! lateral diffusive
-         CASE ( jptra_trd_bbl )   ;   CALL trd_mld_zint( ptrdx, ptrdy, jpmld_bbl, '3D' )   ! bottom boundary layer
+         CASE ( jptra_trd_xad )        ;   CALL trd_mld_zint( ptrdx, ptrdy, jpmld_xad, '3D' )   ! zonal    advection
+         CASE ( jptra_trd_yad )        ;   CALL trd_mld_zint( ptrdx, ptrdy, jpmld_yad, '3D' )   ! merid.   advection
+         CASE ( jptra_trd_zad )        ;   CALL trd_mld_zint( ptrdx, ptrdy, jpmld_zad, '3D' )   ! vertical advection
+         CASE ( jptra_trd_ldf )        ;   CALL trd_mld_zint( ptrdx, ptrdy, jpmld_ldf, '3D' )   ! lateral  diffusion
+         CASE ( jptra_trd_bbl )        ;   CALL trd_mld_zint( ptrdx, ptrdy, jpmld_bbl, '3D' )   ! bottom boundary layer
          CASE ( jptra_trd_zdf )
-            IF( ln_traldf_iso )   THEN
-               CALL trd_mld_zint( ptrdx, ptrdy, jpmld_ldf, '3D' )   ! vertical diffusion (K_z)
-            ELSE
-               CALL trd_mld_zint( ptrdx, ptrdy, jpmld_zdf, '3D' )   ! vertical diffusion (K_z)
+            IF( ln_traldf_iso ) THEN   ;   CALL trd_mld_zint( ptrdx, ptrdy, jpmld_ldf, '3D' )   ! lateral  diffusion (K_z)
+            ELSE                       ;   CALL trd_mld_zint( ptrdx, ptrdy, jpmld_zdf, '3D' )   ! vertical diffusion (K_z)
             ENDIF
-         CASE ( jptra_trd_dmp )   ;   CALL trd_mld_zint( ptrdx, ptrdy, jpmld_dmp, '3D' )   ! internal 3D restoring (tradmp)
-         CASE ( jptra_trd_qsr )   ;   CALL trd_mld_zint( ptrdx, ptrdy, jpmld_for, '3D' )   ! air-sea : penetrative sol radiat
+         CASE ( jptra_trd_dmp )        ;   CALL trd_mld_zint( ptrdx, ptrdy, jpmld_dmp, '3D' )   ! internal 3D restoring (tradmp)
+         CASE ( jptra_trd_qsr )        ;   CALL trd_mld_zint( ptrdx, ptrdy, jpmld_for, '3D' )   ! air-sea : penetrative sol radiat
          CASE ( jptra_trd_nsr )
-            ptrdx(:,:,2:jpk) = 0.e0   ;   ptrdy(:,:,2:jpk) = 0.e0
-            CALL trd_mld_zint( ptrdx, ptrdy, jpmld_for, '2D' )                             ! air-sea : non penetr sol radiat
-         CASE ( jptra_trd_bbc )   ;   CALL trd_mld_zint( ptrdx, ptrdy, jpmld_bbc, '3D' )   ! bottom bound cond (geoth flux)
-         CASE ( jptra_trd_atf )   ;   CALL trd_mld_zint( ptrdx, ptrdy, jpmld_atf, '3D' )   ! asselin numerical
-         CASE ( jptra_trd_npc )   ;   CALL trd_mld_zint( ptrdx, ptrdy, jpmld_npc, '3D' )   ! non penetr convect adjustment
+            ptrdx(:,:,2:jpk) = 0._wp   ;   ptrdy(:,:,2:jpk) = 0._wp
+            CALL trd_mld_zint( ptrdx, ptrdy, jpmld_for, '2D' )                                  ! air-sea : non penetr sol radiat
+         CASE ( jptra_trd_bbc )        ;   CALL trd_mld_zint( ptrdx, ptrdy, jpmld_bbc, '3D' )   ! bottom bound cond (geoth flux)
+         CASE ( jptra_trd_atf )        ;   CALL trd_mld_zint( ptrdx, ptrdy, jpmld_atf, '3D' )   ! asselin numerical
+         CASE ( jptra_trd_npc )        ;   CALL trd_mld_zint( ptrdx, ptrdy, jpmld_npc, '3D' )   ! non penetr convect adjustment
          END SELECT
-
-      ENDIF
-      !
-      CALL wrk_dealloc( jpi, jpj, ztswu, ztswv, ztbfu, ztbfv, z2dx, z2dy )
+         !
+      ENDIF
+      !
+      CALL wrk_dealloc( jpi, jpj, ztswu, ztswv )
       !
    END SUBROUTINE trd_mod
 
+
+   SUBROUTINE trd_budget( ptrdx, ptrdy, ktrd, ctype, kt )
+      !!---------------------------------------------------------------------
+      !!                  ***  ROUTINE trd_budget  ***
+      !! 
+      !! ** Purpose : integral constraint diagnostics for momentum and/or tracer trends
+      !!              
+      !!----------------------------------------------------------------------
+      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   ptrdx   ! Temperature or U trend 
+      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   ptrdy   ! Salinity    or V trend
+      INTEGER                   , INTENT(in   ) ::   ktrd    ! tracer trend index
+      CHARACTER(len=3)          , INTENT(in   ) ::   ctype   ! momentum or tracers trends type 'DYN'/'TRA'
+      INTEGER                   , INTENT(in   ) ::   kt      ! time step
+      !!
+      INTEGER ::   ji, jj   ! dummy loop indices
+      REAL(wp), POINTER, DIMENSION(:,:)  :: ztswu, ztswv, z2dx, z2dy   ! 2D workspace 
+      !!----------------------------------------------------------------------
+
+      CALL wrk_alloc( jpi, jpj, ztswu, ztswv, z2dx, z2dy )
+
+      IF( MOD(kt,nn_trd) == 0 .OR. kt == nit000 .OR. kt == nitend ) THEN
+         !
+         IF( lk_trdtra .AND. ctype == 'TRA' ) THEN       ! active tracer trends
+            SELECT CASE ( ktrd )
+            CASE ( jptra_trd_ldf )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpt_ldf, ctype )   ! lateral diff
+            CASE ( jptra_trd_zdf )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpt_zdf, ctype )   ! vertical diff (Kz)
+            CASE ( jptra_trd_bbc )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpt_bbc, ctype )   ! bottom boundary cond
+            CASE ( jptra_trd_bbl )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpt_bbl, ctype )   ! bottom boundary layer
+            CASE ( jptra_trd_npc )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpt_npc, ctype )   ! static instability mixing
+            CASE ( jptra_trd_dmp )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpt_dmp, ctype )   ! damping
+            CASE ( jptra_trd_qsr )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpt_qsr, ctype )   ! penetrative solar radiat.
+            CASE ( jptra_trd_nsr )   ;   z2dx(:,:) = ptrdx(:,:,1)                          ! non solar radiation
+                                         z2dy(:,:) = ptrdy(:,:,1)
+                                         CALL trd_icp( z2dx , z2dy , jpicpt_nsr, ctype )
+            CASE ( jptra_trd_xad )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpt_xad, ctype )   ! x- horiz adv
+            CASE ( jptra_trd_yad )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpt_yad, ctype )   ! y- horiz adv
+            CASE ( jptra_trd_zad )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpt_zad, ctype )   ! z- vertical adv 
+                                         !                                                 ! surface flux
+                                         IF( lk_vvl ) THEN                                      ! variable volume = zero
+                                            z2dx(:,:) = 0._wp
+                                            z2dy(:,:) = 0._wp
+                                         ELSE                                                   ! constant volume = wn*tsn/e3t
+                                            z2dx(:,:) = wn(:,:,1) * tsn(:,:,1,jp_tem) / fse3t(:,:,1)
+                                            z2dy(:,:) = wn(:,:,1) * tsn(:,:,1,jp_sal) / fse3t(:,:,1)
+                                         ENDIF
+                                         CALL trd_icp( z2dx , z2dy , jpicpt_zl1, ctype ) 
+            END SELECT
+         ENDIF
+
+         IF( lk_trddyn .AND. ctype == 'DYN' ) THEN       ! momentum trends 
+            !
+            SELECT CASE ( ktrd )
+            CASE( jpdyn_trd_hpg )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpd_hpg, ctype )   ! hydrost. pressure gradient
+            CASE( jpdyn_trd_spg )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpd_spg, ctype )   ! surface pressure grad.
+            CASE( jpdyn_trd_pvo )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpd_pvo, ctype )   ! planetary vorticity
+            CASE( jpdyn_trd_rvo )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpd_rvo, ctype )   ! relative  vorticity or metric term
+            CASE( jpdyn_trd_keg )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpd_keg, ctype )   ! KE gradient         or hor. advection
+            CASE( jpdyn_trd_zad )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpd_zad, ctype )   ! vertical  advection 
+            CASE( jpdyn_trd_ldf )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpd_ldf, ctype )   ! lateral  diffusion 
+            CASE( jpdyn_trd_zdf )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpd_zdf, ctype )   ! vertical diffusion (icluding bfr & tau)
+                                        ztswu(:,:) = ( utau_b(ji,jj) + utau(ji,jj) ) / ( fse3u(:,:,1) * rau0 )
+                                        ztswv(:,:) = ( vtau_b(ji,jj) + vtau(ji,jj) ) / ( fse3v(:,:,1) * rau0 )
+                                        CALL trd_icp( ztswu, ztswv, jpicpd_swf, ctype )   ! wind stress trends
+            CASE( jpdyn_trd_bfr )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpd_bfr, ctype )   ! bottom friction trends
+            END SELECT
+            !
+         ENDIF
+         !
+      ENDIF
+      !
+      CALL wrk_dealloc( jpi, jpj, ztswu, ztswv, z2dx, z2dy )
+      !
+   END SUBROUTINE trd_budget
+
+
+   SUBROUTINE trd_3Diom( ptrdx, ptrdy, ktrd, ctype, kt )
+      !!---------------------------------------------------------------------
+      !!                  ***  ROUTINE trd_3Diom  ***
+      !! 
+      !! ** Purpose :   output 3D trends using IOM
+      !!----------------------------------------------------------------------
+      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   ptrdx   ! Temperature or U trend 
+      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   ptrdy   ! Salinity    or V trend
+      INTEGER                   , INTENT(in   ) ::   ktrd    ! tracer trend index
+      CHARACTER(len=3)          , INTENT(in   ) ::   ctype   ! momentum or tracers trends type 'DYN'/'TRA'
+      INTEGER                   , INTENT(in   ) ::   kt      ! time step
+      !!
+      INTEGER ::   ji, jj, jk   ! dummy loop indices
+      REAL(wp), POINTER, DIMENSION(:,:)   ::   z2dx, z2dy, ztswu, ztswv    ! 2D workspace 
+      REAL(wp), POINTER, DIMENSION(:,:,:) ::   z3dx, z3dy                  ! 3D workspace 
+      !!----------------------------------------------------------------------
+
+       IF( lk_trdtra .AND. ctype == 'TRA' ) THEN       ! active tracer trends
+         !
+!!gm Rq: mask the trends already masked in trd_tra, but lbc_lnk should probably be added
+         !
+         SELECT CASE( ktrd )
+         CASE( jptra_trd_xad )   ;   CALL iom_put( "ttrd_xad", ptrdx )        ! x- horizontal advection
+                                     CALL iom_put( "strd_xad", ptrdy )
+         CASE( jptra_trd_yad )   ;   CALL iom_put( "ttrd_yad", ptrdx )        ! y- horizontal advection
+                                     CALL iom_put( "strd_yad", ptrdy )
+         CASE( jptra_trd_zad )   ;   CALL iom_put( "ttrd_zad", ptrdx )        ! z- vertical   advection
+                                     CALL iom_put( "strd_zad", ptrdy )
+                                     IF( .NOT.lk_vvl ) THEN                   ! cst volume : adv flux through z=0 surface
+                                        z2dx(:,:) = wn(:,:,1) * tsn(:,:,1,jp_tem) / fse3t(:,:,1)
+                                        z2dy(:,:) = wn(:,:,1) * tsn(:,:,1,jp_sal) / fse3t(:,:,1)
+                                        CALL iom_put( "ttrd_sad", z2dx )
+                                        CALL iom_put( "strd_sad", z2dy )
+                                     ENDIF
+         CASE( jptra_trd_ldf )   ;   CALL iom_put( "ttrd_ldf", ptrdx )        ! lateral diffusion
+                                     CALL iom_put( "strd_ldf", ptrdy )
+         CASE( jptra_trd_zdf )   ;   CALL iom_put( "ttrd_zdf", ptrdx )        ! vertical diffusion (including Kz contribution)
+                                     CALL iom_put( "strd_zdf", ptrdy )
+         CASE( jptra_trd_dmp )   ;   CALL iom_put( "ttrd_dmp", ptrdx )        ! internal restoring (damping)
+                                     CALL iom_put( "strd_dmp", ptrdy )
+         CASE( jptra_trd_bbl )   ;   CALL iom_put( "ttrd_bbl", ptrdx )        ! bottom boundary layer
+                                     CALL iom_put( "strd_bbl", ptrdy )
+         CASE( jptra_trd_npc )   ;   CALL iom_put( "ttrd_npc", ptrdx )        ! static instability mixing
+                                     CALL iom_put( "strd_npc", ptrdy )
+         CASE( jptra_trd_qsr )   ;   CALL iom_put( "ttrd_qsr", ptrdx )        ! penetrative solar radiat. (only on temperature)
+         CASE( jptra_trd_nsr )   ;   CALL iom_put( "ttrd_qns", ptrdx(:,:,1) ) ! non-solar     radiation   (only on temperature)
+         CASE( jptra_trd_bbc )   ;   CALL iom_put( "ttrd_bbc", ptrdx )        ! geothermal heating   (only on temperature)
+         END SELECT
+      ENDIF
+
+      IF( lk_trddyn .AND. ctype == 'DYN' ) THEN       ! momentum trends 
+            !
+         ptrdx(:,:,:) = ptrdx(:,:,:) * umask(:,:,:)                       ! mask the trends
+         ptrdy(:,:,:) = ptrdy(:,:,:) * vmask(:,:,:)
+!!gm NB : here a lbc_lnk should probably be added
+         !
+         SELECT CASE( ktrd )
+         CASE( jpdyn_trd_hpg )   ;   CALL iom_put( "utrd_hpg", ptrdx )    ! hydrostatic pressure gradient
+                                     CALL iom_put( "vtrd_hpg", ptrdy )
+         CASE( jpdyn_trd_spg )   ;   CALL iom_put( "utrd_spg", ptrdx )    ! surface pressure gradient
+                                     CALL iom_put( "vtrd_spg", ptrdy )
+         CASE( jpdyn_trd_pvo )   ;   CALL iom_put( "utrd_pvo", ptrdx )    ! planetary vorticity
+                                     CALL iom_put( "vtrd_pvo", ptrdy )
+         CASE( jpdyn_trd_rvo )   ;   CALL iom_put( "utrd_rvo", ptrdx )    ! relative  vorticity     (or metric term)
+                                     CALL iom_put( "vtrd_rvo", ptrdy )
+         CASE( jpdyn_trd_keg )   ;   CALL iom_put( "utrd_keg", ptrdx )    ! Kinetic Energy gradient (or had)
+                                     CALL iom_put( "vtrd_keg", ptrdy )
+            z3dx(:,:,:) = 0._wp                                           ! U.dxU & V.dyV (approximation)
+            z3dy(:,:,:) = 0._wp
+            DO jk = 1, jpkm1														         ! no mask as un,vn are masked
+               DO jj = 2, jpjm1
+                  DO ji = 2, jpim1
+                     z3dx(ji,jj,jk) = un(ji,jj,jk) * ( un(ji+1,jj,jk) - un(ji-1,jj,jk) ) / ( 2._wp * e1u(ji,jj) )
+                     z3dy(ji,jj,jk) = vn(ji,jj,jk) * ( vn(ji,jj+1,jk) - vn(ji,jj-1,jk) ) / ( 2._wp * e2v(ji,jj) )
+                  END DO
+               END DO
+            END DO
+            CALL lbc_lnk( z3dx, 'U', -1. )   ;    CALL lbc_lnk( z3dy, 'V', -1. )
+                                     CALL iom_put( "utrd_udx", z3dx  ) 
+                                     CALL iom_put( "vtrd_vdy", z3dy  )
+         CASE( jpdyn_trd_zad )   ;   CALL iom_put( "utrd_zad", ptrdx )    ! vertical   advection
+                                     CALL iom_put( "vtrd_zad", ptrdy )
+         CASE( jpdyn_trd_ldf )   ;   CALL iom_put( "utrd_ldf", ptrdx )    ! lateral diffusion
+                                     CALL iom_put( "vtrd_ldf", ptrdy )
+         CASE( jpdyn_trd_zdf )   ;   CALL iom_put( "utrd_zdf", ptrdx )    ! vertical diffusion 
+                                     CALL iom_put( "vtrd_zdf", ptrdy )
+                                     !                                    ! wind stress trends
+                                     z2dx(:,:) = ( utau_b(ji,jj) + utau(ji,jj) ) / ( fse3u(:,:,1) * rau0 )
+                                     z2dy(:,:) = ( vtau_b(ji,jj) + vtau(ji,jj) ) / ( fse3v(:,:,1) * rau0 )
+                                     CALL iom_put( "utrd_tau", z2dx )
+                                     CALL iom_put( "vtrd_tau", z2dy )
+         CASE( jpdyn_trd_bfr )   ;   CALL iom_put( "utrd_bfr", ptrdx )    ! bottom friction term
+                                     CALL iom_put( "vtrd_bfr", ptrdy )
+         END SELECT
+         !
+      ENDIF
+      !
+      CALL wrk_dealloc( jpi, jpj     , z2dx, z2dy, ztswu, ztswv )
+      CALL wrk_dealloc( jpi, jpj, jpk, z3dx, z3dy )
+      !
+   END SUBROUTINE trd_3Diom
+
 #else
    !!----------------------------------------------------------------------
-   !!   Default case :                                         Empty module
-   !!----------------------------------------------------------------------
-   USE trdmod_oce      ! ocean variables trends
-   USE trdvor          ! ocean vorticity trends 
-   USE trdicp          ! ocean bassin integral constraints properties
-   USE trdmld          ! ocean active mixed layer tracers trends 
+   !!   Default case :           Empty module          No trend diagnostics
    !!----------------------------------------------------------------------
 CONTAINS
-   SUBROUTINE trd_mod(ptrd3dx, ptrd3dy, ktrd , ctype, kt )   ! Empty routine
-      REAL(wp) ::   ptrd3dx(:,:,:), ptrd3dy(:,:,:)
+   SUBROUTINE trd_mod( ptrdx, ptrdy, ktrd, ctype, kt )   ! Empty routine
+      REAL ::   ptrdx(:,:,:), ptrdy(:,:,:)
       INTEGER  ::   ktrd, kt                            
       CHARACTER(len=3) ::  ctype                  
-      WRITE(*,*) 'trd_3d: You should not have seen this print! error ?', ptrd3dx(1,1,1), ptrd3dy(1,1,1)
-      WRITE(*,*) ' "   ": You should not have seen this print! error ?', ktrd, ctype, kt
+      WRITE(*,*) 'trd_mod: You should not have seen this print! error ?',   &
+         &       ptrdx(1,1,1), ptrdy(1,1,1), ktrd, ctype, kt
    END SUBROUTINE trd_mod
 #endif
@@ -251,6 +354,9 @@
       !!----------------------------------------------------------------------
       USE in_out_manager          ! I/O manager
-      !!    
-      NAMELIST/namtrd/ nn_trd, nn_ctls, cn_trdrst_in, cn_trdrst_out, ln_trdmld_restart, rn_ucf, ln_trdmld_instant
+
+      NAMELIST/namtrd/ ln_3D_trd_d, ln_KE_trd, ln_vor_trd, ln_ML_trd_d,   &
+         &             ln_3D_trd_t, ln_PE_trd, ln_glo_trd, ln_ML_trd_t,   &
+         &             nn_trd , cn_trdrst_in , ln_trdmld_restart,         &
+         &             nn_ctls, cn_trdrst_out, ln_trdmld_instant, rn_ucf
       !!----------------------------------------------------------------------
 
@@ -264,11 +370,25 @@
             WRITE(numout,*) ' ~~~~~~~~~~~~~'
             WRITE(numout,*) '   Namelist namtrd : set trends parameters'
-            WRITE(numout,*) '      frequency of trends diagnostics   nn_trd             = ', nn_trd
-            WRITE(numout,*) '      control surface type              nn_ctls            = ', nn_ctls
-            WRITE(numout,*) '      restart for ML diagnostics        ln_trdmld_restart  = ', ln_trdmld_restart
-            WRITE(numout,*) '      instantaneous or mean ML T/S      ln_trdmld_instant  = ', ln_trdmld_instant
-            WRITE(numout,*) '      unit conversion factor            rn_ucf             = ', rn_ucf
+            WRITE(numout,*) '      U & V trends: 3D output                 ln_3D_trd_d        = ', ln_3D_trd_d
+            WRITE(numout,*) '      T & S trends: 3D output                 ln_3D_trd_t        = ', ln_3D_trd_t
+            WRITE(numout,*) '      Kinetic   Energy trends                 ln_KE_trd          = ', ln_KE_trd
+            WRITE(numout,*) '      Potential Energy trends                 ln_PE_trd          = ', ln_PE_trd
+            WRITE(numout,*) '      Barotropic vorticity trends             ln_vor_trd         = ', ln_vor_trd
+            WRITE(numout,*) '      check domain averaged dyn & tra trends  ln_glo_trd         = ', ln_glo_trd
+            WRITE(numout,*) '      U & V trends: Mixed Layer averaged      ln_ML_trd_d        = ', ln_3D_trd_d
+            WRITE(numout,*) '      T & S trends: Mixed Layer averaged      ln_ML_trd_t        = ', ln_3D_trd_t
+     !
+            WRITE(numout,*) '      frequency of trends diagnostics (glo)   nn_trd             = ', nn_trd
+            WRITE(numout,*) '      control surface type            (mld)   nn_ctls            = ', nn_ctls
+            WRITE(numout,*) '      restart for ML diagnostics              ln_trdmld_restart  = ', ln_trdmld_restart
+            WRITE(numout,*) '      instantaneous or mean ML T/S            ln_trdmld_instant  = ', ln_trdmld_instant
+            WRITE(numout,*) '      unit conversion factor                  rn_ucf             = ', rn_ucf
         ENDIF
       ENDIF
+      !
+      IF( ln_KE_trd .OR. ln_PE_trd .OR. ln_ML_trd_d )   &
+         CALL ctl_stop( 'KE, PE, aur ML on momentum are not yet coded we stop' )
+!!gm  : Potential BUG : 3D output only for vector invariant form!  add a ctl_stop or code the flux form case
+!!gm  : bug/pb for vertical advection of tracer in vvl case: add T.dt[eta] in the output... 
       !
       IF( lk_trddyn .OR. lk_trdtra )    CALL trd_icp_init       ! integral constraints trends
Index: /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/TRD/trdmod_oce.F90
===================================================================
--- /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/TRD/trdmod_oce.F90	(revision 3315)
+++ /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/TRD/trdmod_oce.F90	(revision 3316)
@@ -14,5 +14,13 @@
 
    !                                                     !!* Namelist namtrd:  diagnostics on dynamics/tracer trends *
-   INTEGER , PUBLIC  ::   nn_trd  = 10                    !: time step frequency dynamics and tracers trends
+   LOGICAL , PUBLIC  ::   ln_3D_trd_d   = .FALSE.         !: (T) 3D momentum              trends or (F) not
+   LOGICAL , PUBLIC  ::   ln_3D_trd_t   = .FALSE.         !: (T) 3D tracer               trends or (F) not
+   LOGICAL , PUBLIC  ::   ln_PE_trd     = .FALSE.         !: (T) 3D Potential Energy     trends or (F) not
+   LOGICAL , PUBLIC  ::   ln_KE_trd     = .FALSE.         !: (T) 3D Kinetic   Energy     trends or (F) not
+   LOGICAL , PUBLIC  ::   ln_vor_trd    = .FALSE.         !: (T) 3D barotropic vorticity trends or (F) not
+   LOGICAL , PUBLIC  ::   ln_glo_trd    = .FALSE.         !: (T) global domain averaged diag for T, T^2, KE, and PE
+   LOGICAL , PUBLIC  ::   ln_ML_trd_t   = .FALSE.         !: (T) 2D tracer   trends averaged over the mixed layer 
+   LOGICAL , PUBLIC  ::   ln_ML_trd_d   = .FALSE.         !: (T) 2D momentum trends averaged over the mixed layer 
+   INTEGER , PUBLIC  ::   nn_trd  = 10                    !: time step frequency for ln_glo_trd
    INTEGER , PUBLIC  ::   nn_ctls =  0                    !: control surface type for trends vertical integration
    REAL(wp), PUBLIC  ::   rn_ucf   = 1.                   !: unit conversion factor (for netCDF trends outputs)
@@ -43,32 +51,29 @@
    INTEGER, PUBLIC, PARAMETER ::   jptra_trd_zad =  3     !: z- vertical   advection
    INTEGER, PUBLIC, PARAMETER ::   jptra_trd_ldf =  4     !: lateral       diffusion
-   INTEGER, PUBLIC, PARAMETER ::   jptra_trd_zdf =  5     !: vertical diffusion (Kz)
-   INTEGER, PUBLIC, PARAMETER ::   jptra_trd_bbc =  6     !: Bottom Boundary Condition (geoth. flux) 
-   INTEGER, PUBLIC, PARAMETER ::   jptra_trd_bbl =  7     !: Bottom Boundary Layer (diffusive/convective)
-   INTEGER, PUBLIC, PARAMETER ::   jptra_trd_npc =  8     !: static instability mixing
-   INTEGER, PUBLIC, PARAMETER ::   jptra_trd_dmp =  9     !: damping
+   INTEGER, PUBLIC, PARAMETER ::   jptra_trd_zdf =  5     !: vertical      diffusion
+   INTEGER, PUBLIC, PARAMETER ::   jptra_trd_bbc =  6     !: Bottom Boundary Condition (geoth. heating) 
+   INTEGER, PUBLIC, PARAMETER ::   jptra_trd_bbl =  7     !: Bottom Boundary Layer (diffusive and/or advective)
+   INTEGER, PUBLIC, PARAMETER ::   jptra_trd_npc =  8     !: non-penetrative convection treatment
+   INTEGER, PUBLIC, PARAMETER ::   jptra_trd_dmp =  9     !: internal restoring (damping)
    INTEGER, PUBLIC, PARAMETER ::   jptra_trd_qsr = 10     !: penetrative solar radiation
-   INTEGER, PUBLIC, PARAMETER ::   jptra_trd_nsr = 11     !: non solar radiation
-   INTEGER, PUBLIC, PARAMETER ::   jptra_trd_atf = 12     !: Asselin correction
-#if defined key_top
-   !                                                     !!!* Passive tracers trends indexes
+   INTEGER, PUBLIC, PARAMETER ::   jptra_trd_nsr = 11     !: non         solar radiation
+   INTEGER, PUBLIC, PARAMETER ::   jptra_trd_atf = 12     !: Asselin time filter
+   !
+   !                                                     !!!* Passive tracers trends indexes (use if "key_top" defined)
    INTEGER, PUBLIC, PARAMETER ::   jptra_trd_sms  = 13    !: sources m. sinks
    INTEGER, PUBLIC, PARAMETER ::   jptra_trd_radn = 14    !: corr. trn<0 in trcrad
    INTEGER, PUBLIC, PARAMETER ::   jptra_trd_radb = 15    !: corr. trb<0 in trcrad (like atf)
-#endif
-   
+   !
    !                                                     !!!* Momentum trends indexes
    INTEGER, PUBLIC, PARAMETER ::   jpdyn_trd_hpg =  1     !: hydrostatic pressure gradient 
-   INTEGER, PUBLIC, PARAMETER ::   jpdyn_trd_keg =  2     !: kinetic energy gradient
-   INTEGER, PUBLIC, PARAMETER ::   jpdyn_trd_rvo =  3     !: relative vorticity
-   INTEGER, PUBLIC, PARAMETER ::   jpdyn_trd_pvo =  4     !: planetary vorticity
-   INTEGER, PUBLIC, PARAMETER ::   jpdyn_trd_ldf =  5     !: lateral diffusion
-   INTEGER, PUBLIC, PARAMETER ::   jpdyn_trd_had =  6     !: horizontal advection
-   INTEGER, PUBLIC, PARAMETER ::   jpdyn_trd_zad =  7     !: vertical advection
-   INTEGER, PUBLIC, PARAMETER ::   jpdyn_trd_zdf =  8     !: vertical diffusion
-   INTEGER, PUBLIC, PARAMETER ::   jpdyn_trd_spg =  9     !: surface pressure gradient
-   INTEGER, PUBLIC, PARAMETER ::   jpdyn_trd_dat = 10     !: damping term
-   INTEGER, PUBLIC, PARAMETER ::   jpdyn_trd_swf = 11     !: surface wind forcing
-   INTEGER, PUBLIC, PARAMETER ::   jpdyn_trd_bfr = 12     !: bottom friction 
+   INTEGER, PUBLIC, PARAMETER ::   jpdyn_trd_spg =  2     !: surface     pressure gradient
+   INTEGER, PUBLIC, PARAMETER ::   jpdyn_trd_keg =  3     !: kinetic energy gradient  or horizontal advection
+   INTEGER, PUBLIC, PARAMETER ::   jpdyn_trd_rvo =  4     !: relative  vorticity      or metric term
+   INTEGER, PUBLIC, PARAMETER ::   jpdyn_trd_pvo =  5     !: planetary vorticity
+   INTEGER, PUBLIC, PARAMETER ::   jpdyn_trd_zad =  6     !: vertical advection
+   INTEGER, PUBLIC, PARAMETER ::   jpdyn_trd_ldf =  7     !: horizontal diffusion   
+   INTEGER, PUBLIC, PARAMETER ::   jpdyn_trd_zdf =  8     !: vertical   diffusion
+   INTEGER, PUBLIC, PARAMETER ::   jpdyn_trd_swf =  9     !: surface stress
+   INTEGER, PUBLIC, PARAMETER ::   jpdyn_trd_bfr = 10     !: bottom  stress 
 
    !!----------------------------------------------------------------------
Index: /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/TRD/trdtra.F90
===================================================================
--- /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/TRD/trdtra.F90	(revision 3315)
+++ /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/TRD/trdtra.F90	(revision 3316)
@@ -2,29 +2,28 @@
    !!======================================================================
    !!                       ***  MODULE  trdtra  ***
-   !! Ocean diagnostics:  ocean tracers trends
+   !! Ocean diagnostics:  ocean tracers trends pre-processing
    !!=====================================================================
-   !! History :  1.0  !  2004-08  (C. Talandier) Original code
-   !!            2.0  !  2005-04  (C. Deltel)    Add Asselin trend in the ML budget
-   !!            3.3  !  2010-06  (C. Ethe) merge TRA-TRC 
+   !! History :  3.3  !  2010-06  (C. Ethe) creation for the TRA/TRC merge
+   !!            3.5  !  2012-02  (G. Madec) update the comments 
    !!----------------------------------------------------------------------
 #if  defined key_trdtra || defined key_trdmld || defined key_trdmld_trc 
    !!----------------------------------------------------------------------
-   !!   trd_tra      : Call the trend to be computed
-   !!----------------------------------------------------------------------
-   USE dom_oce          ! ocean domain 
-   USE trdmod_oce       ! ocean active mixed layer tracers trends 
-   USE trdmod           ! ocean active mixed layer tracers trends 
-   USE trdmod_trc       ! ocean passive mixed layer tracers trends 
-   USE in_out_manager   ! I/O manager
-   USE lib_mpp          ! MPP library
-   USE wrk_nemo        ! Memory allocation
-
+   !!   trd_tra       : pre-process the tracer trends and calll trd_mod(_trc)
+   !!   trd_tra_adv   : transform a div(U.T) trend into a U.grad(T) trend
+   !!----------------------------------------------------------------------
+   USE dom_oce        ! ocean domain 
+   USE trdmod_oce     ! ocean active mixed layer tracers trends 
+   USE trdmod         ! ocean active mixed layer tracers trends 
+   USE trdmod_trc     ! ocean passive mixed layer tracers trends 
+   USE in_out_manager ! I/O manager
+   USE lib_mpp        ! MPP library
+   USE wrk_nemo       ! Memory allocation
 
    IMPLICIT NONE
    PRIVATE
 
-   PUBLIC   trd_tra          ! called by all  traXX modules
+   PUBLIC   trd_tra   ! called by all tra_... modules
  
-   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: trdtx, trdty, trdt  !:
+   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   trdtx, trdty, trdt   ! use to store the temperature trends
 
    !! * Substitutions
@@ -32,5 +31,5 @@
 #  include "vectopt_loop_substitute.h90"
    !!----------------------------------------------------------------------
-   !! NEMO/OPA 4.0 , NEMO Consortium (2011)
+   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
    !! $Id$
    !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
@@ -53,86 +52,79 @@
       !!                  ***  ROUTINE trd_tra  ***
       !! 
-      !! ** Purpose : Dispatch all trends computation, e.g. vorticity, mld or 
-      !!              integral constraints
+      !! ** Purpose : pre-process tracer trends
       !!
-      !! ** Method/usage : For the mixed-layer trend, the control surface can be either
-      !!       a mixed layer depth (time varying) or a fixed surface (jk level or bowl). 
-      !!      Choose control surface with nn_ctls in namelist NAMTRD :
-      !!        nn_ctls = 0  : use mixed layer with density criterion 
-      !!        nn_ctls = 1  : read index from file 'ctlsurf_idx'
-      !!        nn_ctls > 1  : use fixed level surface jk = nn_ctls
-      !!----------------------------------------------------------------------
-      !
-      INTEGER                         , INTENT(in)           ::  kt      ! time step
-      CHARACTER(len=3)                , INTENT(in)           ::  ctype   ! tracers trends type 'TRA'/'TRC'
-      INTEGER                         , INTENT(in)           ::  ktra    ! tracer index
-      INTEGER                         , INTENT(in)           ::  ktrd    ! tracer trend index
-      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in)           ::  ptrd    ! tracer trend  or flux
-      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in), OPTIONAL ::  pun     ! velocity 
-      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in), OPTIONAL ::  ptra    ! Tracer variablea
+      !! ** Method  : - mask the trend
+      !!              - advection (ptra present) converte the incoming flux (U.T) 
+      !!              into trend (U.T => -U.grat(T)=div(U.T)-T.div(U)) through a 
+      !!              call to trd_tra_adv
+      !!              - 'TRA' case : regroup T & S trends
+      !!              - send the trends to trd_mod(_trc) for further processing
+      !!----------------------------------------------------------------------
+      INTEGER                         , INTENT(in)           ::   kt      ! time step
+      CHARACTER(len=3)                , INTENT(in)           ::   ctype   ! tracers trends type 'TRA'/'TRC'
+      INTEGER                         , INTENT(in)           ::   ktra    ! tracer index
+      INTEGER                         , INTENT(in)           ::   ktrd    ! tracer trend index
+      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in)           ::   ptrd    ! tracer trend  or flux
+      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in), OPTIONAL ::   pun     ! now velocity 
+      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in), OPTIONAL ::   ptra    ! now tracer variable
       !
       REAL(wp), POINTER, DIMENSION(:,:,:)  ::  ztrds
       !!----------------------------------------------------------------------
-
+      !
       CALL wrk_alloc( jpi, jpj, jpk, ztrds )
-
-      IF( .NOT. ALLOCATED( trdtx ) ) THEN       ! allocate trdtra arrays
+      !
+      IF( .NOT. ALLOCATED( trdtx ) ) THEN      ! allocate trdtra arrays
          IF( trd_tra_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'trd_tra : unable to allocate arrays' )
       ENDIF
       
-      ! Control of optional arguments
-      IF( ctype == 'TRA' .AND. ktra == jp_tem ) THEN 
-         IF( PRESENT( ptra ) ) THEN    
-            SELECT CASE( ktrd )            ! shift depending on the direction
-            CASE( jptra_trd_xad )  ;  CALL trd_tra_adv( ptrd, pun, ptra, 'X', trdtx ) 
-            CASE( jptra_trd_yad )  ;  CALL trd_tra_adv( ptrd, pun, ptra, 'Y', trdty ) 
-            CASE( jptra_trd_zad )  ;  CALL trd_tra_adv( ptrd, pun, ptra, 'Z', trdt  ) 
+      IF( ctype == 'TRA' .AND. ktra == jp_tem ) THEN   !==  Temperature trend  ==!
+         !
+         IF( PRESENT( ptra ) ) THEN                       ! advection: transform flux into trend
+            SELECT CASE( ktrd )     
+            CASE( jptra_trd_xad )   ;   CALL trd_tra_adv( ptrd, pun, ptra, 'X', trdtx ) 
+            CASE( jptra_trd_yad )   ;   CALL trd_tra_adv( ptrd, pun, ptra, 'Y', trdty ) 
+            CASE( jptra_trd_zad )   ;   CALL trd_tra_adv( ptrd, pun, ptra, 'Z', trdt  ) 
             END SELECT
-         ELSE
-            trdt(:,:,:) = ptrd(:,:,:)
-            IF( ktrd == jptra_trd_bbc .OR. ktrd == jptra_trd_qsr ) THEN
-               ztrds(:,:,:) = 0.
-               CALL trd_mod( trdt, ztrds, ktrd, ctype, kt )
-            END IF
+         ELSE                                             ! other trends:
+            trdt(:,:,:) = ptrd(:,:,:) * tmask(:,:,:)                      ! mask & store
+            IF( ktrd == jptra_trd_bbc .OR. ktrd == jptra_trd_qsr ) THEN   ! qsr, bbc: on temperature only
+               ztrds(:,:,:) = 0._wp
+               CALL trd_mod( trdt, ztrds, ktrd, ctype, kt )                  ! send to trd_mod
+            ENDIF
+         ENDIF
+         !
+      ENDIF
+
+      IF( ctype == 'TRA' .AND. ktra == jp_sal ) THEN      !==  Salinity trends  ==!
+         !
+         IF( PRESENT( ptra ) ) THEN      ! advection: transform the advective flux into a trend
+            SELECT CASE( ktrd )          !            and send T & S trends to trd_mod
+            CASE( jptra_trd_xad )   ;   CALL trd_tra_adv( ptrd , pun  , ptra, 'X'  , ztrds ) 
+                                        CALL trd_mod    ( trdtx, ztrds, ktrd, ctype, kt    )
+            CASE( jptra_trd_yad )   ;   CALL trd_tra_adv( ptrd , pun  , ptra, 'Y'  , ztrds ) 
+                                    ;   CALL trd_mod    ( trdty, ztrds, ktrd, ctype, kt    )
+            CASE( jptra_trd_zad )   ;   CALL trd_tra_adv( ptrd , pun  , ptra, 'Z'  , ztrds ) 
+                                        CALL trd_mod    ( trdt , ztrds, ktrd, ctype, kt    )
+            END SELECT
+         ELSE                            ! other trends: mask and send T & S trends to trd_mod
+            ztrds(:,:,:) = ptrd(:,:,:) * tmask(:,:,:)
+            CALL trd_mod( trdt, ztrds, ktrd, ctype, kt )  
+         ENDIF
+         !
+      ENDIF
+
+      IF( ctype == 'TRC' ) THEN                           !==  passive tracer trend  ==!
+         !
+         IF( PRESENT( ptra ) ) THEN                          ! advection: transform flux into a trend
+            SELECT CASE( ktrd )
+            CASE( jptra_trd_xad )   ;   CALL trd_tra_adv( ptrd , pun , ptra, 'X', ztrds ) 
+            CASE( jptra_trd_yad )   ;   CALL trd_tra_adv( ptrd , pun , ptra, 'Y', ztrds ) 
+            CASE( jptra_trd_zad )   ;   CALL trd_tra_adv( ptrd , pun , ptra, 'Z', ztrds ) 
+            END SELECT
+         ELSE                                                ! other trends: mask
+            ztrds(:,:,:) = ptrd(:,:,:) * tmask(:,:,:)
          END IF
-      END IF
-
-      IF( ctype == 'TRA' .AND. ktra == jp_sal ) THEN 
-         IF( PRESENT( ptra ) ) THEN    
-            SELECT CASE( ktrd )            ! shift depending on the direction
-            CASE( jptra_trd_xad )  
-                                CALL trd_tra_adv( ptrd, pun, ptra, 'X', ztrds ) 
-                                CALL trd_mod( trdtx, ztrds, ktrd, ctype, kt   )
-            CASE( jptra_trd_yad )  
-                                CALL trd_tra_adv( ptrd, pun, ptra, 'Y', ztrds ) 
-                                CALL trd_mod( trdty, ztrds, ktrd, ctype, kt   )
-            CASE( jptra_trd_zad )    
-                                CALL trd_tra_adv( ptrd, pun, ptra, 'Z', ztrds ) 
-                                CALL trd_mod( trdt , ztrds, ktrd, ctype, kt   )
-            END SELECT
-         ELSE
-            ztrds(:,:,:) = ptrd(:,:,:)
-            CALL trd_mod( trdt, ztrds, ktrd, ctype, kt )  
-         END IF
-      END IF
-
-      IF( ctype == 'TRC' ) THEN
-         !
-         IF( PRESENT( ptra ) ) THEN  
-            SELECT CASE( ktrd )            ! shift depending on the direction
-            CASE( jptra_trd_xad )  
-                                CALL trd_tra_adv( ptrd, pun, ptra, 'X', ztrds ) 
-                                CALL trd_mod_trc( ztrds, ktra, ktrd, kt       )
-            CASE( jptra_trd_yad )  
-                                CALL trd_tra_adv( ptrd, pun, ptra, 'Y', ztrds ) 
-                                CALL trd_mod_trc( ztrds, ktra, ktrd, kt       )
-            CASE( jptra_trd_zad )    
-                                CALL trd_tra_adv( ptrd, pun, ptra, 'Z', ztrds ) 
-                                CALL trd_mod_trc( ztrds, ktra, ktrd, kt       )
-            END SELECT
-         ELSE
-            ztrds(:,:,:) = ptrd(:,:,:)
-            CALL trd_mod_trc( ztrds, ktra, ktrd, kt       )  
-         END IF
+         !                                 
+         CALL trd_mod_trc( ztrds, ktra, ktrd, kt )           ! send trend to trd_mod_trc
          !
       ENDIF
@@ -147,39 +139,39 @@
       !!                  ***  ROUTINE trd_tra_adv  ***
       !! 
-      !! ** Purpose :   transformed the i-, j- or k-advective flux into thes
-      !!              i-, j- or k-advective trends, resp.
-      !! ** Method  :   i-advective trends = -un. di-1[T] = -( di-1[fi] - tn di-1[un] )
-      !!                k-advective trends = -un. di-1[T] = -( dj-1[fi] - tn dj-1[un] )
-      !!                k-advective trends = -un. di+1[T] = -( dk+1[fi] - tn dk+1[un] )
-      !!----------------------------------------------------------------------
-      REAL(wp)        , INTENT(in ), DIMENSION(jpi,jpj,jpk) ::   pf      ! advective flux in one direction
-      REAL(wp)        , INTENT(in ), DIMENSION(jpi,jpj,jpk) ::   pun     ! now velocity  in one direction
-      REAL(wp)        , INTENT(in ), DIMENSION(jpi,jpj,jpk) ::   ptn     ! now or before tracer 
-      CHARACTER(len=1), INTENT(in )                         ::   cdir    ! X/Y/Z direction
-      REAL(wp)        , INTENT(out), DIMENSION(jpi,jpj,jpk) ::   ptrd    ! advective trend in one direction
+      !! ** Purpose :   transformed a advective flux into a masked advective trends
+      !!
+      !! ** Method  :   use the following transformation: -div(U.T) = - U grad(T) + T.div(U)
+      !!       i-advective trends = -un. di-1[T] = -( di-1[fi] - tn di-1[un] )
+      !!       j-advective trends = -un. di-1[T] = -( dj-1[fi] - tn dj-1[un] )
+      !!       k-advective trends = -un. di+1[T] = -( dk+1[fi] - tn dk+1[un] )
+      !!                where fi is the incoming advective flux.
+      !!----------------------------------------------------------------------
+      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in   ) ::   pf      ! advective flux in one direction
+      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in   ) ::   pun     ! now velocity   in one direction
+      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in   ) ::   ptn     ! now or before tracer 
+      CHARACTER(len=1)                , INTENT(in   ) ::   cdir    ! X/Y/Z direction
+      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(  out) ::   ptrd    ! advective trend in one direction
       !
       INTEGER  ::   ji, jj, jk   ! dummy loop indices
-      INTEGER  ::   ii, ij, ik   ! index shift function of the direction
-      REAL(wp) ::   zbtr         ! local scalar
-      !!----------------------------------------------------------------------
-
-      SELECT CASE( cdir )            ! shift depending on the direction
-      CASE( 'X' )   ;   ii = 1   ; ij = 0   ;   ik = 0      ! i-advective trend
-      CASE( 'Y' )   ;   ii = 0   ; ij = 1   ;   ik = 0      ! j-advective trend
-      CASE( 'Z' )   ;   ii = 0   ; ij = 0   ;   ik =-1      ! k-advective trend
+      INTEGER  ::   ii, ij, ik   ! index shift as function of the direction
+      !!----------------------------------------------------------------------
+      !
+      SELECT CASE( cdir )      ! shift depending on the direction
+      CASE( 'X' )   ;   ii = 1   ;   ij = 0   ;   ik = 0      ! i-trend
+      CASE( 'Y' )   ;   ii = 0   ;   ij = 1   ;   ik = 0      ! j-trend
+      CASE( 'Z' )   ;   ii = 0   ;   ij = 0   ;   ik =-1      ! k-trend
       END SELECT
-
-      !                              ! set to zero uncomputed values
-      ptrd(jpi,:,:) = 0.e0   ;   ptrd(1,:,:) = 0.e0
-      ptrd(:,jpj,:) = 0.e0   ;   ptrd(:,1,:) = 0.e0
-      ptrd(:,:,jpk) = 0.e0
-      !
-      !
-      DO jk = 1, jpkm1
+      !
+      !                        ! set to zero uncomputed values
+      ptrd(jpi,:,:) = 0._wp   ;   ptrd(1,:,:) = 0._wp
+      ptrd(:,jpj,:) = 0._wp   ;   ptrd(:,1,:) = 0._wp
+      ptrd(:,:,jpk) = 0._wp
+      !
+      DO jk = 1, jpkm1         ! advective trend
          DO jj = 2, jpjm1
             DO ji = fs_2, fs_jpim1   ! vector opt.
-               zbtr    = 1.e0/ ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
-               ptrd(ji,jj,jk) = - zbtr * (      pf (ji,jj,jk) - pf (ji-ii,jj-ij,jk-ik)                    &
-                 &                          - ( pun(ji,jj,jk) - pun(ji-ii,jj-ij,jk-ik) ) * ptn(ji,jj,jk)  )
+               ptrd(ji,jj,jk) = - (     pf (ji,jj,jk) - pf (ji-ii,jj-ij,jk-ik)                        &
+                 &                  - ( pun(ji,jj,jk) - pun(ji-ii,jj-ij,jk-ik) ) * ptn(ji,jj,jk)  )   &
+                 &              / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )  * tmask(ji,jj,jk)
             END DO
          END DO
@@ -188,5 +180,5 @@
    END SUBROUTINE trd_tra_adv
 
-#   else
+#else
    !!----------------------------------------------------------------------
    !!   Default case :          Dummy module           No trend diagnostics
@@ -196,15 +188,13 @@
    SUBROUTINE trd_tra( kt, ctype, ktra, ktrd, ptrd, pu, ptra )
       !!----------------------------------------------------------------------
-      INTEGER                         , INTENT(in)           ::  kt      ! time step
-      CHARACTER(len=3)                , INTENT(in)           ::  ctype   ! tracers trends type 'TRA'/'TRC'
-      INTEGER                         , INTENT(in)           ::  ktra    ! tracer index
-      INTEGER                         , INTENT(in)           ::  ktrd    ! tracer trend index
-      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in)           ::  ptrd    ! tracer trend 
-      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in), OPTIONAL ::  pu      ! velocity 
-      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in), OPTIONAL ::  ptra    ! Tracer variable 
-      WRITE(*,*) 'trd_3d: You should not have seen this print! error ?', ptrd(1,1,1), ptra(1,1,1), pu(1,1,1),   &
-         &                                                               ktrd, ktra, ctype, kt
+      CHARACTER(len=3)                , INTENT(in)           ::  ctype   
+      INTEGER                         , INTENT(in)           ::  kt, ktra, ktrd
+      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in)           ::  ptrd    
+      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in), OPTIONAL ::  pu, ptra    ! Tracer variable 
+      WRITE(*,*) 'trd_tra: You should not have seen this print! error ?',   &
+         &   ptrd(1,1,1), ptra(1,1,1), pu(1,1,1), ktrd, ktra, ctype, kt
    END SUBROUTINE trd_tra
-#   endif
+#endif
+
    !!======================================================================
 END MODULE trdtra
Index: /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/TRD/trdvor.F90
===================================================================
--- /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/TRD/trdvor.F90	(revision 3315)
+++ /branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/TRD/trdvor.F90	(revision 3316)
@@ -4,6 +4,7 @@
    !! Ocean diagnostics:  momentum trends
    !!=====================================================================
-   !! History :  1.0  !  04-2006  (L. Brunier, A-M. Treguier) Original code 
-   !!            2.0  !  04-2008  (C. Talandier) New trends organization
+   !! History :  1.0  !  2006-01  (L. Brunier, A-M. Treguier) Original code 
+   !!            2.0  !  2008-04  (C. Talandier) New trends organization
+   !!            3.5  !  2012-02  (G. Madec) regroup beta.V computation with pvo trend
    !!----------------------------------------------------------------------
 #if defined key_trdvor   ||   defined key_esopa
@@ -109,5 +110,4 @@
       !!      trends output in netCDF format using ioipsl
       !!----------------------------------------------------------------------
-      !
       INTEGER                     , INTENT(in   ) ::   ktrd       ! ocean trend index
       REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) ::   putrdvor   ! u vorticity trend 
@@ -229,15 +229,23 @@
       END DO
 
-      ! Save Beta.V term to avoid average before Curl
-      ! Beta.V : intergration, no average
-      IF( ktrd == jpvor_bev ) THEN 
+      ! Planetary vorticity: 2nd computation (Beta.V term) store the vertical sum
+      ! as Beta.V term need intergration, not average
+      IF( ktrd == jpvor_pvo ) THEN 
          zubet(:,:) = zudpvor(:,:)
          zvbet(:,:) = zvdpvor(:,:)
-      ENDIF
-
-      ! Average except for Beta.V
+         DO ji = 1, jpim1
+            DO jj = 1, jpjm1
+               vortrd(ji,jj,jpvor_bev) = (    zvbet(ji+1,jj) - zvbet(ji,jj)     &
+                  &                       - ( zubet(ji,jj+1) - zubet(ji,jj) ) ) / ( e1f(ji,jj) * e2f(ji,jj) )
+            END DO
+         END DO
+         ! Average of the Curl and Surface mask
+         vortrd(:,:,jpvor_bev) = vortrd(:,:,jpvor_bev) * hur(:,:) * fmask(:,:,1)
+      ENDIF
+      !
+      ! Average 
       zudpvor(:,:) = zudpvor(:,:) * hur(:,:)
       zvdpvor(:,:) = zvdpvor(:,:) * hvr(:,:)
-   
+      !
       ! Curl
       DO ji=1,jpim1
@@ -247,24 +255,6 @@
          END DO
       END DO
-
       ! Surface mask
       vortrd(:,:,ktrd) = vortrd(:,:,ktrd) * fmask(:,:,1)
-
-      ! Special treatement for the Beta.V term
-      ! Compute the Curl of the Beta.V term which is not averaged
-      IF( ktrd == jpvor_bev ) THEN
-         DO ji=1,jpim1
-            DO jj=1,jpjm1
-               vortrd(ji,jj,jpvor_bev) = (    zvbet(ji+1,jj) - zvbet(ji,jj)     &
-                  &                       - ( zubet(ji,jj+1) - zubet(ji,jj) ) ) / ( e1f(ji,jj) * e2f(ji,jj) )
-            END DO
-         END DO
-
-         ! Average on the Curl
-         vortrd(:,:,jpvor_bev) = vortrd(:,:,jpvor_bev) * hur(:,:)
-
-         ! Surface mask
-         vortrd(:,:,jpvor_bev) = vortrd(:,:,jpvor_bev) * fmask(:,:,1)
-      ENDIF
    
       IF( ndebug /= 0 ) THEN
@@ -330,6 +320,6 @@
 
       ! Curl
-      DO ji=1,jpim1
-         DO jj=1,jpjm1
+      DO ji = 1, jpim1
+         DO jj = 1, jpjm1
             vor_avr(ji,jj) = (  ( zvn(ji+1,jj) - zvn(ji,jj) )    &
                &              - ( zun(ji,jj+1) - zun(ji,jj) ) ) / ( e1f(ji,jj) * e2f(ji,jj) ) * fmask(ji,jj,1)
