diff --git a/docs/TODO b/docs/TODO
index b549518c07d4afd3e1a029d937a3a62033a80423..e849c953adf25cfa78a33393e65644a9e60cc38d 100644
--- a/docs/TODO
+++ b/docs/TODO
@@ -22,6 +22,9 @@ Merge pb:
        Ryad a fait des tests pour regarder impact des allocatable sur CPU => temps * 2
        Code à nettoyer quelque soit l'option retenue
        Dernier code de Ryad: /home/gmap/mrpm/khatib/public/modset/mods_ice4_nucleation_wrapper.tgz et/ou /home/gmap/mrpm/khatib/public/modset/ice4_nucleation_wrapper.f90
+- shallow_mf:
+       Dans Méso-NH: shallow_mf doit être appelé avec PDX=XDXHAT(1) et PDY=XDYHAT(1)
+       Dans AROME: où trouver la taille de maille?
 
 - rain_ice_red: le cas test MesoNH n'est pas bit repro (diffs > 1% sur rapports de melange)
                 sur la modif src/mesonh/rain_ice_red au commit bdd10dd (First rain_ice new/red merge)
@@ -45,6 +48,7 @@ Pb identifiés à corriger plus tard:
        => à corriger
 - seules les options oper ont été testées, il manque des test pour sedim_after, nmaxiter, xmrstep, xtstep, autoconv, rainfr
 - arome/ini_cmfshall devrait s'appeler ini_param_mfshall
+- th_r_from_thl_rt appelée partout, il faudrait limiter à OTEST
 
 Répertoire arome/ext et mesonh/ext contient les codes non PHYEX qu'il faut modifier dans le pack pour qu'il puisse être compilé.
 Ce répertoire devra être vidé à la fin du phasage, les modifications nécessaires ayadevront avoir été fournies par ailleurs
diff --git a/src/arome/ext/aro_shallow_mf.F90 b/src/arome/ext/aro_shallow_mf.F90
new file mode 100644
index 0000000000000000000000000000000000000000..2a6a3bf6548bc35254837c64c934037147d2b3c9
--- /dev/null
+++ b/src/arome/ext/aro_shallow_mf.F90
@@ -0,0 +1,250 @@
+!     ######spl
+      SUBROUTINE  ARO_SHALLOW_MF(KKL, KLON,KLEV, KRR, KRRL, KRRI,KSV,     &
+                HMF_UPDRAFT, HMF_CLOUD, HFRAC_ICE, OMIXUV,            &
+                ONOMIXLG,KSV_LGBEG,KSV_LGEND,                         &
+                KTCOUNT, PTSTEP,                                      &
+                PZZ, PZZF, PDZZF,                                            &
+                PRHODJ, PRHODREF,                                     &
+                PPABSM, PEXNM,                                        &
+                PSFTH,PSFRV,                                          &
+                PTHM,PRM,                                             &
+                PUM,PVM,PTKEM,PSVM,                                   &
+                PDUDT_MF,PDVDT_MF,                                    &
+                PDTHLDT_MF,PDRTDT_MF,PDSVDT_MF,                       &
+                PSIGMF,PRC_MF,PRI_MF,PCF_MF,PFLXZTHVMF,                      &
+                PTHL_UP,PRT_UP,PRV_UP,PRC_UP,PRI_UP,                  &
+                PU_UP, PV_UP, PTHV_UP, PW_UP, PFRAC_UP, PEMF)
+
+      USE PARKIND1, ONLY : JPRB
+      USE YOMHOOK , ONLY : LHOOK, DR_HOOK
+!     ##########################################################################
+!
+!!****  * -  interface to call SHALLOW_MF :
+!!             computation of turbulence "mass flux" fluxes and their divergence
+!!
+!!
+!!
+!!    PURPOSE
+!!    -------
+!!
+!!
+!!
+!!
+!!**  METHOD
+!!    ------
+!!
+!!
+!!
+!!    EXTERNAL
+!!    --------
+!!      Subroutine SHALLOW_MF (routine de MesoNH)
+!!
+!!    IMPLICIT ARGUMENTS
+!!    ------------------
+!!
+!!    REFERENCE
+!!    ---------
+!!
+!!      Documentation AROME
+!!
+!!    AUTHOR
+!!    ------
+!!    S.Malardel
+!!
+!!    MODIFICATIONS
+!!    -------------
+!!      Original    07/2006
+!!      Y. Seity : new arguments for EDMF scheme 04/2009
+!!      S. Riette 18 May 2010: aro_shallow_mf and shallow_mf interfaces changed
+!!      S. Riette Jan 2012: support for both order of vertical levels
+!!
+!-------------------------------------------------------------------------------
+!
+!*       0.    DECLARATIONS
+!              ------------
+!
+USE MODD_PARAMETERS, ONLY: JPVEXT, JPHEXT
+!
+USE MODI_SHALLOW_MF
+!
+IMPLICIT NONE
+!
+!*       0.1   Declarations of dummy arguments :
+!
+!
+!
+INTEGER,                  INTENT(IN)   :: KKL      ! +1 if grid goes from ground to
+                                                   ! atmosphere top, -1 otherwise
+INTEGER,                  INTENT(IN)   :: KLON     !NPROMA under CPG
+INTEGER,                  INTENT(IN)   :: KLEV     !Number of vertical levels
+INTEGER,                  INTENT(IN)   :: KRR      ! Number of moist variables
+INTEGER,                  INTENT(IN)   :: KRRL     ! Number of liquide water variables
+INTEGER,                  INTENT(IN)   :: KRRI     ! Number of ice variables
+INTEGER,                  INTENT(IN)   :: KSV      ! Number of passive scalar variables
+!
+CHARACTER (LEN=4), INTENT(IN)   :: HMF_UPDRAFT  ! Type of Mass Flux Scheme
+CHARACTER (LEN=4), INTENT(IN)   :: HMF_CLOUD    ! Type of statistical cloud scheme
+CHARACTER*1,       INTENT(IN)   :: HFRAC_ICE    ! partition liquid/ice scheme
+LOGICAL,                        INTENT(IN)   :: OMIXUV    ! True if mixing of momentum
+!
+LOGICAL,                INTENT(IN)   :: ONOMIXLG  ! False if mixing of lagrangian tracer
+INTEGER,                INTENT(IN)   :: KSV_LGBEG ! first index of lag. tracer
+INTEGER,                INTENT(IN)   :: KSV_LGEND ! last  index of lag. tracer
+
+INTEGER,                  INTENT(IN)   :: KTCOUNT  ! Temporal loop counter
+REAL,                     INTENT(IN)   :: PTSTEP   ! Time step
+!
+REAL, DIMENSION(KLON,KLEV),   INTENT(IN)   :: PZZ     ! Height of layer boundaries
+REAL, DIMENSION(KLON,KLEV),   INTENT(IN)   :: PZZF    ! Height of level
+REAL, DIMENSION(KLON,KLEV),   INTENT(IN)   :: PDZZF    !thikness between layers
+
+REAL, DIMENSION(KLON,KLEV),   INTENT(IN)   :: PRHODJ  ! Dry density * Jacobian
+REAL, DIMENSION(KLON,KLEV),   INTENT(IN)   :: PRHODREF  ! Dry density
+!
+REAL, DIMENSION(KLON,KLEV),   INTENT(IN)   :: PPABSM      ! Pressure at time t-1
+REAL, DIMENSION(KLON,KLEV),   INTENT(IN)   :: PEXNM   ! Exner function
+!
+! normal surface fluxes of theta and Rv
+REAL, DIMENSION(KLON),          INTENT(IN)   ::  PSFTH,PSFRV
+!   prognostic variables at t- deltat
+!
+!   thermodynamical variables which are transformed in conservative var.
+REAL, DIMENSION(KLON,KLEV),   INTENT(IN)   ::  PTHM       ! pot. temp.
+REAL, DIMENSION(KLON,KLEV,KRR), INTENT(IN) ::  PRM         ! mixing ratio
+REAL, DIMENSION(KLON,KLEV),   INTENT(IN)   ::  PUM,PVM       ! momentum
+REAL, DIMENSION(KLON,KLEV),   INTENT(IN)   ::  PTKEM
+REAL, DIMENSION(KLON,KLEV,KSV), INTENT(IN) ::  PSVM         ! passive scalar
+                                                             ! variables for EDMF scheme
+REAL, DIMENSION(KLON,KLEV),   INTENT(OUT)::  PDUDT_MF     ! tendency of U   by massflux scheme
+REAL, DIMENSION(KLON,KLEV),   INTENT(OUT)::  PDVDT_MF     ! tendency of V   by massflux scheme
+REAL, DIMENSION(KLON,KLEV),   INTENT(OUT)::  PDTHLDT_MF   ! tendency of thl by massflux scheme
+REAL, DIMENSION(KLON,KLEV),   INTENT(OUT)::  PDRTDT_MF    ! tendency of rt  by massflux scheme
+REAL, DIMENSION(KLON,KLEV,KSV), INTENT(OUT)::  PDSVDT_MF    ! tendency of Sv  by massflux scheme
+
+REAL, DIMENSION(KLON,KLEV), INTENT(OUT)   ::  PSIGMF,PRC_MF,PRI_MF,PCF_MF ! cloud info for the cloud scheme
+REAL, DIMENSION(KLON,KLEV), INTENT(OUT)   ::  PFLXZTHVMF           ! Thermal production for TKE scheme
+REAL, DIMENSION(KLON,KLEV), INTENT(INOUT) ::  PTHL_UP   ! Thl updraft characteristics
+REAL, DIMENSION(KLON,KLEV), INTENT(INOUT) ::  PRT_UP    ! Rt  updraft characteristics
+REAL, DIMENSION(KLON,KLEV), INTENT(INOUT) ::  PRV_UP    ! Vapor updraft characteristics
+REAL, DIMENSION(KLON,KLEV), INTENT(INOUT) ::  PU_UP     ! U wind updraft characteristics
+REAL, DIMENSION(KLON,KLEV), INTENT(INOUT) ::  PV_UP     ! V wind updraft characteristics
+REAL, DIMENSION(KLON,KLEV), INTENT(INOUT) ::  PRC_UP    ! cloud content updraft characteristics
+REAL, DIMENSION(KLON,KLEV), INTENT(INOUT) ::  PRI_UP    ! ice content   updraft characteristics
+REAL, DIMENSION(KLON,KLEV), INTENT(INOUT) ::  PTHV_UP   ! Thv   updraft characteristics
+REAL, DIMENSION(KLON,KLEV), INTENT(INOUT) ::  PW_UP     ! vertical speed updraft characteristics
+REAL, DIMENSION(KLON,KLEV), INTENT(INOUT) ::  PFRAC_UP  ! updraft fraction
+REAL, DIMENSION(KLON,KLEV), INTENT(INOUT) ::  PEMF      ! updraft mass flux
+!
+!
+!*       0.2   Declarations of local variables :
+!
+INTEGER :: JRR           ! Loop index for the moist
+INTEGER :: IIB           ! Define the physical domain
+INTEGER :: IIE           !
+INTEGER :: IJB           !
+INTEGER :: IJE           !
+INTEGER :: IKB           !
+INTEGER :: IKE           !
+INTEGER :: IKA, IKU
+INTEGER :: JI, JJ, JL, JK !
+INTEGER ::II
+INTEGER, DIMENSION(size(PRHODJ,1)) :: IKLCL,IKETL,IKCTL
+REAL,DIMENSION(size(PRHODJ,1),size(PRHODJ,2)) :: ZFLXZTHMF,ZFLXZRMF,ZFLXZUMF,ZFLXZVMF
+REAL,DIMENSION(size(PRHODJ,1),size(PRHODJ,2)) :: ZDETR,ZENTR
+!
+!
+
+REAL          ::  ZIMPL        ! degree of implicitness
+!
+!
+!
+!------------------------------------------------------------------------------
+!
+!*       1.     PRELIMINARY COMPUTATIONS
+!               ------------------------
+!
+REAL(KIND=JPRB) :: ZHOOK_HANDLE
+IF (LHOOK) CALL DR_HOOK('ARO_SHALLOW_MF',0,ZHOOK_HANDLE)
+
+
+IIB=1+JPHEXT
+IIE=SIZE(PZZ,1) - JPHEXT
+IJB=1+JPHEXT
+IJE=1 - JPHEXT
+IF(KKL==1)THEN
+  IKA=1
+  IKU=SIZE(PZZ,2)
+ELSE
+  IKA=SIZE(PZZ,2)
+  IKU=1
+ENDIF
+IKB=IKA+KKL*JPVEXT
+IKE=IKU-KKL*JPVEXT
+!
+!
+!------------------------------------------------------------------------------
+!
+!*       2.   INITIALISATION
+!
+!             ---------------
+
+
+ZIMPL=1.
+!ZIMPL=0.
+! tableau a recalculer a chaque pas de temps
+! attention, ZDZZ est l'altitude entre deux niveaux (et pas l'�paisseur de la couche)
+
+!DO JL = IIB,IIE
+!   DO JK = 2, SIZE(PZZF,2)-1
+!      ZDZZ(JL,JK)=PZZF(JL,JK)-PZZF(JL,JK-KKL)
+!   ENDDO
+!   ZDZZ(JL,IKA)=PZZF(JL,IKA)-(1.5*PZZ(JL,IKA)-0.5*PZZ(JL,IKA+KKL)) ! must work with JPVEXT=0 or 1
+!   ZDZZ(JL,IKU)=PZZF(JL,IKU)-PZZF(JL,IKU-KKL) ! excluded from the loop because depending on KKL, IKU can be 1 or SIZE()
+!ENDDO
+!
+!
+!------------------------------------------------------------------------------
+!
+!
+!*       3.   MULTIPLICATION PAR RHODJ
+!             POUR OBTENIR LES TERMES SOURCES DE MESONH
+!
+!         -----------------------------------------------
+
+!
+!------------------------------------------------------------------------------
+!
+!
+!*       4.   APPEL DE LA TURBULENCE MESONH
+!
+!         ---------------------------------
+!
+  CALL SHALLOW_MF(KKA=IKA,KKU=IKU,KKL=KKL,KRR=KRR,KRRL=KRRL,KRRI=KRRI,                    &
+     &HMF_UPDRAFT=HMF_UPDRAFT, HMF_CLOUD=HMF_CLOUD,HFRAC_ICE=HFRAC_ICE,OMIXUV=OMIXUV,     &
+     &ONOMIXLG=ONOMIXLG,KSV_LGBEG=KSV_LGBEG,KSV_LGEND=KSV_LGEND,                          &
+     &PIMPL_MF=ZIMPL, PTSTEP=PTSTEP, PTSTEP_MET=PTSTEP, PTSTEP_SV=PTSTEP,                 &
+     &PDZZ=PDZZF,PZZ=PZZ,                                                                  &
+     &PRHODJ=PRHODJ,PRHODREF=PRHODREF,                                                    &
+     &PPABSM=PPABSM,PEXNM=PEXNM,                                                          &
+     &PSFTH=PSFTH,PSFRV=PSFRV,                                                            &
+     &PTHM=PTHM,PRM=PRM,PUM=PUM,PVM=PVM,PTKEM=PTKEM,PSVM=PSVM,                            &
+     &PDUDT_MF=PDUDT_MF,PDVDT_MF=PDVDT_MF,                                                &
+     &PDTHLDT_MF=PDTHLDT_MF,PDRTDT_MF=PDRTDT_MF,PDSVDT_MF=PDSVDT_MF,                      &
+     &PSIGMF=PSIGMF,PRC_MF=PRC_MF,PRI_MF=PRI_MF,PCF_MF=PCF_MF,PFLXZTHVMF=PFLXZTHVMF,      &
+     &PFLXZTHMF=ZFLXZTHMF,PFLXZRMF=ZFLXZRMF,PFLXZUMF=ZFLXZUMF,PFLXZVMF=ZFLXZVMF,          &
+     &PTHL_UP=PTHL_UP,PRT_UP=PRT_UP,PRV_UP=PRV_UP,PRC_UP=PRC_UP,PRI_UP=PRI_UP,            &
+     &PU_UP=PU_UP, PV_UP=PV_UP, PTHV_UP=PTHV_UP, PW_UP=PW_UP,                             &
+     &PFRAC_UP=PFRAC_UP,PEMF=PEMF,PDETR=ZDETR,PENTR=ZENTR,                                &
+     &KKLCL=IKLCL,KKETL=IKETL,KKCTL=IKCTL                                                 )
+!
+!
+!------------------------------------------------------------------------------
+!
+!
+!*       5.   DIVISION PAR RHODJ DES TERMES SOURCES DE MESONH
+!             (on obtient des termes homog�nes � des tendances)
+!
+!             -----------------------------------------------
+!
+IF (LHOOK) CALL DR_HOOK('ARO_SHALLOW_MF',1,ZHOOK_HANDLE)
+END SUBROUTINE ARO_SHALLOW_MF
diff --git a/src/arome/ext/aro_shallow_mf.h b/src/arome/ext/aro_shallow_mf.h
new file mode 100644
index 0000000000000000000000000000000000000000..07cdb7a7d38e75a76c8c45716c095ddcd22a0e1c
--- /dev/null
+++ b/src/arome/ext/aro_shallow_mf.h
@@ -0,0 +1,66 @@
+INTERFACE
+ SUBROUTINE ARO_SHALLOW_MF(KKL, KLON,KLEV, KRR, KRRL, KRRI,KSV,&
+ & HMF_UPDRAFT, HMF_CLOUD, HFRAC_ICE, OMIXUV,&
+ & ONOMIXLG,KSV_LGBEG,KSV_LGEND,&
+ & KTCOUNT, PTSTEP,&
+ & PZZ, PZZF,PDZZF,&
+ & PRHODJ, PRHODREF,&
+ & PPABSM, PEXNM,&
+ & PSFTH,PSFRV,&
+ & PTHM,PRM,&
+ & PUM,PVM,PTKEM,PSVM,&
+ & PDUDT_MF,PDVDT_MF,&
+ & PDTHLDT_MF,PDRTDT_MF,PDSVDT_MF,&
+ & PSIGMF,PRC_MF,PRI_MF,PCF_MF,PFLXZTHVMF,&
+ & PTHL_UP,PRT_UP,PRV_UP,PRC_UP,PRI_UP,&
+ & PU_UP, PV_UP, PTHV_UP, PW_UP, PFRAC_UP, PEMF) 
+USE PARKIND1  ,ONLY : JPIM     ,JPRB
+INTEGER(KIND=JPIM), INTENT(IN) :: KKL
+INTEGER(KIND=JPIM), INTENT(IN) :: KLON
+INTEGER(KIND=JPIM), INTENT(IN) :: KLEV
+INTEGER(KIND=JPIM), INTENT(IN) :: KRR
+INTEGER(KIND=JPIM), INTENT(IN) :: KRRL
+INTEGER(KIND=JPIM), INTENT(IN) :: KRRI
+INTEGER(KIND=JPIM), INTENT(IN) :: KSV
+CHARACTER (LEN=4), INTENT(IN) :: HMF_UPDRAFT
+CHARACTER (LEN=4), INTENT(IN) :: HMF_CLOUD
+CHARACTER*1, INTENT(IN) :: HFRAC_ICE
+LOGICAL, INTENT(IN) :: OMIXUV
+LOGICAL, INTENT(IN) :: ONOMIXLG
+INTEGER(KIND=JPIM), INTENT(IN) :: KSV_LGBEG
+INTEGER(KIND=JPIM), INTENT(IN) :: KSV_LGEND
+INTEGER(KIND=JPIM), INTENT(IN) :: KTCOUNT
+REAL(KIND=JPRB), INTENT(IN) :: PTSTEP
+REAL(KIND=JPRB), DIMENSION(KLON,KLEV), INTENT(IN) :: PZZ
+REAL(KIND=JPRB), DIMENSION(KLON,KLEV), INTENT(IN) :: PZZF
+REAL(KIND=JPRB), DIMENSION(KLON,KLEV), INTENT(IN) :: PDZZF
+REAL(KIND=JPRB), DIMENSION(KLON,KLEV), INTENT(IN) :: PRHODJ
+REAL(KIND=JPRB), DIMENSION(KLON,KLEV), INTENT(IN) :: PRHODREF
+REAL(KIND=JPRB), DIMENSION(KLON,KLEV), INTENT(IN) :: PPABSM
+REAL(KIND=JPRB), DIMENSION(KLON,KLEV), INTENT(IN) :: PEXNM
+REAL(KIND=JPRB), DIMENSION(KLON), INTENT(IN) :: PSFTH,PSFRV
+REAL(KIND=JPRB), DIMENSION(KLON,KLEV), INTENT(IN) :: PTHM
+REAL(KIND=JPRB), DIMENSION(KLON,KLEV,KRR), INTENT(IN) :: PRM
+REAL(KIND=JPRB), DIMENSION(KLON,KLEV), INTENT(IN) :: PUM,PVM
+REAL(KIND=JPRB), DIMENSION(KLON,KLEV), INTENT(IN) :: PTKEM
+REAL(KIND=JPRB), DIMENSION(KLON,KLEV,KSV), INTENT(IN) :: PSVM
+REAL(KIND=JPRB), DIMENSION(KLON,KLEV), INTENT(OUT):: PDUDT_MF
+REAL(KIND=JPRB), DIMENSION(KLON,KLEV), INTENT(OUT):: PDVDT_MF
+REAL(KIND=JPRB), DIMENSION(KLON,KLEV), INTENT(OUT):: PDTHLDT_MF
+REAL(KIND=JPRB), DIMENSION(KLON,KLEV), INTENT(OUT):: PDRTDT_MF
+REAL(KIND=JPRB), DIMENSION(KLON,KLEV,KSV), INTENT(OUT):: PDSVDT_MF
+REAL(KIND=JPRB), DIMENSION(KLON,KLEV), INTENT(OUT) :: PSIGMF,PRC_MF,PRI_MF,PCF_MF
+REAL(KIND=JPRB), DIMENSION(KLON,KLEV), INTENT(OUT) :: PFLXZTHVMF
+REAL(KIND=JPRB), DIMENSION(KLON,KLEV), INTENT(INOUT) :: PTHL_UP
+REAL(KIND=JPRB), DIMENSION(KLON,KLEV), INTENT(INOUT) :: PRT_UP
+REAL(KIND=JPRB), DIMENSION(KLON,KLEV), INTENT(INOUT) :: PRV_UP
+REAL(KIND=JPRB), DIMENSION(KLON,KLEV), INTENT(INOUT) :: PU_UP
+REAL(KIND=JPRB), DIMENSION(KLON,KLEV), INTENT(INOUT) :: PV_UP
+REAL(KIND=JPRB), DIMENSION(KLON,KLEV), INTENT(INOUT) :: PRC_UP
+REAL(KIND=JPRB), DIMENSION(KLON,KLEV), INTENT(INOUT) :: PRI_UP
+REAL(KIND=JPRB), DIMENSION(KLON,KLEV), INTENT(INOUT) :: PTHV_UP
+REAL(KIND=JPRB), DIMENSION(KLON,KLEV), INTENT(INOUT) :: PW_UP
+REAL(KIND=JPRB), DIMENSION(KLON,KLEV), INTENT(INOUT) :: PFRAC_UP
+REAL(KIND=JPRB), DIMENSION(KLON,KLEV), INTENT(INOUT) :: PEMF
+END SUBROUTINE ARO_SHALLOW_MF
+END INTERFACE
diff --git a/src/arome/gmkpack_ignored_files b/src/arome/gmkpack_ignored_files
index d6a401c5fc5882ea03c26f90c5e295b7043c5580..167df4f8ee9b06289c49a59b218c8a175842d2af 100644
--- a/src/arome/gmkpack_ignored_files
+++ b/src/arome/gmkpack_ignored_files
@@ -149,3 +149,9 @@ phyex/turb/compute_bl89_ml.F90
 phyex/turb/modi_compute_bl89_ml.F90
 phyex/turb/compute_entr_detr.F90
 phyex/turb/modi_compute_entr_detr.F90
+phyex/turb/compute_updraft_raha.F90
+phyex/turb/modi_compute_updraft_raha.F90
+phyex/turb/modi_compute_updraft_rhcj10.F90
+phyex/turb/compute_updraft_rhcj10.F90
+phyex/turb/modi_compute_updraft.F90
+phyex/turb/compute_updraft.F90
diff --git a/src/arome/turb/compute_updraft_rhcj10.F90 b/src/arome/turb/mode_compute_updraft_rhcj10.F90
similarity index 68%
rename from src/arome/turb/compute_updraft_rhcj10.F90
rename to src/arome/turb/mode_compute_updraft_rhcj10.F90
index a120b45ba043f4e4013513e61c4778ceec70015b..b88864cc20b2af9f5e0b26cf36ed004429c16e5c 100644
--- a/src/arome/turb/compute_updraft_rhcj10.F90
+++ b/src/arome/turb/mode_compute_updraft_rhcj10.F90
@@ -1,5 +1,16 @@
+!MNH_LIC Copyright 2012-2019 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence
+!MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt
+!MNH_LIC for details. version 1.
+!-----------------------------------------------------------------
 !     ######spl
-      SUBROUTINE COMPUTE_UPDRAFT_RHCJ10(KKA,KKB,KKE,KKU,KKL,HFRAC_ICE,       &
+     MODULE MODE_COMPUTE_UPDRAFT_RHCJ10
+!    ###########################
+!
+IMPLICIT NONE
+CONTAINS
+!
+SUBROUTINE COMPUTE_UPDRAFT_RHCJ10(KKA,KKB,KKE,KKU,KKL,HFRAC_ICE,       &
                                  OENTR_DETR,OMIXUV,               &
                                  ONOMIXLG,KSV_LGBEG,KSV_LGEND,    &
                                  PZZ,PDZZ,                        &
@@ -13,9 +24,6 @@
                                  PEMF,PDETR,PENTR,                &
                                  PBUO_INTEG,KKLCL,KKETL,KKCTL,    &
                                  PDEPTH     )
-
-      USE PARKIND1, ONLY : JPRB
-      USE YOMHOOK , ONLY : LHOOK, DR_HOOK
 !     #################################################################
 !!
 !!****  *COMPUTE_UPDRAFT_RHCJ10* - calculates caracteristics of the updraft 
@@ -42,11 +50,13 @@
 !!     AUTHOR
 !!     ------
 !!     Y. Bouteloup (2012)
-!!     R. Honert Janv 2013 ==> corection of some coding bugs
+!!     R. Honert Janv 2013 ==> corection of some bugs
 !!     R. El Khatib 15-Oct-2014 Optimization
 !!     Q.Rodier  01/2019 : support RM17 mixing length
 !! --------------------------------------------------------------------------
-!
+
+! WARNING ==>  This updraft is not yet ready to use scalar variables 
+
 !*      0. DECLARATIONS
 !          ------------
 !
@@ -57,6 +67,8 @@ USE MODE_TH_R_FROM_THL_RT_1D, ONLY: TH_R_FROM_THL_RT_1D
 USE MODI_SHUMAN_MF, ONLY: MZF_MF, MZM_MF, GZ_M_W_MF
 
 USE MODE_COMPUTE_BL89_ML, ONLY: COMPUTE_BL89_ML
+USE PARKIND1, ONLY : JPRB
+USE YOMHOOK , ONLY : LHOOK, DR_HOOK
 
 
 IMPLICIT NONE
@@ -64,13 +76,12 @@ IMPLICIT NONE
 !*                    1.1  Declaration of Arguments
 !
 !
-!
 INTEGER,                INTENT(IN)   :: KKA          ! near ground array index
 INTEGER,                INTENT(IN)   :: KKB          ! near ground physical index
 INTEGER,                INTENT(IN)   :: KKE          ! uppest atmosphere physical index
 INTEGER,                INTENT(IN)   :: KKU          ! uppest atmosphere array index
 INTEGER,                INTENT(IN)   :: KKL          ! +1 if grid goes from ground to atmosphere top, -1 otherwise
-CHARACTER*1,            INTENT(IN)   :: HFRAC_ICE    ! partition liquid/ice scheme
+CHARACTER(LEN=1),       INTENT(IN)   :: HFRAC_ICE    ! partition liquid/ice scheme
 LOGICAL,                INTENT(IN) :: OENTR_DETR! flag to recompute entrainment, detrainment and mass flux
 LOGICAL,                INTENT(IN) :: OMIXUV    ! True if mixing of momentum
 LOGICAL,                INTENT(IN)   :: ONOMIXLG  ! False if mixing of lagrangian tracer
@@ -78,8 +89,6 @@ INTEGER,                INTENT(IN)   :: KSV_LGBEG ! first index of lag. tracer
 INTEGER,                INTENT(IN)   :: KSV_LGEND ! last  index of lag. tracer
 REAL, DIMENSION(:,:), INTENT(IN)   :: PZZ       !  Height at the flux point
 REAL, DIMENSION(:,:), INTENT(IN)   :: PDZZ      !  Metrics coefficient
-!REAL, DIMENSION(:,:), INTENT(IN)   :: PEXNMH    !  Exner on flux level
-!REAL, DIMENSION(:,:), INTENT(INOUT) :: PCPTPPHY_UP  ! CpT+PHY of updraft
  
 REAL, DIMENSION(:),   INTENT(IN)   ::  PSFTH,PSFRV
 ! normal surface fluxes of theta,rv,(u,v) parallel to the orography
@@ -90,7 +99,7 @@ REAL, DIMENSION(:,:),   INTENT(IN) ::  PRHODREF   ! dry density of the
 REAL, DIMENSION(:,:),   INTENT(IN) ::  PUM        ! u mean wind
 REAL, DIMENSION(:,:),   INTENT(IN) ::  PVM        ! v mean wind
 REAL, DIMENSION(:,:),   INTENT(IN) ::  PTKEM      ! TKE at t-dt
-
+!
 REAL, DIMENSION(:,:),   INTENT(IN)   ::  PEXNM       ! Exner function at t-dt
 REAL, DIMENSION(:,:),   INTENT(IN)   ::  PTHM           ! pot. temp. at t-dt
 REAL, DIMENSION(:,:),   INTENT(IN)   ::  PRVM           ! vapor mixing ratio at t-dt
@@ -116,7 +125,6 @@ INTEGER, DIMENSION(:),  INTENT(INOUT)::  KKLCL,KKETL,KKCTL! LCL, ETL, CTL
 REAL, DIMENSION(:),     INTENT(OUT)   :: PDEPTH           ! Deepness of cloud
 !                       1.2  Declaration of local variables
 !
-!
 ! Mean environment variables at t-dt at flux point
 REAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) ::    ZTHM_F,ZRVM_F    ! Theta,rv of
                                                                   ! updraft environnement
