Index: /trunk/LMDZ.GENERIC/libf/phystd/callcorrk.F90
===================================================================
--- /trunk/LMDZ.GENERIC/libf/phystd/callcorrk.F90	(revision 1132)
+++ /trunk/LMDZ.GENERIC/libf/phystd/callcorrk.F90	(revision 1133)
@@ -697,5 +697,5 @@
             else
                do nw=1,L_NSPECTV
-                  stel_fract(nw)= stel(nw) * fract(ig)
+                  stel_fract(nw)= stel(nw) * fract(ig)*(1.-eclipse(ig))
                end do
             endif
Index: /trunk/LMDZ.GENERIC/libf/phystd/callkeys.h
===================================================================
--- /trunk/LMDZ.GENERIC/libf/phystd/callkeys.h	(revision 1132)
+++ /trunk/LMDZ.GENERIC/libf/phystd/callkeys.h	(revision 1133)
@@ -30,5 +30,6 @@
      &   , co2supsat                                                    &
      &   , cloudlvl                                                     &
-     &   , pceil                                                    
+     &   , pceil                                                        &
+     &   , rings_shadow                                         
 
       logical callrad,corrk,calldifv,UseTurbDiff                        &
@@ -36,5 +37,6 @@
      &   , season,diurnal,tlocked,lwrite                                &
      &   , callstats,calleofdump                                        &
-     &   , callgasvis,continuum,H2Ocont_simple,graybody
+     &   , callgasvis,continuum,H2Ocont_simple,graybody                 &
+     &   , rings_shadow
 
       logical enertest
Index: /trunk/LMDZ.GENERIC/libf/phystd/inifis.F
===================================================================
--- /trunk/LMDZ.GENERIC/libf/phystd/inifis.F	(revision 1132)
+++ /trunk/LMDZ.GENERIC/libf/phystd/inifis.F	(revision 1133)
@@ -153,4 +153,9 @@
          write(*,*) "tlocked = ",tlocked
 
+         write(*,*) "Saturn ring shadowing ?"
+         rings_shadow = .false.
+         call getin("rings_shadow", rings_shadow)
+         write(*,*) "rings_shadow = ", rings_shadow
+         
 ! Test of incompatibility:
 ! if tlocked, then diurnal should be false
Index: /trunk/LMDZ.GENERIC/libf/phystd/physiq.F90
===================================================================
--- /trunk/LMDZ.GENERIC/libf/phystd/physiq.F90	(revision 1132)
+++ /trunk/LMDZ.GENERIC/libf/phystd/physiq.F90	(revision 1133)
@@ -11,5 +11,5 @@
       use watercommon_h, only : RLVTT, Psat_water,epsi
       use gases_h
-      use radcommon_h, only: sigma
+      use radcommon_h, only: sigma, eclipse
       use radii_mod, only: h2o_reffrad, co2_reffrad, h2o_cloudrad
       use aerosol_mod
@@ -112,5 +112,4 @@
 !           Loops converted to F90 matrix format: J. Leconte (2012)
 !           No more ngridmx/nqmx, F90 commons and adaptation to parallel: A. Spiga (2012)
-!     
 !==================================================================
 
@@ -460,4 +459,7 @@
          ALLOCATE(fract(ngrid))
 
+         !! this is defined in radcommon_h
+         ALLOCATE(eclipse(ngrid))
+
 !        variables set to 0
 !        ~~~~~~~~~~~~~~~~~~
@@ -700,4 +702,17 @@
                ! WARNING: this function appears not to work in 1D
 
+           endif
+
+           !! Eclipse incoming sunlight (e.g. Saturn ring shadowing)
+           if (rings_shadow) then
+             if (diurnal) then
+               call rings(ngrid, declin, ptime, rad, 0., eclipse)
+               call writediagfi(ngrid,"shad","rings"," ", 2, 1.-eclipse)
+               write(*,*) "Ring shadow activated."
+             else
+               write(*,*) "Ring shadow + diurnal=F not supported yet."
+             endif
+           else
+               eclipse(:) = 0.0
            endif
 
@@ -1891,5 +1906,4 @@
         IF ( ALLOCATED(sinlon)) DEALLOCATE(sinlon)
         IF ( ALLOCATED(coslon)) DEALLOCATE(coslon)
-
       endif
 
