Towards an improved treatment of cloud-radiation interaction in weather and climate models: exploring the potential of the Tripleclouds method for various cloud types using libRadtran 2.0.4

. Although the (cid:58)(cid:58)(cid:58) The (cid:58) representation of unresolved clouds in radiation schemes of coarse-resolution weather and climate models has progressed noticeably over the past years. (cid:58)(cid:58)(cid:58)(cid:58)(cid:58)(cid:58)(cid:58)(cid:58)(cid:58)(cid:58)(cid:58) Nevertheless, a lot of room remains for improvement, as the current picture is by no means complete. The main objective of the present study is to advance the cloud-radiation interaction parameterization, focusing on the issues related to model misrepresentation of cloud horizontal inhomogeneity. This subject is addressed with the state-of-the-art Tripleclouds radiative solver, the fundamental feature of which is the inclusion of the opti-5 cally thicker and thinner cloud fraction, where the thicker is associated with the presence of convective updraft elements. The research challenge is to optimally set the pair of cloud condensates characterizing the two cloudy regions and the corresponding geometrical split of layer cloudiness. A diverse cloud ﬁeld data set was collected for the analysis, comprising case studies of stratocumulus, cirrus and cumulonimbus. The primary goal is to assess the validity of global cloud variability estimate along with various condensate distribution assumptions. More sophisticated parameterizations are subsequently explored, optimiz-10 ing the treatment of overcast as well as extremely heterogeneous cloudiness. The radiative diagnostics including atmospheric heating rate and net surface ﬂux are consistently studied using the Tripleclouds method, evaluated against a three-dimensional radiation computation. The performance of Tripleclouds mostly signiﬁcantly surpasses the conventional calculation on horizontally homogeneous cloudiness. The effect of horizontal photon transport is further quantiﬁed. The overall conclusions are intrinsically different for each particular cloud type, encouraging endeavors to enhance the use of cloud regime dependent 15 methodologies in next-generation atmospheric models. This study highlighting the Tripleclouds potential for three essential cloud types signiﬁes the need for more research examining a broader spectrum of cloud morphologies

Authors' reply to reviewer Robin Hogan of the paper: "Towards an improved treatment of cloud-radiation interaction in weather and climate models: exploring the potential of the Tripleclouds method for various cloud types using libRadtran 2.0.4",submitted for publication in GMD by Nina Črnivec and Bernhard Mayer.We thank the reviewer Robin Hogan for good assessment of our work and many valuable comments and suggestions, which helped improving the quality of our original manuscript.Below please find reviewer's comments in blue and our reply in black.Changes in the revised paper are marked with quotation marks and additional indent.The line numbers refer to those in the original manuscript.This paper presents an interesting evaluation of the Tripleclouds (TC) method for representing cloud structure in a radiation scheme using benchmark Monte Carlo calculations and three contrasting, realistic 3D cloud scenes.While the McICA method is a more common way to treat cloud structure in operational models, Tripleclouds is attractive because it is free from stochastic noise.The paper goes into some depth in testing and optimizing the various ways in which Tripleclouds can be configured to reduce the errors.However, I have some major concerns that should be addressed before this paper can be published.The paper is generally well written otherwise.

