## Abstract

To overcome the complexities associated with combining or comparing multisensor data, a statistical gridding algorithm is introduced for projecting data from their unique instrument domain to a uniform space–time domain. The algorithm has two components: 1) a spatial gridding phase in which geophysical properties are filtered on the basis of a set of criteria (e.g., time of day or viewing angle) and then aggregated into nearest-neighbor clusters as defined by equal-angle grid cells and 2) a temporal gridding phase in which daily statistics are calculated per grid cell from which longer time-aggregate statistics are derived. The sensitivity of the gridding algorithm is demonstrated using a month (1–31 August 2009) of level 2 *Aqua*/Moderate Resolution Imaging Spectroradiometer (MODIS) cloud-top pressure (CTP) retrievals as an example. Algorithm sensitivity is tested for grid size, number of days in the definition of a time average, viewing angle, and minimum number of observations per grid cell per day. The average CTP for high-level clouds from a number of different polar-orbiting instruments are compared on a 1° × 1° global grid. With the data projected onto a single grid, differences in CTP retrieval algorithms are highlighted. The authors conclude that this gridding algorithm greatly facilitates the intercomparison of CTP (or any other geophysical parameter) and algorithms in a dynamic environment. Its simplicity lends transparency to understanding the behavior of a given parameter and makes it useful for both research and operational use.

## 1. Introduction

Geophysical properties, such as atmospheric and surface temperatures, water vapor, aerosols, and clouds, have been retrieved operationally from satellite measurements in both geostationary (regional) and polar (global) orbits for more than three decades. While this potentially allows large-scale, long-term analysis, issues arise when comparing measurements from different time periods and sensors. Data combination and intercomparison efforts become difficult because of differences in technology, horizontal and vertical resolution, satellite overpass time, sampling frequency, and retrieval algorithm design strategies, among others. These difficulties are especially acute in climate studies that rely on decades of continuous global data records from a variety of satellite sensors.

A number of approaches to enhance data continuity exist. One is the creation of instrument-specific gridded aggregates. Historically, georeferenced calibrated sensor measurements at native instrument resolution are referred to as level 1B (L1B), the retrieved geophysical properties at instrument resolution as level 2 (L2), and the geophysical properties gridded and aggregated over a time period as level 3 (L3). The latter is useful, not only because of lower data volumes and improved data accessibility, but also as a statistical description of the retrieved properties. The Moderate Resolution Imaging Spectroradiometer (MODIS) L3 product, for example, contains hundreds of statistical data records on a 1° × 1° grid of L2 properties, such as aerosols, water vapor, clouds, and atmospheric profiles, in daily, 8-day, and monthly time aggregates (King et al. 2003). With close to 1000 records, however, the monthly product can easily become too cumbersome to work with. In addition, it may still not meet all research needs because of the static definition of grid size, time, and type of statistics, for example, mean instead of mode. There is also the possibility that a reprocessing of the entire product series becomes necessary as new knowledge on the statistical nature of a parameter (e.g., nonnormal behavior) is gained. This is both costly and complicated. This said, the MODIS L3 product is perhaps one of the most comprehensive to date.

Clouds play a major role in weather and climate patterns and are widely studied as indicators of change and trends (Ramanathan et al. 1989; Liou et al. 1993). The development of globally gridded cloud climatologies to date is either multi- or single sensor cloud properties sampled to a uniform grid. Cloud climatologies include, but are not limited to, the International Satellite Cloud Climatology Project (ISCCP; Schiffer and Rossow 1983), the University of Wisconsin High Resolution Infrared Radiometer Sounder (UW HIRS; Wylie and Menzel 1999; Wylie et al. 2005), the Pathfinder Atmospheres-Extended (PATMOS-x; Heidinger and Pavolonis 2009), the Television and Infrared Observation Satellite (TIROS)-N Operational Vertical Sounder (TOVS) Path-B measurements (Stubenrauch et al. 2006), and another based on TOVS (Susskind et al. 1997). While there may be uniformity of data in a single climatology, the same issues arise when comparing or combining different climatologies. Limitations in temporal scope, lack of transparency as well as consensus among research teams on gridding requirements (pertaining to grid size, type of statistic, time period, etc.) make data comparison efforts problematic.

There are examples in the literature where these stated problems are partly overcome by resampling different gridded cloud products to a common geographic grid using interpolation or zonal statistics (e.g., Hou et al. 1993; Wang and Rossow 1995; Jin et al. 1996; Stubenrauch et al. 1999 a,b,c; Hahn et al. 2001; Pavolonis and Key 2003; Thomas et al. 2004; Zhang et al. 2005). However, resampling a statistical dataset to a different space–time configuration introduces errors that are often untraceable in the final solution. Moreover, unknown errors may incorrectly depict trends.

In this paper, we propose a different approach. Instead of focusing on the development of a static global gridded data product, our goal is to develop a method with which any georeferenced dataset (independent of source instrument) of land, ocean, or atmosphere, at any level of processing (L1B/2), can be projected from its nonuniform instrument domain to a uniform space–time domain. As such, the gridded output is tailored to specific research needs, created for a user-defined (not product-defined) length of time, from any suite of instruments relevant to the study. To be clear, we do not attempt in any way to criticize the strengths and weaknesses of existing climatologies or propose the replacement of L3 products in this paper. Our intent is simply to offer a different approach to a long-standing problem.

