Retrieval of CO2 Column Abundances from Near- Infrared Spectroscopic Measurements

A Candidacy Report by Vijay Natraj

Abstract

A retrieval algorithm was developed to obtain the column averaged dry air mixing ratio of CO2(XCO2) from spectroscopic measurements of absorption in three near-infrared (NIR) bands. The radiative transfer computations were done using the discrete ordinate method, and optimal estimation theory was used for the inversion. To test the algorithm, column O2 was retrieved from measurements of absorption in the O2 A-band over sea surface. After correcting for instrument calibration effects, excellent fit was obtained with the measured spectra. The column O2 was retrieved with an error of around 1%. With spectral averaging using multiple soundings, it is possible to reduce the error to 0.1%. This suggests the feasibility of retrievingXCO2 with precisions better than 0.3%.

Table of Contents

1Introduction1

1.1Problem Description4

1.2Fundamentals of Atmospheric Radiative Transfer5

1.3Spectroscopy and the Underlying Physics of the Measurements10

2Retrieval Strategy14

2.1Discrete Ordinate Method16

2.2Convolution19

2.3Inverse Method20

3Case Study: Retrieval of Aircraft Measurements23

4Future Work31

5Conclusion33

Acknowledgments34

References34

  1. Introduction

Atmospheric CO2 is an efficient greenhouse gas. The CO2 concentration has increased from 280 to 370 parts per million (ppm) since the beginning of the industrial era [IPCC, 1996; Schnell, 2001]. There is growing apprehension that this will adversely alter the global climate [Cicerone, 2001; IPCC, 1996].Measurements from a global network of surface stations[CDIAC, 2001; Schnell, 2001] indicate that the biosphere and oceans have absorbed almost half of the carbon emitted during the past 20 years.However, the nature and geographic distribution of these CO2 sinks are not adequately understood, precluding accurate predictions of their response to future climate or land use changes[Cicerone, 2001; IPCC, 1996].One major concern is that these sinks may saturate in the future, accelerating the build-up of atmospheric CO2[Cox, 2000; Friedlingstein, 2001]. Existing models and measurements also fail to explain why the atmospheric CO2 concentration increases vary from 1 to 6 gigatons of carbon (GtC) per year in response to steadily rising emission rates (see Fig. 1) [Conway, 1994; Frolking, 1996; Houghton, 2000; Keeling, 1995; Le Quere, 2000; Lee, 1998; Randerson, 1999; Randerson, 1997].

The principal shortcoming of the ground-based network is its sparse spatial sampling. Accurate, time-dependent, spatially-resolved, global maps of the column-averaged CO2 dry-air mole fraction (XCO2) will dramatically improve our understanding of its surface sources and sinks. Modelling studies [Rayner, 2001] confirm that source-sink inversion algorithms employing global, space-based measurements of XCO2will outperform those using the existing GLOBALVIEW-CO2(GV-CO2) data if the space-based measurements have accuracies better than 2.5 ppm. Figs. 2 and 3 illustrate the resulting improvements in retrieved carbon flux errors on global and continental (~2107 km2) scales, respectively, for XCO2 errors of 1ppm (0.3%).

Fig. 1: (A) 40-year history of atmospheric CO2 buildup. (B) Observed variations in annual atmospheric CO2 accumulation (CO2) compared with fossil fuel emissions[CDIAC, 2001; Schnell, 2001]. Significant changes in carbon sequestration occur on annual time scales.

Fig. 2: Error in retrieved global carbon flux (GtC/yr) vs. space-based XCO2 measurement accuracy [Rayner, 2001]. The OCO 1-ppm accuracy (circle) is needed to outperform the GV-CO2 network (dashed–line).

Fig. 3: Global maps of carbon flux errors for 26 continent/ocean-basin-sized zones retrieved from inversion studies. A) Studies using data from the 56 GV-CO2 stations produce flux residuals that exceed 1 GtC/yr in some zones. B) Inversion tests using global XCO2 pseudo-data with 1 ppm accuracy reduce the flux errors to <0.5 GtC/yr/zone [Rayner, 2001].

