## Abstract

Ice crystal terminal velocities govern the lifetime of radiatively complex, climatologically important, low-latitude tropopause cirrus clouds. To better understand cloud lifetimes, the terminal velocities of low-latitude tropopause cirrus cloud particles have been estimated using data from aircraft field campaigns. Data used in this study were collected during the Cirrus Regional Study of Tropical Anvils and Cirrus Layers–Florida Area Cirrus Experiment (CRYSTAL-FACE) and the Pre-Aura Validation Experiment (Pre-AVE). Particle properties were measured with the NCAR video ice particle sampler (VIPS) probe, thus providing information about particles in a poorly understood size range. Data used in this study were limited to high-altitude nonconvective thin clouds with temperatures between −56° and −86°C.

Realistic particle terminal velocity estimates require accurate values of particle projected area and mass. Exponential functions were used to predict the dimensional properties of ice particles smaller than 200 microns and were found to predict ice water content measurements well when compared to power-law representations. The shapes of the particle size distributions were found to be monomodal and were well represented by exponential or gamma functions. Incorporating these findings into terminal velocity calculations led to lower values of mass-weighted terminal velocities for particle populations than are currently predicted for low-latitude ice clouds. New parameterizations for individual particle properties as well as particle size distribution properties are presented and compared to commonly used parameterizations. Results from this study are appropriate for use in estimating the properties of low-latitude thin and subvisible cirrus at temperatures lower than −56°C.

## 1. Introduction

Low-latitude in situ generated cirrus and anvil cirrus near the tropopause have a significant effect on the radiation budget of the planet (McFarquhar et al. 2000). Accurate knowledge of the microphysical properties of these clouds is important to better understand their impact on climate. Recently, low-latitude tropopause cirrus clouds have been investigated during several field studies. The 2002 Cirrus Regional Study of Tropical Anvils and Cirrus Layers–Florida Area Cirrus Experiment (CRYSTAL-FACE) and the 2004 Pre-Aura Validation Experiment (Pre-AVE) projects both investigated upper troposphere cirrus with the NASA WB-57 aircraft. During both of these experiments, the National Center for Atmospheric Research (NCAR) video ice particle sampler (VIPS) probe (McFarquhar and Heymsfield 1997, hereafter MH97) was used to provide particle size distribution and particle projected area information for particles from 10 to 350 microns in maximum dimension. The microphysical characteristics of low-latitude tropopause cirrus particles as small as 10 microns have rarely been measured reliably owing to the limitations of electro–optical probes (Korolev et al. 1998).

Observations show that low-latitude tropopause cirrus clouds are common (Hartmann et al. 1992). Nee et al. (1998) observed thin or subvisible cirrus clouds over Taiwan in 50% of their lidar observations between 1993 and 1995. McFarquhar et al. (2000) reported that thin cirrus layers were present above 15 km in 29% of lidar observations during the Central Equatorial Pacific Experiment (CEPEX) and the Tropical Ocean Global Atmosphere Coupled Ocean–Atmosphere Response Experiment (TOGA COARE). McFarquhar et al. also showed, using a delta-four stream radiative transfer model, that subvisible cirrus layers have heating rates of up to 1.0 K day^{−1} with a radiative forcing of 1.2 W m^{−2}, equivalent to about 0.7°C change in surface temperature globally. They point out that, due to the extent and frequency of occurrence of the clouds, the radiative effect of subvisible cirrus in the tropics should not be ignored.

Climate models and cloud-resolving models (CRMs) would benefit from improved parameterizations of cirrus cloud properties (Fu 1996). CRMs and general circulation models (GCMs) remove ice particles from cloud layers either by parameterized terminal velocity values or by autoconversion schemes that convert “cloud ice” to “snow.” Different autoconversion rates lead to different rates of removal of cloud ice from layers, indirectly implying a terminal velocity for the cloud ice. Models that utilize an explicit terminal velocity parameterization generally represent the falling of an ensemble of cloud particles with a single terminal velocity value. Variance from this value adds complexity that may be computationally prohibitive.

The NCAR Community Atmosphere Model (CAM) uses an effective radius *r _{e}*, predicted from temperature as described in Kristjansson et al. (2000), to predict a terminal velocity. To investigate small particle effects on climate models, Mitchell et al. (2008) used particle size distributions (PSDs) with varying concentrations of small particles to represent low-latitude cirrus. Resulting model runs using CAM differed significantly, showing that GCM results are highly dependent on the realism of the microphysical parameterizations. Brown and Heymsfield (2001) showed significantly different results when ice particle density was changed in a CRM. They found that discrepancies in size for the largest particles between the model and observations for high-altitude anvil cirrus were likely due to excessive autoconversion and aggregation rates for cloud ice. Values of mass-weighted terminal velocity

*V*in the Colorado State University (CSU) System for Atmospheric Modeling (SAM) are parameterized in terms of ice water content (IWC). The SAM

_{m}*V*values for low-latitude cirrus are parameterized from observations taken during the Tropical Rainfall Measurement Mission (TRMM) Kwajalein Experiment (KWAJEX; Heymsfield 2003). The observations shown in Heymsfield (2003) do not extend to temperatures lower than −50°C and are from convectively active regions. Heymsfield et al. (2007) give a summary of commonly used