The implementation of a gridding algorithm, as opposed to relying on static gridded products, has distinct advantages: (i) any geophysical property of interest at any level of processing can be projected to the space–time grid; (ii) all research relevant data are subject to the same set of filters and projected onto a common geographic grid; (iii) there are no predefined limitations on the type of parameter or statistic, grid size, time period, or number of gridded properties—these are all user-defined; (iv) the algorithm can be run as research needs arise, thereby reducing the need to aggregate and store L3 products; and (v) the algorithm is straightforward and transparent for application in both research and operational environments.

The gridding algorithm proposed here has two components. First is the “space” gridding phase where geophysical properties are filtered based on a set of criteria (e.g., time of day, viewing angle) and then aggregated into nearest-neighbor clusters as defined by equal-angle grid cells. Second is the “time” gridding phase where daily statistics are built per grid cell from which longer time-aggregate statistics are derived.

The data used in this study are discussed in section 2. The space–time gridding algorithm is described in section 3. For clarification, the gridding of L2 MODIS cloud-top pressure (CTP) retrievals is used as an algorithm application example throughout. However, the reader should bear in mind that CTP can be substituted with any L1B/2 geophysical parameter. Sections 4a and 4b detail the approach adopted in this paper to characterize the space–time algorithm and test its sensitivities to sample size. The latter is determined by gridcell size, aggregate time period, and viewing angle limitations. Discussion of results follows in section 4c and a comparison is made of CTP monthly gridded averages from four polar-orbiting instruments in section 5. Conclusions are reached in section 6 followed by recommendations for future research.

## 2. Data

All data processed in this paper are for the time period 1–31 August 2009 and from instruments on polar-orbiting platforms that obtain near global coverage in approximately two days.

### a. Multispectral imagers

The Advanced Very High Resolution Radiometer (AVHRR) on *Meteorological Operation-A* (*MetOp-A*) senses in three visible (VIS) to near-infrared (NIR) (0.5–1.64 *μ*m) and three IR spectral bands (3.55–12.5 *μ*m). The *MetOp-A* orbit has a local morning crossing time of 0930. AVHRR cloud properties are retrieved with the PATMOS-x cloud retrieval algorithm (Heidinger and Pavolonis 2009). Globally gridded PATMOS-x level 2 retrieval products are available on a 0.1° grid from the Cooperative Institute for Meteorological Satellite Studies (CIMSS; online at cimss.ssec.wisc.edu/patmosx). The PATMOS-x gridding algorithm employs nearest-neighbor logic with the result that there is only one CTP retrieval per 0.1° × 0.1° grid cell with no statistical data abstraction.

The MODIS on *Aqua* is sensing in 36 spectral bands ranging from 0.6 to 14.2 *μ*m. The *Aqua* platform has a local afternoon crossing time of 1330 and is part of the A-Train constellation of sensors, which include the *Cloud-Aerosol Lidar and Infrared Pathfinder Satellite Observations* (*CALIPSO*) and *CloudSat*. The operational 5-km level 2 cloud product (MYD06) is described in the literature (King et al. 2003; Platnick et al. 2003; Menzel et al. 2008). MODIS has a swath dimension of 2330 km with a viewing angle limit of 65.3°. This results in 406 level 2 pixels across track per scan line. The characterization and sensitivity studies in this paper are based on the MODIS level 2 cloud retrieval product (version 5). Data are available from the National Aeronautics and Space Administration (NASA) Goddard Space Flight Center’s level 1 and Atmosphere Archive and Distribution System (LAADS; online at ladsweb.nascom.nasa.gov).

### b. Hyperspectral sounder

The Atmospheric Infrared Sounder (AIRS) is mounted on *Aqua* and measures emitted infrared radiation in 2378 spectral channels spanning the range from 3.7 to 15 *μ*m (Aumann et al. 2003; Chahine et al. 2006). AIRS and MODIS are thus coincident in time. The level 2 retrieval algorithm is described by Susskind et al. (2003). AIRS radiance observations as well as the standard level 2 retrieval products (that include CTP retrievals along with surface properties and the atmospheric profiles of temperature and gases) are found at the NASA Goddard Earth Sciences Data and Information Services Center (GES DISC, available online at daac.gsfc.nasa.gov). AIRS CTP retrievals are available at 45 km × 45 km horizontal resolution, or roughly 0.4° × 0.4°. Unlike MODIS, AIRS has a viewing angle limit of 48.95°. This together with the reduced footprint size leaves 90 fields of view (FOV) across track per scan line.

### c. Lidar

