Index: LMDZ6/trunk/libf/phylmdiso/ice_sursat_mod.F90
===================================================================
--- LMDZ6/trunk/libf/phylmdiso/ice_sursat_mod.F90	(revision 4065)
+++ LMDZ6/trunk/libf/phylmdiso/ice_sursat_mod.F90	(revision 4065)
@@ -0,0 +1,1 @@
+link ../phylmd/ice_sursat_mod.F90
Index: LMDZ6/trunk/libf/phylmdiso/phys_local_var_mod.F90
===================================================================
--- LMDZ6/trunk/libf/phylmdiso/phys_local_var_mod.F90	(revision 4064)
+++ LMDZ6/trunk/libf/phylmdiso/phys_local_var_mod.F90	(revision 4065)
@@ -16,5 +16,9 @@
       REAL, SAVE, ALLOCATABLE :: u_seri(:,:), v_seri(:,:)
       !$OMP THREADPRIVATE(u_seri, v_seri)
-      REAL, SAVE, ALLOCATABLE :: l_mixmin(:,:,:),l_mix(:,:,:),tke_dissip(:,:,:), wprime(:,:,:)
+      REAL, SAVE, ALLOCATABLE :: rneb_seri(:,:)
+      !$OMP THREADPRIVATE(rneb_seri)
+      REAL, SAVE, ALLOCATABLE :: d_rneb_dyn(:,:)
+      !$OMP THREADPRIVATE(d_rneb_dyn)
+      REAL, SAVE, ALLOCATABLE :: l_mixmin(:,:,:),l_mix(:,:,:),tke_dissip(:,:,:),wprime(:,:,:)
       !$OMP THREADPRIVATE(l_mixmin, l_mix, tke_dissip,wprime)
       REAL, SAVE, ALLOCATABLE :: tr_seri(:,:,:)
@@ -602,4 +606,45 @@
 !$OMP THREADPRIVATE(zn2mout)
 
+      REAL, SAVE, ALLOCATABLE :: qclr(:,:)
+      !$OMP THREADPRIVATE(qclr)
+      REAL, SAVE, ALLOCATABLE :: qcld(:,:)
+      !$OMP THREADPRIVATE(qcld)
+      REAL, SAVE, ALLOCATABLE :: qss(:,:)
+      !$OMP THREADPRIVATE(qss)
+      REAL, SAVE, ALLOCATABLE :: qvc(:,:)
+      !$OMP THREADPRIVATE(qvc)
+      REAL, SAVE, ALLOCATABLE :: rnebclr(:,:)
+      !$OMP THREADPRIVATE(rnebclr)
+      REAL, SAVE, ALLOCATABLE :: rnebss(:,:)
+      !$OMP THREADPRIVATE(rnebss)
+      REAL, SAVE, ALLOCATABLE :: gamma_ss(:,:)
+      !$OMP THREADPRIVATE(gamma_ss)
+      REAL, SAVE, ALLOCATABLE :: N1_ss(:,:)
+      !$OMP THREADPRIVATE(N1_ss)
+      REAL, SAVE, ALLOCATABLE :: N2_ss(:,:)
+      !$OMP THREADPRIVATE(N2_ss)
+      REAL, SAVE, ALLOCATABLE :: drneb_sub(:,:)
+      !$OMP THREADPRIVATE(drneb_sub)
+      REAL, SAVE, ALLOCATABLE :: drneb_con(:,:)
+      !$OMP THREADPRIVATE(drneb_con)
+      REAL, SAVE, ALLOCATABLE :: drneb_tur(:,:)
+      !$OMP THREADPRIVATE(drneb_tur)
+      REAL, SAVE, ALLOCATABLE :: drneb_avi(:,:)
+      !$OMP THREADPRIVATE(drneb_avi)
+      REAL, SAVE, ALLOCATABLE :: zqsatl(:,:)
+      !$OMP THREADPRIVATE(zqsatl)
+      REAL, SAVE, ALLOCATABLE :: zqsats(:,:)
+      !$OMP THREADPRIVATE(zqsats)
+      REAL, SAVE, ALLOCATABLE :: Tcontr(:,:)
+      !$OMP THREADPRIVATE(Tcontr)
+      REAL, SAVE, ALLOCATABLE :: qcontr(:,:)
+      !$OMP THREADPRIVATE(qcontr)
+      REAL, SAVE, ALLOCATABLE :: qcontr2(:,:)
+      !$OMP THREADPRIVATE(qcontr2)
+      REAL, SAVE, ALLOCATABLE :: fcontrN(:,:)
+      !$OMP THREADPRIVATE(fcontrN)
+      REAL, SAVE, ALLOCATABLE :: fcontrP(:,:)
+      !$OMP THREADPRIVATE(fcontrP)
+
 #ifdef CPP_StratAer
 !