_{m}*V*parameterizations based on cloud IWC.

_{m}Particle terminal velocities have been the subject of investigations for decades. Jayaweera and Cottis (1969) measured the terminal velocities of plate- and column-shaped objects in liquids. The object size and shape and the viscosity of the liquid were chosen to match the Reynolds number values common for atmospheric ice particles. Fukuta (1969) and Kajikawa (1973) measured the terminal velocities of plate- and column-shaped crystals photographically in a laboratory chamber. Locatelli and Hobbs (1974) experimentally measured the terminal velocities of precipitation-sized particles.

In addition to experimental work, there have been numerous theoretical studies to estimate ice crystal terminal velocities. Using drag coefficients and particle properties, Heymsfield (1972) calculated terminal velocities that agreed reasonably well with experimental values. Heymsfield and Kajikawa (1987) calculated theoretical terminal velocities for precipitation particles and compared them to observations. Mitchell (1996, hereafter M96) used power-law relationships between the Reynolds number and the Best number to estimate particle terminal velocities. M96 defined four power-law relationships between the Reynolds number and Best number to cover the range of particle sizes present in the atmosphere. Khovorostyanov and Curry (2002) found a continuous formulation, which was further refined by Mitchell and Heymsfield (2005) to more accurately treat aggregates. The advantage of the Reynolds number to Best number approach for calculating particle terminal velocities is that a continuous function can be determined for all cloud particle sizes present in the atmosphere with knowledge of their projected area and mass (Khovorostyanov and Curry 2002).

A major challenge in calculating accurate terminal velocities is that the morphology of small particles in low-latitude tropopause cirrus is currently poorly characterized. Heymsfield (1986) showed ice particle replicas collected by an aircraftborne Formvar replicator in the Marshall Island region including pristine vapor-grown plates, trigonal plates, and hollow columns. In contrast, Lawson et al. (2007) showed significant concentrations of irregularly shaped particles as well as pristine plates in subvisible cirrus from the Costa Rica Aura Validation Experiment (CR-AVE) with lower resolution electronic imaging.

This study reports on the microphysical properties of low-latitude tropopause cirrus. The measurements analyzed are from time periods when the aircraft was sampling decaying anvil cirrus or in situ generated cirrus near the tropopause. The results of this study are applicable to low-latitude and subvisible cirrus clouds at temperatures lower than −56°C. For the conditions sampled, it is shown that commonly used parameterizations poorly estimate particle properties. This is partially because the parameterizations were developed from measurements taken in convectively active or thick stratiform clouds at warmer temperatures and then extrapolated to colder temperatures where particle properties are different. The thin cloud measurements used in this study were sampled in nonconvective regions and are composed mostly of particles smaller than 200 *μ*m. Dataset properties are discussed in section 2. Microphysical properties of the cloud particles are described and terminal velocities are calculated in section 3. In section 4, the properties of low-latitude cirrus PSDs are discussed and parameterizations for fits are presented. Section 5 shows ensemble properties and parameterizations for the modeling and remote sensing communities. Conclusions and recommendations are presented in section 6.

## 2. Low-latitude tropopause cirrus dataset

### a. Instruments and data processing

The VIPS probe was used to measure ice particle sizes from 10 to 350 *μ*m for both CRYSTAL-FACE and Pre-AVE on the NASA WB-57 aircraft. The VIPS works by allowing ice crystals to impact a transparent moving belt, which is then imaged and recorded by a video microscope. The images are digitized from the video and analyzed with image analysis software (MH97). Particle breakup, that is, the breakup of large particles on probe inlets, is thought to be insignificant for the time periods analyzed owing to the low to nonexistent concentrations of large particles (*D* > 250 *μ*m). Derived properties for each particle include their projected area and maximum dimension (defined as the maximum dimension measured by the particle probe), which are binned into size distributions. Results are generally averaged for 5 s, equivalent to 800 m of flight at typical WB-57 flight speeds (∼160 m s^{−1}). The physical dimensions of the opening for crystals in the VIPS as well as the aircraft speed determine the sample volume of the instrument. The usable sample volume of the VIPS is approximately 0.84 L s^{−1} at 160 m s^{−1}. The VIPS was optimized for thin clouds, and thick clouds quickly covered the field of view of the instrument. Particle collection efficiency is reduced for particles smaller than 10 microns, although calculations show that the collection efficiency is 94% or greater for larger particles (Ranz and Wong 1952).

For CRYSTAL-FACE, the WB-57 was equipped with a Cloud Aerosol Precipitation Spectrometer (CAPS) probe designed by Droplet Measurement Technology (DMT). The CAPS probe is composed of several detectors, including the cloud imaging probe (CIP), a 2D optical array probe (OAP) similar to the Particle Measuring Systems (PMS) two-dimensional cloud probe (2DC) but with faster electronics. The CIP was used to measure the PSD from 75 microns to several millimeters. The CIP, with a sample volume of 15 L s^{−1} for 200-micron particles at 160 m s^{−1}, is well suited to measure low concentrations of 200-*μ*m and larger particles. Data from the CIP were processed as described in Heymsfield et al. (2002a) with additional processing to remove potential artifacts from particle shattering (Field et al. 2006). In the overlapping size range (100 to 200 microns), the CIP and VIPS probes were found to agree well. For the Pre-AVE data used in this study, the VIPS was the only available instrument for PSD measurement. Given that the Pre-AVE flight was the coldest (between −76° and −86°C), it is believed that the VIPS accurately measured the entire size range present in the cloud as there were few particles larger than 150 *μ*m. Additional measurements used in this study were made by the University of Colorado closed-path tunable diode laser hygrometer (CLH) and the Harvard University Lyman-*α* total water photofragment-fluorescence hygrometer (HT). Both instruments measure ice water content; the results were combined with PSD measurements to estimate individual particle masses.

