source: trunk/LMDZ.VENUS/libf/phyvenus/calc_rayleigh.F90 @ 3556

Last change on this file since 3556 was 2560, checked in by slebonnois, 3 years ago

SL: Implementation of SW computation based on generic model. Switch between this new SW module or old module that reads R. Haus tables implemented with a key (solarchoice)

File size: 2.7 KB
Line 
1      subroutine calc_rayleigh
2
3!==================================================================
4!     
5!     Purpose
6!     -------
7!     Average the Rayleigh scattering in each band, weighting the
8!     average by the blackbody function at temperature tstellar.
9!     Works for an arbitrary mix of gases (currently CO2, N2 and
10!     H2 are possible).
11!     
12!     Authors
13!     -------
14!     Robin Wordsworth (2010)
15!     Jeremy Leconte (2012): Added option for variable gas. Improved water rayleigh (Bucholtz 1995).
16!
17!     Called by
18!     ---------
19!     setspv.F
20!     
21!     Calls
22!     -----
23!     none
24!     
25!==================================================================
26
27      use radinc_h, only: L_NSPECTV
28      use radcommon_h, only: WAVEV, BWNV, DWNV, tstellar, tauray, scalep
29      use gases_h
30
31      implicit none
32#include "YOMCST.h"
33
34      real*8 wl
35      integer N,Nfine,ifine,igas
36      parameter(Nfine=500.0)
37      real*8 :: P0 = 9.423D+6   ! Rayleigh scattering reference pressure in pascals.   
38
39      real*8 tauvar,tausum,tausumvar,tauwei,tauweivar,bwidth,bstart
40      double precision df
41
42      real*8 tauconsti(ngasmx)
43      real*8 tauvari(ngasmx)
44
45      integer icantbewrong
46
47      ! tau0/p0=tau/p (Hansen 1974)
48      ! Then calculate tau0 = (8*pi^3*p_1atm)/(3*N0^2) * 4*delta^2/(RG*RMD*lambda^4)
49      ! where delta=n-1 and N0 is an amagat
50
51!! ADAPTED FOR FIXED VENUS ATMOSPHERE: 96.5% CO2, 3.5% N2
52!! using CO2 = 1 ; N2 = 2
53      tauconsti(1) = (8.7/RG)*1.527*scalep/P0
54      tauconsti(2) = (9.81/RG)*8.569E-3*scalep/(P0/93.0)
55 
56      do N=1,L_NSPECTV
57         
58         tausum = 0.0
59         tauwei = 0.0
60         tausumvar = 0.0
61         tauweivar = 0.0
62         bstart = 10000.0/BWNV(N+1)
63         bwidth = (10000.0/BWNV(N)) - (10000.0/BWNV(N+1))
64         do ifine=1,Nfine
65            wl=bstart+dble(ifine)*bwidth/Nfine
66
67            tauvar=0.0
68            tauvari(1) = (1.0+0.013/wl**2)/wl**4
69            tauvari(2) = (1.0+0.0113/wl**2+0.00013/wl**4)/wl**4
70            tauvar=tauconsti(1)*tauvari(1)*0.965 + tauconsti(2)*tauvari(2)*0.035
71
72            call blackl(dble(wl*1e-6),dble(tstellar),df)
73            df=df*bwidth/Nfine
74            tauwei=tauwei+df
75            tausum=tausum+tauvar*df
76         
77         enddo
78         TAURAY(N)=tausum/tauwei
79         ! we multiply by scalep here (100) because plev, which is used in optcv,
80         ! is in units of mBar, so we need to convert.
81
82      end do
83
84      IF (L_NSPECTV > 6) THEN
85        icantbewrong = L_NSPECTV-6
86        print*,'At 1 atm and lambda = ',WAVEV(icantbewrong),' um'
87        print*,'tau_R = ',TAURAY(icantbewrong)*1013.25
88        print*,'sig_R = ',TAURAY(icantbewrong)*RG*RMD*1.67e-27*100, &
89               'cm^2 molecule^-1'
90      ENDIF
91
92    end subroutine calc_rayleigh
Note: See TracBrowser for help on using the repository browser.