The Cloud-Aerosol Lidar with Orthogonal Polarization (CALIOP; Winker et al. 2004, 2007), onboard the CALIPSO, is also an A-Train instrument. Its local crossing time lags 1–2 min behind that of *Aqua*. The level 2 cloud-layer product constitutes vertical profiles of clouds retrieved at 5-km resolution. However, unlike the passive sensors, there is only a single footprint across track. A strong advantage of CALIOP is that it has a very high vertical resolution of 60 m (Vaughan et al. 2004; available online at www-calipso.larc.nasa.gov/). This makes it ideal for validation of retrievals from passive sensors [see Stein et al. (2011) for a comparison of MODIS ice cloud properties with *CloudSat*].

## 3. Space–time algorithm description

From here on we refer to any space-based measurement of a geophysical parameter, be it radiances, retrieved clouds, surface, or layered thermodynamic properties, simply as “measurement.” Moreover, we make the distinction between “measurement” and “observation” since a geophysical parameter may not always be measured at each observation opportunity in native instrument resolution (or FOV).

The space–time algorithm described here allows daily measurements from any space-based instrument (polar-orbiting or geostationary) to be gridded to an *n*° × *n*° equal-angle grid (*n* > 0). A daily statistic is calculated for each grid cell from which a time aggregate can be generated to neutralize sampling differences among instrument measurements. The algorithm consists of two phases, the first of which is the physical data resampling phase, or “space gridding,” where measurements are projected from their instrument observation grid to a uniform equal-angle grid. All measurements are retained in this phase and clustered into nearest-neighbor sets according to the specified gridcell size. (Note that this is unlike the single value that is retained per grid cell in more traditional gridding routines.) Second, the statistical data reduction phase, or “time gridding,” is where parameter-specific daily statistics are calculated for each grid cell, or nearest-neighbor set. As soon as the desired number of days is accumulated, a time-aggregate statistical grid is calculated.

By dividing the gridding routine into two distinct phases, opportunity is created for statistical data exploration, calculation of compound statistics (i.e., combining more than one measurement in calculation of descriptive statistics), or higher-order statistical moments. In this sense the space–time grid is first and foremost an analysis space with the time-aggregate grid as the final output (or data product).

Below we discuss the space–time algorithm in more detail and illustrate aspects of it by using MODIS CTP retrievals (section 2a) as an algorithm application. Application to other geophysical parameters will be the subject of future work.

### a. Space gridding

One of the main advantages of gridding is data reduction. However, premature data reduction can severely limit or misdirect research. Thus, an important first step in any gridding routine is to carefully assess the research requirements and define rules for filtering the measurements. For example, in this paper we demonstrate space–time gridding for daytime high clouds retrieved from near-nadir observations. The corresponding filtering criteria are sun zenith angle ≤ 84° (daytime), CTP retrieval ≤440 hPa (high clouds), and satellite zenith angle ≤ 32° (near-nadir observations). Note that these criteria can be provided from any coincident ancillary datasets; they do not have to depend exclusively on the measurements. Once the filters have been applied, only the measurements retained are gridded. This is an important consideration: from this step onward any opportunity for inquiry that relies on measurements beyond those that have been filtered by the criteria will be lost. However, if a complete set of criteria is defined, then statistical analysis of the data in gridded space will be more descriptive and meaningful.

Once filtered, the measurements are projected from their native instrument resolution to a uniform gridded resolution. The method developed for the space–time algorithm is referred to as “snap-to-grid.” Each measurement is mapped to its nearest-neighbor (NN) grid cell (sample size is at least one measurement). This differs from the standard approach where each grid cell is assigned its NN measurement (sample size is one). The set of measurements accumulated per grid cell makes statistical data exploration possible. This is useful for gaining insight that will ultimately inform phase two of the algorithm (section 3b). Beyond the standard descriptive measures—such as mean, variance, mode, and spread—one can gain insight into probable (common) or even possible (rare) scenarios at this stage. Moreover, knowledge of the sample distribution plus degree of linearity and entropy per grid cell can greatly enhance gridded analysis.

The snap-to-grid routine is performed in three distinct steps. The geographic coordinate (*x*, *y*) in latitude–longitude degrees (with {*x* ∈ *R*: –90 < *x* 90} and {*y* ∈ *R*: –180 < *y* < 180}) of each measurement is 1) converted to a nonnegative real number set, *r*_{1} = *x* + 90, {*r*_{1} ∈ *R*: 0 < *r*_{1} < 180} and *r*_{2} = *y* + 180, {*r*_{2} ∈ *R*: 0 < *r*_{2} < 360}; 2) adjusted according to the gridcell size (gsize) in units (degrees) with {gsize ∈ *R*: 0 < gsize}, as follows = *r*_{1}/gsize and = *r*_{2}/gsize; and 3) rounded to the nearest integer such that { ∈ *Z*: 0 < < 180/gsize} and { ∈ *Z*: 0 < < 360/gsize} (*R* and *Z* are the real and integer number sets, respectively). Note that the output gridded coordinate point may not be unique. Several geographic coordinates may map onto the same gridded coordinate. As such, a set of nearest-neighbor measurements is accumulated for each grid cell without any distance calculation. This greatly simplifies the algorithm.

