Parameterizing cloud top effective radii from satellite retrieved values , accounting for vertical photon transport : quantification and correction of the resulting bias in droplet concentration and liquid water path retrievals

Droplet concentration (Nd) and liquid water path (LWP) retrievals from passive satellite retrievals of cloud optical depth (τ ) and effective radius (re) usually assume the model of an idealized cloud in which the liquid water content (LWC) increases linearly between cloud base and cloud top (i.e. at a fixed fraction of the adiabatic LWC). Generally it is assumed that the retrieved re value is that at the top of the cloud. In reality, barring re retrieval biases due to cloud heterogeneity, the retrieved re is representative of smaller values that occur lower down in the cloud due to the vertical penetration of photons at the shortwave-infrared wavelengths used to retrieve re. This inconsistency will cause an overestimate of Nd and an underestimate of LWP (referred to here as the “penetration depth bias”), which this paper quantifies via a parameterization of the cloud top re as a function of the retrieved re and τ . Here we estimate the relative re underestimate for a range of idealized modelled adiabatic clouds using bispectral retrievals and plane-parallel radiative transfer. We find a tight relationship between gre = r cloud top e /r retrieved e and τ and that a 1-D relationship approximates the modelled data well. Using this relationship we find that gre values and hence Nd and LWP biases are higher for the 2.1 μm channel re retrieval (re2.1) compared to the 3.7 μm one (re3.7). The theoretical bias in the retrieved Nd is very large for optically thin clouds, but rapidly reduces as cloud thickness increases. However, it remains above 20 % for τ < 19.8 and τ < 7.7 for re2.1 and re3.7, respectively. We also provide a parameterization of penetration depth in terms of the optical depth below cloud top (dτ ) for which the retrieved re is likely to be representative. The magnitude of the Nd and LWP biases for climatological data sets is estimated globally using 1 year of daily MODIS (MODerate Imaging Spectroradiometer) data. Screening criteria are applied that are consistent with those required to help ensure accurate Nd and LWP retrievals. The results show that the SE Atlantic, SE Pacific and Californian stratocumulus regions produce fairly large overestimates due to the penetration depth bias with mean biases of 32–35 % for re2.1 and 15–17 % for re3.7. For the other stratocumulus regions examined the errors are smaller (24–28 % for re2.1 and 10–12 % for re3.7). Significant time variability in the percentage errors is also found with regional mean standard deviations of 19–37 % of the regional mean percentage error for re2.1 and 32–56 % for re3.7. This shows that it is important to apply a daily correction to Nd for the penetration depth error rather than a time–mean correction when examining daily data. We also examine the seasonal variation of the bias and find that the biases in the SE Atlantic, SE Pacific and Californian stratocumulus regions exhibit the most seasonality, with the largest errors occurring in the December, January and February (DJF) season. LWP biases are smaller in magnitude than those for Nd (−8 to −11 % for re2.1 and −3.6 to −6.1 % for re3.7). Published by Copernicus Publications on behalf of the European Geosciences Union. 4274 D. P. Grosvenor et al.: Quantifying and correcting the effect of vertical penetration assumptions In reality, and especially for more heterogeneous clouds, the vertical penetration error will be combined with a number of other errors that affect both the re and τ , which are potentially larger and may compensate or enhance the bias due to vertical penetration depth. Therefore caution is required when applying the bias corrections; we suggest that they are only used for more homogeneous clouds.

Abstract.Droplet concentration (N d ) and liquid water path (LWP) retrievals from passive satellite retrievals of cloud optical depth (τ ) and effective radius (r e ) usually assume the model of an idealized cloud in which the liquid water content (LWC) increases linearly between cloud base and cloud top (i.e. at a fixed fraction of the adiabatic LWC).Generally it is assumed that the retrieved r e value is that at the top of the cloud.In reality, barring r e retrieval biases due to cloud heterogeneity, the retrieved r e is representative of smaller values that occur lower down in the cloud due to the vertical penetration of photons at the shortwave-infrared wavelengths used to retrieve r e .This inconsistency will cause an overestimate of N d and an underestimate of LWP (referred to here as the "penetration depth bias"), which this paper quantifies via a parameterization of the cloud top r e as a function of the retrieved r e and τ .Here we estimate the relative r e underestimate for a range of idealized modelled adiabatic clouds using bispectral retrievals and plane-parallel radiative transfer.We find a tight relationship between g re = r cloud top e /r retrieved e and τ and that a 1-D relationship approximates the modelled data well.Using this relationship we find that g re values and hence N d and LWP biases are higher for the 2.1 µm channel r e retrieval (r e2.1 ) compared to the 3.7 µm one (r e3.7 ).The theoretical bias in the retrieved N d is very large for optically thin clouds, but rapidly reduces as cloud thickness increases.However, it remains above 20 % for τ < 19.8 and τ < 7.7 for r e2.1 and r e3.7 , respectively.We also provide a parameteriza-tion of penetration depth in terms of the optical depth below cloud top (dτ ) for which the retrieved r e is likely to be representative.
The magnitude of the N d and LWP biases for climatological data sets is estimated globally using 1 year of daily MODIS (MODerate Imaging Spectroradiometer) data.Screening criteria are applied that are consistent with those required to help ensure accurate N d and LWP retrievals.The results show that the SE Atlantic, SE Pacific and Californian stratocumulus regions produce fairly large overestimates due to the penetration depth bias with mean biases of 32-35 % for r e2.1 and 15-17 % for r e3.7 .For the other stratocumulus regions examined the errors are smaller (24-28 % for r e2.1 and 10-12 % for r e3.7 ).Significant time variability in the percentage errors is also found with regional mean standard deviations of 19-37 % of the regional mean percentage error for r e2.1 and 32-56 % for r e3.7 .This shows that it is important to apply a daily correction to N d for the penetration depth error rather than a time-mean correction when examining daily data.We also examine the seasonal variation of the bias and find that the biases in the SE Atlantic, SE Pacific and Californian stratocumulus regions exhibit the most seasonality, with the largest errors occurring in the December, January and February (DJF) season.LWP biases are smaller in magnitude than those for N d (−8 to −11 % for r e2.1 and −3.6 to −6.1 % for r e3.7 ).
Published by Copernicus Publications on behalf of the European Geosciences Union.D. P. Grosvenor et al.: Quantifying and correcting the effect of vertical penetration assumptions In reality, and especially for more heterogeneous clouds, the vertical penetration error will be combined with a number of other errors that affect both the r e and τ , which are potentially larger and may compensate or enhance the bias due to vertical penetration depth.Therefore caution is required when applying the bias corrections; we suggest that they are only used for more homogeneous clouds.

