Validation of a modiﬁed AVHRR aerosol optical depth retrieval algorithm over Central Europe

. The Advanced Very High Resolution Radiometer (AVHRR) carried on board the National Oceanic and Atmospheric Administration (NOAA) and the Meteorological Operational Satellite (MetOp) polar orbiting satellites is the only instrument offering more than 25 years of satellite data to analyse aerosols on a daily basis. The present study assessed a modiﬁed AVHRR aerosol optical depth τ a retrieval over land for Europe. The algorithm might also be applied to other parts of the world with similar surface characteristics like Europe, only the aerosol properties would have to be adapted to a new region. The initial approach used a relationship between Sun photometer measurements from the Aerosol Robotic Network (AERONET) and the satellite data to post-process the retrieved τ a . Herein a quasi-stand-alone procedure, which is more suitable for the pre-AERONET era, is presented. In addition, the estimation of surface reﬂectance, the aerosol model, and other processing steps have been adapted. The method’s cross-platform ap-plicability was tested by validating τ a from NOAA-17 and NOAA-18 AVHRR at 15 AERONET sites in Central Europe


Introduction
Aerosols do not only affect the climate (IPPC, 2007), but also have a major influence on visibility (Horvath, 1995) and as air pollutants can become hazardous for human health (Samet et al., 2000). Spaceborne techniques have helped to increase our knowledge of the spatial and temporal distribution of aerosols and to learn more about their physical and chemical properties over the past years. Nowadays, various satellite sensors are used to retrieve different aerosol characteristics (e.g., Yu et al., 2006;Kokhanovsky et al., 2007;Kokhanovsky and De Leeuw, 2009).
Remote sensing of aerosol properties with the Advanced Very High Resolution Radiometer (AVHRR) has a long tradition (Knapp and Stowe, 2002), even though aerosol detection was not the initial aim of this instrument. Nonetheless, the AVHRR, carried on board the National Oceanic and Atmospheric Administration (NOAA) and the Meteorological Operational Satellite (MetOp) polar orbiting satellites, is the only instrument offering the opportunity to analyse more than 25 years of satellite data in moderate Published by Copernicus Publications on behalf of the European Geosciences Union. 1256 M. Riffler et al.: Modified AVHRR AOD retrieval over land spatial resolution on a daily basis. Most studies of aerosol properties from AVHRR focused on areas over the oceans (Mishchenko et al., 1999;Geogdzhayev et al., 2002;Ignatov et al., 2004;Mishchenko and Geogdzhayev, 2007), where the spectral properties of the underlying surface are generally well-known and the surface reflectance is low. One of the first AVHRR retrievals over land (Soufflet et al., 1997) used the dark dense vegetation (DDV) method suggest by Kaufman and Sendra (1988) to derive the aerosol optical depth τ a over boreal forest for a limited data set. The disadvantage of this method is the low occurrence of DDV; e.g., in Europe less than 1% of the pixels (Borde et al., 2003) contain DDV. The retrieval over brighter and heterogeneous land surfaces is more complicated and prone to introduce larger errors due to reduced aerosol signal and temporally unstable surface reflectance. Knapp and Stowe (2002) derived τ a from the reduced resolution (110 × 110 km 2 ) Pathfinder-Atmosphere (PATMOS) data set using also bright surface targets. Hauser et al. (2005a) applied a similar technique to derive τ a at full resolution (1.1 × 1.1 km 2 ) for Central Europe using NOAA-16 AVHRR and they further qualitatively compared monthly and seasonal means from this product with MODIS collection 004 data (Hauser et al., 2005b). A full resolution τ a climatology from AVHRR covering land surfaces is still missing.
This study is a step towards such a climatology for Europe. In the AVHRR algorithm from Hauser et al. (2005a) a relationship between the satellite data and ground-based Sun photometer measurements from the Aerosol Robotic Network (AERONET) was established to post-process the data resulting in a better retrieval quality. Applied to NOAA-17 and NOAA-18 AVHRR this post-correction did not lead to the same improvement. Therefore, we developed and assessed a quasi-stand-alone algorithm which is not dependent on AERONET data for correction. However, optical aerosol properties like single scattering albedo and scattering phase function are derived from AERONET inversion products, since a single-channel approach allows only to derive one piece of aerosol information independently. After analysing calibration stability of the historical AVHRR time series from our archive dating back to 1985, such a procedure may be more suitable for the pre-AERONET era. Moreover, certain retrieval steps from the original algorithm have been adapted, e.g., assumption on the aerosol model and estimation of the surface reflectance. The cross-platform applicability of the method is tested using two satellites with overlapping time period, namely NOAA-17 and NOAA-18 AVHRR, and the retrieved τ a is validated with Sun photometer measurements of AERONET. Furthermore, the performance of the AVHRR aerosol retrieval is related to products from the newer generation Medium Resolution Imaging Spectrometer (MERIS) on board the Environmental Satellite (ENVISAT) and the widely used Moderate Resolution Imaging Spectroradiometer (MODIS) on board Aqua and Terra. Section 2 summarises the data used in this study, Sect. 3 describes the AVHRR τ a retrieval and the differences between the original and the proposed method in more detail. A comprehensive validation of individual overpasses and monthly aggregated data is shown in Sect. 4. Conclusions are drawn in Sect. 5.

