Index: /branches/UKMO/dev_r5518_v3.4_asm_nemovar_community_bgc_ersem/NEMOGCM/NEMO/OPA_SRC/ASM/asmbal.F90
===================================================================
--- /branches/UKMO/dev_r5518_v3.4_asm_nemovar_community_bgc_ersem/NEMOGCM/NEMO/OPA_SRC/ASM/asmbal.F90	(revision 6983)
+++ /branches/UKMO/dev_r5518_v3.4_asm_nemovar_community_bgc_ersem/NEMOGCM/NEMO/OPA_SRC/ASM/asmbal.F90	(revision 6983)
@@ -0,0 +1,136 @@
+MODULE asmbal
+   !!======================================================================
+   !!                       ***  MODULE asmbal  ***
+   !! Assimilation balancing interface: Write to file the balancing increments
+   !!                                   calculated for biogeochemistry
+   !!======================================================================
+   !!----------------------------------------------------------------------
+   !!   'key_asminc' : Switch on the assimilation increment interface
+   !!----------------------------------------------------------------------
+   !!   asm_bal_wri  : Write out the background state
+   !!----------------------------------------------------------------------
+   !! * Modules used
+   USE par_kind, ONLY: &   ! Precision variables
+      & wp
+   USE iom                 ! I/O module
+   USE asminc              ! Main assimilation increments module
+   USE asmpar              ! Parameters for the assimilation interface
+
+   IMPLICIT NONE
+
+   !! * Routine accessibility
+   PRIVATE
+   PUBLIC asm_bal_wri   !: Write out the background state
+
+CONTAINS
+
+   SUBROUTINE asm_bal_wri( kt )
+      !!-----------------------------------------------------------------------
+      !!
+      !!                  ***  ROUTINE asm_bal_wri ***
+      !!
+      !! ** Purpose : Write to file the balancing increments
+      !!              calculated for biogeochemistry
+      !!
+      !! ** Method  : Write to file the balancing increments
+      !!              calculated for biogeochemistry
+      !!
+      !! ** Action  :
+      !!                   
+      !! References :
+      !!
+      !! History    :
+      !!        ! 2014-08 (D. Ford)  Initial version, based on asm_bkg_wri
+      !!-----------------------------------------------------------------------
+      !! * Arguments
+      INTEGER, INTENT( IN ) :: kt        ! Current time-step
+      !! * Local declarations
+      CHARACTER(LEN=50) :: cl_asmbal     ! Filename (with extension)
+      LOGICAL           :: llok          ! Check if file exists
+      INTEGER           :: inum          ! File unit number
+      REAL(wp)          :: zdate         ! Date
+      !!-----------------------------------------------------------------------
+     
+      ! Set things up
+      zdate = REAL( ndastp )
+      WRITE(cl_asmbal, FMT='(A,".nc")' ) TRIM( c_asmbal )
+
+      ! Check if file exists
+      INQUIRE( FILE = TRIM( cl_asmbal ), EXIST = llok )
+      IF( .NOT. llok ) THEN
+         IF(lwp) WRITE(numout,*) ' Setting up assimilation balancing increments file '// &
+            &                    TRIM( c_asmbal ) // ' at timestep = ', kt
+
+         ! Define the output file       
+         CALL iom_open( c_asmbal, inum, ldwrt = .TRUE., kiolib = jprstlib)
+
+         ! Write the information
+         CALL iom_rstput( kt, kt, inum, 'rdastp' , zdate   )
+
+         IF ( ln_logchltotinc .OR. ln_logchlpftinc ) THEN
+#if defined key_fabm
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_ersem_chl1', logchl_balinc_ersem_chl1 )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_ersem_chl2', logchl_balinc_ersem_chl2 )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_ersem_chl3', logchl_balinc_ersem_chl3 )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_ersem_chl4', logchl_balinc_ersem_chl4 )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_ersem_p1c',  logchl_balinc_ersem_p1c  )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_ersem_p1n',  logchl_balinc_ersem_p1n  )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_ersem_p1p',  logchl_balinc_ersem_p1p  )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_ersem_p1s',  logchl_balinc_ersem_p1s  )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_ersem_p2c',  logchl_balinc_ersem_p2c  )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_ersem_p2n',  logchl_balinc_ersem_p2n  )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_ersem_p2p',  logchl_balinc_ersem_p2p  )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_ersem_p3c',  logchl_balinc_ersem_p3c  )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_ersem_p3n',  logchl_balinc_ersem_p3n  )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_ersem_p3p',  logchl_balinc_ersem_p3p  )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_ersem_p4c',  logchl_balinc_ersem_p4c  )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_ersem_p4n',  logchl_balinc_ersem_p4n  )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_ersem_p4p',  logchl_balinc_ersem_p4p  )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_ersem_z4c',  logchl_balinc_ersem_z4c  )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_ersem_z5c',  logchl_balinc_ersem_z5c  )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_ersem_z5n',  logchl_balinc_ersem_z5n  )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_ersem_z5p',  logchl_balinc_ersem_z5p  )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_ersem_z6c',  logchl_balinc_ersem_z6c  )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_ersem_z6n',  logchl_balinc_ersem_z6n  )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_ersem_z6p',  logchl_balinc_ersem_z6p  )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_ersem_n1p',  logchl_balinc_ersem_n1p  )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_ersem_n3n',  logchl_balinc_ersem_n3n  )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_ersem_n4n',  logchl_balinc_ersem_n4n  )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_ersem_n5s',  logchl_balinc_ersem_n5s  )
+#elif defined key_medusa && defined key_foam_medusa
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_medusa_chn', logchl_balinc_medusa_chn )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_medusa_chd', logchl_balinc_medusa_chd )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_medusa_phn', logchl_balinc_medusa_phn )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_medusa_phd', logchl_balinc_medusa_phd )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_medusa_pds', logchl_balinc_medusa_pds )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_medusa_zmi', logchl_balinc_medusa_zmi )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_medusa_zme', logchl_balinc_medusa_zme )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_medusa_din', logchl_balinc_medusa_din )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_medusa_sil', logchl_balinc_medusa_sil )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_medusa_fer', logchl_balinc_medusa_fer )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_medusa_det', logchl_balinc_medusa_det )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_medusa_dtc', logchl_balinc_medusa_dtc )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_medusa_dic', logchl_balinc_medusa_dic )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_medusa_alk', logchl_balinc_medusa_alk )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_medusa_oxy', logchl_balinc_medusa_oxy )
+#elif defined key_hadocc
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_hadocc_nut', logchl_balinc_hadocc_nut )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_hadocc_phy', logchl_balinc_hadocc_phy )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_hadocc_zoo', logchl_balinc_hadocc_zoo )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_hadocc_det', logchl_balinc_hadocc_det )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_hadocc_dic', logchl_balinc_hadocc_dic )
+            CALL iom_rstput( kt, kt, inum, 'logchl_balinc_hadocc_alk', logchl_balinc_hadocc_alk )
+#endif
+         ENDIF
+
+         CALL iom_close( inum )
+      ELSE
+         CALL ctl_warn( TRIM( cl_asmbal ) // ' already exists ', &
+            &           ' Therefore not writing out balancing increments at this timestep', &
+            &           ' Something has probably gone wrong somewhere' )
+         IF(lwp) WRITE(numout,*) ' Failed to set up assimilation balancing increments file '// &
+            &                    TRIM( c_asmbal ) // ' at timestep = ', kt
+      ENDIF
+                                 
+   END SUBROUTINE asm_bal_wri
+END MODULE asmbal
Index: /branches/UKMO/dev_r5518_v3.4_asm_nemovar_community_bgc_ersem/NEMOGCM/NEMO/OPA_SRC/ASM/asminc.F90
===================================================================
--- /branches/UKMO/dev_r5518_v3.4_asm_nemovar_community_bgc_ersem/NEMOGCM/NEMO/OPA_SRC/ASM/asminc.F90	(revision 6982)
+++ /branches/UKMO/dev_r5518_v3.4_asm_nemovar_community_bgc_ersem/NEMOGCM/NEMO/OPA_SRC/ASM/asminc.F90	(revision 6983)
@@ -22,4 +22,5 @@
    !!   ssh_asm_inc    : Apply the SSH increment
    !!   seaice_asm_inc : Apply the seaice increment
+   !!   logchl_asm_inc : Apply the logchl increment
    !!----------------------------------------------------------------------
    USE wrk_nemo         ! Memory Allocation
@@ -52,4 +53,26 @@
    USE bdy_oce, ONLY: bdytmask  
 #endif  
+#if defined key_fabm
+   USE asmlogchlbal_ersem, ONLY: &
+      & asm_logchl_bal_ersem
+   USE trc, ONLY: &
+      & trn,      &
+      & trb
+   USE par_fabm
+#elif defined key_medusa && defined key_foam_medusa
+   USE asmlogchlbal_medusa, ONLY: &
+      & asm_logchl_bal_medusa
+   USE trc, ONLY: &
+      & trn,      &
+      & trb
+   USE par_medusa
+#elif defined key_hadocc
+   USE asmlogchlbal_hadocc, ONLY: &
+      & asm_logchl_bal_hadocc
+   USE trc, ONLY: &
+      & trn,      &
+      & trb
+   USE par_hadocc
+#endif
 
    IMPLICIT NONE
@@ -62,4 +85,5 @@
    PUBLIC   ssh_asm_inc    !: Apply the SSH increment
    PUBLIC   seaice_asm_inc !: Apply the seaice increment
+   PUBLIC   logchl_asm_inc !: Apply the seaice increment
 
 #if defined key_asminc
@@ -69,4 +93,5 @@
 #endif
    LOGICAL, PUBLIC :: ln_bkgwri = .FALSE.      !: No output of the background state fields
+   LOGICAL, PUBLIC :: ln_balwri = .FALSE.      !: No output of the assimilation balancing increments
    LOGICAL, PUBLIC :: ln_asmiau = .FALSE.      !: No applying forcing with an assimilation increment
    LOGICAL, PUBLIC :: ln_asmdin = .FALSE.      !: No direct initialization
@@ -74,5 +99,7 @@
    LOGICAL, PUBLIC :: ln_dyninc = .FALSE.      !: No dynamics (u and v) assimilation increments
    LOGICAL, PUBLIC :: ln_sshinc = .FALSE.      !: No sea surface height assimilation increment
-   LOGICAL, PUBLIC :: ln_seaiceinc             !: No sea ice concentration increment
+   LOGICAL, PUBLIC :: ln_seaiceinc = .FALSE.   !: No sea ice concentration increment
+   LOGICAL, PUBLIC :: ln_logchltotinc = .FALSE. !: No total log10(chlorophyll) increment
+   LOGICAL, PUBLIC :: ln_logchlpftinc = .FALSE. !: No PFT   log10(chlorophyll) increment
    LOGICAL, PUBLIC :: ln_salfix = .FALSE.      !: Apply minimum salinity check
    LOGICAL, PUBLIC :: ln_temnofreeze = .FALSE. !: Don't allow the temperature to drop below freezing
