Enterprise LST Algorithm Development and Its Evaluation with NOAA 20 Data

: Satellite land surface temperatures (LSTs) have been routinely produced for decades from a variety of polar-orbiting and geostationary satellites, which makes it possible to generate LST climate data globally. However, consistency of the satellite LSTs from di ﬀ erent satellite missions is a concern for such purpose; an enterprise satellite LST algorithm is desired for the LST production through di ﬀ erent satellite missions, or at the least, through series satellites of a satellite mission. The enterprise LST algorithm employs the split window technique and uses the emissivity explicitly in its formula. This research focuses on the enterprise LST algorithm design, development and its evaluations with the National Oceanic and Atmospheric Administration’s (NOAA) 20 (N20) Visible Infrared Imaging Radiometer Suite (VIIRS) data available since 5 January 2018. In this study, the enterprise LST algorithm was evaluated using simulation dataset consisting of over 2000 proﬁles from SeeBor collection and the results show a bias of 0.19 K and 0.34 K and standard deviation of 0.48 K and 0.69 K for nighttime and daytime, respectively. The in situ observations from seven NOAA Surface Radiation budget (SURFRAD) sites and two Baseline Surface Radiation Network (BSRN) sites were used for LST validation. The results indicate a bias of − 0.3 K and a root mean square error (RMSE) of 2.06 K for SURFRAD stations and a bias of 0.2 K and a RMSE of ~2 K for BSRN sites. Further, the cross-satellite analysis presents a bias of 0.7 K and an RMSE of 1.9 K for comparisons with AQUA MODIS LST (MYD11_L2, Collection 6). The enterprise N20 VIIRS LST product reached the provisional maturity in February 2019 and is ready for users to use in their applications. − and


Introduction
Land surface temperature (LST) is one of the most important parameters in the weather and climate system controlling surface heat and water exchange with the atmosphere [1]. It has been produced as baseline products for those weather service satellite missions such as Earth Observation Satellites (EOS) and Joint Polar-orbiting Satellite System (JPSS) [2] and has been listed as one of the essential variables (ECV) in Global Climate Observing System (GCOS) [3]. However, consistency of the satellite LSTs from different satellite missions is a concern for generating the LST climate data record and for cross-satellite comparisons; an enterprise satellite LST algorithm is desired for the LST production through different satellite missions, or at the least, through series satellites of a satellite mission.
Many factors may contribute to the LST difference among multiple sensors such as algorithm difference, sensor difference, field of viewing geometry difference, cloud shadow impact along the path etc. Strategies should be developed specifically to eliminate the cross-sensor difference. A primary objective of the enterprise LST development is to provide a state-of-the-art LST algorithm that is applicable to multiple sensors and capable to produce consistent LST products. The enterprise LST algorithm should be simple, robust and flexible so that it can be applied to sensors onboard both geostationary orbit (GEO) satellite missions (e.g., ABI on GOES-R) and Low Earth orbit (LEO) satellite missions (e.g., Visible Infrared Imaging Radiometer Suite (VIIRS) on JPSS series). The enterprise algorithm development employs the same LST derivation algorithm over different sensors and performs the same regression procedure for the related look-up table (LUT) generations. In addition, the enterprise LST development also takes into account the consistency of the related input data. For example, the sensor data is in a harmonized format; the total column precipitable water (TPW) is from numerical weather prediction (NWP). In addition, it used the same land sea mask from NASA and both cloud mask and aerosol optical depth (AOD) input were derived from enterprise algorithm. Besides, concurrent with the enterprise LST development, the land surface emissivity (LSE) product, a very important input for LST retrieval, is also developed to provide the thermal bands emissivity, emissivity uncertainty as well as broadband emissivity of the spectral range between 8-14 µm [4]. All these are put together to ensure the consistency of the enterprise satellite LST products.
Over the past several decades, many algorithms have been proposed to treat the characteristics of various sensors onboard different satellites with different assumptions and approximations for the radiative transfer equation and LSEs. These algorithms can be roughly grouped into three categories: single-channel methods, multi-channel methods, and multi-angle methods, provided that the LSEs are known a priori. If the LSEs are not known, then the algorithms can be categorized into three types: stepwise retrieval method, simultaneous retrieval of LSEs and LST with known atmospheric information, and simultaneous retrieval of LST, LSEs, and atmospheric information [5]. The most widely used approach is the split window (SW) technique, in which the atmospheric effects are compensated using two or more adjacent thermal infrared (TIR) channels (typically at 10-12.5µm). The SW approach was first proposed by McMillin (1975) [6] and successfully applied to sea surface temperature (SST) retrievals. This method is computationally efficient and does not require accurate atmospheric profiles. Encouraged by its success in satellite SST retrieval, many SW approaches have been proposed for LST retrieval [7][8][9][10][11]) and widely used for producing the operational LST products [12][13][14][15]. Only SW algorithms will thus be considered in the development of the enterprise algorithm.
In this study, we focus on the development of an enterprise LST algorithm for the JPSS mission. The N20 satellite, was successfully launched in November, 2017 and the LST data has been available since 5 January, 2018. N20 LST is currently generated using the same derivation algorithm as the Suomi National Polar-Orbiting Partnership (S-NPP) satellite, i.e., surface type dependent split window algorithm. The enterprise LST algorithm will be firstly applied over N20 satellite data and its application to GOES16 satellite is in the plan.
The paper is organized as follows: the enterprise LST algorithm design and development is described in the second section, followed with the algorithm evaluation using simulation data, ground data and cross satellite LST data. The fourth section is discussions while conclusion is given in the final section.