Data
The aerosol optical depth products from three different satellite sensors (AVHRR, MERIS, and MODIS) flown on board five satellites (NOAA-17, NOAA-18, ENVISAT, Aqua, and Terra), which cover mid-morning and afternoon orbits, are validated for the geographical area of Central Europe (40.5 • N to 50 • N, 0 • E to 17 • E) between August 2005 and December 2007. This time frame was selected because comprehensive data from all satellites were available for this period.

AERONET
A total of 15 Sun photometer measurement sites from AERONET (Holben et al., 1998) providing level 2.0 data (cloud screened and quality assured) are used for the validation in the investigated area (cf. Table 1). They represent a wide range of aerosol compositions (e.g., urban, rural) and topographic situations (e.g., flat or complex terrain). Additionally, aerosol micro-physical properties from level 1.5 inversion products  are used to determine the aerosol properties in the investigated area.

MERIS
MERIS is flown on ENVISAT in a morning orbit constellation, it crosses the equator at 10:00 a.m. local solar time (equator crossing time, EXT). We use the standard level 2 aerosol product (Santer et al., 1999) from the European Space Agency (ESA), which is based on the detection of DDV (Kaufman and Sendra, 1988) selected with the atmospheric resistant vegetation index, as the reflectance properties of dark targets are known well enough to retrieve the aerosol signal accurately. The main disadvantage of this technique is the rare occurrence of DDV pixels over Central Europe. In order to get better spatial coverage, the method has been extended to brighter surface types (Santer et al., 2007), which in turn has introduced more errors in the τ a and Angstrom wavelength exponent (α) retrievals (Vidot et al., 2008). We use τ a at 443 nm and α at 1 × 1 km 2 resolution to calculate τ a at 550 nm.

MODIS
For MODIS, Terra (MOD04L2, morning orbit with EXT∼10:30 a.m.) and Aqua (MYD04L2, afternoon orbit with EXT∼01:30 p.m.) collection 005 data from the second generation algorithm (Levy et al., 2007a) are used. For the initial MODIS aerosol retrieval, Kaufman et al. (1997a) presented a method to identify DDV pixels by making use of the correlation between ρ SFC at 2.12 µm (corrected for atmospheric effects) and visible wavelengths (0.47 and 0.66 µm). In Remer et al. (2005) the product was extended to land surfaces with ρ SFC (0.66 µm) ≤0.125. Compared to these prior MODIS products, the current surface reflectance estimation is based on an improved parametrization (ratios between visible and 2.12 µm channels, VISvs2.12), now considering variations by surface type (based on the short-wave infrared normalised difference vegetation index) and angular variability in the VISvs2.12 reflectance relationship. Additionally, the assumed aerosol properties have been updated according to Levy et al. (2007b) and new look-up tables have been computed considering polarisation. Finally, three channels are simultaneously inverted and the 2.12 µm channel is no longer assumed to be totally transparent for aerosols. According to Remer et al. (2005) and Levy et al. (2007a) the accuracy of the MODIS aerosol product over land is τ a =±0.05±0.15τ a .

AVHRR
We use daily full resolution (1.1 × 1.1 km 2 ) NOAA-17 mid-morning (EXT∼10:25 a.m.) and NOAA-18 afternoon (EXT∼01:45 p.m.) overpasses. Calibration, geocoding, orthorectification, cloud and cloud shadow detection are part of the pre-processing of the data and are described in more detail in Hauser et al. (2005a). Meteorological fields of total column ozone, vertically integrated water vapour, and sea level pressure from the European Centre for Medium-Range Weather Forecasts (ECMWF) between 09:00 and 15:00 UTC are interpolated to the time of the satellite overpass and are utilised for the atmospheric correction of water vapour, ozone influence, and Rayleigh optical depth. In order to estimate τ a from AVHRR over land surfaces, a modified version of the single-channel (∼630 nm), multitemporal approach of Hauser et al. (2005a) is applied. The general principle and modifications are described in the following section.

