Index: branches/UKMO/dev_r6501_GO6_package_trunk/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90
===================================================================
--- branches/UKMO/dev_r6501_GO6_package_trunk/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90	(revision 6503)
+++ branches/UKMO/dev_r6501_GO6_package_trunk/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90	(revision 6507)
@@ -106,5 +106,8 @@
    INTEGER, PARAMETER ::   jpr_e3t1st = 41   ! first T level thickness 
    INTEGER, PARAMETER ::   jpr_fraqsr = 42   ! fraction of solar net radiation absorbed in the first ocean level
-   INTEGER, PARAMETER ::   jprcv      = 42   ! total number of fields received
+   INTEGER, PARAMETER ::   jpr_ts_ice = 43   ! skin temperature of sea-ice (used for melt-ponds)
+   INTEGER, PARAMETER ::   jpr_grnm   = 44   ! Greenland ice mass
+   INTEGER, PARAMETER ::   jpr_antm   = 45   ! Antarctic ice mass
+   INTEGER, PARAMETER ::   jprcv      = 45   ! total number of fields received
 
    INTEGER, PARAMETER ::   jps_fice   =  1   ! ice fraction sent to the atmosphere
@@ -136,5 +139,10 @@
    INTEGER, PARAMETER ::   jps_e3t1st = 27   ! first level depth (vvl)
    INTEGER, PARAMETER ::   jps_fraqsr = 28   ! fraction of solar net radiation absorbed in the first ocean level
-   INTEGER, PARAMETER ::   jpsnd      = 28   ! total number of fields sended
+   INTEGER, PARAMETER ::   jps_a_p    = 29   ! meltpond fraction  
+   INTEGER, PARAMETER ::   jps_ht_p   = 30   ! meltpond depth (m) 
+   INTEGER, PARAMETER ::   jps_kice   = 31   ! ice surface layer thermal conductivity
+   INTEGER, PARAMETER ::   jps_sstfrz = 32   ! sea-surface freezing temperature
+   INTEGER, PARAMETER ::   jps_fice1  = 33   ! first-order ice concentration (for time-travelling ice coupling)
+   INTEGER, PARAMETER ::   jpsnd      = 33   ! total number of fields sended
 
    !                                  !!** namelist namsbc_cpl **
@@ -147,8 +155,9 @@
    END TYPE FLD_C
    !                                   ! Send to the atmosphere  
-   TYPE(FLD_C) ::   sn_snd_temp, sn_snd_alb, sn_snd_thick, sn_snd_crt, sn_snd_co2                        
+   TYPE(FLD_C) ::   sn_snd_temp, sn_snd_alb, sn_snd_thick, sn_snd_crt, sn_snd_co2, sn_snd_cond, sn_snd_mpnd, sn_snd_sstfrz, sn_snd_thick1
+
    !                                   ! Received from the atmosphere
    TYPE(FLD_C) ::   sn_rcv_w10m, sn_rcv_taumod, sn_rcv_tau, sn_rcv_dqnsdt, sn_rcv_qsr, sn_rcv_qns, sn_rcv_emp, sn_rcv_rnf
-   TYPE(FLD_C) ::   sn_rcv_cal, sn_rcv_iceflx, sn_rcv_co2                        
+   TYPE(FLD_C) ::   sn_rcv_cal, sn_rcv_iceflx, sn_rcv_co2, sn_rcv_ts_ice, sn_rcv_grnm, sn_rcv_antm
    !                                   ! Other namelist parameters
    INTEGER     ::   nn_cplmodel           ! Maximum number of models to/from which NEMO is potentialy sending/receiving data
@@ -214,8 +223,11 @@
       REAL(wp), POINTER, DIMENSION(:,:) ::   zacs, zaos
       !!
