# Entrapment: An Important Mechanism to Explain the Shortwave 3D Radiative Effect of Clouds

Robin Hogan, Mark Fielding, Howard Barker, Najda Villefranque, Sophia

Schäfer

► To cite this version:

Robin Hogan, Mark Fielding, Howard Barker, Najda Villefranque, Sophia Schäfer. Entrapment: An Important Mechanism to Explain the Shortwave 3D Radiative Effect of Clouds. Journal of the Atmospheric Sciences, 2019, 76 (7), pp.2123-2141. 10.1175/JAS-D-18-0366.1. hal-04790704

HAL Id: hal-04790704

https://ut3-toulouseinp.hal.science/hal-04790704v1

Submitted on 21 Nov 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Entrapment: An Important Mechanism to Explain the Shortwave 3D Radiative Effect of Clouds

ROBIN J. HOGAN AND MARK D. FIELDING

European Centre for Medium-Range Weather Forecasts, Reading, United Kingdom

HOWARD W. BARKER

Environment and Climate Change Canada, Toronto, Canada

NAJDA VILLEFRANQUE

CNRM, Météo-France/CNRS, and LAPLACE, Université de Toulouse, Toulouse, France

SOPHIA A. K. SCHÄFER

Meteorological Institute, Ludwig-Maximilians-Universität München, Munich, Germany

(Manuscript received 19 December 2018, in final form 1 May 2019)

## ABSTRACT

Several mechanisms have previously been proposed to explain differences between the shortwave reflectance of realistic cloud scenes computed using the 1D independent column approximation (ICA) and 3D solutions of the radiative transfer equation. When the sun is low in the sky, interception of sunlight by cloud sides tends to increase reflectance relative to ICA estimates that neglect this effect. When the sun is high, 3D radiative transfer tends to make clouds less reflective, which we argue is explained by the mechanism of “entrapment” whereby horizontal transport of radiation beneath a cloud layer increases the chances, relative to the ICA, of light being absorbed by cloud or the surface. It is especially important for multilayered cloud scenes. We describe modifications to the previously described Speedy Algorithm for Radiative Transfer through Cloud Sides (SPARTACUS) to represent different entrapment assumptions, and test their impact on 65 contrasting scenes from a cloud-resolving model. When entrapment is represented explicitly via a calculation of the mean horizontal distance traveled by reflected light, SPARTACUS predicts a mean “3D radiative effect” (the difference in top-of-atmosphere irradiances between 3D and ICA calculations) of 8.1 W m⁻² for overhead sun. This is within 2% of broadband Monte Carlo calculations on the same scenes. The importance of entrapment is highlighted by the finding that the extreme assumptions in SPARTACUS of “zero entrapment” and “maximum entrapment” lead to corresponding mean 3D radiative effects of 1.7 and 19.6 W m⁻², respectively.

## 1. Introduction

A key challenge in atmospheric modeling for both weather and climate prediction is to improve the interaction of clouds with solar and thermal infrared radiation. While the representation of subgrid cloud horizontal structure and vertical overlap is often now quite sophisticated (e.g., Pincus et al. 2003; Hill et al. 2015; Di Giuseppe and Tompkins 2015), a process missing from all operational models is the horizontal transport of radiation within grid boxes. This was characterized by Hogan and Shonk (2013) as entirely

associated with flow of radiation through cloud sides, and led to the development of the Speedy Algorithm for Radiative Transfer through Cloud Sides (SPARTACUS; Hogan et al. 2016; Schäfer et al. 2016). This solver is now available as an option in the radiation scheme of the forecast model used by the European Centre for Medium-Range Weather Forecasts (ECMWF; Hogan and Bozzo 2018).

In the shortwave, the main effect of transport through cloud sides is “side illumination”: the enhanced interception of direct sunlight when the sun is low in the sky, which increases the reflectance of the scene and equivalently the magnitude of the cloud radiative effect (CRE). However, Barker et al. (2015) reported Monte

Corresponding author: Robin J. Hogan, r.j.hogan@ecmwf.int
<figure><img src="image_2.png" /></figure>

FIG. 1. Schematic illustrating SPARTACUS's three possible treatments of entrapment of solar radiation underneath clouds in the case of two randomly overlapped cloud layers each of cloud fraction 1/2: (a) zero entrapment, (b) explicit entrapment, and (c) maximum entrapment. In each panel the degree to which the downwelling irradiance has been attenuated by the clouds above is indicated by the darkness of the shading. The black-headed arrows depict representative light paths discussed in the text. The double-headed arrow in (b) indicates the horizontal distance traveled by a single light ray reflected below half-level 3.5; the corresponding *mean* horizontal distances are computed in section 4.

Carlo calculations showing that in realistic cloud scenes, the effect of introducing 3D transport was more typically to *reduce* the magnitude of the CRE, particularly when the sun is high in the sky.

Várnai and Davies (1999) characterized 3D solar radiative effects in single-layer cloud scenes in terms of four mechanisms, two of which reduce the reflectance of a cloudy scene and are therefore candidates to explain this behavior. Their “downward escape” mechanism explains how forward-scattered sunlight inside a cloud has a chance to escape through the side of a cloud and reach the surface, whereas in the corresponding independent column approximation (ICA) case it would remain within the cloud and have more chance of being scattered back to space. Welch and Wielicki (1984), Hogan and Shonk (2013), and Barker et al. (2016) argued that this escape mechanism explains why 3D effects reduce the reflectance of cumulus, stratocumulus, aircraft contrails and stochastically generated cloud fields, for high-sun conditions. This process is represented by SPARTACUS, but in this paper we present evidence to show that it is not significant enough to explain the results of Barker et al. (2015), which were for a wide range of realistic and often multilayered cloud scenes.

The second candidate mechanism from Várnai and Davies (1999) is “upward trapping,” which incorporates all light rays that (i) are reflected back to space in the ICA case but not when 3D transport is included and (ii) have a longer pathlength in 3D than ICA. Their diagram to explain how this mechanism typically acts in single-layer cloud scenes depicted an upward-traveling light ray passing horizontally through the side of a cloud above, a process that is already represented by SPARTACUS. However, it is also possible, particularly in multilayered

scenes, for trapping to occur without any transport through cloud sides, but rather as a consequence of horizontal transport entirely within a clear or cloudy region and the upward-reflected ray then intercepting the base of a cloud above. This process is not explicitly handled by SPARTACUS. It was alluded to by Barker and Davies (1992) who considered idealized single-layer clouds over a reflective surface, but not studied in detail for realistic multilayered cloud scenes.

In this paper, we seek to quantify the importance of this mechanism, which we refer to as “entrapment.” In section 2, we describe it in more detail and present a simple mathematical example to illustrate how it reduces the scene reflectance. In section 3, we describe how the limits of zero and maximum entrapment may be represented in SPARTACUS. This is followed by two sections on the more complex “explicit” entrapment calculations: section 4 describes how we estimate the horizontal distance traveled by reflected radiation, with validation against monochromatic Monte Carlo simulations, while section 5 describes how the distance traveled is used to compute how much entrapment occurs, accounting for the fractal nature of clouds. Readers uninterested in the internal workings of SPARTACUS may wish to skip sections 3–5. Then in section 6, estimates of the broadband shortwave 3D radiative effect by the new SPARTACUS solver are evaluated by comparing to Monte Carlo calculations performed on 65 diverse high-resolution scenes from a cloud-resolving model.

## 2. The concept of entrapment

The schematic in Fig. 1 illustrates how entrapment can change the reflectance of a cloud scene. Figure 1a
depicts the behavior assumed in the ICA, in which horizontal transport is ignored: incoming solar radiation scattered upward by the first cloud layer it encounters is likely to escape to space since it passes back through the same clear-sky atmosphere (similar to the “opposition effect” in vegetation; e.g., Hapke et al. 1996). Figure 1b illustrates the process of entrapment by clouds when 3D transport is permitted: radiation passing down through a clear-sky (or less optically thick) part of the atmosphere may be reflected back upward at a slantwise angle and encounter the base of a cloud due to horizontal transport within either the clear-sky or cloudy region. The depiction of upward trapping by Várnai and Davies (1999) was similar except that the two cloud layers were part of the same cloud, and the reflected ray was intercepted by the edge rather than the base of the upper layer. Since the area presented by the base of a cloud is usually much larger than its edge, we would expect the impact of trapping by the base to be greater, on average. Note that entrapment can also occur over reflective surfaces where the upward reflection is by the surface rather than a cloud.

The interception of radiation by the upper cloud layer reduces the reflectance of the scene, but the magnitude of this effect depends on how far the radiation migrates horizontally in the gap between the two cloud layers relative to the size of the clouds in the upper layer. Figure 1c depicts the extreme case in which radiation is completely horizontally homogenized in clear-sky layers. This “maximum entrapment” is actually the behavior of the original shortwave implementation of SPARTACUS described by Hogan et al. (2016), as well as other solvers such as the three-region solver in the original Edwards and Slingo (1996) radiation scheme that was adapted by Shonk and Hogan (2008) to become the “Tripleclouds” solver. Shonk and Hogan (2008) described this radiative homogenization as “anomalous horizontal transport,” which is not really accurate, as at least some of this transport occurs in reality. Nonetheless, their method to remove it and thus to move from maximum entrapment (Fig. 1c) to zero entrapment (Fig. 1a) provides the starting point for representing more realistic explicit entrapment (Fig. 1b) in SPARTACUS.

To demonstrate the importance of the difference between zero and maximum entrapment, we can make some idealizations and derive an analytic expression for scene reflectance in both cases. As in Fig. 1, we consider the two cloud layers each to have a cloud fraction of 1/2, a total cloud cover of 3/4, and to be in a vacuum over a surface with an albedo of zero. Both cloud layers have the same reflectance $R$, and scatter conservatively so that their transmittance is $T = 1 - R$. In the zero-entrapment case we apply the ICA, splitting the scene

