Index: LMDZ6/trunk/DefLists/field_def_lmdz.xml
===================================================================
--- LMDZ6/trunk/DefLists/field_def_lmdz.xml	(revision 5654)
+++ LMDZ6/trunk/DefLists/field_def_lmdz.xml	(revision 5655)
@@ -165,4 +165,6 @@
 	<field id="cloud_cover_sw"    long_name="Ecrad SW cloud cover"    unit="-" />
         <field id="cloud_cover_sw_s2"    long_name="Ecrad SW cloud cover 2e call"    unit="-" />
+	<field id="ZFLUX_DIR_SUN" long_name="ZFLUX_DIR_SUN" unit="W/m2" />
+        <field id="ZFLUX_DIR_SUN_s2" long_name="ZFLUX_DIR_SUN_s2" unit="W/m2" />
         <field id="sols_s2"    long_name="Solar rad. at surf. 2e call"    unit="W/m2" />
         <field id="sols0_s2"    long_name="Solar rad. at surf."    unit="W/m2" />
@@ -1026,4 +1028,8 @@
         <field id="rlu_s2"    long_name="LW upward radiation"    unit="W m-2" />
         <field id="rld_s2"    long_name="LW downward radiation"    unit="W m-2" />	
+	<field id="ZFLUX_DIR" long_name="Direct-beam shortwave radiation" unit="W m-2" />
+        <field id="ZFLUX_DIR_CLEAR" long_name="Direct-beam CS shortwave radiation" unit="W m-2" />
+        <field id="ZFLUX_DIR_s2" long_name="Direct-beam shortwave radiation" unit="W m-2" />
+        <field id="ZFLUX_DIR_CLEAR_s2" long_name="Direct-beam CS shortwave radiation" unit="W m-2" />
         <field id="kz"    long_name="Kz melange"    unit="m2/s" enabled="FALSE" />
         <field id="evu"    long_name="Eddy viscosity coefficient for Momentum Variables"    unit="m2 s-1" />
Index: LMDZ6/trunk/libf/phylmd/ecrad/lmdz/radiation_scheme_mod.f90
===================================================================
--- LMDZ6/trunk/libf/phylmd/ecrad/lmdz/radiation_scheme_mod.f90	(revision 5654)
+++ LMDZ6/trunk/libf/phylmd/ecrad/lmdz/radiation_scheme_mod.f90	(revision 5655)
@@ -160,6 +160,6 @@
 
 ! Direct component of surface flux into horizontal plane
-REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_DIR(KLON)
-REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_DIR_CLEAR(KLON)
+REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_DIR(KLON,KLEV+1)
+REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_DIR_CLEAR(KLON,KLEV+1)
 ! As PFLUX_DIR but into a plane perpendicular to the sun
 REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_DIR_INTO_SUN(KLON)
@@ -585,9 +585,11 @@
 !PFLUX_SW_UP_CLEAR_SURF(KIDIA:KFDIA) = flux%sw_up_clear(KIDIA:KFDIA,KLEV+1)
 !PFLUX_LW_UP_CLEAR_SURF(KIDIA:KFDIA) = flux%lw_up_clear(KIDIA:KFDIA,KLEV+1)
-PFLUX_DIR(KIDIA:KFDIA) = flux%sw_dn_direct(KIDIA:KFDIA,KLEV+1)
-PFLUX_DIR_CLEAR(KIDIA:KFDIA) = flux%sw_dn_direct_clear(KIDIA:KFDIA,KLEV+1)
+! Direct component of flux into horizontal plane
+PFLUX_DIR(KIDIA:KFDIA,:) = flux%sw_dn_direct(KIDIA:KFDIA,:)
+PFLUX_DIR_CLEAR(KIDIA:KFDIA,:) = flux%sw_dn_direct_clear(KIDIA:KFDIA,:)
 PFLUX_DIR_INTO_SUN(KIDIA:KFDIA) = 0.0_JPRB
 WHERE (PMU0(KIDIA:KFDIA) > EPSILON(1.0_JPRB))
-  PFLUX_DIR_INTO_SUN(KIDIA:KFDIA) = PFLUX_DIR(KIDIA:KFDIA) / PMU0(KIDIA:KFDIA)
+! Direct Surface component of flux into a plane perpendicular to the sun        
+  PFLUX_DIR_INTO_SUN(KIDIA:KFDIA) = PFLUX_DIR(KIDIA:KFDIA,KLEV+1) / PMU0(KIDIA:KFDIA)
 END WHERE
 ! Top-of-atmosphere downwelling flux
