## Abstract

In this paper, an empirical equation is presented that can be used to estimate shortwave cloud optical thickness from measurements and analysis of shortwave broadband irradiances. When applied to a time series of broadband observations, this method can predict cloud optical thickness distributions that are very similar to those obtained using the algorithm of Min and Harrison (henceforth the Min algorithm). For a given site, medians of the Min algorithm–derived and empirically derived distributions differ by less than 10%. This level of agreement holds over a wide geographical range of sites. The equation is designed for fully overcast skies, surface albedos less than 0.3, and a cosine of the solar zenith angle greater than 0.15.

## Introduction

Considerable effort has been expended, and continues to be expended, to characterize the spatial and temporal distribution of the shortwave cloud optical thickness over the globe. This effort relies on satellite and ground-based measurements of the reflection or transmission of solar radiation from which one can infer cloud optical thickness. Transmission-based algorithms use spectral irradiances (Min and Harrison 1996; henceforth referred to as the “Min algorithm”) or broadband irradiances (Leontyeva and Stamnes 1994; Dong et al. 1997; Boers 1997; Barker et al. 1998). These algorithms assume a plane-parallel homogenous cloud layer, and as noted by Leontyeva and Stamnes (1994), they calculate an “effective” cloud optical depth *τ*_{c} of an imaginary plane-parallel, homogenous cloud layer that produces the same radiation field at the surface as the actual cloud layer.

The Min algorithm runs relatively quickly on a desktop personal computer (PC) but requires the total and diffuse components of the spectral irradiance at one wavelength, usually 415 nm. Unfortunately, this information is often not readily available, thus, limiting the opportunity of using the Min algorithm at a wide variety of sites. On the other hand, shortwave broadband irradiances are commonly measured, but the algorithms that use these measurements are more difficult to use because they require information about the atmospheric state (e.g., vertical profiles of water vapor), and the algorithms are often computationally very slow. These algorithms are therefore impractical to use for routine calculations of *τ*_{c} over many sites and over long periods of time. With this limitation in mind, it is desirable to find another method of calculating *τ*_{c} that is computationally quick while retaining the flexibility of using only broadband shortwave irradiances as input. To this end, we present here a simple empirical relationship for estimating *τ*_{c} from broadband irradiances.

This relationship can be thought of as a summarization of a more complicated algorithm, somewhat in the vein of a lookup table or a parameterization based on the output of an algorithm. For example, Dong et al. (1998) and Dong and Mace (2003) present simple parameterizations for the cloud droplet effective radius *r*_{e} and *τ*_{c} valid for midlatitude boundary layer stratus clouds or arctic stratus clouds. These parameterizations were developed by running a broadband algorithm (Dong et al. 1997) to derive *r*_{e}. Then the algorithm's output was used to derive a simple equation for *r*_{e} based on liquid water path (LWP), a transmission ratio (the ratio of the measured total downwelling shortwave irradiance to the clear-sky downwelling irradiance *C,* described below), and the cosine of the solar zenith angle *μ*_{0}. Cloud optical depth follows immediately using the relationship *τ*_{c} = (3LWP)/(2*r*_{e}).

This approach requires measurements of the LWP, which are often not available, thus limiting the use of these parameterizations to sites where LWP is known. In this paper, we take an even more simple approach by deriving an equation for *τ*_{c} that does not require LWP; the only inputs required by our equation are surface shortwave broadband irradiances, thus opening up the possibility of estimating *τ*_{c} at measurement sites where shortwave broadband irradiances are available but other, ancillary measurements (e.g., LWP) are not.

The potential user of the technique presented here should be aware of several caveats. First, we must emphasize that cloud optical thickness obtained from our simple equation should be considered an estimate only, and for some applications one might want to consider the use of a full-blown algorithm (although to what extent such an algorithm can accurately calculate the optical thickness is unclear). To aid potential users in deciding whether our method is appropriate for their needs, we evaluate its accuracy in section 4. Second, our approach is based on algorithms that invoke the plane-parallel, homogenous cloud assumption that is most closely satisfied for fully overcast skies. Therefore, we cannot consider cases of broken cloudiness here. Last, the measured surface irradiances are influenced by the transmission through the cloud, as well as the transmission through the atmospheric layers above and below the cloud. Without detailed information of the vertical distribution of the atmosphere's optical properties, it is not possible to distinguish between the total atmospheric column optical depth and the optical depth of the clouds in the column. This difference is often ignored in cloud optical depth algorithms and we ignore it here.