### b. Observations

The dataset used in this study comprises data from three flights. The first part is comprised of measurements from the 23 July 2002 CRYSTAL-FACE flight. During this flight, the NASA WB-57 investigated an area of convection over the Florida peninsula. The WB-57 repeatedly sampled anvil cirrus in varying stages of development. The temperature during the passes was −56°C to −65°C at an altitude between 12 and 13 km. Initial passes through the cloud field intercepted decaying anvil cloud from convection on previous days. Later in the flight, more recently generated anvil clouds were intercepted as well as fresh convection.

The second part of the dataset is from WB-57 measurements on 26 July 2002, also during CRYSTAL-FACE. The WB-57 flew south over the Caribbean Sea and sampled tropopause cirrus near the Honduran coast. Satellite data showed gravity waves containing cirrus that radiated away from Honduran convection. The WB-57 descended from 15 to 12 km through the cloud from a temperature of −75°C to −60°C. Airborne laser infrared absorption spectrometer (ALIAS) HDO (a heavy isotope of H_{2}O) data suggest that the cloud particles were formed from local vapor rather than having been uplifted from lower altitudes (Webster and Heymsfield, 2003), supporting the gravity wave formation conjecture.

During the CRYSTAL-FACE project, it was necessary to make some modifications to the VIPS that resulted in the majority of particles being slightly out of focus. The modifications caused the film strip to be partially out of the plane of focus on one side of the field of view while the other side of the field of view remained in focus. To correct the maximum dimension measurements, a subset of the data was analyzed at a much higher sample rate, allowing individual particles to be imaged at different points as they moved across the field of view. The “in focus” size was compared to the “out of focus” sizes at different locations in the field of view. It was found that the focus affected the particle size predictably, and a parameterization was developed to correct particle sizes based on the location of the particle within the field of view. The final processing of the CRYSTAL-FACE data was limited to areas of the film that had a correction factor of 15% or less for particle length.

The third set of VIPS measurements is from a Pre-AVE flight when the WB-57 was ferrying from Houston, Texas, to San Jose, Costa Rica, on 24 January 2004. The WB-57 intercepted a tropopause cirrus layer between 15 and 17 km in altitude. After 10 minutes of sampling at the 15-km level the aircraft began a slow ascent to the top of the cloud layer near 17 km. The temperature of this layer was between −76° and −86°C. For Pre-AVE, the VIPS had been modified to correct the focus problem described in the CRYSTAL-FACE section.

## 3. Terminal velocities of low-latitude tropopause cirrus particles

This section will focus on estimating the values of projected area *A* and mass *m* from the dataset as area and mass are necessary values for the Reynolds number–Best number *V _{t}* calculation (Re–

*X*). Emphasis will be placed on sub-200-

*μ*m particles since they make up the bulk of the particles measured in the dataset. The Re–

*X*approach will then be applied and the developed

*V*relationship will be compared to common

_{t}*V*parameterizations. A brief explanation of the Re–

_{t}*X*approach is given in the appendix.

### a. Particle projected area

Area to maximum dimension (or area ratio) relationships are common in the literature. The area ratio is the measured projected area of the particle as imaged by the probe divided by the area of the smallest circle that would completely cover the particle image. M96 summarizes several different maximum dimension to area parameterizations. For the dataset, one area to maximum dimension relationship was determined from measurements and used for the remainder of the calculations. Figure 1 shows the datasets used for calculation of the maximum dimension to area ratio relationship. Figures 1a and 1b show the relationship between particle area ratio and maximum dimension for particles observed by the cloud particle imager (CPI) (Fig. 1a) and CIP (Fig. 1b) during the 23 July CRYSTAL-FACE time periods. The bold lines in each panel represent the parameterized values used in the remainder of the study. The lighter lines represent the median values for 10-*μ*m size bins as well as the median ±1 std dev. In further calculations, an uncertainty of ±20% was used because one standard deviation was generally 20% different from the median values for the bulk of the particles measured in the dataset. The CIP area ratio measurements are untrustworthy for sizes smaller than 200 *μ*m owing to the resolution. At larger particle sizes the results from the two probes begin to converge. The CPI area ratio measurements are reasonable for sizes as small as 50 *μ*m, where focus starts to decrease the image quality. The area ratio values measured from VIPS data for Pre-AVE (not shown) were similar in values and variability to the CPI data but were restricted to the small particle sizes observed and therefore not of use for particles larger than 100 *μ*m. For the VIPS data the median area ratio value was 0.78 for sizes smaller than 100 *μ*m, which is in agreement with the bold line parameterized value. Embryonic vapor-grown particles are thought to be droxtal-like in shape (Zhang et al. 2004). Droxtals have area ratio values of about 0.87, depending on the specific formation, and as such the area ratio values used in this research effort were chosen to reach 1.0 for particle embryos. Figure 1c shows the parameterization developed from this dataset along with several other parameterizations from the literature. The shaded area represents the uncertainty in area used in the remainder of the calculations.

