Scheme for calculation of multi-layer cloudiness and precipitation for climate models of intermediate complexity

In this study we present a scheme for calculating the characteristics of multi-layer cloudiness and precipitation for Earth system models of intermediate complexity (EMICs). This scheme considers three-layer stratiform cloudiness and single-column convective clouds. It distinguishes between ice and droplet clouds as well. Precipitation is calculated by using cloud lifetime, which depends on cloud type and phase as well as on statistics of synoptic and convective disturbances. The scheme is tuned to observations by using an ensemble simulation forced by the ERA-40-derived climatology for 1979–2001. Upon calibration, the scheme realistically reproduces basic features of fields of cloud fractions, cloud water path, and precipitation. The simulated globally and annually averaged total cloud fraction is 0.59, and the simulated globally averaged annual precipitation is 100 cm yr−1. Both values agree with empirically derived values. The simulated cloud water path is too small, probably because the simulated vertical extent of stratiform clouds is too small. Geographical distribution and seasonal changes of calculated cloud fraction and precipitation are broadly realistic as well. However, some important regional biases still remain in the scheme, e.g. too little precipitation in the tropics. We discuss possibilities for future improvements in the scheme.


Introduction
Clouds are an important part of the climate system, linking hydrological processes with radiative transfer and atmospheric dynamics.Since the mid-1990s, climate models include prognostic cloud schemes calculating cloud fractions (i.e. the fractional areal coverage by clouds) and cloud water content (Solomon et al., 2007;Zhang et al., 2005;Williams and Tselioudis, 2007).While such schemes are quite elaborate in the state-of-the-art models, some unresolved problems remain (Stephens, 2005;Williams and Tselioudis, 2007;Cesana and Chepfer, 2012).In particular, there is ample evidence that uncertainty in cloud response to external, e.g.anthropogenic forcing, constitutes the largest part of the overall uncertainty in the response of global climate models (Stephens, 2005;Bony et al., 2006;Dufresne and Bony, 2008;Soden and Vecchi, 2011).
For Earth system models of intermediate complexity (EMICs) (Claussen et al., 2002;Petoukhov et al., 2005;Zickfeld et al., 2013;Eby et al., 2013) this problem is even more actual.Most models of this type contain quite simplified cloudiness schemes, frequently accounting only for effective single-layer clouds (see, e.g.Table of EMICs at http://www.pik-potsdam.de/emics/toe_05-06-07.pdf).Such an approach obviously precludes to resolve dominant influence of upper-level clouds on long-wave radiative transfer in the atmosphere, and low-level clouds on the respective shortwave transfer (Stephens, 1978;Liou, 2002).In addition, from simulations with general circulation models it is expected that global warming is accompanied by smaller (larger) cloud fractions in the lower (upper) troposphere (e.g.Solomon et al., 2007).When accounting only for single-layer clouds, simplified climate models cannot reproduce these changes in cloud fractions.Further, one-layer cloud schemes may provide only limited representation of aerosol-cloud interaction (first and second aerosol indirect effects related to changes Table 1.List of symbols used throughout the paper.Long dash in the first column indicates that the corresponding variable is nondimensional.Variable modifiers: j indicates cloud type (= sl, sm, sh, co), and k stands for cloud phase (= drop, ice).

variable and units description
H b,j [m] height of cloud base H t,j [m] height of cloud top H EBL [m] height of the equivalent barotropic level H PBL [m] height of the top of the planetary boundary layer H trop [m] height of the tropopause h j [m] cloud thickness c j [-] cloud fraction c tot [-] total cloud fraction total precipitation q v [kg(H 2 O) kg(air) −1 ] specific humidity q v,0 [kg(H 2 O) kg(air) −1 ] specific humidity at the surface cloud water/ice content per unit volume vertically integrated cloud water/ice content per unit area w e [m s −1 ]  effective vertical velocity (see Eq. 2) w ls [m s −1 ] large-scale vertical velocity w syn [m s −1 ] synoptic-scale standard deviation of vertical velocity in cloud albedo and lifetime respectively; both effects results from an impact of hydroscopic aerosols on the size of clouds droplets and ice crystals, e.g.Charlson et al., 1992;Solomon et al., 2007).Among EMICs which currently have an effective singlelayer cloudiness scheme are the models developed at the Potsdam Institute for Climate Impacts Research (Climber-2, Petoukhov et al., 2000;Ganopolski et al., 2001, andClimber-3α, Montoya et al., 2005) and at the A. M. Obukhov Institute of Atmospheric Physics, Russian Academy of Sciences IAP RAS CM; see Mokhov and Eliseev (2012).Currently, both institutes are developing new versions of the EMICs (Coumou et al., 2011;Eliseev et al., 2011).As a part of this programme, we are working out a new cloud-precipitation scheme.This scheme describes three-layer stratiform clouds and one effective type of convection clouds.
In the present paper, the current version of the scheme is described and tested offline for the present-day climate.

Governing equations
The developed scheme considers four cloud types within a given grid cell.The first three cloud types describe lowlevel, mid-level, and upper-level stratiform clouds (thereafter denoted with the subscripts sl, sm, and sh respectively).This distinction corresponds to observational experience at large horizontal scales (Tian and Curry, 1989;Mazin and Khrgian, 1989).The fourth cloud type is denoted by subscript co and represents convective (cumulus) clouds.
The scheme is designed for use in Earth system models of intermediate complexity.This is the reason why we tried to keep all equations as simple as possible.The latter precludes usage of more elaborated approaches which are implemented in the state-of-the-art global circulation models.In the present scheme, some equations are just derived heuristically and tuned to observations (see below).
Values of basic variables are listed in Table 1.

