Index: /branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/lib_fortran.F90
===================================================================
--- /branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/lib_fortran.F90	(revision 2306)
+++ /branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/lib_fortran.F90	(revision 2307)
@@ -4,15 +4,18 @@
    !! Fortran utilities:  includes some low levels fortran functionality
    !!======================================================================
-   !! History : 3.2  !  2010-05  Michael Dunphy, Rachid BENSHILA Original code
-   !!----------------------------------------------------------------------
-   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
-   !! $Id$ 
-   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
-   !!----------------------------------------------------------------------
-   USE par_oce 
-   USE par_kind
-   USE lib_mpp         ! distributed memory computing
-   USE dom_oce 
-   USE in_out_manager
+   !! History :  3.2  !  2010-05  (M. Dunphy, R. Benshila)  Original code
+   !!----------------------------------------------------------------------
+
+   !!----------------------------------------------------------------------
+   !!   glob_sum    : generic interface for global masked summation over 
+   !!                 the interior domain for 1 or 2 2D or 3D arrays
+   !!                 it works only for T points   
+   !!   SIGN        : generic interface for SIGN to overwrite f95 behaviour
+   !!                 of intrinsinc sign function
+   !!----------------------------------------------------------------------
+   USE par_oce          ! Ocean parameter
+   USE lib_mpp          ! distributed memory computing
+   USE dom_oce          ! ocean domain
+   USE in_out_manager   ! I/O manager
 
    IMPLICIT NONE
@@ -25,122 +28,222 @@
 
    INTERFACE glob_sum
-#if defined key_mpp_rep
-      MODULE PROCEDURE mpp_sum_cmpx
-#else
       MODULE PROCEDURE glob_sum_2d, glob_sum_3d,glob_sum_2d_a, glob_sum_3d_a 
-#endif
    END INTERFACE
 
 #if defined key_nosignedzeo   
    INTERFACE SIGN
-      MODULE PROCEDURE SIGN_SCALAR, SIGN_ARRAY_1D, SIGN_ARRAY_2D, SIGN_ARRAY_3D, &
-                       SIGN_ARRAY_1D_A, SIGN_ARRAY_2D_A, SIGN_ARRAY_3D_A,  & 
-                       SIGN_ARRAY_1D_B, SIGN_ARRAY_2D_B, SIGN_ARRAY_3D_B 
+      MODULE PROCEDURE SIGN_SCALAR, SIGN_ARRAY_1D, SIGN_ARRAY_2D, SIGN_ARRAY_3D,   &
+         &             SIGN_ARRAY_1D_A, SIGN_ARRAY_2D_A, SIGN_ARRAY_3D_A,          & 
+         &             SIGN_ARRAY_1D_B, SIGN_ARRAY_2D_B, SIGN_ARRAY_3D_B 
    END INTERFACE
 #endif
 
+   !!----------------------------------------------------------------------
+   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
+   !! $Id$ 
+   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
+   !!----------------------------------------------------------------------
 CONTAINS 
 
-   FUNCTION glob_sum_2d( ptab )
+#if ! defined key_mpp_rep
+   FUNCTION glob_sum_2d( ptab ) RESULT( glob_sum )
       !!-----------------------------------------------------------------------
       !!                  ***  FUNCTION  glob_sum_2D  ***
       !!
-      !! ** Purpose : perform a sum on the global domain of a 2D array
-      !!-----------------------------------------------------------------------
-      REAL(wp), DIMENSION(:,:),INTENT(in) :: ptab
-      REAL(wp) :: glob_sum_2d
-      !!-----------------------------------------------------------------------
-
-      glob_sum_2d = SUM( ptab(:,:)*tmask_i(:,:) )
-      IF( lk_mpp )   CALL mpp_sum( glob_sum_2d )
-       
+      !! ** Purpose : perform a masked sum on the inner global domain of a 2D array
+      !!-----------------------------------------------------------------------
+      REAL(wp), INTENT(in), DIMENSION(:,:) ::   ptab       ! input 2D array
+      REAL(wp)                             ::   glob_sum   ! global masked sum
+      !!-----------------------------------------------------------------------
+      !
+      glob_sum = SUM( ptab(:,:)*tmask_i(:,:) )
+      IF( lk_mpp )   CALL mpp_sum( glob_sum )
+      !
    END FUNCTION glob_sum_2d
