AN IMPROVED TEMPERATURE SPATIAL INTERPOLATION METHOD FOR SPACEBORNE LIDAR ATMOSPHERIC CORRECTION

As one of the most important meteorological elements, temperature is an indispensable meteorological parameter for the atmospheric correction of spaceborne LiDAR ranging. Given a limited number of surface meteorological observation stations, the temperature values for all region of LiDAR observation need to be interpolated using appropriate spatial interpolation methods. In this paper, based on the monthly surface observation values in individual years (1981-2010) of Sichuan province observation stations, we firstly analyze the effects of three common interpolation methods, including inverse distance weighting (IDW), ordinary kriging (OK) and gradient plus inverse distance squared (GIDS). To solve the problem of low interpolation accuracy in severely undulating terrain area, an improved gradient distance inverse square method based on the adiabatic lapse rate (GIDS-ALR) is proposed. The experimental results show that the GIDS-ALR has an obvious improvement in the effect of severely undulating terrain, where the absolute error has been improved by more than 43% in average. Additionally, the temperature-interpolated MAE is reduced by more than 30%. The effectiveness and applicability of the proposed method is verified in this paper. * Corresponding author


INTRODUCTION
With the development of remote sensing technology, spaceborne LiDAR has becoming one of the fastest and most accurate remote sensing methods for obtaining surface elevation information (Winker et al., 2003;Zhang et al., 2017). The optical path difference generated by the laser pulse passing through the atmosphere is one of the error sources that cannot be ignored in LiDAR ranging (Schutz et al., 2002;Tang et al., 2016). To perform atmospheric correction on the LiDAR ranging value, it is necessary to obtain the atmospheric parameters such as temperature, humidity, and pressure at the laser foot point (Niell et al., 1996;Xin et al., 2011). The temperature is one of the indispensable meteorological parameters for atmospheric correction of spaceborne LiDAR ranging. However, sometimes we cannot directly obtain the atmospheric parameters in the region of LiDAR observation, so the method of spatial interpolation needs to be used.
In view the fact that elevation and topography have the most significant influence on the spatial interpolation of temperature, scholars have done much research on the influence of elevation and terrain on air temperature interpolation (Willmott et al., 1995;Ma et al., 2015;Stewart et al., 2017). Taking into account the complexity of the algorithm and the accuracy of the interpolation results, the Inverse Distance Weighting (IDW), Ordinary Kriging (OK) and Gradient plus Inverse Distance Squared (GIDS) are the three most commonly used methods Li et al., 2005;Li et al., 2011;Li et al., 2014;Mardikis et al., 2005). This paper makes a detailed comparison and analysis of the above three method effects using the monthly surface observation values in individual years  of Sichuan province observation stations. Among these three methods, GIDS performed better than the other two. In order to solve the problem of relatively poor interpolation effect in greatly undulate terrain, we proposed an improved approach--Gradient plus Inverse Distance Squared based on Adiabatic Lapse Rate (GIDS-ALR), in which the multiple linear regression model and temperature adiabatic lapse rate are comprehensively considered and employed to acquire the regression coefficient. The experimental results show an obvious improvement, with the absolute error and the MAE of temperature-interpolated, respectively reduced by more than 43% and 30% in the severely undulating terrain area.