-      NAMELIST/namsbc_cpl/  sn_snd_temp, sn_snd_alb   , sn_snd_thick, sn_snd_crt   , sn_snd_co2,      &
-         &                  sn_rcv_w10m, sn_rcv_taumod, sn_rcv_tau  , sn_rcv_dqnsdt, sn_rcv_qsr,      &
-         &                  sn_rcv_qns , sn_rcv_emp   , sn_rcv_rnf  , sn_rcv_cal   , sn_rcv_iceflx,   &
-         &                  sn_rcv_co2 , nn_cplmodel  , ln_usecplmask
+      NAMELIST/namsbc_cpl/  sn_snd_temp, sn_snd_alb   , sn_snd_thick , sn_snd_crt   , sn_snd_co2,     &
+         &                  sn_snd_cond, sn_snd_mpnd  , sn_snd_sstfrz, sn_snd_thick1,                 &
+         &                  sn_rcv_w10m, sn_rcv_taumod, sn_rcv_tau   , sn_rcv_dqnsdt, sn_rcv_qsr,     &
+         &                  sn_rcv_qns , sn_rcv_emp   , sn_rcv_rnf   , sn_rcv_cal   , sn_rcv_iceflx,  &
+         &                  sn_rcv_co2 , sn_rcv_grnm  , sn_rcv_antm  , sn_rcv_ts_ice, nn_cplmodel  ,  &
+         &                  ln_usecplmask, ln_coupled_iceshelf_fluxes, rn_greenland_calving_fraction, &
+         &                  rn_antarctica_calving_fraction, rn_iceshelf_fluxes_tolerance
       !!---------------------------------------------------------------------
       !
@@ -256,4 +268,6 @@
          WRITE(numout,*)'      runoffs                         = ', TRIM(sn_rcv_rnf%cldes   ), ' (', TRIM(sn_rcv_rnf%clcat   ), ')'
          WRITE(numout,*)'      calving                         = ', TRIM(sn_rcv_cal%cldes   ), ' (', TRIM(sn_rcv_cal%clcat   ), ')'
+         WRITE(numout,*)'      Greenland ice mass              = ', TRIM(sn_rcv_grnm%cldes  ), ' (', TRIM(sn_rcv_grnm%clcat  ), ')'
+         WRITE(numout,*)'      Antarctica ice mass             = ', TRIM(sn_rcv_antm%cldes  ), ' (', TRIM(sn_rcv_antm%clcat  ), ')'
          WRITE(numout,*)'      sea ice heat fluxes             = ', TRIM(sn_rcv_iceflx%cldes), ' (', TRIM(sn_rcv_iceflx%clcat), ')'
          WRITE(numout,*)'      atm co2                         = ', TRIM(sn_rcv_co2%cldes   ), ' (', TRIM(sn_rcv_co2%clcat   ), ')'
@@ -267,6 +281,14 @@
          WRITE(numout,*)'                      - mesh          = ', sn_snd_crt%clvgrd
          WRITE(numout,*)'      oce co2 flux                    = ', TRIM(sn_snd_co2%cldes   ), ' (', TRIM(sn_snd_co2%clcat   ), ')'
+         WRITE(numout,*)'      ice effective conductivity      = ', TRIM(sn_snd_cond%cldes   ), ' (', TRIM(sn_snd_cond%clcat   ), ')'
+         WRITE(numout,*)'      meltponds fraction & depth      = ', TRIM(sn_snd_mpnd%cldes  ), ' (', TRIM(sn_snd_mpnd%clcat   ), ')'
+         WRITE(numout,*)'      sea surface freezing temp       = ', TRIM(sn_snd_sstfrz%cldes   ), ' (', TRIM(sn_snd_sstfrz%clcat   ), ')'
+
          WRITE(numout,*)'  nn_cplmodel                         = ', nn_cplmodel
          WRITE(numout,*)'  ln_usecplmask                       = ', ln_usecplmask
+         WRITE(numout,*)'  ln_coupled_iceshelf_fluxes          = ', ln_coupled_iceshelf_fluxes
+         WRITE(numout,*)'  rn_greenland_calving_fraction       = ', rn_greenland_calving_fraction
+         WRITE(numout,*)'  rn_antarctica_calving_fraction      = ', rn_antarctica_calving_fraction
+         WRITE(numout,*)'  rn_iceshelf_fluxes_tolerance        = ', rn_iceshelf_fluxes_tolerance
       ENDIF
 
@@ -381,5 +403,5 @@
       srcv(jpr_snow)%clname = 'OTotSnow'      ! Snow = solid precipitation
       srcv(jpr_tevp)%clname = 'OTotEvap'      ! total evaporation (over oce + ice sublimation)