An exponential fit was deemed more suitable for the area ratio to maximum dimension relationship for particles 200 *μ*m and smaller. This choice was made due to the trend observed in the CPI imagery area ratio values from 50 to 200 *μ*m, along with the expectation that the smallest particles are droxtal-like and not spherical. The exponential relationship fit to the CPI data is coincident with the area ratio expected for droxtals in the size range expected (20–30 *μ*m). For particles larger than 200 *μ*m a power-law relationship was fit to the convergence of the CIP and CPI particle measurements. The equations for the fit to the dataset are given below:

where Ar is the area ratio and *D* is the particle maximum dimension (in cm).

### b. Particle mass

The second critical component of individual particle terminal velocity calculations is the individual ice particle mass. Particle mass is often parameterized through mass–dimensional (m-D) relationships, which relate the particle mass to the maximum dimension and can be derived from direct measurements of particle mass and dimension. Mitchell et al. (1990) developed several relationships for snow collected at the surface in the Sierra Nevada mountain range in California. This work was furthered by Baker and Lawson (2006, hereafter BL06), who included additional dimensions such as particle width and perimeter into their improved m-D relationship. The BL06 power-law relationship containing perimeter could not be applied to this dataset due to the poor focus for part of the dataset rendered the perimeter measurement impossible.

Indirect estimates of mass–dimensional relationships are sometimes necessary because it is not possible to measure individual particle mass from two-dimensional images. Typical m-D relationships are power-laws of the form *m* = *aD ^{b}*, where

*m*is the mass,

*a*and

*b*are constants, and

*D*is the maximum measured dimension. A commonly used m-D relationship was developed by Brown and Francis (1995, hereafter BF95). The BF95 relationship is typical in that smaller particle masses (

*D*< 100

*μ*m) must be limited to the mass of spherical particles. Since OAPs and scattering probes cannot provide detailed images of particles this small, potential errors are often ignored. Spherical particles would necessarily have area ratio values of 1.0, which Fig. 1 demonstrates is not the case. While the sphericity assumption does not significantly affect integrated properties for broad PSDs, it may lead to significant errors for integrated properties of narrow PSDs.

Individual particle mass is often expressed in terms of density. To calculate particle density, the mass of an ice particle is divided by the volume of a sphere with the same maximum dimension. By measuring volume in cubic centimeters, the units are expressed as grams per cubic centimeter. The particle population-mean effective bulk density () (Black 1990; Heymsfield et al. 2004) can also be calculated to deduce a reasonable range for particle mass. Here is calculated similarly to individual particle density. The mass of all of particles in the PSD, as measured by a total water instrument, is divided by the summed volume of equal maximum dimension spheres for all of particles in the PSD. The directly measured IWC measurements from the CLH or HT probes are used in conjunction with the volume calculated from the PSD measurements. For the CRYSTAL-FACE time periods, the CLH measurements were used for the calculations. Even though the CLH has been reported to have a higher uncertainty (Davis et al. 2007) than the HT (Weinstock et al. 2006), it was noted that the relationship between the CLH measured IWC and the IWC calculated from the forward scattering spectrometer probe (FSSP) was more consistent than for the HT measurements at very low values. As the CLH instrument was not on the aircraft for Pre-AVE, the HT measurement was used for that portion of the dataset. In general, the HT measurements were higher than for the CLH.

Table 1 shows the values calculated directly from the PSDs and the direct IWC measurements. The median values of calculated from the PSD measurements and the direct IWC measurements is between 0.185 and 0.39 g cm^{−3} for the different measurement periods in the dataset. For comparison, the for a period when the WB57 entered a fresh anvil composed of thin plate-shaped particles and aggregates at −55°C was 0.11 g cm^{−3}. The standard deviations are also shown in Table 1 as well as the values calculated when assuming systematic +30% or −15% uncertainties. These uncertainty values were chosen because they produced values that were plausible for all flights. Although the values seem high compared to other data (Heymsfield et al. 2004), they are low compared to values that would be calculated by assuming the BF95 particle mass estimates. The bottom two rows show the percent difference between the measured IWC values and the BF95 calculated as well as the BL06 mass to maximum dimension relationship. This highlights the problems caused by assuming that all particles smaller than 100 microns are solid ice spheres, which is common when power-law m-D relationships are developed using primarily large particles. The disagreement with BL06 was expected as the dimensional properties of the measured BL06 particles (larger than 200 microns measured in snowstorms at the surface between 0° and −5°C) are unlikely to be similar to this dataset (generally smaller than 200 microns, observed at −56° to −86°C).

To calculate m-D relationships from aircraft data, it is necessary to have good PSD and IWC measurements. One then chooses a functional form for the m-D relationship and iteratively determines which values in the m-D relationship lead to the best estimates of the measured IWC when integrated over the PSD (BF95; Heymsfield et al. 2004). For example, if the power-law form of the m-D relationship is to be used (*m* = *aD ^{b}*), one would make an estimate of the likely value of the

*b*parameter and then optimize the