General principle
The physical processes of reflection, scattering and absorption of solar radiation in the earth-atmosphere system can be explained with the reflection function of a cloud-free and vertically homogeneous atmosphere (e.g., King et al., 1999;Liou, 2002) where θ 0 is the solar zenith angle, θ v is the view or satellite zenith angle, φ is the relative azimuth angle between sun and satellite, ρ atm is the atmospheric path reflectance including all atmospheric scatterers, T d (θ 0 ) and T u (θ v ) are the transmittance of the downward and upward path, respectively, s is the spherical albedo, and ρ sfc is the reflectance of the Earth's surface. The dimensionless reflectance can be defined as ρ=π I /cosθ 0 E 0 with I the measured intensity (or radiance) and E 0 the solar irradiance at the top-of-atmosphere (TOA). Usually, the value of s is small, and if a dark surface (with low ρ sfc ) is considered, the product sρ sfc can be neglected. Then it follows that where the atmospheric part is divided into the molecular and aerosol contribution. The physics behind the aerosol contribution over a dark target can be explained using the single scattering approximation (e.g., King et al., 1999) ρ a (θ 0 ,θ v ,φ) = ω 0 τ a P ( ) 4 cos θ 0 cos θ v with ω 0 the single scattering albedo and P ( ) the aerosol scattering phase function of the scattering angle . From Eqs. (2) and (3) one can see that the retrieval of τ a requires the decomposition of the satellite signal into the part reflected from the earth's surface and the part contributed from the atmosphere. The determination of ρ SFC , whose accuracy is a critical parameter in the derivation of τ a , is done by using a multitemporal technique which employs a time series of 45 consecutive days. Within this period a particular pixel is observed under various θ v . This is based on the assumption that the radiometric properties of the underlying land surface are almost stable and observations with low aerosol concentrations (background conditions) exist. If the TOA reflectance ρ TOA is considered as a function of θ v , ρ SFC can be estimated as the resulting minimum reflectance. A continuous bidirectional reflectance function is obtained by connecting the minimum values using the concept of convex hull (cf. Fig. 1).
The single scattering approximation in Eq.
(3) points out that additional information about the optical aerosol properties ω 0 and P ( ) is necessary in order to derive τ a . Since a single-channel approach allows only to retrieve one quantity of aerosol information independently, in our case τ a , the other aerosol properties for the investigated region need to be chosen in advance. This is usually done by using pre-defined aerosol models describing the physical aerosol properties size distribution and refractive index in the investigated region (e.g., d' Almeida et al., 1991) from which ω 0 and P ( ) can be calculated. Finally, ρ TOA and ρ SFC together with the pre-defined aerosol characteristics and the atmospheric conditions in terms of integrated water vapour content and ozone column from numerical weather prediction data are implemented into the Simplified Method for Atmospheric Correction (SMAC) semi-empirical radiative transfer scheme (Rahman and Dedieu, 1994), updated to 6S (Vermote et al., 1997), to calculate τ a .
To enhance the quality of the retrieved τ a , a spatial filtering procedure within a moving window of 15 × 15 pixels is applied to each pixel, which ensures that small-scale structures of the τ a variability are retained. In order to take into account undetected clouds or cloud shadow, only values within the 10 to 40 percentile of the filter window are used to calculate the mean τ a , which is comparable to the method of the MODIS algorithm in pre-selecting particular pixels for processing (Kaufman et al., 1997b;Levy et al., 2007a).

Differences to the original approach
The original NOAA AVHRR aerosol retrieval from Hauser et al. (2005a) uses a post-processing procedure based on a relationship between Sun photometer (AERONET) measurements and the satellite data (NOAA-16 AVHRR) to enhance the quality of the aerosol product. Employed on NOAA-17 and NOAA-18 AVHRR, the post-correction with the original approach did not lead to the same increase of the retrieval accuracy as for NOAA-16. For this reason and having in mind the derivation of a long-term τ a climatology dating back to 1985, we were seeking for a stand-alone algorithm which may also be more suitable for the pre-AERONET era. In the modified approach we adapted several steps of the original processing. In the next paragraphs we will shortly examine the changes introduced with the new procedure.

Surface reflectance
When using the original approach of Hauser et al. (2005a) to estimate ρ SFC , the correlation between the satellite-retrieved and Sun photometer-measured τ a was found to depend on the mean surface reflectance ρ SFC , with R=−0.71, meaning that over bright surfaces the resulting error is larger than over dark surface types. This also held true for the root mean square error RMSE. With the help of ESA's GlobCover product (Bicheron et al., 2008), areas mainly containing cropland (∼39% of the land surface pixels in the investigated area) were identified to result in the largest errors with the original convex hull scheme used to detect the minimum reflectance. Therefore, we adapted the scheme by dividing θ v into four bins (<−30 • /−30 • −0 • /0 • −30 • /30 • <) and searching the convex hull in each bin separately. The differences between the old and new scheme, indicated in Fig. 1, are especially obvious in the backward scattering region. Kaufman and Tanré (1996) showed that an error in ρ SFC of 0.01 translates into τ a ∼0.1. With this adaptation, the correlation for the resulting τ a was slightly improved and the RMSE could substantially be reduced in regions mainly dominated by cropland, e.g., at Avignon from 0.15 (0.12) to 0.09 (0.08) for NOAA-17 (NOAA-18). After applying this change to the ρ SFC estimation, the dependency between R and ρ SFC almost vanished (R=−0.10). But still, the darkest surfaces (Ispra, Venise, Villefranche; cf. Table 3) exhibit the best correlations and lowest errors.

Aerosol model
As shown in Sect. 3.1, the single scattering albedo ω 0 and the scattering phase function P ( ) are two important parameters in the derivation of τ a . They depend on the chosen physical aerosol properties (size distribution and refractive index). The analysis of these micro-physical properties from the AERONET inversion products turned out that the continental aerosol model used so far was too absorbing for the investigated area. In Table 1 we display the average and standard deviation of ω 0 at each AERONET site for the evaluated period. The continental model (original approach) in 6S results in ω 0 =0.89. Ten AERONET sites in Table 1 exhibit a less absorbing aerosol type than the continental model of 6S, whereas at three sites a more absorbing type is predominant. The exceptionally low value in Davos is expected to be an error caused by usually low τ a at this high altitude leading to substantial errors in the retrieval of ω 0 . Calculating the average for all AERONET sites, with the exclusion of Davos, reveals ω 0 =0.91. Zhang et al. (2001) have shown that a ω 0 =0.03 at τ a =0.5 corresponds to an error of 10% in the retrieval with increasing error for larger τ a . The standard deviation of ω 0 at some sites shows substantial variations causing temporarily large deviations from the assumed properties of the continental model which will lead to additional errors.
Thus, we derived a new set of aerosol models from AERONET inversion data splitting the aerosol properties, the average refractive index and size distribution, similar to Levy et al. (2007b) into the four seasons. The highest accuracy of these data can be expected from measurements with τ a (440 nm) ≥ 0.5 . Applying this threshold to the AERONET sites in the study region for most of them not enough points were available to calculate the seasonal means. Therefore, we changed the threshold for τ a (440 nm) from 0.5 to 0.2 (cf.  and subsequently calculated the averages from the level 1.5 data. According to , the quality of the derived size distribution does not depend on this threshold, whereas for the refractive index errors of up to 50% may be possible. The resulting ω 0 for the entire area and the various seasons now is 0.85 for winter, 0.90 for spring and fall, and 0.91 for summer. A comparison between the different scattering phase functions is shown in Fig. 2a. Although the differences between the curves are only minor, calculating the reciprocal of the product ω 0 P ( ) from the single scattering approximation in Eq. (3) for the most frequent range of between 90 • and 160 • results in 2-30% higher values for the seasonal aerosol properties compared to the continental model of 6S with even larger deviations for >160 • , as is demonstrated in Fig. 2b. The differences in the resulting τ a will be discussed in more detail in Sect. 4.1.

Additional changes
Compared to Hauser et al. (2005a) we reduced the size of the filter window from 25 × 25 to 15 × 15 pixels, with the aim not to smooth too much of the small-scale aerosol structures. Despite the reduction of 64% in the filter area, this procedure would still mix values from low and high altitudes in mountainous regions. Matthias et al. (2004) demonstrated with measurements from the European Aerosol Research Lidar Network (EARLINET) that a 500 m height difference corresponds to τ a ∼0.05-0.15 in the lower planetary boundary layer for several locations in Central Europe. Therefore, we apply a 1 km digital elevation model (HYDRO1k, US Geological Survey) and restrict the pixels of the filter window to be within ± 250 m of the central pixel. Finally, owing to quality reasons we apply two criteria: (1) We restrict the retrieval to pixels with ρ SFC ≤0.07 and (2) the regional standard deviation of τ a within the 15 × 15 region must not exceed 0.1. These thresholds were adjusted by comparing the retrieved τ a with AERONET measurements and requiring that approximately 90% of the data should be retained after the filtering procedure. Even though the brightest surface types are excluded, we do not get considerably less observations than in the original approach. Only areas mainly dominated by bright surfaces are excluded in the processing.

New versus old approach
In order to demonstrate the performance of the new method compared to the original one, we applied the modified algorithm to the same data set as in Hauser et al. (2005a). Therefore, NOAA-16 AVHRR data were evaluated for the same set of AERONET sites (nine locations) and the same time period (from May 2001 to December 2002) as in the above mentioned publication. Table 2 shows the results from two exemplary AERONET sites, Avignon and Ispra, and the summary of all sites (cf. Hauser et al., 2005a). The surrounding of Avignon is typically dominated by cropland which resulted in larger uncertainties in the ρ SFC estimation with the original approach, as has been explained before. Ispra, on the other hand, is situated in an area covered by mostly dark surface types like forests and lake water. Focusing on the linear correlation coefficient R one can see that the proposed method clearly outperforms the original approach for Avignon and for the summary of nine sites. At Ispra the values of R resemble each other. The same holds true for the standard deviation σ and the linear regression equation y=A+Bx.
Although areas widely covered by bright land cover types (ρ SFC >0.07) are excluded in the τ a retrieval, the better performance over large areas together with the quasiindependent design (stand-alone) are obvious advantages of the improved method.  (Hauser et al., 2005a) and the new method to derive τ a from AVHRR over land surfaces. The differences between the various methods are explained in Sect. 3. R denotes the linear correlation coefficient, σ is the standard deviation, A and B are the coefficients of the linear regression equation y=A+Bx.

Site
R

Results and discussion
For the comparison of satellite-retrieved and Sun photometer-measured τ a , former validation studies (e.g., Chu et al., 2002) have recommended to take both spatial and temporal variability of the aerosol distribution into account. The first requirement is achieved by building average τ a over an area of 50 × 50 km 2 . The box size corresponds to a maximum number of 25 (2500) MODIS (AVHRR, MERIS) τ a values per box, where at least 20% of the pixels had to be valid to calculate the average. The temporal variability was treated as in the above-mentioned study using at least two out of five possible AERONET measurements within ± 30 min of the satellite overpass to calculate the time average. At the AERONET sites τ a at 550 nm is computed using a second-order polynomial fit in a lognormal space. Figures 3 and 4 give a first impression on the validation results between AERONET and both NOAA spacecraft. They display scatter-density plots and the linear regression equation for the summary of all AERONET sites listed in Table 1 and five additional locations as examples for different land cover types (GlobCover) and aerosol conditions (AERONET). Avignon is mainly dominated by cropland with weakly absorbing and small-sized particles (indicated by the size dependentÅngstrom exponent α in Table 1); Fontainebleau is a mixture of cropland and forest in a region  Table 3.   Hess et al. (1998); Ispra is surrounded by forest, lake water, and is influenced by the highly industrialised and polluted Po Valley; Laegeren is a hilly region with mixed land cover and mostly non-absorbing aerosols with smaller sizes than the continental model of Hess et al. (1998); Villefranche is situated on a peninsula at the Mediterranean adjacent to a coastal region with highly variable topography and varying aerosol conditions.

Validation of AVHRR
Detailed statistical results are presented in Table 3. A total of 2369 match-ups N for NOAA-17 and 2334 for NOAA-18 have been found during the investigated period. The correlation coefficient R, which expresses the overall ability to retrieve the aerosol signal, equals R=0.69 for NOAA-17 and R=0.71 for NOAA-18. The root mean square error RMSE (∼0.1) and the standard deviation σ (∼0.07) represent the magnitude of random errors which are proportional to radiometric noise, variations in ρ SFC , and subpixel cloud contamination.
For most of the sites, the differences between NOAA-17 and NOAA-18 are small and could be due to undetected subpixel clouds or uncertainties in calibration (Li et al., 2009). The time-of-day should not largely affect the retrieval, which is also demonstrated by the similar performance of both MODIS spacecraft, despite their different equator crossing times. Small differences could be introduced by evolving cloudiness, like small-sized cumulus clouds, during the day, especially in mountainous areas. Bidirectional reflectance distribution function (BRDF) effects may slightly influence the results, but as has been shown in Fig. 1 in the estimation of ρ SFC , we account for variations of the viewing geometry. Some uncertainty will remain due to changes in sun elevation during the 45-day period used for the retrieval of ρ SFC and due to intra-daily variations caused by varying overflight times. Thus, two different BRDF models (Rahman et al., 1993;Wu et al., 1995) were tested in order to pre-correct the data for BRDF effects, but the resulting ρ SFC turned out to be less accurate (similar or higher RMSE, reduced correlation) than without the correction. In addition, we also tried to shorten the 45-day period for the estimation of ρ SFC , which showed some improvements at locations with low cloud coverage. On the other hand, in regions or during periods with high cloud coverage the chance to get cloud-free observations with low aerosol concentrations is reduced. Interestingly enough, at the offshore platform of Venise the performance of NOAA-17 is clearly worse than that of NOAA-18. At the same time Villefranche, a location on a peninsula containing both water and land, shows good correlations for both satellites. A possible explanation could be the influence of different illumination conditions when viewing the shallow coastal waters and water sediments around Venise which is also supported by the different values of ρ SFC in Table 3. A comparison of monthly averages does not reveal a systematic difference between the two satellites and the mean relative difference of ∼2% is similar to the one between Terra and Aqua (∼2% for the investigated period).
Considering the linear regression equation y=A+Bx between AVHRR and AERONET, both satellites have a similar offset A and slope B. A nonzero offset indicates a biased retrieval at low τ a which may be associated with sensor calibration errors or erroneous ρ SFC estimation (Zhang et al., 2001). A slight overestimation at low τ a can be found for most of the sites. Deviations of B from unity are associated with incorrect assumptions of the aerosol micro-physical properties (Zhang et al., 2001). Before we introduced the new aerosol model, the continental aerosol model of the former aerosol processing led to a slope B of 0.79 (0.82) for . With the introduction of the new model this could be clearly improved, now resulting in B=0.94 for both satellites. This means that especially high aerosol concentrations are detected more accurately with the new set of aerosol models. The offset did not change and the other statistical parameters remained almost similar compared to the continental model.
Aside from the aerosol model, a second effect may cause slopes of less than unity. The accuracy of the calculations with the utilised semi-empirical radiative transfer code SMAC decreases for τ a >0.8 (Rahman and Dedieu, 1994); therefore, we limited the retrieval to τ a ≤0.8 and higher aerosol loads than this upper limit are not present in the comparison between AERONET and AVHRR. However, less than 1% of the AERONET observations were larger than 0.8 during the investigated period and the impact on the statistics is assumed to be almost negligible. Nonetheless, this limitation may especially affect the monthly statistics during periods with high aerosol loads in certain regions.

Comparison between AVHRR, MERIS, and MODIS
This section compares the performance of the AVHRR τ a retrieval with two other medium resolution sensors, MERIS and MODIS. Instead of relating the results from AVHRR to other validation studies, we want to perform an accurate comparison using the same set of AERONET sites for all satellites during the same time period with the approach explained at the beginning of this section. The results are summarised in Table 3 and Figs. 3-5.
Looking at the scatter-density plots with the summary of all AERONET sites  it is obvious that MODIS shows the best performance, with slightly higher R found for Aqua than for Terra. This is not surprising keeping in mind that, in contrast to AVHRR, MODIS has a much better spectral coverage making use of seven spectral channels for the characterisation and retrieval of aerosol properties. However, the correlations of AVHRR and MERIS resemble each other, although Fig. 5 demonstrates that τ a from the latter is biased high. Considering σ and RMSE, AVHRR mostly reveals lower values than MERIS. The results from MERIS are consistent with Vidot et al. (2008), who found similar errors and explain a part of the high RMSE with the presence of thin clouds, which the current MERIS algorithm is not able to flag. Due to the lack of infra-red channels, the MERIS level 2 cloud screening method in general experiences difficulties to detect clouds accurately (Gómez-Chova et al., 2007). Apparently, MODIS exhibits the lowest error rates (RMSE=0.04), comparable to other MODIS validation studies (Chu et al., 2002;Levy et al., 2007a;Remer et al., 2008). Analogous to AVHRR, Aqua and Terra perform differently at Venise, with a higher RMSE of the morning crossing platform Terra.
Analysing A and B of the linear regression equation (y=A+Bx), at the majority of locations and in the summary of all sites every satellite shows A>0 and B <1 (cf. Table 3) meaning that at low aerosol concentrations τ a is overestimated and at high concentrations underestimated. This effect is most pronounced for MERIS. Compared to the work of Vidot et al. (2008), who found a very low offset (A∼−0.02-0.03) in Europe by using three AERONET sites (Ispra, Lille, Minsk), we obtain A=0.11 for MERIS using 15 European AERONET sites. A part of the discrepancy can be explained with biased α values (Kokhanovsky et al., 2007;Vidot et al., 2008) used for the calculation of τ a at 550 nm and similar to the RMSE, with the occurrence of some extreme outliers. Vidot et al. (2008) apply a smaller region of 12 × 12 km 2 around each AERONET site. We also investigated the influence of the box size around each AERONET site using averaging areas between 10 × 10 km 2 and 50 × 50 km 2 , but the effect on the statistical results turned out to be minor and cannot explain the discrepancies between our results and the one found by Vidot et al. (2008). Aqua and Terra show the lowest offset (0.01) which may also be caused, at least partly, by allowing small negative values of τ a (Levy et al., 2007a) in order to improve the statistics of aggregated (spatial and temporal) products.
As described in the previous section, underestimation of τ a at high aerosol burden can be caused by underestimated light absorption in the aerosol models. In the past, the AVHRR retrieval was based on the same aerosol properties for whole Europe for the entire year. As demonstrated before, a new set of aerosol models helped to improve the slope of the AVHRR aerosol product substantially, now predicting high aerosol loads more accurately and resulting in a slope similar to the MODIS aerosol product. Regarding the MODIS retrieval, a nonabsorbing model (ω 0 ∼0.95) with dynamic size parameters depending on the optical depth is used for the investigated area (Levy et al., 2007b), which results in a slope close to unity. Considerable deviations are only found for the mountain site Davos, which may, to a certain extent, be explained with the small sample number. The MERIS retrieval applies 12 different aerosol models without absorption resulting in B=0.67. Beside the aerosol model, such a large deviation from unity could also be explained with a relatively high overestimation of low aerosol concentrations and as mentioned before, with errors introduced with biased α values.
Finally, looking at the expected error of the MODIS aerosol product over land τ a =±0.05±0.15τ a (e.g., Levy et al., 2007a) it turned out that 86% of Aqua and Terra retrievals are found within these boundaries in the investigated area. For both NOAA-17 and NOAA-18 63% of τ a and for MERIS 56% of τ a can be found within these limits.

Spatial distribution
First, in order to show spatial differences between the aerosol products, we qualitatively compare MODIS as a reference with AVHRR (MERIS). Figure 6 is an example of the monthly mean τ a during April 2007. This month was chosen since stable weather conditions due to high surface pressure yielded to a high observation frequency and also led to a high aerosol burden in the Po River Valley. To calculate the average at a particular pixel at least 10% of the days had to provide valid values.
The overall pattern between MODIS and AVHRR looks comparable with respect to the geographical pattern, but the magnitude of τ a especially in the Po River Valley is clearly underestimated by the latter. This is in accordance with the results of the regression analysis (cf. Table 4) in which the slope of NOAA's AVHRR is clearly below unity meaning that high aerosol loads are underestimated. A more detailed analysis of the data from Ispra revealed that mainly the last third of this particular month was characterised by very high aerosol burden (τ a ). During the first two thirds of April 2007 the underestimation of τ a with AVHRR is only about 13% compared to the AERONET average, while for the last part the proposed method underestimates the measured values by approximately 35%. Considering the entire month both AVHRRs underestimate the AERONET value (τ a =0.44) by 25%. Even the monthly averages from Aqua and Terra MODIS are lower by 16% and 23%, respectively, than the AERONET mean. In addition to this situation, the central parts of the Po Valley are dominated by bright surface types which more often exceed the ρ SFC limit of 0.07. Therefore, poor sampling in these areas leads to large deviations between the satellite-derived and Sun photometer-measured monthly τ a during this particular month. Looking at the entire area processed for this study, approximately 11% of the pixels result in ρ SFC brighter than 0.07, depending on the season. In the vicinity of the Rhone River delta (∼43.5 • N, 4.5 • E) AVHRR also displays lower τ a than MODIS for the same reasons as for the Po Valley. Levy et al. (2009) discussed the impact of different averaging and weighting techniques when aggregating satellite data for climate studies. They found out that the choice of method can have a large effect on the resulting averages. Thus, for climate studies this needs further consideration and investigation.
Focusing on the Alpine region (stretch from 44 • N/7 • E to 47.5 • N/15 • E), AVHRR shows more details compared to MODIS, with clearly emerging valleys which could mean that the former aerosol product seems more appealing for climate or air pollution studies in such regions. Not only does MERIS capture the high concentrations in the Po River Valley well, but is also capable of detecting τ a in the Alpine valleys. Nonetheless, MERIS exhibits a substantial overestimation of τ a in most regions compared to MODIS. Furthermore, the MERIS image contains more areas with data gaps (gray) which may be due to the smaller swath width of MERIS leading to fewer observations per geographic location than for the other instruments.

AERONET versus satellites
Finally, this comparison demonstrates the capability of the polar orbiting platforms to represent the monthly average conditions at ground level which is of interest for climate studies. Therefore, at each location shown in Table 1 the monthly average τ a was calculated for AERONET using all available Sun photometer measurements and for the satellites using all retrieved τ a during a particular month. The spatial averaging of the satellite data is done as described at the beginning of this section. For both AERONET and satellite at  Table 4.  Table 4 and Fig. 7. In contrast to the linear regression of the individual collocations, all instruments show increasing offsets and, except for MODIS, decreasing slopes. Myhre et al. (2005) found a similar behavior for aerosol retrievals over oceans and argued that sampling issues are the main reason. The study of Anderson et al. (2003) supports this findings by investigating mesoscale variations of tropospheric aerosols. They found that the concentration of the airborne particles is barely homogeneous over time scales of more than 12 h. This is supported by the fact that the largest offset is found for MERIS, having a lower repeat cycle at a certain geographical region compared to the other two sensor types. At high aerosol loads, τ a is most clearly underestimated by AVHRR, probably caused by the cutoff at τ a >0.8. As explained before, the cutoff is due to limitations of SMAC (Rahman and Dedieu, 1994). Nonetheless, with the implementation of the new aerosol model the results could be improved. The continental aerosol model led to a slope B of 0.69 (0.62) for NOAA-17 (NOAA-18), whereas now the adapted model leads to B=0.79 (B=0.69). A more detailed inspection of the data revealed that some outliers can either be explained with low frequency of AERONET or satellite observations for certain months. According to Levy et al. (2009), the use of different averaging techniques changes the resulting monthly mean τ a . After filtering the data with a criterion of at least 20% of days with valid measurements per month, the statistics were slightly improved.
Apparently, a reduction of R can be found for most sensors with regard to the individual overpasses. Again, the update of the aerosol properties improved the relationship between the monthly means from AERONET and both AVHRRs, for NOAA-17 (NOAA-18) from R=0.60 (R=0.56) applying the continental model to R=0.67 (R=0.63) by using the new aerosol model. The RMSE, however, is almost similar to the individual match-ups, solely the MERIS RMSE significantly increased. As explained above, high temporal variability in τ a and sampling issues are the major cause for this behaviour. Table 4 supports this finding since the strongest decrease of R and increase of RMSE is found for MERIS. With a repeat cycle of approximately three days, the capability of reproducing the "true" τ a with MERIS is limited.
In Fig. 8 we compare the aggregated monthly averages of all sites together for each of the 29 months between August 2005 and December 2007. Most of the time the reference curve of AERONET is at the lower end of the monthly values meaning that the satellite sensors usually overestimate them. Nonetheless, the yearly cycle is altogether represented well by all of them. As before, most of the differences may be explained with sampling issues which explains that they are largest for MERIS. Considering MODIS and AVHRR, R  (RMSE) of the aggregated monthly averages is significantly higher (lower) than for the averages of the single sites shown in Table 4 with R of 0.88, 0.82, 0.72 and 0.72 for Aqua, Terra, NOAA-17 and NOAA-18, respectively; the RMSE can be found in the range of 0.041 (Terra) to 0.053 (NOAA-18).
Only for MERIS the statistics did not improve a lot with R=0.51 and RMSE=0.21.

Intersatellite comparison
Since we are strongly interested in the compilation of an AVHRR long-term climatology, the robustness of the proposed method is an important issue. Hence, we shortly examine the intersatellite differences and compare the consistency of the AVHRR retrieval with the one of MODIS. The scatter plots in Fig. 9 and Table 5 reveal that the cross-satellite performance of NOAA's AVHRR is similar to MODIS with almost the same R and somewhat lower RMSE. Moving from small-scale (individual averages, Fig. 9a) to large-scale (aggregated averages, Fig. 9b) areas, for both instruments the  Table 1 separately. For the aggregated values (b), all sites are averaged for a particular month. Statistical results can be found in Table 5. statistics significantly improve. From this we can conclude that despite local differences the overall performance is good. This is in accordance with the study of Kaufman et al. (2000) who found out that intra-daily variations of τ a are usually low and, therefore, do not largely influence the retrieval of aerosol properties. These are important findings with regard to the compilation of an AVHRR τ a climatology and they demonstrate that the proposed method works platformindependently.
We assessed the performance of a modified version of an AVHRR aerosol optical depth retrieval algorithm for land surfaces. In contrast to the original approach, it is a quasistand-alone procedure now no longer needing AERONET data to post-correct the retrieved τ a . Nonetheless, since a single-channel approach allows only to derive on piece of aerosol information independently, AERONET inversion products are needed to derive optical aerosol properties like single scattering albedo and scattering phase function in the investigated region. The scheme to derive ρ SFC has been adjusted which turned out to substantially reduce the RMSE over areas mainly dominated by cropland (∼39% of land surface pixels). Moreover, a detailed analysis of AERONET inversion products turned out that the aerosol model used for the AVHRR retrieval so far was absorbing too much in the investigated area. Based on aerosol size distribution and refractive index data from AERONET a new set of aerosol models was derived and the retrieval at high aerosol concentrations could be improved, which in turn has a positive effect on the statistics when deriving monthly averaged τ a . With the aim to compile a long-term τ a climatology for Europe making use of our AVHRR archive dating back to 1985, the method should be robust and cross-platform applicable. Therefore, two AVHRR sensors with overlapping time period flown on board NOAA-17 and NOAA-18 were used in this study. The results were validated against τ a from AERONET Sun photometer measurements and were compared with two other medium resolution satellite sensors, MERIS and MODIS. Related to the newer generation sensors, AVHRR reveals good results and also the cross-platform differences are similar to those between Terra and Aqua MODIS. The low sensitivity of the τ a retrieval to the time-of-day is in accordance with results from other studies. This is an important step towards a climatology. Finally, the linear regression analysis of monthly average τ a between each AERONET site and the various satellites results in lower correlations with higher intercepts and lower slopes than for the individual match-ups. Simultaneously, when aggregating all AERONET sites to monthly means and comparing these with the satellites, the overall capability of them to reproduce the true atmospheric conditions is mostly good. This means that variations at small-scale are harder to capture than those of regional or even larger scale. An explanation for this may mainly be sampling issues which is also supported by the fact that AVHRR and MODIS, having a larger swath width than MERIS and, therefore, a higher repeat cycle for a specific geographic location, show the smallest increase in the error statistics. This is of importance when selecting satellite data for long-term studies.
Owing to sampling issues differences between timeaggregated ground-based and satellite measurements will always occur. Therefore, considering climate studies, not only should we focus on the absolute values of a comparison, but also investigate whether the trend of the aerosol concentration can be captured right.
In order to detect trends in the atmospheric aerosol load, accurate calibration is of primary importance. Initially, AVHRR was not foreseen to provide information for quantitative remote sensing, but rather to contribute to geostationary satellite observations at high latitudes. Hence, no onboard calibration unit for the visible and near-infrared channels is included in the sensor's design and the instrument was observed to degrade with time. Therefore, for the use of historical AVHRR data sets relying on data from several different satellites calibration is an important factor.
Hence, based on the robust method introduced and validated in this study and making use of the large AVHRR data archive at our institute, the proposed method is currently tested on prior AVHRR sensors and the effects of different calibration schemes are investigated in a trend analysis. This shall finally lead to the compilation of an aerosol long-term climatology covering Central Europe.