Rayleigh Scattering routine in Visual Basic

(c) 2013 by Barton Paul Levenson

Language:
Basic

Dialect:
Microsoft Visual Basic Express 2010

Discussion:
This routine finds the reflectivity in a given wavelength band.

Input:

lo: low band boundary in microns
hi: high band boundary in microns
m: mass of your layer, kg
h: height of your layer, m
P: local pressure, Pa
T: local temperature, K
Pw: local water vapor partial pressure, Pa

Output:

R: the reflectance of the layer.

I posted an incorrect version of this in C# in 2013. It gave bogus results because I was mixing cgs and SI units. Apologies! I will take that page down. Use this one. Translate between languages as needed.

    Public Const amu As Double = 1.66053892E-27#        ' kg
    Public Const area As Double = 510066000000000.0#    ' m^2
    Public Const mAtm As Double = 5.148E+18#            ' kg
    Public Const MW As Double = 28.94# * amu		' kg



    ' Rayleigh is after Houghton 1977 p. 73, but with SI units.
    Function Rayleigh(lo As Double, hi As Double, m As Double, _
    h As Double, P As Double, T As Double, Pw As Double)
        Const k As Double = 13.61434846#

        Dim lambda As Double = 0.5# * (lo + hi)
        Dim volume As Double = h * area
        Dim rho As Double = m / volume
        Dim Nprime As Double = rho / MW
        Dim n As Double = refrac(lambda, P, T, Pw)

        Dim term1 As Double = k / (Nprime * lambda ^ 4)
        Dim term2 As Double = (n - 1.0#) ^ 2
        Dim sigma As Double = term1 * term2
        Dim path As Double = Nprime * h
        Dim tau As Double = 1.66# * sigma * path

        Return (1.0# - Exp(-tau))
    End Function ' Rayleigh4



    ' refrac returns the refractive index for air, Edlen 1966.
    Function refrac(lambda As Double, P As Double, T As Double, _
    Pw As Double) As Double
        Const a0 As Double = 83.42#
        Const a1 As Double = 185.08#
        Const a2 As Double = 4.11#
        Const b1 As Double = 114000.0#
        Const b2 As Double = 62400.0#
        Const c0 As Double = 43.49#
        Const c1 As Double = 17000.0#

        Const P0 As Double = 101325.0#
        Const T0 As Double = 288.15#

        Dim nu As Double = 10000.0# / lambda    ' Find wavenumber.
        Dim Ns As Double = a0 + a1 / (1.0# - (nu / b1) ^ 2)
        Ns += a2 / (1.0# - (nu / b2) ^ 2)
        Ns *= (P / P0) * (T0 / T)
        Ns -= (c0 - (nu / c1) ^ 2) * Pw / P0

        Return (0.000001# * Ns + 1.0#)
    End Function ' refrac