Vicarious Calibration of Beijing-1 Multispectral Imagers

Abstract

For on-orbit calibration of the Beijing-1 multispectral imagers (Beijing-1/MS), a field calibration campaign was performed at the Dunhuang calibration site during September and October of 2008. Based on the in situ data and images from Beijing-1 and Terra/Moderate Resolution Imaging Spectroradiometer (MODIS), three vicarious calibration methods (i.e., reflectance-based, irradiance-based, and cross-calibration) were used to calculate the top-of-atmosphere (TOA) radiance of Beijing-1. An analysis was then performed to determine or identify systematic and accidental errors, and the overall uncertainty was assessed for each individual method. The findings show that the reflectance-based method has an uncertainty of more than 10% if the aerosol optical depth (AOD) exceeds 0.2. The cross-calibration method is able to reach an error level within 7% if the images are selected carefully. The final calibration coefficients were derived from the irradiance-based data for 6 September 2008, with an uncertainty estimated to be less than 5%.

1. Introduction

The Beijing-1 micro-satellite was launched on 27 October 2005 and has a life expectancy of more than five years. It is China’s first micro-satellite for applied earth observation and forms part of the Disaster Monitoring Constellation. Beijing-1 weighs 166 kg and carries two imagers: One scans high-resolution (4 m) panchromatic images with a 24 km swath and the other provides medium-resolution (32 m) multispectral images with an ultra-wide 600 km imaging swath. The multispectral imagers comprise six cameras that are divided into two groups: bank0 and bank1 (Figure 1). Each of the two banks is equipped with three cameras of the same type, which acquire images in the red, green, and near-infrared bands, respectively. The fields of view of the two camera banks overlap by approximately 20 km (Figure 1). Given the off-pointing (±30° from nadir) imaging capability and large imaging swath, Beijing-1 can revisit the same location within 2 or 3 days [1,2]. Supported by a variety of real-time and stored-data operating modes and large storage capacity, Beijing-1 is a highly flexible and wide-ranging mission. Since its launch, Beijing-1 has acquired considerable amounts of imaging data, which have been widely used in agriculture, land use, and other applications [3].

Limited by its size and weight, Beijing-1 has no on-board calibration systems. Therefore, the radiometric performance of Beijing-1 is evaluated by using an approach of vicarious calibration. In the late 1980s, the Remote Sensing Group at the University of Arizona developed three vicarious techniques of absolute calibration that rely on in situ measurements. These methods are referred to as the reflectance-, irradiance-, and radiance-based techniques [4,5] and have been used successfully, e.g., for the SPOT HRV [6], Landsat TM/ETM [5,7,8], Daedalus scanner [9], Airborne Visible and Infrared Spectrometer [10], EO-1 Hyperion [11] FY [12], CBERS [13], and HJ-1 [14]. A variation of the radiance-based approach is cross-calibration, in which the calibration of one satellite sensor is transferred to other uncalibrated systems. This method is especially useful for sensors with large footprints and is more difficult to apply to the reflectance- and irradiance-based approaches [15]. Many sensors, such as Landsat TM/ETM [16,17], EO-1 Hyperion [18], Beijing-1/MS [19], SeaWIFS [20], and CBERS [21] have been calibrated via the cross-calibration method. Due to its high accuracy of calibration and excellent radiometric performance [22], the Moderate Resolution Imaging Spectroradiometer (MODIS) is often used as the reference to calibrate other sensors [20,21,23].

To achieve in-orbit calibration of the Beijing-1/MS, a field calibration campaign was performed at the China Dunhuang Calibration Site during September and October of 2008. Using the in situ calibration data, Beijing-1/MS and Terra/MODIS images, three vicarious calibration methods (i.e., reflectance-based, irradiance-based, and cross-calibration) were used independently to calculate the top-of-atmosphere (TOA) radiance (L_TOA). The results of the vicarious methods were then applied to an analysis to obtain the final calibration coefficients of the Beijing-1/MS.

2. Calibration Site and Data Sets

2.1. Calibration Site

The Dunhuang calibration site (40.2°N, 94.3°E) is located in the Gobi Desert in northwestern China, about 35 km west of the city of Dunhuang, Gansu Province. The calibration area is approximately 1,160 m above sea level. The entire target area (30 km × 30 km) for vicarious calibration is situated on a stabilized alluvial fan (see Figure 2). The area used for the vicarious calibration measurements of the high- or medium-spatial resolution sensors is approximately 400 m × 400 m (shown as area A in Figure 2) located in the center of the alluvial fan with the surface covered by cemented gravels. The local atmosphere is dry and typically has low levels of aerosol loading, which is beneficial for the calibration experiments. The atmospheric aerosol characteristics at the site are typical of a rural continental location, although some larger particles have been observed, possibly originating from sand dunes located to the northwest [24–26].