@@ -849,6 +851,6 @@
 
 ! Direct component of surface flux into horizontal plane
-REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_DIR(KLON)
-REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_DIR_CLEAR(KLON)
+REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_DIR(KLON,KLEV+1)
+REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_DIR_CLEAR(KLON,KLEV+1)
 ! As PFLUX_DIR but into a plane perpendicular to the sun
 REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_DIR_INTO_SUN(KLON)
@@ -1283,9 +1285,9 @@
 !PFLUX_SW_UP_CLEAR_SURF(KIDIA:KFDIA) = flux%sw_up_clear(KIDIA:KFDIA,KLEV+1)
 !PFLUX_LW_UP_CLEAR_SURF(KIDIA:KFDIA) = flux%lw_up_clear(KIDIA:KFDIA,KLEV+1)
-PFLUX_DIR(KIDIA:KFDIA) = flux%sw_dn_direct(KIDIA:KFDIA,KLEV+1)
-PFLUX_DIR_CLEAR(KIDIA:KFDIA) = flux%sw_dn_direct_clear(KIDIA:KFDIA,KLEV+1)
+PFLUX_DIR(KIDIA:KFDIA,:) = flux%sw_dn_direct(KIDIA:KFDIA,:)
+PFLUX_DIR_CLEAR(KIDIA:KFDIA,:) = flux%sw_dn_direct_clear(KIDIA:KFDIA,:)
 PFLUX_DIR_INTO_SUN(KIDIA:KFDIA) = 0.0_JPRB
 WHERE (PMU0(KIDIA:KFDIA) > EPSILON(1.0_JPRB))
-  PFLUX_DIR_INTO_SUN(KIDIA:KFDIA) = PFLUX_DIR(KIDIA:KFDIA) / PMU0(KIDIA:KFDIA)
+  PFLUX_DIR_INTO_SUN(KIDIA:KFDIA) = PFLUX_DIR(KIDIA:KFDIA,KLEV+1) / PMU0(KIDIA:KFDIA)
 END WHERE
 ! Top-of-atmosphere downwelling flux
Index: LMDZ6/trunk/libf/phylmd/phys_output_ctrlout_mod.F90
===================================================================
--- LMDZ6/trunk/libf/phylmd/phys_output_ctrlout_mod.F90	(revision 5654)
+++ LMDZ6/trunk/libf/phylmd/phys_output_ctrlout_mod.F90	(revision 5655)
@@ -826,4 +826,16 @@
 !AI Ecrad 3Deffect
 #ifdef CPP_ECRAD
+  TYPE(ctrl_out), SAVE :: o_ZFLUX_DIR_SUN = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
+    'ZFLUX_DIR_SUN', 'Direct radiation into a plane perpendicular to the sun', 'W m-2', (/ ('', i=1, 10) /))
+  TYPE(ctrl_out), SAVE :: o_ZFLUX_DIR_SUN_s2 = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
+    'ZFLUX_DIR_SUN_s2', 'Direct radiation into a plane perpendicular to the sun', 'W m-2', (/ ('', i=1, 10) /))
+  TYPE(ctrl_out), SAVE :: o_ZFLUX_DIR = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
+    'ZFLUX_DIR', 'Direct radiation', 'W m-2', (/ ('', i=1, 10) /))
+  TYPE(ctrl_out), SAVE :: o_ZFLUX_DIR_s2 = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
+    'ZFLUX_DIR_s2', 'Direct radiation', 'W m-2', (/ ('', i=1, 10) /))
+  TYPE(ctrl_out), SAVE :: o_ZFLUX_DIR_CLEAR = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
+    'ZFLUX_DIR_CLEAR', 'Direct radiation', 'W m-2', (/ ('', i=1, 10) /))
+  TYPE(ctrl_out), SAVE :: o_ZFLUX_DIR_CLEAR_s2 = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
+    'ZFLUX_DIR_CLEAR_s2', 'Direct radiation', 'W m-2', (/ ('', i=1, 10) /))
   TYPE(ctrl_out), SAVE :: o_cloud_cover_sw = ctrl_out((/ 11, 11, 10, 10, 10, 10, 11, 11, 11, 11/), &
     'cloud_cover_sw', 'Ecrad SW cloud cover', '-', (/ ('', i=1, 10) /))