Introduction
Clouds have a major impact on Earth's radiative balance (Hartmann et al., 1992) and small changes in their properties are predicted to have large radiative impacts (e.g.Latham et al., 2008).The amount of shortwave flux reflected by fully overcast warm (liquid water) clouds for a given sun and scattering angle, or the reflectance of a cloud, is primarily determined by the cloud optical depth (τ ), which in turn can often be characterized by the liquid water path (LWP; the vertical integral of liquid water content) and the cloud droplet number concentration (N d ).For a given cloud updraft, N d is determined by the number concentration and physicochemical properties of aerosols.Thus, couching cloud reflectance in terms of N d links the cloud albedo to aerosol and microphysical effects via the Twomey (1974) effect, making N d a very useful quantity to determine observationally.N d can also influence cloud macrophysical feedbacks via its control on rain formation (Albrecht, 1989;Stevens et al., 1998;Ackerman et al., 2004;Berner et al., 2013;Feingold et al., 2015) and stratocumulus cloud top entrainment (Ackerman et al., 2004;Bretherton et al., 2007;Hill et al., 2009).
Satellite observations of clouds and N d are immensely useful for studying clouds, cloud-aerosol interactions and for model evaluation since they afford large spatial and temporal coverage.A method to obtain N d from passive satellite observations (e.g. from the MODerate Imaging Spectroradiometer, MODIS; Salomonson et al., 1998) of τ and the cloud droplet effective radius (r e ) for stratiform liquid clouds has been previously demonstrated (Han et al., 1998;Brenguier et al., 2000;Nakajima et al., 2001;Szczodrak et al., 2001;Boers et al., 2006;Quaas et al., 2006;Bennartz, 2007;Grosvenor and Wood, 2014;Bennartz and Rausch, 2017) and is described further below.For more details see the Grosvenor et al. (2018) review paper on this technique, which also describes the known sources of error.In cloudy environments, aerosol optical depth cannot be retrieved from satellites, making cloud property observations such as N d and the cloud droplet effective radius (r e ) the only useful indicator of the influence of aerosol on clouds.An advantage of using N d rather than r e to study cloud-aerosol interactions is that r e is also determined by the cloud water content and thus is a function of cloud macrophysical properties.N d , in contrast, is only weakly controlled by cloud macrophysics, allowing some separation of microphysical and macrophysical effects.
However, retrievals of N d from space are still somewhat experimental and there is a lack of comprehensive validation of the retrievals and the assumptions required.There is a need to characterize and quantify the associated errors; in this paper we focus on doing this for one source of N d error using a 1-year N d data set for stratocumulus clouds from MODIS.
2 The adiabatic N d and LWP retrieval model and the vertical penetration depth bias N d and LWP are retrieved from passive satellite retrievals of r e and τ using an adiabatic cloud model that is described below.However, as shown in Platnick (2000) and Bennartz and Rausch (2017), for a retrieval free from other error sources (e.g.those due to cloud heterogeneity), the retrieved r e is representative of the r e value lower down in the cloud due to the vertical penetration of photons at the shortwave-infrared (SWIR) wavelengths used to retrieve r e .In contrast, the retrieved τ is comprised of contributions from the extinction coefficient β ext (h), where h represents height from cloud base, throughout the whole cloud profile: Here h = 0 represents cloud base and h = H is cloud top.β ext (h) is defined as where r is the droplet radius and n(r) is the droplet size number distribution within a cloud unit volume such that r) represents the ratio between the extinction and the geometric cross section of a given droplet and can be approximated by its asymptotic value of 2 (van de Hulst, 1957) since droplet radii are generally much larger than the wavelength of light concerned (typically 0.6 to 0.85 µm) such that the geometric optics limit is almost reached.
r e and liquid water content (LWC) at a given height are, respectively, defined as and where ρ w is the density of liquid water.Combining Eqs. ( 3) and ( 4) and inserting into Eq.( 2) gives To determine the form of r e (h) in the above equation in terms of L(h) and N d (h) we can utilize the fact that the "k" value, which is a measure of the width of the droplet size distribution (lower values indicate wider distributions), has been shown to be approximately constant in stratocumulus clouds (Martin et al., 1994;Pawlowska et al., 2006;Painemal and Zuidema, 2011).In this study we adopt a value of k = 0.72, which is the value assumed by the MODIS retrieval (Zhang, 2013).r v is the volume radius, defined as where we have used Eq. ( 4) to insert LWC and Eq. ( 6) to write r v as a function of k and r e .Now we utilize the assumptions that N d (h) is constant with height and that LWC(h) is a constant fraction, f ad , of the adiabatic LWC.The latter equates to where c w is the rate of increase of LWC with height (dLWC/dz, with units kg m −4 ) for a moist adiabatic ascent and is referred to as the "condensation rate" in Brenguier et al. (2000), or the "water content lapse rate" in Painemal and Zuidema (2011).See Ahmad et al. (2013) for a derivation.c w is a constant for a given temperature and pressure.Allowing these assumptions, using Eq. ( 7) to substitute for r e in Eq. ( 5) and combining with Eqs. ( 1) and ( 8) we can write At this stage, H * is any arbitrary height above cloud base and τ * is thus the optical depth between the cloud base and that height.H * can be expressed as a function of r e (H * ), k, N d and some constants by using Eqs.( 7) and (8).Then, given r e (H * ) and τ * , N d can be calculated as follows: Generally, when retrieving N d it is assumed that the r e obtained from satellite is representative of that from cloud top, i.e. r e (H * ) = r e (H ) (e.g.Bennartz, 2007;Painemal and Zuidema, 2011).This would then mean that τ * is the full cloud optical depth (τ ) as retrieved by the satellite and thus could be used in Eq. ( 10) above to obtain N d .However, since the r e obtained by satellite is actually equal to r e (H * ) then τ * < τ and thus τ * should be used in Eq. ( 10) instead of the retrieved τ ; the problem lies in the fact that τ * is unknown.However, in this paper we fit a simple function for τ * as a function of τ based on radiative transfer (RT) modelling of a variety of idealized clouds.
Alternatively, Eq. ( 10) can be formulated using the retrieved τ over the full cloud depth (setting τ * = τ ) and the cloud top r e (setting r e (H * ) = r e (H )).The problem then becomes one of estimating r e (H ) from the retrieved r e (H * ).
Here we formulate a parameterization of r e (H )/r e (H * ) as a function of τ .Note that either the τ or r e corrections should be applied to correct N d , but not both together.
Then we estimate the error introduced in N d retrievals for 1 year of MODIS data due to the usual assumptions of r e (H * ) = r e (H ) and τ * = τ , on the assumption that there are no other biases affecting the r e retrieval.We label this bias the "vertical penetration bias".
The method of correcting r e has the advantage over the τ correction since it also allows a correction to the retrieval of LWP.LWP can be estimated (see e.g.Szczodrak et al., 2001) using LWP = 5 9 ρ w r e (H )τ. (11) For a corrected LWP the cloud top r e and the retrieved (total) τ values should be used.Since the retrieved r e (H * ) is likely to be underestimated due to the vertical penetration depth bias, LWP would otherwise be underestimated and the correct value can be obtained by using the parameterized r e (H ) instead.
3 Data and methods

