Geometric Accuracy Improvement Method for High-Resolution Optical Satellite Remote Sensing Imagery Combining Multi-Temporal SAR Imagery and GLAS Data

Abstract

With the widespread availability of satellite data, a single region can be described using multi-source and multi-temporal remote sensing data, such as high-resolution (HR) optical imagery, synthetic aperture radar (SAR) imagery, and space-borne laser altimetry data. These have become the main source of data for geopositioning. However, due to the limitation of the direct geometric accuracy of HR optical imagery and the effect of the small intersection angle of HR optical imagery in stereo pair orientation, the geometric accuracy of HR optical imagery cannot meet the requirements for geopositioning without ground control points (GCPs), especially in uninhabited areas, such as forests, plateaus, or deserts. Without satellite attitude error, SAR usually provides higher geometric accuracy than optical satellites. Space-borne laser altimetry technology can collect global laser footprints with high altitude accuracy. Therefore, this paper presents a geometric accuracy improvement method for HR optical satellite remote sensing imagery combining multi-temporal SAR Imagery and GLAS data without GCPs. Based on the imaging mechanism, the differences in the weight matrix determination of the HR optical imagery and SAR imagery were analyzed. The laser altimetry data with high altitude accuracy were selected and applied as height control point in combined geopositioning. To validate the combined geopositioning approach, GaoFen2 (GF2) optical imagery, GaoFen6 (GF6) optical imagery, GaoFen3 (GF3) SAR imagery, and the Geoscience Laser Altimeter System (GLAS) footprint were tested. The experimental results show that the proposed model can be effectively applied to combined geopositioning to improve the geometric accuracy of HR optical imagery. Moreover, we found that the distribution and weight matrix determination of SAR images and the distribution of GLAS footprints are the crucial factors influencing geometric accuracy. Combined geopositioning using multi-source remote sensing data can achieve a plane accuracy of 1.587 m and an altitude accuracy of 1.985 m, which is similar to the geometric accuracy of geopositioning of GF2 with GCPs.

1. Introduction

High-resolution (HR) optical satellite imagery has become an important data source for geopositioning due to the wide coverage, short return period, and good image interpretation [1]. Due to the limitation of the direct geometric accuracy of HR optical imagery, some ground control points (GCPs) are usually adopted to compensate for geometric errors in HR optical imagery to meet the geometric accuracy requirements of geopositioning [2,3,4,5,6]. However, the lack of access to some uninhabited areas, such as forests, plateaus, or deserts, prevents the acquisition of GCPs, so the plane accuracy cannot meet the requirements for geopositioning. HR optical satellite images are usually acquired with nadir or near-nadir imaging geometry. Therefore, affected by the small intersection angle of HR optical imagery in stereo pair orientation, the altitude accuracy also cannot meet the geopositioning requirements.

Given the widespread availability of satellite data, a single region can be described by multi-source and multi-temporal remote sensing data, such as SAR imagery and laser altimetry data. Without satellite attitude error, SAR usually provides higher geometric accuracy than optical satellites [7]. However, SAR satellite imagery has more noise than HR optical imagery [8]. SAR satellite images usually have a small coverage area and cannot meet the needs of large-area geopositioning. Space-borne laser altimetry technology can collect global laser footprints with high altitude accuracy [9]. However, the laser footprints of current space-borne laser altimetry technology are more discrete and cannot express the detailed information of the ground surface. Therefore, combined geopositioning using multi-source and multi-temporal remote sensing data could be applied to improve the geometric accuracy of optical imagery without GCPs.

Several studies have focused on geometric accuracy improvement methods for HR optical satellite remote sensing imagery using multi-source and multi-temporal remote sensing data. Jeong and Kim studied the geometric accuracy and imaging geometry of dual-sensor stereo pairs in comparison with those of conventional single-sensor stereo pairs. Then, they investigated the influences of convergence angle, bisector elevation angle, and asymmetry angle on geometric accuracy for dual-sensor stereo pair orientation [10]. Jeong et al. investigated the positioning accuracy of multiple-satellite images. Using IKONOS, QuickBird, and KOMPSAT-2 data for Daejeon, Korea, they revealed the potential, limitations, and important considerations for geopositioning applications using images from multiple satellites [11]. Based on the rigorous combined geopositioning model for optical and SAR imagery, Xing et al. discussed the construction of an observation equation, computation of the weight matrix, and the frame of computation process, and further analyzed geometric accuracy with GCPs [12]. Zhang et al. validated the application of the universal rational function model (RFM) for optical and SAR stereo pair combined orientation [13]. Cheng et al. analyzed the influences of the number of SAR images, orbit direction of the SAR satellite, and the distribution of SAR images on geometric accuracy using the RFM model [14]. The first Earth observation laser altimetry satellite, Ice Cloud and Land Elevation Satellite (ICESat-1), carrying the GLAS, was launched in 2003 [15,16,17,18]. With the advantages of high altitude accuracy, global distribution, and a large number of collected footprints, GLAS footprints are widely used for height reference in the calibration and validation of digital elevation models (DEMs) [19]. Some researchers have employed the GLAS footprint as height control points in optical imagery block adjustment, which significantly improved altitude accuracy [20]. Li et al. added two GLAS height control points for every ZiYuan-3 image in the block area, and the vertical accuracy satisfied the requirements doe 1:50,000 scale mapping [21]. However, few studies have been dedicated to the study of geometric accuracy improvement of HR optical images combined with SAR images and space-borne laser altimetry footprints, which could play an important role in geopositioning without GCPs. Therefore, analyzing the key factors influencing the geometric accuracy of combined geopositioning using multi-source remote sensing data without GCPs is required, which is important for high accuracy geopositioning in uninhabited areas.