2.2. Beijing-1/MS and Terra/MODIS Images

Beijing-1/MS and Terra/MODIS images of the Dunhuang calibration site were acquired on 3, 6, 8, 13 September and 22 October 2008 (shown in

Table 1. Acquisition information of Beijing-1/MS and Terra/ Moderate Resolution Imaging Spectroradiometer (MODIS) images of Dunhuang.

	Beijing-1/MS					Terra/MODIS
Date	Imaging Cameras	Pass time (UTC)	Solar Zenith(°)	Viewing Zenith(°)	Relative Azimuth(°)	Pass time (UTC)	Solar Zenith(°)	Viewing Zenith(°)	Relative Azimuth(°)
	3 September	Bank0 & 1	3:35	43.2	12.3	41.5	4:02	39.6	55.04	43.7
6 September	Bank0	3:57	41.0	27.6	121.0	4:33	37.2	15.41	49.9
8 September	Bank0 & 1	3:39	43.9	4.6	55.6	4:17	31.6	34.78	37.6
13 September	Bank0 & 1	3:43	44.7	3.6	124.7	4:37	39.4	6.8	55.3
22 October	Bank0 & 1	3:36	57.3	5.85	69.2	4:42	52.3	2.22	114.7

). Figure 2 is the image of the Dunhuang site obtained on 13 September 2008. This image is a subset of the image mosaic generated from the imagery of bank0 and bank1. It can be seen from Figure 2 that there is a small difference between the images of bank0 and bank1. In addition, two dark images were acquired over the Pacific Ocean on the nights of 18 September in order to reduce the dark current values before vicarious calibration for Beijing-1/MS.

Terra/MODIS images were acquired on the same days as the Beijing-1/MS images and were used as references for cross-calibration. The relevant information from the Terra/MODIS images is shown in Table 1. Channels 4, 1, and 2 of Terra/MODIS were selected for reference with the Beijing-1 multispectral green, red, and near-infrared (NIR) channels. The main characteristics of the Beijing-1/MS and Terra/MODIS bands are listed in

Table 2. Main characteristics of Beijing-1/MS and Terra/MODIS bands.

Beijing-1/MS Bands	Spectral Range (nm)	Center Wavelength (nm)	Resolution (m)	MODIS Bands	Spectral Range (nm)	Center Wavelength (nm)	Resolution (m)
Green	523–605	564	32	4	545–565	554	500
Red	629–690	659	32	1	620–670	646	250
Near infra-red (NIR)	774–900	837	32	2	841–876	856	250

. From Table 1, it is obvious that the observation geometries of Beiing-1/MS and Terra/MODIS on 13 September and 22 October are similar to those collected on other days. This means that the data sets acquired on 13 September and 22 October are more suitable for cross-calibration than those collected on the other days.

2.3. Surface Reflectance

In situ measurements of the surface reflectance were taken from three areas, marked as A, B and C in Figure 2. The measurements in Area A were taken in the 2–3-hour period before and after the Beijing-1 overpass on 3, 6, and 13 September. The GPS trace of the measurements of the surface spectral reflectance on 13 September is shown in the lower right corner of Figure 2. The bi-directional reflectance of the Dunhuang site was measured in Area B on 8 September. The measurements in Area C were taken on 18 October 2008.

All of the surface reflectance measurements were carefully processed. Figure 3a shows the average reflectance of Areas A and C. The three average reflectance curves of Area A are uniform across the spectrum and show very little variation. The reflectance values measured on 3 and 6 September are very similar, but are slightly greater than on 13 September. The average relative difference between 13 September and those of 3 September or 6 September is about 5%. The decrease in the average spectral reflectance on 13 September may be caused by surface humidity at the calibration site, because light rain occurred during the night of 7 September, which caused surface humidity to increase and the reflection ratio to decrease on 13 September. The reflectance of Area C is similar in shape, but significantly larger in value than that of Area A. The differences between the three kinds of targets benefit the results of the calibration based on the data set.

Measurements of bi-directional reflectance in Area B were acquired from 4:08 UTC to 5:19 UTC on 8 September, covering the viewing zenith angles from 0° to 70° at steps of 14°, and relative azimuth angles from 0° to 150° at steps of 30°. The bi-directional reflectance of the Dunhuang site was convolved using the spectral response function of Beijing-1/MS or Terra/MODIS to produce the band-weight reflectance. The ratios were calculated for the band-weight reflectance at 0° zenith angle relative to those of other zenith angles. The relative bi-directional reflectance factor (BRF, which does not account for the reflectance of the reference panel) was calculated and the results for Beijing-1/MS are shown in Figure 3b. The results suggest a generally increasing trend towards the backward direction in all the bands, which is a unique characteristic of large particles over the Gobi.