The Orbiting Carbon Observatory (OCO) mission was thus set up to make the first space-based measurements of atmospheric CO2 with the accuracy, precision, resolution, and coverage needed to characterize the geographic distribution of CO2 sources and sinks and quantify their variability. The OCO satellite will carry three high-resolution, grating spectrometers, designed to measure CO2 and O2 near-infrared (NIR) absorptions from reflected sunlight. Remote sensing retrieval algorithms will process these data to yield estimates of XCO2 with accuracies near 0.3%. Chemical transport models will use OCO XCO2 data and other measurements to retrieve the spatial distribution of CO2 sources and sinks on regional scales over two annual cycles.

1.1Problem Description

The aim of this work is to develop and test algorithms for the retrieval of XCO2 from spectrometric measurements in three NIR bands. Such a retrieval is a highly nontrivial problem and is the subject of extensive research. Any retrieval problem can be broadly divided into two main components, viz., a forward model and an inverse method.

The forward model is an approximate scheme to describe the radiative transfer in the atmosphere. It involves computation of the absorption coefficients, and using these along with key parameters such as the scattering and absorption optical depths, the single scattering albedo and the surface reflectance to compute the radiances at specific solar zenith and viewing angles. These terms will be defined in section 1.2.The numerical scheme used for the forward model will be described in section 2.1. The computed radiances have to pass through a convolution function to simulate the instrument response. This will be discussed in further detail in section 2.2.

The idea behind retrieval is to compare the measured spectrum with the computed spectrum, and iteratively improve the computed spectrum to best match the observed spectrum. The standard method used is the Optimal Estimation Theory, which was propounded by Rodgers [Rodgers, 2000] and will be considered at length in section 2.3.

1.2Fundamentals of Atmospheric Radiative Transfer[Goody, 1989; Liou, 2002]

The specific intensity of radiation, or radiance, denoted by I, is the main quantitative characteristic of a radiation field. It depends on the frequency of radiation υ, the coordinates r of the point under consideration, and the direction n of the ray. Physically, it refers to the flux of energy at point r in direction n per solid angle per frequency. The dependence of the radiance upon these values is usually denoted as Iυ(r,n), where the subscript is usually dropped for economy of notation.

In atmospheric transfer problems, it is convenient to use a dimensionless quantity called the optical depth, rather than the actual depth. The optical depth τ can be defined as

(1)

where σ is the cross-section of an individual molecule, n is the number density of molecules and l is the path length. The product of σ and n is called the extinction coefficientk. The optical depth measures the amount of extinction a beam of light experiences travelling between two points. For historical reasons, the optical depth is defined to be 0 at the top of the atmosphere and increases as we go towards the surface.

The fundamental equation of radiative transfer, which was first postulated by Schwarzschild, is as follows (when the atmosphere is approximated as being vertically stratified and horizontally homogeneous; commonly called plane-parallel):

(2)

where μ is the cosine of the zenith angle, and J is the source function. The first term on the right hand side represents attenuation due to absorption and scattering of a radiance stream as it propagates through the atmosphere, and the source function represents the strengthening of the radiance stream. For solar radiation it arises from photons scattered in the path from all other directions. The presence of this scattering source term ensures that the radiation field is no longer merely a function of local sources and sinks, but of the entire atmospheric radiation field and of its transport over large distances. In practice this makes the solution much more difficult to obtain. If there is no source function, then the above equation reduces to the familiar Beer’s law.

Solution of equation (2) has been the subject of intense research over the years, starting from Chandrasekhar’s classical solution [Chandrasekhar, 1950] to various numerical methods[van de Hulst, 1980]. The technique adopted in this work will be explained in section 2.1. However, a discussion of the fundamental quantities involved in any radiative transfer problem is in order.

The four main quantities influencing the radiation field are the optical depth, the surface reflectance, the single scattering albedo and the phase function. The surface reflectance, denoted by a, is the ratio of the intensity reflected from the surface to that incident on it. It is a function of wavelength. A further complication arises because it can also depend on the zenith and azimuthal angles. This dependence is called the bidirectional reflectance distribution function (BRDF), and will be discussed further in section 3.