Calculation of τ and r e corrections
In order to calculate and we have performed r e retrievals on idealized clouds using a similar algorithm to that used for MODIS retrievals.We produced idealized clouds that span a large range of stratocumulus-like clouds as represented by combinations of N d and LWP.We chose 41 values between N d = 10 and 1000 cm −3 that were equally spaced in log space and 91 values between LWP = 20 and 200 g m −2 spaced equally in linear space.All of the possible combinations from this sampling were used to sample the 2-D (N d , LWP) phase space.
For each combination, discretized adiabatic model profiles following the form of those described in Sect. 1 (i.e. with a vertically constant N d and LWC that increases linearly with height) were generated using c w = 1.81 × 10 −6 kg m −4 , f ad = 0.8 and a vertical spacing of 1 m.The droplet size distributions at each height were represented by a modified gamma distribution with a k value of 0.72, i.e. representative of an effective variance of 0.1.One-dimensional RT calculations, assuming plane-parallel clouds, were performed on these profiles using DISORT (Discrete Ordinates Radiative Transfer Program; Stamnes et al., 1988) radiation code in order to simulate reflectances at wavelengths of 0.86, 2.1 and 3.7 µm, matching those measured by MODIS to retrieve τ and r e over an ocean surface.Note that MODIS provides r e retrievals using both 2.1 and 3.7 µm wavelengths, which are hereafter referred to as r e2.1 and r e3.7 , respectively.The MODIS r e3.7 retrieval requires a correction to account for the contribution to the observed radiance from thermal emission, which is based on the observed 11 µm radiance (Platnick and Valero, 1995;King et al., 2015;Platnick et al., 2017).We account for this in our retrievals by removing the thermal contribution during the RT calculation instead of via the 11 µm radiance, which should produce a consistent end result.The RT calculations were performed assuming a black surface, a clear atmosphere (i.e.gaseous absorption is neglected), using a solar zenith angle (SZA) of 20 • and a nadir viewing angle.These reflectances were then used to retrieve τ and r e values using the Nakajima and King (1990) bispectral method, as operationally used by MODIS.To do so, a lookup table was built from reflectances calculated for a range of clouds that were assumed to be plane-parallel in nature, as assumed for the operational MODIS retrievals; i.e. these clouds were uniform in the vertical and horizontal with infinite horizontal extent.Again, a black surface and a k value of 0.72 were assumed along with the same viewing geometry as for the RT calculations on the adiabatic clouds.A fixed depth of 1 km was assumed with cloud base at an altitude of 1 km and cloud top at 2 km, although the cloud depth has no major effect on the reflectances generated for a given τ and r e .g re was then calculated using the retrieved and model top r e values.dτ was calculated by choosing the value from the model profile of τ , as measured from cloud top downwards, that corresponded to the location where the model profile r e matched the retrieved r e .
Figure 1a shows a 2-D histogram of g re values as a function of τ for the 2.1 µm retrieval.It shows that when plotted in this way g re forms a fairly tight relationship with τ so that for a given τ only a small range of g re values are possible.This suggests that the relationship can be parameterized based upon a 1-D relationship fitted to these data with little loss of accuracy.The median value of each τ bin is also plot- The black line is the median g re in each τ bin after smoothing over τ interval windows of 0.2.The white line is the fit to the mean curve using Eq. ( 14).ted (after smoothing over τ windows of 0.2) and this is the relationship used in this paper.g re is seen to decrease with τ with a gradient that decreases with τ .Similarly, Fig. 2 shows dτ vs. τ , which also shows a tight relationship that is suited to a 1-D parameterization.Fourth-order polynomial curves can be fitted (using the least squares method) to the median value relationships that take the form and The coefficients of these fits are given in Tables 1 and 2 along with the maximum errors for the fit (relative to the mean or median line) for the range shown.The curves (white lines in Figs.1a and 2a) fit the mean data well with maximum absolute differences of 0.001 and 0.09, respectively, for the g re and dτ curves.However, there will be some error when using this relationship (or the mean value relationship) due  Figures 1b and 2b show the same results for the 3.7 µm retrieval.Again tight 1-D relationships are suggested.Here, though, g re and dτ values are lower for a given τ and the curves are steeper at lower τ values, but flatten off much more rapidly.By τ = 7.5 there is little dependence of dτ on τ and dτ saturates at a mean value of ∼ 2.6.The fit estimates for the curves (Eqs.14 and 15 and Tables 1 and 2) again match the actual curves closely with a maximum absolute error in g re and dτ of 0.003 and 0.14, respectively.

MODIS data
For the MODIS data we use 1 year (2008) of MODIS Aqua data and follow a similar methodology to that used in Grosvenor and Wood (2014) in order to create a data set akin to the MODIS Level-3 (L3) product (King et al., 1997;Oreopoulos, 2005).We processed MODIS Collection 5.1 joint Level-2 (L2) swaths into 1 • × 1 • grid boxes.Joint L2 swaths are subsampled versions of the full L2 swaths (sampling every fifth 1 km pixel) that also contains fewer parameters.We process the data from L2 to L3 in order to allow the filtering out of data at high SZAs and to provide both r e2.1 and r e3.7 retrievals.
For this work we relax the screening methodology slightly from that used in Grosvenor and Wood (2014) since here we are interested in the effects of the vertical penetration N d bias upon a more general global data set.We applied the following restrictions to each 1 • ×1 • sample that goes into the daily average (since multiple overpasses per day are possible) in order to attempt to remove some artifacts that may cause biases: 1.At least 50 joint L2 1 km resolution pixels from the MODIS swath that did not suffer from sunglint were required to have been sampled within each grid box.2. At least 80 % of the available (non-sunglint) pixels were required to be of liquid phase based upon the "primary cloud retrieval phase flag".Analysis was only performed on these pixels.A high cloud fraction helps to ensure that the clouds are not broken, since broken clouds are known to cause biases in retrieved optical properties due to photon scattering through the sides of clouds.Often retrievals of N d are restricted to high cloud fraction fields for this reason (Bennartz, 2007;Painemal and Zuidema, 2011) and so we focus on such data points here.
3. The only pixels used were those remaining after (2) for which the "cloud mask status" indicated that the cloud mask could be determined, the "cloud mask cloudiness flag" was set to "confident cloudy", successful simultaneous retrievals of both τ and r e for the 2.1µm channel were performed and the cloud water path confidence from the MODIS L2 quality flags was designated as "very good confidence" (the highest level possible).This is a little different from the official MODIS L3 product where a set of cloud products are provided that are weighted using the quality assurance (QA) flags.Rather than weighting our L3-like product with the QA flags we have simply restricted our analysis to pixels with the highest confidence for water path.