To alleviate the bi-directional effects on radiometric calibration, the relative BRF was applied to convert the nadir reflectance of the calibration site into the viewing geometry of Beijing-1/MS and MODIS, as in Equation (1).

$$\rho \left(\lambda ,{\theta }_v,\phi \right)=\rho (\lambda ,0,0)\times r\left({\theta }_v,\phi \right)$$

(1)

Here, θ_v and φ view the zenith and relative azimuth angles, respectively; r is the relative BRF; λ is the wavelength of the sensor.

2.4. Atmospheric Data

Atmospheric data were collected at the same time as the surface reflectance measurements. A CE-318 sunphotometer was used to measure the intensity of the sun over the Dunhuang calibration site before and after the Beijing-1 overpass. These measurements were used to derive the total and estimated instantaneous atmospheric optical depths and total columnar water vapor. The aerosol optical depth (AOD) in each channel was calculated using the Beer’s attenuation law [27], and the influence of the spectral response function was taken into account to improve the retrieved accuracy [28]. Then, the aerosol optical depth within the 550 nm channel can be calculated via logarithmic interpolation of those of the 440 nm and 675 nm channels (see Figure 4(Left)).

In addition to the CE-318 data, the vertical profiles of atmospheric pressure, temperature, humidity, and ozone were also measured concurrently with the satellite overpass on 3, 6, and 13 September. The data were processed using interpolation and extrapolation methods, and were input into the relative computation module according to the format of the Second Simulation of a Satellite Signal in the Solar Spectrum Vector (6SV) model [29]; from the ground to 100 km altitude, the atmospheric pressure, temperature, moisture content, and ozone content are recorded at 34 levels. The data show almost no variation in temperature and atmospheric profiles on 3, 6, and 13 September. A slight variation in the moisture profile, indicating lower moisture content on 6 September than on the other days, is in accordance with the moisture content observed throughout the entire layer of the atmosphere.

To acquire the diffuse-to-global irradiance ratio data, an OL754 illuminance meter was used to measure the total solar radiation and sky diffuse radiation on 3, 6, and 13 September. Every group of data was measured three times successively. First, the global solar irradiance L₁ was measured. Subsequently, the sky diffuse radiation L₂ was assessed with a light barrier. Finally, the global solar irradiance L₃ was determined. The calculation of the diffuse-to-global irradiance ratio was divided into groups and executed by Equation (2). Figure 4(Right) shows the diffuse-to-global irradiance ratios of the Dunhuang experiments in September 2008. The band-weight diffuse-to-global irradiance ratio of each channel of Beijing-1 can be calculated through convolution of the response function of the Beijing-1/MS and the calculated diffuse-to-global irradiance ratio. Correspondingly, the diffuse-to-global irradiance ratio of the sun observation direction α_s and that of the sensor observation direction α_v can be calculated using interpolation or extrapolation approaches based on measurements over viewing angles

$$a=2L_2/\left(L_1+L_3\right)$$

(2)

3. Methodology

As mentioned in the introduction section, there are many factors particular to the Beijing-1/MS that need to be considered in the calibration. The first is that Beijing-1 has no on-board calibration systems. This means that the vicarious method is the only way to calibrate Beijing-1/MS. The second is that each channel of the multispectral imagers comprises two independent groups of similar cameras, and the integration time of the cameras changes automatically in response to the incoming scene radiation intensity. The third is that Beijing-1/MS has an ultra-wide 600 km imaging swath and a short revisit period. This capability provides the opportunity to calibrate Beijing-1/MS several times via independent methods over a short period of time.

Considering these special factors, an analysis is performed, and the overall workflow is illustrated in Figure 5. Based on images from Beijing-1/MS and MODIS, and the in situ data, three vicarious calibration methods (reflectance-based, irradiance-based, and cross-calibration) were used to calculate the L_TOA of Beijing-1/MS. The L_TOA data sets were then comprehensively analyzed and suspect results were identified and removed. The final Beijing-1/MS calibration was determined by comparing the results of the three independent calibration methods.

3.1. Normalization of Beijing-1/MS Images

(a) Dark current subtraction from Beijing-1/MS level 0 images

