Index: /branches/UKMO/dev_r5518_GSI7_GSI8_landice_bitcomp_medusa/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90
===================================================================
--- /branches/UKMO/dev_r5518_GSI7_GSI8_landice_bitcomp_medusa/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90	(revision 6667)
+++ /branches/UKMO/dev_r5518_GSI7_GSI8_landice_bitcomp_medusa/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90	(revision 6668)
@@ -50,4 +50,9 @@
 #endif
    USE lib_fortran, ONLY: glob_sum
+#if defined key_medusa
+!   USE trc
+   USE trcbio_medusa, ONLY: f_co2flux2d, dms_surf2d
+   USE sms_medusa, ONLY: f_pco2a, dust
+#endif
 
    IMPLICIT NONE
@@ -106,5 +111,7 @@
    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 ::   jpr_atm_pco2 = 46          ! Incoming atm CO2 flux
+   INTEGER, PARAMETER ::   jpr_atm_dust = 47          ! Incoming atm aggregate dust 
+   INTEGER, PARAMETER ::   jprcv      = 47            ! total number of fields received
 
    INTEGER, PARAMETER ::   jps_fice   =  1            ! ice fraction sent to the atmosphere
@@ -142,4 +149,10 @@
    INTEGER, PARAMETER ::   jps_fice1  = 33            ! first-order ice concentration (for time-travelling ice coupling)
    INTEGER, PARAMETER ::   jpsnd      = 33            ! total number of fields sended
+   INTEGER, PARAMETER ::   jps_bio_co2 = 34           ! MEDUSA air-sea CO2 flux in
+   INTEGER, PARAMETER ::   jps_bio_dms = 35           ! MEDUSA DMS surface concentration in
+   INTEGER, PARAMETER ::   jpsnd      = 35            ! total number of fields sent
+
+   REAL(wp), PARAMETER :: dms_unit_conv = 1.0e+6      ! Coversion factor to get outgong DMS in standard units for coupling
+	                                              ! i.e. specifically nmol/L (= umol/m3)
 
    !                                                         !!** namelist namsbc_cpl **
@@ -153,8 +166,11 @@
    ! Send to the atmosphere                           !
    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
+   TYPE(FLD_C) ::   sn_snd_bio_co2, sn_snd_bio_dms                       
 
    ! 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, sn_rcv_ts_ice, sn_rcv_grnm, sn_rcv_antm
+   TYPE(FLD_C) ::   sn_rcv_atm_pco2, sn_rcv_atm_dust                         
+
    ! Other namelist parameters                        !
    INTEGER     ::   nn_cplmodel            ! Maximum number of models to/from which NEMO is potentialy sending/receiving data
@@ -231,4 +247,11 @@
          &                  rn_antarctica_calving_fraction, rn_iceshelf_fluxes_tolerance
       !!---------------------------------------------------------------------
+
+      ! Add MEDUSA related fields to namelist
+      NAMELIST/namsbc_cpl/  sn_snd_bio_co2, sn_snd_bio_dms,                                           &
+         &                  sn_rcv_atm_pco2, sn_rcv_atm_dust
+
+      !!---------------------------------------------------------------------
+
       !
       IF( nn_timing == 1 )  CALL timing_start('sbc_cpl_init')
@@ -255,5 +278,5 @@
       ENDIF
       IF( lwp .AND. ln_cpl ) THEN                        ! control print
-         WRITE(numout,*)'  received fields (mutiple ice categogies)'
+         WRITE(numout,*)'  received fields (mutiple ice categories)'
          WRITE(numout,*)'      10m wind module                 = ', TRIM(sn_rcv_w10m%cldes  ), ' (', TRIM(sn_rcv_w10m%clcat  ), ')'
          WRITE(numout,*)'      stress module                   = ', TRIM(sn_rcv_taumod%cldes), ' (', TRIM(sn_rcv_taumod%clcat), ')'
@@ -272,4 +295,6 @@
          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   ), ')'
+         WRITE(numout,*)'      atm pco2                        = ', TRIM(sn_rcv_atm_pco2%cldes), ' (', TRIM(sn_rcv_atm_pco2%clcat), ')'
+         WRITE(numout,*)'      atm dust                        = ', TRIM(sn_rcv_atm_dust%cldes), ' (', TRIM(sn_rcv_atm_dust%clcat), ')'
          WRITE(numout,*)'  sent fields (multiple ice categories)'
          WRITE(numout,*)'      surface temperature             = ', TRIM(sn_snd_temp%cldes  ), ' (', TRIM(sn_snd_temp%clcat  ), ')'
