Index: LMDZ4/trunk/libf/phylmd/cloudth.F90
===================================================================
--- LMDZ4/trunk/libf/phylmd/cloudth.F90	(revision 1410)
+++ LMDZ4/trunk/libf/phylmd/cloudth.F90	(revision 1411)
@@ -142,19 +142,7 @@
 ! Calcul des écart-types pour s
 !-----------------------------------------------------------------------------------------------------------------
-!      alpha=0.5*(fraca(ind1,ind2)+fraca(ind1,ind2-1))
-!      sigma1s=(Max(1.2*fraca(ind1,ind2)*(sth-senv)**2,(senv/100)**2))**0.5
-!      sigma2s=((0.021/(fraca(ind1,ind2)+0.0)**0.5)*(sth-senv)**2+(sth/100)**2)**0.5 
-!      sigma1s=(1.5**0.5)*(fraca(ind1,ind2)**0.65)*(sth-senv)+0.000045
-!      sigma2s=0.1265*(sth-senv)/(fraca(ind1,ind2)+0.0)**0.35+0.000045     
-
-      sigma1s=(1.1**0.5)*(fraca(ind1,ind2)**0.6)/(1-fraca(ind1,ind2))*((sth-senv)**2)**0.5+0.00003
-      sigma2s=0.11*((sth-senv)**2)**0.5/(fraca(ind1,ind2)+0.02)**0.4+0.00003  
-
-!      sigma1s=(1.5**0.5)*(alpha**0.65)*(sth-senv)+0.000045
-!      sigma2s=0.1265*(sth-senv)/alpha**0.35+0.000045
-      
-!      sigma1s=(2.8**0.5)*(0.1**0.7)*(sth-senv)+0.00002
-!      sigma2s=((0.126/(0.1)**0.3)*(sth-senv))+0.00002
-
+
+      sigma1s=(1.1**0.5)*(fraca(ind1,ind2)**0.6)/(1-fraca(ind1,ind2))*((sth-senv)**2)**0.5+ratqs(ind1,ind2)*po(ind1)
+      sigma2s=0.11*((sth-senv)**2)**0.5/(fraca(ind1,ind2)+0.02)**0.4+0.002*zqta(ind1,ind2)  
 
  
@@ -214,15 +202,7 @@
       senv=aenv*(po(ind1)-zqsatenv(ind1,ind2)) 
       
-!      if (zqenv(ind1).gt.0.005) then
-      sigma1s=0.005*zqenv(ind1)
-!      else 
-!      sigma1s=0.0005
-!      endif
-
-!      sigma1s=ratqs(ind1,ind2)*zqenv(ind1)
-
-
-
-!      sigma1s=0.00003
+
+      sigma1s=ratqs(ind1,ind2)*zqenv(ind1)
+
       sqrt2pi=sqrt(2.*pi)
       xenv=senv/(sqrt(2.)*sigma1s)
Index: LMDZ4/trunk/libf/phylmd/fisrtilp.F
===================================================================
--- LMDZ4/trunk/libf/phylmd/fisrtilp.F	(revision 1410)
+++ LMDZ4/trunk/libf/phylmd/fisrtilp.F	(revision 1411)
@@ -50,4 +50,6 @@
       REAL zthl(klon,klev)
 
+      logical lognormale(klon)
+
 cAA
 c Coeffients de fraction lessivee : pour OFF-LINE
@@ -140,5 +142,5 @@
       zdelq=0.0
       
-!      print*,'CLOUDTH4 A. JAM'
+      print*,'NUAGES4 A. JAM'
       IF (appel1er) THEN
 c
@@ -289,5 +291,5 @@
      .                            /RG/dtime
 
-c pour la glace, on r�vapore toute la pr�ip dans la couche du dessous
+c pour la glace, on rï¿½vapore toute la prï¿½ip dans la couche du dessous
 c la glace venant de la couche du dessus est simplement dans la couche
 c du dessous.
@@ -343,5 +345,5 @@
 c   zqn   : eau totale dans le nuage
 c   zcond : eau condensee moyenne dans la maille.
-c           on prend en compte le r�hauffement qui diminue la partie condensee
+c           on prend en compte le rï¿½hauffement qui diminue la partie condensee
 c
 c   Version avec les raqts
@@ -365,5 +367,5 @@
            enddo 
 
