The primary purpose of this study is to assess the performance of 1D solar radiative transfer codes that are used currently both for research and in weather and climate models. Emphasis is on interpretation and handling of unresolved clouds. Answers are sought to the following questions: (i) How well do 1D solar codes interpret and handle columns of information pertaining to partly cloudy atmospheres? (ii) Regardless of the adequacy of their assumptions about unresolved clouds, do 1D solar codes perform as intended?
One clear-sky and two plane-parallel, homogeneous (PPH) overcast cloud cases serve to elucidate 1D model differences due to varying treatments of gaseous transmittances, cloud optical properties, and basic radiative transfer. The remaining four cases involve 3D distributions of cloud water and water vapor as simulated by cloud-resolving models. Results for 25 1D codes, which included two line-by-line (LBL) models (clear and overcast only) and four 3D Monte Carlo (MC) photon transport algorithms, were submitted by 22 groups. Benchmark, domain-averaged irradiance profiles were computed by the MC codes. For the clear and overcast cases, all MC estimates of top-of-atmosphere albedo, atmospheric absorptance, and surface absorptance agree with one of the LBL codes to within ±2%. Most 1D codes underestimate atmospheric absorptance by typically 15–25 W m–2 at overhead sun for the standard tropical atmosphere regardless of clouds.
Depending on assumptions about unresolved clouds, the 1D codes were partitioned into four genres: (i) horizontal variability, (ii) exact overlap of PPH clouds, (iii) maximum/random overlap of PPH clouds, and (iv) random overlap of PPH clouds. A single MC code was used to establish conditional benchmarks applicable to each genre, and all MC codes were used to establish the full 3D benchmarks. There is a tendency for 1D codes to cluster near their respective conditional benchmarks, though intragenre variances typically exceed those for the clear and overcast cases. The majority of 1D codes fall into the extreme category of maximum/random overlap of PPH clouds and thus generally disagree with full 3D benchmark values. Given the fairly limited scope of these tests and the inability of any one code to perform extremely well for all cases begs the question that a paradigm shift is due for modeling 1D solar fluxes for cloudy atmospheres.
1. Introduction and objectives
Much of the rich structure of the earth's climate can be traced to interactions between radiation and the four-dimensional distribution of the three phases of water. While representation in large-scale atmospheric models (LSAMs) of all such interactions deserves attention and work, those involving clouds are recognized generally as notoriously difficult and in need of much work (e.g., Houghton et al. 1995). This recognition stems from the obvious limitations of modeling unresolved clouds and radiative transfer in LSAMs along with numerical experiments that show that seemingly small systematic changes in cloud properties have significant impacts on simulated regional and global climate (e.g., Senior 1999). Moreover, studies that intercompared cloud radiative feedbacks in LSAMs (Cess et al. 1996, 1997) came to the conclusion that different representations of cloud-related processes in LSAMs may account for much of the uncertainty associated with estimates of climate sensitivity.
Preceding and overlapping the Cess et al. studies, and other LSAM intercomparisons (Lambert and Boer 2001), was the Intercomparison of Radiation Codes in Climate Models (ICRCCM) program (World Climate Research Programme 1984; Ellingson and Fouquart 1991; Fouquart et al. 1991). Fouquart et al. (1991) demonstrated that when several 1D atmospheric solar radiative transfer codes operated on the same atmospheric profile, the range of estimated irradiances often exceeded 20% with root-mean-square differences relative to medians of typically 4%–10%; single-layer homogeneous overcast clouds and aerosols showed the largest disparities.
ICRCCM's main message to climate modeling groups was that their solar codes may be imparting fictitious radiative forcings that, in all likelihood, have adverse effects on simulated climate. A logical question leading from these results is as follows: since different radiation codes operating on identical (and simple) atmospheric columns yield significantly different irradiance profiles, how much of the disparity in estimates of climate sensitivity is due to different treatments of radiative transfer? To this day, this question remains unexplored.
Since ICRCCM, most codes have been modified to better resolve the solar spectrum, spectroscopic databases have been updated [high-resolution transmission molecular absorption database (HITRAN) upgrades], and new gaseous transmittance and cloud optical property parametrizations have become active (e.g., Hu and Stamnes 1993; Fu and Liou 1993; Edwards and Slingo 1996). Yet, most codes are built on the extreme, idealistic, systematic assumption of maximum/random overlap of homogeneous clouds (Geleyn and Hollingsworth 1979). Over the past two decades, however, several studies have demonstrated that different overlap configurations have major impacts on simulated radiation budgets (e.g., Morcrette and Fouquart 1986; Stubenrauch et al. 1997; Barker et al. 1999). Likewise, when horizontal fluctuations in cloud are neglected, estimates of domain-averaged irradiances can be in serious error (e.g., Cahalan et al. 1994a; Barker et al. 1996, 1999; Pincus et al. 1999). There are far fewer, and inconclusive, assessments of ensuing dynamical impacts (Tiedtke 1996; Liang and Wang 1997; Morcrette 1993). This is because multilayer 1D solar codes that handle cloud horizontal and vertical variability have only begun to be tested (e.g., Barker and Fu 2000), and also because the additional input they require is unavailable as yet. Clearly, there is a need to develop both parameterizations for characterizing unresolved clouds and 1D radiative transfer algorithms that utilize them.
Given developments in modeling radiative transfer since ICRCCM and given no formal indication of how 1D solar codes respond to the general scenario of multiple layers of nonovercast clouds, it seems timely to initiate a follow-on to ICRCCM. Thus, the primary objective of this study is to assess how well operational, and experimental, 1D solar radiative transfer codes predict broadband irradiances and heating rate profiles when they operate on columns of information pertaining to partially cloudy atmospheres generated by 3D cloud-resolving models.
Before this objective can be realized, however, it is essential that clear-sky (i.e., cloud and aerosol free) and plane-parallel, homogeneous (PPH) overcast clouds be considered first. This is because variances and biases for simple cases carry over to complicated cloud cases. As it was never the intention of this study to repeat ICRCCM, the number of clear-sky and homogeneous overcast cloud cases considered here is minimal. Hence, this study is more an extension of ICRCCM rather than a repeat.
The second section outlines some details of this intercomparison. The third and fourth sections present the model atmospheres, benchmark calculations, and intercomparison results for simple and complex atmospheres. The fifth section contains conclusions.
2. Participation and general conditions
Tables 1–4 list participants and properties of their model(s). Results were submitted for 25 1D solar codes and four 3D Monte Carlo (MC) codes. Most participants completed calculations for 12 atmospheric profiles, though results for only 7 cases are presented here. For all cases, calculations were performed for solar zenith angle cosines µ0 from 0.05 to 1.0, in increments of 0.05, with spectrally invariant (Lambertian) surface albedo αs of 0.2. In contrast, ICRCCM used µ0 of 0.26 and 0.866 with αs of 0.2 and 0.8. Each participant was asked to provide top-of-atmosphere (TOA) albedo αp, total atmospheric absorptance aatm, and surface absorptance asfc for each µ0 integrated spectrally across wavelength intervals 0.2–0.7, 0.7–5.0, and 0.2–5.0 µm, neglecting terrestrial emission. Upwelling and downwelling irradiances were also requested at each model level for the spectral range 0.2–5.0 µm at µ0 = 0.25, 0.5, and 1.0.
In all cases, cloud optical properties were modeled as pure liquid water spheres with effective radius re = 10 µm, even when it was obvious that cloud water should have been either ice or precipitation. This avoided ambiguity in specifying ice crystal structure, precipitation size distributions, and associated spectral optical properties. Another intercomparison might inform participants about the presence of clouds and have them set the optical properties. This may work with codes used in operational LSAMs, but many participating codes lacked such automation. On the other hand, had specific optical properties been assigned, numerous sets would have been required to facilitate all combinations of spectral intervals (see Table 3). Finally, while all cloud particles were assumed to be liquid spheres with re = 10 µm, individuals were free to generate corresponding spectral optical properties but requested to use the Henyey–Greenstein (1941) phase function when applicable.
3. Clear-sky and homogeneous overcast clouds
The purpose of this section is twofold. First, it is essential to demonstrate that the 3D Monte Carlo algorithms used to set benchmarks for the complex 3D cloud fields agree both among themselves and with a detailed line-by-line (LBL) model that was designated as the benchmark for the simple atmospheres. As yet, LBL models cannot perform 3D radiative transfer calculations for complicated atmospheres to benchmark standards (Partain et al. 2000). Therefore, 3D radiative transfer codes must be compared to LBLs for simple conditions. If the 3D codes perform well for simple conditions, this bolsters confidence in them as providers of benchmarks for complex cloud cases. The second point of this section is to establish ranges of 1D code results for simple conditions.
a. Description of cases
For the single clear-sky (denoted as CLEAR) and two homogeneous overcast clouds, the McClatchey et al. (1972) tropical (TRO) atmosphere was used. This atmospheric profile consists of 64 layers, extends from the surface to ∼100 km, and has prescribed values of pressure, temperature, density, water vapor, and ozone mixing ratios. Carbon dioxide mixing ratio was set to 360 ppm and it was suggested that oxygen be included too.
The first cloudy case, denoted as CLOUD A, used the TRO clear-sky atmosphere, but the layer between 3.5 and 4 km was filled with a uniform overcast cloud with qc = 0.159 g kg–1. For pure liquid spherical droplets with re = 10 µm, visible optical depth was ∼10 (Slingo 1989). The second cloudy case, CLOUD B, was a uniform overcast cloud positioned between 10.5 and 11 km with qc = 0.034 g kg–1, which gives a visible optical depth of ∼1. Thus, CLOUDs A and B resemble closely cases 49 and 46 of Fouquart et al. (1991).
b. Benchmark calculations
It is essential for this study to demonstrate that results designated as benchmarks are indeed such. Mlawer et al. (2000) compared the LBL developed at Atmospheric Environmental Research, Inc. (AER), known as the Code for High-Resolution Accelerated Radiative Transfer (CHARTS), to high-quality spectral surface irradiance measurements made at the U.S. Department of Energy Atmospheric Radiation Measurement (ARM) Program Southern Great Plains (SGP) site (e.g., Stokes and Schwartz 1994) by the Rotating Shadowband Spectroradiometer (RSS; Harrison et al. 1999) and the Absolute Solar Transmittance Interferometer (ASTI; Murcray et al. 1996; Mlawer et al. 2000). In these comparisons, atmospheric profiles were based on coincident soundings (cf. Brown et al. 1999; Clough et al. 1989). While aerosol optical properties were unknown, values of single-scattering albedo and asymmetry parameter required by CHARTS to yield errors in spectral diffuse irradiance always less than 5% were consistent with recognized aerosol models (Kato et al. 1997).
Figure 1 shows an example of how well CHARTS performed for a small µ0 and heavy water vapor burden. While these experiments had adequate radiometric measurements available only at the surface, photons contributing to surface irradiance typically experience the longest paths and most scattering events. Thus, it is difficult to imagine being able to predict spectral surface irradiance as well as in Fig. 1 and have poor corresponding estimates of αp and aatm. It may, therefore, be safe to say that CHARTS represents the modeling standard for clear skies; identification of such a standard was not possible during ICRCCM.
For the current study, CHARTS employed optical properties for homogeneous overcast clouds based on Mie scattering calculations (Wiscombe 1979, 1980) using refractive indices from Segelstein (1981) integrated spectrally over 5-nm intervals for a gamma distribution of droplets with re = 10 µm and effective variance of 0.1.
Due to time and resource constraints, CHARTS completed only a limited set: µ0 = 0.25, 0.5, and 1.0 for the simple cases. Figure 2a shows broadband values of αp, aatm, and asfc as functions of µ0 for the CLEAR experiment from CHARTS and the four MC codes. This shows that the MC codes agree extremely well among themselves and that relative deviations from CHARTS rarely exceed ±2%. These results are corroborated in Table 5a, which lists broadband, visible, and near-IR values of αp, aatm, and asfc for µ0 = 0.5.
Establishing cloudy-sky benchmarks was more difficult. First, since three of the four MC codes have moderate spectral resolution, their cloud optical properties depend on spectral weighting (e.g., Li et al. 1997). Second, underlying assumed droplet size distributions varied slightly though all used re = 10 µm. Despite these uncertainties, Figs. 2b,c show that the MC codes and CHARTS agree almost as well as in the CLEAR experiment. For CLOUD A, estimates of αp and aatm differ slightly and stem from different droplet optical properties and spectral weightings. The mean number of photon scattering events (averaged over all µ0 and wavelengths) for CLOUD B is ∼3, while for CLOUD A it is ∼17. This explains why minor differences in cloud optical properties are apparent in results for CLOUD A.
Tables 5b,c list broadband, visible, and near-IR values for CLOUDs A and B at µ0 = 0.5 for the MC codes and CHARTS. While broadband MC results differ from CHARTS by always <4%, differences climb to ∼9% for the near-IR band. Table 6 lists heating rates in the cloudy layers of CLOUDs A and B for the MC codes and CHARTS. Overall, these results suggest that the MC codes can be used to assess the 1D codes.
c. Contrast with ICRCCM
Figure 3 shows root-mean-square (rms) differences relative to the median for all 1D results, as well as the range of broadband values of αp, aatm, and asfc. Available corresponding values for aatm and asfc from ICRCCM are also indicated. For the CLEAR case, both rms and range are remarkably similar to ICRCCM values. Note, however, that most of the codes documented here were low spectral resolution while in ICRCCM about half were classified as high resolution.
For CLOUD A, estimates of aatm are much more consistent than they were in ICRCCM. This may be due to the current widespread use of a small number of cloud optical property parameterizations (e.g., Slingo 1989; see Table 2), though it is interesting to note that participants here were only told to use re = 10 µm, while in ICRCCM they were given exact size distributions. The cohesion of αp seen for CLOUD A degrades for CLOUD B and the corresponding range increases by a factor of ∼2. It is not clear why this is the case as both CLOUDs A and B were both specified to be liquid phase with re = 10 µm. Save for aatm for CLOUD A, rms's and ranges for aatm and asfc for CLOUDs A and B are up slightly from ICRCCM levels.
d. Results: Clear sky
Figures 4a–c show median values and interquartile ranges (the range containing 25% of all results either side of the median) for αp, aatm, and asfc as functions of µ0 for the 1D models in the CLEAR experiment. It also shows results from the Pennsylvania State University (PSU) MC as it could be partitioned at ∼0.7 µm. As Table 5 shows, the PSU model represents the other MCs well. Broadband αp for the 1D codes agree very well with the benchmarks though they show a tendency to underestimate aatm and overestimate asfc. Compared to the models of 10–15 years ago, most 1D codes in use today predict slightly more aatm and slightly less asfc and αp [see Table 9 in Fouquart et al. (1991)]. The overestimates of asfc shown here are, however, similar in magnitude to those shown by Kinne et al. (1998). The particularly good agreement in the visible portion of the spectrum probably stems from better parameterizations of Rayleigh scattering [cf. Fouquart et al.'s (1991) Fig. 4]. For the near IR, however, there is a marked underestimation of aatm in step with an equally severe overestimation of asfc. This makes sense for when surface albedos are low (as in this case); transmitted photons tend to have longer optical pathlengths than reflected photons. Hence, weak atmospheric absorptance is expected to show more in the transmitted field than the reflected field.
The difference between benchmark and 1D median values of broadband aatm can be described by
which has a (spherical) mean value of
This implies that, when integrated over the sun-up period, a typical 1D code in use today underestimates aatm by ∼11 W m–2 for clear-sky tropical conditions. Since the linear relation Δaatm = bµ0 is so well defined, this suggests that the magnitude of b will apply for somewhat smaller water vapor burdens too. Thus, the underestimations shown here tend to echo those reported by Arking (1999) and are almost certainly due to lack of a water vapor continuum (Mlawer et al. 2000). Several models using gaseous transmittance parameterizations based on old spectroscopic databases and no water vapor continuum follow closely Δaatm ≈ –0.02µ0, which translates into underestimates of diurnal mean atmospheric absorption of about 18 W m–2. For overhead sun, some models underestimate atmospheric absorption by 40 W m–2.
Figure 4d shows the three quartiles for 1D heating rates at µ0 = 0.5 as well as the mean and standard deviation for the benchmark codes. Below ∼350 mb, benchmark values are systematically larger than the vast majority of 1D estimates. This points to water vapor lines and continuum not included, or represented poorly, in old spectroscopic databases upon which many 1D models are based. The 3D codes used here employ parameterizations based on newer databases that include the continuum.
e. Results: Planar, overcast clouds
Figures 5 and 6 show quartiles of αp, aatm, and asfc as functions of µ0 for the 1D models for CLOUDs A and B. For both clouds, αp predicted by most 1D codes are greater than the benchmarks for high sun and less for low sun regardless of spectral range. The overestimation at large µ0 is particularly evident in the near IR and is likely tied to too little aatm, which carries over from CLEAR. It may also be associated with widespread use of Slingo's (1989) parameterized asymmetry parameter, which, as discussed later, is too small across most of the solar spectrum. Also, about 1–2 W m–2 of this underestimation likely stems from extensive use of two-stream approximations (Räisänen 2002). Likewise, underestimation of αp at small µ0 resembles errors inherent in delta-two-stream approximations (King and Harshvardhan 1986).
Regarding aatm, the 1D codes do fairly well in the visible range, but in the near IR and integrated over the entire spectrum, they underestimate systematically relative to both the 3D codes and CHARTS by amounts similar to those for CLEAR, the exception being CLOUD B at small µ0 where cloud transmittances (i.e., asfc) are too large due again to the predominance of delta two streams.
Table 6 lists heating rates for the cloud layer for all 3D models and CHARTS, as well as the quartile ranges for the 1D codes. Somewhat contrary to aatm, as shown in Figs. 5 and 6, heating rate estimates for clouds by 1D models are quite good. Clear relations could not be identified between spectral resolution, optical property parameterization, and cloud heating rate. It is possible that empirical alterations to droplet single-scattering albedos in low-resolution codes have countered excessive absorption that is known to plague these models (Slingo 1989).
A comment on cloud optical property parameterization is warranted at this stage. Figure 7 shows high-resolution values of single-scattering albedo ω0 and asymmetry parameter g computed by Mie theory (Wiscombe 1979, 1980) using refractive indices from Segelstein (1981) for a gamma droplet size distribution with re of 10 µm and effective variance of 0.1. Also shown are estimates from Slingo's (1989) and Hu and Stamnes's (1993) parameterizations at re = 10 µm. Though differences are fairly minor for ω0, Slingo's values for g are systematically low by ∼0.02 across the entire visible and near IR, which will yield clouds that are too reflective. A difference of 0.02 may sound small, but re = 6 µm is required to bring Hu and Stamnes's estimates of g into agreement with Slingo's at re = 10 µm. A difference of 4 µm in re is larger than what is expected to have occurred in the most polluted regions over the industrial era due to increased sulphate aerosol concentrations (e.g., Lohmann and Feichter 1997).
As Table 2 shows, Slingo's parameterization is especially popular despite Hu and Stamnes's being applicable over a wider range of re. This is likely due to Hu and Stamnes's use of noninteger exponents, which degrade computational efficiency to a point that is, apparently, unacceptable for GCMs and numerical weather prediction (NWP) models. There exists, however, an unpublished third-order polynomial version of Hu and Stamnes's parameterization (Y. Hu 1996, personal communication) that is efficient when Horner's algorithm is applied. Finally, Hu and Stamnes's parameterization exhibits minor discontinuities at certain re where fitted coefficients change.
4. Realistic 3D cloud fields
This section is the apex of this study: to assess how well 1D solar codes interpret and handle unresolved nonovercast clouds. Cloud fields produced by 3D cloud-resolving models (CRMs) were used to represent complicated, realistic cloudy atmospheric columns. From the 3D CRM fields, 1D columns of cloud properties were produced to represent fields generated by a hypothetical LSAM. As the LSAM was assumed to have done a perfect job at estimating unresolved cloud, all 1D models were provided with perfect profiles of whatever they needed to operate. Most required profiles of cloud fraction, mean cloud water mixing ratio qc, and mean water vapor mixing ratio inside and outside of clouds. Some models required profiles of mean logarithmic qc and cumulative upward and downward cloud fractions. The former gives an indication of horizontal variability of cloud water, while the latter provides information on cloud overlap.
Having established that for simple atmospheres all four Monte Carlo codes agree extremely well among themselves and with CHARTS, this bolsters confidence in them being used as benchmarks for the complex cloud cases. A novelty of this study is its use of conditional benchmarks. These are discussed in detail in the following section, which is followed by descriptions of the cloud fields and then by results.
a. Benchmark calculations
It is worth noting that an infinite number of 3D cloud fields can all give rise to the same set of 1D profiles of cloud information. Thus, the full 3D benchmark irradiances produced by the MC models are, in principle, only a single sample from a potentially diverse population (Barker et al. 1999). Moreover, all 1D codes make explicit simplifying assumptions about unresolved clouds and thereby, knowingly, omit aspects of unresolved cloud structure. So to provide 1D modelers with only full 3D MC results would not be too informative. For example, a sound a priori assumption is that 1D codes using PPH clouds with maximum/random overlap and 3D MC codes will produce results that generally disagree. So, while full 3D results are informative as representatives of the truth, it would be equally informative to have benchmarks for maximum/random overlap of PPH clouds too.
It is a simple task to create benchmark irradiances for PPH clouds following maximum/random overlap (and other idealized scenarios) with an MC code, though it is likely a luxury that most 1D modelers do not have. As such, one MC code was used to establish conditional benchmarks for the 3D fields. The 1D codes were then sorted according to how they handled unresolved clouds and compared to their respective conditional benchmarks. The following series of benchmark calculations was performed (cf. Barker et al. 1999).
1) Full 3D
In this case, the actual CRM fields with inherent horizontal grid spacings Δx were used by all MC algorithms along with cyclic horizontal boundary conditions. These experiments provided benchmark estimates of broadband, domain-averaged αp, aatm, asfc, and heating rate profiles. All results were averaged around solar azimuth ϕ by injecting each photon at ϕ = 2πR, where R is a uniform random number between 0 and 1.
Four independent 3D MC algorithms were used to establish these benchmarks (see Table 4). This redundancy was primarily to convince all that the benchmark results are trustworthy and that only a single MC was needed to compute conditional benchmarks.
2) Independent column approximation
These cases used the exact 3D CRM fields except that Δx was infinite. In so doing, all information about spatial distributions of cloud was available but there was no indication that clouds had finite horizontal extent, including sides. The independent column approximation (ICA) has often been shown to be a very good approximation to the full 3D solution for both single-layer, planar clouds (Cahalan et al. 1994b) and multilayer, towering clouds (Barker et al. 1998, 1999). The largest differences between ICA and the full 3D solution are at small µ0.
3) Exact overlap
With this model the exact positions of cloudy cells, as dictated by the CRM, are retained, but Δx is infinite and all clouds in a layer have the same mean cloud mixing ratio defined as
and qc(i, j, k) is the mixing ratio in the (i, j)th cell of the kth layer. Therefore, this model represents the ideal 1D PPH model that is capable of handling a continuum of cloud overlap rates (e.g., Barker et al. 1999).
4) Maximum/random overlap
For this model, when clouds occur in adjacent layers they are maximally overlapped and any excess cloud is positioned at random across its layer. If a clear layer separates two cloudy layers, they are assumed to overlap randomly. Also, Δx is infinite. This model was advocated in the observational study of Tian and Curry (1989) and proposed initially by Geleyn and Hollingsworth (1979).
These fields are created by first randomly distributing the Nc(k) cloudy cells, all of which have qc(k), across the uppermost cloud layer. Moving down, if Nc(k – 1) < Nc(k), randomly position cloudy cells in the (k – 1)th layer beneath cloudy cells in the kth layer. If Nc(k – 1) > Nc(k), place Nc(k) cells in layer k – 1 beneath the cloudy cells in layer k and randomly distribute the remaining Nc(k – 1) – Nc(k) cloudy cells in layer k – 1 among the remaining N – Nc(k) cloudless cells in layer k – 1, where N is the total number of cells per layer. If a cloudless layer is encountered, position cloudy cells of the next cloud layer that is encountered randomly across the layer, just as was done for the uppermost cloud layer.
Limited use was made of this category's close relative: maximum overlap. With maximum overlap, cloudy cells are always placed into the minimum number of vertical columns in a domain. Total cloud fractions for maximum overlap are always less than or equal to those for maximum/random overlap.
5) Random overlap
These fields were created by simply taking the Nc(k) cloudy cells containing qc(k) and redistributing them at random across the layer regardless of cloud state in adjacent layers. Again, Δx is infinite. Though LSAMs do not apply this model globally, it is an interesting extreme case that isolates the random component of the maximum/random model, and one that a few groups implemented here for diagnostic reasons. Table 7 summarizes the models just presented.
b. Description of cases
The fields used here were simulated by several very different CRMs. When selecting these fields, the intention was to capture a range of demanding tests, not to represent clouds with the largest globally averaged radiative impacts.
In all cases, water vapor was averaged horizontally within a layer into clear- and cloudy-sky portions. Both 1D and 3D codes used these fields. Though this averaging can have significant impacts on local heating rates, it had very little impact on domain averages. Above model domains up to ∼100 km, skies were cloudless, temperature profiles followed the TRO atmosphere, and water vapor was set almost to zero, which tests the ability of codes to represent sharp vertical moisture gradients such as off the west coasts of midlatitude continents. Ozone followed the TRO profile through the entire atmosphere.
1) Marine stratocumulus
This marine boundary layer cloud was provided by B. Stevens (1998, personal communication) and was generated with conditions during the Atlantic Trade Wind Experiment (ATEX). The cloud extends from about 0.7 to 1.6 km above the surface with the top of the model domain at ∼3.2 km. Vertical grid spacing ranges from 0.02 to 0.04 km, domain size is (6.8 km)2, and Δx = 0.1 km. At the time, it was necessary to consider such a small domain in order to address small-scale variability.
Figure 8a shows the field and summarizes key radiative properties. At about 1.5 km, where most of the cloud sits in a layer about 150 m thick, mean visible extinction coefficient 〈β〉 is ∼100 km–1. Cloud lower down is likely precipitation, as fractional amounts Ac are typically less than 0.05 and 〈β〉 are huge at ∼200 km–1. Also shown is
where σβ is standard deviation of β across a layer. The more inhomogeneous the cloud, the smaller νβ. Through use of a single value for re, νβ is identical to νqc. For this case νβ ≈ 1, though for the entire cloud ντ, for visible optical depth τ, is only 0.45. Had the clouds in this case overlapped maximally, ντ would have been closer to 1 but as Fig. 8a shows, overlap falls between maximum and random; total cloud amount is 0.565 with the cloudiest layer Ac ≈ 0.3. Also, had maximum overlap occurred, mean visible optical depth 〈τ〉, which is 14.4, would have been larger.
2) Open cells
This field represents a cold-air outbreak over warm water (Anderson et al. 1997). As Fig. 8b shows, this field consists of vigorous open cellular convection with clouds from about 0.8 to 7.5 km, which is the top of the model domain. Vertical grid spacing is 0.15 km, domain size is (50 km)2, and Δx = 0.39 km. This domain size is similar to those in numerical weather prediction models. The amount of water vapor in this case is roughly one third that of the others. Though much of the cloud in this simulation should be mixed phase, cloud water was treated as liquid droplets with re = 10 µm. In the main body of the cloud between 4 and 7 km, 〈β〉 is ∼30 km–1. Below the main body 〈β〉 increases to about 70 km–1 but Ac are small. As a result, 〈τ〉 for the entire cloud is 52.9. Again, due to cloud overlap being between random and maximum, ντ is only 0.53 while νβ > 1 for the cloudiest layers. Moreover, total cloud fraction is 0.9 though the cloudiest layer has Ac = 0.7.
3) Deep tropical convection
Figures 8c,d summarize two fields extracted from Grabowski et al.'s (1998) simulation of phase III of the Global Atmospheric Research Programme Atlantic Tropical Experiment (GATE). GATE A (Fig. 8c) consists of nonsquall clusters of organized convection with the anvil removed so as to mimic towering, or developing, clouds. Only liquid phase hydrometeors were considered. GATE B (Fig. 8d) is a squall line and both liquid and ice phases were considered. The ice crystals (referred to as ice A by Grabowski et al.) were, however, treated as though they were liquid spheres. The fields contain deep convective clouds reaching up to 16 km, domains are (400 km)2 and rival those in global climate models, and Δx = 2 km. There are 35 layers of unequal geometric thickness extending up to ∼20 km. Domains this large represent the coarsest GCMs. Rain, snow, and graupel were neglected with little radiative consequence.
For GATE A, 〈β〉 are typically ∼60 km–1 and νβ ≈ 1. Owing again to overlap of inhomogeneous layers, ντ is quite small at only 0.25. Thus, this field is extremely variable and represents a demanding test for 1D codes. Note that clouds overlap almost maximally above 6 km but acquire a substantial random component at lower altitudes. Total cloud fraction is 0.46 and 〈τ〉 ≈ 59.
GATE B's profiles are similar to those of GATE A at altitudes less than 4 km. GATE B is, however, dominated by near-overcast clouds above 8 km where 〈β〉 are less than 10 km–1 and νβ ≈ 0.3. These small values of νβ are due to cohabitation of anvil and dense core regions. Thus, with 〈τ〉 ≈ 85 and ντ ≈ 0.3, this case is as demanding as GATE A.
5. Results: 3D cloud fields
Figures 9–12 show αp as a function of µ0 and heating rate profiles for µ0 = 0.5 for the four 3D cloud fields. For each field, values are shown for full 3D, exact overlap, maximum/random overlap, and random overlap as documented in section 4. The 1D models were partitioned into their respective class and plotted anonymously. Benchmark values for the specific genre of cloud treatment are shown on each plot. For reference, full 3D benchmarks are shown on each plot too (means ± standard deviations for the four MC codes). In the following sections, results for each genre are discussed. By referring to Table 3, one can deduce the genre into which each model falls. The focus here is on αp and heating rate profiles.
The 1D codes included in the ICA genre are those that attempt to deal with horizontal variations in cloud beyond cloud fraction Ac. As Figs. 9–12 reveal, the ICA benchmarks are consistently closest to the full 3D values. In each case though, characteristic ICA biases for αp are apparent: too reflective at large µ0 and too transmissive at small µ0. Clearly, this class of 1D code is still very much in the experimentation phase as the variance of results from case to case is large. For example, for GATE A they agree among themselves fairly well and follow the conditional benchmark, although obvious problems exist for heating in the lower reaches of cloud. In fact, the outlier in this case is the European Centre for Medium-Range Weather Forecasts (ECMWF) model, which uses maximum/random overlap and represents horizontal variability by simply reducing cloud optical depth everywhere all the time by a factor of 0.7 (Cahalan et al. 1994a; Tiedtke 1996). Results for this model are apparent in the other fields for it is the one that resembles corresponding maximum/random benchmarks.
GATE B exemplifies the importance of addressing horizontal variability of cloud, as the αp benchmarks for exact, maximum/random, and random overlap are all similar and much larger than the full 3D values, which are tracked well by the ICA. The same can be said for heating rates. This is because GATE B is dominated by a thick, yet highly variable, near-overcast anvil that renders choice of overlap scheme irrelevant. Contrary to GATE A's results, consistency of αp and heating rate estimates by 1D ICA-like models are terrible.
The OPEN CELLS case proved to be a challenge for both the 1D codes and the ICA benchmark. For full 3D transfer at large µ0, substantial amounts of radiation leak out cloud sides in a downward direction (Welch and Wielicki 1984). It is unlikely that 1D codes will ever address this successfully. The 1D codes tend to perform well for ATEX with the ECMWF code still an outlier. Had the cloud in ATEX been represented by one or two layers at typical LSAM vertical resolution (i.e., hundreds of meters rather than tens of meters), the ECMWF code, and the other 1D codes in this class, would have performed well. The rest of the codes that do not account for horizontal fluctuations would have done worse, for they would have been blind to implicit horizontal variability inherent in the overlap of many homogeneous, partially cloudy layers.
Despite its crude, empirical treatment of overlap, the Meteorological Service of Canada (MSC1) code (Oreopoulos and Barker 1999) performs fairly well. It applies a downward reduction to cloud τ, thereby recognizing that irradiances are variable in the horizontal and that they tend to be small in dense cloud beneath other dense clouds. The MSC1 code requires an estimate of horizontal variability [i.e., Eq. (3)] and techniques that provide it are beginning to emerge (Considine et al. 1997; Jeffery and Austin 2003). Nevertheless, this code, like other 1D codes, is rigid in the sense that altering assumptions about unresolved cloud is difficult, and that to include systematic (known) relations such as horizontal variations in droplet size may be next to impossible.
b. Exact overlap
The 1D codes in this category were those that did not address horizontal variability of cloud but did utilize profiles of accumulated cloud fraction. Thus, they attempted to capture the main features of cloud overlap (i.e., a profile of cloud fraction presented to the direct beam). Figures 9–12 show that, as expected, exact overlap benchmarks of αp are systematically greater than those for full 3D and ICA on account of homogenized clouds that cover the sky to the same extent as those in the full 3D and ICA.
As the 1D codes in this class agree well for CLEAR and CLOUDs A and B, most of the variance for the 3D clouds arises from different handling of cloud overlap. For GATE B and OPEN CELLS these codes agree fairly well, for cloud was quite extensive and details about overlap are minor. For ATEX and GATE A, however, overlap is important (see Figs. 8a,c) and disagreement among 1D codes is large for both αp and heating rate.
Like the ICA-style models, 1D models that require, or at least those that can utilize, information about exact overlap are requesting information that is currently unavailable in LSAMs. For a discussion on the overlapping structure of real clouds derived from cloud-profiling radar, see Hogan and Illingworth (2000).
The large spike in heating rate at 350 mb for OPEN CELLS is due to the sudden encounter of downwelling irradiance with water vapor in the top layer of the cloud model domain. This is associated with correlated-k distribution models that use relatively few k values and were never expected to be applied to such extreme conditions. This feature was also present for ATEX, but above 800 mb, so it does not appear on the plots in Fig. 9.
c. Maximum/random overlap
All the 1D codes in this class use horizontally homogeneous clouds in conjunction with the maximum/random overlap assumption. It is the most populated class, which signifies the growing popularity of this scheme. Nevertheless, this scheme is risky for it is an extreme approximation that, when applied systematically, can be expected to underestimate αp and overestimate asfc (see Fig. 8). Figures 9–12 show that use of homogeneous clouds often fortuitously counter this bias and result in overall radiative responses that are in better agreement with the full 3D benchmarks than are the exact overlap benchmarks. It is important to note, however, that derivatives of radiative fluxes with respect to cloud properties are almost minimized for this model; surely this will impact estimates of climate sensitivity.
For the most part, 1D codes that employ maximum/random overlap track their conditional benchmark values reasonably well in magnitude, but more so in dαp/dµ0. What is interesting, however, is that all of these codes, save for those that underestimate water vapor absorptance severely, underestimate αp for GATE A and ATEX, the cases where the maximum/random overlap assumption has a large impact. Figure 13a is a replot of αp as a function of µ0 for maximum/random overlap of GATE A (i.e., Fig. 11), but it also shows the conditional benchmark for maximum overlap as defined in section 4a(4). While the benchmark maximum overlap values of αp are slightly less than those for maximum/random, they are still greater than most 1D model values.
Figure 13b shows αp for CLOUDs A and B as predicted by all 1D codes that use maximum/random overlap and corresponding values for a MC code. Note that for homogeneous clouds the 1D models tend to reflect more than the benchmark while the opposite is true for the 3D cloud cases. To answer why this occurs is beyond the scope of this report. It is an important issue, however, for if the 1D codes were to begin addressing horizontal variability, their estimates of αp would be suppressed even more than those in Figs. 9–12. Furthermore, note that the relative range of 1D results is roughly double that for the homogeneous clouds. This suggests differences, or errors, in the interpretation of maximum/random overlap by 1D codes.
There is a tendency for the 1D models to underestimate heating rates beneath layers with maximum Ac yet rally back respectably, often overestimating, near the surface. The explanation for this is that low-level clouds are shielded too much by clouds aloft. A similar problem afflicts the 1D codes in the ICA class.
d. Random overlap
Although operational 1D codes do not invoke random overlap everywhere all the time, it is invoked when cloud layers, or blocks of cloudy layers, are separated by at least one clear layer. Thus, when a 1D code reverts to random overlap for the cloud fields used here, it should agree with the conditional benchmarks. As Figs. 9–12 show, this is not always the case. Particularly puzzling are estimates of αp for GATE A. First, the curve that cuts diagonally across the plot (across the plot for ATEX too) utilizes Briegleb's (1992) scheme, which uses overcast clouds with optical depths scaled by a factor of Ac3/2. This scheme does, however, produce maximum heating near the correct altitude. The other two models underestimate benchmark values of αp by significant amounts, not withstanding the fact that they underabsorb in the CLEAR and overcast tests. It is not clear whether the 1D codes are trapping too much upwelling radiation via multiple internal reflections or whether they are not presenting enough cloud to the direct beam. Judging by the heating profiles, photons appear to be getting absorbed at the proper altitude, despite one of the codes seriously underestimating absorption by droplets.
It is unlikely that the benchmarks are incorrect or unrepresentative (cf. Barker et al. 1999). Overall, these results are a bit disconcerting for they suggest that the seemingly routine treatment of random overlap in 1D codes may often be incorrect. This may not bode well for models using maximum/random overlap.
e. Climatological context
To put these results into a climatological context, consider cloud radiative forcing (CRF) at the TOA, which is defined as
where αclrp is TOA albedo from the CLEAR experiment, αcldp is the corresponding value from one of the cloudy experiments; and Sȯ is normal TOA irradiance, which depends on latitude θ and day number J. Defining CRF3D(µ0) as CRF for any of the 3D codes, let the error for a particular 1D code be defined as
The diurnal-mean error in CRF for a 1D model is, therefore,
where Δt is time between sunrise (or sunset) and solar noon, and ΔCRF in the integrand of (4c) was fit with cubic splines for easy numerical integration.
The left column of Fig. 14 shows annual marches of CRF3D for the four 3D cloud fields. The other three columns show contour plots of 〈ΔCRF〉θ,J for three 1D models that address clouds differently; one is from the ICA genre that attempts to capture both horizontal variability and overlap, while the two others represent classes using homogeneous clouds with exact overlap and maximum/random overlap. The three 1D codes used here tracked their respective conditional benchmarks best. While the results in Fig. 14 are highly idealized, they nevertheless give indications of the annual march of errors in CRF to be expected from different assumptions about unresolved clouds.
The ICA genre performs best with 〈ΔCRF〉θ,J errors rarely exceeding 10%; for GATE A and ATEX they are less than 5%, or 2 W m–2. While the 1D code representing exact overlap performed well relative to its benchmarks, it is clear that use of such a code in an LSAM would result in a systematic, and rather severe, enhancement of CRF by as much as 20–50 W m–2 (assuming clouds were predicted accurately). The 1D code representing the maximum/random overlap genre does well for ATEX, though it would have done poorly had there not been so much vertical resolution working in its favor, and poor for GATE B due to omission of horizontal variability. When cloud overlap has a fairly strong maximum component, as in ATEX, or when the 1D bias as a function of µ0 is both positive and negative as in OPEN CELLS, errors in CRF can be small and comparable to ICA errors.
6. Summary and conclusions
The primary objective of this study was to assess how well current 1D radiative transfer codes interpret and handle unresolved clouds. The philosophy adopted here was to use 3D cloud-resolving model (CRM) fields as surrogates for fields that a large-scale atmospheric model (LSAM) would have produced had it sufficient spatial resolution. Applying 3D Monte Carlo (MC) photon transport algorithms to the CRM fields yielded profiles of benchmark irradiances that were used to asses 1D codes. In addition, conditional MC benchmarks were generated according to each method, or genre, of handling cloud horizontal variability and cloud overlap by 1D codes. Conditional benchmarks allowed 1D codes to be evaluated against results within the same genre rather than simply against full 3D results.
To achieve the objective stated above, it was necessary to consider a small number of clear-sky and homogeneous overcast cloud cases. The cases considered revealed that most 1D codes used for research and by weather and climate models underestimate atmospheric absorption of solar radiation. For overhead sun and the standard tropical atmosphere, this underestimation was ∼20 W m–2 relative to the benchmark models, which included an LBL code that has been compared extensively to detailed observations. The majority of this bias almost certainly results from the lack of water vapor continuum absorption in 1D models. Clearly, these errors carry over to, and complicate assessments involving, the 3D cloud fields.
The most common class of 1D model in this study assumes maximum/random overlap of homogeneous clouds. The adoption by many large-scale modeling groups of codes hard-wired with this scheme suggests the emergence of a maximum/random overlap paradigm. Results presented here imply that not all 1D codes of this variety are performing exactly as intended for they tend to systematically underestimate albedo relative to their conditional benchmark. Moreover, maximum/random overlap of homogeneous clouds is simply an extreme and incorrect set of assumptions (e.g., Hogan and Illingworth 2000; Mace and Benson-Troth 2002) that can result in substantial biases even when the code is known to perform perfectly. Forcing a code based on maximum/random overlap of cloud to acknowledge horizontal fluctuations should not be considered a solution either, for it will inevitably result in even greater underestimates of albedo due to increased transmittance of too few clouds exposed to direct beam.
There was a marked paucity of avant-garde 1D models that were able to use the independent column approximation as their conditional benchmark. This would still be the case if this study were initiated today. While those that did participate often performed best, by no means do they agree among themselves, and often they exhibit nonnegligible biases. While the ideal 1D code would account for horizontal variable cloud amid a generalized overlap scheme, there is no guarantee that biases would be reduced to satisfactory levels and that it would be flexible enough to deal with increasingly detailed portrayal of unresolved optical properties (be they cloud or other constituents including the surface). Indeed, it is difficult to imagine how 1D codes that resemble those tested here could even begin to address issues such as horizontally variable droplet sizes, vertical correlations in condensate, and cloud-surface type correlations. Furthermore, it can be expected that increasingly complex codes will demand more computation time. This is an issue that 1D modelers must accept if they intend their codes to be used operationally.
This study has generalized King and Harshvardhan's (1986) assessment of single-layer, monochromatic, two-stream approximations: no single multilayered, broadband, 1D solar code performs well for all conditions. This is partly due to gaseous transmittance parameterizations, cloud optical property parameterizations, and the two-stream approximation (as employed by most), but mostly because of inappropriate cloud overlap assumptions, incorrect application of overlap assumptions, neglect of horizontal variability of cloud, and inappropriate assumptions about horizontal variability. It may be that there are far too many nonlinear aspects to the problem to be able to satisfy all needs with a single application of a code that operates on tangible information and uses a justifiable amount of computation time. As pessimistic as this sounds, it should be seen as a challenge by begging the question that the paradigm of constructing 1D radiative transfer codes that directly incorporate assumptions about unresolved cloud structure is inadequate. As such, the nature of subgrid-scale parameterization should be reconsidered and new methodologies invented for computing radiative heating in large-scale models.
HWB thanks J. Abraham (Director MRB-MSC) for financial support and B. A. Wielicki for computing resources at NASA Langley. While participants appear alphabetically in the author list, additional acknowledgement goes to those at AER for performing additional benchmark testing, S. Cusack for work in the initial stages that helped smooth the procedure, J.-J. Morcrette, P. Räisänen, and B. Bonnel for performing additional computations, and E. E. Clothiaux for extensive proofreading of this manuscript. Thanks also to G. L. Potter for initial motivation, R. G. Ellingson, and W. B. Rossow for presenting our case, W. J. Wiscombe for encouragement, R. Pincus for discussions, and an anonymous reviewer for comments that improved the readability of this manuscript. Providers of cloud-resolving model data are also acknowledged: W. Grabowski, V. Grubišić, M. Moncreiff, B. Stevens, and X. Wu.
Corresponding author address: Dr. Howard Barker, Cloud Physics Research Division (ARMP), Meteorological Service of Canada, 4905 Dufferin St., Downsview, ON M3H 5T4, Canada. Email: firstname.lastname@example.org