Index: NEMO/trunk/src/OCE/ICB/icb_oce.F90
===================================================================
--- NEMO/trunk/src/OCE/ICB/icb_oce.F90	(revision 10665)
+++ NEMO/trunk/src/OCE/ICB/icb_oce.F90	(revision 10691)
@@ -89,4 +89,5 @@
    REAL(wp), PUBLIC, DIMENSION(:,:), ALLOCATABLE ::   ua_e, va_e
    REAL(wp), PUBLIC, DIMENSION(:,:), ALLOCATABLE ::   ssh_e
+   REAL(wp), PUBLIC, DIMENSION(:,:), ALLOCATABLE ::   tmask_e, umask_e, vmask_e
 #if defined key_si3 || defined key_cice
    REAL(wp), PUBLIC, DIMENSION(:,:), ALLOCATABLE ::   ui_e, vi_e
@@ -169,6 +170,6 @@
       !
       ! expanded arrays for bilinear interpolation
-      ALLOCATE( uo_e(0:jpi+1,0:jpj+1) , ua_e(0:jpi+1,0:jpj+1) ,   &
-         &      vo_e(0:jpi+1,0:jpj+1) , va_e(0:jpi+1,0:jpj+1) ,   &
+      ALLOCATE( uo_e(0:jpi+1,0:jpj+1) , ua_e(0:jpi+1,0:jpj+1) ,    &
+         &      vo_e(0:jpi+1,0:jpj+1) , va_e(0:jpi+1,0:jpj+1) ,    &
 #if defined key_si3 || defined key_cice
          &      ui_e(0:jpi+1,0:jpj+1) ,                            &
@@ -183,4 +184,8 @@
       icb_alloc = icb_alloc + ill
 
+      ALLOCATE( tmask_e(0:jpi+1,0:jpj+1), umask_e(0:jpi+1,0:jpj+1), vmask_e(0:jpi+1,0:jpj+1), &
+         &      STAT=ill)
+      icb_alloc = icb_alloc + ill
+
       ALLOCATE( nicbfldpts(jpi) , nicbflddest(jpi) , nicbfldproc(jpni) , &
          &      nicbfldnsend(jpni), nicbfldexpect(jpni) , nicbfldreq(jpni), STAT=ill)
Index: NEMO/trunk/src/OCE/ICB/icbini.F90
===================================================================
--- NEMO/trunk/src/OCE/ICB/icbini.F90	(revision 10665)
+++ NEMO/trunk/src/OCE/ICB/icbini.F90	(revision 10691)
@@ -224,4 +224,12 @@
       src_calving_hflx(:,:) = 0._wp
 
+      ! definition of extended surface masked needed by icb_bilin_h
+      tmask_e(:,:) = 0._wp   ;   tmask_e(1:jpi,1:jpj) = tmask(:,:,1)
+      umask_e(:,:) = 0._wp   ;   umask_e(1:jpi,1:jpj) = umask(:,:,1)
+      vmask_e(:,:) = 0._wp   ;   vmask_e(1:jpi,1:jpj) = vmask(:,:,1)
+      CALL lbc_lnk_icb( 'icbini', tmask_e, 'T', +1._wp, 1, 1 )
+      CALL lbc_lnk_icb( 'icbini', umask_e, 'T', +1._wp, 1, 1 )
+      CALL lbc_lnk_icb( 'icbini', vmask_e, 'T', +1._wp, 1, 1 )
+      !
       ! assign each new iceberg with a unique number constructed from the processor number
       ! and incremented by the total number of processors
Index: NEMO/trunk/src/OCE/ICB/icbutl.F90
===================================================================
--- NEMO/trunk/src/OCE/ICB/icbutl.F90	(revision 10665)
+++ NEMO/trunk/src/OCE/ICB/icbutl.F90	(revision 10691)
@@ -131,4 +131,6 @@
       !!             is half the off shore value, wile the normal-to-the-coast value is zero.
       !!             This is OK as a starting point.
+      !!       !!pm  HARD CODED: - rho_air now computed in sbcblk (what are the effect ?)
+      !!                         - drag coefficient (should it be namelist parameter ?)
       !!
       !!----------------------------------------------------------------------
@@ -142,16 +144,16 @@
       !!----------------------------------------------------------------------
 
-      pe1 = icb_utl_bilin_e( e1t, e1u, e1v, e1f, pi, pj )     ! scale factors
+      pe1 = icb_utl_bilin_e( e1t, e1u, e1v, e1f, pi, pj )      ! scale factors
       pe2 = icb_utl_bilin_e( e2t, e2u, e2v, e2f, pi, pj )
       !
-      puo  = icb_utl_bilin_h( uo_e, pi, pj, 'U' )             ! ocean velocities
-      pvo  = icb_utl_bilin_h( vo_e, pi, pj, 'V' )
-      psst = icb_utl_bilin_h( tt_e, pi, pj, 'T' )             ! SST
-      pcn  = icb_utl_bilin_h( fr_e , pi, pj, 'T' )            ! ice concentration
-      pff  = icb_utl_bilin_h( ff_e , pi, pj, 'F' )            ! Coriolis parameter
-      !
-      pua  = icb_utl_bilin_h( ua_e , pi, pj, 'U' )            ! 10m wind
-      pva  = icb_utl_bilin_h( va_e , pi, pj, 'V' )            ! here (ua,va) are stress => rough conversion from stress to speed
-      zcd  = 1.22_wp * 1.5e-3_wp                              ! air density * drag coefficient
+      puo  = icb_utl_bilin_h( uo_e, pi, pj, 'U', .false.  )    ! ocean velocities
+      pvo  = icb_utl_bilin_h( vo_e, pi, pj, 'V', .false.  )
+      psst = icb_utl_bilin_h( tt_e, pi, pj, 'T', .true.   )    ! SST
+      pcn  = icb_utl_bilin_h( fr_e, pi, pj, 'T', .true.   )    ! ice concentration
+      pff  = icb_utl_bilin_h( ff_e, pi, pj, 'F', .false.  )    ! Coriolis parameter
+      !
+      pua  = icb_utl_bilin_h( ua_e, pi, pj, 'U', .true.   )    ! 10m wind
+      pva  = icb_utl_bilin_h( va_e, pi, pj, 'V', .true.   )    ! here (ua,va) are stress => rough conversion from stress to speed
+      zcd  = 1.22_wp * 1.5e-3_wp                               ! air density * drag coefficient 
       zmod = 1._wp / MAX(  1.e-20, SQRT(  zcd * SQRT( pua*pua + pva*pva)  )  )
       pua  = pua * zmod                                       ! note: stress module=0 necessarly implies ua=va=0
@@ -159,7 +161,7 @@
 
 #if defined key_si3
-      pui = icb_utl_bilin_h( ui_e , pi, pj, 'U' )              ! sea-ice velocities
-      pvi = icb_utl_bilin_h( vi_e , pi, pj, 'V' )
-      phi = icb_utl_bilin_h( hicth, pi, pj, 'T' )              ! ice thickness
+      pui = icb_utl_bilin_h( ui_e , pi, pj, 'U', .false. )    ! sea-ice velocities
+      pvi = icb_utl_bilin_h( vi_e , pi, pj, 'V', .false. )
+      phi = icb_utl_bilin_h( hicth, pi, pj, 'T', .true.  )    ! ice thickness
 #else
       pui = 0._wp
@@ -169,13 +171,13 @@
 
       ! Estimate SSH gradient in i- and j-direction (centred evaluation)
-      pssh_i = ( icb_utl_bilin_h( ssh_e, pi+0.1_wp, pj, 'T' ) -   &
-         &       icb_utl_bilin_h( ssh_e, pi-0.1_wp, pj, 'T' )  ) / ( 0.2_wp * pe1 )
-      pssh_j = ( icb_utl_bilin_h( ssh_e, pi, pj+0.1_wp, 'T' ) -   &
-         &       icb_utl_bilin_h( ssh_e, pi, pj-0.1_wp, 'T' )  ) / ( 0.2_wp * pe2 )
+      pssh_i = ( icb_utl_bilin_h( ssh_e, pi+0.1_wp, pj, 'T', .true. ) -   &
+         &       icb_utl_bilin_h( ssh_e, pi-0.1_wp, pj, 'T', .true. )  ) / ( 0.2_wp * pe1 )
+      pssh_j = ( icb_utl_bilin_h( ssh_e, pi, pj+0.1_wp, 'T', .true. ) -   &
+         &       icb_utl_bilin_h( ssh_e, pi, pj-0.1_wp, 'T', .true. )  ) / ( 0.2_wp * pe2 )
       !
    END SUBROUTINE icb_utl_interp
 
 
-   REAL(wp) FUNCTION icb_utl_bilin_h( pfld, pi, pj, cd_type )
+   REAL(wp) FUNCTION icb_utl_bilin_h( pfld, pi, pj, cd_type, plmask )
       !!----------------------------------------------------------------------
       !!                  ***  FUNCTION icb_utl_bilin  ***
@@ -191,7 +193,10 @@
       REAL(wp)                            , INTENT(in) ::   pi, pj    ! targeted coordinates in (i,j) referential
       CHARACTER(len=1)                    , INTENT(in) ::   cd_type   ! type of pfld array grid-points: = T , U , V or F points
+      LOGICAL                             , INTENT(in) ::   plmask    ! special treatment of mask point
       !
       INTEGER  ::   ii, ij   ! local integer
       REAL(wp) ::   zi, zj   ! local real
+      REAL(wp) :: zw1, zw2, zw3, zw4
+      REAL(wp), DIMENSION(4) :: zmask
       !!----------------------------------------------------------------------
       !
@@ -224,18 +229,43 @@
       ! find position in this processor. Prevent near edge problems (see #1389)
       !
-      IF    ( ii < mig( 1 ) ) THEN   ;   ii = 1
-      ELSEIF( ii > mig(jpi) ) THEN   ;   ii = jpi
+      IF    ( ii < mig( 1 ) ) THEN   ;   ii = 1   ; 
+      ELSEIF( ii > mig(jpi) ) THEN   ;   ii = jpi ;
       ELSE                           ;   ii = mi1(ii)
       ENDIF
-      IF    ( ij < mjg( 1 ) ) THEN   ;   ij = 1
-      ELSEIF( ij > mjg(jpj) ) THEN   ;   ij = jpj
+      IF    ( ij < mjg( 1 ) ) THEN   ;   ij = 1   ;
+      ELSEIF( ij > mjg(jpj) ) THEN   ;   ij = jpj ;
       ELSE                           ;   ij  = mj1(ij)
       ENDIF
       !
-      IF( ii == jpi )   ii = ii-1      
-      IF( ij == jpj )   ij = ij-1
-      !
-      icb_utl_bilin_h = ( pfld(ii,ij  ) * (1.-zi) + pfld(ii+1,ij  ) * zi ) * (1.-zj)   &
-         &            + ( pfld(ii,ij+1) * (1.-zi) + pfld(ii+1,ij+1) * zi ) *     zj
+      IF( ii == jpi ) ii = ii-1
+      IF( ij == jpj ) ij = ij-1
+      !
+      ! define mask array 
+      IF (plmask) THEN
+         ! land value is not used in the interpolation
+         SELECT CASE ( cd_type )
+         CASE ( 'T' )
+            zmask = (/tmask_e(ii,ij), tmask_e(ii+1,ij), tmask_e(ii,ij+1), tmask_e(ii+1,ij+1)/)
+         CASE ( 'U' )
+            zmask = (/umask_e(ii,ij), umask_e(ii+1,ij), umask_e(ii,ij+1), umask_e(ii+1,ij+1)/)
+         CASE ( 'V' )
+            zmask = (/vmask_e(ii,ij), vmask_e(ii+1,ij), vmask_e(ii,ij+1), vmask_e(ii+1,ij+1)/)
+         CASE ( 'F' )
+            ! F case only used for coriolis, ff_f is not mask so zmask = 1
+            zmask = 1.
+         END SELECT
+      ELSE
+         ! land value is used during interpolation
+         zmask = 1.
+      END iF
+      !
+      ! compute weight
+      zw1 = zmask(1) * (1._wp-zi) * (1._wp-zj)
+      zw2 = zmask(2) *        zi  * (1._wp-zj)
+      zw3 = zmask(3) * (1._wp-zi) *        zj
+      zw4 = zmask(4) *        zi  *        zj
+      !
+      ! compute interpolated value
+      icb_utl_bilin_h = ( pfld(ii,ij)*zw1 + pfld(ii+1,ij)*zw2 + pfld(ii,ij+1)*zw3 + pfld(ii+1,ij+1)*zw4 ) / MAX(1.e-20, zw1+zw2+zw3+zw4) 
       !
    END FUNCTION icb_utl_bilin_h