MAJOR COMMENTS
1.The use of maximum-random overlap is well behind the state-of-the-art and will have led to appreciable biases in your GCM and TC schemes.The cirrus case is the most obvious: Hogan and Kew (2005, Table 2) found that for this exact case, a GCM-type radiation scheme with maximum-random overlap underestimated the TOA cloud radiative forcing by 42% in the longwave and 48% in the shortwave, compared to the same GCM-type scheme but with the "true" overlap.I expect significant overlap biases to be present in the other two scenes, and indeed I can only explain your TC flux biases with respect to ICA as being due to overlap being too maximal.To test this, simply compute the total cloud cover predicted from the cloud-fraction profile and your maximum-random overlap assumption, and compare it to the actual total cloud cover of the scene.Note that even the overcast stratocumulus case can have an overlap bias because of the overlapping of subgrid cloud structures, which ought to be represented by non-maximum overlap of the high-LWC and low-LWC cloudy regions (some schemes use a vertical decorrelation scale for cloud heterogeneities of half that used for cloud boundaries, i.e. around 1 km).I don't believe it would be too complicated to incorporate realistic overlap: this can be done using overlap matrices introduced by Shonk et al. (2010), and an example implementation is in the open-source ecRad package -see the radiation_overlap.F90 source file on GitHub.The capability has been in older TC implementations for longer -for example, Shonk and Hogan (2010) evaluated the impact of horizontal heterogeneity and vertical overlap using an implementation of Tripleclouds incorporating exponential-random overlap in the Met Office radiation scheme, and stressed the importance of making both horizontal heterogeneity and vertical overlap as realistic as possible.
We agree that the maximum-random overlap is not the state-of-the-art approach, thank you for pointing this out.We have thus removed the term "state-of-the-art" when describing our Tripleclouds radiative solver (line 5, caption of Fig. 1, line 91).Unfortunately, we did not manage to incorporate more realistic overlap within the scope of this study, but we plan to generalize the overlap rules in the future.We have extended the final section "Summary and conclusions" as follows: "Furthermore, in order to properly consider clouds in sheared environments, the vertical overlap rules in the present TC implementation have to be generalized."In order to eliminate the effects of vertical wind shear in the present study, however, we have intentionally selected cloud cases, where the realistic overlap should be close to the assumed maximum-random overlap.
Firstly, we would like to emphasize that we did not use the exact same cirrus as Hogan and Kew (2005), although it is true that both cloud fields relate to the same case study (June, 24 th , 1999).In particular, Hogan and Kew (2005) used coarser horizontal grid spacing of 1.56 km.We adopted the cirrus data presented in doctoral dissertation of Sophia Schäfer, where the simulation of Hogan and Kew (2005) was rerun with a higher resolution (50 m horizontal grid spacing) to enable proficient geometry analysis.Finally, we have additionally smeared the latter case onto grid with 100 m horizontal grid spacing to facilitate the radiation simulations.We have additionally clarified these issues in Section 2.1.2(footnote): "It should be noted, however, that the studies of Hogan and Kew (2005) and Zhong et al. (2008) use coarser horizontal grid spacing (1.56 km).We adopted the cirrus data from Schäfer (2016), where the simulation of Hogan and Kew (2005) had been rerun with higher resolution (horizontal grid spacing of 50 m), whereby we eventually smeared the data onto the grid with horizontal grid spacing of 100 m to facilitate the Monte Carlo radiative transfer simulations."Furthermore, we probably also use different vertical resolution, which should also affect the difference in vertical overlap errors.We use constant vertical grid spacing of 109 m, whereas Hogan and Kew (2005) state: "This yields a horizontal resolution of 1.56 km and between 45 and 110 m in the vertical (depending on the case)..."(?).Finally, the extent of the calculation domains in Hogan and Kew (2005) and Schäfer (2016) differs as well.In summary, the two cirrus cases have different 3-D geometry, so that the overlap errors do not necessarily match.
In order to gain further insight into vertical overlap issues (as you advised), we have computed the actual total cloud cover of the 3-D scene and compared it to the maximum of the cloud-fraction vertical profile.The latter namely implies the total cloud cover corresponding to the maximum overlap assumption (if the cloud layer is vertically continuous without cloud-free regions, as is the case throughout our study, then the maximumrandom overlap simplifies to the maximum overlap).For our cirrus case, for example, we found the following: There is indeed discrepancy between the actual total cloud cover of the 3-D scene (which is about 76 %) and maximum layer cloud fraction (which is about 47 %; see Fig. 3).However, it should be emphasized that in our maximum-random overlap implementation, we employ pairwise overlap (in this manner the matrix problem is faster to solve).This means that the total cloud cover predicted by our overlap scheme is generally larger than the maximal layer cloud fraction and is therefore closer to the actual total cloud cover!In other words, the pairwise overlap "luckily" relaxes the stiffness of the maximum assumption, acting to reduce the error due to true overlap not being maximal.
Regarding your comment: "..., and indeed I can only explain your TC flux biases with respect to ICA as being due to overlap being too maximal."-We don't see this so clearly (?).In general, the TC biases with respect to the ICA are due to several reasons: -Overlap being too maximal in the TC compared to realistic overlap, which is taken into account in the ICA (but hopefully the error is small in our case); -Inhomogeneity being misrepresented in the TC (either underestimated or overestimated) compared to the full realistic inhomogeneity, which is taken into account in the ICA; -Effective (partial) treatment of 3-D radiative effects in the TC, which are entirely neglected in the ICA.Please see also our answer to comment 18, which further elucidates these issues.
2. When comparing TC results to the benchmark, a number of errors are compounded and there is insufficient effort to separate them out for the reader: (1) the overlap parameterization is biased (see above); (2) the split percentile of 50% cannot adequately represent large FSDs; (3) a value of FSD=0.75 is used at various points which is different (often much lower) than the "truth"; and (4) TC is benchmarked against at 3D model when it makes no attempt to represent 3D effects.This leads to tuning of one parameterization to fix a problem caused by another (see comments 18 and 19).Wouldn't it be more satisfactory to evaluate TC against ICA when all the inputs are correct, in order to identify the intrinsic error in TC, then the impact of not knowing the exact overlap or FSD, or not representing 3D effects, can be quantified separately?Indeed, this is the approach taken by Hogan et al. (2019): their Fig.7a-c shows that when the correct overlap and FSD is used, the shortwave bias against ICA over 65 scenes is less than 3 W m-2.This was then a firm foundation for them to look at the representation of 3D effects, whereas in your case it would be a firm foundation to look primarily at errors due to parameterizing FSD (although you also have some interesting 3D results).
We agree that when we compare the TC results to the benchmark, a number of errors are compounded.The ecRad package (and the study of Hogan et al., 2019) contains the Tripleclouds radiative solver (Shonk and Hogan, 2008) merged with SPARTACUS (Schäfer et al., 2016), where you have the possibility to consider: horizontal cloud inhomogeneity + general (true) vertical overlap + proper 3-D effects (SPARTACUS).Our TC scheme, on the other hand, can only account for horizontal inhomogeneity (and to a minor extent for 3-D effects), therefore we can not easily separate all different effects out for the reader.We indeed tried to discuss the contribution of the various effects to the extent that we could (horizontal inhomogeneity, 3-D effects).
Please see our answer to comments 18 and 19, where we explain how we have extended the discussion by mentioning the vertical overlap issues, the problem of compensating errors etc. (in order to better explain the various error sources).
Nevertheless, our primary aim was to explore various TC configurations, which can in practice be used in GCMs.To that end, the 3-D experiment is the proper benchmark, because it is the best proximity to the real world (and not the ICA).Please see also our answer to comment 13, which discusses the issues related to FSD parameterization.
3. The proposed solution to the high-FSD problem in section 4 should be compared to the solution to the same problem presented by Hogan et al. (2019, appendix).Their solution is derived from theoretical distributions, rather than by trying to minimize an error against a benchmark calculation in which several other errors (notably in overlap) are also present.The structure of the present manuscript is a little frustrating -section 3 contains some puzzling errors that are only addressed, or even properly mentioned, when the reader gets to section 4. Why not flag up the need to address the problem sooner in the paper?Thank you for this suggestion.We have reconstructed the manuscript as you advised.We have thus introduced the various Tripleclouds experiments (baseline as well as optimizations for both overcast and highly heterogeneous cloudiness) already in Section 2. We have simultaneously reconstructed the "Results and discussion" section.Unfortunately we can not test the full solution to the high-FSD problem presented in Hogan et al. (2019), because in our Tripleclouds solver the cloud fraction scaling factor (which determines the geometrical splitting of cloud fraction in two parts) is implemented as a constant (i.e., not height-dependent).Hogan et al. (2019) presents the solution, where the "cloud fraction scaling factor" is a function of FSD (with the true high FSD being height-dependent).Once we manage to generalize the vertical overlap rules as you suggested in the first place (for arbitrary height-dependent cloud fraction scaling factors), we will be able to test the solution of Hogan et al. (2019) for high FSDs on our cloud data (which would be interesting to do).SPECIFIC COMMENTS 4. Abstract line 6: while the optically thicker part could be used to represent convection, most clouds are not convection and the use of a thicker and thinner parts are simply a first-order approximation to the horizontal distribution of optical depth that is found in stratiform clouds.
We agree.We have thus shortened the sentence to: "This subject is addressed with the Tripleclouds radiative solver, the fundamental feature of which is the inclusion of the optically thicker and thinner cloud fraction." 5. Figure 1: Caption should stress that this is a schematic; the vertical resolution shown here is coarser than any operational model.
Thank you for pointing this out.We have stressed these issues: "Note that the schematics are illustrative and that operational models employ finer vertical resolution." 6. Introduction: This seems unnecessarily long.The primary purpose of the paper concerns testing the Tripleclouds scheme for representing horizontal cloud heterogeneity in a 1D radiative transfer context, for which the appropriate benchmark is the Independent Column Approximation.Many of the references and discussion concern 3D radiative effects, which seems not so central to the topic of this paper; a little shortening would therefore seem in order.The introduction should mention the McICA scheme of Pincus et al. (2003), which is much more commonly used to represent cloud structure than Tripleclouds in current global models.Figure 1, panel 2, could just as easily be used to illustrate the McICA scheme as a cloud-resolving model.
We have shortened the Introduction as you advised.We have omitted a few sentences, redundant words and a few references (concerning 3-D radiative effects and some other).Please see the manuscript file with marked changes.We have also mentioned the McICA scheme of Pincus et al. (2003) -thank you for pointing this out.The added sentence is as follows: "While an alternative technique known as the McICA (Pincus et al., 2003) is currently operationally employed in the majority of coarse-resolution models, the TC scheme is attractive because it is free from stochastic noise."In addition, we would like to let you know that we have thoroughly described/compared the McICA and the Tripleclouds scheme (in terms of computational efficiency, stochastic noise and related issues) in the Introduction of the parent paper Črnivec and Mayer (2020).Please see also Fig. 1(e) of Črnivec and Mayer (2020), which nicely illustrates the McICA algorithm.We didn't want to exactly repeat all these contents in the present paper.Finally, we have retained the 3-D experiment as the benchmark, since it is the best proximity to the real world.7. Line 120: Define FSD here, saying particularly that it is the standard deviation divided by the mean, in both cases considering only the non-zero water content values in the horizontal LWC distribution.
Thank you for pointing this out.The FSD parameter was originally defined in Section 2, but we agree that it is better to introduce the definition earlier within Introduction.New formulation: "The latter is conveniently defined by the fractional standard deviation (FSD) of cloud condensate as well as the shape of condensate distribution.The parameter FSD (e.g., of LWC) is defined as the standard deviation (σLWC) divided by the mean (LWC), whereby only the non-zero values in the horizontal LWC distribution are considered." 8. Figure 2: the linear colour scale is not really suitable for the cirrus cloud since it is entirely white for optical depths up to 5, yet a significant fraction of the radiative impact of this cloud in the longwave will be from optical depths less than 5.
Thank you for this remark.We have improved the figure by applying the logarithmic colour scale.9. Line 155: Please say how the distributions were fitted (least-squares or fitting three of the moments of the distribution?)and I'm not sure what use a Gaussian distribution is, except as an excuse to use the 16th percentile, since it is unbounded on the lower end and so negative water contents are predicted.I'm also curious as to whether you can say that either lognormal or gamma are really better; Hogan and Illingworth (JAS 2003, Figs. 4-5) found that there was little to choose between them when comparing to real data, and often the gamma and lognormal were much closer to each other than either were to the noisy distributions of individual scenes.
We have extended the sentence as follows: "To gain further insight about the subgrid cloud variability, the theoretical distributions (Gaussian, gamma, lognormal) were fitted to the actual LWC distribution in each vertical layer, so that they have the same mean and standard deviation as the actual data."Thank you for pointing out some interesting results by Hogan and Illingworth (2003).It is true also in our case that lognormal and gamma distributions performed similarly well when fitted to the actual data in several layers.We have extended the discussion as follows: "The investigation revealed that the actual LWC distribution throughout the majority of the upper portion of the cloud, where radiative effect is maximal, is best approximated with the lognormal distribution (best fit in 5/8 of top layers), followed by the gamma distribution which performs similarly well."Nevertheless, there are also layers where one distribution clearly outperforms the other: see for example Fig. 14 of Črnivec and Mayer (2020) where the same analysis was applied for the cumulus cloud field (and the gamma distribution is clearly closer to the real distribution compared to the lognormal distribution).
The reason why we chose the three distributions (including the Gaussian distribution) is given later in the text (lines 223-225), where we also state that the lognormal and gamma distributions are more realistic.In addition, the comparison of the TC(FSD) experiment assuming Gaussian distribution with the TC(LP) experiment gives insight into the error due to FSD parameterization (e.g., global constant).10.Line 202: \bar{LWC} as defined here is not the layer-mean LWC, but the in-cloud mean LWC, i.e. ignoring the clear region.
Yes, we are aware of this (we thought it is clear that only cloudy pixels in the layer are considered when computing the average LWC).Nevertheless, we agree that the following sentence is better: "… is determined under conservation constraints of the in-cloud mean LWC (denoted as LWC):" 11.Line 234: a short further discussion is required -like the bottom row of your Fig. 3, Hogan et al. (2019, Fig. 10a-d) also found much larger FSD values than reported by Shonk et al. (2010).This is not possible to represent if the two cloudy regions are assumed to have the same area, but in the appendix to that paper they showed how an improved representation of large-FSD distributions could be achieved by making the optically thinner region occupy a larger area.
Thank you for this remark.We agree that the baseline TC configuration (where the two cloudy regions have the same area) is not best suitable for cloud scenes with large FSD.We have thus extended Section 2.2, by adding an extra subsection 2.2.3 entitled "Optimization for highly heterogeneous cloud scenarios" (dealing with large FSDs in the way that you described: making optically thinner region occupy a larger area).In our initial paper these issues were addressed later in Section 4 ("Parameter optimizations"), where we have also pointed out the study of Hogan et al. (2019).Nevertheless, we agree that it is better to raise these issues already in the Methodology section (Section 2).
12. Eq. 8: surely to be a bias, x and y should be averaged over more than one event?Otherwise it is an instantaneous error.Also it should be clarified whether x and y represent horizontal averages (e.g. of heating rate)... but are they also vertical averages of heating rate through the cloud layer?"Cloud-layer RMSE" is ambiguous -does the "layer" refer to model layer or the entire layer of cloud?
Correct, it is an instantaneous error.We have changed the terms "absolute bias" and "relative bias" to "absolute error" and "relative error" throughout the entire paper as well as on the figures.
In the case of 3-D and ICA experiments, x and y indeed represent horizontal averages.This is already stated in the text: "The 3-D and ICA experiment were both performed on the high-resolution cloud grid, with the result subsequently horizontally averaged across the domain."The GCM and TC experiments, on the other hand, are single-column experiments.We always show vertical profiles of heating rate, therefore we did not perform any vertical averaging.
The "cloud-layer RMSE" refers to the entire layer of cloud, as it is already explained in the text (line 280): "The cloud-layer RMSE denotes the RMSE evaluated throughout the vertical extent of the cloud layer of each particular cloud field case study." 13. Table 1 and text: I understand that the TC(FSD) method uses FSD=0.75, but it is unclear from either this table or the text whether the TC(LP) method uses the 16th percentile of the *true* in-cloud LWC distribution for each cases, or for an idealized (e.g.lognormal/gamma) distribution with an FSD of 0.75.If the former, then surely the main difference between the two methods is whether the true distribution is used or not, information that is not stated (at least not clearly).It would seem much more satisfactory in all cases to use the observed FSD in order that we are evaluating the intrinsic TC method, not the rather old and simple parameterization of Shonk et al. (2010).
The TC(LP) method uses the 16 th percentile of the true in-cloud LWC distribution.This is written at the beginning of Section 2.2 (lines 198-204), where we have summarized the original TC method of Shonk and Hogan (2008) which we apply in our study (term the "observed distribution" is used to denote the "true distribution").In order to clarify that we are employing this particular LP methodology as well as to clarify the difference between the LP and the FSD method (that you pointed out), we have reformulated the text in Section 2.3: "The Tripleclouds solver was employed in conjunction with the LP method based on the observed condensate distribution and the FSD method utilizing the distribution assumption in various configurations as outlined in Sect.2.2." To explain the TC(FSD) method once more: The TC(FSD) method in the baseline configuration utilizes an idealized distribution assumption (Gaussian/gamma/lognormal) together with FSD of 0.75.We intentionally didn't use the observed (actual) FSD based on high-resolution 3-D cloud data, because we assume that this exact value is not known in a GCM (in fact, this is our main motivation to test various FSD parameterizations).In other words -our primary aim was to answer the question: How to properly configure Tripleclouds for use in GCMs, where FSD parameterizations have to be used: e.g., global FSD estimate of Shonk et al. (2010), which is certainly rather old and simple; but we have also tested more recent and more sophisticated FSD parameterizations of Boutle et al. (2014) and Hill et al. (2015).The latter were initially introduced later in Section 4, whereas in a revised paper we have introduced them already in Section 2.
We could extend the discussion by comparing TC(FSD=parametrized) and TC(FSD=actual) in order to evaluate the TC error related to the FSD parameterization.We have performed these experiments (not shown in the paper), although the configuration of these experiments is actually not so straightforward -e.g., for the cumulonimbus.In the mixed-phase region of the cumulonimbus, namely, the actual FSD of liquid phase is very large (implying that simultaneously the asymmetrical cloud fraction splitting should be applied), whereas the actual FSD of ice is close to the global estimate (implying that the split percentile of 50 is adequate) -how can one properly deal with these issues (i.e., which split percentile should be used?-because the cloud fraction of liquid and ice phase are not separate parameters for TC solver...) is left for the future study.
The figure below shows the results of the "test experiments" for cumulonimbus, where the split percentile was simply held at 50: shown is the heating rate error in the baseline FSD experiment using global constant ("G") and using the actual FSD ("A") with various distribution assumptions.Whereas the heating rate error within the anvil is reduced when the actual FSD is used (and one could indeed estimate the error related to the FSD parameterization), the heating rate within the stratiform mixed-phase region is only slightly changed (but this is not a fair evaluation of the error related to the FSD parameterization due to SP of 50).

