Skip to content
Snippets Groups Projects
rain_ice.F90 84.2 KiB
Newer Older
  • Learn to ignore specific revisions
  •               IF (IITER(JL)<INB_ITER_MAX .AND. ZMAXB(JL)>PARAMI%XMRSTEP) THEN
    
                  ENDIF
                ENDDO
              ENDIF ! LL_ANY_ITER
            ENDIF ! XMRSTEP/=0.
            !
            !***       4.3 New values of variables for next iteration
            !
    
                ZVART(JL, JV)=ZVART(JL, JV)+ZA(JL, JV)*ZMAXTIME(JL)+ZB(JL, JV)
    
    #ifdef REPRO55
              ZCIT(JL)=ZCIT(JL) * MAX(0., -SIGN(1., -ZVART(JL,IRI)))
    #else
              IF (ZVART(JL,IRI)<=0.) ZCIT(JL) = 0.
    #endif
    
              ZTIME(JL)=ZTIME(JL)+ZMAXTIME(JL)
            ENDDO
    
            !
            !***       4.4 Mixing ratio change due to each process
            !
    
                ZTOT_RVHENI (JMICRO+JL-1)=ZTOT_RVHENI (JMICRO+JL-1)+ZRVHENI_MR(JL)
                ZTOT_RCHONI (JMICRO+JL-1)=ZTOT_RCHONI (JMICRO+JL-1)+ZRCHONI   (JL)*ZMAXTIME(JL)
                ZTOT_RRHONG (JMICRO+JL-1)=ZTOT_RRHONG (JMICRO+JL-1)+ZRRHONG_MR(JL)
                ZTOT_RVDEPS (JMICRO+JL-1)=ZTOT_RVDEPS (JMICRO+JL-1)+ZRVDEPS   (JL)*ZMAXTIME(JL)
                ZTOT_RIAGGS (JMICRO+JL-1)=ZTOT_RIAGGS (JMICRO+JL-1)+ZRIAGGS   (JL)*ZMAXTIME(JL)
                ZTOT_RIAUTS (JMICRO+JL-1)=ZTOT_RIAUTS (JMICRO+JL-1)+ZRIAUTS   (JL)*ZMAXTIME(JL)
                ZTOT_RVDEPG (JMICRO+JL-1)=ZTOT_RVDEPG (JMICRO+JL-1)+ZRVDEPG   (JL)*ZMAXTIME(JL)
                ZTOT_RCAUTR (JMICRO+JL-1)=ZTOT_RCAUTR (JMICRO+JL-1)+ZRCAUTR   (JL)*ZMAXTIME(JL)
                ZTOT_RCACCR (JMICRO+JL-1)=ZTOT_RCACCR (JMICRO+JL-1)+ZRCACCR   (JL)*ZMAXTIME(JL)
                ZTOT_RREVAV (JMICRO+JL-1)=ZTOT_RREVAV (JMICRO+JL-1)+ZRREVAV   (JL)*ZMAXTIME(JL)
                ZTOT_RCRIMSS(JMICRO+JL-1)=ZTOT_RCRIMSS(JMICRO+JL-1)+ZRCRIMSS  (JL)*ZMAXTIME(JL)
                ZTOT_RCRIMSG(JMICRO+JL-1)=ZTOT_RCRIMSG(JMICRO+JL-1)+ZRCRIMSG  (JL)*ZMAXTIME(JL)
                ZTOT_RSRIMCG(JMICRO+JL-1)=ZTOT_RSRIMCG(JMICRO+JL-1)+ZRSRIMCG  (JL)*ZMAXTIME(JL)+ZRSRIMCG_MR(JL)
                ZTOT_RRACCSS(JMICRO+JL-1)=ZTOT_RRACCSS(JMICRO+JL-1)+ZRRACCSS  (JL)*ZMAXTIME(JL)
                ZTOT_RRACCSG(JMICRO+JL-1)=ZTOT_RRACCSG(JMICRO+JL-1)+ZRRACCSG  (JL)*ZMAXTIME(JL)
                ZTOT_RSACCRG(JMICRO+JL-1)=ZTOT_RSACCRG(JMICRO+JL-1)+ZRSACCRG  (JL)*ZMAXTIME(JL)
                ZTOT_RSMLTG (JMICRO+JL-1)=ZTOT_RSMLTG (JMICRO+JL-1)+ZRSMLTG   (JL)*ZMAXTIME(JL)
                ZTOT_RCMLTSR(JMICRO+JL-1)=ZTOT_RCMLTSR(JMICRO+JL-1)+ZRCMLTSR  (JL)*ZMAXTIME(JL)
                ZTOT_RICFRRG(JMICRO+JL-1)=ZTOT_RICFRRG(JMICRO+JL-1)+ZRICFRRG  (JL)*ZMAXTIME(JL)
                ZTOT_RRCFRIG(JMICRO+JL-1)=ZTOT_RRCFRIG(JMICRO+JL-1)+ZRRCFRIG  (JL)*ZMAXTIME(JL)
                ZTOT_RICFRR (JMICRO+JL-1)=ZTOT_RICFRR (JMICRO+JL-1)+ZRICFRR   (JL)*ZMAXTIME(JL)
                ZTOT_RCWETG (JMICRO+JL-1)=ZTOT_RCWETG (JMICRO+JL-1)+ZRCWETG   (JL)*ZMAXTIME(JL)
                ZTOT_RIWETG (JMICRO+JL-1)=ZTOT_RIWETG (JMICRO+JL-1)+ZRIWETG   (JL)*ZMAXTIME(JL)
                ZTOT_RRWETG (JMICRO+JL-1)=ZTOT_RRWETG (JMICRO+JL-1)+ZRRWETG   (JL)*ZMAXTIME(JL)
                ZTOT_RSWETG (JMICRO+JL-1)=ZTOT_RSWETG (JMICRO+JL-1)+ZRSWETG   (JL)*ZMAXTIME(JL)
                ZTOT_RWETGH (JMICRO+JL-1)=ZTOT_RWETGH (JMICRO+JL-1)+ZRWETGH   (JL)*ZMAXTIME(JL)+ZRWETGH_MR(JL)
                ZTOT_RCDRYG (JMICRO+JL-1)=ZTOT_RCDRYG (JMICRO+JL-1)+ZRCDRYG   (JL)*ZMAXTIME(JL)
                ZTOT_RIDRYG (JMICRO+JL-1)=ZTOT_RIDRYG (JMICRO+JL-1)+ZRIDRYG   (JL)*ZMAXTIME(JL)
                ZTOT_RRDRYG (JMICRO+JL-1)=ZTOT_RRDRYG (JMICRO+JL-1)+ZRRDRYG   (JL)*ZMAXTIME(JL)
                ZTOT_RSDRYG (JMICRO+JL-1)=ZTOT_RSDRYG (JMICRO+JL-1)+ZRSDRYG   (JL)*ZMAXTIME(JL)
                ZTOT_RGMLTR (JMICRO+JL-1)=ZTOT_RGMLTR (JMICRO+JL-1)+ZRGMLTR   (JL)*ZMAXTIME(JL)
                ZTOT_RCWETH (JMICRO+JL-1)=ZTOT_RCWETH (JMICRO+JL-1)+ZRCWETH   (JL)*ZMAXTIME(JL)
                ZTOT_RIWETH (JMICRO+JL-1)=ZTOT_RIWETH (JMICRO+JL-1)+ZRIWETH   (JL)*ZMAXTIME(JL)
                ZTOT_RSWETH (JMICRO+JL-1)=ZTOT_RSWETH (JMICRO+JL-1)+ZRSWETH   (JL)*ZMAXTIME(JL)
                ZTOT_RGWETH (JMICRO+JL-1)=ZTOT_RGWETH (JMICRO+JL-1)+ZRGWETH   (JL)*ZMAXTIME(JL)
                ZTOT_RRWETH (JMICRO+JL-1)=ZTOT_RRWETH (JMICRO+JL-1)+ZRRWETH   (JL)*ZMAXTIME(JL)
                ZTOT_RCDRYH (JMICRO+JL-1)=ZTOT_RCDRYH (JMICRO+JL-1)+ZRCDRYH   (JL)*ZMAXTIME(JL)
                ZTOT_RIDRYH (JMICRO+JL-1)=ZTOT_RIDRYH (JMICRO+JL-1)+ZRIDRYH   (JL)*ZMAXTIME(JL)
                ZTOT_RSDRYH (JMICRO+JL-1)=ZTOT_RSDRYH (JMICRO+JL-1)+ZRSDRYH   (JL)*ZMAXTIME(JL)
                ZTOT_RRDRYH (JMICRO+JL-1)=ZTOT_RRDRYH (JMICRO+JL-1)+ZRRDRYH   (JL)*ZMAXTIME(JL)
                ZTOT_RGDRYH (JMICRO+JL-1)=ZTOT_RGDRYH (JMICRO+JL-1)+ZRGDRYH   (JL)*ZMAXTIME(JL)
                ZTOT_RDRYHG (JMICRO+JL-1)=ZTOT_RDRYHG (JMICRO+JL-1)+ZRDRYHG   (JL)*ZMAXTIME(JL)
                ZTOT_RHMLTR (JMICRO+JL-1)=ZTOT_RHMLTR (JMICRO+JL-1)+ZRHMLTR   (JL)*ZMAXTIME(JL)
                ZTOT_RIMLTC (JMICRO+JL-1)=ZTOT_RIMLTC (JMICRO+JL-1)+ZRIMLTC_MR(JL)
                ZTOT_RCBERI (JMICRO+JL-1)=ZTOT_RCBERI (JMICRO+JL-1)+ZRCBERI   (JL)*ZMAXTIME(JL)
    
              ENDDO
            ENDIF
            !
            !***       4.5 Next loop
            !
            LSOFT=.TRUE. ! We try to adjust tendencies (inner while loop)
          ENDDO
        ENDDO
    
    
          !Z..T variables contain the external tendency, we substract it
    
              ZVART(JL, JV) = ZVART(JL, JV) - ZEXTPK(JL, JV) * PTSTEP
    
            ENDDO
          ENDDO
        ENDIF
    
    !-------------------------------------------------------------------------------
    !
    !*       5.     UNPACKING DIAGNOSTICS
    !               ---------------------
    
          ZWR(I1(JL),I2(JL),IRV)=ZVART(JL, IRV)
          ZWR(I1(JL),I2(JL),IRC)=ZVART(JL, IRC)
          ZWR(I1(JL),I2(JL),IRR)=ZVART(JL, IRR)
          ZWR(I1(JL),I2(JL),IRI)=ZVART(JL, IRI)
          ZWR(I1(JL),I2(JL),IRS)=ZVART(JL, IRS)
          ZWR(I1(JL),I2(JL),IRG)=ZVART(JL, IRG)
    
    
    !==========================================================================================================
    
    
    
    !*       6.     COMPUTES THE SLOW COLD PROCESS SOURCES OUTSIDE OF ODMICRO POINTS
    
    !               ----------------------------------------------------------------
    !
    
        IF (.NOT. ODMICRO(JIJ, JK)) THEN
          LLW3D(JIJ, JK)=.TRUE.
          ZW3D(JIJ, JK)=ZZ_LSFACT(JIJ, JK)/PEXN(JIJ, JK)
        ELSE
          LLW3D(JIJ, JK)=.FALSE.
        ENDIF
    
    CALL ICE4_NUCLEATION(CST, PARAMI, ICEP, ICED, D%NIJT*D%NKT, LLW3D(:,:), &
                         PTHT(:, :), PPABST(:, :), PRHODREF(:, :), &                                       
                         PEXN(:, :), ZW3D(:, :), ZT(:, :), &                                                           
                         PRVT(:, :), &                                                                                 
                         PCIT(:, :), ZZ_RVHENI_MR(:, :))
    
    !
    !-------------------------------------------------------------------------------
    !
    !*       7.     TOTAL TENDENCIES
    !               ----------------
    !
    !
    !***     7.1    total tendencies limited by available species
    !
    
    DO JK = IKTB, IKTE
      DO CONCURRENT (JIJ=IIJB:IIJE)
    
        !LV/LS
        ZZ_LSFACT(JIJ,JK)=ZZ_LSFACT(JIJ,JK)/PEXNREF(JIJ,JK)
        ZZ_LVFACT(JIJ,JK)=ZZ_LVFACT(JIJ,JK)/PEXNREF(JIJ,JK)
    
        !Tendency dure to nucleation on non ODMICRO points
        ZZ_RVHENI(JIJ,JK) = MIN(PRVS(JIJ,JK), ZZ_RVHENI_MR(JIJ,JK)/PTSTEP)
    
        !Hydrometeor tendencies is the difference between old state and new state (can be negative)
        ZWR(JIJ,JK,IRV)=(ZWR(JIJ,JK,IRV)-PRVT(JIJ,JK))*ZINV_TSTEP
        ZWR(JIJ,JK,IRC)=(ZWR(JIJ,JK,IRC)-PRCT(JIJ,JK))*ZINV_TSTEP
        ZWR(JIJ,JK,IRR)=(ZWR(JIJ,JK,IRR)-PRRT(JIJ,JK))*ZINV_TSTEP
        ZWR(JIJ,JK,IRI)=(ZWR(JIJ,JK,IRI)-PRIT(JIJ,JK))*ZINV_TSTEP
        ZWR(JIJ,JK,IRS)=(ZWR(JIJ,JK,IRS)-PRST(JIJ,JK))*ZINV_TSTEP
        ZWR(JIJ,JK,IRG)=(ZWR(JIJ,JK,IRG)-PRGT(JIJ,JK))*ZINV_TSTEP
        IF(KRR==7) THEN
          ZWR(JIJ,JK,IRH)=(ZWR(JIJ,JK,IRH)-PRHT(JIJ,JK))*ZINV_TSTEP
        ENDIF
    
        !Theta tendency computed from hydrometeors tendencies
        ZWR(JIJ,JK, ITH) = (ZWR(JIJ,JK,IRC)+ZWR(JIJ,JK,IRR))*ZZ_LVFACT(JIJ,JK)+ &
                           & (ZWR(JIJ,JK,IRI)+ZWR(JIJ,JK,IRS)+ZWR(JIJ,JK,IRG)+ &
                           &  ZWR(JIJ,JK,IRH))*ZZ_LSFACT(JIJ,JK)
    
        !We apply these tendencies to the S variables
        !including the nucleation part
        PTHS(JIJ,JK) = PTHS(JIJ,JK) + ZWR(JIJ,JK,ITH)+ZZ_RVHENI(JIJ,JK)*ZZ_LSFACT(JIJ,JK)
        PRVS(JIJ,JK) = PRVS(JIJ,JK) + ZWR(JIJ,JK,IRV)-ZZ_RVHENI(JIJ,JK)
        PRCS(JIJ,JK) = PRCS(JIJ,JK) + ZWR(JIJ,JK,IRC)
        PRRS(JIJ,JK) = PRRS(JIJ,JK) + ZWR(JIJ,JK,IRR)
        PRIS(JIJ,JK) = PRIS(JIJ,JK) + ZWR(JIJ,JK,IRI)+ZZ_RVHENI(JIJ,JK)
        PRSS(JIJ,JK) = PRSS(JIJ,JK) + ZWR(JIJ,JK,IRS)
        PRGS(JIJ,JK) = PRGS(JIJ,JK) + ZWR(JIJ,JK,IRG)
        IF (KRR==7) THEN
          PRHS(JIJ,JK) = PRHS(JIJ,JK) + ZWR(JIJ,JK,IRH)
        ENDIF
    
    IF(BUCONF%LBU_ENABLE) THEN
      IF (BUCONF%LBUDGET_TH) THEN
    
            ZZ_DIFF(JIJ, JK) = ZZ_LSFACT(JIJ, JK) - ZZ_LVFACT(JIJ, JK)
    
        ZW(I1TOT(JL), I2TOT(JL)) = ZTOT_RVHENI(JL) * ZINV_TSTEP
    
      IF (BUCONF%LBUDGET_TH) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_TH), 'HENU',  ZW(:, :)*ZZ_LSFACT(:, :)*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RV) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RV), 'HENU', -ZW(:, :)                *PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RI) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RI), 'HENU',  ZW(:, :)                *PRHODJ(:, :))
    
    #else
      IF (BUCONF%LBUDGET_TH) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_TH), 'HIN',  ZW(:, :)*ZZ_LSFACT(:, :)*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RV) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RV), 'HIN', -ZW(:, :)                *PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RI) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RI), 'HIN',  ZW(:, :)                *PRHODJ(:, :))
    #endif
    
        ZW(I1TOT(JL), I2TOT(JL)) = ZTOT_RCHONI(JL) * ZINV_TSTEP
    
      IF (BUCONF%LBUDGET_TH) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_TH), 'HON',  ZW(:, :)*ZZ_DIFF(:, :)*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RC) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RC), 'HON', -ZW(:, :)              *PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RI) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RI), 'HON',  ZW(:, :)              *PRHODJ(:, :))
    
        ZW(I1TOT(JL), I2TOT(JL)) = ZTOT_RRHONG(JL) * ZINV_TSTEP
    
      IF (BUCONF%LBUDGET_TH) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_TH), 'SFR',  ZW(:, :)*ZZ_DIFF(:, :)*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RR) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RR), 'SFR', -ZW(:, :)              *PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RG) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RG), 'SFR',  ZW(:, :)              *PRHODJ(:, :))
    
        ZW(I1TOT(JL), I2TOT(JL)) = ZTOT_RVDEPS(JL) * ZINV_TSTEP
    
      IF (BUCONF%LBUDGET_TH) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_TH), 'DEPS',  ZW(:, :)*ZZ_LSFACT(:, :)*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RV) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RV), 'DEPS', -ZW(:, :)                *PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RS) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RS), 'DEPS',  ZW(:, :)                *PRHODJ(:, :))
    
        ZW(I1TOT(JL), I2TOT(JL)) = ZTOT_RIAGGS(JL) * ZINV_TSTEP
    
      IF (BUCONF%LBUDGET_RI) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RI), 'AGGS', -ZW(:, :)*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RS) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RS), 'AGGS',  ZW(:, :)*PRHODJ(:, :))
    
        ZW(I1TOT(JL), I2TOT(JL)) = ZTOT_RIAUTS(JL) * ZINV_TSTEP
    
      IF (BUCONF%LBUDGET_RI) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RI), 'AUTS', -ZW(:, :)*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RS) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RS), 'AUTS',  ZW(:, :)*PRHODJ(:, :))
    
        ZW(I1TOT(JL), I2TOT(JL)) = ZTOT_RVDEPG(JL) * ZINV_TSTEP
    
      IF (BUCONF%LBUDGET_TH) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_TH), 'DEPG',  ZW(:, :)*ZZ_LSFACT(:, :)*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RV) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RV), 'DEPG', -ZW(:, :)                *PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RG) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RG), 'DEPG',  ZW(:, :)                *PRHODJ(:, :))
    
          ZW(I1TOT(JL), I2TOT(JL)) = ZTOT_RCAUTR(JL) * ZINV_TSTEP
    
        IF (BUCONF%LBUDGET_RC) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RC), 'AUTO', -ZW(:, :)*PRHODJ(:, :))
        IF (BUCONF%LBUDGET_RR) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RR), 'AUTO',  ZW(:, :)*PRHODJ(:, :))
    
          ZW(I1TOT(JL), I2TOT(JL)) = ZTOT_RCACCR(JL) * ZINV_TSTEP
    
        IF (BUCONF%LBUDGET_RC) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RC), 'ACCR', -ZW(:, :)*PRHODJ(:, :))
        IF (BUCONF%LBUDGET_RR) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RR), 'ACCR',  ZW(:, :)*PRHODJ(:, :))
    
          ZW(I1TOT(JL), I2TOT(JL)) = ZTOT_RREVAV(JL) * ZINV_TSTEP
    
        IF (BUCONF%LBUDGET_TH) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_TH), 'REVA', -ZW(:, :)*ZZ_LVFACT(:, :)*PRHODJ(:, :))
        IF (BUCONF%LBUDGET_RV) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RV), 'REVA',  ZW(:, :)                *PRHODJ(:, :))
        IF (BUCONF%LBUDGET_RR) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RR), 'REVA', -ZW(:, :)                *PRHODJ(:, :))
    
        ZW1(I1TOT(JL), I2TOT(JL)) = ZTOT_RCRIMSS(JL) * ZINV_TSTEP
    
        ZW2(I1TOT(JL), I2TOT(JL)) = ZTOT_RCRIMSG(JL) * ZINV_TSTEP
    
        ZW3(I1TOT(JL), I2TOT(JL)) = ZTOT_RSRIMCG(JL) * ZINV_TSTEP
    
        CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_TH), 'RIM', (ZW1(:, :)+ZW2(:, :))*ZZ_DIFF(:, :)*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RC) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RC), 'RIM', (-ZW1(:, :)-ZW2(:, :))*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RS) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RS), 'RIM', ( ZW1(:, :)-ZW3(:, :))*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RG) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RG), 'RIM', ( ZW2(:, :)+ZW3(:, :))*PRHODJ(:, :))
    
        ZW1(I1TOT(JL), I2TOT(JL)) = ZTOT_RRACCSS(JL) * ZINV_TSTEP
    
        ZW2(I1TOT(JL), I2TOT(JL)) = ZTOT_RRACCSG(JL) * ZINV_TSTEP
    
        ZW3(I1TOT(JL), I2TOT(JL)) = ZTOT_RSACCRG(JL) * ZINV_TSTEP
    
        CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_TH), 'ACC', (ZW1(:, :)+ZW2(:, :) )*ZZ_DIFF(:, :)*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RR) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RR), 'ACC', (-ZW1(:, :)-ZW2(:, :))*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RS) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RS), 'ACC', ( ZW1(:, :)-ZW3(:, :))*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RG) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RG), 'ACC', ( ZW2(:, :)+ZW3(:, :))*PRHODJ(:, :))
    
        ZW(I1TOT(JL), I2TOT(JL)) = ZTOT_RSMLTG(JL) * ZINV_TSTEP
    
      IF (BUCONF%LBUDGET_RS) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RS), 'CMEL', -ZW(:, :)*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RG) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RG), 'CMEL',  ZW(:, :)*PRHODJ(:, :))
      ZW(:,:) = 0.
    
        ZW(I1TOT(JL), I2TOT(JL)) = ZTOT_RCMLTSR(JL) * ZINV_TSTEP
    
      IF (BUCONF%LBUDGET_RC) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RC), 'CMEL', -ZW(:, :)*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RR) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RR), 'CMEL',  ZW(:, :)*PRHODJ(:, :))
    
        ZW1(I1TOT(JL), I2TOT(JL)) = ZTOT_RICFRRG(JL) * ZINV_TSTEP
    
        ZW2(I1TOT(JL), I2TOT(JL)) = ZTOT_RRCFRIG(JL) * ZINV_TSTEP
    
        ZW3(I1TOT(JL), I2TOT(JL)) = ZTOT_RICFRR(JL) * ZINV_TSTEP
    
        CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_TH), 'CFRZ', ZW2(:, :)*ZZ_DIFF(:, :)*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RR) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RR), 'CFRZ', (-ZW2(:, :)+ZW3(:, :))*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RI) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RI), 'CFRZ', (-ZW1(:, :)-ZW3(:, :))*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RG) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RG), 'CFRZ', ( ZW1(:, :)+ZW2(:, :))*PRHODJ(:, :))
    
        ZW1(I1TOT(JL), I2TOT(JL)) = ZTOT_RCWETG(JL) * ZINV_TSTEP
    
        ZW2(I1TOT(JL), I2TOT(JL)) = ZTOT_RRWETG(JL) * ZINV_TSTEP
    
        ZW3(I1TOT(JL), I2TOT(JL)) = ZTOT_RIWETG(JL) * ZINV_TSTEP
    
        ZW4(I1TOT(JL), I2TOT(JL)) = ZTOT_RSWETG(JL) * ZINV_TSTEP
    
        CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_TH), 'WETG', (ZW1(:, :)+ZW2(:, :))*ZZ_DIFF(:, :)*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RC) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RC), 'WETG', -ZW1(:, :)    *PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RR) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RR), 'WETG', -ZW2(:, :)    *PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RI) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RI), 'WETG', -ZW3(:, :)    *PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RS) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RS), 'WETG', -ZW4(:, :)    *PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RG) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RG), 'WETG', (ZW1(:, :)+ZW2(:, :)+ZW3(:, :)+ZW4(:, :)) &
                                                                          &                             *PRHODJ(:, :))
    
          ZW(I1TOT(JL), I2TOT(JL)) = ZTOT_RWETGH(JL) * ZINV_TSTEP
    
        IF (BUCONF%LBUDGET_RG) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RG), 'GHCV', -ZW(:, :)*PRHODJ(:, :))
        IF (BUCONF%LBUDGET_RH) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RH), 'GHCV',  ZW(:, :)*PRHODJ(:, :))
    
        ZW1(I1TOT(JL), I2TOT(JL)) = ZTOT_RCDRYG(JL) * ZINV_TSTEP
    
        ZW2(I1TOT(JL), I2TOT(JL)) = ZTOT_RRDRYG(JL) * ZINV_TSTEP
    
        ZW3(I1TOT(JL), I2TOT(JL)) = ZTOT_RIDRYG(JL) * ZINV_TSTEP
    
        ZW4(I1TOT(JL), I2TOT(JL)) = ZTOT_RSDRYG(JL) * ZINV_TSTEP
    
        CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_TH), 'DRYG', (ZW1(:, :)+ZW2(:, :) )*ZZ_DIFF(:, :)*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RC) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RC), 'DRYG', -ZW1(:, :)     *PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RR) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RR), 'DRYG', -ZW2(:, :)     *PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RI) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RI), 'DRYG', -ZW3(:, :)     *PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RS) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RS), 'DRYG', -ZW4(:, :)     *PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RG) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RG), 'DRYG', (ZW1(:, :)+ZW2(:, :)+ZW3(:, :)+ZW4(:, :)) &
                                                                          &                              *PRHODJ(:, :))
    
      ZW(:,:) = 0.
    
        ZW(I1TOT(JL), I2TOT(JL)) = ZTOT_RGMLTR(JL) * ZINV_TSTEP
    
      IF (BUCONF%LBUDGET_TH) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_TH), 'GMLT', -ZW(:, :)*ZZ_DIFF(:, :)*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RR) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RR), 'GMLT',  ZW(:, :)              *PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RG) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RG), 'GMLT', -ZW(:, :)              *PRHODJ(:, :))
    
          ZW1(I1TOT(JL), I2TOT(JL)) = ZTOT_RCWETH(JL) * ZINV_TSTEP
    
          ZW2(I1TOT(JL), I2TOT(JL)) = ZTOT_RRWETH(JL) * ZINV_TSTEP
    
          ZW3(I1TOT(JL), I2TOT(JL)) = ZTOT_RIWETH(JL) * ZINV_TSTEP
    
          ZW4(I1TOT(JL), I2TOT(JL)) = ZTOT_RSWETH(JL) * ZINV_TSTEP
    
          ZW5(I1TOT(JL), I2TOT(JL)) = ZTOT_RGWETH(JL) * ZINV_TSTEP
    
          CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_TH), 'WETH', (ZW1(:, :)+ZW2(:, :))*ZZ_DIFF(:, :)*PRHODJ(:, :))
        IF (BUCONF%LBUDGET_RC) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RC), 'WETH', -ZW1(:, :)    *PRHODJ(:, :))
        IF (BUCONF%LBUDGET_RR) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RR), 'WETH', -ZW2(:, :)    *PRHODJ(:, :))
        IF (BUCONF%LBUDGET_RI) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RI), 'WETH', -ZW3(:, :)    *PRHODJ(:, :))
        IF (BUCONF%LBUDGET_RS) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RS), 'WETH', -ZW4(:, :)    *PRHODJ(:, :))
    
        IF (BUCONF%LBUDGET_RG) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RG), 'WETH', -ZW5(:, :)    *PRHODJ(:, :))
    
        IF (BUCONF%LBUDGET_RH) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RH), 'WETH', (ZW1(:, :)+ZW2(:, :)+ZW3(:, :)+ &
                                                                            &ZW4(:, :)+ZW5(:, : ))  *PRHODJ(:, :))
    
    #if defined(REPRO48) || defined(REPRO55)
    
          ZW(I1TOT(JL), I2TOT(JL)) = ZTOT_RGWETH(JL) * ZINV_TSTEP
    
        IF (BUCONF%LBUDGET_RG) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RG), 'HGCV', (-ZW5(:, :)-ZW(:, :))*PRHODJ(:, :))
    
        IF (BUCONF%LBUDGET_RG) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RG), 'HGCV', -ZW(:, :)*PRHODJ(:, :))
    
    #endif
    #if defined(REPRO48) || defined(REPRO55)
    
        IF (BUCONF%LBUDGET_RH) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RH), 'HGCV',  ZW(:, :)*PRHODJ(:, :))
    
          ZW1(I1TOT(JL), I2TOT(JL)) = ZTOT_RCDRYH(JL) * ZINV_TSTEP
    
          ZW2(I1TOT(JL), I2TOT(JL)) = ZTOT_RRDRYH(JL) * ZINV_TSTEP
    
          ZW3(I1TOT(JL), I2TOT(JL)) = ZTOT_RIDRYH(JL) * ZINV_TSTEP
    
          ZW4(I1TOT(JL), I2TOT(JL)) = ZTOT_RSDRYH(JL) * ZINV_TSTEP
    
          ZW5(I1TOT(JL), I2TOT(JL)) = ZTOT_RGDRYH(JL) * ZINV_TSTEP
    
    #if defined(REPRO48) || defined(REPRO55)
    
        !ZW6 must be removed when REPRO* will be suppressed
    
          ZW6(I1TOT(JL), I2TOT(JL)) = ZTOT_RDRYHG(JL) * ZINV_TSTEP
    
          CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_TH), 'DRYH', (ZW1(:, :)+ZW2(:, :))*ZZ_DIFF(:, :)*PRHODJ(:, :))
        IF (BUCONF%LBUDGET_RC) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RC), 'DRYH', -ZW1(:, :)                 *PRHODJ(:, :))
        IF (BUCONF%LBUDGET_RR) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RR), 'DRYH', -ZW2(:, :)                 *PRHODJ(:, :))
        IF (BUCONF%LBUDGET_RI) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RI), 'DRYH', -ZW3(:, :)                 *PRHODJ(:, :))
        IF (BUCONF%LBUDGET_RS) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RS), 'DRYH', -ZW4(:, :)                 *PRHODJ(:, :))
        IF (BUCONF%LBUDGET_RG) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RG), 'DRYH', (-ZW5(:, :)+ZW6(:, :))     *PRHODJ(:, :))
        IF (BUCONF%LBUDGET_RH) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RH), 'DRYH', (ZW1(:, :)+ZW2(:, :)+ZW3(:, :)+   &
                                                                            &ZW4(:, :)+ZW5(:, :)-ZW6(:, :)) &
                                                                            &                             *PRHODJ(:, :))
    
    #if defined(REPRO48) || defined(REPRO55)
    #else
        !When REPRO48 will be suppressed, ZW6 must be removed
    
          ZW(I1TOT(JL), I2TOT(JL)) = ZTOT_RDRYHG(JL) * ZINV_TSTEP
    
        IF (BUCONF%LBUDGET_RG) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RG), 'HGCV', -ZW(:, :)*PRHODJ(:, :))
        IF (BUCONF%LBUDGET_RH) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RH), 'HGCV',  ZW(:, :)*PRHODJ(:, :))
    
          ZW(I1TOT(JL), I2TOT(JL)) = ZTOT_RHMLTR(JL) * ZINV_TSTEP
    
        IF (BUCONF%LBUDGET_TH) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_TH), 'HMLT', -ZW(:, :)*ZZ_DIFF(:, :)*PRHODJ(:, :))
        IF (BUCONF%LBUDGET_RR) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RR), 'HMLT',  ZW(:, :)              *PRHODJ(:, :))
        IF (BUCONF%LBUDGET_RH) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RH), 'HMLT', -ZW(:, :)              *PRHODJ(:, :))
    
        ZW(I1TOT(JL), I2TOT(JL)) = ZTOT_RIMLTC(JL) * ZINV_TSTEP
    
      IF (BUCONF%LBUDGET_TH) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_TH), 'IMLT', -ZW(:, :)*ZZ_DIFF(:, :)*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RC) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RC), 'IMLT',  ZW(:, :)              *PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RI) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RI), 'IMLT', -ZW(:, :)              *PRHODJ(:, :))
    
        ZW(I1TOT(JL), I2TOT(JL)) = ZTOT_RCBERI(JL) * ZINV_TSTEP
    
      IF (BUCONF%LBUDGET_TH) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_TH), 'BERFI',  ZW(:, :)*ZZ_DIFF(:, :)*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RC) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RC), 'BERFI', -ZW(:, :)              *PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RI) CALL BUDGET_STORE_ADD_PHY(D, TBUDGETS(NBUDGET_RI), 'BERFI',  ZW(:, :)              *PRHODJ(:, :))
    
      IF (BUCONF%LBUDGET_TH) CALL BUDGET_STORE_INIT_PHY(D, TBUDGETS(NBUDGET_TH), 'CORR', PTHS(:, :)*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RV) CALL BUDGET_STORE_INIT_PHY(D, TBUDGETS(NBUDGET_RV), 'CORR', PRVS(:, :)*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RC) CALL BUDGET_STORE_INIT_PHY(D, TBUDGETS(NBUDGET_RC), 'CORR', PRCS(:, :)*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RR) CALL BUDGET_STORE_INIT_PHY(D, TBUDGETS(NBUDGET_RR), 'CORR', PRSS(:, :)*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RI) CALL BUDGET_STORE_INIT_PHY(D, TBUDGETS(NBUDGET_RI), 'CORR', PRIS(:, :)*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RS) CALL BUDGET_STORE_INIT_PHY(D, TBUDGETS(NBUDGET_RS), 'CORR', PRSS(:, :)*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RG) CALL BUDGET_STORE_INIT_PHY(D, TBUDGETS(NBUDGET_RG), 'CORR', PRGS(:, :)*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RH .AND. KRR==7) CALL BUDGET_STORE_INIT_PHY(D, TBUDGETS(NBUDGET_RH), 'CORR', PRHS(:, :)*PRHODJ(:, :))
    
    !NOTE:
    !  This call cannot be moved before the preeceding budget calls because,
    !  with AROME, the BUDGET_STORE_INIT does nothing. The equivalent is done only
    !  once before the physics call and copies of the S variables evolve automatically
    !  internally to the budget (DDH) machinery at each BUDGET_STORE_ADD and
    !  BUDGET_STORE_END calls. Thus, the difference between the DDH internal version
    !  of the S variables and the S variables used in the folowing BUDGET_STORE_END
    !  call must only be due to the correction of negativities.
    !
    !We correct negativities with conservation
    
    CALL CORRECT_NEGATIVITIES(D, KRR, PRVS, PRCS, PRRS, &
    
                             &PRIS, PRSS, PRGS, &
                             &PTHS, ZZ_LVFACT, ZZ_LSFACT, PRHS)
    
    
      IF (BUCONF%LBUDGET_TH) CALL BUDGET_STORE_END_PHY(D, TBUDGETS(NBUDGET_TH), 'CORR', PTHS(:, :)*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RV) CALL BUDGET_STORE_END_PHY(D, TBUDGETS(NBUDGET_RV), 'CORR', PRVS(:, :)*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RC) CALL BUDGET_STORE_END_PHY(D, TBUDGETS(NBUDGET_RC), 'CORR', PRCS(:, :)*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RR) CALL BUDGET_STORE_END_PHY(D, TBUDGETS(NBUDGET_RR), 'CORR', PRRS(:, :)*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RI) CALL BUDGET_STORE_END_PHY(D, TBUDGETS(NBUDGET_RI), 'CORR', PRIS(:, :)*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RS) CALL BUDGET_STORE_END_PHY(D, TBUDGETS(NBUDGET_RS), 'CORR', PRSS(:, :)*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RG) CALL BUDGET_STORE_END_PHY(D, TBUDGETS(NBUDGET_RG), 'CORR', PRGS(:, :)*PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RH .AND. KRR==7) CALL BUDGET_STORE_END_PHY(D, TBUDGETS(NBUDGET_RH), 'CORR', PRHS(:, :)*PRHODJ(:, :))
    
    !
    !-------------------------------------------------------------------------------
    !
    !*       8.     COMPUTE THE SEDIMENTATION (RS) SOURCE
    !               -------------------------------------
    !
    
      IF (BUCONF%LBUDGET_RC .AND. PARAMI%LSEDIC) CALL BUDGET_STORE_INIT_PHY(D, TBUDGETS(NBUDGET_RC), 'SEDI', PRCS(:, :) * PRHODJ(:, :))
    
      IF (BUCONF%LBUDGET_RR)              CALL BUDGET_STORE_INIT_PHY(D, TBUDGETS(NBUDGET_RR), 'SEDI', PRRS(:, :) * PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RI)              CALL BUDGET_STORE_INIT_PHY(D, TBUDGETS(NBUDGET_RI), 'SEDI', PRIS(:, :) * PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RS)              CALL BUDGET_STORE_INIT_PHY(D, TBUDGETS(NBUDGET_RS), 'SEDI', PRSS(:, :) * PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RG)              CALL BUDGET_STORE_INIT_PHY(D, TBUDGETS(NBUDGET_RG), 'SEDI', PRGS(:, :) * PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RH .AND. KRR==7) CALL BUDGET_STORE_INIT_PHY(D, TBUDGETS(NBUDGET_RH), 'SEDI', PRHS(:, :) * PRHODJ(:, :))
    
              ZRCT(JIJ,JK)=PRCS(JIJ,JK)*PTSTEP
              ZRRT(JIJ,JK)=PRRS(JIJ,JK)*PTSTEP
              ZRIT(JIJ,JK)=PRIS(JIJ,JK)*PTSTEP
              ZRST(JIJ,JK)=PRSS(JIJ,JK)*PTSTEP
              ZRGT(JIJ,JK)=PRGS(JIJ,JK)*PTSTEP
              ZRHT(JIJ,JK)=PRHS(JIJ,JK)*PTSTEP
    
          CALL ICE4_SEDIMENTATION_STAT(D, CST, ICEP, ICED, &
    
                                      &PTSTEP, KRR, OSEDIC, PDZZ, &
    
                                      &PRCS, ZRCT, PRRS, ZRRT, PRIS, ZRIT,&
                                      &PRSS, ZRST, PRGS, ZRGT,&
                                      &PINPRC, PINPRR, ZINPRI, PINPRS, PINPRG, &
                                      &PSEA=PSEA, PTOWN=PTOWN, &
                                      &PINPRH=PINPRH, PRHT=ZRHT, PRHS=PRHS, PFPR=PFPR)
        ELSE
    
              ZRCT(JIJ,JK)=PRCS(JIJ,JK)*PTSTEP
              ZRRT(JIJ,JK)=PRRS(JIJ,JK)*PTSTEP
              ZRIT(JIJ,JK)=PRIS(JIJ,JK)*PTSTEP
              ZRST(JIJ,JK)=PRSS(JIJ,JK)*PTSTEP
              ZRGT(JIJ,JK)=PRGS(JIJ,JK)*PTSTEP
    
          CALL ICE4_SEDIMENTATION_STAT(D, CST, ICEP, ICED, &
    
                                      &PTSTEP, KRR, OSEDIC, PDZZ, &
    
                                      &PRCS, ZRCT, PRRS, ZRRT, PRIS, ZRIT,&
                                      &PRSS, ZRST, PRGS, ZRGT,&
                                      &PINPRC, PINPRR, ZINPRI, PINPRS, PINPRG, &
                                      &PSEA=PSEA, PTOWN=PTOWN, &
                                      &PFPR=PFPR)
        ENDIF
    
        PINPRS(IIJB:IIJE) = PINPRS(IIJB:IIJE) + ZINPRI(IIJB:IIJE)
    
        !No negativity correction here as we apply sedimentation on PR.S*PTSTEP variables
    
        !SR: It *seems* that we must have two separate calls for ifort
        IF(KRR==7) THEN
    
          CALL ICE4_SEDIMENTATION_SPLIT(D, CST, ICEP, ICED, PARAMI, &
    
                                       &PRHODREF, PPABST, PTHT, ZT, PRHODJ, &
    
                                       &PRCS, PRCT, PRRS, PRRT, PRIS, PRIT, PRSS, PRST, PRGS, PRGT,&
                                       &PINPRC, PINPRR, ZINPRI, PINPRS, PINPRG, &
                                       &PSEA=PSEA, PTOWN=PTOWN, &
                                       &PINPRH=PINPRH, PRHT=PRHT, PRHS=PRHS, PFPR=PFPR)
        ELSE
    
          CALL ICE4_SEDIMENTATION_SPLIT(D, CST, ICEP, ICED, PARAMI, &
    
                                       &PRHODREF, PPABST, PTHT, ZT, PRHODJ, &
    
                                       &PRCS, PRCT, PRRS, PRRT, PRIS, PRIT, PRSS, PRST, PRGS, PRGT,&
                                       &PINPRC, PINPRR, ZINPRI, PINPRS, PINPRG, &
                                       &PSEA=PSEA, PTOWN=PTOWN, &
                                       &PFPR=PFPR)
        ENDIF
    
        PINPRS(IIJB:IIJE) = PINPRS(IIJB:IIJE) + ZINPRI(IIJB:IIJE)
    
        !We correct negativities with conservation
        !SPLI algorith uses a time-splitting. Inside the loop a temporary m.r. is used.
        !   It is initialized with the m.r. at T and is modified by two tendencies:
        !   sedimentation tendency and an external tendency which represents all other
        !   processes (mainly advection and microphysical processes). If both tendencies
    
        !   are negative, sedimentation can remove a species at a given sub-timestep. From
    
        !   this point sedimentation stops for the remaining sub-timesteps but the other tendency
        !   will be still active and will lead to negative values.
    
        !   We could prevent the algorithm to not consume too much a species, instead we apply
    
        CALL CORRECT_NEGATIVITIES(D, KRR, PRVS, PRCS, PRRS, &
    
                                 &PRIS, PRSS, PRGS, &
                                 &PTHS, ZZ_LVFACT, ZZ_LSFACT, PRHS)
      ELSE
    
        CALL PRINT_MSG(NVERB_FATAL, 'GEN', 'RAIN_ICE', 'no sedimentation scheme for PARAMI%CSEDIM='//PARAMI%CSEDIM)
    
      IF (BUCONF%LBUDGET_RC .AND. PARAMI%LSEDIC)   CALL BUDGET_STORE_END_PHY(D, TBUDGETS(NBUDGET_RC), 'SEDI', PRCS(:, :) * PRHODJ(:, :))
    
      IF (BUCONF%LBUDGET_RR)                CALL BUDGET_STORE_END_PHY(D, TBUDGETS(NBUDGET_RR), 'SEDI', PRRS(:, :) * PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RI)                CALL BUDGET_STORE_END_PHY(D, TBUDGETS(NBUDGET_RI), 'SEDI', PRIS(:, :) * PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RS)                CALL BUDGET_STORE_END_PHY(D, TBUDGETS(NBUDGET_RS), 'SEDI', PRSS(:, :) * PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RG)                CALL BUDGET_STORE_END_PHY(D, TBUDGETS(NBUDGET_RG), 'SEDI', PRGS(:, :) * PRHODJ(:, :))
      IF (BUCONF%LBUDGET_RH .AND. KRR==7)   CALL BUDGET_STORE_END_PHY(D, TBUDGETS(NBUDGET_RH), 'SEDI', PRHS(:, :) * PRHODJ(:, :))
    
    
      !"sedimentation" of rain fraction
      IF (PRESENT(PRHS)) THEN
    
        CALL ICE4_RAINFR_VERT(D, ICED, PRAINFR, PRRS(:,:)*PTSTEP, &
                           &PRSS(:,:)*PTSTEP, PRGS(:,:)*PTSTEP, PRHS(:,:)*PTSTEP)
    
        CALL ICE4_RAINFR_VERT(D, ICED, PRAINFR, PRRS(:,:)*PTSTEP, &
                           &PRSS(:,:)*PTSTEP, PRGS(:,:)*PTSTEP)
    
      ENDIF
    ENDIF
    !
    !-------------------------------------------------------------------------------
    !
    !*       9.     COMPUTE THE FOG DEPOSITION TERM 
    !               -------------------------------------
    !
    
    IF (PARAMI%LDEPOSC) THEN !cloud water deposition on vegetation
      IF (BUCONF%LBU_ENABLE .AND. BUCONF%LBUDGET_RC) &
    
         & CALL BUDGET_STORE_INIT_PHY(D, TBUDGETS(NBUDGET_RC), 'DEPO', PRCS(:, :)*PRHODJ(:, :))
    
      DO JIJ = IIJB, IIJE
        PINDEP(JIJ) = PARAMI%XVDEPOSC * PRCT(JIJ, IKB) * PRHODREF(JIJ, IKB) / CST%XRHOLW
        PRCS(JIJ, IKB) = PRCS(JIJ, IKB) - PARAMI%XVDEPOSC * PRCT(JIJ, IKB) / PDZZ(JIJ, IKB)
    
      IF (BUCONF%LBU_ENABLE .AND. BUCONF%LBUDGET_RC) &
    
         & CALL BUDGET_STORE_END_PHY(D, TBUDGETS(NBUDGET_RC), 'DEPO', PRCS(:, :)*PRHODJ(:, :))
    
    ENDIF
    
    IF (LHOOK) CALL DR_HOOK('RAIN_ICE', 1, ZHOOK_HANDLE)
    !
    CONTAINS
      !
    
      SUBROUTINE CORRECT_NEGATIVITIES(D, KRR, PRV, PRC, PRR, &
    
                                     &PRI, PRS, PRG, &
                                     &PTH, PLVFACT, PLSFACT, PRH)
      !
      IMPLICIT NONE
      !
    
      TYPE(DIMPHYEX_t),                     INTENT(IN)    :: D
      INTEGER,                              INTENT(IN)    :: KRR
    
      REAL, DIMENSION(D%NIJT, D%NKT), INTENT(INOUT) :: PRV, PRC, PRR, PRI, PRS, PRG, PTH
      REAL, DIMENSION(D%NIJT, D%NKT), INTENT(IN)    :: PLVFACT, PLSFACT
      REAL, DIMENSION(D%NIJT, D%NKT), OPTIONAL, INTENT(INOUT) :: PRH
    
      REAL(KIND=JPRB) :: ZHOOK_HANDLE
      !
      IF (LHOOK) CALL DR_HOOK('RAIN_ICE:CORRECT_NEGATIVITIES', 0, ZHOOK_HANDLE)
      !
      !We correct negativities with conservation
    
          ! 1) deal with negative values for mixing ratio, except for vapor
          ZW         =PRC(JIJ,JK)-MAX(PRC(JIJ,JK), 0.)
          PRV(JIJ,JK)=PRV(JIJ,JK)+ZW
          PTH(JIJ,JK)=PTH(JIJ,JK)-ZW*PLVFACT(JIJ,JK)
          PRC(JIJ,JK)=PRC(JIJ,JK)-ZW
    
          ZW         =PRR(JIJ,JK)-MAX(PRR(JIJ,JK), 0.)
          PRV(JIJ,JK)=PRV(JIJ,JK)+ZW
          PTH(JIJ,JK)=PTH(JIJ,JK)-ZW*PLVFACT(JIJ,JK)
          PRR(JIJ,JK)=PRR(JIJ,JK)-ZW
    
          ZW         =PRI(JIJ,JK)-MAX(PRI(JIJ,JK), 0.)
          PRV(JIJ,JK)=PRV(JIJ,JK)+ZW
          PTH(JIJ,JK)=PTH(JIJ,JK)-ZW*PLSFACT(JIJ,JK)
          PRI(JIJ,JK)=PRI(JIJ,JK)-ZW
    
          ZW         =PRS(JIJ,JK)-MAX(PRS(JIJ,JK), 0.)
          PRV(JIJ,JK)=PRV(JIJ,JK)+ZW
          PTH(JIJ,JK)=PTH(JIJ,JK)-ZW*PLSFACT(JIJ,JK)
          PRS(JIJ,JK)=PRS(JIJ,JK)-ZW
    
          ZW         =PRG(JIJ,JK)-MAX(PRG(JIJ,JK), 0.)
          PRV(JIJ,JK)=PRV(JIJ,JK)+ZW
          PTH(JIJ,JK)=PTH(JIJ,JK)-ZW*PLSFACT(JIJ,JK)
          PRG(JIJ,JK)=PRG(JIJ,JK)-ZW
    
          IF(KRR==7) THEN
            ZW         =PRH(JIJ,JK)-MAX(PRH(JIJ,JK), 0.)
            PRV(JIJ,JK)=PRV(JIJ,JK)+ZW
            PTH(JIJ,JK)=PTH(JIJ,JK)-ZW*PLSFACT(JIJ,JK)
            PRH(JIJ,JK)=PRH(JIJ,JK)-ZW
          ENDIF
    
          ! 2) deal with negative vapor mixing ratio
    
          ! for rc and ri, we keep ice fraction constant
          ZW=MIN(1., MAX(ICED%XRTMIN(1)-PRV(JIJ,JK), 0.) / &
                    &MAX(PRC(JIJ,JK)+PRI(JIJ,JK), 1.E-20)) ! Proportion of rc+ri to convert into rv
          PTH(JIJ,JK)=PTH(JIJ,JK)-ZW* &
                       &(PRC(JIJ,JK)*PLVFACT(JIJ,JK)+PRI(JIJ,JK)*PLSFACT(JIJ,JK))
          PRV(JIJ,JK)=PRV(JIJ,JK)+ZW*(PRC(JIJ,JK)+PRI(JIJ,JK))
          PRC(JIJ,JK)=(1.-ZW)*PRC(JIJ,JK)
          PRI(JIJ,JK)=(1.-ZW)*PRI(JIJ,JK)
    
          ZW=MIN(MAX(PRR(JIJ,JK), 0.), &
                &MAX(ICED%XRTMIN(1)-PRV(JIJ,JK), 0.)) ! Quantity of rr to convert into rv
          PRV(JIJ,JK)=PRV(JIJ,JK)+ZW
          PRR(JIJ,JK)=PRR(JIJ,JK)-ZW
          PTH(JIJ,JK)=PTH(JIJ,JK)-ZW*PLVFACT(JIJ,JK)
    
          ZW=MIN(MAX(PRS(JIJ,JK), 0.), &
                &MAX(ICED%XRTMIN(1)-PRV(JIJ,JK), 0.)) ! Quantity of rs to convert into rv
          PRV(JIJ,JK)=PRV(JIJ,JK)+ZW
          PRS(JIJ,JK)=PRS(JIJ,JK)-ZW
          PTH(JIJ,JK)=PTH(JIJ,JK)-ZW*PLSFACT(JIJ,JK)
    
          ZW=MIN(MAX(PRG(JIJ,JK), 0.), &
                &MAX(ICED%XRTMIN(1)-PRV(JIJ,JK), 0.)) ! Quantity of rg to convert into rv
          PRV(JIJ,JK)=PRV(JIJ,JK)+ZW
          PRG(JIJ,JK)=PRG(JIJ,JK)-ZW
          PTH(JIJ,JK)=PTH(JIJ,JK)-ZW*PLSFACT(JIJ,JK)
    
          IF(KRR==7) THEN
            ZW=MIN(MAX(PRH(JIJ,JK), 0.), &
                  &MAX(ICED%XRTMIN(1)-PRV(JIJ,JK), 0.)) ! Quantity of rh to convert into rv
            PRV(JIJ,JK)=PRV(JIJ,JK)+ZW
            PRH(JIJ,JK)=PRH(JIJ,JK)-ZW
            PTH(JIJ,JK)=PTH(JIJ,JK)-ZW*PLSFACT(JIJ,JK)
          ENDIF
    
        ENDDO
      ENDDO
      !
      IF (LHOOK) CALL DR_HOOK('RAIN_ICE:CORRECT_NEGATIVITIES', 1, ZHOOK_HANDLE)
      !
      END SUBROUTINE CORRECT_NEGATIVITIES