Cloud vertical boundaries and extent
Height of convective clouds base H b,co is related to the planetary boundary layer height H PBL : where C H,co is a constant.In addition, effective vertical velocity w e = w ls + a w e ,1 w syn + a w e ,2 w oro + a w e ,3 w conv (2) is checked to be positive at this level.Otherwise, it is assumed that no convection occurs at a given geographic location.Eq. ( 2) was introduced by (Petoukhov et al., 2000, their Eq. 36) in order to represent well-known relations of cloud fraction with intensity of large-scale ascent w ls (in terms of the Climber and IAP RAS models, "large-scale" dynamics are dynamics with horizontal spatial scales larger than the Rossby deformation radius and with temporal scales larger than several days), with intensity of synoptic-scale stirring (expressed via w syn ) and with the orographically forced disturbances (for which w oro = ∇ H • u(0) serve as a substitute; here ∇ H is a horizontal gradient operator, and u(0) is nearsurface horizontal wind).Effective vertical velocity in Eq. ( 2) is calculated similar to Eq. ( 36) in Petoukhov et al. (2000), but with coefficients a w e ,1 and a w e ,2 depending on cloud type.An additional modification with respect to Eq. ( 36) in Petoukhov et al. (2000) is due to convective stirring: the term a w E ,5 w conv is introduced with Here q v (0) is near-surface specific humidity, w conv,0 = 0.01 m s −1 , and q v,0 = 0.01 kg(H 2 O) kg(air) −1 .In the scheme, a w E ,3 is set to zero for stratiform clouds.Thus, the last term in Eq. ( 2) is applied only to convective clouds.
It should be noted that H b,co may depend on atmospheric moisture content as well.However, this dependence is rather weak (Mazin and Khrgian, 1989, p. 173, Eq. 1) and simply ignored in our work.
In the current setup, heights of stratiform cloud bases are related either to the height of the planetary boundary layer H PBL or to the height of the equivalent barotropic level H EBL (which is defined as a level at which motions are equivalent to the vertical average of motions in the corresponding column; this level is close to the 500 hPa isobaric level (Holton, 2004, p. 450)) or to the height of the tropopause H trop (Petoukhov et al., 1998(Petoukhov et al., , 2003)).All H PBL , H EBL , and H trop are external parameters of the scheme.The relations read where C H,sl , C H,sm , and C H,sh are parameters.This roughly corresponds to the observational evidence summarised in pp.162-175 of the book by Mazin and Khrgian (1989).
In particular, low-level stratiform cloud bases are typically located close to H PBL .Mid-level cloud bases are located, as a whole, slightly below equivalent barotropic level; their heights only weakly depend on season.Bases of the upperlevel stratiform clouds are shallower in the higher latitudes than in the lower latitudes; the same is true for the tropopause height (Hoinka, 1998).
Calculation of geometric thickness of stratiform clouds is similar to that used in Petoukhov et al. (2000): where c j is cloud fraction, j ∈ {sl, sm, sh} stands for cloud type, parameter h j,0 depends on this type, l h is constant, and the dependence on temperature is Here C h and T h,0 are constants.Cloud temperature T j is assigned to the respective value at cloud base: Finally, heights of the stratiform cloud tops are computed according to H t,j = H b,j + h j .Heights of convective cloud tops are related to the height of the tropopause Geometric thickness of convective clouds is calculated as H t,co − H b,co .In Eq. ( 7), C t,co is a function of specific humidity (via vertical velocity due to convective stirring w conv , see Eq. 3): with an additional constraint that C t,co is smaller than the prescribed value C t,co,max .In Eq. ( 8), C t,co,1 and C t,co,2 are constants.
Thus calculated heights are associated with the nearest vertical level corresponding to input variables.

Cloud fraction
For stratiform clouds, cloud fractions (i.e. the fraction of the area covered by clouds of type j ) are calculated similar to Eq. ( 35) in Petoukhov et al. (2000): Here RH(H b,j ) is relative humidity at cloud bases, l RH,j is a constant, and In Eq. ( 10), C c,1,j , C c,2,j , and w e,0 are constants.Similar to other diagnostic cloud schemes used in global climate models, our scheme assumes positive correlation between relative humidity and cloud fraction.In turn, F c,w e ,j represents the impact of synoptic-scale, convective, and orographic stirring on cloud formation.Basically, this stirring enhances cloud fraction, but saturates at high w e .
Convective clouds are allowed to develop only if w e is positive.If this condition is fulfilled, convective cloud fraction is computed according to Eq. ( 38) in Petoukhov et al. (2000): Here c co,0 , w e,co , and q v,co are constants.Similar to Eqs. ( 9) and (10), Eq. ( 11) was derived heuristically assuming that c co should increase if either atmospheric moisture content or synoptic-scale, convective, and orographic stirring is increased.Again, both dependences should saturate at large q v (0) and w e respectively.
Because stratiform and convective clouds may coexist at a given layer within a given grid cell, it is checked that in every layer, c co + c j ≤ 1, where c j is either c sl or c sm or c sh .