*a*value so that the mean IWC calculated by integrating the m-D relationship over a series of PSD measurements would correctly predict the mean value of the directly measured IWC. The standard deviation between the measured IWC and the values calculated from the PSDs and the m-D relationship is then calculated. A new

*b*value, and subsequent

*a*value, are then calculated in an effort to find a better fit to the data, which is determined by the

*a, b*pair with the smallest standard deviation. For highly variable data, it is best to calculate the standard deviation in logarithmic space so as not to overemphasize the measurements with the highest IWC values.

Although the exact morphology of sub-100-*μ*m particles is poorly understood, small particle imagery clearly indicates that the assumption of ice spheres is a poor choice. Reasonable density values for pristine crystals range from 0.75 g cm^{−3} for droxtal-shaped ice crystals (Zhang et al. 2004) to less 0.1 g cm^{−3} for the pristine hollow column-shaped particles shown in Heymsfield (1986). Since the CPI and VIPS show quasi-spherical irregular particles in the 20- to 50-*μ*m range and the values for the PSDs are as high as 0.4 g cm^{−3}, it seems reasonable to characterize the particles as more droxtal-like. For these reasons, an exponential representation of particle density for particles smaller than 200 microns was chosen, as an exponential is expected to provide an improved representation of particle mass values for small particles. For particles larger than 200 microns, a power-law relation is used. The power-law parameters were chosen so that the combined function was continuous at the junction at 200 *μ*m.

The m-D relationship was calculated by varying the parameters of the exponential and power-law components of the m-D relationships while retaining agreement at 200 microns and then comparing the resulting total mass integrated across the PSD to the directly measured IWC. The combined function with the smallest error between measured and calculated IWC for the dataset was selected. The m-D relationship is shown in Eq. (2) and is the bold line plotted in Fig. 2:

where *m* is the particle mass in grams and *D* is the particle maximum dimension centimeters. The *b* value of 1.9 in the power-law portion of Eq. (2) was chosen because it is thought to be a representative value of the fractal properties of aggregates (Heymsfield et al. 2007). As total mass in the dataset was dominated by *D* < 200 *μ*m particles, the *D* > 200 *μ*m portion of Eq. (2) by itself should be considered to have high uncertainty. It will be shown later that particles with *D* > 200 *μ*m have very little influence on the results of this work. Figure 2 also shows several parameterizations from the literature as well as a power-law fit to the dataset (thin solid line). The best power-law fit to the dataset (*m* = 0.001*D*^{1.9}) had a similar error to that found for Eq. (2), although the power-law–only approach leads to particles smaller than 37 *μ*m being considered ice spheres, which contradicts the Heymsfield (1986) images. The shaded region in Fig. 2 represents the range of density values that could be expected with reasonable estimates of uncertainty in the IWC measurement and values. The upper limit is calculated using a possible uncertainty of +30% for the IWC whereas the lower limit is calculated using an estimate of −15%. A 30% uncertainty is within the published uncertainty range for the CLH for low IWC values, but it is higher than the uncertainty published for the HT instrument (−18%).

The choice of exponential functions to represent mass and area properties of particles may seem unusual, although it does have advantages. As shown in Figs. 1c and 2, parameterizations using power-law relationships for mass and area can intersect with the maximum reasonable values of 1.0 for area ratio and 0.91 g cm^{−3} for density. When calculating particle terminal velocities, it is important to have area and mass relationships that are reasonable in a relative sense (i.e., a particle with the mass of a sphere of the same maximum dimension should have an area ratio of 1.0). Single power-law representations led to two modes in particle terminal velocity calculations. These modes corresponded to the regions where spherical and nonspherical particles were assumed. Attempts to fit the data with two power-law relationships resulted in higher uncertainties when the PSD-integrated values were compared to the measured values.

### c. Particle terminal velocities

The *V _{t}* values were calculated for the observed particle range using the Mitchell and Heymsfield (2005) Re–

*X*technique along with the particle area and mass values reported earlier in this section. The appendix gives a brief description of the Re–

*X*technique as well as a theoretically based

*V*parameterization. All terminal velocities were normalized to 150 hPa and −70°C, typical values for the dataset. Figure 3 shows

_{t}*V*plotted with respect to particle size using the Re–

_{t}*X*method. The dashed line represents

*V*calculated using BF95 particle mass values, and the dotted line indicates

_{t}*V*from a parameterization for unrimed plane dendrites developed by Locatelli and Hobbs (1974). The Locatelli and Hobbs parameterizations are still used today (Garvert et al. 2005), although the error encountered when extrapolated to such small particle sizes is obvious. The asterisks in Fig. 3b represent the laboratory measurements of the mean terminal velocities of pristine plate- and column-shaped particles from the tables in Kajikawa (1973). The error bars surrounding the bold line represent the potential uncertainty in terminal velocity when the area and mass uncertainties discussed in the earlier subsections are considered. The most obvious discrepancy between the theoretical curves and the experimental results are the two Kajikawa data points (nearly overlain) for 20-

_{t}*μ*m particles, which show terminal velocity values near 5 cm s

^{−1}. This uncertainty will also be considered when relationships for mass-weighted terminal velocity are discussed. The reason for this discrepancy is unknown. The 20-

*μ*m Kajikawa particles have significantly higher terminal velocities than would be expected for spherical ice particles. The BF95 curve represents the maximum terminal velocity estimated for solid spherical ice particles smaller than 90