@@ -139,7 +147,7 @@ REAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) :: ZTHS_UP,ZTHSM
 
 REAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) ::  ZCOEF  ! diminution coefficient for too high clouds 
                         
-REAL, DIMENSION(SIZE(PSFTH,1) )            ::  ZWTHVSURF  ! Surface w'thetav'
+REAL                                       ::  ZWTHVSURF  ! Surface w'thetav'
 
 REAL  :: ZRVORD       ! RV/RD
 
@@ -159,7 +167,6 @@ LOGICAL                          ::  GLMIX
 LOGICAL, DIMENSION(SIZE(PTHM,1))              :: GWORK1
 LOGICAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) :: GWORK2
 
-
 INTEGER  :: ITEST
 
 REAL, DIMENSION(SIZE(PTHM,1)) :: ZRC_UP, ZRI_UP, ZRV_UP, ZRSATW, ZRSATI
@@ -173,7 +180,6 @@ REAL, DIMENSION(SIZE(PTHM,1))              ::  ZW_MAX               ! w**2  max
 REAL, DIMENSION(SIZE(PTHM,1))              ::  ZZTOP                ! Top of the updraft
 REAL, DIMENSION(SIZE(PTHM,1))              ::  ZQTM,ZQT_UP
 
-
 REAL  :: ZDEPTH_MAX1, ZDEPTH_MAX2 ! control auto-extinction process
 
 REAL  :: ZTMAX,ZRMAX, ZEPS  ! control value
@@ -182,7 +188,6 @@ REAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) :: ZSHEAR,ZDUDZ,ZDVDZ ! vertical wind
 REAL(KIND=JPRB) :: ZHOOK_HANDLE
 IF (LHOOK) CALL DR_HOOK('COMPUTE_UPDRAFT_RHCJ10',0,ZHOOK_HANDLE)
 
-
 ! Thresholds for the  perturbation of
 ! theta_l and r_t at the first level of the updraft
 
@@ -195,6 +200,7 @@ ZEPS=1.E-15
 ! Initialisation of the constants   
 ZRVORD   = (XRV / XRD) 
 
+! depth are different in compute_updraft (3000. and 4000.) ==> impact is small
 ZDEPTH_MAX1=4500. ! clouds with depth infeRIOr to this value are keeped untouched
 ZDEPTH_MAX2=5000. ! clouds with depth superior to this value are suppressed
 
@@ -231,13 +237,11 @@ PTHV_UP(:,:)=0.
 PBUO_INTEG=0.
 ZBUO      =0.
 
-!no ice cloud coded yet 
+!no ice cloud coded yet
 PRI_UP(:,:)=0.
 PFRAC_ICE_UP(:,:)=0.
 PRSAT_UP(:,:)=PRVM(:,:) ! should be initialised correctly but is (normaly) not used
 
-
-
 ! Initialisation of environment variables at t-dt
 
 ! variables at flux level
@@ -247,11 +251,15 @@ ZUM_F  (:,:) = MZM_MF(PUM(:,:), KKA, KKU, KKL)
 ZVM_F  (:,:) = MZM_MF(PVM(:,:), KKA, KKU, KKL)
 ZTKEM_F(:,:) = MZM_MF(PTKEM(:,:), KKA, KKU, KKL)
 
-!DO JSV=1,ISV 
-! IF (ONOMIXLG .AND. JSV >= KSV_LGBEG .AND. JSV<= KSV_LGEND) CYCLE
+! This updraft is not yet ready to use scalar variables
+!DO JSV=1,ISV
+!  IF (ONOMIXLG .AND. JSV >= KSV_LGBEG .AND. JSV<= KSV_LGEND) CYCLE
+! *** SR merge AROME/Méso-nh: following two lines come from the AROME version
 !   ZSVM_F(:,KKB:IKU,JSV) = 0.5*(PSVM(:,KKB:IKU,JSV)+PSVM(:,1:IKU-1,JSV))
 !   ZSVM_F(:,1,JSV)       = ZSVM_F(:,KKB,JSV) 
-!END DO  
+! *** the following single line comes from the Meso-NH version
+!  ZSVM_F(:,:,JSV) = MZM_MF(KKA,KKU,KKL,PSVM(:,:,JSV))
+!END DO
 
 !          Initialisation of updraft characteristics 
 PTHL_UP(:,:)=ZTHLM_F(:,:)
@@ -259,6 +267,7 @@ PRT_UP(:,:)=ZRTM_F(:,:)
 PU_UP(:,:)=ZUM_F(:,:)
 PV_UP(:,:)=ZVM_F(:,:)
 PSV_UP(:,:,:)=0.
+! This updraft is not yet ready to use scalar variables
 !IF (ONOMIXLG .AND. JSV >= KSV_LGBEG .AND. JSV<= KSV_LGEND) then
 !    PSV_UP(:,:,:)=ZSVM_F(:,:,:)
 !ENDIF
@@ -267,11 +276,11 @@ PSV_UP(:,:,:)=0.
 ! thetal_up,rt_up,thetaV_up, w,Buoyancy term and mass flux (PEMF)
 
 DO JI=1,IIJU
-PTHL_UP(JI,KKB)= ZTHLM_F(JI,KKB)+MAX(0.,MIN(ZTMAX,(PSFTH(JI)/SQRT(ZTKEM_F(JI,KKB)))*XALP_PERT))
-PRT_UP(JI,KKB) = ZRTM_F(JI,KKB)+MAX(0.,MIN(ZRMAX,(PSFRV(JI)/SQRT(ZTKEM_F(JI,KKB)))*XALP_PERT)) 
-
-ZQT_UP(JI) = PRT_UP(JI,KKB)/(1.+PRT_UP(JI,KKB))
-ZTHS_UP(JI,KKB)=PTHL_UP(JI,KKB)*(1.+XLAMBDA_MF*ZQT_UP(JI))
+  PTHL_UP(JI,KKB)= ZTHLM_F(JI,KKB)+MAX(0.,MIN(ZTMAX,(PSFTH(JI)/SQRT(ZTKEM_F(JI,KKB)))*XALP_PERT))
+  PRT_UP(JI,KKB) = ZRTM_F(JI,KKB)+MAX(0.,MIN(ZRMAX,(PSFRV(JI)/SQRT(ZTKEM_F(JI,KKB)))*XALP_PERT)) 
+  
+  ZQT_UP(JI) = PRT_UP(JI,KKB)/(1.+PRT_UP(JI,KKB))
+  ZTHS_UP(JI,KKB)=PTHL_UP(JI,KKB)*(1.+XLAMBDA_MF*ZQT_UP(JI))
 ENDDO
 
 ZTHM_F (:,:) = MZM_MF(PTHM (:,:), KKA, KKU, KKL)
@@ -281,9 +290,9 @@ ZRVM_F (:,:) = MZM_MF(PRVM(:,:), KKA, KKU, KKL)
 
 ! thetav at mass and flux levels 
 DO JK=1,IKU
-DO JI=1,IIJU
-ZTHVM_F(JI,JK)=ZTHM_F(JI,JK)*((1.+ZRVORD*ZRVM_F(JI,JK))/(1.+ZRTM_F(JI,JK)))
-ENDDO
+  DO JI=1,IIJU
+    ZTHVM_F(JI,JK)=ZTHM_F(JI,JK)*((1.+ZRVORD*ZRVM_F(JI,JK))/(1.+ZRTM_F(JI,JK)))
+  ENDDO
 ENDDO
 
 PTHV_UP(:,:)= ZTHVM_F(:,:)
@@ -302,13 +311,13 @@ CALL TH_R_FROM_THL_RT_1D(HFRAC_ICE,PFRAC_ICE_UP(:,KKB),ZPRES_F(:,KKB), &
              PRV_UP(:,KKB),PRC_UP(:,KKB),PRI_UP(:,KKB),ZRSATW(:),ZRSATI(:))
 
 DO JI=1,IIJU
-! compute updraft thevav and buoyancy term at KKB level             
-PTHV_UP(JI,KKB) = ZTH_UP(JI,KKB)*((1+ZRVORD*PRV_UP(JI,KKB))/(1+PRT_UP(JI,KKB))) 
-! compute mean rsat in updraft
-PRSAT_UP(JI,KKB) = ZRSATW(JI)*(1-PFRAC_ICE_UP(JI,KKB)) + ZRSATI(JI)*PFRAC_ICE_UP(JI,KKB)
+  ! compute updraft thevav and buoyancy term at KKB level             
+  PTHV_UP(JI,KKB) = ZTH_UP(JI,KKB)*((1+ZRVORD*PRV_UP(JI,KKB))/(1+PRT_UP(JI,KKB))) 
+  ! compute mean rsat in updraft
+  PRSAT_UP(JI,KKB) = ZRSATW(JI)*(1-PFRAC_ICE_UP(JI,KKB)) + ZRSATI(JI)*PFRAC_ICE_UP(JI,KKB)
 ENDDO
 
-!Tout est commente pour tester dans un premier temps la s�paration en deux de la 
+!Tout est commente pour tester dans un premier temps la separation en deux de la 
 !  boucle verticale, une pour w et une pour PEMF                                                            
 
 ZG_O_THVREF=XG/ZTHVM_F
@@ -316,9 +325,9 @@ ZG_O_THVREF=XG/ZTHVM_F
 ! Calcul de la fermeture de Julien Pergaut comme limite max de PHY
 
 DO JK=KKB,KKE-KKL,KKL   !  Vertical loop
-DO JI=1,IIJU
-  ZZDZ(JI,JK)   = MAX(ZEPS,PZZ(JI,JK+KKL)-PZZ(JI,JK))  ! <== Delta Z between two flux level
-ENDDO
+  DO JI=1,IIJU
+    ZZDZ(JI,JK)   = MAX(ZEPS,PZZ(JI,JK+KKL)-PZZ(JI,JK))  ! <== Delta Z between two flux level
+  ENDDO
 ENDDO
 
 ! compute L_up
@@ -341,21 +350,21 @@ CALL COMPUTE_BL89_ML(KKA,KKB,KKE,KKU,KKL,PDZZ,ZTKEM_F(:,KKB),ZG_O_THVREF(:,KKB),
 ZLUP(:)=MAX(ZLUP(:),1.E-10)
 
 DO JI=1,IIJU
-! Compute Buoyancy flux at the ground
-ZWTHVSURF(JI) = (ZTHVM_F(JI,KKB)/ZTHM_F(JI,KKB))*PSFTH(JI)+     &
-                (0.61*ZTHM_F(JI,KKB))*PSFRV(JI)
-! Mass flux at KKB level (updraft triggered if PSFTH>0.)
-
-IF (ZWTHVSURF(JI)>0.010) THEN
-  PEMF(JI,KKB) = XCMF * ZRHO_F(JI,KKB) * ((ZG_O_THVREF(JI,KKB))*ZWTHVSURF(JI)*ZLUP(JI))**(1./3.)
-  PFRAC_UP(JI,KKB)=MIN(PEMF(JI,KKB)/(SQRT(ZW_UP2(JI,KKB))*ZRHO_F(JI,KKB)),XFRAC_UP_MAX)
-  PEMF(JI,KKB) = ZRHO_F(JI,KKB)*PFRAC_UP(JI,KKB)*SQRT(ZW_UP2(JI,KKB))
-!  ZW_UP2(JI,KKB)=(PEMF(JI,KKB)/(PFRAC_UP(JI,KKB)*ZRHO_F(JI,KKB)))**2
-  GTEST(JI)=.TRUE.
-ELSE
-  PEMF(JI,KKB) =0.
-  GTEST(JI)=.FALSE.
-ENDIF
+  ! Compute Buoyancy flux at the ground
+  ZWTHVSURF = (ZTHVM_F(JI,KKB)/ZTHM_F(JI,KKB))*PSFTH(JI)+     &
+              (0.61*ZTHM_F(JI,KKB))*PSFRV(JI)
+
+  ! Mass flux at KKB level (updraft triggered if PSFTH>0.)
+  IF (ZWTHVSURF>0.010) THEN ! <==  Not 0 Important to have stratocumulus !!!!!
+    PEMF(JI,KKB) = XCMF * ZRHO_F(JI,KKB) * ((ZG_O_THVREF(JI,KKB))*ZWTHVSURF*ZLUP(JI))**(1./3.)
+    PFRAC_UP(JI,KKB)=MIN(PEMF(JI,KKB)/(SQRT(ZW_UP2(JI,KKB))*ZRHO_F(JI,KKB)),XFRAC_UP_MAX)
+    PEMF(JI,KKB) = ZRHO_F(JI,KKB)*PFRAC_UP(JI,KKB)*SQRT(ZW_UP2(JI,KKB))
+  !  ZW_UP2(JI,KKB)=(PEMF(JI,KKB)/(PFRAC_UP(JI,KKB)*ZRHO_F(JI,KKB)))**2
+    GTEST(JI)=.TRUE.
+  ELSE
+    PEMF(JI,KKB) =0.
+    GTEST(JI)=.FALSE.
+  ENDIF
 ENDDO
 
   
@@ -375,14 +384,12 @@ GTESTLCL(:)=.FALSE.
 ZW_MAX(:)      = 0.
 ZZTOP(:)       = 0.
 
-!GTEST(:) = (ZW_UP2(:,KKB)>ZEPS)
-
 DO JK=KKB,KKE-KKL,KKL
 
 ! IF the updraft top is reached for all column, stop the loop on levels
 
-!  ITEST=COUNT(GTEST)
-!  IF (ITEST==0) CYCLE
+  !ITEST=COUNT(GTEST)
+  !IF (ITEST==0) CYCLE
 
 !       Computation of entrainment and detrainment with KF90
 !       parameterization in clouds and LR01 in subcloud layer
@@ -391,10 +398,10 @@ DO JK=KKB,KKE-KKL,KKL
 ! to find the LCL (check if JK is LCL or not)
 
   DO JI=1,IIJU
-  IF ((PRC_UP(JI,JK)+PRI_UP(JI,JK)>0.).AND.(.NOT.(GTESTLCL(JI)))) THEN
+    IF ((PRC_UP(JI,JK)+PRI_UP(JI,JK)>0.).AND.(.NOT.(GTESTLCL(JI)))) THEN
       KKLCL(JI) = JK           
       GTESTLCL(JI)=.TRUE.
-  ENDIF
+    ENDIF
   ENDDO
 
 
@@ -421,46 +428,43 @@ DO JK=KKB,KKE-KKL,KKL
       ZDZ(JI)   = MAX(ZEPS,PZZ(JI,JK+KKL)-PZZ(JI,JK))
       ZTEST(JI) = XA1*ZBUO(JI,JK) -  XB*ZW_UP2(JI,JK)
 
-!  Ancien calcul de la vitesse
-    
+      !  Ancien calcul de la vitesse
       ZCOE(JI)      = ZDZ(JI)
       IF (ZTEST(JI)>0.) THEN
         ZCOE(JI)    = ZDZ(JI)/(1.+ XBETA1)
       ENDIF
 
-!  Calcul de la vitesse
-
+      !  Convective Vertical speed computation
       ZWCOE(JI)         = (1.-XB*ZCOE(JI))/(1.+XB*ZCOE(JI))
       ZBUCOE(JI)        =  2.*ZCOE(JI)/(1.+XB*ZCOE(JI))
 
-! Second Rachel bug correction (XA1 has been forgotten)
+      ! Second Rachel bug correction (XA1 has been forgotten)
      
       ZW_UP2(JI,JK+KKL) = MAX(ZEPS,ZW_UP2(JI,JK)*ZWCOE(JI) + XA1*ZBUO(JI,JK)*ZBUCOE(JI) )  
       ZW_MAX(JI) = MAX(ZW_MAX(JI), SQRT(ZW_UP2(JI,JK+KKL)))
       ZWUP_MEAN(JI)     = MAX(ZEPS,0.5*(ZW_UP2(JI,JK+KKL)+ZW_UP2(JI,JK)))
  
-!  Entrainement et detrainement
+      !  Entrainement and detrainement
 
-! First Rachel bug correction (Parenthesis around 1+beta1 ==> impact is small)   
-     PENTR(JI,JK)  = MAX(0.,(XBETA1/(1.+XBETA1))*(XA1*ZBUO(JI,JK)/ZWUP_MEAN(JI)-XB))
-     ZDETR_BUO(JI) = MAX(0., -(XBETA1/(1.+XBETA1))*XA1*ZBUO(JI,JK)/ZWUP_MEAN(JI))
-     ZDETR_RT(JI)  = XC*SQRT(MAX(0.,(PRT_UP(JI,JK) - ZRTM_F(JI,JK))) / MAX(ZEPS,ZRTM_F(JI,JK)) / ZWUP_MEAN(JI))
-     PDETR(JI,JK)  = ZDETR_RT(JI)+ZDETR_BUO(JI)
+      ! First Rachel bug correction (Parenthesis around 1+beta1 ==> impact is small)   
+      PENTR(JI,JK)  = MAX(0.,(XBETA1/(1.+XBETA1))*(XA1*ZBUO(JI,JK)/ZWUP_MEAN(JI)-XB))
+      ZDETR_BUO(JI) = MAX(0., -(XBETA1/(1.+XBETA1))*XA1*ZBUO(JI,JK)/ZWUP_MEAN(JI))
+      ZDETR_RT(JI)  = XC*SQRT(MAX(0.,(PRT_UP(JI,JK) - ZRTM_F(JI,JK))) / MAX(ZEPS,ZRTM_F(JI,JK)) / ZWUP_MEAN(JI))
+      PDETR(JI,JK)  = ZDETR_RT(JI)+ZDETR_BUO(JI)
    
-! If the updraft did not stop, compute cons updraft characteritics at jk+1
-
-     ZZTOP(JI) = MAX(ZZTOP(JI),PZZ(JI,JK+KKL))
-     ZMIX2(JI) = (PZZ(JI,JK+KKL)-PZZ(JI,JK))*PENTR(JI,JK) !&
+      ! If the updraft did not stop, compute cons updraft characteritics at jk+1
+      ZZTOP(JI) = MAX(ZZTOP(JI),PZZ(JI,JK+KKL))
+      ZMIX2(JI) = (PZZ(JI,JK+KKL)-PZZ(JI,JK))*PENTR(JI,JK) !&
                 
-!  Utilisation de thetaS
-     ZQTM(JI) = PRTM(JI,JK)/(1.+PRTM(JI,JK))            
-     ZTHSM(JI,JK) = PTHLM(JI,JK)*(1.+XLAMBDA_MF*ZQTM(JI))
-     ZTHS_UP(JI,JK+KKL)=(ZTHS_UP(JI,JK)*(1.-0.5*ZMIX2(JI)) + ZTHSM(JI,JK)*ZMIX2(JI)) &
-                          /(1.+0.5*ZMIX2(JI))                   
-     PRT_UP(JI,JK+KKL) =(PRT_UP (JI,JK)*(1.-0.5*ZMIX2(JI)) + PRTM(JI,JK)*ZMIX2(JI))  &
-                          /(1.+0.5*ZMIX2(JI))
-     ZQT_UP(JI) = PRT_UP(JI,JK+KKL)/(1.+PRT_UP(JI,JK+KKL))
-     PTHL_UP(JI,JK+KKL)=ZTHS_UP(JI,JK+KKL)/(1.+XLAMBDA_MF*ZQT_UP(JI))                      
+      !  Utilisation de thetaS
+      ZQTM(JI) = PRTM(JI,JK)/(1.+PRTM(JI,JK))            
+      ZTHSM(JI,JK) = PTHLM(JI,JK)*(1.+XLAMBDA_MF*ZQTM(JI))
+      ZTHS_UP(JI,JK+KKL)=(ZTHS_UP(JI,JK)*(1.-0.5*ZMIX2(JI)) + ZTHSM(JI,JK)*ZMIX2(JI)) &
+                           /(1.+0.5*ZMIX2(JI))                   
+      PRT_UP(JI,JK+KKL) =(PRT_UP (JI,JK)*(1.-0.5*ZMIX2(JI)) + PRTM(JI,JK)*ZMIX2(JI))  &
+                           /(1.+0.5*ZMIX2(JI))
+      ZQT_UP(JI) = PRT_UP(JI,JK+KKL)/(1.+PRT_UP(JI,JK+KKL))
+      PTHL_UP(JI,JK+KKL)=ZTHS_UP(JI,JK+KKL)/(1.+XLAMBDA_MF*ZQT_UP(JI))                      
     ENDIF  ! GTEST
   ENDDO
   
@@ -468,34 +472,36 @@ DO JK=KKB,KKE-KKL,KKL
   IF(OMIXUV) THEN
     IF(JK/=KKB) THEN
       DO JI=1,IIJU
-      IF(GTEST(JI)) THEN
-        PU_UP(JI,JK+KKL) = (PU_UP (JI,JK)*(1-0.5*ZMIX2(JI)) + PUM(JI,JK)*ZMIX2(JI)+ &
-                          0.5*XPRES_UV*(PZZ(JI,JK+KKL)-PZZ(JI,JK))*&
-                          ((PUM(JI,JK+KKL)-PUM(JI,JK))/PDZZ(JI,JK+KKL)+&
-                           (PUM(JI,JK)-PUM(JI,JK-KKL))/PDZZ(JI,JK))        )   &
-                          /(1+0.5*ZMIX2(JI))
-        PV_UP(JI,JK+KKL) = (PV_UP (JI,JK)*(1-0.5*ZMIX2(JI)) + PVM(JI,JK)*ZMIX2(JI)+ &
-                          0.5*XPRES_UV*(PZZ(JI,JK+KKL)-PZZ(JI,JK))*&
-                          ((PVM(JI,JK+KKL)-PVM(JI,JK))/PDZZ(JI,JK+KKL)+&
-                           (PVM(JI,JK)-PVM(JI,JK-KKL))/PDZZ(JI,JK))    )   &
-                          /(1+0.5*ZMIX2(JI))
-      ENDIF
+        IF(GTEST(JI)) THEN
+          PU_UP(JI,JK+KKL) = (PU_UP (JI,JK)*(1-0.5*ZMIX2(JI)) + PUM(JI,JK)*ZMIX2(JI)+ &
+                            0.5*XPRES_UV*(PZZ(JI,JK+KKL)-PZZ(JI,JK))*&
+                            ((PUM(JI,JK+KKL)-PUM(JI,JK))/PDZZ(JI,JK+KKL)+&
+                             (PUM(JI,JK)-PUM(JI,JK-KKL))/PDZZ(JI,JK))        )   &
+                            /(1+0.5*ZMIX2(JI))
+          PV_UP(JI,JK+KKL) = (PV_UP (JI,JK)*(1-0.5*ZMIX2(JI)) + PVM(JI,JK)*ZMIX2(JI)+ &
+                            0.5*XPRES_UV*(PZZ(JI,JK+KKL)-PZZ(JI,JK))*&
+                            ((PVM(JI,JK+KKL)-PVM(JI,JK))/PDZZ(JI,JK+KKL)+&
+                             (PVM(JI,JK)-PVM(JI,JK-KKL))/PDZZ(JI,JK))    )   &
+                            /(1+0.5*ZMIX2(JI))
+        ENDIF
       ENDDO
     ELSE
       DO JI=1,IIJU
-      IF(GTEST(JI)) THEN
-        PU_UP(JI,JK+KKL) = (PU_UP (JI,JK)*(1-0.5*ZMIX2(JI)) + PUM(JI,JK)*ZMIX2(JI)+ &
-                          0.5*XPRES_UV*(PZZ(JI,JK+KKL)-PZZ(JI,JK))*&
-                          ((PUM(JI,JK+KKL)-PUM(JI,JK))/PDZZ(JI,JK+KKL))        ) &
-                          /(1+0.5*ZMIX2(JI))
-        PV_UP(JI,JK+KKL) = (PV_UP (JI,JK)*(1-0.5*ZMIX2(JI)) + PVM(JI,JK)*ZMIX2(JI)+ &
-                          0.5*XPRES_UV*(PZZ(JI,JK+KKL)-PZZ(JI,JK))*&
-                          ((PVM(JI,JK+KKL)-PVM(JI,JK))/PDZZ(JI,JK+KKL))    )   &
-                          /(1+0.5*ZMIX2(JI))
-      ENDIF
+        IF(GTEST(JI)) THEN
+          PU_UP(JI,JK+KKL) = (PU_UP (JI,JK)*(1-0.5*ZMIX2(JI)) + PUM(JI,JK)*ZMIX2(JI)+ &
+                            0.5*XPRES_UV*(PZZ(JI,JK+KKL)-PZZ(JI,JK))*&
+                            ((PUM(JI,JK+KKL)-PUM(JI,JK))/PDZZ(JI,JK+KKL))        ) &
+                            /(1+0.5*ZMIX2(JI))
+          PV_UP(JI,JK+KKL) = (PV_UP (JI,JK)*(1-0.5*ZMIX2(JI)) + PVM(JI,JK)*ZMIX2(JI)+ &
+                            0.5*XPRES_UV*(PZZ(JI,JK+KKL)-PZZ(JI,JK))*&
+                            ((PVM(JI,JK+KKL)-PVM(JI,JK))/PDZZ(JI,JK+KKL))    )   &
+                            /(1+0.5*ZMIX2(JI))
+        ENDIF
       ENDDO
     ENDIF
   ENDIF
+  
+! This updraft is not yet ready to use scalar variables
 !  DO JSV=1,ISV 
 !     IF (ONOMIXLG .AND. JSV >= KSV_LGBEG .AND. JSV<= KSV_LGEND) CYCLE
 !      WHERE(GTEST) 
@@ -514,50 +520,51 @@ DO JK=KKB,KKE-KKL,KKL
           ZRV_UP(:),ZRC_UP(:),ZRI_UP(:),ZRSATW(:),ZRSATI(:))
 
   DO JI=1,IIJU
-  IF(GTEST(JI)) THEN
-    ZT_UP(JI) = ZTH_UP(JI,JK+KKL)*PEXNM(JI,JK+KKL)
-    ZCP(JI) = XCPD + XCL * ZRC_UP(JI)
-    ZLVOCPEXN(JI)=(XLVTT + (XCPV-XCL) *  (ZT_UP(JI)-XTT) ) / ZCP(JI) / PEXNM(JI,JK+KKL)
-    PRC_UP(JI,JK+KKL)=MIN(0.5E-3,ZRC_UP(JI))  ! On ne peut depasser 0.5 g/kg (autoconversion donc elimination !)
-    PTHL_UP(JI,JK+KKL) = PTHL_UP(JI,JK+KKL)+ZLVOCPEXN(JI)*(ZRC_UP(JI)-PRC_UP(JI,JK+KKL))
-    PRV_UP(JI,JK+KKL)=ZRV_UP(JI)
-    PRI_UP(JI,JK+KKL)=ZRI_UP(JI)
-    PRT_UP(JI,JK+KKL)  = PRC_UP(JI,JK+KKL) + PRV_UP(JI,JK+KKL)
-    PRSAT_UP(JI,JK+KKL) = ZRSATW(JI)*(1-PFRAC_ICE_UP(JI,JK+KKL)) + ZRSATI(JI)*PFRAC_ICE_UP(JI,JK+KKL)
-! Compute the updraft theta_v, buoyancy and w**2 for level JK+1   
-!      PTHV_UP(:,JK+KKL) = PTH_UP(:,JK+KKL)*((1+ZRVORD*PRV_UP(:,JK+KKL))/(1+PRT_UP(:,JK+KKL)))
+    IF(GTEST(JI)) THEN
+      ZT_UP(JI) = ZTH_UP(JI,JK+KKL)*PEXNM(JI,JK+KKL)
+      ZCP(JI) = XCPD + XCL * ZRC_UP(JI)
+      ZLVOCPEXN(JI)=(XLVTT + (XCPV-XCL) *  (ZT_UP(JI)-XTT) ) / ZCP(JI) / PEXNM(JI,JK+KKL)
+      PRC_UP(JI,JK+KKL)=MIN(0.5E-3,ZRC_UP(JI))  ! On ne peut depasser 0.5 g/kg (autoconversion donc elimination !)
+      PTHL_UP(JI,JK+KKL) = PTHL_UP(JI,JK+KKL)+ZLVOCPEXN(JI)*(ZRC_UP(JI)-PRC_UP(JI,JK+KKL))
+      PRV_UP(JI,JK+KKL)=ZRV_UP(JI)
+      PRI_UP(JI,JK+KKL)=ZRI_UP(JI)
+      PRT_UP(JI,JK+KKL)  = PRC_UP(JI,JK+KKL) + PRV_UP(JI,JK+KKL)
+      PRSAT_UP(JI,JK+KKL) = ZRSATW(JI)*(1-PFRAC_ICE_UP(JI,JK+KKL)) + ZRSATI(JI)*PFRAC_ICE_UP(JI,JK+KKL)
+
+      ! Compute the updraft theta_v, buoyancy and w**2 for level JK+1   
+      !PTHV_UP(:,JK+KKL) = PTH_UP(:,JK+KKL)*((1+ZRVORD*PRV_UP(:,JK+KKL))/(1+PRT_UP(:,JK+KKL)))
       PTHV_UP(JI,JK+KKL) = ZTH_UP(JI,JK+KKL)*(1.+0.608*PRV_UP(JI,JK+KKL) - PRC_UP(JI,JK+KKL))