<figure><img src="image_3.png" /><figcaption>FIG. 2. Reflectance of the idealized scene discussed in section 2, composed of two randomly overlapped nonabsorbing cloud layers, each with a cloud fraction of 1/2, over a black surface in vacuum. The solid line depicts (2) and the dashed line depicts (3).</figcaption></figure>

into four columns of equal width. One is clear sky with a reflectance of zero, two contain a single cloud layer with reflectance $R$, and the final column consists of two cloud layers, which the Adding Method (Lacis and Hansen 1974) predicts to have a reflectance of

$$R^* = R + T^2 R/(1 - R^2) = 2R/(1 + R). \quad (1)$$

Thus, the scene reflectance in the zero-entrapment limit is the weighted sum of the reflectance of the four columns:

$$R_{\text{zero}} = R/2 + R^*/4 = R(1 + R/2)/(1 + R). \quad (2)$$

This is depicted by the solid line in Fig. 2. In the limit of a perfectly reflective cloud ($R = 1$), the scene reflectance becomes $R_{\text{zero}} = 3/4$.

Now consider the other extreme: maximum entrapment as depicted in Fig. 1c. Since radiation passing through a clear-sky layer retains no memory of its horizontal location with respect to the clouds it has passed through in other layers, we consider the entire domain as a single column. Thus the reflectance of a layer with a cloud fraction of 1/2 is $R/2$, while its transmittance is $1 - R/2$. Applying the Adding Method to obtain the reflectance of the scene simply involves replacing $R$ by $R/2$ in (1), yielding

$$R_{\text{max}} = 2R/(2 + R). \quad (3)$$

This is depicted by the dashed line in Fig. 2, and in the limit of a perfectly reflecting cloud the scene reflectance becomes $R_{\text{max}} = 2/3$. This is significantly less than
the ICA value of 3/4, suggesting that entrapment is an important process to treat when modeling the 3D shortwave radiative effects of clouds.

### 3. Representing entrapment scenarios in SPARTACUS

Here we explain how SPARTACUS may be modified to represent zero and maximum entrapment, illustrated in Fig. 1, as well as the first step in representing explicit entrapment. The symbols used in more than one equation in sections 3-5 are defined in appendix B. SPARTACUS uses the Tripleclouds approach of splitting each cloudy layer into three regions, one clear (denoted a) and two cloudy (denoted b and c) with different optical depths. The radiation problem can then be written in terms of vectors and matrices; for example, $\mathbf{u} = (u^a u^b u^c)^T$ is a column vector containing the upwelling diffuse irradiances at a particular height in each of the three regions. Note that the irradiance component $u^i$ is defined to be power in region j per unit area of the entire grid box, not per unit area of region j.

To represent the full range of entrapment scenarios, two aspects of the SPARTACUS implementation of Hogan et al. (2016) require modification, both in the upward pass of the adding method [their Eqs. (24)-(30)]. The first change is to describe the direct solar beam in terms of $D_{i-1/2}$, the albedo to direct radiation of the entire atmosphere and surface below half level $i - 1/2$. This accompanies $A_{i-1/2}$, which is the corresponding albedo to downwelling diffuse radiation. (As shown in Fig. 1a, we index full atmospheric layers by i, counting down from the highest layer $i = 1$, and half level $i - 1/2$ refers to the interface between layers $i - 1$ and i.) This change mirrors the application by Hogan et al. (2018) of SPARTACUS to vegetation. Both of these albedos are matrices of the form

$$ \mathbf{A} = \begin{pmatrix} A^{aa} & A^{ba} & A^{ca} \\ A^{ab} & A^{bb} & A^{cb} \\ A^{ac} & A^{bc} & A^{cc} \end{pmatrix}, \quad (4) $$

where $A^{ik}$ is the fraction of diffuse downwelling radiation in region j that is reflected up in region k. These definitions ensure that $\mathbf{u}$ at any given height is equal to the sum of reflection of the downward diffuse irradiance $\mathbf{v}$ and the downward direct irradiance $\mathbf{s}$ [see (40) of Hogan et al. 2018]: $\mathbf{u} = \mathbf{A}\mathbf{v} + \mathbf{D}\mathbf{s}$.

The second change needed to represent entrapment concerns how these albedo matrices are translated from the top of one layer to the base of the layer above, accounting for the overlap of the clouds and associated regions in the two layers. We follow Hogan et al. (2016)

and define $A_{\text{below } i-1/2}$ as the albedo of the atmosphere just below half-level $i - 1/2$, so using the region definitions of layer $i$. Likewise, $A_{\text{above } i-1/2}$ is the albedo just above this half level, using the region definitions of layer $i - 1$. Equation (30) of Hogan et al. (2016) relates the two according to the maximum-entrapment assumption:

$$ A_{\text{above } i-1/2} = U_{i-1/2} A_{\text{below } i-1/2} V_{i-1/2}, \quad (5) $$

where $\mathbf{U}$ and $\mathbf{V}$ are the upward and downward overlap matrices. They are defined such that $\mathbf{u}_{\text{above}} = \mathbf{U}\mathbf{u}_{\text{below}}$ and $\mathbf{v}_{\text{below}} = \mathbf{V}\mathbf{v}_{\text{above}}$; that is, they map irradiances passing through a half-level on to the regions of the layer the radiation is entering. To see how (5) leads to maximum entrapment, consider what happens at half-level 3.5 in the two-region example of Fig. 1, where the matrices would be

$$ \begin{align} \mathbf{U}_{3.5} &= \begin{pmatrix} 1/2 & 0 \\ 1/2 & 0 \end{pmatrix}, \notag \\ \mathbf{A}_{\text{below } 3.5} &= \begin{pmatrix} A & 0 \\ 0 & 0 \end{pmatrix}, \notag \\ \mathbf{V}_{3.5} &= \begin{pmatrix} 1 & 1 \\ 0 & 0 \end{pmatrix}, \tag{6} \end{align} $$

where A is the albedo of the atmosphere below half level 3.5. Applying (5) yields

$$ \mathbf{A}_{\text{above } 3.5} = \begin{pmatrix} A/2 & A/2 \\ A/2 & A/2 \end{pmatrix}. \quad (7) $$

This confirms that radiation exiting the base of either the clear or cloudy regions in layer 3 has an equal probability of being reflected back up into either of these two regions.

Shonk and Hogan (2008) described how to eliminate this horizontal transport in their Tripleclouds solver and thereby achieve zero entrapment; their solution may be written as

$$ \mathbf{a}_{\text{above } i-1/2} = \mathbf{V}_{i-1/2}^{\text{T}} \mathbf{a}_{\text{below } i-1/2}, \quad (8) $$

where $\mathbf{a}$ is a column vector containing the reflectances of each region with the assumption that light is always reflected up from the same region it enters. Since Tripleclouds neglects lateral radiation flows between regions, $\mathbf{A}$ is diagonal and $\mathbf{a}$ simply contains its diagonal elements. To apply the zero-entrapment assumption to a SPARTACUS simulation that includes lateral flows between regions, $\mathbf{a}_{\text{below } i-1/2}$ is defined such that its jth element contains the sum of the jth column of $A_{\text{below } i-1/2}$. Physically this means that flows represented by the white-headed arrows in Fig. 3, which involve reflection up into a different region of the upper layer, are forced
<figure><img src="image_4.png" /><figcaption>(b) Diagonal elements</figcaption></figure>

<figure><img src="image_5.png" /></figure>

Layer:

i-1

i

FIG. 3. Schematic illustrating the radiation paths considered when computing entrapment at half-level i − 1 / 2. (a) The four paths that contribute to A ab below i − 1 / 2 , which is the off-diagonal element of A below i − 1 / 2 that represents the fraction of radiation downwelling into region a (clear sky) of layer i that is scattered back up into region b (cloud) of the same layer. (b) The four paths that contribute to A ba below i − 1 / 2 , which is the diagonal element of A below i − 1 / 2 that represents the fraction of radiation downwelling into region a of layer i that is scattered back up in the same region. Note that in practice we include radiation paths that pass through the base of layer i and are reflected back by the lower layers. The white-headed arrows represent radiation flows involving reflection back up into a different region of layer i − 1.

to be reflected up into the same region. After applying (8), A above i − 1 / 2 is defined to be a diagonal matrix with its diagonal elements given by a above i − 1 / 2 and the rest of SPARTACUS is unchanged. Applying zero entrapment to the matrices in (6) yields

$$
\mathbf{A}_{\text{above } 3.5} = \begin{pmatrix} A & 0 \\ 0 & A \end{pmatrix}, \qquad (9)
$$

which indicates that, as required, radiation exiting the
base of either region of layer 3 is reflected back up into
the same region.

To represent explicit entrapment we first assume that
radiation that travels far enough horizontally to pass
through a cloud edge (represented by the off-diagonal
elements of A below i − 1 / 2 and shown in Fig. 3a) has lost
memory of the region it originated from in layer i − 1, so
can be treated by maximum entrapment. By contrast, the
destination of reflected radiation that does not pass
through a cloud edge (the diagonal elements of A below i − 1 / 2
illustrated in Fig. 3b) is computed taking account of the
mean horizontal distance traveled. Mathematically this is
achieved (in the two-region case) by

$$
\mathbf{A}_{\text{above } i-1/2} = \mathbf{U}_{i-1/2} \begin{pmatrix} 0 & A_{\text{below } i-1/2}^{\text{ba}} \\ A_{\text{below } i-1/2}^{\text{ab}} & 0 \end{pmatrix} \mathbf{V}_{i-1/2} \\
+ \sum_{j=a,b} A_{\text{below } i-1/2}^{jj} \mathbf{L}_{i-1/2}^{j} \begin{pmatrix} V_{i-1/2}^{aj} & 0 \\ 0 & V_{i-1/2}^{bj} \end{pmatrix}, \tag{10}
$$