Index: LMDZ6/trunk/libf/phylmd/phys_output_var_mod.f90
===================================================================
--- LMDZ6/trunk/libf/phylmd/phys_output_var_mod.f90	(revision 5654)
+++ LMDZ6/trunk/libf/phylmd/phys_output_var_mod.f90	(revision 5655)
@@ -179,4 +179,11 @@
   !$OMP THREADPRIVATE(cloud_cover_sw, cloud_cover_sw_s2)
 
+  ! Direct-beam shortwave into a horizontal plane
+  REAL, SAVE, ALLOCATABLE :: ZFLUX_DIR(:,:), ZFLUX_DIR_CLEAR(:,:), &
+         &     ZFLUX_DIR_s2(:,:), ZFLUX_DIR_CLEAR_s2(:,:)
+  !$OMP THREADPRIVATE(ZFLUX_DIR, ZFLUX_DIR_CLEAR, ZFLUX_DIR_s2, ZFLUX_DIR_CLEAR_s2)
+  REAL, SAVE, ALLOCATABLE :: ZFLUX_DIR_SUN(:), ZFLUX_DIR_SUN_s2(:)
+  !$OMP THREADPRIVATE(ZFLUX_DIR_SUN, ZFLUX_DIR_SUN_s2)
+
 CONTAINS
 
@@ -252,7 +259,13 @@
     ALLOCATE(icc3dstra(klon, klev))
 
-    ! cloud_cover_sw, cloud_cover_sw_s2 from Ecrad
+! cloud_cover_sw, cloud_cover_sw_s2 from Ecrad
     ALLOCATE(cloud_cover_sw(klon))
     ALLOCATE(cloud_cover_sw_s2(klon))
+    ALLOCATE(ZFLUX_DIR(klon,klev+1))
+    ALLOCATE(ZFLUX_DIR_CLEAR(klon,klev+1))
+    ALLOCATE(ZFLUX_DIR_s2(klon,klev+1))
+    ALLOCATE(ZFLUX_DIR_CLEAR_s2(klon,klev+1))
+    ALLOCATE(ZFLUX_DIR_SUN(klon))
+    ALLOCATE(ZFLUX_DIR_SUN_s2(klon))
 
   END SUBROUTINE phys_output_var_init
@@ -308,4 +321,6 @@
     !AI cloud_cover_sw, cloud_cover_sw_s2 from Ecrad
     DEALLOCATE(cloud_cover_sw, cloud_cover_sw_s2)
+    DEALLOCATE(ZFLUX_DIR, ZFLUX_DIR_CLEAR, ZFLUX_DIR_s2, ZFLUX_DIR_CLEAR_s2)
+    DEALLOCATE(ZFLUX_DIR_SUN, ZFLUX_DIR_SUN_s2)
 
   END SUBROUTINE phys_output_var_end
Index: LMDZ6/trunk/libf/phylmd/phys_output_write_mod.F90
===================================================================
--- LMDZ6/trunk/libf/phylmd/phys_output_write_mod.F90	(revision 5654)
+++ LMDZ6/trunk/libf/phylmd/phys_output_write_mod.F90	(revision 5655)
@@ -257,5 +257,7 @@
          o_LWupSFC_s2, o_LWdnSFC_s2, o_rlu_s2, o_rld_s2, &
          o_rlucs_s2, o_rldcs_s2, o_rsu_s2, o_rsd_s2, &
-         o_rsucs_s2, o_rsdcs_s2   
+         o_rsucs_s2, o_rsdcs_s2, o_ZFLUX_DIR, o_ZFLUX_DIR_s2, &
+         o_ZFLUX_DIR_CLEAR, o_ZFLUX_DIR_CLEAR_s2, &
+         o_ZFLUX_DIR_SUN, o_ZFLUX_DIR_SUN_s2   
 #endif
 
@@ -324,5 +326,6 @@
         sollwdown_s2,lwdn0_s2,lwup_s2,lwdn_s2,lwup0_s2,lwdn0_s2
     USE phys_output_var_mod, ONLY: cloud_cover_sw, &
