Index: branches/UKMO/dev_1d_bugfixes_tocommit/NEMOGCM/NEMO/TOP_SRC/MEDUSA/bio_medusa_init.F90
===================================================================
--- branches/UKMO/dev_1d_bugfixes_tocommit/NEMOGCM/NEMO/TOP_SRC/MEDUSA/bio_medusa_init.F90	(revision 13190)
+++ branches/UKMO/dev_1d_bugfixes_tocommit/NEMOGCM/NEMO/TOP_SRC/MEDUSA/bio_medusa_init.F90	(revision 13354)
@@ -35,5 +35,6 @@
       USE bio_medusa_mod
       USE par_oce,           ONLY: jpi, jpj, jpk
-      USE sms_medusa,        ONLY: jdms, pgrow_avg, ploss_avg, phyt_avg, mld_max
+      USE sms_medusa,        ONLY: jdms, pgrow_avg, ploss_avg, phyt_avg, mld_max, &
+         &                         pgrow_avg_3d, ploss_avg_3d, phyt_avg_3d
       USE trc,               ONLY: ln_diatrc, med_diag, nittrc000 
       USE in_out_manager,    ONLY: lwp, numout
@@ -199,4 +200,7 @@
       ploss_avg(:,:) = 0.0
       phyt_avg(:,:)  = 0.0
+      pgrow_avg_3d(:,:,:) = 0.0
+      ploss_avg_3d(:,:,:) = 0.0
+      phyt_avg_3d(:,:,:)  = 0.0
       IF( kt == nittrc000 ) THEN
          mld_max(:,:) = 0.0
Index: branches/UKMO/dev_1d_bugfixes_tocommit/NEMOGCM/NEMO/TOP_SRC/MEDUSA/plankton.F90
===================================================================
--- branches/UKMO/dev_1d_bugfixes_tocommit/NEMOGCM/NEMO/TOP_SRC/MEDUSA/plankton.F90	(revision 13190)
+++ branches/UKMO/dev_1d_bugfixes_tocommit/NEMOGCM/NEMO/TOP_SRC/MEDUSA/plankton.F90	(revision 13354)
@@ -46,4 +46,5 @@
                                    ln_foam_medusa,                         &
                                    pgrow_avg, ploss_avg, phyt_avg,         &
+                                   pgrow_avg_3d, ploss_avg_3d, phyt_avg_3d, &
                                    xkphd, xkphn, xkzme, xkzmi,             &
                                    xmetapd, xmetapn, xmetazme, xmetazmi,   &
@@ -229,4 +230,12 @@
                                      ((zphn(ji,jj) + zphd(ji,jj)) *       &
                                       fse3t(ji,jj,jk) * fq0)
+                  !!
+                  pgrow_avg_3d(ji,jj,jk) = (fprn(ji,jj) * zphn(ji,jj)) +  &
+                                           (fprd(ji,jj) * zphd(ji,jj))
+                  ploss_avg_3d(ji,jj,jk) = fgmepd(ji,jj) + fdpd(ji,jj) +  &
+                                           fdpd2(ji,jj)                +  &
+                                           fgmepn(ji,jj) + fdpn(ji,jj) +  &
+                                           fdpn2(ji,jj)  + fgmipn(ji,jj)
+                  phyt_avg_3d(ji,jj,jk)  = zphn(ji,jj) + zphd(ji,jj)
                ENDIF
             ENDDO
Index: branches/UKMO/dev_1d_bugfixes_tocommit/NEMOGCM/NEMO/TOP_SRC/MEDUSA/sms_medusa.F90
===================================================================
--- branches/UKMO/dev_1d_bugfixes_tocommit/NEMOGCM/NEMO/TOP_SRC/MEDUSA/sms_medusa.F90	(revision 13190)
+++ branches/UKMO/dev_1d_bugfixes_tocommit/NEMOGCM/NEMO/TOP_SRC/MEDUSA/sms_medusa.F90	(revision 13354)
@@ -362,4 +362,7 @@
    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: ploss_avg  !: Mixed layer average phytoplankton loss
    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: phyt_avg   !: Mixed layer average phytoplankton
+   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: pgrow_avg_3d  !: Mixed layer average phytoplankton growth
+   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: ploss_avg_3d  !: Mixed layer average phytoplankton loss
+   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: phyt_avg_3d   !: Mixed layer average phytoplankton
    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: mld_max    !: Maximum mixed layer depth
 !!
@@ -438,4 +441,6 @@
       !* Fields for ocean colour data assimilation
       ALLOCATE( pgrow_avg(jpi,jpj)   , ploss_avg(jpi,jpj)   ,       &
+         &      pgrow_avg_3d(jpi,jpj,jpk) , ploss_avg_3d(jpi,jpj,jpk) , &
+         &      phyt_avg_3d(jpi,jpj,jpk)  ,                         &
          &      phyt_avg(jpi,jpj)    , mld_max(jpi,jpj)     ,    STAT=ierr(9) )
 #endif
Index: branches/UKMO/dev_1d_bugfixes_tocommit/NEMOGCM/NEMO/TOP_SRC/trcrst.F90
===================================================================
--- branches/UKMO/dev_1d_bugfixes_tocommit/NEMOGCM/NEMO/TOP_SRC/trcrst.F90	(revision 13190)
+++ branches/UKMO/dev_1d_bugfixes_tocommit/NEMOGCM/NEMO/TOP_SRC/trcrst.F90	(revision 13354)
@@ -368,4 +368,15 @@
             mld_max(:,:)   = 0.0
          ENDIF
+         IF( iom_varid( numrtr, 'pgrow_avg_3d', ldstop = .FALSE. ) > 0 ) THEN
+            IF(lwp) WRITE(numout,*) ' MEDUSA pgrow_avg present - reading in ...'
+            CALL iom_get( numrtr, jpdom_autoglo, 'pgrow_avg_3d',  pgrow_avg_3d(:,:,:)  )
+            CALL iom_get( numrtr, jpdom_autoglo, 'ploss_avg_3d',  ploss_avg_3d(:,:,:)  )
+            CALL iom_get( numrtr, jpdom_autoglo, 'phyt_avg_3d',   phyt_avg_3d(:,:,:)   )
+         ELSE
+            IF(lwp) WRITE(numout,*) ' MEDUSA pgrow_avg_3d absent - setting to zero ...'
+            pgrow_avg_3d(:,:,:) = 0.0
+            ploss_avg_3d(:,:,:) = 0.0
+            phyt_avg_3d(:,:,:)  = 0.0
+         ENDIF
       ENDIF
 
@@ -553,4 +564,7 @@
          CALL iom_rstput( kt, nitrst, numrtw, 'phyt_avg',  phyt_avg(:,:)  )
          CALL iom_rstput( kt, nitrst, numrtw, 'mld_max',   mld_max(:,:)   )
+         CALL iom_rstput( kt, nitrst, numrtw, 'pgrow_avg_3d', pgrow_avg_3d(:,:,:) )
+         CALL iom_rstput( kt, nitrst, numrtw, 'ploss_avg_3d', ploss_avg_3d(:,:,:) )
+         CALL iom_rstput( kt, nitrst, numrtw, 'phyt_avg_3d',  phyt_avg_3d(:,:,:)  )
       ENDIF
 !!
