Index: LMDZ5/branches/testing/libf/phylmd/physiq.F
===================================================================
--- LMDZ5/branches/testing/libf/phylmd/physiq.F	(revision 1795)
+++ LMDZ5/branches/testing/libf/phylmd/physiq.F	(revision 1864)
@@ -29,4 +29,5 @@
       USE phys_state_var_mod ! Variables sauvegardees de la physique
       USE phys_output_var_mod ! Variables pour les ecritures des sorties
+      USE phys_output_write_mod
       USE fonte_neige_mod, ONLY  : fonte_neige_get_vars
       USE phys_output_mod
@@ -44,4 +45,5 @@
       use conf_phys_m, only: conf_phys
       use radlwsw_m, only: radlwsw
+      use phyaqua_mod, only: zenang_an
       USE control_mod
 #ifdef REPROBUS
@@ -49,4 +51,5 @@
 #endif
       USE indice_sol_mod
+      USE phytrac_mod, ONLY : phytrac
 
 !IM stations CFMIP
@@ -173,5 +176,4 @@
       PARAMETER (ok_region=.FALSE.)
 c======================================================================
-      real weak_inversion(klon),dthmin(klon)
       real seuil_inversion
       save seuil_inversion
@@ -184,5 +186,4 @@
       REAL zz,znum,zden
       REAL wmax_th(klon)
-      REAL zmax_th(klon)
       REAL tau_overturning_th(klon)
 
@@ -220,5 +221,5 @@
       REAL u(klon,klev)
       REAL v(klon,klev)
-      REAL t(klon,klev),theta(klon,klev),thetal(klon,klev)
+      REAL t(klon,klev),thetal(klon,klev)
 c thetal: ligne suivante a decommenter si vous avez les fichiers     MPL 20130625
 c fth_fonctions.F90 et parkind1.F90
@@ -227,5 +228,4 @@
       REAL qx(klon,klev,nqtot)
       REAL flxmass_w(klon,klev)
-      REAL omega(klon,klev) ! vitesse verticale en Pa/s
       REAL d_u(klon,klev)
       REAL d_v(klon,klev)
@@ -242,5 +242,4 @@
       real clw(klon,klev),elij(klon,klev,klev)
       real epmlmMm(klon,klev,klev),eplaMm(klon,klev)
-      real wdtrainA(klon,klev),wdtrainM(klon,klev)
 ! RomP <<<
 !IM definition dynamique o_trac dans phys_output_open
@@ -281,6 +280,4 @@
 c
 c
-c prw: precipitable water
-      real prw(klon)
 
       REAL convliq(klon,klev)  ! eau liquide nuageuse convective
@@ -295,6 +292,4 @@
 c flwp, fiwp = Liquid Water Path & Ice Water Path (kg/m2)
 c flwc, fiwc = Liquid Water Content & Ice Water Content (kg/kg)
-      REAL flwp(klon), fiwp(klon)
-      REAL flwc(klon,klev), fiwc(klon,klev)
       REAL flwp_c(klon), fiwp_c(klon)
       REAL flwc_c(klon,klev), fiwc_c(klon,klev)
@@ -302,5 +297,4 @@
       REAL flwc_s(klon,klev), fiwc_s(klon,klev)
 
-      REAL evap_pot(klon,nbsrf)
 
 cIM ISCCP simulator v3.4
@@ -528,10 +522,4 @@
 c$OMP THREADPRIVATE(itap)
 c
-      real slp(klon) ! sea level pressure
-c
-      REAL fevap(klon,nbsrf)
-      REAL fluxlat(klon,nbsrf)
-c
-      REAL qsol(klon)
       REAL,save ::  solarlong0
 c$OMP THREADPRIVATE(solarlong0)
@@ -544,6 +532,4 @@
 c
       INTEGER igwd,idx(klon),itest(klon)
-c
-      REAL agesno(klon,nbsrf)
 c
 c      REAL,allocatable,save :: run_off_lic_0(:)
@@ -565,5 +551,4 @@
 
       REAL mip(klon,klev)  ! mass flux shed by the adiab ascent at each level
-      REAL Vprecip(klon,klev+1)   ! precipitation vertical profile
 c
       REAL wape_prescr, fip_prescr
@@ -589,6 +574,4 @@
 c$OMP THREADPRIVATE(wake_s_min_lsp)
 
-      real ale_wake(klon)
-      real alp_wake(klon)
 
       real ok_wk_lsp(klon)