The mean 1
is restricted to values lower than 3.2 km.This is done to avoid both deeper clouds for which N d retrievals are likely to be problematic due to the increased likelihood of a breakdown of the assumptions required to estimate N d , such as a constant fraction of the adiabatic value for LWC and vertically constant N d , and increased retrieval issues due to cloud heterogeneity.CTH is calculated from the MODIS 1 • × 1 • mean cloud top temperature (CTT) and the sea surface temperature (SST) using the method of Zuidema et al. (2009).SST data were obtained from version 2 of the NOAA Optimum Interpolation (OI) Sea Surface Temperature data set (NOAA_OI_SST_V2) that provides weekly SST data at 1 • × 1 • resolution.This was interpolated to daily data on the assumption that SST does not vary significantly over sub-weekly timescales.
5. The mean 1 • ×1 • SZA was restricted to ≤ 65 • following the identification of biases in the retrieved τ , r e and N d at high SZAs (Grosvenor and Wood, 2014).

Only 1 •
× 1 • grid points with mean τ > 5 were considered for the N d data set due to larger uncertainties from instrument error and other sources of reflectance error for τ and r e retrievals at low τ (Zhang and Plantnick, 2011;Sourdeval et al., 2016).
Following this screening, the 1 • × 1 • grid boxes associated with each MODIS Aqua overpass were averaged into daily mean values for ocean covered surfaces only.Figure 3a shows the number of days from the year of data examined in this study (year 2008) that fulfilled the above criteria and thus are likely to produce a good N d retrieval.Regions with high numbers of days where useful N d retrievals can be made have been selected for closer examination in this study; they are listed in Table 3 along with information on the mean and maximum numbers of days of good data.The permanent marine stratocumulus decks are among those selected, namely those in the SE Pacific off the western coast of S. America (Region no.1), in the SE Atlantic off the western coast of southern Africa (Region no. 2), off the coast of California and the Baja Peninsula (Region no.3), in the Bering Sea off the SW coast of Alaska (Region no. 6) and in the Barents Sea to the north of Scandinavia (Region no.8).These regions are where the highest numbers of selected days occur with values ranging up to a maximum of 141 days (for the Bering Sea region).The Barents Sea region has the lowest maximum number of days out of this group, reflecting the fact that N d retrievals cannot be made during a lot of the winter season in this region due to a lack of sunlight.The Southern Ocean (Region no. 5) and the NW Atlantic (Region no.7) regions frequently produce stratocumulus, although it is often associated with the cold sectors of cyclones and so its location from day to day is more transient.These regions are also affected by high SZAs in the winter seasons, which also restricts the number of retrievals possible there.The East China Sea region (Region no.4) produces the lowest mean and maximum numbers of days since the stratocumulus areas are mostly restricted to near the coast and occur mostly in the winter season.
N d was calculated for both r e2.1 and r e3.7 using Eq. ( 10) from the 1 • × 1 • daily mean τ , r e and CTT.This was done by using the retrieved τ value in Eq. ( 10) along with both the retrieved r e value (i.e.assuming that r e (H * ) = r e (H ) as is often assumed for N d retrievals) and by estimating r e (H * ) using the retrieved r e along with the g re values that were calculated as described above.This therefore gives N d data sets for the "standard" method and a corrected method, allowing the differences between the two to be examined.A similar process was applied for the LWP retrieval.

