2025 Volume 103 Issue 1 Pages 87-109
Himawari-8/9 is a next-generation Japanese geostationary earth orbit (GEO) meteorological satellite with an onboard sensor – the Advanced Himawari Imager (AHI). Because Himawari-8/9 AHI observe the Earth’s hemispheres every 10 min with multiple spectral bands, AHI providing an unprecedented opportunity to facilitate its observation datasets are expected to be a new data source for terrestrial monitoring in terms of mitigating cloud contaminations. Estimation of land surface reflectance (LSR) is crucial in quantitative terrestrial monitoring. In this study, we aimed to develop a method for estimating the LSR and angular-adjusted LSR of the Himawari-8/9 AHI using the look-up table (LUT) based Second Simulation of a Satellite Signal in the Solar Spectrum Vector (6SV) Radiative Transfer Model (RTM) and kernel-driven semi-empirical bidirectional reflectance distribution function (BRDF) model. The estimated LSR underwent evaluation and inter-comparison through two distinct approaches: ray-matching and estimating angular-adjusted LSR. Ray-matching of the obtained data pairs with the MODIS LSR product shows that the correlation coefficients (r) for all bands are greater than 0.86 at low latitudes. Angular-adjusted LSRs estimated using AHI time-series data at mid-latitudes also show good agreement with MODIS (r > 0.5), particularly the red and near-infrared bands (r > 0.9). The results obtained by our method are in high agreement with those calculated using the reference aerosol optical thickness (AOT) (r > 0.98). Our findings highlight the potential application of our methodology to other GEO satellites for high-frequency terrestrial monitoring at continental to global scales.
Land surface reflectance (LSR) is indispensable in terrestrial monitoring, and quantifies the fraction of solar radiation reflected off Earth’s surface, which is intrinsically linked to surface properties, as well as the geometry of illumination and observation (Lee et al. 2020a, 2022). Traditionally, sensors from low earth orbit (LEO) satellites, such as the Terra/Aqua and Suomi national polar-orbiting partnership (NPP), have been widely employed as the primary sources of LSR data (Liang et al. 2002; Vermote et al. 2014). However, these data are often limited by cloud cover and infrequent observations (Fensholt et al. 2011). Particularly in the tropical regions, a month of missing data may be observed (Nagai et al. 2014). In contrast, third-generation geostationary earth orbit (GEO) satellites, including GOES Advanced Baseline Imager (ABI) (Schmit et al. 2017), Himawari-8 Advanced Himawari Imager (AHI) (Bessho et al. 2016), FY-4 Advanced Geostationary Radiation Imager (AGRI) (Yang et al. 2018), GK2A Advanced Meteorological Imager (AMI) (Lee et al. 2017), and MTG-I Flexible Combined Imager (FCI) (Holmlund et al. 2021), have emerged as promising alternatives owing to their high temporal observation frequency and multiple solar reflective bands (Miura et al. 2019; Wang et al. 2020).
The AHI onboard Himawari-8/9, developed by the Japan Meteorological Agency (JMA), offers improvements in sensor capabilities and spatiotemporal resolution compared with its predecessor, the Multifunctional Transport Satellite (MTSAT)-2 Imager (Bessho et al. 2016). Although the primary design focus of AHI is weather observation and forecasting, it has proven beneficial for other applications such as disaster detection (Higuchi 2021; Miura and Nagai 2020), vegetation monitoring (Zhang et al. 2021; Yamamoto et al. 2023) and snow cover estimation (Wang et al. 2019).
Nonetheless, the environmental monitoring potential of its Himawari-8/9 AHI remains underutilized owing to the lack of rigorously evaluated and validated publicly accessible LSR datasets. Therefore, reliable methods for LSR retrievals are essential for meaningful AHI research. Radiative transfer models (RTM) are the most widely used method for estimating the LSR from GEO satellites because of their ability to simulate solar radiation transmission in the atmosphere (Vermote et al. 1997). Initially, a simplified method for atmospheric correction (SMAC) was employed to estimate GEO-based LSR for the Meteosat Second Generation (MSG) Spinning Enhanced Visible and InfraRed Imager (SEVIRI) (Proud et al. 2010). While effective in numerous situations, this method was occasionally limited in accurately estimating LSR under complex atmospheric conditions and at medium-to-high viewing angles. Peng (2020) designed the LSR and albedo estimation algorithms for GOES-R using a look-up table (LUT) based on the Second Simulation of the Satellite Signal in the Solar Spectrum Vector (6SV) RTM. Lee et al. (2020a) designed a GK-2A Land Surface Albedo estimation algorithm using a 6SV-based LUT method. Li et al. (2019) and Wang et al. (2020) describe that the NASA GeoNEX group provides AHI LSRs estimated by the Multi-Angle Implementation of Atmospheric Correction (MAIAC) algorithm.
Currently, there are gaps in knowledge and limitations in the evaluation of GEO-based LSR. The evaluation of GEO LSR products involves two methods: ray-matching with LEO satellites (Yu and Wu 2016; Qin and McVicar 2018) and inter-comparison using angular-adjusted LSR (Li et al. 2022). The applicable regions of the ray-matching method are limited to low-latitude regions due to the differences in the observation condition (sun-target-sensor geometry) of the GEO sensors and the LEO sensors (Qin and McVicar 2018), leaving a gap in mid-latitude evaluations. To overcome those gaps, certain studies have used Moderate Resolution Imaging Spectroradiometer (MODIS) Bidirectional Reflectance Distribution Function (BRDF) products to evaluate GEO LSR (Tran et al. 2020; Li et al. 2022). Meanwhile, Zhang et al. (2022) suggested the feasibility of BRDF estimated by time-series data from GEO satellites, which could become a potential method for evaluating GEO satellite products.
This study aims to refine the existing AHI LSR estimation method and to provide a reliable AHI LSR dataset. Moreover, we aim to release an open-source framework for estimating LSR of GEO satellites, promoting research on land monitoring by GEO satellites. Building on the work of these researchers, we refined the workflow for estimating Himawari/AHI LSR and its evaluation by applying methods from existing studies and improving data selection and processing. By utilizing a ray-matching and a kernel-driven semiempirical BRDF model to simulate LSR, we conducted an inter-comparison with MODIS (Terra/Aqua) LSR products, evaluating the estimated AHI LSR at low to mid-latitudes. In addition, we explored the feasibility of applying the methodologies and materials used in this study to other GEO satellites.
This study focuses on regions within the observation coverage of the Himawari-8/9 AHI, primarily encompassing East and Southeast Asia and Oceania (Fig. 1), according to the gridded dataset provided by the Center for Environmental Remote Sensing (CEReS) at Chiba University. This region experienced various terrestrial surface changes, such as land use changes due to deforestation and plantation in Southeast Asia (Vadrevu et al. 2019), and is extensive fire and subsequent ecological changes in Australia (Gibson et al. 2020; Abram et al. 2021).
Observation region for CEReS Himawari/AHI gridded data. The red circles in the figure represent the AHI view zenith angles, and the central star symbol represents the Himawari-8 AHI sub-satellite point. The background land cover data is from MCD12Q1. AHI, Advanced Himawari Imager; CEReS, Center for Environmental Remote Sensing.
The Himawari-8 satellite was launched on October 7, 2014, and its data services became operational on July 7, 2015. The AHI has 16 spectral bands, from visible to thermal infrared, to capture various meteorological phenomena and atmospheric components. Offering two distinct observation intervals, the AHI provides imagery every 10 min for the Full Disk (FD) region encompassing the entire hemisphere, and every 2.5 min for the Japanese region or the target area (Bessho et al. 2016). Himawari-9 was launched on November 2, 2016. Himawari-9 is a satellite of the same type with Himawari-8 and serves as a backup for Himawari-8. Observation began on December 13, 2022, replacing Himawari-8.
In this study, we used the CEReS gridded Himawari-8/9 AHI dataset, which includes additional precise geometric correction and reprojection to the latitude/longitude coordinates from Himawari standard data (HSD) distributed by the JMA (Takenaka et al. 2020). The accuracy of the geometric correction was evaluated, and the dataset can be utilized for studies requiring high geometric precision (Yamamoto et al. 2020), such as land monitoring (Zhao et al. 2022; Yamamoto et al. 2023). Furthermore, annually updated calibration coefficients were applied to the dataset in order to minimize the effects of sensor degradation (https://www.cr.chiba-u.jp/databases/GEO/H8_9/FD/index_en_V20190123.html). The specifications of the CEReS gridded Himawari-8/9 AHI dataset are presented in Table 1. In this study, Bands 01 to 06 of the AHI were used for the LSR estimation. The spatial resolution of Band 03 was reduced from 0.005° to 0.01° by taking a mean to align with that of Band 04.
We utilized cloud mask data based on an algorithm reported by Yamamoto et al. (2018), which sets thresholds for reflectance and brightness temperature in visible and infrared band data to identify clouds. The algorithm provides confidence of clear sky at 0.02° spatial resolution, and all pixels with confidence below 0.95 are identified as cloudy in this study. Because the cloud mask data are provided with 0.02° spatial resolution, we used the nearest neighbor (NN) method to convert the data to 0.01° spatial resolution to match the spatial resolution of solar reflective bands (0.01°).
Moreover, we generated angular data for the AHI, including the solar zenith angle (SZA), solar azimuth angle (SAA), view zenith angle (VZA), and view azimuth angle (VAA). SAA and SZA were calculated by the Variations Séculaires des Orbites Planétaire (VSOP) 87 theory (Bretagnon and Francou 1988), a mathematical and analytical theory developed to accurately predict the orbital positions of the planets in the solar system over time. VZA and VAA were calculated using the geometric relationship between the satellites and each grid. To save computational time, we generated angular data with spatial and temporal resolutions of 0.04° and 10 minutes, respectively, and used NN interpolation to match the spatial resolution to that of the AHI.
b. Auxiliary data for LSR estimationIn this study, the Copernicus Atmosphere Monitoring Service Reanalysis (CAMSRA) ECMWF Atmospheric Composition Reanalysis 4 (EAC4) atmospheric composition dataset was utilized to retrieve key parameters for atmospheric correction, total column ozone, total column water vapor, total aerosol optical thickness (AOT) at 550 nm, and aerosol components. This dataset encompasses aerosol and atmospheric chemistry information every 3 h at 0.75° spatial resolution and is derived from the assimilation of satellite inversion data using the integrated forecasting system of ECMWF (Inness et al. 2019; Koffi and Bergamaschi 2018). To ensure compatibility with the temporal and spatial resolutions of the AHI data (10 min, 0.01°, and 0.02°), linear interpolation was used for both temporal and spatial interpolation.
The multi-error-removed improved-terrain (MERIT) digital elevation model (DEM) (Yamazaki et al. 2017) was used for topographic data. This high-resolution representation of Earth’s surface is an invaluable resource for geoscience research (Yamazaki et al. 2017; Uuemaa et al. 2020). Similarly, we converted the dataset from a spatial resolution of 3” to match the AHI spatial resolutions (0.05°, 0.01°, and 0.02°) by averaging the grids within each AHI grid.
c. Terra and Aqua MODIS LSR productsWe used August 2018 Terra/Aqua MODIS daily LSR grid datasets [Terra: MOD09GA, Aqua: MYD09GA (Vermote et al. 2002)] as reference datasets to quantitatively evaluate the estimated AHI LSR. These datasets include Bands 01 through 07 of MODIS. By referring to the Spectral Response Functions (SRFs) of AHI and MODIS (Fig. 2), we used MODIS Bands 01–07 (except 05), which are close in central wavelength to the AHI bands (Bands 01–06). The accuracy of the latest Collection 6 LSR product is 0.005 + 0.05 × LSR or more (Vermote and Kotchenova 2008). We applied mean value resampling to match spatial resolution of MODIS to AHI, i.e. MODIS data with 500 m spatial resolution using 2 × 2 pixels corresponds to one AHI 1 km spatial resolution pixel and 4 × 4 corresponds to AHI 2 km spatial resolution pixels.
Spectral response functions of (solid line) AHI and (dashed line) MODIS, (green line) green vegetation spectral curve. AHI, Advanced Himawari Imager; MODIS, Moderate Resolution Imaging Spectroradiometer; SRF, Spectral Response Functions.
Additionally, we utilized the MODIS BRDF parameter product (MCD43A1) (Schaaf et al. 2002) to estimate the MODIS LSR at the AHI observation angle. MCD43A1 is a 16-day synthesized product with a spatial resolution of 500 m and includes BRDF parameters in 7 spectral bands.
d. Aerosol Robotic Network and SKYNETUncertainties can arise in the estimation of LSR owing to the accuracy of the input parameters related to the atmosphere and aerosols. To address this issue, our study evaluated interpolated data from the CAMSRAEAC4 with in-situ (23 sites) AOT from the Aerosol Robotic Network (AERONET) and SKYNET for the years 2018–2019.
AERONET, a global network of ground-based aerosol monitoring stations, provides long-term continuous datasets of aerosol optical and microphysical properties. These data have been used extensively in aerosol research, atmospheric correction of satellite data, and air quality monitoring (O’Neill et al. 2003). Similar to AERONET, SKYNET operates as a ground-based observation network (Aoki and Fujiyoshi 2003; Irie et al. 2017), thereby providing valuable data for validating satellite-based observations (Damiani et al. 2018; Hori et al. 2018). Given the extensive presence of AERONET stations in Europe and North America, we incorporated data from SKYNET to improve coverage in East Asia. AOT550 was obtained using Ångström exponent (AE) after calculating the AOT at 500 nm (Ångström 1929). To further mitigate cloud contamination, we used an AHI cloud mask (Yamamoto et al. 2023) for additional screening.
We further used two sites, Beijing (BEJ) and Miyako (MYK), to assess the uncertainties introduced by the input AOT using observations. These two sites were chosen for their contrasting aerosol characteristics. We estimated the local LSR at noon using the in-situ AOT and the AOT interpolated from the CAMSRA-EAC4 dataset in 2018.
2.3 LSR estimation and BRDF modelingThe estimation of the AHI LSR and BRDF kernel model parameters consisted of two primary steps (Fig. 3). The first step involved retrieving the AHI LSR by utilizing the 6SV RTM LUT in conjunction with auxiliary data. In the second step, the estimated LSR within the framework of the kernel driven BRDF model was used to estimate the BRDF kernel model parameters and subsequently compute the angular adjusted LSRs. The associated subsections in Fig. 3 offer detailed explanations of each step of the algorithm.
Flowchart of the study. AHI, Advanced Himawari Imager; BRDF, Bidirectional Reflectance Distribution Function; CEReS, Center for Environmental Remote Sensing; LSR, Land Surface Reflectance; LUT, Look-up Table; RTM, Radiative transfer models.
We used the vector version of 6S, commonly referred to as 6SV, to estimate the LSR (Vermote et al. 2006). The 6SV employs the radiative transfer theory and successive order of scattering method to accurately simulate the atmospheric effects from the sun to the target and sensor (Vermote et al. 1997, 2006; Kotchenova et al. 2006; Kotchenova and Vermote 2007). Its implementation has become widespread across various studies, including those focusing on the retrieval of AOT (Xie et al. 2022), the estimation of land surface albedo, and LSR (Kotchenova et al. 2008; Lee et al. 2020). We first assume that the land surface is Lambertian and conducted atmospheric corrections to the AHI data. 6SV provides three coefficients (Xap, Xb and Xc in Eq. 1) as outputs under the input conditions. For Lambertian surfaces and above the sea level, the expression for the LSR within the 6SV is expressed by Eq. (1),
![]() |
with
![]() |
where ρt is target surface reflectance, ρTOA is the top of atmosphere (TOA) reflectance, ρe is total atmospheric reflectance due to aerosol and molecular scatterings, θs is SZA, θv is VZA, φs is SAA, φv is VAA, φs − φv is relative azimuth angle (RAA), zt is target altitude, S is the spherical albedo of the atmosphere, T↓ is the total downward transmittance, T↑ is the total upward transmittance, and Tg is the gaseous transmission of atmospheric gases, including H2O, O2, O3, CO2, N2O and CH4.
We employed an LUT-based method commonly used for atmospheric correction (Seong et al. 2020; Peng 2020; Kim et al. 2022) to reduce the processing time. Table 2 presents the design details of the LUT used in this study, which were based on previous studies (Seong et al. 2020). The range of input parameters used in the 6SV model was limited to 0° to 80° for SZA and VZA with 5° incremental steps. The RAA ranged from 0° to 180° in 10° increments. The atmospheric input conditions consisted of irregularly spaced values, including the 12 AOT values. The surface elevation range considered ranged from 0 to 8 km with 2 km increments, the ozone range was 0.20 atm-cm to 0.40 atm-cm with 0.05 atm-cm increments, and the TPW range was 0 to 7 g cm−2 with 1 g cm−2 increments, based on statistics from CAMSRA-EAC4 between 2015 and 2020 within the FD region. Additionally, only maritime and continental aerosol types were employed because of the difficulty in classifying the continental and urban types and dominance of islands and continents over urban areas in the Himawari-8/AHI observational region.
The classification of aerosol types was based on matching the total columns of the five aerosol types (dust, sea salt, organic aerosol, black carbon, and sulphate) in the CAMSRA-EAC4 dataset with the aerosol types in the 6SV RTM, as aerosol-type data are not readily available (Shen et al. 2019). The classification was based on the percentage contribution of each aerosol component. We defined the grids in which the maximum component is the marine component (Oceanic in 6SV RTM) as maritime aerosol and the rest as continental aerosol. Table 3 presents the correspondence between the aerosol types used in the CAMSRA-EAC4 and those used in the 6SV. In 6SV RTM, the aerosol profile is assumed to be exponential with a scale height of 2 km.
In this study, a kernel-driven semi-empirical BRDF model was employed to estimate BRDF information (Matsuoka et al. 2016). This approach is valuable for estimating the reflectance of challenging-to-measure surfaces and modeling the reflectance under various illumination and observation conditions (Lucht et al. 2000; Schaaf et al. 2002). Similar to MODIS, we assumed the land surface to be isotropic and performed BRDF information estimation. The angular-adjusted LSR is defined by Eq. (2).
![]() |
where ρ is the angular-adjusted LSR; Kvol and Kgeo are the volume and geometric kernel values, respectively; and fiso, fvol, and fgeo refer to the kernel model parameters.
The Ross-Thick (RTK) kernel was chosen as the volume kernel, whereas the Li-Sparse-Reciprocal (LSR) kernel was selected as the geometric kernel, same as the MODIS BRDF product (Schaaf et al. 2002). The RTK-LSR combination model has been widely recognized for its effectiveness in inverting GEO satellite BRDF kernel model parameters (Matsuoka et al. 2016; Zhang et al. 2022). The calculation of each kernel is described in Eqs. (A1) and (A2) in the Appendix.
We performed multiple linear regressions on a pixel-by-pixel basis using the time series to estimate the BRDF kernel model parameters. Furthermore, we retrieved angular-adjusted LSR by selecting the center day of a rolling three-day window. We restricted this synthesis period to the local timeframe of 10:00–17:00 local time, ensuring that the data collected would be representative of daytime conditions. A time series of consecutive cloud-free periods was selected to minimize the effect of clouds on the multiple linear regression.
2.4 EvaluationThe following comparisons were made to evaluate the retrieved LSR and assess their accuracy. (1) Intercomparison between Spectral Band Adjustment Factors (SBAF)-adjusted AHI LSR vs. MOD09 LSR using ray-matched pairs obtained over a tropical Asia region in August 2018. (2) Inter-comparison between SBAF-adjusted, angle-adjusted AHI LSR vs. MOD09 LSR for a 0.3°-by-0.3° area located at the center of Australia at January 3, 2018. where angular-adjusted LSR was made using a synthesis period of January 2–4, 2018. (3) Inter-comparison between SBAF-adjusted AHI LSR vs. angle-adjusted MODIS LSR (MCD43A1) for the same 0.3°-by-0.3° area located at the center of Australia on January 3, 2018. (4) Inter-comparison between the spatially interpolated CAMSRA-EAC4 AOT vs. AERONET/SKYNET AOT over the Himawari FD region in 2018, and (5) Inter-comparison between AHI local-noon LSR estimated with CAMSRA-EAC4 AOT vs. the same LSR but estimated with in-situ AOT at two SKYNET sites in 2018.
a. Spectral band adjustment factorTo account for the differences in the SRFs between AHI and MODIS, which could lead to divergent results when observing identical targets (Chander et al. 2013; Li et al. 2019; Okuyama et al. 2018), we employed the SBAF tool (Scarino et al. 2016). This tool has gained extensive acceptance in sensor cross-calibration studies owing to its efficacy in rectifying spectral discrepancies (Kim et al. 2021; Yu and Wu 2016). Utilizing Eq. (3), along with the SBAFs supplied, we adjusted the MODIS LSR to better align with AHI measurements for a more accurate intercomparison.
![]() |
where ρMODIS is the MODIS LSR and ρAHI, MODIS is MODIS LSR after employing SBAF.
However, Band 06 (2.3 µm) of AHI falls outside the spectral coverage of the Scanning Imaging Absorption Spectrometer for Atmospheric Chartography (SCIAMACHY) used by the SBAF. Consequently, in this study, we focused on adjusting bands 01 to 05 of AHI for 11 landcovers.
b. Ray-matchingThe ray-matching method was employed to generate data pairs for comparison between the LEO and GEO sensors. Ray-matching screening data pairs with similar observation and illumination conditions by constraining the sun, sensor, and target geometry. This method has been extensively applied to numerous LEO-GEO inter-calibration studies including AHI-MODIS (Qin and McVicar 2018), AHI-VIIRS (Yu and Wu 2016), and ABI-VIIRS (Uprety et al. 2020; Jing et al. 2020).
We performed ray-matching between the AHI data and daily LSR products from MODIS (Terra: MOD09 GA, Aqua: MYD09GA), setting a matching interval of 0.02° in August 2018. We further evaluated the quality of this match using a scatter plot of the filtered pixels. This process was accomplished using the AHI-MODIS screening criteria established by Qin and McVicar (2018), where the difference of the VZA and VAA was less than 1° and 10°, respectively (see Table 4). Furthermore, to reduce the uncertainties introduced by factors, such as cloud cover, resampling, and aerosols, we introduced four additional matching conditions (See Table 4). The spatial distribution of the matching results (Fig. 4) shows only very limited regions in low latitudes were selected.
To address the limitation of the ray-matching method, which can only derive data pairs in low-latitude regions (Fig. 4), we employed the estimated BRDF kernel model parameters to estimate the angular-adjusted LSR of AHI at the MODIS observation angle. Because of the fixed observation angle of AHI, we selected a flat area covering a size of 0.3° × 0.3° as the evaluation area in this study (Fig. 5). The topography of the area is relatively simple, with a southeast to northwest aspect (Fig. S1). The land cover type is open shrubland based on MCD12Q1. Through visual inspection, we selected January 2 to 4, 2018, as the synthetic period of continuous clear skies. The angular-adjusted LSR was estimated by calculating the kernel value of the observation angle of MODIS/Terra and Aqua and intercomparison with MOD/MYD09GA.
Spatial distribution of areas where AHI-MODIS sensor observation conditions match. The light blue and orange regions indicate the matching regions of AHI-MODIS/Terra and AHI-MODIS/Aqua, respectively. AHI, Advanced Himawari Imager; MODIS, Moderate Resolution Imaging Spectroradiometer.
Location and image of Angular-adjusted LSR estimation area. LSR, Land Surface Reflectance.
The inputted AOT (i.e., CAMSRA-EAC4 data) were consistent with in-situ observations (Fig. 6, Table S1). Most sites exhibited low RMSE values (RMSE < 0.15), except urban locations such as Beijing, Taipei_CWB, and Mandalay_MTU. In contrast, the Lauder site in New Zealand (Liley and Forgan 2009) exhibited minimal AOT and the lowest RMSE.
Spatial distribution of RMSE of observation sites for AOT used in the validation of CAMSRA-EAC4 data. The values above each site name indicate the RMSE calculated from the comparison between CAMSRA-EAC4 data and site observation data over 2018 and 2019. AERONET, Aerosol Robotic Network; RMSE, root mean square error.
Regardless of the differences in AOT between in-situ observation and our inputs, the impact of the error in AOT on LSR was small (Fig. 7). In the results of the MYK site, where the AOT error is small, the correlation coefficient r of all the bands is greater than 0.99, the RMSE is below 0.005, and all the data are closely distributed on both sides of the 1:1 line. Conversely, at the BEJ site, where the AOT error is larger, the r is greater than 0.94 for all six bands. However, the bias of AHI band 01 to 04 is larger than that of band 05 and 06. The majority of the points in these bands are concentrated around the 1:1 line, but certain values that are over- or under-estimated.
Scatter plot of LSR estimated using (x-axis) in-situ AOTand (y-axis) CAMSRA-EAC4 AOT in Beijing and Miyako. The red dashed line and black soil line represent the regression line and 1:1 line, respectively. r, RMSE, and bias are the correlation coefficient, root-mean-square error, and bias, respectively. The density color bar, ranging from blue to red, denotes the concentration of data points. Bias is calculated as AHI LSR – MODIS LSR. AOT: aerosol optical thickness; CAMSRA-EAC4: Copernicus Atmosphere Monitoring Service Reanalysis ECMWF Atmospheric Composition Reanalysis; LSR: Land Surface Reflectance.
As a result of intercomparing ray-matching conditions over low latitude regions (Fig. 4), the estimated AHI LSR (post-SBAF adjustment) exhibits strong consistency with those of MODIS, as indicated by r exceeding 0.82, and most data points clustered around the one-to-one line (Fig. 8). Table 6 presents statistics on land cover types for the data pairs used in the AHI-MODIS intercomparison in Fig. 8. Unlike the results for AHI-MODIS/Terra, the AHI-MODIS/Aqua comparison shows two distinct clusters of scatter points across all bands except for Band 01. However, these clusters are distributed on both sides of the 1:1 line. Thus, AHI LSRs estimated in this study are consistent with those of MODIS over the geographic regions covered by the ray-matching.
Scatter plots of (x-axis) MODIS LSR and (y-axis) AHI estimated LSR. All matching points were obtained in August 2018. The red dashed line represents the regression line, and the black solid line represents the 1:1 line. r, RMSE, and bias are the correlation coefficient, root-mean-square error, and bias, respectively. The density color bar, ranging from blue to red, denotes the concentration of data points. Bias is calculated as AHI LSR – MODIS LSR. AHI, Advanced Himawari Imager; LSR, Land Surface Reflectance; MODIS, Moderate Resolution Imaging Spectroradiometer; RMSE; root mean square error.
The relationships between AHI and MODIS LSR varied for each band. For Band 03 (0.64 µm) and Band 04 (0.86 µm), which are of importance for vegetation monitoring, the regression slopes are close to one, ranging from 0.889 to 1.011. These bands also maintain high r between 0.923 to 0.978, demonstrating their strong linear relationship. The biases for these bands are minimal, within ±0.003. Similarly, the results for Band 01 (0.47 µm) and Band 05 (1.6 µm) also indicate a good agreement, with regression slopes ranging from 0.76 to 1.084 and r between 0.82 and 0.945. The biases for these bands are slightly positive, between 0.003 and 0.02, indicating a tendency for the AHI LSR to be slightly higher than the MODIS LSR. Conversely, the regression slopes for Band 06 for both Terra and Aqua are lower at approximately 0.719. Particularly, Band 02 (0.51 µm) exhibits the least consistency among all the bands. Furthermore, Band 02’s linear regression slope for Terra is low at 0.549, and inconsistent between Terra and Aqua.
c. Angular-adjusted LSRWithout matching the observation condition (sun-target-sensor geometry), the AHI and MODIS LSR were not as consistent as when those were matched (Fig. 9). Direct comparison of MODIS LSR data with AHI LSR estimates—conducted without matching and BRDF adjustments, even if reveals a strong linear relationship in AHI Bands 03 and 04 (r > 0.78), and the linear regression coefficients are close to 1 (1.012 to 1.117), but there are twice as many RMSEs and bias as in the corresponding results of Fig. 8. Similarly in the results for the rest of the bands, this increased bias remains even though r > 0.5 and the slope of the linear regression near 1 (e.g., AHI-MODIS/Terra Band 05). The evaluation indexes of AHI-MODIS/Terra and AHI-MODIS/Aqua in AHI Bands 01 and 02 are close to each other. Conversely, in the other spectral bands, the performance of AHI-MODIS/Aqua is somewhat inferior to that of AHI-MODIS/Terra, particularly in Band 06, which exhibits the highest RMSE at 0.072.
The angular-adjusted LSR estimated using the BRDF information shows good linear relationship with MODIS LSR (Fig. 10), for all bands, with r exceeded 0.5. In the case of Bands 03 and 04, the r was improving to over 0.89, while the regression slopes are close to one (0.904 to 1.169), demonstrating their strong consistency. Also, the bias for these bands are minimal, 0.007 for Band 03, and 0.017 for Band 04. However, in Bands 01 and 02, the dynamic range is narrower due to the centralized data distribution and lower value domain. Even though r is greater than 0.733 for all bands except AHI-MODIS/Terra Band 01, the linear regression slope is low, ranging from 0.443 to 0.649. In Band 02 (green), even though the domain of values is small, a large bias (0.015) still occurs. Linear regression slopes for AHI Bands 05 and 06 improved to over 0.56, and r exceeded 0.62, indicating a good consistency. In Bands 03-06, the angular-adjusted LSR computed under the Terra observational condition agrees better with the MODIS/Terra LSR than the corresponding MODIS/Aqua observational condition.
3.2 Estimated AHI LSRThe RGB image composited using LSR provides a better representation of the true color of the ground (Fig. 11). We generated cloud-free RGB images using 7 days of data from UTC 02:00 to 06:00 with a gamma value of 2.2. The LSR composite image (Fig. 11b) effectively removed the atmospheric haze and more clearly showed the ground surface than the TOA composite image (Fig. 11a). The LSR composite displays these regions relatively vividly, which is likely closer to the actual color of the terrain, particularly in the case of Australia’s distinct red soils and desert landscapes. However, in areas with high VZA and SZA, the LSR composite image is reddish in color compared to the TOA composite image.
Time series of normalized difference vegetation index (NDVI) show clear seasonal pattern with overall larger NDVIs in LSR compared with TOA values (Fig. 12). In the Deciduous Needleleaf Forest site (FHK), there is a clear seasonal pattern with NDVI values peaking in the summer and declining towards the winter. The LSR NDVI shows a high amplitude in these seasonal peaks and troughs compared to the TOA NDVI. The deciduous broadleaf forest site (TKY) follows a similar seasonal trend as the needleleaf forest, with NDVI values rising during the growing season and falling during the deciduous season. The LSR NDVI peaks are higher than those from the TOA, and the seasonal changes are marked. For the open forest savanna site (AU_Dry), the NDVI values show little variation throughout the two years, with both the LSR and TOA reflectance NDVI relatively constant and close together. The seasonal fluctuations in NDVI values for evergreen broadleaf forests (YNF) were much more subdued. While LSR NDVI remained at a high value (0.8) throughout the year, TOA NDVI was lower and had slight fluctuations.
We employed the CAMSRA-EAC4 dataset as the input for the 6SV model. Unlike the daily CAMS near-real-time dataset used by Lee et al. (2020a) and Seong et al. (2020), the CAMSRA-EAC4 provides improved temporal resolution from daily to 3 hourly. Consistency between CAMSRA-EAC4 and in-situ AOT data has also improved compared with the daily CAMS near-real-time dataset (Lee et al. 2022). In addition, the interpolated CAMSRA-EAC4 based water vapor, ozone and AOT550 were also consistent those by AERONET measurements (Fig. S2). The r between CAMSRA-EAC4 and AERONET AOT exceeded 0.82. For water vapor and ozone, this coefficient was even higher, surpassing 0.9. The global coverage of CAMSRA-EAC4 allows it to be used as an input for atmospheric corrections from other GEO satellites. Thus, our current selection of CAMSRA-EAC4 data as inputs of atmospheric parameters is one of the best available datasets if we aim to apply it to other GEO satellites.
Scatter plots of (x-axis) MODIS LSR and (y-axis) AHI LSR directly inter-comparison. The red dashed line and black soil line represent the regression line and 1:1 line, respectively. r, RMSE, and bias are the correlation coefficient, root-mean-square error, and bias, respectively. The density color bar, ranging from blue to red, denotes the concentration of data points. Bias is calculated as AHI LSR – MODIS LSR. Bias is calculated by AHI LSR – MODIS LSR.
We developed an aerosol type map that categorizes aerosols into maritime and continental types using the CAMSRA-EAC4 dataset. This categorization offers new insights, helping bridge data gaps in current 6SV atmospheric correction research. In the 6SV RTM, continental aerosol types are classified into three components: dust-like (70 %), water-soluble (29 %), and soot (1 %). Similarly, urban aerosols are defined with the following components: 17 % dust-like, 61 % water-soluble, and 22 % soot (Vermote et al. 2006). We classified black carbon and organic matter of CAMSRA-EAC4 as the soot and water-soluble components, respectively (Table 3). Since 84 % of organic matter is water-soluble in CAMSRA-EAC4 (Inness et al. 2019), the remaining particles create discrepancies in classifying between continental and urban aerosol types.
Scatter plots of (x-axis) MODIS LSR and (y-axis) AHI angular-adjusted LSR. The synthesis period was January 2–4, 2018 and the MODIS data acquisition date was January 3, 2018. The red dashed line and black soil line represent the regression line and 1:1 line, respectively. r, RMSE, and bias are the correlation coefficient, root-mean-square error, and bias, respectively. The density color bar, ranging from blue to red, denotes the concentration of data points. Bias is calculated as AHI LSR – MODIS LSR. Bias is calculated by AHI LSR – MODIS LSR.
The ray-matching results revealed a strong linear relationship (all of r > 0.78 in ray-matching with MODIS; Fig. 8) underscores the strong consistency between the estimated AHI LSR and the MODIS LSR product. In particular, r is greater than 0.9 in both Band 03 and 04, which is close to the evaluation results reported by Li et al. (2019). Compared to the results in the previous study on AHI-LEO sensor ray-matching (Yu and Wu 2016; Qin and McVicar 2018), even though we performed rigorous cloud screening, our result exhibited outliers (Fig. 8). Ray-matching is a common method used in sensor cross-calibration studies; therefore, its matching area encompasses land, oceans, and clouds. When applied ray-matching to land surface product evaluation, the ocean and clouds become interferences. Moreover, because the homogeneity of the land is worse than that of the sea, it is more likely to be affected by geo-location errors.
Comparison of Himawari-8 AHI RGB compositing images using (a) top-of- atmosphere reflectance (taken at 0300 UTC on May 1, 2019) and (b) Land Surface Reflectance (May 1–7, 2019). AHI, Advanced Himawari Imager.
NDVI time series at four sites in 2018 and 2019: TKY (36.15°N, 137.42°E); FHK (35.44°N, 138.76°E); YNF (26.75°N, 128.21°E); AU_Dry (15.26°S, 132.37°E). Orange and green represent the TOA reflectance and LSR, respectively. NDVI, Normalized Difference Vegetation Index; LSR, Land Surface Reflectance; EVI2, TOA, top of atmosphere.
The differences between the Terra and Aqua LSRs stem from the different land covers of the corresponding matching areas. In the matching results of August 2018, the available matching area of AHI-MODIS/Terra mainly constitutes forests, while the matching area of AHI-MODIS/Aqua includes cropland and natural vegetation mosaics in addition to forests (Table 6). Owing to the spectral properties of vegetation, this land cover difference was prominent in AHI Bands 02, 03, and 04 (Fig. 8). In AHI Band 02, the center wavelength of AHI is shorter than that of the corresponding band of MODIS (See Fig. 2), which is manifested by the lower LSR of AHI than that of MODIS in the results of AHI-MODIS/Terra. In Band 03, the highest density value in AHI-MODIS/Terra is lower than that of AHI-Aqua, and similarly in Band 04 AHI-MODIS/Terra is higher than that of AHI-MODIS/Aqua.
Furthermore, the band and land cover dependency ray-matching results showed LSR discrepancies were found only in Band 02 over evergreen broadleaf and woody savannas among different land cover types (Fig. S4).
c. Compensating the effect of different wavelength ranges of AHI and MODISThe application of SBAF effectively reduced differences between AHI and MODIS bands (Table 5). In AHI Bands 01, 03, and 04, the incorporation of SBAF enhances both the slope of the linear regression and the r, but bias did not change significantly. AHI band 05 has a slight improvement in r (from 0.879 to 0.892), but bias decreases from 0.001 to −0.012, and the LSR of AHI is higher than the adjusted MODIS LSR. Owing to the AHI’s band 02 center wavelength of which near towards the blue, the application of SBAF exacerbates the discrepancy between the two datasets, from r from 0.836 to 0.801. This result aligns with the findings presented in Qin and McVicar (2018).
d. BRDF estimation and angular-adjusted LSRThe significance of the BRDF effect becomes particularly pronounced when comparing the AHI with the MODIS LSR data, particularly only simultaneously without BRDF correction (see Section 3.1b). Notably, our BRDF-correction approach the r and linear regression slope for the Angular-adjusted LSR under MODIS observation conditions, while concurrently reducing bias. However, the inherent limitations of AHI as a single-angle sensor, unable to capture multi-angular surface observations, necessitate a nuanced consideration of hillslope effects on AHI’s BRDF estimations, as underscored by Matsuoka et al. (2016). This limitation is evident in the results presented in Fig. 10, where the observation angles of Terra (VAA ≈ 280°) and Aqua (VAA ≈ 85°) lead to differences in the results. In the case of our validation area (Fig. 5), the slope extends from northwest to southeast (Fig. S1). Terra and Aqua observe the slope from different directions, while AHI and MODIS/Terra observed the slope from similar directions. Therefore, LSR of AHI-MODIS/Terra achieves better consistency compared with that of AHI-MODIS/Aqua.
Due to the difference in observational modes between LEO and GEO, there is a difference between the LEO-based BRDF information and the GEO-based BRDF information. We used MODIS BRDF product to evaluate our BRDF parameter, as detailed in Fig. S5. By employing MODIS BRDF parameters from the MCD43A1 product, we computed the angular-adjusted LSR corresponding to AHI observation angles. The resulting scatter distributions and evaluation metrics across various bands closely align with those observed in Fig. 10.
e. Different approaches to retrieve LSR of GEO satellite dataThere are two primary methods for estimation of LSR of GEO satellite data. The first method is the traditional RTM approach, which requires additional atmospheric data as inputs, adopted in this study. The second involves simultaneously estimating AOT and LSR without relying on external atmospheric data [e.g., MAIAC; Li et al. (2019) and a Coupled RTM (Ma et al. 2021)]. Traditional RTM approaches have been used for decades (MODIS, VIIRS) due to simplicity and computational efficiency. On the other hand, methods based on simultaneously estimating AOT and LSR have been applied widely in recent years (Lyapustin et al. 2018; Li et al. 2019). Further intercomparison of outputs based on these two approaches is needed to quantify the uncertainties caused by different approaches.
4.2 Applicability to other 3rd generation GEO sensorsOur methodology for calculating the LSR is designed to be globally applicable and can be adapted for use with other GEO satellites. By simply modifying the SRFs and observation angles to match those of different satellites, this method offers a versatile framework for extending the LSR calculations to a broader range of satellite systems. Upon the release of this method and its corresponding code, creating a hyper-temporal and high-spatial resolution global dataset that integrates data from Himawari-8/9 and other geostationary satellites would become feasible. This facilitated the development of a homogeneous global dataset using a uniform processing algorithm.
We showed some practical evidence to support the potential applicability of our algorithm for application to other GEO satellites. First, evaluation of CAMSRA-EAC4 data using in-situ AOT revealed remarkable consistency. Therefore, CAMSRA-EAC4 data are reasonable for global application. Second, essential atmospheric input parameters obtained by CAMSRA-EAC4 data can produce similar performance compared with MODIS products. We performed LSR estimation for MOD02 (MODIS TOA reflectance product) using CAMSRA-EAC4 data and compared the results with those obtained using MOD09 (MODIS LSR product) in Fig. S6, the estimated LSR by MODIS using our algorithm with CAMSRA-EAC4 data produced similar values with MODIS LSR product. Third, SBAF effectively reduced inter-sensor variability caused by response functions (Table 5 and see discussion Section 4.1a). Lastly, leveraging globally observable sensors, such as MODIS and VIIRS, as intermediaries facilitates the cross-calibration of various geostationary satellites, thereby enriching the quality and comparability of the data collected.
4.3 Limitations and FutureFirst, although the AOT of CAMSRA-EAC4 closely matches to the in-situ AOT, AOT errors continue to persist owing to high pollution levels and biomass burning (Hoque et al. 2020). Moreover, because only two aerosol models were used in this study, limitations exist in urban or desert aerosols (Shen et al. 2019). Additionally, the spatial resolution of CAMSRA-EAC4 (0.75°) is much coarser than that of AHI (about 0.01°). Meanwhile, the current twice-yearly update frequency of the CAMSRA-EAC4 data presents a timeliness challenge that needs to be addressed. The Himawari-8 hourly AOT data provided by JAXA and the ECMWF Reanalysis v5 (ERA5) dataset can be used as alternative data owing to its real-time capability.
Second, evaluation by angular-adjusted LSR in mid-latitude region is still not a direct evaluation. A potential method involves exploring ray-matching using LEO satellites capable of capturing a higher VZA. Current off-nadir sensors, such as GCOM-C/SGLI (Imaoka et al. 2010) and Terra/MISR (Diner et al. 1998), diverge from nadir sensors (e.g., MODIS, VIIRS). Their unique design, featuring both forward and backward tilt angles, allows them to capture different observation angles compared to conventional nadir observations.
Third, this study fails to adequately address topographic correction, particularly in challenging terrains, such as mountainous regions. These areas are often affected by geolocation errors and require precise orthorectification methods for maintaining high geometric accuracy, as mentioned by Matsuoka and Yoshioka (2023).
Lastly, challenges persist in the edge regions where the VZA and SZA exceed 70–80°. As evidenced in the results section and supported by studies such as those by Kim et al. (2022), these edge regions exhibit a “reddening effect” and are susceptible to atmospheric over-correction. In addition, in these regions, the method used for interpolating LUTs affects the accuracy of LSR estimation (Lee et al. 2015; Lee et al. 2020b).
In this study, we formulated and implemented an algorithm for estimating the LSR and angular-adjusted LSR from the Himawari-8/9 AHI data. The evaluation of the proposed method was inter-compared using LSR products from the MODIS sensors onboard the Terra and Aqua satellites. This algorithm encompasses the estimation of LSR using the 6SV RTM with CAMSRA-EAC4 data, and the derivation of angular-adjusted LSR based on the BRDF parameters estimated using a kernel-driven BRDF model.
During the evaluation process, we conducted a comprehensive evaluation using ray-matching, angular-adjusted AHI LSR, and angular-adjusted MODIS LSR. Our results show that the LSR values estimated using our proposed algorithm maintain a high level of agreement at both low and mid-latitudes, thus providing researchers with a high-frequency AHI LSR product. In addition to this, we evaluated the interpolated CAMSRA-EAC4 data using in-situ AOT and estimated the LSR using these two datasets to assess the uncertainty in the estimates. The uncertainty introduced into the LSR estimation process was lower than that observed in previous studies. This finding underscores the potential of CAMSRA-EAC4’s high temporal resolution for use in GEO satellite LSR estimation studies.
In conclusion, our study not only successfully estimated the Himawari-8/9 AHI LSR but also presents a promising algorithm that can potentially be adapted for LSR estimation studies involving other GEO satellites, such as FY-4A AGRI and GOES ABI. Finally, we have publicly availed the code and data from this study; our data contributes to global-scale terrestrial monitoring at higher time scales.
The datasets generated and/or analyzed in this study are publicly available and can be accessed at ftp://modis.cr.chiba-u.ac.jp/ichii/SEND_NEW/H8AHI_SR/. The code used for the analysis is available on GitHub at https://github.com/Lw46/.
The supplementary material includes Fig. S1 to Fig. S6, and Table S1. Figure S1 is similar to Fig. 5 but provides the DEM of the ROI along with the directions of observation for AHI, MODIS/Terra, and MODIS/Aqua. Figure S2 shows scatter plots comparing in-situ data from all AERONET sites (x-axis) with CAMSRA-EAC4 (y-axis), including (a) AOT550, (b) water vapor, and (c) ozone. Figure S3 presents scatter plots of AHI NDVI vs. MODIS NDVI for three cases: ray-matching (Fig. 8), angular-adjusted (Fig. 9), and direct comparison (Fig. 10). Figure S4 provides scatter plots of MODIS LSR (x-axis) and AHI estimated LSR (y-axis) for four dominant land cover types shown in Fig. 8. Figure S5 displays scatter plots comparing MODIS angular-corrected LSR (MCD43A1, x-axis) and AHI LSR (y-axis) in the ROI (Fig. 5). Figure S6 shows scatter plots of MODIS LSR from MOD09 (x-axis) and LSR estimated from MOD02 with CAMSRA-EAC4 data input and the 6SV RTM (y-axis) in the ROI (Fig. 5). Table S1 summarizes the evaluation indicators for the AERONET and SKYNET sites used in this study.
This study was supported by JSPS Core-to-Core Program (Grant Number: JPJSCCA20220008), JSPS KAKENHI (Grant Number: JP20H04320, JP21K 12227, JP22H03727, JP22H05004, JP20K20487, JP21 K05669, JP23KJ0304), JST SPRING (Grant Number JPMJSP2109), MEXT Virtual Laboratory (VL) project, the Environment Research and Technology Development Fund (JPMEERF20215005) of the Environmental Restoration and Conservation Agency of Japan and the JAXA 3rd research announcement on the Earth Observations (grant number 19RT000351).
Ross-Thick:
![]() |
where,
![]() |
Li-Sparse- Reciprocal:
![]() |
where,
![]() |