@@ -100,4 +127,59 @@
    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::   seaice_bkginc         ! Increment to the background sea ice conc
 
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_bkginc !: Increment to background logchl
+#if defined key_fabm
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_ersem_chl1  !: Increment to ERSEM diatom chl   from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_ersem_chl2  !: Increment to ERSEM nanophy chl  from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_ersem_chl3  !: Increment to ERSEM picophy chl  from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_ersem_chl4  !: Increment to ERSEM microphy chl from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_ersem_p1c   !: Increment to ERSEM diatom c     from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_ersem_p1n   !: Increment to ERSEM diatom n     from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_ersem_p1p   !: Increment to ERSEM diatom p     from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_ersem_p1s   !: Increment to ERSEM diatom s     from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_ersem_p2c   !: Increment to ERSEM nanophy c    from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_ersem_p2n   !: Increment to ERSEM nanophy n    from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_ersem_p2p   !: Increment to ERSEM nanophy p    from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_ersem_p3c   !: Increment to ERSEM picophy c    from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_ersem_p3n   !: Increment to ERSEM picophy n    from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_ersem_p3p   !: Increment to ERSEM picophy p    from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_ersem_p4c   !: Increment to ERSEM microphy c   from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_ersem_p4n   !: Increment to ERSEM microphy n   from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_ersem_p4p   !: Increment to ERSEM microphy p   from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_ersem_z4c   !: Increment to ERSEM mesozoo c    from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_ersem_z5c   !: Increment to ERSEM microzoo c   from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_ersem_z5n   !: Increment to ERSEM microzoo n   from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_ersem_z5p   !: Increment to ERSEM microzoo p   from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_ersem_z6c   !: Increment to ERSEM het flag c   from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_ersem_z6n   !: Increment to ERSEM het flag n   from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_ersem_z6p   !: Increment to ERSEM het flag p   from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_ersem_n1p   !: Increment to ERSEM phosphate    from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_ersem_n3n   !: Increment to ERSEM nitrate      from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_ersem_n4n   !: Increment to ERSEM ammonium     from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_ersem_n5s   !: Increment to ERSEM silicate     from logchl
+#elif defined key_medusa && defined key_foam_medusa
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_medusa_chn  !: Increment to MEDUSA nondiatom chl from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_medusa_chd  !: Increment to MEDUSA diatom chl    from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_medusa_phn  !: Increment to MEDUSA nondiatom n   from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_medusa_phd  !: Increment to MEDUSA diatom n      from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_medusa_pds  !: Increment to MEDUSA diatom s      from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_medusa_zmi  !: Increment to MEDUSA microzoop n   from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_medusa_zme  !: Increment to MEDUSA mesozoop n    from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_medusa_din  !: Increment to MEDUSA nitrate       from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_medusa_sil  !: Increment to MEDUSA silicate      from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_medusa_fer  !: Increment to MEDUSA iron          from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_medusa_det  !: Increment to MEDUSA detritus n    from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_medusa_dtc  !: Increment to MEDUSA detritus c    from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_medusa_dic  !: Increment to MEDUSA dic           from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_medusa_alk  !: Increment to MEDUSA alkalinity    from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_medusa_oxy  !: Increment to MEDUSA oxygen        from logchl
+#elif defined key_hadocc
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_hadocc_nut  !: Increment to HadOCC nutrient      from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_hadocc_phy  !: Increment to HadOCC phytoplankton from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_hadocc_zoo  !: Increment to HadOCC zooplankton   from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_hadocc_det  !: Increment to HadOCC detritus      from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_hadocc_dic  !: Increment to HadOCC DIC           from logchl
+   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_balinc_hadocc_alk  !: Increment to HadOCC alkalinity    from logchl
+#endif
+
    INTEGER :: mld_choice        = 4   !: choice of mld criteria to use for physics assimilation
                                       !: 1) hmld      - Turbocline/mixing depth                           [W points]
@@ -106,4 +188,12 @@
                                       !: 4) hmld_tref - Temperature criterion (0.2 K change from surface) [T points]
 
+   INTEGER :: mld_choice_bgc    = 5   !: choice of mld criteria to use for physics assimilation
+                                      !: 1) hmld      - Turbocline/mixing depth                           [W points]
+                                      !: 2) hmlp      - Density criterion (0.01 kg/m^3 change from 10m)   [W points]
+                                      !: 3) hmld_kara - Kara MLD                                          [Interpolated]
+                                      !: 4) hmld_tref - Temperature criterion (0.2 K change from surface) [T points]
+                                      !: 5) hmlpt     - Density criterion (0.01 kg/m^3 change from 10m)   [T points]
+
+   INTEGER :: nn_asmpfts        = 0   !: number of logchl PFTs assimilated
 
    !! * Substitutions
@@ -157,10 +247,14 @@
       LOGICAL :: lk_surft      ! Logical: T => Increments file contains surft variable 
                                !               so only apply surft increments.
+      !
+      CHARACTER(LEN=2) :: cl_pftstr
       !!
-      NAMELIST/nam_asminc/ ln_bkgwri,                                      &
+      NAMELIST/nam_asminc/ ln_bkgwri, ln_balwri,                           &
          &                 ln_trainc, ln_dyninc, ln_sshinc,                &
+         &                 ln_logchltotinc, ln_logchlpftinc,               &
          &                 ln_asmdin, ln_asmiau,                           &
          &                 nitbkg, nitdin, nitiaustr, nitiaufin, niaufn,   &
-         &                 ln_salfix, salfixmin, nn_divdmp, mld_choice
+         &                 ln_salfix, salfixmin, nn_divdmp, mld_choice,    &
+         &                 mld_choice_bgc
       !!----------------------------------------------------------------------
 
@@ -171,8 +265,11 @@
       ! Set default values
       ln_bkgwri = .FALSE.
+      ln_balwri = .FALSE.
       ln_trainc = .FALSE.
       ln_dyninc = .FALSE.
       ln_sshinc = .FALSE.
       ln_seaiceinc = .FALSE.
+      ln_logchltotinc = .FALSE.
+      ln_logchlpftinc = .FALSE.
       ln_asmdin = .FALSE.
       ln_asmiau = .TRUE.
@@ -185,4 +282,13 @@
       nitiaufin = 150
       niaufn    = 0
+#if defined key_fabm
+      nn_asmpfts = 4
+#elif defined key_medusa && defined key_foam_medusa
+      nn_asmpfts = 2
+#elif defined key_hadocc
+      nn_asmpfts = 1
+#else
+      nn_asmpfts = 0
+#endif
 
       REWIND( numnam_ref )              ! Namelist nam_asminc in reference namelist : Assimilation increment
@@ -194,4 +300,10 @@
 902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nam_asminc in configuration namelist', lwp )
       IF(lwm) WRITE ( numond, nam_asminc )
+
+      IF ( ln_logchltotinc ) THEN
+         nn_asmpfts = 1
+      ELSE IF ( .NOT.( ln_logchlpftinc ) ) THEN
+         nn_asmpfts = 0
+      ENDIF
 
       ! Control print
@@ -202,4 +314,5 @@
          WRITE(numout,*) '   Namelist namasm : set assimilation increment parameters'
          WRITE(numout,*) '      Logical switch for writing out background state          ln_bkgwri = ', ln_bkgwri
+         WRITE(numout,*) '      Logical switch for writing out balancing increments      ln_balwri = ', ln_balwri
          WRITE(numout,*) '      Logical switch for applying tracer increments            ln_trainc = ', ln_trainc
          WRITE(numout,*) '      Logical switch for applying velocity increments          ln_dyninc = ', ln_dyninc
@@ -207,4 +320,7 @@
          WRITE(numout,*) '      Logical switch for Direct Initialization (DI)            ln_asmdin = ', ln_asmdin
          WRITE(numout,*) '      Logical switch for applying sea ice increments        ln_seaiceinc = ', ln_seaiceinc
+         WRITE(numout,*) '      Logical switch for applying total logchl incs      ln_logchltotinc = ', ln_logchltotinc
+         WRITE(numout,*) '      Logical switch for applying PFT   logchl incs      ln_logchlpftinc = ', ln_logchlpftinc
+         WRITE(numout,*) '      Number of logchl PFTs assimilated                       nn_asmpfts = ', nn_asmpfts
          WRITE(numout,*) '      Logical switch for Incremental Analysis Updating (IAU)   ln_asmiau = ', ln_asmiau
          WRITE(numout,*) '      Timestep of background in [0,nitend-nit000-1]            nitbkg    = ', nitbkg
@@ -216,4 +332,5 @@
          WRITE(numout,*) '      Minimum salinity after applying the increments           salfixmin = ', salfixmin
          WRITE(numout,*) '      Choice of MLD for physics assimilation                  mld_choice = ', mld_choice
+         WRITE(numout,*) '      Choice of MLD for BGC assimilation                  mld_choice_bgc = ', mld_choice_bgc
       ENDIF
 
@@ -263,6 +380,8 @@
 
       IF (      ( ( .NOT. ln_asmdin ).AND.( .NOT. ln_asmiau ) ) &
-           .AND.( ( ln_trainc ).OR.( ln_dyninc ).OR.( ln_sshinc ) .OR. ( ln_seaiceinc) )) &
-         & CALL ctl_stop( ' One or more of ln_trainc, ln_dyninc, ln_sshinc and ln_seaiceinc is set to .true.', &
+         & .AND.( ( ln_trainc ).OR.( ln_dyninc ).OR.( ln_sshinc ).OR.( ln_seaiceinc ).OR. &
+         &        ( ln_logchltotinc ).OR.( ln_logchlpftinc ) )) &
+         & CALL ctl_stop( ' One or more of ln_trainc, ln_dyninc, ln_sshinc, ln_seaiceinc,', &
+         &                ' ln_logchltotinc, and ln_logchlpftinc is set to .true.', &
          &                ' but ln_asmdin and ln_asmiau are both set to .false. :', &
          &                ' Inconsistent options')
@@ -273,6 +392,7 @@
 
       IF ( ( .NOT. ln_trainc ).AND.( .NOT. ln_dyninc ).AND.( .NOT. ln_sshinc ).AND.( .NOT. ln_seaiceinc ) &
-         &                     )  &
-         & CALL ctl_warn( ' ln_trainc, ln_dyninc, ln_sshinc and ln_seaiceinc are set to .false. :', &
+         & .AND.( .NOT. ln_logchltotinc ).AND.( .NOT. ln_logchlpftinc ) )  &
+         & CALL ctl_warn( ' ln_trainc, ln_dyninc, ln_sshinc, ln_seaiceinc,', &
+         &                ' ln_logchltotinc, and ln_logchlpftinc are set to .false. :', &
          &                ' The assimilation increments are not applied')
 
@@ -294,4 +414,15 @@
          &                ' Background time step for Direct Initialization is outside', &
          &                ' the cycle interval')
+
+      IF ( ( ln_logchltotinc ).AND.( ln_logchlpftinc ) ) THEN
+         CALL ctl_stop( ' ln_logchltotinc and ln_logchlpftinc both set:', &
+            &           ' These options are not compatible')
+      ENDIF
+
+      IF ( ( ln_balwri ).AND.( .NOT. ( ( ln_logchltotinc ).OR.( ln_logchlpftinc ) ) ) ) THEN
+         CALL ctl_warn( ' Balancing increments are only calculated for logchl', &
+            &           ' Not assimilating logchl, so ln_balwri will be set to .false.')
+         ln_balwri = .FALSE.
+      ENDIF
 
       IF ( nstop > 0 ) RETURN       ! if there are any errors then go no further
@@ -396,5 +527,116 @@
       ssh_iau(:,:)    = 0.0
 #endif