*μ*m using the Re–

*X*method.

## 4. PSD properties of low-latitude tropopause cirrus

In contrast to parameterizations developed in recent publications (Ivanova et al. 2001; MH97), the measured PSDs for this dataset were monomodal. In previous studies, PSDs have typically been measured with a two-dimensional optical array probe for particles larger than 50 microns and a forward scattering-type probe for particles smaller than 50 microns. Results from forward scattering probes have been shown to be contaminated by the shattering of large particles on the leading edges of the probes (Field et al. 2003; Heymsfield et al. 2006). This can lead to artificially high concentrations of small particles in the presence of large ice particles, giving the appearance of bimodality (Larsen et al. 1998). Ryan (2000) points out that the commonly observed transition from the small mode to large mode in a measured bimodal size distribution takes place at the cutoff point between different probes. In this study, the VIPS measured the PSD from 10 to 350 *μ*m, including typical mode transition points, and no signs of bimodality were noted.

Figure 4 shows average PSDs from four different regions of low-latitude tropopause cirrus sampled in the dataset. Each of the PSDs is an average of at least 10 min of data. The values for the time periods shown in Figs. 4a,c,d are those shown in columns 1, 3, and 4 of Table 1. In the case of Figs. 4c,d, the averaged PSDs include data from descents or ascents through cloud layers 2 and 1.7 km thick, respectively. This demonstrates that, even with the vertical averaging necessary to account for model grid points with significant vertical separation, monomodal PSDs are still appropriate for this cloud type. The PSD shown in Fig. 4b is from the top of a fresh anvil. The data from that pass was not highly reliable owing to the VIPS being mostly covered with overlapping particles due to high concentrations for a significant portion of the pass. As such, this portion of the dataset was not used for further calculations, although the monomodal properties of the PSD are of interest. For all of the cloud types sampled, the size distribution is smooth and continuous between 20 and 200 *μ*m.

The results from MH97 are often used to represent low-latitude tropopause cirrus properties. Data used in developing the MH97 parameterizations were discussed in Heymsfield and McFarquhar (1996). MH97 parameterized PSDs using measurements from an earlier version of the VIPS and the PMS 2DC probe for sizing particles 100 *μ*m and larger for flights between −20° and −70°C in convective systems. Only 70 representative VIPS time periods were analyzed owing to the laborious nature of the process, compared to the 480 used in this study. Five seconds of VIPS data at the beginning of each 1-min time period was analyzed and joined with each 1 min of data from the other probes (MH97). Because the research aircraft was not equipped with an IWC instrument, total IWC was estimated using “habit specific” power-law–type m-D relationships. The MH97 PSD parameterization poorly represents the properties of the thin clouds characterized in this study. The moderate bimodality apparent in the MH97 parameterization can be attributed to the presence of higher concentrations of small particles in active systems or could possibly be due to the shattering of large particles on the forward surfaces of the 2DC probe (Field et al. 2006).

Figure 5 shows a representative PSD with parameterized fits using the MH97 distribution, a gamma distribution, and an exponential distribution. The gamma distribution is based on the function *N*(*D*) = *N _{o}D^{μ}e*

^{−λΓD}, where

*N*is the intercept,

_{o}*λ*

_{Γ}is the slope, and

*μ*is the dispersion. The

*μ*parameter of the gamma representation allows the fit curve to bend downward for smaller particles, whereas the exponential fit intersects the concentration axis perpendicularly. Exponential fits are of the same form as the gamma distribution except that

*μ*is 0:

*N*(

*D*) =

*N*

_{o}e^{−λD}. The exponential fit was calculated using the techniques described by Zhang et al. (2007), which takes into account the different moments of the size distribution. For a given second moment (extinction) and third moment (IWC), two families of

*λ*and

*N*values can be calculated. The best exponential fit was taken to be where the two

_{o}*λ*–

*N*curves crossed. This leads to an exponential size distribution parameterization that predicts both extinction and IWC reasonably. For the two moments, the

_{o}*N*to

_{o}*λ*relationships are

and

where *A* is the integrated projected area of the PSD. The values of 94 and 38 come from integrating the exponential parts of Eqs. (1) and (2), respectively.

The second row of plots in Fig. 5 shows the percentage relative error for each representation compared to the measured PSD. This is calculated by taking the difference between the observed and fitted values for each PSD size bin for the entire dataset. This concentration difference and Eq. (2) are then used to calculate the difference in mass between the measured PSD and the parameterized PSD, which is then normalized to the total mass in the PSD. The results show that the MH97 parameterization has up to an order of magnitude more potential error than the exponential or gamma fits, which have similar discrepancies.

Figure 6 shows the parameters for exponential and gamma fits plotted versus temperature. The first and second rows of plots show the relationships between different parameters in the gamma function. The third row shows the parameters for the exponential distribution for the dataset. Where reasonable trends are apparent, trend lines are shown to the median values. Error bars indicate the 25th and 75th percentile values in each grouping. Where fit lines are drawn, equations for the fits to the median values as well as fits to the 25th and 75th percentiles are given in Table 2. These relationships are significantly different from those derived in other publications (Heymsfield et al. 2002a; Heymsfield 2003), primarily because the clouds measured in this study were much thinner than those measured in previous studies.