and similarly for D above i − 1 / 2. The first term on the right-
hand side is the maximum-entrapment formula, (5),
applied just to the off-diagonal elements. The second
term on the right-hand side contains A j j below i − 1 / 2 , the j th

diagonal of A below i − 1 / 2 , as well as elements from the j th
row of V i − 1 / 2 .

Matrix L j i − 1 / 2 is central to SPARTACUS’s represen-
tation of explicit entrapment: it expresses how much
lateral transport occurs for reflected radiation within
region j of layer i , but also accounting for radiation
passing down through the layers below. Its elements
quantify the weight of each of the arrows in Fig. 3b.
Since it repartitions radiation between regions without
changing the total energy, its columns sum to 1 (i.e., it
is a left stochastic matrix). If we wished for the diagonal
elements of A below i − 1 / 2 to be treated with zero entrap-
ment (thereby eliminating the white-headed arrows in
Fig. 3b), then we could simply set L j i − 1 / 2 to the identity
matrix leading to the second term on the right-hand side
behaving exactly as (8). Otherwise, L j i − 1 / 2 is computed
in a two-step process. Section 4 describes and validates
the first step: computing the mean horizontal distance
traveled by reflected radiation. Section 5 describes
and validates the second step: using the mean horizontal
distance traveled to compute the elements of matrix
L j i − 1 / 2 and hence how much radiation is “trapped” by
passing across the dashed line in Fig. 3b.

4. Explicit entrapment: Horizontal distance
traveled by reflected radiation

The white-headed arrow in Fig. 1b illustrates the
horizontal distance traveled by a single light ray re-
flected below half-level 3.5, and includes the horizontal
distance associated with both the downward and upward
parts of the journey. This section deals with the task of
estimating the *mean* horizontal distance traveled by re-
flected radiation below a particular half-level *i* − 1/2,
considering all possible light paths, including those that
penetrate down and up through multiple layers. However,
we exclude light paths that pass laterally through cloud
boundaries, as these are treated by maximum entrapment embodied in the first term on the right-hand side of (10). Since the downward part of the journey depends on whether the incoming radiation at half-level i − 1/2 is diffuse or direct, we compute separate diffuse and direct mean horizontal distances, denoted $x_{i-1/2}$ and $y_{i-1/2}$, respectively. But note that light rays contributing to $y_{i-1/2}$ can only strictly be considered to be “direct” sunlight until their first scattering event.

a. Method

We consider a plane-parallel atmosphere (i.e., containing one region in each layer), and later adapt it to multiple regions. The horizontal distance traveled by direct radiation passing through a single layer i is

$$ \Delta y_i = \Delta z_i \tan \theta_0, \quad (11) $$

where $\Delta z_i$ is the physical thickness of the layer and $\theta_0$ is the solar zenith angle. It was shown by Schäfer et al. (2016) that the equivalent expression for isotropic diffuse radiation is

$$ \Delta x_i = \Delta z_i \pi / 2. \quad (12) $$

We use these as building blocks in our estimate of distance traveled during reflection from multiple layers. Consider how the two-stream Adding Method is applied for diffuse radiation in a single layer i: given the layer reflectance $R_i$ and transmittance $T_i$, the scene albedo at the top of the layer is given by

$$ A_{i-1/2} = R_i + T_i^2 A_{i+1/2} [1 + R_i A_{i+1/2} + (R_i A_{i+1/2})^2 + \dots], \quad (13) $$

where $A_{i+1/2}$ is the scene albedo at the bottom of the layer and the terms in the square brackets represent multiple reflections between layer i and the atmosphere below. This is a geometric series that reduces to

$$ A_{i-1/2} = R_i + T_i^2 A_{i+1/2} / (1 - R_i A_{i+1/2}). \quad (14) $$

It can be adapted to compute $x_{i-1/2}$, the mean horizontal distance traveled by reflected radiation below the top of the layer, as the weighted average of the distances associated with each order of scattering in (13):

$$ \begin{aligned} A_{i-1/2} x_{i-1/2} = & R_i \hat{x}_i + T_i^2 A_{i+1/2} \\ & \times \{2 \hat{x}_i + x_{i+1/2} \\ & + R_i A_{i+1/2} [\hat{x}_i + \sqrt{2} (\hat{x}_i + x_{i+1/2})] \\ & + (R_i A_{i+1/2})^2 [\hat{x}_i + \sqrt{3} (\hat{x}_i + x_{i+1/2})] \\ & + \dots \}, \end{aligned} \quad (15) $$

where $x_{i+1/2}$ is the mean distance traveled by reflected radiation below the base of layer i. The first term on the right-hand side contains $\hat{x}_i$, the mean horizontal distance traveled by radiation that is reflected by layer i, rather than penetrating the layer and being reflected by the layers below. We assume that, on average, such radiation penetrates to the center of the layer before being reflected back out (hence traveling a distance $\Delta x_i/2$ on each of the downward and upward paths), and that the azimuthal scattering angle at the point of reflection is random (so these distances should be added in quadrature) leading to

$$ \hat{x}_i = [2(\Delta x_i/2)^2]^{1/2} = \Delta x_i / \sqrt{2}. \quad (16) $$

This neglects the additional distance associated with multiple scattering entirely within a layer, but the good performance reported in section 4b suggests that this is a small effect in practice.

The first line in the curly brackets in (15) represents radiation that passes down through the entire layer i and back up again, so the horizontal distance associated with transiting the layer is twice that of radiation reflected by the layer (the $2\hat{x}_i$ term), and is added to the distance associated with reflection by the layers below half level i + 1/2 (the $x_{i+1/2}$ term). The subsequent lines in the curly brackets in (15) include multiple reflections between layer i and the layers below half level i + 1/2. Each reflection adds $\hat{x}_i + x_{i+1/2}$ to the distance traveled, but due to the random azimuthal scattering angle, they again should be added in quadrature leading to the square-root weighting term.

Equation (15) may be rearranged to obtain

$$ \begin{aligned} x_{i-1/2} = & \hat{x}_i + \frac{T_i^2 A_{i+1/2}}{A_{i-1/2}} (\hat{x}_i + x_{i+1/2}) \\ & \times \sum_{j=0}^{\infty} \sqrt{j+1} (R_i A_{i+1/2})^j. \end{aligned} \quad (17) $$

Can we reduce this infinite series to a closed-form expression as before? Equation (14) exploited the fact that $\sum_{j=0}^{\infty} a^j = (1-a)^{-1}$, and by differentiating this expression we find that $\sum_{j=0}^{\infty} (j+1)a^j = (1-a)^{-2}$. The infinite series in (17) lies, in some sense, between these two series, and we find empirically that to a good approximation $\sum_{j=0}^{\infty} \sqrt{j+1} a^j \approx (1-a)^{-1.5}$, which has errors of less than 10% for $a < 0.9$. Thus we approximate (17) by

$$ x_{i-1/2} \approx \hat{x}_i + \frac{T_i^2 A_{i+1/2} (\hat{x}_i + x_{i+1/2})}{A_{i-1/2} (1 - R_i A_{i+1/2})^{1.5}}. \quad (18) $$

This equation may be applied sequentially from the surface up through the atmosphere to obtain a profile
of $x_{i-1/2}$. The surface value is $x_{n+1/2} = 0$, since reflection right at the surface is not associated with any horizontal transport.

Next we seek an equivalent expression for $y_{i-1/2}$, the horizontal distance traveled by reflected direct radiation. The equivalent expression to (14) for the albedo to direct radiation is

$$D_{i-1/2} = S_i^+ + T_i \frac{S_i^- A_{i+1/2} + E_i D_{i+1/2}}{1 - R_i A_{i+1/2}}, \quad (19)$$

where three new layer properties have been introduced: $E_i$ is the fraction of direct radiation that penetrates the layer without being scattered, $S_i^-$ is the fraction that penetrates the layer but is scattered on the way (so emerging into the diffuse stream at the base of the layer), while $S_i^+$ is the reflectance of the layer to direct radiation.

The mean horizontal transport associated with reflection by the $S_i^+$ term we denote as $\hat{y}_i$, and is assumed to be associated with direct radiation that penetrates to the center of the layer followed by a scattering event and a diffuse path back to the top of the layer. Therefore, we add $\Delta y_i/2$ and $\Delta x_i/2$ in quadrature:

$$\hat{y}_i = \frac{1}{2}(\Delta y_i^2 + \Delta x_i^2)^{1/2}. \quad (20)$$

The two terms on the numerator of the right-hand side of (19) represent the two ways that direct radiation can penetrate the layer: with and without being scattered on the way. When scattering occurs, the diffuse albedo $A_{i+1/2}$ is involved, which is associated with diffuse horizontal distance $x_{i+1/2}$, whereas when scattering does not occur, the direct albedo $D_{i+1/2}$ is involved, which is associated with direct horizontal distance $y_{i+1/2}$. Subsequent reflections, governed by the denominator of (19), all involve the addition of the diffuse horizontal distance $x_{i+1/2}$. Applying the logic of (15) leads to terms involving the addition, in quadrature, of $y_{i+1/2}$ (associated with a direct reflection) and one or more of $x_{i+1/2}$ (associated with one or more internal diffuse reflections). Unfortunately, these terms do not reduce conveniently to a closed-form expression like (18). Therefore, we make the approximation that in any term involving the combination of $y_{i+1/2}$ and $x_{i+1/2}$, $y_{i+1/2}$ can be replaced by $x_{i+1/2}$, leading to