-      IF ( ( ln_trainc ).OR.( ln_dyninc ).OR.( ln_sshinc ).OR.( ln_seaiceinc ) ) THEN
+      IF ( ( ln_logchltotinc ).OR.( ln_logchlpftinc ) ) THEN
+         IF ( ln_logchltotinc ) THEN
+            ALLOCATE( logchl_bkginc(jpi,jpj,1))
+         ELSE IF ( ln_logchlpftinc ) THEN
+            ALLOCATE( logchl_bkginc(jpi,jpj,nn_asmpfts))
+         ENDIF
+         logchl_bkginc(:,:,:) = 0.0
+#if defined key_fabm
+         ALLOCATE( logchl_balinc_ersem_chl1(jpi,jpj,jpk) )
+         ALLOCATE( logchl_balinc_ersem_chl2(jpi,jpj,jpk) )
+         ALLOCATE( logchl_balinc_ersem_chl3(jpi,jpj,jpk) )
+         ALLOCATE( logchl_balinc_ersem_chl4(jpi,jpj,jpk) )
+         ALLOCATE( logchl_balinc_ersem_p1c(jpi,jpj,jpk)  )
+         ALLOCATE( logchl_balinc_ersem_p1n(jpi,jpj,jpk)  )
+         ALLOCATE( logchl_balinc_ersem_p1p(jpi,jpj,jpk)  )
+         ALLOCATE( logchl_balinc_ersem_p1s(jpi,jpj,jpk)  )
+         ALLOCATE( logchl_balinc_ersem_p2c(jpi,jpj,jpk)  )
+         ALLOCATE( logchl_balinc_ersem_p2n(jpi,jpj,jpk)  )
+         ALLOCATE( logchl_balinc_ersem_p2p(jpi,jpj,jpk)  )
+         ALLOCATE( logchl_balinc_ersem_p3c(jpi,jpj,jpk)  )
+         ALLOCATE( logchl_balinc_ersem_p3n(jpi,jpj,jpk)  )
+         ALLOCATE( logchl_balinc_ersem_p3p(jpi,jpj,jpk)  )
+         ALLOCATE( logchl_balinc_ersem_p4c(jpi,jpj,jpk)  )
+         ALLOCATE( logchl_balinc_ersem_p4n(jpi,jpj,jpk)  )
+         ALLOCATE( logchl_balinc_ersem_p4p(jpi,jpj,jpk)  )
+         ALLOCATE( logchl_balinc_ersem_z4c(jpi,jpj,jpk)  )
+         ALLOCATE( logchl_balinc_ersem_z5c(jpi,jpj,jpk)  )
+         ALLOCATE( logchl_balinc_ersem_z5n(jpi,jpj,jpk)  )
+         ALLOCATE( logchl_balinc_ersem_z5p(jpi,jpj,jpk)  )
+         ALLOCATE( logchl_balinc_ersem_z6c(jpi,jpj,jpk)  )
+         ALLOCATE( logchl_balinc_ersem_z6n(jpi,jpj,jpk)  )
+         ALLOCATE( logchl_balinc_ersem_z6p(jpi,jpj,jpk)  )
+         ALLOCATE( logchl_balinc_ersem_n1p(jpi,jpj,jpk)  )
+         ALLOCATE( logchl_balinc_ersem_n3n(jpi,jpj,jpk)  )
+         ALLOCATE( logchl_balinc_ersem_n4n(jpi,jpj,jpk)  )
+         ALLOCATE( logchl_balinc_ersem_n5s(jpi,jpj,jpk)  )
+         logchl_balinc_ersem_chl1(:,:,:) = 0.0
+         logchl_balinc_ersem_chl2(:,:,:) = 0.0
+         logchl_balinc_ersem_chl3(:,:,:) = 0.0
+         logchl_balinc_ersem_chl4(:,:,:) = 0.0
+         logchl_balinc_ersem_p1c(:,:,:)  = 0.0
+         logchl_balinc_ersem_p1n(:,:,:)  = 0.0
+         logchl_balinc_ersem_p1p(:,:,:)  = 0.0
+         logchl_balinc_ersem_p1s(:,:,:)  = 0.0
+         logchl_balinc_ersem_p2c(:,:,:)  = 0.0
+         logchl_balinc_ersem_p2n(:,:,:)  = 0.0
+         logchl_balinc_ersem_p2p(:,:,:)  = 0.0
+         logchl_balinc_ersem_p3c(:,:,:)  = 0.0
+         logchl_balinc_ersem_p3n(:,:,:)  = 0.0
+         logchl_balinc_ersem_p3p(:,:,:)  = 0.0
+         logchl_balinc_ersem_p4c(:,:,:)  = 0.0
+         logchl_balinc_ersem_p4n(:,:,:)  = 0.0
+         logchl_balinc_ersem_p4p(:,:,:)  = 0.0
+         logchl_balinc_ersem_z4c(:,:,:)  = 0.0
+         logchl_balinc_ersem_z5c(:,:,:)  = 0.0
+         logchl_balinc_ersem_z5n(:,:,:)  = 0.0
+         logchl_balinc_ersem_z5p(:,:,:)  = 0.0
+         logchl_balinc_ersem_z6c(:,:,:)  = 0.0
+         logchl_balinc_ersem_z6n(:,:,:)  = 0.0
+         logchl_balinc_ersem_z6p(:,:,:)  = 0.0
+         logchl_balinc_ersem_n1p(:,:,:)  = 0.0
+         logchl_balinc_ersem_n3n(:,:,:)  = 0.0
+         logchl_balinc_ersem_n4n(:,:,:)  = 0.0
+         logchl_balinc_ersem_n5s(:,:,:)  = 0.0
+#elif defined key_medusa && defined key_foam_medusa
+         ALLOCATE( logchl_balinc_medusa_chn(jpi,jpj,jpk) )
+         ALLOCATE( logchl_balinc_medusa_chd(jpi,jpj,jpk) )
+         ALLOCATE( logchl_balinc_medusa_phn(jpi,jpj,jpk) )
+         ALLOCATE( logchl_balinc_medusa_phd(jpi,jpj,jpk) )
+         ALLOCATE( logchl_balinc_medusa_pds(jpi,jpj,jpk) )
+         ALLOCATE( logchl_balinc_medusa_zmi(jpi,jpj,jpk) )
+         ALLOCATE( logchl_balinc_medusa_zme(jpi,jpj,jpk) )
+         ALLOCATE( logchl_balinc_medusa_din(jpi,jpj,jpk) )
+         ALLOCATE( logchl_balinc_medusa_sil(jpi,jpj,jpk) )
+         ALLOCATE( logchl_balinc_medusa_fer(jpi,jpj,jpk) )
+         ALLOCATE( logchl_balinc_medusa_det(jpi,jpj,jpk) )
+         ALLOCATE( logchl_balinc_medusa_dtc(jpi,jpj,jpk) )
+         ALLOCATE( logchl_balinc_medusa_dic(jpi,jpj,jpk) )
+         ALLOCATE( logchl_balinc_medusa_alk(jpi,jpj,jpk) )
+         ALLOCATE( logchl_balinc_medusa_oxy(jpi,jpj,jpk) )
+         logchl_balinc_medusa_chn(:,:,:) = 0.0
+         logchl_balinc_medusa_chd(:,:,:) = 0.0
+         logchl_balinc_medusa_phn(:,:,:) = 0.0
+         logchl_balinc_medusa_phd(:,:,:) = 0.0
+         logchl_balinc_medusa_pds(:,:,:) = 0.0
+         logchl_balinc_medusa_zmi(:,:,:) = 0.0
+         logchl_balinc_medusa_zme(:,:,:) = 0.0
+         logchl_balinc_medusa_din(:,:,:) = 0.0
+         logchl_balinc_medusa_sil(:,:,:) = 0.0
+         logchl_balinc_medusa_fer(:,:,:) = 0.0
+         logchl_balinc_medusa_det(:,:,:) = 0.0
+         logchl_balinc_medusa_dtc(:,:,:) = 0.0
+         logchl_balinc_medusa_dic(:,:,:) = 0.0
+         logchl_balinc_medusa_alk(:,:,:) = 0.0
+         logchl_balinc_medusa_oxy(:,:,:) = 0.0
+#elif defined key_hadocc
+         ALLOCATE( logchl_balinc_hadocc_nut(jpi,jpj,jpk) )
+         ALLOCATE( logchl_balinc_hadocc_phy(jpi,jpj,jpk) )
+         ALLOCATE( logchl_balinc_hadocc_zoo(jpi,jpj,jpk) )
+         ALLOCATE( logchl_balinc_hadocc_det(jpi,jpj,jpk) )
+         ALLOCATE( logchl_balinc_hadocc_dic(jpi,jpj,jpk) )
+         ALLOCATE( logchl_balinc_hadocc_alk(jpi,jpj,jpk) )
+         logchl_balinc_hadocc_nut(:,:,:) = 0.0
+         logchl_balinc_hadocc_phy(:,:,:) = 0.0
+         logchl_balinc_hadocc_zoo(:,:,:) = 0.0
+         logchl_balinc_hadocc_det(:,:,:) = 0.0
+         logchl_balinc_hadocc_dic(:,:,:) = 0.0
+         logchl_balinc_hadocc_alk(:,:,:) = 0.0
+#endif
+      ENDIF
+      IF ( ( ln_trainc ).OR.( ln_dyninc ).OR.( ln_sshinc ).OR.( ln_seaiceinc ) &
+         &  .OR.( ln_logchltotinc ).OR.( ln_logchlpftinc ) ) THEN
 
          !--------------------------------------------------------------------
@@ -527,4 +769,25 @@
             ! to allow for differences in masks
             WHERE( ABS( seaice_bkginc(:,:) ) > 1.0e+10 ) seaice_bkginc(:,:) = 0.0
