Changeset 9431


Ignore:
Timestamp:
2018-03-26T16:49:51+02:00 (2 years ago)
Author:
dford
Message:

Allow assimilation of non-log chlorophyll, and prepare for assimilation of PFTs.

Location:
branches/UKMO/dev_r5518_GO6_package_asm_3d_bgc/NEMOGCM/NEMO/OPA_SRC/ASM
Files:
3 edited

Legend:

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

    r9382 r9431  
    6161   USE asmpco2bal, ONLY:    & ! pCO2 balancing for MEDUSA 
    6262      & asm_pco2_bal 
     63   USE sms_medusa             ! MEDUSA parameters 
    6364   USE par_medusa             ! MEDUSA parameters 
    6465   USE mocsy_wrapper, ONLY: & ! MEDUSA carbonate system 
     
    819820   !!=========================================================================== 
    820821 
     822   SUBROUTINE asm_bgc_unlog_2d( pbkg, pinc_log, pinc_nonlog ) 
     823      !!------------------------------------------------------------------------ 
     824      !!                    ***  ROUTINE asm_bgc_init_incs  *** 
     825      !! 
     826      !! ** Purpose :   convert log increments to non-log 
     827      !! 
     828      !! ** Method  :   need to account for model background, 
     829      !!                cannot simply do 10^log_inc. Need to: 
     830      !!                1) Add log_inc to log10(background) to get log10(analysis) 
     831      !!                2) Take 10^log10(analysis) to get analysis 
     832      !!                3) Subtract background from analysis to get non-log incs 
     833      !! 
     834      !! ** Action  :   return non-log increments 
     835      !! 
     836      !! References :    
     837      !!------------------------------------------------------------------------ 
     838      !! 
     839      REAL(wp), INTENT(in   ), DIMENSION(jpi,jpj) :: pbkg        ! Background 
     840      REAL(wp), INTENT(in   ), DIMENSION(jpi,jpj) :: pinc_log    ! Log incs 
     841      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) :: pinc_nonlog ! Non-log incs 
     842      ! 
     843      INTEGER                                     :: ji, jj      ! Loop counters 
     844      !! 
     845      !!------------------------------------------------------------------------ 
     846 
     847      DO jj = 1, jpj 
     848         DO ji = 1, jpi 
     849            IF ( pbkg(ji,jj) > 0.0 ) THEN 
     850               pinc_nonlog(ji,jj) = 10**( LOG10( pbkg(ji,jj) ) + & 
     851                  &                       pinc_log(ji,jj) )      & 
     852                  &                 - pbkg(ji,jj) 
     853            ELSE 
     854               pinc_nonlog(ji,jj) = 0.0 
     855            ENDIF 
     856         END DO 
     857      END DO 
     858 
     859   END SUBROUTINE asm_bgc_unlog_2d 
     860 
     861   !!=========================================================================== 
     862   !!=========================================================================== 
     863   !!=========================================================================== 
     864 
    821865   SUBROUTINE phyto2d_asm_inc( kt, ll_asmdin, ll_asmiau, kcycper, pwgtiau ) 
    822866      !!------------------------------------------------------------------------ 
     
    837881      REAL(wp), DIMENSION(kcycper), INTENT(IN) :: pwgtiau   ! IAU weights 
    838882      ! 
    839       INTEGER  :: jk              ! Loop counter 
    840       INTEGER  :: it              ! Index 
    841       REAL(wp) :: zincwgt         ! IAU weight for current time step 
    842       REAL(wp) :: zincper         ! IAU interval in seconds 
    843       !!------------------------------------------------------------------------ 
    844  
    845       IF ( ln_slchldiainc .OR. ln_slchlnoninc .OR. & 
    846          & ln_schltotinc  .OR. ln_slphytotinc .OR. ln_slphydiainc .OR. & 
    847          & ln_slphynoninc ) THEN 
    848          CALL ctl_stop( ' No PFT assimilation quite yet' ) 
    849       ENDIF 
     883      INTEGER                      :: jk          ! Loop counter 
     884      INTEGER                      :: it          ! Index 
     885      REAL(wp)                     :: zincwgt     ! IAU weight for current time step 
     886      REAL(wp)                     :: zincper     ! IAU interval in seconds 
     887      REAL(wp), DIMENSION(jpi,jpj) :: zmld        ! Mixed layer depth 
     888      REAL(wp), DIMENSION(jpi,jpj) :: zinc_chltot ! Local chltot incs 
     889      REAL(wp), DIMENSION(jpi,jpj) :: zbkg_chltot ! Local chltot bkg 
     890      REAL(wp), DIMENSION(jpi,jpj) :: zinc_phytot ! Local phytot incs 
     891      REAL(wp), DIMENSION(jpi,jpj) :: zbkg_phytot ! Local phytot bkg 
     892#if defined key_medusa && defined key_foam_medusa 
     893      REAL(wp), DIMENSION(jpi,jpj) :: zinc_chldia ! Local chldia incs 
     894      REAL(wp), DIMENSION(jpi,jpj) :: zbkg_chldia ! Local chldia bkg 
     895      REAL(wp), DIMENSION(jpi,jpj) :: zinc_chlnon ! Local chlnon incs 
     896      REAL(wp), DIMENSION(jpi,jpj) :: zbkg_chlnon ! Local chlnon bkg 
     897      REAL(wp), DIMENSION(jpi,jpj) :: zinc_phydia ! Local phydia incs 
     898      REAL(wp), DIMENSION(jpi,jpj) :: zbkg_phydia ! Local phydia bkg 
     899      REAL(wp), DIMENSION(jpi,jpj) :: zinc_phynon ! Local phynon incs 
     900      REAL(wp), DIMENSION(jpi,jpj) :: zbkg_phynon ! Local phynon bkg 
     901#endif 
     902      !!------------------------------------------------------------------------ 
    850903       
    851904      IF ( kt <= nit000 ) THEN 
    852905 
     906         ! Un-log any log increments for passing to balancing routines 
     907         ! Total chlorophyll 
     908         IF ( ln_slchltotinc ) THEN 
     909#if defined key_medusa && defined key_foam_medusa 
     910            zbkg_chltot(:,:) = tracer_bkg(:,:,1,jpchn) + tracer_bkg(:,:,1,jpchd) 
     911#elif defined key_hadocc 
     912            zbkg_chltot(:,:) = chl_bkg(:,:,1) 
     913#endif 
     914            CALL asm_bgc_unlog_2d( zbkg_chltot, slchltot_bkginc, zinc_chltot ) 
     915         ELSE IF ( ln_schltotinc ) THEN 
     916            zinc_chltot(:,:) = schltot_bkginc(:,:) 
     917         ELSE 
     918            zinc_chltot(:,:) = 0.0 
     919         ENDIF 
     920 
     921#if defined key_medusa && defined key_foam_medusa 
     922         ! Diatom chlorophyll 
     923         IF ( ln_slchldiainc ) THEN 
     924            zbkg_chldia(:,:) = tracer_bkg(:,:,1,jpchd) 
     925            CALL asm_bgc_unlog_2d( zbkg_chldia, slchldia_bkginc, zinc_chldia ) 
     926         ELSE 
     927            zinc_chldia(:,:) = 0.0 
     928         ENDIF 
     929#endif 
     930 
     931#if defined key_medusa && defined key_foam_medusa 
     932         ! Non-diatom chlorophyll 
     933         IF ( ln_slchlnoninc ) THEN 
     934            zbkg_chlnon(:,:) = tracer_bkg(:,:,1,jpchn) 
     935            CALL asm_bgc_unlog_2d( zbkg_chlnon, slchlnon_bkginc, zinc_chlnon ) 
     936         ELSE 
     937            zinc_chlnon(:,:) = 0.0 
     938         ENDIF 
     939#endif 
     940 
     941         ! Total phytoplankton carbon 
     942         IF ( ln_slphytotinc ) THEN 
     943#if defined key_medusa && defined key_foam_medusa 
     944            zbkg_phytot(:,:) = (trn(:,:,1,jpphn) * xthetapn) + (trn(:,:,1,jpphd) * xthetapd) 
     945#elif defined key_hadocc 
     946            zbkg_phytot(:,:) = trn(:,:,1,jp_had_phy) * c2n_p 
     947#endif 
     948            CALL asm_bgc_unlog_2d( zbkg_phytot, slphytot_bkginc, zinc_phytot ) 
     949         ELSE 
     950            zinc_phytot(:,:) = 0.0 
     951         ENDIF 
     952 
     953#if defined key_medusa && defined key_foam_medusa 
     954         ! Diatom phytoplankton carbon 
     955         IF ( ln_slphydiainc ) THEN 
     956            zbkg_phydia(:,:) = trn(:,:,1,jpphd) * xthetapd 
     957            CALL asm_bgc_unlog_2d( zbkg_phydia, slphydia_bkginc, zinc_phydia ) 
     958         ELSE 
     959            zinc_phydia(:,:) = 0.0 
     960         ENDIF 
     961#endif 
     962 
     963#if defined key_medusa && defined key_foam_medusa 
     964         ! Non-diatom phytoplankton carbon 
     965         IF ( ln_slphynoninc ) THEN 
     966            zbkg_phynon(:,:) = trn(:,:,1,jpphn) * xthetapn 
     967            CALL asm_bgc_unlog_2d( zbkg_phynon, slphynon_bkginc, zinc_phynon ) 
     968         ELSE 
     969            zinc_phynon(:,:) = 0.0 
     970         ENDIF 
     971#endif 
     972 
     973         ! Select mixed layer 
     974         IF ( ll_asmdin ) THEN 
     975            CALL ctl_warn( ' Doing direct initialisation with ocean colour assimilation', & 
     976               &           ' Mixed layer depth taken to be background maximum mld_max_bkg' ) 
     977            zmld(:,:) = mld_max_bkg(:,:) 
     978         ELSE 
     979            SELECT CASE( mld_choice_bgc ) 
     980            CASE ( 1 )                   ! Turbocline/mixing depth [W points] 
     981               zmld(:,:) = hmld(:,:) 
     982            CASE ( 2 )                   ! Density criterion (0.01 kg/m^3 change from 10m) [W points] 
     983               zmld(:,:) = hmlp(:,:) 
     984            CASE ( 3 )                   ! Kara MLD [Interpolated] 
     985#if defined key_karaml 
     986               IF ( ln_kara ) THEN 
     987                  zmld(:,:) = hmld_kara(:,:) 
     988               ELSE 
     989                  CALL ctl_stop( ' Kara mixed layer requested for phyto2d assimilation,', & 
     990                     &           ' but ln_kara=.false.' ) 
     991               ENDIF 
     992#else 
     993               CALL ctl_stop( ' Kara mixed layer requested for phyto2d assimilation,', & 
     994                  &           ' but is not defined' ) 
     995#endif 
     996            CASE ( 4 )                   ! Temperature criterion (0.2 K change from surface) [T points] 
     997               !zmld(:,:) = hmld_tref(:,:) 
     998               CALL ctl_stop( ' hmld_tref mixed layer requested for LogChl assimilation,', & 
     999                  &           ' but is not available in this version' ) 
     1000            CASE ( 5 )                   ! Density criterion (0.01 kg/m^3 change from 10m) [T points] 
     1001               zmld(:,:) = hmlpt(:,:) 
     1002            END SELECT 
     1003         ENDIF 
     1004 
    8531005         zincper = (nitiaufin_r - nitiaustr_r + 1) * rdt 
    8541006 
    8551007#if defined key_medusa && defined key_foam_medusa 
    856          CALL asm_logchl_bal_medusa( slchltot_bkginc, zincper, mld_choice_bgc, & 
    857             &                        rn_maxchlinc, ln_phytobal, ll_asmdin,  & 
    858             &                        pgrow_avg_bkg, ploss_avg_bkg,           & 
    859             &                        phyt_avg_bkg, mld_max_bkg,              & 
     1008         CALL asm_logchl_bal_medusa( (ln_slchltotinc .OR. ln_schltotinc), & 
     1009            &                        zinc_chltot,                         & 
     1010            &                        ln_slchldiainc,                      & 
     1011            &                        zinc_chldia,                         & 
     1012            &                        ln_slchlnoninc,                      & 
     1013            &                        zinc_chlnon,                         & 
     1014            &                        ln_slphytotinc,                      & 
     1015            &                        zinc_phytot,                         & 
     1016            &                        ln_slphydiainc,                      & 
     1017            &                        zinc_phydia,                         & 
     1018            &                        ln_slphynoninc,                      & 
     1019            &                        zinc_phynon,                         & 
     1020            &                        zincper,                             & 
     1021            &                        rn_maxchlinc, ln_phytobal, zmld,     & 
     1022            &                        pgrow_avg_bkg, ploss_avg_bkg,        & 
     1023            &                        phyt_avg_bkg, mld_max_bkg,           & 
    8601024            &                        tracer_bkg, phyto2d_balinc ) 
    8611025#elif defined key_hadocc 
    862          CALL asm_logchl_bal_hadocc( slchltot_bkginc, zincper, mld_choice_bgc, & 
    863             &                        rn_maxchlinc, ln_phytobal, ll_asmdin,  & 
    864             &                        pgrow_avg_bkg, ploss_avg_bkg,           & 
    865             &                        phyt_avg_bkg, mld_max_bkg,              & 
    866             &                        chl_bkg(:,:,1), cchl_p_bkg(:,:,1),      & 
     1026         CALL asm_logchl_bal_hadocc( (ln_slchltotinc .OR. ln_schltotinc), & 
     1027            &                        zinc_chltot,                         & 
     1028            &                        ln_slphytotinc,                      & 
     1029            &                        zinc_phytot,                         & 
     1030            &                        zincper,                             & 
     1031            &                        rn_maxchlinc, ln_phytobal, zmld,     & 
     1032            &                        pgrow_avg_bkg, ploss_avg_bkg,        & 
     1033            &                        phyt_avg_bkg, mld_max_bkg,           & 
     1034            &                        cchl_p_bkg(:,:,1),                   & 
    8671035            &                        tracer_bkg, phyto2d_balinc ) 
    8681036#else 
    869          CALL ctl_stop( 'Attempting to assimilate slchltot, ', & 
     1037         CALL ctl_stop( 'Attempting to assimilate phyto2d data, ', & 
    8701038            &           'but not defined a biogeochemical model' ) 
    8711039#endif 
     
    8871055            IF(lwp) THEN 
    8881056               WRITE(numout,*)  
    889                WRITE(numout,*) 'logchl_asm_inc : logchl IAU at time step = ', & 
     1057               WRITE(numout,*) 'phyto2d_asm_inc : phyto2d IAU at time step = ', & 
    8901058                  &  kt,' with IAU weight = ', pwgtiau(it) 
    8911059               WRITE(numout,*) '~~~~~~~~~~~~' 
  • branches/UKMO/dev_r5518_GO6_package_asm_3d_bgc/NEMOGCM/NEMO/OPA_SRC/ASM/asmlogchlbal_hadocc.F90

    r9292 r9431  
    2222   USE par_oce,       ONLY: jpi, jpj, jpk  ! domain array sizes 
    2323   USE dom_oce,       ONLY: gdepw_n        ! domain information 
    24    USE zdfmxl                              ! mixed layer depth 
    2524   USE iom                                 ! i/o 
    2625   USE par_hadocc                          ! HadOCC parameters 
     
    6867CONTAINS 
    6968 
    70    SUBROUTINE asm_logchl_bal_hadocc( logchl_bkginc, aincper, mld_choice_bgc, & 
    71       &                              k_maxchlinc, ld_logchlbal, ld_asmdin,   & 
    72       &                              pgrow_avg_bkg, ploss_avg_bkg,           & 
    73       &                              phyt_avg_bkg, mld_max_bkg,              & 
    74       &                              chl_bkg, cchl_p_bkg,                    & 
    75       &                              tracer_bkg, logchl_balinc ) 
     69   SUBROUTINE asm_logchl_bal_hadocc( ld_chltot,                      & 
     70      &                              pinc_chltot,                    & 
     71      &                              ld_phytot,                      & 
     72      &                              pinc_phytot,                    & 
     73      &                              pincper,                        & 
     74      &                              p_maxchlinc, ld_phytobal, pmld, & 
     75      &                              pgrow_avg_bkg, ploss_avg_bkg,   & 
     76      &                              phyt_avg_bkg, mld_max_bkg,      & 
     77      &                              cchl_p_bkg,                     & 
     78      &                              tracer_bkg, phyto2d_balinc ) 
    7679      !!--------------------------------------------------------------------------- 
    7780      !!                    ***  ROUTINE asm_logchl_bal_hadocc  *** 
    7881      !! 
    79       !! ** Purpose :   calculate increments to HadOCC from logchl increments 
    80       !! 
    81       !! ** Method  :   convert logchl increments to chl increments 
    82       !!                call nitrogen balancing scheme 
    83       !! 
    84       !! ** Action  :   populate logchl_balinc 
     82      !! ** Purpose :   calculate increments to HadOCC from 2d phytoplankton increments 
     83      !! 
     84      !! ** Method  :   call nitrogen balancing scheme 
     85      !! 
     86      !! ** Action  :   populate phyto2d_balinc 
    8587      !! 
    8688      !! References :   Hemmings et al., 2008, J. Mar. Res. 
     
    8890      !!--------------------------------------------------------------------------- 
    8991      !! 
    90       REAL(wp), INTENT(in   ), DIMENSION(jpi,jpj)           :: logchl_bkginc  ! logchl increments 
    91       REAL(wp), INTENT(in   )                               :: aincper        ! Assimilation period 
    92       INTEGER,  INTENT(in   )                               :: mld_choice_bgc ! MLD criterion 
    93       REAL(wp), INTENT(in   )                               :: k_maxchlinc    ! Max chl increment 
    94       LOGICAL,  INTENT(in   )                               :: ld_logchlbal   ! Balancing y/n 
    95       LOGICAL,  INTENT(in   )                               :: ld_asmdin      ! Direct initialisation y/n 
     92      LOGICAL,  INTENT(in   )                               :: ld_chltot      ! Assim chltot y/n 
     93      REAL(wp), INTENT(inout), DIMENSION(jpi,jpj)           :: pinc_chltot    ! chltot increments 
     94      LOGICAL,  INTENT(in   )                               :: ld_phytot      ! Assim phytot y/n 
     95      REAL(wp), INTENT(inout), DIMENSION(jpi,jpj)           :: pinc_phytot    ! phytot increments 
     96      REAL(wp), INTENT(in   )                               :: pincper        ! Assimilation period 
     97      REAL(wp), INTENT(in   )                               :: p_maxchlinc    ! Max chl increment 
     98      LOGICAL,  INTENT(in   )                               :: ld_phytobal    ! Balancing y/n 
     99      REAL(wp), INTENT(inout), DIMENSION(jpi,jpj)           :: pmld           ! Mixed layer depth 
    96100      REAL(wp), INTENT(in   ), DIMENSION(jpi,jpj)           :: pgrow_avg_bkg  ! Avg phyto growth 
    97101      REAL(wp), INTENT(in   ), DIMENSION(jpi,jpj)           :: ploss_avg_bkg  ! Avg phyto loss 
    98102      REAL(wp), INTENT(in   ), DIMENSION(jpi,jpj)           :: phyt_avg_bkg   ! Avg phyto 
    99103      REAL(wp), INTENT(in   ), DIMENSION(jpi,jpj)           :: mld_max_bkg    ! Max MLD 
    100       REAL(wp), INTENT(in   ), DIMENSION(jpi,jpj)           :: chl_bkg        ! Surface chlorophyll 
    101104      REAL(wp), INTENT(in   ), DIMENSION(jpi,jpj)           :: cchl_p_bkg     ! Surface C:Chl 
    102105      REAL(wp), INTENT(in   ), DIMENSION(jpi,jpj,jpk,jptra) :: tracer_bkg     ! State variables 
    103       REAL(wp), INTENT(  out), DIMENSION(jpi,jpj,jpk,jptra) :: logchl_balinc ! Balancing increments 
     106      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj,jpk,jptra) :: phyto2d_balinc ! Balancing increments 
    104107      !! 
    105108      INTEGER                                               :: ji, jj, jk, jn ! Loop counters 
    106109      INTEGER                                               :: jkmax          ! Loop index 
    107110      INTEGER,                 DIMENSION(6)                 :: i_tracer       ! Tracer indices 
    108       REAL(wp),                DIMENSION(jpi,jpj)           :: chl_inc        ! Chlorophyll increments 
    109       REAL(wp),                DIMENSION(jpi,jpj)           :: zmld           ! Mixed layer depth 
    110111      REAL(wp),                DIMENSION(16)                :: modparm        ! Model parameters 
    111112      REAL(wp),                DIMENSION(20)                :: assimparm      ! Assimilation parameters 
     
    116117      !!--------------------------------------------------------------------------- 
    117118 
    118       ! Convert log10(chlorophyll) increment back to a chlorophyll increment 
    119       ! In order to transform logchl incs to chl incs, need to account for model 
    120       ! background, cannot simply do 10^logchl_bkginc. Need to: 
    121       ! 1) Add logchl inc to log10(background) to get log10(analysis) 
    122       ! 2) Take 10^log10(analysis) to get analysis 
    123       ! 3) Subtract background from analysis to get chl incs 
    124       ! If k_maxchlinc > 0 then cap total absolute chlorophyll increment at that value 
    125       DO jj = 1, jpj 
    126          DO ji = 1, jpi 
    127             IF ( chl_bkg(ji,jj) > 0.0 ) THEN 
    128                chl_inc(ji,jj) = 10**( LOG10( chl_bkg(ji,jj) ) + logchl_bkginc(ji,jj) ) - chl_bkg(ji,jj) 
    129                IF ( k_maxchlinc > 0.0 ) THEN 
    130                   chl_inc(ji,jj) = MAX( -1.0 * k_maxchlinc, MIN( chl_inc(ji,jj), k_maxchlinc ) ) 
     119      IF ( ld_chltot ) THEN 
     120         ! If p_maxchlinc > 0 then cap total absolute chlorophyll increment at that value 
     121         DO jj = 1, jpj 
     122            DO ji = 1, jpi 
     123               IF ( p_maxchlinc > 0.0 ) THEN 
     124                  pinc_chltot(ji,jj) = MAX( -1.0 * p_maxchlinc, MIN( pinc_chltot(ji,jj), p_maxchlinc ) ) 
    131125               ENDIF 
    132             ELSE 
    133                chl_inc(ji,jj) = 0.0 
    134             ENDIF 
     126            END DO 
    135127         END DO 
    136       END DO 
    137        
    138       ! Select mixed layer 
    139       IF ( ld_asmdin ) THEN 
    140          CALL ctl_warn( ' Doing direct initialisation of HadOCC with chlorophyll assimilation', & 
    141             &           ' Mixed layer depth taken to be background maximum mld_max_bkg' ) 
    142          zmld(:,:) = mld_max_bkg(:,:) 
    143128      ELSE 
    144          SELECT CASE( mld_choice_bgc ) 
    145          CASE ( 1 )                   ! Turbocline/mixing depth [W points] 
    146             zmld(:,:) = hmld(:,:) 
    147          CASE ( 2 )                   ! Density criterion (0.01 kg/m^3 change from 10m) [W points] 
    148             zmld(:,:) = hmlp(:,:) 
    149          CASE ( 3 )                   ! Kara MLD [Interpolated] 
    150 #if defined key_karaml 
    151             IF ( ln_kara ) THEN 
    152                zmld(:,:) = hmld_kara(:,:) 
    153             ELSE 
    154                CALL ctl_stop( ' Kara mixed layer requested for LogChl assimilation,', & 
    155                   &           ' but ln_kara=.false.' ) 
    156             ENDIF 
    157 #else 
    158             CALL ctl_stop( ' Kara mixed layer requested for LogChl assimilation,', & 
    159                &           ' but is not defined' ) 
    160 #endif 
    161          CASE ( 4 )                   ! Temperature criterion (0.2 K change from surface) [T points] 
    162             !zmld(:,:) = hmld_tref(:,:) 
    163             CALL ctl_stop( ' hmld_tref mixed layer requested for LogChl assimilation,', & 
    164                &           ' but is not available in this version' ) 
    165          CASE ( 5 )                   ! Density criterion (0.01 kg/m^3 change from 10m) [T points] 
    166             zmld(:,:) = hmlpt(:,:) 
    167          END SELECT 
     129         CALL ctl_stop( ' No PFT assimilation quite yet' ) 
    168130      ENDIF 
    169131       
    170       IF ( ld_logchlbal ) THEN   ! Nitrogen balancing 
     132      IF ( ld_phytobal ) THEN   ! Nitrogen balancing 
    171133 
    172134         ! Set up model parameters to be passed into Hemmings balancing routine 
     
    229191         CALL bio_analysis( jpi, jpj, jpk, ZDZ(:,:,:), i_tracer, modparm,               & 
    230192            &               n2be_p, n2be_z, n2be_d, assimparm,                          & 
    231             &               INT(aincper), 1, kmt(:,:), tmask(:,:,:),                    & 
    232             &               zmld(:,:), mld_max_bkg(:,:), chl_inc(:,:), cchl_p_bkg(:,:), & 
     193            &               INT(pincper), 1, kmt(:,:), tmask(:,:,:),                    & 
     194            &               pmld(:,:), mld_max_bkg(:,:), pinc_chltot(:,:), cchl_p_bkg(:,:), & 
    233195            &               nbal_active, phyt_avg_bkg(:,:),                             & 
    234196            &               gl_active, pgrow_avg_bkg(:,:), ploss_avg_bkg(:,:),          & 
     
    240202 
    241203         ! Save balancing increments 
    242          logchl_balinc(:,:,:,jp_had_nut) = outincs(:,:,:,i_tracer(1)) 
    243          logchl_balinc(:,:,:,jp_had_phy) = outincs(:,:,:,i_tracer(2)) 
    244          logchl_balinc(:,:,:,jp_had_zoo) = outincs(:,:,:,i_tracer(3)) 
    245          logchl_balinc(:,:,:,jp_had_pdn) = outincs(:,:,:,i_tracer(4)) 
    246          logchl_balinc(:,:,:,jp_had_dic) = outincs(:,:,:,i_tracer(5)) 
    247          logchl_balinc(:,:,:,jp_had_alk) = outincs(:,:,:,i_tracer(6)) 
     204         phyto2d_balinc(:,:,:,jp_had_nut) = outincs(:,:,:,i_tracer(1)) 
     205         phyto2d_balinc(:,:,:,jp_had_phy) = outincs(:,:,:,i_tracer(2)) 
     206         phyto2d_balinc(:,:,:,jp_had_zoo) = outincs(:,:,:,i_tracer(3)) 
     207         phyto2d_balinc(:,:,:,jp_had_pdn) = outincs(:,:,:,i_tracer(4)) 
     208         phyto2d_balinc(:,:,:,jp_had_dic) = outincs(:,:,:,i_tracer(5)) 
     209         phyto2d_balinc(:,:,:,jp_had_alk) = outincs(:,:,:,i_tracer(6)) 
    248210       
    249211      ELSE   ! No nitrogen balancing 
    250212       
    251213         ! Initialise phytoplankton increment to zero 
    252          logchl_balinc(:,:,:,jp_had_phy) = 0.0 
     214         phyto2d_balinc(:,:,:,jp_had_phy) = 0.0 
    253215          
    254216         ! Convert surface chlorophyll increment to phytoplankton nitrogen 
    255          logchl_balinc(:,:,1,jp_had_phy) = ( cchl_p_bkg(:,:) / (mw_carbon * c2n_p) ) * chl_inc(:,:) 
     217         phyto2d_balinc(:,:,1,jp_had_phy) = ( cchl_p_bkg(:,:) / (mw_carbon * c2n_p) ) * pinc_chltot(:,:) 
    256218          
    257219         ! Propagate through mixed layer 
     
    261223               jkmax = jpk-1 
    262224               DO jk = jpk-1, 1, -1 
    263                   IF ( ( zmld(ji,jj) >  gdepw_n(ji,jj,jk)   ) .AND. & 
    264                      & ( zmld(ji,jj) <= gdepw_n(ji,jj,jk+1) ) ) THEN 
    265                      zmld(ji,jj) = gdepw_n(ji,jj,jk+1) 
     225                  IF ( ( pmld(ji,jj) >  gdepw_n(ji,jj,jk)   ) .AND. & 
     226                     & ( pmld(ji,jj) <= gdepw_n(ji,jj,jk+1) ) ) THEN 
     227                     pmld(ji,jj) = gdepw_n(ji,jj,jk+1) 
    266228                     jkmax = jk 
    267229                  ENDIF 
     
    269231               ! 
    270232               DO jk = 2, jkmax 
    271                   logchl_balinc(ji,jj,jk,jp_had_phy) = logchl_balinc(ji,jj,1,jp_had_phy) 
     233                  phyto2d_balinc(ji,jj,jk,jp_had_phy) = phyto2d_balinc(ji,jj,1,jp_had_phy) 
    272234               END DO 
    273235               ! 
     
    276238 
    277239         ! Set other balancing increments to zero 
    278          logchl_balinc(:,:,:,jp_had_nut) = 0.0 
    279          logchl_balinc(:,:,:,jp_had_zoo) = 0.0 
    280          logchl_balinc(:,:,:,jp_had_pdn) = 0.0 
    281          logchl_balinc(:,:,:,jp_had_dic) = 0.0 
    282          logchl_balinc(:,:,:,jp_had_alk) = 0.0 
     240         phyto2d_balinc(:,:,:,jp_had_nut) = 0.0 
     241         phyto2d_balinc(:,:,:,jp_had_zoo) = 0.0 
     242         phyto2d_balinc(:,:,:,jp_had_pdn) = 0.0 
     243         phyto2d_balinc(:,:,:,jp_had_dic) = 0.0 
     244         phyto2d_balinc(:,:,:,jp_had_alk) = 0.0 
    283245       
    284246      ENDIF 
     
    291253   !!---------------------------------------------------------------------- 
    292254CONTAINS 
    293    SUBROUTINE asm_logchl_bal_hadocc( logchl_bkginc, aincper, mld_choice_bgc, & 
    294       &                              k_maxchlinc, ld_logchlbal,              & 
    295       &                              pgrow_avg_bkg, ploss_avg_bkg,           & 
    296       &                              phyt_avg_bkg, mld_max_bkg,              & 
    297       &                              chl_bkg, cchl_p_bkg,                    & 
    298       &                              tracer_bkg, logchl_balinc ) 
    299       REAL    :: logchl_bkginc(:,:) 
    300       REAL    :: aincper 
    301       INTEGER :: mld_choice_bgc 
    302       REAL    :: k_maxchlinc 
    303       LOGICAL :: ld_logchlbal 
     255   SUBROUTINE asm_logchl_bal_hadocc( ld_chltot,                      & 
     256      &                              pinc_chltot,                    & 
     257      &                              ld_phytot,                      & 
     258      &                              pinc_phytot,                    & 
     259      &                              pincper,                        & 
     260      &                              p_maxchlinc, ld_phytobal, pmld, & 
     261      &                              pgrow_avg_bkg, ploss_avg_bkg,   & 
     262      &                              phyt_avg_bkg, mld_max_bkg,      & 
     263      &                              cchl_p_bkg,                     & 
     264      &                              tracer_bkg, phyto2d_balinc ) 
     265      LOGICAL :: ld_chltot 
     266      REAL    :: pinc_chltot(:,:) 
     267      LOGICAL :: ld_phytot 
     268      REAL    :: pinc_phytot(:,:) 
     269      REAL    :: pincper 
     270      REAL    :: p_maxchlinc 
     271      LOGICAL :: ld_phytobal 
     272      REAL    :: pmld(:,:) 
    304273      REAL    :: pgrow_avg_bkg(:,:) 
    305274      REAL    :: ploss_avg_bkg(:,:) 
    306275      REAL    :: phyt_avg_bkg(:,:) 
    307276      REAL    :: mld_max_bkg(:,:) 
    308       REAL    :: chl_bkg(:,:) 
    309277      REAL    :: cchl_p_bkg(:,:) 
    310278      REAL    :: tracer_bkg(:,:,:,:) 
    311       REAL    :: logchl_balinc(:,:,:,:) 
     279      REAL    :: phyto2d_balinc(:,:,:,:) 
    312280      WRITE(*,*) 'asm_logchl_bal_hadocc: You should not have seen this print! error?' 
    313281   END SUBROUTINE asm_logchl_bal_hadocc 
  • branches/UKMO/dev_r5518_GO6_package_asm_3d_bgc/NEMOGCM/NEMO/OPA_SRC/ASM/asmlogchlbal_medusa.F90

    r9292 r9431  
    2323   USE par_oce,       ONLY: jpi, jpj, jpk  ! domain array sizes 
    2424   USE dom_oce,       ONLY: gdepw_n        ! domain information 
    25    USE zdfmxl                              ! mixed layer depth 
    2625   USE zdftmx,        ONLY: ln_tmx_itf, &  ! Indonesian Throughflow 
    2726      &                     mask_itf       ! tidal mixing mask 
     
    7069CONTAINS 
    7170 
    72    SUBROUTINE asm_logchl_bal_medusa( logchl_bkginc, aincper, mld_choice_bgc, & 
    73       &                              k_maxchlinc, ld_logchlbal, ld_asmdin,   & 
    74       &                              pgrow_avg_bkg, ploss_avg_bkg,           & 
    75       &                              phyt_avg_bkg, mld_max_bkg,              & 
    76       &                              tracer_bkg, logchl_balinc ) 
     71   SUBROUTINE asm_logchl_bal_medusa( ld_chltot,                      & 
     72      &                              pinc_chltot,                    & 
     73      &                              ld_chldia,                      & 
     74      &                              pinc_chldia,                    & 
     75      &                              ld_chlnon,                      & 
     76      &                              pinc_chlnon,                    & 
     77      &                              ld_phytot,                      & 
     78      &                              pinc_phytot,                    & 
     79      &                              ld_phydia,                      & 
     80      &                              pinc_phydia,                    & 
     81      &                              ld_phynon,                      & 
     82      &                              pinc_phynon,                    & 
     83      &                              pincper,                        & 
     84      &                              p_maxchlinc, ld_phytobal, pmld, & 
     85      &                              pgrow_avg_bkg, ploss_avg_bkg,   & 
     86      &                              phyt_avg_bkg, mld_max_bkg,      & 
     87      &                              tracer_bkg, phyto2d_balinc ) 
    7788      !!--------------------------------------------------------------------------- 
    7889      !!                    ***  ROUTINE asm_logchl_bal_medusa  *** 
    7990      !! 
    80       !! ** Purpose :   calculate increments to MEDUSA from logchl increments 
    81       !! 
    82       !! ** Method  :   convert logchl increments to chl increments 
    83       !!                average up MEDUSA to look like HadOCC 
     91      !! ** Purpose :   calculate increments to MEDUSA from 2d phytoplankton increments 
     92      !! 
     93      !! ** Method  :   average up MEDUSA to look like HadOCC 
    8494      !!                call nitrogen balancing scheme 
    8595      !!                separate back out to MEDUSA 
    8696      !! 
    87       !! ** Action  :   populate logchl_balinc 
     97      !! ** Action  :   populate phyto2d_balinc 
    8898      !! 
    8999      !! References :   Hemmings et al., 2008, J. Mar. Res. 
     
    91101      !!--------------------------------------------------------------------------- 
    92102      !! 
    93       REAL(wp), INTENT(in   ), DIMENSION(jpi,jpj)           :: logchl_bkginc  ! logchl increments 
    94       REAL(wp), INTENT(in   )                               :: aincper        ! Assimilation period 
    95       INTEGER,  INTENT(in   )                               :: mld_choice_bgc ! MLD criterion 
    96       REAL(wp), INTENT(in   )                               :: k_maxchlinc    ! Max chl increment 
    97       LOGICAL,  INTENT(in   )                               :: ld_logchlbal   ! Balancing y/n 
    98       LOGICAL,  INTENT(in   )                               :: ld_asmdin      ! Direct initialisation y/n 
     103      LOGICAL,  INTENT(in   )                               :: ld_chltot      ! Assim chltot y/n 
     104      REAL(wp), INTENT(inout), DIMENSION(jpi,jpj)           :: pinc_chltot    ! chltot increments 
     105      LOGICAL,  INTENT(in   )                               :: ld_chldia      ! Assim chldia y/n 
     106      REAL(wp), INTENT(inout), DIMENSION(jpi,jpj)           :: pinc_chldia    ! chldia increments 
     107      LOGICAL,  INTENT(in   )                               :: ld_chlnon      ! Assim chlnon y/n 
     108      REAL(wp), INTENT(inout), DIMENSION(jpi,jpj)           :: pinc_chlnon    ! chlnon increments 
     109      LOGICAL,  INTENT(in   )                               :: ld_phytot      ! Assim phytot y/n 
     110      REAL(wp), INTENT(inout), DIMENSION(jpi,jpj)           :: pinc_phytot    ! phytot increments 
     111      LOGICAL,  INTENT(in   )                               :: ld_phydia      ! Assim phydia y/n 
     112      REAL(wp), INTENT(inout), DIMENSION(jpi,jpj)           :: pinc_phydia    ! phydia increments 
     113      LOGICAL,  INTENT(in   )                               :: ld_phynon      ! Assim phynon y/n 
     114      REAL(wp), INTENT(inout), DIMENSION(jpi,jpj)           :: pinc_phynon    ! phynon increments 
     115      REAL(wp), INTENT(in   )                               :: pincper        ! Assimilation period 
     116      REAL(wp), INTENT(in   )                               :: p_maxchlinc    ! Max chl increment 
     117      LOGICAL,  INTENT(in   )                               :: ld_phytobal    ! Balancing y/n 
     118      REAL(wp), INTENT(inout), DIMENSION(jpi,jpj)           :: pmld           ! Mixed layer depth 
    99119      REAL(wp), INTENT(in   ), DIMENSION(jpi,jpj)           :: pgrow_avg_bkg  ! Avg phyto growth 
    100120      REAL(wp), INTENT(in   ), DIMENSION(jpi,jpj)           :: ploss_avg_bkg  ! Avg phyto loss 
     
    102122      REAL(wp), INTENT(in   ), DIMENSION(jpi,jpj)           :: mld_max_bkg    ! Max MLD 
    103123      REAL(wp), INTENT(in   ), DIMENSION(jpi,jpj,jpk,jptra) :: tracer_bkg     ! State variables 
    104       REAL(wp), INTENT(  out), DIMENSION(jpi,jpj,jpk,jptra) :: logchl_balinc ! Balancing increments 
     124      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj,jpk,jptra) :: phyto2d_balinc ! Balancing increments 
    105125      !! 
    106126      INTEGER                                               :: ji, jj, jk, jn ! Loop counters 
     
    120140      REAL(wp)                                              :: zrat_chn_phn   ! Ratio of jpchn:jpphn 
    121141      REAL(wp)                                              :: zrat_dtc_det   ! Ratio of jpdtc:jpdet 
    122       REAL(wp),                DIMENSION(jpi,jpj)           :: chl_inc        ! Chlorophyll increments 
    123142      REAL(wp),                DIMENSION(jpi,jpj)           :: medusa_chl     ! MEDUSA total chlorophyll 
    124143      REAL(wp),                DIMENSION(jpi,jpj)           :: cchl_p         ! C:Chl for total phy 
    125       REAL(wp),                DIMENSION(jpi,jpj)           :: zmld           ! Mixed layer depth 
    126144      REAL(wp),                DIMENSION(16)                :: modparm        ! Model parameters 
    127145      REAL(wp),                DIMENSION(20)                :: assimparm      ! Assimilation parameters 
     
    132150      !!--------------------------------------------------------------------------- 
    133151 
    134       ! Convert log10(chlorophyll) increment back to a chlorophyll increment 
    135       ! In order to transform logchl incs to chl incs, need to account for model 
    136       ! background, cannot simply do 10^logchl_bkginc. Need to: 
    137       ! 1) Add logchl inc to log10(background) to get log10(analysis) 
    138       ! 2) Take 10^log10(analysis) to get analysis 
    139       ! 3) Subtract background from analysis to get chl incs 
    140       ! If k_maxchlinc > 0 then cap total absolute chlorophyll increment at that value 
    141       medusa_chl(:,:) = tracer_bkg(:,:,1,jpchn) + tracer_bkg(:,:,1,jpchd) 
    142       DO jj = 1, jpj 
    143          DO ji = 1, jpi 
    144             IF ( medusa_chl(ji,jj) > 0.0 ) THEN 
    145                chl_inc(ji,jj) = 10**( LOG10( medusa_chl(ji,jj) ) + logchl_bkginc(ji,jj) ) - medusa_chl(ji,jj) 
    146                IF ( k_maxchlinc > 0.0 ) THEN 
    147                   chl_inc(ji,jj) = MAX( -1.0 * k_maxchlinc, MIN( chl_inc(ji,jj), k_maxchlinc ) ) 
     152      IF ( ld_chltot ) THEN 
     153         ! If p_maxchlinc > 0 then cap total absolute chlorophyll increment at that value 
     154         DO jj = 1, jpj 
     155            DO ji = 1, jpi 
     156               IF ( p_maxchlinc > 0.0 ) THEN 
     157                  pinc_chltot(ji,jj) = MAX( -1.0 * p_maxchlinc, MIN( pinc_chltot(ji,jj), p_maxchlinc ) ) 
    148158               ENDIF 
    149             ELSE 
    150                chl_inc(ji,jj) = 0.0 
    151             ENDIF 
    152          END DO 
    153       END DO 
    154        
    155       ! Select mixed layer 
    156       IF ( ld_asmdin ) THEN 
    157          CALL ctl_warn( ' Doing direct initialisation of MEDUSA with chlorophyll assimilation', & 
    158             &           ' Mixed layer depth taken to be background maximum mld_max_bkg' ) 
    159          zmld(:,:) = mld_max_bkg(:,:) 
     159            END DO 
     160         END DO 
    160161      ELSE 
    161          SELECT CASE( mld_choice_bgc ) 
    162          CASE ( 1 )                   ! Turbocline/mixing depth [W points] 
    163             zmld(:,:) = hmld(:,:) 
    164          CASE ( 2 )                   ! Density criterion (0.01 kg/m^3 change from 10m) [W points] 
    165             zmld(:,:) = hmlp(:,:) 
    166          CASE ( 3 )                   ! Kara MLD [Interpolated] 
    167 #if defined key_karaml 
    168             IF ( ln_kara ) THEN 
    169                zmld(:,:) = hmld_kara(:,:) 
    170             ELSE 
    171                CALL ctl_stop( ' Kara mixed layer requested for LogChl assimilation,', & 
    172                   &           ' but ln_kara=.false.' ) 
    173             ENDIF 
    174 #else 
    175             CALL ctl_stop( ' Kara mixed layer requested for LogChl assimilation,', & 
    176                &           ' but is not defined' ) 
    177 #endif 
    178          CASE ( 4 )                   ! Temperature criterion (0.2 K change from surface) [T points] 
    179             !zmld(:,:) = hmld_tref(:,:) 
    180             CALL ctl_stop( ' hmld_tref mixed layer requested for LogChl assimilation,', & 
    181                &           ' but is not available in this version' ) 
    182          CASE ( 5 )                   ! Density criterion (0.01 kg/m^3 change from 10m) [T points] 
    183             zmld(:,:) = hmlpt(:,:) 
    184          END SELECT 
     162         CALL ctl_stop( ' No PFT assimilation quite yet' ) 
    185163      ENDIF 
    186164       
    187       IF ( ld_logchlbal ) THEN   ! Nitrogen balancing 
     165      IF ( ld_phytobal ) THEN   ! Nitrogen balancing 
    188166 
    189167         ! Set up model parameters to be passed into Hemmings balancing routine. 
     
    269247         CALL bio_analysis( jpi, jpj, jpk, gdepw_n(:,:,2:jpk), i_tracer, modparm,   & 
    270248            &               n2be_p, n2be_z, n2be_d, assimparm,                      & 
    271             &               INT(aincper), 1, INT(SUM(tmask,3)), tmask(:,:,:),       & 
    272             &               zmld(:,:), mld_max_bkg(:,:), chl_inc(:,:), cchl_p(:,:), & 
     249            &               INT(pincper), 1, INT(SUM(tmask,3)), tmask(:,:,:),       & 
     250            &               pmld(:,:), mld_max_bkg(:,:), pinc_chltot(:,:), cchl_p(:,:), & 
    273251            &               nbal_active, phyt_avg_bkg(:,:),                         & 
    274252            &               gl_active, pgrow_avg_bkg(:,:), ploss_avg_bkg(:,:),      & 
     
    280258          
    281259         ! Loop over each grid point partioning the increments 
    282          logchl_balinc(:,:,:,:) = 0.0 
     260         phyto2d_balinc(:,:,:,:) = 0.0 
    283261         DO jk = 1, jpk 
    284262            DO jj = 1, jpj 
     
    290268                     zfrac_phd = 1.0 - zfrac_phn 
    291269                     zrat_pds_phd = tracer_bkg(ji,jj,jk,jppds) / tracer_bkg(ji,jj,jk,jpphd) 
    292                      logchl_balinc(ji,jj,jk,jpphn) = outincs(ji,jj,jk,i_tracer(2)) * zfrac_phn 
    293                      logchl_balinc(ji,jj,jk,jpphd) = outincs(ji,jj,jk,i_tracer(2)) * zfrac_phd 
    294                      logchl_balinc(ji,jj,jk,jppds) = logchl_balinc(ji,jj,jk,jpphd) * zrat_pds_phd 
     270                     phyto2d_balinc(ji,jj,jk,jpphn) = outincs(ji,jj,jk,i_tracer(2)) * zfrac_phn 
     271                     phyto2d_balinc(ji,jj,jk,jpphd) = outincs(ji,jj,jk,i_tracer(2)) * zfrac_phd 
     272                     phyto2d_balinc(ji,jj,jk,jppds) = phyto2d_balinc(ji,jj,jk,jpphd) * zrat_pds_phd 
    295273 
    296274                     ! Chlorophyll split up based on existing ratios to phytoplankton nitrogen 
    297                      ! Not using chl_inc directly as it's only 2D 
    298                      ! This method should give same results at surface as splitting chl_inc would 
     275                     ! Not using pinc_chltot directly as it's only 2D 
     276                     ! This method should give same results at surface as splitting pinc_chltot would 
    299277                     zrat_chn_phn = tracer_bkg(ji,jj,jk,jpchn) / tracer_bkg(ji,jj,jk,jpphn) 
    300278                     zrat_chd_phd = tracer_bkg(ji,jj,jk,jpchd) / tracer_bkg(ji,jj,jk,jpphd) 
    301                      logchl_balinc(ji,jj,jk,jpchn) = logchl_balinc(ji,jj,jk,jpphn) * zrat_chn_phn 
    302                      logchl_balinc(ji,jj,jk,jpchd) = logchl_balinc(ji,jj,jk,jpphd) * zrat_chd_phd 
     279                     phyto2d_balinc(ji,jj,jk,jpchn) = phyto2d_balinc(ji,jj,jk,jpphn) * zrat_chn_phn 
     280                     phyto2d_balinc(ji,jj,jk,jpchd) = phyto2d_balinc(ji,jj,jk,jpphd) * zrat_chd_phd 
    303281                  ENDIF 
    304282 
     
    307285                     zfrac_zmi = tracer_bkg(ji,jj,jk,jpzmi) / (tracer_bkg(ji,jj,jk,jpzmi) + tracer_bkg(ji,jj,jk,jpzme)) 
    308286                     zfrac_zme = 1.0 - zfrac_zmi 
    309                      logchl_balinc(ji,jj,jk,jpzmi) = outincs(ji,jj,jk,i_tracer(3)) * zfrac_zmi 
    310                      logchl_balinc(ji,jj,jk,jpzme) = outincs(ji,jj,jk,i_tracer(3)) * zfrac_zme 
     287                     phyto2d_balinc(ji,jj,jk,jpzmi) = outincs(ji,jj,jk,i_tracer(3)) * zfrac_zmi 
     288                     phyto2d_balinc(ji,jj,jk,jpzme) = outincs(ji,jj,jk,i_tracer(3)) * zfrac_zme 
    311289                  ENDIF 
    312290 
    313291                  ! Nitrogen nutrient straight from balancing scheme 
    314                   logchl_balinc(ji,jj,jk,jpdin) = outincs(ji,jj,jk,i_tracer(1)) 
     292                  phyto2d_balinc(ji,jj,jk,jpdin) = outincs(ji,jj,jk,i_tracer(1)) 
    315293 
    316294                  ! Nitrogen detritus straight from balancing scheme 
    317                   logchl_balinc(ji,jj,jk,jpdet) = outincs(ji,jj,jk,i_tracer(4)) 
     295                  phyto2d_balinc(ji,jj,jk,jpdet) = outincs(ji,jj,jk,i_tracer(4)) 
    318296 
    319297                  ! DIC straight from balancing scheme 
    320                   logchl_balinc(ji,jj,jk,jpdic) = outincs(ji,jj,jk,i_tracer(5)) 
     298                  phyto2d_balinc(ji,jj,jk,jpdic) = outincs(ji,jj,jk,i_tracer(5)) 
    321299 
    322300                  ! Alkalinity straight from balancing scheme 
    323                   logchl_balinc(ji,jj,jk,jpalk) = outincs(ji,jj,jk,i_tracer(6)) 
     301                  phyto2d_balinc(ji,jj,jk,jpalk) = outincs(ji,jj,jk,i_tracer(6)) 
    324302 
    325303                  ! Remove diatom silicate increment from nutrient silicate to conserve mass 
    326                   IF ( ( tracer_bkg(ji,jj,jk,jpsil) - logchl_balinc(ji,jj,jk,jppds) ) > 0.0 ) THEN 
    327                      logchl_balinc(ji,jj,jk,jpsil) = logchl_balinc(ji,jj,jk,jppds) * (-1.0) 
     304                  IF ( ( tracer_bkg(ji,jj,jk,jpsil) - phyto2d_balinc(ji,jj,jk,jppds) ) > 0.0 ) THEN 
     305                     phyto2d_balinc(ji,jj,jk,jpsil) = phyto2d_balinc(ji,jj,jk,jppds) * (-1.0) 
    328306                  ENDIF 
    329307 
     
    331309                     ! Carbon detritus based on existing ratios 
    332310                     zrat_dtc_det = tracer_bkg(ji,jj,jk,jpdtc) / tracer_bkg(ji,jj,jk,jpdet) 
    333                      logchl_balinc(ji,jj,jk,jpdtc) = logchl_balinc(ji,jj,jk,jpdet) * zrat_dtc_det 
     311                     phyto2d_balinc(ji,jj,jk,jpdtc) = phyto2d_balinc(ji,jj,jk,jpdet) * zrat_dtc_det 
    334312                  ENDIF 
    335313 
    336314                  ! Do nothing with iron or oxygen for the time being 
    337                   logchl_balinc(ji,jj,jk,jpfer) = 0.0 
    338                   logchl_balinc(ji,jj,jk,jpoxy) = 0.0 
     315                  phyto2d_balinc(ji,jj,jk,jpfer) = 0.0 
     316                  phyto2d_balinc(ji,jj,jk,jpoxy) = 0.0 
    339317                   
    340318               END DO 
     
    345323       
    346324         ! Initialise individual chlorophyll increments to zero 
    347          logchl_balinc(:,:,:,jpchn) = 0.0 
    348          logchl_balinc(:,:,:,jpchd) = 0.0 
     325         phyto2d_balinc(:,:,:,jpchn) = 0.0 
     326         phyto2d_balinc(:,:,:,jpchd) = 0.0 
    349327          
    350328         ! Split up total surface chlorophyll increments 
     
    354332                  zfrac_chn = tracer_bkg(ji,jj,1,jpchn) / medusa_chl(ji,jj) 
    355333                  zfrac_chd = 1.0 - zfrac_chn 
    356                   logchl_balinc(ji,jj,1,jpchn) = chl_inc(ji,jj) * zfrac_chn 
    357                   logchl_balinc(ji,jj,1,jpchd) = chl_inc(ji,jj) * zfrac_chd 
     334                  phyto2d_balinc(ji,jj,1,jpchn) = pinc_chltot(ji,jj) * zfrac_chn 
     335                  phyto2d_balinc(ji,jj,1,jpchd) = pinc_chltot(ji,jj) * zfrac_chd 
    358336               ENDIF 
    359337            END DO 
     
    366344               jkmax = jpk-1 
    367345               DO jk = jpk-1, 1, -1 
    368                   IF ( ( zmld(ji,jj) >  gdepw_n(ji,jj,jk)   ) .AND. & 
    369                      & ( zmld(ji,jj) <= gdepw_n(ji,jj,jk+1) ) ) THEN 
    370                      zmld(ji,jj) = gdepw_n(ji,jj,jk+1) 
     346                  IF ( ( pmld(ji,jj) >  gdepw_n(ji,jj,jk)   ) .AND. & 
     347                     & ( pmld(ji,jj) <= gdepw_n(ji,jj,jk+1) ) ) THEN 
     348                     pmld(ji,jj) = gdepw_n(ji,jj,jk+1) 
    371349                     jkmax = jk 
    372350                  ENDIF 
     
    374352               ! 
    375353               DO jk = 2, jkmax 
    376                   logchl_balinc(ji,jj,jk,jpchn) = logchl_balinc(ji,jj,1,jpchn) 
    377                   logchl_balinc(ji,jj,jk,jpchd) = logchl_balinc(ji,jj,1,jpchd) 
     354                  phyto2d_balinc(ji,jj,jk,jpchn) = phyto2d_balinc(ji,jj,1,jpchn) 
     355                  phyto2d_balinc(ji,jj,jk,jpchd) = phyto2d_balinc(ji,jj,1,jpchd) 
    378356               END DO 
    379357               ! 
     
    382360 
    383361         ! Set other balancing increments to zero 
    384          logchl_balinc(:,:,:,jpphn) = 0.0 
    385          logchl_balinc(:,:,:,jpphd) = 0.0 
    386          logchl_balinc(:,:,:,jppds) = 0.0 
    387          logchl_balinc(:,:,:,jpzmi) = 0.0 
    388          logchl_balinc(:,:,:,jpzme) = 0.0 
    389          logchl_balinc(:,:,:,jpdin) = 0.0 
    390          logchl_balinc(:,:,:,jpsil) = 0.0 
    391          logchl_balinc(:,:,:,jpfer) = 0.0 
    392          logchl_balinc(:,:,:,jpdet) = 0.0 
    393          logchl_balinc(:,:,:,jpdtc) = 0.0 
    394          logchl_balinc(:,:,:,jpdic) = 0.0 
    395          logchl_balinc(:,:,:,jpalk) = 0.0 
    396          logchl_balinc(:,:,:,jpoxy) = 0.0 
     362         phyto2d_balinc(:,:,:,jpphn) = 0.0 
     363         phyto2d_balinc(:,:,:,jpphd) = 0.0 
     364         phyto2d_balinc(:,:,:,jppds) = 0.0 
     365         phyto2d_balinc(:,:,:,jpzmi) = 0.0 
     366         phyto2d_balinc(:,:,:,jpzme) = 0.0 
     367         phyto2d_balinc(:,:,:,jpdin) = 0.0 
     368         phyto2d_balinc(:,:,:,jpsil) = 0.0 
     369         phyto2d_balinc(:,:,:,jpfer) = 0.0 
     370         phyto2d_balinc(:,:,:,jpdet) = 0.0 
     371         phyto2d_balinc(:,:,:,jpdtc) = 0.0 
     372         phyto2d_balinc(:,:,:,jpdic) = 0.0 
     373         phyto2d_balinc(:,:,:,jpalk) = 0.0 
     374         phyto2d_balinc(:,:,:,jpoxy) = 0.0 
    397375 
    398376      ENDIF 
     
    404382         DO jn = 1, jptra 
    405383            DO jk = 1, jpk 
    406                logchl_balinc(:,:,jk,jn) = logchl_balinc(:,:,jk,jn) * ( 1.0 - mask_itf(:,:) ) 
     384               phyto2d_balinc(:,:,jk,jn) = phyto2d_balinc(:,:,jk,jn) * ( 1.0 - mask_itf(:,:) ) 
    407385            END DO 
    408386         END DO 
     
    416394   !!---------------------------------------------------------------------- 
    417395CONTAINS 
    418    SUBROUTINE asm_logchl_bal_medusa( logchl_bkginc, aincper, mld_choice_bgc, & 
    419       &                              k_maxchlinc, ld_logchlbal,              & 
    420       &                              pgrow_avg_bkg, ploss_avg_bkg,           & 
    421       &                              phyt_avg_bkg, mld_max_bkg,              & 
    422       &                              tracer_bkg, logchl_balinc ) 
    423       REAL    :: logchl_bkginc(:,:) 
    424       REAL    :: aincper 
    425       INTEGER :: mld_choice_bgc 
    426       REAL    :: k_maxchlinc 
    427       LOGICAL :: ld_logchlbal 
     396   SUBROUTINE asm_logchl_bal_medusa( ld_chltot,                      & 
     397      &                              pinc_chltot,                    & 
     398      &                              ld_chldia,                      & 
     399      &                              pinc_chldia,                    & 
     400      &                              ld_chlnon,                      & 
     401      &                              pinc_chlnon,                    & 
     402      &                              ld_phytot,                      & 
     403      &                              pinc_phytot,                    & 
     404      &                              ld_phydia,                      & 
     405      &                              pinc_phydia,                    & 
     406      &                              ld_phynon,                      & 
     407      &                              pinc_phynon,                    & 
     408      &                              pincper,                        & 
     409      &                              p_maxchlinc, ld_phytobal, pmld, & 
     410      &                              pgrow_avg_bkg, ploss_avg_bkg,   & 
     411      &                              phyt_avg_bkg, mld_max_bkg,      & 
     412      &                              tracer_bkg, phyto2d_balinc ) 
     413      LOGICAL :: ld_chltot 
     414      REAL    :: pinc_chltot(:,:) 
     415      LOGICAL :: ld_chldia 
     416      REAL    :: pinc_chldia(:,:) 
     417      LOGICAL :: ld_chlnon 
     418      REAL    :: pinc_chlnon(:,:) 
     419      LOGICAL :: ld_phytot 
     420      REAL    :: pinc_phytot(:,:) 
     421      LOGICAL :: ld_phydia 
     422      REAL    :: pinc_phydia(:,:) 
     423      LOGICAL :: ld_phynon 
     424      REAL    :: pinc_phynon(:,:) 
     425      REAL    :: pincper 
     426      REAL    :: p_maxchlinc 
     427      LOGICAL :: ld_phytobal 
     428      REAL    :: pmld(:,:) 
    428429      REAL    :: pgrow_avg_bkg(:,:) 
    429430      REAL    :: ploss_avg_bkg(:,:) 
     
    431432      REAL    :: mld_max_bkg(:,:) 
    432433      REAL    :: tracer_bkg(:,:,:,:) 
    433       REAL    :: logchl_balinc(:,:,:,:) 
     434      REAL    :: phyto2d_balinc(:,:,:,:) 
    434435      WRITE(*,*) 'asm_logchl_bal_medusa: You should not have seen this print! error?' 
    435436   END SUBROUTINE asm_logchl_bal_medusa 
Note: See TracChangeset for help on using the changeset viewer.