A. V. Eliseev et al.: Scheme for cloudiness and precipitation for EMICs
If this condition is not met, convective cloud fraction is reduced to c co = 1 − c j .In other words, if both stratiform and convective clouds coexist in a given grid cell, the former is considered to be favoured.
Total cloud fractions are computed by overlapping clouds at different levels.Convective clouds are always considered as a single column with maximum overlap between individual computational layers.For stratiform clouds, a random overlap between low-, mid-, and upper-level clouds is always used.However, if, for example, H t,co > H b,sm , then the area covered by cumulus clouds is removed from the latter random overlap for low-and mid-level stratiform clouds.A similar approach, but extended to the upper-level stratiform clouds as well, is used if H t,co > H b,sh .

Cloud water and ice content
For stratiform clouds, the cloud water path (commonly defined as the column amount of liquid and frozen water in clouds) is calculated after Eq. (2) on p. 332 in Mazin and Khrgian (1989): where T f = 273.16K, and α W and r MK are constants; here the subscript "MK" indicates that this equation is adapted from Mazin and Khrgian (1989).Cloud water content is then distributed vertically, assuming that lateral boundaries of stratiform clouds are vertical and W j profile is homogeneous within the cloud.For convective clouds, total cloud water path W co is calculated by integrating the respective vertical profile over the cloud depth Here Q co (z) is volumetric cloud water/ice content which is computed using Eq.(1) on p. 337 in Mazin and Khrgian (1989): where In turn, maximum volumetric water/ice content in convective clouds, Q co,max is approximated based on the results reported in Fig. 2 on the same page in Mazin and Khrgian (1989): where b 1,MK , b 2,MK , and b 3,MK are constants.In applying Eq. ( 16), an additional check that Q co,max ≥ 0 is performed.For all cloud types, ice and droplet clouds are distinguished.The molar fraction of frozen and non-frozen water molecules, f ice and f drop respectively, at a given height z is calculated according to Rotstayn (1997): The values of T m,1 and T m,2 are assumed to be independent of cloud type.Total cloud water path (per grid cell) W tot is calculated as a weighted mean of W j (j = sl, sm, sh, co) assuming the same overlap as for total cloud fraction.

