Ignore:
Timestamp:
2019-02-01T17:27:20+01:00 (18 months ago)
Author:
dford
Message:

Implement biogeochemistry assimilation for FABM-ERSEM.

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/UKMO/AMM15_v3_6_STABLE_package_collate_BGC_DA/NEMOGCM/NEMO/OPA_SRC/ASM/asmbgc.F90

    r10574 r10622  
    1111   !! 'key_medusa'          : MEDUSA model 
    1212   !! 'key_roam'            : MEDUSA extras for carbonate cycle 
     13   !! 'key_fabm'            : ERSEM model coupled via FABM 
    1314   !! 'key_karaml'          : Kara mixed layer depth 
    1415   !!--------------------------------------------------------------------------- 
     
    3839      & ln_kara,            & 
    3940#endif    
    40       & hmld,               &  
     41      & hmld,               & 
     42      & hmld_tref,          & 
    4143      & hmlp,               & 
    4244      & hmlpt  
     
    4547      & c_asmbal,           & 
    4648      & nitbkg_r,           & 
     49      & nitavgbkg_r,        & 
    4750      & nitdin_r,           & 
    4851      & nitiaustr_r,        & 
    4952      & nitiaufin_r 
     53   USE asmphyto2dbal_medusa, ONLY: & ! phyto2d balancing for MEDUSA 
     54      & asm_phyto2d_bal_medusa 
     55   USE asmphyto2dbal_hadocc, ONLY: & ! phyto2d balancing for HadOCC 
     56      & asm_phyto2d_bal_hadocc 
     57   USE asmphyto2dbal_ersem,  ONLY: & ! phyto2d balancing for ERSEM 
     58      & asm_phyto2d_bal_ersem 
     59   USE asmpco2bal, ONLY:    & ! pCO2 balancing 
     60      & asm_pco2_bal 
    5061#if defined key_top 
    5162   USE trc, ONLY:           & ! passive tracer variables 
    5263      & trn,                & 
    53       & trb 
     64      & trb,                & 
     65      & nittrc000 
    5466   USE par_trc, ONLY:       & ! passive tracer parameters 
    5567      & jptra 
    5668#endif 
    5769#if defined key_medusa 
    58    USE asmphyto2dbal_medusa, ONLY: & ! phyto2d balancing for MEDUSA 
    59       & asm_phyto2d_bal_medusa 
    60    USE asmpco2bal, ONLY:    & ! pCO2 balancing for MEDUSA 
    61       & asm_pco2_bal 
    6270   USE sms_medusa             ! MEDUSA parameters 
    6371   USE par_medusa             ! MEDUSA parameters 
     
    7280      & mld_max 
    7381#elif defined key_hadocc 
    74    USE asmphyto2dbal_hadocc, ONLY: & ! phyto2d balancing for HadOCC 
    75       & asm_phyto2d_bal_hadocc 
    76    USE asmpco2bal, ONLY:    & ! pCO2 balancing for HadOCC 
    77       & asm_pco2_bal 
    7882   USE par_hadocc             ! HadOCC parameters 
    7983   USE had_bgc_const          ! HadOCC parameters 
     
    8690   USE had_bgc_const, ONLY: & ! HadOCC C:Chl ratio 
    8791      & cchl_p 
     92#elif defined key_fabm 
     93   USE par_fabm               ! FABM-ERSEM parameters 
     94   USE trc, ONLY:           & ! FABM-ERSEM diagnostic variables 
     95      & pgrow_avg,          & 
     96      & ploss_avg,          & 
     97      & phyt_avg,           & 
     98      & mld_max 
    8899#endif 
    89100 
     
    97108   PUBLIC  asm_bgc_init_bkg       ! called by asm_inc_init in asminc.F90 
    98109   PUBLIC  asm_bgc_bal_wri        ! called by nemo_gcm in nemogcm.F90 
     110   PUBLIC  asm_bgc_bkg_alloc      ! called by asm_bkg_wri in asmbkg.F90 
     111   PUBLIC  asm_bgc_bkg_tavg       ! called by asm_bkg_wri in asmbkg.F90 
    99112   PUBLIC  asm_bgc_bkg_wri        ! called by asm_bkg_wri in asmbkg.F90 
    100113   PUBLIC  phyto2d_asm_inc        ! called by bgc_asm_inc in asminc.F90 
     
    109122   LOGICAL, PUBLIC :: ln_slchldiainc = .FALSE. !: No surface diatom     log10(chlorophyll) increment 
    110123   LOGICAL, PUBLIC :: ln_slchlnoninc = .FALSE. !: No surface non-diatom log10(chlorophyll) increment 
     124   LOGICAL, PUBLIC :: ln_slchlnaninc = .FALSE. !: No surface nanophyto  log10(chlorophyll) increment 
     125   LOGICAL, PUBLIC :: ln_slchlpicinc = .FALSE. !: No surface picophyto  log10(chlorophyll) increment 
     126   LOGICAL, PUBLIC :: ln_slchldininc = .FALSE. !: No surface dinoflag   log10(chlorophyll) increment 
    111127   LOGICAL, PUBLIC :: ln_schltotinc  = .FALSE. !: No surface total      chlorophyll        increment 
    112128   LOGICAL, PUBLIC :: ln_slphytotinc = .FALSE. !: No surface total      log10(phyto C)     increment 
     
    146162   REAL(wp), DIMENSION(:,:),     ALLOCATABLE :: slchldia_bkginc ! slchldia inc 
    147163   REAL(wp), DIMENSION(:,:),     ALLOCATABLE :: slchlnon_bkginc ! slchlnon inc 
     164   REAL(wp), DIMENSION(:,:),     ALLOCATABLE :: slchlnan_bkginc ! slchlnan inc 
     165   REAL(wp), DIMENSION(:,:),     ALLOCATABLE :: slchlpic_bkginc ! slchlpic inc 
     166   REAL(wp), DIMENSION(:,:),     ALLOCATABLE :: slchldin_bkginc ! slchldin inc 
    148167   REAL(wp), DIMENSION(:,:),     ALLOCATABLE :: schltot_bkginc  ! schltot inc 
    149168   REAL(wp), DIMENSION(:,:),     ALLOCATABLE :: slphytot_bkginc ! slphytot inc 
     
    172191   REAL(wp), DIMENSION(:,:,:),   ALLOCATABLE :: chl_bkg         ! Background chl 
    173192   REAL(wp), DIMENSION(:,:,:),   ALLOCATABLE :: cchl_p_bkg      ! Background c:chl 
     193#if defined key_fabm 
     194   REAL(wp), DIMENSION(:,:,:),   ALLOCATABLE :: totalk_bkg      ! Background total alkalinity 
     195#endif 
     196 
     197   REAL(wp), DIMENSION(:,:),     SAVE, ALLOCATABLE :: pgrow_avg_tavg  ! Time average of pgrow_avg 
     198   REAL(wp), DIMENSION(:,:),     SAVE, ALLOCATABLE :: ploss_avg_tavg  ! Time average of ploss_avg 
     199   REAL(wp), DIMENSION(:,:),     SAVE, ALLOCATABLE :: phyt_avg_tavg   ! Time average of phyt_avg 
     200   REAL(wp), DIMENSION(:,:,:,:), SAVE, ALLOCATABLE :: trn_tavg        ! Time average of trn 
     201#if defined key_hadocc 
     202   REAL(wp), DIMENSION(:,:,:),   SAVE, ALLOCATABLE :: HADOCC_CHL_tavg ! Time average of HADOCC_CHL 
     203   REAL(wp), DIMENSION(:,:,:),   SAVE, ALLOCATABLE :: cchl_p_tavg     ! Time average of cchl_p 
     204#elif defined key_fabm 
     205   REAL(wp), DIMENSION(:,:,:),   SAVE, ALLOCATABLE :: totalk_tavg     ! Time average of O3_TA 
     206#endif 
    174207 
    175208#  include "domzgr_substitute.h90" 
     
    190223      !!------------------------------------------------------------------------ 
    191224 
    192 #if ! defined key_top || ( ! defined key_hadocc && ! defined key_medusa ) 
     225#if ! defined key_top || ( ! defined key_hadocc && ! defined key_medusa && ! defined key_fabm ) 
    193226      CALL ctl_stop( ' Attempting to assimilate biogeochemical observations', & 
    194227         &           ' but no compatible biogeochemical model is available' ) 
     
    208241#endif 
    209242 
     243#if defined key_fabm 
     244      IF ( ln_pphinc ) THEN 
     245         CALL ctl_stop( ' Cannot currently assimilate pH into FABM-ERSEM' ) 
     246      ENDIF 
     247#endif 
     248 
    210249      IF ( ( ln_phytobal ).AND.                                       & 
    211250         & ( .NOT. ln_slchltotinc ).AND.( .NOT. ln_slchldiainc ).AND. & 
    212          & ( .NOT. ln_slchlnoninc ).AND.( .NOT. ln_schltotinc  ).AND. & 
    213          & ( .NOT. ln_slphytotinc ).AND.( .NOT. ln_slphydiainc ).AND. & 
    214          & ( .NOT. ln_slphynoninc ) ) THEN 
     251         & ( .NOT. ln_slchlnoninc ).AND.( .NOT. ln_slchlnaninc ).AND. & 
     252         & ( .NOT. ln_slchlpicinc ).AND.( .NOT. ln_slchldininc ).AND. & 
     253         & ( .NOT. ln_schltotinc  ).AND.( .NOT. ln_slphytotinc ).AND. & 
     254         & ( .NOT. ln_slphydiainc ).AND.( .NOT. ln_slphynoninc ) ) THEN 
    215255         CALL ctl_warn( ' Cannot calculate phytoplankton balancing increments', & 
    216256            &           ' if not assimilating ocean colour,',                   & 
     
    221261      IF ( ( ln_balwri ).AND.                                         & 
    222262         & ( .NOT. ln_slchltotinc ).AND.( .NOT. ln_slchldiainc ).AND. & 
    223          & ( .NOT. ln_slchlnoninc ).AND.( .NOT. ln_schltotinc  ).AND. & 
     263         & ( .NOT. ln_slchlnoninc ).AND.( .NOT. ln_slchlnaninc ).AND. & 
     264         & ( .NOT. ln_slchlpicinc ).AND.( .NOT. ln_slchldininc ).AND. & 
     265         & ( .NOT. ln_schltotinc  ).AND.                              & 
    224266         & ( .NOT. ln_slphytotinc ).AND.( .NOT. ln_slphydiainc ).AND. & 
    225267         & ( .NOT. ln_slphynoninc ).AND.( .NOT. ln_plchltotinc ).AND. & 
     
    244286 
    245287      IF ( ( ln_slchltotinc .OR. ln_schltotinc  ) .AND. & 
    246          & ( ln_slchldiainc .OR. ln_slchlnoninc ) ) THEN 
     288         & ( ln_slchldiainc .OR. ln_slchlnoninc .OR. ln_slchlnaninc .OR. & 
     289         &   ln_slchlpicinc .OR. ln_slchldininc ) ) THEN 
    247290         CALL ctl_stop( ' Can only assimilate total or PFT surface chlorophyll, not both' ) 
    248291      ENDIF 
     
    250293      IF ( ln_phytobal .AND.                                      & 
    251294         & ( ( ln_slchlnoninc .AND.( .NOT. ln_slchldiainc ) ).OR. & 
     295#if defined key_fabm 
     296         &   (  ( ln_slchldiainc .OR. ln_slchlnaninc .OR. ln_slchlpicinc .OR. ln_slchldininc ) .AND. & 
     297         &      ( ( .NOT. ln_slchldiainc ) .OR. ( .NOT. ln_slchlnaninc ) .OR. & 
     298         &        ( .NOT. ln_slchlpicinc ) .OR. ( .NOT. ln_slchldininc ) ) ) .OR. & 
     299#endif 
    252300         &   ( ln_slchldiainc .AND.( .NOT. ln_slchlnoninc ) ) ) ) THEN 
    253301         CALL ctl_stop( ' Cannot calculate phytoplankton balancing increments', & 
     
    306354      ENDIF 
    307355       
     356      IF ( ln_slchlnaninc ) THEN 
     357         ALLOCATE( slchlnan_bkginc(jpi,jpj) ) 
     358         CALL asm_bgc_read_incs_2d( knum, 'bckinslchlnan', slchlnan_bkginc ) 
     359      ENDIF 
     360       
     361      IF ( ln_slchlpicinc ) THEN 
     362         ALLOCATE( slchlpic_bkginc(jpi,jpj) ) 
     363         CALL asm_bgc_read_incs_2d( knum, 'bckinslchlpic', slchlpic_bkginc ) 
     364      ENDIF 
     365       
     366      IF ( ln_slchldininc ) THEN 
     367         ALLOCATE( slchldin_bkginc(jpi,jpj) ) 
     368         CALL asm_bgc_read_incs_2d( knum, 'bckinslchldin', slchldin_bkginc ) 
     369      ENDIF 
     370       
    308371      IF ( ln_schltotinc ) THEN 
    309372         ALLOCATE( schltot_bkginc(jpi,jpj) ) 
     
    379442       
    380443      IF ( ln_slchltotinc .OR. ln_slchldiainc .OR. ln_slchlnoninc .OR. & 
     444         & ln_slchlnaninc .OR. ln_slchlpicinc .OR. ln_slchldininc .OR. & 
    381445         & ln_schltotinc  .OR. ln_slphytotinc .OR. ln_slphydiainc .OR. & 
    382446         & ln_slphynoninc ) THEN 
     
    516580      !!------------------------------------------------------------------------ 
    517581 
    518 #if defined key_top && ( defined key_hadocc || defined key_medusa ) 
     582#if defined key_top && ( defined key_hadocc || defined key_medusa || defined key_fabm ) 
    519583      IF ( ln_slchltotinc .OR. ln_slchldiainc .OR. ln_slchlnoninc .OR. & 
     584         & ln_slchlnaninc .OR. ln_slchlpicinc .OR. ln_slchldininc .OR. & 
    520585         & ln_schltotinc  .OR. ln_slphytotinc .OR. ln_slphydiainc .OR. & 
    521586         & ln_slphynoninc .OR. ln_plchltotinc .OR. ln_pchltotinc ) THEN 
     
    557622         CALL iom_get( inum, jpdom_autoglo, 'medusa_phd', tracer_bkg(:,:,:,jpphd) ) 
    558623         CALL iom_get( inum, jpdom_autoglo, 'medusa_pds', tracer_bkg(:,:,:,jppds) ) 
     624#elif defined key_fabm 
     625         CALL iom_get( inum, jpdom_autoglo, 'ersem_chl1', tracer_bkg(:,:,:,jp_fabm_chl1) ) 
     626         CALL iom_get( inum, jpdom_autoglo, 'ersem_chl2', tracer_bkg(:,:,:,jp_fabm_chl2) ) 
     627         CALL iom_get( inum, jpdom_autoglo, 'ersem_chl3', tracer_bkg(:,:,:,jp_fabm_chl3) ) 
     628         CALL iom_get( inum, jpdom_autoglo, 'ersem_chl4', tracer_bkg(:,:,:,jp_fabm_chl4) ) 
     629         CALL iom_get( inum, jpdom_autoglo, 'ersem_p1c',  tracer_bkg(:,:,:,jp_fabm_p1c)  ) 
     630         CALL iom_get( inum, jpdom_autoglo, 'ersem_p1n',  tracer_bkg(:,:,:,jp_fabm_p1n)  ) 
     631         CALL iom_get( inum, jpdom_autoglo, 'ersem_p1p',  tracer_bkg(:,:,:,jp_fabm_p1p)  ) 
     632         CALL iom_get( inum, jpdom_autoglo, 'ersem_p1s',  tracer_bkg(:,:,:,jp_fabm_p1s)  ) 
     633         CALL iom_get( inum, jpdom_autoglo, 'ersem_p2c',  tracer_bkg(:,:,:,jp_fabm_p2c)  ) 
     634         CALL iom_get( inum, jpdom_autoglo, 'ersem_p2n',  tracer_bkg(:,:,:,jp_fabm_p2n)  ) 
     635         CALL iom_get( inum, jpdom_autoglo, 'ersem_p2p',  tracer_bkg(:,:,:,jp_fabm_p2p)  ) 
     636         CALL iom_get( inum, jpdom_autoglo, 'ersem_p3c',  tracer_bkg(:,:,:,jp_fabm_p3c)  ) 
     637         CALL iom_get( inum, jpdom_autoglo, 'ersem_p3n',  tracer_bkg(:,:,:,jp_fabm_p3n)  ) 
     638         CALL iom_get( inum, jpdom_autoglo, 'ersem_p3p',  tracer_bkg(:,:,:,jp_fabm_p3p)  ) 
     639         CALL iom_get( inum, jpdom_autoglo, 'ersem_p4c',  tracer_bkg(:,:,:,jp_fabm_p4c)  ) 
     640         CALL iom_get( inum, jpdom_autoglo, 'ersem_p4n',  tracer_bkg(:,:,:,jp_fabm_p4n)  ) 
     641         CALL iom_get( inum, jpdom_autoglo, 'ersem_p4p',  tracer_bkg(:,:,:,jp_fabm_p4p)  ) 
    559642#endif 
    560643          
     
    588671            CALL iom_get( inum, jpdom_autoglo, 'medusa_alk', tracer_bkg(:,:,:,jpalk) ) 
    589672            CALL iom_get( inum, jpdom_autoglo, 'medusa_oxy', tracer_bkg(:,:,:,jpoxy) ) 
     673#elif defined key_fabm 
     674            CALL iom_get( inum, jpdom_autoglo, 'ersem_z4c',  tracer_bkg(:,:,:,jp_fabm_z4c)  ) 
     675            CALL iom_get( inum, jpdom_autoglo, 'ersem_z5c',  tracer_bkg(:,:,:,jp_fabm_z5c)  ) 
     676            CALL iom_get( inum, jpdom_autoglo, 'ersem_z5n',  tracer_bkg(:,:,:,jp_fabm_z5n)  ) 
     677            CALL iom_get( inum, jpdom_autoglo, 'ersem_z5p',  tracer_bkg(:,:,:,jp_fabm_z5p)  ) 
     678            CALL iom_get( inum, jpdom_autoglo, 'ersem_z6c',  tracer_bkg(:,:,:,jp_fabm_z6c)  ) 
     679            CALL iom_get( inum, jpdom_autoglo, 'ersem_z6n',  tracer_bkg(:,:,:,jp_fabm_z6n)  ) 
     680            CALL iom_get( inum, jpdom_autoglo, 'ersem_z6p',  tracer_bkg(:,:,:,jp_fabm_z6p)  ) 
     681            CALL iom_get( inum, jpdom_autoglo, 'ersem_n1p',  tracer_bkg(:,:,:,jp_fabm_n1p)  ) 
     682            CALL iom_get( inum, jpdom_autoglo, 'ersem_n3n',  tracer_bkg(:,:,:,jp_fabm_n3n)  ) 
     683            CALL iom_get( inum, jpdom_autoglo, 'ersem_n4n',  tracer_bkg(:,:,:,jp_fabm_n4n)  ) 
     684            CALL iom_get( inum, jpdom_autoglo, 'ersem_n5s',  tracer_bkg(:,:,:,jp_fabm_n5s)  ) 
     685            CALL iom_get( inum, jpdom_autoglo, 'ersem_o2o',  tracer_bkg(:,:,:,jp_fabm_o2o)  ) 
     686            CALL iom_get( inum, jpdom_autoglo, 'ersem_o3c',  tracer_bkg(:,:,:,jp_fabm_o3c)  ) 
     687            CALL iom_get( inum, jpdom_autoglo, 'ersem_o3ba', tracer_bkg(:,:,:,jp_fabm_o3ba) ) 
     688            CALL iom_get( inum, jpdom_autoglo, 'ersem_o3ta', totalk_bkg(:,:,:)              ) 
     689            totalk_bkg(:,:,:) = totalk_bkg(:,:,:) * tmask(:,:,:) 
    590690#endif 
    591691         ELSE IF ( ln_spco2inc .OR. ln_sfco2inc .OR. ln_pphinc ) THEN 
     
    596696            CALL iom_get( inum, jpdom_autoglo, 'medusa_dic', tracer_bkg(:,:,:,jpdic) ) 
    597697            CALL iom_get( inum, jpdom_autoglo, 'medusa_alk', tracer_bkg(:,:,:,jpalk) ) 
     698#elif defined key_fabm 
     699            CALL iom_get( inum, jpdom_autoglo, 'ersem_o3c',  tracer_bkg(:,:,:,jp_fabm_o3c)  ) 
     700            CALL iom_get( inum, jpdom_autoglo, 'ersem_o3ba', tracer_bkg(:,:,:,jp_fabm_o3ba) ) 
     701            CALL iom_get( inum, jpdom_autoglo, 'ersem_o3ta', totalk_bkg(:,:,:)              ) 
     702            totalk_bkg(:,:,:) = totalk_bkg(:,:,:) * tmask(:,:,:) 
    598703#endif 
    599704            CALL iom_get( inum, jpdom_autoglo, 'mld_max',   mld_max_bkg   ) 
     
    622727         CALL iom_get( inum, jpdom_autoglo, 'medusa_dic', tracer_bkg(:,:,:,jpdic) ) 
    623728         CALL iom_get( inum, jpdom_autoglo, 'medusa_alk', tracer_bkg(:,:,:,jpalk) ) 
     729#elif defined key_fabm 
     730         CALL iom_get( inum, jpdom_autoglo, 'ersem_o3c',  tracer_bkg(:,:,:,jp_fabm_o3c)  ) 
     731         CALL iom_get( inum, jpdom_autoglo, 'ersem_o3ba', tracer_bkg(:,:,:,jp_fabm_o3ba) ) 
     732         CALL iom_get( inum, jpdom_autoglo, 'ersem_o3ta', totalk_bkg(:,:,:)              ) 
     733         totalk_bkg(:,:,:) = totalk_bkg(:,:,:) * tmask(:,:,:) 
    624734#endif 
    625735         CALL iom_get( inum, jpdom_autoglo, 'mld_max',   mld_max_bkg   ) 
     
    687797 
    688798         IF ( ln_slchltotinc .OR. ln_slchldiainc .OR. ln_slchlnoninc .OR. & 
     799            & ln_slchlnaninc .OR. ln_slchlpicinc .OR. ln_slchldininc .OR. & 
    689800            & ln_schltotinc  .OR. ln_slphytotinc .OR. ln_slphydiainc .OR. & 
    690801            & ln_slphynoninc ) THEN 
     
    716827               CALL iom_rstput( kt, kt, inum, 'phy2d_alk', phyto2d_balinc(:,:,:,jp_had_alk) ) 
    717828            ENDIF 
     829#elif defined key_fabm 
     830            CALL iom_rstput( kt, kt, inum, 'phy2d_chl1', phyto2d_balinc(:,:,:,jp_fabm_chl1) ) 
     831            CALL iom_rstput( kt, kt, inum, 'phy2d_chl2', phyto2d_balinc(:,:,:,jp_fabm_chl2) ) 
     832            CALL iom_rstput( kt, kt, inum, 'phy2d_chl3', phyto2d_balinc(:,:,:,jp_fabm_chl3) ) 
     833            CALL iom_rstput( kt, kt, inum, 'phy2d_chl4', phyto2d_balinc(:,:,:,jp_fabm_chl4) ) 
     834            CALL iom_rstput( kt, kt, inum, 'phy2d_p1c',  phyto2d_balinc(:,:,:,jp_fabm_p1c)  ) 
     835            CALL iom_rstput( kt, kt, inum, 'phy2d_p1n',  phyto2d_balinc(:,:,:,jp_fabm_p1n)  ) 
     836            CALL iom_rstput( kt, kt, inum, 'phy2d_p1p',  phyto2d_balinc(:,:,:,jp_fabm_p1p)  ) 
     837            CALL iom_rstput( kt, kt, inum, 'phy2d_p1s',  phyto2d_balinc(:,:,:,jp_fabm_p1s)  ) 
     838            CALL iom_rstput( kt, kt, inum, 'phy2d_p2c',  phyto2d_balinc(:,:,:,jp_fabm_p2c)  ) 
     839            CALL iom_rstput( kt, kt, inum, 'phy2d_p2n',  phyto2d_balinc(:,:,:,jp_fabm_p2n)  ) 
     840            CALL iom_rstput( kt, kt, inum, 'phy2d_p2p',  phyto2d_balinc(:,:,:,jp_fabm_p2p)  ) 
     841            CALL iom_rstput( kt, kt, inum, 'phy2d_p3c',  phyto2d_balinc(:,:,:,jp_fabm_p3c)  ) 
     842            CALL iom_rstput( kt, kt, inum, 'phy2d_p3n',  phyto2d_balinc(:,:,:,jp_fabm_p3n)  ) 
     843            CALL iom_rstput( kt, kt, inum, 'phy2d_p3p',  phyto2d_balinc(:,:,:,jp_fabm_p3p)  ) 
     844            CALL iom_rstput( kt, kt, inum, 'phy2d_p4c',  phyto2d_balinc(:,:,:,jp_fabm_p4c)  ) 
     845            CALL iom_rstput( kt, kt, inum, 'phy2d_p4n',  phyto2d_balinc(:,:,:,jp_fabm_p4n)  ) 
     846            CALL iom_rstput( kt, kt, inum, 'phy2d_p4p',  phyto2d_balinc(:,:,:,jp_fabm_p4p)  ) 
     847            IF ( ln_phytobal ) THEN 
     848               CALL iom_rstput( kt, kt, inum, 'phy2d_z4c',   phyto2d_balinc(:,:,:,jp_fabm_z4c)  ) 
     849               CALL iom_rstput( kt, kt, inum, 'phy2d_z5c',   phyto2d_balinc(:,:,:,jp_fabm_z5c)  ) 
     850               CALL iom_rstput( kt, kt, inum, 'phy2d_z5n',   phyto2d_balinc(:,:,:,jp_fabm_z5n)  ) 
     851               CALL iom_rstput( kt, kt, inum, 'phy2d_z5p',   phyto2d_balinc(:,:,:,jp_fabm_z5p)  ) 
     852               CALL iom_rstput( kt, kt, inum, 'phy2d_z6c',   phyto2d_balinc(:,:,:,jp_fabm_z6c)  ) 
     853               CALL iom_rstput( kt, kt, inum, 'phy2d_z6n',   phyto2d_balinc(:,:,:,jp_fabm_z6n)  ) 
     854               CALL iom_rstput( kt, kt, inum, 'phy2d_z6p',   phyto2d_balinc(:,:,:,jp_fabm_z6p)  ) 
     855               CALL iom_rstput( kt, kt, inum, 'phy2d_n1p',   phyto2d_balinc(:,:,:,jp_fabm_n1p)  ) 
     856               CALL iom_rstput( kt, kt, inum, 'phy2d_n3n',   phyto2d_balinc(:,:,:,jp_fabm_n3n)  ) 
     857               CALL iom_rstput( kt, kt, inum, 'phy2d_n4n',   phyto2d_balinc(:,:,:,jp_fabm_n4n)  ) 
     858               CALL iom_rstput( kt, kt, inum, 'phy2d_n5s',   phyto2d_balinc(:,:,:,jp_fabm_n5s)  ) 
     859               CALL iom_rstput( kt, kt, inum, 'phy2d_o2o',   phyto2d_balinc(:,:,:,jp_fabm_o2o)  ) 
     860               CALL iom_rstput( kt, kt, inum, 'phy2d_o3c',   phyto2d_balinc(:,:,:,jp_fabm_o3c)  ) 
     861               CALL iom_rstput( kt, kt, inum, 'phy2d_o3ba',  phyto2d_balinc(:,:,:,jp_fabm_o3ba) ) 
     862            ENDIF 
    718863#endif 
    719864         ENDIF 
     
    728873#elif defined key_hadocc 
    729874            CALL iom_rstput( kt, kt, inum, 'phy3d_phy', phyto3d_balinc(:,:,:,jp_had_phy) ) 
     875#elif defined key_fabm 
     876            CALL iom_rstput( kt, kt, inum, 'phy3d_chl1', phyto3d_balinc(:,:,:,jp_fabm_chl1) ) 
     877            CALL iom_rstput( kt, kt, inum, 'phy3d_chl2', phyto3d_balinc(:,:,:,jp_fabm_chl2) ) 
     878            CALL iom_rstput( kt, kt, inum, 'phy3d_chl3', phyto3d_balinc(:,:,:,jp_fabm_chl3) ) 
     879            CALL iom_rstput( kt, kt, inum, 'phy3d_chl4', phyto3d_balinc(:,:,:,jp_fabm_chl4) ) 
     880            CALL iom_rstput( kt, kt, inum, 'phy3d_p1c',  phyto3d_balinc(:,:,:,jp_fabm_p1c)  ) 
     881            CALL iom_rstput( kt, kt, inum, 'phy3d_p1n',  phyto3d_balinc(:,:,:,jp_fabm_p1n)  ) 
     882            CALL iom_rstput( kt, kt, inum, 'phy3d_p1p',  phyto3d_balinc(:,:,:,jp_fabm_p1p)  ) 
     883            CALL iom_rstput( kt, kt, inum, 'phy3d_p1s',  phyto3d_balinc(:,:,:,jp_fabm_p1s)  ) 
     884            CALL iom_rstput( kt, kt, inum, 'phy3d_p2c',  phyto3d_balinc(:,:,:,jp_fabm_p2c)  ) 
     885            CALL iom_rstput( kt, kt, inum, 'phy3d_p2n',  phyto3d_balinc(:,:,:,jp_fabm_p2n)  ) 
     886            CALL iom_rstput( kt, kt, inum, 'phy3d_p2p',  phyto3d_balinc(:,:,:,jp_fabm_p2p)  ) 
     887            CALL iom_rstput( kt, kt, inum, 'phy3d_p3c',  phyto3d_balinc(:,:,:,jp_fabm_p3c)  ) 
     888            CALL iom_rstput( kt, kt, inum, 'phy3d_p3n',  phyto3d_balinc(:,:,:,jp_fabm_p3n)  ) 
     889            CALL iom_rstput( kt, kt, inum, 'phy3d_p3p',  phyto3d_balinc(:,:,:,jp_fabm_p3p)  ) 
     890            CALL iom_rstput( kt, kt, inum, 'phy3d_p4c',  phyto3d_balinc(:,:,:,jp_fabm_p4c)  ) 
     891            CALL iom_rstput( kt, kt, inum, 'phy3d_p4n',  phyto3d_balinc(:,:,:,jp_fabm_p4n)  ) 
     892            CALL iom_rstput( kt, kt, inum, 'phy3d_p4p',  phyto3d_balinc(:,:,:,jp_fabm_p4p)  ) 
    730893#endif 
    731894         ENDIF 
     
    738901            CALL iom_rstput( kt, kt, inum, 'pco2_dic', pco2_balinc(:,:,:,jp_had_dic) ) 
    739902            CALL iom_rstput( kt, kt, inum, 'pco2_alk', pco2_balinc(:,:,:,jp_had_alk) ) 
     903#elif defined key_fabm 
     904            CALL iom_rstput( kt, kt, inum, 'pco2_o3c',  pco2_balinc(:,:,:,jp_fabm_o3c)  ) 
     905            CALL iom_rstput( kt, kt, inum, 'pco2_o3ba', pco2_balinc(:,:,:,jp_fabm_o3ba) ) 
    740906#endif 
    741907         ELSE IF ( ln_sfco2inc ) THEN 
     
    746912            CALL iom_rstput( kt, kt, inum, 'fco2_dic', pco2_balinc(:,:,:,jp_had_dic) ) 
    747913            CALL iom_rstput( kt, kt, inum, 'fco2_alk', pco2_balinc(:,:,:,jp_had_alk) ) 
     914#elif defined key_fabm 
     915            CALL iom_rstput( kt, kt, inum, 'fco2_o3c',  pco2_balinc(:,:,:,jp_fabm_o3c)  ) 
     916            CALL iom_rstput( kt, kt, inum, 'fco2_o3ba', pco2_balinc(:,:,:,jp_fabm_o3ba) ) 
    748917#endif 
    749918         ENDIF 
     
    756925            CALL iom_rstput( kt, kt, inum, 'ph_dic', ph_balinc(:,:,:,jp_had_dic) ) 
    757926            CALL iom_rstput( kt, kt, inum, 'ph_alk', ph_balinc(:,:,:,jp_had_alk) ) 
     927#elif defined key_fabm 
     928            CALL iom_rstput( kt, kt, inum, 'ph_o3c',  ph_balinc(:,:,:,jp_fabm_o3c)  ) 
     929            CALL iom_rstput( kt, kt, inum, 'ph_o3ba', ph_balinc(:,:,:,jp_fabm_o3ba) ) 
    758930#endif 
    759931         ENDIF 
     
    774946   !!=========================================================================== 
    775947 
    776    SUBROUTINE asm_bgc_bkg_wri( kt, knum ) 
     948   SUBROUTINE asm_bgc_bkg_alloc 
     949      !!------------------------------------------------------------------------ 
     950      !!                    ***  ROUTINE asm_bgc_bkg_alloc  *** 
     951      !! 
     952      !! ** Purpose :   allocate time-average arrays for background 
     953      !! 
     954      !! ** Method  :   allocate time-average arrays for background 
     955      !! 
     956      !! ** Action  :   allocate time-average arrays for background 
     957      !! 
     958      !! References :   asm_bkg_wri 
     959      !!------------------------------------------------------------------------ 
     960      !! 
     961      INTEGER :: ierror 
     962      !! 
     963      !!------------------------------------------------------------------------ 
     964       
     965      ALLOCATE( pgrow_avg_tavg(jpi,jpj), STAT=ierror ) 
     966      IF( ierror > 0 ) THEN 
     967         CALL ctl_stop( 'asm_bgc_bkg_alloc: unable to allocate pgrow_avg_tavg' ) 
     968      ENDIF 
     969      pgrow_avg_tavg(:,:) = 0.0 
     970       
     971      ALLOCATE( ploss_avg_tavg(jpi,jpj), STAT=ierror ) 
     972      IF( ierror > 0 ) THEN 
     973         CALL ctl_stop( 'asm_bgc_bkg_alloc: unable to allocate ploss_avg_tavg' ) 
     974      ENDIF 
     975      ploss_avg_tavg(:,:) = 0.0 
     976       
     977      ALLOCATE( phyt_avg_tavg(jpi,jpj), STAT=ierror ) 
     978      IF( ierror > 0 ) THEN 
     979         CALL ctl_stop( 'asm_bgc_bkg_alloc: unable to allocate phyt_avg_tavg' ) 
     980      ENDIF 
     981      phyt_avg_tavg(:,:) = 0.0 
     982       
     983      ALLOCATE( trn_tavg(jpi,jpj,jpk,jptra), STAT=ierror ) 
     984      IF( ierror > 0 ) THEN 
     985         CALL ctl_stop( 'asm_bgc_bkg_alloc: unable to allocate trn_tavg' ) 
     986      ENDIF 
     987      trn_tavg(:,:,:,:) = 0.0 
     988 
     989#if defined key_hadocc 
     990      ALLOCATE( HADOCC_CHL_tavg(jpi,jpj,jpk), STAT=ierror ) 
     991      IF( ierror > 0 ) THEN 
     992         CALL ctl_stop( 'asm_bgc_bkg_alloc: unable to allocate HADOCC_CHL_tavg' ) 
     993      ENDIF 
     994      HADOCC_CHL_tavg(:,:,:) = 0.0 
     995 
     996      ALLOCATE( cchl_p_tavg(jpi,jpj,jpk), STAT=ierror ) 
     997      IF( ierror > 0 ) THEN 
     998         CALL ctl_stop( 'asm_bgc_bkg_alloc: unable to allocate cchl_p_tavg' ) 
     999      ENDIF 
     1000      cchl_p_tavg(:,:,:) = 0.0 
     1001 
     1002#elif defined key_fabm 
     1003      ALLOCATE( totalk_tavg(jpi,jpj,jpk), STAT=ierror ) 
     1004      IF( ierror > 0 ) THEN 
     1005         CALL ctl_stop( 'asm_bgc_bkg_alloc: unable to allocate totalk_tavg' ) 
     1006      ENDIF 
     1007      totalk_tavg(:,:,:) = 0.0 
     1008#endif 
     1009                                  
     1010   END SUBROUTINE asm_bgc_bkg_alloc 
     1011 
     1012   !!=========================================================================== 
     1013   !!=========================================================================== 
     1014   !!=========================================================================== 
     1015 
     1016   SUBROUTINE asm_bgc_bkg_tavg(kt, pnumtimes_tavg) 
    7771017      !!------------------------------------------------------------------------ 
    7781018      !!                    ***  ROUTINE asm_bgc_bkg_wri  *** 
     
    7871027      !!------------------------------------------------------------------------ 
    7881028      !! 
    789       INTEGER, INTENT(in   ) :: kt     ! Current time-step 
    790       INTEGER, INTENT(in   ) :: knum   ! i/o unit of increments file 
    791       !! 
    792       !!------------------------------------------------------------------------ 
    793  
     1029      INTEGER,  INTENT(in   ) :: kt               ! Current time-step 
     1030      !! 
     1031      REAL(wp), INTENT(in   ) :: pnumtimes_tavg   ! No of times to average over 
     1032      !! 
     1033      !!------------------------------------------------------------------------ 
     1034       
     1035      IF (kt == nittrc000) THEN 
     1036         pgrow_avg_tavg(:,:) = 0.0 
     1037         ploss_avg_tavg(:,:) = 0.0 
     1038         phyt_avg_tavg(:,:)  = 0.0 
     1039         trn_tavg(:,:,:,:)   = 0.0 
    7941040#if defined key_hadocc 
    795       CALL iom_rstput( kt, nitbkg_r, knum, 'pgrow_avg'   , pgrow_avg             ) 
    796       CALL iom_rstput( kt, nitbkg_r, knum, 'ploss_avg'   , ploss_avg             ) 
    797       CALL iom_rstput( kt, nitbkg_r, knum, 'phyt_avg'    , phyt_avg              ) 
    798       CALL iom_rstput( kt, nitbkg_r, knum, 'mld_max'     , mld_max               ) 
    799       CALL iom_rstput( kt, nitbkg_r, knum, 'hadocc_nut'  , trn(:,:,:,jp_had_nut) ) 
    800       CALL iom_rstput( kt, nitbkg_r, knum, 'hadocc_phy'  , trn(:,:,:,jp_had_phy) ) 
    801       CALL iom_rstput( kt, nitbkg_r, knum, 'hadocc_zoo'  , trn(:,:,:,jp_had_zoo) ) 
    802       CALL iom_rstput( kt, nitbkg_r, knum, 'hadocc_pdn'  , trn(:,:,:,jp_had_pdn) ) 
    803       CALL iom_rstput( kt, nitbkg_r, knum, 'hadocc_dic'  , trn(:,:,:,jp_had_dic) ) 
    804       CALL iom_rstput( kt, nitbkg_r, knum, 'hadocc_alk'  , trn(:,:,:,jp_had_alk) ) 
    805       CALL iom_rstput( kt, nitbkg_r, knum, 'hadocc_chl'  , HADOCC_CHL(:,:,:)     ) 
    806       CALL iom_rstput( kt, nitbkg_r, knum, 'hadocc_cchl' , cchl_p(:,:,:)         ) 
     1041         HADOCC_CHL_tavg(:,:,:) = 0.0 
     1042         cchl_p_tavg(:,:,:)     = 0.0 
     1043#elif defined key_fabm 
     1044         totalk_tavg(:,:,:)  = 0.0 
     1045#endif 
     1046      ENDIF 
     1047       
     1048      pgrow_avg_tavg(:,:) = pgrow_avg_tavg(:,:) + pgrow_avg(:,:) / pnumtimes_tavg 
     1049      ploss_avg_tavg(:,:) = ploss_avg_tavg(:,:) + ploss_avg(:,:) / pnumtimes_tavg 
     1050      phyt_avg_tavg(:,:)  = phyt_avg_tavg(:,:)  + phyt_avg(:,:)  / pnumtimes_tavg 
     1051      trn_tavg(:,:,:,:)   = trn_tavg(:,:,:,:)   + trn(:,:,:,:)   / pnumtimes_tavg 
     1052#if defined key_hadocc 
     1053      HADOCC_CHL_tavg(:,:,:) = HADOCC_CHL_tavg(:,:,:) + HADOCC_CHL(:,:,:) / pnumtimes_tavg 
     1054      cchl_p_tavg(:,:,:)     = cchl_p_tavg(:,:,:)     + cchl_p(:,:,:)     / pnumtimes_tavg 
     1055#elif defined key_fabm 
     1056      totalk_tavg(:,:,:)  = totalk_tavg(:,:,:) + & 
     1057         &                  fabm_get_interior_diagnostic_data(model, jp_fabm_o3ta) / pnumtimes_tavg 
     1058      totalk_tavg(:,:,:)  = totalk_tavg(:,:,:) * tmask(:,:,:) 
     1059#endif 
     1060                                  
     1061   END SUBROUTINE asm_bgc_bkg_tavg 
     1062 
     1063   !!=========================================================================== 
     1064   !!=========================================================================== 
     1065   !!=========================================================================== 
     1066 
     1067   SUBROUTINE asm_bgc_bkg_wri( kt, knum, ld_avgbkg ) 
     1068      !!------------------------------------------------------------------------ 
     1069      !!                    ***  ROUTINE asm_bgc_bkg_wri  *** 
     1070      !! 
     1071      !! ** Purpose :   write out bgc background 
     1072      !! 
     1073      !! ** Method  :   write out bgc background 
     1074      !! 
     1075      !! ** Action  :   write out bgc background 
     1076      !! 
     1077      !! References :   asm_bkg_wri 
     1078      !!------------------------------------------------------------------------ 
     1079      !! 
     1080      INTEGER, INTENT(in   ) :: kt          ! Current time-step 
     1081      INTEGER, INTENT(in   ) :: knum        ! i/o unit of increments file 
     1082      LOGICAL, INTENT(in   ) :: ld_avgbkg   ! Averaged background? 
     1083      !! 
     1084      INTEGER                :: nitbgcbkg_r ! Period referenced to nit000 
     1085      !! 
     1086      !!------------------------------------------------------------------------ 
     1087       
     1088      IF (ld_avgbkg) THEN 
     1089         nitbgcbkg_r = nitavgbkg_r 
     1090      ELSE 
     1091         nitbgcbkg_r = nitbkg_r 
     1092         pgrow_avg_tavg(:,:) = pgrow_avg(:,:) 
     1093         ploss_avg_tavg(:,:) = ploss_avg(:,:) 
     1094         phyt_avg_tavg(:,:)  = phyt_avg(:,:) 
     1095         trn_tavg(:,:,:,:)   = trn(:,:,:,:) 
     1096#if defined key_hadocc 
     1097         HADOCC_CHL_tavg(:,:,:) = HADOCC_CHL(:,:,:) 
     1098         cchl_p_tavg(:,:,:)     = cchl_p(:,:,:) 
     1099#elif defined key_fabm 
     1100         totalk_tavg(:,:,:)  = fabm_get_interior_diagnostic_data(model, jp_fabm_o3ta) 
     1101         totalk_tavg(:,:,:)  = totalk_tavg(:,:,:) * tmask(:,:,:) 
     1102#endif 
     1103      ENDIF 
     1104 
     1105#if defined key_hadocc 
     1106      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'pgrow_avg'   , pgrow_avg_tavg             ) 
     1107      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ploss_avg'   , ploss_avg_tavg             ) 
     1108      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'phyt_avg'    , phyt_avg_tavg              ) 
     1109      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'mld_max'     , mld_max                    ) 
     1110      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'hadocc_nut'  , trn_tavg(:,:,:,jp_had_nut) ) 
     1111      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'hadocc_phy'  , trn_tavg(:,:,:,jp_had_phy) ) 
     1112      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'hadocc_zoo'  , trn_tavg(:,:,:,jp_had_zoo) ) 
     1113      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'hadocc_pdn'  , trn_tavg(:,:,:,jp_had_pdn) ) 
     1114      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'hadocc_dic'  , trn_tavg(:,:,:,jp_had_dic) ) 
     1115      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'hadocc_alk'  , trn_tavg(:,:,:,jp_had_alk) ) 
     1116      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'hadocc_chl'  , HADOCC_CHL(:,:,:)     ) 
     1117      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'hadocc_cchl' , cchl_p(:,:,:)         ) 
    8071118#elif defined key_medusa 
    808       CALL iom_rstput( kt, nitbkg_r, knum, 'pgrow_avg'   , pgrow_avg        ) 
    809       CALL iom_rstput( kt, nitbkg_r, knum, 'ploss_avg'   , ploss_avg        ) 
    810       CALL iom_rstput( kt, nitbkg_r, knum, 'phyt_avg'    , phyt_avg         ) 
    811       CALL iom_rstput( kt, nitbkg_r, knum, 'mld_max'     , mld_max          ) 
    812       CALL iom_rstput( kt, nitbkg_r, knum, 'medusa_chn'  , trn(:,:,:,jpchn) ) 
    813       CALL iom_rstput( kt, nitbkg_r, knum, 'medusa_chd'  , trn(:,:,:,jpchd) ) 
    814       CALL iom_rstput( kt, nitbkg_r, knum, 'medusa_phn'  , trn(:,:,:,jpphn) ) 
    815       CALL iom_rstput( kt, nitbkg_r, knum, 'medusa_phd'  , trn(:,:,:,jpphd) ) 
    816       CALL iom_rstput( kt, nitbkg_r, knum, 'medusa_pds'  , trn(:,:,:,jppds) ) 
    817       CALL iom_rstput( kt, nitbkg_r, knum, 'medusa_zmi'  , trn(:,:,:,jpzmi) ) 
    818       CALL iom_rstput( kt, nitbkg_r, knum, 'medusa_zme'  , trn(:,:,:,jpzme) ) 
    819       CALL iom_rstput( kt, nitbkg_r, knum, 'medusa_din'  , trn(:,:,:,jpdin) ) 
    820       CALL iom_rstput( kt, nitbkg_r, knum, 'medusa_sil'  , trn(:,:,:,jpsil) ) 
    821       CALL iom_rstput( kt, nitbkg_r, knum, 'medusa_fer'  , trn(:,:,:,jpfer) ) 
    822       CALL iom_rstput( kt, nitbkg_r, knum, 'medusa_det'  , trn(:,:,:,jpdet) ) 
    823       CALL iom_rstput( kt, nitbkg_r, knum, 'medusa_dtc'  , trn(:,:,:,jpdtc) ) 
    824       CALL iom_rstput( kt, nitbkg_r, knum, 'medusa_dic'  , trn(:,:,:,jpdic) ) 
    825       CALL iom_rstput( kt, nitbkg_r, knum, 'medusa_alk'  , trn(:,:,:,jpalk) ) 
    826       CALL iom_rstput( kt, nitbkg_r, knum, 'medusa_oxy'  , trn(:,:,:,jpoxy) ) 
     1119      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'pgrow_avg'   , pgrow_avg_tavg        ) 
     1120      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ploss_avg'   , ploss_avg_tavg        ) 
     1121      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'phyt_avg'    , phyt_avg_tavg         ) 
     1122      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'mld_max'     , mld_max               ) 
     1123      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'medusa_chn'  , trn_tavg(:,:,:,jpchn) ) 
     1124      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'medusa_chd'  , trn_tavg(:,:,:,jpchd) ) 
     1125      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'medusa_phn'  , trn_tavg(:,:,:,jpphn) ) 
     1126      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'medusa_phd'  , trn_tavg(:,:,:,jpphd) ) 
     1127      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'medusa_pds'  , trn_tavg(:,:,:,jppds) ) 
     1128      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'medusa_zmi'  , trn_tavg(:,:,:,jpzmi) ) 
     1129      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'medusa_zme'  , trn_tavg(:,:,:,jpzme) ) 
     1130      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'medusa_din'  , trn_tavg(:,:,:,jpdin) ) 
     1131      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'medusa_sil'  , trn_tavg(:,:,:,jpsil) ) 
     1132      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'medusa_fer'  , trn_tavg(:,:,:,jpfer) ) 
     1133      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'medusa_det'  , trn_tavg(:,:,:,jpdet) ) 
     1134      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'medusa_dtc'  , trn_tavg(:,:,:,jpdtc) ) 
     1135      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'medusa_dic'  , trn_tavg(:,:,:,jpdic) ) 
     1136      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'medusa_alk'  , trn_tavg(:,:,:,jpalk) ) 
     1137      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'medusa_oxy'  , trn_tavg(:,:,:,jpoxy) ) 
     1138#elif defined key_fabm 
     1139      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'pgrow_avg'   , pgrow_avg_tavg               ) 
     1140      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ploss_avg'   , ploss_avg_tavg               ) 
     1141      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'phyt_avg'    , phyt_avg_tavg                ) 
     1142      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'mld_max'     , mld_max                      ) 
     1143      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ersem_chl1'  , trn_tavg(:,:,:,jp_fabm_chl1) ) 
     1144      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ersem_chl2'  , trn_tavg(:,:,:,jp_fabm_chl2) ) 
     1145      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ersem_chl3'  , trn_tavg(:,:,:,jp_fabm_chl3) ) 
     1146      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ersem_chl4'  , trn_tavg(:,:,:,jp_fabm_chl4) ) 
     1147      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ersem_p1c'   , trn_tavg(:,:,:,jp_fabm_p1c)  ) 
     1148      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ersem_p1n'   , trn_tavg(:,:,:,jp_fabm_p1n)  ) 
     1149      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ersem_p1p'   , trn_tavg(:,:,:,jp_fabm_p1p)  ) 
     1150      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ersem_p1s'   , trn_tavg(:,:,:,jp_fabm_p1s)  ) 
     1151      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ersem_p2c'   , trn_tavg(:,:,:,jp_fabm_p2c)  ) 
     1152      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ersem_p2n'   , trn_tavg(:,:,:,jp_fabm_p2n)  ) 
     1153      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ersem_p2p'   , trn_tavg(:,:,:,jp_fabm_p2p)  ) 
     1154      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ersem_p3c'   , trn_tavg(:,:,:,jp_fabm_p3c)  ) 
     1155      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ersem_p3n'   , trn_tavg(:,:,:,jp_fabm_p3n)  ) 
     1156      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ersem_p3p'   , trn_tavg(:,:,:,jp_fabm_p3p)  ) 
     1157      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ersem_p4c'   , trn_tavg(:,:,:,jp_fabm_p4c)  ) 
     1158      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ersem_p4n'   , trn_tavg(:,:,:,jp_fabm_p4n)  ) 
     1159      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ersem_p4p'   , trn_tavg(:,:,:,jp_fabm_p4p)  ) 
     1160      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ersem_z4c'   , trn_tavg(:,:,:,jp_fabm_z4c)  ) 
     1161      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ersem_z5c'   , trn_tavg(:,:,:,jp_fabm_z5c)  ) 
     1162      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ersem_z5n'   , trn_tavg(:,:,:,jp_fabm_z5n)  ) 
     1163      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ersem_z5p'   , trn_tavg(:,:,:,jp_fabm_z5p)  ) 
     1164      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ersem_z6c'   , trn_tavg(:,:,:,jp_fabm_z6c)  ) 
     1165      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ersem_z6n'   , trn_tavg(:,:,:,jp_fabm_z6n)  ) 
     1166      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ersem_z6p'   , trn_tavg(:,:,:,jp_fabm_z6p)  ) 
     1167      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ersem_n1p'   , trn_tavg(:,:,:,jp_fabm_n1p)  ) 
     1168      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ersem_n3n'   , trn_tavg(:,:,:,jp_fabm_n3n)  ) 
     1169      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ersem_n4n'   , trn_tavg(:,:,:,jp_fabm_n4n)  ) 
     1170      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ersem_n5s'   , trn_tavg(:,:,:,jp_fabm_n5s)  ) 
     1171      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ersem_o2o'   , trn_tavg(:,:,:,jp_fabm_o2o)  ) 
     1172      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ersem_o3c'   , trn_tavg(:,:,:,jp_fabm_o3c)  ) 
     1173      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ersem_o3ba'  , trn_tavg(:,:,:,jp_fabm_o3ba) ) 
     1174      CALL iom_rstput( kt, nitbgcbkg_r, knum, 'ersem_o3ta'  , totalk_tavg                  ) 
    8271175#endif 
    8281176 
     
    9111259      REAL(wp), DIMENSION(jpi,jpj) :: zinc_phynon ! Local phynon incs 
    9121260      REAL(wp), DIMENSION(jpi,jpj) :: zbkg_phynon ! Local phynon bkg 
     1261#elif defined key_fabm 
     1262      REAL(wp), DIMENSION(jpi,jpj) :: zinc_chldia ! Local chldia incs 
     1263      REAL(wp), DIMENSION(jpi,jpj) :: zbkg_chldia ! Local chldia bkg 
     1264      REAL(wp), DIMENSION(jpi,jpj) :: zinc_chldin ! Local chldin incs 
     1265      REAL(wp), DIMENSION(jpi,jpj) :: zbkg_chldin ! Local chldin bkg 
     1266      REAL(wp), DIMENSION(jpi,jpj) :: zinc_chlnan ! Local chlnan incs 
     1267      REAL(wp), DIMENSION(jpi,jpj) :: zbkg_chlnan ! Local chlnan bkg 
     1268      REAL(wp), DIMENSION(jpi,jpj) :: zinc_chlpic ! Local chlpic incs 
     1269      REAL(wp), DIMENSION(jpi,jpj) :: zbkg_chlpic ! Local chlpic bkg 
    9131270#endif 
    9141271      !!------------------------------------------------------------------------ 
     
    9261283#elif defined key_hadocc 
    9271284            zbkg_chltot(:,:) = chl_bkg(:,:,1) 
     1285#elif defined key_fabm 
     1286            zbkg_chltot(:,:) = tracer_bkg(:,:,1,jp_fabm_chl1) + & 
     1287               &               tracer_bkg(:,:,1,jp_fabm_chl2) + & 
     1288               &               tracer_bkg(:,:,1,jp_fabm_chl3) + & 
     1289               &               tracer_bkg(:,:,1,jp_fabm_chl4) 
    9281290#endif 
    9291291            CALL asm_bgc_unlog_2d( zbkg_chltot, slchltot_bkginc, zinc_chltot ) 
     
    9341296         ENDIF 
    9351297 
    936 #if defined key_medusa 
     1298#if defined key_medusa || defined key_fabm 
    9371299         ! Diatom chlorophyll 
    9381300         IF ( ln_slchldiainc ) THEN 
     1301#if defined key_medusa 
    9391302            zbkg_chldia(:,:) = tracer_bkg(:,:,1,jpchd) 
     1303#elif defined key_fabm 
     1304            zbkg_chldia(:,:) = tracer_bkg(:,:,1,jp_fabm_chl1) 
     1305#endif 
    9401306            CALL asm_bgc_unlog_2d( zbkg_chldia, slchldia_bkginc, zinc_chldia ) 
    9411307         ELSE 
     
    9541320#endif 
    9551321 
     1322#if defined key_fabm 
     1323         ! Nanophytoplankton chlorophyll 
     1324         IF ( ln_slchlnaninc ) THEN 
     1325            zbkg_chlnan(:,:) = tracer_bkg(:,:,1,jp_fabm_chl2) 
     1326            CALL asm_bgc_unlog_2d( zbkg_chlnan, slchlnan_bkginc, zinc_chlnan ) 
     1327         ELSE 
     1328            zinc_chlnan(:,:) = 0.0 
     1329         ENDIF 
     1330 
     1331         ! Picophytoplankton chlorophyll 
     1332         IF ( ln_slchlpicinc ) THEN 
     1333            zbkg_chlpic(:,:) = tracer_bkg(:,:,1,jp_fabm_chl3) 
     1334            CALL asm_bgc_unlog_2d( zbkg_chlpic, slchlpic_bkginc, zinc_chlpic ) 
     1335         ELSE 
     1336            zinc_chlpic(:,:) = 0.0 
     1337         ENDIF 
     1338 
     1339         ! Dinoflagellate chlorophyll 
     1340         IF ( ln_slchldininc ) THEN 
     1341            zbkg_chldin(:,:) = tracer_bkg(:,:,1,jp_fabm_chl4) 
     1342            CALL asm_bgc_unlog_2d( zbkg_chldin, slchldin_bkginc, zinc_chldin ) 
     1343         ELSE 
     1344            zinc_chldin(:,:) = 0.0 
     1345         ENDIF 
     1346#endif 
     1347 
    9561348         ! Total phytoplankton carbon 
    9571349         IF ( ln_slphytotinc ) THEN 
     
    9881380         ! Select mixed layer 
    9891381         IF ( ll_asmdin ) THEN 
    990 #if defined key_top && ( defined key_hadocc || defined key_medusa ) 
     1382#if defined key_top && ( defined key_hadocc || defined key_medusa || defined key_fabm ) 
    9911383            CALL ctl_warn( ' Doing direct initialisation with ocean colour assimilation', & 
    9921384               &           ' Mixed layer depth taken to be background maximum mld_max_bkg' ) 
     
    10141406#endif 
    10151407            CASE ( 4 )                   ! Temperature criterion (0.2 K change from surface) [T points] 
    1016                !zmld(:,:) = hmld_tref(:,:) 
    1017                CALL ctl_stop( ' hmld_tref mixed layer requested for phyto2d assimilation,', & 
    1018                   &           ' but is not available in this version' ) 
     1408               zmld(:,:) = hmld_tref(:,:) 
    10191409            CASE ( 5 )                   ! Density criterion (0.01 kg/m^3 change from 10m) [T points] 
    10201410               zmld(:,:) = hmlpt(:,:) 
     
    10531443            &                         cchl_p_bkg(:,:,1),                   & 
    10541444            &                         tracer_bkg, phyto2d_balinc ) 
     1445#elif defined key_fabm 
     1446         CALL asm_phyto2d_bal_ersem(  (ln_slchltotinc .OR. ln_schltotinc), & 
     1447            &                         zinc_chltot,                         & 
     1448            &                         ln_slchldiainc,                      & 
     1449            &                         zinc_chldia,                         & 
     1450            &                         ln_slchlnaninc,                      & 
     1451            &                         zinc_chlnan,                         & 
     1452            &                         ln_slchlpicinc,                      & 
     1453            &                         zinc_chlpic,                         & 
     1454            &                         ln_slchldininc,                      & 
     1455            &                         zinc_chldin,                         & 
     1456            &                         zincper,                             & 
     1457            &                         rn_maxchlinc, ln_phytobal, zmld,     & 
     1458            &                         pgrow_avg_bkg, ploss_avg_bkg,        & 
     1459            &                         phyt_avg_bkg, mld_max_bkg,           & 
     1460            &                         totalk_bkg,                          & 
     1461            &                         tracer_bkg, phyto2d_balinc ) 
    10551462#else 
    10561463         CALL ctl_stop( 'Attempting to assimilate phyto2d data, ', & 
     
    10981505                  &                         phyto2d_balinc(:,:,:,jp_had0:jp_had1) * zincwgt 
    10991506            END WHERE 
     1507#elif defined key_fabm 
     1508            WHERE( phyto2d_balinc(:,:,:,jp_fabm0:jp_fabm1) > 0.0_wp .OR. & 
     1509               &   trn(:,:,:,jp_fabm0:jp_fabm1) + phyto2d_balinc(:,:,:,jp_fabm0:jp_fabm1) * zincwgt > 0.0_wp ) 
     1510               trn(:,:,:,jp_fabm0:jp_fabm1) = trn(:,:,:,jp_fabm0:jp_fabm1) + & 
     1511                  &                           phyto2d_balinc(:,:,:,jp_fabm0:jp_fabm1) * zincwgt 
     1512               trb(:,:,:,jp_fabm0:jp_fabm1) = trb(:,:,:,jp_fabm0:jp_fabm1) + & 
     1513                  &                           phyto2d_balinc(:,:,:,jp_fabm0:jp_fabm1) * zincwgt 
     1514            END WHERE 
    11001515#endif 
    11011516 
     
    11321547               trb(:,:,:,jp_had0:jp_had1) = trn(:,:,:,jp_had0:jp_had1) 
    11331548            END WHERE 
     1549#elif defined key_fabm 
     1550            WHERE( phyto2d_balinc(:,:,:,jp_fabm0:jp_fabm1) > 0.0_wp .OR. & 
     1551               &   trn(:,:,:,jp_fabm0:jp_fabm1) + phyto2d_balinc(:,:,:,jp_fabm0:jp_fabm1) > 0.0_wp ) 
     1552               trn(:,:,:,jp_fabm0:jp_fabm1) = trn(:,:,:,jp_fabm0:jp_fabm1) + & 
     1553                  &                           phyto2d_balinc(:,:,:,jp_fabm0:jp_fabm1) 
     1554               trb(:,:,:,jp_fabm0:jp_fabm1) = trn(:,:,:,jp_fabm0:jp_fabm1) 
     1555            END WHERE 
    11341556#endif 
    11351557  
     
    11631585      REAL(wp), DIMENSION(kcycper), INTENT(IN) :: pwgtiau   ! IAU weights 
    11641586      ! 
    1165       INTEGER                          :: ji, jj, jk   ! Loop counters 
    1166       INTEGER                          :: it           ! Index 
    1167       REAL(wp)                         :: zincwgt      ! IAU weight for timestep 
    1168       REAL(wp)                         :: zfrac_chn    ! Fraction of jpchn 
    1169       REAL(wp)                         :: zfrac_chd    ! Fraction of jpchd 
    1170       REAL(wp)                         :: zrat_phn_chn ! jpphn:jpchn ratio 
    1171       REAL(wp)                         :: zrat_phd_chd ! jpphd:jpchd ratio 
    1172       REAL(wp)                         :: zrat_pds_chd ! jppds:jpchd ratio 
    1173       REAL(wp), DIMENSION(jpi,jpj,jpk) :: chl_inc      ! Chlorophyll increments 
    1174       REAL(wp), DIMENSION(jpi,jpj,jpk) :: bkg_chl      ! Chlorophyll background 
     1587      INTEGER                          :: ji, jj, jk    ! Loop counters 
     1588      INTEGER                          :: it            ! Index 
     1589      REAL(wp)                         :: zincwgt       ! IAU weight for timestep 
     1590      REAL(wp)                         :: zfrac_chn     ! Fraction of jpchn 
     1591      REAL(wp)                         :: zfrac_chd     ! Fraction of jpchd 
     1592      REAL(wp)                         :: zfrac_chl1    ! Fraction of jp_fabm_chl1 
     1593      REAL(wp)                         :: zfrac_chl2    ! Fraction of jp_fabm_chl2 
     1594      REAL(wp)                         :: zfrac_chl3    ! Fraction of jp_fabm_chl3 
     1595      REAL(wp)                         :: zfrac_chl4    ! Fraction of jp_fabm_chl4 
     1596      REAL(wp)                         :: zrat_phn_chn  ! jpphn:jpchn ratio 
     1597      REAL(wp)                         :: zrat_phd_chd  ! jpphd:jpchd ratio 
     1598      REAL(wp)                         :: zrat_pds_chd  ! jppds:jpchd ratio 
     1599      REAL(wp)                         :: zrat_p1c_chl1 ! jp_fabm_p1c:jp_fabm_chl1 ratio 
     1600      REAL(wp)                         :: zrat_p1n_chl1 ! jp_fabm_p1n:jp_fabm_chl1 ratio 
     1601      REAL(wp)                         :: zrat_p1p_chl1 ! jp_fabm_p1p:jp_fabm_chl1 ratio 
     1602      REAL(wp)                         :: zrat_p1s_chl1 ! jp_fabm_p1s:jp_fabm_chl1 ratio 
     1603      REAL(wp)                         :: zrat_p2c_chl2 ! jp_fabm_p2c:jp_fabm_chl2 ratio 
     1604      REAL(wp)                         :: zrat_p2n_chl2 ! jp_fabm_p2n:jp_fabm_chl2 ratio 
     1605      REAL(wp)                         :: zrat_p2p_chl2 ! jp_fabm_p2p:jp_fabm_chl2 ratio 
     1606      REAL(wp)                         :: zrat_p3c_chl3 ! jp_fabm_p3c:jp_fabm_chl3 ratio 
     1607      REAL(wp)                         :: zrat_p3n_chl3 ! jp_fabm_p3n:jp_fabm_chl3 ratio 
     1608      REAL(wp)                         :: zrat_p3p_chl3 ! jp_fabm_p3p:jp_fabm_chl3 ratio 
     1609      REAL(wp)                         :: zrat_p4c_chl4 ! jp_fabm_p4c:jp_fabm_chl4 ratio 
     1610      REAL(wp)                         :: zrat_p4n_chl4 ! jp_fabm_p4n:jp_fabm_chl4 ratio 
     1611      REAL(wp)                         :: zrat_p4p_chl4 ! jp_fabm_p4p:jp_fabm_chl4 ratio 
     1612      REAL(wp), DIMENSION(jpi,jpj,jpk) :: chl_inc       ! Chlorophyll increments 
     1613      REAL(wp), DIMENSION(jpi,jpj,jpk) :: bkg_chl       ! Chlorophyll background 
    11751614      !!------------------------------------------------------------------------ 
    11761615 
     
    11891628#elif defined key_hadocc 
    11901629            bkg_chl(:,:,:) = chl_bkg(:,:,:) 
     1630#elif defined key_fabm 
     1631            bkg_chl(:,:,:) = tracer_bkg(:,:,:,jp_fabm_chl1) + & 
     1632               &             tracer_bkg(:,:,:,jp_fabm_chl2) + & 
     1633               &             tracer_bkg(:,:,:,jp_fabm_chl3) + & 
     1634               &             tracer_bkg(:,:,:,jp_fabm_chl4) 
    11911635#endif 
    11921636            DO jk = 1, jpk 
     
    12401684#elif defined key_hadocc 
    12411685         phyto3d_balinc(:,:,:,jp_had_phy) = ( cchl_p_bkg(:,:,:) / (mw_carbon * c2n_p) ) * chl_inc(:,:,:) 
     1686#elif defined key_fabm 
     1687         ! Loop over each grid point partioning the increments based on existing ratios 
     1688         DO jk = 1, jpk 
     1689            DO jj = 1, jpj 
     1690               DO ji = 1, jpi 
     1691                  IF ( ( tracer_bkg(ji,jj,jk,jp_fabm_chl1) > 0.0 ) .AND. & 
     1692                     & ( tracer_bkg(ji,jj,jk,jp_fabm_chl2) > 0.0 ) .AND. & 
     1693                     & ( tracer_bkg(ji,jj,jk,jp_fabm_chl3) > 0.0 ) .AND. & 
     1694                     & ( tracer_bkg(ji,jj,jk,jp_fabm_chl4) > 0.0 ) ) THEN 
     1695                     zfrac_chl1 = tracer_bkg(ji,jj,jk,jp_fabm_chl1) / bkg_chl(ji,jj,jk) 
     1696                     zfrac_chl2 = tracer_bkg(ji,jj,jk,jp_fabm_chl2) / bkg_chl(ji,jj,jk) 
     1697                     zfrac_chl3 = tracer_bkg(ji,jj,jk,jp_fabm_chl3) / bkg_chl(ji,jj,jk) 
     1698                     zfrac_chl4 = tracer_bkg(ji,jj,jk,jp_fabm_chl4) / bkg_chl(ji,jj,jk) 
     1699                     phyto3d_balinc(ji,jj,jk,jp_fabm_chl1) = chl_inc(ji,jj,jk) * zfrac_chl1 
     1700                     phyto3d_balinc(ji,jj,jk,jp_fabm_chl2) = chl_inc(ji,jj,jk) * zfrac_chl2 
     1701                     phyto3d_balinc(ji,jj,jk,jp_fabm_chl3) = chl_inc(ji,jj,jk) * zfrac_chl3 
     1702                     phyto3d_balinc(ji,jj,jk,jp_fabm_chl4) = chl_inc(ji,jj,jk) * zfrac_chl4 
     1703                     zrat_p1c_chl1 = tracer_bkg(ji,jj,jk,jp_fabm_p1c) / tracer_bkg(ji,jj,jk,jp_fabm_chl1) 
     1704                     zrat_p1n_chl1 = tracer_bkg(ji,jj,jk,jp_fabm_p1n) / tracer_bkg(ji,jj,jk,jp_fabm_chl1) 
     1705                     zrat_p1p_chl1 = tracer_bkg(ji,jj,jk,jp_fabm_p1p) / tracer_bkg(ji,jj,jk,jp_fabm_chl1) 
     1706                     zrat_p1s_chl1 = tracer_bkg(ji,jj,jk,jp_fabm_p1s) / tracer_bkg(ji,jj,jk,jp_fabm_chl1) 
     1707                     zrat_p2c_chl2 = tracer_bkg(ji,jj,jk,jp_fabm_p2c) / tracer_bkg(ji,jj,jk,jp_fabm_chl2) 
     1708                     zrat_p2n_chl2 = tracer_bkg(ji,jj,jk,jp_fabm_p2n) / tracer_bkg(ji,jj,jk,jp_fabm_chl2) 
     1709                     zrat_p2p_chl2 = tracer_bkg(ji,jj,jk,jp_fabm_p2p) / tracer_bkg(ji,jj,jk,jp_fabm_chl2) 
     1710                     zrat_p3c_chl3 = tracer_bkg(ji,jj,jk,jp_fabm_p3c) / tracer_bkg(ji,jj,jk,jp_fabm_chl3) 
     1711                     zrat_p3n_chl3 = tracer_bkg(ji,jj,jk,jp_fabm_p3n) / tracer_bkg(ji,jj,jk,jp_fabm_chl3) 
     1712                     zrat_p3p_chl3 = tracer_bkg(ji,jj,jk,jp_fabm_p3p) / tracer_bkg(ji,jj,jk,jp_fabm_chl3) 
     1713                     zrat_p4c_chl4 = tracer_bkg(ji,jj,jk,jp_fabm_p4c) / tracer_bkg(ji,jj,jk,jp_fabm_chl4) 
     1714                     zrat_p4n_chl4 = tracer_bkg(ji,jj,jk,jp_fabm_p4n) / tracer_bkg(ji,jj,jk,jp_fabm_chl4) 
     1715                     zrat_p4p_chl4 = tracer_bkg(ji,jj,jk,jp_fabm_p4p) / tracer_bkg(ji,jj,jk,jp_fabm_chl4) 
     1716                     phyto3d_balinc(ji,jj,jk,jp_fabm_p1c) = phyto3d_balinc(ji,jj,jk,jp_fabm_chl1) * zrat_p1c_chl1 
     1717                     phyto3d_balinc(ji,jj,jk,jp_fabm_p1n) = phyto3d_balinc(ji,jj,jk,jp_fabm_chl1) * zrat_p1n_chl1 
     1718                     phyto3d_balinc(ji,jj,jk,jp_fabm_p1p) = phyto3d_balinc(ji,jj,jk,jp_fabm_chl1) * zrat_p1p_chl1 
     1719                     phyto3d_balinc(ji,jj,jk,jp_fabm_p1s) = phyto3d_balinc(ji,jj,jk,jp_fabm_chl1) * zrat_p1s_chl1 
     1720                     phyto3d_balinc(ji,jj,jk,jp_fabm_p2c) = phyto3d_balinc(ji,jj,jk,jp_fabm_chl2) * zrat_p2c_chl2 
     1721                     phyto3d_balinc(ji,jj,jk,jp_fabm_p2n) = phyto3d_balinc(ji,jj,jk,jp_fabm_chl2) * zrat_p2n_chl2 
     1722                     phyto3d_balinc(ji,jj,jk,jp_fabm_p2p) = phyto3d_balinc(ji,jj,jk,jp_fabm_chl2) * zrat_p2p_chl2 
     1723                     phyto3d_balinc(ji,jj,jk,jp_fabm_p3c) = phyto3d_balinc(ji,jj,jk,jp_fabm_chl3) * zrat_p3c_chl3 
     1724                     phyto3d_balinc(ji,jj,jk,jp_fabm_p3n) = phyto3d_balinc(ji,jj,jk,jp_fabm_chl3) * zrat_p3n_chl3 
     1725                     phyto3d_balinc(ji,jj,jk,jp_fabm_p3p) = phyto3d_balinc(ji,jj,jk,jp_fabm_chl3) * zrat_p3p_chl3 
     1726                     phyto3d_balinc(ji,jj,jk,jp_fabm_p4c) = phyto3d_balinc(ji,jj,jk,jp_fabm_chl4) * zrat_p4c_chl4 
     1727                     phyto3d_balinc(ji,jj,jk,jp_fabm_p4n) = phyto3d_balinc(ji,jj,jk,jp_fabm_chl4) * zrat_p4n_chl4 
     1728                     phyto3d_balinc(ji,jj,jk,jp_fabm_p4p) = phyto3d_balinc(ji,jj,jk,jp_fabm_chl4) * zrat_p4p_chl4 
     1729                  ENDIF 
     1730               END DO 
     1731            END DO 
     1732         END DO 
    12421733#else 
    12431734         CALL ctl_stop( 'Attempting to assimilate p(l)chltot, ', & 
     
    12851776                  &                         phyto3d_balinc(:,:,:,jp_had0:jp_had1) * zincwgt 
    12861777            END WHERE 
     1778#elif defined key_fabm 
     1779            WHERE( phyto3d_balinc(:,:,:,jp_fabm0:jp_fabm1) > 0.0_wp .OR. & 
     1780               &   trn(:,:,:,jp_fabm0:jp_fabm1) + phyto3d_balinc(:,:,:,jp_fabm0:jp_fabm1) * zincwgt > 0.0_wp ) 
     1781               trn(:,:,:,jp_fabm0:jp_fabm1) = trn(:,:,:,jp_fabm0:jp_fabm1) + & 
     1782                  &                           phyto3d_balinc(:,:,:,jp_fabm0:jp_fabm1) * zincwgt 
     1783               trb(:,:,:,jp_fabm0:jp_fabm1) = trb(:,:,:,jp_fabm0:jp_fabm1) + & 
     1784                  &                           phyto3d_balinc(:,:,:,jp_fabm0:jp_fabm1) * zincwgt 
     1785            END WHERE 
    12871786#endif 
    12881787 
     
    13181817                  &                         phyto3d_balinc(:,:,:,jp_had0:jp_had1) 
    13191818               trb(:,:,:,jp_had0:jp_had1) = trn(:,:,:,jp_had0:jp_had1) 
     1819            END WHERE 
     1820#elif defined key_fabm 
     1821            WHERE( phyto3d_balinc(:,:,:,jp_fabm0:jp_fabm1) > 0.0_wp .OR. & 
     1822               &   trn(:,:,:,jp_fabm0:jp_fabm1) + phyto3d_balinc(:,:,:,jp_fabm0:jp_fabm1) > 0.0_wp ) 
     1823               trn(:,:,:,jp_fabm0:jp_fabm1) = trn(:,:,:,jp_fabm0:jp_fabm1) + & 
     1824                  &                           phyto3d_balinc(:,:,:,jp_fabm0:jp_fabm1) 
     1825               trb(:,:,:,jp_fabm0:jp_fabm1) = trn(:,:,:,jp_fabm0:jp_fabm1) 
    13201826            END WHERE 
    13211827#endif 
     
    14471953            &               pco2_balinc(:,:,1,jp_had_dic), pco2_balinc(:,:,1,jp_had_alk) ) 
    14481954 
     1955#elif defined key_fabm 
     1956         ! Account for phytoplankton balancing if required 
     1957         IF ( ln_phytobal ) THEN 
     1958            dic_bkg_temp(:,:) = tracer_bkg(:,:,1,jp_fabm_o3c) + phyto2d_balinc(:,:,1,jp_fabm_o3c) 
     1959            alk_bkg_temp(:,:) = totalk_bkg(:,:,1)             + phyto2d_balinc(:,:,1,jp_fabm_o3ba) 
     1960         ELSE 
     1961            dic_bkg_temp(:,:) = tracer_bkg(:,:,1,jp_fabm_o3c) 
     1962            alk_bkg_temp(:,:) = totalk_bkg(:,:,1) 
     1963         ENDIF 
     1964 
     1965         CALL asm_pco2_bal( pco2_bkginc_temp(:,:), dic_bkg_temp(:,:), alk_bkg_temp(:,:), & 
     1966            &               tem_bkg_temp(:,:), sal_bkg_temp(:,:),                        & 
     1967            &               pco2_balinc(:,:,1,jp_fabm_o3c), pco2_balinc(:,:,1,jp_fabm_o3ba) ) 
     1968 
    14491969#else 
    14501970         CALL ctl_stop( 'Attempting to assimilate pCO2/fCO2, ', & 
     
    14541974         ! Select mixed layer 
    14551975         IF ( ll_asmdin ) THEN 
    1456 #if defined key_hadocc || defined key_medusa 
     1976#if defined key_hadocc || defined key_medusa || defined key_fabm 
    14571977            CALL ctl_warn( ' Doing direct initialisation with pCO2 assimilation', & 
    14581978               &           ' Mixed layer depth taken to be background maximum mld_max_bkg' ) 
     
    14802000#endif 
    14812001            CASE ( 4 )                   ! Temperature criterion (0.2 K change from surface) [T points] 
    1482                !zmld(:,:) = hmld_tref(:,:) 
    1483                CALL ctl_stop( ' hmld_tref mixed layer requested for pCO2 assimilation,', & 
    1484                   &           ' but is not available in this version' ) 
     2002               zmld(:,:) = hmld_tref(:,:) 
    14852003            CASE ( 5 )                   ! Density criterion (0.01 kg/m^3 change from 10m) [T points] 
    14862004               zmld(:,:) = hmlpt(:,:) 
     
    15552073                  &                         pco2_balinc(:,:,:,jp_had0:jp_had1) * zincwgt 
    15562074            END WHERE 
     2075#elif defined key_fabm 
     2076            WHERE( pco2_balinc(:,:,:,jp_fabm0:jp_fabm1) > 0.0_wp .OR. & 
     2077               &   trn(:,:,:,jp_fabm0:jp_fabm1) + pco2_balinc(:,:,:,jp_fabm0:jp_fabm1) * zincwgt > 0.0_wp ) 
     2078               trn(:,:,:,jp_fabm0:jp_fabm1) = trn(:,:,:,jp_fabm0:jp_fabm1) + & 
     2079                  &                           pco2_balinc(:,:,:,jp_fabm0:jp_fabm1) * zincwgt 
     2080               trb(:,:,:,jp_fabm0:jp_fabm1) = trb(:,:,:,jp_fabm0:jp_fabm1) + & 
     2081                  &                           pco2_balinc(:,:,:,jp_fabm0:jp_fabm1) * zincwgt 
     2082            END WHERE 
    15572083#endif 
    15582084 
     
    15892115                  &                         pco2_balinc(:,:,:,jp_had0:jp_had1) 
    15902116               trb(:,:,:,jp_had0:jp_had1) = trn(:,:,:,jp_had0:jp_had1) 
     2117            END WHERE 
     2118#elif defined key_fabm 
     2119            WHERE( pco2_balinc(:,:,:,jp_fabm0:jp_fabm1) > 0.0_wp .OR. & 
     2120               &   trn(:,:,:,jp_fabm0:jp_fabm1) + pco2_balinc(:,:,:,jp_fabm0:jp_fabm1) > 0.0_wp ) 
     2121               trn(:,:,:,jp_fabm0:jp_fabm1) = trn(:,:,:,jp_fabm0:jp_fabm1) + & 
     2122                  &                           pco2_balinc(:,:,:,jp_fabm0:jp_fabm1) 
     2123               trb(:,:,:,jp_fabm0:jp_fabm1) = trn(:,:,:,jp_fabm0:jp_fabm1) 
    15912124            END WHERE 
    15922125#endif 
     
    18542387 
    18552388         IF ( ln_pno3inc ) THEN 
    1856 #if defined key_hadocc || defined key_medusa 
     2389#if defined key_hadocc || defined key_medusa || defined key_fabm 
    18572390#if defined key_hadocc 
    18582391            it = jp_had_nut 
    18592392#elif defined key_medusa 
    18602393            it = jpdin 
     2394#elif defined key_fabm 
     2395            it = jp_fabm_n3n 
    18612396#endif 
    18622397            IF ( ln_phytobal ) THEN 
     
    18782413 
    18792414         IF ( ln_psi4inc ) THEN 
     2415#if defined key_medusa || defined key_fabm 
    18802416#if defined key_medusa 
    18812417            it = jpsil 
     2418#elif defined key_fabm 
     2419            it = jp_fabm_n5s 
     2420#endif 
    18822421            IF ( ln_phytobal ) THEN 
    18832422               psi4_bkginc(:,:,:) = psi4_bkginc(:,:,:) - phyto2d_balinc(:,:,:,it) 
     
    18982437 
    18992438         IF ( ln_pdicinc ) THEN 
    1900 #if defined key_hadocc || defined key_medusa 
     2439#if defined key_hadocc || defined key_medusa || defined key_fabm 
    19012440#if defined key_hadocc 
    19022441            it = jp_had_dic 
    19032442#elif defined key_medusa 
    19042443            it = jpdic 
     2444#elif defined key_fabm 
     2445            it = jp_fabm_o3c 
    19052446#endif 
    19062447            IF ( ln_phytobal ) THEN 
     
    19222463 
    19232464         IF ( ln_palkinc ) THEN 
    1924 #if defined key_hadocc || defined key_medusa 
     2465#if defined key_hadocc || defined key_medusa || defined key_fabm 
    19252466#if defined key_hadocc 
    19262467            it = jp_had_alk 
    19272468#elif defined key_medusa 
    19282469            it = jpalk 
     2470#elif defined key_fabm 
     2471            it = jp_fabm_o3ba 
    19292472#endif 
    19302473            IF ( ln_phytobal ) THEN 
     
    19462489 
    19472490         IF ( ln_po2inc ) THEN 
     2491#if defined key_medusa || defined key_fabm 
    19482492#if defined key_medusa 
    19492493            it = jpoxy 
     2494#elif defined key_fabm 
     2495            it = jp_fabm_o2o 
     2496#endif 
    19502497            IF ( ln_phytobal ) THEN 
    19512498               po2_bkginc(:,:,:) = po2_bkginc(:,:,:) - phyto2d_balinc(:,:,:,it) 
     
    20042551                  trb(:,:,:,jpdin) = trb(:,:,:,jpdin) + pno3_bkginc(:,:,:) * zincwgt 
    20052552               END WHERE 
     2553#elif defined key_fabm 
     2554               WHERE( pno3_bkginc(:,:,:) > 0.0_wp .OR. & 
     2555                  &   trn(:,:,:,jp_fabm_n3n) + pno3_bkginc(:,:,:) * zincwgt > 0.0_wp ) 
     2556                  trn(:,:,:,jp_fabm_n3n) = trn(:,:,:,jp_fabm_n3n) + pno3_bkginc(:,:,:) * zincwgt 
     2557                  trb(:,:,:,jp_fabm_n3n) = trb(:,:,:,jp_fabm_n3n) + pno3_bkginc(:,:,:) * zincwgt 
     2558               END WHERE 
    20062559#else 
    20072560               CALL ctl_stop ( ' bgc3d_asm_inc: no compatible BGC model defined' ) 
     
    20152568                  trn(:,:,:,jpsil) = trn(:,:,:,jpsil) + psi4_bkginc(:,:,:) * zincwgt 
    20162569                  trb(:,:,:,jpsil) = trb(:,:,:,jpsil) + psi4_bkginc(:,:,:) * zincwgt 
     2570               END WHERE 
     2571#elif defined key_fabm 
     2572               WHERE( psi4_bkginc(:,:,:) > 0.0_wp .OR. & 
     2573                  &   trn(:,:,:,jp_fabm_n5s) + psi4_bkginc(:,:,:) * zincwgt > 0.0_wp ) 
     2574                  trn(:,:,:,jp_fabm_n5s) = trn(:,:,:,jp_fabm_n5s) + psi4_bkginc(:,:,:) * zincwgt 
     2575                  trb(:,:,:,jp_fabm_n5s) = trb(:,:,:,jp_fabm_n5s) + psi4_bkginc(:,:,:) * zincwgt 
    20172576               END WHERE 
    20182577#else 
     
    20342593                  trb(:,:,:,jpdic) = trb(:,:,:,jpdic) + pdic_bkginc(:,:,:) * zincwgt 
    20352594               END WHERE 
     2595#elif defined key_fabm 
     2596               WHERE( pdic_bkginc(:,:,:) > 0.0_wp .OR. & 
     2597                  &   trn(:,:,:,jp_fabm_o3c) + pdic_bkginc(:,:,:) * zincwgt > 0.0_wp ) 
     2598                  trn(:,:,:,jp_fabm_o3c) = trn(:,:,:,jp_fabm_o3c) + pdic_bkginc(:,:,:) * zincwgt 
     2599                  trb(:,:,:,jp_fabm_o3c) = trb(:,:,:,jp_fabm_o3c) + pdic_bkginc(:,:,:) * zincwgt 
     2600               END WHERE 
    20362601#else 
    20372602               CALL ctl_stop ( ' bgc3d_asm_inc: no compatible BGC model defined' ) 
     
    20522617                  trb(:,:,:,jpalk) = trb(:,:,:,jpalk) + palk_bkginc(:,:,:) * zincwgt 
    20532618               END WHERE 
     2619#elif defined key_fabm 
     2620               WHERE( palk_bkginc(:,:,:) > 0.0_wp .OR. & 
     2621                  &   trn(:,:,:,jp_fabm_o3ba) + palk_bkginc(:,:,:) * zincwgt > 0.0_wp ) 
     2622                  trn(:,:,:,jp_fabm_o3ba) = trn(:,:,:,jp_fabm_o3ba) + palk_bkginc(:,:,:) * zincwgt 
     2623                  trb(:,:,:,jp_fabm_o3ba) = trb(:,:,:,jp_fabm_o3ba) + palk_bkginc(:,:,:) * zincwgt 
     2624               END WHERE 
    20542625#else 
    20552626               CALL ctl_stop ( ' bgc3d_asm_inc: no compatible BGC model defined' ) 
     
    20632634                  trn(:,:,:,jpoxy) = trn(:,:,:,jpoxy) + po2_bkginc(:,:,:) * zincwgt 
    20642635                  trb(:,:,:,jpoxy) = trb(:,:,:,jpoxy) + po2_bkginc(:,:,:) * zincwgt 
     2636               END WHERE 
     2637#elif defined key_fabm 
     2638               WHERE( po2_bkginc(:,:,:) > 0.0_wp .OR. & 
     2639                  &   trn(:,:,:,jp_fabm_o2o) + po2_bkginc(:,:,:) * zincwgt > 0.0_wp ) 
     2640                  trn(:,:,:,jp_fabm_o2o) = trn(:,:,:,jp_fabm_o2o) + po2_bkginc(:,:,:) * zincwgt 
     2641                  trb(:,:,:,jp_fabm_o2o) = trb(:,:,:,jp_fabm_o2o) + po2_bkginc(:,:,:) * zincwgt 
    20652642               END WHERE 
    20662643#else 
     
    20912668            ! Initialize the now fields with the background + increment 
    20922669            ! Background currently is what the model is initialised with 
    2093 #if defined key_hadocc 
    2094             CALL ctl_warn( ' Doing direct initialisation of HadOCC with 3D BGC assimilation', & 
     2670            CALL ctl_warn( ' Doing direct initialisation with 3D BGC assimilation', & 
    20952671               &           ' Background state is taken from model rather than background file' ) 
    2096 #elif defined key_medusa 
    2097             CALL ctl_warn( ' Doing direct initialisation of MEDUSA with 3D BGC assimilation', & 
    2098                &           ' Background state is taken from model rather than background file' ) 
    2099 #endif 
    21002672 
    21012673            IF ( ln_pno3inc ) THEN 
     
    21122684                  trb(:,:,:,jpdin) = trn(:,:,:,jpdin) 
    21132685               END WHERE 
     2686#elif defined key_fabm 
     2687               WHERE( pno3_bkginc(:,:,:) > 0.0_wp .OR. & 
     2688                  &   trn(:,:,:,jp_fabm_n3n) + pno3_bkginc(:,:,:) > 0.0_wp ) 
     2689                  trn(:,:,:,jp_fabm_n3n) = trn(:,:,:,jp_fabm_n3n) + pno3_bkginc(:,:,:) 
     2690                  trb(:,:,:,jp_fabm_n3n) = trn(:,:,:,jp_fabm_n3n) 
     2691               END WHERE 
    21142692#else 
    21152693               CALL ctl_stop ( ' bgc3d_asm_inc: no compatible BGC model defined' ) 
     
    21232701                  trn(:,:,:,jpsil) = trn(:,:,:,jpsil) + psi4_bkginc(:,:,:) 
    21242702                  trb(:,:,:,jpsil) = trn(:,:,:,jpsil) 
     2703               END WHERE 
     2704#elif defined key_fabm 
     2705               WHERE( psi4_bkginc(:,:,:) > 0.0_wp .OR. & 
     2706                  &   trn(:,:,:,jp_fabm_n5s) + psi4_bkginc(:,:,:) > 0.0_wp ) 
     2707                  trn(:,:,:,jp_fabm_n5s) = trn(:,:,:,jp_fabm_n5s) + psi4_bkginc(:,:,:) 
     2708                  trb(:,:,:,jp_fabm_n5s) = trn(:,:,:,jp_fabm_n5s) 
    21252709               END WHERE 
    21262710#else 
     
    21422726                  trb(:,:,:,jpdic) = trn(:,:,:,jpdic) 
    21432727               END WHERE 
     2728#elif defined key_fabm 
     2729               WHERE( pdic_bkginc(:,:,:) > 0.0_wp .OR. & 
     2730                  &   trn(:,:,:,jp_fabm_o3c) + pdic_bkginc(:,:,:) > 0.0_wp ) 
     2731                  trn(:,:,:,jp_fabm_o3c) = trn(:,:,:,jp_fabm_o3c) + pdic_bkginc(:,:,:) 
     2732                  trb(:,:,:,jp_fabm_o3c) = trn(:,:,:,jp_fabm_o3c) 
     2733               END WHERE 
    21442734#else 
    21452735               CALL ctl_stop ( ' bgc3d_asm_inc: no compatible BGC model defined' ) 
     
    21602750                  trb(:,:,:,jpalk) = trn(:,:,:,jpalk) 
    21612751               END WHERE 
     2752#elif defined key_fabm 
     2753               WHERE( palk_bkginc(:,:,:) > 0.0_wp .OR. & 
     2754                  &   trn(:,:,:,jp_fabm_o3ba) + palk_bkginc(:,:,:) > 0.0_wp ) 
     2755                  trn(:,:,:,jp_fabm_o3ba) = trn(:,:,:,jp_fabm_o3ba) + palk_bkginc(:,:,:) 
     2756                  trb(:,:,:,jp_fabm_o3ba) = trn(:,:,:,jp_fabm_o3ba) 
     2757               END WHERE 
    21622758#else 
    21632759               CALL ctl_stop ( ' bgc3d_asm_inc: no compatible BGC model defined' ) 
     
    21712767                  trn(:,:,:,jpoxy) = trn(:,:,:,jpoxy) + po2_bkginc(:,:,:) 
    21722768                  trb(:,:,:,jpoxy) = trn(:,:,:,jpoxy) 
     2769               END WHERE 
     2770#elif defined key_fabm 
     2771               WHERE( po2_bkginc(:,:,:) > 0.0_wp .OR. & 
     2772                  &   trn(:,:,:,jp_fabm_o2o) + po2_bkginc(:,:,:) > 0.0_wp ) 
     2773                  trn(:,:,:,jp_fabm_o2o) = trn(:,:,:,jp_fabm_o2o) + po2_bkginc(:,:,:) 
     2774                  trb(:,:,:,jp_fabm_o2o) = trn(:,:,:,jp_fabm_o2o) 
    21732775               END WHERE 
    21742776#else 
Note: See TracChangeset for help on using the changeset viewer.