$$y_{i-1/2} \approx \hat{y}_i + \frac{T_i}{D_{i-1/2}} \times \{[S_i^- A_{i+1/2} \xi + E_i D_{i+1/2}(\xi-1)](\hat{x}_i + x_{i+1/2}) + E_i D_{i+1/2}(\hat{y}_i + y_{i+1/2})\}, \quad (21)$$

where

$$\xi = (1 - R_i A_{i+1/2})^{-1.5}. \quad (22)$$

In a plane-parallel atmosphere, (14), (18), (19), and (21) may be applied sequentially from the surface to top of atmosphere (TOA), to obtain profiles of $x_{i-1/2}$ and $y_{i-1/2}$. We stress that even though the calculation proceeds in a single upward pass through the atmosphere, the computed mean horizontal distances include both the downward and upward parts of the journey of reflected light rays. We need the horizontal distance traveled in partially cloudy profiles, so seek vectors $\mathbf{x}$ and $\mathbf{y}$ whose $j$th elements contain the horizontal distances associated with region $j$. As illustrated in Fig. 3b, each region is considered independently, so we may still use these four equations to step the elements of $\mathbf{x}$ and $\mathbf{y}$ from the base of the layer to the top. The other inputs to these equations, $R_i$, $T_i$, $E_i$, $S_i^+$, $A_{i+1/2}$, and $D_{i+1/2}$ are taken as the diagonal elements to the corresponding matrices available in the SPARTACUS computation. Physically, the diagonal elements are used because we are interested in horizontal transport that remains within a region in layer $i$; radiation that passes laterally between regions in layer $i$ was dealt with by (10). The final aspect to deal with partially cloudy profiles is to translate $\mathbf{x}$ and $\mathbf{y}$ from the regions below half-level $i-1/2$ to the regions above. We use the relevant overlap matrix similarly to the operation in (8):

$$\mathbf{x}_{\text{above } i-1/2} = \mathbf{V}_{i-1/2}^{\text{T}} \mathbf{x}_{\text{below } i-1/2}, \quad (23)$$

and likewise for $\mathbf{y}$.

## b. Evaluation

Here we evaluate the estimates of mean horizontal distance traveled by reflected radiation as a function of height (the values of $\mathbf{x}$ and $\mathbf{y}$ above), using Monte Carlo calculations by the model of Villefranque et al. (2019, manuscript submitted to J. Adv. Model. Earth Syst.), which implements ray-tracing techniques from computer graphics and permits the paths of individual photons to be tracked. The results are shown in Fig. 4 for four cloud scenes and three solar zenith angles in simulations using periodic boundary conditions in the horizontal. All are at a single wavelength in vacuum with idealized cloud optical properties over a Lambertian surface with an albedo of 0.2. The first profile (Figs. 4a–d) consists of a plane-parallel cloud layer containing isotropic scatterers with an optical depth of 1 and a single-scattering albedo of 0.999999. Beneath the cloud, all reflection is from the surface so the mean horizontal distance traveled increases linearly with height above
<figure><img src="image_6.png" /><figcaption>FIG. 4. The mean horizontal distances traveled by reflected light for diffuse and direct downwelling radiation (x and y in the text, respectively), as computed from Monte Carlo simulations (solid lines) and estimated by SPARTACUS (dashed lines) using the method described in section 4. Each row of panels corresponds to a different atmospheric profile described in the text, with (a), (e), (i), (m) the cloud fraction and subsequent columns showing results for solar zenith angles (SZAs) of (b), (f), (j), (n) 0°, (c), (g), (k), (o) 44°, and (d), (h), (l), (p) 70°. Panels (m)–(p) correspond to the cloud shown in Fig. 5c. The error bars show the error on the mean from the Monte Carlo calculations; note that the Monte Carlo model can only compute mean distances given a sufficient number of photons, so it does not show x above cloud top or y deep in an optically thick cloud.</figcaption></figure>

the surface z. The direct mean horizontal distance y increases with θ₀ due to (11), but is not zero for overhead sun since it includes the return journey from the surface to height z as diffuse radiation. Within the cloud a fraction of downwelling radiation is reflected by the cloud, rather than the surface, and so the mean

horizontal distance is reduced. We see that the SPARTACUS estimates using the method described above are accurate to around 10% for x and 3% for y. The second profile (Figs. 4e–h) is the same but with an optical depth of 5. The SPARTACUS errors are somewhat larger at around 18%.
The last two profiles contain more realistic clouds. Both assume an asymmetry factor of 0.86, appropriate for liquid clouds in the midvisible. SPARTACUS then performs the usual delta-Eddington scaling, treating some of the forward scattered light as if it had not been scattered at all. To achieve a fair comparison in terms of the definition of direct and “diffuse” radiation, the Monte Carlo model takes the delta-Eddington-scaled extinction coefficient, and assumes a Henyey-Greenstein scattering phase function using the delta-Eddington-scaled asymmetry factor value of 0.462. Figures 4i–l show the results for a 6.4 km × 6.4 km large-eddy simulation of cumulus clouds from Brown et al. (2002), which was also used by Hogan et al. (2016) and is based on an observed case from the Atmospheric Radiation Measurement (ARM) program. Figures 4m–p show the results for a 100 km × 100 km scene from a 250-m simulation by the Canadian Global Environmental Multiscale (GEM) model of a multilayer liquid cloud (Pacific scene 16). The GEM scenes are described in detail in section 6. In both the ARM and the GEM cases, the typical SPARTACUS errors are 25% for x and 6% for y. Given the simplifications involved, SPARTACUS performs very well in estimating horizontal distance traveled, and should be adequate to feed into the final step for computing entrapment.

## 5. Explicit entrapment: How much radiation is trapped?

### a. Method

Here we use the estimated mean horizontal distances from the previous section to compute the matrix $L_{i-1/2}^j$ in (10). The meaning of $L_{i-1/2}^j$ in the two-region case can be explained by considering how it acts on $w_{i-1/2}^j(x) = (w^{ja} w^{jb})^T$, where $w^{jk}$ is the fraction of radiative energy in region j of layer i that lies beneath region k in layer i − 1. It is denoted as a function of the mean horizontal distance traveled, x, since the radiation entered the layer. Matrix $L_{i-1/2}^j$ is then defined such that

$$w_{i-1/2}^j(x_{i-1/2}^j) = L_{i-1/2}^j w_{i-1/2}^j(0), \quad (24)$$

where $w^j(0)$ represents the radiation partitioning at the instant the radiation enters the layer from above (when the horizontal distance traveled is zero) and $w_{i-1/2}^j(x_{i-1/2}^j)$ is the radiation partitioning at the instant it leaves the layer after being reflected (when the mean distance traveled is $x_{i-1/2}^j$). Hence, the diagonal elements of $L_{i-1/2}^j$ are represented by the black-headed arrows in Fig. 3b and the off-diagonals by the white-headed arrows. Note that (10) includes the albedo $A_{\text{below } i-1/2}^{jj}$ so $L_{i-1/2}^j$ describes only the redistribution of

energy, and therefore the sum of the elements in $w_{i-1/2}^j(x)$ is constant at 1 for any x.

How can we compute $L_{i-1/2}^j$, accounting for the possibility of radiation passing across the dashed line in Fig. 3b more than once? This can occur via multiple scattering events in the layers below $i-1/2$, as well as via radiation passing beneath small clouds and emerging on the other side. The problem is analogous to the original SPARTACUS problem of working out the net exchange of radiation through cloud sides, and so we can use the same method. Essentially, we wish to solve a system of coupled differential equations of the form

$$\frac{d w_{i-1/2}^j}{dx} = \Gamma w_{i-1/2}^j, \quad (25)$$

where the matrix $\Gamma$ contains the rates of radiation exchange between the “subregions” of region j (with subregions defined by the regions of layer $i-1$ that they are beneath, illustrated by the dashed line in Fig. 3b), and

$$\Gamma = \begin{pmatrix} -f^{ab} & +f^{ba} \\ +f^{ab} & -f^{ba} \end{pmatrix}, \quad (26)$$

where $f^{kl}$ is the rate at which radiation (direct or diffuse) is transported from subregion k to l, per unit increase in horizontal distance traveled x. The solution to (25) is (24) but with $L_{i-1/2}^j$ written as a matrix exponential:

$$L_{i-1/2}^j = \exp(\Gamma x_{i-1/2}^j). \quad (27)$$

The repeated elements in (26) mean that the matrix exponential may be computed efficiently for both the two- and three-region cases as described in the appendix of Hogan et al. (2018). The same method is used to compute the contribution to the direct albedo matrix, but using the direct horizontal distance traveled $y_{i-1/2}^j$.

### b. Representing fractal behavior

To test the validity of this approach, we use the contrasting binary cloud scenes shown in Figs. 5a–d, which have been generated by applying an optical-depth threshold to four of the GEM simulations described in section 6. A scalar field is defined containing a value of 1 in the clear (black) areas and 0 in the cloudy (white) areas, which can be thought of as solar radiation that has passed through the gaps between the clouds. Gaussian smoothing is then applied to the field with varying smoothing scales x, representing horizontal radiation transport beneath the cloud. Previous studies of the interaction of radiation and clouds have found a Gaussian to be reasonably good at describing the horizontal distribution of diffuse radiation originating from a point
<figure><img src="image_7.png" /><figcaption>FIG. 5. Evaluation of the part of the algorithm that computes the fraction of radiation trapped due to horizontal transport. (a)–(d) Binary cloud masks obtained by applying a visible-optical-depth threshold to four GEM scenes. (e)–(h) The corresponding fraction of downwelling clear-sky radiation that is trapped beneath a cloud as a function of the mean horizontal distance traveled $x$; the black line shows the result of a Gaussian smoothing of the cloud masks (treated as truth), while the gray lines show two candidate methods for SPARTACUS that take as input the cloud scale $S$ [values indicated above (e)–(h)] and cloud cover.</figcaption></figure>

