Widespread persistent contrails over the western Great Lakes during 9 October 2000 were examined using commercial flight data, coincident meteorological data, and satellite remote sensing data from several platforms. The data were analyzed to determine the atmospheric conditions under which the contrails formed and to measure several physical properties of the contrails, including areal coverage, spreading rates, fall speeds, and optical properties. Most of the contrails were located between 10.6 and 11.8 km in atmospheric conditions consistent with a modified form of the Appleman contrail formation theory. However, the Rapid Update Cycle-2 analyses have a dry bias in the upper-tropospheric relative humidity with respect to ice (RHI), as indicated by persistent contrail generation during the outbreak where RHI ≥ 85%. The model analyses show that synoptic-scale vertical velocities affect the formation of persistent contrails. Areal coverage by linear contrails peaked at 30 000 km2, but the maximum contrail-generated cirrus coverage was over twice as large. Contrail spreading rates averaged around 2.7 km h−1, and the contrails were visible in the 4-km Geostationary Operational Environmental Satellite (GOES) imagery approximately 1 h after formation. Contrail fall speed estimates were between 0.00 and 0.045 m s−1 based on observed contrail advection rates. Optical depth measurements ranged from 0.1 to 0.6, with consistent differences between remote sensing methods. Contrail formation density was roughly correlated with air traffic density after the effects of competing cloud coverage, humidity, and vertical velocity were considered. Improved tropospheric humidity measurements are needed for realistic simulations of contrail and cirrus development.
Contrails are a source of anthropogenic cloudiness. Similar in physical properties to natural cirrus, contrails affect the atmospheric radiation budget and may influence climate. Because air traffic is expected to grow by 2%–5% annually for the next 50 years, contrail coverage will also increase and may produce a significant amount of radiative forcing by 2050 (Minnis et al. 1999). Unfortunately, current estimates of global contrail radiative forcing are extremely uncertain, and range over almost two orders of magnitude (Minnis et al. 1999; Meyer et al. 2002) as a result of large discrepancies in global estimates of contrail coverage, contrail altitude, contrail optical properties, and the co-occurrence of contrails over other clouds.
The relation of contrail properties to meteorological conditions is still unknown. Recent theoretical estimates of global contrail coverage (Sausen et al. 1998; Ponater et al. 2002) are tuned to early measurements of linear contrail coverage determined visually from infrared satellite imagery. The estimates differ based on the parameterization used to diagnose contrails and the numerical weather analyses employed to determine the ambient conditions. Estimates of contrail optical depth and particle size are typically based on aircraft measurements (Strauss et al. 1997; Schröder et al. 2000) or satellite retrievals (Minnis et al. 1998a,b; Duda et al. 1998), and have not been related to the surrounding atmospheric conditions.
Better understanding and more accurate prediction of contrail climatic effects will require better estimates of global contrail properties and their relation to the atmospheric state. Of particular interest are predictions of widespread outbreaks of persistent contrails in otherwise clear skies; these contrails are likely to have the largest impact on the global radiative energy budget. Determining the frequency of occurrence and the magnitude of such outbreaks would help refine estimates of global contrail radiative forcing. Additionally, development of reliable methods for diagnosing persistent contrails and their physical and radiative properties from numerical weather analyses is essential for predicting future contrail climate impacts.
Toward these goals, this paper studies widespread persistent contrails over the western Great Lakes during 9 October 2000 by using a combination of actual flight data, coincident meteorological data, and satellite remote sensing. This framework will enable the study of the physical properties of contrails in relation to the meteorological environment. To date, no previous study has used all three data sources to examine regional properties of persistent contrails. Walters et al. (2000) used a combination of radiosondes, surface observations of contrail formation, and flight data from Federal Aviation Administration (FAA) flight logs and air traffic control data to verify contrail forecast algorithms near Dayton, Ohio. Jackson et al. (2001) used a similar combination of data over eastern Massachusetts to develop a contrail forecasting model. Neither study discriminated between persistent and nonpersistent contrails, estimated contrail optical properties and areal coverage, nor evaluated time-dependent changes in the contrails. They only evaluated contrail occurrence at a single location. This research is intended to determine under which conditions a persistent contrail outbreak can occur, to measure the spreading of the contrail amount over time, to observe the development of contrail properties during the outbreak, and to determine empirical humidity thresholds required for contrail formation and persistence. The results of the analyses should be valuable for improving the parameterization of contrails in climate models.
Satellite imagery, flight track data, and meteorological data from numerical weather prediction model analyses were used to match individual contrails with flight tracks, and to determine the meteorological conditions under which the contrails formed. A brief description of each dataset is provided below.
a. Satellite data
A large outbreak of persistent contrails developed over the Upper Peninsula of Michigan and the surrounding region on 9 October 2000. Many of the contrails were visible in the eighth Geostationary Operational Environmental Satellite (GOES-8) imager 4-km channel-4 (10.7 μm) infrared brightness temperature (T10.7) images and in the channel-4 minus channel-5 (12.0 μm) split-window (T12.0) brightness temperature difference (BTD) imagery for at least 3–4 h before dissipating. Data from the National Oceanic and Atmospheric Administration's (NOAA's) Advanced Very High Resolution Radiometer (AVHRR) 1-km imager on NOAA-14 taken at 2050 UTC were also used (see Fig. 1), as were multispectral 1-km data from the Moderate Resolution Imaging Spectroradiometer (MODIS) on the Terra satellite.
Contrails were identified manually as bright linear features in half-hourly GOES-8 BTD images from 1015 to 2245 UTC 9 October 2000 (data from 1715 UTC were missing). The contrails were divided into linear segments, and their coordinates were logged into a database.
b. Air traffic data
Commercial air traffic data from the FlyteTrax product (FT; FlyteComm, Inc., San Jose, California) were used to determine flight tracks over the Great Lakes area on 9 October 2000. The database consists of 2- or 5-min readings of aircraft position (latitude, longitude, altitude) for every nonmilitary flight (identified by flight number) over the United States and Canada. A series of these readings for a particular aircraft flight is referred to here as a flight track, and each reading is called a flight track point. Although the FT database does not include military flights, it contains most of the air traffic over the Great Lakes region. All flight track points above 7.6 km (25 000 ft), between 36.5° and 51°N and between 95° and 79°W, and for the times between 1000 and 2200 UTC, were included in the analysis.
c. Meteorological data
Atmospheric profiles of height, temperature, humidity, and horizontal and vertical wind speeds were derived from the 40-km resolution, 1-h Rapid Update Cycle-2 (RUC) analyses (Benjamin et al. 2004a,b) in 25- hPa intervals from 450 hPa (approximately 6.4 km over the Great Lakes) to 175 hPa (12.7 km). The RUC data were linearly interpolated at each pressure level to a 1° × 1° grid at 15-min intervals from 1000 to 2300 UTC 9 October 2000. To minimize the processing time, only a portion of the RUC horizontal domain centered over the Great Lakes (from 36° to 54°N and from 97° to 75°W) was used.
3. Flight track analysis
The satellite, flight track, and meteorological data were used together to determine the altitude of contrail formation and the temperature and humidity conditions under which the contrails formed. As a first step in this process, the flight track data were advected to the times of the satellite observations by using the RUC model winds. The advected flight tracks were then compared with the satellite observations of contrails to determine which aircraft produced contrails. From the RUC meteorological data and contrail position data, the atmospheric conditions in the contrail formation regions could be determined.
a. Advection of flight tracks
The RUC horizontal winds were used to advect the aircraft flight tracks for each 15-min interval throughout the period of available data (1000–2200 UTC, hereafter called the analysis period). To simplify the advection scheme, the flight track altitudes indicated in the FT database were converted to pressure levels by using a domain-mean pressure–height profile that was computed hourly from the RUC data. The model winds were interpolated to these levels and used to advect the flight tracks for each time interval. After adjusting the model winds to true north–south (u) and east–west (υ) wind speed vectors, the wind speed and bearing of the winds were computed for each 15-min time step, and great circle navigation formulas were used to compute the new latitude and longitude at each time step.
The flight tracks were also advected vertically using the RUC-computed vertical velocities. The RUC vertical velocities were converted from units of hectopascals per second to meters per second using the hydrostatic approximation. The fall speed of the contrail particles was usually set to zero, although some simulations using prescribed contrail fall speeds were run to study the effects of particle sedimentation on contrail advection. Values ranging from −0.03 to 0.15 m s−1 were tested, and the results are discussed in section 4e.
b. Matching flight tracks with contrails
The flight tracks were matched using a contrail-centered coordinate system assuming a spherical earth (see Fig. 2). The great circle through each linear contrail segment serves as the &ldquo=uator” of the coordinate system. The “prime meridian” of the each segment is the great circle placed perpendicular to the first great circle and centered on the western end of the contrail. The length of the contrail segment therefore represents the longitudinal range of the cloud.
With this coordinate system, the proximity of each flight track to a contrail is evaluated. For each time step, the shortest distance (Ds) of each flight track point (including current flights and flights advected from an earlier time) from the equator was computed using a series of geometric vector operations. If the flight track did not extend across the entire longitudinal range of the contrail, the flight track was not considered. For all other flight tracks, the mean Ds was computed for each flight track segment across the longitudinal range of the contrail. The flight track with the minimum average Ds for all applicable time steps throughout the simulation was chosen as the flight track matching the contrail.
Table 1 shows the flights matched to 11 of the contrails during the analysis period. These contrails were arbitrarily chosen, but are expected to be representative of the entire contrail population. To reduce the possibility of mismatches, most of the contrails selected for Table 1 were relatively isolated from nearby contrails. For most of the flights, the selected flight track was obvious as the mean Ds of the selected track was usually less than 12 km and much less than the values from other flights. We note that the selected flight track usually was the best match for all times and for a variety of fall speeds, even in cases like ACA163 where other flight tracks had comparable (but larger) values of Ds. For two of the flights, ACA113 and AAL55, our confidence of a match is less due to the larger value of Ds during the analysis period. For all of the selected flights, the altitude of the aircraft ranges between 10.6 km (35 000 ft) and 11.8 km (39 000 ft, roughly between 240 and 200 hPa), and usually 1 h is necessary for the contrail to become wide enough to become visible in the GOES-8 4-km imagery. The typical width for a contrail was 8 km approximately 2 h after it was first seen in the GOES-8 imagery.
4. Retrieved contrail properties
Several physical properties of the contrails were determined from the available datasets. The areal coverage of the contrails was determined both subjectively by manual inspection of the BTD imagery and objectively by the automated method of Mannstein et al. (1999). The RUC wind analyses were used to estimate both the rate of contrail spreading seen in the satellite imagery and the fall speeds of the contrails. Multiple satellite remote sensing retrievals of contrail height, optical thickness, particle size, and other properties were made for selected groups of contrails. The numerical weather model analyses were also used to determine the conditions necessary for contrail formation. Finally, the FT database was used to relate air traffic density to the observed contrail formation density.
a. Areal coverage of contrails
The contrails formed nearly continuously over the western Great Lakes area throughout the analysis period. Once formed, the contrails moved to the southwest into Wisconsin, Indiana, and Illinois where they massed into nearly continuous decks of cirrus clouds before dissipating by 0200 UTC 10 October 2000. Automated, objective estimates of linear contrail coverage were made from the GOES-8 data in an area covering Wisconsin, southern Ontario, the Upper Peninsula of Michigan, and Lakes Michigan and Superior using the image processing methodology of Mannstein et al. (1999). The areas of analysis were determined by manual inspection of the GOES-8 imagery and circumscribed by analysis boxes. As the contrails advected, the analysis boxes were moved (and expanded) to track not only new linear contrails, but also the older contrails until they dissipated. An example of the analysis regions is shown in Fig. 3. A subjective analysis of the same regions was also performed by manually inspecting the T10.7 minus T12.0 images. The contrail coverage in each analysis box was based on the fraction of pixels in which both the magnitude of the mean BTD and the standard deviation in BTD of the 3 × 3 pixel grid centered on the pixel were larger than arbitrary thresholds. Although the manual inspection was subjective, it allowed for the tracking of contrails that become too diffuse to be identified in the automated method.
The Mannstein et al. (1999) algorithm was designed to work with 1-km AVHRR data, while the imagery from GOES-8 has a nominal resolution of 4 km. The viewing zenith angle in this case is around 53° resulting in a pixel size close to 7 km over the area of interest. In addition, the GOES-8 T10.7 and T12.0 images sometimes have across-track line errors that can produce noticeable striping. The objective algorithm can occasionally misclassify this artifact as a contrail (Meyer et al. 2002). To alleviate this problem, a weighted five-pixel average was used to smooth the differences between scan lines. For each pixel in each scan line, half of the weight is assigned to the line being smoothed, while a total of 40% is divided between the corresponding pixels in the two adjacent lines. The remaining 10% is divided between the corresponding pixels two scan lines away.
Figure 4 shows the time series of contrail coverage determined from both the objective technique and the manual inspection. The plot also shows the areal coverage determined by the objective algorithm at 2050 UTC from the 1-km AVHRR data. The close correspondence between the GOES- and AVHRR-derived coverage suggests that, at least in this example, the Mannstein et al. (1999) algorithm may be applicable to 4-km GOES-8 data. The areal coverages determined from the objective and subjective analyses match each other well between 1015 and 1645 UTC, but the coverage from the subjective method quickly becomes more than twice that retrieved from the objective algorithm after 1745 UTC, when many of the contrails in the analysis region became too diffuse and too closely clustered together for the algorithm to resolve. At 1745 UTC, the spreading linear contrails covered an area of nearly 25 000 km2 over the western Great Lakes region. By 2045 UTC, the linear contrails that could still be identified as such covered nearly 30 000 km2. This coverage does not include the cirrus clouds in northern Illinois and Wisconsin that were formed from older contrails, which more than double the contrail coverage. After 2045 UTC, both the linear contrail coverage and the coverage by contrail-generated cirrus slowly decrease at the end of the analysis period.
b. Contrail spreading
The increases in contrail areal coverage presented in the previous section are not only caused by increases in the number of contrails, but also by the spreading of individual contrails. Jensen et al. (1998) note that the vertical shear of the horizontal wind is the dominant process driving contrail spreading. From modeling studies they show that both particle sedimentation and the vertical motions within the contrail resulting from longwave radiative heating allow wind shear to spread the contrail. Several factors affect the magnitude and depth of the vertical motions in the contrail, including atmospheric stability, the magnitude of the radiative heating, the contrail ice water path, and the cloud microphysics. Jensen et al. (1998) simulate only one case of contrail growth; the relation between the contrail spreading rates and these factors is not completely known. The conditions for contrail spreading in this analysis were less favorable than those in the Jensen et al. (1998) study since the magnitude of the wind shear between 225 and 250 hPa (where most of the aircraft flights were located) usually ranged between 2 and 4 m s−1 km−1, compared to 6 m s−1 km−1 in Jensen et al. (1998). The actual shear affecting the contrails was even less since the direction of the vertical shear was usually not perpendicular to the flight path/contrail direction. Based on the predominant contrail direction (NW to SE), the mean vertical shear of the horizontal wind over the Upper Peninsula of Michigan and eastern Lake Superior (hereafter called upMI/eLS) during the analysis period was only 1.5 m s−1 km−1, while it was 2.8 m s−1 km−1 over eastern Wisconsin (hereafter called eWI). In addition, the atmosphere was slightly more stable in this study. The lapse rate between 225 and 250 hPa was between 7.6 and 8.5 K km−1 throughout the analysis period over both upMI/eLS and eWI.
Estimates of contrail spreading from GOES-8 observations are difficult to perform accurately due to the coarse resolution of the imagery. Observation of groups of trails in the 4-km GOES-8 imagery suggests that for flight times between 1400 and 1900 UTC, the trails are about 6 km wide 2.25 h after flight time, and about 10 km wide 3.75 h after flight time. An analysis of the RUC wind field over the western Great Lakes shows that the trails formed in a region of light shear (1–3 m s−1 km−1) that extends to a depth of 1 or 2 km below contrail depth, depending on the time of formation. Figure 5 shows the spreading of a contrail based on the RUC wind analyses. The spreading was computed from the wind field by comparing the advection of the contrails assuming a zero sedimentation (vertical deepening rate) and the advection of the contrail that is falling (deepening) by a specified constant rate. Figure 5 shows that once the contrail particles reach the level of larger shear a couple kilometers below the cloud, they are quickly pulled apart. Thus, only contrails with small, slowly falling particles would remain intact. Figure 5 also shows the contrail spreading rates estimated from the satellite observations. The observed contrail spreading rates are most consistent with contrails deepening at a rate around 0.060–0.125 m s−1.
c. Remote sensing analyses
Several remote sensing analyses were performed on the evolving contrails using GOES, MODIS, and AVHRR imagery. From the AVHRR data, an optical depth retrieval was performed using a fixed cloud temperature (FCT) approach similar to that of Meyer et al. (2002). First, the effective emissivity (ɛ) of each contrail pixel (identified by the objective algorithm) was computed from the observed T10.7, the observed background temperature Tb, and the contrail temperature Tc using
where B is the Planck function at 10.7 μm. For this study, Tc was assumed to be 215 K based on the regional mean temperature at 225 hPa. Here, Tb is calculated as the average temperature of all pixels at a distance of three pixels horizontally, vertically, or diagonally from a contrail pixel that are not adjacent to any other contrail pixel. In addition, Tb must be greater than Tc to ensure that the background pixels measure points below the contrail. If no pixels meeting these criteria are found, the mean background temperature calculated for all other contrail pixels within the local 10′ grid box is used. The visible optical depth (τ) of the contrail pixels is then computed from the effective emissivity by inverting the parameterization of Minnis et al. (1993b),
which accounts for the infrared scattering. In (2), μ is the cosine of the viewing zenith angle, and the coefficients a = −0.458 and b = 1.033 are for an axisymmetrical 20-μm hexagonal ice column. Differences in contrail particle size are expected to have little impact on the computed optical depth; calculations based on the much larger cirrus uncinus (CU) distribution (De = 123 μm) of Minnis et al. (1993a) show only a 4% relative change in τ. The average visible optical depth of the linear contrails in the 2050 UTC AVHRR image is 0.14 based on the FCT approach.
The longwave radiative forcing of the contrails was computed directly from the observed radiances. The unit contrail longwave forcing is given as
where Qc and Qb are the longwave fluxes for the contrails and background, respectively. The broadband fluxes were determined from the narrowband radiances observed by the satellite using the method of Minnis et al. (1991). The mean longwave forcing of the linear contrails in the AVHRR image (Fig. 1) was 10.4 W m−2.
The retrieval was also performed on two series of GOES-8 images. The first (loop A) followed several contrails over the Upper Peninsula of Michigan and the surrounding area between 1745 and 1915 UTC, while the second (loop B) tracked developing contrails over the eastern end of Lake Superior that drifted over the Upper Peninsula of Michigan and Lake Michigan from 1845 UTC until 2045 UTC. In both cases, most of the contrails in the analysis boxes had appeared in the GOES-8 imagery 1 h before the start of the time series. The results of the retrievals are shown in Fig. 6a. The mean retrieved optical depths in loop A increase from 0.10 at 1745 UTC to 0.17 at 1845 UTC, while the contrails in loop B start at 0.16 and increase to 0.31 one hour later before decreasing to 0.20 at 2045 UTC. In both time series the contrails reach their maximum in optical depth about 2 h after appearing in the satellite imagery, or 3 h after their initial formation. In this analysis, the contrail longwave radiative forcing is nearly exactly proportional to the cloud optical depth (Fig. 6b), and varies from 7.8 to 13.1 W m−2 in loop A and from 12.2 to 21.5 W m−2 in loop B.
The visible-infrared solar-infrared split-window technique (VISST; Minnis et al. 1995) was also used to infer the physical properties of the evolving contrails from the GOES-8 and MODIS data. VISST simultaneously determines effective particle size, visible optical depth (τυ), cloud temperature (Tc), and effective cloud height (zc) from the radiances observed from several GOES channels. Minnis et al. (1995) describe the algorithm in detail. Using clear-sky radiance estimates for each channel and the sun–satellite geometry, VISST computes the spectral radiances expected for both liquid droplet and ice water clouds for a range of particle sizes and optical depths. The effective diameters (De) for the hexagonal ice column model clouds vary from 6 to 135 μm. The effective cloud height is the altitude or pressure from the nearest vertical temperature profile that corresponds to Tc. The model cloud radiances are computed using the cloud emittance and reflectance parameterizations of Minnis et al. (1998a) and a visible channel surface– atmosphere–cloud reflectance parameterization (Arduini et al. 2002). VISST determines the cloud properties by matching the observed visible (0.65 μm), solar infrared (3.7 μm), and infrared (10.7 μm) radiances to the same quantities computed with the parameterizations and corrected to the top of the atmosphere. The process is iterative and can compute results for both ice and liquid water clouds. In the retrievals presented here, the model assumed all contrail clouds were ice.
VISST retrievals of ice-cloud properties using GOES and other multispectral satellite data have been compared to both in situ aircraft measurements (Young et al. 1998; Duda et al. 2002) and passive and active radiometric measurements at surface sites (Mace et al. 1998; Dong et al. 2002). The VISST retrievals show good agreement with the measurements for cloud height and cloud particle size. The derived heights for optically thin cirrus clouds are generally located somewhere between the bottom and center of the cloud (X. Dong 2003, personal communication; Minnis et al. 2002b).
Figure 6a includes the optical depths retrieved by VISST for loops A and B. Also shown in Fig. 6a is a VISST analysis using Terra MODIS data from the 1801 UTC overpass. The VISST optical depths are larger than the corresponding retrievals using the FCT method, but the VISST optical depths from GOES and from MODIS are in approximate agreement.
Another section of the contrail outbreak was analyzed from 1645 to 2145 UTC with the VISST retrieval. During this period, several contrails starting over the Upper Peninsula of Michigan at 1645 UTC spread into a diffuse cloud mass over Wisconsin that could not be identified as linear contrails by the automated algorithm. The white box over Wisconsin in Fig. 3 shows the location of the contrail-generated cirrus at 1945 UTC. Two soundings were used to initialize the retrievals throughout the period (see Fig. 7a). The 0000 UTC sounding from Gaylord, Michigan (APX), was taken in a region of midlevel cloudiness. The cloudiness appears in the APX sounding as a moist layer (RH as large as 80%) from 2.6 to 4.8 km. The Green Bay, Wisconsin (GRB), sounding at 0000 UTC shows two layers from 4.2 to 5.1 km and from 5.8 to 7.0 km with humidities near 50%, but no midlevel cloudiness is apparent at GRB between 1645 and 2145 UTC.
Figure 8a shows a time series of VISST-retrieved ice crystal sizes for the analyzed scene. Except for the large particle diameters (60 μm) measured at 1645 UTC, the mean particle diameters remain around 27 μm throughout the period. These particle sizes are somewhat smaller than those reported by Duda et al. (2001) for contrails near Hawaii. The smaller sizes are likely the result of competition for atmospheric moisture between the numerous contrails. The large particle sizes retrieved at 1645 UTC may be caused by differences between the backscattering properties of the actual contrail particles compared to the properties used in VISST. The scattering angle for the observations at 1645 UTC was about 179°. Another source of uncertainty in De is the occurrence of some low-level cloudiness in the box at that time.
The cloud heights retrieved from VISST are shown in Fig. 8b. The maximum heights from the VISST retrievals are 2–3 km higher than the mean values, but show similar trends during the day. Contrail heights were derived independently from stereography using combinations of GOES-8 with MODIS, AVHRR, and simultaneous GOES-10 (at 135°W) imagery following the method of Fujita (1982). The stereographic results within the VISST analysis region are also shown in Fig. 8b. Since the heights had to be determined from easily recognizable, bright features in the imagery, the stereographic estimates correspond to the highest portions of the contrail cirrus and, thus, are expected to be higher than the VISST results. They are generally within 1 km of the maximum VISST heights from GOES-8. The stereographic height estimate was verified in one contrail by using a concurrent estimate from cloud shadows; both methods derived a height of 10.7 km. The heights of freshly formed contrails (1–1.5 h after formation) ranged from 10.7 to 10.85 km.
d. Contrail formation conditions
From the classical Appleman (1953) criteria for persistent contrail formation, persistent contrails form when the atmospheric temperature reaches a critical value, and the relative humidity with respect to ice (RHI) is 100% or greater. Widespread cirrus formation is expected to occur above the humidity level required for homogeneous nucleation. Such cirrus would obscure and deplete much of the moisture needed for contrail formation. From the theoretical model for homogeneous nucleation developed by Sassen and Dodd (1989), the nucleation humidity threshold at the mean flight level of the matched flight tracks would be 154%. Thus, we would expect the contrail outbreak area to occur within the RHI range from 100% to 154%. Using the method of Schrader (1997) and an appropriate engine efficiency factor for commercial aircraft (Busen and Schumann 1995), the regions of persistent contrail formation at 225 hPa were computed using the RUC model data. The engine efficiency factor (η) is related to the ratio of heat to water vapor emitted by the aircraft, and Sausen et al. (1998) report that η = 0.3 ± 0.05 for the present commercial fleet. Although the engine efficiency factor affects the critical temperature at which contrails can form, the RHI must still be 100% or larger for contrail persistence. Assuming ice-supersaturated conditions, a change in η from 0.2 to 0.4 would result in a lowering of the mean persistent contrail formation level from approximately 295 hPa (9.36 km, 30 700 ft) to 315 hPa (8.92 km, 29 300 ft). Even when accounting for variations in the engine efficiency factor, the contrail coverage on 9 October 2000 exceeded the area predicted by classical theory, indicating that contrails formed in regions with RUC-analyzed RHI < 100.
The RUC model, like all numerical prediction models, computes a grid-averaged value of relative humidity that can be less than the RHI at individual locations within the grid. Thus, persistent contrails can exist even when the true grid RHI is less than 100%. In addition, upper-tropospheric humidity measurements (from radiosondes and other platforms) used in the RUC model analyses often underestimate RHI (Elliot and Gaffen 1991; Miloshevich et al. 2001), and a dry bias in the model is likely. Miloshevich et al. (2001) computed humidity corrections for radiosondes at low temperatures based on comparisons with frost-point hygrometer data, and found that radiosondes underestimate upper- tropospheric humidity as much as 60% at typical contrail formation temperatures (−50°C).
To determine an empirical humidity threshold where persistent contrails are possible, the locations of the persistent contrails identified from the GOES-8 imagery were compared to the RUC analyses of RHI. Figure 9 shows the paths of all flight tracks at 1745 UTC that occur in the area where the RHI was estimated to be 85% or greater. Nearly all of the contrails fall within this threshold. The largest RHI found in the area where the contrails appear in Fig. 9 was 105%. Figure 1 shows that cirrus and other clouds covered the remaining portion of the high RHI area south and east of the contrails, where RHI reached as high as 110%.
Since the temperature at the mean pressure level of the matched flight tracks was 215 K, the corrected RHI threshold for the contrails in this case would be 126%. The largest corrected RHI in the region containing contrails would be 146%. This corrected RHI range falls within the limits imposed by classical contrail formation conditions and the Sassen–Dodd criteria for homogeneous nucleation of cirrus. Figure 7b shows the RHI profile for the APX and GRB soundings adjusted by Eqs. (4) and (5) for temperatures below 243 K. The modified soundings have humidities similar to the RUC analyses. The APX sounding has adjusted RHIs over 100% between 11.0 and 12.7 km, while the GRB sounding has ice supersaturation between 11.5 and 12.6 km. The height of the supersaturated layer in the soundings is slightly higher than the corresponding RUC analyses, and is likely due to the effect of cold temperatures on the radiosonde sensor time constant (Miloshevitch et al. 2001).
e. Contrail fall speeds
Figure 10 shows the matched flight tracks and contrail segments identified from the GOES-8 imagery at 1745 UTC. Although several of the contrails are well matched to individual flights, the origins of some trails are indeterminate. Some may be due to the presence of military flights, or result from limitations in the flight- matching scheme. Good matches between the flight tracks and all of the contrails are also unlikely due to uncertainties and errors arising from several sources. These include the initial flight track positions, the estimation of contrail location from the 4-km GOES-8 imagery, and the advection of the flight tracks using the RUC model winds. The advection of the flight tracks is probably subject to the most uncertainty since the simulated tracks remained at nearly the same pressure level throughout the analysis period (the RUC vertical velocities were small), but the contrails are actually falling continually and are influenced by winds at different levels. The mean altitude (pressure level) of the flight tracks that most closely matched the contrails was 10.9 km (36 000 ft, 225 hPa), a value that corresponds closely to the results of the stereographic analysis in Fig. 8b.
Reports of contrail fall speed measurements are rare. Konrad and Howard (1974) measured the fall speed of contrail-generated streamers using radar. They measured fall speeds between 0.4 and 1.4 m s−1 for the precipitation streamers, but no measurements of contrails were mentioned. Schumann (1994) and Sussmann (1999) report lidar measurements of contrail height in the wake of commercial aircraft at cruise altitude. These measurements occurred during the first minute after contrail formation when a downward motion in the vortex structure was induced by the lift of the aircraft. They report mean contrail fall speeds of from 2.4 to 3.1 m s−1. After the wake vortices break up, Brunt–Väisälä oscillations return the exhaust back near its original altitude. Once the effects of the wake motions dissipate within the first few tens of minutes after contrail formation, the dispersion of the contrail is dominated by its interaction with the ambient atmosphere (Lewellen and Lewellen 2001), and the effects of particle sedimentation become important. The measurement of contrail fall speed beyond that at its initial formation has to our knowledge not yet been treated in the literature.
The effects of particle sedimentation on the advection of the contrails were investigated by running the flight advection (matching scheme) several times using different contrail fall speed values. The selected speeds ranged from −0.03 m s−1 (i.e., slowly rising) to 0.15 m s−1. Estimates of ice crystal fall speeds for maximum dimensions matching the retrieved particle sizes (from 22 to 60 μm) range from about 0.01 to 0.06 m s−1 (Pruppacher and Klett 1980; Heymsfield and Iaquinta 1999). A mean contrail fall speed is assumed here because the expected range of fall speeds for any particular contrail is small (since the variation in particle size is also small).
As mentioned earlier, the selected flight tracks for the 11 flights shown in Table 1 usually were the best matches for all times, even when different fall speeds were used. For some flights, the mean Ds dropped noticeably when a fall speed was included. Figure 11 shows the mean Ds computed for the flight tracks as a function of the contrail fall speed. The lowest values usually occur when the fall speed is from 0.0 to 0.045 m s−1. For most flights, the Ds values become larger when the fall speed exceeds 0.06 m s−1. Figure 11 also shows the mean Ds computed for all 11 flights as a function of contrail fall speed (squares). The smallest errors occur when the fall speeds are between 0.015 and 0.03 m s−1. This result suggests that the fall speeds of most of the contrail particles are consistent with the particle sizes retrieved from the VISST technique.
The RHI fields from the RUC analyses confirm the conclusion that the contrail fall speeds must be small. The mean RHI over upMI/eLS usually remains above the 85% level throughout the analysis period at both 200 and 225 hPa, while at 250 and 275 hPa the mean RHI exceeds 85% until approximately 1800 UTC. At lower levels the mean RHI over upMI/eLS is usually less than 85% throughout the analysis period. Over eWI, the mean RHI is above the 85% threshold at both 200 and 225 hPa throughout most of the analysis period, near the threshold at 250 hPa, and below the threshold at lower levels. Thus, the RUC humidity analyses suggest that the contrail particles would remain in ice-saturated air as long as they did not fall below approximately 260 hPa (10.1 km, or 33 100 ft). Since most flight tracks are around 11 km, and many contrails appear to persist for at least 4 h, the maximum net contrail fall speeds based on the RUC RHI analyses would be 0.06 m s−1. Although the results from Fig. 5 suggest a contrail-deepening rate as large as 0.125 m s−1, the net contrail fall speeds may be less since vertical motions in the contrail would loft many of the cloud particles to higher altitudes. In the Jensen et al. (1998) study, the simulated contrail top appears to rise at a rate around 0.10 m s−1, which would be consistent with both the contrail fall speeds presented in Fig. 11 and the contrail deepening rates in Fig. 5.
The fall speed of an individual contrail was measured by using stereography and flight track data. Flight CDN902 was one of the flights identified from the GOES-8 imagery in Table 1. Estimates of the contrail height from both stereography and cloud shadows were 10.7 km at 1945 UTC. The flight track data show that the flight passed over the region about 60 min earlier. Since the flight track altitude of CDN902 was 37 000 ft (11.28 km), the estimated contrail fall speed is 0.16 m s−1. This fall rate is somewhat larger than the rate derived from contrail advection and spreading.
f. Contrail density versus air traffic density
To correlate contrail density with air traffic data, the GOES-8 imagery was manually inspected over a 12-h period (1000–2200 UTC) to determine when and where contrails initially formed. Since most contrails appeared in the images approximately 1 h after their probable flight track of origin, the contrails were advected backward in time by the mean wind at 225 hPa to estimate their initial location. A time step of 30–50 min back in time produced the best overall agreement between air traffic and contrail density. The contrail densities within 1° × 1° grid boxes from 42° to 50°N and 92° to 80°W were computed for each hour, and averaged over the 12-h period. Similarly, the air traffic density over the same period and domain was computed. To improve the match between the air traffic and contrail density data, several contrail formation and visibility factors were considered to determine those regions where contrails would not form nor be seen. First, the GOES-8 imagery was inspected to find where and when noncontrail overcast existed that would obscure contrail detection. Second, those areas with 225-hPa RHI that were below the 85% threshold were selected to account for the contrail formation conditions (the temperatures throughout the region were cold enough to satisfy the modified Appleman criteria). An inspection of the satellite imagery and the RUC vertical velocity field at 225 hPa (based on 1-h forecasts) showed that no contrails formed in regions where the RUC vertical velocity was less than −0.015 m s−1, even when the RHI was above the 85% threshold. Since the RUC uses isentropic coordinates at 225 hPa, the vertical velocity analyses are expected to be able to distinguish between rising and falling motion at this level (Benjamin et al. 2004b; S. G. Benjamin 2003, personal communication). Thus, all areas where the vertical velocity was less than −0.015 m s−1 (sufficiently strong negative motion) were selected, and all regions meeting these cloudiness, humidity, or vertical velocity criteria were excluded from the air traffic and contrail density calculations. Figure 12 shows a plot of the mean contrail density versus air traffic density for the areas of the grid not affected by the formation and visibility factors. The contrails were advected backward in time by 30 min in the figure. The contrail density appears to be roughly correlated with air traffic density, even though the contrails from only 1 in 10 flights became widespread enough to be detected in the 4-km GOES-8 imagery.
Several physical properties were retrieved from an outbreak of contrails over the western Great Lakes on 9 October 2000 by using a combination of satellite, meteorological, and air traffic datasets. Contrail-generated cirrus more than doubled the areal coverage compared to linear contrail coverage estimates derived using the objective analysis. The rate of contrail spreading observed from the GOES-8 imagery was 2.7 km h−1, which allowed the contrails to be visible in the GOES imagery approximately 1 h after formation. This rate was about half the rate reported in the simulation by Jensen et al. (1998), but their study has larger wind shear and slightly less stable conditions. The estimates of spreading based on the RUC wind analyses suggest that the contrails deepened at a rate around 0.060–0.125 m s−1, while the net fall speed of the contrails were less than 0.045 m s−1. The difference between these rates implies that at least part of the contrail deepening was the result of lofting by vertical motions within the contrail. In contrast, the fall speed derived for contrails from stereography and flight track data suggests that contrail fall speeds may have been slightly larger (0.12–0.16 m s−1).
The two remote sensing analyses gave divergent results of some cloud properties. The stereographic results and the minimum level of saturated air in the soundings suggest that the mean VISST cloud heights are far too low and optical depths are too large while τc from the FCT method is more accurate. The FCT results indicate that the mean contrail τc was between 0.1 and 0.2, and the maximum contrail optical depths occurred approximately 3 h after formation. Optical depths from VISST are typically about twice as large as those retrieved from the FCT method. For a nebulous cloud like contrails or cirrus, the VISST-retrieved value of Tc should correspond to some altitude between the cloud top and bottom. Thus, the retrieved temperatures should be at some level below the cloud top. For spreading contrails, the difference between cloud top and the effective radiating altitude should be less than 1 km. The thinness of the contrails makes retrievals of contrail optical depth more difficult for visible-channel-based retrievals. VISST optical depths corresponding to the maximum heights are probably the most accurate values because they are closest to and follow the trends in the stereographically determined cloud-top heights, but they do not represent the entire cloud field. Although the mean VISST-derived optical depths are too large, the relative trends in height and optical depth are reasonable. Furthermore, the longwave radiative forcing computed using the VISST-retrieved cloud properties should be comparable to that from the FCT approach (10.4 W m−2) because the greater cloud temperatures would offset the greater optical depths.
Figure 9 demonstrates that the air traffic density in the region where contrails formed was relatively uniform. The mean air traffic density per 1° × 1° grid cell for flights with altitudes greater than or equal to 30 000 ft was usually only four per hour in the region where contrails formed, which may account for the overall uniformity and clarity of the identified contrail segments. As suggested by the nonlinear shape of the air traffic–contrail density relationship in Fig. 12, observed contrail densities probably reach an asymptotic limit in the GOES imagery at higher air traffic densities. This limit would be due to both the smearing of multiple contrails in satellite imagery and the production of fewer contrails growing wide enough to be seen in the imagery because of competition for water vapor. The formation of the cirrus deck seen in Fig. 1 is an example of smearing of the individual contrails in areas with dense air traffic. The matching of contrail density data with air traffic data showed that approximately 10% of the flights produced contrails wide enough to be seen in GOES-8 imagery.
This study has provided the first detailed analysis of contrail formation and motion using actual flight data integrated into numerical weather analyses. A flight track advection and matching scheme was developed that allows for the determination of the atmospheric conditions associated with individual contrails seen in satellite imagery. The contrail outbreak analyzed here reveals a clear case of relative humidities that exceed ice saturation at flight level over an area where no natural cirrus were forming over the course of the day. These results also reinforce the need for adjustment of the RHI fields from rawinsondes and numerical weather prediction models, even in high-resolution models such as the RUC, which also have sophisticated cloud/moisture parameterizations. Grid-averaging effects and the difficulty in measuring upper-tropospheric moisture lead to a dry bias in model analyses. More studies are necessary to conclusively determine the magnitude of this adjustment. Applying an RHI adjustment, however, should bring theory and observations into closer agreement.
The analysis results also highlight the difficulties in retrieving microphysical properties of optically thin ice clouds using reflected solar radiation. Use of radiances at other wavelengths such as 1.38, 3.9, 8.5, and 12 μm may yield more accurate objective results. Contrail fall speed estimates derived from several sources show that to first order the contrail height can be assumed to be constant throughout its lifetime and, therefore, the contrail temperature can be specified to effect a simple infrared retrieval like that used by Meyer et al. (2002). The fall speeds reported here (0–0.045 m s−1) are the first published estimates of fall speeds in mature contrails, and are consistent with ice crystal fall rates for the particle sizes retrieved from VISST.
Contrail formation density was roughly correlated with air traffic density after the effects of competing cloud coverage, humidity, and vertical velocity were considered. The RUC vertical velocity analyses showed that vertical velocity has an important effect on contrail formation, as areas with negative vertical velocities over 0.015 m s−1 in magnitude produced no contrails even when the RHI was above the threshold where contrails persisted in other regions.
A surprising result from this study is that the linear contrail detection algorithm developed by Mannstein et al. (1999) for application to 1-km infrared satellite imagery produced nearly identical amounts of contrail coverage using 4-km GOES-8 and 1-km AVHRR data. The single case study, however, is insufficient to conclude that comparable results will always be obtained. Many more comparisons are needed to reach any firm conclusions. The results also confirmed earlier studies that the linear-contrail coverage approach underestimates the contrail-induced cirrus coverage. Although contrail spreading has been examined in several previous studies, there have been no direct comparisons between the automated retrievals and contrail cirrus coverage determined subjectively from satellite loops. The factor of 2 difference found here for the aged contrail cirrus versus the automated linear contrail coverage is smaller than estimated in the previous studies but is still a significant factor.
The analysis procedures developed here provide the framework for realistic simulations of contrail life cycles and interactions. Despite the difficulties in determining upper-tropospheric humidity, the locations of the contrails in this example confirm that the growth and longevity of contrails are explained by the general Appleman temperature (relative humidity) criteria. The frequency of atmospheric conditions favorable for the development of widespread persistent contrails could be estimated by examining numerical weather model analyses like the RUC if the fields are adjusted for RHI and the dry bias. Much additional research is required, however, to parameterize contrail spreading, precipitation, and dissipation. Although it is currently possible to simulate contrails and cirrus to some degree, it will be difficult to realistically simulate contrails and cirrus clouds using the current relative humidity datasets without adjustments because the mechanisms for naturally producing cirrus clouds require values of RHI that do not appear in current operational measuring systems. Without the proper humidities, the physics of the formation process must be altered to produce cirrus in detailed models. Improved measurements of moisture in the upper troposphere, therefore, remain a critical need for accurate modeling of upper-level clouds and contrails.
This research was supported by Grants NASA NAG-1-2135 and NAG-1-02044 and by the NASA Pathfinder program. The authors wish to thank Sunny Sun-Mack and Kirk Ayers for the VISST retrievals, Donald Garber for his help on the stereographic analyses, and Hermann Mannstein for providing a version of his analysis program.
Corresponding author address: David P. Duda, NASA Langley Research Center, MS 420, Hampton, VA 23681-2199. Email: email@example.com