-      srcv(jpr_ievp)%clname = 'OIceEvap'      ! evaporation over ice = sublimation
+      srcv(jpr_ievp)%clname = 'OIceEvp'      ! evaporation over ice = sublimation
       srcv(jpr_sbpr)%clname = 'OSubMPre'      ! sublimation - liquid precipitation - solid precipitation 
       srcv(jpr_semp)%clname = 'OISubMSn'      ! ice solid water budget = sublimation - solid precipitation
@@ -395,4 +417,7 @@
       END SELECT
       !
+      !Set the number of categories for coupling of sublimation
+      IF ( TRIM( sn_rcv_emp%clcat ) == 'yes' ) srcv(jpr_ievp)%nct = jpl
+      !
       !                                                      ! ------------------------- !
       !                                                      !     Runoffs & Calving     !   
@@ -408,4 +433,6 @@
       !
       srcv(jpr_cal   )%clname = 'OCalving'   ;   IF( TRIM( sn_rcv_cal%cldes ) == 'coupled' )   srcv(jpr_cal)%laction = .TRUE.
+      srcv(jpr_grnm  )%clname = 'OGrnmass'   ;   IF( TRIM( sn_rcv_grnm%cldes ) == 'coupled' )   srcv(jpr_grnm)%laction = .TRUE.
+      srcv(jpr_antm  )%clname = 'OAntmass'   ;   IF( TRIM( sn_rcv_antm%cldes ) == 'coupled' )   srcv(jpr_antm)%laction = .TRUE.
       !
       !                                                      ! ------------------------- !
@@ -481,4 +508,17 @@
          srcv(jpr_topm:jpr_botm)%laction = .TRUE.
       ENDIF
+      
+#if defined key_cice && ! defined key_cice4
+      !                                                      ! ----------------------------- !
+      !                                                      !  sea-ice skin temperature     !   
+      !                                                      !  used in meltpond scheme      !
+      !                                                      !  May be calculated in Atm     !
+      !                                                      ! ----------------------------- !
+      srcv(jpr_ts_ice)%clname = 'OTsfIce'
+      IF ( TRIM( sn_rcv_ts_ice%cldes ) == 'ice' ) srcv(jpr_ts_ice)%laction = .TRUE.
+      IF ( TRIM( sn_rcv_ts_ice%clcat ) == 'yes' ) srcv(jpr_ts_ice)%nct = jpl
+      !TODO: Should there be a consistency check here?
+#endif
+
       !                                                      ! ------------------------------- !
       !                                                      !   OPA-SAS coupling - rcv by opa !   
@@ -598,10 +638,10 @@
       !                                                      ! ------------------------- !
       ssnd(jps_toce)%clname = 'O_SSTSST'
-      ssnd(jps_tice)%clname = 'O_TepIce'
+      ssnd(jps_tice)%clname = 'OTepIce'
       ssnd(jps_tmix)%clname = 'O_TepMix'
       SELECT CASE( TRIM( sn_snd_temp%cldes ) )
       CASE( 'none'                                 )       ! nothing to do
       CASE( 'oce only'                             )   ;   ssnd( jps_toce )%laction = .TRUE.
-      CASE( 'oce and ice' , 'weighted oce and ice' )
+      CASE( 'oce and ice' , 'weighted oce and ice' , 'oce and weighted ice')
          ssnd( (/jps_toce, jps_tice/) )%laction = .TRUE.
          IF ( TRIM( sn_snd_temp%clcat ) == 'yes' )  ssnd(jps_tice)%nct = jpl
@@ -637,6 +677,11 @@
       ssnd(jps_hice)%clname = 'OIceTck'
       ssnd(jps_hsnw)%clname = 'OSnwTck'
+      ssnd(jps_a_p)%clname  = 'OPndFrc'
+      ssnd(jps_ht_p)%clname = 'OPndTck'
+      ssnd(jps_fice1)%clname = 'OIceFrd'
       IF( k_ice /= 0 ) THEN
          ssnd(jps_fice)%laction = .TRUE.                  ! if ice treated in the ocean (even in climato case)
+         ssnd(jps_fice1)%laction = .TRUE.                 ! First-order regridded ice concentration, to be used
+	                                                  ! in producing atmos-to-ice fluxes
 ! Currently no namelist entry to determine sending of multi-category ice fraction so use the thickness entry for now
          IF ( TRIM( sn_snd_thick%clcat ) == 'yes' ) ssnd(jps_fice)%nct = jpl