source (e.g., Hogan and Battaglia 2008; Wissmeier et al. 2013). The fraction of the total scalar field that is then in the cloudy parts of the domain is the “trapped fraction,” and is shown by the black lines in Figs. 5e–h. The dotted lines show the cloud cover, which corresponds to the trapped fraction one would expect if the radiative energy were completely homogenized horizontally (maximum entrapment).

To apply the matrix-exponential method described above to estimate the trapped fraction we need to define the lateral exchange rates in (26). From the geometric arguments of Hogan and Shonk (2013), we would expect

$$f^{kl} = L^{kl}/\pi c^k, \qquad (28)$$

where $L^{kl}$ is the length of the interface between regions $k$ and $l$ (i.e., the perimeter length of the clouds) per unit area of the domain, and $c^k$ is the fractional area of region $k$, both of which can be obtained by analyzing the binary cloud fields. We then apply (26) and (27) to obtain $L$, the bottom-left element of which is the fraction of radiative energy beneath the cloud shown by the solid gray lines in Figs. 5e–h. It can be seen that this is a good prediction of the trapped fraction for low values of $x$, but for larger $x$ it overestimates entrapment significantly.

This can be explained by the fractal nature of clouds. The Hogan and Shonk (2013) definition of lateral exchange

rates imposes a length scale on the cloud field based on the perimeter length. For example, we could define the effective cloud scale $S$ [see (20) of Hogan et al. 2018] such that normalized perimeter length is

$$L^{ab} = 4c^a(1 - c^a)/S, \qquad (29)$$

where $S$ can be thought of as the size that equally sized squares would need to have if they were placed randomly on a grid and their fractional cover and total perimeter length were equal to the values for the actual cloud field. The $S$ values for the scenes in Fig. 5 are shown above panels e–h. Substitution of (29) into (28) gives $f^{ab} = 4(1 - c^a)/\pi S$, indicating that the rate of exchange between regions is inversely proportional to $S$. Thus, if all clouds indeed had a diameter of around $S$ then we would expect the trapped fraction to quickly approach the asymptotic value of the cloud cover for $x \gg S$.

In reality the clouds span a wide range of scales, and the presence of very large clouds reduces the trapped fraction for larger $x$. Another way of looking at this is to recognize that since clouds are fractal, the effective perimeter length $\hat{L}^{ab}$ ought to be a power-law function of the length scale at which the cloud field is being probed. This is usually written as $\hat{L}^{ab} \propto x^{1-D}$, where $D$ is the fractal dimension. Many studies have estimated the
fractal dimension of clouds, with D = 1.5 being a reasonable representative value (e.g., Cahalan and Joseph 1989; Gotoh and Fujii 1998; Wood and Field 2011), implying $\hat{L}^{ab} \propto x^{-1/2}$.

The dashed gray lines in Figs. 5e-h show the result of making effective perimeter length depend on $x$ in this way, with the formula leading to the best fit given by

$$ \hat{L}^{ab} = L^{ab} \min(1, \sqrt{0.4S/x}). \quad (30) $$

This formula caps the effective perimeter length $\hat{L}^{ab}$ to be no larger than the measured value $L^{ab}$ for $x < 0.4S$. The much better agreement with the curves computed from actual cloud fields gives us confidence that this formula is appropriate to use in computing entrapment within SPARTACUS. Note that we have not found it necessary to use (30) also for computing lateral exchange across cloud boundaries, because the relevant length scale would be the horizontal distance traveled by radiation as it passes through an individual model layer, which is much less than $x$.

### c. Treatment of overhanging clouds

Applying these findings in SPARTACUS presents one further issue to resolve, since as shown in Fig. 3b we are not dealing with radiative exchange between regions, but exchange between the subregions of region $j$ in layer $i$, defined according to the regions above them in layer $i-1$. Unfortunately the perimeter length of the interface between these subregions is not completely defined by the variables available to SPARTACUS. Consider the case of two layers, each with a cloud fraction of 0.5, an overlap parameter of $\alpha_{i-1/2} = 0.5$ (i.e., halfway between maximum and random overlap) such that the combined cloud cover of the two layers is 0.625, and a particular value for the effective normalized cloud perimeter length in the top layer $\hat{L}_{i-1}^{ab}$. Figure 6a depicts the way in which the clouds could be configured that maximizes the perimeter length between the subregions of region $a$ in layer $i$, which in this 2D diagram is illustrated by the number of overhanging clouds shown by the dashed lines.

Figure 6b depicts an alternative configuration at the other extreme: in the left half of the grid box the clouds are overlapped maximally and in the right half they are overlapped randomly. Since the overhangs are only associated with the randomly overlapped half of the grid box, there are fewer overhangs. This treatment of cloud overlap was explored by Shonk et al. (2010), who showed that the area of the grid box associated with maximum overlap of region $j$ could be written as $\alpha_{i-1/2} \min(c_{i-1}^j, c_i^j)$. We are concerned with removing the fraction of region $j$ in layer $i$

<figure><img src="image_8.png" /><figcaption>(b) Fewest overhangs: $\zeta=0$</figcaption></figure>

<figure><img src="image_9.png" /></figure>

FIG. 6. Illustration of two model cloud configurations with the same cloud fraction and overlap parameter: (a) the configuration that maximizes the number of “overhangs” where explicit entrapment must be calculated (shown by the dashed vertical lines with curved arrows across them) and (b) the configuration that minimizes the number of overhangs. In SPARTACUS this is controlled by the $\zeta$ factor shown above each panel.

that is maximally overlapped, so divide through by $c_i^j$ to obtain the fraction of region $j$ that is randomly overlapped:

$$ C_{i-1/2}^j = 1 - \alpha_{i-1/2} \min(c_{i-1}^j, c_i^j) / c_i^j. \quad (31) $$

Thus, the most overhangs (Fig. 6a) is obtained by using $\hat{L}_{i-1}^{kl}$ for the effective perimeter length between regions $k$ and $l$, while the fewest overhangs (Fig. 6b) is obtained by using $\hat{L}_{i-1}^{kl} C_{i-1/2}^j$. There is no theoretical or observational evidence to support which is the most likely, so we introduce a user-defined “overhang factor” $\zeta$ that varies the effective perimeter length linearly between most overhangs ($\zeta=1$) and fewest overhangs ($\zeta=0$). This factor is left as a parameter to be determined in section 6 according to which leads to SPARTACUS predicting 3D radiative effects most accurately. It turns out that this property is of second-order importance compared to whether entrapment is treated at all.

## 6. Results

In this section we evaluate the shortwave 3D radiative effect predicted by the new SPARTACUS implementation in the ecRad radiation scheme (Hogan and Bozzo 2018), and investigate the impact of various different treatments of entrapment. We have used 65 scenes generated from simulations by Environment and Climate Change Canada’s GEM model (Girard et al. 2014), using the configuration described by Leroyer et al. (2014) with the Milbrandt and Yau (2005) double-moment bulk microphysics cloud scheme.

Each scene measures 100 km × 100 km, has a horizontal resolution of 250 m and employs 56 vertical
layers. The simulations were originally performed to generate synthetic satellite data from two swaths: an Atlantic swath on 7 December 2014 from Greenland to the Dominican Republic, from which 39 scenes were extracted, and a Pacific swath on 24 June 2015 from Hawaii to Tonga, from which a further 26 scenes have been extracted. Thus, the scenes span a wide range of cloud conditions.

Both ICA and 3D Monte Carlo shortwave radiative transfer calculations have been performed on these scenes using the model of Barker et al. (2003), which tracks photons through sequences of scattering events until they are either absorbed by a particle, molecule, or the surface, or exit the domain's top. Calculations were performed for solar zenith angles at 5° intervals between 0° and 85°, but with random solar azimuth angle, and assuming a periodic domain. The Rapid Radiative Transfer Model for GCMs (RRTM-G) of Iacono et al. (2008) was used to represent gas absorption, the Yi et al. (2013) scheme for ice optical properties and Mie theory for liquid droplets. Scattering by air molecules and cloud particles were handled by the Rayleigh and Henyey-Greenstein phase functions, respectively. To simplify the comparison with 1D radiation schemes, all calculations assumed a Lambertian surface with an albedo of 0.05.

We first compare the TOA cloud radiative effect between Monte Carlo calculations run in an ICA mode (Monte Carlo ICA) and Tripleclouds (the SPARTACUS control) for the same scenes, that is, in the absence of 3D radiative transfer. In addition to cloud fraction and gridbox-mean liquid and ice mixing ratio, Tripleclouds takes as input the overlap parameter at each half-level and the fractional standard deviation of in-cloud water content, FSD, in each layer. It was found that the original implementation of Tripleclouds was not capable of accurately representing the effect of horizontal heterogeneities for FSD > 1.5, which occurs in many of these scenes. Appendix A describes an improvement to Tripleclouds and SPARTACUS that has overcome this problem. Figures 7a-c reveal that the resulting root-mean-square error (RMSE) in CRE predicted by Tripleclouds is around 10% and its bias is only 1%–2%.

The differences between 3D and ICA Monte Carlo calculations of CRE for the 65 scenes are summarized by the black box-and-whisker plots in Figs. 8a and 8b, and the mean by the thick black line. We see that for $\theta_0$ in the range 0°–75°, 3D effects increase both TOA and surface CRE (i.e., make them less negative), and therefore act in the sense of warming the Earth system by making clouds less reflective. In individual cases, 3D effects can act to cool the Earth system by up to $7 \text{ W m}^{-2}$ at $\theta_0 = 75°$, via interception of direct solar radiation by cloud sides, but