-! A corriger pour utiliser q et non r !!!!      
-    ZMIX1(JI)=ZZDZ(JI,JK)*(PENTR(JI,JK)-PDETR(JI,JK))
-  ENDIF
+      ! A corriger pour utiliser q et non r !!!!      
+      ZMIX1(JI)=ZZDZ(JI,JK)*(PENTR(JI,JK)-PDETR(JI,JK))
+    ENDIF
   ENDDO
 
   DO JI=1,IIJU
-  IF(GTEST(JI)) THEN
-    PEMF(JI,JK+KKL)=PEMF(JI,JK)*EXP(ZMIX1(JI))
-  ENDIF
+    IF(GTEST(JI)) THEN
+      PEMF(JI,JK+KKL)=PEMF(JI,JK)*EXP(ZMIX1(JI))
+    ENDIF
   ENDDO
 
   DO JI=1,IIJU
-  IF(GTEST(JI)) THEN
-! Updraft fraction must be smaller than XFRAC_UP_MAX
-    PFRAC_UP(JI,JK+KKL)=MIN(XFRAC_UP_MAX,PEMF(JI,JK+KKL)/(SQRT(ZW_UP2(JI,JK+KKL))*ZRHO_F(JI,JK+KKL)))
-    PEMF(JI,JK+KKL) = ZRHO_F(JI,JK+KKL)*PFRAC_UP(JI,JK+KKL)*SQRT(ZW_UP2(JI,JK+KKL))
-  ENDIF
+    IF(GTEST(JI)) THEN
+      ! Updraft fraction must be smaller than XFRAC_UP_MAX
+      PFRAC_UP(JI,JK+KKL)=MIN(XFRAC_UP_MAX, &
+                             &PEMF(JI,JK+KKL)/(SQRT(ZW_UP2(JI,JK+KKL))*ZRHO_F(JI,JK+KKL)))
+      PEMF(JI,JK+KKL) = ZRHO_F(JI,JK+KKL)*PFRAC_UP(JI,JK+KKL)*SQRT(ZW_UP2(JI,JK+KKL))
+    ENDIF
   ENDDO
 
 ! Test if the updraft has reach the ETL
   DO JI=1,IIJU
-  IF (GTEST(JI).AND.(PBUO_INTEG(JI,JK)<=0.)) THEN
+    IF (GTEST(JI) .AND. (PBUO_INTEG(JI,JK)<=0.)) THEN
       KKETL(JI) = JK+KKL
-  ENDIF
+    ENDIF
   ENDDO
 
 
 ! Test is we have reached the top of the updraft
-
   DO JI=1,IIJU
-  IF (GTEST(JI).AND.((ZW_UP2(JI,JK+KKL)<=ZEPS).OR.(PEMF(JI,JK+KKL)<=ZEPS))) THEN
+    IF (GTEST(JI) .AND. ((ZW_UP2(JI,JK+KKL)<=ZEPS).OR.(PEMF(JI,JK+KKL)<=ZEPS))) THEN
       ZW_UP2   (JI,JK+KKL)=ZEPS
       PEMF     (JI,JK+KKL)=0.
       GTEST    (JI)       =.FALSE.
@@ -569,13 +576,11 @@ DO JK=KKB,KKE-KKL,KKL
       PTHV_UP  (JI,JK+KKL)=ZTHVM_F(JI,JK+KKL)
       PFRAC_UP (JI,JK+KKL)=0.
       KKCTL    (JI)       =JK+KKL
-
-  ENDIF
+    ENDIF
   ENDDO
 
 ENDDO   ! Fin de la boucle verticale 
 
-
 PW_UP(:,:)=SQRT(ZW_UP2(:,:))
 PEMF(:,KKB) =0.
 
@@ -594,14 +599,15 @@ GWORK2(:,:) = SPREAD( GWORK1(:), DIM=2, NCOPIES=IKU )
 ZCOEF(:,:) = SPREAD( (1.-(PDEPTH(:)-ZDEPTH_MAX1)/(ZDEPTH_MAX2-ZDEPTH_MAX1)), DIM=2, NCOPIES=IKU)
 ZCOEF(:,:)=MIN(MAX(ZCOEF(:,:),0.),1.)
 DO JK=1, IKU
-DO JI=1,IIJU
-IF (GWORK2(JI,JK)) THEN
-   PEMF(JI,JK)     = PEMF(JI,JK)     * ZCOEF(JI,JK) 
-   PFRAC_UP(JI,JK) = PFRAC_UP(JI,JK) * ZCOEF(JI,JK) 
-ENDIF
-ENDDO
+  DO JI=1,IIJU
+    IF (GWORK2(JI,JK)) THEN
+      PEMF(JI,JK)     = PEMF(JI,JK)     * ZCOEF(JI,JK) 
+      PFRAC_UP(JI,JK) = PFRAC_UP(JI,JK) * ZCOEF(JI,JK) 
+    ENDIF
+  ENDDO
 ENDDO
 
 IF (LHOOK) CALL DR_HOOK('COMPUTE_UPDRAFT_RHCJ10',1,ZHOOK_HANDLE)
 
 END SUBROUTINE COMPUTE_UPDRAFT_RHCJ10
+END MODULE MODE_COMPUTE_UPDRAFT_RHCJ10
diff --git a/src/arome/turb/modi_compute_updraft.F90 b/src/arome/turb/modi_compute_updraft.F90
deleted file mode 100644
index 81f2363391efe31d55677362d3bc344e90f2b82e..0000000000000000000000000000000000000000
--- a/src/arome/turb/modi_compute_updraft.F90
+++ /dev/null
@@ -1,79 +0,0 @@
-!     ######spl
-     MODULE MODI_COMPUTE_UPDRAFT
-!    ###########################
-!
-INTERFACE
-!
-!     #################################################################
-      SUBROUTINE COMPUTE_UPDRAFT(KKA,KKB,KKE,KKU,KKL, HFRAC_ICE,  &
-                                 OENTR_DETR,OMIXUV,               &
-                                 ONOMIXLG,KSV_LGBEG,KSV_LGEND,    &
-                                 PZZ,PDZZ,                        &
-                                 PSFTH,PSFRV,                     &
-                                 PPABSM,PRHODREF,PUM,PVM,PTKEM,   &
-                                 PTHM,PRVM,PTHLM,PRTM,            &
-                                 PSVM,PTHL_UP,PRT_UP,             &
-                                 PRV_UP,PRC_UP,PRI_UP,PTHV_UP,    &
-                                 PW_UP,PU_UP, PV_UP, PSV_UP,      &
-                                 PFRAC_UP,PFRAC_ICE_UP,PRSAT_UP,  &
-                                 PEMF,PDETR,PENTR,                &
-                                 PBUO_INTEG,KKLCL,KKETL,KKCTL,    &
-                                 PDEPTH)
-!     #################################################################
-!
-!*                    1.1  Declaration of Arguments
-!
-!
-!
-INTEGER,                INTENT(IN)   :: KKA          ! near ground array index
-INTEGER,                INTENT(IN)   :: KKB          ! near ground physical index
-INTEGER,                INTENT(IN)   :: KKE          ! uppest atmosphere physical index
-INTEGER,                INTENT(IN)   :: KKU          ! uppest atmosphere array index
-INTEGER,                INTENT(IN)   :: KKL          ! +1 if grid goes from ground to atmosphere top, -1 otherwise
-CHARACTER*1,            INTENT(IN)   :: HFRAC_ICE    ! partition liquid/ice scheme
-LOGICAL,                INTENT(IN) :: OENTR_DETR! flag to recompute entrainment, detrainment and mass flux
-LOGICAL,                INTENT(IN) :: OMIXUV    ! True if mixing of momentum
-LOGICAL,                INTENT(IN)   :: ONOMIXLG  ! False if mixing of lagrangian tracer
-INTEGER,                INTENT(IN)   :: KSV_LGBEG ! first index of lag. tracer
-INTEGER,                INTENT(IN)   :: KSV_LGEND ! last  index of lag. tracer
-REAL, DIMENSION(:,:), INTENT(IN)   :: PZZ       !  Height at the flux point
-REAL, DIMENSION(:,:), INTENT(IN)   :: PDZZ      !  Metrics coefficient
- 
-REAL, DIMENSION(:),   INTENT(IN)   ::  PSFTH,PSFRV
-! normal surface fluxes of theta,rv,(u,v) parallel to the orography
-!
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PPABSM     ! Pressure at t-dt
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PRHODREF   ! dry density of the
-                                                  ! reference state
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PUM        ! u mean wind
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PVM        ! v mean wind
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PTKEM      ! TKE at t-dt
-!
-REAL, DIMENSION(:,:),   INTENT(IN)   ::  PTHM           ! liquid pot. temp. at t-dt
-REAL, DIMENSION(:,:),   INTENT(IN)   ::  PRVM           ! vapor mixing ratio at t-dt
-REAL, DIMENSION(:,:),   INTENT(IN)   ::  PTHLM,PRTM     ! cons. var. at t-dt
-
-REAL, DIMENSION(:,:,:), INTENT(IN)   ::  PSVM           ! scalar var. at t-dt
-
-REAL, DIMENSION(:,:),   INTENT(OUT)  ::  PTHL_UP,PRT_UP   ! updraft properties
-REAL, DIMENSION(:,:),   INTENT(OUT)  ::  PU_UP, PV_UP     ! updraft wind components
-REAL, DIMENSION(:,:),   INTENT(INOUT)::  PRV_UP,PRC_UP, & ! updraft rv, rc
-                                         PRI_UP,PTHV_UP,& ! updraft ri, THv
-                                         PW_UP,PFRAC_UP,& ! updraft w, fraction
-                                         PFRAC_ICE_UP,&   ! liquid/solid fraction in updraft
-                                         PRSAT_UP         ! Rsat
-
-REAL, DIMENSION(:,:,:), INTENT(OUT)  ::  PSV_UP           ! updraft scalar var. 
-                                         
-REAL, DIMENSION(:,:),   INTENT(INOUT)::  PEMF,PDETR,PENTR ! Mass_flux,
-                                                          ! entrainment, detrainment
-REAL, DIMENSION(:,:),   INTENT(INOUT) :: PBUO_INTEG       ! Integrated Buoyancy 
-INTEGER, DIMENSION(:),  INTENT(INOUT)::  KKLCL,KKETL,KKCTL! LCL, ETL, CTL                                           
-REAL, DIMENSION(:),     INTENT(OUT)   :: PDEPTH           ! Deepness of cloud
-
-
-END SUBROUTINE COMPUTE_UPDRAFT
-
-END INTERFACE
-!
-END MODULE MODI_COMPUTE_UPDRAFT
diff --git a/src/arome/turb/modi_compute_updraft_raha.F90 b/src/arome/turb/modi_compute_updraft_raha.F90
deleted file mode 100644
index be726b2c80e027078c050676199474803e5f12f6..0000000000000000000000000000000000000000
--- a/src/arome/turb/modi_compute_updraft_raha.F90
+++ /dev/null
@@ -1,80 +0,0 @@
-!     ######spl
-     MODULE MODI_COMPUTE_UPDRAFT_RAHA
-!    ###########################
-!
-INTERFACE
-!
-      SUBROUTINE COMPUTE_UPDRAFT_RAHA(KKA,KKB,KKE,KKU,KKL,HFRAC_ICE,       &
-                                 OENTR_DETR,OMIXUV,               &
-                                 ONOMIXLG,KSV_LGBEG,KSV_LGEND,    &
-                                 PZZ,PDZZ,                        &
-                                 PSFTH,PSFRV,                     &
-                                 PPABSM,PRHODREF,PUM,PVM, PTKEM,  &
-                                 PEXNM,PTHM,PRVM,PTHLM,PRTM,      &
-                                 PSVM,PTHL_UP,PRT_UP,             &
-                                 PRV_UP,PRC_UP,PRI_UP,PTHV_UP,    &
-                                 PW_UP,PU_UP, PV_UP, PSV_UP,      &
-                                 PFRAC_UP,PFRAC_ICE_UP,PRSAT_UP,  &
-                                 PEMF,PDETR,PENTR,                &
-                                 PBUO_INTEG,KKLCL,KKETL,KKCTL,    &
-                                 PDEPTH     )
-!     #################################################################
-!!
-
-!*                    1.1  Declaration of Arguments
-!
-!
-!
-INTEGER,                INTENT(IN)   :: KKA          ! near ground array index
-INTEGER,                INTENT(IN)   :: KKB          ! near ground physical index
-INTEGER,                INTENT(IN)   :: KKE          ! uppest atmosphere physical index
-INTEGER,                INTENT(IN)   :: KKU          ! uppest atmosphere array index
-INTEGER,                INTENT(IN)   :: KKL          ! +1 if grid goes from ground to atmosphere top, -1 otherwise
-CHARACTER*1,            INTENT(IN)   :: HFRAC_ICE    ! partition liquid/ice scheme
-LOGICAL,                INTENT(IN) :: OENTR_DETR! flag to recompute entrainment, detrainment and mass flux
-LOGICAL,                INTENT(IN) :: OMIXUV    ! True if mixing of momentum
-LOGICAL,                INTENT(IN)   :: ONOMIXLG  ! False if mixing of lagrangian tracer
-INTEGER,                INTENT(IN)   :: KSV_LGBEG ! first index of lag. tracer
-INTEGER,                INTENT(IN)   :: KSV_LGEND ! last  index of lag. tracer
-REAL, DIMENSION(:,:), INTENT(IN)   :: PZZ       !  Height at the flux point
-REAL, DIMENSION(:,:), INTENT(IN)   :: PDZZ      !  Metrics coefficient
- 
-REAL, DIMENSION(:),   INTENT(IN)   ::  PSFTH,PSFRV
-! normal surface fluxes of theta,rv,(u,v) parallel to the orography
-!
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PPABSM     ! Pressure at t-dt
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PRHODREF   ! dry density of the
-                                                  ! reference state
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PUM        ! u mean wind
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PVM        ! v mean wind
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PTKEM      ! TKE at t-dt
-
-REAL, DIMENSION(:,:),   INTENT(IN)   ::  PEXNM       ! Exner function at t-dt
-REAL, DIMENSION(:,:),   INTENT(IN)   ::  PTHM           ! liquid pot. temp. at t-dt
-REAL, DIMENSION(:,:),   INTENT(IN)   ::  PRVM           ! vapor mixing ratio at t-dt
-REAL, DIMENSION(:,:),   INTENT(IN)   ::  PTHLM,PRTM     ! cons. var. at t-dt
-
-REAL, DIMENSION(:,:,:), INTENT(IN)   ::  PSVM           ! scalar var. at t-dt
-
-REAL, DIMENSION(:,:),   INTENT(OUT)  ::  PTHL_UP,PRT_UP   ! updraft properties
-REAL, DIMENSION(:,:),   INTENT(OUT)  ::  PU_UP, PV_UP     ! updraft wind components
-REAL, DIMENSION(:,:),   INTENT(INOUT)::  PRV_UP,PRC_UP, & ! updraft rv, rc
-                                         PRI_UP,PTHV_UP,& ! updraft ri, THv
-                                         PW_UP,PFRAC_UP,& ! updraft w, fraction
-                                         PFRAC_ICE_UP,&   ! liquid/solid fraction in updraft
-                                         PRSAT_UP         ! Rsat
-
-REAL, DIMENSION(:,:,:), INTENT(OUT)  ::  PSV_UP           ! updraft scalar var. 
-                                         
-REAL, DIMENSION(:,:),   INTENT(INOUT)::  PEMF,PDETR,PENTR ! Mass_flux,
-                                                          ! detrainment,entrainment
-REAL, DIMENSION(:,:),   INTENT(INOUT) :: PBUO_INTEG       ! Integrated Buoyancy 
-INTEGER, DIMENSION(:),  INTENT(INOUT)::  KKLCL,KKETL,KKCTL! LCL, ETL, CTL                                     
-REAL, DIMENSION(:),     INTENT(OUT)   :: PDEPTH           ! Deepness of cloud
-
-
-END SUBROUTINE COMPUTE_UPDRAFT_RAHA
-
-END INTERFACE
-!
-END MODULE MODI_COMPUTE_UPDRAFT_RAHA
diff --git a/src/arome/turb/modi_compute_updraft_rhcj10.F90 b/src/arome/turb/modi_compute_updraft_rhcj10.F90
deleted file mode 100644
index 74cfccd50036695b603181e1742c53ec670b5ad8..0000000000000000000000000000000000000000
--- a/src/arome/turb/modi_compute_updraft_rhcj10.F90
+++ /dev/null
@@ -1,80 +0,0 @@
-!     ######spl
-     MODULE MODI_COMPUTE_UPDRAFT_RHCJ10
-!    ###########################
-!
-INTERFACE
-!
-      SUBROUTINE COMPUTE_UPDRAFT_RHCJ10(KKA,KKB,KKE,KKU,KKL,HFRAC_ICE,       &
-                                 OENTR_DETR,OMIXUV,               &
-                                 ONOMIXLG,KSV_LGBEG,KSV_LGEND,    &
-                                 PZZ,PDZZ,                        &
-                                 PSFTH,PSFRV,                     &
-                                 PPABSM,PRHODREF,PUM,PVM, PTKEM,  &
-                                 PEXNM,PTHM,PRVM,PTHLM,PRTM,      &
-                                 PSVM,PTHL_UP,PRT_UP,             &
-                                 PRV_UP,PRC_UP,PRI_UP,PTHV_UP,    &
-                                 PW_UP,PU_UP, PV_UP, PSV_UP,      &
-                                 PFRAC_UP,PFRAC_ICE_UP,PRSAT_UP,  &
-                                 PEMF,PDETR,PENTR,                &
-                                 PBUO_INTEG,KKLCL,KKETL,KKCTL,    &
-                                 PDEPTH     )
-!     #################################################################
-!!
-
-!*                    1.1  Declaration of Arguments
-!
-!
-!
-INTEGER,                INTENT(IN)   :: KKA          ! near ground array index
-INTEGER,                INTENT(IN)   :: KKB          ! near ground physical index
-INTEGER,                INTENT(IN)   :: KKE          ! uppest atmosphere physical index
-INTEGER,                INTENT(IN)   :: KKU          ! uppest atmosphere array index
-INTEGER,                INTENT(IN)   :: KKL          ! +1 if grid goes from ground to atmosphere top, -1 otherwise
-CHARACTER*1,            INTENT(IN)   :: HFRAC_ICE    ! partition liquid/ice scheme
-LOGICAL,                INTENT(IN) :: OENTR_DETR! flag to recompute entrainment, detrainment and mass flux
-LOGICAL,                INTENT(IN) :: OMIXUV    ! True if mixing of momentum
-LOGICAL,                INTENT(IN)   :: ONOMIXLG  ! False if mixing of lagrangian tracer
-INTEGER,                INTENT(IN)   :: KSV_LGBEG ! first index of lag. tracer
-INTEGER,                INTENT(IN)   :: KSV_LGEND ! last  index of lag. tracer
-REAL, DIMENSION(:,:), INTENT(IN)   :: PZZ       !  Height at the flux point
-REAL, DIMENSION(:,:), INTENT(IN)   :: PDZZ      !  Metrics coefficient
- 
-REAL, DIMENSION(:),   INTENT(IN)   ::  PSFTH,PSFRV
-! normal surface fluxes of theta,rv,(u,v) parallel to the orography
-!
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PPABSM     ! Pressure at t-dt
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PRHODREF   ! dry density of the
-                                                  ! reference state
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PUM        ! u mean wind
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PVM        ! v mean wind
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PTKEM      ! TKE at t-dt
-
-REAL, DIMENSION(:,:),   INTENT(IN)   ::  PEXNM       ! Exner function at t-dt
-REAL, DIMENSION(:,:),   INTENT(IN)   ::  PTHM           ! liquid pot. temp. at t-dt
-REAL, DIMENSION(:,:),   INTENT(IN)   ::  PRVM           ! vapor mixing ratio at t-dt
-REAL, DIMENSION(:,:),   INTENT(IN)   ::  PTHLM,PRTM     ! cons. var. at t-dt
-
-REAL, DIMENSION(:,:,:), INTENT(IN)   ::  PSVM           ! scalar var. at t-dt
-
-REAL, DIMENSION(:,:),   INTENT(OUT)  ::  PTHL_UP,PRT_UP   ! updraft properties
-REAL, DIMENSION(:,:),   INTENT(OUT)  ::  PU_UP, PV_UP     ! updraft wind components
-REAL, DIMENSION(:,:),   INTENT(INOUT)::  PRV_UP,PRC_UP, & ! updraft rv, rc
-                                         PRI_UP,PTHV_UP,& ! updraft ri, THv
-                                         PW_UP,PFRAC_UP,& ! updraft w, fraction
-                                         PFRAC_ICE_UP,&   ! liquid/solid fraction in updraft
-                                         PRSAT_UP         ! Rsat
-
-REAL, DIMENSION(:,:,:), INTENT(OUT)  ::  PSV_UP           ! updraft scalar var. 
-                                         
-REAL, DIMENSION(:,:),   INTENT(INOUT)::  PEMF,PDETR,PENTR ! Mass_flux,
-                                                          ! detrainment,entrainment
-REAL, DIMENSION(:,:),   INTENT(INOUT) :: PBUO_INTEG       ! Integrated Buoyancy 
-INTEGER, DIMENSION(:),  INTENT(INOUT)::  KKLCL,KKETL,KKCTL! LCL, ETL, CTL                                     
-REAL, DIMENSION(:),     INTENT(OUT)   :: PDEPTH           ! Deepness of cloud
-
-
-END SUBROUTINE COMPUTE_UPDRAFT_RHCJ10
-
-END INTERFACE
-!
-END MODULE MODI_COMPUTE_UPDRAFT_RHCJ10
diff --git a/src/arome/turb/shallow_mf.F90 b/src/arome/turb/shallow_mf.F90
index 0e2ceb3adea8159d0e29d4529c13e088dc217c80..e5c5c59b075bc132cea1eacedf172110933a0d06 100644
--- a/src/arome/turb/shallow_mf.F90
+++ b/src/arome/turb/shallow_mf.F90
@@ -64,9 +64,9 @@ USE MODD_PARAMETERS, ONLY: JPVEXT
 USE MODD_PARAM_MFSHALL_n
 
 USE MODE_THL_RT_FROM_TH_R_MF, ONLY: THL_RT_FROM_TH_R_MF
-USE MODI_COMPUTE_UPDRAFT
-USE MODI_COMPUTE_UPDRAFT_RHCJ10
-USE MODI_COMPUTE_UPDRAFT_RAHA
+USE MODE_COMPUTE_UPDRAFT, ONLY: COMPUTE_UPDRAFT
+USE MODE_COMPUTE_UPDRAFT_RHCJ10, ONLY: COMPUTE_UPDRAFT_RHCJ10
+USE MODE_COMPUTE_UPDRAFT_RAHA, ONLY: COMPUTE_UPDRAFT_RAHA
 USE MODE_MF_TURB, ONLY: MF_TURB
 USE MODE_MF_TURB_EXPL, ONLY: MF_TURB_EXPL
 USE MODE_COMPUTE_MF_CLOUD, ONLY: COMPUTE_MF_CLOUD
@@ -167,6 +167,10 @@ INTEGER, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) :: IERR
 !!! 1. Initialisation
 
 REAL(KIND=JPRB) :: ZHOOK_HANDLE
+
+REAL :: PDX=0., PDY=0.
+
+
 IF (LHOOK) CALL DR_HOOK('SHALLOW_MF',0,ZHOOK_HANDLE)
 
 ! vertical boundaries
@@ -211,7 +215,8 @@ IF (HMF_UPDRAFT == 'EDKF') THEN
                        PTHL_UP,PRT_UP,PRV_UP,PRC_UP,PRI_UP,      &
                        PTHV_UP, PW_UP, PU_UP, PV_UP, ZSV_UP,     &
                        PFRAC_UP,ZFRAC_ICE_UP,ZRSAT_UP,PEMF,PDETR,&