Study Area and Data
To thoroughly compare and analyze the experimental results of interpolation methods, Sichuan Province is preferred as an experimental area, taking the distribution of surface meteorological stations and terrain fluctuations into account. Sichuan Province locates in the transitional zone between the Qinghai-Tibet Plateau and the Yangtze River Plain, where its hypsography has great disparity and topography is complex and diverse. The elevation of western plateau and mountainous areas is more than 3,000 meters, while that of the eastern basin and hills is mostly between 500 and 2000 meters.
The experimental data includes: the 3D coordinates, temperature values and DEM data of the 219 surface meteorological stations in Sichuan Province and its neighbouring regions. The temperature data, China's Monthly Surface Observation Values in Individual Years , is downloaded from the China Meteorological Data Network (http://data.cma.cn/site/index.html).
The DEM data, SRTMDEMUTM 90M resolution digital elevation data products in Sichuan Province and its neighbouring regions, is sourced from Geospatial Data Cloud (http://www.gscloud.cn/). The distribution of meteorological stations is shown in Figure 1, which is coloured by elevation, with circles representing sample stations and triangles test stations.

Inverse Distance Weighting:
IDW is one of the most commonly used and simplest methods of spatial interpolation which is often used for temperature interpolation (Lu et al.,2008;Chen et al., 2012). It uses a weighted average of the temperature values from nearby sample stations to estimate the magnitude of that attribute at target locations. The weight of a specific sample station is assigned based on the distances between the target location and the station. For a specific target location 0 p , we have: , i d means the distance between 0 p and i p , k is the power of distance which generally ranges from 1~4( k is set to 2 in this paper), and n is the number of sample stations.

Ordinary Kriging:
OK is the most frequently used method in the group of Kriging interpolation techniques (Holdaway et al.,1996;Jarvis et al., 2001). Based on the original data of regionalized variables and the structural characteristics of the variation function, we get the optimal linear unbiased estimation of regionalized variable of target locations by two steps: firstly, variation function model is established through analyzing the structure of space random field thoroughly; secondly, the Kriging equations are built and solved on the basis of the variation function model. The estimated value can be calculated by the following equitation: where  i refers to the weight of the value of sample station i p , with the condition that the weights must add up to 1: which aims at making the mathematical expectation of the estimated variable equal to the estimated parameter and minimizing the variance of the variable. . It is a generalized IDW method, additionally considering the relations between the meteorological element and latitude, longitude, elevation. The element's estimated value can be calculated by:

Gradient plus Inverse Distance Squared: GIDS was first proposed by
where ( , , ) x y z denote the latitude, longitude and elevation of the target position, and ( , , ) is established and the least squares estimation (LSM) method is used to solved these three coefficients. The regression model is as follows: Tp means the estimated value of meteorological element of sample station, and 0 C is regression constant.

Validation Method
In this paper, the "actual verification" method is adopted to conduct performance confirmation (Shu-mm et al., 2005). We uniformly selected 24 representative stations as test stations from a total of 219 meteorological stations, while the rest as sample stations. Mean Absolute Error (MAE) and Root Mean Squared Error (RMSE) were used as evaluation criteria to analyze and evaluate the accuracy of different temperature interpolation methods in different terrain conditions. The MAE and RMSE can be expressed as following:

Comparison of Interpolation Results
In order to compare the accuracy of these three temperature interpolation methods including IDW, GIDS, and OK, 24 meteorological stations are chosen for test in Sichuan Province according to different terrain characteristic. The terrain types of the 24 test stations cover the plateau in western Sichuan, the mountain in the middle of Sichuan, the basin in East Sichuan and so on. Based on the validation method described above, the monthly surface observation values in individual years  of test stations are estimated, and the MAE and the RMSE are calculated as well.

Figure 2. Temperature interpolation results
The interpolated temperature curve of each station is shown in Figure 2, in which January, April, July and October are selected as the representative months of the four seasons respectively. It is obvious that the interpolation results of the GIDS method are much closer to the actual temperature, implying a better interpolation effect, compared with the IDW and OK method. Most of the interpolation errors of GIDS are equivalent or less than that of IDW and OK, except for that of Songpan and Muli stations in January. The MAE and RMSE of the temperature using the three interpolation methods are given in Table 1

Improvement of GIDS Method
In this paper, the absolute error curves of 24 test stations under different terrain environments in 12 months are studied. It is found that the 12-month absolute error trend is similar. In order to fully explain the experimental results of the GIDS method with the best interpolation effect, this paper takes the errors curve of February with the largest error as an example for analysis.
The interpolation temperature absolute error curve of IDW, OK and GIDS in February are shown in Fig.3. As can be seen from the figure, in the stations located in the Sichuan basin region with less topographic relief, the interpolation results are of high accuracy with the interpolation errors controlled within 1.5°C. However, the best method GIDS (See Section 3.2) still has a relatively large error from Shiqu to Miyi stations, which are located in severely undulating terrain areas. The main reason can be explained by the fact that the spatial interpolation methods above applied the Tobler's First Law of Geography, that is, the meteorological elements are spatial autocorrelation.
However, topography will inevitably reduce the correlation when it comes to the mountainous regions. In addition, the GIDS method treats the effects of x, y and z on the temperature equally, conflict with the fact that the variation of temperature in the horizontal and vertical directions is different.

Figure 3 Absolute errors of IDW, OK and GIDS
Considering the contradiction and the irrationality of interpolation, the GIDS based on the Adiabatic Lapse Rate (GIDS-ALR) is presented in this paper. This improved method considers the effects of elevation and latitude and longitude separately, and accurately calculates the partial regression coefficients x C , y C , and z C , that describe the relationship between the temperature and the latitude, longitude, and elevation at the point to be interpolated. The detailed steps of this method are as follows: Step 1: Determine z C . Taking into account the influence law of elevation on temperature that the partial regression coefficient of temperature and elevation z C is close to the adiabatic lapse rate g T in the same region, the empirical value of g T in corresponding position is used as z C in each sample station in this paper, i.e. = z g CT , where the adiabatic lapse rate g T is commonly taken to be a constant. In this paper, empirical values 0.0065 and 0.0060℃/m are tested separately (Jiang et al., 2016;Sun Yat-Sen University et al., 1978).
Step 2: Solve x C , y C . Since z C has been considered in step 1, a multivariate linear regression model for calculating x C and y C was established: where ˆ( ) mi Tp denotes the estimated value of modified temperature by g T to eliminate the impact of elevation, concerning the model concentrate on the effect of x C and y C . The calculation formula of () mi Tpis given: Step 3: Put x C , y C and z C obtained above into equation (4) to get the estimated value we need.
To compares the interpolation results of GIDS and GIDS-ALR methods, experimental analysis is conducted based on different terrain features. The absolute error curve of the two methods is shown in Figure 4.
It can be seen that as to Songpan station with the largest error by GIDS, after applying the GIDS-ALR method proposed in this paper, the error is reduced from 6.5720°C to 1.9407°C, and the absolute error improved by 70.47%. In addition, at the test stations Shiqu, Litang, Daocheng and Miyi, the errors by GIDS-ALR are also significantly smaller than those by GIDS, with the errors decreased by 43.62% and 46.68% when the adiabatic lapse rate g T valuing -0.0060 and -0.0065, respectively. The reason for the error reduction can be analysed from the elevation difference between the test station and the corresponding sample stations. The elevation of the Songpan station is 2852.1 meters, and the average elevation of the sample stations is 870.5455 meters which is less than that of the Songpan station. The differences of the elevation between the other test stations and the average elevation of their corresponding sample stations are between 300 meters and 1112 meters. In these regions where the topography is undulating, the GIDS-ALR method can achieve higher accuracy, by considering the influence of elevation and latitude, longitude separately, and taking into account empirical adiabatic lapse rate simultaneously. For the 15 test stations from Shiqu to Jinyang with relatively violent topographic fluctuations, when the adiabatic lapse rate g T values -0.0060 and -0.0065, the MAE of the interpolated temperature obtained by the GIDS-ALR method is reduced by 33.37% and 31% compared with GIDS. Summing up, the GIDS-ALR method proposed in this paper can overcome the problem of larger interpolation error in the severely undulating terrain effectively. And results in the gentle terrain area can also be slightly improved. Due to the same experiential value of adiabatic lapse rate applied to different sample stations, the temperature values of sample stations cannot be corrected very well, which has certain effect on the interpolation results, and is also to be improved.

CONCLUSIONS
Facing the demand of spatial interpolation of temperature needed during the correction of LiDAR ranging, this paper firstly compares and analyzes three interpolation methods of IDW, OK and GIDS. The results show the MAE and RMSE of GIDS are reduced by over 50% than those of IDW and OK, showing that GIDS method works best. In addition, the interpolated temperatures using the three methods in severely undulating terrain areas are all deviated greatly from the real value. According to the results of comparative analysis above, aiming at the problem of poor effect in the dramatic terrain fluctuations, the GIDS-ALR method is proposed. This improved method considers the effects of elevation and latitude and longitude separately making it calculating the partial regression coefficients more accurately. Compared with GIDS, the MAE of interpolated temperature calculated by the GIDS-ALR is much lower, which effectively verify the validity of the proposed method. Subsequent to this method, the algorithm can also be optimized in the parameter settings of adiabatic lapse rate to improve performance. In the future, more validation experiments should be taken to verify the effectiveness and applicability of the GIDS-ALR. We have so far do preliminary test in some other areas like Xinjiang Province and Hainan Province, in which the results are consistence with the conclusion in this paper. Furthermore, relevant depth experiments and analysis need to done to improve and perfect the method.