@@ -280,4 +305,6 @@
          WRITE(numout,*)'                      - orientation   = ', sn_snd_crt%clvor
          WRITE(numout,*)'                      - mesh          = ', sn_snd_crt%clvgrd
+         WRITE(numout,*)'      bio co2 flux                    = ', TRIM(sn_snd_bio_co2%cldes), ' (', TRIM(sn_snd_bio_co2%clcat), ')'
+         WRITE(numout,*)'      bio dms flux                    = ', TRIM(sn_snd_bio_dms%cldes), ' (', TRIM(sn_snd_bio_dms%clcat), ')'
          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   ), ')'
@@ -495,4 +522,20 @@
       !                                                      ! ------------------------- !
       srcv(jpr_co2 )%clname = 'O_AtmCO2'   ;   IF( TRIM(sn_rcv_co2%cldes   ) == 'coupled' )    srcv(jpr_co2 )%laction = .TRUE.
+
+
+      !                                                      ! --------------------------------------- !    
+      !                                                      ! Incoming CO2 and DUST fluxes for MEDUSA !
+      !                                                      ! --------------------------------------- !  
+      srcv(jpr_atm_pco2)%clname = 'OATMPCO2'
+
+      IF (TRIM(sn_rcv_atm_pco2%cldes) == 'coupled') THEN
+        srcv(jpr_atm_pco2)%laction = .TRUE.
+      END IF
+               
+      srcv(jpr_atm_dust)%clname = 'OATMDUST'   
+      IF (TRIM(sn_rcv_atm_dust%cldes) == 'coupled')  THEN
+        srcv(jpr_atm_dust)%laction = .TRUE.
+      END IF
+	 
       !                                                      ! ------------------------- !
       !                                                      !   topmelt and botmelt     !   
@@ -991,4 +1034,10 @@
       REAL(wp), POINTER, DIMENSION(:,:) ::   ztx, zty, zmsk, zemp, zqns, zqsr
       !!----------------------------------------------------------------------
+
+      ! RSRH temporary arrays for testing, just to recieve incoming MEDUSA related fields
+      ! until we know where they need to go.
+      REAL(wp), ALLOCATABLE :: atm_pco2(:,:)
+      REAL(wp), ALLOCATABLE :: atm_dust(:,:)
+
       !
       IF( nn_timing == 1 )  CALL timing_start('sbc_cpl_rcv')
@@ -1128,4 +1177,17 @@
          !  
       ENDIF
+
+#if defined key_medusa
+      ! RSRH Allocate temporary arrays to receive incoming fields during testing
+      ALLOCATE(atm_pco2(jpi,jpj))
+      ALLOCATE(atm_dust(jpi,jpj))
+
+      IF( srcv(jpr_atm_pco2)%laction) atm_pco2(:,:) = frcv(jpr_atm_pco2)%z3(:,:,1)
+      IF( srcv(jpr_atm_dust)%laction) atm_dust(:,:) = frcv(jpr_atm_dust)%z3(:,:,1)
+     
+      ! RSRH Deallocate temporary arrays.
+      DEALLOCATE(atm_pco2)
+      DEALLOCATE(atm_dust)
+#endif
 
 #if defined key_cpl_carbon_cycle
@@ -2168,4 +2230,18 @@
       !
 #endif
+
+
+
+#if defined key_medusa
+      IF( ssnd(jps_bio_co2)%laction ) CALL cpl_prism_snd( jps_bio_co2, isec, f_co2flux2d, info )
+
+      IF( ssnd(jps_bio_dms)%laction )  THEN
+          ! We need to multiply DMS by a conversion factor to get values in the standard units expected in
+          ! the coupling space.
+          ztmp1(:,: ) = dms_surf2d(:,:) * dms_unit_conv
+         CALL cpl_prism_snd( jps_bio_dms, isec, RESHAPE ( ztmp1, (/jpi,jpj,1/) ), info )
+      ENDIF
+#endif
+
       !                                                      ! ------------------------- !
       IF( ssnd(jps_ocx1)%laction ) THEN                      !      Surface current      !
