Index: trunk/NEMO/TOP_SRC/TRP/trcadv_muscl2.F90
===================================================================
--- trunk/NEMO/TOP_SRC/TRP/trcadv_muscl2.F90	(revision 186)
+++ trunk/NEMO/TOP_SRC/TRP/trcadv_muscl2.F90	(revision 202)
@@ -10,8 +10,7 @@
    !!----------------------------------------------------------------------
    !! * Modules used
-   USE oce_trc         ! ocean dynamics and active tracers
-   USE trc             ! ocean passive tracers
-   USE trcbbl          ! tracers: bottom boudary layer
-   USE lbclnk
+   USE oce_trc         ! ocean dynamics and active tracers variables
+   USE trc             ! ocean passive tracers variables
+   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
 
    IMPLICIT NONE
@@ -125,18 +124,9 @@
          DO jk = 1, jpkm1
             DO jj = 2, jpj
-               DO ji = fs_2, fs_jpim1   ! vector opt.
-                  z0u = zt1(ji,jj,jk) * zt1(ji-1,jj,jk)
-                  IF( z0u > 0. ) THEN
-                     ztp1(ji,jj,jk) = 0.5 * ( zt1(ji,jj,jk)+zt1(ji-1,jj,jk) )
-                  ELSE
-                     ztp1(ji,jj,jk) = 0.e0
-                  ENDIF
-
-                  z0v = zt2(ji,jj,jk) * zt2(ji,jj-1,jk)
-                  IF( z0v > 0. ) THEN
-                     ztp2(ji,jj,jk) = 0.5 * ( zt2(ji,jj,jk)+zt2(ji,jj-1,jk) )
-                  ELSE
-                     ztp2(ji,jj,jk) = 0.e0
-                  ENDIF
+               DO ji = fs_2, jpi   ! vector opt.
+                  ztp1(ji,jj,jk) =                    ( zt1(ji,jj,jk) + zt1(ji-1,jj  ,jk) )   &
+                     &           * ( 0.25 + SIGN( 0.25, zt1(ji,jj,jk) * zt1(ji-1,jj  ,jk) ) )
+                  ztp2(ji,jj,jk) =                    ( zt2(ji,jj,jk) + zt2(ji  ,jj-1,jk) )   &
+                     &           * ( 0.25 + SIGN( 0.25, zt2(ji,jj,jk) * zt2(ji  ,jj-1,jk) ) )
                END DO
             END DO
@@ -287,6 +277,9 @@
          END DO
 
-
-
+         IF(l_ctl) THEN         ! print mean trends (used for debugging)
+            ztra = SUM( tra(2:nictl,2:njctl,1:jpkm1,jn) * tmask(2:nictl,2:njctl,1:jpkm1) )
+            WRITE(numout,*) ' trc/had  - ',ctrcnm(jn),' : ', ztra-tra_ctl(jn), ' muscl' 
+            tra_ctl(jn) = ztra 
+         ENDIF
 
          ! II. Vertical advective fluxes
@@ -306,10 +299,6 @@
             DO jj = 1, jpj
                DO ji = 1, jpi
-                  z0w = zt1(ji,jj,jk) * zt1(ji,jj,jk+1) 
-                  IF( z0w > 0. ) THEN
-                     ztp1(ji,jj,jk) = 0.5 * ( zt1(ji,jj,jk) + zt1(ji,jj,jk+1) )
-                  ELSE
-                     ztp1(ji,jj,jk) = 0.e0
-                  ENDIF
+                  ztp1(ji,jj,jk) =                    ( zt1(ji,jj,jk) + zt1(ji,jj,jk+1) )   &
+                     &           * ( 0.25 + SIGN( 0.25, zt1(ji,jj,jk) * zt1(ji,jj,jk+1) ) )
                END DO
             END DO
@@ -360,11 +349,9 @@
 
          ! surface values
-#if defined key_dynspg_fsc
-         ! free surface-constant volume
-         zt1(:,:, 1 ) = zwn(:,:,1) * trb(:,:,1,jn)
-#else
-         ! rigid lid : flux set to zero
-         zt1(:,:, 1 ) = 0.e0
-#endif
+         IF( lk_dynspg_fsc .OR. lk_dynspg_fsc_tsk ) THEN        ! free surface-constant volume
+            zt1(:,:, 1 ) = zwn(:,:,1) * trb(:,:,1,jn)
+         ELSE                                                   ! rigid lid : flux set to zero
+            zt1(:,:, 1 ) = 0.e0
+         ENDIF
 
          ! bottom values
@@ -389,11 +376,11 @@
          END DO
 
-         IF( l_ctl .AND. lwp ) THEN         ! print mean trends (used for debugging)
-            ztra = SUM( tra(2:jpim1,2:jpjm1,1:jpkm1,jn) * tmask(2:jpim1,2:jpjm1,1:jpkm1) )
+         IF(l_ctl) THEN         ! print mean trends (used for debugging)
+            ztra = SUM( tra(2:nictl,2:njctl,1:jpkm1,jn) * tmask(2:nictl,2:njctl,1:jpkm1) )
             WRITE(numout,*) ' trc/zad  - ',ctrcnm(jn),' : ', ztra-tra_ctl(jn), ' muscl2' 
             tra_ctl(jn) = ztra 
          ENDIF
 
-      ENDDO
+      END DO
 
    END SUBROUTINE trc_adv_muscl2