plenty of other scenes have a warming effect even for large $\theta_0$.

The various SPARTACUS simulations, which all use Tripleclouds as their 1D control, enable us to elucidate the role of entrapment in explaining this behavior. SPARTACUS has been run taking as input the observed cloud perimeter length in each layer $L^{ab}$ and the length of the contour separating the optically thinner and thicker parts of the in-cloud extinction distribution $L^{bc}$. In practice the former is characterized by the cloud effective scale defined by (29), and the latter by a cloud heterogeneity scale $S^{\text{het}}$ given analogously by

$$L^{bc} = 4c^c(1 - c^c)/S^{\text{het}}. \quad (32)$$

The blue lines in Figs. 8a and 8b show the zero-entrapment scenario, which was illustrated in Fig. 1a. The 3D effect is much weaker overall; at TOA, the greatest mean cooling is $1.9 \text{ W m}^{-2}$ at $\theta_0 = 75°$, and the greatest mean warming is $1.7 \text{ W m}^{-2}$ at $\theta_0 = 0°$. This general pattern can be explained by the mechanisms of shortwave side illumination and downward escape discussed by Hogan and Shonk (2013), but downward escape is clearly too weak a mechanism to explain the strong 3D effect found in the Monte Carlo simulations for overhead sun. At the other extreme, the green lines show the maximum entrapment scenario, which was illustrated in Fig. 1c and involves complete horizontal homogenization of radiation in each clear or cloudy region. This has a strong warming effect, reaching $19.6 \text{ W m}^{-2}$ at TOA for overhead sun, but is over twice as strong as the reference Monte Carlo calculations.

The other two SPARTACUS simulations in Figs. 8a and 8b are much closer to the reference calculations: the red and pink lines show results using explicit entrapment described in sections 4 and 5, with the two treatments of cloud overhangs illustrated in Fig. 6. It is clear that the least-overhang scenario ($\zeta = 0$) agrees with the Monte Carlo calculations best at TOA; there is still a slight overestimate of the warming effect of 3D radiative transfer, but it is less than 8% for $\theta_0 \le 55°$. Therefore, $\zeta = 0$ is our preferred SPARTACUS configuration for the remainder of the paper. Figures 7d-f compare individual 3D radiative effects between Monte Carlo and this SPARTACUS configuration, and while there is some scatter, the correlation coefficients of 0.79–0.86 confirm that there is skill in predicting 3D effects for individual cases.

Figure 8c shows the change to total atmospheric absorption when 3D effects are included. The Monte Carlo calculations show an increase in absorption by around $1 \text{ W m}^{-2}$ at most solar zenith angles, which is comparable to the findings of Barker et al. (2016). Both the main
<figure><img src="image_10.png" /><figcaption>FIG. 7. (a)–(c) Comparison of TOA shortwave CRE between reference ICA calculations using the Monte Carlo (MC) model on the horizontal axis and Tripleclouds (the SPARTACUS control) on the vertical axis for each of the 65 GEM scenes at three different solar zenith angles of 20°, 50°, and 80°, respectively. (d)–(f) Corresponding comparison of the 3D radiative effect, that is, the difference in CRE between 3D and 1D calculations. The SPARTACUS calculations use explicit entrapment with an overhang factor of $\zeta = 0$, also shown by the red lines in Fig. 8. Each panel also states the mean of the MC calculations, and the bias and RMSE of the SPARTACUS/Tripleclouds calculations with respect to MC.</figcaption></figure>

3D mechanisms contribute to this effect: side illumination at large solar zenith angle enhances the interception and hence absorption by clouds, while Fig. 3b shows that entrapment increases the pathlength of radiation in clear skies beneath cloud, enhancing water vapor absorption. SPARTACUS with explicit entrapment leads to around $2 \text{ W m}^{-2}$ greater atmospheric absorption than Tripleclouds, on average, which is twice the 3D effect in the Monte Carlo simulations. This is related to the presence of a handful of outliers among the SPARTACUS simulations (shown by red dots in Fig. 8c); indeed, if we were to look at the median rather than the mean of the 65 cases then it would suggest instead that SPARTACUS tends to underestimate the 3D effect on atmospheric absorption.

To investigate the factors that influence the nature of 3D radiative transfer in individual cases, and the fidelity with which they are captured by SPARTACUS, we analyze the radiation fields for the four contrasting GEM scenes depicted in Fig. 9. Vertical profiles of the four main inputs to SPARTACUS are shown in

Figs. 10a–d. Atlantic case 6 consists of cumulus clouds with some vertical development; the small effective cloud scale of $S \approx 1$ km, and hence large cloud-side area, leads to significant shortwave side illumination, with Fig. 10e showing a $7 \text{ W m}^{-2}$ increase in the reflectance of the scene at $\theta_0 = 75°$ (a cooling effect). The dependence on $\theta_0$ is well captured by SPARTACUS, including the change in sign to a $7 \text{ W m}^{-2}$ warming at $\theta_0 = 0°$, or a –9% change to CRE (see the caption of Fig. 10 for the total CRE values for overhead sun), but the similarity between the various entrapment configurations highlights that in this case the 3D effect for high sun is mainly due to downward escape.

The three remaining scenes, by contrast, appear to be dominated by entrapment. Atlantic case 14 contains deep frontal cloud with considerable small-scale structure. The zero-entrapment simulation in Fig. 10f shows the significant cooling effect of side illumination, but the explicit-entrapment simulation shows that this is overwhelmed by entrapment and indeed the net warming by 3D effects is up to $39 \text{ W m}^{-2}$ (a –8% change to CRE).
<figure><img src="image_11.png" /><figcaption>FIG. 8. Difference between simulations with and without 3D effects in (a) TOA CRE, (b) surface CRE, and (c) atmospheric absorption. The thick black and red lines show the average of the 65 GEM scenes for the Monte Carlo model and SPARTACUS in its preferred configuration of explicit entrapment with $\zeta = 0$. The box plots represent the corresponding median and interquartile range of the 65 scenes, with the whiskers representing the 5th and 95th percentiles. The other thick lines show averages of different SPARTACUS configurations.</figcaption></figure>

A key factor is the large vertical extent of the cloud, which means that radiation passing down through the gaps in the clouds can travel a large distance horizontally before being reflected back up to its original level, increasing the trapping. Atlantic case 32 contains much more homogeneous and overcast boundary layer cloud. The zero-entrapment simulation has a 3D radiative effect of less than $1 \text{ W m}^{-2}$, confirming that cloudside effects are weak. With entrapment included, the 3D effect is a warming of up to $6-7 \text{ W m}^{-2}$, with good agreement between SPARTACUS and Monte Carlo. In absolute terms this effect is significant, but this scene is the most reflective of the four and in relative terms it is only a -1% change to CRE. Additional Tripleclouds

and SPARTACUS calculations in which the in-cloud heterogeneity is removed (i.e., setting FSD = 0) lead to the 3D effect almost entirely disappearing, which suggest that it is due to trapping associated with cloud heterogeneity, similar to one of the mechanisms proposed by Várnai and Davies (1999). Finally, Pacific case 25 consists of remnants of deep convection including anvils with $S \approx 10$ km. Again the entrapment mechanism appears to dominate.

## 7. Conclusions

Cloud scenes have varied and complex structures, and consequently it can be very challenging to understand

<figure><img src="image_12.png" /><figcaption>FIG. 9. Three-dimensional visualizations of four contrasting GEM scenes, where the grid axes are marked in kilometers. The results of 3D radiation calculations on these scenes are shown in Fig. 10.</figcaption></figure>
<figure><img src="image_13.png" /><figcaption>FIG. 10. (a)–(d) The four main cloud geometry parameters used by SPARTACUS for the four GEM scenes shown in Fig. 9, where cloud fraction in gray varies between 0 and 1 across the horizontal axis, and the cloud and heterogeneity scales are defined by (29) and (32). (e)–(h) The corresponding calculations of the difference in TOA shortwave CRE between 3D and 1D calculations, where the black lines show Monte Carlo results and the colored lines show various configurations of SPARTACUS. For reference, the overhead-sun CREs for the four scenes according to the 3D Monte Carlo model were −79, −476, −832, and −184 W m⁻², respectively.</figcaption></figure>

the magnitude and even sign of the differences between radiation calculations with and without horizontal transport. The simplest mechanism to understand shortwave 3D radiative transfer is side illumination, which enhances cloud reflectance. This has led many previous studies to focus on cloud types with a relatively large cloud-side area such as cumulus (Benner and Evans 2001; Pincus et al. 2005) and aircraft contrails (Gounou and Hogan 2007). However, Barker et al. (2015) analyzed a much more varied and representative set of cloud fields and found that shortwave 3D transport tends to reduce the reflectance of clouds overall and hence has a warming effect on the Earth system. In this paper we propose the mechanism of entrapment to explain this behavior. Entrapment is similar to one of the mechanisms suggested by Várnai and Davies (1999) for single-layer cloud scenes, but an important insight is that it need not involve transport through cloud sides. It tends to be strongest in deep, multilayer scenes, which are common in reality but have tended to be ignored in previous case studies, presumably due to them being regarded as too complex to interpret.

We have described modifications to the shortwave SPARTACUS solver of Hogan et al. (2016) to incorporate

an explicit calculation of entrapment, making use of the effective cloud scale variable already provided as input to SPARTACUS. This involves a novel method to estimate the mean horizontal distance traveled by reflected radiation, something that could be useful in other contexts, for example in determining when the radiation scheme of a cloud-resolving model ought to represent lateral exchange of radiation between grid boxes. We have also found it necessary to explicitly represent the fractal dimension of cloud perimeters.