- 
-   FUNCTION glob_sum_3d( ptab )
+   
+   
+   FUNCTION glob_sum_3d( ptab ) RESULT( glob_sum )
       !!-----------------------------------------------------------------------
       !!                  ***  FUNCTION  glob_sum_3D  ***
       !!
-      !! ** Purpose : perform a sum on the global domain of a 3D array
-      !!-----------------------------------------------------------------------
-      REAL(wp), DIMENSION(:,:,:) :: ptab
-      REAL(wp) :: glob_sum_3d
-      !
+      !! ** Purpose : perform a masked sum on the inner global domain of a 3D array
+      !!-----------------------------------------------------------------------
+      REAL(wp), INTENT(in), DIMENSION(:,:,:) ::   ptab       ! input 3D array
+      REAL(wp)                               ::   glob_sum   ! global masked sum
+      !!
       INTEGER :: jk
       !!-----------------------------------------------------------------------
-       
-      GLOB_SUM_3D = 0.e0
+      !
+      glob_sum = 0.e0
       DO jk = 1, jpk
-         glob_sum_3d = glob_sum_3d + SUM( ptab(:,:,jk)*tmask_i(:,:) )
-      END DO
-      IF( lk_mpp )   CALL mpp_sum( glob_sum_3d )
-       
+         glob_sum = glob_sum + SUM( ptab(:,:,jk)*tmask_i(:,:) )
+      END DO
+      IF( lk_mpp )   CALL mpp_sum( glob_sum )
+      !
    END FUNCTION glob_sum_3d
 
-   FUNCTION glob_sum_2d_a( ptab1, ptab2 )
+
+   FUNCTION glob_sum_2d_a( ptab1, ptab2 ) RESULT( glob_sum )
       !!-----------------------------------------------------------------------
       !!                  ***  FUNCTION  glob_sum_2D _a ***
       !!
-      !! ** Purpose : perform a sum on the global domain of two 2D array
-      !!-----------------------------------------------------------------------
-      REAL(wp), DIMENSION(:,:) :: ptab1, ptab2
-      REAL(wp), DIMENSION(2)   :: glob_sum_2d_a
-      !!-----------------------------------------------------------------------
-                    
-      glob_sum_2d_a(1) = SUM( ptab1(:,:)*tmask_i(:,:) )
-      glob_sum_2d_a(2) = SUM( ptab2(:,:)*tmask_i(:,:) )
-      IF( lk_mpp )   CALL mpp_sum( glob_sum_2d_a,2 )
-       
+      !! ** Purpose : perform a masked sum on the inner global domain of two 2D array
+      !!-----------------------------------------------------------------------
+      REAL(wp), INTENT(in), DIMENSION(:,:) ::   ptab1, ptab2   ! input 2D array
+      REAL(wp)            , DIMENSION(2)   ::   glob_sum       ! global masked sum
+      !!-----------------------------------------------------------------------
+      !             
+      glob_sum(1) = SUM( ptab1(:,:)*tmask_i(:,:) )
+      glob_sum(2) = SUM( ptab2(:,:)*tmask_i(:,:) )
+      IF( lk_mpp )   CALL mpp_sum( glob_sum, 2 )
+      !
    END FUNCTION glob_sum_2d_a
  
-   FUNCTION glob_sum_3d_a( ptab1, ptab2 )
+ 
+   FUNCTION glob_sum_3d_a( ptab1, ptab2 ) RESULT( glob_sum )
       !!-----------------------------------------------------------------------
       !!                  ***  FUNCTION  glob_sum_3D_a ***
       !!
