[2560] | 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 |
---|