Results
Following Eq. ( 10), the ratio between the uncorrected and corrected N d values can be shown to be Figure 4 shows how the relative N d bias varies as a function of retrieved τ when using an r e that has been corrected using the g re from Fig. 1 (mean curve, black line).At τ = 5 the relative error is 46 % for the r e2.1 retrieval and 28 % for the r e3.7 retrieval.At higher τ the errors reduce rapidly but remain above 10 % for the r e2.1 retrieval over the τ range shown.For the r e3.7 retrieval the relative error drops below 10 % for τ >∼ 13.Thus, the overall degree of error due to this effect will be determined by the distribution of τ for the regions of interest, which we take into consideration here using MODIS data for a representative N d data set.
Alternatively, if the correction is formulated in terms of a correction to τ we obtain The equation shows that, for a constant dτ , the relative N d bias due to an uncorrected τ value would increase with decreasing τ as τ − dτ approaches zero.e in Eq. 10) to those from the corrected calculation (using the corrected r e for r * e as calculated from the retrieved value and g re ; the g re values used are those shown by the black line in Fig. 1) vs. retrieved τ .
Figure 3b shows the time-mean τ for the data set as filtered by criteria 1-7 above, i.e. to replicate the type of filtering that would likely be performed for N d retrievals.Table 3 lists the regional means of these time-mean values along with the regional means of the standard deviations of τ over time.It shows that the mean τ values of the tropical and subtropical regions are generally lower than those at higher latitudes.The East China Sea, Barents Sea, NW Atlantic and Southern Ocean regions exhibit the highest mean τ values out of those examined and so should be expected to show the lowest N d biases due to the vertical penetration effect.The SE Atlantic region (and the region to the west of Africa in general) show low τ and can be expected to give high N d biases.Table 3 also lists the fraction of days for which τ ≤ 10 (f τ ≤10 ).τ = 10 is the value above which N d biases drop below 31 % for the 2.1 µm retrieval and below 14 % for r e3.7 according to Fig. 4.
Thus f τ ≤10 indicates the fraction of days for which daily N d biases will be greater than 31 % for that channel.The values in the table indicate that even in the least affected region (Barents Sea) this will occur for 21 % of the days.For the SE Atlantic and SE Pacific region the percentages rise to 69 % and 53 % of the days, respectively.Thus, the vertical penetration depth N d bias is prevalent in all regions for which N d data sets are likely to be used and particularly so in the subtropical stratocumulus regions where N d retrievals have been widely used and studied.
The overall bias is now estimated using 1 year of actual MODIS data in order to obtain a realistic distribution of τ values.However, it should be noted that the data set used is deliberately filtered in order to only retain data points that are likely to give useful N d data, namely low liquid clouds www.atmos-meas-tech.net/11/4273/2018/Atmos.Meas.Tech., 11, 4273-4289, 2018 Table 3. Regional statistics for the various marine stratocumulus regions shown in Fig. 3. Shown are the mean and maximum number of days that fulfill the screening criteria in order to be considered as useful N d retrievals; the regional means and standard deviations (σ ) of the time-averaged optical depths (τ ) for the screened data set; and the regional mean of the fraction of days for which τ ≤ 10 (f τ ≤10 ), which is calculated using only data from grid points for which the number of days with N d data was ≥ 15.
No with extensive 1 • × 1 • cloud fractions, i.e. predominately stratocumulus.This is done in order to assess biases for the types of clouds that N d data sets will typically be used to study.
Figure 5 shows a map of the mean percentage biases and Table 4 gives the regional means of the values in the map.Considering firstly the biases for the r e2.1 retrieval, the biases are highest in the tropics and subtropics.The regional mean bias is 34.5 % for the SE Atlantic region (Region no. 2), which is the stratocumulus region that seems to suffer the most.The biases are a little lower for the other major stratocumulus regions; e.g. for the SE Pacific region (Region no. 1) and the Californian region (Region no. 3) the mean biases are 32 %, although the biases increase further west, where the dominant cloud regime tends to shift towards trade cumulus clouds.The remaining regions all have mean biases of 24-28 %.The Barents Sea region (Region no.8) has a value of only 23.7 %, representing the stratocumulus region with lowest mean bias.These results indicate higher τ values for the clouds in the East China Sea, Southern Ocean, Bering Sea, NW Atlantic and Barents Sea regions relative to the Californian and SE Pacific stratocumulus regions, with the SE Atlantic region exhibiting the lowest τ values.This is confirmed by the mean τ values shown in Table 3.The biases for the r e3.7 retrieval display the same spatial patterns as for r e2.1 , but are significantly lower; the mean value in the region with the maximum bias (SE Atlantic, Region no. 2) is 17 % and in the region with the lowest bias (Barents Sea, Region no.8) it is 10 %.
The regional mean LWP biases are also listed in Table 4.They are negative since an r e underestimation from the vertical penetration effect leads to an LWP underestimate (see Eq. 11).The biases are also smaller in magnitude than for N d due to the smaller sensitivity of LWP to r e inherent in the latter equation.They are anticorrelated with the N d biases such that the region with largest N d bias (SE Atlantic) has the largest negative LWP bias of −11.1 %.The smallest magnitude bias occurs in that Barents Sea region (−8 %).
It is also useful to know how variable the biases are from day to day for a given point in space since this will determine how useful the application of a single offset bias correction might be for correcting N d biases for daily data.Figure 6 shows the time variability of the bias in the form of the relative standard deviations (over time) of the percentage N d biases.It reveals that the percentage bias in N d generally has a larger relative standard deviation at latitudes above around  40 • with values typically ranging up to around 30-50 % (of the mean percentage N d bias) for the 2.1 µm retrieval.Relative variability is greater for the 3.7 µm retrieval, perhaps due to the much lower mean percentage errors.Some of the selected regions show more variability than others, in particular the Barents Sea and East China Sea regions.
Table 4 gives the regional means of the relative standard deviations revealing values that range from approximately 20 to 40 % of the mean percentage biases for the 2.1 µm retrieval and 30-60 % for the 3.7 µm one.This shows that the application of a single annual mean offset bias correction is likely to lead to fairly large biases for the N d estimates for individual days for regions where the mean N d errors are significant.If daily data are used to determine relationships between cloud properties and N d without correcting for the biases examined here then significant variability in N d might be introduced that may affect those relationships via non-linear effects.
Figure 7 shows how the percentage N d biases change with season for the r e2.1 retrieval only.Interestingly, the highest biases tend to occur in the DJF season for the SE Pacific and SE Atlantic stratocumulus regions, indicating that τ values are lower in DJF for those seasons.The September-October-November (SON) season also generally produces higher biases than March-April-May (MAM) and June-July-August (JJA) for those regions, particularly for SE Pacific.For the East China Sea region the biases are lower in SON and DJF seasons than in the other seasons.We note that there are little data in this region for JJA since there are few low-altitude clouds with large regional liquid cloud fractions there in this season.Either the other regions do not show a large amount of seasonal variability or N d data are only available for part of the year due to a lack of sunlight in the winter months.

