Index: LMDZ4/trunk/libf/phylmd/newmicro.F
===================================================================
--- LMDZ4/trunk/libf/phylmd/newmicro.F	(revision 1334)
+++ LMDZ4/trunk/libf/phylmd/newmicro.F	(revision 1337)
@@ -11,4 +11,8 @@
 
       USE dimphy
+      USE phys_local_var_mod, only: scdnc,cldncl,reffclwtop,lcc,
+     .                              reffclws,reffclwc,cldnvi,lcc3d,
+     .                              lcc3dcon,lcc3dstra
+      USE phys_state_var_mod, only: rnebcon,clwcon
       IMPLICIT none
 c======================================================================
@@ -46,4 +50,19 @@
 #include "radepsi.h"
 #include "radopt.h"
+c choix de l'hypothese de recouvrememnt nuageuse
+      LOGICAL RANDOM,MAXIMUM_RANDOM,MAXIMUM,FIRST
+      parameter (RANDOM=.FALSE., MAXIMUM_RANDOM=.TRUE., MAXIMUM=.FALSE.)
+c Hypoyhese de recouvrement : MAXIMUM_RANDOM
+      INTEGER flag_max
+      REAL phase3d(klon, klev),dh(klon, klev),pdel(klon, klev),
+     .     zrho(klon, klev)
+      REAL tcc(klon), ftmp(klon), lcc_integrat(klon), height(klon)
+      REAL thres_tau,thres_neb
+      PARAMETER (thres_tau=0.3, thres_neb=0.001)
+      REAL t_tmp
+      REAL gravit
+      PARAMETER (gravit=9.80616)  !m/s2
+      REAL pqlwpcon(klon, klev), pqlwpstra(klon, klev) 
+c
       REAL paprs(klon,klev+1), pplay(klon,klev)
       REAL t(klon,klev)
@@ -475,5 +494,183 @@
          pcl(i)=1.-pcl(i)
       ENDDO