The single-scattering albedo ω0 represents the fraction of energy scattered to that removed from the radiance stream under consideration. For conservative scattering, ω0 = 1, and for pure absorption, ω0 = 0. In an inhomogeneous atmosphere, the single scattering albedo is a function of optical depth.

When light strikes a particle with an index of refraction different from its environment, the light is refracted. The angle at which the light is bent is a function of the size and shape of the particle as well as the wavelength of the incident light and the incidence angle of the light. In general, each particle will have a different scattering profile. This scattering profile is called the phase function.For simplicity an average phase function which adequately describes the most important features of the scattering process is used. This average phase function is further constrained by assuming that the probability of scattering from one direction into another is a function only of the angle between the two directions. The phase function Pthen describes the amount of light scattered from the incident direction into the scattered direction.

Using spherical trigonometry cosine laws, the scattering angle Θcan be obtained in terms of the directions of the incident and scattered radiation:

(3)

where μ and μ’ are the cosines of the incident and scattered zenith angles, and and are the incident and scattered azimuthal angles respectively. The convention is adopted that μ > 0 refers to upward radiance streams and μ < 0 refers to downward radiance streams.

There are a number of ways in which the phase function may be normalised, but the most natural is that used by the astrophysicists. They treat the phase function as a probability distribution; consequently, their normalisation condition requires the integral of the phase function over all angles to equal unity.

(4)

The simplest phase function is the isotropic phase function

(5)

The factor of 1/4πresults from the normalisation condition and the fact that there are 4πsteradians in a sphere.

Particles that are very small in comparison to the wavelength of the radiation (such as oxygen and nitrogen molecules) exhibit a type of scattering called Rayleigh scattering. The Rayleigh phase function is given by the following relation:

(6)

If the phase function is not isotropic, then a parameter called the asymmetry parameter is used to describe the degree of anisotropy of the phase function. This parameter is often denoted by g and is defined as the integral over all angles of the phase function multiplied by the cosine of the angle. In other words, it is the first moment of the phase function.

A single-parameter analytic form for an anisotropic, scattering phase function is the Henyey-Greenstein phase function:

(7)

1.3Spectroscopy and the Underlying Physics of the Measurements

High-resolution (R=~20,000) NIR spectrometry of complete rotation-vibration bands was chosen to obtain XCO2[Kuang, 2002]. This technique has the advantage that, at NIR wavelengths, thermal emission from the surface, the atmosphere, and the instrument are all negligible compared to reflected sunlight, simplifying radiative transfer and instrument calibration. Also, absorption of sunlight by the NIR CO2 bands is most sensitive to the CO2 concentration in the boundary layer, where the effects of sources and sinks are most readily detected. Thermal infrared soundings can yield observations of CO2 from space, but these measurements have limited sensitivity to CO2 concentrations near the surface (see Fig. 4) and their accuracy (~1%) is insufficient to meet the 0.3% requirement described earlier. Solar occultations also lack the ability to sense the boundary layer. LIDAR has been advocated for monitoring CO2 from space, but even the best ground-based trace gas LIDARS have accuracies around 5%.

Fig. 4: A comparison of the averaging kernels for column CO2 soundings using NIR absorption of reflected sunlight and thermal IR emission. NIR measurements are much more sensitive to surface phenomena.

The CO2 column abundance will be retrieved from the absorption band centred at 1.61m (see Fig. 5). This spectral region (called the weak CO2 band in this report) is virtually free of absorption by other gases, simplifying CO2 retrievals. Also,most of the spectral lines are not saturated, such that their absorption increases almost linearly with the CO2 abundance and path length.

However, the absorption also depends on factors unrelated to the CO2 abundance (e.g. clouds, aerosols, atmospheric temperature profile, surface pressure). Measurements of a reference gas, whose concentration is uniform, constant, and well known, are needed to convert the observed CO2 column abundance to the dry air mole fraction, XCO2. Molecular oxygen is the best available candidate. XCO2 can be derived from soundings of the total column abundances of CO2 and O2 as follows:

(8)

where 0.2095 is the O2 mole fraction.