@@ -609,5 +592,4 @@
       REAL wake_dtPBL(klon,klev)      ! Wake : differential heating (wake - unpertubed) PBL
       REAL wake_dqPBL(klon,klev)      ! Wake : differential moistening (wake - unpertubed) PBL
-      REAL wake_omg(klon,klev)        ! Wake : velocity difference (wake - unpertubed)
       REAL wake_ddeltat(klon,klev),wake_ddeltaq(klon,klev)
       REAL wake_dp_deltomg(klon,klev) ! Wake : gradient vertical de wake_omg
@@ -615,5 +597,4 @@
 c
 cpourquoi y'a pas de save??
-      REAL wake_h(klon)               ! Wake : hauteur de la poche froide
 c
       INTEGER wake_k(klon)            ! Wake sommet
@@ -645,7 +626,5 @@
 cRR:fin declarations poches froides
 c=======================================================================================================
-
-      REAL zw2(klon,klev+1)
-      REAL fraca(klon,klev+1)        
+        
       REAL ztv(klon,klev),ztva(klon,klev)
       REAL zpspsk(klon,klev)
@@ -665,12 +644,6 @@
       real env_tke_max0(klon)                                ! TKE dans l'environnement au LCL 
 
-c---Spectre de thermiques de type 2 au LCL
-      real n2(klon),s2(klon)
-      real ale_bl_stat(klon)
-
 c---D\'eclenchement stochastique
       integer :: tau_trig(klon)
-      real proba_notrig(klon)
-      real random_notrig(klon)
 
 c--------Statistical Boundary Layer Closure: ALP_BL--------
@@ -680,10 +653,4 @@
       real env_tke_max(klon,klev)                            ! Profil de TKE dans l'environnement
 
-c---Fermeture statistique
-      real alp_bl_det(klon)                                     ! ALP d\'terministe du thermique unique
-      real alp_bl_fluct_m(klon)                                 ! ALP li\'ee aux fluctuations de flux de masse sous-nuageux
-      real alp_bl_fluct_tke(klon)                               ! ALP li\'ee aux fluctuations d'\'energie cin\'etique sous-nuageuse 
-      real alp_bl_conv(klon)                                    ! ALP li\'ee \`a grande \'echelle
-      real alp_bl_stat(klon)                                    ! ALP totale
 
 ccc fin nrlmd le 10/04/2012
@@ -697,12 +664,8 @@
 c Variables locales:
 c
-      REAL cdragh(klon) ! drag coefficient pour T and Q
-      REAL cdragm(klon) ! drag coefficient pour vent
 cAA
 cAA  Pour phytrac 
       REAL u1(klon)             ! vents dans la premiere couche U
       REAL v1(klon)             ! vents dans la premiere couche V
-
-      REAL zxffonte(klon), zxfqcalving(klon),zxfqfonte(klon)
 
 c@$$      LOGICAL offline           ! Controle du stockage ds "physique"
@@ -713,5 +676,4 @@
 ! RomP >>> 
       REAL beta_prec_fisrt(klon,klev) ! taux de conv de l'eau cond (fisrt)
-      REAL beta_prec(klon,klev)       ! taux de conv de l'eau cond (utilise)
 ! RomP <<<
       INTEGER       :: iii
@@ -722,23 +684,7 @@
 c
 cIM 050204 END
-      REAL evap(klon), devap(klon) ! evaporation et sa derivee
-      REAL sens(klon), dsens(klon) ! chaleur sensible et sa derivee
-
-      REAL bils(klon) ! bilan de chaleur au sol
-
-      REAL wfbilo(klon,nbsrf) ! bilan d'eau, pour chaque
-C                             ! type de sous-surface et pondere par la fraction
-      REAL wfbils(klon,nbsrf) ! bilan de chaleur au sol, pour chaque
-C                             ! type de sous-surface et pondere par la fraction
-      REAL slab_wfbils(klon)  ! bilan de chaleur au sol pour le cas de slab, sur les points d'ocean
-
-      REAL fder(klon)         
-      REAL ve(klon) ! integr. verticale du transport meri. de l'energie
-      REAL vq(klon) ! integr. verticale du transport meri. de l'eau
-      REAL ue(klon) ! integr. verticale du transport zonal de l'energie
-      REAL uq(klon) ! integr. verticale du transport zonal de l'eau
-c
-      REAL frugs(klon,nbsrf)
-      REAL zxrugs(klon) ! longueur de rugosite
+      REAL devap(klon) ! evaporation et sa derivee
+      REAL dsens(klon) ! chaleur sensible et sa derivee
+
 c
 c Conditions aux limites