-      
+
+c ========================================================
+! DIAGNOSTICS CALCULATION FOR CMIP5 PROTOCOL
+c ========================================================
+!! change by Nicolas Yan (LSCE) 
+!! Cloud Droplet Number Concentration (CDNC) : 3D variable
+!! Fractionnal cover by liquid water cloud (LCC3D) : 3D variable
+!! Cloud Droplet Number Concentration at top of cloud (CLDNCL) : 2D variable
+!! Droplet effective radius at top of cloud (REFFCLWTOP) : 2D variable
+!! Fractionnal cover by liquid water at top of clouds (LCC) : 2D variable
+      IF (ok_newmicro) THEN
+         IF (ok_aie) THEN
+            DO k = 1, klev
+               DO i = 1, klon
+                  phase3d(i,k)=1-zfice2(i,k)
+                  IF (pclc(i,k) .LE. seuil_neb) THEN
+                     lcc3d(i,k)=seuil_neb*phase3d(i,k)
+                  ELSE
+                     lcc3d(i,k)=pclc(i,k)*phase3d(i,k)
+                  ENDIF
+                  scdnc(i,k)=lcc3d(i,k)*cdnc(i,k) ! m-3
+               ENDDO
+            ENDDO
+
+            DO i=1,klon
+               lcc(i)=0.
+               reffclwtop(i)=0.
+               cldncl(i)=0.
+               IF(RANDOM .OR. MAXIMUM_RANDOM) tcc(i) = 1.
+               IF(MAXIMUM) tcc(i) = 0.
+            ENDDO
+     
+            FIRST=.TRUE.
+
+            DO i=1,klon
+               DO k=klev-1,1,-1 !From TOA down
+
+
+            ! Test, if the cloud optical depth exceeds the necessary
+            ! threshold:
+
+             IF (pcltau(i,k).GT.thres_tau .AND. pclc(i,k).GT.thres_neb)
+     .                                                             THEN
+               ! To calculate the right Temperature at cloud top,
+               ! interpolate it between layers:
+                  t_tmp = t(i,k) +
+     .              (paprs(i,k+1)-pplay(i,k))/(pplay(i,k+1)-pplay(i,k))
+     .              * ( t(i,k+1) - t(i,k) )
+
+                  IF(MAXIMUM) THEN
+                    IF(FIRST) THEN
+                       write(*,*)'Hypothese de recouvrement: MAXIMUM'
+                       FIRST=.FALSE.
+                    ENDIF
+                    flag_max= -1.
+                    ftmp(i) = MAX(tcc(i),pclc(i,k))
+                  ENDIF
+
+                  IF(RANDOM) THEN
+                    IF(FIRST) THEN
+                       write(*,*)'Hypothese de recouvrement: RANDOM'
+                       FIRST=.FALSE.
+                    ENDIF
+                    flag_max= 1.
+                    ftmp(i) = tcc(i) * (1-pclc(i,k))
+                  ENDIF
+
+                  IF(MAXIMUM_RANDOM) THEN
+                    IF(FIRST) THEN
+                       write(*,*)'Hypothese de recouvrement: MAXIMUM_
+     .                         RANDOM'
+                       FIRST=.FALSE.
+                    ENDIF
+                    flag_max= 1.
+                    ftmp(i) = tcc(i) *
+     .              (1. - MAX(pclc(i,k),pclc(i,k+1))) /
+     .              (1. - MIN(pclc(i,k+1),1.-thres_neb))
+                  ENDIF
+c Effective radius of cloud droplet at top of cloud (m)
+                  reffclwtop(i) = reffclwtop(i) + rad_chaud_tab(i,k) * 
+     .           1.0E-06 * phase3d(i,k) * ( tcc(i) - ftmp(i))*flag_max
+c CDNC at top of cloud (m-3)
+                  cldncl(i) = cldncl(i) + cdnc(i,k) * phase3d(i,k) * 
+     .                 (tcc(i) - ftmp(i))*flag_max
+c Liquid Cloud Content at top of cloud
+                  lcc(i) = lcc(i) + phase3d(i,k) * (tcc(i)-ftmp(i))*
+     .                    flag_max
+c Total Cloud Content at top of cloud
+                  tcc(i)=ftmp(i)
+              
+          ENDIF ! is there a visible, not-too-small cloud?  
+          ENDDO ! loop over k
+
+          IF(RANDOM .OR. MAXIMUM_RANDOM) tcc(i)=1.-tcc(i)
+         ENDDO ! loop over i
+
+!! Convective and Stratiform Cloud Droplet Effective Radius (REFFCLWC  REFFCLWS)
+            DO i = 1, klon
+               DO k = 1, klev
+                  pqlwpcon(i,k)=rnebcon(i,k)*clwcon(i,k) ! fraction eau liquide convective
+                  pqlwpstra(i,k)=pclc(i,k)*phase3d(i,k)-pqlwpcon(i,k) ! fraction eau liquide stratiforme
+                  IF (pqlwpstra(i,k) .LE. 0.0) pqlwpstra(i,k)=0.0
+! Convective Cloud Droplet Effective Radius (REFFCLWC) : variable 3D
+                  reffclwc(i,k)=1.1
+     &                 *((pqlwpcon(i,k)*pplay(i,k)/(RD * T(i,k)))
+     &                 /(4./3*RPI*1000.*cdnc(i,k)) )**(1./3.)
+                  reffclwc(i,k) = MAX(reffclwc(i,k) * 1e6, 5.)
+
+! Stratiform Cloud Droplet Effective Radius (REFFCLWS) : variable 3D
+                  IF ((pclc(i,k)-rnebcon(i,k)) .LE. seuil_neb) THEN ! tout sous la forme convective
+                     reffclws(i,k)=0.0
+                     lcc3dstra(i,k)= 0.0
+                  ELSE
+                     reffclws(i,k) = (pclc(i,k)*phase3d(i,k)*
+     &                               rad_chaud_tab(i,k)-
+     &                            pqlwpcon(i,k)*reffclwc(i,k))
+                     IF(reffclws(i,k) .LE. 0.0) reffclws(i,k)=0.0
+                     lcc3dstra(i,k)=pqlwpstra(i,k)
+                 ENDIF
+!Convertion from um to m
+                  IF(rnebcon(i,k). LE. seuil_neb) THEN
+                    reffclwc(i,k) = reffclwc(i,k)*seuil_neb*clwcon(i,k)
+     &                              *1.0E-06
+                    lcc3dcon(i,k)= seuil_neb*clwcon(i,k)
+                  ELSE
+                    reffclwc(i,k) = reffclwc(i,k)*pqlwpcon(i,k)
+     &                              *1.0E-06
+                    lcc3dcon(i,k) = pqlwpcon(i,k)
+                  ENDIF
+
+                  reffclws(i,k) = reffclws(i,k)*1.0E-06
+
+               ENDDO !klev
+            ENDDO !klon 
+
+!! Column Integrated Cloud Droplet Number (CLDNVI) : variable 2D
+            DO k = 1, klev
+               DO i = 1, klon
+                   pdel(i,k) = paprs(i,k)-paprs(i,k+1)
+                   zrho(i,k)=pplay(i,k)/t(i,k)/RD                  ! kg/m3
+                   dh(i,k)=pdel(i,k)/(gravit*zrho(i,k)) ! hauteur de chaque boite (m)
+               ENDDO
+            ENDDO
+c
+            DO i = 1, klon
+               cldnvi(i)=0.
+               lcc_integrat(i)=0.
+               height(i)=0.
+               DO k = 1, klev
+                  cldnvi(i)=cldnvi(i)+cdnc(i,k)*lcc3d(i,k)*dh(i,k)
+                  lcc_integrat(i)=lcc_integrat(i)+lcc3d(i,k)*dh(i,k)
+                  height(i)=height(i)+dh(i,k)
+               ENDDO ! klev
+               lcc_integrat(i)=lcc_integrat(i)/height(i)
+               IF (lcc_integrat(i) .LE. 1.0E-03) THEN
+                  cldnvi(i)=cldnvi(i)*lcc(i)/seuil_neb
+               ELSE
+                  cldnvi(i)=cldnvi(i)*lcc(i)/lcc_integrat(i)
+               ENDIF
+            ENDDO ! klon
+            
+            DO i = 1, klon
+               DO k = 1, klev
+                IF (scdnc(i,k) .LE. 0.0) scdnc(i,k)=0.0
+                IF (reffclws(i,k) .LE. 0.0) reffclws(i,k)=0.0
+                IF (reffclwc(i,k) .LE. 0.0) reffclwc(i,k)=0.0
+                IF (lcc3d(i,k) .LE. 0.0) lcc3d(i,k)=0.0
+                IF (lcc3dcon(i,k) .LE. 0.0) lcc3dcon(i,k)=0.0
+                IF (lcc3dstra(i,k) .LE. 0.0) lcc3dstra(i,k)=0.0
+               ENDDO
+               IF (reffclwtop(i) .LE. 0.0) reffclwtop(i)=0.0
+               IF (cldncl(i) .LE. 0.0) cldncl(i)=0.0
+               IF (cldnvi(i) .LE. 0.0) cldnvi(i)=0.0
+               IF (lcc(i) .LE. 0.0) lcc(i)=0.0
+            ENDDO
+
+         ENDIF !ok_aie
+      ENDIF !ok newmicro
+c
 C
       RETURN