VIIRS Sensor Characteristics
VIIRS onboard S-NPP, N20 and future JPSS series, a scanning radiometer, collects visible and infrared imagery and radiometric measurements of the land, atmosphere, cryosphere, and oceans. VIIRS extends and improves upon a series of measurements initiated by the AVHRR and MODIS. The advanced design of VIIRS will provide users with finer spatial resolution both at nadir and scan edges. VIIRS has 22 spectral bands covering wavelengths from 0.4 to 12.5 µm, providing data for the production of more than 20 Environmental Data Records (EDRs) including LST EDR. The M bands includes two split window channels used for the LST retrieval as shown in Table 1. Also note that Advanced Baseline Imager (ABI), the primary instrument on the GOES-R satellite series for imaging Earth's weather, oceans and environment, includes two split window channels. Therefore the split window LST retrieval algorithm can be applied to both ABI and VIIRS sensors.
Although the N20 satellite is designed to be the same as S-NPP satellite, the sensor response function (SRF) shows some minor difference particularly at the spectral of band 15 (M15) as shown in Figure 1. The spectral resolution of N20 VIIRS SRF is about 50 nm, which is much coarser than the VIIRS SRF for S-NPP satellite with the spectral resolution of 5 nm. It is found that the coarse resolution makes it sensitive to the spectral selection for the BT calculation. Therefore the linear interpolation from 50 nm to 5 nm was applied for the N20 VIIRS SRF.  Also note that Advanced Baseline Imager (ABI), the primary instrument on the GOES-R satellite series for imaging Earth's weather, oceans and environment, includes two split window channels. Therefore the split window LST retrieval algorithm can be applied to both ABI and VIIRS sensors.
Although the N20 satellite is designed to be the same as S-NPP satellite, the sensor response function (SRF) shows some minor difference particularly at the spectral of band 15 (M15) as shown in Figure 1. The spectral resolution of N20 VIIRS SRF is about 50 nm, which is much coarser than the VIIRS SRF for S-NPP satellite with the spectral resolution of 5 nm. It is found that the coarse resolution makes it sensitive to the spectral selection for the BT calculation. Therefore the linear interpolation from 50 nm to 5 nm was applied for the N20 VIIRS SRF.

Simulations Analysis
Considering that the polar orbiting satellite is of global coverage, the global representativeness of the profiles is essential. The selected profile should have a wide coverage of different climate and atmospheric conditions such as low and high temperature and wet and dry atmosphere. We used SEEBOR profile collection, which in its latest version Seebor 5.0, consists of 15704 global profiles of temperature, moisture, and ozone at 101 pressure levels for clear sky conditions. The profiles are taken from NOAA-88, an ECMWF 60L training set, TIGR-3, ozonesondes from 8 NOAA Climate Monitoring and Diagnostics Laboratory (CMDL) sites, and radiosondes of 2004 in the Sahara desert [16]. Each situation is described, from the surface to the top of the atmosphere, by the values of the temperature, TPW and ozone concentrations on a given pressure grid. Quality checks were applied to all the profiles along with the following saturation criteria: The relative humidity (RH) value of the profiles must be less than 99% at each level below the 250 hPa pressure level for clear sky conditions. In this study additional cloud screening has been applied in order to further reduce the impact from

Simulations Analysis
Considering that the polar orbiting satellite is of global coverage, the global representativeness of the profiles is essential. The selected profile should have a wide coverage of different climate and atmospheric conditions such as low and high temperature and wet and dry atmosphere. We used SEEBOR profile collection, which in its latest version Seebor 5.0, consists of 15704 global profiles of temperature, moisture, and ozone at 101 pressure levels for clear sky conditions. The profiles are taken from NOAA-88, an ECMWF 60L training set, TIGR-3, ozonesondes from 8 NOAA Climate Monitoring and Diagnostics Laboratory (CMDL) sites, and radiosondes of 2004 in the Sahara desert [16]. Each situation is described, from the surface to the top of the atmosphere, by the values of the temperature, TPW and ozone concentrations on a given pressure grid. Quality checks were applied to all the profiles along with the following saturation criteria: The relative humidity (RH) value of the profiles must be less than 99% at each level below the 250 hPa pressure level for clear sky conditions. In this study additional cloud screening has been applied in order to further reduce the impact from cloud residue using the method described in Galve (2008) [17]. In this method, profiles with relative humidity (RH) Remote Sens. 2019, 11, 2003 4 of 18 larger than 90% or two consecutive layers with RH larger than 85% or RH larger than 80% within the first two kilometers of vertical profiles are considered as cloudy or foggy condition so as excluded in the regression procedure. Table 2 summarizes the SEEBOR profile collection with cloud screened using above method with nearly 50% profiles are excluded. MODTRAN radiative transfer model (RTM) (version 5.2), one of the most accurate RTM models, was used to simulate the top of atmosphere radiance with the vertical atmospheric parameters such as temperature, ozone, and relative humidity in profiles. Since the MODTRAN simulation is time consuming, specific strategies were used in running the simulation. As there are still a considerable number of profiles after cloud screening, we vary the prescribed LST in five ranges LST-5, LST -2.5,LST, LST +2.5, and LST + 5 and emissivity ranging between 0.90 and 0.9999 at a step of 0.025. Then linear interpolation is performed for other Ts at 1 K step and emissivity at 0.00125 step. Considering the viewing angle effect has nonlinear feature, we iterated the sensor zenith angle from 0 to 70 degree with 5-degree increment. The mean sensor channel radiance was determined by integrating over the SRF and then converted into corresponding brightness temperatures using the Planck function. This step constructs the comprehensive database for generating the coefficients LUT of the LST algorithms. When the enterprise LST algorithm is applied to another sensor, the radiance integration and conversion step will be repeated with corresponding sensor SRF.
Above SEEBOR simulation database was then split into two subsets-one subset for the generation of the LST LUT and the other one for independent theoretical LUT evaluation. For LUT generation, it is very important to choose profiles that are cloud clear, with global distribution, and representativeness of regional climate. Besides, factors such as the seasonal effect and elevation should also be considered in the profile selection. In this study, the method described in Zhou et. al [18] was used for the profile selection as shown in Equation (1). In this method latitude and longitude distribution, elevation, TPW and seasonal distribution were considered.
where x 0 and y 0 are the longitude and latitude of the selected profile, respectively. x i and y i represent the longitude and latitude of the profiles to be selected. The unit for x and y is degree. z sur f 0 is the altitude of the selected profile and z sur f i is the altitude of the profiles to be selected. w 0 and w i , with unit of gcm −2 , are the TPW for the selected profile and the profile to be selected, respectively. m 0 and m i are months of the selected profile and the profile to be selected, respectively. Threshold were set as: ∆x is set as 60 degrees for high latitude and 30 degrees for mid/low latitude, respectively; ∆y is set as 15 degrees for high latitude and 10 degrees for mid/low latitude, respectively; ∆z is set as 1000 m; ∆w is set as 0.5 gcm −2 ; ∆m is set as 2 months. Considering the day night difference, the profiles are separated for daytime and nighttime. After the cloud screening, only 29 profiles are left with TPW greater than 5.0 gcm −2 , therefore they are all retained to better represent the humid atmospheric conditions. The surface temperature spans from 222 K to 330 K, which TPW reaches up to 6.3 gcm −2 and the difference between surface and air temperature ranges from −15 K to 29 K for daytime; for nighttime profiles, the surface temperature spans from 229 K to 330 K, while TPW reaches up to 6.1 gcm −2 and the difference between surface and air temperature varies from −15 K to 27 K. Figure 2 shows the global distribution of the profiles selected for the regression analysis for daytime and nighttime, respectively.
Remote Sens. 2019, 11, x FOR PEER REVIEW 5 of 18 and are months of the selected profile and the profile to be selected, respectively. Threshold were set as: Δx is set as 60 degrees for high latitude and 30 degrees for mid/low latitude, respectively; Δy is set as 15 degrees for high latitude and 10 degrees for mid/low latitude, respectively; Δz is set as 1000 m; Δw is set as 0.5 gcm −2 ; Δm is set as 2 months.
Considering the day night difference, the profiles are separated for daytime and nighttime. After the cloud screening, only 29 profiles are left with TPW greater than 5.0 gcm -2 , therefore they are all retained to better represent the humid atmospheric conditions. The surface temperature spans from 222 K to 330 K, which TPW reaches up to 6.3 gcm -2 and the difference between surface and air temperature ranges from -15 K to 29 K for daytime; for nighttime profiles, the surface temperature spans from 229 K to 330 K, while TPW reaches up to 6.1 gcm -2 and the difference between surface and air temperature varies from -15 K to 27 K. Figure 2 shows the global distribution of the profiles selected for the regression analysis for daytime and nighttime, respectively. The regression procedure is composed of two functional components. The first component is the data selection. Based on the simulation database mentioned above, this component selects all data required for the algorithm LUT generation. It takes into account the sensor noise, emissivity uncertainty and its global distribution. The sensor noise is set to NeDT level for the two split window channels. One year of global emissivity data was analyzed and frequency of emissivity data pairs was obtained and used in the data selection. Difference between surface temperature and air temperature was set to (-15,20) and (-15,15) for daytime and nighttime, respectively. Followed the data selection, regression of the LST algorithm will be performed to get the algorithm coefficients LUT.