-      !! ** Purpose : perform a sum on the global domain of two 3D array
-      !!-----------------------------------------------------------------------
-      REAL(wp), DIMENSION(:,:,:) :: ptab1, ptab2
-      REAL(wp), DIMENSION(2)     :: glob_sum_3d_a
-      !
+      !! ** Purpose : perform a masked sum on the inner global domain of two 3D array
+      !!-----------------------------------------------------------------------
+      REAL(wp), INTENT(in), DIMENSION(:,:,:) ::   ptab1, ptab2   ! input 3D array
+      REAL(wp)            , DIMENSION(2)     ::   glob_sum       ! global masked sum
+      !!
       INTEGER :: jk
       !!-----------------------------------------------------------------------
-       
-      glob_sum_3d_a(:) = 0.e0
+      !
+      glob_sum(:) = 0.e0
       DO jk = 1, jpk
-         glob_sum_3d_a(1) = glob_sum_3d_a(1) + SUM( ptab1(:,:,jk)*tmask_i(:,:) )
-         glob_sum_3d_a(2) = glob_sum_3d_a(2) + SUM( ptab2(:,:,jk)*tmask_i(:,:) )
-      END DO
-      IF( lk_mpp )   CALL mpp_sum( glob_sum_3d_a,2 )
-       
+         glob_sum(1) = glob_sum(1) + SUM( ptab1(:,:,jk)*tmask_i(:,:) )
+         glob_sum(2) = glob_sum(2) + SUM( ptab2(:,:,jk)*tmask_i(:,:) )
+      END DO
+      IF( lk_mpp )   CALL mpp_sum( glob_sum, 2 )
+      !
    END FUNCTION glob_sum_3d_a
 
-#if defined key_mpp_rep  
-   FUNCTION mpp_sum_cmpx( pval )
-      !!----------------------------------------------------------------------
-      !!                  ***  FUNCTION  mpp_sum_cmpx ***
+#else  
+   !!----------------------------------------------------------------------
+   !!   'key_mpp_rep'                                   MPP reproducibility
+   !!----------------------------------------------------------------------
+   
+   FUNCTION glob_sum_2d( ptab ) RESULT( glob_sum )
+      !!----------------------------------------------------------------------
+      !!                  ***  FUNCTION  glob_sum_2d ***
       !!
       !! ** Purpose : perform a sum in calling DDPDD routine
-      !!
-      !!----------------------------------------------------------------------
-      REAL(wp) :: mpp_sum_cmpx
-      !
-      REAL(wp), DIMENSION(jpi,jpj), INTENT(IN) :: &
-         & pval
-      COMPLEX(wp):: ctmp
-      REAL(wp) ::ztmp
-      INTEGER :: ji,jj
-      !!-----------------------------------------------------------------------
-     
+      !!----------------------------------------------------------------------
+      REAL(wp), INTENT(in), DIMENSION(jpi,jpj) ::   ptab
+      REAL(wp)                                 ::   glob_sum   ! global masked sum
+      !!
+      COMPLEX(wp)::   ctmp
+      REAL(wp)   ::   ztmp
+      INTEGER    ::   ji, jj   ! dummy loop indices
+      !!-----------------------------------------------------------------------
+      !
       ztmp = 0.e0
-      ctmp = CMPLX(0.e0,0.e0,wp)
-      DO jj = 1,jpj
+      ctmp = CMPLX( 0.e0, 0.e0, wp )
+      DO jj = 1, jpj
          DO ji =1, jpi
-         ztmp =  pval(ji,jj) * tmask_i(ji,jj)
-         CALL DDPDD(CMPLX(ztmp,0.e0,wp),ctmp)
+         ztmp =  ptab(ji,jj) * tmask_i(ji,jj)
+         CALL DDPDD( CMPLX( ztmp, 0.e0, wp ), ctmp )
          END DO
       END DO
       IF( lk_mpp )   CALL mpp_sum( ctmp )   ! sum over the global domain