## Procedure

The goal of our empirical approach is to find some approximate functional relationship between the solar broadband irradiances and cloud optical thickness,

where *r,* to be determined, is some combination of broadband irradiances and the cosine of the solar zenith angle *μ*_{0}, and *τ*^{e}_{c} is the estimated cloud optical thickness (explicitly denoted by the superscript “*e*”). The variable *r* and the function *f* are chosen to minimize the difference between *τ*^{e}_{c} and a set of effective cloud optical thicknesses *τ*_{c}, derived from some algorithm; this set forms our standard of comparison. To construct this set, we must choose an algorithm to calculate the *τ*_{c} that compose this set. Two possible choices are the narrowband Min algorithm, which provides an optical thickness at 415 nm, and a broadband code [briefly described in Barnard et al. (2001)] based on the Santa Barbara Discrete Ordinates (DISORT) Atmospheric Radiative Transfer (SBDART) model (Ricchiazzi et al. 1998). The broadband code outputs an optical thickness at 550 nm, while for the code's internal use, cloud optical thicknesses at other wavelengths are found using the relationship *τ*_{λ} = *τ*_{550nm}(*Q*^{e}_{λ}/*Q*^{e}_{550nm}), where *τ*_{λ} is the cloud optical thickness at the wavelength *λ,* and *Q*^{e}_{λ} is the extinction efficiency.