Algorithms Determination
Under the clear sky condition, the sensor measured radiance within the TIR spectral range (8-14 µm) for the narrow band of a sensor channel includes three radiance contributions from surface outgoing emission, atmospheric path emission and atmospheric downwelling sky irradiance reflected by the surface. And the earth surface emitted radiance is a function of temperature and emissivity and is attenuated along the atmospheric path to the sensor. For LST derivation, a method has to be developed to decouple the emissivity and the atmospheric absorption effects from the satellite received radiance. Most commonly used methods are to linearize the process and combine multi-channel for optimal atmospheric correction. The Split Window (SW) algorithm gained its popularity for generating the operational LST product for various sensors (e.g., MODIS, VIIRS, Spinning Enhanced Visible and InfraRed Imager (SEVIRI) etc.) due to its advantage of simplicity, robustness, computational efficiency and not requiring the accurate atmospheric profile. The regression procedure is composed of two functional components. The first component is the data selection. Based on the simulation database mentioned above, this component selects all data required for the algorithm LUT generation. It takes into account the sensor noise, emissivity uncertainty and its global distribution. The sensor noise is set to NeDT level for the two split window channels. One year of global emissivity data was analyzed and frequency of emissivity data pairs was obtained and used in the data selection. Difference between surface temperature and air temperature was set to (−15,20) and (−15,15) for daytime and nighttime, respectively. Followed the data selection, regression of the LST algorithm will be performed to get the algorithm coefficients LUT.