-      mpp_sum_cmpx= REAL(ctmp,wp)
-       
-   END FUNCTION mpp_sum_cmpx
+      glob_sum = REAL(ctmp,wp)
+      !
+   END FUNCTION glob_sum_2d   
+
+
+   FUNCTION glob_sum_3d( ptab ) RESULT( glob_sum )
+      !!----------------------------------------------------------------------
+      !!                  ***  FUNCTION  glob_sum_3d ***
+      !!
+      !! ** Purpose : perform a sum on a 3D array in calling DDPDD routine
+      !!----------------------------------------------------------------------
+      REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk) ::   ptab
+      REAL(wp)                                     ::   glob_sum   ! global masked sum
+      !!
+      COMPLEX(wp)::   ctmp
+      REAL(wp)   ::   ztmp
+      INTEGER    ::   ji, jj, jk   ! dummy loop indices
+      !!-----------------------------------------------------------------------
+      !
+      ztmp = 0.e0
+      ctmp = CMPLX( 0.e0, 0.e0, wp )
+      DO jk = 1, jpk
+         DO jj = 1, jpj
+            DO ji =1, jpi
+            ztmp =  ptab(ji,jj,jk) * tmask_i(ji,jj)
+            CALL DDPDD( CMPLX( ztmp, 0.e0, wp ), ctmp )
+            END DO
+         END DO    
+      END DO
+      IF( lk_mpp )   CALL mpp_sum( ctmp )   ! sum over the global domain
+      glob_sum = REAL(ctmp,wp)
+      !
+   END FUNCTION glob_sum_3d   
+
+
+   FUNCTION glob_sum_2d_a( ptab1, ptab2 ) RESULT( glob_sum )
+      !!----------------------------------------------------------------------
+      !!                  ***  FUNCTION  glob_sum_2d_a ***
+      !!
+      !! ** Purpose : perform a sum on two 2D arrays in calling DDPDD routine
+      !!----------------------------------------------------------------------
+      REAL(wp), INTENT(in), DIMENSION(jpi,jpj) ::   ptab1, ptab2
+      REAL(wp)                                 ::   glob_sum   ! global masked sum
+      !!
+      COMPLEX(wp)::   ctmp
+      REAL(wp)   ::   ztmp
+      INTEGER    ::   ji, jj   ! dummy loop indices
+      !!-----------------------------------------------------------------------
+      !
+      ztmp = 0.e0
+      ctmp = CMPLX( 0.e0, 0.e0, wp )
+      DO jj = 1, jpj
+         DO ji =1, jpi
+         ztmp =  ptab1(ji,jj) * tmask_i(ji,jj)
+         CALL DDPDD( CMPLX( ztmp, 0.e0, wp ), ctmp )
+         ztmp =  ptab2(ji,jj) * tmask_i(ji,jj)
+         CALL DDPDD( CMPLX( ztmp, 0.e0, wp ), ctmp )
+         END DO
+      END DO
+      IF( lk_mpp )   CALL mpp_sum( ctmp )   ! sum over the global domain
+      glob_sum = REAL(ctmp,wp)
+      !
+   END FUNCTION glob_sum_2d_a   
+
+
+   FUNCTION glob_sum_3d_a( ptab1, ptab2 ) RESULT( glob_sum )
+      !!----------------------------------------------------------------------
+      !!                  ***  FUNCTION  glob_sum_3d_a ***
+      !!
+      !! ** Purpose : perform a sum on two 3D array in calling DDPDD routine
+      !!----------------------------------------------------------------------
+      REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk) ::   ptab1, ptab2
+      REAL(wp)                                     ::   glob_sum   ! global masked sum
+      !!
+      COMPLEX(wp)::   ctmp
+      REAL(wp)   ::   ztmp
+      INTEGER    ::   ji, jj, jk   ! dummy loop indices
+      !!-----------------------------------------------------------------------
+      !
+      ztmp = 0.e0
+      ctmp = CMPLX( 0.e0, 0.e0, wp )
+      DO jk = 1, jpk
+         DO jj = 1, jpj
+            DO ji =1, jpi
+            ztmp =  ptab1(ji,jj,jk) * tmask_i(ji,jj)
+            CALL DDPDD( CMPLX( ztmp, 0.e0, wp ), ctmp )
+            ztmp =  ptab2(ji,jj,jk) * tmask_i(ji,jj)
+            CALL DDPDD( CMPLX( ztmp, 0.e0, wp ), ctmp )
+            END DO
+         END DO    
+      END DO
+      IF( lk_mpp )   CALL mpp_sum( ctmp )   ! sum over the global domain
+      glob_sum = REAL(ctmp,wp)
+      !
+   END FUNCTION glob_sum_3d_a   
+
 
    SUBROUTINE DDPDD( ydda, yddb )