@@ -655,4 +700,41 @@
       CASE default   ;   CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_snd_thick%cldes' )
       END SELECT
+
+      !                                                      ! ------------------------- !
+      !                                                      ! Ice Meltponds             !
+      !                                                      ! ------------------------- !
+#if defined key_cice && ! defined key_cice4
+      ! Meltponds only CICE5 
+      ssnd(jps_a_p)%clname = 'OPndFrc'   
+      ssnd(jps_ht_p)%clname = 'OPndTck'   
+      SELECT CASE ( TRIM( sn_snd_mpnd%cldes ) )
+      CASE ( 'none' )
+         ssnd(jps_a_p)%laction = .FALSE.
+         ssnd(jps_ht_p)%laction = .FALSE.
+      CASE ( 'ice only' ) 
+         ssnd(jps_a_p)%laction = .TRUE.
+         ssnd(jps_ht_p)%laction = .TRUE.
+         IF ( TRIM( sn_snd_mpnd%clcat ) == 'yes' ) THEN
+            ssnd(jps_a_p)%nct = jpl
+            ssnd(jps_ht_p)%nct = jpl
+         ELSE
+            IF ( jpl > 1 ) THEN
+               CALL ctl_stop( 'sbc_cpl_init: use weighted ice option for sn_snd_mpnd%cldes if not exchanging category fields' )
+            ENDIF
+         ENDIF
+      CASE ( 'weighted ice' ) 
+         ssnd(jps_a_p)%laction = .TRUE.
+         ssnd(jps_ht_p)%laction = .TRUE.
+         IF ( TRIM( sn_snd_mpnd%clcat ) == 'yes' ) THEN
+            ssnd(jps_a_p)%nct = jpl 
+            ssnd(jps_ht_p)%nct = jpl 
+         ENDIF
+      CASE default   ;   CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_snd_mpnd%cldes' )
+      END SELECT
+#else
+      IF( TRIM( sn_snd_mpnd%cldes /= 'none' ) THEN
+         CALL ctl_stop('Meltponds can only be used with CICEv5')
+      ENDIF
+#endif
 
       !                                                      ! ------------------------- !
@@ -687,4 +769,36 @@
       !                                                      ! ------------------------- !
       ssnd(jps_co2)%clname = 'O_CO2FLX' ;  IF( TRIM(sn_snd_co2%cldes) == 'coupled' )    ssnd(jps_co2 )%laction = .TRUE.
+      !
+      
+      !                                                      ! ------------------------- !
+      !                                                      ! Sea surface freezing temp !
+      !                                                      ! ------------------------- !
+      ssnd(jps_sstfrz)%clname = 'O_SSTFrz' ; IF( TRIM(sn_snd_sstfrz%cldes) == 'coupled' )  ssnd(jps_sstfrz)%laction = .TRUE.
+      !
+      !                                                      ! ------------------------- !
+      !                                                      !    Ice conductivity       !
+      !                                                      ! ------------------------- !
+      ! Note that ultimately we will move to passing an ocean effective conductivity as well so there
+      ! will be some changes to the parts of the code which currently relate only to ice conductivity
+      ssnd(jps_kice )%clname = 'OIceKn'
+      SELECT CASE ( TRIM( sn_snd_cond%cldes ) )
+      CASE ( 'none' )
+         ssnd(jps_kice)%laction = .FALSE.
+      CASE ( 'ice only' )
+         ssnd(jps_kice)%laction = .TRUE.
+         IF ( TRIM( sn_snd_cond%clcat ) == 'yes' ) THEN
+            ssnd(jps_kice)%nct = jpl
+         ELSE
+            IF ( jpl > 1 ) THEN
+               CALL ctl_stop( 'sbc_cpl_init: use weighted ice option for sn_snd_cond%cldes if not exchanging category fields' )
+            ENDIF
+         ENDIF
+      CASE ( 'weighted ice' )
+         ssnd(jps_kice)%laction = .TRUE.
+         IF ( TRIM( sn_snd_cond%clcat ) == 'yes' ) ssnd(jps_kice)%nct = jpl
+      CASE default   ;   CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_snd_cond%cldes' )
+      END SELECT
+      !
+      
 
       !                                                      ! ------------------------------- !
@@ -783,5 +897,28 @@
       ncpl_qsr_freq = 86400 / ncpl_qsr_freq
 
-      CALL wrk_dealloc( jpi,jpj,   zacs, zaos )
+      IF( ln_coupled_iceshelf_fluxes ) THEN
+          ! Crude masks to separate the Antarctic and Greenland icesheets. Obviously something
+          ! more complicated could be done if required.
+          greenland_icesheet_mask = 0.0
+          WHERE( gphit >= 0.0 ) greenland_icesheet_mask = 1.0
+          antarctica_icesheet_mask = 0.0
+          WHERE( gphit < 0.0 ) antarctica_icesheet_mask = 1.0
+
+          ! initialise other variables
+          greenland_icesheet_mass_array(:,:) = 0.0
+          antarctica_icesheet_mass_array(:,:) = 0.0
+
+          IF( .not. ln_rstart ) THEN
+             greenland_icesheet_mass = 0.0 
+             greenland_icesheet_mass_rate_of_change = 0.0 
+             greenland_icesheet_timelapsed = 0.0
+             antarctica_icesheet_mass = 0.0 
+             antarctica_icesheet_mass_rate_of_change = 0.0 
+             antarctica_icesheet_timelapsed = 0.0
+          ENDIF
+
+      ENDIF
+
+      CALL wrk_dealloc( jpi,jpj, zacs, zaos )
       !
       IF( nn_timing == 1 )   CALL timing_stop('sbc_cpl_init')
@@ -841,7 +978,10 @@
       !!
       LOGICAL  ::   llnewtx, llnewtau      ! update wind stress components and module??
-      INTEGER  ::   ji, jj, jn             ! dummy loop indices
+      INTEGER  ::   ji, jj, jl, jn         ! dummy loop indices
       INTEGER  ::   isec                   ! number of seconds since nit000 (assuming rdt did not change since nit000)
       REAL(wp) ::   zcumulneg, zcumulpos   ! temporary scalars     
+      REAL(wp) ::   zgreenland_icesheet_mass_in, zantarctica_icesheet_mass_in
+      REAL(wp) ::   zgreenland_icesheet_mass_b, zantarctica_icesheet_mass_b
+      REAL(wp) ::   zmask_sum, zepsilon      
       REAL(wp) ::   zcoef                  ! temporary scalar
       REAL(wp) ::   zrhoa  = 1.22          ! Air density kg/m3
@@ -894,5 +1034,5 @@
                IF( srcv(jpr_otx2)%laction ) THEN
                   CALL rot_rep( frcv(jpr_otx2)%z3(:,:,1), frcv(jpr_oty2)%z3(:,:,1), srcv(jpr_otx2)%clgrid, 'en->j', zty )   
-               ELSE	
+               ELSE
                   CALL rot_rep( frcv(jpr_otx1)%z3(:,:,1), frcv(jpr_oty1)%z3(:,:,1), srcv(jpr_otx1)%clgrid, 'en->j', zty )  
                ENDIF
@@ -991,4 +1131,23 @@
 #endif
 
+#if defined key_cice && ! defined key_cice4
+      !  ! Sea ice surface skin temp:
+      IF( srcv(jpr_ts_ice)%laction ) THEN
+        DO jl = 1, jpl
+          DO jj = 1, jpj
+            DO ji = 1, jpi
+              IF (frcv(jpr_ts_ice)%z3(ji,jj,jl) > 0.0) THEN
+                tsfc_ice(ji,jj,jl) = 0.0
+              ELSE IF (frcv(jpr_ts_ice)%z3(ji,jj,jl) < -60.0) THEN
+                tsfc_ice(ji,jj,jl) = -60.0
+              ELSE
+                tsfc_ice(ji,jj,jl) = frcv(jpr_ts_ice)%z3(ji,jj,jl)
+              ENDIF
+            END DO
+          END DO
+        END DO
+      ENDIF
+#endif
+
       !  Fields received by SAS when OASIS coupling
       !  (arrays no more filled at sbcssm stage)
@@ -1105,4 +1264,63 @@
          !
       ENDIF
+      
+      !                                                        ! land ice masses : Greenland
+      zepsilon = rn_iceshelf_fluxes_tolerance
+
+      IF( srcv(jpr_grnm)%laction ) THEN
+         greenland_icesheet_mass_array(:,:) = frcv(jpr_grnm)%z3(:,:,1)
+         ! take average over ocean points of input array to avoid cumulative error over time
+         zgreenland_icesheet_mass_in = SUM( greenland_icesheet_mass_array(:,:) * tmask(:,:,1) )
+         IF(lk_mpp) CALL mpp_sum( zgreenland_icesheet_mass_in )
+         zmask_sum = SUM( tmask(:,:,1) )
+         IF(lk_mpp) CALL mpp_sum( zmask_sum ) 
+         zgreenland_icesheet_mass_in = zgreenland_icesheet_mass_in / zmask_sum
+         greenland_icesheet_timelapsed = greenland_icesheet_timelapsed + rdt         
+         IF( ABS( zgreenland_icesheet_mass_in - greenland_icesheet_mass ) > zepsilon ) THEN
+            zgreenland_icesheet_mass_b = greenland_icesheet_mass
+            
+            ! Only update the mass if it has increased
+            IF ( (zgreenland_icesheet_mass_in - greenland_icesheet_mass) > 0.0 ) THEN
+               greenland_icesheet_mass = zgreenland_icesheet_mass_in
+            ENDIF
+            
+            IF( zgreenland_icesheet_mass_b /= 0.0 ) &
+           &     greenland_icesheet_mass_rate_of_change = ( greenland_icesheet_mass - zgreenland_icesheet_mass_b ) / greenland_icesheet_timelapsed 
+            greenland_icesheet_timelapsed = 0.0_wp       
+         ENDIF
+         IF(lwp) WRITE(numout,*) 'Greenland icesheet mass (kg) read in is ', zgreenland_icesheet_mass_in
+         IF(lwp) WRITE(numout,*) 'Greenland icesheet mass (kg) used is    ', greenland_icesheet_mass
+         IF(lwp) WRITE(numout,*) 'Greenland icesheet mass rate of change (kg/s) is ', greenland_icesheet_mass_rate_of_change
+         IF(lwp) WRITE(numout,*) 'Greenland icesheet seconds lapsed since last change is ', greenland_icesheet_timelapsed
+      ENDIF
+
+      !                                                        ! land ice masses : Antarctica
+      IF( srcv(jpr_antm)%laction ) THEN
+         antarctica_icesheet_mass_array(:,:) = frcv(jpr_antm)%z3(:,:,1)
+         ! take average over ocean points of input array to avoid cumulative error from rounding errors over time
+         zantarctica_icesheet_mass_in = SUM( antarctica_icesheet_mass_array(:,:) * tmask(:,:,1) )
+         IF(lk_mpp) CALL mpp_sum( zantarctica_icesheet_mass_in )
+         zmask_sum = SUM( tmask(:,:,1) )
+         IF(lk_mpp) CALL mpp_sum( zmask_sum ) 
+         zantarctica_icesheet_mass_in = zantarctica_icesheet_mass_in / zmask_sum
+         antarctica_icesheet_timelapsed = antarctica_icesheet_timelapsed + rdt         
+         IF( ABS( zantarctica_icesheet_mass_in - antarctica_icesheet_mass ) > zepsilon ) THEN
+            zantarctica_icesheet_mass_b = antarctica_icesheet_mass
+            
+            ! Only update the mass if it has increased
+            IF ( (zantarctica_icesheet_mass_in - antarctica_icesheet_mass) > 0.0 ) THEN
+               antarctica_icesheet_mass = zantarctica_icesheet_mass_in
+            END IF
+            
+            IF( zantarctica_icesheet_mass_b /= 0.0 ) &
+          &      antarctica_icesheet_mass_rate_of_change = ( antarctica_icesheet_mass - zantarctica_icesheet_mass_b ) / antarctica_icesheet_timelapsed 
+            antarctica_icesheet_timelapsed = 0.0_wp       
+         ENDIF
+         IF(lwp) WRITE(numout,*) 'Antarctica icesheet mass (kg) read in is ', zantarctica_icesheet_mass_in
+         IF(lwp) WRITE(numout,*) 'Antarctica icesheet mass (kg) used is    ', antarctica_icesheet_mass
+         IF(lwp) WRITE(numout,*) 'Antarctica icesheet mass rate of change (kg/s) is ', antarctica_icesheet_mass_rate_of_change
+         IF(lwp) WRITE(numout,*) 'Antarctica icesheet seconds lapsed since last change is ', antarctica_icesheet_timelapsed
+      ENDIF
+
       !
       CALL wrk_dealloc( jpi,jpj,   ztx, zty, zmsk, zemp, zqns, zqsr )
@@ -1400,5 +1618,32 @@
          ztprecip(:,:) = frcv(jpr_rain)%z3(:,:,1) + zsprecip(:,:)  ! May need to ensure positive here
          zemp_tot(:,:) = frcv(jpr_tevp)%z3(:,:,1) - ztprecip(:,:)
+#if defined key_cice
+         IF ( TRIM(sn_rcv_emp%clcat) == 'yes' ) THEN
+            ! zemp_ice is the sum of frcv(jpr_ievp)%z3(:,:,1) over all layers - snow
+            zemp_ice(:,:) = - frcv(jpr_snow)%z3(:,:,1)
+            DO jl=1,jpl
+               zemp_ice(:,:   ) = zemp_ice(:,:) + frcv(jpr_ievp)%z3(:,:,jl)
+            ENDDO
+            ! latent heat coupled for each category in CICE
+            qla_ice(:,:,1:jpl) = - frcv(jpr_ievp)%z3(:,:,1:jpl) * lsub
+         ELSE
+            ! If CICE has multicategories it still expects coupling fields for
+            ! each even if we treat as a single field
+            ! The latent heat flux is split between the ice categories according
+            ! to the fraction of the ice in each category
+            zemp_ice(:,:) = frcv(jpr_ievp)%z3(:,:,1) - frcv(jpr_snow)%z3(:,:,1)
+            WHERE ( zicefr(:,:) /= 0._wp ) 
+               ztmp(:,:) = 1./zicefr(:,:)
+            ELSEWHERE 
+               ztmp(:,:) = 0.e0
+            END WHERE  
+            DO jl=1,jpl
+               qla_ice(:,:,jl) = - a_i(:,:,jl) * ztmp(:,:) * frcv(jpr_ievp)%z3(:,:,1) * lsub 
+            END DO
+            WHERE ( zicefr(:,:) == 0._wp )  qla_ice(:,:,1) = -frcv(jpr_ievp)%z3(:,:,1) * lsub 
+         ENDIF
+#else         
          zemp_ice(:,:) = frcv(jpr_ievp)%z3(:,:,1) - frcv(jpr_snow)%z3(:,:,1)
+#endif                  
             CALL iom_put( 'rain'         , frcv(jpr_rain)%z3(:,:,1)              )   ! liquid precipitation 
          IF( iom_use('hflx_rain_cea') )   &
@@ -1795,4 +2040,15 @@
                CASE default                  ;   CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_temp%clcat' )
                END SELECT
+            CASE( 'oce and weighted ice' )   ;   ztmp1(:,:) =   tsn(:,:,1,jp_tem) + rt0 
+               SELECT CASE( sn_snd_temp%clcat )
+               CASE( 'yes' )   
+                  ztmp3(:,:,1:jpl) = tn_ice(:,:,1:jpl) * a_i(:,:,1:jpl)
+               CASE( 'no' )
+                  ztmp3(:,:,:) = 0.0
+                  DO jl=1,jpl
+                     ztmp3(:,:,1) = ztmp3(:,:,1) + tn_ice(:,:,jl) * a_i(:,:,jl)
+                  ENDDO
+               CASE default                  ;   CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_temp%clcat' )
+               END SELECT
             CASE( 'mixed oce-ice'        )   
                ztmp1(:,:) = ( ztmp1(:,:) + rt0 ) * zfr_l(:,:) 
@@ -1865,4 +2121,14 @@
          END SELECT
          IF( ssnd(jps_fice)%laction )   CALL cpl_snd( jps_fice, isec, ztmp3, info )
+      ENDIF
+      
+      ! Send ice fraction field (first order interpolation), for weighting UM fluxes to be passed to NEMO
+      IF (ssnd(jps_fice1)%laction) THEN
+         SELECT CASE (sn_snd_thick1%clcat)
+         CASE( 'yes' )   ;   ztmp3(:,:,1:jpl) =  a_i(:,:,1:jpl)
+         CASE( 'no' )    ;   ztmp3(:,:,1) = fr_i(:,:)
+         CASE default    ;   CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_thick1%clcat' )
+      END SELECT
+         CALL cpl_snd (jps_fice1, isec, ztmp3, info)
       ENDIF
       
@@ -1910,4 +2176,48 @@
          IF( ssnd(jps_hsnw)%laction )   CALL cpl_snd( jps_hsnw, isec, ztmp4, info )
       ENDIF
+      !
+      ! Send meltpond fields 
+      IF( ssnd(jps_a_p)%laction .OR. ssnd(jps_ht_p)%laction ) THEN
+         SELECT CASE( sn_snd_mpnd%cldes) 
+         CASE( 'weighted ice' ) 
+            SELECT CASE( sn_snd_mpnd%clcat ) 
+            CASE( 'yes' ) 
+               ztmp3(:,:,1:jpl) =  a_p(:,:,1:jpl) * a_i(:,:,1:jpl) 
+               ztmp4(:,:,1:jpl) =  ht_p(:,:,1:jpl) * a_i(:,:,1:jpl) 
+            CASE( 'no' ) 
+               ztmp3(:,:,:) = 0.0 
+               ztmp4(:,:,:) = 0.0 
+               DO jl=1,jpl 
+                 ztmp3(:,:,1) = ztmp3(:,:,1) + a_p(:,:,jpl) * a_i(:,:,jpl) 
+                 ztmp4(:,:,1) = ztmp4(:,:,1) + ht_p(:,:,jpl) * a_i(:,:,jpl) 
+               ENDDO 
+            CASE default    ;   CALL ctl_stop( 'sbc_cpl_mpd: wrong definition of sn_snd_mpnd%clcat' ) 
+            END SELECT 
+         CASE( 'ice only' )    
+            ztmp3(:,:,1:jpl) = a_p(:,:,1:jpl) 
+            ztmp4(:,:,1:jpl) = ht_p(:,:,1:jpl) 
+         END SELECT 
+         IF( ssnd(jps_a_p)%laction )   CALL cpl_snd( jps_a_p, isec, ztmp3, info )    
+         IF( ssnd(jps_ht_p)%laction )   CALL cpl_snd( jps_ht_p, isec, ztmp4, info )    
+         !
+         ! Send ice effective conductivity
+         SELECT CASE( sn_snd_cond%cldes)
+         CASE( 'weighted ice' )   
+            SELECT CASE( sn_snd_cond%clcat )
+            CASE( 'yes' )   
+               ztmp3(:,:,1:jpl) =  kn_ice(:,:,1:jpl) * a_i(:,:,1:jpl)
+            CASE( 'no' )
+               ztmp3(:,:,:) = 0.0
+               DO jl=1,jpl
+                 ztmp3(:,:,1) = ztmp3(:,:,1) + kn_ice(:,:,jl) * a_i(:,:,jl)
+               ENDDO
+            CASE default                  ;   CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_cond%clcat' )
+            END SELECT
+         CASE( 'ice only' )   
+           ztmp3(:,:,1:jpl) = kn_ice(:,:,1:jpl)
+         END SELECT
+         IF( ssnd(jps_kice)%laction )   CALL cpl_snd( jps_kice, isec, ztmp3, info )
+      ENDIF
+      !
       !
 #if defined key_cpl_carbon_cycle
@@ -2089,5 +2399,8 @@
       IF( ssnd(jps_rnf   )%laction )  CALL cpl_snd( jps_rnf   , isec, RESHAPE ( rnf , (/jpi,jpj,1/) ), info )
       IF( ssnd(jps_taum  )%laction )  CALL cpl_snd( jps_taum  , isec, RESHAPE ( taum, (/jpi,jpj,1/) ), info )
-
+      
+      ztmp1(:,:) = sstfrz(:,:) + rt0
+      IF( ssnd(jps_sstfrz)%laction )  CALL cpl_snd( jps_sstfrz, isec, RESHAPE ( ztmp1, (/jpi,jpj,1/) ), info )
+      !
       CALL wrk_dealloc( jpi,jpj,       zfr_l, ztmp1, ztmp2, zotx1, zoty1, zotz1, zitx1, zity1, zitz1 )
       CALL wrk_dealloc( jpi,jpj,jpl,   ztmp3, ztmp4 )
