Scientific impact of MODIS C 5 calibration degradation and C 6 + improvements

The Collection 6 (C6) MODIS (Moderate Resolution Imaging Spectroradiometer) land and atmosphere data sets are scheduled for release in 2014. C6 contains significant revisions of the calibration approach to account for sensor aging. This analysis documents the presence of systematic temporal trends in the visible and near-infrared (500 m) bands of the Collection 5 (C5) MODIS Terra and, to lesser extent, in MODIS Aqua geophysical data sets. Sensor degradation is largest in the blue band (B3) of the MODIS sensor on Terra and decreases with wavelength. Calibration degradation causes negative global trends in multiple MODIS C5 products including the dark target algorithm’s aerosol optical depth over land and Ångström exponent over the ocean, global liquid water and ice cloud optical thickness, as well as surface reflectance and vegetation indices, including the normalized difference vegetation index (NDVI) and enhanced vegetation index (EVI). As the C5 production will be maintained for another year in parallel with C6, one objective of this paper is to raise awareness of the calibration-related trends for the broad MODIS user community. The new C6 calibration approach removes major calibrations trends in the Level 1B (L1B) data. This paper also introduces an enhanced C6+ calibration of the MODIS data set which includes an additional polarization correction (PC) to compensate for the increased polarization sensitivity of MODIS Terra since about 2007, as well as detrending and Terra– Aqua cross-calibration over quasi-stable desert calibration sites. The PC algorithm, developed by the MODIS ocean biology processing group (OBPG), removes residual scan angle, mirror side and seasonal biases from aerosol and surface reflectance (SR) records along with spectral distortions of SR. Using the multiangle implementation of atmospheric correction (MAIAC) algorithm over deserts, we have also developed a detrending and cross-calibration method which removes residual decadal trends on the order of several tenths of 1 % of the top-of-atmosphere (TOA) reflectance in the visible and near-infrared MODIS bands B1–B4, and provides a good consistency between the two MODIS sensors. MAIAC analysis over the southern USA shows that the C6+ approach removed an additional negative decadal trend of Terra1NDVI∼ 0.01 as compared to Aqua data. This change is particularly important for analysis of vegetation dynamics and trends in the tropics, e.g., Amazon rainforest, where the morning orbit of Terra provides considerably more cloudfree observations compared to the afternoon Aqua measurements.

