Vicarious calibration of S-NPP VIIRS reflective solar M-bands against MODIS Aqua over dark water scenes

The Visible Infrared Imaging Radiometer Suite (VIIRS) is being used to continue the record of Earth Science observations and data products produced routinely from National Aeronautics and Space Administration (NASA) Moderate Resolution Imaging Spectroradiometer (MODIS) measurements. However, the absolute calibration of VIIRS’s reflected solar bands is thought to be biased, leading to offsets in derived data products such as aerosol optical depth (AOD) when similar algorithms are applied to the different sensors. This study presents a vicarious calibration of these VIIRS bands against MODIS 5 Aqua over dark water scenes, finding corrections to VIIRS between approximately +2 % and -7 % (dependent on band) are needed to bring the two into alignment, and indications of relative trending of up to ∼0.45 % per year in some bands. The derived vicarious gains are also applied in an AOD retrieval, and are shown to decrease the bias and total error in AOD across the midvisible spectral region compared to the standard VIIRS NASA calibration. The resulting bias characteristics are similar to those of NASA MODIS AOD data products, which is encouraging in terms of multisensor data continuity. 10


Introduction
Launched in late 2011, the Suomi National Polar-Orbiting Partnership (S-NPP) satellite is a precursor to the Joint Polar Satellite System (JPSS), which represents the next generation of the USA's operational Earth observation satellites.One of the instruments aboard S-NPP (and the JPSS series) is the Visible Infrared Imaging Radiometer Suite (VIIRS; Cao et al., 2013Cao et al., , 2014)), designed to continue the types of observations made by the Defence Meteorological Satellite Program (DMSP) Advanced Very High Resolution Radiometers (AVHRR, data record 1978 onwards) and Earth Observing System (EOS) sensors such as the Sea-viewing Wide Field-of-view Sensor (SeaWiFS, 1997-2010) and Terra/Aqua Moderate Resolution Imaging Spectroradiometers (MODIS, 2000 onwards).All of these instruments are passive broad-swath imaging radiometers, measuring top of atmosphere (TOA) radiance in a set of reflective solar bands (RSBs) and (except SeaWiFS) thermal emissive bands (TEBs).
Data from DMSP and EOS-era instruments have been used for a broad variety of Earth science applications, including the study of tropospheric aerosols, and a number of algorithms have been developed to create aerosol optical depth (AOD) data methodology applied.Section 4 illustrates the results of the analysis, and Section 5 shows the improvement in retrieved AOD resulting from the calibration exercise by applying the SOAR algorithm to VIIRS scenes passing over AERONET sites.
2 Sensor characteristics MODIS (Barnes et al., 1998, Toller et al., 2013) and VIIRS (Cao et al., 2013, Xiong et al., 2016) are both spaceborne broadswath single-viewing multispectral passive imaging radiometers.VIIRS records data in 22 moderate-resolution bands (Mbands) across the visible and thermal infrared spectral regions with a nominal pixel size of 750 m at the centre of the swath.
MODIS has a total of 36 bands covering the same spectral region, with nominal pixel sizes of 250 m-1 km at the centre of the swath (dependent on band).Each of these VIIRS M-bands has a central wavelength close to one more MODIS bands.
Table 1 shows the band pairs used in this analysis, although MODIS has additional bands, including some others across the visible spectral region of interest here.Note however that some of the MODIS bands designed for ocean colour applications saturate at radiances found over land or cloudy scenes; some of the VIIRS RSBs bands are dual-gain and so do not saturate in many of these cases.In this analysis (and also in MODIS routine aerosol algorithm processing) the relevant so-called MODIS 'land bands' (MODIS B1-B7), where a close match is available, are used instead of these ocean colour bands, even when the latter have a closer central wavelength.For simplicity, wavelengths of MODIS/VIIRS band pairs will be referred to using the notation given in the right column of Table 1, although full sensor relative spectral response functions (RSRs) were used for all radiative transfer calculations presented in this work.The RSRs for the bands used are shown in Figure 1.
VIIRS additionally has five imagery-resolution bands (I-bands) with a nominal pixel size of 375 m and band centres close to some M-band positions, and a Day-Night Band (DNB), which is an enhanced follow-on to the DMSP Operational Line Scanner (OLS) (Lee et al., 2006).Neither the I-bands nor DNB are used in the present DB or SOAR algorithms so will not be discussed further.Likewise, the VIIRS and MODIS TEBs will not be discussed further.
As mentioned previously, stability of MODIS/VIIRS RSBs is monitored and maintained using the SD, SDSM, and lunar rolls (Sun et al., 2007, Wu et al., 2013, Xiong et al., 2016).As a result the RSB absolute calibration for each band is tied to the measured reflectance (ρ i ) rather than radiance, where . (1) In the above L is the spectral radiance passing into the satellite field of view, E the downwelling solar spectral irradiance at TOA, and Φ i the sensor RSR for band i, all functions of wavelength λ.The factor D is the Earth-Sun distance in astronomical units (variable throughout the year), and µ 0 the cosine of the solar zenith angle, which affect the total solar radiation received.
L and consequently ρ i depend on the surface/atmospheric states and observation geometry, omitted here for simplicity of notation.Equation 1 is often simplified by considering the total radiance observed by a band and precomputing the spectrally-integrated solar irradiance across the band, which then only varies as a function of day of year, leading to the more common form The radiative transfer codes used in the DB and SOAR algorithms operate in units of Sun-normalised radiance, I/F 0 , to minimise numerical instabilities at large solar zenith angles (i.e. as 1/µ 0 tends to infinity).As a result the discussion in the present study also uses units of Sun-normalised radiance, which does not affect the adoption of the results presented herein for other applications.Working in Sun-normalised rather than total radiance also has the advantages of accounting for the effects of the sensors' different RSRs on F 0 for each band, and slightly different solar zenith angles at the times of MODIS and VIIRS 10 observations.
Both sensors suffer from a 'bow-tie distortion' which affects the size, shape, and overlap of pixels from nadir to scan edge (Wolfe et al., 2012(Wolfe et al., , 2013)).Essentially, as the detector scans across-track pixels become broader and elongated, and pixels  from consecutive scans overlap, which has consequences for retrieval characteristics as a function of scan angle, and can affect statistics of AOD retrievals (Sayer et al., 2015a).VIIRS incorporates several design features to reduce this distortion.The VIIRS native pixel size is actually smaller than the nominal M-band size in the across-track direction.The scan is divided into three regions (in both directions).From nadir out to a scan angle of 31.72 • , three pixels are aggregated across-track; from 31.72 • -44.86 • two pixels are aggregated, and from 44.86 • to the edge of scan (56.28 • , corresponding to a view zenith angle around 75 • ) no aggregation is performed.This limits across-track distortion at the end of each aggregation zone to a factor of two, compared to a factor of about six without this oversampling and aggregation.Additionally, at the outer two aggregation zones, two and four pixels respectively are deleted from the edge of scan (so-called 'bowtie deletion') to minimise the degree to which consecutive scans overlap (although not all overlap is removed by this).S-NPP and Aqua are both in Sun-synchronous orbits with daytime Equatorial local solar crossing times at center of swath is around 13:30 UTC.However, orbit altitudes (averages of 705 km for Aqua and 839 km for S-NPP) and inclination are different which mean that the sensors are not always observing the same location (near-simultaneous observations within 10 minutes occur roughly every other day, covering only a segment of the world each time).Both satellites' orbital repeat cycles are 16-day.
The MODIS swath width is approximately 2,330 km, providing near-global daily daytime coverage (there are gaps between consecutive orbits at low latitudes) while VIIRS has a sufficiently broad swath (3,040 km) that consecutive orbits overlap, even at the Equator.Both sensors have some degree of overlap between consecutive orbits at mid-and high latitudes.