The dark sea images were averaged in the row direction to generate the dark current of each CCD. The dark current was found to be independent of integration time [30]. The dark current was firstly applied to all the Beijing-1/MS data from 3 September to 22 October before performing radiometric calibration.

(b) Normalization of Integration Time

The integration time of the Beijing-1/MS was designed to change automatically in response to incoming scene radiation intensity in order to obtain a better image. Laboratory tests prior to launch showed that the DN scales up with integration time in a linear manner [30] described by Equation (3):

$$L=\left({\mathit{DN}}_0-{\mathit{DN}}_{\mathit{dark}}\right)\times A\times \left(I_s/I_t\right)$$

(3)

where L is the TOA spectral radiance (in W/m²/sr/μm); I_s is the standard integration time; I_t is the target image integration time; A is the gain value of the calibration coefficient; DN_dark is the dark current; and DN₀ is the original output digital number from the level 0 product of Beijing-1/MS. Note that the dark current DN_dark was read from the dark images acquired on 18 September 2008.

According to Equation (3), taking the imaging integration time on 3 September as the standard time, the integration time normalization coefficients of the calibration data I_s/I_t were calculated (listed in

Table 3. Normalization coefficients of integration time of Dunhuang images.

Date	Integration Time (μs)	Is/It
3 September	650	1
6 September	643	1.010886
8 September	658	0.987842
13 September	668	0.973054
22 October	880	0.738636

) and applied to all the Beijing-1/MS data before performing radiometric calibration.

(c) Normalization of DN Values of Bank0 and Bank1

Each channel of the Beijing-1/MS overlaps by approximately 20 km. Because of the very similar design, the same atmospheric path, and the same incident field of the two banks of cameras, the disparity in DN values within the overlapping zones represents the systematic difference between the two banks of cameras. Therefore, the systematic difference can be eliminated by normalizing the DN values in the overlapping zone of the two banks of cameras. Linear regression (Equation (4)) of DN values in the overlapping zone was used to calculate the relationship between the bank0 and bank1 cameras. The linear regression coefficients (listed in

Table 4. Linear correlation coefficients of bank0 and bank1.

Channel	a	b	Goodness of Fit/R2
G	0.9686	−1.8872	0.9972
R	0.9586	−1.1938	0.9968
NIR	0.9403	−1.4212	0.9969

) were then used to transfer the Beijing-1/MS data from bank1 to bank0. After this transformation, only the data of bank0 was calibrated. The calibration coefficients of bank1 can be derived from the results of bank0 by the inversion of Equation (4).

$$D{N}_{\mathit{bank}0}=a\times D{N}_{\mathit{bank}1}+b$$

(4)

here, DN_bank0 and DN_bank1 are the respective DN values of images from the cameras of bank0 and bank1, a is the proportionality coefficient of DN_bank0 and DN_bank1, and b is the constant of this expression, listed in Table 4.

3.2. Reflectance-Based and Irradiance-Based Methods

The reflectance-based method requires accurate measurement of the spectral reflectance of the ground target and measurement of spectral extinction depths and other meteorological parameters. Scattering and absorption in the atmosphere are computed using approximate radiative transfer models and codes, such as 6S, 6SV [29,31]. The code output is a TOA radiance value for a given ground reflectance. This radiance is compared with the average DN from the image of the ground area, to give a calibration coefficient in units of counts per unit radiance. The irradiance-based method is also called the improved reflectance-based method. It uses all the measured data from the reflectance-based method together with measurements of the ratio of diffuse-to-global spectral irradiance at ground level. This additional measurement helps reduce the uncertainties in the aerosol model used for the scattering computations [32]. The principle used in the reflectance- and irradiance-based methods to calculate the TOA spectral reflectance is shown by Equation (5) [31] and Equation (6) [12], respectively. Both methods use Equation (7) to transform the TOA spectral reflectance into the TOA radiance.

$$\rho *\left({\theta }_s,{\theta }_v,{\varphi }_{v-s}\right)={\rho }_a\left({\theta }_s,{\theta }_v,{\varphi }_{v-s}\right)+\frac{{\rho }_t}{1-s\cdot {\rho }_t}\cdot T\left({\theta }_s\right)\cdot T\left({\theta }_v\right)$$

(5)

$$\rho *\left({\theta }_s,{\theta }_v,{\varphi }_{v-s}\right)=T_g\left[{\rho }_a\left({\theta }_s,{\theta }_v,{\varphi }_{v-s}\right)+\frac{e^{-\tau /{\mu }_s}}{1-a_s}\cdot \rho \left(1-\rho s\right)\cdot \frac{e^{-\tau /{\mu }_v}}{1-a_v}\right]$$

(6)