-        cloud_cover_sw_s2
+        cloud_cover_sw_s2, ZFLUX_DIR, ZFLUX_DIR_s2, &
+        ZFLUX_DIR_CLEAR, ZFLUX_DIR_CLEAR_s2, ZFLUX_DIR_SUN, ZFLUX_DIR_SUN_s2
 #endif
 
@@ -1282,4 +1285,10 @@
        CALL histwrite_phy(o_rlucs_s2, lwup0_s2)
        CALL histwrite_phy(o_rldcs_s2, lwdn0_s2)
+       CALL histwrite_phy(o_ZFLUX_DIR, ZFLUX_DIR)
+       CALL histwrite_phy(o_ZFLUX_DIR_s2, ZFLUX_DIR_s2)
+       CALL histwrite_phy(o_ZFLUX_DIR_CLEAR, ZFLUX_DIR_CLEAR)
+       CALL histwrite_phy(o_ZFLUX_DIR_CLEAR_s2, ZFLUX_DIR_CLEAR_s2)
+       CALL histwrite_phy(o_ZFLUX_DIR_SUN, ZFLUX_DIR_SUN)
+       CALL histwrite_phy(o_ZFLUX_DIR_SUN_s2, ZFLUX_DIR_SUN_s2)
     ENDIF !ok_3Deffect
 #endif       
Index: LMDZ6/trunk/libf/phylmd/physiq_mod.F90
===================================================================
--- LMDZ6/trunk/libf/phylmd/physiq_mod.F90	(revision 5654)
+++ LMDZ6/trunk/libf/phylmd/physiq_mod.F90	(revision 5655)
@@ -94,5 +94,6 @@
     USE lmdz_cloud_optics_prop_ini, ONLY : cloud_optics_prop_ini
     USE phys_output_var_mod, ONLY :      cloudth_sth,cloudth_senv,cloudth_sigmath,cloudth_sigmaenv
-    USE phys_output_var_mod, ONLY : cloud_cover_sw, cloud_cover_sw_s2
+    USE phys_output_var_mod, ONLY : cloud_cover_sw, cloud_cover_sw_s2, ZFLUX_DIR, ZFLUX_DIR_CLEAR, &
+  &      ZFLUX_DIR_s2, ZFLUX_DIR_CLEAR_s2, ZFLUX_DIR_SUN, ZFLUX_DIR_SUN_s2
 
 
@@ -4648,4 +4649,5 @@
                ZLWFT0_i, ZFLDN0, ZFLUP0, &
                ZSWFT0_i, ZFSDN0, ZFSUP0, &
+               ZFLUX_DIR, ZFLUX_DIR_CLEAR, ZFLUX_DIR_SUN, &
                cloud_cover_sw)
 
@@ -4726,4 +4728,5 @@
                      ZLWFT0_i, ZFLDN0, ZFLUP0, &
                      ZSWFT0_i, ZFSDN0, ZFSUP0, &
+                     ZFLUX_DIR, ZFLUX_DIR_CLEAR, ZFLUX_DIR_SUN, &
                      cloud_cover_sw)
           ENDIF !ok_4xCO2atm
@@ -4772,4 +4775,5 @@
                      ZLWFT0_i, ZFLDN0, ZFLUP0, &
                      ZSWFT0_i, ZFSDN0, ZFSUP0, &
+                     ZFLUX_DIR, ZFLUX_DIR_CLEAR, ZFLUX_DIR_SUN, &
                      cloud_cover_sw_s2)
           ENDIF ! ok_3Deffect
Index: LMDZ6/trunk/libf/phylmd/radlwsw_m.F90
===================================================================
--- LMDZ6/trunk/libf/phylmd/radlwsw_m.F90	(revision 5654)
+++ LMDZ6/trunk/libf/phylmd/radlwsw_m.F90	(revision 5655)
@@ -45,4 +45,5 @@
        ZLWFT0_i, ZFLDN0, ZFLUP0, &
        ZSWFT0_i, ZFSDN0, ZFSUP0, &
+       ZFLUX_DIR, ZFLUX_DIR_CLEAR, ZFLUX_DIR_INTO_SUN, &
        cloud_cover_sw)
 
@@ -333,17 +334,20 @@
     REAL(KIND=8) ZAEROSOL_OLD(KLON,6,KLEV), &  ! 
          ZAEROSOL(KLON,KLEV,naero_spc) !