+         ENDIF
+
+         IF ( ln_logchltotinc ) THEN
+            CALL iom_get( inum, jpdom_autoglo, 'bckinlogchl', logchl_bkginc(:,:,1), 1 )
+            ! Apply the masks
+            logchl_bkginc(:,:,1) = logchl_bkginc(:,:,1) * tmask(:,:,1)
+            ! Set missing increments to 0.0 rather than 1e+20
+            ! to allow for differences in masks
+            WHERE( ABS( logchl_bkginc(:,:,1) ) > 1.0e+10 ) logchl_bkginc(:,:,1) = 0.0
+         ENDIF
+
+         IF ( ln_logchlpftinc ) THEN
+            DO jt = 1, nn_asmpfts
+               WRITE(cl_pftstr,'(I2.2)') jt
+               CALL iom_get( inum, jpdom_autoglo, 'bckinlogchl'//cl_pftstr, logchl_bkginc(:,:,jt), 1 )
+               ! Apply the masks
+               logchl_bkginc(:,:,jt) = logchl_bkginc(:,:,jt) * tmask(:,:,1)
+               ! Set missing increments to 0.0 rather than 1e+20
+               ! to allow for differences in masks
+               WHERE( ABS( logchl_bkginc(:,:,jt) ) > 1.0e+10 ) logchl_bkginc(:,:,jt) = 0.0
+            END DO
          ENDIF
 
@@ -1256,4 +1519,424 @@
 
    END SUBROUTINE seaice_asm_inc
+
+   SUBROUTINE logchl_asm_inc( kt )
+      !!----------------------------------------------------------------------
+      !!                    ***  ROUTINE logchl_asm_inc  ***
+      !!          
+      !! ** Purpose : Apply the chlorophyll assimilation increments.
+      !!
+      !! ** Method  : Calculate increments to state variables using nitrogen
+      !!              balancing.
+      !!              Direct initialization or Incremental Analysis Updating.
+      !!
+      !! ** Action  : 
+      !!----------------------------------------------------------------------
+      INTEGER, INTENT(IN) :: kt   ! Current time step
+      !
+      INTEGER  :: jk              ! Loop counter
+      INTEGER  :: it              ! Index
+      REAL(wp) :: zincwgt         ! IAU weight for current time step
+      REAL(wp) :: zincper         ! IAU interval in seconds
+      !!----------------------------------------------------------------------
+
+      IF ( kt <= nit000 ) THEN
+
+         zincper = (nitiaufin_r - nitiaustr_r + 1) * rdt
+
+#if defined key_fabm
+         CALL asm_logchl_bal_ersem( ln_logchlpftinc,          nn_asmpfts,               &
+            &                       mld_choice_bgc,           logchl_bkginc,            &
+            &                       logchl_balinc_ersem_chl1, logchl_balinc_ersem_chl2, &
+            &                       logchl_balinc_ersem_chl3, logchl_balinc_ersem_chl4, &
+            &                       logchl_balinc_ersem_p1c,  logchl_balinc_ersem_p1n,  &
+            &                       logchl_balinc_ersem_p1p,  logchl_balinc_ersem_p1s,  &
+            &                       logchl_balinc_ersem_p2c,  logchl_balinc_ersem_p2n,  &
+            &                       logchl_balinc_ersem_p2p,  logchl_balinc_ersem_p3c,  &
+            &                       logchl_balinc_ersem_p3n,  logchl_balinc_ersem_p3p,  &
+            &                       logchl_balinc_ersem_p4c,  logchl_balinc_ersem_p4n,  &
+            &                       logchl_balinc_ersem_p4p,  logchl_balinc_ersem_z4c,  &
+            &                       logchl_balinc_ersem_z5c,  logchl_balinc_ersem_z5n,  &
+            &                       logchl_balinc_ersem_z5p,  logchl_balinc_ersem_z6c,  &
+            &                       logchl_balinc_ersem_z6n,  logchl_balinc_ersem_z6p,  &
+            &                       logchl_balinc_ersem_n1p,  logchl_balinc_ersem_n3n,  &
+            &                       logchl_balinc_ersem_n4n,  logchl_balinc_ersem_n5s )
+#elif defined key_medusa && defined key_foam_medusa
+         !CALL asm_logchl_bal_medusa()
+         CALL ctl_stop( 'Attempting to assimilate logchl into MEDUSA, ', &
+            &           'but not fully implemented yet' )
+#elif defined key_hadocc
+         !CALL asm_logchl_bal_hadocc()
+         CALL ctl_stop( 'Attempting to assimilate logchl into HadOCC, ', &
+            &           'but not fully implemented yet' )
+#else
+         CALL ctl_stop( 'Attempting to assimilate logchl, ', &
+            &           'but not defined a biogeochemical model' )
+#endif
+
+      ENDIF
+
+      IF ( ln_asmiau ) THEN
+
+         !--------------------------------------------------------------------
+         ! Incremental Analysis Updating
+         !--------------------------------------------------------------------
+
+         IF ( ( kt >= nitiaustr_r ).AND.( kt <= nitiaufin_r ) ) THEN
+
+            it = kt - nit000 + 1
+            zincwgt = wgtiau(it)      ! IAU weight for the current time step
+            ! note this is not a tendency so should not be divided by rdt
+
+            IF(lwp) THEN
+               WRITE(numout,*) 
+               WRITE(numout,*) 'logchl_asm_inc : logchl IAU at time step = ', &
+                  &  kt,' with IAU weight = ', wgtiau(it)
+               WRITE(numout,*) '~~~~~~~~~~~~'
+            ENDIF
+
+            ! Update the biogeochemical tendencies
+#if defined key_fabm
+            DO jk = 1, jpkm1
+               ! Add directly to trn and trb, rather than to tra, as not a tendency
+               trn(:,:,jk,jp_fabm_m1+jp_fabm_chl1) = trn(:,:,jk,jp_fabm_m1+jp_fabm_chl1) + &
+                  &                                  logchl_balinc_ersem_chl1(:,:,jk) * zincwgt
+               trn(:,:,jk,jp_fabm_m1+jp_fabm_chl2) = trn(:,:,jk,jp_fabm_m1+jp_fabm_chl2) + &
+                  &                                  logchl_balinc_ersem_chl2(:,:,jk) * zincwgt
+               trn(:,:,jk,jp_fabm_m1+jp_fabm_chl3) = trn(:,:,jk,jp_fabm_m1+jp_fabm_chl3) + &
+                  &                                  logchl_balinc_ersem_chl3(:,:,jk) * zincwgt
+               trn(:,:,jk,jp_fabm_m1+jp_fabm_chl4) = trn(:,:,jk,jp_fabm_m1+jp_fabm_chl4) + &
+                  &                                  logchl_balinc_ersem_chl4(:,:,jk) * zincwgt
+               trn(:,:,jk,jp_fabm_m1+jp_fabm_p1c)  = trn(:,:,jk,jp_fabm_m1+jp_fabm_p1c) + &
+                  &                                  logchl_balinc_ersem_p1c(:,:,jk)  * zincwgt
+               trn(:,:,jk,jp_fabm_m1+jp_fabm_p1n)  = trn(:,:,jk,jp_fabm_m1+jp_fabm_p1n) + &
+                  &                                  logchl_balinc_ersem_p1n(:,:,jk)  * zincwgt
+               trn(:,:,jk,jp_fabm_m1+jp_fabm_p1p)  = trn(:,:,jk,jp_fabm_m1+jp_fabm_p1p) + &
+                  &                                  logchl_balinc_ersem_p1p(:,:,jk)  * zincwgt
+               trn(:,:,jk,jp_fabm_m1+jp_fabm_p1s)  = trn(:,:,jk,jp_fabm_m1+jp_fabm_p1s) + &
+                  &                                  logchl_balinc_ersem_p1s(:,:,jk)  * zincwgt
+               trn(:,:,jk,jp_fabm_m1+jp_fabm_p2c)  = trn(:,:,jk,jp_fabm_m1+jp_fabm_p2c) + &
+                  &                                  logchl_balinc_ersem_p2c(:,:,jk)  * zincwgt
+               trn(:,:,jk,jp_fabm_m1+jp_fabm_p2n)  = trn(:,:,jk,jp_fabm_m1+jp_fabm_p2n) + &
+                  &                                  logchl_balinc_ersem_p2n(:,:,jk)  * zincwgt
+               trn(:,:,jk,jp_fabm_m1+jp_fabm_p2p)  = trn(:,:,jk,jp_fabm_m1+jp_fabm_p2p) + &
+                  &                                  logchl_balinc_ersem_p2p(:,:,jk)  * zincwgt
+               trn(:,:,jk,jp_fabm_m1+jp_fabm_p3c)  = trn(:,:,jk,jp_fabm_m1+jp_fabm_p3c) + &
+                  &                                  logchl_balinc_ersem_p3c(:,:,jk)  * zincwgt
+               trn(:,:,jk,jp_fabm_m1+jp_fabm_p3n)  = trn(:,:,jk,jp_fabm_m1+jp_fabm_p3n) + &
+                  &                                  logchl_balinc_ersem_p3n(:,:,jk)  * zincwgt
+               trn(:,:,jk,jp_fabm_m1+jp_fabm_p3p)  = trn(:,:,jk,jp_fabm_m1+jp_fabm_p3p) + &
+                  &                                  logchl_balinc_ersem_p3p(:,:,jk)  * zincwgt
+               trn(:,:,jk,jp_fabm_m1+jp_fabm_p4c)  = trn(:,:,jk,jp_fabm_m1+jp_fabm_p4c) + &
+                  &                                  logchl_balinc_ersem_p4c(:,:,jk)  * zincwgt
+               trn(:,:,jk,jp_fabm_m1+jp_fabm_p4n)  = trn(:,:,jk,jp_fabm_m1+jp_fabm_p4n) + &
+                  &                                  logchl_balinc_ersem_p4n(:,:,jk)  * zincwgt
+               trn(:,:,jk,jp_fabm_m1+jp_fabm_p4p)  = trn(:,:,jk,jp_fabm_m1+jp_fabm_p4p) + &
+                  &                                  logchl_balinc_ersem_p4p(:,:,jk)  * zincwgt
+               trn(:,:,jk,jp_fabm_m1+jp_fabm_z4c)  = trn(:,:,jk,jp_fabm_m1+jp_fabm_z4c) + &
+                  &                                  logchl_balinc_ersem_z4c(:,:,jk)  * zincwgt
+               trn(:,:,jk,jp_fabm_m1+jp_fabm_z5c)  = trn(:,:,jk,jp_fabm_m1+jp_fabm_z5c) + &
+                  &                                  logchl_balinc_ersem_z5c(:,:,jk)  * zincwgt
+               trn(:,:,jk,jp_fabm_m1+jp_fabm_z5n)  = trn(:,:,jk,jp_fabm_m1+jp_fabm_z5n) + &
+                  &                                  logchl_balinc_ersem_z5n(:,:,jk)  * zincwgt
+               trn(:,:,jk,jp_fabm_m1+jp_fabm_z5p)  = trn(:,:,jk,jp_fabm_m1+jp_fabm_z5p) + &
+                  &                                  logchl_balinc_ersem_z5p(:,:,jk)  * zincwgt
+               trn(:,:,jk,jp_fabm_m1+jp_fabm_z6c)  = trn(:,:,jk,jp_fabm_m1+jp_fabm_z6c) + &
+                  &                                  logchl_balinc_ersem_z6c(:,:,jk)  * zincwgt
+               trn(:,:,jk,jp_fabm_m1+jp_fabm_z6n)  = trn(:,:,jk,jp_fabm_m1+jp_fabm_z6n) + &
+                  &                                  logchl_balinc_ersem_z6n(:,:,jk)  * zincwgt
+               trn(:,:,jk,jp_fabm_m1+jp_fabm_z6p)  = trn(:,:,jk,jp_fabm_m1+jp_fabm_z6p) + &
+                  &                                  logchl_balinc_ersem_z6p(:,:,jk)  * zincwgt
+               trn(:,:,jk,jp_fabm_m1+jp_fabm_n1p)  = trn(:,:,jk,jp_fabm_m1+jp_fabm_n1p) + &
+                  &                                  logchl_balinc_ersem_n1p(:,:,jk)  * zincwgt
+               trn(:,:,jk,jp_fabm_m1+jp_fabm_n3n)  = trn(:,:,jk,jp_fabm_m1+jp_fabm_n3n) + &
+                  &                                  logchl_balinc_ersem_n3n(:,:,jk)  * zincwgt
+               trn(:,:,jk,jp_fabm_m1+jp_fabm_n4n)  = trn(:,:,jk,jp_fabm_m1+jp_fabm_n4n) + &
+                  &                                  logchl_balinc_ersem_n4n(:,:,jk)  * zincwgt
+               trn(:,:,jk,jp_fabm_m1+jp_fabm_n5s)  = trn(:,:,jk,jp_fabm_m1+jp_fabm_n5s) + &
+                  &                                  logchl_balinc_ersem_n5s(:,:,jk)  * zincwgt
+
+               trb(:,:,jk,jp_fabm_m1+jp_fabm_chl1) = trb(:,:,jk,jp_fabm_m1+jp_fabm_chl1) + &
+                  &                                  logchl_balinc_ersem_chl1(:,:,jk) * zincwgt
+               trb(:,:,jk,jp_fabm_m1+jp_fabm_chl2) = trb(:,:,jk,jp_fabm_m1+jp_fabm_chl2) + &
+                  &                                  logchl_balinc_ersem_chl2(:,:,jk) * zincwgt
+               trb(:,:,jk,jp_fabm_m1+jp_fabm_chl3) = trb(:,:,jk,jp_fabm_m1+jp_fabm_chl3) + &
+                  &                                  logchl_balinc_ersem_chl3(:,:,jk) * zincwgt
+               trb(:,:,jk,jp_fabm_m1+jp_fabm_chl4) = trb(:,:,jk,jp_fabm_m1+jp_fabm_chl4) + &
+                  &                                  logchl_balinc_ersem_chl4(:,:,jk) * zincwgt
+               trb(:,:,jk,jp_fabm_m1+jp_fabm_p1c)  = trb(:,:,jk,jp_fabm_m1+jp_fabm_p1c) + &
+                  &                                  logchl_balinc_ersem_p1c(:,:,jk)  * zincwgt
+               trb(:,:,jk,jp_fabm_m1+jp_fabm_p1n)  = trb(:,:,jk,jp_fabm_m1+jp_fabm_p1n) + &
+                  &                                  logchl_balinc_ersem_p1n(:,:,jk)  * zincwgt
+               trb(:,:,jk,jp_fabm_m1+jp_fabm_p1p)  = trb(:,:,jk,jp_fabm_m1+jp_fabm_p1p) + &
+                  &                                  logchl_balinc_ersem_p1p(:,:,jk)  * zincwgt
+               trb(:,:,jk,jp_fabm_m1+jp_fabm_p1s)  = trb(:,:,jk,jp_fabm_m1+jp_fabm_p1s) + &
+                  &                                  logchl_balinc_ersem_p1s(:,:,jk)  * zincwgt
+               trb(:,:,jk,jp_fabm_m1+jp_fabm_p2c)  = trb(:,:,jk,jp_fabm_m1+jp_fabm_p2c) + &
+                  &                                  logchl_balinc_ersem_p2c(:,:,jk)  * zincwgt
+               trb(:,:,jk,jp_fabm_m1+jp_fabm_p2n)  = trb(:,:,jk,jp_fabm_m1+jp_fabm_p2n) + &
+                  &                                  logchl_balinc_ersem_p2n(:,:,jk)  * zincwgt
+               trb(:,:,jk,jp_fabm_m1+jp_fabm_p2p)  = trb(:,:,jk,jp_fabm_m1+jp_fabm_p2p) + &
+                  &                                  logchl_balinc_ersem_p2p(:,:,jk)  * zincwgt
+               trb(:,:,jk,jp_fabm_m1+jp_fabm_p3c)  = trb(:,:,jk,jp_fabm_m1+jp_fabm_p3c) + &
+                  &                                  logchl_balinc_ersem_p3c(:,:,jk)  * zincwgt
+               trb(:,:,jk,jp_fabm_m1+jp_fabm_p3n)  = trb(:,:,jk,jp_fabm_m1+jp_fabm_p3n) + &
+                  &                                  logchl_balinc_ersem_p3n(:,:,jk)  * zincwgt
+               trb(:,:,jk,jp_fabm_m1+jp_fabm_p3p)  = trb(:,:,jk,jp_fabm_m1+jp_fabm_p3p) + &
+                  &                                  logchl_balinc_ersem_p3p(:,:,jk)  * zincwgt
+               trb(:,:,jk,jp_fabm_m1+jp_fabm_p4c)  = trb(:,:,jk,jp_fabm_m1+jp_fabm_p4c) + &
+                  &                                  logchl_balinc_ersem_p4c(:,:,jk)  * zincwgt
+               trb(:,:,jk,jp_fabm_m1+jp_fabm_p4n)  = trb(:,:,jk,jp_fabm_m1+jp_fabm_p4n) + &
+                  &                                  logchl_balinc_ersem_p4n(:,:,jk)  * zincwgt
+               trb(:,:,jk,jp_fabm_m1+jp_fabm_p4p)  = trb(:,:,jk,jp_fabm_m1+jp_fabm_p4p) + &
+                  &                                  logchl_balinc_ersem_p4p(:,:,jk)  * zincwgt
+               trb(:,:,jk,jp_fabm_m1+jp_fabm_z4c)  = trb(:,:,jk,jp_fabm_m1+jp_fabm_z4c) + &
+                  &                                  logchl_balinc_ersem_z4c(:,:,jk)  * zincwgt
+               trb(:,:,jk,jp_fabm_m1+jp_fabm_z5c)  = trb(:,:,jk,jp_fabm_m1+jp_fabm_z5c) + &
+                  &                                  logchl_balinc_ersem_z5c(:,:,jk)  * zincwgt
+               trb(:,:,jk,jp_fabm_m1+jp_fabm_z5n)  = trb(:,:,jk,jp_fabm_m1+jp_fabm_z5n) + &
+                  &                                  logchl_balinc_ersem_z5n(:,:,jk)  * zincwgt
+               trb(:,:,jk,jp_fabm_m1+jp_fabm_z5p)  = trb(:,:,jk,jp_fabm_m1+jp_fabm_z5p) + &
+                  &                                  logchl_balinc_ersem_z5p(:,:,jk)  * zincwgt
+               trb(:,:,jk,jp_fabm_m1+jp_fabm_z6c)  = trb(:,:,jk,jp_fabm_m1+jp_fabm_z6c) + &
+                  &                                  logchl_balinc_ersem_z6c(:,:,jk)  * zincwgt
+               trb(:,:,jk,jp_fabm_m1+jp_fabm_z6n)  = trb(:,:,jk,jp_fabm_m1+jp_fabm_z6n) + &
+                  &                                  logchl_balinc_ersem_z6n(:,:,jk)  * zincwgt
+               trb(:,:,jk,jp_fabm_m1+jp_fabm_z6p)  = trb(:,:,jk,jp_fabm_m1+jp_fabm_z6p) + &
+                  &                                  logchl_balinc_ersem_z6p(:,:,jk)  * zincwgt
+               trb(:,:,jk,jp_fabm_m1+jp_fabm_n1p)  = trb(:,:,jk,jp_fabm_m1+jp_fabm_n1p) + &
+                  &                                  logchl_balinc_ersem_n1p(:,:,jk)  * zincwgt
+               trb(:,:,jk,jp_fabm_m1+jp_fabm_n3n)  = trb(:,:,jk,jp_fabm_m1+jp_fabm_n3n) + &
+                  &                                  logchl_balinc_ersem_n3n(:,:,jk)  * zincwgt
+               trb(:,:,jk,jp_fabm_m1+jp_fabm_n4n)  = trb(:,:,jk,jp_fabm_m1+jp_fabm_n4n) + &
+                  &                                  logchl_balinc_ersem_n4n(:,:,jk)  * zincwgt
+               trb(:,:,jk,jp_fabm_m1+jp_fabm_n5s)  = trb(:,:,jk,jp_fabm_m1+jp_fabm_n5s) + &
+                  &                                  logchl_balinc_ersem_n5s(:,:,jk)  * zincwgt
+            END DO
+
+#elif defined key_medusa && defined key_foam_medusa
+            DO jk = 1, jpkm1
+               ! Add directly to trn and trb, rather than to tra, as not a tendency
+               trn(:,:,jk,jpchn) = trn(:,:,jk,jpchn) + logchl_balinc_medusa_chn(:,:,jk) * zincwgt
+               trn(:,:,jk,jpchd) = trn(:,:,jk,jpchd) + logchl_balinc_medusa_chd(:,:,jk) * zincwgt
+               trn(:,:,jk,jpphn) = trn(:,:,jk,jpphn) + logchl_balinc_medusa_phn(:,:,jk) * zincwgt
+               trn(:,:,jk,jpphd) = trn(:,:,jk,jpphd) + logchl_balinc_medusa_phd(:,:,jk) * zincwgt
+               trn(:,:,jk,jpzmi) = trn(:,:,jk,jpzmi) + logchl_balinc_medusa_zmi(:,:,jk) * zincwgt
+               trn(:,:,jk,jpzme) = trn(:,:,jk,jpzme) + logchl_balinc_medusa_zme(:,:,jk) * zincwgt
+               trn(:,:,jk,jpdin) = trn(:,:,jk,jpdin) + logchl_balinc_medusa_din(:,:,jk) * zincwgt
+               trn(:,:,jk,jpsil) = trn(:,:,jk,jpsil) + logchl_balinc_medusa_sil(:,:,jk) * zincwgt
+               trn(:,:,jk,jpfer) = trn(:,:,jk,jpfer) + logchl_balinc_medusa_fer(:,:,jk) * zincwgt
+               trn(:,:,jk,jpdet) = trn(:,:,jk,jpdet) + logchl_balinc_medusa_det(:,:,jk) * zincwgt
+               trn(:,:,jk,jppds) = trn(:,:,jk,jppds) + logchl_balinc_medusa_pds(:,:,jk) * zincwgt
+               trn(:,:,jk,jpdtc) = trn(:,:,jk,jpdtc) + logchl_balinc_medusa_dtc(:,:,jk) * zincwgt
+               trn(:,:,jk,jpdic) = trn(:,:,jk,jpdic) + logchl_balinc_medusa_dic(:,:,jk) * zincwgt
+               trn(:,:,jk,jpalk) = trn(:,:,jk,jpalk) + logchl_balinc_medusa_alk(:,:,jk) * zincwgt
+               trn(:,:,jk,jpoxy) = trn(:,:,jk,jpoxy) + logchl_balinc_medusa_oxy(:,:,jk) * zincwgt
+
+               trb(:,:,jk,jpchn) = trb(:,:,jk,jpchn) + logchl_balinc_medusa_chn(:,:,jk) * zincwgt
+               trb(:,:,jk,jpchd) = trb(:,:,jk,jpchd) + logchl_balinc_medusa_chd(:,:,jk) * zincwgt
+               trb(:,:,jk,jpphn) = trb(:,:,jk,jpphn) + logchl_balinc_medusa_phn(:,:,jk) * zincwgt
+               trb(:,:,jk,jpphd) = trb(:,:,jk,jpphd) + logchl_balinc_medusa_phd(:,:,jk) * zincwgt
+               trb(:,:,jk,jpzmi) = trb(:,:,jk,jpzmi) + logchl_balinc_medusa_zmi(:,:,jk) * zincwgt
+               trb(:,:,jk,jpzme) = trb(:,:,jk,jpzme) + logchl_balinc_medusa_zme(:,:,jk) * zincwgt
+               trb(:,:,jk,jpdin) = trb(:,:,jk,jpdin) + logchl_balinc_medusa_din(:,:,jk) * zincwgt
+               trb(:,:,jk,jpsil) = trb(:,:,jk,jpsil) + logchl_balinc_medusa_sil(:,:,jk) * zincwgt
+               trb(:,:,jk,jpfer) = trb(:,:,jk,jpfer) + logchl_balinc_medusa_fer(:,:,jk) * zincwgt
+               trb(:,:,jk,jpdet) = trb(:,:,jk,jpdet) + logchl_balinc_medusa_det(:,:,jk) * zincwgt
+               trb(:,:,jk,jppds) = trb(:,:,jk,jppds) + logchl_balinc_medusa_pds(:,:,jk) * zincwgt
+               trb(:,:,jk,jpdtc) = trb(:,:,jk,jpdtc) + logchl_balinc_medusa_dtc(:,:,jk) * zincwgt
+               trb(:,:,jk,jpdic) = trb(:,:,jk,jpdic) + logchl_balinc_medusa_dic(:,:,jk) * zincwgt
+               trb(:,:,jk,jpalk) = trb(:,:,jk,jpalk) + logchl_balinc_medusa_alk(:,:,jk) * zincwgt
+               trb(:,:,jk,jpoxy) = trb(:,:,jk,jpoxy) + logchl_balinc_medusa_oxy(:,:,jk) * zincwgt
+            END DO
+
+#elif defined key_hadocc
+            DO jk = 1, jpkm1
+               ! Add directly to trn and trb, rather than to tra, as not a tendency
+               trn(:,:,jk,jp_had_nut) = trn(:,:,jk,jp_had_nut) + logchl_balinc_hadocc_nut(:,:,jk) * zincwgt
+               trn(:,:,jk,jp_had_phy) = trn(:,:,jk,jp_had_phy) + logchl_balinc_hadocc_phy(:,:,jk) * zincwgt
+               trn(:,:,jk,jp_had_zoo) = trn(:,:,jk,jp_had_zoo) + logchl_balinc_hadocc_zoo(:,:,jk) * zincwgt
+               trn(:,:,jk,jp_had_pdn) = trn(:,:,jk,jp_had_pdn) + logchl_balinc_hadocc_det(:,:,jk) * zincwgt
+               trn(:,:,jk,jp_had_dic) = trn(:,:,jk,jp_had_dic) + logchl_balinc_hadocc_dic(:,:,jk) * zincwgt
+               trn(:,:,jk,jp_had_alk) = trn(:,:,jk,jp_had_alk) + logchl_balinc_hadocc_alk(:,:,jk) * zincwgt
+            
+               trb(:,:,jk,jp_had_nut) = trb(:,:,jk,jp_had_nut) + logchl_balinc_hadocc_nut(:,:,jk) * zincwgt
+               trb(:,:,jk,jp_had_phy) = trb(:,:,jk,jp_had_phy) + logchl_balinc_hadocc_phy(:,:,jk) * zincwgt
+               trb(:,:,jk,jp_had_zoo) = trb(:,:,jk,jp_had_zoo) + logchl_balinc_hadocc_zoo(:,:,jk) * zincwgt
+               trb(:,:,jk,jp_had_pdn) = trb(:,:,jk,jp_had_pdn) + logchl_balinc_hadocc_det(:,:,jk) * zincwgt
+               trb(:,:,jk,jp_had_dic) = trb(:,:,jk,jp_had_dic) + logchl_balinc_hadocc_dic(:,:,jk) * zincwgt
+               trb(:,:,jk,jp_had_alk) = trb(:,:,jk,jp_had_alk) + logchl_balinc_hadocc_alk(:,:,jk) * zincwgt
+            END DO
+#endif
+           
+            ! Do not deallocate arrays - needed by asm_bal_wri
+            ! which is called at end of model run
+
+         ENDIF
+
+      ELSEIF ( ln_asmdin ) THEN 
+
+         !--------------------------------------------------------------------
+         ! Direct Initialization
+         !--------------------------------------------------------------------
+         
+         IF ( kt == nitdin_r ) THEN
+
+            neuler = 0                    ! Force Euler forward step
+
+#if defined key_fabm
+            ! Initialize the now fields with the background + increment
+            ! Background currently is what the model is initialised with
+            CALL ctl_warn( ' Doing direct initialisation of ERSEM with chlorophyll assimilation', &
+               &           ' Background state is taken from model rather than background file' )
+            trn(:,:,:,jp_fabm_m1+jp_fabm_chl1) = trn(:,:,:,jp_fabm_m1+jp_fabm_chl1) + &
+               &                                 logchl_balinc_ersem_chl1(:,:,:)
+            trn(:,:,:,jp_fabm_m1+jp_fabm_chl2) = trn(:,:,:,jp_fabm_m1+jp_fabm_chl2) + &
+               &                                 logchl_balinc_ersem_chl2(:,:,:)
+            trn(:,:,:,jp_fabm_m1+jp_fabm_chl3) = trn(:,:,:,jp_fabm_m1+jp_fabm_chl3) + &
+               &                                 logchl_balinc_ersem_chl3(:,:,:)
+            trn(:,:,:,jp_fabm_m1+jp_fabm_chl4) = trn(:,:,:,jp_fabm_m1+jp_fabm_chl4) + &
+               &                                 logchl_balinc_ersem_chl4(:,:,:)
+            trn(:,:,:,jp_fabm_m1+jp_fabm_p1c)  = trn(:,:,:,jp_fabm_m1+jp_fabm_p1c)  + &
+               &                                 logchl_balinc_ersem_p1c(:,:,:)
+            trn(:,:,:,jp_fabm_m1+jp_fabm_p1n)  = trn(:,:,:,jp_fabm_m1+jp_fabm_p1n)  + &
+               &                                 logchl_balinc_ersem_p1n(:,:,:)
+            trn(:,:,:,jp_fabm_m1+jp_fabm_p1p)  = trn(:,:,:,jp_fabm_m1+jp_fabm_p1p)  + &
+               &                                 logchl_balinc_ersem_p1p(:,:,:)
+            trn(:,:,:,jp_fabm_m1+jp_fabm_p1s)  = trn(:,:,:,jp_fabm_m1+jp_fabm_p1s)  + &
+               &                                 logchl_balinc_ersem_p1s(:,:,:)
+            trn(:,:,:,jp_fabm_m1+jp_fabm_p2c)  = trn(:,:,:,jp_fabm_m1+jp_fabm_p2c)  + &
+               &                                 logchl_balinc_ersem_p2c(:,:,:)
+            trn(:,:,:,jp_fabm_m1+jp_fabm_p2n)  = trn(:,:,:,jp_fabm_m1+jp_fabm_p2n)  + &
+               &                                 logchl_balinc_ersem_p2n(:,:,:)
+            trn(:,:,:,jp_fabm_m1+jp_fabm_p2p)  = trn(:,:,:,jp_fabm_m1+jp_fabm_p2p)  + &
+               &                                 logchl_balinc_ersem_p2p(:,:,:)
+            trn(:,:,:,jp_fabm_m1+jp_fabm_p3c)  = trn(:,:,:,jp_fabm_m1+jp_fabm_p3c)  + &
+               &                                 logchl_balinc_ersem_p3c(:,:,:)
+            trn(:,:,:,jp_fabm_m1+jp_fabm_p3n)  = trn(:,:,:,jp_fabm_m1+jp_fabm_p3n)  + &
+               &                                 logchl_balinc_ersem_p3n(:,:,:)
+            trn(:,:,:,jp_fabm_m1+jp_fabm_p3p)  = trn(:,:,:,jp_fabm_m1+jp_fabm_p3p)  + &
+               &                                 logchl_balinc_ersem_p3p(:,:,:)
+            trn(:,:,:,jp_fabm_m1+jp_fabm_p4c)  = trn(:,:,:,jp_fabm_m1+jp_fabm_p4c)  + &
+               &                                 logchl_balinc_ersem_p4c(:,:,:)
+            trn(:,:,:,jp_fabm_m1+jp_fabm_p4n)  = trn(:,:,:,jp_fabm_m1+jp_fabm_p4n)  + &
+               &                                 logchl_balinc_ersem_p4n(:,:,:)
+            trn(:,:,:,jp_fabm_m1+jp_fabm_p4p)  = trn(:,:,:,jp_fabm_m1+jp_fabm_p4p)  + &
+               &                                 logchl_balinc_ersem_p4p(:,:,:)
+            trn(:,:,:,jp_fabm_m1+jp_fabm_z4c)  = trn(:,:,:,jp_fabm_m1+jp_fabm_z4c)  + &
+               &                                 logchl_balinc_ersem_z4c(:,:,:)
+            trn(:,:,:,jp_fabm_m1+jp_fabm_z5c)  = trn(:,:,:,jp_fabm_m1+jp_fabm_z5c)  + &
+               &                                 logchl_balinc_ersem_z5c(:,:,:)
+            trn(:,:,:,jp_fabm_m1+jp_fabm_z5n)  = trn(:,:,:,jp_fabm_m1+jp_fabm_z5n)  + &
+               &                                 logchl_balinc_ersem_z5n(:,:,:)
+            trn(:,:,:,jp_fabm_m1+jp_fabm_z5p)  = trn(:,:,:,jp_fabm_m1+jp_fabm_z5p)  + &
+               &                                 logchl_balinc_ersem_z5p(:,:,:)
+            trn(:,:,:,jp_fabm_m1+jp_fabm_z6c)  = trn(:,:,:,jp_fabm_m1+jp_fabm_z6c)  + &
+               &                                 logchl_balinc_ersem_z6c(:,:,:)
+            trn(:,:,:,jp_fabm_m1+jp_fabm_z6n)  = trn(:,:,:,jp_fabm_m1+jp_fabm_z6n)  + &
+               &                                 logchl_balinc_ersem_z6n(:,:,:)
+            trn(:,:,:,jp_fabm_m1+jp_fabm_z6p)  = trn(:,:,:,jp_fabm_m1+jp_fabm_z6p)  + &
+               &                                 logchl_balinc_ersem_z6p(:,:,:)
+            trn(:,:,:,jp_fabm_m1+jp_fabm_n1p)  = trn(:,:,:,jp_fabm_m1+jp_fabm_n1p)  + &
+               &                                 logchl_balinc_ersem_n1p(:,:,:)
+            trn(:,:,:,jp_fabm_m1+jp_fabm_n3n)  = trn(:,:,:,jp_fabm_m1+jp_fabm_n3n)  + &
+               &                                 logchl_balinc_ersem_n3n(:,:,:)
+            trn(:,:,:,jp_fabm_m1+jp_fabm_n4n)  = trn(:,:,:,jp_fabm_m1+jp_fabm_n4n)  + &
+               &                                 logchl_balinc_ersem_n4n(:,:,:)
+            trn(:,:,:,jp_fabm_m1+jp_fabm_n5s)  = trn(:,:,:,jp_fabm_m1+jp_fabm_n5s)  + &
+               &                                 logchl_balinc_ersem_n5s(:,:,:)
+
+            trb(:,:,:,jp_fabm_m1+jp_fabm_chl1) = trn(:,:,:,jp_fabm_m1+jp_fabm_chl1)
+            trb(:,:,:,jp_fabm_m1+jp_fabm_chl2) = trn(:,:,:,jp_fabm_m1+jp_fabm_chl2)
+            trb(:,:,:,jp_fabm_m1+jp_fabm_chl3) = trn(:,:,:,jp_fabm_m1+jp_fabm_chl3)
+            trb(:,:,:,jp_fabm_m1+jp_fabm_chl4) = trn(:,:,:,jp_fabm_m1+jp_fabm_chl4)
+            trb(:,:,:,jp_fabm_m1+jp_fabm_p1c)  = trn(:,:,:,jp_fabm_m1+jp_fabm_p1c)
+            trb(:,:,:,jp_fabm_m1+jp_fabm_p1n)  = trn(:,:,:,jp_fabm_m1+jp_fabm_p1n)
+            trb(:,:,:,jp_fabm_m1+jp_fabm_p1p)  = trn(:,:,:,jp_fabm_m1+jp_fabm_p1p)
+            trb(:,:,:,jp_fabm_m1+jp_fabm_p1s)  = trn(:,:,:,jp_fabm_m1+jp_fabm_p1s)
+            trb(:,:,:,jp_fabm_m1+jp_fabm_p2c)  = trn(:,:,:,jp_fabm_m1+jp_fabm_p2c)
+            trb(:,:,:,jp_fabm_m1+jp_fabm_p2n)  = trn(:,:,:,jp_fabm_m1+jp_fabm_p2n)
+            trb(:,:,:,jp_fabm_m1+jp_fabm_p2p)  = trn(:,:,:,jp_fabm_m1+jp_fabm_p2p)
+            trb(:,:,:,jp_fabm_m1+jp_fabm_p3c)  = trn(:,:,:,jp_fabm_m1+jp_fabm_p3c)
+            trb(:,:,:,jp_fabm_m1+jp_fabm_p3n)  = trn(:,:,:,jp_fabm_m1+jp_fabm_p3n)
+            trb(:,:,:,jp_fabm_m1+jp_fabm_p3p)  = trn(:,:,:,jp_fabm_m1+jp_fabm_p3p)
+            trb(:,:,:,jp_fabm_m1+jp_fabm_p4c)  = trn(:,:,:,jp_fabm_m1+jp_fabm_p4c)
+            trb(:,:,:,jp_fabm_m1+jp_fabm_p4n)  = trn(:,:,:,jp_fabm_m1+jp_fabm_p4n)
+            trb(:,:,:,jp_fabm_m1+jp_fabm_p4p)  = trn(:,:,:,jp_fabm_m1+jp_fabm_p4p)
+            trb(:,:,:,jp_fabm_m1+jp_fabm_z4c)  = trn(:,:,:,jp_fabm_m1+jp_fabm_z4c)
+            trb(:,:,:,jp_fabm_m1+jp_fabm_z5c)  = trn(:,:,:,jp_fabm_m1+jp_fabm_z5c)
+            trb(:,:,:,jp_fabm_m1+jp_fabm_z5n)  = trn(:,:,:,jp_fabm_m1+jp_fabm_z5n)
+            trb(:,:,:,jp_fabm_m1+jp_fabm_z5p)  = trn(:,:,:,jp_fabm_m1+jp_fabm_z5p)
+            trb(:,:,:,jp_fabm_m1+jp_fabm_z6c)  = trn(:,:,:,jp_fabm_m1+jp_fabm_z6c)
+            trb(:,:,:,jp_fabm_m1+jp_fabm_z6n)  = trn(:,:,:,jp_fabm_m1+jp_fabm_z6n)
+            trb(:,:,:,jp_fabm_m1+jp_fabm_z6p)  = trn(:,:,:,jp_fabm_m1+jp_fabm_z6p)
+            trb(:,:,:,jp_fabm_m1+jp_fabm_n1p)  = trn(:,:,:,jp_fabm_m1+jp_fabm_n1p)
+            trb(:,:,:,jp_fabm_m1+jp_fabm_n3n)  = trn(:,:,:,jp_fabm_m1+jp_fabm_n3n)
+            trb(:,:,:,jp_fabm_m1+jp_fabm_n4n)  = trn(:,:,:,jp_fabm_m1+jp_fabm_n4n)
+            trb(:,:,:,jp_fabm_m1+jp_fabm_n5s)  = trn(:,:,:,jp_fabm_m1+jp_fabm_n5s)
+
+#elif defined key_medusa && defined key_foam_medusa
+            ! Initialize the now fields with the background + increment
+            ! Background currently is what the model is initialised with
+            CALL ctl_warn( ' Doing direct initialisation of MEDUSA with chlorophyll assimilation', &
+               &           ' Background state is taken from model rather than background file' )
+            trn(:,:,:,jpchn) = trn(:,:,:,jpchn) + logchl_balinc_medusa_chn(:,:,:)
+            trn(:,:,:,jpchd) = trn(:,:,:,jpchd) + logchl_balinc_medusa_chd(:,:,:)
+            trn(:,:,:,jpphn) = trn(:,:,:,jpphn) + logchl_balinc_medusa_phn(:,:,:)
+            trn(:,:,:,jpphd) = trn(:,:,:,jpphd) + logchl_balinc_medusa_phd(:,:,:)
+            trn(:,:,:,jpzmi) = trn(:,:,:,jpzmi) + logchl_balinc_medusa_zmi(:,:,:)
+            trn(:,:,:,jpzme) = trn(:,:,:,jpzme) + logchl_balinc_medusa_zme(:,:,:)
+            trn(:,:,:,jpdin) = trn(:,:,:,jpdin) + logchl_balinc_medusa_din(:,:,:)
+            trn(:,:,:,jpsil) = trn(:,:,:,jpsil) + logchl_balinc_medusa_sil(:,:,:)
+            trn(:,:,:,jpfer) = trn(:,:,:,jpfer) + logchl_balinc_medusa_fer(:,:,:)
+            trn(:,:,:,jpdet) = trn(:,:,:,jpdet) + logchl_balinc_medusa_det(:,:,:)
+            trn(:,:,:,jppds) = trn(:,:,:,jppds) + logchl_balinc_medusa_pds(:,:,:)
+            trn(:,:,:,jpdtc) = trn(:,:,:,jpdtc) + logchl_balinc_medusa_dtc(:,:,:)
+            trn(:,:,:,jpdic) = trn(:,:,:,jpdic) + logchl_balinc_medusa_dic(:,:,:)
+            trn(:,:,:,jpalk) = trn(:,:,:,jpalk) + logchl_balinc_medusa_alk(:,:,:)
+            trn(:,:,:,jpoxy) = trn(:,:,:,jpoxy) + logchl_balinc_medusa_oxy(:,:,:)
+
+            trb(:,:,:,jpchn) = trn(:,:,:,jpchn)
+            trb(:,:,:,jpchd) = trn(:,:,:,jpchd)
+            trb(:,:,:,jpphn) = trn(:,:,:,jpphn)
+            trb(:,:,:,jpphd) = trn(:,:,:,jpphd)
+            trb(:,:,:,jpzmi) = trn(:,:,:,jpzmi)
+            trb(:,:,:,jpzme) = trn(:,:,:,jpzme)
+            trb(:,:,:,jpdin) = trn(:,:,:,jpdin)
+            trb(:,:,:,jpsil) = trn(:,:,:,jpsil)
+            trb(:,:,:,jpfer) = trn(:,:,:,jpfer)
+            trb(:,:,:,jpdet) = trn(:,:,:,jpdet)
+            trb(:,:,:,jppds) = trn(:,:,:,jppds)
+            trb(:,:,:,jpdtc) = trn(:,:,:,jpdtc)
+            trb(:,:,:,jpdic) = trn(:,:,:,jpdic)
+            trb(:,:,:,jpalk) = trn(:,:,:,jpalk)
+            trb(:,:,:,jpoxy) = trn(:,:,:,jpoxy)
+
+#elif defined key_hadocc
+            ! Initialize the now fields with the background + increment
+            ! Background currently is what the model is initialised with
+            CALL ctl_warn( ' Doing direct initialisation of HadOCC with chlorophyll assimilation', &
+               &           ' Background state is taken from model rather than background file' )
+            trn(:,:,:,jp_had_nut) = trn(:,:,:,jp_had_nut) + logchl_balinc_hadocc_nut(:,:,:)
+            trn(:,:,:,jp_had_phy) = trn(:,:,:,jp_had_phy) + logchl_balinc_hadocc_phy(:,:,:)
+            trn(:,:,:,jp_had_zoo) = trn(:,:,:,jp_had_zoo) + logchl_balinc_hadocc_zoo(:,:,:)
+            trn(:,:,:,jp_had_pdn) = trn(:,:,:,jp_had_pdn) + logchl_balinc_hadocc_det(:,:,:)
+            trn(:,:,:,jp_had_dic) = trn(:,:,:,jp_had_dic) + logchl_balinc_hadocc_dic(:,:,:)
+            trn(:,:,:,jp_had_alk) = trn(:,:,:,jp_had_alk) + logchl_balinc_hadocc_alk(:,:,:)
+
+            trb(:,:,:,jp_had_nut) = trn(:,:,:,jp_had_nut)
+            trb(:,:,:,jp_had_phy) = trn(:,:,:,jp_had_phy)
+            trb(:,:,:,jp_had_zoo) = trn(:,:,:,jp_had_zoo)
+            trb(:,:,:,jp_had_pdn) = trn(:,:,:,jp_had_pdn)
+            trb(:,:,:,jp_had_dic) = trn(:,:,:,jp_had_dic)
+            trb(:,:,:,jp_had_alk) = trn(:,:,:,jp_had_alk)
+#endif
+ 
+            ! Do not deallocate arrays - needed by asm_bal_wri
+            ! which is called at end of model run
+         ENDIF
+         !
+      ENDIF
+      !
+   END SUBROUTINE logchl_asm_inc
    
    !!======================================================================
Index: /branches/UKMO/dev_r5518_v3.4_asm_nemovar_community_bgc_ersem/NEMOGCM/NEMO/OPA_SRC/ASM/asmlogchlbal_ersem.F90
===================================================================
--- /branches/UKMO/dev_r5518_v3.4_asm_nemovar_community_bgc_ersem/NEMOGCM/NEMO/OPA_SRC/ASM/asmlogchlbal_ersem.F90	(revision 6983)
+++ /branches/UKMO/dev_r5518_v3.4_asm_nemovar_community_bgc_ersem/NEMOGCM/NEMO/OPA_SRC/ASM/asmlogchlbal_ersem.F90	(revision 6983)
@@ -0,0 +1,308 @@
+MODULE asmlogchlbal_ersem
+   !!======================================================================
+   !!                       ***  MODULE asmlogchlbal_ersem  ***
+   !! Calculate increments to ERSEM based on surface logchl increments
+   !!======================================================================
+   !! History :  3.6  ! 2016-09  (D. Ford)     Original code
+   !!----------------------------------------------------------------------
+#if defined key_asminc && defined key_fabm
+   !!----------------------------------------------------------------------
+   !!   'key_asminc'         : assimilation increment interface
+   !!   'key_fabm'           : FABM-ERSEM coupling
+   !!----------------------------------------------------------------------
+   !!   asm_logchl_bal_ersem : routine to calculate increments to ERSEM
+   !!----------------------------------------------------------------------
+   USE par_kind,    ONLY: wp             ! kind parameters
+   USE par_oce,     ONLY: jpi, jpj, jpk  ! domain array sizes
+   USE dom_oce,     ONLY: gdepw_n        ! domain information
+   USE zdfmxl                            ! mixed layer depth
+   USE iom                               ! i/o
+   USE trc,         ONLY: trn            ! ERSEM state variables
+   USE par_fabm                          ! FABM parameters
+
+   IMPLICIT NONE
+   PRIVATE                   
+
+   PUBLIC asm_logchl_bal_ersem
+
+   REAL(wp), PARAMETER :: tiny = 1e-20
+
+CONTAINS
+
+   SUBROUTINE asm_logchl_bal_ersem( ln_logchlpftinc,          npfts,                    &
+      &                             mld_choice_bgc,           logchl_bkginc,            &      
+      &                             logchl_balinc_ersem_chl1, logchl_balinc_ersem_chl2, &
+      &                             logchl_balinc_ersem_chl3, logchl_balinc_ersem_chl4, &
+      &                             logchl_balinc_ersem_p1c,  logchl_balinc_ersem_p1n,  &
+      &                             logchl_balinc_ersem_p1p,  logchl_balinc_ersem_p1s,  &
+      &                             logchl_balinc_ersem_p2c,  logchl_balinc_ersem_p2n,  &
+      &                             logchl_balinc_ersem_p2p,  logchl_balinc_ersem_p3c,  &
+      &                             logchl_balinc_ersem_p3n,  logchl_balinc_ersem_p3p,  &
+      &                             logchl_balinc_ersem_p4c,  logchl_balinc_ersem_p4n,  &
+      &                             logchl_balinc_ersem_p4p,  logchl_balinc_ersem_z4c,  &
+      &                             logchl_balinc_ersem_z5c,  logchl_balinc_ersem_z5n,  &
+      &                             logchl_balinc_ersem_z5p,  logchl_balinc_ersem_z6c,  &
+      &                             logchl_balinc_ersem_z6n,  logchl_balinc_ersem_z6p,  &
+      &                             logchl_balinc_ersem_n1p,  logchl_balinc_ersem_n3n,  &
+      &                             logchl_balinc_ersem_n4n,  logchl_balinc_ersem_n5s )
+      !!---------------------------------------------------------------------------
+      !!                    ***  ROUTINE asm_logchl_bal_ersem  ***
+      !!
+      !! ** Purpose :   calculate increments to ERSEM from logchl increments
+      !!
+      !! ** Method  :   forthcoming...
+      !!
+      !! ** Action  :   forthcoming...
+      !!
+      !! References :   forthcoming...
+      !!---------------------------------------------------------------------------
+      !!
+      LOGICAL,  INTENT(in   )                           :: ln_logchlpftinc
+      INTEGER,  INTENT(in   )                           :: npfts
+      INTEGER,  INTENT(in   )                           :: mld_choice_bgc
+      REAL(wp), INTENT(in   ), DIMENSION(jpi,jpj,npfts) :: logchl_bkginc
+      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj,jpk)   :: logchl_balinc_ersem_chl1
+      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj,jpk)   :: logchl_balinc_ersem_chl2
+      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj,jpk)   :: logchl_balinc_ersem_chl3
+      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj,jpk)   :: logchl_balinc_ersem_chl4
+      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj,jpk)   :: logchl_balinc_ersem_p1c
+      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj,jpk)   :: logchl_balinc_ersem_p1n
+      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj,jpk)   :: logchl_balinc_ersem_p1p
+      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj,jpk)   :: logchl_balinc_ersem_p1s
+      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj,jpk)   :: logchl_balinc_ersem_p2c
+      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj,jpk)   :: logchl_balinc_ersem_p2n
+      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj,jpk)   :: logchl_balinc_ersem_p2p
+      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj,jpk)   :: logchl_balinc_ersem_p3c
+      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj,jpk)   :: logchl_balinc_ersem_p3n
+      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj,jpk)   :: logchl_balinc_ersem_p3p
+      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj,jpk)   :: logchl_balinc_ersem_p4c
+      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj,jpk)   :: logchl_balinc_ersem_p4n
+      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj,jpk)   :: logchl_balinc_ersem_p4p
+      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj,jpk)   :: logchl_balinc_ersem_z4c
+      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj,jpk)   :: logchl_balinc_ersem_z5c
+      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj,jpk)   :: logchl_balinc_ersem_z5n
+      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj,jpk)   :: logchl_balinc_ersem_z5p
+      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj,jpk)   :: logchl_balinc_ersem_z6c
+      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj,jpk)   :: logchl_balinc_ersem_z6n
+      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj,jpk)   :: logchl_balinc_ersem_z6p
+      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj,jpk)   :: logchl_balinc_ersem_n1p
+      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj,jpk)   :: logchl_balinc_ersem_n3n
+      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj,jpk)   :: logchl_balinc_ersem_n4n
+      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj,jpk)   :: logchl_balinc_ersem_n5s
+      !!
+      INTEGER                      :: ji, jj, jk
+      INTEGER                      :: jkmax
+      REAL(wp)                     :: chl_tot, chl_inc
+      REAL(wp), DIMENSION(jpi,jpj) :: zmld
+      !!---------------------------------------------------------------------------
+
+      ! Split surface logchl incs into surface Chl1-4 incs
+      !
+      IF ( ln_logchlpftinc ) THEN
+         !
+         ! Assimilating separate PFTs, so separately transform each from LogChl to Chl
+         !
+         IF ( npfts /= 4 ) THEN
+            CALL ctl_stop( 'If assimilating PFTs into ERSEM, nn_asmpfts must be 4' )
+         ENDIF
+         DO jj = 1, jpj
+            DO ji = 1, jpi
+               IF ( trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl1) > 0.0 ) THEN
+                  logchl_balinc_ersem_chl1(ji,jj,1) = 10**( LOG10( trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl1) + &
+                     &                                             tiny ) + logchl_bkginc(ji,jj,1) ) -    &
+                     &                                trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl1) - tiny
+               ENDIF
+               IF ( trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl2) > 0.0 ) THEN
+                  logchl_balinc_ersem_chl2(ji,jj,1) = 10**( LOG10( trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl2) + &
+                     &                                             tiny ) + logchl_bkginc(ji,jj,2) ) -    &
+                     &                                trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl2) - tiny
+               ENDIF
+               IF ( trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl3) > 0.0 ) THEN
+                  logchl_balinc_ersem_chl3(ji,jj,1) = 10**( LOG10( trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl3) + &
+                     &                                             tiny ) + logchl_bkginc(ji,jj,3) ) -    &
+                     &                                trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl3) - tiny
+               ENDIF
+               IF ( trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl4) > 0.0 ) THEN
+                  logchl_balinc_ersem_chl4(ji,jj,1) = 10**( LOG10( trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl4) + &
+                     &                                             tiny ) + logchl_bkginc(ji,jj,4) ) -    &
+                     &                                trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl4) - tiny
+               ENDIF
+            END DO
+         END DO
+      ELSE
+         !
+         ! Assimilating total Chl, so transform total from LogChl to Chl
+         ! and split between PFTs according to the existing background ratios
+         !
+         DO jj = 1, jpj
+            DO ji = 1, jpi
+               IF ( ( trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl1) > 0.0 ) .AND. &
+                  & ( trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl2) > 0.0 ) .AND. &
+                  & ( trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl3) > 0.0 ) .AND. &
+                  & ( trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl4) > 0.0 ) ) THEN
+                  chl_tot = trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl1) + trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl2) + &
+                     &      trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl3) + trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl4)
+                  chl_inc = 10**( LOG10( chl_tot + tiny ) + logchl_bkginc(ji,jj,1) ) - chl_tot - tiny
+                  logchl_balinc_ersem_chl1(ji,jj,1) = chl_inc * trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl1) / chl_tot
+                  logchl_balinc_ersem_chl2(ji,jj,1) = chl_inc * trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl2) / chl_tot
+                  logchl_balinc_ersem_chl3(ji,jj,1) = chl_inc * trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl3) / chl_tot
+                  logchl_balinc_ersem_chl4(ji,jj,1) = chl_inc * trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl4) / chl_tot
+               ENDIF
+            END DO
+         END DO
+      ENDIF
+
+      ! Propagate surface Chl1-4 incs through mixed layer
+      ! First, choose mixed layer definition
+      !
+      SELECT CASE( mld_choice_bgc )
+      CASE ( 1 )                   ! Turbocline/mixing depth [W points]
+         zmld(:,:) = hmld(:,:)
+      CASE ( 2 )                   ! Density criterion (0.01 kg/m^3 change from 10m) [W points]
+         zmld(:,:) = hmlp(:,:)
+      CASE ( 3 )                   ! Kara MLD [Interpolated]
+#if defined key_karaml
+         IF ( ln_kara ) THEN
+            zmld(:,:) = hmld_kara(:,:)
+         ELSE
+            CALL ctl_stop( ' Kara mixed layer requested for LogChl assimilation,', &
+               &           ' but ln_kara=.false.' )
+         ENDIF
+#else
+         CALL ctl_stop( ' Kara mixed layer requested for LogChl assimilation,', &
+            &           ' but is not defined' )
+#endif
+      CASE ( 4 )                   ! Temperature criterion (0.2 K change from surface) [T points]
+         zmld(:,:) = hmld_tref(:,:)
+      CASE ( 5 )                   ! Density criterion (0.01 kg/m^3 change from 10m) [T points]
+         zmld(:,:) = hmlpt(:,:)
+      END SELECT
+      !
+      ! Now set MLD to bottom of a level and propagate incs equally through mixed layer
+      !
+      DO jj = 1, jpj
+         DO ji = 1, jpi
+            !
+            jkmax = jpk-1
+            DO jk = jpk-1, 1, -1
+               IF ( ( zmld(ji,jj) >  gdepw_n(ji,jj,jk)   ) .AND. &
+                  & ( zmld(ji,jj) <= gdepw_n(ji,jj,jk+1) ) ) THEN
+                  zmld(ji,jj) = gdepw_n(ji,jj,jk+1)
+                  jkmax = jk
+               ENDIF
+            END DO
+            !
+            DO jk = 2, jkmax
+               logchl_balinc_ersem_chl1(ji,jj,jk) = logchl_balinc_ersem_chl1(ji,jj,1)
+               logchl_balinc_ersem_chl2(ji,jj,jk) = logchl_balinc_ersem_chl2(ji,jj,1)
+               logchl_balinc_ersem_chl3(ji,jj,jk) = logchl_balinc_ersem_chl3(ji,jj,1)
+               logchl_balinc_ersem_chl4(ji,jj,jk) = logchl_balinc_ersem_chl4(ji,jj,1)
+            END DO
+            !
+         END DO
+      END DO
+
+      ! Balance P[1-4][c,n,p,s]
+      ! Maintain existing background ratios
+      !
+      DO jk = 1, jpk
+         DO jj = 1, jpj
+            DO ji = 1, jpi
+               IF ( trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl1) > 0.0 ) THEN
+                  logchl_balinc_ersem_p1c(ji,jj,jk) = logchl_balinc_ersem_chl1(ji,jj,jk)   * &
+                     &                                trn(ji,jj,jk,jp_fabm_m1+jp_fabm_p1c) / &
+                     &                                trn(ji,jj,jk,jp_fabm_m1+jp_fabm_chl1)
+                  logchl_balinc_ersem_p1n(ji,jj,jk) = logchl_balinc_ersem_chl1(ji,jj,jk)   * &
+                     &                                trn(ji,jj,jk,jp_fabm_m1+jp_fabm_p1n) / &
+                     &                                trn(ji,jj,jk,jp_fabm_m1+jp_fabm_chl1)
+                  logchl_balinc_ersem_p1p(ji,jj,jk) = logchl_balinc_ersem_chl1(ji,jj,jk)   * &
+                     &                                trn(ji,jj,jk,jp_fabm_m1+jp_fabm_p1p) / &
+                     &                                trn(ji,jj,jk,jp_fabm_m1+jp_fabm_chl1)
+                  logchl_balinc_ersem_p1s(ji,jj,jk) = logchl_balinc_ersem_chl1(ji,jj,jk)   * &
+                     &                                trn(ji,jj,jk,jp_fabm_m1+jp_fabm_p1s) / &
+                     &                                trn(ji,jj,jk,jp_fabm_m1+jp_fabm_chl1)
+               ENDIF
+               IF ( trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl2) > 0.0 ) THEN
+                  logchl_balinc_ersem_p2c(ji,jj,jk) = logchl_balinc_ersem_chl2(ji,jj,jk)   * &
+                     &                                trn(ji,jj,jk,jp_fabm_m1+jp_fabm_p2c) / &
+                     &                                trn(ji,jj,jk,jp_fabm_m1+jp_fabm_chl2)
+                  logchl_balinc_ersem_p2n(ji,jj,jk) = logchl_balinc_ersem_chl2(ji,jj,jk)   * &
+                     &                                trn(ji,jj,jk,jp_fabm_m1+jp_fabm_p2n) / &
+                     &                                trn(ji,jj,jk,jp_fabm_m1+jp_fabm_chl2)
+                  logchl_balinc_ersem_p2p(ji,jj,jk) = logchl_balinc_ersem_chl2(ji,jj,jk)   * &
+                     &                                trn(ji,jj,jk,jp_fabm_m1+jp_fabm_p2p) / &
+                     &                                trn(ji,jj,jk,jp_fabm_m1+jp_fabm_chl2)
+               ENDIF
+               IF ( trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl3) > 0.0 ) THEN
+                  logchl_balinc_ersem_p3c(ji,jj,jk) = logchl_balinc_ersem_chl3(ji,jj,jk)   * &
+                     &                                trn(ji,jj,jk,jp_fabm_m1+jp_fabm_p3c) / &
+                     &                                trn(ji,jj,jk,jp_fabm_m1+jp_fabm_chl3)
+                  logchl_balinc_ersem_p3n(ji,jj,jk) = logchl_balinc_ersem_chl3(ji,jj,jk)   * &
+                     &                                trn(ji,jj,jk,jp_fabm_m1+jp_fabm_p3n) / &
+                     &                                trn(ji,jj,jk,jp_fabm_m1+jp_fabm_chl3)
+                  logchl_balinc_ersem_p3p(ji,jj,jk) = logchl_balinc_ersem_chl3(ji,jj,jk)   * &
+                     &                                trn(ji,jj,jk,jp_fabm_m1+jp_fabm_p3p) / &
+                     &                                trn(ji,jj,jk,jp_fabm_m1+jp_fabm_chl3)
+               ENDIF
+               IF ( trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl4) > 0.0 ) THEN
+                  logchl_balinc_ersem_p4c(ji,jj,jk) = logchl_balinc_ersem_chl4(ji,jj,jk)   * &
+                     &                                trn(ji,jj,jk,jp_fabm_m1+jp_fabm_p4c) / &
+                     &                                trn(ji,jj,jk,jp_fabm_m1+jp_fabm_chl4)
+                  logchl_balinc_ersem_p4n(ji,jj,jk) = logchl_balinc_ersem_chl4(ji,jj,jk)   * &
+                     &                                trn(ji,jj,jk,jp_fabm_m1+jp_fabm_p4n) / &
+                     &                                trn(ji,jj,jk,jp_fabm_m1+jp_fabm_chl4)
+                  logchl_balinc_ersem_p4p(ji,jj,jk) = logchl_balinc_ersem_chl4(ji,jj,jk)   * &
+                     &                                trn(ji,jj,jk,jp_fabm_m1+jp_fabm_p4p) / &
+                     &                                trn(ji,jj,jk,jp_fabm_m1+jp_fabm_chl4)
+               ENDIF
+            END DO
+         END DO
+      END DO
+
+   END SUBROUTINE asm_logchl_bal_ersem
+
+#else
+   !!----------------------------------------------------------------------
+   !!   Default option : Empty routine
+   !!----------------------------------------------------------------------
+CONTAINS
+   SUBROUTINE asm_logchl_bal_ersem( ln_logchlpftinc,          npfts,                    &
+      &                             mld_choice_bgc,           logchl_bkginc,            &      
+      &                             logchl_balinc_ersem_chl1, logchl_balinc_ersem_chl2, &
+      &                             logchl_balinc_ersem_chl3, logchl_balinc_ersem_chl4, &
+      &                             logchl_balinc_ersem_p1c,  logchl_balinc_ersem_p1n,  &
+      &                             logchl_balinc_ersem_p1p,  logchl_balinc_ersem_p1s,  &
+      &                             logchl_balinc_ersem_p2c,  logchl_balinc_ersem_p2n,  &
+      &                             logchl_balinc_ersem_p2p,  logchl_balinc_ersem_p3c,  &
+      &                             logchl_balinc_ersem_p3n,  logchl_balinc_ersem_p3p,  &
+      &                             logchl_balinc_ersem_p4c,  logchl_balinc_ersem_p4n,  &
+      &                             logchl_balinc_ersem_p4p,  logchl_balinc_ersem_z4c,  &
+      &                             logchl_balinc_ersem_z5c,  logchl_balinc_ersem_z5n,  &
+      &                             logchl_balinc_ersem_z5p,  logchl_balinc_ersem_z6c,  &
+      &                             logchl_balinc_ersem_z6n,  logchl_balinc_ersem_z6p,  &
+      &                             logchl_balinc_ersem_n1p,  logchl_balinc_ersem_n3n,  &
+      &                             logchl_balinc_ersem_n4n,  logchl_balinc_ersem_n5s )
+      LOGICAL :: ln_logchlpftinc
+      INTEGER :: npfts
+      INTEGER :: mld_choice_bgc
+      REAL    :: logchl_bkginc(:,:,:)
+      REAL    :: logchl_balinc_ersem_chl1(:,:,:), logchl_balinc_ersem_chl2(:,:,:)
+      REAL    :: logchl_balinc_ersem_chl3(:,:,:), logchl_balinc_ersem_chl4(:,:,:)
+      REAL    :: logchl_balinc_ersem_p1c(:,:,:),  logchl_balinc_ersem_p1n(:,:,:)
+      REAL    :: logchl_balinc_ersem_p1p(:,:,:),  logchl_balinc_ersem_p1s(:,:,:)
+      REAL    :: logchl_balinc_ersem_p2c(:,:,:),  logchl_balinc_ersem_p2n(:,:,:)
+      REAL    :: logchl_balinc_ersem_p2p(:,:,:),  logchl_balinc_ersem_p3c(:,:,:)
+      REAL    :: logchl_balinc_ersem_p3n(:,:,:),  logchl_balinc_ersem_p3p(:,:,:)
+      REAL    :: logchl_balinc_ersem_p4c(:,:,:),  logchl_balinc_ersem_p4n(:,:,:)
+      REAL    :: logchl_balinc_ersem_p4p(:,:,:),  logchl_balinc_ersem_z4c(:,:,:)
+      REAL    :: logchl_balinc_ersem_z5c(:,:,:),  logchl_balinc_ersem_z5n(:,:,:)
+      REAL    :: logchl_balinc_ersem_z5p(:,:,:),  logchl_balinc_ersem_z6c(:,:,:)
+      REAL    :: logchl_balinc_ersem_z6n(:,:,:),  logchl_balinc_ersem_z6p(:,:,:)
+      REAL    :: logchl_balinc_ersem_n1p(:,:,:),  logchl_balinc_ersem_n3n(:,:,:)
+      REAL    :: logchl_balinc_ersem_n4n(:,:,:),  logchl_balinc_ersem_n5s )
+      WRITE(*,*) 'asm_logchl_bal_ersem: You should not have seen this print! error?'
+   END SUBROUTINE asm_logchl_bal_ersem
+#endif
+
+   !!======================================================================
+END MODULE asmlogchlbal_ersem
Index: /branches/UKMO/dev_r5518_v3.4_asm_nemovar_community_bgc_ersem/NEMOGCM/NEMO/OPA_SRC/ASM/asmpar.F90
===================================================================
--- /branches/UKMO/dev_r5518_v3.4_asm_nemovar_community_bgc_ersem/NEMOGCM/NEMO/OPA_SRC/ASM/asmpar.F90	(revision 6982)
+++ /branches/UKMO/dev_r5518_v3.4_asm_nemovar_community_bgc_ersem/NEMOGCM/NEMO/OPA_SRC/ASM/asmpar.F90	(revision 6983)
@@ -20,7 +20,10 @@
       & c_asmtrj = 'assim_trj',                  & !: Filename for storing the 
                                                    !: reference trajectory