Index: /trunk/LMDZ.GENERIC/libf/phystd/radcommon_h.F90
===================================================================
--- /trunk/LMDZ.GENERIC/libf/phystd/radcommon_h.F90	(revision 1132)
+++ /trunk/LMDZ.GENERIC/libf/phystd/radcommon_h.F90	(revision 1133)
@@ -135,4 +135,6 @@
       save Cmk
 
+      ! extinction of incoming sunlight (Saturn's rings, eclipses, etc...)
+      REAL, DIMENSION(:), ALLOCATABLE :: eclipse
 
       end module radcommon_h
Index: /trunk/LMDZ.GENERIC/libf/phystd/rings.F90
===================================================================
--- /trunk/LMDZ.GENERIC/libf/phystd/rings.F90	(revision 1133)
+++ /trunk/LMDZ.GENERIC/libf/phystd/rings.F90	(revision 1133)
@@ -0,0 +1,156 @@
+SUBROUTINE rings(ngrid, declin, ptime, rad, flat, eclipse)
+! Calculates Saturn's rings shadowing
+! Authors: M. Sylvestre, M. Capderou, S. Guerlet, A. Spiga
+
+    use comdiurn_h
+    use comgeomfi_h
+ 
+    implicit none   
+
+    INTEGER, INTENT(IN) :: ngrid  ! horizontal grid dimension
+    REAL, INTENT(IN) :: declin    ! latitude of the subsolar point
+    REAL, INTENT(IN) :: ptime     ! UTC time in sol fraction : ptime=0.5 at noon
+    REAL, INTENT(IN) :: rad       ! equatorial radius of the planet
+    REAL, INTENT(IN) :: flat      ! flattening of the planet 
+    REAL, DIMENSION(ngrid), INTENT(OUT) :: eclipse ! absorption of the light by the rings    
+    
+    REAL :: rpol   ! polar radius of the planet
+    REAL :: e      ! shape excentricity of the planet : (1-e*e) = (1-f)*(1-f)    
+    INTEGER, PARAMETER :: nb_a=1 ! number of subdivisions of the A ring
+    INTEGER, PARAMETER :: nb_b=3 ! number of subdivisions of the B ring
+    INTEGER, PARAMETER :: nb_c=1 ! number of subdivisions of the C ring
+    INTEGER :: i
+
+    ! arrays for the rings. TBD: dynamical?
+    REAL, DIMENSION(nb_a) :: A_Rint ! internal radii of the subdivisions of the A ring 
+    REAL, DIMENSION(nb_a) :: A_Rext ! external radii of the subdivisions of the A ring
+    REAL, DIMENSION(nb_b) :: B_Rint ! internal radii of the subdivisions of the B ring
+    REAL, DIMENSION(nb_b) :: B_Rext ! external radii of the subdivisions of the B ring 
+    REAL, DIMENSION(nb_c) :: C_Rint ! internal radii of the subdivisions of the C ring
+    REAL, DIMENSION(nb_c) :: C_Rext ! external radii of the subdivisions of the C ring 
+                
+    ! Opacities of the rings : for each one we can give different opacities for each part
+    REAL, DIMENSION(nb_a) :: tau_A ! opacity of the A ring
+    REAL, DIMENSION(nb_b) :: tau_B ! opacity of the B ring
+    REAL, DIMENSION(nb_c) :: tau_C ! opacity of the C ring
+                
+    ! Parameters used to calculate if a point is under a ring subdivision's shadow
+    REAL :: phi_S                             ! subsolar point longitude
+    REAL, PARAMETER :: pi=acos(-1.0)    
+    REAL, DIMENSION(:), ALLOCATABLE:: x, y, z ! cartesian coordinates of the points on the planet
+    REAL :: xs, ys, zs                        ! cartesian coordinates of the points of the subsolar point
+    REAL, DIMENSION(:), ALLOCATABLE :: k      ! parameter (?)
+    REAL, DIMENSION(:), ALLOCATABLE :: N      ! parameter to compute cartesian coordinates on a ellipsoidal planet
+    REAL, DIMENSION(:), ALLOCATABLE :: r      ! distance at which the incident ray of sun crosses the equatorial plane
+                                              ! measured from the center of the planet   
+    REAL :: Ns                                ! (same for the subsolar point)
+   
+    ! equinox --> no shadow (AS: why is this needed?)
+    if(declin .eq. 0.) then
+        eclipse = 0.
+        return 
+    endif 
+
+! 1) INITIALIZATION
+
+    ! Generic
+    rpol = (1.- flat)*rad
+    e = sqrt(2*flat - flat**2)
+    ALLOCATE(x(ngrid))
+    ALLOCATE(y(ngrid))
+    ALLOCATE(z(ngrid))
+    ALLOCATE(k(ngrid))
+    ALLOCATE(N(ngrid))
+    ALLOCATE(r(ngrid))
+    eclipse = 2000.
+
+    ! Size of the rings
+    A_Rint(1) = 2.01*rad
+    A_Rext(1) = 2.26*rad
+    B_Rint(1) = 1.55*rad
+    B_Rext(1) = 1.67*rad
+    B_Rint(2) = 1.67*rad
+    B_Rext(2) = 1.83*rad
+    B_Rint(3) = 1.83*rad
+    B_Rext(3) = 1.92*rad
+    C_Rint(1) = 1.21*rad
+    C_Rext(1) = 1.53*rad
+    
+    ! Opacities of the rings (TBD: update with most recent values?)
+    tau_A(1) = 0.4
+    tau_B(1) = 0.8
+    tau_B(2) = 2.
+    tau_B(3) = 1.4
+    tau_C(1) = 0.08
+
+    ! Convert to cartesian coordinates
+    N = rad/sqrt(1-(e**2)*sinlat**2)
+    x = N*coslat*coslon
+    y = N*coslat*sinlon
+    z = N*(1-e**2)*sinlat
+
+! 2) LOCATION OF THE SUBSOLAR POINT 
+ 
+    ! subsolar longitude is deduced from time fraction ptime
+    ! SG: the minus sign is important! ... otherwise subsolar point adopts a reverse rotation
+    phi_S = -(ptime - 0.5)*2.*pi 
+    write(*,*) 'subsol point coords : ', declin*180./pi, phi_S*180./pi
+
+    ! subsolar latitude is declin (declination of the sun)
+    ! now convert in cartesian coordinates : 
+    Ns = rad/sqrt(1-(e**2)*sin(declin)**2)
+    xs = Ns*cos(declin)*cos(phi_S)
+    ys = Ns*cos(declin)*sin(phi_S)
+    zs = Ns*(1-e**2)*sin(declin)
+
+! 3) WHERE DOES THE INCIDENT RAY OF SUN CROSS THE EQUATORIAL PLAN ?
+
+    k = -z/zs
+    r = (k*xs + x)**2 + (k*ys + y)**2 
+    r = sqrt(r)
+
+! 4) SO WHERE ARE THE SHADOW OF THESE RINGS ?
+
+    ! Summer hemisphere is not under the shadow of the rings
+    where(lati*declin .gt. 0.)
+       eclipse = 1000.
+    end where
+
+    ! No shadow of the rings by night
+    where(x*xs + y*ys + z*zs .lt. 0.)
+       eclipse = 1000.
+    end where
+
+    ! if the incident rays of sun cross a ring, there is a shadow
+    do i=1, nb_A 
+        where(r .ge. A_Rint(i) .and. r .le. A_Rext(i) .and. eclipse .ne. 1000.)
+            eclipse = 1. - exp(-tau_A(i)/cos(declin))
+        end where
+    end do 
+
+    do i=1, nb_B 
+        where(r .ge. B_Rint(i) .and. r .le. B_Rext(i) .and. eclipse .ne. 1000.)
+            eclipse = 1. - exp(-tau_B(i)/cos(declin))
+        end where
+    enddo
+    
+    do i=1, nb_C 
+        where(r .ge. C_Rint(i) .and. r .le. C_Rext(i) .and. eclipse .ne. 1000.)
+            eclipse = 1. - exp(-tau_C(i)/cos(declin))
+        end where
+    enddo
+
+    ! At the other places and the excluded ones, eclipse is 0. 
+    where(eclipse .eq. 2000. .or. eclipse .eq. 1000.)
+        eclipse = 0. 
+    end where 
+
+! 5) CLEAN THE PLACE
+    DEALLOCATE(x)
+    DEALLOCATE(y)
+    DEALLOCATE(z)
+    DEALLOCATE(k)
+    DEALLOCATE(N)
+    DEALLOCATE(r)
+
+END SUBROUTINE rings
Index: /trunk/LMDZ.UNIVERSAL/libf/phygeneric/rings.F90
===================================================================
--- /trunk/LMDZ.UNIVERSAL/libf/phygeneric/rings.F90	(revision 1133)
+++ /trunk/LMDZ.UNIVERSAL/libf/phygeneric/rings.F90	(revision 1133)
@@ -0,0 +1,1 @@
+link ../../../LMDZ.GENERIC/libf/phystd/rings.F90