+    ! Interm
+    REAL(KIND=8), dimension(klon,klev+1) :: ZFLUX_DIR_i, &           ! Direct compt of surf flux into horizontal plane
+                                             ZFLUX_DIR_CLEAR_i        ! CS Direct
     ! OUTPUTS
-    REAL(KIND=8) ZFLUX_DIR(klon), &           ! Direct compt of surf flux into horizontal plane
-         ZFLUX_DIR_CLEAR(klon), &     ! CS Direct 
-         ZFLUX_DIR_INTO_SUN(klon), &  ! 
-         ZFLUX_UV(klon), &            ! UV flux
-         ZFLUX_PAR(klon), &           ! photosynthetically active radiation similarly
-         ZFLUX_PAR_CLEAR(klon), &     ! CS photosynthetically 
-         ZFLUX_SW_DN_TOA(klon), &     ! DN SW flux at TOA
-         ZEMIS_OUT(klon)              ! effective broadband emissivity
-
-    REAL(KIND=8) ZLWDERIVATIVE(klon,klev+1)   ! LW derivatives
-    REAL(KIND=8) ZSWDIFFUSEBAND(klon,NSW), &  ! SW DN flux in diffuse albedo band 
-         ZSWDIRECTBAND(klon,NSW)      ! SW DN flux in direct albedo band
+    REAL(KIND=8), dimension(klon,klev+1) :: ZFLUX_DIR, &           ! Direct compt of surf flux into horizontal plane
+                                                         ZFLUX_DIR_CLEAR        ! CS Direct 
+    REAL(KIND=8), dimension(klon)        ::  ZFLUX_DIR_INTO_SUN, &  ! 
+                                                 ZFLUX_UV, &            ! UV flux
+                                                 ZFLUX_PAR, &           ! photosynthetically active radiation similarly
+                                                 ZFLUX_PAR_CLEAR, &     ! CS photosynthetically 
+                                                 ZFLUX_SW_DN_TOA(klon), &     ! DN SW flux at TOA
+                                                 ZEMIS_OUT              ! effective broadband emissivity
+    REAL(KIND=8), dimension(klon,klev+1) ::  ZLWDERIVATIVE   ! LW derivatives
+
+    ! REAL(KIND=8) ZSWDIFFUSEBAND(klon,NSW), &  ! SW DN flux in diffuse albedo band 
+    !     ZSWDIRECTBAND(klon,NSW)      ! SW DN flux in direct albedo band
     REAL(KIND=8) SOLARIRAD
     REAL(KIND=8) seuilmach
@@ -1392,5 +1396,5 @@
                   & ZFSUP_i, ZFLUX_i(:,1,:), ZFCUP_i, ZFLUC_i(:,1,:), &
                   !       Surf Direct flux : ATTENTION
-                  & ZFLUX_DIR, ZFLUX_DIR_CLEAR, ZFLUX_DIR_INTO_SUN, &
+                  & ZFLUX_DIR_i, ZFLUX_DIR_CLEAR_i, ZFLUX_DIR_INTO_SUN, &
                   !       UV and para flux
                   & ZFLUX_UV, ZFLUX_PAR, ZFLUX_PAR_CLEAR, &
@@ -1436,5 +1440,5 @@
                   & ZFSUP_i, ZFLUX_i(:,1,:), ZFCUP_i, ZFLUC_i(:,1,:), &
                   !       Surf Direct flux : ATTENTION
-                  & ZFLUX_DIR, ZFLUX_DIR_CLEAR, ZFLUX_DIR_INTO_SUN, &
+                  & ZFLUX_DIR_i, ZFLUX_DIR_CLEAR_i, ZFLUX_DIR_INTO_SUN, &
                   !       UV and para flux
                   & ZFLUX_UV, ZFLUX_PAR, ZFLUX_PAR_CLEAR, &
@@ -1496,4 +1500,7 @@
                 ZFLDNC0(i,k+1)= -1.*ZFLCCDWN_i(i,klev+1-k)
                 ZFLUPC0(i,k+1)= ZFLCCUP_i(i,klev+1-k)
+                ! Direct flux
+                ZFLUX_DIR(i,k+1) = ZFLUX_DIR_i(i,klev+1-k)
+                ZFLUX_DIR_CLEAR  = ZFLUX_DIR_CLEAR_i(i,klev+1-k)
                 IF (ok_volcan) THEN
                    ZSWADAERO(i,k+1)=ZSWADAERO(i,klev+1-k)*fract(i) !--NL