Precipitation
Precipitation rate is computed as a sum of large-scale (stratiform) and convective precipitation: Large-scale precipitation is calculated by summing the contributions from all stratiform clouds in a given grid cell: P ls = P ls,sl + P ls,sm + P ls,sh , with P ls,j = c j • f drop P ls,j,drop + f ice P ls,j,ice .
In turn, where j indicates cloud type, k stands for cloud phase (either droplet or ice), Q j is volumetric water content in clouds, and τ j,k is the lifetime of cloud type j in phase k.
Convective precipitation is attributed to cumulus clouds.It is calculated by integrating precipitation in the vertical direction where p co represents the contribution to P co from the infinitesimally thin vertical layer.The latter is and t he contribution from convective clouds in phase k reads For all cloud types, lifetime is calculated similar to that in Petoukhov et al. (2000) where j ∈ {sl, sm, sh, co}, k ∈ {drop, ice}, F c,w e ,j is the same as in Eq. ( 10), and a τ is a constant.We assume that lifetimes for liquid and frozen parts of clouds of all types in a given grid cell are linearly related to each other: where k τ,ice is a constant.For liquid stratiform clouds (j ∈ {sl, sm, sh, } we assume that lifetime is a weak function of cloud fraction with the constant τ 0 .Thus, clouds which are more horizontally extensive exist longer than smaller (presumably broken) clouds.A similar assumption for convective clouds is not needed because these clouds basically exists as systems of localised towers.Therefore, for liquid convective clouds, and k τ,conv is a constant.Note that the partition between ice and liquid cloud particles may be changed during their fall to the ground.As a result, it is impractical to use f ice or f drop to calculate rain-or snowfall rate at the surface.It is assumed to be calculated by the model's land surface scheme based on surface temperature.

An approach
First of all, the scheme was tuned manually to arrive at the parameter values listed in Table 2.This was done in order to set a reasonable starting point for the automated calibration procedure described below.This parameter set as well as the simulations with this set are referred to as initial.
In the latter automated calibration, governing parameters of the scheme were sampled by using the Latin hypercube sampling (McKay et al., 1979;Stein, 1987).We chose only to sample the parameters which are either most uncertain or those which modify the results of calculations with the Table 2. List of the standard values of the governing parameters of the scheme.Long dash in the first column indicates that the corresponding variable is non-dimensional, and in the last column it shows that a specific parameter is not applied to cumulus clouds.variable and units value 0.990 0.990 0.990 0.998 scheme most strongly.In addition, some parameters are redundant in the scheme (e.g.any change of w conv,0 may be compensated by an opposite relative change in the value of a w E ,3 ), and for some it is unclear how to prescribe their prior ranges without a loss of consistency with observations (e.g.all parameters adapted from Mazin and Khrgian, 1989, and denoted by subscript MK).The parameters which are varied in the presented simulations are listed in Table 3.This table also contains the ranges in which these parameters are varied.For all parameters, uniform (non-informative) priors were chosen.Total sample size in parameter space was 5000.For comparison with observations, only such variables are chosen for which relatively reliable data sets exist.Those variables are total cloud fraction c tot , total (vertically integrated over the whole atmospheric depth) cloud water and Table 3. List of the perturbed parameters of the scheme together with their priory ranges.Long dash in the first column indicates that the corresponding variable is non-dimensional.The symbols "SL", "SM", "SH", and "CO" indicate particular cloud types according to classification used in the scheme.In the last column, Bayesian mean and standard deviation are shown.
ice content W tot , and total precipitation rate P tot .In addition, to assess partition between stratiform and convective clouds, a contribution to P tot from large-scale and convective precipitation is assessed as well.Total score for the scheme is calculated by multiplying the individual skills for cloud fractions S c , cloud water path S W , and precipitation S P : The goal of the optimisation procedure is Skill score for cloud fraction is constructed from its globally and annually averaged value, and fields for annual mean, January, and July cloud fractions: where is N (X; X m , σ X ) is a normal distribution function of variable X with mean X m and standard deviation σ X .In turn, c tot,g,ann is the globally and annually averaged total cloud fraction.Here and below, indices M and O stand for modelled and observed fields, respectively.Skills S c,ann , S c,Jan , and S c,Jul are computed as in Taylor (2001): where X stands for any of "c,ann", "c,Jan", and "c,Jul", and function In Eq. ( 31) r X is the coefficient of the spatial correlation between area-weighted modelled and observed fields of X, and A X is the so-called relative spatial variation calculated according to where A 2 X,M is the spatial average of X M − X M,g 2 , and X M,g is a globally (but not necessarily annually) averaged value of the modelled field X M .In turn, A X,O is defined similar to A X,M but for the observed field.Skill score for cloud water content is calculated by using an equation similar to Eq. ( 28): The meaning of terms on the right-hand side of Eq. ( 33) is analogous to that in Eq. ( 28).This is only applied for total (vertically integrated) cloud water path W tot .The procedure to calculate terms on the right-hand side of Eq. ( 33) is again similar to Eqs. ( 29) and (30).
Precipitation skill score is S P = S P,g S P,ann S P,Jan S P,Jul .
Because it is important to distinguish between large-scale and convective precipitation, P ls and P co respectively, individual terms in Eq. ( 34) are calculated differently from their counterparts in Eqs. ( 28) and (33).In particular, S P,g = S P,tot,g S P,rat,g , where S P,tot,g = N P tot,g,ann,M ; P tot,g,ann,O , σ P tot,g,ann,O , S P,rat,g = N p rat,g,ann,M ; p rat,g,ann,O , σ p rat,g,ann,O .Here P tot = P ls + P co , p rat = P co /P ls .Further, S P,ann = S P,tot,ann S p,rat,ann . (37) Here S P,tot,ann = T P,tot,ann , S p,rat,ann = T p,rat,ann .
The terms S P,Jan and S P,Jul are calculated by using equations similar to Eqs. ( 37)-(39) but with respective monthly mean fields in place of annual mean ones.
After that, sampled parameters were subjected to Bayesian averaging (Kass and Raftery, 1995;Hoeting et al., 1999) using total scores S as weights.The ensemble means for all sampled parameters obtained in this way were considered as a calibrated parameter set thereafter in this paper (Table 3), and their standard deviations were considered as a measure of respective allowable range width.We checked different procedures to obtain this optimal parameter set.In particular, we have tried to zero weights if S's were smaller than the half of their maximum.In this approach, ensemble mean values were basically unchanged but their standard deviations were smaller.In addition, we have tried to manually select a best-performing sample and use its parameters as optimal.However, in the latter approach no parameter sample was superior with respect to their Bayesian means.

Forcing data and observational data sets
The simulations were forced by the monthly mean ERA-40 reanalysis (Simmons and Gibson, 2000) climatology for 1979-2001.Synoptic-scale standard deviations of vertical velocity were calculated by using the 2.5-6 day Murakami filter identically to that used by Petoukhov et al. (2008), and were converted to z coordinates assuming geostrophy.Height of the planetary boundary layer was set equal to 1.5 km, and the value 5.5 km was used for the height of the equivalent barotropic level (Charney and Eliassen, 1949;Hoskins and Karoly, 1981) -The MODIS Science Team (MODIS-ST) data set (Frey et al., 2008).Instead of the CERES algorithm, 14 of 36 spectral channels of MODIS instruments (with the central wavelengths from 0.66 to 13.94 µm) are used in the MODIS-ST cloud mask algorithm to discriminate cloud pixels from clear sky.
-ERA-40 reanalysis data (Simmons and Gibson, 2000).This data set is affected by imperfections of the forecast model.This is especially true for cloud-related variables belonging to the so-called class "C".However, because our simulations will be forced by the ERA-40 data, it is instructive to compare simulation output with that reanalysis data.
Basically, satellite retrievals reliably detect total cloud fraction.However, because of the "satellite view" of cloud layers (upper cloud layers may mask lower ones) mid-and lowerlevel cloud fractions detection is not straightforward.This is the basic reason why only total cloud fractions were used for calibration, and not cloud fractions in different layers.
Another reason is the above-mentioned (see Sect. 2) difference between the definition of the cloud layers in the present scheme and that used in common cloud products.An extensive intercomparison between these data sets was reported by Chernokulsky andMokhov (2010, 2012).We set σ c tot,g,ann,O to 0.1, which is a typical value for interannual standard deviation of globally averaged total cloud fractions as estimated by using the ISCCP data.Cloud water path W tot was evaluated against the CERES retrievals (Minnis et al., 2011).In this data set, the cloud water path is computed as function of cloud optical depth and appropriate effective particle size.For S W,g , σ c W tot ,g,ann,O is set ad hoc to 0.1 × c W tot ,g,ann,O .
Total precipitation is compared with the GPCP-2.2data set (Global Precipitation Climatology Project, version 2.2, an update from Huffman et al., 2009).Lacking purely empirical data about the subdivision of total precipitation into large-scale and convective precipitation, we have calibrated the scheme by using the p rat calculated based on ERA-40 data.Note that while global annual precipitation fractions differ by 29 % (Table 4), the spatial pattern of precipitation rate in ERA-40 is close to that in GPCP data.For the GPCP-2.2data, σ P tot,g,ann,O = 1.5 mm mo −1 and σ p rat,g,ann,O = 0.1.
We arbitrarily divided these data into training and comparison sets.The training set consists of ISCCP data for cloud fraction, CERES data for cloud water path, GPCP data for total precipitation, and ERA-40 data for fraction of large-scale precipitation in a total one.All other data were used only for comparison.
For the above-mentioned data, multi-year monthly means were constructed for 2001-2006.This period formally differs from that for the forcing data.However, this is not a crucial point for our calibration because the scope of this paper is to determine climatological means.

Results of calibration
Basically, the scheme with calibrated parameters agrees better with observations relative to its counterpart with the initial parameter set.This is evident even at the global scale, with most marked improvement for cloud water path W tot (Table 4).Slight deterioration is visible for the fraction of convective precipitation in total precipitation.

Cloud fractions
At the global scale, cloud fractions simulated by the scheme with the calibrated parameter set equal 0.59, which is slightly below the observational range 0.60-0.67(Table 4); more extensive comparison of different empirical data sets leads to the value 0.66 ± 0.02 (Chernokulsky and Mokhov, 2010).The simulated value for the scheme with the calibrated parameters set is very close to that for the version with the initial set of parameters.
When averaged over the Northern Hemisphere, total cloud fractions for each calendar month stay within the uncertainty range calculated from different empirical data sets (Fig. 1a).This is true even if reanalysis data are discarded and comparison is limited only to satellite data.The agreement is worse for the Southern Hemisphere, where total cloud fractions are underestimated throughout the year.For both hemispheres, our scheme correctly simulates minimum (maximum) cloud  fraction during the cold (warm) part of the year.However, the amplitude of the annual cycle for modelled c tot is greater than the satellite-derived one, especially in the Southern Hemisphere.
The scheme broadly reproduces the geographical pattern of cloud fractions.Similar to observations, annual mean total cloud fraction, c tot , attains maxima in northern and southern mid-latitudes, where c tot is typically between 0.7 and 0.9 (Fig. 2a and b).This is in general agreement with empirically derived values over oceans (Fig. 2c-f).However, over land our scheme with the initial parameter set overestimates total cloud fraction in this latitudes, since satellite-based data show smaller cloud fractions (from 0.5 to 0.7 from ISCCP and MODIS, and even from 0.3 to 0.7 from CERES).This bias is slightly diminished upon calibration.This is accompanied by reduced total cloud fraction over mid-latitudinal oceans, which worsens the agreement with observations.In the subtropics, the simulated total cloud fractions range from 0.1 to 0.5, which is too small in comparison to observations.Note that too-deep subtropical minima of c tot become shallower upon calibration.The fraction of convective clouds over the Indo-Pacific warm pool and over the Amazon Basin in our scheme (0.7 and larger) generally agrees with observations.
Basic conclusions made for the performance of the scheme for annual mean total cloud fractions may be translated to c tot   fields for individual months (Figs. 3 and 4).For all months, the scheme realistically reproduces total cloud fractions over mid-latitudinal oceans, but overestimates c tot over land at the same latitudes.That overestimate is more marked in winter than in summer, which is consistent with the overestimated amplitude of the annual cycle of c tot .Subtropical minima are too deep throughout the year.However, the scheme correctly places abundant convective clouds near the Equator in the winter hemisphere.
Comparison of the simulated cloud fractions in different layers with observations is not straightforward.The first reason for that is due to differences in the classification of cloud layers between the proposed scheme on the one hand and common satellite cloud products on the other.In our scheme, clouds belong to a particular layer depending on the height of cloud bases (see Sect. 2.2).As a result, convective clouds always belong to the lower layer in our scheme.This is in contrast with satellite retrievals, which classify clouds based on their tops.There, convective clouds may be classified to either low-, mid-, or upper-level clouds depending on the vertical extent of convective cloud ensemble.Another reason leading to difficulties in comparison of cloud fractions in individual layers is the above-mentioned "satellite   view" of cloud layers in common cloud satellite products (see Sect. 3.2).
However, some comparison may be performed with the results reported by Mace et al. (2009), who used the same classification scheme as we do for the merged lidar and radar observations from CALIPSO and CloudSat satellites.For this comparison, however, we have to keep in mind that Mace et al. (2009) reported only one year of measurements (from July 2006 to June 2007), which is quite different from the long-term climatology which we attempt to simulate here.For the reader's convenience, Fig. 5 is redrawn in the Supplement of the present paper (Fig. S1) in a fashion compatible with relevant figures from Mace et al. (2009).In turn, the latter figures are reproduced in Fig. S2 of the Supplement with permission of Wiley and Sons Inc.
In particular, our annual mean low-level cloud fraction c l = c sl +c co may be compared to their Fig.10a.Our scheme with the initial parameter set simulates c l between 0.6 and 0.7 over the mid-latitudes and in the areas of tropical convection, and from 0.1 to 0.5 in the subtropics (Fig. 5a).The largest c l (above 0.7) is simulated over the Antarctic.Upon calibration, c l in the northern and southern mid-latitudes is from 0.4 to 0.6, and over the Antarctic it is increased to values of 0.8 and larger.(Fig. 5b).In the middle latitudes, the calibration improves the agreement with the retrievals by Mace et al. (2009).In the subtropics, c l is somewhat increased, which again improves the agreement.Maxima of c l over the Indo-Pacific warm pool and over Amazonia become broader, which deteriorates our simulations.
One limitation of the present scheme is the lack of stratocumulus (Sc) decks over the eastern parts of the oceans.Annual mean stratocumulus cloud fraction in these regions fractions may be as large as 0.6 (Wood, 2012), and yields about 80-90 % of all low-level cloud fraction here.Our scheme produces low-level cloud fractions in these regions smaller than 0.2, which underestimate markedly the observed one.It is likely that this underestimate is due to neglecting the impact of atmospheric inversions on cloud formation.Such inversions suppress moisture fluxes from the planetary boundary layer to the free troposphere.In turn, under these conditions, vertical profile of specific (and relative) humidity may deviate strongly from the respective monthly averaged profile.An implementation of this impact may be one future improvement for our scheme.Note, however, that ERA-40 data underestimate the satellite-derived cloud fraction in these regions as well.This is an example of a problem that most contemporary cloudiness schemes in global climate models (GCMs) have in representing stratocumulus decks.In particular, Lauer and Hamilton (2013) reported that the latest generation of these models, the CMIP5 (Coupled models Intercomparison Project, phase 5) GCMs, underestimate the amount of subtropical stratocumulus decks by 30-50 %.
Mid-level cloud fractions c m ≡ c sm may be compared with Fig. 11 from Mace et al. (2009).In the version with the initial parameter set, mid-level cloud fractions range from 0.3 to 0.5 in mid-latitudes and the convective regions in the tropics (Fig. 5c).In other tropical and subtropical regions, c m is below 0.2 everywhere.Upon calibration, everywhere in the tropics and subtropics c m < 0.1, and in higher latitudes c m is between 0.1 and 0.2 (Fig. 5d).This drastically improves agreement with the hydrometeor fractions with bases from 3 to 6 km reported in Fig. 11 of Mace et al. (2009).
The modelled upper-level cloud fractions c h ≡ c sh markedly increase during calibration.In the version with the initial parameter set, annual mean c h is below 0.2 everywhere over the globe (Fig. 5e).Upon calibration, c h increases to 0.2-0.4 in the middle and high latitudes of the Northern and Southern Hemisphere as well as over convective regions in the tropics (Fig. 5f).As compared to Fig. 12a      the scheme's performance.In particular, extratropical upperlevel cloud fractions become broadly realistic, while there is an underestimation of c h in the areas of tropical convection by a factor of 2.

Cloud water path
Cloud water path (per model grid cell) W tot is markedly increased during calibration.In the initial version, globally and annually averaged W tot is equal to 66 g(H 2 O) m −2 , which is about one-half the respective value derived from CERES, 125 g(H 2 O) m −2 (Table 4).After calibration, modelled W tot increases to 82 g(H 2 O) m −2 , which is again too small in comparison to observations, but the agreement is better.However, we note that W tot in the CERES data is uncertain by a factor of 2. Further, the state-of-the-art general circulation models exhibit quite diverse simulation of this variable as well.In particular, the present-day globally averaged ensemble mean W tot for CMIP5 GCMs is 87 g(H 2 O) m −2 (which is close to our value), and the respective intermodel range is 37-167 g(H 2 O) m −2 (Lauer and Hamilton, 2013).The modelled cloud water path averaged over the Northern and Southern Hemisphere show maxima in summer (Fig. 6).
Calibration slightly decreases W tot in the extratropics throughout the year and markedly increases it in the tropics.Annual mean cloud water path in both versions of the scheme is from 20 to 80 g(H 2 O) m −2 (Fig. 7a).Over land it broadly agrees with the CERES data, while the cloud water path is too small over the storm-track-affected regions.Over oceans, it is a strong underestimate (Fig. 7c).In the tropics, calibration increases annual mean W tot by 20-50 %.As a result, the calibrated values of W tot in the tropics agree slightly better with the CERES data than the initial ones.In addition, W tot is too small in comparison to the CERES data.However, in these regions the CERES data suffer from large uncertainty (Minnis et al., 2011).
In winter, the cloud water path is severely underestimated, especially over land (Figs. 8 and 9).While one has to bear in mind that there is large uncertainty in the CERES retrievals in the high latitudes, an underestimate is clear in the middle latitudes.In summer, the mid-latitudinal cloud water path is somewhat small in comparison to the CERES data, but is reasonable as a whole.In contrast, in the tropics, W tot is somewhat too high, but the latter bias is markedly smaller than that in the middle latitudes in winter.
The largest contribution to W tot comes from low-level stratiform clouds during all seasons, and from mid-level stratiform clouds during the warm part of the year (not shown).In the tropics, the contribution from convective clouds is also valuable.
In the tropics, the underestimation of W tot is likely at least partly related to the above-mentioned lack of stratocumulus decks in the model.In the storm tracks, the respective underestimate is likely due to combination of the processes which are neglected in our scheme.First, geometric thickness of stratiform clouds is likely too small in our model.In particular, the typical thickness of low-level stratiform clouds h sl in middle latitudes is from 150 to 300 m.The latter is markedly smaller than (very limited) observational data summarised by Mazin and Khrgian (1989, p. 188), for which h sl ≥ 300 m.We note that low-level stratiform clouds are major contributors to W tot in the middle latitudes.This is similar for upperlevel stratiform clouds: in our calculations, h sh in middle latitudes is slightly larger than 100 m, while according to observations these clouds could be as thick as 1 km (Mazin and Khrgian, 1989, p. 189).Thus, the scheme might be improved by revising Eq. ( 5).Additional source of error in W tot is due to underestimated cloud fraction in the storm tracks (recall that our W tot is per grid cell rather than per cloudy part of the cell).Finally, the current version of the scheme completely lacks cloud-aerosol interaction, which increases clouds' lifetimes and therefore enhances their water content (Twomey, 1974;Albrecht, 1989;Hobbs, 1993;Lohmann and Feichter, 2005).
General underestimation of cloud water path in the tropics is common for the contemporary global climate models as well, e.g. for the CMIP5 GCMs (Jiang et al., 2012;Lauer and Hamilton, 2013).In the storm tracks, the same models display a diverse behaviour, with some models underestimating W tot , and others overestimating it.As a result, biases in our simulations are within the characteristic range of contemporary climate models.

Precipitation
Annual global precipitation changes insignificantly during calibration.In the version with the initial parameter set it is equal to 101 cm yr −1 , and in the calibrated version it is 100 cm yr −1 .Both values are within the range set by the GPCP data (88 cm yr −1 ) and by the ERA-40 (113 cm yr −1 ) (Table 4).The fraction of large-scale precipitation in the initial version is 0.48, which is an underestimate relative to the ERA-40 data (0.53).It becomes even smaller (0.45) after calibration (Table 4).
For monthly precipitation averaged over the Northern and Southern Hemisphere, both initial and calibrated versions reasonably agree with empirical climatologies (Fig. 10).The calibration enhances precipitation in the storm tracks, and suppresses it elsewhere.In the calibrated version, monthly precipitation in the Northern (Southern) Hemisphere changes from 7 cm mo −1 (6 cm mo −1 ) in winter to 11 cm mo −1 (14 cm mo −1 ) in summer.
Upon calibration, annual precipitation slightly decreases in the middle latitudes and in the monsoon-affected region and markedly increases in the tropics (Fig. 11a and b).In the calibrated version, precipitation P tot ranges from 90 to 180 cm yr −1 in the middle latitudes.This is a decrease by about one-fourth from the initial version.In turn, in the moist tropics, the calibrated precipitation ranges from 180 to 300 cm yr −1 , which is a respective increase by a factor 1.5.In dry subtropics, precipitation does not change markedly during calibration, being below 60 cm yr −1 .In most regions, the calibrated annual precipitation values agree better with the GPCP and ERA-40 climatologies than the initial ones.
One observes the marked decrease in the calibrated values relative to initial ones in the middle latitudes of the winter hemisphere (Figs.12a, b and 13a, b).In January, precipitation over eastern Eurasia is diminished from 2-5 cm mo −1 in the initial version to 1-3 cm mo −1 .The latter much better agrees with the empirical data in comparison to the former one (Fig. 12c and d).Over northern mid-latitudinal oceans, winter precipitation decreases by a factor of 2 or 3 during calibration.In the calibrated version it ranges from 4 to 16 cm mo −1 over northern mid-latitudinal oceans in January, and from 6 to 20 cm mo −1 over southern mid-latitudinal oceans in July.Both ranges are in agreement with empirical data (Figs. 12c,d and 13c,d).It should be noted that the above-mentioned severe underestimate of the fraction of stratocumulus decks in the subtropics should not severely affect simulation of precipitation because these clouds are non-precipitating ones (Houze, 1994).However, because our precipitation is somewhat too high in middle latitudes, and the cloud water path is too small there, it is likely that the calibrated lifetimes for stratiform clouds are too small, probably by a factor of 2.

Conclusions
This paper presents a scheme for calculation of the characteristics of multi-layer cloudiness and associated precipitation designed for climate models of intermediate complexity (EMICs).In contrast to the schemes previously used in such models, the scheme considers three-layer stratiform cloudiness and single-column convective clouds.It distinguishes between ice and droplet clouds as well.All main cloudiness characteristics (cloud fraction, cloud water path) are calculated interactively.Precipitation is calculated by using cloud lifetime, which depends on cloud type and phase as well as on statistics of synoptic and convective disturbances.
A novel approach for tuning this scheme was used.This approach was based on sampling of major governing parameters of the scheme.The corresponding cost function was constructed based on total cloud fraction, cloud water path, and precipitation, taking into account global mean values and annual mean, January, and July spatial distributions for these variables.Bayesian averaging was used to calculate the optimal parameters set.
After calibration, the scheme realistically reproduces the main characteristics of cloudiness and precipitation.The simulated globally and annually averaged total cloud fraction is 0.59, and the simulated globally averaged annual precipitation is 100 cm yr −1 .Both values agree with empirically derived values.
The scheme agrees with observations for total cloud fractions over mid-latitudinal oceans, but overestimates c tot over land at the same latitudes.The latter overestimate is more marked in winter than in summer.Subtropical minima are too deep throughout the year.The scheme correctly places abundant convective clouds near the Equator in the winter hemisphere, while it underestimates upper-level cloud fraction in the regions with strong convection.
Cloud water path is severely underestimated by the scheme.In particular, major storm tracks contain too little water.Cloud water path of tropical convective clouds is reproduced reasonably.However, we note that satellite retrievals for this variable are very uncertain, and the state-ofthe-art general circulation models exhibit quite a diverse simulation of W tot .
Calibration improves the simulation of total precipitation as well as the simulation of fraction of large-scale precipitation in total precipitation.However, important regional biases still remain in the scheme, e.g.too little precipitation in the tropical area.
Note that our calibration is not a simple "fitting exercise".In particular, cloud fractions at different layers were not trained explicitly during calibration.Nevertheless, they agree with available (rather limited) empirical data.This gives confidence in the physical basis of our scheme.Our scheme, while exhibiting some important biases, is a substantial improvement for Earth system models of intermediate complexity.However, regional and seasonal biases still present in the calibrated version show that there is a room for improvement of the scheme.One line of improvement may be the implementation of stratiform clouds originating from convective anvils in the upper troposphere in the presence of strong winds (Mazin and Khrgian, 1989;Houze, 1994).This may ameliorate too small upper-level cloud fractions in the tropical convective regions in the current version of the scheme.Another major improvement should be an implementation of stratocumulus decks representation, which should improve cloudiness over the eastern parts of the oceans.This implementation has to take into account inversion layers trapping convection and limiting vertical development of convective clouds (Wood, 2012).In addition, this version of the scheme lacks aerosol-cloud interaction (Twomey, 1974;Albrecht, 1989;Hobbs, 1993;Lohmann and Feichter, 2005).We note that one approach to include the latter in climate models of intermediate complexity has been developed by Bauer et al. (2008).An updated version of their scheme is planned to be implemented in our scheme in the future.

Fig. 1 .
Figure 1: 1 Fig. 1.Total cloud fractions (fraction) averaged over the Northern and Southern Hemispheres (a and b respectively) for the model w guess and calibrated parameter sets (gray and black lines correspondingly) as well as for the ISCCP, MODIS, CERES, and ERA-4 (red, yellow, green, and blue curves correspondingly).

Fig. 1 .
Fig.1.Total cloud fractions (fraction) averaged over the Northern and Southern Hemisphere (a and b, respectively) for the model with initial guess and calibrated parameter sets (grey and black lines respectively) as well as for the ISCCP, MODIS, CERES, and ERA-40 data sets (red, yellow, green, and blue curves respectively).

Fig. 2 .
Fig. 2. Annual mean modelled total cloud fraction (fractions) for initial and calibrated parameter sets (a and b, respectively) in comparison to the ISCCP D2, CERES, MODIS, and ERA-40 climatologies (c, d, e, and f, respectively).

Fig. 5 .
Fig. 5. Annual mean low-level cloud fraction c l ≡ c sl + c co (a and b), mid-level cloud fraction c m ≡ c sm (c and d), and upper-level cloud fraction c h ≡ c sh (e and f) for the model versions with initial (a, c and e) and calibrated (b, d and f) parameter sets.The units are fractions.
Figure 6: 6 water path (g(H2O)m −2 ) averaged over the Northern and Southern Hemispheres (a and b respectively) for the model with ibrated parameter sets (gray and black lines correspondingly) as well as for the CERES data set (green).

Fig. 6 .
Fig.6.Cloud water path (g(H 2 O) m −2 ) averaged over the Northern and Southern Hemisphere (a and b respectively) for the model with initial and calibrated parameter sets (grey and black lines respectively) as well as for the CERES data set (green).

Fig. 7 .
Fig. 7. Annual mean modelled total cloud water path ( g(H 2 O) m −2 ) for initial and calibrated parameter sets (a and b, respectively) in comparison to the CERES climatology (c).

Fig. 10 .
Fig.10.Total precipitation (cm mo −1 ) averaged over the Northern and Southern Hemisphere (a and b respectively) for the model with standard and calibrated parameter sets (grey and black lines respectively) as well as for the GPCP and ERA-40 data sets (magenta and blue curves respectively).
For globally and annually averaged cloud fraction S c,g = N c tot,g,ann,M ; c tot,g,ann,O , σ c tot,g,ann,O , ) Geosci.Model Dev., 6, 1745-1765, 2013 www.geosci-model-dev.net/6/1745/2013/ The lowermost level was located at the Earth's surface, and the next one was at H PBL .Other levels were equally spaced in height up to the tropopause.The latter was diagnosed from the monthly mean ERA-40 data using the conventional definition for thermal tropopause.For total cloud fractions, the following monthly climatologies were used: (Minnis et al., 2011)ection, 21 discrete computational levels were Geosci.Model Dev., 6, 1745-1765, 2013 www.geosci-model-dev.net/6/1745/2013/used.-TheCloudsand the Earth's Radiant Energy System (CERES)(Minnis et al., 2011).This data set was created by simultaneous retrievals of cloud properties and broadband radiative fluxes from the instruments on two LEO Terra and Aqua satellites from the Earth Observing System.The data from the Terra

Table 4 .
Globally and annually averaged values as calculated by the proposed scheme with two parameter sets in comparison with the available observational data.