@@ -1037,4 +1082,14 @@
 
       ALLOCATE(zn2mout(klon,6))
+
+! Supersaturation
+      ALLOCATE(rneb_seri(klon,klev))
+      ALLOCATE(d_rneb_dyn(klon,klev))
+      ALLOCATE(qclr(klon,klev), qcld(klon,klev), qss(klon,klev), qvc(klon,klev))
+      ALLOCATE(rnebclr(klon,klev), rnebss(klon,klev), gamma_ss(klon,klev))
+      ALLOCATE(N1_ss(klon,klev), N2_ss(klon,klev))
+      ALLOCATE(drneb_sub(klon,klev), drneb_con(klon,klev), drneb_tur(klon,klev), drneb_avi(klon,klev))
+      ALLOCATE(zqsatl(klon,klev), zqsats(klon,klev))
+      ALLOCATE(Tcontr(klon,klev), qcontr(klon,klev), qcontr2(klon,klev), fcontrN(klon,klev), fcontrP(klon,klev))
 
 #ifdef CPP_StratAer
@@ -1390,4 +1445,14 @@
       DEALLOCATE(zn2mout)
 
+! Supersaturation
+      DEALLOCATE(rneb_seri)
+      DEALLOCATE(d_rneb_dyn)
+      DEALLOCATE(qclr, qcld, qss, qvc)
+      DEALLOCATE(rnebclr, rnebss, gamma_ss)
+      DEALLOCATE(N1_ss, N2_ss)
+      DEALLOCATE(drneb_sub, drneb_con, drneb_tur, drneb_avi)
+      DEALLOCATE(zqsatl, zqsats)
+      DEALLOCATE(Tcontr, qcontr, qcontr2, fcontrN, fcontrP)
+
 #ifdef CPP_StratAer
 ! variables for strat. aerosol CK
Index: LMDZ6/trunk/libf/phylmdiso/phys_output_ctrlout_mod.F90
===================================================================
--- LMDZ6/trunk/libf/phylmdiso/phys_output_ctrlout_mod.F90	(revision 4064)
+++ LMDZ6/trunk/libf/phylmdiso/phys_output_ctrlout_mod.F90	(revision 4065)
@@ -1940,4 +1940,52 @@
     'albslw3', 'Surface albedo LW3', '-', (/ ('', i=1, 10) /))
 