Measurement biases produce systematic offsets in XCO2. These biases can be introduced by the instrument (e.g. uncertainties in the zero offset and instrument lineshape, stray light contamination, detector nonlinearity),surface and atmosphericproperties (uncertainties in cloud and aerosol opacity, temperature, surface pressure, surface reflectance, topography, etc.), and/or incompleteness and errors in the existing linelist database. While the instrument errors can be minimised by proper calibration, the other systematic errors need to be accounted for.

Fortunately, XCO2 can be measured to much higher precision and accuracy than column CO2 alone, because many of the systematic errors affect the CO2 and O2 columns similarly (e.g., cloud, aerosol, surface pressure), and therefore cancel in the ratio of these quantities.

Aircraft measurements [O'Brien, 1992]have shown that O2 A-Band measurements can provide column O2 (i.e., surface pressure) estimates with accuracies of 1 mbar (0.1%).Airborne particles can absorb or scatter sunlight back to space before it traverses the full atmospheric column, precluding true column CO2 measurements in regions occupied by opaque clouds. The O2 A-band is sensitive to clouds and small aerosol particles[Heidinger, 2000; Stephens, 2000]. The fact that the band contains both weak and strong lines(see Fig. 5) provides additional information on the vertical distribution of clouds and aerosols.O2 A-band measurements are therefore used to characterize the clouds and aerosols in each sounding, so that those with too much scattering can be rejected.For less opaque soundings (optical depths 0.030.3), cloud and aerosol scattering can introduce errors in XCO2 by adding uncertainty to the photon path length. However, O2A-band observations alone are not adequate for characterizing the scattering by water ice clouds and aerosols in the weak CO2 band because their optical properties (optical depth, single scattering albedo, and phase function) can vary substantially between the two bands.

Aerosol effects on the photon path length can be characterised by combining simultaneous spectroscopic observations from the O2 A-band and the CO2 band near 2.06 m (called the strong CO2 band in this report). Unlike the weak CO2 band, the CO2 absorption near 2.06 m is produced by strongly saturated lines that have a relatively weak (square root) dependence on the CO2 concentration and enhanced sensitivity to clouds and aerosols (see Fig. 5). Measurements of water vapour absorption lines and CO2 hot bands within the 2.06 m region also provide explicit observational constraints on these properties.

Fig. 5: Simulated atmospheric transmission (solar zenith angle=40) for the three OCO spectrometers (instrument effects included), including all absorbing species in each interval (O2-green, CO2-black, H2O-blue)

The spectral range for each band includes the complete band as well as some continuum at both ends. By using the entire band, biases due to uncertainties in atmospheric temperatures are minimised. The continuum at the band edges provides additional information about the wavelength dependent optical properties of the surface reflectance and airborne particles.

  1. Retrieval Strategy

The principal characteristics and flow of the XCO2 retrieval algorithm are presented schematically in Fig. 6. For each sounding, the retrieval process begins with an assumed environmental state, defined by the surface pressure ps, surface reflectancea, vertical temperature profileT(z), mixing ratios of CO2, water vapor, and other trace gases, [X(z)], and cloud and aerosol optical depth distributions (c and a, respectively). These parameters can be initialized from known climatology, or from adjacent retrievals. This information is combined with pre-tabulated, wavelength-dependent gas, aerosol, and cloud optical properties. The gas absorption cross-sections, (p,T), in the three spectral regions are derived and tabulated as functions of p and T using a line-by-line model [Meadows, 1996] and spectral line databases such as HITRAN. For clouds and aerosols, the wavelength-dependent optical properties (absorption and scattering cross sections, and phase functions) for liquid water, water ice crystals, and a wide variety of common aerosol types have been derived for a range of common size distributions using Mie scattering (liquid water, sulfates), T-Matrix (small dust), and geometric optics (cirrus cloud) codes. This tabulated data is combined with the atmospheric state parameters and information about the viewing geometry and solar zenith angle, and used in a multi-layer, spectrum resolving (line-by-line), multiple-scattering model to generate angle-dependent radiance spectra for the three bands[Crisp, 1997]. These synthetic spectra are then processed with a model that simulates the instrument’s spectral response to the incident radiation, and produces results that can be compared directly to the calibrated spectra.