Discussion
There are some caveats to the results that we presented here that we now discuss.We have shown that, theoretically, the  Platnick (2000) showed that dτ has some dependence on viewing geometry and so the consideration of a wider range of view and SZAs should ideally be made.
In addition, a liquid water condensation rate (c w ) value of 1.81 × 10 −6 kg m −4 was assumed for the model adiabatic clouds, which corresponds to a cloud temperature of 278 K at a pressure of 850 hPa.In reality, cloud temperatures and hence c w will vary, mainly as a function of cloud temperature.We have performed sensitivity tests using a value of 1.0 × 10 −6 kg m −4 , which corresponds to a cloud temperature of 262 K.This is likely to be close to the coldest temperature attained by boundary layer clouds over the oceans that are coupled to the surface, which are generally the types of clouds for which the droplet concentration retrievals are applied.The results (not shown) reveal mean differences (across the τ values tested) in the mean dτ line (i.e. the white line in Fig. 2a) of 4.2 % (the maximum difference was 10 %) for the 2.1 µm retrieval and 3.6 % (maximum of 7.6 %) for the 3.7 µm one.For g re the differences were of the opposite sign and much smaller, with mean differences in the median g re line of −0.45 % (the most negative difference was −0.6 %) for the 2.1 µm retrieval and −0.3 % (−0.66 % was the most negative difference) for the 3.7 µm one.Therefore, the effect of c w changes is relatively minor.These results would also apply for equivalent changes in the cloud adiabaticity (i.e. the value of f ad ).
The modelling of the idealized clouds and the correction rests on the assumption that r e increases monotonically with height within the cloud (following the adiabatic assumption), but there is some suggestion that the development of precipitation-sized droplets might lead to larger droplets being preferentially found below cloud top (Chang and Li, 2002;Nakajima et al., 2010a, b;Suzuki et al., 2010).However, Zhang et al. (2012) found that MODIS retrievals of r e performed on model-generated clouds were not significantly affected by the presence of precipitation.Also, during the VOCALS field campaign in the SE Pacific region, aircraft observations showed that r e generally did increase with height up to cloud top (Painemal and Zuidema, 2011), indicating that this is not a problem at least for the near-coastal clouds tested.Further offshore the likelihood of precipitation increases as clouds become more cumulus-like and so for those clouds the issue may be greater and hence more caution should be exercised when interpreting the results presented here for such regions.
Evaporation effects related to entrainment also have the potential to reduce r e , N d and LWC near cloud top and hence negate some of the assumptions upon which the N d retrievals rest.However, we argue that the entrainment effect upon r e is likely to be minimal for two reasons: firstly, the evidence suggests that for stratocumulus clouds extreme inhomogeneous mixing occurs at cloud top, which reduces the LWC and N d , but does not change r e (Burnet and Brenguier, 2007;Brenguier et al., 2011;Painemal and Zuidema, 2011).Secondly, the results of Painemal and Zuidema (2011) indicate that entrainment occurs within approximately the first 0.5 optical depths from cloud top on average; the penetration depths calculated here are considerably larger than this for reasonably thick clouds (Fig. 2).The effect of the reduced N d and LWC within the entrainment zone is not so clear-cut; this would negate the assumption of a vertically constant N d and monotonically increasing LWC used to formulate the total τ .However, given the likely small τ contribution from the entrainment region relative to the total τ , this effect is likely to be small.
It is also clear that the suggested correction for the vertical penetration effect should only be applied to the retrievals of N d with consideration of other bias sources.These other potential error sources are numerous and include r e biases due to subpixel heterogeneity (Zhang and Plantnick, 2011;Zhang et al., 2012Zhang et al., , 2016)), 3-D radiative effects (Marshak et al., 2006), assumptions regarding the degree of cloud adiabaticity (f ad in Eq. 10; Janssen et al., 2011;Merk et al., 2016), the choice of k value (assumed constant; Brenguier et al., 2011;Merk et al., 2016), the assumption of a vertically uniform N d , the assumed droplet size distribution shape and width (Zhang, 2013), viewing geometry effects (Várnai and Davies, 1999;Horváth, 2004;Várnai and Marshak, 2007;Kato and Marshak, 2009;Liang et al., 2009Liang et al., , 2015;;Di Girolamo et al., 2010;Maddux et al., 2010;Liang and Girolamo, 2013;Grosvenor and Wood, 2014;Bennartz and Rausch, 2017) and upper-level cloud and aerosol layers (Haywood et al., 2004;Bennartz and Harshvardhan, 2007;Davis et al., 2009;Meyer et al., 2013;Adebiyi et al., 2015;Sourdeval et al., 2013Sourdeval et al., , 2016)).These errors have the potential to bias N d in a way that opposes the positive bias expected from the vertical penetration effect such that the overall biases may cancel out.Indeed, the largest source of error in N d is likely that from r e biases given the sensitivity of N d to r e in Eq. ( 10).MODIS r e has generally been shown to be biased positively compared to aircraft observations (Painemal and Zuidema, 2011;King et al., 2013), which would lead to a negative N d error when taken alone.Thus, the application of the correction described in this paper in isolation has the potential to enhance any negative bias in N d caused by a positive r e bias.
Our paper quantifies the vertical penetration bias in isolation to the other effects mentioned above.It should be questioned, though, whether the presence of cloud heterogeneity and other effects somehow prevent the effects of the vertical stratification from influencing the retrieved r e , making it irrelevant.This could be a potential explanation for why it is often observed that r e2.1 is larger than r e3.7 (Zhang and Plantnick, 2011) in contrast to the direction expected from adiabatic clouds given the vertical penetration effect, since it is known that subpixel heterogeneity effects tend to cause a positive r e2.1 bias relative to r e3.7 (Zhang et al., 2012).We argue, though, that the vertical stratification effect occurs in addition to other effects (e.g.heterogeneity) with the latter cancelling out and often exceeding the former such that the positive difference between r e2.1 and r e3.7 would be even larger without the vertical stratification effect.The cancellation of biases may also explain why VOCALS aircraft measurements (Painemal and Zuidema, 2011) tended to show that r e2.1 and r e3.7 were very similar.
We also note that there are many situations when the expected result due to vertical stratification of r e does occur (i.e.r e3.7 > r e2.1 ), as demonstrated in Painemal et al. (2013) and Fig. 8.This shows ratios between r e3.7 and r e2.1 for an example MODIS scene in the SE Pacific stratocumulus region.Ratios using the uncorrected MODIS r e values are shown, which shows that the ratio exceeds 1 for most of the stratocumulus cloud region (the clouds that adjoin the coast) with ratios ranging from around 1.1 to 1.2.In the more broken clouds the ratio is less than 1, which is likely a result of cloud heterogeneity.However, it would be expected that N d retrievals would not be applied to such clouds.The figure also shows the ratios calculated using r e3.7 and r e2.1 values that have been corrected using the g re factors.If the differences between r e3.7 and r e2.1 were caused by vertical stratification alone and our parameterization were correctly predicting the cloud top r e for both MODIS channels, then this ratio should be equal to 1.This is the case for the clouds close to the coast, indicating that our parameterization is working well for these clouds.The ratios are a little higher than 1 further north and west (around 1.05-1.08)indicating that either our parameterization is not working correctly for these clouds or other factors are causing relative differences between r e3.7 and r e2.1 (e.g.subadiabaticity, cloud heterogeneity).Figure 9 shows the percentage of pixels for which r e3.7 > r e2.1 for 90 days of MODIS SE Pacific observations divided into four different heterogeneity bins.Heterogeneity is characterized by the H σ parameter (Liang et al., 2009), which is the standard deviation of the 250 m resolution 0.86 µm reflectance (R 0.86 ) divided by the mean R 0.86 .It is clear that for many regions relative r e values that are consistent with an adiabatic profile occur more than 50 % of the time, particularly when the cloud heterogeneity is low.This suggests that it may be possible to use H σ to determine the situations in which the bias correction is more applicable.However, it is hard to definitively prove our argument within the scope of this study, particularly for more heterogeneous regions, since it would likely require computationally expensive 3-D RT modelling of known cloud fields (e.g. from LES models), followed by r e and τ retrievals.Painemal and Zuidema (2011) actually demonstrated that MODIS N d agreed rather well with N d from aircraft for the SE Pacific region despite a fairly large positive r e bias; this was thought to be due to the fortuitous cancellation of (for N d ) the r e bias with biases in the k parameter and f ad .How- using uncorrected re values.Right: using re values that have been corrected using the parameterizations for gre (for both re3.7 and re2.1).A ratio of one is expected for the plot on the right if the relative differences between re3.7 and re2.1 are caused by vertical stratification alone and if the parameterization is correctly predicting the relative differences.
region#1) for the 2.1 µm retrieval, assuming perfect initial agreement.This indicates that another N d bias may have been operating in order to give the good observed agreement.
The MODIS retrieval uses reflectances from both a visible and a shortwave infra-red (SWIR) wavelength channel with the former being primarily determined by τ and the latter by r e .However, a bi-spectral retrieval is used and so there is also some 10 sensitivity of the retrieved τ to the SWIR reflectance, which will be representative of the r e below cloud top due to the vertical penetration effect.This, combined with the fact that the MODIS forward retrieval model assumes a vertically uniform cloud, will result in the retrieved τ being biased relative to the real value (assuming the real cloud has an adiabatic profile).Figure 10 shows the difference between the retrieved and model profile τ ; the bias is negative and smaller in magnitude than 5 % for the 3.7 µm retrieval.They are slightly larger for the 2.1 µm retrieval, but still lower in magnitude than 5 %, except at r e 7 µm.
Although it should be noted that some of this bias may be due to other causes related to the inconsistencies between the vertically uniform and adiabatic models, rather than the r e vertical penetration bias.Since the retrieved N d is proportional to the square root of τ , this will lead to small N d biases.Biases in LWP will be similar to those in τ since LWP is proportional to 5 τ , but r e biases are still likely to dominate (e.g., see Fig. 1).Thus, we have not pursued this further.
In this paper we have only considered retrievals over the ocean, although retrievals over land for τ and r e are available for MODIS.MODIS surface albedo uncertainties are likely to be much higher over land than over the oceans (King et al., 2004;Rosenfeld et al., 2004;Bréon and Doutriaux-Boucher, 2005) since the surface albedo is much more variable over land.In Using r e values that have been corrected using the parameterizations for g re (for both r e3.7 and r e2.1 ).A ratio of 1 is expected for the plot on the right if the relative differences between r e3.7 and r e2.1 are caused by vertical stratification alone and if the parameterization is correctly predicting the relative differences.ever, the agreement between aircraft and MODIS N d seen in Painemal and Zuidema (2011) would deteriorate if a correction for the N d bias due to the penetration depth effect discussed here were also applied.Table 4 indicates that the result would be a MODIS N d underestimate of around 32 % (average for SE Pacific, Region no. 1) for the 2.1 µm retrieval, assuming perfect initial agreement.This indicates that another N d bias may have been operating in order to give the good observed agreement.