sites. The PC algorithm, developed by the MODIS ocean biology processing group (OBPG), removes residual scan angle, mirror side and seasonal biases from aerosol and surface reflectance (SR) records along with spectral distortions of SR. Using the multiangle implementation of atmospheric correction (MAIAC) algorithm over deserts, we have also developed a detrending and cross-calibration method which removes residual decadal trends on the order of several tenths of 1 % of the top-of-atmosphere (TOA) reflectance in the visible and near-infrared MODIS bands B1-B4, and provides a good consistency between the two MODIS sensors. MA-IAC analysis over the southern USA shows that the C6+ approach removed an additional negative decadal trend of Terra NDVI ∼ 0.01 as compared to Aqua data. This change is particularly important for analysis of vegetation dynamics and trends in the tropics, e.g., Amazon rainforest, where the morning orbit of Terra provides considerably more cloudfree observations compared to the afternoon Aqua measurements. Barnes, 2006). The SD degradation is caused by exposure to solar radiation, which reduces its reflectivity over time particularly at shorter wavelengths. In the normal Earth view (EV) mode, the solar diffuser (Spectralon plate) is shielded from solar radiation by the SD door, which opens only during the SD calibration cycle. A SD attenuation screen with 7.8 % transmittance is also part of the MODIS calibration system, which is used to calibrate the high gain bands. The detailed documentation describing MODIS Terra onboard calibration is contained on the website of the MODIS Characterization Support Team (MCST: http://mcst.gsfc.nasa.gov/calibration/ information).
Two MODIS Terra events had a significant, long-lasting effect on the sensor performance and its degradation over time. The first event occurred during prelaunch thermal vacuum testing, when a portion of the nadir aperture door was overheated. As a result, a strip of door paint (epoxy) evaporated and coated part of the optics and the scanning mirror. The affected parts were visually cleaned, but either some contamination remained or protective coating was damaged resulting in differences between the two sides of the scan mirror as well as affecting their response versus scan angle (RVS). The RVS was not subsequently recharacterized and its exact state at launch was unknown. This is important, because once in orbit, the shape of the MODIS prelaunch RVS is being used to calibrate the instrument.
The second event was related to the SD door operation anomaly that occurred in May 2003. It led to the decision to keep the SD door permanently open and fix the SD screen in the closed state. From that time on, the SD plate has been degrading at a faster pace due to increased exposure to solar radiation. The decreasing reflectivity of the SD plate provides a smaller signal for the SD stability monitor (SDSM) and potentially reduces its capability to accurately track sensor response change over time.
For these reasons, trends in MODIS Terra calibration were not obvious in land and atmosphere products from the pure calibration prospective up until about 2009-2010. In general, detection of a sensor calibration trend and its distinction from physical changes in the Earth system is a very challenging and often unproductive task. It requires accumulation of time series of sufficient length, comprehensive statistical analyses, and often a confluence of independent "trend" reports from various geophysical products in the Earth science disciplines; this has been the case with MODIS Terra.
The earliest reports on MODIS Terra calibration degradation including change of polarization sensitivity came from the MODIS ocean biology processing group (OBPG) (Franz et el., 2008). This discipline is very sensitive to the sensor calibration and its temporal changes because the signal used in biomonitoring of the world ocean ("ocean color") is small. The significant temporal drifts of the MODIS Terra ocean color record from the more stable MODIS Aqua and, particularly, the benchmark SeaWiFS (Sea-viewing Wide Fieldof-view Sensor) record were documented in Kwiatkowska et al. (2008). These studies also showed that the calibration degradation effect is largest in the shortwave spectrum and rapidly decreases with wavelength.
The ocean color analysis uses the set of MODIS ocean bands (B8-B16) that have a spatial resolution of 1 km and cover a spectral range of 0.41-0.87 µm. Except for channels B8 and B9, they have a low dynamic range and a low saturation threshold required for accurate analysis of the weak ocean signal. The standard MODIS dark target (DT) aerosol algorithm (Levy et al., 2007a) over land, cloud optical properties over land and ocean (Platnick et al., 2003;King et al., 2013), and land atmospheric correction (Vermote and Kotchenova, 2008) algorithms rely mostly on MODIS "land" channels B1-B7 (0.47-2.13 µm) with 250-500 m nadir resolution. With longer time series from MODIS Terra, the evidence of systematic trends, caused by calibration degradation, was found in the MODIS C5 NDVI (normalized difference vegetation index; Wang et al., 2012), aerosol (Levy et al., 2010) and cloud optical thickness data sets. Our current analysis shows the evidence of trends in the surface reflectance and vegetation index C5 records from MODIS Terra.
The upcoming MODIS Collection 6 is based on a new calibration approach that removes major nonpolarimetric calibration trends from the MODIS Terra and Aqua data. This paper pursues several main goals: (1) to provide a physical level explanation of the previous (C5) and current (C6) MCST calibration approaches for the broad user community (with references to the detailed technical descriptions); (2) give an overview of the long-term calibration trends in aerosol, cloud, and surface reflectance C5 products from MODIS Terra and Aqua; and (3) describe further C6+ enhancement of MODIS Terra calibration developed in this work, which include the following: -an additional polarization correction of MODIS Terra data based on the OBPG model, and -detrending and cross-calibration adjustment of Terra and Aqua based on analysis of quasi-stable desert calibration sites.

Effect of MODIS calibration trends on C5 atmospheric and land products
The calibration degradation of MODIS Terra has affected the aerosol and cloud optical products (e.g., MOD04_L2 and MOD06_L2) of the MODIS atmosphere discipline group. While detailed trend analysis is beyond the scope of this paper, Figs. 1 and 2 provide an illustration using the time series of aerosol and cloud products, reported in the joint MODIS atmosphere team C5 Level 3 monthly product for the period 2002-2013. In all plots, the mean value for each month is an area-weighted average of monthly mean "grid points" (1 • × 1 • ) between ±60 • latitude. Each monthly grid value is aggregated from the daily grid means (e.g., Levy et al., 2009). However the aerosol time series in Fig. 1 gives equal weighting to each day's contribution, whereas the cloud optical thickness time series in Fig. 2 weights each day by the number of retrieved cloudy pixels.
Changes in MODIS Terra calibration were first noticed by the "dark-target" (DT) aerosol algorithm (e.g., Levy et al., 2007a, b;Remer et al., 2008) group from an apparent divergence between Terra and Aqua AOD (aerosol optical depth; 0.55 µm) over land (Levy et al., 2010). This divergence, which was not supported by the ground-based AERONET (Aerosol Robotic Network; Holben et al., 1998) sun photometry, becomes more obvious with the longer time series (Fig. 1, top). While the Aqua record showed no change over time, the Terra AOD record exhibited a significant downward linear trend of about −27 % decade −1 . A detailed Levy et al.'s (2010) comparison of the long-term Terra-C5 AOD record with AERONET (Holben et al., 1998) data over land revealed a gradual change of the retrieval bias from positive (pre-2004) to negative, while the Aqua-C5 AOD record showed a good consistency with AERONET over time. These results imply changes in MODIS Terra sensor characteristics, mainly in band B3 (0.47 µm) which is the primary band for the DT aerosol retrievals over land.
A similar analysis over the ocean showed no significant AOD trends in both sensors and a small positive Terra AOD bias with respect to Aqua (Fig. 1, middle). The DT ocean algorithm uses a different set of bands: while the total AOD is derived from the B2 (0.86 µm) reflectance, it relies on other bands, mainly (B4, B1) to determine spectral dependence of AOD, or Ångström exponent (AE). The AE time series over ocean ( Fig. 1, bottom) reveals a negative trend of MODIS Terra with respect to Aqua indicating calibration changes in bands B4, B1 and B2 over time.
These analyses from the DT algorithm group, supported by the other teams' findings, stimulated research on trend characterization and correction which eventually evolved into the C6 calibration approach described in the next section. Levy et al. (2013) tested the C6 calibration updates with the DT aerosol algorithm and concluded that Terra and Aqua will converge on their global trends.
The MOD06 cloud optical and microphysical property product includes cloud optical thickness, effective radius, and derived water path data sets (Platnick et al., 2003;Hubanks et al., 2008). The information content for cloud optical thickness primarily comes from a spectral band that is practically nonabsorbing for cloud particles, specifically, MODIS bands B1 (0.66 µm) and B2 (0.86 µm), over nonsnow/ice land and ocean surfaces, respectively. Similar to the discussion of band B3 above, both of these spectral bands saw significant degradation in C5 L1B data. Figure 2 shows mean cloud optical thickness (COT) for liquid-phase clouds, and for MODIS Terra and Aqua for the common pe- riod through July 2013 for land and ocean retrievals (averaged over ±60 • latitude). A trend fit of about −5 and −4 % decade −1 is found for MODIS Aqua for land and ocean retrievals, respectively; a fit of −17 and −13 % decade −1 is found for MODIS Terra for land and ocean, respectively. A similar result is found for the ice-phase COT time series. The significantly larger trend for MODIS Terra is a result of the reflectance calibration trending artifact in MODIS bands B1 and B2.
The analyses from the DT algorithm group stimulated other studies, e.g., by Wang et al. (2012), to understand how calibration drifts in B3 could lead to an apparent drift in NDVI (B1, B2). To illustrate the MODIS Terra C5 calibration artifact on NDVI, we selected a 500 km tile in Georgia, USA, and conducted processing of the common Terra/Aqua time period using the multiangle implementation of atmospheric correction (MAIAC) algorithm (Lyapustin et al., 2011a(Lyapustin et al., , b, 2012. MAIAC uses the time series analysis and processing of groups of pixels for cloud detection and simultaneous retrievals of aerosol properties and surface bidirectional reflectance. Its provisional status products will become publicly available shortly after completion of MODIS Land C6 reprocessing in 2014. Figure 3a gives an example of MODIS Terra top-of-atmosphere (TOA) RGB (red, green, blue) reflectance over selected area and results of MAIAC processing. The monthly mean area-average NDVI are presented in Fig. 3b showing a decreasing trend of about −0.01 decade −1 in Terra C5 (red) as compared to Aqua C6 (blue) data.
The described multidisciplinary analysis clearly shows sensor degradation artifacts present in MODIS C5 science products. While the observed Aqua trends may contain both natural and calibration-related components (for instance, a small negative C6 Aqua NDVI trend may be related to an urban development of the region, which includes the Atlanta metropolitan area, as well as to the natural variations during 11.6 years of analysis from October 2002 to April 2013), the significant divergence of C5 Terra geophysical records is a clear indication of the calibration trend. The MODIS Collection 6 will significantly reduce these artifacts in science products. At the same time, the MODIS C5 record will be maintained in parallel for one more year. One of the goals of this paper is to raise awareness of the broad MODIS user community to the calibration-related artifacts in MODIS C5 products which for Terra can be assessed as global "decreasing" decadal trends of ∼ 27 % AOD,∼ 17 % COT, and ∼ 0.01 NDVI over land.

Overview of radiometric calibration of MODIS reflective bands
The MODIS scanning and calibration geometry is schematically shown in Fig. 4 with angles of incidence (AOI) defined relative to the normal to the scan mirror. During one rotation of the double-sided scan mirror, it observes Earth in the AOI range of 10.5-65.5 • which corresponds to 1354 pixels along the scan line. It also views the solar diffuser at AOI = 50.25 • , equivalent to the EV pixel 978, and the moon via the space view (SV) port at 11.25 • AOI, EV pixel 17 (Sun et al., 2003). The moon observations usually require a spacecraft roll maneuver that is conducted about nine times a year. Thus, the standard MODIS calibration protocol is limited to only two angles (at SD and moon AOIs) out of the full RVS function. Prior to C6, the standard approach included trending of these two points with linear interpolation/extrapolation if necessary to estimate the RVS at other angles. By the C6 time frame, enough evidence had accumulated to indicate that the MODIS Terra RVS change is nonlinear and that the previous (C5) approach was not sufficient (Sun et al., 2012).
To track RVS change, the MCST C6 calibration algorithm introduced the EV monitoring of stable desert calibration sites recommended by the Committee on Earth Observation Satellites (CEOS) (http://calval.cr.usgs.gov/sites_ catalog_map.php). In this case, all AOIs can be characterized independently via the surface BRDF (bidirectional reflectance distribution function). The current C6 approach (Sun et al., 2012;Toller et al., 2013) includes (1)   trending at SV (moon view) AOI; (3) obtaining the normalization coefficient (moon/EV) 11.25 • at the moon AOI; and (4) obtaining full RVS (all AOIs) by multiplying EV RVS by the normalization coefficient. The new routine achieves a good agreement between the EV trending at different desert sites and the moon trending at the moon AOI. Obviously, accuracy of this vicarious calibration approach is not as high as the one based on the direct moon view. However, this "selfcalibration" approach has provided a much needed improvement in the full RVS characterization.
Examples of C6 L1B calibration improvement are presented in Fig. 5. It shows the time series of near-nadir TOA reflectance in band B3 over the desert site Libya 4. While Terra C5 L1B data exhibit an obvious spectrally dependent TOA reflectance trend, the C6 data set indicates that major long-term calibration trends were removed. The bottom plot shows that the respective Aqua B3 record is much more stable and C5 to C6 change is minimal.

Polarization correction of MODIS Terra data
MODIS has no onboard capability to track changes in its polarization sensitivity (Sun and Xiong, 2007). It was optimized during the design phase and was characterized before launch. Prelaunch measurements showed polarization amplitudes increasing toward the higher mirror AOIs and adding up to 2 % for most ocean-color bands, except for band B8 (412 nm), where the amplitude was ∼ 5 %; it was below 2 % for the visible land bands (B3, B1, B4) (Meister et al., 2005).
The early evidence of changing MODIS Terra calibration were obtained by the ocean color team . The comparison of the ocean color products between MODIS Terra and SeaWiFS/MODIS Aqua revealed systematic seasonal and latitudinal differences and trends which can be explained by polarization sensitivity of MODIS Terra and calibration change over time. In the following, OBPG has developed a MODIS Terra-Aqua cross-calibration approach to assess Terra polarization sensitivity as a function of time. This method was originally prototyped and applied to MODIS ocean bands . Recently, OBPG has extended their analysis to MODIS land bands, which makes our present analysis possible.
The cross-calibration approach relies on Level 3 ocean color product (water leaving radiance) and aerosol optical depth produced by the ocean color algorithm from MODIS Aqua, and NCEP (National Centers for Environmental Prediction) wind speed. The MODIS sensor on the Aqua platform had much less and more predictable change in orbit and its calibration stability was well monitored (MCST: http://mcst.gsfc.nasa.gov/calibration/information).
Given atmosphere-ocean properties and assuming lack of diurnal variation over open ocean, one can evaluate "theoretical" TOA radiance under clear skies for the Terra overpass and view geometry (L t ). The measured MODIS Terra TOA radiance (L m ) is related to the theoretical value via elements of the Mueller matrix: where M 11 is an additional gain (RVS) factor and polarization sensitivity coefficients are normalized by M 11 . The components of the Stokes vector, Q and U , are computed for the Rayleigh atmosphere with contribution from ocean glint. To minimize errors related to the glint model, only geometries with glint reflectance < 0.0001 are considered. This approach, applied over the global ocean, produced monthly correction coefficients (M 11 , m 12 , m 13 ) starting from 2002. The polarization correction (PC) coefficients are functions of wavelength, scan mirror side (MS), angle of incidence (or scan angle related to the pixel number along the scan, 1-1354), and band detector. To limit the data set for the 500 m channels (20 detectors), the M coefficients are generated for the neighboring detector pairs, in other words, at 1 km resolution. For the near-nadir geometries, where azimuthal angle changes rapidly, one should also consider rotation angle α, which is the angle between the sensor reference (scan) plane, with respect to which the Mueller matrix is defined, and plane of observation, defined by the normal to the surface and view direction, with respect to which the Q and U components are calculated (for detail, see Meister et al., 2012;Kwiatkowska et al., 2008;Gordon et al., 1997). Once PC coefficients are produced, Eq. (1) is used to obtain the corrected radiance L t from the MODIS Terra measurements L m . For land applications, we account for additional surface height/pressure variations and disregard the effect of glint in Q and U computations. For this work, surface height (pressure) correction is based on linear interpolation between two lookup tables of Q and U generated at relative pressures of P = 0.7 and 1. The radiative transfer was conducted with the atmospheric polarization computations (APC) vector code (Korkin et al., 2013).
The described polarization correction was previously implemented and tested in the MODIS Deep Blue aerosol algorithm (Jeong et al., 2011) which reported an improved agreement with AERONET (Holben et al., 1998) measurements and with the MISR (Multiangle Imaging SpectroRadiometer; e.g., Kahn et al., 2010) global aerosol product over land. Figure 6 illustrates physical effects of MODIS Terra calibration degradation for the mid-Atlantic USA, east coast, DOY 349 (day of year) 2012 (left), and application of polarization correction (right). Each column shows RGB TOA reflectance and MAIAC RGB BRF and AOT 0.47 (aerosol optical thickness). The calibration artifacts appear in aerosol retrievals and atmospheric correction as 10 km stripes representing the residual scan mirror-side difference. The surface BRF image also shows spectral distortions as an artifi- cial purple color. In this case, surface reflectance in the blue band is biased high and exceeds the respective values in the green and red bands, which does not occur in nature. These artifacts are shown near the edge of scan where the magnitude of distortions is largest, however they disappear in the corrected image (right).
Furthermore, PC effect was quantified using MODIS Terra data for the previously mentioned 500 km tile (Fig. 3a) for which 3 L1B data sets were developed and processed by the MAIAC algorithm: (1) original uncorrected C6 L1B; (2) gain (M 11 ) correction only; and (3) full correction. The monthly averaged values over the clear-sky pixels were obtained for the TOA reflectance and for MAIAC AOT and spectral BRF. Analysis of such data set allows us to quantify PC effect on TOA radiance and on MAIAC products, as well as to com-Atmos. Meas. Tech., 7, 4353-4365 pare the relative influence of gain (M 11 ) and polarization (m 12 , m 13 ) factors. Figure 7a shows the effect of full polarization correction (corrected-uncorrected) on MODIS Terra B3 TOA reflectance. The red and black colors refer to mirror sides 1 and 2, respectively. The black curve (MS 2) shows a relatively smooth increase of bias and lack of seasonal effects until ∼ 2006 which indicates that polarization sensitivity is low and gain (M 11 ) correction is the major factor. From 2007, polarization sensitivity rapidly develops over the next 4 years, which is expressed in a strong increase of seasonality. The polarization sensitivity effect is maximal during the winter months when the view geometries with scattering angle near 90 • , where Rayleigh polarization reaches maximum, are most frequent because of the low sun elevation over the horizon. The polarization sensitivity of mirror side 1 is about a factor of 2 lower, and it develops later in ∼ 2009.
We should mention that the presented results for the TOA reflectance, based on monthly cloud-free averages over a large area, are modulated by variable cloud cover, aerosols and heterogeneity of the surface. Thus, it should correctly represent the overall trends in calibration drift and seasonal effect in general, while adding the short-term (monthly) noise from our sampling strategy.
The relative impact of gain (RVS, or M 11 )-only correction vs. polarization correction for B3 can be seen in Fig. 7b built for November of 2012. Here, the black line shows RVS-only correction and the red line shows the full correction for the TOA reflectance. The horizontal axis shows view zenith angle (VZA) from the beginning of scan (BOS) on the left to the end of scan (EOS) on the right. The total correction changes sign from BOS to EOS and is larger by magnitude at EOS. There is also a striking difference in polarization sensitivity between the two mirror sides: while magnitude of polarization correction is almost constant across the scan for MS1, it is negligible at BOS and grows almost linearly towards EOS for MS2.
The red band (B1) analysis (not shown) reveals little difference between the two mirror sides, as well as the practically negligible effect from the polarization correction. The low PC effect in comparison to the blue band is explained by both the reduced magnitude of the m 12 coefficient and significantly lower Rayleigh signal (and polarization) at 0.65 µm.
Finally, Fig. 8 illustrates the PC effect on the MAIAC AOT and BRF in bands B3 and B1. Uncorrected and corrected results from MODIS Terra are shown in red and black colors. For reference, the blue color shows the respective results from C6 L1B MODIS Aqua data. The introduced correction slightly increases aerosol optical depth without a noticeable trend, while the entire MAIAC trend is allocated in the surface BRF. This analysis reveals additional decadal trends of bands, respectively, which was not captured in the current C6 L1B calibration. Thus, while the MCST C6 L1B calibration removed most sensor degradation trends, the OBPG's cross-calibration (polarization) analysis over the ocean has captured more subtle decadal trends at the level of several tenths of 1 % in reflectance units.

MODIS detrending and Terra-Aqua
cross-calibration over desert sites

Calibration analysis over desert sites
While the above analysis showed a clear improvement of Terra C6 L1B calibration after polarization correction, it nevertheless cannot be regarded as final. To yield a more detailed picture of the accuracy of C6 calibration, we applied MA-IAC processing to quasi-stable CEOS-recommended desert calibration sites. As before, we used the full time records from C6 MODIS Terra, with and without polarization correction, and from C6 MODIS Aqua. The 50 × 50 km 2 subsets of MODIS data were provided by the MODIS Adaptive Processing System (MODAPS). To limit the effect of the view geometry variability, MAIAC BRF data were normalized to the fixed nadir view (VZA = 0 • ) and 45 • solar zenith angle (SZA) using the retrieved BRDF model (BRF n ) for each 1 km pixel. Figure 9 shows the time series of clear-sky monthly average BRF n over the subset area in five MODIS bands: B3, B8 (0.412 µm), B2 (0.87 µm), B1, and B4 (0.55 µm) for the site Libya 4. We use the same color scheme as before (e.g., Fig. 8). The top two plots show that application of polar-ization correction improves agreement of Terra BRF n with Aqua in B3 and B8. PC is not applied in the near infrared band (B2) where Terra and Aqua show noticeable trends of the opposite sign. Similarly, it is easy to see that both MODIS Terra and Aqua display nonzero trends in all five plotted bands.
Next, while geometric normalization significantly reduces variability of the surface reflectance, some residual variability of 0.015-0.02, related to seasonality (SZA), remains. Careful inspection of data for B1 and B4 shows that this residual seasonality is more similar between C6 Aqua and uncorrected C6 Terra, while Terra with polarization correction displays "out-of-phase" variations. Moreover, the green band (B4) data show an amplification of the variability with time, the effect being more obvious for some other desert sites. We could not trace the "out-of-phase" variations and the "variability amplification" to any particular cause. As both effects are not desirable, and keeping in mind that the OBPG's correction approach is inherently less accurate at longer wavelengths where the ocean is dark (e.g., B1), a decision was made to limit polarization correction of MODIS Terra to the shortwave bands (B3, B8-B10) only.

Terra-Aqua detrending and cross-calibration
The surface BRF n data cannot be used to remove residual trends shown in Fig. 9, as this procedure should be applied to the top-of-atmosphere reflectance. To solve this problem, we recomputed the expected TOA reflectance (R TOA n ) for the normalized view geometry (VZA = 0 • , SZA = 45 • ) using MAIAC-retrieved parameters, including cloud mask, column  water vapor, aerosol properties, and spectral surface BRDF. The time series of daily area-average R TOA n for bands 3 and 2 is shown in Fig. 10 where the plots on the left show C6 Terra (with PC for B3) and the plots on the right show C6 Aqua data. To avoid sampling bias (due to variable cloudiness, aerosol variability etc.), daily instead of monthly values are used in this case. The TOA normalized reflectance R TOA n provides the required detrending (slope) coefficients for each band as shown in the plots.
This procedure has been applied to seven CEOS desert sites independently. As a result of this analysis, we selected four sites (Libya 1, Libya 2, Libya 4, Egypt 1) which gave relatively similar trends (within a factor of 2-3 difference). Three other sites were excluded for different reasons: Niger showed a factor of 3-5 stronger seasonality resulting in an unreliable trend, while Sudan 1 and Mali 1 produced much larger and opposing trends. The final trend (slope per unit of reflectance per decade) was obtained as an average of the four sites selected for Terra and Aqua separately.
In the next steps, we (1) used this value to obtain the new "detrended" L1B data set, (2) repeated MAIAC processing, and (3) generated the new geometry-normalized R TOA n data set. The "detrended" R TOA n data set is shown in Figure 11 for the same bands 3 and 2. Comparing Figs. 10 and 11, one can see that this procedure effectively reduced the residual spectral trends by a factor of 5-10. With residual trends effectively removed, the R TOA n data, by virtue of geometric normalization, can now be used to obtain a cross-calibration coefficient between two MODIS sensors on Terra and Aqua. This additional gain factor is applied to MODIS Terra since MODIS Aqua has been a better characterized and more stable instrument. As before, the gain coefficient is obtained as an average of the four selected sites. The summary of detrending and cross-calibration gain coefficients is provided in Tables 1 and 2, respectively.

Final analysis with Georgia tile
As a result of the conducted work, the calibration procedure has been extended beyond the C6 L1B level to include OBPG's polarization correction in bands B8-B10 and B3 (for MODIS Terra), detrending (both sensors), and crosscalibration gain adjustment (MODIS Terra), which we further call C6+. As a final test, we applied MAIAC processing to the Georgia (USA) tile to compare different versions of calibration including C5, C6 and C6+ for Terra, and C6, C6+ for Aqua. Figure 12 shows the clear-sky monthly surface BRF n time series for B1-B3 and B8, using the following color scheme: Terra C5 (brown), C6 (red), C6+ (dashed black); and Aqua C6 (blue), C6+ (dashed blue). The difference between the Aqua C6 and C6+ versions is subtle, and can be better distinguished in trend lines except for B3 where the BRF n difference is more apparent. For MODIS Terra, the C5 → C6 difference is large in bands B2, B3, and B8, while the C6 → C6+ change is significant in bands B1, B3 and B8. While the BRF n metric may not be an ideal indicator of success over a selected tile which changes over time from urban development and agriculture, and is affected by the shortterm climate variability, our ultimate criterion is the change of trend lines in the right direction towards the reduction of trends in MODIS Terra and closure between Terra and Aqua C6+ versions for all bands. Plots for bands B3 and B8 show the overcorrection in the C6 L1B version from the introduced RVS trending over the desert sites. In B8, this procedure results in unstable growth of the BRF n seasonal amplitude over time, which is then canceled by the OBPG's polarization correction. The main reason for this instability is the lack of accounting for the sensor's polarization sensitivity during the C6 RVS trending. This emphasizes the need for further improvement of the MCST calibration routine which should simultaneously account for the changes in RVS and in polarization sensitivity of the sensor. It is not clear how this could be achieved, but latest investigations of the Climate Absolute Radiance and Refractivity Observatory (CLARREO; Wielicki et al., 2013) team in collaboration with MCST and OBPG hold promise. Figure 13 shows two more products widely used in the land-applied analysis and modeling -the NDVI (Tucker, 1979) and EVI (enhanced vegetation index; Liu and Huete, 1985). As before, the largest difference for MODIS Terra appears between C5 and C6 versions, while the C6 → C6+ change is smaller. The final C6+ trends from Terra and Aqua (dashed black and blue lines) are practically indistinguishable, which validates the developed detrending and crosscalibration procedure.
Finally, Table 3 shows assessments of decadal changes in NDVI and EVI for different versions of MODIS calibration including C5 Terra and C6+ Terra and Aqua for the Georgia (USA) tile. The total decadal change NDVI from Terra C5 to C6+ is close to 0.01, which is equivalent to the global gross primary production (GPP) change of 1 Pg C (annually) and has significant implications for the global carbon modeling. The Terra-Aqua difference in decadal NDVI changed has reduced by about a factor of 3 in C6+ version.

Conclusions
Aging of Earth observing sensors begins as soon as they start in-orbit operations. This happens for a number of reasons, the main one being exposure to the solar and cosmic radiation. MODIS on Terra has had a more rapid in-orbit degradation accompanied with changes in the response vs. scan angle (RVS) and increased polarization sensitivity. Until ∼ 2007, these changes were not detected through MODIS calibration, and they were not obvious in MODIS Terra science products. This work provides the latest quantitative characterization of trends in different MODIS C5 Terra and Aqua products including DT AOD (over land), global COT, surface reflectance and NDVI/EVI. Due to a longer record, MODIS Terra data (with stronger calibration-related trends) are often used to uncover long-term changes in the Earth system. One of the goals of this paper is to provide a disclaimer for the geophysical trend studies based on the MODIS C5 data set.
The new C6 calibration approach removes major calibrations trends in MODIS Level 1B data. At the same time, analysis by the ocean biology processing group detected changes in the MODIS Terra polarization sensitivity and developed a polarization correction method through Terra-Aqua crosscalibration over clear-sky ocean scenes. Based on MAIAC analysis, we show that the OBPG PC removes the residual scan angle, mirror side and seasonal errors from aerosol and surface reflectance records along with spectral distortions of SR (surface reflectance). Our further MAIAC-based analysis over CEOS desert calibration sites revealed residual decadal trends on the order of several tenths of 1 % in the TOA reflectance in the visible and near-infrared MODIS bands B1-B4 as well as a systematic Terra-Aqua bias. To remove these artifacts, we introduced a MODIS Terra and Aqua detrending and cross-calibration method. Effectively, this very extensive analysis has led to the new C6+ MODIS data set which augments C6 calibration with PC for MODIS Terra bands B8-B10 and B3, followed by detrending of both sensors and by an additional gain adjustment for MODIS Terra to match the Aqua TOA record.
MAIAC science analysis over the southern USA shows that the C6+ version will provide the most reliable MODIS record with the best consistency between Terra and Aqua measurements. The latter will significantly benefit multiple algorithms which rely on the time series analysis and which use or may use the combined MODIS Terra-Aqua record, such as BRDF/albedo algorithm (Schaaf et al., 2002), change detection , MAIAC, etc. The removal of additional negative decadal trend artifacts from Terra NDVI ∼ 0.01 ( EVI ∼ 0.02) has implications for the global carbon modeling and analysis of vegetation dynamics, especially over the tropics (Hilker et al., 2012(Hilker et al., , 2014. Specifically, this result may explain some recently reported trends in gross and net primary productivity or vegetation greenness (e.g., Zhao and Running, 2010). As a result, implementation of the C6+ calibration may help address the problem of a "missing carbon sink" (e.g., Myneni et al., 2001;Pan et al., 2011).