@@ -804,13 +750,7 @@
       REAL diafra(klon,klev)  ! fraction nuageuse
       REAL cldliq(klon,klev)  ! eau liquide nuageuse
-      REAL cldfra(klon,klev)  ! fraction nuageuse
-      REAL cldtau(klon,klev)  ! epaisseur optique
-      REAL cldemi(klon,klev)  ! emissivite infrarouge
 c
 CXXX PB 
       REAL fluxq(klon,klev, nbsrf)   ! flux turbulent d'humidite
-      REAL fluxt(klon,klev, nbsrf)   ! flux turbulent de chaleur
-      REAL fluxu(klon,klev, nbsrf)   ! flux turbulent de vitesse u
-      REAL fluxv(klon,klev, nbsrf)   ! flux turbulent de vitesse v
 c
       REAL zxfluxt(klon, klev)
@@ -818,8 +758,5 @@
       REAL zxfluxu(klon, klev)
       REAL zxfluxv(klon, klev)
-CXXX
-c
-      REAL fsollw(klon, nbsrf)   ! bilan flux IR pour chaque sous surface
-      REAL fsolsw(klon, nbsrf)   ! flux solaire absorb. pour chaque sous surface
+
 c Le rayonnement n'est pas calcule tous les pas, il faut donc
 c                      sauvegarder les sorties du rayonnement
@@ -834,9 +771,7 @@
       REAL conv_q(klon,klev) ! convergence de l'humidite (kg/kg/s)
       REAL conv_t(klon,klev) ! convergence de la temperature(K/s)
-c
-      REAL cldl(klon),cldm(klon),cldh(klon) !nuages bas, moyen et haut
-      REAL cldt(klon),cldq(klon) !nuage total, eau liquide integree
-c
-      REAL zxtsol(klon), zxqsurf(klon), zxsnow(klon), zxfluxlat(klon)
+
+c
+      REAL zxsnow(klon)
       REAL zxsnow_dummy(klon)
 c
@@ -857,6 +792,5 @@
       REAL t_coup
       PARAMETER (t_coup=234.0)
-c
-      REAL zphi(klon,klev)
+
 cym A voir plus tard !!
 cym      REAL zx_relief(iim,jjmp1)
@@ -864,18 +798,12 @@
 c
 c Grandeurs de sorties
-      REAL s_pblh(klon), s_lcl(klon), s_capCL(klon)
-      REAL s_oliqCL(klon), s_cteiCL(klon), s_pblt(klon)
-      REAL s_therm(klon), s_trmb1(klon), s_trmb2(klon)
+      REAL s_capCL(klon)
+      REAL s_oliqCL(klon), s_cteiCL(klon)
+      REAL s_trmb1(klon), s_trmb2(klon)
       REAL s_trmb3(klon)
 cKE43
 c Variables locales pour la convection de K. Emanuel (sb):
-c
-      REAL upwd(klon,klev)      ! saturated updraft mass flux
-      REAL dnwd(klon,klev)      ! saturated downdraft mass flux
-      REAL dnwd0(klon,klev)     ! unsaturated downdraft mass flux
+
       REAL tvp(klon,klev)       ! virtual temp of lifted parcel
-      REAL plcl(klon)           ! Lifting Condensation Level
-      REAL plfc(klon)           ! Level of Free Convection
-      REAL wbeff(klon)          ! saturated updraft velocity at LFC
       CHARACTER*40 capemaxcels  !max(CAPE)
 
@@ -902,5 +830,4 @@
 c eva: evaporation de l'eau liquide nuageuse
 c vdf: couche limite (Vertical DiFfusion)
-      REAL rneb(klon,klev)
 
 ! tendance nulles
@@ -914,13 +841,7 @@
 cIM 081204 END
 c
-      REAL pmfu(klon,klev), pmfd(klon,klev)
       REAL pen_u(klon,klev), pen_d(klon,klev)
       REAL pde_u(klon,klev), pde_d(klon,klev)
       INTEGER kcbot(klon), kctop(klon), kdtop(klon)
-      REAL pmflxr(klon,klev+1), pmflxs(klon,klev+1)
-      REAL prfl(klon,klev+1), psfl(klon,klev+1)
-c
-      REAL rain_lsc(klon)
-      REAL snow_lsc(klon)
 c
       REAL ratqsc(klon,klev)