@@ -159,27 +262,30 @@
       !! References : Using Acurate Arithmetics to Improve Numerical
       !!              Reproducibility and Sability in Parallel Applications
-      !!              Yun HE and Chris H. Q. DING, Journal of Supercomputing
-      !!                                            18, 259-277, 2001 
-      !!----------------------------------------------------------------------
-
-      COMPLEX(wp), INTENT(in)     :: ydda
-      COMPLEX(wp), INTENT(inout)  :: yddb
-      
+      !!              Yun HE and Chris H. Q. DING, Journal of Supercomputing 18, 259-277, 2001 
+      !!----------------------------------------------------------------------
+      COMPLEX(wp), INTENT(in   ) ::   ydda
+      COMPLEX(wp), INTENT(inout) ::   yddb
+      !
       REAL(wp) :: zerr, zt1, zt2  ! local work variables
-
+      !!-----------------------------------------------------------------------
+      !
       ! Compute ydda + yddb using Knuth's trick.
-      zt1  = real(ydda) + real(yddb)
-      zerr = zt1 - real(ydda)
-      zt2  = ((real(yddb) - zerr) + (real(ydda) - (zt1 - zerr))) &
-            + aimag(ydda) + aimag(yddb)
-
+      zt1  = REAL(ydda) + REAL(yddb)
+      zerr = zt1 - REAL(ydda)
+      zt2  = ( (REAL(yddb) - zerr) + (REAL(ydda) - (zt1 - zerr)) )   &
+         &   + AIMAG(ydda)         + AIMAG(yddb)
+      !
       ! The result is t1 + t2, after normalization.
-      yddb = cmplx ( zt1 + zt2, zt2 - ((zt1 + zt2) - zt1),wp )
-      
+      yddb = CMPLX( zt1 + zt2, zt2 - ((zt1 + zt2) - zt1), wp )
+      !
    END SUBROUTINE DDPDD
 #endif
 
 #if defined key_nosignedzero
-   FUNCTION SIGN_SCALAR(pa,pb)
+   !!----------------------------------------------------------------------
+   !!   'key_nosignedzero'                                         F90 SIGN
+   !!----------------------------------------------------------------------
+   
+   FUNCTION SIGN_SCALAR( pa, pb )
       !!-----------------------------------------------------------------------
       !!                  ***  FUNCTION SIGN_SCALAR  ***
@@ -188,14 +294,13 @@
       !!-----------------------------------------------------------------------
       REAL(wp) :: pa,pb          ! input
-      REAL(wp) :: SIGN_SCALAR  ! result
-      IF ( pb >= 0.e0) THEN
-         SIGN_SCALAR = ABS(pa)
-      ELSE
-         SIGN_SCALAR =-ABS(pa)
+      REAL(wp) :: SIGN_SCALAR    ! result
+      !!-----------------------------------------------------------------------
+      IF ( pb >= 0.e0) THEN   ;   SIGN_SCALAR = ABS(pa)
+      ELSE                    ;   SIGN_SCALAR =-ABS(pa)
       ENDIF
-
    END FUNCTION SIGN_SCALAR
 
-   FUNCTION SIGN_ARRAY_1D(pa,pb) 
+
+   FUNCTION SIGN_ARRAY_1D( pa, pb ) 
       !!-----------------------------------------------------------------------
       !!                  ***  FUNCTION SIGN_ARRAY_1D  ***
@@ -203,13 +308,12 @@
       !! ** Purpose : overwrite f95 behaviour of intrinsinc sign function
       !!-----------------------------------------------------------------------