This paper presents a geometric accuracy improvement method of HR optical satellite remote sensing imagery combining multi-temporal SAR Imagery and GLAS data. Based on error analysis of HR optical and SAR images, we prove that the affine transformation model in image space can effectively compensate for the unmodelled errors in the combined geopositioning model. Considering the imaging mechanisms of HR optical and SAR satellites, the weight determination of the combined geopositioning model is proposed in this paper. Based on the characteristics of GLAS products, the footprints located in flat areas with high altitude accuracy were selected as height control points for HR optical imagery combined geopositioning. Then, the observation equation was constructed using the GLAS control points. In our experiments, GF2 HR optical images, GF6 HR optical images, GF3 SAR images, and the GLAS footprint covering Dongying area in China were tested. The corresponding digital orthophotograph model (DOM) and DEM were applied as reference data to evaluate plane accuracy and altitude accuracy, respectively. The experimental results showed that the proposed model can be adopted for combined multi-source and multi-temporal remote sensing data geopositioning to improve geometric accuracy of HR optical satellite remote sensing imagery without GCPs. The distribution and weight matrix settings for SAR images and the number of GLAS footprints were found to be the crucial factors influencing the geometric accuracy. For GF2, multi-source combined geopositioning achieved a plane accuracy of 1.587 m and an altitude accuracy of 1.985 m, which is similar to the geometric accuracy of geopositioning of GF2 with GCPs. The geometric accuracies for GF6 multi-source combined geopositioning were a plane accuracy of 2.908 m and an altitude accuracy of 7.669 m, which is similar to the geometric accuracy of geopositioning of GF6 with GCPs.

The remainder of the paper is organized as follows: Section 2 demonstrates the methodology used to improve the geometric accuracy of HR optical satellite remote sensing imagery with the combined geopositioning method. Section 3 describes the experimental area and data and presents the experimental results and discussion. Section 4 summarizes our work and presents the study conclusions and research perspectives.

2. Methodology

2.1. Rational Polynomial Coefficients Model for Combined Geopositioning

The RFM model, in the form of rational polynomial coefficients (RPCs), connects the object coordinate $(B,L,H)$ with the image coordinate $(l,s)$ [22,23,24,25,26,27,28]. $B,L,$ and H are the longitude, latitude, and ellipsoid height, respectively; l and s are the line and sample coordinates, respectively. This can be expressed by the following equation:

$$\left\{\begin{array}{c}{\displaystyle l_n=\frac{Nu{m}_l(U,V,W)}{De{n}_l(U,V,W)}}\hfill \\ {\displaystyle s_n=\frac{Nu{m}_s(U,V,W)}{De{n}_s(U,V,W)}}\hfill \end{array}\right.$$

(1)

$$\begin{array}{c}\begin{array}{c}\hfill Nu{m}_l(U,V,W)=a_1+a_{2}V+a_{3}U+a_{4}W+a_{5}VU+a_{6}VW+a_{7}UW+a_8{V}^2+a_9{U}^2+a_{10}{W}^2+a_{11}UVW\\ \hfill +a_{12}{V}^3+a_{13}V{U}^2+a_{14}V{W}^2+a_{15}{V}^{2}U+a_{16}{U}^3+a_{17}U{W}^2+a_{18}{V}^{2}W+a_{19}{U}^{2}W+a_{20}{W}^3\end{array}\hfill \\ \begin{array}{c}\hfill De{n}_l(U,V,W)=b_1+b_{2}V+b_{3}U+b_{4}W+b_{5}VU+b_{6}VW+b_{7}UW+b_8{V}^2+b_9{U}^2+b_{10}{W}^2+b_{11}UVW\\ \hfill +b_{12}{V}^3+b_{13}V{U}^2+b_{14}V{W}^2+b_{15}{V}^{2}U+b_{16}{U}^3+b_{17}U{W}^2+b_{18}{V}^{2}W+b_{19}{U}^{2}W+b_{20}{W}^3\end{array}\hfill \\ \begin{array}{c}\hfill Nu{m}_s(U,V,W)=c_1+c_{2}V+c_{3}U+c_{4}W+c_{5}VU+c_{6}VW+c_{7}UW+c_8{V}^2+c_9{U}^2+c_{10}{W}^2+c_{11}UVW\\ \hfill +c_{12}{V}^3+c_{13}V{U}^2+c_{14}V{W}^2+c_{15}{V}^{2}U+c_{16}{U}^3+c_{17}U{W}^2+c_{18}{V}^{2}W+c_{19}{U}^{2}W+c_{20}{W}^3\end{array}\hfill \\ \begin{array}{c}\hfill De{n}_s(U,V,W)=d_1+d_{2}V+d_{3}U+d_{4}W+d_{5}VU+d_{6}VW+d_{7}UW+d_8{V}^2+d_9{U}^2+d_{10}{W}^2+d_{11}UVW\\ \hfill +d_{12}{V}^3+d_{13}V{U}^2+d_{14}V{W}^2+d_{15}{V}^{2}U+d_{16}{U}^3+d_{17}U{W}^2+d_{18}{V}^{2}W+d_{19}{U}^{2}W+d_{20}{W}^3\end{array}\hfill \end{array}$$

(2)

where $a_i,b_i,c_i,d_i(i=1,2,\dots ,20)$ are the rational polynomial coefficients, and the values of $b_1,d_1$ are usually set to 1. $(l_n,s_n)$, $(U,V,W)$ are the normalized coordinates of the image coordinate $(l,s)$ and the object coordinate $(Lat,Lon,Height)$, respectively, which ensures the stability of calculation of the RFM model. According to the following equations, the image coordinate and the object coordinate are normalized between –1 and 1.

$$\left\{\begin{array}{c}{\displaystyle l_n=\frac{l-LineOff}{LineScale}}\hfill \\ {\displaystyle s_n=\frac{s-SampleOff}{SampleScale}}\hfill \end{array}\right.$$

(3)

where $LineOff$ and $SampleOff$ are the translation value of the image coordinate. $LineScale$ and $SampleScale$ are the scale value of the image coordinate.

$$\left\{\begin{array}{c}{\displaystyle U=\frac{Lon-LonOff}{LonScale}}\hfill \\ {\displaystyle V=\frac{Lat-LatOff}{LatScale}}\hfill \\ {\displaystyle W=\frac{Height-HeiOff}{HeiScale}}\hfill \end{array}\right.$$

(4)

where $LonOff$, $LatOff$, and $HeiOff$ are the translation value of the object coordinate. $LonScale$, $LatScale$, and $HeiScale$ are the scale value of the object coordinate.

With the advantage of sensor independence, the RPC model can be used for both HR optical imagery and SAR imagery. Therefore, the RPC model was adopted here for combined geopositioning.

2.2. Orientation Error

Some errors occur during satellite data acquisition and processing, such as sensor errors, platform ephemeris errors, attitude measurement errors, and so on, causing offset errors, scale errors, and rotation errors in remote sensing imagery, as listed in

Table 1. Geopositioning error sources of optical satellite and synthetic aperture radar (SAR) satellite.

Error	Error Source of Optical Satellite	Error Source of SAR Satellite
Along-track offset error	Time error, along-track position error, pitch angle error, along-track detector errors	Drift error of spacecraft time, along-track position error, radial position error, velocity error
Along-track linear or scale error	Lens distortion, focal length error, position error along the sub-satellite direction	Relative error of stable local oscillator, along-track velocity error
Cross-track offset error	Cross-track position error, roll angle error, cross-track detector errors, atmosphere refraction error	Delay error between signal transmitter and signal receiver, cross-track position error, radial position error, atmospheric delay
Cross-track linear or scale error	Cross-track lens distortion, focal length error, position error along the sub-satellite direction	Incident angle error, radial position error, atmospheric delay

.

The error sources of optical satellite and SAR satellite cause offset errors and linear or scale errors along or across the track direction. Therefore, to achieve perfect space intersection, the affine transformation model in image space is applied to compensate for existing geometric errors. The RPC model for combined geopositioning is as follows:

$$\left\{\begin{array}{c}l+a_0+a_1\times l+a_2\times s=F_l(Lat,Lon,Height)\hfill \\ s+b_0+b_1\times l+b_2\times s=F_s(Lat,Lon,Height)\hfill \end{array}\right.$$

(5)

where $a_0$ and $b_0$ compensate for the offset error and $a_1,a_2,b_1$, and $b_2$ compensate for the scale and rotation error.

Without GCPs, a data defect exists in the combined geopositioning model, which results in equation divergence [29]. Therefore, a pseudo-observation equation of the image affine transformation parameters was employed here to achieve combined geopositioning without GCPs. The integrated model for combined geopositioning using multi-source remote sensing data can be represented in matrix form as:

$$\begin{array}{cc}{V}_1=B_G{X}_G+B_A{X}_A-L_1\hfill & P_1\hfill \\ V_2=I{X}_A-L_2\hfill & P_2\hfill \end{array}$$

(6)

$$\begin{array}{c}{B}_G=\left[\begin{array}{ccc}{\displaystyle \frac{\partial l}{\partial Lat}}\hfill & {\displaystyle \frac{\partial l}{\partial Lon}}\hfill & {\displaystyle \frac{\partial l}{\partial Height}}\hfill \\ {\displaystyle \frac{\partial s}{\partial Lat}}\hfill & {\displaystyle \frac{\partial s}{\partial Lon}}\hfill & {\displaystyle \frac{\partial s}{\partial Height}}\hfill \end{array}\right],X_G=\left[\begin{array}{c}▵Lat\hfill \\ ▵Lon\hfill \\ ▵Height\hfill \end{array}\right]\hfill \end{array}$$

(7)

$$\begin{array}{c}{B}_A=\left[\begin{array}{cccccc}{\displaystyle \frac{\partial l}{\partial a_0}}\hfill & {\displaystyle \frac{\partial l}{\partial a_1}}\hfill & {\displaystyle \frac{\partial l}{\partial a_2}}\hfill & {\displaystyle \frac{\partial l}{\partial b_0}}\hfill & {\displaystyle \frac{\partial l}{\partial b_1}}\hfill & {\displaystyle \frac{\partial l}{\partial b_2}}\hfill \\ {\displaystyle \frac{\partial s}{\partial a_0}}\hfill & {\displaystyle \frac{\partial s}{\partial a_1}}\hfill & {\displaystyle \frac{\partial s}{\partial a_2}}\hfill & {\displaystyle \frac{\partial s}{\partial b_0}}\hfill & {\displaystyle \frac{\partial s}{\partial b_1}}\hfill & {\displaystyle \frac{\partial s}{\partial b_2}}\hfill \end{array}\right],X_A=\left[\begin{array}{c}▵a_0\hfill \\ ▵a_1\hfill \\ ▵a_2\hfill \\ ▵b_1\hfill \\ ▵b_2\hfill \\ ▵b_3\hfill \end{array}\right]\hfill \end{array}$$

(8)

where the first observation equation is relevant to the tie point and the second observation equation is the pseudo-observation equation. $X_G$ and $X_A$ are the object coordinates of tie points and the affine transformation parameters of images, respectively; $B_G$ and $B_A$ are coefficient matrixes derived from the derivatives of $X_G$ and $X_A$; I is the unit coefficient matrix; $L_1$ and $L_2$ are observations in the observation equations; and $P_1$ and $P_2$ are the weight matrixes. The weight $P_1$ relating to the tie points is set to one pixel because the measurement accuracy of image coordinates is usually achieved at one pixel. The weight $P_2$ relating to the affine transformation parameters is set according to the physical meaning and priori error of each parameter, expressed as:

$$\begin{array}{c}\left\{\begin{array}{c}{\displaystyle P_{a_0}=\frac{1}{{{\sigma }_{a_0}}^2}}\hfill \\ {\displaystyle P_{b_0}=\frac{1}{{{\sigma }_{b_0}}^2}}\hfill \end{array}\right.\left\{\begin{array}{c}{\displaystyle P_{a_1}=P_{b_1}=\frac{1}{{{\sigma }_1}^2}}\hfill \\ {\displaystyle P_{a_2}=P_{b_2}=\frac{1}{{{\sigma }_2}^2}}\hfill \end{array}\right.\hfill \end{array}$$

(9)

where ${\sigma }_{a_0}$ and ${\sigma }_{b_0}$ are the offset errors along the line direction (row direction) and sample direction (column direction), respectively; ${\sigma }_1$ and ${\sigma }_2$ are the scale and rotation errors along the line direction and sample direction, respectively.

2.3. Combined Geopositioning with SAR

The relationships between object distance and image distance are different between HR optical imagery and SAR imagery due to the discrepancy in the imaging mechanisms. HR optical imagery applies central projection, as shown in Figure 1a. Therefore, the relationship between the object distance $△D$ and image distance $△d$ can be represented as follows:

$$△d=\frac{△D}{R_{opt}}$$

(10)

where $R_{opt}$ is the spatial resolution of HR optical imagery.

However, SAR imagery applies slant range projection in the cross-track direction. The sample distance in SAR imagery reflects the slant range distance to the object, not the cross-track distance, as shown in Figure 1b. Therefore, the relationship between the cross-track distance $△Y$ and sample distance $△s$ in the image can be expressed as follows:

$$△s=\frac{△Y\times sin{\alpha }_{inc}}{R_{SAR\_sample}}$$

(11)

where ${\alpha }_{inc}$ is the incident angle and $R_{SAR\_sample}$ is the resolution along the sample direction in SAR images.

Therefore, the weight determination for the offset error compensation parameters differs. Based on the a priori geometric accuracy, offset error along the line direction and sample direction for HR optical imagery and SAR imagery can be calculated as follows:

$$\begin{array}{c}\left\{{\begin{array}{c}{\displaystyle {\sigma }_{a_0}=\frac{{\sigma }_X}{R_{opt}}}\hfill \\ {\displaystyle {\sigma }_{b_0}=\frac{{\sigma }_Y}{R_{opt}}}\hfill \end{array}}_{opt}\right.\left\{{\begin{array}{c}{\displaystyle {\sigma }_{a_0}=\frac{{\sigma }_X}{R_{SAR\_line}}}\hfill \\ {\displaystyle {\sigma }_{b_0}=\frac{{\sigma }_Y\times sin{\alpha }_{inc}}{R_{SAR\_sample}}}\hfill \end{array}}_{SAR}\right.\end{array}$$

(12)

where ${\sigma }_X$ and ${\sigma }_Y$ are the along-track and cross-track geometric accuracy, respectively, and are usually calculated by the priori geometric accuracy $\sigma$, i.e., ${\sigma }_X={\sigma }_Y=\sigma /\sqrt{2}$. $R_{SAR\_line}$ is the resolution along the line direction in SAR images.

The weight determination for the scale and rotation error compensation parameters is similar for both HR optical imagery and SAR imagery. They are usually determined by the maximum pixels M resulting from the scale, rotation, and size of the image [14], as follows:

$$\left\{\begin{array}{c}{\displaystyle {\sigma }_1=\frac{M}{H}}\hfill \\ {\displaystyle {\sigma }_2=\frac{M}{W}}\hfill \end{array}\right.$$

(13)

where H and W are the height and width of imagery, respectively.

2.4. Combined Geopositioning with GLAS

In flat areas, the altitude accuracy of the GLAS footprint can be as high as 15 cm [30,31,32,33,34,35]. However, due to the impact of outliers, clouds, slope, and vegetation, the altitude accuracy of the GLAS footprint is significantly lower, and cannot meet the altitude accuracy needs of geopositioning with high accuracy [36,37]. Therefore, to acquire GLAS footprints with high altitude accuracy, an extraction criterion including multiple constraints was applied. The second version of ASTER GDEM (GDEM2), attribute characteristic in the GLA14 product and waveform characteristic in the GLA01 product were used to select footprints with the thresholds defined in

Table 2. The criterion for Geoscience Laser Altimeter System (GLAS) footprint extraction.

Characteristic		Condition
Altitude difference compared with GDEM2		<20
GLA01 product	Number of peaks of received waveform	=1
	Sigma of Gaussian of received waveform	≤3.2
GLA14 product	Attitude quality indicator	=0
	Surface reflectance parameter	<1
	Saturation correction flag	=0

. Using GDEM2, the footprints with mistaken height values can be discarded. The GLA01 product is an altimetry data product that contains the received waveforms information required to produce higher accuracy range and elevation products. The number of peaks and the width of the waveform can reflect the characteristics of the topography and features. Given the constraints of waveforms, the footprints located in flat areas can be extracted. The GLA14 product is the land products file that contains the parameters of observation conditions and processing quality of the footprint. With the constraints of GLA14, footprints with high altitude accuracy can be extracted.

Based on the object coordinate of the GLAS footprint and RPC model, the image point of the footprint can be calculated by applying the back-projection principle. However, due to geometric error in HR optical imagery, a gap exists in the footprint location between the reference DOM and the HR optical imagery as shown by the green cross in Figure 2a,b. These image points within the same footprint are not homologous points located in the same texture area in HR optical and SAR images, as shown by the green cross in Figure 2b,c. Therefore, due to the high geometric accuracy of SAR, the image point of the footprint is first back-projected to the SAR imagery using the RPC model to ensure correction of image coordinates. The corresponding image point in HR optical imagery is then obtained manually. Unfortunately, many footprints are located in poor texture areas, resulting in difficulties for corresponding image point selection. The extracted footprint located in flat terrain and reflect the average height of a circle with a diameter of 70 m. Therefore, an obvious feature point, such as road intersection, within a range of 35 m of the footprint in SAR image is selected as the final image point of the footprint, as indicated by the red cross in Figure 2c.

The GLAS footprints are used as the height control point, so that only the plane coordinate should be solved in the combined geopositioning model, as shown in the following equation:

$$\begin{array}{c}{B}_G=\left[\begin{array}{cc}{\displaystyle \frac{\partial l}{\partial Lat}}\hfill & {\displaystyle \frac{\partial l}{\partial Lon}}\hfill \\ {\displaystyle \frac{\partial s}{\partial Lat}}\hfill & {\displaystyle \frac{\partial s}{\partial Lon}}\hfill \end{array}\right],X_G=\left[\begin{array}{c}▵Lat\hfill \\ ▵Lon\hfill \end{array}\right]\hfill \end{array}$$

(14)

3. Experimental Results and Evaluation

3.1. Experimental Data

To validate the application of the combined geopositioning model, three GF2 HR optical images, two GF6 HR optical images, two GF3 SAR images, and nine GLAS footprints covering the Dongying region in China were acquired. The GF2 satellite is the Chinese high-resolution optical remote sensing satellite that was launched on 19 August 2014. It is equipped with a double-camera push broom imaging system with a high resolution of 0.81 m and swath width of 45 km [38]. The GF6 satellite is characterized by high resolution and wide coverage, and was launched on 2 June 2018 [39]. The GF3 satellite is the first C-band multi-polarization SAR satellite remote sensing satellite in China, which was launched on 10 August 2016 [40]. It can observe the Earth in 12 operation modes with the resolution from 1 to 500 m and a swath width from 5 to 650 km to meet diverse application needs. The priori geometric accuracies of GF2 imagery, GF6 imagery, and GF3 imagery were approximately 30 m, 16 m, and 5 m, respectively. Detailed information about the GF2 HR optical satellite, GF6 HR optical satellite and GF3 SAR satellite is provided in

Table 3. Satellite parameters of GF2, GF6, and GF3.

	GF2	GF6	GF3
Orbit altitude	631 km	644.5 km	Orbit altitude	506 km
Orbit type	SunSync	SunSync	Orbit type	SunSync
Life	5–8 years	8 years	Life	8 years
	B1: 450–520 nm	B1: 450–520 nm
	B2: 520–590 nm	B2: 520–600 nm
Spectral range	B3: 630–690 nm	B3: 630–690 nm	Band	C
	B4: 770–890 nm	B4: 760–900 nm
	Pan: 450–900 nm	Pan: 450–900 nm
Field of view (${}^{\circ }$)	2.1	8.6	Incident angle (${}^{\circ }$)	10–60
			Scanning angle (${}^{\circ }$)	±20
Resolution (m)	0.81	2	Imaging mode	Spotlight mode (1, 10 × 10)
Swath width (km)	45	90	(Resolution (m),	Hyperfine strip model (3, 30)
			swath width (km))	Fine strip mode 1 (5, 50)
				Fine strip mode 2 (10, 100)
				Standard strip model (25, 130)
				Narrow scan mode (50, 300)
				Wide scan mode (100, 500)
				Full polarization strip model 1 (8, 30)
				Full polarization strip model 2 (25, 40)
				Wave imaging mode (10, 5 × 5)
				Global observation model (500, 650)
				Extended incident angle model (25, 130 or 80)

.

When the resolution difference between the optical imagery and GF3 SAR imagery is large, tie points selection is difficult, which influences the geometric accuracy with combined data. Therefore, the GF3 SAR imagery in spotlight mode was applied in this study.

Table 4. Description of experimental data for combined geopositioning.

Image	Acquisition Date	Pass	Incident Angle (${}^{\circ }$)	Size (pixels)	Resolution (m)
GF2-1	2 February 2015	DEC	15.502	29,200 × 27,620	0.81 × 0.81
GF2-2	21 September 2015	DEC	−8.403	29,200 × 27,620	0.81 × 0.81
GF2-3	26 August 2016	DEC	−12.444	29,200 × 27,620	0.81 × 0.81
GF6-1	13 March 2019	DEC	−0.004	44,500 × 48,312	2 × 2
GF6-2	17 March 2019	DEC	−0.003	44,500 × 48,312	2 × 2
GF3-1	9 April 2017	DEC	43.796–44.342	30,642 × 13,544	0.34 × 0.56
GF3-2	16 April 2017	ASC	39.224–39.881	33,760 × 13,382	0.33 × 0.56
Reference data	Resolution (m)	Accuracy (m)	Range
DOM	0.271	0.251	Lat: 37.3408–37.5170
DEM	1.806	0.273	Lon: 118.3691–118.8056

lists the detailed information of the experimental data. Figure 3 shows the coverage area of the experimental data. Figure 4 depicts the overview images of the experimental data.

In the overlap area, tie points were selected manually. The measurement accuracy of the tie points and image points of the GLAS footprints was achieved at one pixel. To analyze the geometric accuracy, a DOM with a resolution of 0.271 m and a DEM with a resolution of 1.806 m were used to evaluate the plane accuracy and altitude accuracy, respectively, as shown in Figure 5. The plane accuracy of the DOM and the altitude accuracy of the DEM were higher than 0.3 m. The true object coordinates of the tie points were obtained from the reference data and were used as check points (CKPs) to analyze geometric accuracy.

3.2. Combined Geopositioning of Optical and SAR Imagery

To analyze the influence of the distribution of SAR on geometric accuracy in combined geopositioning, three experiments were conducted using three GF2 images without GF3 images, with one GF3 image, and with two GF3 images. Another experiment using two GF3 images was conducted. The root mean squared error (RMSE) and maximum error (MAX) of plane and altitude were calculated to evaluate geometric accuracy, as listed in part A of

Table 5. Geometric accuracy of GF2, GF3, and GLAS combined geopositioning.

	Experimental Data	Average	RMSE (m)			MAX (m)
		Convergent Angle (${}^{\circ }$)	Plane	Altitude		Plane	Altitude
A	GF2	18.630	13.085	75.857		17.489	80.069
	GF3	73.140	1.532	3.272		2.684	5.597
	GF2 + GF3-1	27.644	11.977	16.095		13.763	20.732
	GF2 + GF3-1 + GF3-2	34.845	1.724	3.784		3.279	8.903
	GF2 + GF3-1 + GF3-2 (with the same weight)	34.845	8.107	6.061		11.191	8.952
B	GF2 + one GLAS footprint	-	14.129	3.208		17.507	6.404
	GF2 + four GLAS footprints	-	13.890	2.865		17.410	5.842
	GF2 + nine GLAS footprints	-	13.945	2.882		17.447	5.560
	GF2 + four height control points	-	13.974	1.805		17.264	4.367
C	GF2 + GF3 + four GLAS footprints	-	1.587	1.985		3.455	7.424
	GF2 + four control points	-	1.945	1.803		4.885	3.789

. The plane error and altitude error of CKs are shown in Figure 6. The results of these experiments reveal that:

• Due to the small convergent angle, the geometric accuracy of GF2 geopositioning was very low, with a plane RMSE of 13.085 m and an altitude RMSE of 75.857 m and a plane MAX of 17.489 m and an altitude MAX of 80.069 m at the convergent angle of 18.630${}^{\circ }$. As the number of GF3 images increased in the third and fourth rows of part A, the convergent angle increased gradually, which further improved the combined geopositioning geometric accuracy. The convergent angle of GF2 + GF3-1 geopositioning was 27.644${}^{\circ }$. The geometric accuracy of GF2 + GF3-1 geopositioning was a plane RMSE of 11.977 m and an altitude RMSE of 16.095 m and a plane MAX of 13.763 m and an altitude MAX of 20.732 m. The geometric accuracy of GF2 + GF3-1 + GF3-2 geopositioning was a plane RMSE of 1.724 m and an altitude RMSE of 3.784 m and a plane MAX of 3.279 m and an altitude MAX of 8.903 m in the convergent angle of 34.845${}^{\circ }$. Therefore, the convergent angle plays an important role in the improvement of geometric accuracy.
• The geometric accuracy of GF2 and GF3 combined geopositioning (plane RMSE of 1.724 m, altitude RMSE of 3.784 m, plane MAX of 3.279 m, and altitude MAX of 8.903 m) was lower than the accuracy of GF3 geopositioning (plane RMSE of 1.532 m, altitude RMSE of 3.272 m, plane MAX of 2.684 m, and altitude MAX of 5.597 m). The lower geometric accuracy of the GF2 imagery limited the geometric accuracy of combined geopositioning. Thus, the geometric accuracy for combined geopositioning is lower than that attained using SAR geopositioning.
• The geometric accuracy of GF2 and GF3 combined geopositioning with the same weight setting (plane RMSE of 8.107 m, altitude RMSE of 6.061 m, plane MAX of 11.191 m, and altitude MAX of 8.952 m) is lower than the accuracy of GF2 and GF3 combined geopositioning with the different weight setting. The weight setting is an important factor influencing the geometric accuracy for combined geopositioning.
• The geometric error in the overlap area of the GF2 and GF3 images was usually lower than that in the area not covered by GF3 images as shown in Figure 6d. This result indicates that the performance of SAR images in combined geopositioning is similar to GCPs, eliminating system errors within a certain range of control.

Similarly, we conducted three experiments using two GF6 images without GF3 images, with one GF3 image, and with two GF3 images as listed in part A of

Table 6. Geometric accuracy of GF6, GF3, and GLAS combined geopositioning.

	Experimental Data	Average	RMSE (m)			MAX (m)
		Convergent Angle (${}^{\circ }$)	Plane	Altitude		Plane	Altitude
A	GF6	5.16	12.878	35.901		14.643	53.549
	GF6 + GF3-1	26.20	10.079	20.524		11.654	34.517
	GF6 + GF3-1 + GF3-2	37.43	2.823	9.186		4.131	19.543
	GF6 + GF3-1 + GF3-2 (with the same weight)	37.43	5.280	9.528		7.143	20.989
B	GF6 + one GLAS footprint	-	11.570	11.542		13.400	21.641
	GF6 + four GLAS footprints	-	11.558	11.519		13.388	21.710
	GF6 + nine GLAS footprints	-	11.512	11.418		17.447	5.560
	GF6 + four height control points	-	11.620	12.490		13.414	24.435
C	GF6 + GF3 + four GLAS footprints	-	2.908	7.669		4.343	15.877
	GF6 + four control points	-	1.737	9.290		3.761	20.144

. The plane error and altitude error of CKs are shown in Figure 7. By analyzing the results, we found:

• The geometric accuracy of GF6 geopositioning was a plane RMSE of 12.878 m and an altitude RMSE of 35.901 m at the convergent angle of 5.16${}^{\circ }$, which is better than the geometric accuracy of GF2 geopositioning. The priori geometric accuracies without GCPs of GF2 imagery and GF6 imagery were about 30 and 16 m, respectively. The geometric accuracy without GCPs of HR imagery plays an important role in the combined geopositioning.
• The geometric accuracy of GF6 and GF3-1 image combined geopositioning (plane RMSE of 10.079 m, altitude RMSE of 20.524 m) was lower than the accuracy of GF6 and two GF3 images combined geopositioning (plane RMSE of 2.823 m, altitude RMSE of 9.186 m). The convergent angles of GF6 + GF3-1 combined geopositioning and GF6 + GF3-1 + GF3-2 combined geopositioning were 26.20${}^{\circ }$ and 37.43, respectively. Therefore, the convergent angle is a key factor in the improvement of geometric accuracy, which is similar to the results of GF2 combined geopositioning.
• The geometric accuracy with the same weight setting (plane RMSE of 5.280 m, altitude RMSE of 9.528 m, plane MAX of 7.143 m, and altitude MAX of 20.989 m) was lower than the accuracy of GF6 and GF3 combined geopositioning with the different weight settings, which is similar with the results of GF2 combined geopositioning.
• The geometric error in the overlap area of the GF6 and GF3 images was usually lower than that in the area not covered by GF3 images, as shown in Figure 7c. This result indicates that the performance of SAR images in combined geopositioning is similar to GCPs, eliminating system errors within a certain range of control.

In combined geopositioning without GCPs, the weight matrix becomes a crucial factor influencing the geometric accuracy. Therefore, combined geopositioning with different weight matrixes was studied to analyze the influence of prior accuracy on geometric accuracy. The priori geometric accuracy of GF2 imagery was approximately 30 m and the geometric accuracy of the GF3 imagery was about 5 m. Errors resulting from scale and rotation were usually less than 10 pixels for both data sources. Therefore, the GF2 offset error was set to 10, 30, 50, 100, 150, and 200 m and that of GF3 was set to 1, 3, 5, 10, and 20 m. The GF2 scale and rotation error was set to 1, 3, 5, 10, 15, and 20 pixels and that of GF3 was set to 0.5, 1, 2, 5, and 10 pixels.

Table 7. Geometric accuracy with different offset parameters for GF2 and GF3 combined geopositioning. The first column relates to GF2 and the first row relates to GF3.

	1		3		5		10		20
	Plane	Altitude	Plane	Altitude	Plane	Altitude	Plane	Altitude	Plane	Altitude
10	1.802	3.720	2.135	4.035	3.277	4.699	6.599	7.952	9.304	19.376
30	1.823	3.732	1.743	3.747	1.724	3.784	2.232	3.970	5.432	4.854
50	1.839	3.735	1.766	3.727	1.720	3.716	1.718	3.670	3.727	3.538
100	1.847	3.736	1.780	3.718	1.736	3.686	1.616	3.547	2.669	3.018
150	1.848	3.737	1.783	3.717	1.740	3.682	1.611	3.525	2.488	2.926
200	1.849	3.737	1.784	3.716	1.742	3.681	1.611	3.517	2.330	2.895

lists the geometric accuracies with different offset errors and

Table 8. Geometric accuracy with different scale and rotation parameters for GF2 and GF3 combined geopositioning.

	0.5		1		2		5		10
	Plane	Altitude	Plane	Altitude	Plane	Altitude	Plane	Altitude	Plane	Altitude
1	1.588	4.550	1.718	4.210	2.144	4.042	3.396	4.199	4.115	4.346
3	1.607	4.289	1.643	4.009	1.862	3.825	2.872	3.892	3.951	4.113
5	1.659	4.232	1.661	3.975	1.748	3.786	2.415	3.772	3.632	3.997
10	1.707	4.206	1.707	3.968	1.724	3.784	1.976	3.690	3.037	3.841
15	1.719	4.202	1.723	3.968	1.738	3.790	1.856	3.658	2.750	3.749
20	1.723	4.200	1.730	3.969	1.746	3.793	1.808	3.637	2.600	3.669

lists the geometric accuracies with different scale and rotation error.

To further analyze the influence of weight matrix on geometric positioning accuracy, different weight matrixes with various offset errors and scale and rotation errors were selected, as shown in Figure 8. Figure 8a,b show the plane geometric accuracy analysis and altitude geometric accuracy analysis with different offset parameters, respectively. The black points are the result of geometric accuracy and the colored surfaces are the fitting surfaces of the black points. Figure 8c,d show the plane geometric accuracy analysis and altitude geometric accuracy analysis with different scale and rotation parameters, respectively. By analyzing the geometric accuracy with different weight matrixes, we found the following:

• When the weight of the GF3 offset error compensation parameter was determined, geometric accuracy improved as GF2 offset error increased as shown in each column of Table 7. When the GF3 offset error compensation parameter was set to 1 m, the plane geometric accuracies with the GF2 offset error compensation parameters of 10, 30, 50, 100, 150, and 200 m were 1.802, 1.823, 1.839, 1.847, 1.848, and 1.849 m, respectively. The altitude geometric accuracies with the GF2 offset error compensation parameters of 10, 30, 50, 100, 150, and 200 m were 3.720, 3.732, 3.735, 3.736, 3.737, and 3.737 m, respectively. Similarly, when the offset error of GF2 image was lower and the offset error of GF3 image was larger, the geometric accuracy obviously decreased, as shown in the blue areas of Figure 8a,b. This result indicates that the geometric accuracy of HR optical image should be set lower than the a priori accuracy, which will enhance the credibility of SAR images and further improve the accuracy of HR optical imagery and SAR imagery combined geopositioning.
• The best accuracy, plane RMSE of 1.611 m and altitude RMSE of 3.517 m, was achieved with the geometric accuracy of 10 m for GF3, which is slightly lower than the a priori accuracy of GF3 imagery. Similarly, the fitting surface was concave around the offset error of GF3 fitting of 10 m, which means the best accuracy was obtained in this area, as shown in the red areas of Figure 8a. The a priori accuracy was obtained without the interference of altitude error. However, due to geometric errors, altitude errors exist in the GF3 imagery, which further decreases the actual plane accuracy. Therefore, to achieve a higher geometric accuracy without GCPs, the weight of the offset error compensation parameters should be set slightly lower than the given a priori accuracy of SAR imagery.
• Similarly, when the weights of the GF3 scale and rotation error compensation parameter were determined, the geometric accuracy improved as GF2 scale and rotation error increased, as shown in each column of Table 8. When the GF3 scale and rotation error was set to five pixels, the plane geometric accuracies with the GF2 scale and rotation errors of 1, 3, 5, 10, 15, and 20 pixels were 3.396, 2.872, 2.415, 1.976, 1.856, and 1.808 m, respectively. The altitude geometric accuracies with the GF2 scale and rotation error of 1, 3, 5, 10, 15, and 20 pixels were 4.199, 3.892, 3.772, 3.690, 3.658, and 3.637 m, respectively. The geometric accuracy was relatively stable when the GF3 scale and rotation error was set between 0.5 and 5 pixels as shown in Figure 8c,d. However, geometric accuracy decreased rapidly when the gap between the scale and rotation error setting and the given a priori knowledge was large. The weight matrix setting with a priori knowledge of geometric accuracy of the remote sensing images can provide the highest geometric accuracy for combined geopositioning.

3.3. GLAS Footprint Extraction Evaluation

Based on the reference DEM data in Dongying, Zibo, and Jimo regions of Shandong Province in China, the accuracy of the GLAS footprint extraction with the proposed criteria was evaluated. Some GLAS footprints and the reference DEM data are shown in Figure 9.

The number of original data in the test areas was 5312. The mean and RMSE of the original GLAS footprints were –4.490 and 16.060 m, respectively. Using the criterion of GDEM2 <20, the number of footprints was 5064. The mean and RMSE of the original GLAS footprints were −1.253 and 2.799 m, respectively, as listed in

Table 9. GLAS footprint extraction with the proposed criteria.

Criterion	Number of Data	Mean (m)	RMSE (m)
Original data	5312	−4.490	16.060
With GDEM2 <20	5064	−1.253	2.799
With GDEM2 <20 and criterions in GLA14 product	3257	−1.339	2.615
With GDEM2 <20 and criterions in GLA01 product	3419	−0.865	2.049
With GDEM2 <20 and criterions in GLA14 and GLA01 products	2272	−0.919	1.969

. The footprints with mistaken height value were discard effectively. With the criterion of GDEM2 <20 and the characteristic in the GLA14 product, 2055 footprints were discarded and the altitude accuracy was improved to 2.615 m. With the criterion of GDEM2 <20 and the characteristic in the GLA01 product, only 3419 footprints were retained and altitude accuracy was improved to 2.049 m. Using the GDEM2 data and the criterions in GLA01 and GLA14 products, the mean and RMSE of the extracted footprints were −0.919 m and 1.969 m, respectively. A total of 2272 footprints retained with high altitude accuracy. The proposed method can effectively extract the GLAS footprint with high accuracy and can be used for combined geopositioning to improve the altitude accuracy of HR optical satellite remote sensing images. Nine GLAS footprints distributed uniformly in the experimental region were further selected manually. The altitude accuracy of the GLAS footprints was evaluated using the reference DEM. The mean error and RMSE were −0.514 and 0.654 m, respectively.

3.4. Combined Geopositioning with GLAS

To analyze the influence of the distribution of GLAS footprints on geometric accuracy in combined geopositioning, we performed four experiments using three GF2 images with one GLAS footprint at the center, with four GLAS footprints at the boundaries, with nine evenly distributed footprints, and with four height control points in the corners. The geometric accuracies are listed in part B of Table 5. The geometric errors of the CKPs for combined geopositioning using GLAS footprints and height control points are shown in Figure 10. These results show that:

• The altitude accuracy gradually improved as the number of GLAS footprints increased from one to four while plane accuracy was invariant. The RMSE and maximum altitude errors with one GLAS footprint were 3.208 and 6.404 m, respectively. The RMSE and MAX with four GLAS footprints improved to 2.865 and 5.842 m, respectively. However, when the other five footprints were added, the altitude error did not significantly reduce. The RMSE and MAX altitude errors with nine GLAS footprints were 2.882 and 5.560 m, respectively. The difference between the four GLAS footprints and nine GLAS footprints was 0.017 m. In practice, combined geopositioning with four GLAS footprints is the optimal choice because it achieves high altitude accuracy and saves resources for GLAS footprint acquisition.
• With four GLAS footprints, altitude accuracy significantly improved from 75.857 to 2.865 m, as shown in the second row of part B in Table 5. With four height control prints, the RMSE and MAX were 1.085 and 4.367 m, respectively. The difference between four GLAS footprints and four height control prints was 1.780 m. The altitude accuracy in combined geopositioning was slightly lower than that with four height control points. This was caused by the altitude error in the GLAS footprints.

We also performed four experiments using two GF6 images with one GLAS footprint at the center, with four GLAS footprints at the boundaries, with nine evenly distributed footprints, and with four height control points in the corners. The geometric accuracy is listed in part B of Table 6. The geometric error of CKPs for combined geopositioning using GLAS footprints and height control points is shown in Figure 11. Using one GLAS footprint, the altitude accuracy improved from 35.901 to 11.542 m. However, the altitude accuracy did not improve significantly as the number of GLAS footprints increased from one to nine while plane accuracy was invariant. The altitude accuracy with four height control prints was similar to the results obtained using GLAS footprints because the convergent angle of the GF6 images was too small. Thus, the the convergent angle is the key factor influencing the altitude accuracy.

3.5. Combined Geopositioning with SAR and GLAS

To compare combined geopositioning with the traditional geopositioning method using GCPs, four tie points in the corners were applied as GCPs for HR geopositioning. The plane accuracy and altitude accuracy for GF2 geopositioning were 1.945 and 1.803 m, respectively, as listed in the second row of part C of Table 5. In combined geopositioning with GF2, GF3, and four GLAS footprints, the plane accuracy and altitude accuracy improved to 1.587 and 1.985 m, respectively, as listed in the first row of part C of Table 5. System errors were effectively eliminated by combined geopositioning, as shown in Figure 12a. Similarly, using GCPs, the the plane accuracy and altitude accuracy for GF6 geopositioning were 1.737 and 9.290 m, respectively. For combined geopositioning with GF6, GF3, and four GLAS footprints, the plane accuracy and altitude accuracy were 2.908 and 7.669 m, respectively, as listed in the first row of part C of Table 6 and as shown in Figure 13a. The geometric accuracy of combined geopositioning was similar to that of geopositioning with GCPs for both GF2 geopositioning and GF6 geopositioning, demonstrating that combined geopositioning using multi-source remote sensing data can be applied to large-scale mapping without GCPs.

4. Conclusions

In this paper, we presented a geometric accuracy improvement method for HR optical satellite remote sensing imagery combining multi-temporal SAR Imagery and GLAS data. Based on the imaging mechanisms, we analyzed the weight determination for HR optical imagery and SAR imagery in the combined geopositioning model. We selected GLAS footprints with high altitude accuracy using the proposed criteria, which were applied as height control points in combined geopositioning. Using the GF2, GF6, GF3, and GLAS footprints, we showed that the proposed model can be effectively applied in combined geopositioning using multi-source remote sensing data to significantly improve the geometric accuracy of HR optical satellite remote sensing imagery without GCPs. The influences of the distribution and weight matrix of SAR imagery and the distribution of GLAS footprints on geometric accuracy were synthetically analyzed. This paper provides a comprehensive analysis of geometric accuracy improvement for HR optical satellite remote sensing imagery using combined geopositioning with multi-source and multi-temporal remote sensing data.