Evaluation against Monte Carlo calculations on 65 contrasting scenes from a cloud-resolving model reveals the new SPARTACUS scheme to be capable of predicting the “3D effect,” that is, the difference between cloud radiative effect computed with and without horizontal radiative transport, with a TOA bias of no more than 0.3 W m⁻² for all solar zenith angles, and skill in predicting the dependence of the 3D effect on solar zenith angle in individual scenes. On average, 3D radiative effects tend to make these scenes less reflective (similar to the findings of Barker et al. 2015), implying that entrapment is a more important mechanism than side illumination. However, this result is highly dependent
on the realism of the clouds simulated by the cloud-resolving model; if real clouds were smaller, on average, than those used here then the side-illumination mechanism would be relatively more important.

The modified SPARTACUS is now an option in the ecRad radiation scheme (Hogan and Bozzo 2018) used in the ECMWF model. Hogan and Bozzo (2018) reported the original SPARTACUS with maximum entrapment to be 3.5 times slower than Tripleclouds, and we find that explicit entrapment increases this to around 4.5. While too costly to use operationally, it is fast enough to use for research purposes. The next step will be to use this validated tool to estimate the global impact of 3D radiative transfer, not just in the shortwave but also in the longwave (Schäfer et al. 2016).

*Acknowledgments.* We thank Zhipeng Qu for performing the GEM model simulations. NV acknowledges support from the Agence Nationale de la Recherche (Grant ANR-16-CE01-0010), and from the French Ministry of Higher Education, Research and Innovation through the doctoral school SDU2E of Université de Toulouse.

# APPENDIX A

## Improving Tripleclouds for Very Heterogeneous Scenes

SPARTACUS can be thought of as a 3D extension of the 1D Tripleclouds solver of Shonk and Hogan (2008), so before evaluating SPARTACUS against 3D Monte Carlo calculations we need to be sure that Tripleclouds agrees well with the corresponding ICA calculations. Tripleclouds takes as input a profile of fractional standard deviation (FSD), which is the standard deviation of in-cloud water content or extinction coefficient divided by the in-cloud mean value. Tripleclouds divides the cloud in each layer into two regions (denoted *b* and *c*) of different extinction coefficient. Shonk and Hogan (2008) reported that for FSD up to 2, predicted irradiances agreed best with ICA if the two cloudy regions had equal area, region *b* used the 16th percentile of the full extinction distribution, and the extinction of region *c* was chosen so as to conserve the layer-mean extinction. The implementation of Tripleclouds in ecRad (Hogan and Bozzo 2018) includes the option to represent either a lognormal distribution of optical depth, in which case the ratio of the 16th percentile to the in-cloud mean is given approximately by (44) of Hogan et al. (2016), or a gamma distribution (e.g., Barker et al. 1996) for which this ratio is approximated by

$$r^b = \exp(-\text{FSD} - \text{FSD}^2/2 - \text{FSD}^3/4). \quad (A1)$$

Conservation of mean extinction coefficient requires the ratio of the extinction of region *c* to the in-cloud mean to be $r^c = 2 - r^b$.

Comparison of Tripleclouds to ICA calculations on the scenes described in section 6 revealed the gamma distribution to perform best, but even then Tripleclouds tended to overestimate scene reflectance for the more heterogeneous scenes, some of which have FSD values up to 4. In Fig. A1a we have repeated the analysis of Shonk and Hogan (2008) but for a gamma rather than a lognormal distribution, and considered larger values of FSD. A substantial albedo overestimate is apparent for FSD > 2. The problem arises because for large FSD, $r^b$ tends to 0 and $r^c$ to 2. Since the two cloudy regions have equal area, the actual fractional standard deviation of the Tripleclouds representation of the in-cloud extinction distribution tends to 1, which may be much less than the FSD of the gamma distribution being approximated.

The solution we propose to overcome this problem is twofold. First, a new $r^b$ is defined with a lower limit of 0.025:

$$r_{\text{new}}^b = 0.975r^b + 0.025. \quad (A2)$$

Second, for large values of FSD we increase the fractional area of region *b* and correspondingly reduce that of region *c*: for FSD in the range 1.5–3.75, the fraction of the cloudy area occupied by region *b* increases linearly from 0.5 to 0.9, while outside this range it is capped at 0.5 or 0.9. The extinction of region *c* is still chosen to conserve the layer-mean value. Figure A1b shows that these changes virtually eliminate the albedo bias up to an FSD of 4. This solution has been implemented in both the Tripleclouds and SPARTACUS solvers of ecRad, and is used in section 6.

# APPENDIX B

## List of Symbols

The following list includes symbols used in more than one equation in sections 3–5, and “PP” indicates a variable from section 4a where a plane-parallel atmosphere has been assumed.

* **A** above $i-1/2$ Diffuse albedo of entire atmosphere and surface below interface *i* − 1/2, with matrix elements configured for regions in the layer *above* the interface (region *i* − 1)

* **A** above $i-1/2$ Diffuse albedo of entire atmosphere and surface below interface *i* − 1/2 (PP)

* **c**<sup>*k*</sup> Fraction of layer occupied by region *k*
<figure><img src="image_14.png" /><figcaption>FIG. A1. Albedo bias of the Tripleclouds method for a gamma distribution of cloud optical depths with increasing fractional standard deviation indicated on the vertical axis: (a) equal-area cloudy regions with the optical depth of the first cloudy region given by the 16th percentile of the gamma distribution and (b) the new method in which the area of the two cloudy regions can be different. Following Shonk and Hogan (2008), this analysis approximates the relationship between albedo $\alpha$ and optical depth $\delta$ as $\alpha = 0.2 + 0.525\delta/(\delta + 3.5)$.</figcaption></figure>

<table><tr><td><b>D</b><sub>above i-1/2</sub></td><td>Same as <b>A</b><sub>above i-1/2</sub>, but for direct radiation</td><td><b>T</b><sub>i</sub></td><td>Transmittance of layer <i>i</i> to diffuse radiation (PP)</td></tr><tr><td><b>D</b><sub>i-1/2</sub></td><td>Direct albedo of entire atmosphere and surface below interface <i>i</i> - 1/2 (PP)</td><td><b>u</b></td><td>Vector of upwelling irradiances in each region at a particular height</td></tr><tr><td><b>E</b><sub>i</sub></td><td>Fraction of direct radiation penetrating layer <i>i</i> without being scattered (PP)</td><td><b>U</b><sub>i-1/2</sub></td><td>Upward overlap matrix expressing how upwelling irradiances in each region just below interface <i>i</i> - 1/2 are transported into the regions just above</td></tr><tr><td><b>f</b><sup>kl</sup></td><td>Rate at which radiation passes from sub region <i>k</i> to <i>l</i>, per unit increase in horizontal distance traveled</td><td><b>v</b></td><td>Vector of downwelling diffuse irradiances in each region at a particular height</td></tr><tr><td><b>L</b><sup>kl</sup></td><td>Length of interface between regions <i>k</i> and <i>l</i> normalized by the area of the domain</td><td><b>V</b><sub>i-1/2</sub></td><td>Downward overlap matrix expressing how downwelling irradiances in each region just above interface <i>i</i> - 1/2 are transported into the regions just below</td></tr><tr><td><b>L̂</b><sup>kl</sup></td><td>Effective normalized interface length, accounting for the fractal nature of clouds</td><td><b>w</b><sup>j</sup><sub>i-1/2</sub></td><td>Vector expressing the fraction of radiation in region <i>j</i> of layer <i>i</i> that is beneath each region of layer <i>i</i> - 1</td></tr><tr><td><b>L</b><sup>j</sup><sub>i-1/2</sub></td><td>Matrix expressing how much radiation entering region <i>j</i> of layer <i>i</i> from the top and subsequently reflected back up, is exchanged between the various regions in layer <i>i</i> - 1</td><td><b>x</b><sub>i-1/2</sub></td><td>Mean horizontal distance traveled by reflected diffuse radiation below interface <i>i</i> - 1/2 (PP)</td></tr><tr><td><b>R</b><sub>i</sub></td><td>Diffuse reflectance of layer <i>i</i> (PP)</td><td><b>ŷ</b><sub>i</sub></td><td>Mean horizontal distance traveled by diffuse radiation reflected by layer <i>i</i> (PP)</td></tr><tr><td><b>s</b></td><td>Vector of downwelling direct irradiances in each region at a particular height</td><td><b>y</b><sub>i-1/2</sub></td><td>Mean horizontal distance traveled by reflected direct radiation below interface <i>i</i> - 1/2 (PP)</td></tr><tr><td><b>S</b><sub>i</sub><sup>-</sup></td><td>Reflectance of layer <i>i</i> to direct radiation (PP)</td><td></td><td></td></tr><tr><td><b>S</b><sub>i</sub><sup>+</sup></td><td>Fraction of direct radiation that penetrates layer <i>i</i> and is scattered on the way (PP)</td><td></td><td></td></tr><tr><td><b>S</b>, <b>S</b><sup>het</sup></td><td>Effective cloud scale, cloud heterogeneity scale</td><td></td><td></td></tr></table>
$\hat{y}_i$ Mean horizontal distance traveled by direct radiation reflected by layer *i* (PP)

$\Gamma$ Matrix expressing the rates of radiation exchange between the subregions of a region

$\Delta x_i$ Horizontal distance traveled by diffuse radiation passing through layer *i*

$\Delta y_i$ Horizontal distance traveled by direct radiation passing through layer *i*

$\Delta z_i$ Thickness of layer *i*

$\zeta$ Overhang factor

$\theta_0$ Solar zenith angle

REFERENCES