-                       PENTR,ZBUO_INTEG,KKLCL,KKETL,KKCTL,ZDEPTH )
+                       PENTR,ZBUO_INTEG,KKLCL,KKETL,KKCTL,ZDEPTH,&
+                       PDX,PDY)
 ELSEIF (HMF_UPDRAFT == 'RHCJ') THEN
    CALL COMPUTE_UPDRAFT_RHCJ10(KKA,IKB,IKE,KKU,KKL,HFRAC_ICE,GENTR_DETR,OMIXUV,&
                        ONOMIXLG,KSV_LGBEG,KSV_LGEND,             &
diff --git a/src/common/turb/mode_compute_entr_detr.F90 b/src/common/turb/mode_compute_entr_detr.F90
index dd96237eeff3aebebba6565cb9755c464a25e281..35bf18b95a01ef1e0284eaf4adef0471d7a427c6 100644
--- a/src/common/turb/mode_compute_entr_detr.F90
+++ b/src/common/turb/mode_compute_entr_detr.F90
@@ -180,8 +180,6 @@ REAL(KIND=JPRB) :: ZHOOK_HANDLE
   ZFRAC_ICE(:)=PFRAC_ICE(:) ! to not modify fraction of ice
  
   ZPRE(:)=PPRE_MINUS_HALF(:)
-  ZMIXTHL(:)=0.1
-  ZMIXRT(:)=0.1
 
 !                1.4 Estimation of PPART_DRY
   DO JLOOP=1,SIZE(OTEST)
@@ -350,6 +348,9 @@ DO JLOOP=1,SIZE(OTEST)
     ZMIXRT(JLOOP)  = ZKIC_INIT * &
                (PRTM(JLOOP,KK)-ZDZ*(PRTM(JLOOP,KK)-PRTM(JLOOP,JI))/PDZZ(JLOOP,KK)) + &
                (1. - ZKIC_INIT)*PRT_UP(JLOOP)
+  ELSE
+    ZMIXTHL(JLOOP) = 300.
+    ZMIXRT(JLOOP) = 0.1
   ENDIF
 ENDDO
   CALL TH_R_FROM_THL_RT_1D(HFRAC_ICE,ZFRAC_ICE,&
diff --git a/src/arome/turb/compute_updraft.F90 b/src/common/turb/mode_compute_updraft.F90
similarity index 90%
rename from src/arome/turb/compute_updraft.F90
rename to src/common/turb/mode_compute_updraft.F90
index b60af83a454713702cbbc588d6cbb7ee37d4cc6b..cc5c9672a6c586a2c08289d6384b899286fe210a 100644
--- a/src/arome/turb/compute_updraft.F90
+++ b/src/common/turb/mode_compute_updraft.F90
@@ -1,4 +1,14 @@
+!MNH_LIC Copyright 2004-2019 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence
+!MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt
+!MNH_LIC for details. version 1.
+!-----------------------------------------------------------------
 !     ######spl
+     MODULE MODE_COMPUTE_UPDRAFT
+!    ###########################
+!
+IMPLICIT NONE
+CONTAINS
       SUBROUTINE COMPUTE_UPDRAFT(KKA,KKB,KKE,KKU,KKL,HFRAC_ICE,   &
                                  OENTR_DETR,OMIXUV,               &
                                  ONOMIXLG,KSV_LGBEG,KSV_LGEND,    &
@@ -12,10 +22,8 @@
                                  PFRAC_UP,PFRAC_ICE_UP,PRSAT_UP,  &
                                  PEMF,PDETR,PENTR,                &
                                  PBUO_INTEG,KKLCL,KKETL,KKCTL,    &
-                                 PDEPTH     )
+                                 PDEPTH, PDX, PDY     )
 
-      USE PARKIND1, ONLY : JPRB
-      USE YOMHOOK , ONLY : LHOOK, DR_HOOK
 !     #################################################################
 !!
 !!****  *COMPUTE_UPDRAFT* - calculates caracteristics of the updraft 
@@ -49,14 +57,17 @@
 !!     S. Riette Jan 2012: support for both order of vertical levels
 !!     V.Masson, C.Lac : 02/2011 : SV_UP initialized by a non-zero value
 !!     S. Riette Apr 2013: improvement of continuity at the condensation level
+!!     R.Honnert Oct 2016 : Add ZSURF and Update with AROME
 !!     Q.Rodier  01/2019 : support RM17 mixing length
+!!     R.Honnert 01/2019 : add LGZ (reduction of the mass-flux surface closure with the resolution)
 !! --------------------------------------------------------------------------
 !
 !*      0. DECLARATIONS
 !          ------------
 !
-USE MODD_CST
-USE MODD_PARAM_MFSHALL_n
+USE MODD_CST, ONLY: XG, XRV, XRD
+USE MODD_PARAM_MFSHALL_n, ONLY: LGZ, XALP_PERT, XCMF, XPRES_UV, XFRAC_UP_MAX, &
+                                XABUO, XBENTR, XENTR_DRY, XBDETR, XGZ
 USE MODD_TURB_n, ONLY : CTURBLEN
 
 USE MODE_COMPUTE_ENTR_DETR, ONLY: COMPUTE_ENTR_DETR
@@ -64,7 +75,8 @@ USE MODE_TH_R_FROM_THL_RT_1D, ONLY: TH_R_FROM_THL_RT_1D
 USE MODI_SHUMAN_MF, ONLY: MZM_MF, MZF_MF, GZ_M_W_MF
 
 USE MODE_COMPUTE_BL89_ML, ONLY: COMPUTE_BL89_ML
-
+USE PARKIND1, ONLY : JPRB
+USE YOMHOOK , ONLY : LHOOK, DR_HOOK
 
 IMPLICIT NONE
 
@@ -77,7 +89,7 @@ INTEGER,                INTENT(IN)   :: KKB          ! near ground physical inde
 INTEGER,                INTENT(IN)   :: KKE          ! uppest atmosphere physical index
 INTEGER,                INTENT(IN)   :: KKU          ! uppest atmosphere array index
 INTEGER,                INTENT(IN)   :: KKL          ! +1 if grid goes from ground to atmosphere top, -1 otherwise
-CHARACTER*1,            INTENT(IN)   :: HFRAC_ICE    ! partition liquid/ice scheme
+CHARACTER(LEN=1),       INTENT(IN)   :: HFRAC_ICE    ! partition liquid/ice scheme
 LOGICAL,                INTENT(IN) :: OENTR_DETR! flag to recompute entrainment, detrainment and mass flux
 LOGICAL,                INTENT(IN) :: OMIXUV    ! True if mixing of momentum
 LOGICAL,                INTENT(IN)   :: ONOMIXLG  ! False if mixing of lagrangian tracer
@@ -117,6 +129,7 @@ REAL, DIMENSION(:,:),   INTENT(INOUT)::  PEMF,PDETR,PENTR ! Mass_flux,
 REAL, DIMENSION(:,:),   INTENT(INOUT) :: PBUO_INTEG       ! Integrated Buoyancy 
 INTEGER, DIMENSION(:),  INTENT(INOUT) :: KKLCL,KKETL,KKCTL! LCL, ETL, CTL
 REAL, DIMENSION(:),     INTENT(OUT)   :: PDEPTH           ! Deepness of cloud
+REAL,                   INTENT(IN)    :: PDX, PDY
 !                       1.2  Declaration of local variables
 !
 !
@@ -163,7 +176,7 @@ LOGICAL                          ::  GLMIX
 LOGICAL, DIMENSION(SIZE(PTHM,1))              :: GWORK1
 LOGICAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) :: GWORK2
 
-INTEGER  :: ITEST
+INTEGER  :: ITEST, JLOOP
 
 REAL, DIMENSION(SIZE(PTHM,1)) :: ZRC_UP, ZRI_UP, ZRV_UP,&
                                  ZRSATW, ZRSATI,&
@@ -172,8 +185,11 @@ REAL, DIMENSION(SIZE(PTHM,1)) :: ZRC_UP, ZRI_UP, ZRV_UP,&
 REAL  :: ZDEPTH_MAX1, ZDEPTH_MAX2 ! control auto-extinction process
 
 REAL  :: ZTMAX,ZRMAX  ! control value
+
+REAL, DIMENSION(SIZE(PTHM,1)) :: ZSURF
 REAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) :: ZSHEAR,ZDUDZ,ZDVDZ ! vertical wind shear
 REAL(KIND=JPRB) :: ZHOOK_HANDLE
+!
 IF (LHOOK) CALL DR_HOOK('COMPUTE_UPDRAFT',0,ZHOOK_HANDLE)
 
 ! Thresholds for the  perturbation of
@@ -299,13 +315,10 @@ IF (OENTR_DETR) THEN
     ZDUDZ = MZF_MF(GZ_M_W_MF(PUM,PDZZ, KKA, KKU, KKL), KKA, KKU, KKL)
     ZDVDZ = MZF_MF(GZ_M_W_MF(PVM,PDZZ, KKA, KKU, KKL), KKA, KKU, KKL)
     ZSHEAR = SQRT(ZDUDZ*ZDUDZ + ZDVDZ*ZDVDZ)
-    PRINT*, 'phasage bete sans controle'
-    CALL ABORT
-    STOP
   ELSE
     ZSHEAR = 0. !no shear in bl89 mixing length
   END IF
- !
+  !
   CALL COMPUTE_BL89_ML(KKA,KKB,KKE,KKU,KKL,PDZZ,ZTKEM_F(:,KKB),ZG_O_THVREF(:,KKB),ZTHVM,KKB,GLMIX,.TRUE.,ZSHEAR,ZLUP)
   ZLUP(:)=MAX(ZLUP(:),1.E-10)
 
@@ -314,8 +327,14 @@ IF (OENTR_DETR) THEN
                 (0.61*ZTHM_F(:,KKB))*PSFRV(:)
 
   ! Mass flux at KKB level (updraft triggered if PSFTH>0.)
+  IF (LGZ) THEN
+    ZSURF(:)=TANH(XGZ*SQRT(PDX*PDY)/ZLUP)
+  ELSE
+    ZSURF(:)=1.
+  END IF
   WHERE (ZWTHVSURF(:)>0.)
-    PEMF(:,KKB) = XCMF * ZRHO_F(:,KKB) * ((ZG_O_THVREF(:,KKB))*ZWTHVSURF*ZLUP)**(1./3.)
+    PEMF(:,KKB) = XCMF * ZSURF(:) * ZRHO_F(:,KKB) *  &
+            ((ZG_O_THVREF(:,KKB))*ZWTHVSURF*ZLUP)**(1./3.)
     PFRAC_UP(:,KKB)=MIN(PEMF(:,KKB)/(SQRT(ZW_UP2(:,KKB))*ZRHO_F(:,KKB)),XFRAC_UP_MAX)
     ZW_UP2(:,KKB)=(PEMF(:,KKB)/(PFRAC_UP(:,KKB)*ZRHO_F(:,KKB)))**2
     GTEST(:)=.TRUE.
@@ -401,16 +420,18 @@ DO JK=KKB,KKE-KKL,KKL
 
 
 ! If the updraft did not stop, compute cons updraft characteritics at jk+KKL
-  WHERE(GTEST)     
-    ZMIX2(:) = (PZZ(:,JK+KKL)-PZZ(:,JK))*PENTR(:,JK) !&
-    ZMIX3_CLD(:) = (PZZ(:,JK+KKL)-PZZ(:,JK))*(1.-ZPART_DRY(:))*ZDETR_CLD(:,JK) !&                   
-    ZMIX2_CLD(:) = (PZZ(:,JK+KKL)-PZZ(:,JK))*(1.-ZPART_DRY(:))*ZENTR_CLD(:,JK)
-                
-    PTHL_UP(:,JK+KKL)=(PTHL_UP(:,JK)*(1.-0.5*ZMIX2(:)) + PTHLM(:,JK)*ZMIX2(:)) &
-                          /(1.+0.5*ZMIX2(:))   
-    PRT_UP(:,JK+KKL) =(PRT_UP (:,JK)*(1.-0.5*ZMIX2(:)) + PRTM(:,JK)*ZMIX2(:))  &
-                          /(1.+0.5*ZMIX2(:))
-  ENDWHERE
+  DO JLOOP=1,SIZE(GTEST)
+    IF(GTEST(JLOOP)) THEN
+      ZMIX2(JLOOP) = (PZZ(JLOOP,JK+KKL)-PZZ(JLOOP,JK))*PENTR(JLOOP,JK) !&
+      ZMIX3_CLD(JLOOP) = (PZZ(JLOOP,JK+KKL)-PZZ(JLOOP,JK))*(1.-ZPART_DRY(JLOOP))*ZDETR_CLD(JLOOP,JK) !&                   
+      ZMIX2_CLD(JLOOP) = (PZZ(JLOOP,JK+KKL)-PZZ(JLOOP,JK))*(1.-ZPART_DRY(JLOOP))*ZENTR_CLD(JLOOP,JK)
+                  
+      PTHL_UP(JLOOP,JK+KKL)=(PTHL_UP(JLOOP,JK)*(1.-0.5*ZMIX2(JLOOP)) + PTHLM(JLOOP,JK)*ZMIX2(JLOOP)) &
+                            /(1.+0.5*ZMIX2(JLOOP))   
+      PRT_UP(JLOOP,JK+KKL) =(PRT_UP (JLOOP,JK)*(1.-0.5*ZMIX2(JLOOP)) + PRTM(JLOOP,JK)*ZMIX2(JLOOP))  &
+                            /(1.+0.5*ZMIX2(JLOOP))
+    ENDIF
+  ENDDO
   
 
   IF(OMIXUV) THEN
@@ -551,3 +572,4 @@ ENDIF
 
 IF (LHOOK) CALL DR_HOOK('COMPUTE_UPDRAFT',1,ZHOOK_HANDLE)
 END SUBROUTINE COMPUTE_UPDRAFT
+END MODULE MODE_COMPUTE_UPDRAFT
diff --git a/src/arome/turb/compute_updraft_raha.F90 b/src/common/turb/mode_compute_updraft_raha.F90
similarity index 97%
rename from src/arome/turb/compute_updraft_raha.F90
rename to src/common/turb/mode_compute_updraft_raha.F90
index 05fa47509543c5cd9cdc062c70c61fe052e81016..acf02a06f740b25605915b3f4dd1fd819cdd04e7 100644
--- a/src/arome/turb/compute_updraft_raha.F90
+++ b/src/common/turb/mode_compute_updraft_raha.F90
@@ -1,4 +1,14 @@
+!MNH_LIC Copyright 2012-2019 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence
+!MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt
+!MNH_LIC for details. version 1.
+!-----------------------------------------------------------------
 !     ######spl
+     MODULE MODE_COMPUTE_UPDRAFT_RAHA
+!    ###########################
+!
+IMPLICIT NONE
+CONTAINS
       SUBROUTINE COMPUTE_UPDRAFT_RAHA(KKA,KKB,KKE,KKU,KKL,HFRAC_ICE, &
                                  OENTR_DETR,OMIXUV,                  &
                                  ONOMIXLG,KSV_LGBEG,KSV_LGEND,       &
@@ -14,8 +24,6 @@
                                  PBUO_INTEG,KKLCL,KKETL,KKCTL,       &
                                  PDEPTH     )
 
-      USE PARKIND1, ONLY : JPRB
-      USE YOMHOOK , ONLY : LHOOK, DR_HOOK
 !     #################################################################
 !!
 !!****  *COMPUTE_UPDRAF_RAHA* - calculates caracteristics of the updraft 
@@ -44,7 +52,7 @@
 !!     AUTHOR
 !!     ------
 !!     Y. Bouteloup (2012)
-!!     R. Honert Janv 2013 ==> corection of some coding bugs
+!!     R. Honnert Janv 2013 ==> corection of some coding bugs
 !!     Y. Bouteloup Janv 2014 ==> Allow the use of loops in the both direction
 !! --------------------------------------------------------------------------
 !
@@ -56,6 +64,9 @@ USE MODD_PARAM_MFSHALL_n
 
 USE MODE_TH_R_FROM_THL_RT_1D, ONLY: TH_R_FROM_THL_RT_1D
 USE MODI_SHUMAN_MF, ONLY: MZM_MF
+!
+USE PARKIND1, ONLY : JPRB
+USE YOMHOOK , ONLY : LHOOK, DR_HOOK
 
 IMPLICIT NONE
 
@@ -68,7 +79,7 @@ INTEGER,                INTENT(IN)   :: KKB          ! near ground physical inde
 INTEGER,                INTENT(IN)   :: KKE          ! uppest atmosphere physical index
 INTEGER,                INTENT(IN)   :: KKU          ! uppest atmosphere array index
 INTEGER,                INTENT(IN)   :: KKL          ! +1 if grid goes from ground to atmosphere top, -1 otherwise
-CHARACTER*1,            INTENT(IN)   :: HFRAC_ICE    ! partition liquid/ice scheme
+CHARACTER(LEN=1),       INTENT(IN)   :: HFRAC_ICE    ! partition liquid/ice scheme
 LOGICAL,                INTENT(IN) :: OENTR_DETR! flag to recompute entrainment, detrainment and mass flux
 LOGICAL,                INTENT(IN) :: OMIXUV    ! True if mixing of momentum
 LOGICAL,                INTENT(IN)   :: ONOMIXLG  ! False if mixing of lagrangian tracer
@@ -584,3 +595,4 @@ ENDWHERE
 IF (LHOOK) CALL DR_HOOK('COMPUTE_UPDRAF_RAHA',1,ZHOOK_HANDLE)
 
 END SUBROUTINE COMPUTE_UPDRAFT_RAHA