-      & c_asminc = 'assim_background_increments'   !: Filename for storing the 
+      & c_asminc = 'assim_background_increments', & !: Filename for storing the 
                                                    !: increments to the background
                                                    !: state
+      & c_asmbal = 'assim.balincs'                 !: Filename for storing the 
+                                                   !: balancing increments calculated
+                                                   !: for biogeochemistry
 
    INTEGER, PUBLIC :: nitbkg_r      !: Background time step referenced to nit000
Index: /branches/UKMO/dev_r5518_v3.4_asm_nemovar_community_bgc_ersem/NEMOGCM/NEMO/OPA_SRC/nemogcm.F90
===================================================================
--- /branches/UKMO/dev_r5518_v3.4_asm_nemovar_community_bgc_ersem/NEMOGCM/NEMO/OPA_SRC/nemogcm.F90	(revision 6982)
+++ /branches/UKMO/dev_r5518_v3.4_asm_nemovar_community_bgc_ersem/NEMOGCM/NEMO/OPA_SRC/nemogcm.F90	(revision 6983)
@@ -61,4 +61,5 @@
    USE asminc          ! assimilation increments     
    USE asmbkg          ! writing out state trajectory
+   USE asmbal          ! writing out assimilation balancing increments
    USE diaptr          ! poleward transports           (dia_ptr_init routine)
    USE diadct          ! sections transports           (dia_dct_init routine)