## 5. Ensemble properties and derived relationships

The mass-weighted terminal velocity *V _{m}* for each PSD is calculated by summing for each size bin the product of the mass and

*V*(ΣIWC

_{t}_{bin}

*V*

_{t bin}), then dividing by the total mass for the PSD. To better understand the importance of different size particles in the calculation of

*V*, IWC

_{m}_{bin}

*V*has been cumulatively summed, then normalized by the total IWC for each size distribution. This calculation shows the contribution to

_{t}*V*made by different sizes, and a contour plot showing the 10th through 90th percentiles for the dataset is shown in Fig. 7. This demonstrates that mass flux in these clouds is dominated by the critical yet less well-understood 15–60-

_{m}*μ*m range. This plot also demonstrates, as mentioned earlier, that the parameterized mass of particles larger than 200

*μ*m is of little importance when calculating

*V*for the dataset.

_{m}The CSU System for Atmospheric Modeling parameterizes *V _{m}* by IWC. Figure 8a shows the calculated

*V*values and the SAM parameterized values plotted with respect to IWC. The shaded area represents the uncertainty in

_{m}*V*based on the uncertainties in the

_{m}*V*calculation discussed earlier. The results based on the VIPS derived values are significantly lower than those predicted by the SAM parameterization. The SAM parameterization was developed from Heymsfield (2003) with a treatment of particle mass that was similar to the BF95 values in thicker convective and stratiform clouds. By using the

_{t}*V*uncertainty ranges mentioned earlier, upper and lower limits can be placed on

_{t}*V*. Parameterizations for the median values and the uncertainty ranges for this and following relationships are shown in Table 3. The uncertainty implied by the 20-

_{m}*μ*m Kajikawa (1973) measurement was also investigated by arbitrarily setting a minimum

*V*value of 5 cm s

_{t}^{−1}and recalculating the

*V*values. In general, the median values used to make the parameterizations were affected by less than 0.5 cm s

_{m}^{−1}, indicating that the uncertainty for particles that small is insignificant. The use of power-law functions for mass and area dimensional relationships [

*m*= 0.001

*D*

^{1.9}for mass; Fig. 1c (dotted line) in Heymsfield and Miloshevich (2003) for area] increased the scatter in Fig. 8a and led to

*V*values about 5% higher than those calculated using exponential relationships.

_{m}Effective particle sizes are being estimated by the Moderate Resolution Imaging Spectrometer (MODIS) onboard the *Terra* and *Aqua* satellites. Effective particle radius *r _{e}* can be defined as

*r*= 0.75 IWC/

_{e}*Aρ*(Mitchell 2002), where

_{i}*A*is the particle projected area and

*ρ*is the density of ice. As the ratio of mass to area is a fundamental part of the

_{i}*V*calculation, it is not surprising that there is high correlation to an

_{t}*r*to

_{e}*V*relationship (Fig. 8b). The parameterization for the relationship between

_{m}*r*and

_{e}*V*, as well as equations to calculate the boundary of the uncertainty ranges, is also shown in Table 3. With the results from Cloud–Aerosol Lidar and Infrared Pathfinder Satellite Observations (CALIPSO), it is now possible to estimate the global distribution of tenuous clouds such as those sampled for this study. This study is also useful for the development of parameterizations related to CALIPSO measurements for low-latitude cirrus and anvil cirrus clouds near the tropopause. Figures 9a and 9b show the relationship between extinction (

_{m}*σ*) and

*V*and the relationship between

_{m}*σ*and IWC. Relationships to the medians as well as to uncertainty ranges are shown in Table 3.

## 6. Conclusions and recommendations

A detailed analysis of a low-latitude tropopause cirrus cloud microphysical dataset has been completed, but clearly more work is needed to refine the estimates. High-quality microphysical data measured by the VIPS probe during several flights are used to characterize cloud particle properties including particle sizes in the range poorly measured by other probe combinations. Cloud particle properties necessary for accurate estimates of particle terminal velocity have been thoroughly examined. It was found that common parameterizations are inadequate when applied to particle populations composed mainly of particles smaller than 200 *μ*m. Exponential functions were found to provide a more realistic representation of the area ratio and density when high-quality imagery was considered. The use of exponential functions for mass and area led to more continuous *V _{t}* estimates, and the

*V*parameterization showed less uncertainty when exponential representations were used. Mass-weighted terminal velocities were found to be significantly lower than common low-latitude parameterizations, mainly due to the small sizes of the particles observed, as common parameterizations were developed using data from thicker clouds with higher concentrations of large particles.

_{m}Size distributions measured by the VIPS were continuous through the size range where the small and largeparticle modes often intersect in bimodal parameterizations. Particle size distributions were found to be well modeled and parameterized by gamma distributions, whereas parameterizations for exponential representations were not as universal. Particle size distributions and particle terminal velocities were used to calculate mass-weighted terminal velocities for particle populations. Because particles were generally smaller and predicted to have lower masses, *V _{m}* values were significantly lower than parameterizations that were developed from data collected in convective regions. Parameterizations of

*V*and IWC have been developed based on cloud extinction, which will be useful in analyzing CALIPSO data.

_{m}Recommendations are as follows:

Tropopause thin cirrus layers either in situ generated or anvil generated are significantly different from thicker clouds generated by fresh convection and should be represented differently in models.