-      REAL(wp) :: pa,pb(:)      ! input
+      REAL(wp) :: pa,pb(:)                   ! input
       REAL(wp) :: SIGN_ARRAY_1D(SIZE(pb,1))  ! result
-      WHERE ( pb >= 0.e0 )
-         SIGN_ARRAY_1D = ABS(pa)
-      ELSEWHERE
-         SIGN_ARRAY_1D =-ABS(pa)
-      END WHERE
-
+      !!-----------------------------------------------------------------------
+      WHERE ( pb >= 0.e0 )   ;   SIGN_ARRAY_1D = ABS(pa)
+      ELSEWHERE              ;   SIGN_ARRAY_1D =-ABS(pa)
+      END WHERE
    END FUNCTION SIGN_ARRAY_1D
+
 
    FUNCTION SIGN_ARRAY_2D(pa,pb) 
@@ -221,11 +325,8 @@
       REAL(wp) :: pa,pb(:,:)      ! input
       REAL(wp) :: SIGN_ARRAY_2D(SIZE(pb,1),SIZE(pb,2))  ! result
-
-      WHERE ( pb >= 0.e0 )
-         SIGN_ARRAY_2D = ABS(pa)
-      ELSEWHERE
-         SIGN_ARRAY_2D =-ABS(pa)
-      END WHERE
-
+      !!-----------------------------------------------------------------------
+      WHERE ( pb >= 0.e0 )   ;   SIGN_ARRAY_2D = ABS(pa)
+      ELSEWHERE              ;   SIGN_ARRAY_2D =-ABS(pa)
+      END WHERE
    END FUNCTION SIGN_ARRAY_2D
 
@@ -238,12 +339,11 @@
       REAL(wp) :: pa,pb(:,:,:)      ! input
       REAL(wp) :: SIGN_ARRAY_3D(SIZE(pb,1),SIZE(pb,2),SIZE(pb,3))  ! result
-      WHERE ( pb >= 0.e0 )
-         SIGN_ARRAY_3D = ABS(pa)
-      ELSEWHERE
-         SIGN_ARRAY_3D =-ABS(pa)
-      END WHERE
-
+      !!-----------------------------------------------------------------------
+      WHERE ( pb >= 0.e0 )   ;   SIGN_ARRAY_3D = ABS(pa)
+      ELSEWHERE              ;   SIGN_ARRAY_3D =-ABS(pa)
+      END WHERE
    END FUNCTION SIGN_ARRAY_3D
 
+
    FUNCTION SIGN_ARRAY_1D_A(pa,pb) 
       !!-----------------------------------------------------------------------
@@ -253,13 +353,11 @@
       !!-----------------------------------------------------------------------
       REAL(wp) :: pa(:),pb(:)      ! input
-      REAL(wp) :: SIGN_ARRAY_1D_A(SIZE(b,1))  ! result
-
-      WHERE ( pb >= 0.e0 )
-         SIGN_ARRAY_1D_A = ABS(pa)
-      ELSEWHERE
-         SIGN_ARRAY_1D_A =-ABS(pa)
-      END WHERE
-
+      REAL(wp) :: SIGN_ARRAY_1D_A(SIZE(pb,1))  ! result
+      !!-----------------------------------------------------------------------
+      WHERE ( pb >= 0.e0 )   ;   SIGN_ARRAY_1D_A = ABS(pa)
+      ELSEWHERE              ;   SIGN_ARRAY_1D_A =-ABS(pa)
+      END WHERE
    END FUNCTION SIGN_ARRAY_1D_A
+
 
    FUNCTION SIGN_ARRAY_2D_A(pa,pb) 
@@ -271,12 +369,10 @@
       REAL(wp) :: pa(:,:),pb(:,:)      ! input
       REAL(wp) :: SIGN_ARRAY_2D_A(SIZE(pb,1),SIZE(pb,2))  ! result