A single 3D grid is calculated per geophysical parameter per day. The grid has dimensions [nlat × nlon × nobs], where nlat and nlon indicate the latitude–longitude dimensions respectively, and nobs the number of observations. The magnitude of nlat and nlon are determined by the user-specified grid size; for example, gsize = 0.5° results in nlat = 360 and nlon = 720. The magnitude of nobs varies on the one hand with gsize, and on the other with source instrument spatial and temporal resolution. A coarser output grid resolution (larger gsize values) results in a higher nobs value. Alternatively, a coarser input grid resolution (larger FOV areas) results in a lower nobs value. Similarly, a higher frequency of measurements will result in a higher nobs value. The magnitude of nobs can be greater than or equal to the number of measurements (nmes) per grid cell (nmes ≤ nobs). The magnitude of nmes, in turn, is simply the number of nonzero values in the gridcell vector (with length nobs). The reason why nobs may in some instances not be equal to nmes is that a geophysical parameter is not always measured (retrieved) at each instrument observation; for example, a cloud parameter is not inferred for clear-sky conditions or a land surface temperature is not retrieved in the presence of clouds.

A daily grid has a finite set of cells where nobs = 0. This may be the result of the applied criteria (e.g., cells over polar nighttime regions will be empty if a daytime criteria was specified), or simply due to limitations in instrument coverage. A polar-orbiting instrument such as MODIS samples roughly 70% of the globe per day, but this depends on the swath width. The unsampled portion is what leads to empty cells in a daily grid.

Note that no statistical calculations are made at this phase of the space–time algorithm. It serves strictly as a physical data resampling routine. The only data reduction that occurs at this phase is due to the filtering routines. All data that pass the filters are resampled and stored for statistical calculation during the next phase. But all data that are filtered out are discarded, even the ancillary data are discarded. This simplifies data handling routines and speeds up enquiry during the next phase. For the algorithm application we put forward in this paper, 31 grids were computed, one for each day of August 2009 at 1° × 1° resolution.

### b. Time gridding

Statistical analysis is performed during this phase. Where measurements from different times of the day are accumulated per grid cell during phase one (e.g., as a result of two consecutive orbits, or inclusion of measurements from both descending and ascending orbits), they are statistically summarized into a single daily value during phase two. The advantage of having a set of measurements per grid cell should now become apparent; for instance, it allows the ability to calculate any number of daily statistics, for example, mean, mode, and standard deviation. Traditionally, daily gridded level 3 products contain no statistical information; they serve only a physical data reduction purpose. Statistical analysis with daily level 3 products becomes possible only once they are aggregated over time. We argue, however, that premature data reduction has the potential to not only limit scientific inquiry but also introduce small errors that enlarge over time. Careful handling of statistical data description from the outset has the potential to reduce data misrepresentation and enhance understanding of observed trends and patterns.

In calculating the daily statistical grids, the dimension of a space–time grid is initially reduced from [nlat × nlon × nobs] to [nlat × nlon]. However, the daily statistical grids stacked together form a second 3D grid with dimensions [nlat × nlon × ndays] (with ndays as the number of days). Once an aggregate time statistic is calculated per cell, the final output grid has dimension [nlat × nlon]. In the area covered by the measurements, the time statistic will be calculated from at least one daily statistic, that is, ndays ≥ 1. It remains to determine how the magnitude of ndays affects the final result.

## 4. Algorithm application and sensitivity analysis

### a. Weighted time average

The scope and flexibility of the space–time approach is illustrated by calculating a compound weighted time average of CTP. Note, however, that this is just an example since any descriptive statistic, model, or higher-order moment can be calculated for each daily grid cell since nobs *>* 1 for each populated cell. The compound weighted time average is provided by

where CTP_{m} is the weighted time average of CTP and CTP_{d} the daily average of CTP. In this case, the nmes ≤ nobs since a cloud property may not always be inferred for a given FOV (implying clear sky). This equation calculates a weighted time average by combining daily CTP averages with percentage cloud cover as defined by the fraction nmes/nobs. This follows the line of argument that uncertainty in the final output is reduced when the statistic is weighed by the fraction of occurrence. In other words, in the example presented here the time-aggregate CTP average per grid cell is not just an average of averages but reflects the percentage cloud cover such that the end result is dominated by CTP averages from predominantly cloudy days (thus minimizing contribution from predominantly clear days). An alternative approach could be to weight the daily averages by nmes only, which will give grid cells with a high number of measurements (i.e., those in the center of orbit tracks) a stronger weight in the time average. Far from being prescriptive, Eq. (1) is merely an illustration of the algorithm scope given the proposed gridding configuration.

### b. Sensitivity to sample size

Statistical analysis is sensitive to sample size; an average calculated from a small sample is more sensitive (less robust) to outliers than those calculated from a large sample. In gridded space it is important to maintain statistical robustness since the scale of geophysical analysis makes it easy for small, localized errors to go undetected initially but to accumulate over space and time into questionable trends.

We characterize this sensitivity to sample size by evaluating the space–time algorithm response to minimum sample size thresholds. In addition, we demonstrate how sample size, and by implication the space–time product, is affected by definitions of gridcell size, time period of statistical aggregate, and viewing angle limitations.