% of Re37>Re21 pixels
The MODIS retrieval uses reflectances from both a visible and a SWIR wavelength channel, with the former being primarily determined by τ and the latter by r e .However, a bispectral retrieval is used and so there is also some sensitivity of the retrieved τ to the SWIR reflectance, which will be representative of the r e below cloud top due to the vertical penetration effect.This, combined with the fact that the MODIS forward retrieval model assumes a vertically uniform cloud, will result in the retrieved τ being biased relative to the real value (assuming the real cloud has an adiabatic profile).Figure 10 shows the difference between the retrieved and model profile τ ; the bias is negative and smaller in magnitude than 5 % for the 3.7 µm retrieval.They are slightly larger for the 2.1 µm retrieval, but still lower in magnitude than 5 %, except at r e 7 µm.Although it should be noted that some of this bias may be due to other causes related to the inconsistencies between the vertically uniform and adiabatic models rather than the r e vertical penetration bias.Since the retrieved N d is proportional to the square root of τ , this will lead to small N d biases.Biases in LWP will be similar to those in τ since LWP is proportional to τ , but r e biases are still likely to dominate (e.g.see Fig. 1).Thus, we have not pursued this further.
In this paper we have only considered retrievals over the ocean, although retrievals over land for τ and r e are available for MODIS.MODIS surface albedo are likely to be much higher over land than over the oceans (King et al., 2004;Rosenfeld et al., 2004;Bréon and Doutriaux-Boucher, 2005) since the surface albedo is much more variable over land.In addition, cloud masking is more difficult over land, particularly over non-vegetated surfaces, transitional areas between desert and vegetated surfaces and above high-altitude regions (Platnick et al., 2003).We have ignored land regions in order to avoid such complications and also because stratocumulus clouds are more prevalent over ocean regions (Klein and Hartmann, 1993;Wood, 2012).However, the results shown in this paper may still apply over land.The results of Rosenfeld et al. (2004) and Platnick et al. (2017, their Fig. 14) suggest that surface albedo uncertainties are more important at lower optical depths ( 5) and for the 2.1 µm retrieval (relative to the 3.7 µm one).Thus, for thicker clouds and the 3.7 µm retrieval land surface albedo issues may be less problematic.
Finally, we note that the thermal emission correction for the MODIS r e3.7 (see Sect. 3.1) retrieval has some uncertainty that should be considered; the uncertainty for this is included (combined with other uncertainties) in the MODIS Collection 6 pixel level uncertainty products (Platnick et al., 2017).It is possible that effects additional to those included, such as cloud heterogeneity and surface heterogeneity, may further increase the uncertainty beyond that estimated in the MODIS products, but these are currently not well documented.