@@ -935,5 +856,4 @@
       logical ok_newmicro
       save ok_newmicro
-      real ref_liq(klon,klev), ref_ice(klon,klev)
 c$OMP THREADPRIVATE(ok_newmicro)
       save fact_cldcon,facttemps
@@ -973,8 +893,5 @@
 c Variables locales pour effectuer les appels en serie
 c
-      REAL zx_rh(klon,klev)
 cIM RH a 2m (la surface)
-      REAL rh2m(klon), qsat2m(klon)
-      REAL tpot(klon), tpote(klon)
       REAL Lheat
 
@@ -1091,16 +1008,11 @@
       REAL ZRCPD
 c-jld ec_conser
-      REAL t2m(klon,nbsrf)  ! temperature a 2m
       REAL q2m(klon,nbsrf)  ! humidite a 2m
 
 cIM: t2m, q2m, ustar, u10m, v10m et t2mincels, t2maxcels
-      REAL zt2m(klon), zq2m(klon)             !temp., hum. 2m moyenne s/ 1 maille
-      REAL zustar(klon),zu10m(klon), zv10m(klon)  ! u* et vents a 10m moyennes s/1 maille
       CHARACTER*40 t2mincels, t2maxcels       !t2m min., t2m max
       CHARACTER*40 tinst, tave, typeval
       REAL cldtaupi(klon,klev)  ! Cloud optical thickness for pre-industrial (pi) aerosols
 
-      REAL re(klon, klev)       ! Cloud droplet effective radius
-      REAL fl(klon, klev)  ! denominator of re
 
       REAL re_top(klon), fl_top(klon) ! CDR at top of liquid water clouds
@@ -1228,6 +1140,5 @@
 ! Utilise notamment en 1D mais peut etre active egalement en 3D
 ! en imposant la valeur de igout.
-c======================================================================
-
+c======================================================================d
       if (prt_level.ge.1) then
           igout=klon/2+1/klon
@@ -1250,15 +1161,4 @@
 
 c======================================================================
-
-cym => necessaire pour iflag_con != 2    
-      pmfd(:,:) = 0.
-      pen_u(:,:) = 0.
-      pen_d(:,:) = 0.
-      pde_d(:,:) = 0.
-      pde_u(:,:) = 0.
-      aam=0.
-
-      torsfc=0.
-      forall (k=1: llm) zmasse(:, k) = (paprs(:, k)-paprs(:, k+1)) / rg
 
       if (first) then 
@@ -1301,4 +1201,17 @@
 
       endif  ! first
+
+cym => necessaire pour iflag_con != 2    
+      pmfd(:,:) = 0.
+      pen_u(:,:) = 0.
+      pen_d(:,:) = 0.
+      pde_d(:,:) = 0.
+      pde_u(:,:) = 0.
+      aam=0.
+
+      torsfc=0.
+      forall (k=1: llm) zmasse(:, k) = (paprs(:, k)-paprs(:, k+1)) / rg
+
+
 
        modname = 'physiq'
@@ -1600,5 +1513,5 @@
        call phys_output_open(rlon,rlat,nCFMIP,tabijGCM,
      &                       iGCM,jGCM,lonGCM,latGCM,
-     &                       jjmp1,nlevSTD,clevSTD,
+     &                       jjmp1,nlevSTD,clevSTD,rlevSTD,
      &                       nbteta, ctetaSTD, dtime,ok_veget,
      &                       type_ocean,iflag_pbl,ok_mensuel,ok_journe,
@@ -1606,16 +1519,13 @@
      &                       read_climoz, phys_out_filestations,
      &                       new_aod, aerosol_couple,
-     &                       flag_aerosol_strat )
+     &                       flag_aerosol_strat, pdtphys, paprs, pphis, 
+     &                       pplay, lmax_th, ptconv, ptconvth, ivap, 
+     &                       d_t, qx, d_qx, zmasse, ok_sync)
 c$OMP END MASTER
 c$OMP BARRIER
 
+#undef histISCCP
 #ifdef histISCCP
 #include "ini_histISCCP.h"
-#endif
-
-#ifdef histNMC
-#include "ini_histhfNMC.h"
-#include "ini_histdayNMC.h"
-#include "ini_histmthNMC.h"
 #endif
 
@@ -1969,6 +1879,13 @@
       DO i = 1, klon
          zlvdcp=RLVTT/RCPD/(1.0+RVTMP2*q_seri(i,k))