Given a polar-orbiting instrument configuration, grid cells on instrument swath edges (high viewing angles) contain fewer observations than those at swath center (near nadir). This can cause statistical errors that are exaggerated when days are aggregated into months. For example, a grid cell on the edge of an orbit swath for one day could have nobs = 3 and nmes = 1. The same grid cell could be in the center of an orbit swath for a different day, ending up with nobs = 100 and nmes = 12. Given Eq. (1), the daily weighted average calculated for the edge cell will make a larger contribution to the monthly CTP_{m} statistic than the center cell, even though a single measurement determined the average. This is addressed by setting a minimum limit for nobs to minimize the introduction of geographical misrepresentation (i.e., when a large area is assumed to be homogenously defined by a small isolated area). In single instrument studies, repeated experimentation can lead investigators to identify a static threshold for nobs, for example, 10. However, with multi-instrument comparisons a single static threshold can easily become meaningless since nobs varies greatly depending on the value of gsize and source instrument spatial and temporal resolution. Thus we address two questions. What is the minimum threshold for nobs to ensure a statistically robust sample size for a grid cell? And, is there an objective method with which to determine such a threshold?

Four different nobs thresholds are tested. They are the following: 1) nobs > 0 (all grid cells are retained no matter how low their nobs), 2) nobs > 10 (static threshold), 3) nobs > mean(nobs) − SD(nobs) [grid cells omitted with nobs lower than 16%—half of one standard deviation (SD)—of the mean nobs accumulated in a month], and last 4) nobs > mean(nobs) − [1.5 × SD(nobs)] (grid cells omitted with nobs lower than 3.5%—half of 1.5 SD— of the mean nobs accumulated in a month).

### c. Results and discussion

The space–time algorithm is demonstrated in this section using a specific application of CTP, namely daytime (sun zenith angle ≤ 84°) high cloud CTP retrievals (CTP ≤ 440 hPa) from near-nadir (viewing angle ≤ 32°) *Aqua*/MODIS measurements (section 2a) gridded to an equal-angle grid of 1° × 1° resolution and time averaged over the period 1–31 August 2009. Daily CTP averages are calculated that are then weighed based on the percentage cloudiness and averaged for a month according to Eq. (1). Even though this application is limited to daytime measurements that strictly do not constitute a full “day” of satellite observations (polar-orbiting sensors have both day- and nighttime orbits), we maintain reference to “daily grids” for the sake of simplicity. Note that all quantitative results communicated in this section refer to this application, unless otherwise specified.

On a single day, 288 files of 5 × 5 km MODIS measurements of potential CTP are generated. Each file contains 270 × 406 observations. Because of the polar-orbiting configuration of *Aqua*, the number of MODIS observations per equal-angle 1° × 1° grid cell has a strong latitudinal dependence; that is, the higher latitudes (>50°) are sampled more frequently than the lower latitudes (<50°). For the application, a grid cell at (0°N, 0°E) has a total of 3603 daytime MODIS observations for the month of August 2009, while a grid cell totals 4292 at (60°N, 0°E) and 12 180 at (80°N, 0°E). This could raise a concern that differences in sample size may cause higher degrees of smoothing in the poles than in the lower latitudes. We contend that a limit on the maximum sample size can reduce this problem. Such a limit is implemented in the space–time algorithm by restricting calculation of statistics to a single day. As a result, a time statistic then becomes the accumulation of single daily statistics for a single day and therefore limits the exponential accumulation of observations in the poles that can affect the statistical integrity of global results. Other factors that limit the number of observations accumulated per grid cell per day are the filters applied at the time of gridding, for example, time of day and viewing angle restrictions.

In addition to addressing the maximum sample size, limiting the minimum number of observations (nobs) per grid cell per day is also a useful tactic to prevent statistical instability. Three minimum threshold limits were tested for nobs (Table 1), such that a daily average was calculated for a grid cell if and only if the total observations on that day exceeded the threshold. With the average number of MODIS observations per grid cell per day about 280 (maximum = 501) on a 1° × 1° grid, a threshold of 10 is less than 1% of the nobs for most cells. The AIRS operational L2 product, on the other hand, has a coarser spatial resolution and averages 26 observations per cell per day (maximum = 50) on the same grid. A threshold of 10 is about 50% of the average and will thus severely hinder daily global analysis by removing too many cells from consideration. This demonstrates that a static threshold can severely compromise multisensor comparisons. However, a dynamic threshold that depends on a statistical relationship (e.g., tests 3 and 4 from Table 1) is independent of sensor resolution and thus more useful.

The total number of days in a month where nobs passed the threshold test is displayed for the three values in Figs. 1a–d, respectively. The higher the threshold, the fewer the available daily values; test 1 (that includes all cells in calculation) has overall the highest values for ndays, and test 3 (that filters out cells with the lowest 16% of nobs) overall the lowest ndays values. However, variation in ndays is not limited to the threshold magnitude since ndays can vary as much as four days from cell to cell across an orbit track in the low latitudes for a single threshold. This is the result of orbital sampling and shifting from day to day.