Conclusions
We have described and quantified a positive bias in satellite retrievals of cloud droplet concentration (N d ) and liquid water path (LWP) that make use of the adiabatic cloud assumption to estimate these quantities from satellite-observed cloud optical depth (τ ), effective radius (r e ) and cloud top temperature.We term the bias the "vertical penetration bias".It arises due to the well-documented vertical penetration of photons with wavelengths in the shortwave-infrared range into the upper regions of clouds, so that r e retrievals are representative of values some distance below cloud top (Platnick, 2000;Bennartz and Rausch, 2017) rather than being those at cloud top as assumed by the N d and LWP retrievals.Here we quantified the optical depth as measured from cloud top downwards, dτ , at which the retrieved r e equaled the actual r e for adiabatic clouds covering a large range of total cloud optical depths and N d values.We showed that knowledge of dτ allows a corrected N d to be calculated by subtracting dτ from the observed τ and using that in the N d retrieval instead of τ .We characterized dτ as functions of τ for the 2.1 and 3.7 µm r e retrievals (r e2.1 and r e1.6 , respectively) and found that a 1-D relationship approximates the modelled data well.dτ increases with τ and is larger for r e2.1 than for r e3.7 and so the vertical penetration N d bias affects retrievals based on r e2.1 more than those using r e3.7 .Similarly, we also parameterized the true cloud top effective radius (r e (H )) as a function of the retrieved r e and τ , allowing both a corrected N d and LWP to be calculated by using r e (H ) instead of the retrieved value.Both the dτ and r e correction methods give similar results for the N d retrievals suggesting that the latter is preferable since it also allows for a correction to LWP.However, for some applications it may be useful to be able to parameterize dτ .
We quantified the vertical penetration N d bias for 1-year N d and LWP data sets.The corrections presented here suggest that N d and LWP errors will increase as the τ value of the cloud scene gets lower.For many regions that are considered trustworthy for N d and LWP retrievals (typically stratocumulus regions), there are high frequencies of low τ values and the N d biases are significant.For example, for the SE Pacific and SE Atlantic regions clouds with τ ≤ 10 (for which N d errors are expected to be ≥ 31 % for r e2.1 and ≥ 15 % for r e3.7 ) occur, respectively, 53 and 69 % of the time on average.The mean r e2.1 vertical penetration N d biases for these regions were 32 and 35 %, respectively.Out of the stratocumulus regions examined, these two were the worst affected.For r e3.7 the N d biases were much smaller; for example, mean biases for the SE Pacific and SE Atlantic regions were 15 and 17 %, respectively.N d biases were predicted to be worse for the tropical and subtropical regions than for higher latitudes.The time variability of the biases were also examined and were shown to be significant (regional mean standard deviations of 19-37 and 32-56 % for r e2.1 and r e3.7 , respectively).This indicates that long-term averages of the vertical penetration N d bias corrections are not useful for correcting N d data over short timescales (e.g.daily N d data).We also examined the seasonality of the N d biases and showed that, for the stratocumulus regions, generally the DJF season was worst affected, followed by SON.
LWP biases were of a lower magnitude than those for N d and were negative.The largest biases were again for the SE Atlantic region where the mean bias was −11.1 % and the smallest for the Barents Sea region (−8 %).Biases were also lower when using r e3.7 with a maximum (most negative) bias of −6.1 %.
We caution that the correction for the N d and LWP vertical penetration biases presented here should only be considered in combination with corrections for other biases that affect τ and r e .Zhang et al. (2016) suggest a correction for the subpixel heterogeneity bias effect, but corrections may not currently exist for all biases and it is likely that some unidentified biases still exist.Therefore, we recommend that our correction is currently only applied to homogeneous cloud scenes in order to minimize possible entanglements with biases resulting from heterogeneity effects, which are not accounted for.Such conditions can be obtained by limiting retrievals to associated heterogeneity (H σ ) values (available in MODIS MYD06 Collection 6 products) to less than about 0.1.Otherwise N d and LWP biases could be made worse, for example, in situations where the fortuitous cancellation of opposing errors leads to initially small N d errors.The latter was suspected to have occurred for the comparison between MODIS N d retrievals and in situ aircraft observations as presented in Painemal and Zuidema (2011).We showed that the SE Pacific, which is the region examined in that study, had a mean vertical penetration depth error of 32 %, suggesting that another unidentified N d bias may have been operating in order to give good agreement.
Previous studies have shown that r e3.7 is less prone to biases due to subpixel averaging (Zhang and Plantnick, 2011;Zhang et al., 2012Zhang et al., , 2016)).Thus, combined with the work presented here, this supports the conclusion that r e3.7 likely represents a better choice for use in N d and LWP retrievals.
For future work, it is recommended that additional characterization of dτ and g re is performed for a range of viewing geometries in order to ensure that the results presented here are robust for all cloud retrievals.The use of 3-D radiative transfer calculations and simulated retrievals upon known LES model fields would also be useful for investigating how heterogeneity effects might interact with the vertical penetration effects.Further investigation into how the presence of precipitation affects our assumptions and results is also warranted.

Figure 1 .
Figure1.Two-dimensional histogram of g re as a function of τ for a range of clouds (see text) for the 2.1 µm r e retrieval (a) and the 3.7 µm retrieval (b).The black line is the median g re in each τ bin after smoothing over τ interval windows of 0.2.The white line is the fit to the mean curve using Eq.(14).
(Eq. 15)  to estimate the mean dτ value as a function of τ .The maximum absolute error between the fit and the mean line is also shown. in the g re and dτ values seen in the underlying histograms.

Figure 2 .
Figure 2. As for Fig. 1 except for dτ as a function of τ and using Eq.(15) for the white line.

Figure 3 .
Figure 3. (a) Number of days in 2008 that fulfilled the criteria required to be counted as a valid N d retrieval.See the text for details on the criteria.Various regions of interest are also denoted by the boxes and numbers.(b) Mean optical depth for data set with filtering criteria 1-7 applied (see text).

Figure 4 .
Figure4.The ratios of N d values from the standard MODIS calculation (using the retrieved r e for r * e in Eq. 10) to those from the corrected calculation (using the corrected r e for r * e as calculated from the retrieved value and g re ; the g re values used are those shown by the black line in Fig.1) vs. retrieved τ .

Figure 5 .
Figure 5. Maps of the annual mean percentage error for uncorrected N d retrievals using a year (2008) of daily MODIS data that have been filtered to select data points in which N d retrievals are favourable and therefore most likely to be used for N d data sets (see text for details).The left plot shows the results for the r e2.1 retrieval and the right for the r e3.7 retrieval.

Figure 6 .
Figure 6.As for Fig. 5 except showing the relative (as a percentage) standard deviation of the percentage N d bias over time, i.e. σ % bias /% bias.

Figure 7 .
Figure 7. Seasonal mean percentage N d biases for the r e2.1 retrieval only.

Figure 8 .
Figure 8. Ratios of re3.7 to re2.1 for a MODIS snapshot scene from the Southeast Pacific stratocumulus region from 16th June, 2015.Left:

20Figure 8 .
Figure 8. Ratios of r e3.7 to r e2.1 for a MODIS snapshot scene from Southeast Pacific stratocumulus region from 16 June 2015.(a) Using uncorrected r e values.(b)Using r e values that have been corrected using the parameterizations for g re (for both r e3.7 and r e2.1 ).A ratio of 1 is expected for the plot on the right if the relative differences between r e3.7 and r e2.1 are caused by vertical stratification alone and if the parameterization is correctly predicting the relative differences.

Figure 9 .
Figure 9.The percentage of pixels for which r e3.7 > r e2.1 for 90 days (January, February and March of 2008) of 0.1 • resolution MODIS Collection 6 observations for the SE Pacific stratocumulus region.Only single layer liquid clouds are included and data points have been filtered to exclude τ < 5 and partially cloudy pixels.The four panels are for four different bins of the heterogeneity parameter (the standard deviation of the 250 m resolution 0.86 µm reflectance divided by the mean reflectance) with bin ranges labelled in square brackets above the panels and the x and y axes in degrees longitude and latitude, respectively.

Figure 10 .
Figure 10.Percentage τ bias (retrieved vs. actual value from the input model profile) as a function of τ and r e .

Table 1 .
Coefficients for the fitting curve (Eq.14) to estimate the median g re value as a function of τ .The maximum absolute error between the fit and the median line is also shown.

Table 2 .
Coefficients for the fitting curve

Table 4 .
. Region name Mean no.days Max no.days Mean τ σ τ f τ ≤10 Regional means of the predicted time-mean percentage biases in N d and LWP due to the vertical penetration depth error and regional means of the relative (percentage) standard deviations (over time) of the percentage N d and LWP biases (i.e.regional means of the values in Fig.6and the equivalent for LWP).Bias results are shown for both the 2.1 and the 3.7 µm r e retrievals.