Because the cloud optical thickness shows a very weak dependence over solar wavelengths [i.e., (*Q*^{e}_{λ}/*Q*^{e}_{550nm}) ≈ 1; Hu and Stamnes 1993], we would expect that the two algorithms would yield similar results. This is indeed the case except for very low optical thickness (0 < *τ*_{c} < 10), for which the Min algorithm produces optical thickness that are slightly larger than the broadband optical thicknesses. This difference is about Δ*τ*_{c} ≈ 1 for *τ*_{c} equal to 2, and it gradually decreases to zero as *τ*_{c} approaches 10. This difference scarcely affects distributions of *τ*_{c}; for example, when applied to data taken from the Atmospheric Radiation Measurement (ARM) site at Barrow, Alaska, during the summer of 2002, the median values of the distributions of *τ*_{c} are 12.0 and 11.9 for the broadband and narrowband algorithms, respectively. Because computations using the Min algorithm are much quicker than the broadband algorithm, and because we have no reason to prefer one algorithm over the other, we adopt the Min algorithm as our standard of comparison to find the function *f* and the variable *r* of (1).

The variable *r* was determined in a two-step process. First, we know that *τ*_{c} is related to the transmission *T* of diffuse irradiance, where *T* is defined as the ratio of the downwelling diffuse irradiance at the surface *D,* divided by total downwelling irradiance at the top of the atmosphere *I*_{o}. Because *T* is related to the cloud optical thickness, we expect that *r* would depend on *T.* The constant *I*_{o} is the solar constant multiplied by *μ*_{0} and is presumably well known. The accuracy of the diffuse irradiance measurement *D* depends on the calibration of the broadband sensor that measures it, and for a well calibrated instrument the measurement error is thought to be about 5 W m^{−2} for cloudy skies (Dutton et al. 2001). This level of uncertainty can be significant in the presence of thick clouds when the diffuse broadband irradiance can be much less than 50 W m^{−2}. Under such conditions, calculations of *T* could be seriously in error. However, as noted by Dong et al. (1997), this issue can be avoided by calculating a proxy for *T.* This proxy is formed by substituting the clear-sky total irradiance *C* for *I*_{o} in the ratio *D*/*I*_{o}; the modified transmission is therefore *D*/*C.* An estimate of the total downwelling irradiance that would be measured at the earth's surface in the complete absence of clouds *C* is obtained from the surface measurements using the shortwave flux analysis of Long (2001) and Long and Ackerman (2000). Because *C* is derived from the same instrument suite that measures *D,* the error in the ratio *D*/*C* attributable to calibration error is significantly reduced. In addition to this reduction, the use of the ratio *D*/*C* instead of *D*/*I*_{o} has another important advantage that strengthens the accuracy of our empirical method: seasonal fluctuations of unknown atmospheric variables (e.g., water vapor, aerosol optical thickness) do not significantly degrade the method's accuracy. This advantage will be discussed in section 3.

The final step in defining *r* is finding a way to include the effect of the solar zenith angle in its definition. Motivated by the desire for simplicity, we assume that *r* has the form, *D*/(*C**μ*^{α}_{0}), where *α* is a constant, to be found using a trial-and-error process. To find a specific value of *α,* we take data from ARM sites (these data are described below), assume a value of *α,* and then use the Min algorithm to calculate *τ*_{c} for each value of *r.* We denote the Min algorithm–derived optical thickness as *τ*^{M}_{c} to distinguish it from the estimated optical thickness *τ*^{e}_{c}, defined in (1). Taking *τ*^{M}_{c} as truth, we plot *τ*^{M}_{c} versus *r* over a wide range of *r,* and we hope to find a particular value for *α* such that the *τ*^{M}_{c} collapse on a single curve (as well as possible), implying a near-one-to-one relationship between *r* and *τ*^{M}_{c}. We must point out that the very simple form assumed for *r* may slightly degrade the supposed one-to-one relationship because for a fixed *τ*^{M}_{c}, as *μ*_{0} changes *r* may vary slightly. However, the results shown in sections 3 and 4 suggest that this residual dependence of *τ*^{M}_{c} on *μ*_{0} is not very significant (or the results would not be as good as they are).

The degree to which a one-to-one relationship materialized was assessed visually by examining plots of *τ*^{M}_{c} versus *r* for various values of *α,* and from this process we chose a value of *α* equal to 0.25. The final form for *r* is, therefore,

In this expression, it is important to realize that *D* has been corrected for IR loss (Dutton et al. 2001; Long et al. 2003).

That the exponent should be equal to a number close to ¼ has some theoretical justification. To demonstrate this justification requires several steps. First, we use the delta-Eddington approximation (Joseph et al. 1976) to find an expression for *D.* For simplicity, we assume that the cloud is the sole influence on radiative transfer in the atmosphere, and we further assume that the cloud does not absorb. The single-scattering albedo is, therefore, zero (a good approximation for wavelengths less than 1 *μ*m; see Hu and Stamnes 1993). Invoking these assumptions and applying the delta-Eddington approximation, the expression for the spectral diffuse component of the shortwave flux at the surface is

where *D*_{λ} is the spectral diffuse irradiance at the wavelength *λ, **g* is the scattering asymmetry parameter, *A* is the surface albedo, and the spectral solar irradiance at the top of the atmosphere is *πF*_{λ}. For large cloud optical thickness, *τ*_{λ} ≳ 10, (3) is approximated well for all *μ*_{0} (and plausible *g* values) by the expression

As mentioned above *τ*_{λ} is approximately constant over solar wavelengths, so that *τ*_{λ} = *τ*_{c} ≈ constant, and we can, therefore, estimate the broadband diffuse irradiance at the surface as

where *πF*_{0} (≈10^{3} W m^{−2}) is the broadband solar forcing at the top of the atmosphere between wavelengths where cloud absorption is small (0.25 *μ*m < *λ* ≲ 1 *μ*m); the exact value of this quantity is immaterial for the simple argument presented here. In (5), the solar zenith angle dependence is *μ*_{0} (2 + 3*μ*_{0}) and it is approximated very well by 5*μ*^{1.51}_{0}, where the exponent 1.51 is objectively determined by a nonlinear least squares fit. Therefore, (5) becomes

We now use (6) in the ratio *r* [=*D*/(*C**μ*^{α}_{0})], and, furthermore, we substitute for *C* the expression developed by Long and Ackerman (2000) to describe the clear-sky total broadband flux *C* = *f**μ*^{β}_{0}. In this formula *f* and *β* are fitting coefficients, objectively determined, which vary from one clear day to the next. Typical values of *f* are 1100 W m^{−2}, and for *β,* 1.25 (e.g., see Table 4 in Long and Ackerman 2000). With these substitutions, the ratio *D*/(*C**μ*^{α}_{0}) is

Clearly, in this simple analysis, the relationship between *r* and *τ*_{c} will be independent of *μ*_{0} if we chose *α* to be equal to 0.26. Although the exact value determined for *α* depends on the value chosen for *β,* all plausible *β* values yield *α* values that are close to 0.25, thus, supporting the value determined for *α* = (¼) from our visual assessment mentioned above. Using this value, *τ*^{M}_{c} plotted versus *r* should show very little dependence on *μ*_{0}; that is, *τ*^{M}_{c} plotted versus *r* should collapse on a single curve.

The near collapse of *τ*^{M}_{c}, which is illustrated in Fig. 1, shows log_{e}(*τ*^{M}_{c}) plotted versus the value of *r* defined in (2). The spectral data required to calculate *τ*^{M}_{c} was obtained from the 415-nm channel of the Multi-Filter Rotating Shadowband Radiometer (MFRSR; Harrison et al. 1994) while *r* was found using collocated measurements of broadband irradiances. The data used to construct Fig. 1 came from a geographically diverse set of sites in the ARM system: two sites, E13 and E22, located at ARM's southern Great Plains (SGP) facility; one site from the North Slope of Alaska (NSA) facility at Barrow; and one site at Manus Island, Papua New Guinea, which is part of ARM's tropical western Pacific (TWP) observational network. The number of points used in this analysis, representing 5-min averages of *r* and *τ*^{M}_{c} over the year 2000, ranges from about 6700 at Barrow to about 9100 for the two SGP sites. In Fig. 1, the cloud optical thicknesses from each site have been color coded so that the blue, red, and green symbols represent *τ*^{M}_{c} from the NSA, TWP, and SGP sites, respectively. To remove the effect of broken cloudiness that most acutely violate the plane-parallel homogenous assumption, we only considered cases in which the cloud fraction, as determined from the Long algorithm (Long 2001), is greater than 0.99.

The near collapse evident in Fig. 1 is fortuitously consistent from one geographical site to another, suggesting that the relationship between cloud optical thickness and *r* is universal. An empirical relationship between these quantities—the function *f* in (1)—is found by fitting a curve to the points shown in Fig. 1. The solid black line in this figure shows such a fit, in which the function *f* is assumed to have the form of a hyperbolic arctangent

where *a, **b,* and *c* are the coefficients of the fit, determined by using a least squares routine to minimize the difference between *τ*^{M}_{c} and *τ*^{e}_{c}.

## Sensitivity of the results to fluctuations in cloud properties, surface albedo, and atmospheric constituents

The fit shown in Fig. 1 is an “initial” fit, because we subsequently found that the fit should be modified slightly, based on computer simulations designed to test the sensitivity of the fit to variations in surface albedo and a myriad of other factors that influence atmospheric transmission, and, therefore, retrievals of *τ*_{c}. These other factors include the following:

cosine of the solar zenith angle

*μ*_{0},aerosol optical thickness

*τ*_{aer},abundances of water vapor and ozone,

surface pressure,

cloud particle effective radius

*r*_{e},fractions of ice and water in the cloud, and

inhomogeneities in cloud cover and inhomogeneities within the cloud.

Of these factors, we cannot consider factor 6 because there is no way of knowing the mix of water and ice phase at most of the sites where our method may be applied, and many cloud optical thickness algorithms do not handle mixed phase clouds. Cloud inhomogeneities, factor 7, cannot be considered for similar reasons; that is, most surface measurement sites do not possess the means to estimate cloud inhomogeneities, and, to the best of our knowledge, nearly all of the transmission-based cloud optical depth algorithms invoke the standard plane-parallel, homogenous cloud approximation. [The effect of cloud inhomogeneities on transmission-based and satellite cloud optical depth retrievals is an active area of research; see, e.g., Boers et al. (2000) and Várnai and Marshak (2001).] For factors 1–5, our method assumes that the first of these, *μ*_{0}, can be calculated accurately, but assumes no knowledge of factors 2–5. Thus, we must also investigate to what extent fluctuations of these unknown quantities influence our method.

To examine the effect of these fluctuations, as well as variations in surface albedo on our method, we designed a set of sensitivity studies based on the SBDART model. These sensitivity studies followed a strategy of, for example, running the model for a given surface albedo for clear- and cloudy-sky cases. These cases encompassed a range of *μ*_{0} extending from 0.15 to 0.95. For the clear-sky cases, we assumed an aerosol optical thickness of 0.10 at 550 nm, and these runs provided the variable *C* in (2). The cloudy cases were identical to the clear-sky cases, except a single layer of clouds with *τ*_{c} ranging from 1 to 500 was inserted into the simulation. From these cloudy simulations we obtained *D.* These two sets of simulations allowed the generation of plots similar to Fig. 1. For a “standard” case with a surface pressure equal to 1013 hPa, a broadband surface albedo of 0.15, and an effective radius of 8 *μ*m, Fig. 2 shows log_{e}(*τ*_{c}) plotted versus the variable *r,* as well as a curve fit to the simulated points, [*r,* log_{e}(*τ*_{c})]. Also shown in Fig. 2, by the dashed line, is the initial fit to the Min algorithm–derived optical thicknesses, as shown in Fig. 1. The two fits are nearly identical for *r* in the interval [0.0, 0.5] but for *r* > 0.5, the two curves show some disagreement, which is most pronounced for *r* close to one. This disagreement is consistent with the tendency, mentioned above, for the Min algorithm to produce slightly larger optical thicknesses than a broadband algorithm for small optical thicknesses. The agreement between the two curves is sufficiently good, however, to suggest that the computer simulations are a useful means of testing the sensitivity of fit to the variations mentioned above.

We first examine the extent to which variations of water vapor contribute to errors in our method. In parallel, we also examine the advantages of using the clear-sky flux *C* to normalize the diffuse irradiance in the ratio *r* in (2). Figure 3a shows two equations, with form given by (8), fit to simulated data obtained from the SBDART runs. The line labeled with small triangles is the fit obtained by assuming a vapor content of 5 cm (an amount of vapor typical of ARM's tropical Manus Island site, or ARM's SGP site in the summer). For a given *μ*_{0} and *τ*_{c}, the ratio *r* used for this fit is *D*_{5cm}/(*C*_{5cm}*μ*^{1/4}_{0}), where the subscripts on *D* and *C* indicate that these values are appropriate for 5 cm of vapor. Next, we assume that instead of 5 cm of vapor, the actual vapor in the atmosphere is 1 cm (typical of ARM's Barrow site in summer, or SGP in the winter). We further assume that in the calculation of *r, **C* is erroneously set equal to *C*_{5cm}, and that *τ*_{c} remains unchanged. Because of the increased atmospheric transmission due to decreased vapor, the diffuse flux at the surface of the “1 cm” case must increase over that of the “5 cm” case, causing a corresponding increase in the value of *r* to *D*_{1cm}/(*C*_{5cm}*μ*^{1/4}_{0}); note that *D*_{1cm}/(*C*_{5cm}*μ*^{1/4}_{0}) is greater than *D*_{5cm}/(*C*_{5cm}*μ*^{1/4}_{0}), because *D*_{1cm} is greater than *D*_{5cm}. Thus, when the vapor changes but *C* remains fixed at *C*_{5cm}, *τ*_{c} can no longer be associated with a single value of *r.* The potential error of leaving *C* fixed is illustrated by an additional curve in Fig. 3a indicated by the line labeled by circles. This line shows the fit that would be realized if one assumed a vapor abundance of 5 cm when in fact it was 1 cm. Clearly, both curves are different, and the difference between the curves shows the error caused by assuming that *C* remains fixed. This error is about 16% for *τ*_{c} equal to 10.

However, in our scheme, *C* does not remain fixed! Instead, it varies in response to fluctuations in factors 1–4 listed above. If we redo the analysis above, but instead calculate *r* using a *C* that is appropriate for 1 cm of vapor, we find that the two curves shown in Fig. 3a become nearly identical and the error becomes very small. For example, for a given value of *r,* the range in retrieved *τ*_{c} values is only 0.5 for vapor changes from 1 to 5 cm, for *τ*_{c} around 20. In view of the very large magnitude of the vapor fluctuation considered for this example, this error is small.

A larger error can occur in practice, however, because we cannot determine *C* exactly for each cloudy day. This idealization is not achievable because the Long algorithm (Long and Ackerman 2000; Long 2001), from which we obtain *C,* finds this quantity only for clear, or nearly clear days. For other days, an interpolation scheme is employed to determine *C.* Therefore, for a cloudy day of interest, although we would expect that interpolated *C* is a reasonably accurate representation of *C* for that day, we cannot be certain that this value is perfectly representative of the conditions of that day. Thus, some error may ensue.

To evaluate the magnitude of this error, Fig. 3b shows the different fits that would occur if the assumed water vapor is 2 cm, when in fact it is either 1 or 3 cm. This 1-cm fluctuation is an estimate of the day-to-day variation of vapor at ARM's SGP site determined by looking at a plot of the daily averaged vapor abundances for 2000. (Occasionally, the variation is much greater.) For fluctuations of this magnitude, the error in *τ*_{c} retrievals is relatively small; for example, for *τ*_{c} equal to 10 the potential error is about 5%, or in absolute terms 0.5. The percent error decreases dramatically as *τ*_{c} increases. In regard to the error caused by vapor fluctuations, our conclusions are similar to those stated by Leontyeva and Stamnes (1994). When studying the feasibility of calculating cloud optical thickness from surface broadband irradiance measurements, they concluded that the use of monthly mean climatological values of precipitable vapor to retrieve *τ*_{c} is “quite acceptable in view of other uncertainties.”

The uncertainties in *τ*_{c} retrievals using our simple method, brought about by variations in water vapor, aerosol optical thickness (*τ*_{aer}), and ozone, are shown in Table 1 for three values of *τ*_{c}. The first column in this table indicates the variable that is being changed, all of the others remain fixed. For example, if we derive our empirical equation assuming an amount of vapor equal to 2 cm, and we apply this equation to a case when the actual vapor is 2 cm, the retrieved *τ*_{c} is equal *τ*^{′}_{c}, the diffuse irradiance is equal to *D*′, and the clear-sky irradiance is *C*_{2cm}. Next, we apply the same equation to a case where *D*′ does not change and *C* is assumed to be *C*_{2cm}, but the actual vapor is only 1 cm. Because *D* and *C* are fixed, the value of *r* is unchanged, and our equation would lead us to believe that *τ*_{c} is still equal to *τ*^{′}_{c}. But, to keep the *D* constant in the face of reduced vapor absorption, the actual *τ*_{c} must be larger than *τ*^{′}_{c}. The difference between the actual *τ*_{c} and *τ*^{′}_{c} is the error of the retrieval listed in the Table 1. (Table 1 does not show errors attributable to surface pressure fluctuations because these errors are negligible.)

Table 1 also lists errors associated with effective radius variations. To examine the sensitivity of the retrieved *τ*_{c} to *r*_{e}, we generated curves, similar to the curve depicted in Fig. 2, for effective radii in a range of radii from 6 to 14 *μ*m, in increments of 4 *μ*m. These curves are shown in Fig. 4, in which we see that changes in the effective radius result in slightly different fits. Using this figure as a guide, we can estimate the error in *τ*_{c} if we use an equation appropriate for one value of the effective radius, when in fact the actual effective radius is something else. The magnitude of these errors is relatively small and depends upon the value of the optical thickness. For example, suppose that we apply an empirical equation, appropriate for clouds with effective radii of 10 *μ*m, to clouds with an *r*_{e} of 6 *μ*m and an actual optical thickness of 15. The cloud optical thickness derived from the “10 *μ*m” equation is 16, and it differs from the actual optical thickness by 1, or about 7%. Repeating this scheme for an actual *r*_{e} is 14 *μ*m, we find that the error is about 3%. Errors of this magnitude—roughly 5%—are typical for other combinations of *r*_{e} and cloud optical thickness, as long as the cloud optical thickness is greater than about 8. For very small optical thicknesses, the absolute error in *τ*_{c} remains about 1, but the relative error can be much larger than the 5% stated above.

Next, we examine how surface albedo affects the equation fits. Figure 5 shows three curves derived for different values of the broadband surface albedo: 0.0, 0.15, 0.30; clearly, the form of the equation is more sensitive to surface albedo than the other factors considered above. This sensitivity is sufficiently strong to motivate a modification of (8). Referring to Fig. 5, we note that the albedo is additive in the range of albedos considered here; that is, as the albedo increases, the coefficient “*a*” in (8) increases almost in direct proportion to the increase in albedo. This property suggests a simple refinement to the equation to account for surface albedo by adding the albedo to the equation and recomputing the coefficients. The final, refined equation, including specific values of the fitting coefficients, is

where the surface albedo is limited to a range from 0.0 to 0.3.

The surface albedos used in these analyses were derived from broadband radiometers facing downward, if available. From these data we found that, in the absence of snow cover, the SGP sites have a broadband albedo of about 0.2, while the Barrow site's albedo is about 0.16 in the summer. For Manus Island, the effective albedo is influenced by water, which surrounds the island, and the island surface itself. We do not know how to assign an effective albedo to this site, but we can assume it is roughly bracketed to the albedo for the ocean (0.08) and the albedo for a tropical forest (0.12). [These albedo values are taken from Table C-7 in Stull (1991)]. With these limits in mind we assumed that the albedo at this site is 0.10.

Last, we ask whether variations in *μ*_{0} affect our results. Recall that the definition of the variable *r* in (2) was designed to minimize the effect of variations of *μ*_{0} on the form of the Eq. (8), so that a single equation would suffice regardless of the value of *μ*_{0}. To determine to what extent this goal has been realized, we plot in Fig. 6 the average difference between *τ*^{M}_{c} and *τ*^{e}_{c} versus *μ*_{0}. This plot reveals that, on the average, the error between these quantities exhibits a slight bias with respect to *μ*_{0}, and this bias is largest, about 2, for large *μ*_{0}. In light of the approximate and simple nature of our method, and the magnitude of the absolute errors shown in Table 1, we consider this bias acceptable.

## Assessment of the method's accuracy

A user of (9) needs to know how well this equation can represent *τ*_{c} derived from the Min algorithm. We, therefore, ask, what is the discrepancy between *τ*^{M}_{c} and *τ*^{e}_{c}? This discrepancy can be examined in several ways. First, we visually and quantitatively examine the point-by-point error, defined as the difference between *τ*^{M}_{c} and *τ*^{e}_{c} at a specific time. The mean deviation of *τ*^{M}_{c} − *τ*^{e}_{c}, computed over a long time series, then provides an indication of the point-by-point error. (We define the mean deviation as the average of the |*τ*^{M}_{c} − *τ*^{e}_{c}| over all points.) It should be noted, however, that this error is influenced, and usually enhanced, by sampling and timing differences between the MFRSR and the broadband radiometers. Long (1996) has shown that these differences can have a significant effect on error statistics, even for 5- and 15-min averages of the data.

Figure 7 shows a visual assessment of the point-by-point error by plotting time series of both *τ*^{M}_{c} and *τ*^{e}_{c} for a single day at the SGP central facility. The level of agreement seen on this day is typical for other days as well. The difference in the two time series provides a quick, visual portrayal of the point-by-point error. Clearly, optical thicknesess obtained from the empirical algorithm are highly correlated with the Min algorithm–derived optical thicknesses. The degree to which point-by-point agreement is achieved between *τ*^{M}_{c} and *τ*^{e}_{c} is summarized in Table 2 by tabulating the mean deviation for the sites used in the fitting process leading to (9). This statistic is also shown for three additional SGP site facilities, E2, E16, and E18, that were *not* part of the fitting process. On average, the mean deviation is about 7.

The point-by-point error can be evaluated by examining the distribution of the error, *τ*^{M}_{c} − *τ*^{e}_{c}. For facility E16, this distribution is shown by the solid line in Fig. 8. The distribution is Gaussian-like, with a mean of 1.8 and standard deviation of 6.5. The contribution to the distribution's tails comes from errors associated with very large optical thicknesses (*τ*^{M}_{c}, *τ*^{e}_{c} ≳ 150). If the distribution is restricted to optical thicknesses less than 40, then the distribution is much narrower and looks more Gaussian; the new, restricted distribution is indicated by the dashed line in Fig. 8. The standard deviation of this distribution is 4.2.

Because long-term statistics of cloud properties are more relevant to climate studies than point-by-point comparisons at an instant of time, and because of the influence of sampling issues on point-by-point statistics mentioned above, we compared distributions of *τ*^{M}_{c} and *τ*^{e}_{c} over 2000 for all sites listed in Table 2. Figure 9 shows such distributions obtained from our technique and the Min algorithm for facility E16 for 2000. These distributions represent the segregation of optical thickness into bins that are one unit wide, and the ordinate of the plot indicates the number of counts per bin. Distributions of *τ*^{M}_{c} and *τ*^{e}_{c} are plotted by the solid and dashed curves, respectively. This figure reveals that the empirical method is able to closely replicate the distribution of *τ*^{M}_{c}, including the long tail that extends up to an optical thickness of 200.

The agreement between the distributions of *τ*^{M}_{c} and *τ*^{e}_{c} can be roughly assessed by comparing medians of the distributions. These medians are also listed in Table 2, as well as the absolute difference, and percentage difference between the medians. For a given site, the medians differ by less than 10% of the *τ*^{M}_{c} value, and the absolute difference [median(*τ*^{M}_{c}) − median(*τ*^{e}_{c})] is less than about 2. The method appears to underpredict the median optical depths at the SGP site by about 5%, but more data are needed to determine if this bias persists over the site.

## Conclusions

Many algorithms have been developed to retrieve cloud optical thickness from transmission measurements of solar radiation. Algorithms that use broadband irradiances require several ancillary measurements, and it has been our experience that these algorithms require large amounts of computer time. Because of such considerations, these algorithms may be impractical for routine monitoring of cloud optical thickness. To skirt this limitation, we have developed a simple empirical expression that gives cloud optical thickness as a function of *μ*_{0}, the surface albedo, the broadband diffuse irradiance, and the broadband “clear sky” total irradiance. The last quantity is easily obtained from the Long (2001) algorithm and can be quickly computed on a desktop PC.

Cloud optical thicknesses were calculated using the Min algorithm and empirical methods for several geographically diverse sites for data in the ARM Program for 2000. For these sites, the medians of the distributions of cloud optical thickness, obtained from either the Min algorithm or the empirical method described here, are well within 10% of one another, and the shapes of the distributions are very similar.

When using this method, however, one must be cognizant of the following limitations. First, we do not know how well this method works for high-albedo surfaces, such as snow. We, therefore, recommend at present that this method only be applied for surface albedos less than 0.30. If the albedo is not known, we recommend assuming a value of 0.15. Second, the method is only applicable to fully overcast skies with cloud fractions greater than 0.99. The cloud fraction can be determined using the algorithm of Long (2001). Third, this method is not designed for low solar elevation angles, *μ*_{0} < 0.15. Fourth, we expect that substituting model-derived clear-sky total irradiances for the same quantity computed from the Long (2001) algorithm is likely to lead to inferior results because 1) models have known difficulties accurately predicting the diffuse component of the total irradiance (i.e., see Halthore and Schwartz 2000) and 2) model-derived clear-sky irradiances will not led to the cancellation of calibration factors mentioned in section 2. Therefore, we recommend that potential users of this method employ the Long algorithm and avoid the use of models to find *C.* Last, the level of accuracy exhibited by our simple method may not be adequate for all applications. For situations in which more accuracy is needed, one might consider using a complete, physically based algorithm, if one has the necessary input data to run the algorithm.

## Acknowledgments

This research was sponsored by the U.S. Department of Energy's Atmospheric Radiation Measurement Program (ARM) under Contract DE-AC06-76RLO 1830 at Pacific Northwest National Laboratory. The Pacific Northwest National Laboratory is operated for the U.S. Department of Energy by Battelle Memorial Institute. The authors thank Dr. Jim Mather for his comments regarding this work.

## REFERENCES

## Footnotes

*Corresponding author address:* Dr. James C. Barnard, Pacific Northwest National Laboratory, P.O. Box 999, Richland, WA 99352. james.barnard@pnl.gov