When compared to test 1, fewer than 27% of the cells of test 2 register any difference (Figs. 1a,b). However, given test 3, 53% of the cells have a difference of two or more days and 36% of the cells given test 4. Although the use of a minimum threshold reduces the number of days with which a time statistic can be calculated, the alternative of applying no threshold at all has a high probability of introducing geophysical errors in the results. Conceptually, the scenario in which a small spatial measurement (e.g., 5 km × 5 km) determines the grid value for a large uniform area of (e.g., 100 km × 100 km) should be avoided, and statistically, the reasons should become apparent in the discussions below.

The weighted monthly CTP average (Fig. 2a) follows the familiar geographic pattern also detected by Wylie and Menzel (1999) and Wylie et al. (2005). The highest clouds (>200 hPa) are found concentrated in the intertropical convergence zone (ITCZ) with high clouds between 200 and 440 hPa in the midlatitude storm belts of the North Atlantic, North Pacific, and Antarctic Oceans. No high clouds were detected for the month of August 2009 in the subtropical high pressure zones over the oceans and terrestrial subtropical deserts. These areas are typically characterized as having low cloud cover (Wylie et al. 2005).

The monthly average difference between each test 2–4 and test 1 is displayed in Figs. 2b–d. Test 2 has the lowest difference and test 3 the highest. Furthermore, it is clear that applying a minimum threshold filter removes most cells from the edges of orbit tracks (this is especially clear in Figs. 2c,d). To illustrate how a nobs threshold can affect the spatial coverage of a time statistic, daytime CTP averages were calculated for 81.6% of the grid cells with test 1 and 80.3% for test 3. Thus, although a minimum nobs threshold removes cells from consideration, the effect on the time-aggregate spatial coverage appears to be minimal.

Detailed differences in the time average for the four thresholds can be seen for the Indian Ocean (5°N, 40°E) to (15°S, 140°E) in Figs. 2e,f. The difference of test 3 and 4 with test 1 are displayed as either (i) increased (<0), (ii) decreased (>0), or (iii) remained the same (=0) in Figs. 2e,f. The diagonal stripes visible in Fig. 2f may be a result of the shifting orbits from one day to the next. Overall, though, the effect appears to be random and localized. These results are for a month of observations given a specific application. However, results may differ with different time intervals, statistics, geophysical parameters, source instruments, and grid resolutions. Despite that, it remains desirable to retain as many of the observations as possible, while removing outliers, to ensure a robust sample size in the calculation of a time statistic for each grid cell. Test 3 seems to remove too many days from the monthly average calculation (Fig. 1). Test 4 seems to be the most useful since it removes the outliers (cells with the lowest 3.5% of nobs) without severely compromising the total number of daily averages (Fig. 1). In addition, it is independent of instrument spatial resolution unlike test 2.

As was demonstrated in Figs. 2b–d, the monthly cell statistics can be sensitive to nobs thresholds that affect the sample size. For test 3, the difference with test 1 is as high as 288 hPa and globally averages 2 hPa (Table 1). Other factors that affect sample size in the CTP application is gridcell size, viewing angle limitations, and number of days in the time-aggregate statistic. The space–time algorithm sensitivity to variations in these three variables is discussed below.

The *Aqua*/MODIS has a 16-day repeat cycle. Given the Gregorian calendar, a month ranges between 28 and 31 days and therefore does not readily coincide with this cycle. We test here the sensitivity of monthly CTP gridded averages to two time aggregates: 32 days versus a month (31 days in the case of August).

The total number of observations per 1° × 1° grid cell is displayed for 31 (Fig. 3a) and 32 days (Fig. 3b), respectively. In the latter no orbit tracks are visible since the time period coincides exactly with two *Aqua* repeat cycles. This means that all grid cells in a latitude (zonal) band are sampled equally, with a total observations per grid cell that ranges from 4000 to 4500 in the low latitudes. This total number of observations per grid cell stays near 4000 for 31 days (Fig. 3a). However, upon closer inspection, orbital tracks are visible with a difference of up to 1000 (3500–4500) between adjacent grid cells. While it is tempting to want to adapt a multiple of the 16-day cycle for global statistics, the benefit is not apparent when the weighted CTP average is considered in Figs. 3c,d; the differences remain small with a global mean difference −0.22 hPa and standard deviation 17 hPa. These small differences are more evident when the view is zoomed in to a region over the eastern tropical Pacific (15°N, 140°W)–(15°S, 80°W) (Figs. 3e,f). Of the 1891 grid cells in this region, averages are calculated for the same number of cells with only a few cells registering a change in average value. Note that the regional difference visible here will drift from month to month because of satellite orbit changes. For the sake of spatial coverage and perhaps an ease of analysis description, Gregorian months are useful statistical units that do not introduce any obvious artificial features despite the initial uneven sampling.