-
-      WHERE ( pb >= 0.e0 )
-         SIGN_ARRAY_2D_A = ABS(pa)
-      ELSEWHERE
-         SIGN_ARRAY_2D_A =-ABS(pa)
-      END WHERE
-
+      !!-----------------------------------------------------------------------
+      WHERE ( pb >= 0.e0 )   ;   SIGN_ARRAY_2D_A = ABS(pa)
+      ELSEWHERE              ;   SIGN_ARRAY_2D_A =-ABS(pa)
+      END WHERE
    END FUNCTION SIGN_ARRAY_2D_A
+
 
    FUNCTION SIGN_ARRAY_3D_A(pa,pb) 
@@ -288,12 +384,10 @@
       REAL(wp) :: pa(:,:,:),pb(:,:,:)  ! input
       REAL(wp) :: SIGN_ARRAY_3D_A(SIZE(pb,1),SIZE(pb,2),SIZE(pb,3)) ! result
-
-      WHERE ( pb >= 0.e0 )
-         SIGN_ARRAY_3D_A = ABS(pa)
-      ELSEWHERE
-         SIGN_ARRAY_3D_A =-ABS(pa)
-      END WHERE
-
+      !!-----------------------------------------------------------------------
+      WHERE ( pb >= 0.e0 )   ;   SIGN_ARRAY_3D_A = ABS(pa)
+      ELSEWHERE              ;   SIGN_ARRAY_3D_A =-ABS(pa)
+      END WHERE
    END FUNCTION SIGN_ARRAY_3D_A
+
 
    FUNCTION SIGN_ARRAY_1D_B(pa,pb) 
@@ -305,12 +399,10 @@
       REAL(wp) :: pa(:),pb      ! input
       REAL(wp) :: SIGN_ARRAY_1D_B(SIZE(pa,1))  ! result
-
-      IF ( pb >= 0.e0 ) THEN
-         SIGN_ARRAY_1D_B = ABS(pa)
-      ELSE
-         SIGN_ARRAY_1D_B =-ABS(pa)
+      !!-----------------------------------------------------------------------
+      IF( pb >= 0.e0 ) THEN   ;   SIGN_ARRAY_1D_B = ABS(pa)
+      ELSE                    ;   SIGN_ARRAY_1D_B =-ABS(pa)
       ENDIF
-
    END FUNCTION SIGN_ARRAY_1D_B
+
 
    FUNCTION SIGN_ARRAY_2D_B(pa,pb) 
@@ -322,12 +414,10 @@
       REAL(wp) :: pa(:,:),pb      ! input
       REAL(wp) :: SIGN_ARRAY_2D_B(SIZE(pa,1),SIZE(pa,2))  ! result
-
-      IF ( pb >= 0.e0 ) THEN
-         SIGN_ARRAY_2D_B = ABS(pa)
-      ELSE
-         SIGN_ARRAY_2D_B =-ABS(pa)
+      !!-----------------------------------------------------------------------
+      IF( pb >= 0.e0 ) THEN   ;   SIGN_ARRAY_2D_B = ABS(pa)
+      ELSE                    ;   SIGN_ARRAY_2D_B =-ABS(pa)
       ENDIF
-
    END FUNCTION SIGN_ARRAY_2D_B
+
 
    FUNCTION SIGN_ARRAY_3D_B(pa,pb) 
@@ -339,13 +429,11 @@
       REAL(wp) :: pa(:,:,:),pb      ! input
       REAL(wp) :: SIGN_ARRAY_3D_B(SIZE(pa,1),SIZE(pa,2),SIZE(pa,3))  ! result
-
-      IF (pb >= 0.e0 ) THEN
-         SIGN_ARRAY_3D_B = ABS(pa)
-      ELSE
-         SIGN_ARRAY_3D_B =-ABS(pa)
+      !!-----------------------------------------------------------------------
+      IF( pb >= 0.e0 ) THEN   ;   SIGN_ARRAY_3D_B = ABS(pa)
+      ELSE                    ;   SIGN_ARRAY_3D_B =-ABS(pa)
       ENDIF
- 
    END FUNCTION SIGN_ARRAY_3D_B
 #endif
 
+   !!======================================================================
 END MODULE lib_fortran
