vSmartMOM

Appearance

Absorption: Spectral Line Shapes ​

The purpose of this tutorial is to demonstrate how to compute line-shapes using the vSmartMOM Julia package. Here, we focus on rotational-vibrational transition lines, what processes determine line-shapes, and how the line-shapes depend on pressure and temperature.

Basic Tools ​

We are making use of the HITRAN database, see a list of tutorials here.

There are other spectroscopic linelists, but we are focusing on HITRAN only here. In addition, there are more complex line-shapes that we are not treating here, including collisional narrowing and line-mixing, collision-induced absorption (CIA). (The RT pipeline can add CIA and the MT_CKD water-vapour continuum through the absorption block's cia_filesandmtckd_file keys; they are not part of the line-shape calculation shown here.)

Line-shape: ​

If we have a line-strength S (in cm cm/molecule) for a specific transition at , we can compute the cross section as:

where   denotes the line-shape function (in 1/cm), which is normalized to 1:

There are several processes that affect the shape and width of , and we will walk through the most important ones here.

Doppler Broadening ​

Doppler broadening is caused by a simple doppler shift of emitted and absorbed frequencies, caused by the relative velocities of the molecules along the line of sight. A doppler shifted apparent frequency from the centroid frequency can be described as:

where is the relative velocity of the absorbing photon along the line of sight. The doppler shift is then simply

Let us take a simple example with the satellite flying at 7km/s and either staring into the flight direction (technically, it wouldn't see the atmosphere then, but let's ignore this) or looking into the back:

julia
# Speed of light (in m/s)
c = 299792458.0

# Relative velocity (in m/s)
vᵣ = 7000.0

# Center wavenumber (say 1600nm, which is 1e7/1600=6250 cm$^{-1}$)
v₀ = 6250.0

# Doppler shift:
Δ_ν = ( v₀ * vᵣ) / c

# Just writing out doppler shift in wavenumbers and wavelengths
println("Doppler shift = $Δ_ν cm-1")
println("Doppler shift = $(1e7/(v₀-Δ_ν)-1e7/v₀) nm")

Random movements of molecules in gases lead to doppler broadening effects ​

In the one-dimensional case – say along the x-axis –, the speed of molecules is distributed according to the Maxwell-Boltzmann distribution:

where is the particle mass, is the Boltzmann constant, and T is thermodynamic temperature.

We can thus define a Doppler width as

which yields the following line-shape:

where is the 1/e half-width, not the standard deviation. The standard deviation is  , and HWHM is  . See the HITRAN width definitions.

Let us put in some numbers with R=8.3144598 J/K/mol at 6000cm:

T = 220K, 290K

M = 16g/mol (CH) or 44g/mol (CO)

Multiply by about 1.6651 (2) to get the full-width half-maximum (FWHM) of the spectral line.

Natural broadening ​

For an isolated transition with upper-state population lifetime and a negligibly narrow lower state, the natural Lorentz FWHM in ordinary frequency is ; divide by the speed of light in cm/s for cm⁻¹. Natural broadening is usually small for atmospheric infrared lines because radiative lifetimes are long compared with collisional dephasing times. Doppler broadening comes from the velocity distribution, not a shorter lifetime.

Collisional broadening ​

Collisions between molecules quench the excited state and thus reduce the lifetime of the upper state, resulting in a widening of the line width. This behavor gives rise to the so-called Lorentz lineshape.

depends linearly on the number density of the perturbing molecules and the relative speed of the collision partners (thus scales linearly with pressure and with , basically both density and velocity affect the number of collisions per time). We often call this type of broadening pressure broadening, caused by collisions between molecules or atoms, which can supply or remove small amounts of energy during radiative transitions, thereby allowing photons with a broader range of frequencies to produce a particular transition of a molecule.

Voigt lineshape ​

The Voigt line shape is the convolution of Gaussian Doppler and Lorentz profiles. It can be expressed using the complex error (Faddeeva) function; numerical approximations evaluate that function efficiently.

Other more complex lineshapes ​

Once you dig deeper, there are various other more complex line-shapes (and line-mixing effects), which we ignore for now as the Voigt line-shape can provide very reasonable results. See, for instance, here.

The generated docs include a compact Plotly preview of the same line-shape families. The Voigt trace is computed as a numerical convolution of the Doppler and Lorentz profiles on the displayed grid.

Import the required tools:

julia
using CairoMakie
using Pkg.Artifacts
using vSmartMOM
using vSmartMOM.Absorption

Deprecation notice: The functions below (read_hitran, make_hitran_model, absorption_cross_section) belong to the internal vSmartMOM.Absorption module. Its standalone line-by-line model (make_hitran_model, HitranModel) is deprecated in 2.2 and will be removed in 3.0; its first use logs a one-time notice with the corresponding AtmosphericAbsorption.jl calls and their convention differences. AtmosphericAbsorption.jl is the supported absorption dependency going forward; use it for new direct σ(ν, T, p), custom-grid, ABSCO HDF, and GPU workflows. The RT pipeline (model_from_parameters → rt_run) already uses AtmosphericAbsorption.jl. Tables for YAML LUTfiles are still built with make_interpolation_model in 2.x (or saved as AtmosphericAbsorption .absco tables).

Read the HITRAN database for CO, using artifacts in Julia:

julia
co2_par      = vSmartMOM.Absorption.read_hitran(artifact("CO2"), mol=2, iso=1, ν_min=6214.4, ν_max=6214.8);
nothing #hide

Create a Voigt model (all in CPU mode)

julia
line_voigt   = make_hitran_model(co2_par, Voigt(), architecture=CPU())

Create a Doppler model

julia
line_doppler = make_hitran_model(co2_par, Doppler(), architecture=CPU())

Create a Lorentz model

julia
line_lorentz = make_hitran_model(co2_par, Lorentz(), architecture=CPU())
# Specify our wavenumber grid
ν = 6210:0.001:6220;
nothing #hide

Here we compute the cross sections at different pressures (in hPa, 3rd argument) and temperatures (in K, 4th argument)

julia
cs_co2_1atm   = absorption_cross_section(line_voigt, ν, 1013.0     , 296.0);
cs_co2_075atm = absorption_cross_section(line_voigt, ν, 0.75*1013.0, 296.0);
cs_co2_05atm  = absorption_cross_section(line_voigt, ν, 0.5*1013.0 , 296.0);
cs_co2_025atm = absorption_cross_section(line_voigt, ν, 0.25*1013.0, 296.0);
cs_co2_01atm  = absorption_cross_section(line_voigt, ν, 0.1*1013.0 , 296.0);
nothing #hide

Get some more line-shapes just for Doppler and Voigt to better compare the different line shapes

julia
cs_co2_01atm    = absorption_cross_section(line_voigt,   ν, 0.1*1013.0 , 296.0);
cs_co2_doppler  = absorption_cross_section(line_doppler, ν, 0.1*1013.0 , 296.0);
cs_co2_lorentz  = absorption_cross_section(line_lorentz, ν, 0.1*1013.0 , 296.0);
nothing #hide

Using a scaling factor for display purposes:

julia
ff = 1e20;
fig = Figure(size=(700, 450))
ax = Axis(fig[1,1],
    xlabel = "Wavenumber (cm⁻¹)",
    ylabel = "Absorption cross section (10⁻²⁰ cm²/molecule)")
lines!(ax, collect(ν), ff .* cs_co2_1atm,   label="Voigt, 1 atm")
lines!(ax, collect(ν), ff .* cs_co2_075atm, label="Voigt, 0.75 atm")
lines!(ax, collect(ν), ff .* cs_co2_05atm,  label="Voigt, 0.5 atm")
lines!(ax, collect(ν), ff .* cs_co2_025atm, label="Voigt, 0.25 atm")
xlims!(ax, 6214, 6215.2)
axislegend(ax, position=:rt)
fig

From these figures we can see how strongly the line-width is impacted by pressure broadening at ranges from about 250 hPa to 1013 hPa. Especially the line wings are substantially widened, which is typical of Lorentz line-shapes.

julia
fig = Figure(size=(700, 450))
ax = Axis(fig[1,1],
    xlabel = "Wavenumber (cm⁻¹)",
    ylabel = "Normalized absorption cross section")
lines!(ax, collect(ν), cs_co2_01atm   ./ maximum(cs_co2_01atm), label="Voigt, 296K, 0.1 atm")
lines!(ax, collect(ν), cs_co2_doppler ./ maximum(cs_co2_01atm), label="Doppler, 296K")
lines!(ax, collect(ν), cs_co2_lorentz ./ maximum(cs_co2_01atm), label="Lorentz, 296K, 0.1 atm")
xlims!(ax, 6214.4, 6214.8)
axislegend(ax, position=:rt)
fig

We can see the individual impacts of Doppler and Pressure broadening at about 100 hPa, at which Doppler broadening becomes important too. Most of the shape of the wings is still dominated by pressure broadening but the center has substantial contributions from Doppler broadening as well. The Voigt shape shown here is the convolution of the Lorentz and Doppler shape.

From an individual line to a band ​

Here, we will compute an entire band of CO (a few to be precise) and look at some simple behavior, e.g. the re-distribution of individual lines in the P and R branch with changing temperature.

Load a wider spectral range

julia
co2_par_band = vSmartMOM.Absorption.read_hitran(artifact("CO2"), mol=2, iso=1, ν_min=6000.0, ν_max=6400.0);
nothing #hide

Create a Voigt model for that

julia
band_voigt   = make_hitran_model(co2_par_band , Voigt(), architecture=CPU());
nothing #hide

Define the wavenumber grid and compute cross sections at the same pressure but different temperature:

julia
ν_band = 6300:0.01:6400;
σ_co2_Voigt220 = absorption_cross_section(band_voigt, ν_band, 1013.0 , 220.0);
σ_co2_Voigt290 = absorption_cross_section(band_voigt, ν_band, 1013.0 , 290.0);
nothing #hide
julia
fig = Figure(size=(700, 450))
ax = Axis(fig[1,1],
    xlabel = "Wavenumber (cm⁻¹)",
    ylabel = "σ (10⁻²⁰ cm²/molecule)")
lines!(ax, collect(ν_band), ff .* σ_co2_Voigt220, alpha=0.5, label="220K")
lines!(ax, collect(ν_band), ff .* σ_co2_Voigt290, alpha=0.5, label="290K")
xlims!(ax, 6300, 6380)
axislegend(ax, position=:rt)
fig

This is still a bit hard to see as the differences are relatively subtle. However, we can already see that absorptions near the band center are stronger at lower temperatures, while the opposite is true at higher rotational states. This behavior is caused by the distribution of lower rotational states, which only occupy higher Js at higher temperatures. This picture becomes more clear if we look at differences.

julia
fig = Figure(size=(700, 450))
ax = Axis(fig[1,1],
    xlabel = "Wavenumber (cm⁻¹)",
    ylabel = "Δσ (10⁻²⁰ cm²/molecule)")
lines!(ax, collect(ν_band), ff .* (σ_co2_Voigt220 .- σ_co2_Voigt290), label="220K − 290K")
xlims!(ax, 6300, 6380)
axislegend(ax, position=:rt)
fig

The interactive Plotly panel below is a qualitative illustration using five synthetic lines and manually chosen strengths/widths. It is not the HITRAN CO₂ calculation above and must not be used as a numerical reference:

Now we can clearly see how the change in the distribution of the lower rotational state leads to a redistribution of cross sections (a shift from lower to higher ground states with higher T)

Animated views ​

Set ENV["VSMARTMOM_RUN_HEAVY_DOCS"]="true" to generate the optional pressure and temperature animations with CairoMakie.

julia
if get(ENV, "VSMARTMOM_RUN_HEAVY_DOCS", "false") == "true"
T = 290.0
fig = Figure(size=(700, 450))
ax = Axis(fig[1,1],
    xlabel = "Wavenumber (cm⁻¹)",
    ylabel = "σ (10⁻²⁰ cm²/molecule)",
    yscale = log10)
record(fig, joinpath(@__DIR__, "absorption_pressure.gif"), 10:10:1100; framerate=15) do p
    empty!(ax)
    σ = absorption_cross_section(band_voigt, ν_band, Float64(p), T)
    lines!(ax, collect(ν_band), ff .* σ, label="p=$(p) hPa")
    ylims!(ax, 1e-7, 1e-1)
    axislegend(ax, position=:rt)
end

When run locally, this writes absorption_pressure.gif next to the tutorial.

julia
p = 900.0
fig = Figure(size=(700, 450))
ax = Axis(fig[1,1],
    xlabel = "Wavenumber (cm⁻¹)",
    ylabel = "σ (10⁻²⁰ cm²/molecule)",
    yscale = log10)
record(fig, joinpath(@__DIR__, "absorption_temperature.gif"), 10:10:320; framerate=15) do T
    empty!(ax)
    σ = absorption_cross_section(band_voigt, ν_band, p, Float64(T))
    lines!(ax, collect(ν_band), ff .* σ, label="T=$(T) K")
    ylims!(ax, 1e-7, 1e-1)
    axislegend(ax, position=:rt)
end

When run locally, this writes absorption_temperature.gif next to the tutorial.

julia
end

# More extreme case, let's take 10 atmospheres (10,000 hPa)
σ = absorption_cross_section(band_voigt, ν_band, 10000.0, 300.0);
fig = Figure(size=(700, 450))
ax = Axis(fig[1,1],
    xlabel = "Wavenumber (cm⁻¹)",
    ylabel = "σ (10⁻²⁰ cm²/molecule)")
lines!(ax, collect(ν_band), ff .* σ, label="p = 10000 hPa")
axislegend(ax, position=:rt)
fig

We can see that individual lines can be blurred at very high pressures once the broadening becomes wider than the line spacing.

This page was generated using Literate.jl.