-c        zlsdcp=RLSTT/RCPD/(1.0+RVTMP2*q_seri(i,k))
-         zlsdcp=RLVTT/RCPD/(1.0+RVTMP2*q_seri(i,k))
+cjyg<
+c      Attention : Arnaud a propose des formules completement differentes
+c                  A verifier !!!
+         zlsdcp=RLSTT/RCPD/(1.0+RVTMP2*q_seri(i,k))
+         IF (iflag_ice_thermo .EQ. 0) THEN
+           zlsdcp=zlvdcp
+         ENDIF 
+c>jyg
+
          zdelta = MAX(0.,SIGN(1.,RTT-t_seri(i,k)))
          zb = MAX(0.0,ql_seri(i,k))
@@ -2083,5 +2000,5 @@
      e     t_seri,    q_seri,    u_seri,  v_seri,   
      e     pplay,     paprs,     pctsrf,            
-     +     ftsol,     falb1,     falb2,   ustar, u10m,   v10m,
+     +     ftsol,falb1,falb2,ustar,u10m,v10m,wstar,
      s     sollwdown, cdragh,    cdragm,  u1,    v1,
      s     albsol1,   albsol2,   sens,    evap,  
@@ -2354,10 +2271,12 @@
                 nbtr_tmp=nbtr
              END IF
-          CALL concvl (iflag_con,iflag_clos,
+cjyg   iflag_con est dans clesphys
+cc          CALL concvl (iflag_con,iflag_clos,
+          CALL concvl (iflag_clos,
      .        dtime,paprs,pplay,t_undi,q_undi,
      .        t_wake,q_wake,wake_s,
      .        u_seri,v_seri,tr_seri,nbtr_tmp,
      .        ALE,ALP,
-     .        ema_work1,ema_work2,
+     .        sig1,w01,
      .        d_t_con,d_q_con,d_u_con,d_v_con,d_tr,
      .        rain_con, snow_con, ibas_con, itop_con, sigd,
@@ -2392,5 +2311,5 @@
      .        paprs,pplay,t_seri,q_seri,
      .        u_seri,v_seri,tr_seri,ntra,
-     .        ema_work1,ema_work2,
+     .        sig1,w01,
      .        d_t_con,d_q_con,d_u_con,d_v_con,d_tr,
      .        rain_con, snow_con, ibas_con, itop_con,
@@ -2856,5 +2775,6 @@
      .           frac_impa, frac_nucl, beta_prec_fisrt,
      .           prfl, psfl, rhcl, 
-     .           zqasc, fraca,ztv,zpspsk,ztla,zthl,iflag_cldcon )
+     .           zqasc, fraca,ztv,zpspsk,ztla,zthl,iflag_cldcon,
+     .           iflag_ice_thermo)
 
       WHERE (rain_lsc < 0) rain_lsc = 0.
@@ -3729,7 +3649,8 @@
      I     paprs,    pplay,     pmfu,     pmfd,
      I     pen_u,    pde_u,     pen_d,    pde_d,
-     I     cdragh,   coefh(:,:,is_ave),     fm_therm, entr_therm,
+     I     cdragh,   coefh(:,:,is_ave),   fm_therm, entr_therm,
      I     u1,       v1,        ftsol,    pctsrf,
-     I     ustar,     u10m,      v10m,
+     I     zustar,   zu10m,     zv10m,
+     I     wstar(:,is_ave),    ale_bl,         ale_wake,
      I     rlat,     rlon,
      I     frac_impa,frac_nucl, beta_prec_fisrt,beta_prec,
@@ -4031,5 +3952,13 @@
       endif
 
-#include "phys_output_write_new.h"
+
+!On effectue les sorties:
+
+      CALL phys_output_write(itap, pdtphys, paprs, pphis,               &
+     &                  pplay, lmax_th, aerosol_couple,                 &
+     &                  ok_ade, ok_aie, ivap, new_aod, ok_sync,         &
+     &                  ptconv, read_climoz, clevSTD, freq_moyNMC,      &
+     &                  ptconvth, d_t, qx, d_qx, zmasse,                &
+     &                  flag_aerosol_strat)
 
 
@@ -4038,10 +3967,4 @@
 #ifdef histISCCP
 #include "write_histISCCP.h"
-#endif
-
-#ifdef histNMC
-#include "write_histhfNMC.h"
-#include "write_histdayNMC.h"
-#include "write_histmthNMC.h"
 #endif
 