+END MODULE MODE_COMPUTE_UPDRAFT_RAHA
diff --git a/src/mesonh/turb/compute_updraft.f90 b/src/mesonh/turb/compute_updraft.f90
deleted file mode 100644
index 1e04a990bdb31b47569a9d2ab2108cb53979e459..0000000000000000000000000000000000000000
--- a/src/mesonh/turb/compute_updraft.f90
+++ /dev/null
@@ -1,647 +0,0 @@
-!MNH_LIC Copyright 2004-2019 CNRS, Meteo-France and Universite Paul Sabatier
-!MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence
-!MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt
-!MNH_LIC for details. version 1.
-!-----------------------------------------------------------------
-!     ######spl
-     MODULE MODI_COMPUTE_UPDRAFT
-!    ###########################
-!
-INTERFACE
-!
-!     #################################################################
-      SUBROUTINE COMPUTE_UPDRAFT(KKA,KKB,KKE,KKU,KKL, HFRAC_ICE,  &
-                                 OENTR_DETR,OMIXUV,               &
-                                 ONOMIXLG,KSV_LGBEG,KSV_LGEND,    &
-                                 PZZ,PDZZ,                        &
-                                 PSFTH,PSFRV,                     &
-                                 PPABSM,PRHODREF,PUM,PVM,PTKEM,   &
-                                 PTHM,PRVM,PTHLM,PRTM,            &
-                                 PSVM,PTHL_UP,PRT_UP,             &
-                                 PRV_UP,PRC_UP,PRI_UP,PTHV_UP,    &
-                                 PW_UP,PU_UP, PV_UP, PSV_UP,      &
-                                 PFRAC_UP,PFRAC_ICE_UP,PRSAT_UP,  &
-                                 PEMF,PDETR,PENTR,                &
-                                 PBUO_INTEG,KKLCL,KKETL,KKCTL,    &
-                                 PDEPTH)
-!     #################################################################
-!
-!*                    1.1  Declaration of Arguments
-!
-!
-!
-INTEGER,                INTENT(IN)   :: KKA          ! near ground array index
-INTEGER,                INTENT(IN)   :: KKB          ! near ground physical index
-INTEGER,                INTENT(IN)   :: KKE          ! uppest atmosphere physical index
-INTEGER,                INTENT(IN)   :: KKU          ! uppest atmosphere array index
-INTEGER,                INTENT(IN)   :: KKL          ! +1 if grid goes from ground to atmosphere top, -1 otherwise
-CHARACTER(len=1),       INTENT(IN)   :: HFRAC_ICE    ! partition liquid/ice scheme
-LOGICAL,                INTENT(IN) :: OENTR_DETR! flag to recompute entrainment, detrainment and mass flux
-LOGICAL,                INTENT(IN) :: OMIXUV    ! True if mixing of momentum
-LOGICAL,                INTENT(IN)   :: ONOMIXLG  ! False if mixing of lagrangian tracer
-INTEGER,                INTENT(IN)   :: KSV_LGBEG ! first index of lag. tracer
-INTEGER,                INTENT(IN)   :: KSV_LGEND ! last  index of lag. tracer
-REAL, DIMENSION(:,:), INTENT(IN)   :: PZZ       !  Height at the flux point
-REAL, DIMENSION(:,:), INTENT(IN)   :: PDZZ      !  Metrics coefficient
- 
-REAL, DIMENSION(:),   INTENT(IN)   ::  PSFTH,PSFRV
-! normal surface fluxes of theta,rv,(u,v) parallel to the orography
-!
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PPABSM     ! Pressure at t-dt
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PRHODREF   ! dry density of the
-                                                  ! reference state
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PUM        ! u mean wind
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PVM        ! v mean wind
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PTKEM      ! TKE at t-dt
-!
-REAL, DIMENSION(:,:),   INTENT(IN)   ::  PTHM           ! liquid pot. temp. at t-dt
-REAL, DIMENSION(:,:),   INTENT(IN)   ::  PRVM           ! vapor mixing ratio at t-dt
-REAL, DIMENSION(:,:),   INTENT(IN)   ::  PTHLM,PRTM     ! cons. var. at t-dt
-
-REAL, DIMENSION(:,:,:), INTENT(IN)   ::  PSVM           ! scalar var. at t-dt
-
-REAL, DIMENSION(:,:),   INTENT(OUT)  ::  PTHL_UP,PRT_UP   ! updraft properties
-REAL, DIMENSION(:,:),   INTENT(OUT)  ::  PU_UP, PV_UP     ! updraft wind components
-REAL, DIMENSION(:,:),   INTENT(INOUT)::  PRV_UP,PRC_UP, & ! updraft rv, rc
-                                         PRI_UP,PTHV_UP,& ! updraft ri, THv
-                                         PW_UP,PFRAC_UP,& ! updraft w, fraction
-                                         PFRAC_ICE_UP,&   ! liquid/solid fraction in updraft
-                                         PRSAT_UP         ! Rsat
-
-REAL, DIMENSION(:,:,:), INTENT(OUT)  ::  PSV_UP           ! updraft scalar var. 
-                                         
-REAL, DIMENSION(:,:),   INTENT(INOUT)::  PEMF,PDETR,PENTR ! Mass_flux,
-                                                          ! entrainment, detrainment
-REAL, DIMENSION(:,:),   INTENT(INOUT) :: PBUO_INTEG       ! Integrated Buoyancy 
-INTEGER, DIMENSION(:),  INTENT(INOUT)::  KKLCL,KKETL,KKCTL! LCL, ETL, CTL                                           
-REAL, DIMENSION(:),     INTENT(OUT)   :: PDEPTH           ! Deepness of cloud
-
-
-END SUBROUTINE COMPUTE_UPDRAFT
-
-END INTERFACE
-!
-END MODULE MODI_COMPUTE_UPDRAFT
-!     ######spl
-      SUBROUTINE COMPUTE_UPDRAFT(KKA,KKB,KKE,KKU,KKL,HFRAC_ICE,   &
-                                 OENTR_DETR,OMIXUV,               &
-                                 ONOMIXLG,KSV_LGBEG,KSV_LGEND,    &
-                                 PZZ,PDZZ,                        &
-                                 PSFTH,PSFRV,                     &
-                                 PPABSM,PRHODREF,PUM,PVM, PTKEM,  &
-                                 PTHM,PRVM,PTHLM,PRTM,            &
-                                 PSVM,PTHL_UP,PRT_UP,             &
-                                 PRV_UP,PRC_UP,PRI_UP,PTHV_UP,    &
-                                 PW_UP,PU_UP, PV_UP, PSV_UP,      &
-                                 PFRAC_UP,PFRAC_ICE_UP,PRSAT_UP,  &
-                                 PEMF,PDETR,PENTR,                &
-                                 PBUO_INTEG,KKLCL,KKETL,KKCTL,    &
-                                 PDEPTH     )
-
-!     #################################################################
-!!
-!!****  *COMPUTE_UPDRAFT* - calculates caracteristics of the updraft 
-!!                         
-!!
-!!    PURPOSE
-!!    -------
-!!****  The purpose of this routine is to build the updraft model 
-!!
-!
-!!**  METHOD
-!!    ------
-!!
-!!    EXTERNAL
-!!    --------
-!!      
-!!    IMPLICIT ARGUMENTS
-!!    ------------------
-!!
-!!      !!     REFERENCE
-!!     ---------
-!!       Book 1 of Meso-NH documentation (chapter Turbulence)
-!!       Soares et al. 2004 QJ
-!!
-!!     AUTHOR
-!!     ------
-!!     J.Pergaud
-!!     V.Masson : Optimization 07/2010
-!!     S. Riette : 07/2010 : modification for reproducibility  
-!!     S. Riette may 2011: ice added, interface modified
-!!     S. Riette Jan 2012: support for both order of vertical levels
-!!     V.Masson, C.Lac : 02/2011 : SV_UP initialized by a non-zero value
-!!     S. Riette Apr 2013: improvement of continuity at the condensation level
-!!     R.Honnert Oct 2016 : Add ZSURF and Update with AROME
-!!     Q.Rodier  01/2019 : support RM17 mixing length
-!!     R.Honnert 01/2019 : add LGZ (reduction of the mass-flux surface closure with the resolution)
-!! --------------------------------------------------------------------------
-!
-!*      0. DECLARATIONS
-!          ------------
-!
-USE MODD_CST
-USE MODD_PARAM_MFSHALL_n
-USE MODD_TURB_n, ONLY : CTURBLEN
-
-USE MODE_COMPUTE_ENTR_DETR, ONLY: COMPUTE_ENTR_DETR
-USE MODI_TH_R_FROM_THL_RT_1D
-USE MODI_SHUMAN_MF
-
-USE MODE_COMPUTE_BL89_ML, ONLY: COMPUTE_BL89_ML
-USE MODD_GRID_n, ONLY : XDXHAT, XDYHAT
-
-
-IMPLICIT NONE
-
-!*                    1.1  Declaration of Arguments
-!
-!
-!
-INTEGER,                INTENT(IN)   :: KKA          ! near ground array index
-INTEGER,                INTENT(IN)   :: KKB          ! near ground physical index
-INTEGER,                INTENT(IN)   :: KKE          ! uppest atmosphere physical index
-INTEGER,                INTENT(IN)   :: KKU          ! uppest atmosphere array index
-INTEGER,                INTENT(IN)   :: KKL          ! +1 if grid goes from ground to atmosphere top, -1 otherwise
-CHARACTER(len=1),       INTENT(IN)   :: HFRAC_ICE    ! partition liquid/ice scheme
-LOGICAL,                INTENT(IN) :: OENTR_DETR! flag to recompute entrainment, detrainment and mass flux
-LOGICAL,                INTENT(IN) :: OMIXUV    ! True if mixing of momentum
-LOGICAL,                INTENT(IN)   :: ONOMIXLG  ! False if mixing of lagrangian tracer
-INTEGER,                INTENT(IN)   :: KSV_LGBEG ! first index of lag. tracer
-INTEGER,                INTENT(IN)   :: KSV_LGEND ! last  index of lag. tracer
-REAL, DIMENSION(:,:), INTENT(IN)   :: PZZ       !  Height at the flux point
-REAL, DIMENSION(:,:), INTENT(IN)   :: PDZZ      !  Metrics coefficient
- 
-REAL, DIMENSION(:),   INTENT(IN)   ::  PSFTH,PSFRV
-! normal surface fluxes of theta,rv,(u,v) parallel to the orography
-!
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PPABSM     ! Pressure at t-dt
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PRHODREF   ! dry density of the
-                                                  ! reference state
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PUM        ! u mean wind
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PVM        ! v mean wind
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PTKEM      ! TKE at t-dt
-!
-REAL, DIMENSION(:,:),   INTENT(IN)   ::  PTHM           ! liquid pot. temp. at t-dt
-REAL, DIMENSION(:,:),   INTENT(IN)   ::  PRVM           ! vapor mixing ratio at t-dt
-REAL, DIMENSION(:,:),   INTENT(IN)   ::  PTHLM,PRTM     ! cons. var. at t-dt
-
-REAL, DIMENSION(:,:,:), INTENT(IN)   ::  PSVM           ! scalar var. at t-dt
-
-REAL, DIMENSION(:,:),   INTENT(OUT)  ::  PTHL_UP,PRT_UP   ! updraft properties
-REAL, DIMENSION(:,:),   INTENT(OUT)  ::  PU_UP, PV_UP     ! updraft wind components
-REAL, DIMENSION(:,:),   INTENT(INOUT)::  PRV_UP,PRC_UP, & ! updraft rv, rc
-                                         PRI_UP,PTHV_UP,& ! updraft ri, THv
-                                         PW_UP,PFRAC_UP,& ! updraft w, fraction
-                                         PFRAC_ICE_UP,&   ! liquid/solid fraction in updraft
-                                         PRSAT_UP         ! Rsat
-
-REAL, DIMENSION(:,:,:), INTENT(OUT)  ::  PSV_UP           ! updraft scalar var. 
-                                         
-REAL, DIMENSION(:,:),   INTENT(INOUT)::  PEMF,PDETR,PENTR ! Mass_flux,
-                                                          ! detrainment,entrainment
-REAL, DIMENSION(:,:),   INTENT(INOUT) :: PBUO_INTEG       ! Integrated Buoyancy 
-INTEGER, DIMENSION(:),  INTENT(INOUT) :: KKLCL,KKETL,KKCTL! LCL, ETL, CTL
-REAL, DIMENSION(:),     INTENT(OUT)   :: PDEPTH           ! Deepness of cloud
-!                       1.2  Declaration of local variables
-!
-!
-! Mean environment variables at t-dt at flux point
-REAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) ::    &
-                        ZTHM_F,ZRVM_F                 ! Theta,rv of
-                                                      ! updraft environnement
-REAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) ::    &
-                        ZRTM_F, ZTHLM_F, ZTKEM_F,&    ! rt, thetal,TKE,pressure,
-                        ZUM_F,ZVM_F,ZRHO_F,      &    ! density,momentum
-                        ZPRES_F,ZTHVM_F,ZTHVM,   &    ! interpolated at the flux point
-                        ZG_O_THVREF,             &    ! g*ThetaV ref
-                        ZW_UP2,                  &    ! w**2  of the updraft
-                        ZBUO_INTEG_DRY, ZBUO_INTEG_CLD,&! Integrated Buoyancy
-                        ZENTR_CLD,ZDETR_CLD           ! wet entrainment and detrainment
-
-REAL, DIMENSION(SIZE(PSVM,1),SIZE(PTHM,2),SIZE(PSVM,3)) :: &
-                        ZSVM_F ! scalar variables 
-
-                        
-REAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) ::  &
-                        ZTH_UP,                  &    ! updraft THETA 
-                        ZRC_MIX, ZRI_MIX              ! guess of Rc and Ri for KF mixture
-
-REAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) ::  ZCOEF  ! diminution coefficient for too high clouds 
-                        
-REAL, DIMENSION(SIZE(PSFTH,1) )            ::  ZWTHVSURF  ! Surface w'thetav'
-
-REAL  :: ZRDORV       ! RD/RV
-REAL  :: ZRVORD       ! RV/RD
-
-
-REAL, DIMENSION(SIZE(PTHM,1)) :: ZMIX1,ZMIX2,ZMIX3_CLD,ZMIX2_CLD
-
-REAL, DIMENSION(SIZE(PTHM,1)) :: ZLUP         ! Upward Mixing length from the ground
-
-INTEGER  :: ISV                ! Number of scalar variables                               
-INTEGER  :: JK,JI,JSV          ! loop counters
-
-LOGICAL, DIMENSION(SIZE(PTHM,1)) ::  GTEST,GTESTLCL,GTESTETL
-                               ! Test if the ascent continue, if LCL or ETL is reached
-LOGICAL                          ::  GLMIX 
-                               ! To choose upward or downward mixing length
-LOGICAL, DIMENSION(SIZE(PTHM,1))              :: GWORK1
-LOGICAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) :: GWORK2
-
-INTEGER  :: ITEST,JLOOP
-
-REAL, DIMENSION(SIZE(PTHM,1)) :: ZRC_UP, ZRI_UP, ZRV_UP, ZRSATW, ZRSATI,&
-                                 ZPART_DRY
-
-REAL  :: ZDEPTH_MAX1, ZDEPTH_MAX2 ! control auto-extinction process
-
-REAL  :: ZTMAX,ZRMAX  ! control value
-
-REAL, DIMENSION(SIZE(PTHM,1)) :: ZSURF
-REAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) :: ZSHEAR,ZDUDZ,ZDVDZ ! vertical wind shear
-! Thresholds for the  perturbation of
-! theta_l and r_t at the first level of the updraft
-ZTMAX=2.0
-ZRMAX=1.E-3
-!------------------------------------------------------------------------
-
-!                     INITIALISATION
-
-! Initialisation of the constants   
-ZRDORV   = XRD / XRV   !=0.622
-ZRVORD   = (XRV / XRD) 
-
-ZDEPTH_MAX1=3000. ! clouds with depth inferior to this value are keeped untouched
-ZDEPTH_MAX2=4000. ! clouds with depth superior to this value are suppressed
-
-!                 Local variables, internal domain
-!number of scalar variables
-ISV=SIZE(PSVM,3)
-
-IF (OENTR_DETR) THEN
-  ! Initialisation of intersesting Level :LCL,ETL,CTL
-  KKLCL(:)=KKE
-  KKETL(:)=KKE
-  KKCTL(:)=KKE
-
-  !
-  ! Initialisation
-  !* udraft governing variables
-  PEMF(:,:)=0.
-  PDETR(:,:)=0.
-  PENTR(:,:)=0.
-
-  ! Initialisation
-  !* updraft core variables
-  PRV_UP(:,:)=0.
-  PRC_UP(:,:)=0.
-  PRI_UP(:,:)=0.
-  PW_UP(:,:)=0.
-  ZTH_UP(:,:)=0.
-  PFRAC_UP(:,:)=0.
-  PTHV_UP(:,:)=0.
-
-  PBUO_INTEG=0.
-
-  PFRAC_ICE_UP(:,:)=0.
-  PRSAT_UP(:,:)=PRVM(:,:) ! should be initialised correctly but is (normaly) not used
-
-  !cloud/dry air mixture cloud content
-  ZRC_MIX = 0.
-  ZRI_MIX = 0.
-
-END IF
-
-! Initialisation of environment variables at t-dt
-! variables at flux level
-ZTHLM_F(:,:) = MZM_MF(KKA,KKU,KKL,PTHLM(:,:))
-ZRTM_F (:,:) = MZM_MF(KKA,KKU,KKL,PRTM(:,:))
-ZUM_F  (:,:) = MZM_MF(KKA,KKU,KKL,PUM(:,:))
-ZVM_F  (:,:) = MZM_MF(KKA,KKU,KKL,PVM(:,:))
-ZTKEM_F(:,:) = MZM_MF(KKA,KKU,KKL,PTKEM(:,:))
-
-DO JSV=1,ISV
-  IF (ONOMIXLG .AND. JSV >= KSV_LGBEG .AND. JSV<= KSV_LGEND) CYCLE
-  ZSVM_F(:,:,JSV) = MZM_MF(KKA,KKU,KKL,PSVM(:,:,JSV))
-END DO
-!                     
-!          Initialisation of updraft characteristics 
-PTHL_UP(:,:)=ZTHLM_F(:,:)
-PRT_UP(:,:)=ZRTM_F(:,:)
-PU_UP(:,:)=ZUM_F(:,:)
-PV_UP(:,:)=ZVM_F(:,:)
-PSV_UP(:,:,:)=ZSVM_F(:,:,:)
-
-
-! Computation or initialisation of updraft characteristics at the KKB level
-! thetal_up,rt_up,thetaV_up, w2,Buoyancy term and mass flux (PEMF)
-
-PTHL_UP(:,KKB)= ZTHLM_F(:,KKB)+MAX(0.,MIN(ZTMAX,(PSFTH(:)/SQRT(ZTKEM_F(:,KKB)))*XALP_PERT))
-PRT_UP(:,KKB) = ZRTM_F(:,KKB)+MAX(0.,MIN(ZRMAX,(PSFRV(:)/SQRT(ZTKEM_F(:,KKB)))*XALP_PERT)) 
-
-
-IF (OENTR_DETR) THEN
-  ZTHM_F (:,:) = MZM_MF(KKA,KKU,KKL,PTHM (:,:))
-  ZPRES_F(:,:) = MZM_MF(KKA,KKU,KKL,PPABSM(:,:))
-  ZRHO_F (:,:) = MZM_MF(KKA,KKU,KKL,PRHODREF(:,:))
-  ZRVM_F (:,:) = MZM_MF(KKA,KKU,KKL,PRVM(:,:))
-
-  ! thetav at mass and flux levels
-  ZTHVM_F(:,:)=ZTHM_F(:,:)*((1.+ZRVORD*ZRVM_F(:,:))/(1.+ZRTM_F(:,:)))
-  ZTHVM(:,:)=PTHM(:,:)*((1.+ZRVORD*PRVM(:,:))/(1.+PRTM(:,:)))
-
-  PTHV_UP(:,:)=ZTHVM_F(:,:)
-
-  ZW_UP2(:,:)=0.
-  ZW_UP2(:,KKB) = MAX(0.0001,(2./3.)*ZTKEM_F(:,KKB))
-
-
-  ! Computation of non conservative variable for the KKB level of the updraft
-  ! (all or nothing ajustement)
-  PRC_UP(:,KKB)=0.
-  PRI_UP(:,KKB)=0.
-  CALL TH_R_FROM_THL_RT_1D(HFRAC_ICE,PFRAC_ICE_UP(:,KKB),ZPRES_F(:,KKB), &
-             PTHL_UP(:,KKB),PRT_UP(:,KKB),ZTH_UP(:,KKB), &
-             PRV_UP(:,KKB),PRC_UP(:,KKB),PRI_UP(:,KKB),ZRSATW(:),ZRSATI(:))
-
-  ! compute updraft thevav and buoyancy term at KKB level
-  PTHV_UP(:,KKB) = ZTH_UP(:,KKB)*((1+ZRVORD*PRV_UP(:,KKB))/(1+PRT_UP(:,KKB)))
-  ! compute mean rsat in updraft
-  PRSAT_UP(:,KKB) = ZRSATW(:)*(1-PFRAC_ICE_UP(:,KKB)) + ZRSATI(:)*PFRAC_ICE_UP(:,KKB)
-                                                            
-  ! Closure assumption for mass flux at KKB level
-  !
-
-  ZG_O_THVREF=XG/ZTHVM_F
-
-  ! compute L_up
-  GLMIX=.TRUE.
-  ZTKEM_F(:,KKB)=0.
-  !
-  IF(CTURBLEN=='RM17') THEN
-    ZDUDZ = MZF_MF(KKA,KKU,KKL,GZ_M_W_MF(KKA,KKU,KKL,PUM,PDZZ))
-    ZDVDZ = MZF_MF(KKA,KKU,KKL,GZ_M_W_MF(KKA,KKU,KKL,PVM,PDZZ))
-    ZSHEAR = SQRT(ZDUDZ*ZDUDZ + ZDVDZ*ZDVDZ)
-  ELSE
-    ZSHEAR = 0. !no shear in bl89 mixing length
-  END IF
- ! 
-  CALL COMPUTE_BL89_ML(KKA,KKB,KKE,KKU,KKL,PDZZ,ZTKEM_F(:,KKB),ZG_O_THVREF(:,KKB),ZTHVM,KKB,GLMIX,.FALSE.,ZSHEAR,ZLUP)
-  ZLUP(:)=MAX(ZLUP(:),1.E-10)
-
-  ! Compute Buoyancy flux at the ground
-  ZWTHVSURF(:) = (ZTHVM_F(:,KKB)/ZTHM_F(:,KKB))*PSFTH(:)+     &
-                (0.61*ZTHM_F(:,KKB))*PSFRV(:)
-
-  ! Mass flux at KKB level (updraft triggered if PSFTH>0.)
-  IF (LGZ) THEN
-    ZSURF(:)=TANH(XGZ*SQRT(XDXHAT(1)*XDYHAT(1))/ZLUP)
-  ELSE
-    ZSURF(:)=1.
-  END IF
-  WHERE (ZWTHVSURF(:)>0.)
-    PEMF(:,KKB) = XCMF * ZSURF(:) * ZRHO_F(:,KKB) *  &
-            ((ZG_O_THVREF(:,KKB))*ZWTHVSURF*ZLUP)**(1./3.)
-    PFRAC_UP(:,KKB)=MIN(PEMF(:,KKB)/(SQRT(ZW_UP2(:,KKB))*ZRHO_F(:,KKB)),XFRAC_UP_MAX)
-    ZW_UP2(:,KKB)=(PEMF(:,KKB)/(PFRAC_UP(:,KKB)*ZRHO_F(:,KKB)))**2
-    GTEST(:)=.TRUE.
-  ELSEWHERE
-    PEMF(:,KKB) =0.
-    GTEST(:)=.FALSE.
-  ENDWHERE
-ELSE
-  GTEST(:)=PEMF(:,KKB+KKL)>0.
-END IF
-
-!--------------------------------------------------------------------------
-
-!                        3. Vertical ascending loop
-!                           -----------------------
-!
-! If GTEST = T the updraft starts from the KKB level and stops when GTEST becomes F
-!
-!
-GTESTLCL(:)=.FALSE.
-GTESTETL(:)=.FALSE.
-
-!       Loop on vertical level
-DO JK=KKB,KKE-KKL,KKL
-
-! IF the updraft top is reached for all column, stop the loop on levels
-  ITEST=COUNT(GTEST)
-  IF (ITEST==0) CYCLE
-
-!       Computation of entrainment and detrainment with KF90
-!       parameterization in clouds and LR01 in subcloud layer
-
-
-! to find the LCL (check if JK is LCL or not)
-
-  WHERE ((PRC_UP(:,JK)+PRI_UP(:,JK)>0.).AND.(.NOT.(GTESTLCL)))
-      KKLCL(:) = JK           
-      GTESTLCL(:)=.TRUE.
-  ENDWHERE
-
-! COMPUTE PENTR and PDETR at mass level JK
-  IF (OENTR_DETR) THEN
-    IF(JK/=KKB) THEN
-      ZRC_MIX(:,JK) = ZRC_MIX(:,JK-KKL) ! guess of Rc of mixture
-      ZRI_MIX(:,JK) = ZRI_MIX(:,JK-KKL) ! guess of Ri of mixture
-    ENDIF
-    CALL COMPUTE_ENTR_DETR(JK,KKB,KKE,KKL,GTEST,GTESTLCL,HFRAC_ICE,PFRAC_ICE_UP(:,JK),&
-                           PRHODREF(:,JK),ZPRES_F(:,JK),ZPRES_F(:,JK+KKL),&
-                           PZZ(:,:),PDZZ(:,:),ZTHVM(:,:),  &
-                           PTHLM(:,:),PRTM(:,:),ZW_UP2(:,:),ZTH_UP(:,JK),   &
-                           PTHL_UP(:,JK),PRT_UP(:,JK),ZLUP(:),         &
-                           PRC_UP(:,JK),PRI_UP(:,JK),PTHV_UP(:,JK),&
-                           PRSAT_UP(:,JK),ZRC_MIX(:,JK),ZRI_MIX(:,JK),                 &
-                           PENTR(:,JK),PDETR(:,JK),ZENTR_CLD(:,JK),ZDETR_CLD(:,JK),&
-                           ZBUO_INTEG_DRY(:,JK), ZBUO_INTEG_CLD(:,JK), &
-                           ZPART_DRY(:)   )
-    PBUO_INTEG(:,JK)=ZBUO_INTEG_DRY(:,JK)+ZBUO_INTEG_CLD(:,JK)
-
-    IF (JK==KKB) THEN
-       PDETR(:,JK)=0.
-       ZDETR_CLD(:,JK)=0.
-    ENDIF   
- 
-!       Computation of updraft characteristics at level JK+KKL
-    WHERE(GTEST)
-      ZMIX1(:)=0.5*(PZZ(:,JK+KKL)-PZZ(:,JK))*(PENTR(:,JK)-PDETR(:,JK))
-      PEMF(:,JK+KKL)=PEMF(:,JK)*EXP(2*ZMIX1(:))
-    ENDWHERE
-  ELSE
-    GTEST(:) = (PEMF(:,JK+KKL)>0.)
-  END IF 
-  
-
-! stop the updraft if MF becomes negative
-  WHERE (GTEST.AND.(PEMF(:,JK+KKL)<=0.))
-    PEMF(:,JK+KKL)=0.
-    KKCTL(:) = JK+KKL
-    GTEST(:)=.FALSE.
-    PFRAC_ICE_UP(:,JK+KKL)=PFRAC_ICE_UP(:,JK)
-    PRSAT_UP(:,JK+KKL)=PRSAT_UP(:,JK)
-  ENDWHERE
-
-
-! If the updraft did not stop, compute cons updraft characteritics at jk+KKL
-! WHERE(GTEST)     
-  DO JLOOP=1,SIZE(GTEST) 
-   IF (GTEST(JLOOP) ) THEN
-    ZMIX2(JLOOP) = (PZZ(JLOOP,JK+KKL)-PZZ(JLOOP,JK))*PENTR(JLOOP,JK) !&
-    ZMIX3_CLD(JLOOP) = (PZZ(JLOOP,JK+KKL)-PZZ(JLOOP,JK))*(1.-ZPART_DRY(JLOOP))*ZDETR_CLD(JLOOP,JK) !&                   
-    ZMIX2_CLD(JLOOP) = (PZZ(JLOOP,JK+KKL)-PZZ(JLOOP,JK))*(1.-ZPART_DRY(JLOOP))*ZENTR_CLD(JLOOP,JK)
-                
-    !PTHL_UP(JLOOP,JK+KKL)=(PTHL_UP(JLOOP,JK)*(1.-0.5*ZMIX2(JLOOP)) + PTHLM(JLOOP,JK)*ZMIX2(JLOOP)) &
-    !                      /(1.+0.5*ZMIX2(JLOOP))   
-    !PRT_UP(JLOOP,JK+KKL) =(PRT_UP (JLOOP,JK)*(1.-0.5*ZMIX2(JLOOP)) + PRTM(JLOOP,JK)*ZMIX2(JLOOP))  &
-    !                      /(1.+0.5*ZMIX2(JLOOP))
-
-    PTHL_UP(JLOOP,JK+KKL)=PTHL_UP(JLOOP,JK)*EXP(-ZMIX2(JLOOP)) + PTHLM(JLOOP,JK)*(1-EXP(-ZMIX2(JLOOP))) 
-    PRT_UP(JLOOP,JK+KKL) =PRT_UP (JLOOP,JK)*EXP(-ZMIX2(JLOOP)) +  PRTM(JLOOP,JK)*(1-EXP(-ZMIX2(JLOOP)))  
-
-   END IF
-  END DO
-! ENDWHERE
-  
-
-  IF(OMIXUV) THEN
-    IF(JK/=KKB) THEN
-      WHERE(GTEST)
-        PU_UP(:,JK+KKL) = (PU_UP (:,JK)*(1-0.5*ZMIX2(:)) + PUM(:,JK)*ZMIX2(:)+ &
-                          0.5*XPRES_UV*(PZZ(:,JK+KKL)-PZZ(:,JK))*&
-                          ((PUM(:,JK+KKL)-PUM(:,JK))/PDZZ(:,JK+KKL)+&
-                           (PUM(:,JK)-PUM(:,JK-KKL))/PDZZ(:,JK))        )   &
-                          /(1+0.5*ZMIX2(:))
-        PV_UP(:,JK+KKL) = (PV_UP (:,JK)*(1-0.5*ZMIX2(:)) + PVM(:,JK)*ZMIX2(:)+ &
-                          0.5*XPRES_UV*(PZZ(:,JK+KKL)-PZZ(:,JK))*&
-                          ((PVM(:,JK+KKL)-PVM(:,JK))/PDZZ(:,JK+KKL)+&
-                           (PVM(:,JK)-PVM(:,JK-KKL))/PDZZ(:,JK))    )   &
-                          /(1+0.5*ZMIX2(:))
-      ENDWHERE
-    ELSE
-      WHERE(GTEST)
-        PU_UP(:,JK+KKL) = (PU_UP (:,JK)*(1-0.5*ZMIX2(:)) + PUM(:,JK)*ZMIX2(:)+ &
-                          0.5*XPRES_UV*(PZZ(:,JK+KKL)-PZZ(:,JK))*&
-                          ((PUM(:,JK+KKL)-PUM(:,JK))/PDZZ(:,JK+KKL))        )   &
-                          /(1+0.5*ZMIX2(:))
-        PV_UP(:,JK+KKL) = (PV_UP (:,JK)*(1-0.5*ZMIX2(:)) + PVM(:,JK)*ZMIX2(:)+ &
-                          0.5*XPRES_UV*(PZZ(:,JK+KKL)-PZZ(:,JK))*&
-                          ((PVM(:,JK+KKL)-PVM(:,JK))/PDZZ(:,JK+KKL))    )   &
-                          /(1+0.5*ZMIX2(:))
-      ENDWHERE
-
-    ENDIF
-  ENDIF
-  DO JSV=1,ISV 
-     IF (ONOMIXLG .AND. JSV >= KSV_LGBEG .AND. JSV<= KSV_LGEND) CYCLE
-      WHERE(GTEST) 
-           PSV_UP(:,JK+KKL,JSV) = (PSV_UP (:,JK,JSV)*(1-0.5*ZMIX2(:)) + &
-                        PSVM(:,JK,JSV)*ZMIX2(:))  /(1+0.5*ZMIX2(:))
-      ENDWHERE                        
-  END DO  
-  
- IF (OENTR_DETR) THEN
-
-! Compute non cons. var. at level JK+KKL
-  ZRC_UP(:)=PRC_UP(:,JK) ! guess = level just below
-  ZRI_UP(:)=PRI_UP(:,JK) ! guess = level just below
-  CALL TH_R_FROM_THL_RT_1D(HFRAC_ICE,PFRAC_ICE_UP(:,JK+KKL),ZPRES_F(:,JK+KKL), &
-          PTHL_UP(:,JK+KKL),PRT_UP(:,JK+KKL),ZTH_UP(:,JK+KKL),              &
-          ZRV_UP(:),ZRC_UP(:),ZRI_UP(:),ZRSATW(:),ZRSATI(:))
-  WHERE(GTEST)
-    PRC_UP(:,JK+KKL)=ZRC_UP(:)
-    PRV_UP(:,JK+KKL)=ZRV_UP(:)
-    PRI_UP(:,JK+KKL)=ZRI_UP(:)
-    PRSAT_UP(:,JK+KKL) = ZRSATW(:)*(1-PFRAC_ICE_UP(:,JK+KKL)) + ZRSATI(:)*PFRAC_ICE_UP(:,JK+KKL)
-  ENDWHERE
-  
-
-! Compute the updraft theta_v, buoyancy and w**2 for level JK+KKL
-  WHERE(GTEST)
-    PTHV_UP(:,JK+KKL) = ZTH_UP(:,JK+KKL)*((1+ZRVORD*PRV_UP(:,JK+KKL))/(1+PRT_UP(:,JK+KKL)))
-    WHERE (ZBUO_INTEG_DRY(:,JK)>0.)
-      ZW_UP2(:,JK+KKL)  = ZW_UP2(:,JK) + 2.*(XABUO-XBENTR*XENTR_DRY)* ZBUO_INTEG_DRY(:,JK)
-    ELSEWHERE
-      ZW_UP2(:,JK+KKL)  = ZW_UP2(:,JK) + 2.*XABUO* ZBUO_INTEG_DRY(:,JK)
-    ENDWHERE
-    ZW_UP2(:,JK+KKL)  = ZW_UP2(:,JK+KKL)*(1.-(XBDETR*ZMIX3_CLD(:)+XBENTR*ZMIX2_CLD(:)))&
-            /(1.+(XBDETR*ZMIX3_CLD(:)+XBENTR*ZMIX2_CLD(:))) &
-            +2.*(XABUO)*ZBUO_INTEG_CLD(:,JK)/(1.+(XBDETR*ZMIX3_CLD(:)+XBENTR*ZMIX2_CLD(:)))
- ENDWHERE
-
-
-! Test if the updraft has reach the ETL
-  GTESTETL(:)=.FALSE.
-  WHERE (GTEST.AND.(PBUO_INTEG(:,JK)<=0.))
-      KKETL(:) = JK+KKL
-      GTESTETL(:)=.TRUE.
-  ENDWHERE
-
-! Test is we have reached the top of the updraft
-  WHERE (GTEST.AND.((ZW_UP2(:,JK+KKL)<=0.).OR.(PEMF(:,JK+KKL)<=0.)))
-      ZW_UP2(:,JK+KKL)=0.
-      PEMF(:,JK+KKL)=0.
-      GTEST(:)=.FALSE.
-      PTHL_UP(:,JK+KKL)=ZTHLM_F(:,JK+KKL)
-      PRT_UP(:,JK+KKL)=ZRTM_F(:,JK+KKL)
-      PRC_UP(:,JK+KKL)=0.
-      PRI_UP(:,JK+KKL)=0.
-      PRV_UP(:,JK+KKL)=0.
-      PTHV_UP(:,JK+KKL)=ZTHVM_F(:,JK+KKL)
-      PFRAC_UP(:,JK+KKL)=0.
-      KKCTL(:)=JK+KKL
-  ENDWHERE
- 
-! compute frac_up at JK+KKL
-  WHERE (GTEST)
-      PFRAC_UP(:,JK+KKL)=PEMF(:,JK+KKL)/(SQRT(ZW_UP2(:,JK+KKL))*ZRHO_F(:,JK+KKL))
-  ENDWHERE
-
-! Updraft fraction must be smaller than XFRAC_UP_MAX
-  WHERE (GTEST)
-      PFRAC_UP(:,JK+KKL)=MIN(XFRAC_UP_MAX,PFRAC_UP(:,JK+KKL))
-  ENDWHERE
-
-
-! When cloudy and non-buoyant, updraft fraction must decrease
-  
-  WHERE ((GTEST.AND.GTESTETL).AND.GTESTLCL)
-    PFRAC_UP(:,JK+KKL)=MIN(PFRAC_UP(:,JK+KKL),PFRAC_UP(:,JK))
-  ENDWHERE
-
-! Mass flux is updated with the new updraft fraction
-  
-  IF (OENTR_DETR) PEMF(:,JK+KKL)=PFRAC_UP(:,JK+KKL)*SQRT(ZW_UP2(:,JK+KKL))*ZRHO_F(:,JK+KKL)
-
- END IF
-
-ENDDO
-
-IF(OENTR_DETR) THEN
-
-  PW_UP(:,:)=SQRT(ZW_UP2(:,:))
-
-  PEMF(:,KKB) =0.
-
-! Limits the shallow convection scheme when cloud heigth is higher than 3000m.
-! To do this, mass flux is multiplied by a coefficient decreasing linearly
-! from 1 (for clouds of ZDEPTH_MAX1 m of depth) to 0 (for clouds of ZDEPTH_MAX2 m of depth).
-! This way, all MF fluxes are diminished by this amount.
-! Diagnosed cloud fraction is also multiplied by the same coefficient.
-!
-  DO JI=1,SIZE(PTHM,1) 
-     PDEPTH(JI) = MAX(0., PZZ(JI,KKCTL(JI)) -  PZZ(JI,KKLCL(JI)) )
-  END DO
-
-  GWORK1(:)= (GTESTLCL(:) .AND. (PDEPTH(:) > ZDEPTH_MAX1) )
-  GWORK2(:,:) = SPREAD( GWORK1(:), DIM=2, NCOPIES=MAX(KKU,KKA) )
-  ZCOEF(:,:) = SPREAD( (1.-(PDEPTH(:)-ZDEPTH_MAX1)/(ZDEPTH_MAX2-ZDEPTH_MAX1)), DIM=2, NCOPIES=SIZE(ZCOEF,2))
-  ZCOEF=MIN(MAX(ZCOEF,0.),1.)
-  WHERE (GWORK2) 
-    PEMF(:,:)     = PEMF(:,:)     * ZCOEF(:,:)
-    PFRAC_UP(:,:) = PFRAC_UP(:,:) * ZCOEF(:,:)
-  ENDWHERE
-ENDIF   
-END SUBROUTINE COMPUTE_UPDRAFT
diff --git a/src/mesonh/turb/compute_updraft_raha.f90 b/src/mesonh/turb/compute_updraft_raha.f90
deleted file mode 100644
index 1cf8c32b22ea103beb7c41b2f79a11e7abb549d6..0000000000000000000000000000000000000000
--- a/src/mesonh/turb/compute_updraft_raha.f90
+++ /dev/null
@@ -1,666 +0,0 @@
-!MNH_LIC Copyright 2012-2019 CNRS, Meteo-France and Universite Paul Sabatier
-!MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence
-!MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt
-!MNH_LIC for details. version 1.
-!-----------------------------------------------------------------
-!     ######spl
-     MODULE MODI_COMPUTE_UPDRAFT_RAHA
-!    ###########################
-!
-INTERFACE
-!
-!     #################################################################
-      SUBROUTINE COMPUTE_UPDRAFT_RAHA(KKA,KKB,KKE,KKU,KKL,HFRAC_ICE, &
-                                 OENTR_DETR,OMIXUV,                  &
-                                 ONOMIXLG,KSV_LGBEG,KSV_LGEND,       &
-                                 PZZ,PDZZ,                           &
-                                 PSFTH,PSFRV,                        &
-                                 PPABSM,PRHODREF,PUM,PVM, PTKEM,     &
-                                 PEXNM,PTHM,PRVM,PTHLM,PRTM,         &
-                                 PSVM,PTHL_UP,PRT_UP,                &
-                                 PRV_UP,PRC_UP,PRI_UP,PTHV_UP,       &
-                                 PW_UP,PU_UP, PV_UP, PSV_UP,         &
-                                 PFRAC_UP,PFRAC_ICE_UP,PRSAT_UP,     &
-                                 PEMF,PDETR,PENTR,                   &
-                                 PBUO_INTEG,KKLCL,KKETL,KKCTL,       &
-                                 PDEPTH     )
-!     #################################################################
-!
-!*                    1.1  Declaration of Arguments
-!
-!
-INTEGER,                INTENT(IN)   :: KKA          ! near ground array index
-INTEGER,                INTENT(IN)   :: KKB          ! near ground physical index
-INTEGER,                INTENT(IN)   :: KKE          ! uppest atmosphere physical index
-INTEGER,                INTENT(IN)   :: KKU          ! uppest atmosphere array index
-INTEGER,                INTENT(IN)   :: KKL          ! +1 if grid goes from ground to atmosphere top, -1 otherwise
-CHARACTER(len=1),       INTENT(IN)   :: HFRAC_ICE    ! partition liquid/ice scheme
-LOGICAL,                INTENT(IN) :: OENTR_DETR! flag to recompute entrainment, detrainment and mass flux
-LOGICAL,                INTENT(IN) :: OMIXUV    ! True if mixing of momentum
-LOGICAL,                INTENT(IN)   :: ONOMIXLG  ! False if mixing of lagrangian tracer
-INTEGER,                INTENT(IN)   :: KSV_LGBEG ! first index of lag. tracer
-INTEGER,                INTENT(IN)   :: KSV_LGEND ! last  index of lag. tracer
-REAL, DIMENSION(:,:), INTENT(IN)   :: PZZ       !  Height at the flux point
-REAL, DIMENSION(:,:), INTENT(IN)   :: PDZZ      !  Metrics coefficient
- 
-REAL, DIMENSION(:),   INTENT(IN)   ::  PSFTH,PSFRV
-! normal surface fluxes of theta,rv,(u,v) parallel to the orography
-!
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PPABSM     ! Pressure at t-dt
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PRHODREF   ! dry density of the
-                                                  ! reference state
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PUM        ! u mean wind
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PVM        ! v mean wind
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PTKEM      ! TKE at t-dt
-
-REAL, DIMENSION(:,:),   INTENT(IN)   ::  PEXNM       ! Exner function at t-dt
-REAL, DIMENSION(:,:),   INTENT(IN)   ::  PTHM           ! liquid pot. temp. at t-dt
-REAL, DIMENSION(:,:),   INTENT(IN)   ::  PRVM           ! vapor mixing ratio at t-dt
-REAL, DIMENSION(:,:),   INTENT(IN)   ::  PTHLM,PRTM     ! cons. var. at t-dt
-
-REAL, DIMENSION(:,:,:), INTENT(IN)   ::  PSVM           ! scalar var. at t-dt
-
-REAL, DIMENSION(:,:),   INTENT(OUT)  ::  PTHL_UP,PRT_UP   ! updraft properties
-REAL, DIMENSION(:,:),   INTENT(OUT)  ::  PU_UP, PV_UP     ! updraft wind components
-REAL, DIMENSION(:,:),   INTENT(INOUT)::  PRV_UP,PRC_UP    ! updraft rv, rc
-REAL, DIMENSION(:,:),   INTENT(INOUT)::  PRI_UP,PTHV_UP   ! updraft ri, THv
-REAL, DIMENSION(:,:),   INTENT(INOUT)::  PW_UP,PFRAC_UP   ! updraft w, fraction
-REAL, DIMENSION(:,:),   INTENT(INOUT)::  PFRAC_ICE_UP     ! liquid/solid fraction in updraft
-REAL, DIMENSION(:,:),   INTENT(INOUT)::  PRSAT_UP         ! Rsat
-
-REAL, DIMENSION(:,:,:), INTENT(OUT)  ::  PSV_UP           ! updraft scalar var. 
-                                         
-REAL, DIMENSION(:,:),   INTENT(INOUT)::  PEMF,PDETR,PENTR ! Mass_flux,
-                                                          ! detrainment,entrainment
-REAL, DIMENSION(:,:),   INTENT(INOUT) :: PBUO_INTEG       ! Integrated Buoyancy 
-INTEGER, DIMENSION(:),  INTENT(INOUT)::  KKLCL,KKETL,KKCTL! LCL, ETL, CTL                                     
-REAL, DIMENSION(:),     INTENT(OUT)   :: PDEPTH           ! Deepness of cloud
-
-
-END SUBROUTINE COMPUTE_UPDRAFT_RAHA
-
-END INTERFACE
-!
-END MODULE MODI_COMPUTE_UPDRAFT_RAHA
-!
-!     ######spl
-      SUBROUTINE COMPUTE_UPDRAFT_RAHA(KKA,KKB,KKE,KKU,KKL,HFRAC_ICE, &
-                                 OENTR_DETR,OMIXUV,                  &
-                                 ONOMIXLG,KSV_LGBEG,KSV_LGEND,       &
-                                 PZZ,PDZZ,                           &
-                                 PSFTH,PSFRV,                        &
-                                 PPABSM,PRHODREF,PUM,PVM, PTKEM,     &
-                                 PEXNM,PTHM,PRVM,PTHLM,PRTM,         &
-                                 PSVM,PTHL_UP,PRT_UP,                &
-                                 PRV_UP,PRC_UP,PRI_UP,PTHV_UP,       &
-                                 PW_UP,PU_UP, PV_UP, PSV_UP,         &
-                                 PFRAC_UP,PFRAC_ICE_UP,PRSAT_UP,     &
-                                 PEMF,PDETR,PENTR,                   &
-                                 PBUO_INTEG,KKLCL,KKETL,KKCTL,       &
-                                 PDEPTH     )
-
-!     #################################################################
-!!
-!!****  *COMPUTE_UPDRAF_RAHA* - calculates caracteristics of the updraft 
-!!                         
-!!
-!!    PURPOSE
-!!    -------
-!!****  The purpose of this routine is to build the updraft following Rio et al (2010)
-!!      Same as compute_updraft_rhcj10 exept the use of Hourdin et al closure
-!!
-!
-!!**  METHOD
-!!    ------
-!!
-!!    EXTERNAL
-!!    --------
-!!      
-!!    IMPLICIT ARGUMENTS
-!!    ------------------
-!!
-!!      !!     REFERENCE
-!!     ---------
-!!       Rio et al (2010) (Boundary Layer Meteorol 135:469-483)
-!!       Hourdin et al (xxxx)
-!!
-!!     AUTHOR
-!!     ------
-!!     Y. Bouteloup (2012)
-!!     R. Honnert Janv 2013 ==> corection of some coding bugs
-!!     Y. Bouteloup Janv 2014 ==> Allow the use of loops in the both direction
-!! --------------------------------------------------------------------------
-!
-!*      0. DECLARATIONS
-!          ------------
-
-USE MODD_CST
-USE MODD_PARAM_MFSHALL_n
-
-USE MODI_TH_R_FROM_THL_RT_1D
-USE MODI_SHUMAN_MF
-
-IMPLICIT NONE
-
-!*                    1.1  Declaration of Arguments
-!
-!
-!
-INTEGER,                INTENT(IN)   :: KKA          ! near ground array index
-INTEGER,                INTENT(IN)   :: KKB          ! near ground physical index
-INTEGER,                INTENT(IN)   :: KKE          ! uppest atmosphere physical index
-INTEGER,                INTENT(IN)   :: KKU          ! uppest atmosphere array index
-INTEGER,                INTENT(IN)   :: KKL          ! +1 if grid goes from ground to atmosphere top, -1 otherwise
-CHARACTER(len=1),       INTENT(IN)   :: HFRAC_ICE    ! partition liquid/ice scheme
-LOGICAL,                INTENT(IN) :: OENTR_DETR! flag to recompute entrainment, detrainment and mass flux
-LOGICAL,                INTENT(IN) :: OMIXUV    ! True if mixing of momentum
-LOGICAL,                INTENT(IN)   :: ONOMIXLG  ! False if mixing of lagrangian tracer
-INTEGER,                INTENT(IN)   :: KSV_LGBEG ! first index of lag. tracer
-INTEGER,                INTENT(IN)   :: KSV_LGEND ! last  index of lag. tracer
-REAL, DIMENSION(:,:), INTENT(IN)   :: PZZ       !  Height at the flux point
-REAL, DIMENSION(:,:), INTENT(IN)   :: PDZZ      !  Metrics coefficient
- 
-REAL, DIMENSION(:),   INTENT(IN)   ::  PSFTH,PSFRV
-! normal surface fluxes of theta,rv,(u,v) parallel to the orography
-!
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PPABSM     ! Pressure at t-dt
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PRHODREF   ! dry density of the
-                                                  ! reference state
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PUM        ! u mean wind
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PVM        ! v mean wind
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PTKEM      ! TKE at t-dt
-
-REAL, DIMENSION(:,:),   INTENT(IN)   ::  PEXNM       ! Exner function at t-dt
-REAL, DIMENSION(:,:),   INTENT(IN)   ::  PTHM           ! liquid pot. temp. at t-dt
-REAL, DIMENSION(:,:),   INTENT(IN)   ::  PRVM           ! vapor mixing ratio at t-dt
-REAL, DIMENSION(:,:),   INTENT(IN)   ::  PTHLM,PRTM     ! cons. var. at t-dt
-
-REAL, DIMENSION(:,:,:), INTENT(IN)   ::  PSVM           ! scalar var. at t-dt
-
-REAL, DIMENSION(:,:),   INTENT(OUT)  ::  PTHL_UP,PRT_UP   ! updraft properties
-REAL, DIMENSION(:,:),   INTENT(OUT)  ::  PU_UP, PV_UP     ! updraft wind components
-REAL, DIMENSION(:,:),   INTENT(INOUT)::  PRV_UP,PRC_UP    ! updraft rv, rc
-REAL, DIMENSION(:,:),   INTENT(INOUT)::  PRI_UP,PTHV_UP   ! updraft ri, THv
-REAL, DIMENSION(:,:),   INTENT(INOUT)::  PW_UP,PFRAC_UP   ! updraft w, fraction
-REAL, DIMENSION(:,:),   INTENT(INOUT)::  PFRAC_ICE_UP     ! liquid/solid fraction in updraft
-REAL, DIMENSION(:,:),   INTENT(INOUT)::  PRSAT_UP         ! Rsat
-
-REAL, DIMENSION(:,:,:), INTENT(OUT)  ::  PSV_UP           ! updraft scalar var. 
-                                         
-REAL, DIMENSION(:,:),   INTENT(INOUT)::  PEMF,PDETR,PENTR ! Mass_flux,
-                                                          ! detrainment,entrainment
-REAL, DIMENSION(:,:),   INTENT(INOUT) :: PBUO_INTEG       ! Integrated Buoyancy 
-INTEGER, DIMENSION(:),  INTENT(INOUT)::  KKLCL,KKETL,KKCTL! LCL, ETL, CTL                                     
-REAL, DIMENSION(:),     INTENT(OUT)   :: PDEPTH           ! Deepness of cloud
-!                       1.2  Declaration of local variables
-!
-!
-! Mean environment variables at t-dt at flux point
-REAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) ::    ZTHM_F,ZRVM_F,ZRCM_F    ! Theta,rv of
-                                                                  ! updraft environnement
-REAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) :: ZRTM_F, ZTHLM_F, ZTKEM_F      ! rt, thetal,TKE,pressure,
-REAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) :: ZUM_F,ZVM_F,ZRHO_F            ! density,momentum
-REAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) :: ZPRES_F,ZTHVM_F,ZTHVM         ! interpolated at the flux point
-REAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) :: ZG_O_THVREF                   ! g*ThetaV ref
-REAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) :: ZW_UP2                        ! w**2  of the updraft
-
-REAL, DIMENSION(SIZE(PSVM,1),SIZE(PTHM,2),SIZE(PSVM,3)) :: ZSVM_F ! scalar variables 
-                        
-
-                        
-REAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) :: ZTH_UP                        ! updraft THETA 
-REAL, DIMENSION(SIZE(PTHM,1))              :: ZT_UP                         ! updraft T
-REAL, DIMENSION(SIZE(PTHM,1))              :: ZLVOCPEXN                     ! updraft L
-REAL, DIMENSION(SIZE(PTHM,1))              :: ZCP                           ! updraft cp
-REAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) :: ZBUO                          ! Buoyancy 
-REAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) :: ZTHS_UP,ZTHSM
-
-REAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) ::  ZCOEF  ! diminution coefficient for too high clouds 
-                        
-REAL, DIMENSION(SIZE(PSFTH,1) )            ::  ZWTHVSURF  ! Surface w'thetav'
-
-REAL  :: ZRDORV       ! RD/RV
-REAL  :: ZRVORD       ! RV/RD
-
-
-REAL, DIMENSION(SIZE(PTHM,1)) :: ZMIX1,ZMIX2,ZMIX3
-
-REAL, DIMENSION(SIZE(PTHM,1)) :: ZLUP         ! Upward Mixing length from the ground
-
-REAL, DIMENSION(SIZE(PTHM,1)) :: ZDEPTH       ! Deepness limit for cloud
-
-INTEGER  :: ISV                ! Number of scalar variables                               
-INTEGER  :: IKU,IIJU           ! array size in k
-INTEGER  :: JK,JI,JJ,JSV          ! loop counters
-
-LOGICAL, DIMENSION(SIZE(PTHM,1)) ::  GTEST,GTESTLCL,GTESTETL
-                               ! Test if the ascent continue, if LCL or ETL is reached
-LOGICAL                          ::  GLMIX 
-                               ! To choose upward or downward mixing length
-LOGICAL, DIMENSION(SIZE(PTHM,1))              :: GWORK1
-LOGICAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) :: GWORK2
-
-
-INTEGER  :: ITEST
-
-REAL, DIMENSION(SIZE(PTHM,1)) :: ZRC_UP, ZRI_UP, ZRV_UP, ZWP2, ZRSATW, ZRSATI
-
-LOGICAL, DIMENSION(SIZE(PTHM,1)) :: GTEST_FER
-REAL,    DIMENSION(SIZE(PTHM,1)) :: ZPHI,ZALIM_STAR_TOT
-REAL,    DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) :: ZDTHETASDZ,ZALIM_STAR,ZZDZ,ZZZ
-INTEGER, DIMENSION(SIZE(PTHM,1)) :: IALIM
-
-REAL, DIMENSION(SIZE(PTHM,1))              ::  ZTEST,ZDZ,ZWUP_MEAN    ! 
-REAL, DIMENSION(SIZE(PTHM,1))              ::  ZCOE,ZWCOE,ZBUCOE
-REAL, DIMENSION(SIZE(PTHM,1))              ::  ZDETR_BUO, ZDETR_RT
-REAL, DIMENSION(SIZE(PTHM,1))              ::  ZW_MAX               ! w**2  max of the updraft
-REAL, DIMENSION(SIZE(PTHM,1))              ::  ZZTOP                ! Top of the updraft
-REAL, DIMENSION(SIZE(PTHM,1))              ::  ZA,ZB,ZQTM,ZQT_UP
-
-REAL  :: ZDEPTH_MAX1, ZDEPTH_MAX2 ! control auto-extinction process
-
-REAL  :: ZTMAX,ZRMAX, ZEPS  ! control value
-
-
-! Thresholds for the  perturbation of
-! theta_l and r_t at the first level of the updraft
-
-ZTMAX=2.0
-ZRMAX=1.E-3
-ZEPS=1.E-15
-!------------------------------------------------------------------------
-!                     INITIALISATION
-
-! Initialisation of the constants   
-ZRDORV   = XRD / XRV   !=0.622
-ZRVORD   = (XRV / XRD) 
-
-ZDEPTH_MAX1=4500. ! clouds with depth infeRIOr to this value are keeped untouched
-ZDEPTH_MAX2=5000. ! clouds with depth superior to this value are suppressed
-
-!                 Local variables, internal domain
-! Internal Domain
-
-IKU=SIZE(PTHM,2)
-IIJU =SIZE(PTHM,1)
-!number of scalar variables
-ISV=SIZE(PSVM,3)
-
-! Initialisation of intersesting Level :LCL,ETL,CTL
-KKLCL(:)=KKE
-KKETL(:)=KKE
-KKCTL(:)=KKE
-
-! 
-! Initialisation
-!* udraft governing variables
-PEMF(:,:)=0.
-PDETR(:,:)=0.
-PENTR(:,:)=0.
-
-! Initialisation
-!* updraft core variables
-PRV_UP(:,:)=0.
-PRC_UP(:,:)=0.
-
-PW_UP(:,:)=0.
-ZTH_UP(:,:)=0.
-PFRAC_UP(:,:)=0.
-PTHV_UP(:,:)=0.
-
-PBUO_INTEG=0.
-ZBUO      =0.
-
-!no ice cloud coded yet 
-PRI_UP(:,:)=0.
-PFRAC_ICE_UP(:,:)=0.
-PRSAT_UP(:,:)=PRVM(:,:) ! should be initialised correctly but is (normaly) not used
-
-! Initialisation of environment variables at t-dt
-
-! variables at flux level
-ZTHLM_F(:,:) = MZM_MF(KKA,KKU,KKL,PTHLM(:,:))
-ZRTM_F (:,:) = MZM_MF(KKA,KKU,KKL,PRTM(:,:))
-ZUM_F  (:,:) = MZM_MF(KKA,KKU,KKL,PUM(:,:))
-ZVM_F  (:,:) = MZM_MF(KKA,KKU,KKL,PVM(:,:))
-ZTKEM_F(:,:) = MZM_MF(KKA,KKU,KKL,PTKEM(:,:))
-
-!DO JSV=1,ISV 
-! IF (ONOMIXLG .AND. JSV >= KSV_LGBEG .AND. JSV<= KSV_LGEND) CYCLE
-!   ZSVM_F(:,KKB:IKU,JSV) = 0.5*(PSVM(:,KKB:IKU,JSV)+PSVM(:,1:IKU-1,JSV))
-!   ZSVM_F(:,1,JSV)       = ZSVM_F(:,KKB,JSV) 
-!END DO  
-
-!          Initialisation of updraft characteristics 
-PTHL_UP(:,:)=ZTHLM_F(:,:)
-PRT_UP(:,:)=ZRTM_F(:,:)
-PU_UP(:,:)=ZUM_F(:,:)
-PV_UP(:,:)=ZVM_F(:,:)
-PSV_UP(:,:,:)=0.
-!IF (ONOMIXLG .AND. JSV >= KSV_LGBEG .AND. JSV<= KSV_LGEND) then
-!    PSV_UP(:,:,:)=ZSVM_F(:,:,:)
-!ENDIF
-
-! Computation or initialisation of updraft characteristics at the KKB level
-! thetal_up,rt_up,thetaV_up, w�,Buoyancy term and mass flux (PEMF)
-
-PTHL_UP(:,KKB)= ZTHLM_F(:,KKB)+MAX(0.,MIN(ZTMAX,(PSFTH(:)/SQRT(ZTKEM_F(:,KKB)))*XALP_PERT))
-PRT_UP(:,KKB) = ZRTM_F(:,KKB)+MAX(0.,MIN(ZRMAX,(PSFRV(:)/SQRT(ZTKEM_F(:,KKB)))*XALP_PERT)) 
-
-ZQT_UP(:) = PRT_UP(:,KKB)/(1.+PRT_UP(:,KKB))
-ZTHS_UP(:,KKB)=PTHL_UP(:,KKB)*(1.+XLAMBDA_MF*ZQT_UP(:))
-
-ZTHM_F (:,:) = MZM_MF(KKA,KKU,KKL,PTHM (:,:))
-ZPRES_F(:,:) = MZM_MF(KKA,KKU,KKL,PPABSM(:,:))
-ZRHO_F (:,:) = MZM_MF(KKA,KKU,KKL,PRHODREF(:,:))
-ZRVM_F (:,:) = MZM_MF(KKA,KKU,KKL,PRVM(:,:))
-
-! thetav at mass and flux levels 
-ZTHVM_F(:,:)=ZTHM_F(:,:)*((1.+ZRVORD*ZRVM_F(:,:))/(1.+ZRTM_F(:,:)))
-ZTHVM(:,:)=PTHM(:,:)*((1.+ZRVORD*PRVM(:,:))/(1.+PRTM(:,:)))
-
-PTHV_UP(:,:)= ZTHVM_F(:,:)
-PRV_UP (:,:)= ZRVM_F (:,:)
-
-ZW_UP2(:,:)=ZEPS
-ZW_UP2(:,KKB) = MAX(0.0001,(1./6.)*ZTKEM_F(:,KKB))
-GTEST = (ZW_UP2(:,KKB) > ZEPS)  
-
-! Computation of non conservative variable for the KKB level of the updraft
-! (all or nothing ajustement)
-PRC_UP(:,KKB)=0.
-PRI_UP(:,KKB)=0.
-
-CALL TH_R_FROM_THL_RT_1D(HFRAC_ICE,PFRAC_ICE_UP(:,KKB),ZPRES_F(:,KKB), &
-             PTHL_UP(:,KKB),PRT_UP(:,KKB),ZTH_UP(:,KKB), &
-             PRV_UP(:,KKB),PRC_UP(:,KKB),PRI_UP(:,KKB),ZRSATW(:),ZRSATI(:))
-
-! compute updraft thevav and buoyancy term at KKB level             
-PTHV_UP(:,KKB) = ZTH_UP(:,KKB)*((1+ZRVORD*PRV_UP(:,KKB))/(1+PRT_UP(:,KKB))) 
-! compute mean rsat in updraft
-PRSAT_UP(:,KKB) = ZRSATW(:)*(1-PFRAC_ICE_UP(:,KKB)) + ZRSATI(:)*PFRAC_ICE_UP(:,KKB)
-
-!Tout est commente pour tester dans un premier temps la s�paration en deux de la 
-!  boucle verticale, une pour w et une pour PEMF                                                            
-
-ZG_O_THVREF=XG/ZTHVM_F
-
-
-!  Definition de l'alimentation au sens de la fermeture de Hourdin et al
-
-ZALIM_STAR(:,:)   = 0.
-ZALIM_STAR_TOT(:) = 0.    ! <== Normalization of ZALIM_STAR
-IALIM(:)          = KKB   ! <== Top level of the alimentation layer
-
-DO JK=KKB,KKE-KKL,KKL   !  Vertical loop
-  ZZDZ(:,JK)   = MAX(ZEPS,PZZ(:,JK+KKL)-PZZ(:,JK))       ! <== Delta Z between two flux level
-  ZZZ(:,JK)    = MAX(0.,0.5*(PZZ(:,JK+KKL)+PZZ(:,JK)) )  ! <== Hight of mass levels
-  ZDTHETASDZ(:,JK) = (ZTHVM_F(:,JK)-ZTHVM_F(:,JK+KKL))   ! <== Delta theta_v
-  
-  WHERE ((ZTHVM_F(:,JK+KKL)<ZTHVM_F(:,JK)) .AND. (ZTHVM_F(:,KKB)>=ZTHVM_F(:,JK)))
-     ZALIM_STAR(:,JK)  = SQRT(ZZZ(:,JK))*ZDTHETASDZ(:,JK)/ZZDZ(:,JK)
-     ZALIM_STAR_TOT(:) = ZALIM_STAR_TOT(:)+ZALIM_STAR(:,JK)*ZZDZ(:,JK)
-     IALIM(:)          = JK
-  ENDWHERE
-ENDDO
-
-! Normalization of ZALIM_STAR
-DO JK=KKB,KKE-KKL,KKL   !  Vertical loop   
-   WHERE (ZALIM_STAR_TOT > ZEPS)
-      ZALIM_STAR(:,JK)  = ZALIM_STAR(:,JK)/ZALIM_STAR_TOT(:)
-   ENDWHERE
-ENDDO
-ZALIM_STAR_TOT(:) = 0.
-
-
-! --------- END of alimentation calculation  ---------------------------------------      
-  
-
-!--------------------------------------------------------------------------
-
-!                        3. Vertical ascending loop
-!                           -----------------------
-!
-! If GTEST = T the updraft starts from the KKB level and stops when GTEST becomes F
-!
-!
-GTESTLCL(:)=.FALSE.
-GTESTETL(:)=.FALSE.
-
-!       Loop on vertical level to compute W
-
-ZW_MAX(:)      = 0.
-ZZTOP(:)       = 0.
-ZPHI(:) = 0.
-
-
-DO JK=KKB,KKE-KKL,KKL
-
-! IF the updraft top is reached for all column, stop the loop on levels
-
-!  ITEST=COUNT(GTEST)
-!  IF (ITEST==0) CYCLE
-
-!       Computation of entrainment and detrainment with KF90
-!       parameterization in clouds and LR01 in subcloud layer
-
-
-! to find the LCL (check if JK is LCL or not)
-
-  WHERE ((PRC_UP(:,JK)+PRI_UP(:,JK)>0.).AND.(.NOT.(GTESTLCL)))
-      KKLCL(:) = JK           
-      GTESTLCL(:)=.TRUE.
-  ENDWHERE
-
-
-! COMPUTE PENTR and PDETR at mass level JK
-
-    
-!  Buoyancy is computed on "flux" levels where updraft variables are known
-
-  ! Compute theta_v of updraft at flux level JK    
-    
-    ZRC_UP(:)  = PRC_UP(:,JK)
-    ZRI_UP(:)  = PRI_UP(:,JK) ! guess 
-    ZRV_UP(:)  = PRV_UP(:,JK)
-    ZBUO      (:,JK) = ZG_O_THVREF(:,JK)*(PTHV_UP(:,JK) - ZTHVM_F(:,JK))
-    PBUO_INTEG(:,JK) = ZBUO(:,JK)*(PZZ(:,JK+KKL)-PZZ(:,JK))
-    
-    ZDZ(:)   = MAX(ZEPS,PZZ(:,JK+KKL)-PZZ(:,JK))
-    ZTEST(:) = XA1*ZBUO(:,JK) -  XB*ZW_UP2(:,JK)
-
-    ZCOE(:)      = ZDZ(:)
-    WHERE (ZTEST(:)>0.)
-      ZCOE(:)    = ZDZ(:)/(1.+ XBETA1)
-    ENDWHERE   
-
-!  Calcul de la vitesse
-
-    ZWCOE(:)         = (1.-XB*ZCOE(:))/(1.+XB*ZCOE(:))
-    ZBUCOE(:)        =  2.*ZCOE(:)/(1.+XB*ZCOE(:))
-    
-    ZW_UP2(:,JK+KKL) = MAX(ZEPS,ZW_UP2(:,JK)*ZWCOE(:) + XA1*ZBUO(:,JK)*ZBUCOE(:) )  
-    ZW_MAX(:) = MAX(ZW_MAX(:), SQRT(ZW_UP2(:,JK+KKL)))
-    ZWUP_MEAN(:)     = MAX(ZEPS,0.5*(ZW_UP2(:,JK+KKL)+ZW_UP2(:,JK)))
- 
-!  Entrainement et detrainement
-
-   PENTR(:,JK)  = MAX(0.,(XBETA1/(1.+XBETA1))*(XA1*ZBUO(:,JK)/ZWUP_MEAN(:)-XB))
-   
-   ZDETR_BUO(:) = MAX(0., -(XBETA1/(1.+XBETA1))*XA1*ZBUO(:,JK)/ZWUP_MEAN(:))
-   ZDETR_RT(:)  = XC*SQRT(MAX(0.,(PRT_UP(:,JK) - ZRTM_F(:,JK))) / MAX(ZEPS,ZRTM_F(:,JK)) / ZWUP_MEAN(:))
-   PDETR(:,JK)  = ZDETR_RT(:)+ZDETR_BUO(:)
-
-   
-! If the updraft did not stop, compute cons updraft characteritics at jk+1
-  WHERE(GTEST)     
-    ZZTOP(:) = MAX(ZZTOP(:),PZZ(:,JK+KKL))
-    ZMIX2(:) = (PZZ(:,JK+KKL)-PZZ(:,JK))*PENTR(:,JK) !&
-    ZMIX3(:) = (PZZ(:,JK+KKL)-PZZ(:,JK))*PDETR(:,JK) !&                   
-           
-    ZQTM(:) = PRTM(:,JK)/(1.+PRTM(:,JK))            
-    ZTHSM(:,JK) = PTHLM(:,JK)*(1.+XLAMBDA_MF*ZQTM(:))
-    ZTHS_UP(:,JK+KKL)=(ZTHS_UP(:,JK)*(1.-0.5*ZMIX2(:)) + ZTHSM(:,JK)*ZMIX2(:)) &
-                          /(1.+0.5*ZMIX2(:))   
-    PRT_UP(:,JK+KKL)=(PRT_UP (:,JK)*(1.-0.5*ZMIX2(:)) + PRTM(:,JK)*ZMIX2(:))  &
-                          /(1.+0.5*ZMIX2(:))
-    ZQT_UP(:) = PRT_UP(:,JK+KKL)/(1.+PRT_UP(:,JK+KKL))
-    PTHL_UP(:,JK+KKL)=ZTHS_UP(:,JK+KKL)/(1.+XLAMBDA_MF*ZQT_UP(:))                      
-  ENDWHERE
-  
-
-  IF(OMIXUV) THEN
-    IF(JK/=KKB) THEN
-      WHERE(GTEST)
-        PU_UP(:,JK+KKL) = (PU_UP (:,JK)*(1-0.5*ZMIX2(:)) + PUM(:,JK)*ZMIX2(:)+ &
-                          0.5*XPRES_UV*(PZZ(:,JK+KKL)-PZZ(:,JK))*&
-                          ((PUM(:,JK+KKL)-PUM(:,JK))/PDZZ(:,JK+KKL)+&
-                           (PUM(:,JK)-PUM(:,JK-KKL))/PDZZ(:,JK))        )   &
-                          /(1+0.5*ZMIX2(:))
-        PV_UP(:,JK+KKL) = (PV_UP (:,JK)*(1-0.5*ZMIX2(:)) + PVM(:,JK)*ZMIX2(:)+ &
-                          0.5*XPRES_UV*(PZZ(:,JK+KKL)-PZZ(:,JK))*&
-                          ((PVM(:,JK+KKL)-PVM(:,JK))/PDZZ(:,JK+KKL)+&
-                           (PVM(:,JK)-PVM(:,JK-KKL))/PDZZ(:,JK))    )   &
-                          /(1+0.5*ZMIX2(:))
-      ENDWHERE
-    ELSE
-      WHERE(GTEST)
-        PU_UP(:,JK+KKL) = (PU_UP (:,JK)*(1-0.5*ZMIX2(:)) + PUM(:,JK)*ZMIX2(:)+ &
-                          0.5*XPRES_UV*(PZZ(:,JK+KKL)-PZZ(:,JK))*&
-                          ((PUM(:,JK+KKL)-PUM(:,JK))/PDZZ(:,JK+KKL))        )   &
-                          /(1+0.5*ZMIX2(:))
-        PV_UP(:,JK+KKL) = (PV_UP (:,JK)*(1-0.5*ZMIX2(:)) + PVM(:,JK)*ZMIX2(:)+ &
-                          0.5*XPRES_UV*(PZZ(:,JK+KKL)-PZZ(:,JK))*&
-                          ((PVM(:,JK+KKL)-PVM(:,JK))/PDZZ(:,JK+KKL))    )   &
-                          /(1+0.5*ZMIX2(:))
-      ENDWHERE
-
-    ENDIF
-  ENDIF
-!  DO JSV=1,ISV 
-!     IF (ONOMIXLG .AND. JSV >= KSV_LGBEG .AND. JSV<= KSV_LGEND) CYCLE
-!      WHERE(GTEST) 
-!           PSV_UP(:,JK+KKL,JSV) = (PSV_UP (:,JK,JSV)*(1-0.5*ZMIX2(:)) + &
-!                        PSVM(:,JK,JSV)*ZMIX2(:))  /(1+0.5*ZMIX2(:))
-!      ENDWHERE                        
-!  ENDDO  
-  
-  
-! Compute non cons. var. at level JK+KKL
-  ZRC_UP(:)=PRC_UP(:,JK) ! guess = level just below
-  ZRI_UP(:)=PRI_UP(:,JK) ! guess = level just below
-  ZRV_UP(:)=PRV_UP(:,JK)
-  CALL TH_R_FROM_THL_RT_1D(HFRAC_ICE,PFRAC_ICE_UP(:,JK+KKL),ZPRES_F(:,JK+KKL), &
-          PTHL_UP(:,JK+KKL),PRT_UP(:,JK+KKL),ZTH_UP(:,JK+KKL),              &
-          ZRV_UP(:),ZRC_UP(:),ZRI_UP(:),ZRSATW(:),ZRSATI(:))
-  WHERE(GTEST)
-    ZT_UP(:) = ZTH_UP(:,JK+KKL)*PEXNM(:,JK+KKL)
-    ZCP(:) = XCPD + XCL * ZRC_UP(:)
-    ZLVOCPEXN(:)=(XLVTT + (XCPV-XCL) *  (ZT_UP(:)-XTT) ) / ZCP(:) / PEXNM(:,JK+KKL)
-    PRC_UP(:,JK+KKL)=MIN(0.5E-3,ZRC_UP(:))  ! On ne peut depasser 0.5 g/kg (autoconversion donc elimination !)
-    PTHL_UP(:,JK+KKL) = PTHL_UP(:,JK+KKL)+ZLVOCPEXN(:)*(ZRC_UP(:)-PRC_UP(:,JK+KKL))
-    PRV_UP(:,JK+KKL)=ZRV_UP(:)
-    PRI_UP(:,JK+KKL)=ZRI_UP(:)
-    PRT_UP(:,JK+KKL)  = PRC_UP(:,JK+KKL) + PRV_UP(:,JK+KKL)
-    PRSAT_UP(:,JK+KKL) = ZRSATW(:)*(1-PFRAC_ICE_UP(:,JK+KKL)) + ZRSATI(:)*PFRAC_ICE_UP(:,JK+KKL)
-  ENDWHERE
-  
-
-! Compute the updraft theta_v, buoyancy and w**2 for level JK+1   
- WHERE(GTEST)
-!      PTHV_UP(:,JK+KKL) = ZTH_UP(:,JK+KKL)*((1+ZRVORD*PRV_UP(:,JK+KKL))/(1+PRT_UP(:,JK+KKL)))
-      PTHV_UP(:,JK+KKL) = ZTH_UP(:,JK+KKL)*(1.+0.608*PRV_UP(:,JK+KKL) - PRC_UP(:,JK+KKL))
- ENDWHERE
-
-
-! Test if the updraft has reach the ETL
-  GTESTETL(:)=.FALSE.
-  WHERE (GTEST.AND.(PBUO_INTEG(:,JK)<=0.))
-      KKETL(:) = JK+KKL
-      GTESTETL(:)=.TRUE.
-  ENDWHERE
-
-! Test is we have reached the top of the updraft
-
-  WHERE (GTEST.AND.((ZW_UP2(:,JK+KKL)<=ZEPS)))
-      ZW_UP2(:,JK+KKL)=ZEPS
-      GTEST(:)=.FALSE.
-      PTHL_UP(:,JK+KKL)=ZTHLM_F(:,JK+KKL)
-      PRT_UP(:,JK+KKL)=ZRTM_F(:,JK+KKL)
-      PRC_UP(:,JK+KKL)=0.
-      PRI_UP(:,JK+KKL)=0.
-      PRV_UP(:,JK+KKL)=0.
-      PTHV_UP(:,JK+KKL)=ZTHVM_F(:,JK+KKL)
-      PFRAC_UP(:,JK+KKL)=0.
-      KKCTL(:)=JK+KKL
-  ENDWHERE
-
-ENDDO 
-
-! Closure assumption for mass flux at KKB+1 level (Mass flux is supposed to be 0 at KKB level !)                                                 
-!   Hourdin et al 2002 formulation
-
-
-ZZTOP(:) = MAX(ZZTOP(:),ZEPS)
-
-DO JK=KKB+KKL,KKE-KKL,KKL !  Vertical loop
-   WHERE(JK<=IALIM)
-     ZALIM_STAR_TOT(:) = ZALIM_STAR_TOT(:) + ZALIM_STAR(:,JK)*ZALIM_STAR(:,JK)*ZZDZ(:,JK)/PRHODREF(:,JK)
-   ENDWHERE  
-ENDDO   
-
-WHERE (ZALIM_STAR_TOT*ZZTOP > ZEPS)
-   ZPHI(:) =  ZW_MAX(:)/(XR*ZZTOP(:)*ZALIM_STAR_TOT(:))
-ENDWHERE      
-
-GTEST(:) = .TRUE.
-PEMF(:,KKB+KKL) = ZPHI(:)*ZZDZ(:,KKB)*ZALIM_STAR(:,KKB)
-! Updraft fraction must be smaller than XFRAC_UP_MAX
-PFRAC_UP(:,KKB+KKL)=PEMF(:,KKB+KKL)/(SQRT(ZW_UP2(:,KKB+KKL))*ZRHO_F(:,KKB+KKL))
-PFRAC_UP(:,KKB+KKL)=MIN(XFRAC_UP_MAX,PFRAC_UP(:,KKB+KKL))
-PEMF(:,KKB+KKL) = ZRHO_F(:,KKB+KKL)*PFRAC_UP(:,KKB+KKL)*SQRT(ZW_UP2(:,KKB+KKL))
-
-DO JK=KKB+KKL,KKE-KKL,KKL !  Vertical loop
-     
-   GTEST = (ZW_UP2(:,JK) > ZEPS)  
-
-  WHERE (GTEST)
-    WHERE(JK<IALIM)
-       PEMF(:,JK+KKL) = MAX(0.,PEMF(:,JK) + ZPHI(:)*ZZDZ(:,JK)*(PENTR(:,JK) - PDETR(:,JK)))
-    ELSEWHERE
-       ZMIX1(:)=ZZDZ(:,JK)*(PENTR(:,JK)-PDETR(:,JK))
-       PEMF(:,JK+KKL)=PEMF(:,JK)*EXP(ZMIX1(:))
-    ENDWHERE
-
-! Updraft fraction must be smaller than XFRAC_UP_MAX
-    PFRAC_UP(:,JK+KKL)=PEMF(:,JK+KKL)/(SQRT(ZW_UP2(:,JK+KKL))*ZRHO_F(:,JK+KKL))
-    PFRAC_UP(:,JK+KKL)=MIN(XFRAC_UP_MAX,PFRAC_UP(:,JK+KKL))
-    PEMF(:,JK+KKL) = ZRHO_F(:,JK+KKL)*PFRAC_UP(:,JK+KKL)*SQRT(ZW_UP2(:,JK+KKL))
-  ENDWHERE
-
-ENDDO
-
-PW_UP(:,:)=SQRT(ZW_UP2(:,:))
-PEMF(:,KKB) =0.
-
-! Limits the shallow convection scheme when cloud heigth is higher than 3000m.
-! To do this, mass flux is multiplied by a coefficient decreasing linearly
-! from 1 (for clouds of 3000m of depth) to 0 (for clouds of 4000m of depth).
-! This way, all MF fluxes are diminished by this amount.
-! Diagnosed cloud fraction is also multiplied by the same coefficient.
-!
-DO JI=1,SIZE(PTHM,1) 
-   PDEPTH(JI) = MAX(0., PZZ(JI,KKCTL(JI)) -  PZZ(JI,KKLCL(JI)) )
-END DO
-
-GWORK1(:)= (GTESTLCL(:) .AND. (PDEPTH(:) > ZDEPTH_MAX1) )
-GWORK2(:,:) = SPREAD( GWORK1(:), DIM=2, NCOPIES=IKU )
-ZCOEF(:,:) = SPREAD( (1.-(PDEPTH(:)-ZDEPTH_MAX1)/(ZDEPTH_MAX2-ZDEPTH_MAX1)), DIM=2, NCOPIES=IKU)
-ZCOEF=MIN(MAX(ZCOEF,0.),1.)
-WHERE (GWORK2) 
-   PEMF(:,:)     = PEMF(:,:)     * ZCOEF(:,:)
-   PFRAC_UP(:,:) = PFRAC_UP(:,:) * ZCOEF(:,:)
-ENDWHERE
-
-
-END SUBROUTINE COMPUTE_UPDRAFT_RAHA
diff --git a/src/mesonh/turb/compute_updraft_rhcj10.f90 b/src/mesonh/turb/mode_compute_updraft_rhcj10.f90
similarity index 53%
rename from src/mesonh/turb/compute_updraft_rhcj10.f90
rename to src/mesonh/turb/mode_compute_updraft_rhcj10.f90
index 9ea9748840360c03dde0305b5f68b21527c2d769..abcbdc1c09192e02db5e837a06c68f79fbb2e77a 100644
--- a/src/mesonh/turb/compute_updraft_rhcj10.f90
+++ b/src/mesonh/turb/mode_compute_updraft_rhcj10.f90
@@ -4,84 +4,11 @@
 !MNH_LIC for details. version 1.
 !-----------------------------------------------------------------
 !     ######spl