@@ -158,4 +159,5 @@
                 IF( ln_dyninc ) CALL dyn_asm_inc( nit000 - 1 )    ! Dynamics
                 IF( ln_sshinc ) CALL ssh_asm_inc( nit000 - 1 )    ! SSH
+                IF( ln_logchltotinc .OR. ln_logchlpftinc ) CALL logchl_asm_inc( nit000 - 1 )
              ENDIF
           ENDIF
@@ -177,4 +179,6 @@
 
       IF( lk_diaobs   )   CALL dia_obs_wri
+      !
+      IF( ( lk_asminc ).AND.( ln_balwri ) ) CALL asm_bal_wri( nitend )  ! Output balancing increments
       !
       IF( ln_icebergs )   CALL icb_end( nitend )
Index: /branches/UKMO/dev_r5518_v3.4_asm_nemovar_community_bgc_ersem/NEMOGCM/NEMO/OPA_SRC/step.F90
===================================================================
--- /branches/UKMO/dev_r5518_v3.4_asm_nemovar_community_bgc_ersem/NEMOGCM/NEMO/OPA_SRC/step.F90	(revision 6982)
+++ /branches/UKMO/dev_r5518_v3.4_asm_nemovar_community_bgc_ersem/NEMOGCM/NEMO/OPA_SRC/step.F90	(revision 6983)
@@ -242,4 +242,6 @@
       ! Passive Tracer Model
       !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
+      IF( lk_asminc .AND. ln_asmiau .AND. ( ln_logchltotinc .OR. ln_logchlpftinc ) ) &
+         &               CALL logchl_asm_inc( kstp )  ! logchl assimilation
                          CALL trc_stp( kstp )         ! time-stepping
 #endif