The MODIS instrument has a scanning pattern of ±55°. When we restrict our analysis to only those measurements from observations within ±32° of nadir, we exclude a significant number of measurements. The rationale for excluding them is based on the fact that satellite observations at high viewing angles will be biased toward larger cloud fractions; this is caused by the increase in pixel size with viewing angle. Just as in Fig. 3, when the number of observations increases (as they do when the viewing angle restriction is removed), the spatial coverage increases and the time average decreases (because of larger sample sizes). This effect is noticeable in Fig. 4. The weighted time average of CTP is displayed for a region (30°N, 180°W)–(30°S, 80°W) with the near-nadir restriction (Fig. 4a) as well as for all viewing angles (Fig. 4b). The former has zero values for 13.5% of all cells with a spatial average of the remaining cells at 266 hPa. The latter, on the other hand, has 0 values for only 5% with a 250-hPa spatial average of the remaining cells. We conclude that the gridding algorithm is highly sensitive to viewing angle restrictions.

Although the monthly CTP average of all retrievals (with no view angle filtering) produces a smooth result without the enhancement of features due to orbital sampling or gaps (e.g., as displayed in Fig. 4a), it is unrealistic to assume that this could be applied to all instruments as the scanning pattern varies greatly from instrument to instrument (e.g., AIRS has a scanning pattern of ±48.95°). The goal with the space–time gridding algorithm is to produce a tool with which to process retrievals from any polar-orbiting instrument for easy multisensor comparisons. A near-nadir criterion (viewing angle ≤ 32°) is not only scientifically conservative but also applicable to any cross-track scanning instrument. But as demonstrated here, a limit of 32° may be too aggressive for CTP gridding alone (where cloud fraction is irrelevant). A value >32° but less than the smallest instrument viewing angle may serve a multisensor CTP application better.

Finally, the sensitivity of the gridding algorithm is tested to differences in grid size. The time average CTP on a 1° × 1° grid is compared to that on a 2° × 2° grid (Fig. 5). A larger grid cell will have a higher sample size, so a degree of smoothing will occur. This is noticeable when focusing the analysis on an eastern Pacific region (30°N, 180°)–(30°S, 80°W) (Figs. 5c,d). Small regional features disappear when the grid size is enlarged. However, on a global scale many of the large features are insensitive to grid size. The most obvious sensitivity to the grid size is the degree of spatial coverage. Most of the data gaps in Fig. 5a disappear when the grid size is doubled in Fig. 5b. Upon closer inspection, this has the additional effect that new features are introduced on a coarser scale that are not visible on a finer scale. We are in no position to make recommendations as to the accuracy of one over the other, but can merely conclude that globally gridded statistics are very sensitive to grid size. Careful attention should be paid to this when developing a test case or implementing a gridcell size for analysis in real time.

We have demonstrated here some of the spatial and temporal sensitivities of a new type of geophysical space–time gridding algorithm with MODIS CTP retrievals as an example application. Any gridding criteria or variable that causes the sample size to increase per grid cell will result in a lower dynamic range of the gridded measurement and smoother geophysical features in final product. The sensitivity of the algorithm to sample size demonstrates that the application of a general set of rules to the gridding of multiple geophysical parameters may not be appropriate. More meaningful gridded statistical data can be generated when a set of filtering criteria and gridding variables unique to the science application and geophysical parameter of interest is developed. For different teams to compare a given geophysical parameter in their products to investigate similarities and differences, the choices made in filtering the data must be clearly stated. This will ensure that the observed spatial features and trends are less due to introduced errors and more to the actual space–time properties. In the next section we demonstrate how this gridding algorithm can be applied to any number of instruments for a multisensor comparison that is independent of instrument configuration.

## 5. Multisensor gridding application

In the final part of our study, we investigate the suitability of using the space–time algorithm for a multisensor comparison by calculating a gridded monthly CTP average for four different polar-orbiting instruments: MODIS, AVHRR, AIRS, and CALIOP (section 2). This part of the experiment is limited to a visual comparison as the objective is to evaluate the effectiveness of the algorithm for global comparisons. In a future paper we will evaluate the strengths and weaknesses of various cloud retrieval algorithms and different polar-orbiting instruments on a uniform grid.

A weighted average of CTP [Eq. (1)] is calculated for 1–31 August 2009 on a 1° × 1° equal-angle grid. Test 4 was applied to remove grid cells with too few observations. The results in Fig. 6 represent the monthly average of daytime high cloud CTP retrievals from near-nadir observations (viewing angle ≤ 32°).

The intercomparison of cloud height retrievals from different instruments has been the subject of a number of studies (e.g., Wylie et al. 2005, Weisz et al. 2007). The results displayed here reflect the widely known capabilities of the different instruments. The infrared instrument, AIRS, appears to be detecting more high clouds than MODIS Collection 5 especially in the tropics. AVHRR/PATMOS-x replicates the AIRS tropical high cloud results but agrees more with MODIS in the midlatitudes and polar regions. Discrepancies between AIRS and MODIS CTP averages depend not only on the spectral region sampled but also on the spatial resolution of each instrument. The hyperspectral capabilities of AIRS make it possible to capture small vertical features, whereas the high spatial capabilities of MODIS make it possible to detect horizontal variability. The lidar instrument, CALIOP, is very sensitive to the detection of optically thin clouds, but its coverage is limited to a small number of grid cells in one month. CALIOP, however, offers opportunities for validation of cloud results, despite its limited horizontal range.