Gamma functions should be used to represent PSDs. Relationships between state parameters and gamma distribution parameters are more robust than for exponential relationships and were found to be more universal throughout the dataset.

To account for the nonsphericity of small crystals that has been observed in high-quality imaging data of high-altitude, low-latitude cirrus (Heymsfield 1986) exponential representations of particle mass and area to maximum dimension are thought to be an improvement over single power-law representations, which assume that small particles are spherical. Exponential representations of mass–and area–dimensional relationships are consistent with observations, although significant uncertainty remains for particles smaller than 50

*μ*m.More research is needed to better understand the mass and area dimensional relationships for sub-200-

*μ*m tropopause cirrus particles. These relationships are critical for calculation of particle terminal velocity. The current uncertainty leads to significant potential error in low-latitude tropopause cirrus lifetime estimates.

## Acknowledgments

This research was supported by CSU Contract G-3045-9 under the ATM Prime Award 0425247 and NASA Grant NNX07AQ85G. The authors thank Leslee Schmitt for editing the manuscript.

## REFERENCES

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

^{18}O/

^{16}O,

^{17}O/

^{16}O in and out of clouds map dehydration pathways.

**,**

**,**

**,**

### Appendix

#### The Reynolds–Best Number Approach for Calculating Vt

By balancing the gravitational force with the force caused by drag on a falling particle, the terminal velocity *V _{t}* of the particle can be shown to be

where *m* is the particle mass, *g* is the gravitational force, *ρ _{a}* is the density of air,

*A*is the particle projected area in the direction of fall, and

*C*is the coefficient of drag (M96). Since

_{D}*C*is velocity dependent,

_{D}*V*is often calculated using the Reynolds–Best (Re–

_{t}*X*) number approach. A full description of the Re–

*X*approach is omitted for brevity. The reader is encouraged to consult Khovorostyanov and Curry (2002) for a description of the approach and Mitchell and Heymsfield (2005) for the modifications that were included in the calculations used in this study. Briefly, the Re–

*X*approach is based on boundary layer theory and the comparison of empirically derived power-law expressions relating the Reynolds number to the Best number. The Best number is

where *ν* is the kinematic viscosity and *D* is the particle maximum dimension (Heymsfield 1972); *X* is calculated using *D*-dependent mass and projected area relationships. The Mitchell and Heymsfield modifications to the Khovorostyanov and Curry (2002) relationship between Re and *X* are used in this study for the estimation of Re, which is then used to calculate *V _{t}* from Re =

*V*/

_{t}D*ν*. The Re–

*X*technique is applicable in both inertial and viscous-dominated flow regimes. M96 shows that the Re–

*X*approach reduces very nearly to the Stokes

*V*when the area and mass of small solid spheres are considered.

_{t}A theoretical value of *V _{t}* can be developed by noting that M96 showed that

when Re = *αX ^{β}*. By entering the

*m*and

*A*values, developed earlier, Eq. (A3) becomes

Noting that *β* = 0.97 for *X* values smaller than 10, or about *D* < 100 *μ*m [from Eq. (18) in M96], it is reasonable to consider the loss of accuracy in *V _{t}* by assuming

*β*= 1 in exchange for the flexibility gained by eliminating

*ρ*from the equation. (There is an additional

_{a}*ρ*buried in the viscosity term.) By fitting a linear relationship (

_{a}*β*= 1) to the calculation of

*X*and Re values shown in Fig. A1a for the particle range observed in the dataset, it can be shown that

When Eq. (A5) is plotted onto Fig. A1a, it is indistinguishable from the other curves until *X* increases beyond 10, where Eq. (A5) deviates slightly to values that overpredict Re. Combining Eqs. (A4) and (A5) gives the approximation

for −70°C and 150 mb: *V _{t}* is in units of centimeters per second and

*D*is in centimeters. All other constants, or near constants, are combined into the numeral 217 600. The exponent to the

*D*parameter theoretically would be 2.0 when Eqs. (A4) and (A5) are combined, but it was found that an adjustment to 1.9 (with a corresponding adjustment to the constant factor) fit the

*V*values from the Re–

_{t}*X*calculations significantly better. Figure A1c shows the calculated

*V*values (same as bold line in Fig. A1b) along with the parameterized values from Eq. (A6). The uncertainty shown by the error bars in Fig. A1b can be estimated by substituting 137 500 or 320 000 for the 217 600 in Eq. (A6). Figure A1c shows the discrepancy between the Re–

_{t}*X*calculated

*V*and Eq. (A6). Equation (A6) can be adjusted for different pressure levels by

_{t}*V*∝ [

_{t}*ρ*

_{0}/

*ρ*]

_{a}^{0.54}, where

*ρ*

_{0}is the air density at 150 mb and −70°C and

*ρ*is the air density at the level of interest. Equation (A6) should not be used for particles larger than 200

_{a}*μ*m as Eq. (A5) does not apply to larger particles and the area and density parameterizations used in deriving Eq. (A6) were for the sub-200-

*μ*m range only. Equation (A6) should not be used for active convective regions.

## Footnotes

*Corresponding author address:* Carl Schmitt, 3450 Mitchell Lane, P.O. Box 3000, Boulder, CO 80301. Email: schmittc@ucar.edu