Figure:
The heating rate error of various TC experiments for cumulonimbus (evaluated against a 3-D benchmark).Top row: The three baseline TC(FSD) experiments with global FSD estimate ("G") and additional experiments with the actual FSD ("A") together with various distribution assumptions.The bottom row shows only the two experiments using the gamma distribution with different FSDs (for clarity).For simplicity, the split percentile of 50 is used.
Nevertheless, we are primarily interested in the error of TC quantities if the parameterized FSD is used, since this is the experiment, which is actually performed within a GCM.
Finally, we have added the following paragraph at the end of Section 3, exposing general difficulties regarding the TC treatment of mixed-phase region: "In overall summary, the baseline Tripleclouds setup performed well for the apparently most complex deep convective scenario.Nevertheless, improved configurations should be further sought in the future.It would be especially intriguing to contemplate how to better treat the mixed-phase region, where the actual FSD of liquid phase is extremely large.Thus similar optimizations as for the cirrus case study could be introduced, although in the mixed-phase region of the present cumulonimbus, where the actual FSD of ice is close to the global estimate, caution needs to be taken when asymmetrically splitting the cloud fraction."Thank you for this remark, we agree that the "conventional GCM" is not an appropriate label.We have changed the description of this experiment in Table 1 to: "GCM radiation scheme utilizing homogeneous cloudiness".We have also removed the word "conventional" everywhere in the text.We have defined the acronyms LP and FSD in the table.
15. Line 310: Say *why* the "GCM" scheme has his bias: by homogenizing the cloud, the probability of radiation being scattered or emitted is greater.This is even more clearly evident in the cirrus case.
We provided a detailed explanation for the GCM biases (due to homogeneous cloudiness) in the parent study Črnivec and Mayer (2020) -see their Section 5.1.We have thus extended the paragraph as follows: "The physical explanation for the GCM bias arising from homogeneous cloudiness is given by Črnivec and Mayer (2020)." 16. Line 350 and lower panels of Fig. 6: In a fractional sense both the "GCM" and TC errors are very large, especially in the infrared.However, this is not a fair evaluation of the TC method because the fractional standard deviation is extreme in this case (around 3) whereas you are feeding it with a value of 0.75, or using a split percentile of 50% which cannot capture large FSDs.It would be better to discuss your solution to the problem earler.
Although the true FSD of the present cirrus is very large, we aimed to evaluate the TC with existing FSD parameterizations (e.g., global constant of 0.75), since this is an experiment which could be performed in a weather or climate model (where the "true FSD" -i.e., the value based on high-resolution 3-D cloud data -is not known).We therefore think it is fine to evaluate the TC with parametrized FSD (even though the values might be much lower than the "true FSD").We showed in a later section (originally Section 4: "Parameter optimizations") that even using more sophisticated FSD parameterization for ice inhomogeneity of Hill et al., (2015) (which would hopefully result in FSD being closer to reality and thus larger) does not bring any improvements compared to the FSD of 0.75.Simultaneously we have shown that it is crucial to modify the split percentile.In the revised paper, we have exposed these issues earlier (Section 2).17. Figure 8, left panels (and also discussion at lines 431, 437 and elsewhere in section 3.2): It is not the net surface flux that should be shown here, but the cloud radiative effect.This way the true fractional error of the various methods can be worked out.For example, the cirrus case at 60 degrees SZA shows the "GCM" method has a solar bias of -25 W m-2, but the net flux is around 300 W m-2, implying that this is less than a 10% bias.However, it should really be compared to the cloud radiative forcing which Hogan and Kew (2005) estimated to be -39 W m-2 for this case.Thus the error is more like 64%.
Our aim was to consistently analyze the atmospheric heating rate and net surface flux throughout the entire study, because these are the two quantities that are actually being computed in a weather or climate model.In addition, we think that a net surface flux error of 10 % is not negligible and therefore worth showing.
We therefore did not perform the analysis of the top-of-the-atmosphere radiative fluxes (cloud radiative forcing), where the error could perhaps be larger as you suggested.We would like to emphasize again, however, that the cases analyzed in Hogan and Kew (2005) and in the present study are not the same, as they employ different resolution (Hogan and Kew, 2005 used horizontal grid spacing of 1.56 km, whereas we use horizontal grid spacing of 100 m; please see our answer to the first major comment for further details).In summary, the analysis of the low-and high-resolution cirrus might not necessarily bring the same conclusions regarding the magnitude of errors.
18. Section 3.2: I find the discussion of the biases in Fig. 8 rather unsatisfactory because there is inadequate discussion of the role of heterogeneity and overlap, or the problem of compensating errors (although the 3D effect is discussed at length).My interpretation would be as follows: 18a.Stratocumulus: GCM is biased low because it overly homogenizes the cloud, but then the question is why TC is biased high even though the homogeneity is about right.This raises the question as to whether the assumption of maximum overlap of the in-cloud heteorogeneities explains the excessive transmission to the surface.The vertical decorrelation length of in-cloud heterogeneities is typically assumed to be half that of cloud boundaries, so about 1 km, and if this was implemented it would block more of the solar radiation and reduce the positive surface-flux bias.Just an idea.
It is true that the maximum overlap of optically thicker cloudy regions should result in an increased transmission of radiation towards the surface (compared to more realistic overlap considering vertical decorrelation length) -thank you for this idea, which we also included in the discussion.We however later showed (in Section 4.1 of the original manuscript) that the TC is biased high mainly because the inhomogeneity in the baseline experiment is not about right: when using better inhomogeneity estimate of Boutle et al. (2014) the net surface flux bias is significantly reduced.We have extended discussion about the TC net surface flux for the stratocumulus as follows: "When the TC(LP) is applied, the net flux error is mostly slightly reduced, whereas in TC(FSD) baseline experiments the error is increased compared to the GCM.The latter finding is consistent with previous considerations, where it was pointed out that global FSD introduces excessive inhomogeneity to the radiatively important part of the stratocumulus, unrealistically reducing the absorption of solar radiation within the cloud layer.The corresponding increased cloud-layer transmittance, as we demonstrated herein, has important implications for the surface budget, therefore proposed optimization is highlighted in the next section.It should finally be noted that also the assumed maximum overlap of optically thicker cloudy regions could result in a somewhat excessive transmission towards the surface (compared to the situation in a GCM).The vertical decorrelation length of in-cloud heterogeneities is typically assumed to be half that of cloud boundaries (Shonk et al., 2010).If this phenomenon was implemented in the TC scheme it would block more of the solar radiation and reduce the positive surface net flux error." 18b.Cirrus: GCM/TC models all significantly underestimate surface transmission because they all fail to capture the strong cloud heterogeneity.Only later does the reader find Fig. 10 in which a "fix" is presented, but it would be simpler to present the fix at the same time as the original problem.
We agree.In the revised paper we have presented the methodology to fix the treatment of highly heterogeneous scenes already in Section 2.
18c.Cumulonimbus: There is a strong 3D effect for this case so if TC is closer to the 3D benchmark than it is to ICA then it is for the wrong reason: surely you should aim for TC to agree with ICA, then use some other scheme to try to capture the 3D effects on top of TC?The fact that all the TC calculations have a positive bias with respect to ICA is probably due to the incorrect maximum-random overlap assumption in this deep cloud system.
Yes, ideally we would use some other scheme to properly capture the 3-D effects on top of TC (as you do with SPARTACUS in ecRad), but this is unfortunately out of the scope of the present study (we don't yet have the scheme that would take subgrid 3-D cloud-radiative effects into account).So we aim for TC to agree with 3-D rather than with ICA (because the 3-D is the realistic benchmark and we can even capture part of the 3-D effects in our present Tripleclouds implementation; which ICA can not capture at all).See also for example Section 5.1 of Črnivec and Mayer (2020), where we have also shown for the cumulus cloud field that at low Sun the TC bias is smaller than the ICA bias (when both were evaluated against the same 3-D benchmark).This was attributed to the effective treatment of 3-D radiative effects in the TC scheme (which we don't consider as being "wrong" in our case, but rather as "desired").It would of course be better if we had a more advanced 3-D scheme.We have extended the final section "Summary and conclusions" as follows: "Furthermore, in order to properly consider clouds in sheared environments, the vertical overlap rules in the present TC implementation have to be generalized.Finally, if the subgrid horizontal photon transport is to be accounted for in a proficient manner, the two-stream equations need to be extended to include terms representing the in-layer horizontal radiative energy exchange between the cloud and the cloud-free part of the grid box as well as that between the optically thicker and thinner part of the cloud.We currently investigate some of these topics, which will be addressed in a forthcoming study." The fact that all TC calculations have a positive bias with respect to ICA (true in the solar spectral range) can indeed be partly due to overlap being too maximal.This leads to an increased amount of direct (and hence net) solar radiation reaching the ground in the TC calculation compared to that in the ICA.In addition, the aforementioned effective treatment of 3-D effects in the TC scheme (the relevant one in this case: "cloud side escape effect") simultaneously leads to an increased amount of diffuse (and hence net) surface radiation compared to its counterpart in the ICA.We have included these issues in the discussion as follows: "The fact that all TC solar calculations have a positive bias with respect to the ICA can be partly also due to the assumed overlap being too maximal.This leads to an increased amount of direct solar radiation reaching the ground in the TC calculation compared to that in the ICA.It should however be kept in mind that the partial effective treatment of 3-D radiative effects (e.g., cloud side escape) in the TC scheme simultaneously leads to an increased diffuse surface radiation compared to its counterpart in the ICA.Future studies should try to disentangle and quantify these effects." 19. Line 493 and bottom-left panel of Fig. 10: I don't think you can conclude much from this analysis because you have a large error from the maximum-random overlap assumption, so tuning the treatment of cloud horizontal heterogeneity is probably leading to the wrong conclusions.
We don't necessarily have a large error from the maximum-random overlap assumption (see our answer to comment 1).We have however concluded the paragraph with the following remark: "Finally, in more detailed future studies additionally considering vertical overlap issues, it should be kept in mind that the fixing of cloud horizontal inhomogeneity and vertical overlap should be addressed concurrently to avoid the problem of compensating errors (Shonk and Hogan, 2010)." 1 Introduction