These results illustrate the ease with which the uniform space–time gridding algorithm can be applied to any number of instruments. A full gridded cloud retrieval comparison of different instrument algorithms will follow in future research.

## 6. Conclusions

The space–time gridding algorithm proposed in this paper provides a method for combining and comparing data from different space-based instruments on a common geographic coordinate system that is user defined. Gridding requirements for a geophysical parameter strongly depend on the research questions at hand. Static generalized gridded data products, such as level 3, from different sources may be incompatible and limited in scope. Moreover, in the quest for broad application, the creators of level 3 products may oversimplify geophysical data descriptions. With more than three decades of satellite data, the need for data combination and continuity on a large scale has never been greater.

The algorithm proposed here breaks the problem of gridding into two components: spatial data resampling and subsequently temporal data reduction. The requirements for both components are user defined—from the geophysical data filters to the grid size, type of statistic, and length of time. The implementation of a gridding algorithm, as opposed to using static gridded products for analysis, has the advantage that any geophysical property of interest at any level of processing can be gridded with research-specific statistics in mind. This greatly simplifies the gridding problem since key data reduction steps are left up to the user to define. Additionally, this provides the user with the flexibility to treat each geophysical parameter differently.

In addition to the dynamic nature of the algorithm, another novel characteristic is that all measurements retained by the filtering criteria are resampled during phase one. In other words, there is no spatial data reduction. Each grid cell contains a set of measurements from which a daily statistic can be calculated. This is an important consideration: the only type of data reduction in the space–time algorithm is statistically based. Traditionally, daily gridded products are created by selecting a single value per grid cell, whether it is based on a nearest neighbor or random location requirement. However, with such an approach the daily products serve a physical data reduction purpose. Statistical description is achieved only when daily grids are aggregated over time. This space–time algorithm, on the other hand, allows statistical description on a daily basis. The advantage of this approach is illustrated with the calculation of a weighted average. A cloud-top pressure (CTP) average is calculated per grid cell per day and weighed by the percentage cloudiness (number of cloudy versus clear measurements per grid cell). As a result, the time-aggregate average is dominated by the CTP average from predominantly cloudy days. This serves to reduce statistical uncertainty in the final output. Far from being prescriptive, this serves to demonstrate the scope of the space–time algorithm for the calculation of compound statistics and error handling.

Four different threshold values for a minimum number of observations per grid cell are tested. We conclude that a threshold that filters out those cells with fewer than 3.5% (the lowest of a 1.5 standard deviation) of the average number of observations per grid cell is rigorous but not too aggressive to hinder daily global analysis. Moreover it is instrument independent, which is an important consideration for multisensor comparisons.

The sensitivity of the space–time gridding algorithm is investigated for gridding variables that affect sample size, namely grid size, number of days, and viewing angle restrictions. This is approach is illustrated by means of an algorithm application using *Aqua*/MODIS cloud retrievals (MYD06, collection 5). A time average is calculated for daytime high cloud retrievals from near-nadir observations for the period 1–31 August 2009. Based on the results discussed in section 4, a number of conclusions are reached. Regional features in the CTP average are mostly artifacts of the cloud distribution on single days rather than irregular orbital sampling. Thus, any number of days (ndays) can be aggregated for calculation of a time statistic. This is demonstrated by comparing a time statistic calculated from a typical Gregorian month definition (31 days) with one calculated from 32 days (coinciding with two *Aqua* platform repeat cycles). The 31-day average did not introduce obvious artificial features despite strong orbital sampling differences.

The sensitivity of the algorithm to sample size demonstrates that the application of a general set of rules to the gridding of multiple geophysical parameters may not be appropriate. We argue that more meaningful gridded statistical data can be generated when a set of filtering criteria (e.g., viewing angle limits) and gridding variables (e.g., gsize and ndays) are developed that are unique to the science application and geophysical parameter at hand. This will ensure that the observed spatial features and trends are due less to introduced errors and more to the physical space–time properties.

Rather than producing static L3 gridded products, the space–time algorithm described here provides analysis flexibility with an opportunity for data exploration in grid space. This could enhance the interpretation and understanding of observed patterns and trends. Finally, the conceptual and computational ease of the algorithm lends it transparency to dynamically create gridded geophysical properties that are tailored to specific research needs from any suite of space-based measurements.

## Acknowledgments

This work is funded through NASA Grant NNX11AK22G. The views, opinions, and findings contained in this report are those of the author(s) and should not be construed as an official National Oceanic and Atmospheric Administration or U.S. government position, policy, or decision.

## REFERENCES

*CALIPSO*retrieval algorithms and data products.

*Laser Radar Techniques for Atmospheric Sensing,*U. N. Singh, Ed., International Society for Optical Engineering (SPIE Proceedings, Vol. 5575), 16–30.

*Laser Radar Techniques for Atmospheric Sensing,*U. N. Singh, Ed., International Society for Optical Engineering (SPIE Proceedings, Vol. 5575), 8–15.