-     MODULE MODI_COMPUTE_UPDRAFT_RHCJ10
+     MODULE MODE_COMPUTE_UPDRAFT_RHCJ10
 !    ###########################
 !
-INTERFACE
-!
-!     #################################################################
-      SUBROUTINE COMPUTE_UPDRAFT_RHCJ10(KKA,KKB,KKE,KKU,KKL, HFRAC_ICE,  &
-                                 OENTR_DETR,OMIXUV,               &
-                                 ONOMIXLG,KSV_LGBEG,KSV_LGEND,    &
-                                 PZZ,PDZZ,                        &
-                                 PSFTH,PSFRV,                     &
-                                 PPABSM,PRHODREF,PUM,PVM,PTKEM,   &
-                                 PTHM,PRVM,PTHLM,PRTM,            &
-                                 PSVM,PTHL_UP,PRT_UP,             &
-                                 PRV_UP,PRC_UP,PRI_UP,PTHV_UP,    &
-                                 PW_UP,PU_UP, PV_UP, PSV_UP,      &
-                                 PFRAC_UP,PFRAC_ICE_UP,PRSAT_UP,  &
-                                 PEMF,PDETR,PENTR,                &
-                                 PBUO_INTEG,KKLCL,KKETL,KKCTL,    &
-                                 PDEPTH)
-!     #################################################################
-!
-!*                    1.1  Declaration of Arguments
-!
-!
-!
-INTEGER,                INTENT(IN)   :: KKA          ! near ground array index
-INTEGER,                INTENT(IN)   :: KKB          ! near ground physical index
-INTEGER,                INTENT(IN)   :: KKE          ! uppest atmosphere physical index
-INTEGER,                INTENT(IN)   :: KKU          ! uppest atmosphere array index
-INTEGER,                INTENT(IN)   :: KKL          ! +1 if grid goes from ground to atmosphere top, -1 otherwise
-CHARACTER(len=1),       INTENT(IN)   :: HFRAC_ICE    ! partition liquid/ice scheme
-LOGICAL,                INTENT(IN) :: OENTR_DETR! flag to recompute entrainment, detrainment and mass flux
-LOGICAL,                INTENT(IN) :: OMIXUV    ! True if mixing of momentum
-LOGICAL,                INTENT(IN)   :: ONOMIXLG  ! False if mixing of lagrangian tracer
-INTEGER,                INTENT(IN)   :: KSV_LGBEG ! first index of lag. tracer
-INTEGER,                INTENT(IN)   :: KSV_LGEND ! last  index of lag. tracer
-REAL, DIMENSION(:,:), INTENT(IN)   :: PZZ       !  Height at the flux point
-REAL, DIMENSION(:,:), INTENT(IN)   :: PDZZ      !  Metrics coefficient
- 
-REAL, DIMENSION(:),   INTENT(IN)   ::  PSFTH,PSFRV
-! normal surface fluxes of theta,rv,(u,v) parallel to the orography
-!
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PPABSM     ! Pressure at t-dt
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PRHODREF   ! dry density of the
-                                                  ! reference state
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PUM        ! u mean wind
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PVM        ! v mean wind
-REAL, DIMENSION(:,:),   INTENT(IN) ::  PTKEM      ! TKE at t-dt
-!
-REAL, DIMENSION(:,:),   INTENT(IN)   ::  PTHM           ! liquid pot. temp. at t-dt
-REAL, DIMENSION(:,:),   INTENT(IN)   ::  PRVM           ! vapor mixing ratio at t-dt
-REAL, DIMENSION(:,:),   INTENT(IN)   ::  PTHLM,PRTM     ! cons. var. at t-dt
-
-REAL, DIMENSION(:,:,:), INTENT(IN)   ::  PSVM           ! scalar var. at t-dt
-
-REAL, DIMENSION(:,:),   INTENT(OUT)  ::  PTHL_UP,PRT_UP   ! updraft properties
-REAL, DIMENSION(:,:),   INTENT(OUT)  ::  PU_UP, PV_UP     ! updraft wind components
-REAL, DIMENSION(:,:),   INTENT(INOUT)::  PRV_UP,PRC_UP, & ! updraft rv, rc
-                                         PRI_UP,PTHV_UP,& ! updraft ri, THv
-                                         PW_UP,PFRAC_UP,& ! updraft w, fraction
-                                         PFRAC_ICE_UP,&   ! liquid/solid fraction in updraft
-                                         PRSAT_UP         ! Rsat
-
-REAL, DIMENSION(:,:,:), INTENT(OUT)  ::  PSV_UP           ! updraft scalar var. 
-                                         
-REAL, DIMENSION(:,:),   INTENT(INOUT)::  PEMF,PDETR,PENTR ! Mass_flux,
-                                                          ! entrainment, detrainment
-REAL, DIMENSION(:,:),   INTENT(INOUT) :: PBUO_INTEG       ! Integrated Buoyancy 
-INTEGER, DIMENSION(:),  INTENT(INOUT)::  KKLCL,KKETL,KKCTL! LCL, ETL, CTL                                           
-REAL, DIMENSION(:),     INTENT(OUT)   :: PDEPTH           ! Deepness of cloud
-
-
-END SUBROUTINE COMPUTE_UPDRAFT_RHCJ10
-
-END INTERFACE
-!
-END MODULE MODI_COMPUTE_UPDRAFT_RHCJ10
+IMPLICIT NONE
+CONTAINS
 !
 SUBROUTINE COMPUTE_UPDRAFT_RHCJ10(KKA,KKB,KKE,KKU,KKL,HFRAC_ICE,       &
                                  OENTR_DETR,OMIXUV,               &
@@ -99,7 +26,7 @@ SUBROUTINE COMPUTE_UPDRAFT_RHCJ10(KKA,KKB,KKE,KKU,KKL,HFRAC_ICE,       &
                                  PDEPTH     )
 !     #################################################################
 !!
-!!****  *COMPUTE_UPDRAF_RHCJ10* - calculates caracteristics of the updraft 
+!!****  *COMPUTE_UPDRAFT_RHCJ10* - calculates caracteristics of the updraft 
 !!                         
 !!
 !!    PURPOSE
@@ -124,35 +51,37 @@ SUBROUTINE COMPUTE_UPDRAFT_RHCJ10(KKA,KKB,KKE,KKU,KKL,HFRAC_ICE,       &
 !!     ------
 !!     Y. Bouteloup (2012)
 !!     R. Honert Janv 2013 ==> corection of some bugs
-!!     Q.Rodier  01/2019 : support RM17 mixing length 
+!!     R. El Khatib 15-Oct-2014 Optimization
+!!     Q.Rodier  01/2019 : support RM17 mixing length
 !! --------------------------------------------------------------------------
 
 ! WARNING ==>  This updraft is not yet ready to use scalar variables 
 
 !*      0. DECLARATIONS
 !          ------------
-
+!
 USE MODD_CST
 USE MODD_PARAM_MFSHALL_n
 USE MODD_TURB_n, ONLY : CTURBLEN
-USE MODI_COMPUTE_ENTR_DETR
 USE MODI_TH_R_FROM_THL_RT_1D
-USE MODI_SHUMAN_MF
+USE MODI_SHUMAN_MF, ONLY: MZF_MF, MZM_MF, GZ_M_W_MF
 
 USE MODE_COMPUTE_BL89_ML, ONLY: COMPUTE_BL89_ML
+USE PARKIND1, ONLY : JPRB
+USE YOMHOOK , ONLY : LHOOK, DR_HOOK
 
 
 IMPLICIT NONE
 
 !*                    1.1  Declaration of Arguments
-
-
+!
+!
 INTEGER,                INTENT(IN)   :: KKA          ! near ground array index
 INTEGER,                INTENT(IN)   :: KKB          ! near ground physical index
 INTEGER,                INTENT(IN)   :: KKE          ! uppest atmosphere physical index
 INTEGER,                INTENT(IN)   :: KKU          ! uppest atmosphere array index
 INTEGER,                INTENT(IN)   :: KKL          ! +1 if grid goes from ground to atmosphere top, -1 otherwise
-CHARACTER(len=1),       INTENT(IN)   :: HFRAC_ICE    ! partition liquid/ice scheme
+CHARACTER(LEN=1),       INTENT(IN)   :: HFRAC_ICE    ! partition liquid/ice scheme
 LOGICAL,                INTENT(IN) :: OENTR_DETR! flag to recompute entrainment, detrainment and mass flux
 LOGICAL,                INTENT(IN) :: OMIXUV    ! True if mixing of momentum
 LOGICAL,                INTENT(IN)   :: ONOMIXLG  ! False if mixing of lagrangian tracer
@@ -194,13 +123,13 @@ REAL, DIMENSION(:,:),   INTENT(INOUT) :: PBUO_INTEG       ! Integrated Buoyancy
 INTEGER, DIMENSION(:),  INTENT(INOUT)::  KKLCL,KKETL,KKCTL! LCL, ETL, CTL                                     
 REAL, DIMENSION(:),     INTENT(OUT)   :: PDEPTH           ! Deepness of cloud
 !                       1.2  Declaration of local variables
-
+!
 ! Mean environment variables at t-dt at flux point
 REAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) ::    ZTHM_F,ZRVM_F    ! Theta,rv of
                                                                   ! updraft environnement
 REAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) :: ZRTM_F, ZTHLM_F, ZTKEM_F      ! rt, thetal,TKE,pressure,
 REAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) :: ZUM_F,ZVM_F,ZRHO_F            ! density,momentum
-REAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) :: ZPRES_F,ZTHVM_F,ZTHVM         ! interpolated at the flux point
+REAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) :: ZPRES_F,ZTHVM_F               ! interpolated at the flux point
 REAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) :: ZG_O_THVREF                   ! g*ThetaV ref
 REAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) :: ZW_UP2                        ! w**2  of the updraft
 
@@ -213,22 +142,20 @@ REAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) :: ZBUO                          ! Bu
 
 REAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) ::  ZCOEF  ! diminution coefficient for too high clouds 
                         
-REAL, DIMENSION(SIZE(PSFTH,1) )            ::  ZWTHVSURF  ! Surface w'thetav'
+REAL                                       ::  ZWTHVSURF  ! Surface w'thetav'
 
-REAL  :: ZRDORV       ! RD/RV
 REAL  :: ZRVORD       ! RV/RD
 
 
-REAL, DIMENSION(SIZE(PTHM,1)) :: ZMIX1,ZMIX2,ZMIX3
+REAL, DIMENSION(SIZE(PTHM,1)) :: ZMIX1,ZMIX2
 
 REAL, DIMENSION(SIZE(PTHM,1)) :: ZLUP         ! Upward Mixing length from the ground
 
-
 INTEGER  :: ISV                ! Number of scalar variables                               
 INTEGER  :: IKU,IIJU           ! array size in k
-INTEGER  :: JK,JI,JJ,JSV          ! loop counters
+INTEGER  :: JK,JI,JSV          ! loop counters
 
-LOGICAL, DIMENSION(SIZE(PTHM,1)) ::  GTEST,GTESTLCL,GTESTETL
+LOGICAL, DIMENSION(SIZE(PTHM,1)) ::  GTEST,GTESTLCL
                                ! Test if the ascent continue, if LCL or ETL is reached
 LOGICAL                          ::  GLMIX 
                                ! To choose upward or downward mixing length
@@ -239,21 +166,21 @@ INTEGER  :: ITEST
 
 REAL, DIMENSION(SIZE(PTHM,1)) :: ZRC_UP, ZRI_UP, ZRV_UP, ZRSATW, ZRSATI
 
-REAL,    DIMENSION(SIZE(PTHM,1)) :: ZPHI
-REAL,    DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) :: ZZDZ,ZZZ
+REAL,    DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) :: ZZDZ
 
 REAL, DIMENSION(SIZE(PTHM,1))              ::  ZTEST,ZDZ,ZWUP_MEAN    ! 
 REAL, DIMENSION(SIZE(PTHM,1))              ::  ZCOE,ZWCOE,ZBUCOE
 REAL, DIMENSION(SIZE(PTHM,1))              ::  ZDETR_BUO, ZDETR_RT
 REAL, DIMENSION(SIZE(PTHM,1))              ::  ZW_MAX               ! w**2  max of the updraft
 REAL, DIMENSION(SIZE(PTHM,1))              ::  ZZTOP                ! Top of the updraft
-REAL, DIMENSION(SIZE(PTHM,1))              ::  ZBETA1
 
 REAL  :: ZDEPTH_MAX1, ZDEPTH_MAX2 ! control auto-extinction process
 
 REAL  :: ZTMAX,ZRMAX, ZEPS  ! control value
 
 REAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) :: ZSHEAR,ZDUDZ,ZDVDZ ! vertical wind shear
+REAL(KIND=JPRB) :: ZHOOK_HANDLE
+IF (LHOOK) CALL DR_HOOK('COMPUTE_UPDRAFT_RHCJ10',0,ZHOOK_HANDLE)
 
 ! Thresholds for the  perturbation of
 ! theta_l and r_t at the first level of the updraft
@@ -265,7 +192,6 @@ ZEPS=1.E-15
 !                     INITIALISATION
 
 ! Initialisation of the constants   
-ZRDORV   = XRD / XRV   !=0.622
 ZRVORD   = (XRV / XRD) 
 
 ! depth are different in compute_updraft (3000. and 4000.) ==> impact is small
@@ -273,10 +199,6 @@ ZDEPTH_MAX1=4500. ! clouds with depth infeRIOr to this value are keeped untouche
 ZDEPTH_MAX2=5000. ! clouds with depth superior to this value are suppressed
 
 
-!  Initialisation of ZBETA1 ==> I do not remember why I introduced a KLON array for beta1 !
-
-ZBETA1(:) = XBETA1
-   
 !                 Local variables, internal domain
 ! Internal Domain
 
@@ -309,22 +231,27 @@ PTHV_UP(:,:)=0.
 PBUO_INTEG=0.
 ZBUO      =0.
 
+!no ice cloud coded yet
 PRI_UP(:,:)=0.
 PFRAC_ICE_UP(:,:)=0.
 PRSAT_UP(:,:)=PRVM(:,:) ! should be initialised correctly but is (normaly) not used
 
 ! Initialisation of environment variables at t-dt
-! variables at flux level
-ZTHLM_F(:,:) = MZM_MF(KKA,KKU,KKL,PTHLM(:,:))
-ZRTM_F (:,:) = MZM_MF(KKA,KKU,KKL,PRTM(:,:))
-ZUM_F  (:,:) = MZM_MF(KKA,KKU,KKL,PUM(:,:))
-ZVM_F  (:,:) = MZM_MF(KKA,KKU,KKL,PVM(:,:))
-ZTKEM_F(:,:) = MZM_MF(KKA,KKU,KKL,PTKEM(:,:))
 
+! variables at flux level
+ZTHLM_F(:,:) = MZM_MF(PTHLM(:,:), KKA, KKU, KKL)
+ZRTM_F (:,:) = MZM_MF(PRTM(:,:), KKA, KKU, KKL)
+ZUM_F  (:,:) = MZM_MF(PUM(:,:), KKA, KKU, KKL)
+ZVM_F  (:,:) = MZM_MF(PVM(:,:), KKA, KKU, KKL)
+ZTKEM_F(:,:) = MZM_MF(PTKEM(:,:), KKA, KKU, KKL)
 
-! This updraft is not yet ready to use scalar variables 
+! This updraft is not yet ready to use scalar variables
 !DO JSV=1,ISV
 !  IF (ONOMIXLG .AND. JSV >= KSV_LGBEG .AND. JSV<= KSV_LGEND) CYCLE
+! *** SR merge AROME/Méso-nh: following two lines come from the AROME version
+!   ZSVM_F(:,KKB:IKU,JSV) = 0.5*(PSVM(:,KKB:IKU,JSV)+PSVM(:,1:IKU-1,JSV))
+!   ZSVM_F(:,1,JSV)       = ZSVM_F(:,KKB,JSV) 
+! *** the following single line comes from the Meso-NH version
 !  ZSVM_F(:,:,JSV) = MZM_MF(KKA,KKU,KKL,PSVM(:,:,JSV))
 !END DO
 
@@ -334,27 +261,32 @@ PRT_UP(:,:)=ZRTM_F(:,:)
 PU_UP(:,:)=ZUM_F(:,:)
 PV_UP(:,:)=ZVM_F(:,:)
 PSV_UP(:,:,:)=0.
-! This updraft is not yet ready to use scalar variables 
+! This updraft is not yet ready to use scalar variables
 !IF (ONOMIXLG .AND. JSV >= KSV_LGBEG .AND. JSV<= KSV_LGEND) then
 !    PSV_UP(:,:,:)=ZSVM_F(:,:,:)
 !ENDIF
 
 ! Computation or initialisation of updraft characteristics at the KKB level
-! thetal_up,rt_up,thetaV_up, w�,Buoyancy term and mass flux (PEMF)
+! thetal_up,rt_up,thetaV_up, w,Buoyancy term and mass flux (PEMF)
 
-!PTHL_UP(:,KKB)= ZTHLM_F(:,KKB)+MAX(0.,MIN(ZTMAX,(PSFTH(:)/SQRT(ZTKEM_F(:,KKB)))*XALP_PERT))
-!PRT_UP(:,KKB) = ZRTM_F(:,KKB)+MAX(0.,MIN(ZRMAX,(PSFRV(:)/SQRT(ZTKEM_F(:,KKB)))*XALP_PERT)) 
-PTHL_UP(:,KKB)= ZTHLM_F(:,KKB)
-PRT_UP(:,KKB) = ZRTM_F(:,KKB)
+DO JI=1,IIJU
+  !PTHL_UP(JI,KKB)= ZTHLM_F(JI,KKB)+MAX(0.,MIN(ZTMAX,(PSFTH(JI)/SQRT(ZTKEM_F(JI,KKB)))*XALP_PERT))
+  !PRT_UP(JI,KKB) = ZRTM_F(JI,KKB)+MAX(0.,MIN(ZRMAX,(PSFRV(JI)/SQRT(ZTKEM_F(JI,KKB)))*XALP_PERT)) 
+  PTHL_UP(JI,KKB)= ZTHLM_F(JI,KKB)
+  PRT_UP(JI,KKB) = ZRTM_F(JI,KKB)
+ENDDO
 
-ZTHM_F (:,:) = MZM_MF(KKA,KKU,KKL,PTHM (:,:))
-ZPRES_F(:,:) = MZM_MF(KKA,KKU,KKL,PPABSM(:,:))
-ZRHO_F (:,:) = MZM_MF(KKA,KKU,KKL,PRHODREF(:,:))
-ZRVM_F (:,:) = MZM_MF(KKA,KKU,KKL,PRVM(:,:))
+ZTHM_F (:,:) = MZM_MF(PTHM (:,:), KKA, KKU, KKL)
+ZPRES_F(:,:) = MZM_MF(PPABSM(:,:), KKA, KKU, KKL)
+ZRHO_F (:,:) = MZM_MF(PRHODREF(:,:), KKA, KKU, KKL)
+ZRVM_F (:,:) = MZM_MF(PRVM(:,:), KKA, KKU, KKL)
 
 ! thetav at mass and flux levels 
-ZTHVM_F(:,:)=ZTHM_F(:,:)*((1.+ZRVORD*ZRVM_F(:,:))/(1.+ZRTM_F(:,:)))
-ZTHVM(:,:)=PTHM(:,:)*((1.+ZRVORD*PRVM(:,:))/(1.+PRTM(:,:)))
+DO JK=1,IKU
+  DO JI=1,IIJU
+    ZTHVM_F(JI,JK)=ZTHM_F(JI,JK)*((1.+ZRVORD*ZRVM_F(JI,JK))/(1.+ZRTM_F(JI,JK)))
+  ENDDO
+ENDDO
 
 PTHV_UP(:,:)= ZTHVM_F(:,:)
 PRV_UP (:,:)= ZRVM_F (:,:)
@@ -371,12 +303,14 @@ CALL TH_R_FROM_THL_RT_1D(HFRAC_ICE,PFRAC_ICE_UP(:,KKB),ZPRES_F(:,KKB), &
              PTHL_UP(:,KKB),PRT_UP(:,KKB),ZTH_UP(:,KKB), &
              PRV_UP(:,KKB),PRC_UP(:,KKB),PRI_UP(:,KKB),ZRSATW(:),ZRSATI(:))
 