General background
The fundamental role of clouds and their interaction with radiation in weather and climate can hardly be overemphasized (e.g., Boucher et al., 2013;Stevens and Bony, 2013;Bony et al., 2015).Clouds are complex phenomena, since they exhibit an immense variety of shapes and sizes (Randall et al., 2003) and highly variable degrees of inhomogeneity (Shonk et al., 2010;Hill et al., 2012Hill et al., , 2015;;Boutle et al., 2014  thermal radiation, the most common effects are radiatively induced cooling at cloud top and warming at cloud base, which promotes convective instabilities within the cloud (Webster and Stephens, 1980).This radiative destabilization of the cloud layer is impelled primarily by thermal radiation, whereas during daytime solar radiation generally has a stabilizing tendency , albeit the latter strongly depends on solar zenith angle ( Črnivec and Mayer, 2019).In the boundary layer, the atmosphere and thereby clouds are directly affected by the Earth's surface, via the transition of heat, moisture and momentum (Baur et al., 2018).The net (difference between downward and upward) surface radiative flux is a crucial component of surface energy budget (Manabe, 1969).Although solar surface flux is customarily markedly larger than its thermal counterpart and thus mostly dominates during daytime ( Črnivec and Mayer, 2019), the latter is important at nighttime when solar forcing is absent.All in all, radiatively induced temperature changes in clouds and at the surface are firmly linked to a broad range of atmospheric moist thermodynamic, turbulent and microphysical processes, e.g.formation of precipitation (Harrington et al., 2000;Klinger et al., 2019).A skillful representation of these ::::::: coupled processes in numerical models as well as of their mutual interplay still poses many grand challenges to atmospheric scientists across the world (Schneider et al., 2017).
The present study aspires to make progress on the treatment of unresolved cloud-radiation interaction in coarse-resolution weather and climate models, often referred to as the general circulation models (GCMs).Despite significant advances in the last few decades, the radiation schemes of GCMs still only crudely represent the interchange between clouds and radiation, being impaired by the poor representation of model cloudiness itself in the first place (Randall et al., 2003).A handy way to tackle these shortcomings is by means of explicit cloud modeling.Over the past decades, large-eddy simulation (LES) and cloud-resolving models (CRMs) (e.g., Klemp and Wilhelmson, 1978;Tao and Simpson, 1993;Khairoutdinov and Randall, 2003;Stevens et al., 2005) have established themselves as a well-acknowledged tool in cloud physics research (Guichard and Couvreaux, 2017).The idea of the so-called superparameterization (??), where a two-dimensional (2-D) CRM is embedded in individual column of a host GCM, has recently been revisited by applying fully three-dimensional (3-D) LES model as a superparameterization, albeit on massively parallel computers (?).Short-term global predictions using direct LES/CRM simulations extending up to a few months or even years are beginning to be feasible (Bretherton and Khairoutdinov , 2015).
Long-term climate projections utilizing coupled atmosphere-ocean systems in a direct high-resolution mode, however, will not be possible for a next couple of decades even on most powerful supercomputers (Schneider et al., 2017).Similarly, despite remarkable advancements in numerical weather prediction (NWP), which is burdened by the users' demand for real-time forecasts, global NWP at a subkilometer scale will stay challenging in the near future (Bauer et al., 2015).
In this way the radiative transfer was extensively studied for cumulus (Davies, 1978;Kobayashi, 1988;Welch and Wielicki, 1989), stratocumulus (??; Cahalan et al., 1995 ::::::::::::::::::::: Cahalan et al., 1994Cahalan et al., , 1995;;Zuidema and Evans, 1998;Di Giuseppe and Tompkins, 2003a), cirrus (?; Carlin et al., 2002;Hogan and Kew, 2005;Zhong et al., 2008;Fauchez et al., 2014) as well as deep convective and anvil clouds (Barker et al., 1999;?;Di Giuseppe andTompkins, 2003b, 2005 ), although some of the earliest work used either very idealized cuboid, single-layered or 2-D clouds.Moreover, some of these studies assessed solely either the ICA or the subgrid cloud variability bias ::: error : and therefore did not shed any light on their relative contribution.Nevertheless, it was commonly found that classic 3-D radiative effects associated with horizontal photon flow manifest most pronouncedly in areas of notable horizontal gradients of optical properties and are thus regularly related to cloud side boundaries.But also the in-cloud horizontal variations of optical depth were found to impact 3-D radiative transfer and especially the GCM-type approximation.Due to aforementioned reasons, conflicting claims can be found in the litera-ture regarding the magnitude and sign of these biases :::: errors.Better understanding of these effects is required to advance the parameterization of cloud-radiation interaction.

Focus of this study
The present study aims to reinforce earlier studies by establishing 3-D benchmarks and further exploring the validity of ICA for various cloud types.In particular, we intend to assess the ICA suitability when the GCM resolution refines to the meso-scale O(10-100 km).Di Giuseppe and Tompkins (2003b), as an illustration, showed that 3-D radiative effects increase as the GCM resolution approaches the meso-scale.In addition, we strive to investigate more realistic cloud morphologies, as we applied finer horizontal grid spacing in cloud-generating models compared to the previous research (e.g., Di Giuseppe and Tompkins, 2003b; Hogan and Kew, 2005;Zhong et al., 2008).Following the idea of SH08, a second cloudy region has recently been incorporated into the broadly used δ-Eddington twostream method with maximum-random overlap assumption for partial cloudiness by Črnivec and Mayer (2020).The insertion of a second cloudy region in the two-stream framework requires an upgrade ::::::: extension : of vertical overlap rules.This task was accomplished exploiting the core-shell model for convective clouds (Heus and Jonker, 2008;Heiblum et al., 2019), where the convective core characterized by updraft and condensate loading is positioned in the geometrical centre of the cloud, enclosed by the shell related to downdrafts and condensate evaporation.In the terminology of radiative transfer, the maximum-random overlap was thus retained for the entire fractional cloudiness and additionally applied for the optically thicker segment.This extended vertical overlap formulation implicitly places the optically thicker cloudy region towards the cloud interior in the horizontal plane, whereas the other optically thinner region covers the cloud periphery, as depicted in Fig. 1 (fourth panel).
Once two-stream radiative fluxes had been successfully imposed onto a system of three-region atmospheric layersin accordance with vertical overlap rules, the research challenge is to optimally set the pair of liquid/ice water content (LWC/IWC) characterizing the two cloudy regions and the corresponding geometrical split of layer cloudiness.The answer to the posed scientific question is critically dependent on the characteristic of the underlying subgrid cloud horizontal variability.The latter is conveniently defined by the fractional standard deviation (FSD) of cloud condensate as well as the shape of its distribution

Cirrus
The cirrus cloud was generated with the stochastic cloud model Cloudgen of Hogan and Kew (2005), described also in Zhong from radar observations in southern England (Hogan et al., 2003;Hogan and Illingworth, 2003).We chose the cirrus uncinus 170 case study of June, 24th, 1999, which is the first of the three cases discussed by Hogan and Kew (2005) and subsequently also by Zhong et al. (2008) 1 and was adopted herein as it is the case with smallest vertical wind shear.The horizontal domain size is 51.2 x 51.2 km 2 with a grid spacing of 100 m.The vertical extent of the domain is 7 km using constant vertical grid spacing of 109 m. horizontal heterogeneity is largest in the central part of the cloud layer, with FSD exceeding 3.5.The cirrus layer is thus by far not uniform, rather it exhibits cellular structures ("generating cells"), which would in reality be associated with convective motions.The latter produce higher supersaturations (Heymsfield, 1977) and increase cirrus ice crystal residence time (Mitchell, 1994), which leads to an increased IWC within the cells.The layer IWC of the present cirrus is lognormally distributed, as this is the intrinsic Cloudgen property.

Cumulonimbus
The cumulonimbus cloud was simulated with the Goddard Cumulus Ensemble cloud resolving model (GCE-CRM), described in detail by Tao and Simpson (1993) and more recently by Tao et al. (2003).The simulation relates to the convective event observed on 23 February 1999 during the Tropical Rainfall Measuring Mission (Simpson et al., 1988(Simpson et al., , 1996) )  The comparison of theoretical distributions with actual LWC/IWC distributions reveals the following findings: the assumption of gaussianity is void for the present cumulonimbus scenario and in each vertical layer either gamma or lognormal distribution was classified as the best fit.Thus, throughout the majority of the liquid region the actual LWC distribution is best approximated with the gamma distribution (best fit in 58 % of layers).In the mixed-phase region the LWC distribution is best approximated with the lognormal distribution (best fit in 85 % of layers), whereas the IWC distribution is best approximated with the gamma distribution (best fit in 85 % of layers).In the ice region the IWC distribution is best approximated with the gamma distribution within the bottommost 30 % and uppermost 13 % of the region (best fit in all cases), while within the remaining central part the lognormal distribution appears to be the optimal approximation (best fit in 58 % of layers).

Configuring baseline Tripleclouds experiments
2.2.1 ::::::: Baseline ::::::::::: Tripleclouds ::::::::::: experiments To utilize the TC radiative solver, a pair of LWC/IWC defining the two cloudy regions has to be generated in each vertical layer.We refer solely to the liquid phase in the remainder of this section, since analogous considerations apply to the ice phase.
Thus, the LWC values characterizing the optically thin and thick cloudy region are referred to as the LW C cn (Cloud thiN) and LW C ck (Cloud thicK).Based on the analysis of high-resolution cloud radar data, SH08 showed that TC performs well for top-of-the-atmosphere (TOA) radiative fluxes when the LWC cn is chosen to be the value corresponding to the 16th percentile of the observed LWC distribution, whereas the LWC ck is determined under conservation constraints of layer-mean ::: the ::::::: in-cloud :::: This method is referred to as the "lower percentile (LP) method" and utilizes a "split percentile (SP)" of 50, implying that cloudiness in each vertical layer is divided through a median of distribution into two equal parts.
When the TC solver resides in a host GCM model, however, the details about the underlying LWC variability are not known, therefore several assumptions have to be introduced.To obtain the pair (LWC cn , LWC ck ) from LW C, which is indeed available in a GCM, we introduce the so-called LW C-scaling factors for the optically thin and thick cloudy region, termed s cn and s ck respectively, fulfilling .:::::: These ::: are :::::: termed ::: s cn ::: and ::: s ck :::: and ::::: fulfill the following relationships: LW C ck = LW Cs ck .
(3) Different parameters to define the degree of cloud horizontal inhomogeneity are employed in the existing literature and numerical models (e.g., ?:::::::::::::::: Cahalan et al., 1994;Smith and Del Genio, 2001;Carlin et al., 2002;Rossow et al., 2002;Oreopoulos and Cahalan, 2005).A frequently used parameter is the previously mentioned When TC is employed in a host GCM, moreover, an assumption about the shape of LWC distribution has to be made.To this end, we test three assumptions for subgrid cloud condensate distribution: Gaussian distribution, which traditionally prevailed s cn (Gaussian) s ck (Gaussian) s cn (Gamma) s ck (Gamma) s cn (Lognormal) s ck (Lognormal) in many models due to its simplicity, as well as more realistic gamma (supported by Barker et al., 1996;Pincus et al., 1999;Carlin et al., 2002;Rossow et al., 2002) and lognormal distribution (supported by Pincus et al., 1999;Hogan and Illingworth, 2003).For a Gaussian distribution, the 16th percentile (suggested by SH08) is given by: although caution needs to be taken as this expression becomes unphysical for FSD > 1.Similarly, according to Hogan et al.
(2019), for a gamma distribution the 16th percentile is approximated by: Finally, according to Hogan et al. (2016), for a lognormal distribution: For any FSD value, the s cn defined with Eq. ( 6) or ( 7) lies in the range between 0 and 1.The desired conservation of LW C implies s ck = 2 − s cn , where the layer cloudiness is geometrically halved.The approach outlined above, utilizing any of the selected distributional assumptions to generate the LWC pair, is referred to as the "FSD method" in the remainder of this paper.
The resulting LW C-scaling factors for Gaussian, gamma and lognormal distributions as a function of FSD are shown in Fig. 4.

Setup of radiative transfer calculations
The radiative transfer simulations were carried out with the libRadtran software :::::::::::::: (Emde et al., 2016 : ).For consistency, the calculation setup is the same as in Črnivec and Mayer (2020).Here we summarize basic parameter settings and provide additional information regarding the treatment of the ice phase, which was absent in our preceding study.Thus, except for 3-D fields of LWC and IWC as well as the assigned effective radii (parameterized following Bugliaro et al., 2011), all atmospheric conditions are assumed to be horizontally homogeneous and correspond to the US standard atmosphere (Anderson et al., 1986).
The domain extends vertically up to a height of 120 km, which is considered to be the TOA.In the lowest portion of the domain where clouds are located we preserved the original high-resolution vertical grid as inherited from the parent cloud model and interpolated the background standard atmospheric conditions onto this grid.It should be reminded that solely LWC and IWC were used as input for radiation calculations to define the cloud fields, while other hydrometeor categories (i.e., precipitationsized particles, such as rain, snow and graupel) were excluded from the analysis.Optical properties of water droplets (assumed to be spherical) were prescribed following the parameterization of Hu and Stamnes (1993).Optical properties of ice crystals were specified based on the parameterization of Yang et al. (2000), assuming habit of hexagonal columns.Solar zenith angle (SZA) was varied between 0 • (overhead Sun), 30 • and 60 • , whereby downward flux at TOA corresponds to 1365, 1182 and 683 W m −2 .At the surface, the temperature of 288.2 K implies upward flux of 389.5 W m −2 according to the Stefan-Boltzmann law.The surface was assumed to have a constant albedo of 0.25 in the solar spectral range, whereas in the thermal spectral range it was assumed to be black (albedo=0).
Table 1 summarizes the baseline ::::::: radiation experiments carried out in this study.The benchmark experiment was performed using the 3-D Monte Carlo radiative model MYSTIC (Mayer, 2009).The horizontal extent of the domain matched that of each individual cloud field case study.Periodic boundary conditions were applied in the 3-D configuration.In addition, the Monte Carlo experiment in ICA mode was carried out, which is the same as the 3-D experiment, except that periodic boundary conditions are imposed at each grid column.The 3-D and ICA experiment were both performed on the high-resolution cloud grid, with the result subsequently horizontally averaged across the domain.The difference between the ICA and 3-D ::::: results was used to assess the impact of horizontal photon transport on domain-averaged (GCM-scale) radiative quantities.Whereas the exact number of traced photons depends on the particular cloud case, it was held sufficiently high, so that the Monte Carlo noise of domain-averaged quantities was kept below 0.1 %.Moreover, we performed a conventional GCM-type calculation on a layer-averaged fractional cloud .The ::::: using ::: the δ-Eddington two-stream method with maximum-random overlap assumption for partial cloudiness : , ::::: which : was recently implemented into libRadtran ( Črnivec and Mayer, 2019)and was used herein in place of the conventional GCM scheme.The Tripleclouds solver was employed in conjunction with both the LP and FSD method .In order to verify (or discard) the validity of the global FSD estimate for heating rates and surface fluxes, we applied its mean value (0.75) in the initial TC experiments together with all three distributional assumptions ::: the ::: LP :::::: method ::::: based ::: on  For each experiment we diagnosed atmospheric heating rate and net surface flux in the solar and thermal part of the spectrum, as well as their total (integrated) effect.The error is measured by the absolute bias :::: error (Eq.?? : 9), relative bias :::: error (Eq.?? :: 10) and for the heating rate profile additionally by the cloud-layer RMSE (Eq.11): absolute biasabsolute error :::::::::

Results and discussion
We present in Sect.?? ::::: discuss : the atmospheric heating rate ::: and ::: net ::::::: surface ::: flux : for each cloud type, whereas in Sect.?? we discuss the net surface flux in the baseline experiments listed in Table 1.The total heating rate, a physically relevant quantity during daytime, is dominated by thermal cooling.Thus ::: This : persistent cloud-top radiative cooling , : is : a typical feature of marine stratocumulus-topped boundary layers (STBLs; Wood, 2012), : . : It : drives convective instability and controls turbulence within the underlying mixed layer (Randall, 1980;Deardorff, 1981;Stevens et al., 1999), when adequately coupled to a dynamical model.The radiative biases ::::: errors : could therefore importantly affect the evolution of the stratocumulus layer itself.
Radiative heating rate for the stratocumulus cloud.The grey-shaded area denotes the cloud layer.
The examination of radiative biases :::: errors : (Fig. 5, bottom row) reveals that these are maximized in the uppermost part of the stratocumulus layer as well.The disparity between the ICA and 3-D is minor: a maximum difference of −0.2 K day −1 is observed in the solar spectral range for SZA of 60 • (cloud side illumination effect; Jakub andMayer, 2015, 2016).In the thermal spectral range, the ICA underestimates the amount of 3-D cooling by about 2.5 K day −1 in the uppermost grid point of the cloud layer (cloud side cooling effect; Klinger andMayer, 2014, 2016).These comparatively small ICA biases ::::: errors are attributed to the minor cloud top topography (difference between the nearby local height maximum and minimum; Zuidema and Evans, 1998) of the present stratocumulus.The radiative transfer, namely, acts to smooth out structures at spatial scales smaller than the photon mean free path, with the latter corresponding to several hundred meters in STBLs (Marshak et al., 1995).
The 3-D solar effects (Fig. 8, bottom row) are present at all SZAs and maximized at 60 • (cloud side illumination), where the ICA bias :::: error of −0.1 K day −1 is observed throughout the majority of the cloud layer.In the thermal spectral range the ICA bias :::: error : is negligible.Similar results were found by Zhong et al. (2008) (recall that the latter investigated the same midlatitude cirrus, although on coarser grid), who also showed that domain-averaged ICA and 3-D heating rates agree within 0.1 K day −1 in both the longwave and the shortwave.
In the GCM the solar heating rate is overestimated by up to 1.4, 1.2 and 0.7 K day −1 at SZA of 0 • , 30 • and 60 • , respectively.
The height where this maximum bias :::: error is observed corresponds with the height of maximum benchmark heating.In the  thermal spectral range, the GCM bias :::: error : enhances radiatively driven destabilization of the cirrus layer by an overestimation of top cooling by 2.8 K day −1 and a substantial overestimation of base warming (bias :::: error : exceeding 5.6 K day −1 ).The thermal GCM bias :::: error : is in close agreement with that observed by Zhong et al. (2008), whereas the solar GCM bias :::: error is by a factor of 2 to 3 smaller.The latter finding indicates the potential dependence of GCM biases ::::: errors on the initial cloud grid resolution, which could affect the TC experiments as well and has to be more thoroughly examined in the future.The daytime GCM bias :::: error profile closely resembles that of its nighttime counterpart, such that the radiatively driven destabilization of the cirrus layer is persistently substantially escalated.
In the underlying stratiform layer the ICA bias :::: error is comparatively small, but it increases again in the bottommost region of shallow cumuli due to their increased side area, where 3-D radiative effects are maximized ( Črnivec and Mayer, 2019).A similar picture is identified in the thermal spectral range, where the maximal ICA bias ::: error : of 0.1 K day −1 is observed in the anvil and the cumuliform region, whereas in the stratiform layer the 3-D effect is limited.
In both solar and thermal spectral range the GCM reveals large biases :::: errors : within the anvil portion and even larger biases ::::: errors in the stratiform layer underneath (Fig. 11, bottom row).The latter are as expected a manifestation of considerable horizontal inhomogeneity observed in the stratiform region (recall that its actual FSD is two times larger than that of the anvil), which implies that the horizontally homogeneous cloud assumption is violated more in the stratiform region than in the anvil.
If the stratiform layer had not been partially shielded by the anvil, the biases ::::: errors : therein would be even larger.For overhead Sun, for example, we observe an overestimation of solar heating by up to 3.6 K day −1 in the anvil region and 3.5 K day −1 in the stratiform region.Thermal GCM bias :::: error : of cloud top cooling up to −4.2 K day −1 and that of cloud base warming up to 3.9 K day −1 is observed within the anvil.Within the stratiform region, thermal cooling is overestimated with a bias of up to :: an :::: error ::: of −5.2 K day −1 and thermal warming is overestimated by 5.4 K day −1 in the GCM configuration.This indicates a significant need for proper TC usage when treating deep convection.All Tripleclouds :::::: baseline : experiments yield a significant reduction of solar, thermal and total heating rate bias :::: error when compared to conventional :: the : GCM.The TC(LP) experiment performs best, reducing the cloud-layer RMSE two-to threefold.
As an illustration, thermal RMSE of 1.5 K day −1 is reduced to solely 0.6 K day −1 .Although the actual LWC and IWC are better approximated with either lognormal or gamma distribution, the assumption of gaussianity works best in practice.
The reason for this is similar as was for the cirrus case study: the actual FSD of the cumulonimbus is mostly larger than the global estimate.As the assumed gaussianity implies the largest difference between the LWC/IWC pair, it partly accounts for the missing inhomogeneity degree introduced by global FSD.Noteworthy, within the stratiform layer (the liquid phase of which is markedly heterogeneous with FSD similar to that of the cirrus case), TC(FSD) experiments represent a considerable improvement compared to the GCM.This could be partially due to radiatively important effect of ice within the stratiform mixed-phase region: the actual FSD of ice is in close proximity to the global estimate, thus acting to reduce the overall TC error in this region.
3.4 Net surface flux

Stratocumulus
Figure ?? (top row) shows the net surface flux underneath the stratocumulus.The ICA bias is small during daytime and nighttime, maximized at overhead Sun (up to −5 W m −2 or −2 %).This is essentially attributed to the photon cloud side escape effect (Várnai and Davies, 1999), where preferential forward scattering on cloud droplets increases 3-D surface downward radiation (Wissmeier et al., 2013).An increased solar absorption in the homogeneous GCM cloudiness implies reduced transmittance and hence underestimated daytime net surface flux.The bias is largest at overhead Sun (−33 W m −2 or −9 %) and decreases with increasing SZA, whereas at nighttime the GCM bias is minor.
When the TC(LP) is applied, the net flux bias is mostly slightly reduced, whereas in TC(FSD) experiments the bias is increased compared to the conventional GCM, therefore optimizations are tested in Sect.??.
In summary, Tripleclouds in its baseline configurations proved to perform well for the apparently most complex deep convective scenario, where it depleted the majority of GCM biases.In the case of stratocumulus and cirrus a refined TC realization is highly desired.

Optimization for overcast cloud scenarios
It was previously pointed out that in the uppermost overcast part of the stratocumulus, the actual FSD is smaller than the introduced global estimate.This might be partially attributed to the fact that overcast grid boxes do not contain cloud edges, which generally contribute to increased variability.Mixing of cloudy and cloud-free air at the edges of clouds, namely, tends to decrease the mean LWC as well as to increase the spread of LWC, both acting to increase the FSD.A grid box excluding cloud edges will therefore have lower FSD.To that end, we test the parametric FSD relationship proposed by Boutle et al. (2014), denoted as B14, for liquid cloud inhomogeneity, developed based on a rich combination of satellite, in situ and ground-based observations.This parameterization takes into account that variability is generally dependent on grid box size and cloud fraction and exhibits a discontinuity at C=1 capturing the aforementioned cloud edge effect.Thus the FSD of liquid phase for a grid box of horizontal size x kmand cloud fraction C is given by: where: Φ c (x, C) = (xC) 1/3 (0.06xC) 1.5 + 1 −0.17 .
Optimization for stratocumulus (same experiment labeling on all panels). :::: ICA. :::::: Future :::::: studies :::::: should :: try :: to :::::::::: disentangle ::: and ::::::: quantify ::::: these :::::: effects.: The key point worth mentioning when highly heterogeneous scenes as is the cirrus cloud are tackled with the Tripleclouds solver, is that the split percentile (SP) of 50 (geometrically halving layer cloudiness when allocating optically thin and thick portions of the cloud) is not the best choice (Hogan et al., 2019).The examination of IWC distribution in each vertical layer of the present cirrus indeed reveals that these are highly skewed (with a modal value close to zero and a long tail with rarely occurring high IWC).Therefore it seems reasonable to allocate a larger portion of the cloud to : In :::::: overall :::::::: summary, : the optically thinner region.This concurrently implies increasing the weighting of IWC cn and decreasing the weighting of IWC ck , whereby the latter is shifted to a higher value to conserve the layer mean.In order to discern the optimum geometrical partitioning of the cirrus into two parts, we carried out multiple experiments with global FSD, gradually increasing the SP from 50 to 99 (the limit of 100 coincides to the horizontally homogeneous cloud representation).Further, we aim to assess the advantage of more sophisticated FSD parameterizations.We thus evaluate the parameterization for ice cloud inhomogeneity of Hill et al. (2015) H15, developed on the basis of CloudSat (Stephens et al., 2002(Stephens et al., , 2008) ) data products.All TC experiments presented in this section model the subgrid IWC distribution as lognormal.
Optimizations for cirrus (subgrid variability is modeled lognormal in all TC experiments shown).
Figure 10 (top row) shows the cloud-layer RMSE and net surface flux bias using TC with different FSD parameterizations for three selected splitting events, characterized by the SP of 50 (baseline), 75 and 90 (allocating 3/4 and 9/10 of layer cloudiness to the optically thinner region ).For comparison the GCM experiment is shown as well.It is apparent that there is a considerable sensitivity to the choice of geometrical splitting, with the most asymmetrical split (matching the SP of 90) performing best in all cases.Noteworthy, at a given splitting event, the experiments where the FSD is parameterized according to H15 mostly lead to degraded results compared to those with global FSD (in particular at best-split scenario with SP of 90).Thus although the parameterization of H15 incorporates height dependence of horizontal variability (via cloud fraction), it underestimates the actual FSD being even smaller than the global estimate (Fig. 10, middle row, left), which brings the aforementioned radiative degradation.Further research oriented towards advanced retrievals of high cloud inhomogeneity is therefore firmly advocated.
Vertical profiles of solar and thermal heating rate in TC experiments with global FSD for the aforementioned splitting events are further compared with the GCM in Fig. 10.For highly asymmetrical splitting, the cloud-radiative bias throughout the majority of the cirrus layer is significantly reduced.

Summary and conclusions
This study aims to advance the conceptual understanding of radiative transfer in marine stratocumulus, midlatitude cirrus and tropical deep convective clouds.The focus is laid on the issues related to misrepresentation of cloud horizontal inhomogeneity in coarse-resolution weather and climate models, which are tackled with the aid of the Tripleclouds radiative solver.The TC method, primarily introduced by SH08, is an approach accounting for horizontal cloud inhomogeneity by using two regions in each vertical layer to represent the cloud (as opposed to one, which is the convention of traditional cloud models).One of these regions is utilized to represent the optically thicker portion of the cloud, whereas the other region represents the residual optically thinner portion.The challenge is to properly set the pair of liquid/ice water content defining the two cloudy regions and divide the layer cloudiness geometrically in the corresponding two parts.The primary objective of the present work was to evaluate the global FSD estimate together with different assumptions for cloud condensate distribution, which are commonly applied in models.The TC concept was recently integrated in the efficient δ-Eddington two-stream radiation scheme ( Črnivec and Mayer, 2020) and was used herein to answer these questions within the scope of the libRadtran library.For our study we chose three intrinsically contrasting cloud typesgenerated by different cloud models, which should reflect diverse cloud conditions occurring globally.These high-resolution cloud field data allow to gain important insights about small-scale cloud variability and give the opportunity to compare the actual modeled variability with the global estimate or other existing parameterizations.For each cloud type, various TC experiments were evaluated against a 3-D benchmark calculation.These results were compared with the conventional GCM computation utilizing homogeneous layer cloudiness, which can be viewed as the upper limit for the tolerable TC error.Moreover, the ICA approximation was compared with the 3-D benchmark to quantify the bias :::: error related to neglected horizontal photon transport.A systematic investigation of cloud-layer heating rate and net surface flux bias :::: error was provided for each selected cloud case.
It was found that in the majority of applications, the ICA is significantly more accurate than the conventional GCM experiment, indicating :::: GCM ::::::::::: experiment.:::: This :::::::: indicates a large potential for Tripleclouds, which reduces the bias :::: error : related to unresolved cloud structure, but not to horizontal photon transport.Regarding the optimal TC configuration, which aims to minimize the radiative bias :::: error, the exact conclusions drawn depend on each particular cloud typecase study.In general, the simplest TC arrangement using a globally constant FSD and geometrically halving the layer cloudiness worked best for the apparently most complex deep convective scenario.In case of stratocumulus and cirrus, an improved TC performance was highly desired.To that end, the second objective of the present study was to assess recent advanced FSD parameterizations, characterizing systematic departures from global mean cloud variability observed for liquid and ice phase.For the stratocumulus, an optimization in terms of a parametric FSD relationship portraying reduced horizontal variability at overcast conditions lead to a substantially improved TC realization.For extremely heterogeneous cirrus, on the other hand, allocating the greater portion to the optically thinner part of the cloud (e.g., approximately 9/10 of layer cloudiness in the solar part of the spectrum) , proved to be of crucial importance in the TC settings, eliminating the vast majority of GCM biases ::::: errors.These findings are in support of cloud regime dependent approaches, which ought to be further boosted to be used in radiation schemes of next-generation atmospheric models.Whereas current GCMs do not explicitly predict cloud meteorological regimes (e.g., whether model cloudiness appears in the form of cumulus or stratocumulus), they have the ability to diagnose cloud type based on temperature and humidity fields (Norris, 1998).An alternative is to consider the physical processes responsible for cloud formation as imprinted in the parameterization schemes activated to generate cloudiness within a model grid box (e.g., shallow convection as opposed to large-scale saturation).We thereby propose that the TC configuration should be adequately adjusted according to cloud type.

:
figurations refer to distinctive improvements desired in the case of overcast as well as ::for ::::::: overcast ::: and : extremely heterogeneous cloudiness.Together with the parent study ofČrnivec and Mayer (2020), this work presents the first usage of Tripleclouds to consistently study the cloud-layer heating rate and net surface flux.The remainder of this paper is structured as follows: the cloud field data and methodology are presented in Sect. 2. The results of the baseline radiation ::::::: radiative ::::::: transfer experiments are discussed in Sect.3. Building upon this discussion, optimized TC realizations are proposed in Sect.??.Section 4 concludes with a brief summary and outlook.

Figure 2 .
Figure 2. Vertically integrated visible optical thickness of selected cloud field case studies.Note that in case of cumulonimbus, the blue/red shading denotes optical thickness of liquid/ice water content.
Figure 2 (left) visualizes the cloud field in terms of vertically integrated optical thickness.Vertical profiles of averaged LWC, cloud fraction (defined by LWC > 10 −3 g m −3 ) and FSD are shown on Fig. 3 (left).The overcast stratocumulus scene is topped slightly above 1 km height.The FSD, although roughly centered around the global estimate, strongly depends on the position within the cloud layer: it exhibits a maximum (1.2) in the lowest portion of the cloud layer and a minimum (0.3) in its uppermost radiatively important region.

Figure 3 .
Figure 3. Characterization of selected cloud field case studies in terms of averaged LWC/IWC (top row), cloud fraction (middle row) and fractional standard deviation of LWC/IWC (bottom row), whereby the vertical black line shows the mean global FSD estimate and the grey shaded area denotes its uncertainty ::::::::::::: (Shonk et al., 2010 : ).
Large-ScaleBiosphere-Atmosphere (TRMM-LBA) experiment in Amazonia.The horizontal domain size is 64.0 x 64.0 km 2 , with the vertical extent of the domain being 23 km, which is sufficient to allow the growth of tropical cirrus anvil.The grid spacing of 250 m is applied in each horizontal direction and 200 m in the vertical direction.The simulation is described byLang et al. (2007) and briefly byZinner et al. (2008).Due to light environmental winds (Fig.2bofLang et al., 2007), the convection was rather weakly organized. 2Figure2(right) visualizes the cloud field in terms of vertically integrated optical thickness.Vertical profiles of averaged LWC and IWC, the corresponding cloud fraction as well as FSD are shown on Fig.3(right).The case study is characterized by three distinct cloud layers: a liquid phase region extending from 0.8 km to 4.4 km consisting of shallow cumuli, a mixed-phase stratiform region located between 4.4 km and 8.2 km, and an ice phase region extending from 8.2 km to 17.4 km height, encompassing the cumulonimbus deep convective core and the anvil.Remarkably, the stratiform layer is highly heterogeneous, with the maximum FSD of the liquid phase exceeding 3.2.The maximum FSD in the bottommost cumuliform region as well as in the uppermost anvil region reaches approximately half of this value.
:::::::: introduced : fractional standard deviation (FSD) of LWC, which is simply the standard deviation (σ LW C ) divided by the mean (LW C) :: of ::::: LWC.The FSD is a convenient measure, since it accounts for a strong correlation between LW C and σ LW C(Smith and Del Genio, 2001;Carlin et al., 2002;Hill et al., 2012).Based on a comprehensive review of numerous observational studies encompassing diverse cloud data sets,Shonk et al. (2010) converted various variability measures into a single globally applicable FSD parameter, whose .::: Its : mean value and uncertainty are:

Figure 4 .
Figure 4. Scaling factors (s cn , s ck ) of LW C for Gaussian, gamma and lognormal distributions.The black line and the grey area represent mean global FSD and its uncertainty.
RMSE = (y − x) 2 , (11) where x represents the 3-D benchmark and y represents the biased quantity (outcome of ::::::: outcome :: of ::::: either the ICA, GCM or any TC experiment).The cloud-layer RMSE denotes the RMSE evaluated throughout the vertical extent of the cloud layer of each particular cloud field case study.
Figure 5(top row)  shows the radiative heating rate in the benchmark 3-D experiment for the stratocumulus cloud.There is large absorption of solar radiation in the cloud layer, resulting in the maximum heating rate of about 53, 47 and 27 K day −1 at SZA of 0 • , 30 • and 60 • , respectively.The peak heating rates are concentrated in the uppermost part of the cloud layer, since both cloud fraction and LWC increase from cloud base towards cloud top.In the thermal spectral range the cloud layer is subjected to strong cooling, reaching a peak of almost −140 K day −1 at the same height where maximum solar heating is attained.These large heating and cooling rates are partially a manifestation of high vertical resolution(Di Giuseppe and Tompkins, 2003a).

Figure 8 .
Figure 8. Radiative heating rate for the cirrus cloud.The grey-shaded area denotes the cloud layer.

Figure 8 (
Figure8(top row) shows the radiative heating rate in the benchmark 3-D experiment for the cirrus cloud.The solar absorption in the ice layer results in a maximum heating rate of about 3.6, 3.3 and 2.2 K day −1 at SZA of 0 • , 30 • and 60 • , respectively.

Figure 11 .
Figure11.Radiative heating rate for the cumulonimbus cloud.The grey-shaded area denotes the cloud layer.The liquid region is shaded light grey, the mixed-phase region is shaded middle grey and the ice region is shaded dark grey.
14. Table1: "Conventional GCM" is not an appropriate label since many/most operational GCMs now use the McICA method.You could call it "Homogeneous cloud"?For this table to be most easily understood, the acronyms LP and FSD should be defined in the table.

Table 1 .
List of radiative transfer experiments.