Vicarious calibration methodology
This analysis seeks to calibrate VIIRS bands M01-M11 (spectral range 412-2,250 nm) against the corresponding MODIS Aqua bands shown in Table 1.Band M09 is not considered; this band (for both MODIS and VIIRS) is located in a spectral region of strong water vapour absorption, and so is typically used in threshold tests to detect high cloud tops (cirrus or deep convective clouds) rather than in geophysical retrieval algorithms (e.g.Frey et al., 2008), and the quality of its absolute calibration in MODIS is unclear.The steps of this exercise are as follows, and outlined in detail below: 1. Selection of appropriate MODIS/VIIRS pixels to consider.
2. Correction for the effects of absorption by trace gases in the atmosphere.
3. Forward radiative transfer modelling to predict the TOA reflectance which should be observed by VIIRS, given that observed by MODIS.
4. Aggregation of results to a monthly time scale and derivation of vicarious calibration coefficients.

Data description and selection of appropriate pixels
The NASA VIIRS data processing is facilitated by Science Investigator-led Processing Systems (SIPS) for each science discipline.In support of this, the VIIRS Atmospheres SIPS at the University of Wisconsin (http://sips.ssec.wisc.edu)have created 'matchfiles' of collocated VIIRS and MODIS Aqua observations, to more easily compare the two sensors.These matchfiles form the basis of the present analysis, and have been created from the MCST MODIS Aqua Collection 6 data (at 1 km nominal pixel size, the MYD021KM product) and VCST VIIRS Version 1.1 data with Version 1.0.1 calibration (the VL1BM product).
These are the current versions of the L1 data for both sensors used in routine processing of MODIS, and which will be used for the first processing of the VIIRS Deep Blue/SOAR L2 data products.The time period considered is from March 2012 (several months after the S-NPP launch, at which point the VIIRS M-band RSB/TEB data were considered ready for use) to May 2016 (the last date for which the VIIRS Version 1.0.1 calibration has been generated at the time of writing) inclusive.
The matchfiles contain L1 RSB/TEB data, geolocation, and land/sea mask information for MODIS and VIIRS, with VIIRS pixels mapped into MODIS pixels (as, due to a combination of native spatial resolution and bow-tie distortions, VIIRS pixels are typically smaller than MODIS ones).Because of this, the matchfiles contain the mean and standard deviation of VIIRS RSB reflectance/TEB radiance within the area of each MODIS pixel, as well as that corresponding to the nearest VIIRS M-band pixel to the centre of each MODIS pixel.The files also contain the MODIS Collection 6 cloud mask for each MODIS pixel (the MYD35 data product; an updated version of that described by Frey et al., 2008).
Two additional ancillary data sets were used in the analysis.The first is the Goddard Earth Observing System Model, Version 5 (GEOS-5) Forward Processing for Instrument Teams (FPIT), available from http://gmao.gsfc.nasa.gov/GEOS,at Atmos.Meas. Tech. Discuss., doi:10.5194/amt-2016-238, 2016 Manuscript under review for journal Atmos.Meas.Tech.Published: 17 August 2016 c Author(s) 2016.CC-BY 3.0 License.0.5 • latitude, 0.625 • longitude, and 3-hour temporal resolution.Surface winds and O 3 and H 2 O total column abundances were extracted and interpolated linearly to each pixel in the matchfiles.The second is climatologies of oceanic chlorophylla concentration (Chl) derived from the SeaWiFS record, available from http://oceancolor.gsfc.nasa.gov.These are provided at a native spatial resolution of 9 km, and one data set for each of the twelve months of the year; to fill gaps, for the present analysis these were degraded to 0.25 • spatial resolution by taking a median average, and remaining gaps filled using the nearest available month of the year in time (i.e. for a gap in March, data from February and April are used; if also missing then data from January and May, etc).The resulting gap-filled climatology was also interpolated linearly (in units of log 10 Chl) to each pixel in the matchfiles.
Matched MODIS/VIIRS pixels were then selected for further analysis if they passed all of the following criteria: -Difference in observation time under 10 minutes, to minimise changes in the surface/atmospheric state between observations.
-View zenith angle and scattering angle differences both smaller than 3 • , to minimise uncertainties related to the geometric dependence of the scenes viewed (e.g.surface reflectance, atmospheric absorption, and scattering phase matrices).Note that the observation time threshold mentioned above also effectively acts as a threshold on solar zenith angle difference.
-MODIS land mask classifications 'deep inland water', 'moderate or continental ocean', or 'deep ocean', as well as climatological Chl < 1 mgm −3 , to remove cases where the ocean surface reflectance model used (see later) may be less appropriate.
-MODIS cloud mask classification 'confidently clear', and at least 5 km away from any pixel classed as 'confidently cloudy' or 'uncertain'.This distance threshold removes approximately 80 % of 'confidently clear' pixels, although decreases the likelihood of errors resulting from cloud movement between the MODIS/VIIRS overpass times, classification errors, and 3D radiative transfer effects.
-Solar zenith angles smaller than 70 • , to minimise shadow length and parallax effects, and ensure a strong daytime RSB signal.
-Latitude Equatorward of 60 • , as ship-based observations suggest that background oceanic aerosol optical characteristics at polar latitudes can differ from those at lower latitudes (Smirnov et al., 2011).
-Sun glint I/F 0 (from GEOS-5 winds and the model of Cox and Munk, 1954a, b) <0.01 for both sensors, as over-ocean radiative transfer modelling is subject to higher uncertainties in glint hotspots.
-Total column H 2 O less than 3 cm, to decrease uncertainties related to trace gas absorption (because this can have a large absorption in nIR and swIR bands and exhibits fairly large spatiotemporal variability).Note that, although the bulk of matched pixels occur in the mid-and high latitudes, this constraint removes most of the potential matchups in tropical regions (since these often have water vapour amounts in excess of this threshold).
Figure 2 shows an example of the results of applying the above filters to matchfiles created from a pair of MODIS/VIIRS granules, separated by 1 minute in acquisition time.

Correction for trace gas absorption
The next step is to correct the TOA MODIS/VIIRS RSB reflectances for the effects of absorbing trace gases.Note that the corrections described here are also applied in NASA VIIRS DB/SOAR processing.The assumption normally made in AOD retrieval algorithm processing is that the molecular absorption can be decoupled from other contributions to the TOA signal, and so corrected for by applying a 'brightening' factor to the observed TOA reflectances.In the visible spectral region, this is justified because the optical depths of the absorbing species are fairly small, and (particularly in the case of O 3 ) the bulk of the absorption is located higher in the atmosphere than the main other contributors to the signal (Rayleigh scattering, aerosols, and surface reflectance), so corrections for this absorption can be developed with high accuracy (Gueymard, 1995).
The computational advantage of performing such a correction is that it vastly decreases the dimensionality of radiative transfer lookup tables (LUTs) used in the retrieval process, since individual gas species and their variable vertical profiles do not need to be built in to them.
This analysis uses the same approach taken in operational MODIS aerosol processing (Appendix A of Levy et al., 2013), with the additional step that, following Tanré et al. (1992), the effective column H 2 O amount is taken as half the total column H 2 O to better account for the atmospheric vertical structure.In brief, gas absorption is calculated the Line-By-Line Radiative , since their contributions to the total absorption are weak and/or their spatiotemporal variability comparatively low.Note that NO 2 variability, and consequently absorption, can be significant (Ahmad et al., 2007) for sensor bands in the blue spectral region (412 nm, and to a lesser extent, 440-490 nm).However the exclusion of pixels close to land masses (Section 3.1) means that this is not an issue for the data considered here, since background oceanic levels are low due to a lack of strong sources, and the short lifetime of NO 2 means that long-range transport is fairly small.
The gas absorption corrections are calculated for and applied to each pixel and band.After applying these corrections, the MODIS/VIIRS RSB data are effectively that which would be seen by the sensors in the absence of these trace gases in the atmosphere, removing one cause of differences in TOA reflectance between the two instruments.This step is important because, despite the similarity of band central wavelengths (Table 1), RSR shapes can be sufficiently different (Figure 1) that differences in the level of gas absorption for nominally similar bands can be non-negligible in some cases.This is illustrated in Figure 3: while some bands, such as MODIS B4/VIIRS M04 (both centred near 550 nm), are very similar and gas corrections are tightly correlated and close to 1:1, others show more difference in magnitude and/or spread, illustrating the importance of accounting for trace gas absorption accurately when comparing L1 data from the two sensors.

Forward radiative transfer modelling of predicted VIIRS reflectance
The next step is to determine, given the observed MODIS TOA reflectance, what reflectance VIIRS should see for each pixel.
Because of the differences between sensor RSRs (Figure 1), and the slight differences in observation geometry between the two sensors for a given pixel, this requires a radiative transfer forward model, and the results of the analysis will be sensitive to the assumptions made in that forward model.This analysis uses the VLIDORT radiative transfer model (Spurr, 2006), which is the same as is used in the NASA VIIRS SOAR data set, and allows for a detailed description of aerosol properties and surface bidirectional reflectance distribution function (BRDF), a pseudospherical atmosphere, as well as a vector treatment of the atmospheric radiation field, which is important for accurate radiative transfer at short visible wavelengths (e.g.Levy et al., 2004).It is also able to account for the full RSRs of the sensors when performing calculations.This radiative transfer model has some advancements over those used previously by DB/SOAR, has been benchmarked against standard results with good performance (Spurr, 2006), and versions have also been used for other aerosol remote sensing applications (e.g.

Aerosol optical model
As the comparison is restricted to open ocean scenes, it is a reasonable assumption that most of the AOD is contributed by 'clean' (i.e.little continental influence) maritime aerosols (e.g.O'Dowd and de Leeuw, 2007).For this reason, the 'pure marine' aerosol optical model of Sayer et al. (2012b) is used.This model was based on AERONET inversions (Dubovik and King, 2000) from a variety of sites, was applied previously in SeaWiFS SOAR processing (Sayer et al., 2012a), and is also applied in VIIRS SOAR data processing.Real and imaginary aerosol refractive indices were taken from Hess et al. (1998), as there are few measurements of aerosol optical properties across the whole VIIRS spectral range.Specifically, the fine mode uses the 'water soluble' component refractive indices and the coarse mode the 'coarse mode sea salt' component, both for aerosols at 70 % relative humidity.In the radiative transfer simulations the aerosol is assumed to occupy a homogeneous vertical layer from the surface to 1 km altitude; as the aerosol is close to nonabsorbing, and the data are further filtered for low-AOD conditions (discussed further later), the vertical structure has little influence on the modelled TOA signal.
The standard assumption made in the analysis is that the aerosol fine mode fraction (FMF) of optical depth at 550 nm is 0.4, which is a typical value determined from observations in a variety of global oceans (Smirnov et al., 2011).However, to assess the uncertainty resulting from this assumption (discussed later), radiative transfer simulations are also performed for more extreme cases of FMF=0.1 and FMF=0.7.

Surface reflectance model
The ocean surface BRDF is an updated version of the treatment used by Sayer et al. (2012a) for SeaWiFS, and the same model discussed herein is also applied for SOAR VIIRS processing.The BRDF model draws on the widely-used method of Koepke (1984), and includes contributions from oceanic whitecaps, sun glint, and scattering from within the water ('underlight', using the basic formalism of Austin, 1974).Both the whitecap and underlight terms have been updated since the SeaWiFS application.
The wind speed dependence of the whitecap formulation has been updated using the formulation of Callaghan et al. (2008), which tends to slightly decrease the whitecap contribution to the BRDF at most wind speeds, since Callaghan et al. (2008) and other studies suggest that the older formulation used previously (Monahan and Muircheartaigh, 1980) may overestimate the whitecap fraction.
Underlight is calculated using an empirical relationship based on Chl to estimate absorption and scattering from pigments and covarying materials.This relationship was developed for so-called 'Case 1' (largely open-ocean) waters (Morel and Prieur, 1977).Within the underlight component of the reflectance model, several updates have been made to the assumed water absorption/scattering properties (previously taken from Smith and Baker, 1981).Lee et al. (2015) found that prior estimates of the absorption coefficient of water at visible and ultraviolet wavelengths were too high, and so the Lee et al. (2015) coefficients have been adopted instead over their available spectral range (300-550 nm).Pope and Fry (1997) is used for 550-725 nm, which results in smooth continuity with the results of Lee et al. (2015), and Hale and Querry (1973) is used for longer wavelengths (although above 700 nm water absorption is so strong that ocean reflectance depends negligibly on chosen data source).The water scattering coefficient was also updated according to Zhang et al. (2009), and chlorophyll absorption spectrum updated according to Lee et al. (1998) and Vasilkov et al. (2005).At the same time, the spectral range of the parametrisation has been extended to 300-900 nm (from the prior 400-700 nm).Directional aspects of the underlight contribution (so-called f /Q ratio) have also been updated, according to Morel et al. (2002).
The combined effect of these coefficient updates, relative to prior implementations of the same basic model (Sayer et al., 2010(Sayer et al., , 2012a) is an increase of up to a few tens of percent in the underlight contribution to ocean reflectance for the blue and green spectral region (550 nm and shorter wavelengths), which translates to a few percent in TOA reflectance.

LUT creation and application
To use MODIS observations as a predictor for VIIRS, VLIDORT has been used to construct a pair of LUTs of MODIS and VIIRS reflectance for a variety of surface and atmospheric conditions for each band.The node points are shown in Table 2, and their spacing has been chosen such that the linear interpolation error between node points is less than 1 % (relative) error in Sun-normalised radiance, with the average bias across conditions negligible.
The LUT is used by looping over each matched MODIS/VIIRS pixel pair and band, using the measured MODIS TOA Sunnormalised radiance to estimate the AOD at the reference wavelength of 550 nm (based on the ancillary MODIS geometric information, wind speed, and Chl climatology), and then using the derived AOD (together with the ancillary VIIRS geometric information to predict the Sun-normalised radiance VIIRS would be expected to see if its absolute calibration were equal  to that of MODIS.For each pixel, this process is repeated for each band independently, and for each of the three aerosol FMF assumptions (0.1, 0.4, 0.7).This dynamic AOD estimation, rather than assuming e.g. a single AOD across all scenes, ensures that the spectral MODIS TOA Sun-normalised radiance is matched exactly for each pixel and band, and decreases the uncertainties involved in cross-calibrating the two sensors.One important point to note is that, as the calibration of VIIRS is being tied to that of MODIS, what is most important here is not so much the absolute accuracy of the radiative transfer modelling, but the accuracy of the spectral/directional extrapolation between MODIS and VIIRS wavelengths and geometries.
A further filtering step takes place at this stage.Pixels are only retained if, for each band and for each of the FMF assumptions, an exact match to the MODIS reflectance is found within the range of AOD stored in the LUT (i.e.0-0.2, Table 2).This removes residual cases where the forward model might not be appropriate, e.g.residual cases of contamination by clouds or continental (e.g.smoke, dust) aerosols where the aerosol optical model assumption may be significantly in error, or cases where the ancillary data (Chl, wind speeds, trace gas abundances) are significantly in error.This decreases the available data volume, but helps to ensure that the remaining pixels correspond to clean open-ocean cases where the transfer between MODIS and VI-IRS spectral and directional characteristics has been achieved with high fidelity, and retains identical spatiotemporal sampling for all bands.
At the end of this stage, each suitable pixel and band have associated MODIS and VIIRS TOA RSB observations, plus an estimate of the TOA signal which VIIRS would be expected to see, under the assumption that MODIS Aqua's calibration is correct.An example of the results, composited from the points obtained in January 2016, is shown in Figure 4. Figure 4(a) shows the expected spectral dependence of the TOA signal over open ocean -namely, a darkening as wavelength increases, because Rayleigh scattering, oceanic surface reflectance, and aerosol scattering all decline at longer wavelengths.MODIS at 470 nm is notably brighter than VIIRS at the VIIRS band near 490 nm, and to a lesser extent for the 650/670 nm pair, for these reasons.The theoretical predicted VIIRS TOA signal is much closer to the observed VIIRS than observed MODIS for these bands, illustrating again the importance of accounting for the differences in sensor RSRs rather than just comparing the TOA signals directly (i.e.differences resulting from sensor spectral characteristics may be larger than those resulting from sensor calibration errors).The lower panel (Figure 4b) shows the median and variability of the estimated AOD at 550 nm (for the reference FMF=0.4 case) over this month, calculated independently for each band from MODIS observations, as described above.The variability on each point reflects both the real spatiotemporal variability in AOD through the month, as well as the effects of both MODIS calibration errors, and those uncertainties resulting from errors in the forward model assumptions and ancillary data sets used on the analysis.Larger variability is expected for the shorter and longer wavelengths due to increased sensitivity to 5 assumptions about ocean surface reflectance in the former case, and trace gas corrections in the latter case, as well as the fact that assumptions on FMF at 550 nm become more significant at wavelengths further from 550 nm.An attempt at quantification of the effects of these errors on derived vicarious calibration gains is presented in Section 4. In an ideal case, the same average AOD at 550 nm would be returned for all bands; for this particular month the average is 0.071, with a standard deviation between bands of 0.015, similar to the typical level of uncertainty on midvisible AOD obtained from AERONET standard

Aggregation to monthly time scales, and calculation of vicarious calibration correction
The final step in the analysis is to aggregate the pixel-level results to a monthly time scale, and use these to derive vicarious calibration coefficients.The rationale for a monthly time step is that results derived from observations on a single day are likely to have correlated errors (in terms of forward model and ancillary data), as they are drawn from a limited spatial and temporal snapshot of the world.In contrast, averaging to a month should provide sufficient sampling that errors can be averaged out to a large extent.Aggregation also minimises the influence of remaining outliers.At the same time, monthly time scales remain sufficiently short that any longer-term behaviour, such as seasonality or drifts in the relative calibration, can be examined.
For each month and band, the remaining pixels are sorted by the theoretical VIIRS I/F 0 (i.e. the signal VIIRS is expected to report, given the MODIS observations and the spectral/directional differences between the sensors' observations, as modelled through VLIDORT) and divided into 50 equally-populated bins.For each bin, the median theoretical and actual VIIRS I/F 0 10 are recorded; use of medians rather than means decreases the sensitivity to outliers caused by real scene changes or radiative transfer errors.An example of this process for January 2016 is shown in Figure 5.
The vicarious calibration gain correction (to make VIIRS radiatively consistent with MODIS Aqua for a given Earth scene) derived for this month of data is simply the mean ratio between these binned values, and VIIRS can be 'corrected' to be This ratio approach is inherently making the assumption that the vicarious calibration correction is a simple gain scaling factor, that it is linear, and that there is no offset between the two for the darkest scenes; this is generally expected to be the case, and as a result the assumption is common in such vicarious calibration exercises (Franz et al., 2007, Cao et al., 2013, 2014, Uprety et al., 2013, 2014, Lyapustin et al., 2014, Uprety and Cao, 2015, Wang and Cao, 2016).If linear least-squares regression is performed on the binned data (not shown) rather than simply taking the mean of the ratios between the bin-median values, then the effective gain coefficients are similar and offsets are close to zero, although there is a little more month-to-month variation since two free parameters are being determined (offset and gradient), and a deviation in one of these parameters is countered by a deviation of opposite direction in the other.Thus, the numerical effects on corrected TOA reflectance, or I/F 0 , is negligible whether the bin-ratio or linear regression technique is used.
Further, linear coefficients of determination between the binned theoretical and actual VIIRS signals are close to unity (R 2 > 0.995), indicating the assumption of linearity is justified in this case.Note that this does not preclude detector nonlinearities across the whole range of VIIRS brightnesses, because the present analysis is restricted to scenes over ocean, which are (particularly for the swIR bands) fairly dark.

Derived gains, time series, and uncertainty
Performing the steps detailed in Section 3 results in data distributed regionally and seasonally as shown in Figure 6.Some regions are sampled frequently and others never, due to the intersection of the two satellites' 16-day repeat cycles.Seasonal variations are caused predominantly by variations in cloud cover and land (more land in the Northern Hemisphere), and secondarily by changes in solar angle (which affects the Sun glint location, latitudes of daylight, as well as scattering angle differences).
The distribution of these points in time is shown in Figure 7.As well as seasonal variability (caused by the aforementioned factors), the available data volume is larger from 2014 onwards than in 2012 and 2013.This can be explained by the satellites' orbital times.Figure 8 shows the Equatorial local solar crossing times of the ascending (i.e.daytime) nodes for both sensors.
While both are often quoted as a nominal 13:30 UTC crossing time, neither orbits at exactly this time.Aqua's orbit crosses around 13:35 UTC, with small seasonal variability, and is very tightly controlled as it flies as part of the A-Train constellation.
S-NPP was initially closer to 13:25 UTC and this has gradually changed through the mission as a result of spacecraft orbital adjustments, with the two platforms' crossing times within 5 minutes from mid-2014 to mid-2015.As a result, the 10 minute time difference threshold imposed on the analysis (Section 3.1), imposed to minimise the influence of changes in the scene viewed between the two sensors' overpasses is comparatively more restrictive in 2012 and 2013.Note that these changes in crossing time remain within the missions' tolerance requirements.It is also worth noting that these changes in Equator crossing time do not strongly influence the geographic distribution of matched pixels.Compositing the monthly results (e.g. Figure 5) gives the time series shown in Figure 9.The three swIR bands show more noise than the others, which is expected since both ocean and atmosphere are quite dark at these wavelengths, and the resulting ratios of small numbers are inherently less stable than ratios for brighter bands.Several bands show seasonal oscillation and it is unclear at present to what extent this is caused by both sensors' individual radiometric stability, and to what extent this may be related to seasonal changes in geographical sampling (as it is possible that assumptions made in the analysis may be more 5 or less appropriate in different regions, leading to residual geographic error).Nevertheless, these seasonal oscillations, where present, tend to be small (amplitude <0.01, i.e. <1 %, at all wavelengths).In addition, several bands show an apparent trending in the vicarious gain coefficient over the five-year period.For band M01, this is equivalent to VIIRS becoming increasingly relatively brighter than MODIS; in this specific instance, it is thought that VIIRS is the more stable of the two sensors (B.Franz/G.Meister, personal communication, 2016).For the other bands, it is not clear at present which, or both, of the sensors is degrading.For VIIRS bands M07, M08, and M10 the change over the five-year period exceeds 1 % and so the temporal dependence is probably worth accounting for until residual trending of both sensors can be analysed and corrected for by the respective instrument teams.Stability of both sensors is monitored using the on-board solar diffuser and SDSM, as discussed previously; for the bands in question, additional polynomial detrending analyses were performed and implemented for the MODIS C6 reprocessing (Doelling et al., 2015), although the most recent years of MODIS data had not yet been collected at that time, so it is possible that any additional degradation has deviated from these prior models.
The mission-averaged vicarious gains are shown in Table 3, along with (for bands M07, M08, and M10) linear trends.These trends were calculated from least-squares linear regression of the monthly gains with the time ordinate taken as years since the start of January 1, 2010, and the uncertainties presented with these parameters are the standard least-squares linear regression uncertainties.A linear model was used based on visual examination of the data, although there is no particular reason to expect a linear change as opposed to any other specific functional form, or that this behaviour will continue in years to come, so these trends should be interpreted with caution.For the mission-averaged gains, the total uncertainty σ tot was estimated as the quadrature sum of four components: These components were estimated as follows: This section illustrates the results of applying the vicarious calibration gains in Table 3 (including time-dependence, for the relevant bands) to VIIRS data and processing them through the SOAR retrieval algorithm, to illustrate the effects on derived aerosol properties.The VIIRS application of SOAR will be described in a subsequent study, although it is basically an extension of the SeaWiFS application (Sayer et al., 2012a) to incorporate some of the additional features of VIIRS.Like many others, SOAR is a multispectral inversion using LUTs of physically-based radiative transfer results (e.g.Stowe et al., 1997, Tanré et al., 1997, Mishchenko et al., 1999, Thomas et al., 2009, Sayer et al., 2010, Jackson et al., 2013).The algorithm uses VIIRS bands M03, M04, M05, M08, M08, M10, and M11 (i.e.seven of the 10 bands analysed in this study) and provides AOD, FMF, and an indication of best-fit aerosol optical model at a nominal pixel size of 6 × 6 km 2 (8 × 8 M-band pixels).Spectral AOD is determined through the retrieved AOD at 550 nm, together with the retrieved FMF and aerosol optical model.The specific details are of secondary importance here, as the main purpose is to illustrate the effects of the calibration change on the retrieval.To demonstrate these effects, the SOAR algorithm has been applied to VIIRS granules passing over six AERONET sites (Table 4), using the standard NASA L1b products with and without the vicarious gains developed in this study.Comparing the AOD retrievals with AERONET enables a characterisation of both how much the spectral AOD retrieval is affected by the vicarious calibration, and whether these changes have improved the retrieval or not.
Note that SOAR is a multispectral inversion, fitting all bands simultaneously; the underlying radiative transfer is nonlinear in AOD, and the bands are not weighted equally.As a result, changes in an individual band's calibration do not map linearly into retrieved AOD at a given wavelength, and it makes the most sense to analyse the behaviour of the retrieval system as a whole, rather than attempt to assess or infer the effect of changes to the calibration of individual bands.These sites have been chosen based on their locations (coastal/island AERONET sites with typically fairly low AODs, such that aerosol optical model assumptions will be less important contributors to the retrieval error characteristics) and the fact that, for at least part of the VIIRS mission, the sun photometers deployed at these sites included a 1640 nm filter, so that AOD can be validated across the whole spectral range of bands used by SOAR-VIIRS with minimal extrapolation (the majority of AERONET sun photometers lack filters at wavelengths longer than 1020 nm).The analysis protocol is the same as in Sayer Figure 10 provides summary validation statistics for AOD at 550 nm, composited across the AERONET sites.Using the standard L1 calibration the bias is 0.0313, similar in magnitude to that of the NOAA VIIRS ocean AOD products (Huang et al., 2016); applying the vicarious gains removes a little over half of this bias (about 0.016), and also gives a root mean square (RMS) error smaller by about 0.01.It additionally brings another 13 % of matchups in agreement with AERONET to within the data set's expected level of uncertainty.The resulting bias on midvisible AOD is similar to or smaller than preliminary validation results for the current Collection 6 MODIS over-water AOD retrieval algorithm, which is conceptually similar to SOAR (Sayer et al., 2012c, Levy et al., 2013).The correlation coefficient is similar, indicating that the effect is more a shift in the AOD distribution than a change in the scatter.The total data volume changes slightly, as the calibration change affects cloud masking and quality assurance parts of the SOAR algorithm in the two runs (Sayer et al., 2012a).
Figure 11 presents some similar summary statistics, but for spectral AOD rather than AOD at 550 nm.Only the 1649 points where both runs provide a valid matchup are considered in this figure (and these provide the counts listed in Table 4).Similar    4.
accounting for the differences between the sensors' spectral response functions and viewing geometries, to tie the VIIRS calibration to that of MODIS.
The analysis suggests that the standard NASA L1 (version 1.1.0)VIIRS data require scaling by between approximately +2 % and -7 %, depending on band, to bring them into agreement with MODIS Aqua, with indications of relative trending of up to ∼0.45 % per year (over the March 2012-February 2016 time period analysed) in some bands.The relative contribution 5 of the two sensors to these drifts in relative calibration is not yet clear.The precision of these gains is typically of order 0.5 % for visible and nIR bands, and 1 % for swIR bands, due to a combination of trending and increased noise in the latter cases.
Application of the vicarious gains to the SOAR algorithm, to derive spectral AOD over water pixels around six AERONET sites, Atmos.Meas. Tech. Discuss., doi:10.5194/amt-2016-238, 2016 Manuscript under review for journal Atmos.Meas.Tech.Published: 17 August 2016 c Author(s) 2016.CC-BY 3.0 License.illustrates that the vicarious gains do provide an improvement in retrieved AOD in the visible spectral region (470-865 nm), although do not address a bias in AOD retrieved at swIR wavelengths.
Even if MODIS Aqua's calibration is imperfect, as any sensor's necessarily is, this analysis is consistent with prior indications that it is likely to be better than that of VIIRS.Importantly, obtaining radiatively-consistent L1 data increases the likelihood of similar error statistics in downstream L2 data products, which facilitates the creation of long-term data records by combining individual sensors, by minimising discontinuities between products.Additionally, the analysis technique is independent of the AERONET validation data used to evaluate the AOD retrieval algorithms, so the effect of the calibration changes can be evaluated in a way that is not circular, and also has the advantage that the vicarious gains can be applied directly to other retrieval algorithms (i.e.errors in the SOAR algorithm or AERONET data are not aliased into the vicarious gains).Further, the analysis can be easily repeated as additional versions of the source MODIS and VIIRS L1 data become available, to assess whether the offsets, and relative stability, of the sensors has changed.Conceptually it could also be applied to other sensors with frequent orbital overlaps-although, as seen in this analysis, even for sensors with nominally similar orbital overpass times, a shift in a few minutes can have a large effect on the data volume available for analysis.
Several caveats remain.The scenes analysed here (cloud-free oceans) are for brightnesses typical for aerosol retrievals, i.e. fairly dark, and detector nonlinearities may mean that for very bright scenes (e.g.optically-thick clouds or snow) the relative offsets between the two sensors may differ.A second main caveat is that, while similar, VIIRS and MODIS' spectral response functions are different, and these differences are large enough that the calibration exercise must involve the use of radiative transfer models, as the expected differences arising solely from these shifts in spectral response can be larger than the vicarious calibration gains.The effects of this are more severe the more strongly the relevant sensor bands differ.This analysis has accounted for these spectral response differences, which is possible since the spectral dependence of atmospheric and surface scattering and absorption can be accounted for with high accuracy over these cloud-free low-aerosol oceanic scenes, although these differences do contribute to residual error and uncertainty in the derived vicarious gains.The well-calibrated hyperspectral Infrared Atmospheric Sounding Interferometer (IASI) has been used in the thermal IR domain to investigate the calibration and spectral response of other sensors (Sohn et al., 2010, Goldberg et al., 2011), as its hyperspectral bands can be combined to mimic closely the broader bands of MODIS and other sensors; there is no current equivalent, however, for the visible through to swIR spectral domain.The proposed Climate Absolute Radiance and Refractivity Observatory (CLARREO) Reflected Solar Spectrometer (Wu et al., 2015), once available, will add this important capability.

Figure 1 .
Figure 1.Relative spectral response functions for VIIRS (red) and MODIS (blue) bands used in this study (cf.Table 1).

Figure 2 .
Figure 2. Example MODIS/VIIRS match up for two near-coincident granules (beginning one minute apart).The S-NPP VIIRS granule is outlined in red, and MODIS Aqua in blue.Suitable matched pixels are shown in green.
Atmos.Meas.Tech.Discuss., doi:10.5194/amt-2016-238,2016   Manuscript under review for journal Atmos.Meas.Tech.Published: 17 August 2016 c Author(s) 2016.CC-BY 3.0 License.Transfer Model (LBLRTM,Clough et al., 2005), which includes the High Resolution Transmission (HITRAN) gas absorption data base (http://hitran.org), is used in combination with the MODIS and VIIRS RSRs to calculate the absorption by atmospheric trace gases as a function of their amount and vertical profile shape.The effective air mass factors for absorption, which includes the effects of the Earth's curvature and typical gas profiles, are obtained fromGueymard (1995).GEOS-5 data (discussed previously) are used to obtain O 3 and H 2 O column amounts, while climatological abundances are used for the other

Figure 3 .
Figure 3. MODIS/VIIRS spectral trace gas absorption corrections for suitable matched pixels in January 2016.

Figure 4 .
Figure 4. (a) Median spectral observed TOA Sun-normalised radiance for MODIS (blue) and VIIRS (red) bands, and (black dashed) modelled expected VIIRS signal, after the AOD estimation and filtering steps described in Section 3.3.3.(b) Median (points) and central 68 % range (lines) of AOD at 550 nm, estimated from MODIS for the assumed marine aerosol optical model and FMF=0.4.In both panels, data are composites of the 1,482,794 valid points obtained in January 2016.
Atmos.Meas.Tech.Discuss., doi:10.5194/amt-2016-238,2016   Manuscript under review for journal Atmos.Meas.Tech.Published: 17 August 2016 c Author(s) 2016.CC-BY 3.0 License.unbiased with respect to MODIS Aqua by multiplying the TOA signal (whether in reflectance or I/F 0 units) for the band in question by this number.

Figure 6 .Figure 7 .
Figure 6.Seasonal spatial distribution of matched MODIS/VIIRS pixels used for the vicarious calibration exercise, aggregated to a 1 • horizontal grid size.

Figure 8 .
Figure 8. Variation of S-NPP (red) and Aqua (blue) satellite Equatorial local solar crossing times for the ascending (daytime) orbital nodes.The dashed line indicates the nominal 13:30 UTC crossing time often used as shorthand when discussing these platforms.

Figure 9 .
Figure9.Time series of monthly mean vicarious VIIRS gain correction factors derived for each band (c.f.Figure5).Red lines indicate the mean value for a month and the shaded grey area the standard deviation within a month.
et al. (2012a), and is described briefly here.AERONET provides point observations of spectral columnar AOD with a repeat Atmos.Meas.Tech.Discuss., doi:10.5194/amt-2016-238,2016   Manuscript under review for journal Atmos.Meas.Tech.Published: 17 August 2016 c Author(s) 2016.CC-BY 3.0 License.frequency of approximately once per 10-15 minutes (in cloud-free conditions) and an uncertainty of order 0.01-0.02,with the larger uncertainties at shorter wavelengths(Eck et al., 1999).To mitigate the effects of spatiotemporal aerosol variability on the comparison, AERONET Version 2 Level 2 (cloud-screened and quality-assured;Smirnov et al., 2000) observations are averaged over ±30 minutes around the time of the VIIRS overpass, and interpolated to the VIIRS M-band wavelengths using the spectrally-closest AERONET AOD (with the exception that the AERONET 870 nm band is used in preference to 1020 nm due to increased uncertainties in the latter), and the Ångström exponent over the appropriate spectral region.Occasional missing AERONET AOD data at an individual wavelength are gap-filled in the same way.As it lies outside the AERONET spectral range, M11 (2250 nm) AOD is estimated from the AERONET AOD at 1640 nm and an Ångström exponent calculated over the spectral range 870-1640 nm.This spectral interpolation/extrapolation introduces negligible additional uncertainty into the AERONET values.These AERONET AOD data are compared with averaged SOAR-VIIRS retrievals passing algorithm quality assurance checks(Sayer et al., 2012a) whose pixel central locations lie within ±25 km of the AERONET site.

Figure 10 .
Figure 10.Statistics of SOAR-VIIRS 550 nm AOD validation against AERONET (a) without and (b) with the vicarious gains derived in the present study.Statistics are given on the panels: Pearson's linear correlation coefficient R, the median VIIRS-AERONET AOD bias, the root mean square (RMS) error, n the number of matched points, and f the fraction of points in agreement with AERONET within ±(0.03+10 %) Statistics are composites for all sites listed in Table4. .

Figure 11 .
Figure 11.Statistics of SOAR-VIIRS spectral AOD validation against AERONET without (black) and with (red) the vicarious gains derived in the present study.Panels show (a) correlation coefficients, (b) root mean square AOD error, and (c) median AOD bias.Statistics are composites for all sites listed in Table4.

Table 1 .
VIIRS moderate-resolution (M) bands, and band centres of similar MODIS bands used in this study.The final column indicates the shorthand notation adopted for each respective band pair in this study.

Table 2 .
Summary table showing the node points used in the MODIS/VIIRS intercalibration LUT.

Table 3 .
Derived vicarious gains.Mission-averaged gains are presented, as well as linear trends, for those bands where the estimated gain change over the mission to date is larger than 1 %, and statistically distinguishable from zero at the 90 % confidence level.For the trends, times t are defined in terms of years since January 1, 2010.Figures in parentheses indicate the one standard deviation uncertainty estimate.

Table 4 .
Locations of the AERONET sites used, and number of matchups obtained at each.