-! compute updraft thevav and buoyancy term at KKB level             
-PTHV_UP(:,KKB) = ZTH_UP(:,KKB)*((1+ZRVORD*PRV_UP(:,KKB))/(1+PRT_UP(:,KKB))) 
-! compute mean rsat in updraft
-PRSAT_UP(:,KKB) = ZRSATW(:)*(1-PFRAC_ICE_UP(:,KKB)) + ZRSATI(:)*PFRAC_ICE_UP(:,KKB)
+DO JI=1,IIJU
+  ! compute updraft thevav and buoyancy term at KKB level             
+  PTHV_UP(JI,KKB) = ZTH_UP(JI,KKB)*((1+ZRVORD*PRV_UP(JI,KKB))/(1+PRT_UP(JI,KKB))) 
+  ! compute mean rsat in updraft
+  PRSAT_UP(JI,KKB) = ZRSATW(JI)*(1-PFRAC_ICE_UP(JI,KKB)) + ZRSATI(JI)*PFRAC_ICE_UP(JI,KKB)
+ENDDO
 
-!Tout est commente pour tester dans un premier temps la séparation en deux de la 
+!Tout est commente pour tester dans un premier temps la separation en deux de la 
 !  boucle verticale, une pour w et une pour PEMF                                                            
 
 ZG_O_THVREF=XG/ZTHVM_F
@@ -384,8 +318,9 @@ ZG_O_THVREF=XG/ZTHVM_F
 ! Calcul de la fermeture de Julien Pergaut comme limite max de PHY
 
 DO JK=KKB,KKE-KKL,KKL   !  Vertical loop
-  ZZDZ(:,JK)   = MAX(ZEPS,PZZ(:,JK+KKL)-PZZ(:,JK))  ! <== Delta Z between two flux level
-  ZZZ(:,JK)    = 0.5*(PZZ(:,JK+KKL)+PZZ(:,JK))     ! <== Hight of mass levels
+  DO JI=1,IIJU
+    ZZDZ(JI,JK)   = MAX(ZEPS,PZZ(JI,JK+KKL)-PZZ(JI,JK))  ! <== Delta Z between two flux level
+  ENDDO
 ENDDO
 
 ! compute L_up
@@ -393,8 +328,8 @@ GLMIX=.TRUE.
 ZTKEM_F(:,KKB)=0.
 !
 IF(CTURBLEN=='RM17') THEN
-  ZDUDZ = MZF_MF(KKA,KKU,KKL,GZ_M_W_MF(KKA,KKU,KKL,PUM,PDZZ))
-  ZDVDZ = MZF_MF(KKA,KKU,KKL,GZ_M_W_MF(KKA,KKU,KKL,PVM,PDZZ))
+  ZDUDZ = MZF_MF(GZ_M_W_MF(PUM,PDZZ, KKA, KKU, KKL), KKA, KKU, KKL)
+  ZDVDZ = MZF_MF(GZ_M_W_MF(PVM,PDZZ, KKA, KKU, KKL), KKA, KKU, KKL)
   ZSHEAR = SQRT(ZDUDZ*ZDUDZ + ZDVDZ*ZDVDZ)
 ELSE
   ZSHEAR = 0. !no shear in bl89 mixing length
@@ -404,20 +339,23 @@ CALL COMPUTE_BL89_ML(KKA,KKB,KKE,KKU,KKL,PDZZ,ZTKEM_F(:,KKB),ZG_O_THVREF(:,KKB),
                        ZTHVM_F,KKB,GLMIX,.TRUE.,ZSHEAR,ZLUP)
 ZLUP(:)=MAX(ZLUP(:),1.E-10)
 
-! Compute Buoyancy flux at the ground
-ZWTHVSURF(:) = (ZTHVM_F(:,KKB)/ZTHM_F(:,KKB))*PSFTH(:)+     &
-                (0.61*ZTHM_F(:,KKB))*PSFRV(:)
-
-! Mass flux at KKB level (updraft triggered if PSFTH>0.)
-WHERE (ZWTHVSURF(:)>0.010)  ! <==  Not 0 Important to have stratocumulus !!!!!
-  PEMF(:,KKB) = XCMF * ZRHO_F(:,KKB) * ((ZG_O_THVREF(:,KKB))*ZWTHVSURF*ZLUP)**(1./3.)
-  PFRAC_UP(:,KKB)=MIN(PEMF(:,KKB)/(SQRT(ZW_UP2(:,KKB))*ZRHO_F(:,KKB)),XFRAC_UP_MAX)
-  ZW_UP2(:,KKB)=(PEMF(:,KKB)/(PFRAC_UP(:,KKB)*ZRHO_F(:,KKB)))**2
-  GTEST(:)=.TRUE.
-ELSEWHERE
-  PEMF(:,KKB) =0.
-  GTEST(:)=.FALSE.
-ENDWHERE
+DO JI=1,IIJU
+  ! Compute Buoyancy flux at the ground
+  ZWTHVSURF = (ZTHVM_F(JI,KKB)/ZTHM_F(JI,KKB))*PSFTH(JI)+     &
+              (0.61*ZTHM_F(JI,KKB))*PSFRV(JI)
+
+  ! Mass flux at KKB level (updraft triggered if PSFTH>0.)
+  IF (ZWTHVSURF>0.010) THEN ! <==  Not 0 Important to have stratocumulus !!!!!
+    PEMF(JI,KKB) = XCMF * ZRHO_F(JI,KKB) * ((ZG_O_THVREF(JI,KKB))*ZWTHVSURF*ZLUP(JI))**(1./3.)
+    PFRAC_UP(JI,KKB)=MIN(PEMF(JI,KKB)/(SQRT(ZW_UP2(JI,KKB))*ZRHO_F(JI,KKB)),XFRAC_UP_MAX)
+
+    ZW_UP2(JI,KKB)=(PEMF(JI,KKB)/(PFRAC_UP(JI,KKB)*ZRHO_F(JI,KKB)))**2
+    GTEST(JI)=.TRUE.
+  ELSE
+    PEMF(JI,KKB) =0.
+    GTEST(JI)=.FALSE.
+  ENDIF
+ENDDO
 
   
 !--------------------------------------------------------------------------
@@ -429,21 +367,19 @@ ENDWHERE
 !
 !
 GTESTLCL(:)=.FALSE.
-GTESTETL(:)=.FALSE.
 
 
 !       Loop on vertical level to compute W
 
 ZW_MAX(:)      = 0.
 ZZTOP(:)       = 0.
-ZPHI(:) = 0.
 
 DO JK=KKB,KKE-KKL,KKL
 
 ! IF the updraft top is reached for all column, stop the loop on levels
 
-  ITEST=COUNT(GTEST)
-!  IF (ITEST==0) CYCLE  !  <== I do not remember why I removed this ... 
+  !ITEST=COUNT(GTEST)
+  !IF (ITEST==0) CYCLE  !  <== I do not remember why I removed this ...
 
 !       Computation of entrainment and detrainment with KF90
 !       parameterization in clouds and LR01 in subcloud layer
@@ -451,10 +387,12 @@ DO JK=KKB,KKE-KKL,KKL
  
 ! to find the LCL (check if JK is LCL or not)
 
-  WHERE ((PRC_UP(:,JK)+PRI_UP(:,JK)>0.).AND.(.NOT.(GTESTLCL)))
-      KKLCL(:) = JK           
-      GTESTLCL(:)=.TRUE.
-  ENDWHERE
+  DO JI=1,IIJU
+    IF ((PRC_UP(JI,JK)+PRI_UP(JI,JK)>0.).AND.(.NOT.(GTESTLCL(JI)))) THEN
+      KKLCL(JI) = JK           
+      GTESTLCL(JI)=.TRUE.
+    ENDIF
+  ENDDO
 
 
 ! COMPUTE PENTR and PDETR at mass level JK
@@ -471,67 +409,84 @@ DO JK=KKB,KKE-KKL,KKL
                PPABSM(:,JK),PTHL_UP(:,JK),PRT_UP(:,JK),&
                ZTH_UP(:,JK),ZRV_UP,ZRC_UP,ZRI_UP,ZRSATW(:),ZRSATI(:))            
     
-    WHERE (GTEST)
-      PTHV_UP   (:,JK) = ZTH_UP(:,JK)*(1.+ZRVORD*ZRV_UP(:))/(1.+PRT_UP(:,JK))
-      ZBUO      (:,JK) = ZG_O_THVREF(:,JK)*(PTHV_UP(:,JK) - ZTHVM_F(:,JK))    
-      PBUO_INTEG(:,JK) = ZBUO(:,JK)*(PZZ(:,JK+KKL)-PZZ(:,JK))
+  DO JI=1,IIJU
+    IF (GTEST(JI)) THEN
+      PTHV_UP   (JI,JK) = ZTH_UP(JI,JK)*(1.+ZRVORD*ZRV_UP(JI))/(1.+PRT_UP(JI,JK))
+      ZBUO      (JI,JK) = ZG_O_THVREF(JI,JK)*(PTHV_UP(JI,JK) - ZTHVM_F(JI,JK))    
+      PBUO_INTEG(JI,JK) = ZBUO(JI,JK)*(PZZ(JI,JK+KKL)-PZZ(JI,JK))
       
-      ZDZ(:)   = MAX(ZEPS,PZZ(:,JK+KKL)-PZZ(:,JK))
-      ZTEST(:) = XA1*ZBUO(:,JK) -  XB*ZW_UP2(:,JK)
-
-      ZCOE(:)      = ZDZ(:)
-      WHERE (ZTEST(:)>0.)
-        ZCOE(:)    = ZDZ(:)/(1.+ ZBETA1(:))
-      ENDWHERE   
-
-!  Convective Vertical speed computation
-
-      ZWCOE(:)         = (1.-XB*ZCOE(:))/(1.+XB*ZCOE(:))
-      ZBUCOE(:)        =  2.*ZCOE(:)/(1.+XB*ZCOE(:))
-
-! Second Rachel bug correction (XA1 has been forgotten ... not yet tested ...)
-!      ZW_UP2(:,JK+KKL) = MAX(ZEPS,ZW_UP2(:,JK)*ZWCOE(:) + ZBUO(:,JK)*ZBUCOE(:) ) 
-      ZW_UP2(:,JK+KKL) = MAX(ZEPS,ZW_UP2(:,JK)*ZWCOE(:) + XA1*ZBUO(:,JK)*ZBUCOE(:) )  
-      ZW_MAX(:) = MAX(ZW_MAX(:), SQRT(ZW_UP2(:,JK+KKL)))
-      ZWUP_MEAN(:)     = MAX(ZEPS,0.5*(ZW_UP2(:,JK+KKL)+ZW_UP2(:,JK)))
+      ZDZ(JI)   = MAX(ZEPS,PZZ(JI,JK+KKL)-PZZ(JI,JK))
+      ZTEST(JI) = XA1*ZBUO(JI,JK) -  XB*ZW_UP2(JI,JK)
+
+      !  Ancien calcul de la vitesse
+      ZCOE(JI)      = ZDZ(JI)
+      IF (ZTEST(JI)>0.) THEN
+        ZCOE(JI)    = ZDZ(JI)/(1.+ XBETA1)
+      ENDIF
+
+      !  Convective Vertical speed computation
+      ZWCOE(JI)         = (1.-XB*ZCOE(JI))/(1.+XB*ZCOE(JI))
+      ZBUCOE(JI)        =  2.*ZCOE(JI)/(1.+XB*ZCOE(JI))
+
+      ! Second Rachel bug correction (XA1 has been forgotten ... not yet tested ...)
+      !ZW_UP2(JI,JK+KKL) = MAX(ZEPS,ZW_UP2(JI,JK)*ZWCOE(JI) + ZBUO(JI,JK)*ZBUCOE(JI) ) 
+      ZW_UP2(JI,JK+KKL) = MAX(ZEPS,ZW_UP2(JI,JK)*ZWCOE(JI) + XA1*ZBUO(JI,JK)*ZBUCOE(JI) )  
+      ZW_MAX(JI) = MAX(ZW_MAX(JI), SQRT(ZW_UP2(JI,JK+KKL)))
+      ZWUP_MEAN(JI)     = MAX(ZEPS,0.5*(ZW_UP2(JI,JK+KKL)+ZW_UP2(JI,JK)))
  
-!  Entrainement and detrainement
+      !  Entrainement and detrainement
 
-! First Rachel bug correction (Parenthesis around 1+beta1 ==> impact is small)   
-     PENTR(:,JK)  = MAX(0.,(ZBETA1(:)/(1.+ZBETA1(:)))*(XA1*ZBUO(:,JK)/ZWUP_MEAN(:)-XB))
-     ZDETR_BUO(:) = MAX(0., -(ZBETA1(:)/(1.+ZBETA1(:)))*XA1*ZBUO(:,JK)/ZWUP_MEAN(:))
-     ZDETR_RT(:)  = XC*SQRT(MAX(0.,(PRT_UP(:,JK) - ZRTM_F(:,JK))) / MAX(ZEPS,ZRTM_F(:,JK)) / ZWUP_MEAN(:))
-     PDETR(:,JK)  = ZDETR_RT(:)+ZDETR_BUO(:)
+      ! First Rachel bug correction (Parenthesis around 1+beta1 ==> impact is small)   
+      PENTR(JI,JK)  = MAX(0.,(XBETA1/(1.+XBETA1))*(XA1*ZBUO(JI,JK)/ZWUP_MEAN(JI)-XB))
+      ZDETR_BUO(JI) = MAX(0., -(XBETA1/(1.+XBETA1))*XA1*ZBUO(JI,JK)/ZWUP_MEAN(JI))
+      ZDETR_RT(JI)  = XC*SQRT(MAX(0.,(PRT_UP(JI,JK) - ZRTM_F(JI,JK))) / MAX(ZEPS,ZRTM_F(JI,JK)) / ZWUP_MEAN(JI))
+      PDETR(JI,JK)  = ZDETR_RT(JI)+ZDETR_BUO(JI)
    
-! If the updraft did not stop, compute cons updraft characteritics at jk+1
-
-     ZZTOP(:) = MAX(ZZTOP(:),PZZ(:,JK+KKL))
-     ZMIX2(:) = (PZZ(:,JK+KKL)-PZZ(:,JK))*PENTR(:,JK) !&
-     ZMIX3(:) = (PZZ(:,JK+KKL)-PZZ(:,JK))*PDETR(:,JK) !&                   
+      ! If the updraft did not stop, compute cons updraft characteritics at jk+1
+      ZZTOP(JI) = MAX(ZZTOP(JI),PZZ(JI,JK+KKL))
+      ZMIX2(JI) = (PZZ(JI,JK+KKL)-PZZ(JI,JK))*PENTR(JI,JK) !&
                 
-     PTHL_UP(:,JK+KKL)=(PTHL_UP(:,JK)*(1.-0.5*ZMIX2(:)) + PTHLM(:,JK)*ZMIX2(:)) &
-                          /(1.+0.5*ZMIX2(:))   
-     PRT_UP(:,JK+KKL) =(PRT_UP (:,JK)*(1.-0.5*ZMIX2(:)) + PRTM(:,JK)*ZMIX2(:))  &
-                          /(1.+0.5*ZMIX2(:))
-  ENDWHERE  ! GTEST
+      PTHL_UP(JI,JK+KKL)=(PTHL_UP(JI,JK)*(1.-0.5*ZMIX2(JI)) + PTHLM(JI,JK)*ZMIX2(JI)) &
+                           /(1.+0.5*ZMIX2(JI))   
+      PRT_UP(JI,JK+KKL) =(PRT_UP (JI,JK)*(1.-0.5*ZMIX2(JI)) + PRTM(JI,JK)*ZMIX2(JI))  &
+                           /(1.+0.5*ZMIX2(JI))
+    ENDIF  ! GTEST
+  ENDDO
   
 
   IF(OMIXUV) THEN
-    WHERE(GTEST) 
-      PU_UP(:,JK+KKL) = (PU_UP (:,JK)*(1-0.5*ZMIX2(:)) + PUM(:,JK)*ZMIX2(:)+ &
-                        0.5*XPRES_UV*(PZZ(:,JK+KKL)-PZZ(:,JK))*&
-                        ((PUM(:,JK+KKL)-PUM(:,JK))/PDZZ(:,JK+KKL)+&
-                         (PUM(:,JK)-PUM(:,JK-KKL))/PDZZ(:,JK))        )   &
-                        /(1+0.5*ZMIX2(:))
-      PV_UP(:,JK+KKL) = (PV_UP (:,JK)*(1-0.5*ZMIX2(:)) + PVM(:,JK)*ZMIX2(:)+ &
-                        0.5*XPRES_UV*(PZZ(:,JK+KKL)-PZZ(:,JK))*&
-                        ((PVM(:,JK+KKL)-PVM(:,JK))/PDZZ(:,JK+KKL)+&
-                         (PVM(:,JK)-PVM(:,JK-KKL))/PDZZ(:,JK))    )   &
-                        /(1+0.5*ZMIX2(:))
-    ENDWHERE
+    IF(JK/=KKB) THEN
+      DO JI=1,IIJU
+        IF(GTEST(JI)) THEN
+          PU_UP(JI,JK+KKL) = (PU_UP (JI,JK)*(1-0.5*ZMIX2(JI)) + PUM(JI,JK)*ZMIX2(JI)+ &
+                            0.5*XPRES_UV*(PZZ(JI,JK+KKL)-PZZ(JI,JK))*&
+                            ((PUM(JI,JK+KKL)-PUM(JI,JK))/PDZZ(JI,JK+KKL)+&
+                             (PUM(JI,JK)-PUM(JI,JK-KKL))/PDZZ(JI,JK))        )   &
+                            /(1+0.5*ZMIX2(JI))
+          PV_UP(JI,JK+KKL) = (PV_UP (JI,JK)*(1-0.5*ZMIX2(JI)) + PVM(JI,JK)*ZMIX2(JI)+ &
+                            0.5*XPRES_UV*(PZZ(JI,JK+KKL)-PZZ(JI,JK))*&
+                            ((PVM(JI,JK+KKL)-PVM(JI,JK))/PDZZ(JI,JK+KKL)+&
+                             (PVM(JI,JK)-PVM(JI,JK-KKL))/PDZZ(JI,JK))    )   &
+                            /(1+0.5*ZMIX2(JI))
+        ENDIF
+      ENDDO
+    ELSE
+      DO JI=1,IIJU
+        IF(GTEST(JI)) THEN
+          PU_UP(JI,JK+KKL) = (PU_UP (JI,JK)*(1-0.5*ZMIX2(JI)) + PUM(JI,JK)*ZMIX2(JI)+ &
+                            0.5*XPRES_UV*(PZZ(JI,JK+KKL)-PZZ(JI,JK))*&
+                            ((PUM(JI,JK+KKL)-PUM(JI,JK))/PDZZ(JI,JK+KKL))        ) &
+                            /(1+0.5*ZMIX2(JI))
+          PV_UP(JI,JK+KKL) = (PV_UP (JI,JK)*(1-0.5*ZMIX2(JI)) + PVM(JI,JK)*ZMIX2(JI)+ &
+                            0.5*XPRES_UV*(PZZ(JI,JK+KKL)-PZZ(JI,JK))*&
+                            ((PVM(JI,JK+KKL)-PVM(JI,JK))/PDZZ(JI,JK+KKL))    )   &
+                            /(1+0.5*ZMIX2(JI))
+        ENDIF
+      ENDDO
+    ENDIF
   ENDIF
   
-! This updraft is not yet ready to use scalar variables  
+! This updraft is not yet ready to use scalar variables
 !  DO JSV=1,ISV 
 !     IF (ONOMIXLG .AND. JSV >= KSV_LGBEG .AND. JSV<= KSV_LGEND) CYCLE
 !      WHERE(GTEST) 
@@ -548,53 +503,58 @@ DO JK=KKB,KKE-KKL,KKL
   CALL TH_R_FROM_THL_RT_1D(HFRAC_ICE,PFRAC_ICE_UP(:,JK+KKL),ZPRES_F(:,JK+KKL), &
           PTHL_UP(:,JK+KKL),PRT_UP(:,JK+KKL),ZTH_UP(:,JK+KKL),              &
           ZRV_UP(:),ZRC_UP(:),ZRI_UP(:),ZRSATW(:),ZRSATI(:))
-  WHERE(GTEST)
-    PRC_UP(:,JK+KKL)=ZRC_UP(:)
-    PRV_UP(:,JK+KKL)=ZRV_UP(:)
-    PRI_UP(:,JK+KKL)=ZRI_UP(:)
-    PRSAT_UP(:,JK+KKL) = ZRSATW(:)*(1-PFRAC_ICE_UP(:,JK+KKL)) + ZRSATI(:)*PFRAC_ICE_UP(:,JK+KKL)
-  ENDWHERE
-  
-
-! Compute the updraft theta_v, buoyancy and w**2 for level JK+1   
- WHERE(GTEST)
-      PTHV_UP(:,JK+KKL) = ZTH_UP(:,JK+KKL)*((1+ZRVORD*PRV_UP(:,JK+KKL))/(1+PRT_UP(:,JK+KKL)))
- ENDWHERE
-
- WHERE(GTEST) 
-    ZMIX1(:)=ZZDZ(:,JK)*(PENTR(:,JK)-PDETR(:,JK))
-    PEMF(:,JK+KKL)=PEMF(:,JK)*EXP(ZMIX1(:))
 
-! Updraft fraction must be smaller than XFRAC_UP_MAX
-    PFRAC_UP(:,JK+KKL)=PEMF(:,JK+KKL)/(SQRT(ZW_UP2(:,JK+KKL))*ZRHO_F(:,JK+KKL))
-    PFRAC_UP(:,JK+KKL)=MIN(XFRAC_UP_MAX,PFRAC_UP(:,JK+KKL))
-  ENDWHERE  
+  DO JI=1,IIJU
+    IF(GTEST(JI)) THEN
+      PRC_UP(JI,JK+KKL)=ZRC_UP(JI)
+      PRV_UP(JI,JK+KKL)=ZRV_UP(JI)
+      PRI_UP(JI,JK+KKL)=ZRI_UP(JI)
+      PRSAT_UP(JI,JK+KKL) = ZRSATW(JI)*(1-PFRAC_ICE_UP(JI,JK+KKL)) + ZRSATI(JI)*PFRAC_ICE_UP(JI,JK+KKL)
+
+      ! Compute the updraft theta_v, buoyancy and w**2 for level JK+1   
+      PTHV_UP(JI,JK+KKL) = ZTH_UP(JI,JK+KKL)*((1+ZRVORD*PRV_UP(JI,JK+KKL))/(1+PRT_UP(JI,JK+KKL)))
+      ZMIX1(JI)=ZZDZ(JI,JK)*(PENTR(JI,JK)-PDETR(JI,JK))
+    ENDIF
+  ENDDO
+
+  DO JI=1,IIJU
+    IF(GTEST(JI)) THEN
+      PEMF(JI,JK+KKL)=PEMF(JI,JK)*EXP(ZMIX1(JI))
+    ENDIF
+  ENDDO
+
+  DO JI=1,IIJU
+    IF(GTEST(JI)) THEN
+      ! Updraft fraction must be smaller than XFRAC_UP_MAX
+      PFRAC_UP(JI,JK+KKL)=MIN(XFRAC_UP_MAX, &
+                             &PEMF(JI,JK+KKL)/(SQRT(ZW_UP2(JI,JK+KKL))*ZRHO_F(JI,JK+KKL)))
+    ENDIF
+  ENDDO
 
 ! Test if the updraft has reach the ETL
-  GTESTETL(:)=.FALSE.
-  WHERE (GTEST.AND.(PBUO_INTEG(:,JK)<=0.))
-      KKETL(:) = JK+KKL
-      GTESTETL(:)=.TRUE.
-  ENDWHERE
+  DO JI=1,IIJU
+    IF (GTEST(JI) .AND. (PBUO_INTEG(JI,JK)<=0.)) THEN
+      KKETL(JI) = JK+KKL
+    ENDIF
+  ENDDO
 
 
 ! Test is we have reached the top of the updraft
-
-  WHERE (GTEST.AND.((ZW_UP2(:,JK+KKL)<=ZEPS).OR.(PEMF(:,JK+KKL)<=ZEPS)))
-      ZW_UP2   (:,JK+KKL)=ZEPS
-      PEMF     (:,JK+KKL)=0.
-      GTEST    (:)       =.FALSE.
-      PTHL_UP  (:,JK+KKL)=ZTHLM_F(:,JK+KKL)
-      PRT_UP   (:,JK+KKL)=ZRTM_F(:,JK+KKL)
-      PRC_UP   (:,JK+KKL)=0.
-      PRI_UP   (:,JK+KKL)=0.
-      PRV_UP   (:,JK+KKL)=ZRVM_F (:,JK+KKL)
-      PTHV_UP  (:,JK+KKL)=ZTHVM_F(:,JK+KKL)
-      PFRAC_UP (:,JK+KKL)=0.
-      KKCTL    (:)       =JK+KKL
-      
-  ENDWHERE
-
+  DO JI=1,IIJU
+    IF (GTEST(JI) .AND. ((ZW_UP2(JI,JK+KKL)<=ZEPS).OR.(PEMF(JI,JK+KKL)<=ZEPS))) THEN
+      ZW_UP2   (JI,JK+KKL)=ZEPS
+      PEMF     (JI,JK+KKL)=0.
+      GTEST    (JI)       =.FALSE.
+      PTHL_UP  (JI,JK+KKL)=ZTHLM_F(JI,JK+KKL)
+      PRT_UP   (JI,JK+KKL)=ZRTM_F(JI,JK+KKL)
+      PRC_UP   (JI,JK+KKL)=0.
+      PRI_UP   (JI,JK+KKL)=0.
+      PRV_UP   (JI,JK+KKL)=ZRVM_F (JI,JK+KKL)
+      PTHV_UP  (JI,JK+KKL)=ZTHVM_F(JI,JK+KKL)
+      PFRAC_UP (JI,JK+KKL)=0.
+      KKCTL    (JI)       =JK+KKL
+    ENDIF
+  ENDDO
 
 ENDDO   ! Fin de la boucle verticale 
 
@@ -607,19 +567,24 @@ PEMF(:,KKB) =0.
 ! This way, all MF fluxes are diminished by this amount.
 ! Diagnosed cloud fraction is also multiplied by the same coefficient.
 !
-DO JI=1,SIZE(PTHM,1) 
+DO JI=1,IIJU
    PDEPTH(JI) = MAX(0., PZZ(JI,KKCTL(JI)) -  PZZ(JI,KKLCL(JI)) )
-END DO
+ENDDO
 
 GWORK1(:)= (GTESTLCL(:) .AND. (PDEPTH(:) > ZDEPTH_MAX1) )
 GWORK2(:,:) = SPREAD( GWORK1(:), DIM=2, NCOPIES=IKU )
 ZCOEF(:,:) = SPREAD( (1.-(PDEPTH(:)-ZDEPTH_MAX1)/(ZDEPTH_MAX2-ZDEPTH_MAX1)), DIM=2, NCOPIES=IKU)
-ZCOEF=MIN(MAX(ZCOEF,0.),1.)
-WHERE (GWORK2) 
-   PEMF(:,:)     = PEMF(:,:)     * ZCOEF(:,:) 
-   PFRAC_UP(:,:) = PFRAC_UP(:,:) * ZCOEF(:,:) 
-ENDWHERE
-
-END SUBROUTINE COMPUTE_UPDRAFT_RHCJ10
+ZCOEF(:,:)=MIN(MAX(ZCOEF(:,:),0.),1.)
+DO JK=1, IKU
+  DO JI=1,IIJU
+    IF (GWORK2(JI,JK)) THEN
+      PEMF(JI,JK)     = PEMF(JI,JK)     * ZCOEF(JI,JK) 
+      PFRAC_UP(JI,JK) = PFRAC_UP(JI,JK) * ZCOEF(JI,JK) 
+    ENDIF
+  ENDDO
+ENDDO
 
+IF (LHOOK) CALL DR_HOOK('COMPUTE_UPDRAFT_RHCJ10',1,ZHOOK_HANDLE)
 
+END SUBROUTINE COMPUTE_UPDRAFT_RHCJ10
+END MODULE MODE_COMPUTE_UPDRAFT_RHCJ10
diff --git a/src/mesonh/turb/shallow_mf.f90 b/src/mesonh/turb/shallow_mf.f90
index 6d6942708c6d6d0cf3cd1841cb052b8c79a3fcb3..6a98885dd1969802e0b3c99831aca0adee50045d 100644
--- a/src/mesonh/turb/shallow_mf.f90
+++ b/src/mesonh/turb/shallow_mf.f90
@@ -183,7 +183,7 @@ USE MODD_TURB_n, ONLY: CTURBLEN
 USE MODI_THL_RT_FROM_TH_R_MF
 USE MODI_COMPUTE_UPDRAFT
 USE MODI_COMPUTE_UPDRAFT_RHCJ10
-USE MODI_COMPUTE_UPDRAFT_RAHA
+USE MODE_COMPUTE_UPDRAFT_RAHA, ONLY: COMPUTE_UPDRAFT_RAHA
 USE MODE_MF_TURB, ONLY: MF_TURB
 USE MODE_MF_TURB_EXPL, ONLY: MF_TURB_EXPL
 USE MODI_MF_TURB_GREYZONE