Barker, H. W., and J. A. Davies, 1992: Solar radiative fluxes for broken cloud fields above reflecting surfaces. *J. Atmos. Sci.*, **49**, 749–761, https://doi.org/10.1175/1520-0469(1992)049:0749:SRFFBC>2.0.CO;2.

—, B. A. Wielicki, and L. Parker, 1996: A parameterization for computing grid-averaged solar fluxes for inhomogeneous marine boundary layer clouds. Part II: Validation using satellite data. *J. Atmos. Sci.*, **53**, 2304–2316, https://doi.org/10.1175/1520-0469(1996)053:2304:APFCGA>2.0.CO;2.

—, R. K. Goldstein, and D. E. Stevens, 2003: Monte Carlo simulation of solar reflectances for cloudy atmospheres. *J. Atmos. Sci.*, **60**, 1881–1894, https://doi.org/10.1175/1520-0469(2003)060:1881:MCSOSR>2.0.CO;2.

—, J. N. S. Cole, J. Li, B. Yi, and P. Yang, 2015: Estimation of errors in two-stream approximations of the solar radiative transfer equation for cloudy-sky conditions. *J. Atmos. Sci.*, **72**, 4053–4074, https://doi.org/10.1175/JAS-D-15-0033.1.

—, —, —, and K. von Salzen, 2016: A parametrization of 3-D subgrid-scale clouds for conventional GCMs: Assessment using A-Train satellite data and solar radiative transfer characteristics. *J. Adv. Model. Earth Syst.*, **8**, 566–597, https://doi.org/10.1002/2015MS000601.

Benner, T. C., and K. F. Evans, 2001: Three-dimensional solar radiative transfer in small tropical cumulus fields derived from high-resolution imagery. *J. Geophys. Res.*, **106**, 14 975–14 984, https://doi.org/10.1029/2001JD900158.

Brown, A. R., and Coauthors, 2002: Large-eddy simulation of the diurnal cycle of shallow cumulus convection over land. *Quart. J. Roy. Meteor. Soc.*, **128**, 1075–1093, https://doi.org/10.1256/003590002320373210.

Cahalan, R. F., and J. H. Joseph, 1989: Fractal statistics of cloud fields. *Mon. Wea. Rev.*, **117**, 261–272, https://doi.org/10.1175/1520-0493(1989)117:0261:FSOCF>2.0.CO;2.

Di Giuseppe, F., and A. M. Tompkins, 2015: Generalizing cloud overlap treatment to include the effect of wind shear. *J. Atmos. Sci.*, **72**, 2865–2876, https://doi.org/10.1175/JAS-D-14-0277.1.

Edwards, J. M., and A. Slingo, 1996: Studies with a flexible new radiation code: 1. Choosing a configuration for a large scale model. *Quart. J. Roy. Meteor. Soc.*, **122**, 689–719, https://doi.org/10.1002/qj.49712253107.

Girard, C., and Coauthors, 2014: Staggered vertical discretization of the Canadian Environmental Multiscale (GEM) model using a coordinate of the log-hydrostatic-pressure type. *Mon. Wea. Rev.*, **142**, 1183–1196, https://doi.org/10.1175/MWR-D-13-00255.1.

Gotoh, K., and Y. Fujii, 1998: A fractal dimensional analysis on the cloud shape parameters of cumulus over land. *J. Appl. Meteor.*, **37**, 1283–1292, https://doi.org/10.1175/1520-0450(1998)037:1283:AFDAOT>2.0.CO;2.

Gounou, A., and R. J. Hogan, 2007: A sensitivity study of the effect of horizontal photon transport on the radiative forcing of contrails. *J. Atmos. Sci.*, **64**, 1706–1716, https://doi.org/10.1175/JAS3915.1.

Hapke, B., D. DiMucci, R. Nelson, and W. Smythe, 1996: The cause of the hot spot in vegetation canopies in soils: Shadow-hiding versus coherent backscatter. *Remote Sens. Environ.*, **58**, 63–68, https://doi.org/10.1016/0034-4257(95)00257-X.

Hill, P. G., C. J. Morcrette, and I. A. Boutle, 2015: A regime-dependent parametrization of subgrid-scale cloud water content variability. *Quart. J. Roy. Meteor. Soc.*, **141**, 1975–1986, https://doi.org/10.1002/qj.2506.

Hogan, R. J., and A. Battaglia, 2008: Fast lidar and radar multiple-scattering models. Part II: Wide-angle scattering using the time-dependent two-stream approximation. *J. Atmos. Sci.*, **65**, 3636–3651, https://doi.org/10.1175/2008JAS2643.1.

—, and J. K. P. Shonk, 2013: Incorporating the effects of 3D radiative transfer in the presence of clouds into two-stream multilayer radiation schemes. *J. Atmos. Sci.*, **70**, 708–724, https://doi.org/10.1175/JAS-D-12-041.1.

—, and A. Bozzo, 2018: A flexible and efficient radiation scheme for the ECMWF model. *J. Adv. Model. Earth Syst.*, **10**, 1990–2008, https://doi.org/10.1029/2018MS001364.

—, S. A. K. Schäfer, C. Klinger, J.-C. Chiu, and B. Mayer, 2016: Representing 3D cloud-radiation effects in two-stream schemes: 2. Matrix formulation and broadband evaluation. *J. Geophys. Res. Atmos.*, **121**, 8583–8599, https://doi.org/10.1002/2016JD024875.

—, T. Quaife, and R. Braghiere, 2018: Fast matrix treatment of 3-D radiative transfer in vegetation canopies: SPARTACUS-Vegetation 1.1. *Geosci. Model Dev.*, **11**, 339–350, https://doi.org/10.5194/gmd-11-339-2018.

Iacono, M. J., J. S. Delamere, E. J. Mlawer, M. W. Shephard, S. A. Clough, and W. D. Collins, 2008: Radiative forcing by long-lived greenhouse gases: Calculations with the AER radiative transfer models. *J. Geophys. Res.*, **113**, D13103, https://doi.org/10.1029/2008JD009944.

Lacis, A. A., and J. E. Hansen, 1974: A parameterization for the absorption of solar radiation in the Earth's atmosphere. *J. Atmos. Sci.*, **31**, 118–133, https://doi.org/10.1175/1520-0469(1974)031:0118:APFTAO>2.0.CO;2.

Leroyer, S., S. Bélair, S. Husain, and J. Mailhot, 2014: Sub-kilometer numerical weather prediction in an urban coastal area: A case study over the Vancouver metropolitan area. *J. Appl. Meteor. Climatol.*, **53**, 1433–1453, https://doi.org/10.1175/JAMC-D-13-0202.1.

Milbrandt, J. A., and M. K. Yau, 2005: A multimoment bulk microphysics parameterization. Part II: A proposed three-moment closure and scheme description. *J. Atmos. Sci.*, **62**, 3065–3081, https://doi.org/10.1175/JAS3535.1.

Pincus, R., H. W. Barker, and J.-J. Morcrette, 2003: A fast, flexible, approximate technique for computing radiative transfer in inhomogeneous cloud fields. *J. Geophys. Res.*, **108**, 4376, https://doi.org/10.1029/2002JD003322.

—, C. Hannay, and K. F. Evans, 2005: The accuracy of determining three-dimensional radiative transfer effects in cumulus clouds using ground-based profiling instruments. *J. Atmos. Sci.*, **62**, 2284–2293, https://doi.org/10.1175/JAS3464.1.
Schäfer, S. A. K., R. J. Hogan, C. Klinger, and B. Mayer, 2016: Representing 3D effects in two-stream radiation schemes: 1. Longwave considerations and effective cloud edge length. *J. Geophys. Res. Atmos.*, **121**, 8567–8582, https://doi.org/10.1002/2016JD024876.

Shonk, J. K. P., and R. J. Hogan, 2008: Tripleclouds: An efficient method for representing horizontal cloud inhomogeneity in 1D radiation schemes by using three regions at each height. *J. Climate*, **21**, 2352–2370, https://doi.org/10.1175/2007JCLI1940.1.

—, —, J. M. Edwards, and G. G. Mace, 2010: Effect of improving representation of horizontal and vertical cloud structure on the Earth’s global radiation budget. Part I: Review and parametrization. *Quart. J. Roy. Meteor. Soc.*, **136**, 1191–1204, https://doi.org/10.1002/qj.647.

Várnai, T., and R. Davies, 1999: Effects of cloud heterogeneities on shortwave radiation: Comparison of cloud-top variability and internal heterogeneity. *J. Atmos. Sci.*, **56**, 4206–4224,

https://doi.org/10.1175/1520-0469(1999)056<4206:EOCHOS> 2.0.CO;2.

Welch, R. M., and B. A. Wielicki, 1984: Stratocumulus cloud field reflected fluxes: The effect of cloud shape. *J. Atmos. Sci.*, **41**, 3085–3103, https://doi.org/10.1175/1520-0469(1984)041<3085: SCFRFT>2.0.CO;2.

Wissmeier, U., R. Buras, and B. Mayer, 2013: paNTICA: A fast 3D radiative transfer scheme to calculate surface solar irradiance for NWP and LES models. *J. Appl. Meteor. Climatol.*, **52**, 1698–1715, https://doi.org/10.1175/JAMC-D-12-0227.1.

Wood, R., and P. R. Field, 2011: The distribution of cloud horizontal sizes. *J. Climate*, **24**, 4800–4816, https://doi.org/10.1175/2011JCLI4056.1.

Yi, B., P. Yang, B. A. Baum, T. L'Ecuyer, L. Oreopoulos, E. J. Mlawer, A. J. Heymsfield, and K.-K. Liou, 2013: Influence of ice particle surface roughening on the global cloud radiative effect. *J. Atmos. Sci.*, **70**, 2794–2807, https://doi.org/10.1175/JAS-D-13-020.1.