+!--aviation & supersaturation
+  TYPE(ctrl_out), SAVE :: o_oclr = ctrl_out((/ 1, 1, 1, 1, 10, 10, 11, 11, 11, 11/), &
+    'oclr', 'Clear sky total water', 'kg/kg', (/ ('', i=1, 10) /))
+  TYPE(ctrl_out), SAVE :: o_ocld = ctrl_out((/ 1, 1, 1, 1, 10, 10, 11, 11, 11, 11/), &
+    'ocld', 'Cloudy sky total water', 'kg/kg', (/ ('', i=1, 10) /))
+  TYPE(ctrl_out), SAVE :: o_oss = ctrl_out((/ 1, 1, 1, 1, 10, 10, 11, 11, 11, 11/), &
+    'oss', 'ISSR total water', 'kg/kg', (/ ('', i=1, 10) /))
+  TYPE(ctrl_out), SAVE :: o_ovc = ctrl_out((/ 1, 1, 1, 1, 10, 10, 11, 11, 11, 11/), &
+    'ovc', 'In-cloup vapor', 'kg/kg', (/ ('', i=1, 10) /))
+  TYPE(ctrl_out), SAVE :: o_rnebclr = ctrl_out((/ 1, 1, 1, 1, 10, 10, 11, 11, 11, 11/), &
+    'rnebclr', 'Clear sky fraction', '-', (/ ('', i=1, 10) /))
+  TYPE(ctrl_out), SAVE :: o_rnebss = ctrl_out((/ 1, 1, 1, 1, 10, 10, 11, 11, 11, 11/), &
+    'rnebss', 'ISSR fraction', '-', (/ ('', i=1, 10) /))
+  TYPE(ctrl_out), SAVE :: o_rnebseri = ctrl_out((/ 1, 1, 1, 1, 10, 10, 11, 11, 11, 11/), &
+    'rnebseri', 'Cloud fraction', '-', (/ ('', i=1, 10) /))
+  TYPE(ctrl_out), SAVE :: o_gammass = ctrl_out((/ 1, 1, 1, 1, 10, 10, 11, 11, 11, 11/), &
+    'gammass', 'Gamma supersaturation', '', (/ ('', i=1, 10) /))
+  TYPE(ctrl_out), SAVE :: o_N1_ss = ctrl_out((/ 1, 1, 1, 1, 10, 10, 11, 11, 11, 11/), &
+    'N1ss', 'N1', '', (/ ('', i=1, 10) /))
+  TYPE(ctrl_out), SAVE :: o_N2_ss = ctrl_out((/ 1, 1, 1, 1, 10, 10, 11, 11, 11, 11/), &
+    'N2ss', 'N2', '', (/ ('', i=1, 10) /))
+  TYPE(ctrl_out), SAVE :: o_drnebsub = ctrl_out((/ 1, 1, 1, 1, 10, 10, 11, 11, 11, 11/), &
+    'drnebsub', 'Cloud fraction change because of sublimation', 's-1', (/ ('', i=1, 10) /))
+  TYPE(ctrl_out), SAVE :: o_drnebcon = ctrl_out((/ 1, 1, 1, 1, 10, 10, 11, 11, 11, 11/), &
+    'drnebcon', 'Cloud fraction change because of condensation', 's-1', (/ ('', i=1, 10) /))
+  TYPE(ctrl_out), SAVE :: o_drnebtur = ctrl_out((/ 1, 1, 1, 1, 10, 10, 11, 11, 11, 11/), &
+    'drnebtur', 'Cloud fraction change because of turbulence', 's-1', (/ ('', i=1, 10) /))
+  TYPE(ctrl_out), SAVE :: o_drnebavi = ctrl_out((/ 1, 1, 1, 1, 10, 10, 11, 11, 11, 11/), &
+    'drnebavi', 'Cloud fraction change because of aviation', 's-1', (/ ('', i=1, 10) /))
+  TYPE(ctrl_out), SAVE :: o_qsatl = ctrl_out((/ 1, 1, 1, 1, 10, 10, 11, 11, 11, 11/), &
+    'qsatl', 'Saturation with respect to liquid water', '', (/ ('', i=1, 10) /))
+  TYPE(ctrl_out), SAVE :: o_qsats = ctrl_out((/ 1, 1, 1, 1, 10, 10, 11, 11, 11, 11/), &
+    'qsats', 'Saturation with respect to solid water', '', (/ ('', i=1, 10) /))
+  TYPE(ctrl_out), SAVE :: o_flight_m = ctrl_out((/ 1, 1, 1, 1, 10, 10, 11, 11, 11, 11/), &
+    'flightm', 'Flown meters', 'm/s/mesh', (/ ('', i=1, 10) /))
+  TYPE(ctrl_out), SAVE :: o_flight_h2o = ctrl_out((/ 1, 1, 1, 1, 10, 10, 11, 11, 11, 11/), &
+    'flighth2o', 'H2O flight emission', 'kg H2O/s/mesh', (/ ('', i=1, 10) /))
+  TYPE(ctrl_out), SAVE :: o_Tcontr = ctrl_out((/ 1, 1, 1, 1, 11, 11, 11, 11, 11, 11/),&
+    'Tcontr', 'Temperature threshold for contrail formation', 'K', (/ ('',i=1,10) /))
+  TYPE(ctrl_out), SAVE :: o_qcontr = ctrl_out((/ 1, 1, 1, 1, 11, 11, 11, 11, 11, 11/),&
+    'qcontr', 'Specific humidity threshold for contrail formation','Pa', (/ ('', i=1, 10) /))
+  TYPE(ctrl_out), SAVE :: o_qcontr2 = ctrl_out((/ 1, 1, 1, 1, 11, 11, 11, 11, 11, 11/),&
+    'qcontr2', 'Specific humidity threshold for contrail formation','kg/kg', (/ ('', i=1, 10) /))
+  TYPE(ctrl_out), SAVE :: o_fcontrN = ctrl_out((/ 2, 2, 2, 2, 2, 2, 11, 11, 11, 11/),&
+    'fcontrN', 'Fraction with non-persistent contrail in clear-sky', '-', (/ ('', i=1,10)/))
+  TYPE(ctrl_out), SAVE :: o_fcontrP = ctrl_out((/ 2, 2, 2, 2, 2, 2, 11, 11, 11, 11/),&
+    'fcontrP', 'Fraction with persistent contrail in ISSR', '-', (/ ('', i=1,10)/))
+
 !!!!!!!!!!!!! Sorties niveaux standards de pression NMC 
   TYPE(ctrl_out), SAVE :: o_tnondef = ctrl_out((/ 11, 11, 11, 11, 11, 11, 5, 5, 5, 11/), &
Index: LMDZ6/trunk/libf/phylmdiso/physiq_mod.F90
===================================================================
--- LMDZ6/trunk/libf/phylmdiso/physiq_mod.F90	(revision 4064)
+++ LMDZ6/trunk/libf/phylmdiso/physiq_mod.F90	(revision 4065)
@@ -73,4 +73,5 @@
     USE tracinca_mod, ONLY: config_inca
     USE tropopause_m,     ONLY: dyn_tropopause
+    USE ice_sursat_mod,  ONLY: flight_init, airplane
     USE vampir
     USE VERTICAL_LAYERS_MOD, ONLY: aps,bps, ap, bp
@@ -167,5 +168,5 @@
        ! [Variables internes non sauvegardees de la physique]
        ! Variables locales pour effectuer les appels en serie
-       t_seri,q_seri,ql_seri,qs_seri,u_seri,v_seri,tr_seri, &
+       t_seri,q_seri,ql_seri,qs_seri,u_seri,v_seri,tr_seri,rneb_seri, &
        ! Dynamic tendencies (diagnostics)
        d_t_dyn,d_q_dyn,d_ql_dyn,d_qs_dyn,d_u_dyn,d_v_dyn,d_tr_dyn,d_rneb_dyn, &
@@ -514,4 +515,6 @@
     INTEGER isol          ! indice de traceurs pour eau glace
     PARAMETER (isol=3)
+    INTEGER irneb         ! indice de traceurs pour fraction nuageuse LS (optional)
+    PARAMETER (irneb=4)   
     !
     !
@@ -1393,12 +1396,12 @@
        ENDIF
 
-       IF ((iflag_ice_sursat.GT.0).AND.(iflag_ice_thermo.EQ.0)) THEN
-          WRITE (lunout, *) ' iflag_ice_sursat=1 requires iflag_ice_thermo=1 as well'
+       IF (ok_ice_sursat.AND.(iflag_ice_thermo.EQ.0)) THEN
+          WRITE (lunout, *) ' ok_ice_sursat=y requires iflag_ice_thermo=1 as well'
           abort_message='see above'
           CALL abort_physic(modname,abort_message,1)
        ENDIF
 
-       IF ((iflag_ice_sursat.GT.0).AND.(nqo.NE.4)) THEN
-          WRITE (lunout, *) ' iflag_ice_sursat=1 requires 4 H2O tracers ', &
+       IF (ok_ice_sursat.AND.(nqo.NE.4)) THEN
+          WRITE (lunout, *) ' ok_ice_sursat=y requires 4 H2O tracers ', &
                '(H2Ov, H2Ol, H2Oi, rnebi) but nqo=', nqo, '. Might as well stop here.'
           abort_message='see above'
@@ -2403,4 +2406,7 @@
           ELSE IF (nqo.eq.3) THEN
              qs_seri(i,k) = qx(i,k,isol)
+          ELSE IF (nqo.eq.4) THEN
+             qs_seri(i,k) = qx(i,k,isol)
+             rneb_seri(i,k) = qx(i,k,irneb)
           ENDIF
        ENDDO
@@ -4620,5 +4626,5 @@
     IF (ok_new_lscp) THEN
 
-    CALL lscp(phys_tstep,missing_value,paprs,pplay, &
+    CALL lscp(phys_tstep,missing_val,paprs,pplay, &
          t_seri, q_seri,ptconv,ratqs, &
          d_t_lsc, d_q_lsc, d_ql_lsc, d_qi_lsc, rneb, rneb_seri, & 
@@ -4628,5 +4634,5 @@
          prfl, psfl, rhcl,  &
          zqasc, fraca,ztv,zpspsk,ztla,zthl,iflag_cld_th, &
-         iflag_ice_thermo, iflag_ice_sursat)
+         iflag_ice_thermo, ok_ice_sursat)
 
     ELSE