Algorithms Determination
Under the clear sky condition, the sensor measured radiance within the TIR spectral range (8-14 µm) for the narrow band of a sensor channel includes three radiance contributions from surface outgoing emission, atmospheric path emission and atmospheric downwelling sky irradiance reflected by the surface. And the earth surface emitted radiance is a function of temperature and emissivity and is attenuated along the atmospheric path to the sensor. For LST derivation, a method has to be developed to decouple the emissivity and the atmospheric absorption effects from the satellite received radiance. Most commonly used methods are to linearize the process and combine multi-channel for optimal atmospheric correction. The Split Window (SW) algorithm gained its popularity for generating the operational LST product for various sensors (e.g., MODIS, VIIRS, Spinning Enhanced Visible and InfraRed Imager (SEVIRI) etc.) due to its advantage of simplicity, robustness, computational efficiency and not requiring the accurate atmospheric profile.
We studied various SW LST algorithms from the literature [9,[19][20][21][22][23][24][25]. The VIIRS LST algorithm for S-NPP satellite employed the surface type dependent SW approach, where the algorithm coefficients are stratified according to the surface type and day/night conditions [25]. Spectral emissivity is not Remote Sens. 2019, 11,2003 6 of 18 explicitly required in the algorithm. The enterprise LST algorithm will use emissivity explicitly as an input and exclude the term of (T 11 − T 12 ) 2 because: (1) The S-NPP LST algorithm performance presents a strong dependency on the accuracy of the surface type input [26]. (2) The emissivity variability cannot be well characterized by the surface type particularly those with considerable seasonal or annual surface dynamics. In addition, the surface type is generally mixed and heterogeneous at the moderate pixel footprint (e.g.,~1 km) which results in the mis-representation of the actual emissivity at different error levels. (3) The TPW is not explicitly considered which means the wet and dry condition is treated equally in the algorithm for each surface type, which results in the degraded performance in both dry and wet conditions for surface types with wide climate distributions. (4) The (T 11 − T 12 ) 2 term may lead to large LST uncertainty under specific conditions, e.g., when significant brightness temperature (BT) difference between the two split window bands (this can be as high as 16K over Australia) occurs. The simulation database is not sufficient to represent such conditions.
In this study, seven candidate algorithms (Table 3) were selected for enterprise LST production. Yu et al. (2008) [25] showed that, if an algorithm's coefficients are determined for typical TPW amounts, the algorithm accuracy can degrade significantly at large viewing angles unless a corrective term is applied. Therefore, the explicit path length correction term, i.e., (T 11 − T 12 ) (sec θ − 1) is added into the equation of five candidate algorithms including 1, 3, 4, 5 and 6. It means that the equations of the five algorithms are composed of two parts: the base split window algorithm and the path length correction as shown in Table 3. The algorithm 2 and 7 use implicit path correction by stratifying the retrieval with different satellite viewing angles. Table 3. Candidate land surface temperature (LST) retrieval algorithms.
As with most SW algorithms, the candidate algorithms explicitly use LSE values and allows easy incorporation of the new and improved global emissivity products (e.g., [27,28]. To support the enterprise LST development, LSE product was concurrently developed. The LSE product is generated based on the vegetation cover method (VCM), which combines two constant emissivity values from the bare ground and full vegetation coverage situations. Since bare ground emissivity has notable spatial variation, long term historical ASTER and MODIS emissivity products are employed to derive pixel by pixel bare ground emissivity working as background climatology. The daily emissivity is adjusted according to the VIIRS daily green vegetation fraction (GVF) and snow fraction products. The new daily LSE product (NOAA LSE Product thereafter) provides the spectral emissivity for VIIRS split window bands as well as the broadband emissivity at 1 km spatial resolution. In situ validation indicates the LSE product with uncertainty lower than 0.015 [4].
To select a suitable algorithm for the enterprise LST production, we analyzed the theoretical accuracy and sensitivity of the candidate SW algorithms using a comprehensive simulation dataset. As shown in Table 4 that all candidate algorithms except algorithm 4 yield similar accuracy with slight different uncertainty. Algorithm 2 and 7 have nearly the same uncertainty for all combined conditions of day/night, viewing geometry and TPW. The algorithm 7 was selected as the enterprise algorithm because of its consistently good performance over multiple sensors, simplicity and robustness. spatial variation, long term historical ASTER and MODIS emissivity products are employed to derive pixel by pixel bare ground emissivity working as background climatology. The daily emissivity is adjusted according to the VIIRS daily green vegetation fraction (GVF) and snow fraction products.
The new daily LSE product (NOAA LSE Product thereafter) provides the spectral emissivity for VIIRS split window bands as well as the broadband emissivity at 1 km spatial resolution. In situ validation indicates the LSE product with uncertainty lower than 0.015 [4].
To select a suitable algorithm for the enterprise LST production, we analyzed the theoretical accuracy and sensitivity of the candidate SW algorithms using a comprehensive simulation dataset. As shown in Table 4 that all candidate algorithms except algorithm 4 yield similar accuracy with slight different uncertainty. Algorithm 2 and 7 have nearly the same uncertainty for all combined conditions of day/night, viewing geometry and TPW. The algorithm 7 was selected as the enterprise algorithm because of its consistently good performance over multiple sensors, simplicity and robustness.

Algorithms sensitivity analysis
Theoretically all inputs may introduce uncertainty to the final LST product. In this study, two important error sources are considered: sensor noise and the surface emissivity uncertainty. We therefore analyzed the sensitivities of the selected enterprise LST algorithm with respect to the two factors. The simulation dataset described above is used in the sensitivity analysis.

Sensor noise uncertainty
The LST uncertainty δTs due to the sensor data noise can be described as where δT1 and δT2 represent the band uncertainties resulting from the brightness temperature uncertainties at 11 and 12 μm, respectively assuming that the two error sources are uncorrelated. With the selected enterprise algorithm, these two components can be expressed as Therefore, the LST uncertainty with respect to BT uncertainties is where δBT11 and δBT12 are set to be equal to NEdT i.e., 0.070 and 0.072 for the two VRIRS split window bands (Table 1), respectively.

Algorithms Sensitivity Analysis
Theoretically all inputs may introduce uncertainty to the final LST product. In this study, two important error sources are considered: sensor noise and the surface emissivity uncertainty. We therefore analyzed the sensitivities of the selected enterprise LST algorithm with respect to the two factors. The simulation dataset described above is used in the sensitivity analysis.

Sensor noise Uncertainty
The LST uncertainty δT s due to the sensor data noise can be described as where δT 1 and δT 2 represent the band uncertainties resulting from the brightness temperature uncertainties at 11 and 12 µm, respectively assuming that the two error sources are uncorrelated. With the selected enterprise algorithm, these two components can be expressed as Therefore, the LST uncertainty with respect to BT uncertainties is where δBT 11 and δBT 12 are set to be equal to NEdT i.e., 0.070 and 0.072 for the two VIIRS split window bands (Table 1), respectively. The results as shown in Figure 3 indicate that the sensor noise of 0.07 and 0.072 for BT 11 and BT 12 can cause average LST uncertainty of 0.27 K and 0.25 K for daytime and nighttime, respectively. The uncertainty increases with increasing TPW and sensor viewing angle.

Emissivity Uncertainty
Analytically, the LST uncertainty δT s due to the emissivity uncertainty can also be described as Equation (2) where δT 1 and δT 2 represent uncertainties resulting from the uncertainties of the mean emissivity (ε) and emissivity difference (∆ε), respectively. With the selected enterprise LST algorithm, these two components are Therefore, the maximum LST uncertainty for the selected enterprise algorithm is Considering that ε = (ε 11 + ε 12 )/2 and ∆ε = ε 11 − ε 12 , and assuming the emissivity uncertainties in each band are the same, i.e., ∆ε = ∆ε 11 = ∆ε 12 , the maximum uncertainty of the emissivity difference is δ(∆ε) = |δε 11 | + |δε 12 |= 2δε. Thus, the LST uncertainty, δT s , related to the emissivity uncertainty can be calculated using the Equation (6). For illustration purpose, we assumed that the mean emissivity (ε) and emissivity difference (∆ε) are 0.98 and 0.005, respectively. Results presented in Figure 4 show that the LST uncertainty (δT) increases approximately linearly with the emissivity uncertainty and can be as high as 3 K when the emissivity data has a fairly large uncertainty under very moist conditions. Similar sensitivity results were observed for the nighttime cases. Note, however, that the predicted LST uncertainty calculated using Equation (6) represents an extreme situation where δ(∆ε) and δ(ε) emissivity uncertainty worsen the LST retrieval (i.e., the errors are always superimposed on each other rather than weakening each other). In practice, the final LST error may be significantly smaller.

Algorithm Evaluation
Three methods including evaluation from simulation data, ground data evaluation and cross satellite evaluation, are used to evaluate the quality of the enterprise LST algorithm.

Evaluation from Simulation Data
As mentioned before, the simulation data is separated into two groups: one group is training database which is composed of 389 profiles for daytime and 409 profiles for nighttime; all remaining profiles are used for the independent evaluation. The same emissivity distribution is used and statistical analysis is performed over multiple parameters: latitude, TPW, and viewing zenith angles. Table 5 shows the overall evaluation results of the enterprise algorithm with a bias of 0.19 K and 0.34 K and a standard deviation (STD) of 0.48 K and 0.69 K for nighttime and daytime, respectively. In detail, the STD increases with the TPW and viewing angle as shown in Figure 5. When TPW is less than 1.5 gcm −2 , the STD is less than 1.0 K even under large viewing angle at 70 degree. When TPW is between 2.0 to 3.0 gcm −2 , the STD is beyond 1.0 K when the angle is larger than 30 degree. When the angle is beyond 50 degree, the STD tends to have a great increase up to 2.5 K under very moist condition for nighttime. For daytime, the trend is similar but more obvious increase at large angle greater than 60 degree is observed for STD up to 4 K, which suggests a high uncertainty of the LST algorithm under moist condition with large viewing angle. The bias does not show clear trend. It is from -0.5 K to 0.6 K for nighttime and 0.2 K to 0.8 K for daytime. The fitting performance is better at high latitude then followed mid and low latitude as shown in Figure 6.

Ground Data Evaluation
The in-situ LST measurements are often used for the validation of satellite LST products. It is important to note that satellite LST is the average radiative temperature over the satellite footprint at a certain viewing angle. The footprint for VIIRS sensor is 750 m at the nadir and 1.6 km on the scanning edge. The ground LST is usually measured from equipment mounted on the top of a tower, with about 70 m * 70 m field of view. Therefore it requires that the ground site is homogeneous enough to be representative of the satellite pixel. However such site is rarely available in reality. Besides, due to its rapid temporal variation, the in-situ LST measurements must have high temporal frequency.
In this study, the ground measurements from SURFRAD network and BSRN were used to validate the enterprise VIIRS LST. statistical analysis is performed over multiple parameters: latitude, TPW, and viewing zenith angles.   Table 5 shows the overall evaluation results of the enterprise algorithm with a bias of 0.19 K and 0.34 K and a standard deviation (STD) of 0.48 K and 0.69 K for nighttime and daytime, respectively. In detail, the STD increases with the TPW and viewing angle as shown in Figure 5. When TPW is less than 1.5 gcm −2 , the STD is less than 1.0 K even under large viewing angle at 70 degree. When TPW is between 2.0 to 3.0 gcm −2 , the STD is beyond 1.0 K when the angle is larger than 30 degree. When the angle is beyond 50 degree, the STD tends to have a great increase up to 2.5 K under very moist condition for nighttime. For daytime, the trend is similar but more obvious increase at large angle greater than 60 degree is observed for STD up to 4 K, which suggests a high uncertainty of the LST algorithm under moist condition with large viewing angle. The bias does not show clear trend. It is from -0.5 K to 0.6 K for nighttime and 0.2 K to 0.8 K for daytime. The fitting performance is better at high latitude then followed mid and low latitude as shown in Figure 6.

Ground Data Evaluation
The in-situ LST measurements are often used for the validation of satellite LST products. It is important to note that satellite LST is the average radiative temperature over the satellite footprint at a certain viewing angle. The footprint for VIIRS sensor is 750 m at the nadir and 1.6 km on the scanning edge. The ground LST is usually measured from equipment mounted on the top of a tower, with about 70 m * 70 m field of view. Therefore it requires that the ground site is homogeneous enough to be representative of the satellite pixel. However such site is rarely available in reality. Besides, due to its rapid temporal variation, the in-situ LST measurements must have high temporal frequency. In this study, the ground measurements from SURFRAD network and BSRN were used to validate the enterprise VIIRS LST.
The ground LST is derived from its upwelling and downwelling radiance flux by Stefen-Boltzmann law The ground LST is derived from its upwelling and downwelling radiance flux by Stefen-Boltzmann law Remote Sens. 2019, 11, 2003 11 of 18 where F ↑ and F ↓ are upwelling and downwelling longwave radiation flux, respectively; ε is the broadband surface emissivity and σ is the Stefan-Boltzman constant (σ = 5.67051 × 10 −8 Wm −2 k −4 ) and Ts is surface skin temperature. Ts is then obtained by inverting Equation (8): Here, the broadband emissivity is derived from NOAA LSE product as mentioned in Section 2.1. The enterprise LST product is not yet in operation, therefore it was calculated locally using N20 data for the time period from 5 January 2018 to 31 December 2018. Considering the cloud residue effect and ground in-situ data quality control, the following procedures were used for the matchup between the satellite LST and ground in-situ LST measurements: (1) The temporal difference is less than 86 s, which is the typical duration of a single granule (2) Spatially closest pixel is used for the comparison (3) Confidently clear indicated by the cloud mask product (4) To reduce the nearby cloud impact, the 3*3 neighboring pixels are all marked as confidently clear (5) The standard deviation of the band 15 brightness temperature in the neighboring 3 by 3 box is less than the threshold, which is set as 1.5 for all sites. (6) The standard deviation of the 30 min (centered at the matchup time) downwelling radiation from in-situ observations is less than a predetermined threshold, which is set as 1.2 for all sites.

Validation Results of SURFRAD Stations
NOAA SURFRAD network site locations were chosen with the intent of best representing the diverse climates of the United States. Special consideration was given to places where the landform and vegetation are homogeneous over an extended region so that the point measurements would be qualitatively representative of a large area [29]. Besides SURFRAD provides long term high quality downwelling and upwelling infrared radiation, along with other meteorological parameters [30]. It has been widely used for satellite LST validation [1,18,26,[31][32][33].
The uncertainty of the field pyrometer is about 5 W atts /m 2 and the effective foot print over 10 m tower is about 70 m by 70 m. The data is with 3 min interval before 2009 and 1 min thereafter. There are seven SURFRAD stations spatially distributed in the continental US in different climate regions with details shown in Table 6. Besides above procedures, additional temporal filter is applied over the site in Bondville, IL. During late spring and early summer, an obvious daytime discrepancy between the satellite LST and in-situ measurement has been reported from previous studies [26,31,33]. The satellite LST is 6 K to 10 K warmer than the ground in-situ observations. Li et al. (2014) [31] compared 10 years 16-day average NDVI and daily emissivity datasets from the MODIS observations and found that this feature might be caused by anomalous NDVI-emissivity relationship, i.e., emissivity does not change accordingly with the NDVI change during the time period. Guillevic et al. (2012) [33] mentioned that validation results obtained for stations surrounded by croplands present strong seasonal dependency: station observations may be closer/deviate more from the temperature of surrounding fields, according to crop maturity. When crops are short, they will have little shade on ground. Satellite views the scene at a certain angle so as the soil surface emission has a great contribution to the warm satellite LST estimation. Radiometer of the station, however, is always looking down at nadir measuring the surface upwelling radiance from the vegetation canopy. Such observation difference might be the reason leading to lower ground observations than satellite LST estimation. Considering its great impact on the validation, the matchups in May and June were removed from the results over BND site. Table 7 lists the validation results against SURFRAD observations. The overall results show an accuracy of 0.3 K and a precision of 2.04 K. The results vary among the sites. Table 6. Geolocation and surface type of the seven Surface Radiation budget (SURFRAD) stations.

Validation Results of BSRN Sites
The radiometric network BSRN was launched in 1992 by World Climate Research Program (WCRP) to support the research projects of the WCRP and other scientific programs. It provide typically 1-min averaged short-and long-wave surface radiation fluxes of the best possible quality currently available. As of mid-2013, the data import is organized in so-called station-to-archive files, which contain all the data from one station collected during one month. Currently a total of over 7000 station-month datasets from 58 stations are available in the World Radiation Monitoring Center (WRMC) [34]. Though there are 58 stations available in BSRN, only two stations are selected for LST validation. There are several reasons: Firstly, there are only 8 sites with long wave upwelling observations, which are required in the ground LST calculation. Secondly, some sites have retired so that no temporal overlap with VIIRS measurements; thirdly some sites do not satisfy the thermal homogeneity requirements for LST validation particularly some sites are close to water body. The selected two sites are located at Gobabeb, Namibia (GOB) and Cabauw, The Netherlands (CAB).
GOB site is located in the hyper-arid climate of the Namibia desert. Note that for GOB site, the downwelling radiation sensor and upwelling radiation sensor are mounted at different locations with about 10 km distance between them. The downwelling sensor is located near Gobabeb Research and Training Center (latitude 23.56 • S, longitude 15.04 • E) where the nearby Kuiseb river forms a natural boundary between large gravel plains (> 900 km 2 ) and the sand dunes of the Namib Desert [35]. The upwelling sensor is located at latitude 23.519 • S and longitude 15.083 • E, where it is homogeneous gravel plains with very sparsely grass coverage. The matchups between in situ and satellite data is set at the homogeneous upwelling sensor location. The Figure 7 shows the validation results. Some outliers due to sub-pixel clouds can be seen (Station LST higher than VIIRS LST). Green represents nighttime matchups and red represents daytime matchups. VIIRS LST is about 0.5 K colder than ground LST measurements at nighttime but nearly no bias at daytime over Gobabeb, Namibia (GOB) site.
gravel plains with very sparsely grass coverage. The matchups between in situ and satellite data is set at the homogeneous upwelling sensor location. The Figure 7 shows the validation results. Some outliers due to sub-pixel clouds can be seen (Station LST higher than VIIRS LST). Green represents nighttime matchups and red represents daytime matchups. VIIRS LST is about 0.5 K colder than ground LST measurements at nighttime but nearly no bias at daytime over Gobabeb, Namibia (GOB) site. CAB site is mainly covered by cropland with residential housing in the nearby area. The validation result shows an overall bias of 0.2 K and RMSE of 2.03 K. VIIRS LST is about 1.1 K colder and 0.7 K warmer than the ground in-situ LST measurements for daytime and nighttime measurements, respectively.

Cross Satellite Comparison with MODIS LST
In this study, the enterprise VIIRS LST is compared with the operational MODIS LST product. In order to reduce the data noise, firstly the temporal difference is controlled within 10 min of the simultaneous nadir overpass (SNO) scenario. The pixel level matchup requires that both LST products indicate the cloud clear condition; the spatially closest pixel and the viewing angle difference is within 10 degree.  CAB site is mainly covered by cropland with residential housing in the nearby area. The validation result shows an overall bias of 0.2 K and RMSE of 2.03 K. VIIRS LST is about 1.1 K colder and 0.7 K warmer than the ground in-situ LST measurements for daytime and nighttime measurements, respectively.

Cross Satellite Comparison with MODIS LST
In this study, the enterprise VIIRS LST is compared with the operational MODIS LST product. In order to reduce the data noise, firstly the temporal difference is controlled within 10 min of the simultaneous nadir overpass (SNO) scenario. The pixel level matchup requires that both LST products indicate the cloud clear condition; the spatially closest pixel and the viewing angle difference is within 10 degree.
The comparison between N20 and AQUA LST is on regional scale, which is constrained by the SNO availability. The SNO selection starts from February, 2018. Some SNO were not selected due to considerable cloud coverage. The comparison therefore has a temporal and spatial limitation, which only represents the LST difference over North America, Africa, Australia, Mid-North Asia and Greenland area. The bias varies from −1.3 K to 1.74 K and the STD varies from 0.84 K to 3.18 K as shown in Table 8. As shown in Figure 8, the N20 VIIRS LST is on average 0.7 K warmer than AQUA MODIS LST with a STD of 1.8 K. The comparison presents a good coverage of temperature range from about 210 K to 340 K. An excellent agreement is obtained for the temperature below 260 K. The discrepancy increases over 300 K. Some outliers are observed for either VIIRS much higher than MODIS LST or vice versa, which might be associated to cloud contamination. The comparison between N20 and AQUA LST is on regional scale, which is constrained by the SNO availability. The SNO selection starts from February, 2018. Some SNO were not selected due to considerable cloud coverage. The comparison therefore has a temporal and spatial limitation, which only represents the LST difference over North America, Africa, Australia, Mid-North Asia and Greenland area. The bias varies from -1.3 K to 1.74 K and the STD varies from 0.84 K to 3.18 K as shown in Table 8. As shown in Figure 8, the N20 VIIRS LST is on average 0.7 K warmer than AQUA MODIS LST with a STD of 1.8 K. The comparison presents a good coverage of temperature range from about 210 k to 340 K. An excellent agreement is obtained for the temperature below 260K. The discrepancy increases over 300 K. Some outliers are observed for either VIIRS much higher than MODIS LST or vice versa, which might be associated to cloud contamination. The error distributions over sensor zenith angle (STZ) indicate VIISR LST is consistently 0.7 K warmer than MODIS LST when the VIIRS STZ is below 50 degree and then decrease after that. The difference STD generally increases with the angle particularly when it is over 50 degree. The difference also presents characteristics of changes in latitude. As shown in Table 9, the difference of the two sensor LSTs is quite small at high latitude area while it increases at low to mid-latitude area. The bias varies from 0.43 K to 0.92 K and the STD varies from 0.77 K to 2.61 K as shown in Table 9.

Discussion
Regarding the ground data validation, one of the most significant issues is observed over The error distributions over sensor zenith angle (STZ) indicate VIISR LST is consistently 0.7 K warmer than MODIS LST when the VIIRS STZ is below 50 degree and then decrease after that. The difference STD generally increases with the angle particularly when it is over 50 degree. The difference also presents characteristics of changes in latitude. As shown in Table 9, the difference of the two sensor LSTs is quite small at high latitude area while it increases at low to mid-latitude area. The bias varies from 0.43 K to 0.92 K and the STD varies from 0.77 K to 2.61 K as shown in Table 9.

Discussion
Regarding the ground data validation, one of the most significant issues is observed over Goodwin Creek site (GWN) in SURFRAD, which shows overestimation of 2.1 K at nighttime and underestimation of 2.75 K at daytime. The GWN site is located in a rural pasture grassland in the state of Mississippi. As shown in Figure 9, the area surrounding the site consists of a mixture of grass and forest. The spatial variation analysis of a 10 * 10 window of the ASTER GED v3 100 m emissivity (2000-2008) emissivity [27] around the site shows a STD of 0.003 with emissivity ranges from 0.962 to 0.978 around the site.
The lower emissivity appears at right side of the site. Such day night error difference was also reported in Li et al. (2014) [31] for the validation of the MODIS LST product using SURFRAD data, in which the 10-year average daytime LST shows underestimation of up to 7 K particularly during spring and summer, and nighttime LST shows about 2 K overestimation throughout a year. Zhou et al. (2019) [18] reported about 4 K underestimation of AQUA MODIS LST at daytime. Li et al. (2014) [31] pointed out that the sensor zenith angle has a strong correlation with daytime LST difference and wind speed has a strong correlation with nighttime LST difference. The nighttime warmer bias of 1.82 K was reported in AATSR LST validation against GWN observations [36]. It was pointed out that tree crown temperature is usually close to air temperature while open (short grass) areas experience stronger radiative cooling at nighttime. On the contrary, at daytime the short grass area temperature tends to be warmer than the tree crown temperature. In this study, we tried to use average of 5 * 5 satellite pixels LST to compare with the ground observations, but it is found that it did not significantly improve the bias, which indicates the average LST is still a mismatch with the in-situ sensor observation.  [31] pointed out that the sensor zenith angle has a strong correlation with daytime LST difference and wind speed has a strong correlation with nighttime LST difference. The nighttime warmer bias of 1.82 K was reported in AATSR LST validation against GWN observations [36]. It was pointed out that tree crown temperature is usually close to air temperature while open (short grass) areas experience stronger radiative cooling at nighttime. On the contrary, at daytime the short grass area temperature tends to be warmer than the tree crown temperature. In this study, we tried to use average of 5 * 5 satellite pixels LST to compare with the ground observations, but it is found that it did not significantly improve the bias, which indicates the average LST is still a mismatch with the in-situ sensor observation. Another issue is found over the DRA site in SURFRAD, which indicates the satellite LST is underestimated about 2 K for both nighttime and daytime as shown in Table 7. Li et al. (2014) [31] reported that MODIS LST is on average about 3 K lower than the measurements at nighttime and about 2 K lower than the station data at daytime over the DRA site. Zhou et al. (2019) [18] reported 2 K and 4 K underestimation of AQUA MODIS LST for daytime and nighttime, respectively. The DRA site is located about 65 miles northwest of Las Vegas, Nevada. As shown in Figure 10, the site consists of mostly open shrubland with some exposed bare soils. The emissivity map within 1 km of the site using ASTER GED data indicates that it is the most homogeneous site in terms of land cover among the seven sites in the SURFRAD network. The LST underestimation might be caused by the different shrubland cluster so as the shadow effect. Another issue is found over the DRA site in SURFRAD, which indicates the satellite LST is underestimated about 2 K for both nighttime and daytime as shown in Table 7. Li et al. (2014) [31] reported that MODIS LST is on average about 3 K lower than the measurements at nighttime and about 2 K lower than the station data at daytime over the DRA site. Zhou et al. (2019) [18] reported 2 K and 4 K underestimation of AQUA MODIS LST for daytime and nighttime, respectively. The DRA site is located about 65 miles northwest of Las Vegas, Nevada. As shown in Figure 10, the site consists of mostly open shrubland with some exposed bare soils. The emissivity map within 1 km of the site using ASTER GED data indicates that it is the most homogeneous site in terms of land cover among the seven sites in the SURFRAD network. The LST underestimation might be caused by the different shrubland cluster so as the shadow effect.

Conclusions
To produce consistent LST data record from geostationary and polar satellite missions, the enterprise LST algorithm has been proposed. The SW technique is employed for the enterprise LST algorithm derivation due to its robustness, computational efficiency and wide applications in the operational LST products, e.g., MODIS and SEVIRI. And channels availability makes the SW algorithm applicable to both VIIRS and ABI sensor. A comprehensive simulation database has been built for the LUT generation and theoretical analysis, which uses the MODTRAN radiative transfer model for the calculation of the top of atmosphere radiance. Three methods were used to evaluate the quality of the enterprise LST algorithm. Theoretical analysis using over 2000+ SEEBOR profiles show the bias of 0.19 K and 0.34 K; STD of 0.48 K and 0.69 K for nighttime and daytime, respectively.
NOAA 20 was successfully launched on 19 November, 2017 and its data has been available since 5 January, 2018. The ground measurements from seven sites of SURFRAD network and two sites of BSRN network were used for the validation of the enterprise N20 LST. The results show accuracy of -0.32 K and RMSE of 2.17 K over all 7 SURFRAD sites; accuracy of -0.24 K and 0.19 K, RMSE of 1.98 K and 2.03 K over BSRN GOB and CAB site, respectively. Nighttime outperforms daytime with corresponding precision of 1.91 K and 2.45 K for SURFRAD validation. For BSRN site, the nighttime precision is less than 1.5 K whereas the daytime precision is about 2.4 K. And LST discrepancy between the satellite LST and ground LST varies over sites. For SURFRAD stations, validation results suggest significant over/underestimation of satellite LST compared to the ground LST over GWN site and DRA site. Results over other sites indicate a close agreement between the satellite LST estimations and ground measurements with the precision from 1.45 K to 1.78 K. We should bear in mind that the ground validation quality is strongly subject to the spatial heterogeneity of surface temperature within the satellite field of view, which are not well presented by the ground in-situ observations. To well present the satellite LST, the site wide thermal homogeneity should be at a wider spatial area at least 2 km surrounding the site due to the viewing angular effect, which is extremely hard in reality to meet this requirement. Since the heterogeneity is more significant at daytime than at nighttime, many studies suggest using the nighttime only data to minimize the spatial heterogeneity. In this study, the accuracy for nighttime varies from 0 to 2.1 K (GWN site) and the precision from 1.2 K to 1.8 K for all 7 SURFRAD sites. For BSRN sites, the nighttime bias is -0.5 K and 0.78 K and RMSE is 1.34 K and 1.50 K for GOB and CAB site, respectively.
The cross-satellite comparison was performed between the enterprise N20 VIIRS LST and AQUA MODIS LST. The results indicate accuracy of 0.66 K and uncertainty of 1.94 K between N20 VIIRS LST and AQUA MODIS LST. Similar to the ground data validation, the cross comparison are affected by many factors e.g., cloud residue or cloud shadow, difference in field of view between sensors, algorithm difference, regional terrain, time of day, month of year etc. In this study, the comparison between N20 and AQUA has a wider spatial scope including US, Africa, Australia, North mid Asia and polar area and a wider temporal span from February to June.

Conclusions
To produce consistent LST data record from geostationary and polar satellite missions, the enterprise LST algorithm has been proposed. The SW technique is employed for the enterprise LST algorithm derivation due to its robustness, computational efficiency and wide applications in the operational LST products, e.g., MODIS and SEVIRI. And channels availability makes the SW algorithm applicable to both VIIRS and ABI sensor. A comprehensive simulation database has been built for the LUT generation and theoretical analysis, which uses the MODTRAN radiative transfer model for the calculation of the top of atmosphere radiance. Three methods were used to evaluate the quality of the enterprise LST algorithm. Theoretical analysis using over 2000+ SEEBOR profiles show the bias of 0.19 K and 0.34 K; STD of 0.48 K and 0.69 K for nighttime and daytime, respectively.
NOAA 20 was successfully launched on 19 November, 2017 and its data has been available since 5 January, 2018. The ground measurements from seven sites of SURFRAD network and two sites of BSRN network were used for the validation of the enterprise N20 LST. The results show accuracy of −0.32 K and RMSE of 2.17 K over all 7 SURFRAD sites; accuracy of −0.24 K and 0.19 K, RMSE of 1.98 K and 2.03 K over BSRN GOB and CAB site, respectively. Nighttime outperforms daytime with corresponding precision of 1.91 K and 2.45 K for SURFRAD validation. For BSRN site, the nighttime precision is less than 1.5 K whereas the daytime precision is about 2.4 K. And LST discrepancy between the satellite LST and ground LST varies over sites. For SURFRAD stations, validation results suggest significant over/underestimation of satellite LST compared to the ground LST over GWN site and DRA site. Results over other sites indicate a close agreement between the satellite LST estimations and ground measurements with the precision from 1.45 K to 1.78 K. We should bear in mind that the ground validation quality is strongly subject to the spatial heterogeneity of surface temperature within the satellite field of view, which are not well presented by the ground in-situ observations. To well present the satellite LST, the site wide thermal homogeneity should be at a wider spatial area at least 2 km surrounding the site due to the viewing angular effect, which is extremely hard in reality to meet this requirement. Since the heterogeneity is more significant at daytime than at nighttime, many studies suggest using the nighttime only data to minimize the spatial heterogeneity. In this study, the accuracy for nighttime varies from 0 to 2.1 K (GWN site) and the precision from 1.2 K to 1.8 K for all 7 SURFRAD sites. For BSRN sites, the nighttime bias is −0.5 K and 0.78 K and RMSE is 1.34 K and 1.50 K for GOB and CAB site, respectively.
The cross-satellite comparison was performed between the enterprise N20 VIIRS LST and AQUA MODIS LST. The results indicate accuracy of 0.66 K and uncertainty of 1.94 K between N20 VIIRS LST and AQUA MODIS LST. Similar to the ground data validation, the cross comparison are affected by many factors e.g., cloud residue or cloud shadow, difference in field of view between sensors, algorithm difference, regional terrain, time of day, month of year etc. In this study, the comparison between N20 and AQUA has a wider spatial scope including US, Africa, Australia, North mid Asia and polar area and a wider temporal span from February to June.
The enterprise LST algorithm has been running in the operational environment since June, 2019 and the data has been available at NOAA 's the Comprehensive Large Array-data Stewardship System (CLASS). Its application to GOES 16 ABI is ongoing.