$$L=\rho *\cdot \text{cos}\left({\theta }_s\right)\cdot E_0/\left(d^2\cdot \pi \right)$$

(7)

In Equations (5)–(7), θ_s is the sun zenith angle, θ_v is the view zenith angle of the sensor, and ϕ_v−s is the relative azimuth angle between the view azimuth angle and the sun azimuth angle. ρ_t is the measured spectral reflectance of the ground target, and ρ_a is the reflectance that corresponds to the atmosphere path radiance (atmosphere intrinsic reflectance). S is the atmospheric hemisphere reflectance. T(θ_s) and T(θ_v) are the total transmittance of the solar path and the view path, respectively. ρ^* and L are the TOA spectral reflectance and the TOA radiance value of the ground target, respectively. μ_s and μ_v are the values of cosθ_s and cosθ_v, respectively, and a_s and a_v are the diffuse-to-global ratio of the sun direction and the view direction, respectively. E₀ is the TOA value of solar irradiance, d² is the Sun-Earth distance.

Based on the in situ reflectance measurements and the atmospheric parameters (see

Table 5. The atmospheric parameters used for reflectance- and irradiance-based methods.

Date	3 September	6 September	13 September
AOD at 550 nm	0.32	0.20	0.34
Vertical Column of Water Content (g/cm2)	0.74	0.64	0.76
Aerosol Type	Shettle model for background desert aerosol
Atmospheric Model	Mid-latitude summer

), the band-weight TOA radiance of Beijing-1/MS was calculated by the radiative transfer code 6SV [29] according to the reflectance-based and irradiance-based methods. The final results are listed in

Table 6. Calibration coefficients (W/m²/sr/μm/DN) derived from three vicarious calibration methods based on Dunhuang calibration data sets.

Method	Time	Green	Red	NIR
Reflectance-based	3 September	0.8513	1.0148	0.8693
	6 September	0.7959	0.9646	0.8363
	13 September	0.8385	0.9991	0.8570
Irradiance-based	3 September	0.6371	0.7809	0.7161
	6 September	0.6794	0.8401	0.7655
	13 September	0.6504	0.8062	0.7652
Cross-calibration	3 September–22 October	0.6560	0.7783	0.7427

.

3.3. Cross-Calibration Method

Cross-calibration is one of various methods used for post-launch satellite sensor calibration. Here, a given sensor is calibrated against another satellite sensor for which the radiometric calibration is better known, via near-simultaneous imaging of a common ground target [33]. By the cross-calibration method, the band-weight TOA radiance of Beijing-1/MS (L_BJ) is obtained from the band-averaged TOA reflectance of MODIS (ρ_MODIS), according to Equation (8) and (9). The TOA reflectance of MODIS is taken rather than the TOA radiance of MODIS in order to avoid the error caused by different versions of the E₀. When transforming the TOA reflectance into the TOA radiance of MODIS, the E₀ provided by the data is used (see Equation (7)). MODIS uses a combination of E₀ from different resources [34], whereas we use the E₀ in 6SV from Neckel et al. [35].

K_BJ–MODIS is the spectral matching factor, which is the ratio of the TOA reflectance of the uncalibrated spectral band to that of the reference spectral band (see Equation (9)), which is defined to adjust the spectral band difference between the Beijing-1/MS and MODIS. The TOA reflectance of Beijing-1/MS and MODIS were calculated using the 6SV model based on the in situ measurements of the surface and atmospheric data, respectively, and then K_BJ–MODIS was calculated according to Equation (9). Specially, the spectral matching factors for the images acquired on 8 September and 22 October were calculated based on the in situ data sets obtained on 3 September.

$$L_{BJ}=K_{BJ-\mathit{MODIS}}\times \frac{E_{BJ}\times \text{cos}{\theta }_{BJ}}{\pi \times d^2}\times {\rho }_{\mathit{MODIS}}$$

(8)

$$K_{BJ-\mathit{MODIS}}=\frac{{\rho }_{BJ}}{{\rho }_{\mathit{MODIS}}}$$

(9)