-              if (iflag_cldcon.eq.5) then
+              if (iflag_cldcon>=5) then
 
                  call cloudth(klon,klev,k,ztv,
@@ -377,7 +379,14 @@
                  enddo
 
-              else
-
+              endif
+
+! Pour iflag_cldcon<=4, on prend toujours la lognormale
+! Dans le cas iflag_cldcon=5, on prend systÃ©matiquement la bi-gaussienne
+! Dans le cas iflagÃ_cldcon=6, on prend la lognormale en absence des thermiques
+
+            lognormale(:)=
+     .      iflag_cldcon<=4.or.(iflag_cldcon==6.and.fraca(:,k)<1.e-10)
             do i=1,klon
+            if (lognormale(i)) then
             zpdf_sig(i)=ratqs(i,k)*zq(i)
             zpdf_k(i)=-sqrt(log(1.+(zpdf_sig(i)/zq(i))**2))
@@ -391,4 +400,9 @@
             zpdf_e2(i)=sign(min(abs(zpdf_e2(i)),5.),zpdf_e2(i))
             zpdf_e2(i)=1.-erf(zpdf_e2(i))
+            endif
+            enddo
+
+            do i=1,klon
+            if (lognormale(i)) then
             if (zpdf_e1(i).lt.1.e-10) then
                rneb(i,k)=0.
@@ -398,8 +412,8 @@
                zqn(i)=zq(i)*zpdf_e2(i)/zpdf_e1(i)
             endif
+            endif
             
            enddo 
 
-         endif ! iflag_cldcon
 
         endif ! iflag_pdf
Index: LMDZ4/trunk/libf/phylmd/physiq.F
===================================================================
--- LMDZ4/trunk/libf/phylmd/physiq.F	(revision 1410)
+++ LMDZ4/trunk/libf/phylmd/physiq.F	(revision 1411)
@@ -1,3 +1,4 @@
 ! $Id$
+!
 c#define IO_DEBUG
 
@@ -261,4 +262,5 @@
       CHARACTER*4 bb2
       CHARACTER*2 bb3
+c
 
       real twriteSTD(klon,nlevSTD,nfiles)
@@ -863,5 +865,4 @@
       REAL rflag(klon)          ! flag fonctionnement de convect
       INTEGER iflagctrl(klon)          ! flag fonctionnement de convect
-
 c -- convect43:
       INTEGER ntra              ! nb traceurs pour convect4.3
@@ -1232,4 +1233,5 @@
       call phys_state_var_init(read_climoz)
       call phys_output_var_init
+
       print*, '================================================='
 cIM for NMC files
@@ -1251,5 +1253,4 @@
 c         pmflxr=0.
 c         pmflxs=0.
-
         itau_con=0
         first=.false.
@@ -1423,5 +1424,5 @@
            ema_pcb(i)  = 0.
            ema_pct(i)  = 0.
-c          ema_workcbmf(i) = 0.
+           ema_workcbmf(i) = 0.
           ENDDO
 cIM15/11/02 rajout initialisation ibas_con,itop_con cf. SB =>BEG
@@ -1448,5 +1449,5 @@
 c================================================================================
 
-         ENDIF !debut
+         ENDIF
 
            DO i=1,klon
@@ -2170,5 +2171,5 @@
      .        d_t_con,d_q_con,d_u_con,d_v_con,d_tr,
      .        rain_con, snow_con, ibas_con, itop_con, sigd,
-     .        ema_cbmf,upwd,dnwd,dnwd0,
+     .        upwd,dnwd,dnwd0,
      .        Ma,mip,Vprecip,cape,cin,tvp,Tconv,iflagctrl,
      .        pbase,bbase,dtvpdt1,dtvpdq1,dplcldt,dplcldr,qcondc,wd,
@@ -2189,5 +2190,4 @@
 
           ELSE ! ok_cvl
-
 c MAF conema3 ne contient pas les traceurs
           CALL conema3 (dtime,
@@ -2229,7 +2229,8 @@
 
           DO i = 1, klon
-            ema_pcb(i)  = paprs(i,ibas_con(i))
+            ema_pcb(i)  = pbase(i)
           ENDDO
           DO i = 1, klon
+
 ! L'idicage de itop_con peut cacher un pb potentiel
 ! FH sous la dictee de JYG, CR
@@ -2242,5 +2243,8 @@
               endif
             endif
-          ENDDO     
+          ENDDO
+          DO i = 1, klon
+            ema_cbmf(i) = ema_workcbmf(i)
+          ENDDO      
       ELSE IF (iflag_con.eq.0) THEN
           write(lunout,*) 'On n appelle pas la convection'
@@ -2459,5 +2463,5 @@
 c  ==============
 
-! Dans le cas oÃ¹ on active les thermiques, on fait partir l'ajustement
+! Dans le cas où on active les thermiques, on fait partir l'ajustement
 ! a partir du sommet des thermiques.
 ! Dans le cas contraire, on demarre au niveau 1.
@@ -2544,54 +2548,4 @@
          if(prt_level.ge.9) print*,' CLOUDS_GNO OK'
 
-c-----------------------------------------------------------------------
-c   par calcul direct de l'ecart-type
-c-----------------------------------------------------------------------
-
-      else if (iflag_cldcon>=5) then
-         wmax_th(:)=0.
-         zmax_th(:)=0.
-         do k=1,klev
-            do i=1,klon
-               wmax_th(i)=max(wmax_th(i),zw2(i,k))
-               if (detr_therm(i,k).gt.0.) zmax_th(i)=pphi(i,k)/rg
-            enddo
-         enddo
-         tau_overturning_th(:)=zmax_th(:)/max(0.5*wmax_th(:),0.1)
-         if(prt_level.ge.9)
-     &     write(lunout,*)'TAU TH OK ',
-     &     tau_overturning_th(1),detr_therm(1,3)
-
-c On impose que l'air autour de la fraction couverte par le thermique
-c plus son air detraine durant tau_overturning_th soit superieur
-c a 0.1 q_seri
-         zz=0.1
-         do k=1,klev
-            do i=1,klon
-               lambda_th(i,k)=0.5*(fraca(i,k)+fraca(i,k+1))+
-     s         tau_overturning_th(i)*detr_therm(i,k)
-     s         *rg/(paprs(i,k)-paprs(i,k+1))
-               znum=(1.-zz)*q_seri(i,k)
-               zden=zqasc(i,k)-zz*q_seri(i,k)
-               if (znum-lambda_th(i,k)*zden<0.) lambda_th(i,k)=znum/zden
-               lambda_th(i,k)=min(lambda_th(i,k),0.9)
-            enddo
-         enddo
-
-         if(iflag_cldcon==5) then
-            do k=1,klev
-               do i=1,klon
-                  ratqsc(i,k)=sqrt(lambda_th(i,k)/(1.-lambda_th(i,k)))*
-     s            abs((zqasc(i,k)-q_seri(i,k))/q_seri(i,k))
-               enddo
-            enddo
-         else if(iflag_cldcon==6) then
-            do k=1,klev
-               do i=1,klon
-                  ratqsc(i,k)=sqrt(lambda_th(i,k))*
-     s            (zqasc(i,k)-q_seri(i,k))/q_seri(i,k)
-               enddo
-            enddo
-         endif
-
       endif
 
@@ -2657,5 +2611,5 @@
 
       if (iflag_cldcon.eq.1 .or.iflag_cldcon.eq.2
-     s    .or.iflag_cldcon.ge.4) then
+     s    .or.iflag_cldcon.eq.4) then
 
 ! On ajoute une constante au ratqsc*2 pour tenir compte de 
@@ -2692,4 +2646,5 @@
       endif
 
+      print*,'PHSYIQ NUAGES4'
 
 c
@@ -2869,6 +2824,25 @@
 
 c   On prend la somme des fractions nuageuses et des contenus en eau
-      cldfra(:,:)=min(max(cldfra(:,:),rnebcon(:,:)),1.)
-      cldliq(:,:)=cldliq(:,:)+rnebcon(:,:)*clwcon(:,:)
+
+      if (iflag_cldcon>=5) then
+
+! Si on est sur un point touche par la convection profonde et pas
+! par les thermiques, on prend la couverture nuageuse et l'eau nuageuse
+! de la convection profonde.
+
+         print*,'TEST SCHEMA DE NUAGES '
+         do k=1,klev
+            do i=1,klon
+               if (ptconv(i,k).and. .not. ptconvth(i,k)) then
+                   cldfra(i,k)=rnebcon(i,k)
+                   cldliq(i,k)=rnebcon(i,k)*clwcon(i,k)
+               endif
+            enddo
+         enddo
+      else
+! Ancienne version
+         cldfra(:,:)=min(max(cldfra(:,:),rnebcon(:,:)),1.)
+         cldliq(:,:)=cldliq(:,:)+rnebcon(:,:)*clwcon(:,:)
+      endif
 
       ENDIF
@@ -3338,7 +3312,6 @@
 !     s        ref_liq,ref_ice
           call phys_cosp(itap,dtime,freq_cosp,
-     $                   ok_mensuelCOSP,ok_journeCOSP,ok_hfCOSP,
-     $                   ecrit_mth,ecrit_day,ecrit_hf,
-     $                   klon,klev,rlon,rlat,presnivs,overlap,
+     $                 ecrit_mth,ecrit_day,ecrit_hf,overlap,
+     $                   klon,klev,rlon,rlat,presnivs,
      $                   ref_liq,ref_ice,
      $                   pctsrf(:,is_ter)+pctsrf(:,is_lic),
@@ -3349,5 +3322,4 @@
      $                   pmflxr(:,1:klev),pmflxs(:,1:klev),
      $                   mr_ozone,cldtau, cldemi)
-
 !     L          calipso2D,calipso3D,cfadlidar,parasolrefl,atb,betamol,
 !     L          cfaddbze,clcalipso2,dbze,cltlidarradar,
@@ -3497,5 +3469,4 @@
       wwriteSTD(:,:,4)=wlevSTD(:,:)
 c
-cIM initialisation 5eme fichier de sortie 
 cIM ajoute 5eme niveau 170310 BEG
       twriteSTD(:,:,5)=tlevSTD(:,:)
@@ -3610,5 +3581,4 @@
 cIM global posePB#include "write_bilKP_ave.h"
 c
-
 c Sauvegarder les valeurs de t et q a la fin de la physique:
 c
@@ -3669,8 +3639,5 @@
       DO k = 1, klev
       DO i = 1, klon
-cJYG/IM theta en debut du pas de temps
-cJYG/IM       theta(i,k)=t(i,k)*(100000./pplay(i,k))**(RD/RCPD)
-cJYG/IM theta en fin de pas de temps de physique
-        theta(i,k)=t_seri(i,k)*(100000./pplay(i,k))**(RD/RCPD)
+        theta(i,k)=t(i,k)*(100000./pplay(i,k))**(RD/RCPD)
       ENDDO
       ENDDO