Here, L_BJ is the band-weighted TOA radiance of Beijing-1 (in W/(m²μm), ρ_BJ and ρ_MODIS are the band-weighted TOA reflectances of Beijing-1 and MODIS, respectively. E_BJ is the band-averaged exoatmospheric solar irradiance of Beijing-1/MS (in W/(m² μm)); θ is the solar zenith angle, and d² is the Earth–Sun distance in Astronomical Units.

The selected pairs of images for Beijing-1/MS and MODIS were first geometrically corrected to match each other. Then, the small subsets in the uniform portion within the areas of A and C were carefully selected and the mean values were used in cross-calibration. The TOA radiance of Beijing-1/MS was calculated by Equation (8) from the TOA reflectance of MODIS with the corresponding spectral matching factor. The bi-directional effect of the Dunhuang calibration site in Area A accounted for different solar-viewing geometries in Beijing-1/MS and MODIS images acquired on 3, 6, 8, 13 September and 22 October. Pixels in Area C were only used for images acquired on 13 September and 22 October, due to the near-nadir viewing direction and the small difference in the solar-viewing geometries of Beijing-1/MS and MODIS. The final cross-calibration results were derived by the linear regression of the all the five pairs of images for Area A acquired on 3, 6, 8, 13 September and 22 October and two pairs of images for Area C acquired on 13 September and 22 October. The results are shown in Figure 6 and listed in Table 6.

4. Results

The L_TOA calibrated via the three vicarious methods and the corresponding DNs were scattered in the green, red, and NIR channels, respectively (Figure 6). The results for radiometric calibration by the three methods are listed in Table 6. All the calibration coefficients were obtained under the standard integration time (Is = 650 μs) and only apply to bank0 cameras. The coefficients for bank1 cameras could be easily transferred from those for bank0 by Equation (4).

It is interesting to note that: (1) The results of the three vicarious methods are separate whereas the results of one method are gathered; (2) the coefficients derived by the reflectance-based method are obviously larger than those derived from the other two methods; (3) for the irradiance-based method, the coefficients derived on 6 September are larger than those on the other two days; (4) for the reflectance-based and irradiance-based methods, the results on 6 September are more similar than those obtained on the other two days; (5) the linear trends derived from cross-calibration points are sufficiently consistent to provide confidence in the cross-calibration results. The results also indicate that Beijing-1/MS presents a stable radiometric performance within at least two months during a period of three years since the satellite was launched.

During the short duration of the Dunhuang calibration campaigns, it is reasonable to assume that the radiometric performance of Beijing-1/MS is stable and that the results obtained from different methods should be approximately equal. In practice, it is inevitable that various errors in the observed data, data processing, and those inherent in the selected calibration method may decrease the accuracy of the calibration coefficients. One important principle for calibration is that different and independent techniques should be used in order to identify, remove, or account for any errors, especially those, which are systematic in the calibration results [32].

Generally, the calibrated gains of the reflectance-based method are about 9% to 30% greater than those of the irradiance-based and cross-calibration method. This appears to be a systematic error in the reflectance-based method, which may be a result of errors in calculating the aerosol properties. Theoretically, the irradiance-based method is an improved reflectance-based method, which enables greater accuracy [4,36]. It uses the ratio of diffuse-to-global spectral irradiance to minimize the large uncertainty associated with assumptions regarding the aerosol type. Therefore, the results of the reflectance-based method are adjusted to account for the large uncertainty issue, which is inherent to the assumption of aerosol type.

With regard to the results of the irradiance-based method, the gains of 6 September are much greater (the relative error is about 4%–8%) than those of 3 and 13 September. After examining the measurements of aerosol optical depths and ratios of diffuse-to-global irradiance for the three days (see Figures 4), it was determined that: (1) Aerosol optical depth on 6 Septemberis the smallest; (2) the atmosphere is most stable because aerosol optical depth is nearly invariant during the period 08:30 to 13:30; and (3) the measurement of the ratios of diffuse-to-global irradiance on 6 September encompasses the largest range of zenith angle (i.e., from 39° to 73°). Conversely, on the other two days, the atmosphere is less stable, and the measurement of the ratios of diffuse-to-global irradiance encompasses a reduced zenith angle range (i.e., from 33° to 46°, and from 37° to 50° for 3 September and 6 September, respectively). The collective effect of a narrow zenith angle range and unstable atmospheric conditions would introduce large uncertainties when extrapolating the observed ratios of diffuse-to-global irradiance to the values near the satellite viewing direction (e.g., near the nadir observation). Therefore, the results of the irradiance-based method are most accurate on 6 September, and the results of the other two days are not considered in the last results.

In comparison, the gains derived from irradiance-based and cross-calibration methods are more similar to each other, with relative differences of approximately 3% to 8%. This coincidence observed in two different individual methods provides a certain level of confidence for the calibration gains of Beijing-1/MS. The estimation of uncertainty detailed in the following section analysis will further bolster such a conclusion. In addition, taking the uncertainties of the irradiance-based and cross-calibration methods into account, the results derived from using the irradiance-based method are more accurate and reliable than those derived by the cross-calibration method.

In brief, an analysis of the results derived by the three different and independent methods provides the opportunity to determine whether systematic errors exist in one or more of the techniques. Based on the analysis, it is concluded that the results derived on 6 September by the irradiance-based method are chosen as the final calibration coefficients of Beijing-1/MS.

5. Uncertainties Estimation

Different error sources were accounted for in the various methods, and they are reflected in the differences between the final calibration coefficients. The uncertainty in the reflectance-based method may derive from the measurements of ground reflectance, aerosol optical depth, assumption of aerosol type, assumption of atmospheric model, and non-Lambertian reflectance of ground and geometric imaging angles. To minimize the large uncertainty due to the assumptions of aerosol type, the irradiance-based method adds the measurement of the ratios of diffuse-to-global irradiance, which was reported to improve accuracy [4,36]. However, the measurement of the ratios of diffuse-to-global irradiance may also introduce a large uncertainty when the measurement does not cover a wide change of zenith under stable atmospheric conditions. The cross-calibration method has the advantage of independence of synchronous in situ measurements for satellite overpass and can achieve moderate accuracy. The uncertainty in this method may derive from radiometric calibration of the referenced satellite, spectral match functions, and BRDF of the ground cover.

5.1. Uncertainty Estimation for Reflectance- and Irradiance-Based Methods

The calibration uncertainties in the reflectance- and irradiance-based methods have been well discussed by Biggar and other researchers as being 4.9% and 3.5%, respectively [4,36]. However, the contribution due to the assumption of aerosol type (including complex index and size distribution) was estimated under a very low AOD which is less than 0.1 in 550 nm; therefore, that method may not correct for the AOD values as high as 0.2 or 0.3 in our work. During the vicarious calibration campaigns, the assumption of correct aerosol type was of great importance. In practice, it is difficult to accurately determine aerosol properties in field experiments. The behavioral properties of the actual aerosols were often very different from those of the standard aerosol types in the radiative transfer models. These differences may result in systematic uncertainties in the calibration results. Thus, two other aerosol types (i.e., continental model and urban model) were chosen to replace the desert aerosol, and the 6SV model was run again in order to approximate their contributions to the total calibration uncertainty (Figure 7).

It is evident that the reflectance-based method is more sensitive to aerosol type than the other methods. The relative errors due to different aerosol types were more than 10% for this method. The results derived from the continental aerosol type and those derived from the urban aerosol type approached each other. It is interesting that the results obtained using different aerosol types were much more similar on 6 September than on the other two days. This can be attributed to a reduced AOD on 6 September, meaning that the uncertainty due to aerosol properties was not as evident as in the cases with larger aerosol optical depths. In addition, it can be seen that the difference was less within the longer wavelengths than within the shorter wavelengths. This also explains why the results derived via the reflectance-based method and by the irradiance-based method were more similar to each other on 6 September than on the other two days. It was also found that the results derived from the reflectance-based method with the urban aerosol type approached those of the irradiance-based method. However, this may be coincidental, because experience suggests that the urban aerosol type is not suitable for modeling the Gobi environment of Dunhuang. However, the aerosol type of the Dunhuang calibration site should be verified through long-term ground observation.

In contrast, the radiance predicted using the irradiance-based method showed little change under different aerosol types. The relative uncertainties due to different aerosol type are 2%–3%, 1%–2%, and 3%–4% from NIR band to green band for 3 September, 6 September, and 13 September, respectively. Note that the uncertainty increases with AOD. Therefore, the irradiance-based method is believed to provide greater accuracy than the reflectance-based method for cases in which the ratio of diffuse-to-global irradiance has great accuracy. However, the measurement of the ratio of diffuse-to-global irradiance is an additional source of error in the irradiance-based method, and it is unavailable in the mid-latitude regions in the satellite viewing direction because the zenith angle is so small; therefore, extrapolation is necessary to achieve an accurate value in such cases. Furthermore, the measurement should cover a relatively large zenith range, and stable atmospheric conditions are another important requirement. We compared the observations of the ratios of diffuse-to-global irradiance and aerosol optical depths for the three days and found that 6 September had more stable atmospheric conditions, a reduced AOD, and a wider range of observation angles for the measurement of the ratio of diffuse-to-global irradiance. Therefore, the results derived from the irradiance-based method on that day were certain to have a higher accuracy than on the other two days. Referring to the error estimation by Biggar [4,36], we estimated that the contributions of errors of diffuse-to-global ratio measurements to the final calibration uncertainty were 4%, 3%, and 4% referring to 3 September, 6 September, and 13 September, respectively.

Besides the errors associated with aerosol type and measurements of diffuse-to-global ratio, the error due to incomplete measurement of BRDF at the Dunhuang site is estimated to be 2% [26]. This kind of error is present in both the reflectance- and reflectance-based methods. In addition, the accuracy of 6SV is much improved and the relative uncertainty is estimated as 0.4–0.6%. We briefly summarized the estimation of uncertainties discussed above in

Table 7. Estimated uncertainty in the reflectance-based method.

Source	Overall Uncertainty (%)
	3 September	6 September	13 September
Ground reflectance measurement		2.1
Optical depth measurement		1.1
Absorption computation		1.3
Assumption of aerosol type	13.0–16.4	9.9–12.6	15.6–19.4
Vertical distribution		1.0
BRDF error		2.0
Inherent code accuracy		0.6
Uncertainty in the value of μs		0.2
Total uncertainty (root sum of squares)	13.5–16.8	10.5–13.1	16.0–19.7

and

Table 8. Estimated uncertainty in the irradiance-based method.

Source	Overall Uncertainty (%)
	3 September	6 September	13 September
Optical depth measurement		1.1
Diffuse-to-global irradiance measurement	4.0	3.0	4.0
Ground reflectance measurement		2.1
BRDF error		2.0
Assumption of aerosol type	1.6–2.9	1.4–1.9	2.6–4.3
Inherent code accuracy		0.6
Uncertainty in the value of μs and μv		0.1
Total uncertainty (root sum of squares)	5.3–5.8	4.6–4.7	5.7–6.7

, with some uncertainties derived from Biggar [4,36].

5.2. Uncertainty Estimation for Cross-Calibration Method

The accuracy of the cross-calibration method may be reduced by the spectral mismatch of the sensor, BRDF of the ground surface, registration of the data sets, and calibration accuracy of the referenced sensor [37–39]. Among these factors, the first two sources will ultimately result in uncertainty in the spectral matching factor K_BJ–MODIS of Equation (9). Herein, we focus on uncertainty in the spectral matching factor, especially accounting for BRDF error and Lambertian assumption. Such an estimation indicates the approximate magnitude of the overall uncertainty.

Relative errors of K_BJ–MODIS are calculated between the case of Lambertian assumption at the Dunhuang site and that of BRDF-correction on 13 September, ranging from 2.5% in the NIR band to 2.8% in the green band. Similar errors on 22 October may occur due to the similar solar-viewing geometry of the two sensors to that on 13 September. This means that the uncertainty of the non-Lambertian ground characteristics of Area C is no more than 3%. In addition, the values of K_BJ–MODIS for Dunhuang are computed by the other two aerosol types (i.e., continental and urban) from 3 September to 22 October. The relative errors due to aerosol type are approximately to 0.5%–4.2% for different bands and solar-viewing geometries. As the calibrated gains are derived from multiple points, the largest errors should account for each individual error source. The overall uncertainties are listed in

Table 9. Estimated uncertainty in the cross-calibration method.

Source	Overall uncertainty (%)
Uncertainty of MODIS calibration	3.0
Image registration error	1.0
Atmospheric stability	1.0
Assumption of aerosol type	4.2
BRDF error of Dunhuang site (Area A)	2.0
Non-Lambertian ground characteristics (Area C)	3.0
Inherent code accuracy	0.6
Total uncertainty (root sum of squares)	6.4

.

6. Conclusions

In this paper, we proposed a vicarious calibration procedure for Beijing-1/MS in order to improve the robustness of the calibration results. Three vicarious calibration methods (reflectance-based, irradiance-based, and cross-calibration) were used based on five calibration data sets acquired for the Dunhuang calibration site during September and October 2008. An analysis was conducted to identify the methodological and measurement errors in the calibration results. The final calibration gains of Beijing-1/MS were taken to be the results derived via the irradiance-based method on 6 September, within an accuracy of 5%. It is also found that the cross-calibration method provides moderate uncertainty in that less than 7% if the image pairs and cross-calibration sites are chosen with caution. This study clearly shows that the reflectance-based, irradiance-based, and cross-calibration vicarious calibration methods each inevitably introduce inherent uncertainties. The use of several different and independent vicarious calibration methods is essential to determine the systematic errors inherent in one or more methods, and to best identify, remove, or account for such errors in the calibration results.

Future work will focus on the long-term ground observation of atmospheric conditions at selected calibration sites such as Dunhuang, and on measurements of surface reflectance of typical calibration targets. Based on these data sets, the aerosol properties of the calibration site and the BRDF of the calibration targets will be analyzed, and the long-term radiometric performance of Beijing-1 will be assessed by cross-calibration methods.