A surface albedo inversion method based on GF-1 WFV data

Abstract

The present invention proposes a surface albedo inversion method based on GF-1 WFV data, comprising the following steps: preprocessing GF-1 WFV images to obtain a GF-1 surface reflectance dataset for the study area, extracting prior knowledge of surface anisotropy from the MCD43A1 albedo product; forward calculating the prior knowledge of surface anisotropy from the MCD43A1 albedo product to obtain a parametric BRDF for the GF-1 WFV data; integrating the parametric BRDF to obtain the black-and-white sky albedo; and performing a weighted summation to obtain the GF-1 narrowband albedo; analyzing the wideband albedo to screen out the surface feature spectral curve and obtain the GF-1 wideband albedo; and verifying the inversion result of the GF-1 wideband albedo using the wideband albedo calculated from the MCD43A1 albedo product. By comparing the inversion method with the MODIS albedo product, the accuracy of the inversion method meets the required accuracy.

Description

Technical Field

[0001] The present invention relates to the technical field of spatial information remote sensing, and in particular to a surface albedo inversion method based on GF-1 WFV data. Background Art

[0002] Surface albedo is one of the driving factors of the surface radiation energy balance and the interaction between the earth and the atmosphere. Surface albedo is an important parameter that is widely used in surface energy balance, medium- and long-term weather forecasts, and global change research. Currently, there are albedo datasets generated based on satellite remote sensing data, such as GLASS (Global LAnd Surface Satellite) and MODIS (Moderate-resolution imaging spectroradiometer). However, the spatial resolution of these albedo products is low, making it difficult to accurately describe the spatiotemporal variations of surface albedo in areas with high landscape fragmentation, and it is also difficult to observe detailed surface detail feature information. Domestic high-resolution satellite surface albedo inversions are almost all established for Landsat and HJ-1 images, but there are not many surface albedo inversions for domestic high-resolution satellite series.

[0003] There are two methods for inverting surface albedo using remote sensing technology: one is to obtain surface albedo based on a statistical model. Starting from the definition of albedo, the MODTRAN 4.0 model is used to calculate the luminous flux to calculate the surface albedo. However, this method does not take into account surface anisotropy; the other inversion method considers the various anisotropic characteristics of the surface and inverts the surface albedo through the bidirectional reflectance distribution function (BRDF). The BRDF describes the anisotropic characteristics of the surface and is the ratio between the micro-increments of the radiant brightness in the incident direction and the reflected direction. Using the BRDF model to invert the albedo of satellite data requires multi-angle pixel reflectance information, which is more difficult for single-angle satellites.

[0004] The invention patent with application number 202210733214.4 discloses an albedo inversion method based on the classification of surface reflectance anisotropy characteristics. On the basis of the anisotropy smoothing index AFX, an index perpendicular to AFX is proposed, named the vertical anisotropy smoothing index PAFX. The two indices are combined to classify the global surface reflectance anisotropy characteristic information, and the albedo is inverted and verified based on the prior knowledge of the surface reflectance anisotropy characteristics. The invention can more accurately classify the surface reflectance anisotropy characteristics, improve the albedo inversion accuracy when the observation data is insufficient, and is simple and easy to operate, and is convenient to use. However, the above invention does not solve the problem that the GF-1 high-resolution satellite inverted surface albedo lacks multi-angle reflectance data, so the GF-1 surface anisotropy characteristics cannot be obtained, making it difficult to perform surface albedo inversion. Summary of the Invention

[0005] In response to the technical problems of low spatial resolution of surface albedo in existing methods and difficulty in obtaining multi-angle reflectance data in a period of time during the inversion process, the present invention proposes a surface albedo inversion method based on GF-1 WFV data. Based on the advantages of both GF-1 multispectral WFV data and MODIS albedo product (MCD43A1), the GF-1 satellite has a high spatial resolution of 16 meters, which can provide clearer surface features than most current low-resolution surface albedo, while the MODIS albedo product can obtain surface anisotropic characteristics (BRDF).

[0006] In order to achieve the above object, the technical solution of the present invention is implemented as follows: a surface albedo inversion method based on GF-1 WFV data, the steps of which are as follows:

[0007] Step 1: Download the required GF-1 WFV image from the Land Observation Satellite Data Service Platform, and then download the MCD43A1 albedo product of the study area;

[0008] Step 2: Preprocess the GF-1 WFV image and crop the preprocessed data to obtain the GF-1 surface reflectance dataset of the study area, and extract the surface anisotropy prior knowledge of the MCD43A1 albedo product;

[0009] Step 3: Use the GF-1 surface reflectance dataset to perform forward calculations on the surface anisotropy prior knowledge of the MCD43A1 albedo product. Adjust the surface anisotropy prior knowledge to obtain the parametric BRDF of the GF-1 WFV data. Integrate the parametric BRDF to obtain the black and white sky albedo, and then perform a weighted sum to obtain the GF-1 narrowband albedo.

[0010] Step 4: By analyzing the wide-band albedo, the spectral curve of the ground object is screened out. The GF-1 narrow-band albedo obtained in step 3 is analyzed by multiple linear regression to obtain the visible light wide-band conversion model of GF-1 WFV data and the GF-1 wide-band albedo;

[0011] Step 5: Use the wide-band albedo calculated using the MCD43A1 albedo product to verify and analyze the inversion results of the GF-1 wide-band albedo.

[0012] Preferably, the preprocessing of the GF-1 WFV image in step 2 includes radiometric calibration and atmospheric correction, the GF-1 WFV image is subjected to radiometric calibration and atmospheric correction using the 6S model to obtain GF-1 WFV reflectivity data, the GF-1 WFV reflectivity data is reprojected into the WGS-84 coordinate system to obtain a standard geographic coordinate system, and then the corresponding study area is cropped out, and the GF-1 WFV reflectivity data of the study area are combined to obtain the GF-1 surface reflectivity dataset of the study area.

[0013] Preferably, the method for extracting the surface heterogeneity prior knowledge is: using the MRT tool to convert the MCD43A1 albedo product from a sinusoidal projection to a WGS-84 coordinate system and reproject it to obtain a standard geographic coordinate system; the surface heterogeneity prior knowledge of the MCD43A1 albedo product is the surface parameter BRDF', and the surface parameter BRDF' selects the kernel coefficient f of the red, green and blue bands of each pixel in the 7 bands iso 、f vol 、f geo .

[0014] Preferably, the implementation method of step three is:

[0015] Assume that the reflectivity of the GF-1 WFV image is R G , the surface anisotropy prior knowledge of the MCD43A1 albedo product is brought into the kernel-driven model to obtain the MODIS surface reflectance R M , calculate the reflectivity R G Compared with MODIS surface reflectance R M The ratio of:

[0016] By adjusting the reflectivity ratio C λ The surface parameter BRDF' of the surface anisotropy prior knowledge of the MCD43A1 albedo product is adjusted to obtain the parameter BRDF=C of the GF-1 WFV data to be determined. λ ×BRDF′; that is, the kernel coefficient f iso (λ), f vol (λ), f geo (λ);

[0017] Integrate the parameter BRDF of the GF-1 WFV data to obtain the corresponding black sky albedo a bsa and white sky albedo a wsa They are:

[0018]

[0019]

[0020] hi (θ i )=g 0i +g 1i θ 2 +g 2i θ 3 ;

[0021] Among them, θ, g 0i 、g 1i 、g 2i They represent the pixel solar zenith angle and the coefficients of the polynomial fitting expression of each kernel, h i (θ)) is the integral value of the i-th nucleus in the observation hemisphere, H i is the integral value of the nucleus in the incident and observation hemispheres, Respectively represent f iso (λ), f vol (λ), f geo (λ);

[0022] The black sky albedo a bsa (θ) and white sky albedo a wsa Perform weighted combination to obtain the narrow band albedo of band λ: a(θ,λ)=[(1-s(τ))a bsa (θ)]+s(τ)a wsa ;

[0023] Among them, s(τ) is the proportion of sky scattered light, and its size is estimated by the parameters after atmospheric correction.

[0024] Preferably, the kernel-driven model normalizes the MCD43A1 albedo product to the same order of magnitude, and normalizes the kernel coefficients f of the red, green and blue bands of each pixel to iso 、f vol 、f geo Normalization obtains the normalized BRDF' parameters:

[0025]

[0026] Where R is the surface bidirectional reflectivity function, θ is the solar zenith angle, To observe the zenith angle, is the relative azimuth; λ is the wavelength; K vol and K geo They are represented as volume scattering kernel and geometric optics kernel, f iso (λ), f vol (λ), f geo (λ) are constant coefficients related to wavelength.

[0027] Preferably, the volume scattering kernel K voland the geometric optics kernel K geo Select the RossThick kernel function and the LiTransit kernel function respectively;

[0028] The MODIS surface reflectance R M The calculation method is: the BRDF parameter fi of the MCD43A1 albedo product iso (λ), f vol (λ), f geo (λ) is brought into the following kernel-driven model to obtain:

[0029]

[0030] Preferably, the coefficients of the polynomial fitting expressions of the respective kernels are:

[0031]

[0032] Preferably, the implementation method of step 4 is as follows: 190 surface feature spectral curves are screened from the USGS digital spectral library using the SBDART model, including grassland, water, farmland, and urban surface feature types, six atmospheric visibility, eight solar zenith angles, and three atmospheric modes are input, and a multivariate linear regression analysis is established to obtain a visible light wide-band conversion model for GF-1 WFV data:

[0033] α GF-1 =0.443α B +0.317α G +0.240α R ;

[0034] Among them, α B , α G , α R are the narrow-band albedo of the blue light band, green light band, and red light band, respectively, which are obtained through the narrow-band albedo of the band λ; α GF-1 Represents GF-1 broadband albedo.

[0035] Preferably, the broadband albedo is the ratio of the upward and downward radiation fluxes of the surface within a certain wavelength range:

[0036]

[0037] Where α(θ,Λ) is the broadband albedo, Λ is the band range from λ1 to λ2, and F u (θ,Λ) and F d (θ, Λ) are the upward and downward radiation fluxes, respectively, obtained by the SBDART model; α(θ, λ) is the narrow-band albedo of band λ, and θ is the solar apex angle;

[0038] The solar zenith angle ranges from 0 to 80°, with each interval being 10°. The three atmospheric modes include tropical, mid-latitude winter, and mid-latitude summer.

[0039] The method for performing verification analysis in step 5 is:

[0040] Randomly select sample points for root mean square error analysis, and the root mean square error is:

[0041] Where N is the number of sample points, RMSE is less than 0.05; a GF-1 represents the GF-1 wide-band albedo, a MODIS is the MODIS surface albedo.

[0042] Preferably, the MODIS surface albedo a MODIS The calculation method is: the kernel coefficient f of each pixel in the red, green and blue bands of the MCD43A1 albedo product iso 、f vol 、f geo Bring it into the MODIS BRDF model for integration to obtain the black sky albedo and white sky albedo of these three bands, and then substitute the black sky albedo and white sky albedo of these three bands into the narrow band to wide band conversion relationship a of the MODIS surface albedo product itself. MODIS =0.331a R +0.424a B +0.246a G In the MODIS Surface albedo; where a R 、a B 、a G It is the narrow-band albedo of the red, blue, and green light bands of the MCD43A1 albedo product.

[0043] The present invention is based on the kernel-driven model. First, the GF-1 data is radiometrically calibrated and atmospherically corrected to obtain the surface reflectance. Then, the MODIS albedo product (MCD43A1) corresponding to the study area is selected to extract the prior knowledge of surface reflectance anisotropy. The knowledge is fitted with GF-1 satellite data, and finally the narrow-band to wide-band conversion is performed to obtain the GF-1 surface albedo.

[0044] Compared with the existing technology, the present invention has the following advantages: based on the advantages of GF-1 multispectral WFV data and MODIS albedo product (MCD43A1), the kernel-driven model is used to preprocess the GF-1 WFV data to obtain processed images, radiometric calibration and atmospheric correction are performed on the GF-1 data to obtain surface reflectance, and prior knowledge of surface reflectance anisotropy is extracted from the MODIS albedo product MCD43A1 as the underlying surface dataset. This is then fitted with GF-1 WFV image data and converted from narrow band to wide band, thereby obtaining a BRDF dataset under the GF-1 WFV data, and further obtaining a high-resolution surface albedo product based on the GF-1 WFV data. By comparing with the MODIS albedo product, the inversion accuracy can meet the accuracy requirements. BRIEF DESCRIPTION OF THE DRAWINGS

[0045] In order to more clearly illustrate the embodiments of the present invention or the technical solutions in the prior art, the following briefly introduces the drawings required for use in the embodiments or the description of the prior art. Obviously, the drawings described below are only some embodiments of the present invention. For ordinary technicians in this field, other drawings can be obtained based on these drawings without paying any creative work.

[0046] Figure 1 It is a schematic diagram of the process of the present invention.

[0047] Figure 2 This is one of the experimental comparison diagrams of the present invention, wherein (a) is the GF-1 surface albedo and (b) is the MODIS surface albedo.

[0048] Figure 3 is the fitting curve of the sample points of the present invention.

[0049] Figure 4 This is the absolute error diagram of the albedo between GF-1 and MODIS of the present invention. DETAILED DESCRIPTION

[0050] The following will clearly and completely describe the technical solutions in the embodiments of the present invention in conjunction with the accompanying drawings. Obviously, the described embodiments are only part of the embodiments of the present invention, not all of the embodiments. All other embodiments obtained by ordinary technicians in this field based on the embodiments of the present invention without creative work are within the scope of protection of the present invention.

[0051] like Figure 1 As shown in Figure 1, a surface albedo inversion method based on GF-1 WFV data has the following steps:

[0052] Step 1: Download the required GF-1 WFV imagery from the Land Observation Satellite Data Service Platform. The GF-1's four WFV cameras have a swath width of 800 kilometers and a resolution of 16 meters, achieving a perfect combination of wide coverage and high spatial resolution. Then, download the MCD43A1 albedo product for the study area from the EARTHDATA website. The MCD43A1 albedo product is widely used worldwide and is one of the most recognized products in the international remote sensing community.

[0053] Step 2: Preprocess the GF-1 WFV image and crop the preprocessed data to obtain the GF-1 surface reflectance dataset of the study area, and extract the surface anisotropy prior knowledge of the MCD43A1 albedo product.

[0054] The downloaded GF-1 WFV satellite image was preprocessed, including radiometric calibration and atmospheric correction. The 6S model was used to perform radiometric calibration and atmospheric correction on the GF-1 WFV image to obtain the GF-1 WFV reflectance data. The corrected image eliminated the atmospheric effect and became clearer. The corresponding study area was cropped after reprojection to the WGS-84 coordinate system. The WGS coordinate system reprojection can obtain the standard geographic coordinate system, and then the GF-1 surface reflectance dataset of the study area was obtained.

[0055] The MRT tool is used to convert the MCD43A1 albedo product from a sinusoidal projection to the WGS-84 coordinate system and reproject it to obtain a standard geographic coordinate system. The MCD43A1 albedo product provides a 16-day global composite 500m surface parameter BRDF', which includes the kernel coefficient f for each pixel in the seven bands. iso 、f vol 、f geo , where the MRT tool is used to select the kernel coefficient f of the red, green and blue bands iso 、f vol 、f geo .

[0056] Step 3: Use the GF-1 surface reflectance dataset to forward calculate the surface parameter BRDF' of the MCD43A1 albedo product to obtain the MODIS surface reflectance data R M Then, the ratio of the measured reflectivity and the simulated reflectivity in the same observation geometric direction is calculated to obtain C λ , adjust the surface parameter BRDF' to obtain the parameter BRDF of GF-1 WFV data, integrate it to obtain the black and white sky albedo, and then add it in a certain proportion to obtain the GF-1 narrow band albedo.

[0057] The specific implementation method is:

[0058] Step 3.1: Based on the kernel-driven model, extract prior knowledge from the MCD43A1 albedo product. Due to the complex surface anisotropy characteristics, it is necessary to normalize the MCD43A1 albedo product. Normalize the MCD43A1 albedo product to the same order of magnitude and convert f iso (λ), f vol (λ), f geo (λ) The three parts are normalized to obtain the normalized BRDF' parameters:

[0059]

[0060] Where R is the surface bidirectional reflectivity, θ is the solar zenith angle, To observe the zenith angle, is the relative azimuth; λ is the wavelength; K vol and K geo They are represented as volume scattering kernel and geometric optics kernel, f iso (λ), f vol (λ), f geo (λ) are constant coefficients related to wavelength. The BRDF (f iso (λ), f vol (λ), f geo (λ)). In order to avoid negative reflectivity, when choosing K vol and K feo The two kernels are selected as RossThick kernel function and LiTransit kernel function respectively. The LiTransit kernel function can effectively avoid the situation where the fitting reflectivity is negative, and its fitting ability is strong.

[0061] Step 3.2: Assume that the GF-1 WFV reflectivity is R G , by taking the MCD43A1 albedo parameter f iso (λ), f vol (λ), f geo (λ) is brought into the kernel-driven model as in formula (2), and the MODIS surface reflectance R is obtained. M , and then get the reflectivity ratio of the two:

[0062]

[0063]

[0064] By adjusting the reflectivity ratio C λ Adjust the surface parameter BRDF' of the surface anisotropy prior knowledge of the MCD43A1 albedo product, and obtain the BRDF of the GF-1 WFV data to be determined from the MCD43A1 prototype parameters

[0065] BRDF = C λ ×BRDF′ (4).

[0066] Step 3.3: Integrate the parameter BRDF of the GF-1 WFV data to obtain the corresponding black sky albedo a bsa and white sky albedo a wsa They are:

[0067]

[0068]

[0069] h i (θ i )=g 0i +g 1i θ 2 +g 2i θ 3 (7)

[0070] Among them, θ, g 0i 、g 1i 、g 2i They represent the pixel solar zenith angle and the coefficients of the polynomial fitting expression of each kernel, and their values ​​are shown in Table 1 below. The kernel integral h i (θ) is a term calculated in formula (5), H i The term calculated by formula (6) has nothing to do with the observation angle, H i The corresponding integral values ​​in RossThick and LiTransit are 0.189184 and -1.137762 respectively. i (λ) contains f iso (λ), f vol (λ), f geo The three terms (λ) are calculated by formula (4). Multiplying the two together can give the black sky albedo and white sky albedo in the λ band.

[0071] Table 1 Coefficients of polynomial fitting expressions

[0072]

[0073] Step 3.4: In order to obtain the narrowband albedo, the two black and white sky albedos are weighted and combined in a certain proportion:

[0074] a(θ,λ)=[(1-s(τ))a bsa (θ)]+s(τ)a wsa (8)

[0075] Where s(τ) is the fraction of sky scattered light, estimated from the atmospheric correction parameters in step 2. a(θ) is the narrowband albedo of band λ, which will be calculated for the wideband later.

[0076] Step 4: Broadband albedo is defined as the ratio of the upgoing to downgoing radiation flux from the surface over a certain wavelength range:

[0077]

[0078] Where α(θ,Λ) is the broadband albedo, Λ is the band range from λ1 to λ2, and F u (θ,Λ) and F d (θ, Λ) are the upgoing and downgoing radiation fluxes, respectively, obtained through the SBDART model. α(θ, λ) is the narrow-band albedo of band λ, and θ is the solar zenith angle. The SBDART model calculates the upgoing and downgoing radiation fluxes of the Earth. SBDART has a good interface, allowing us to efficiently establish simulated radiation fluxes. Through the SBDART model, 190 typical ground feature spectral curves were screened from the USGS digital spectral library in the United States, including grasslands, water bodies, farmland, cities and other ground feature types. Six atmospheric visibility and eight solar zenith angles (0-80°) with a 10° interval and three atmospheric models were input, including tropical, mid-latitude winter, and mid-latitude summer. A multivariate linear regression analysis was established to obtain a visible light wide-band conversion model suitable for GF-1 WFV data:

[0079] α GF-1 =0.443α B +0.317α G +0.240α R (10)

[0080] Among them, α B , α G , α R are the narrow band albedo of the blue light band, green light band, and red light band, respectively, which are obtained by a(θ) in formula (7), α GF-1 Represents GF-1 broadband albedo.

[0081] Step 5: Use the wide-band albedo calculated by MCD43A1 to verify the GF-1 surface albedo inversion result obtained by formula (10).

[0082] Randomly select sample points and perform root mean square error analysis. The root mean square error represents the square root of the ratio of the square of the deviation between the observed value and the true value to the number of observations n:

[0083]

[0084] Among them, the RMSE accuracy requirement is less than 0.05. GF-1 、a MODIS denote the GF-1 surface albedo and the MODIS surface albedo, respectively. GF-1 It is obtained by formula (10). The kernel coefficient f of each pixel in the red, green and blue bands of MODIS MCD43A1 product is iso 、f vol 、f geo Bring it into the MODIS BRDF model for integration to obtain the black sky albedo and white sky albedo of these three bands. Then substitute the black sky albedo and white sky albedo of these three bands into the narrow band to wide band conversion relationship of the MODIS surface albedo product itself, such as formula (12), to obtain a MODIS Surface albedo; where a R 、a B 、a G It is the narrow-band albedo of the red, blue, and green light bands of the MCD43A1 albedo product.

[0085] a MODIS =0.331a R +0.424a B +0.246a G (12)

[0086] In this example, the GF-1 is a high-temporal satellite that acquires a large dataset. It carries two panchromatic / multispectral (P / MS) cameras and four WFV cameras. The four WFV cameras have a swath width of 800 kilometers and a resolution of 16 meters. The MODIS is a key sensor on the Aqua and Terra satellites. It observes the entire Earth's surface every one to two days and acquires 36 bands, seven of which are used to generate BRDF products.

[0087] In this example, the Zhengzhou Airport Economic Zone in Henan Province, my country, was selected as the study area for surface albedo inversion. First, the GF-1 WFV imagery from March 18, 2020, was preprocessed. The Zhengzhou Airport Economic Zone study area was cropped after preprocessing. BRDF parameters for the three visible light bands were extracted from the MCD43A1 product. These model parameters were then incorporated into their own conversion coefficient equations to obtain the MODIS albedo product. The two albedo products were cross-validated.

[0088] Figure 2 is obtained by bringing in the GF-1 band transfer function, Figure 2It can be seen that the GF-1 albedo and MODIS albedo products are consistent in image distribution. Then 400 sample points were randomly selected in the test area and fitted and analyzed. Figure 3 As shown in the figure, the RMSE obtained can be less than 0.05, which can meet the accuracy requirements of albedo. Figure 4 The mean absolute error of the GF-1 satellite albedo is 0.0185.

[0089] The method first preprocesses GF-1 WFV data to obtain processed images. It then extracts prior knowledge of surface anisotropy and obtains GF-1 surface BRDF parameters through fitting. Black-and-white sky surface albedo is obtained by integrating the kernel-driven model. Finally, a narrow-band to wide-band conversion is performed to obtain the GF-1 surface albedo product. Comparisons with the MODIS albedo product demonstrate that the inversion accuracy meets the required precision.

[0090] This paper proposes using the MODIS albedo product (MCD43A1) to obtain surface anisotropy characteristics, thereby obtaining GF-1 WFV surface anisotropy information, then calculating the GF-1 narrowband albedo, and finally converting the narrowband albedo to a wideband albedo. Using experimental data from the Zhengzhou Airport GF-1 imagery as an example, this paper extracts prior knowledge of surface anisotropy from the MCD43A1 imagery to overcome the GF-1 high-resolution satellite's inability to acquire multi-angle reflectivity data, and thus its inability to obtain surface anisotropy characteristics. This paper provides a solution for inverting surface albedo from GF-1 WFV data, and has significant application value for inverting surface albedo from GF-1 satellites.

[0091] The above description is only a preferred embodiment of the present invention and is not intended to limit the present invention. Any modifications, equivalent substitutions, improvements, etc. made within the spirit and principles of the present invention should be included in the scope of protection of the present invention.

Claims

1. A surface albedo inversion method based on GF-1WFV data, characterized in that: The steps are as follows: Step 1: Download the required GF-1WFV image from the Land Observation Satellite Data Service Platform, and then download the MCD43A1 albedo product of the study area; Step 2: Preprocess the GF-1WFV image and crop the preprocessed data to obtain the GF-1 surface reflectance dataset of the study area, and extract the surface anisotropy prior knowledge of the MCD43A1 albedo product; Step 3: Use the GF-1 surface reflectance dataset to perform forward calculations on the surface anisotropy prior knowledge of the MCD43A1 albedo product. Adjust the surface anisotropy prior knowledge to obtain the parametric BRDF of the GF-1WFV data. Integrate the parametric BRDF to obtain the black and white sky albedo, and then perform a weighted sum to obtain the GF-1 narrowband albedo. Step 4: By analyzing the wide-band albedo, the spectral curve of the ground object is screened out. The GF-1 narrow-band albedo obtained in step 3 is analyzed by multiple linear regression to obtain the visible light wide-band conversion model of GF-1WFV data and the GF-1 wide-band albedo; Step 5: Use the wide-band albedo calculated using the MCD43A1 albedo product to verify and analyze the inversion results of the GF-1 wide-band albedo.

2. The surface albedo inversion method based on GF-1WFV data according to claim 1, characterized in that: In step 2, the GF-1WFV image is preprocessed including radiometric calibration and atmospheric correction. The GF-1WFV image is radiometrically calibrated and atmospherically corrected using the 6S model to obtain the GF-1WFV reflectance data. The GF-1WFV reflectance data is reprojected to the WGS-84 coordinate system to obtain the standard geographic coordinate system, and then the corresponding study area is cropped. The GF-1WFV reflectance data of the study area are combined to obtain the GF-1 surface reflectance dataset of the study area.

3. The surface albedo inversion method based on GF-1WFV data according to claim 1 or 2, characterized in that: The method for extracting the surface heterogeneity prior knowledge is as follows: using the MRT tool to convert the MCD43A1 albedo product from a sinusoidal projection to the WGS-84 coordinate system and reproject it to obtain a standard geographic coordinate system; the surface heterogeneity prior knowledge of the MCD43A1 albedo product is the surface parameter BRDF', which selects the kernel coefficient f of the red, green and blue bands of each pixel in the seven bands. iso 、f vol 、f geo .

4. The surface albedo inversion method based on GF-1WFV data according to claim 3 is characterized in that: The implementation method of step three is: Assume that the reflectivity of GF-1WFV image is R G , the surface anisotropy prior knowledge of the MCD43A1 albedo product is brought into the kernel-driven model to obtain the MODIS surface reflectance R M , calculate the reflectivity R G Compared with MODIS surface reflectance R M The ratio of: By adjusting the reflectivity ratio C λ The surface parameter BRDF' of the surface anisotropy prior knowledge of the MCD43A1 albedo product is adjusted to obtain the parameter BRDF=C of the GF-1 WFV data to be determined. λ ×BRDF′; that is, the kernel coefficient f iso (λ), f vol (λ), f geo (λ); Integrate the parameter BRDF of the GF-1 WFV data to obtain the corresponding black sky albedo a bsa and white sky albedo a wsa They are: h i (i i )-g 0i +g 1i i 2 +g 2i i 3 ; Among them, θ, g 0i 、g 1i 、g 2i They represent the pixel solar zenith angle and the coefficients of the polynomial fitting expression of each kernel, h i (θ)) is the integral value of the i-th nucleus in the observation hemisphere, H i is the integral value of the nucleus in the incident and observation hemispheres, f i (λ) represents the kernel coefficient f iso (λ), f vol (λ), f geo (λ); The black sky albedo a bsa (θ) and white sky albedo a wsa Perform weighted combination to obtain the narrow band albedo of band λ: a(θ,λ)=[(1-s(τ))a bsa (θ)]+s(τ)a wsa ; Among them, s(τ) is the proportion of sky scattered light, and its size is estimated by the parameters after atmospheric correction.

5. The surface albedo inversion method based on GF-1 WFV data according to claim 4, characterized in that: The kernel-driven model normalizes the MCD43A1 albedo product to the same order of magnitude and normalizes the kernel coefficients f of the red, green and blue bands of each pixel. iso 、f vol 、f geo Normalization obtains the normalized BRDF′ parameters: Where R is the surface bidirectional reflectivity function, θ is the solar zenith angle, To observe the zenith angle, is the relative azimuth; λ is the wavelength; K vol and K geo They are represented as volume scattering kernel and geometric optics kernel, f iso (λ), f vol (λ), f geo (λ) are constant coefficients related to wavelength.

6. The surface albedo inversion method based on GF-1 WFV data according to claim 5, characterized in that: The volume scattering kernel K vol and the geometric optics kernel K geo Select the RossThick kernel function and the LiTransit kernel function respectively; The MODIS surface reflectance R M The calculation method is: the BRDF parameter f of the MCD43A1 albedo product iso (λ), f vol (λ), f geo (λ) is brought into the kernel-driven model to obtain:

7. The surface albedo inversion method based on GF-1 WFV data according to claim 6, characterized in that: The coefficients of the polynomial fitting expressions of the respective kernels are:

8. The surface albedo inversion method based on GF-1WFV data according to any one of claims 4 to 7, characterized in that: The implementation method of step 4 is as follows: 190 surface feature spectral curves are screened from the USGS digital spectral library using the SBDART model, including grassland, water, farmland, and urban surface feature types. Six atmospheric visibility levels, eight solar zenith angles, and three atmospheric modes are input, and a multivariate linear regression analysis is established to obtain a visible light broadband conversion model for GF-1WFV data: α GF-1 ＝0.443α B +0.317α G +0.240α R ； Among them, α B , α G , α R are the narrow-band albedo of the blue, green, and red bands of the GF-1WFV data, respectively, obtained through the narrow-band albedo of band λ; α GF-1 Represents GF-1 broadband albedo.

9. The surface albedo inversion method based on GF-1WFV data according to claim 8, characterized in that: The broadband albedo is the ratio of the upgoing to downgoing radiation flux from the Earth's surface within a certain wavelength range: Where α(θ,Λ) is the broadband albedo, Λ is the band range from λ1 to λ2, and F u (θ,Λ) and F d (θ, Λ) are the upward and downward radiation fluxes, respectively, obtained by the SBDART model; α(θ, λ) is the narrow-band albedo of band λ, and θ is the solar apex angle; The solar zenith angle ranges from 0° to 80°, with each interval being 10°. The three atmospheric modes include tropical, mid-latitude winter, and mid-latitude summer. The method for performing verification analysis in step 5 is: Randomly select sample points for root mean square error analysis, and the root mean square error is: Where N is the number of sample points, RMSE is less than 0.05; a GF-1 represents the GF-1 wide-band albedo, a MODIS is the MODIS surface albedo.

10. The surface albedo inversion method based on GF-1WFV data according to claim 9, characterized in that: The MODIS surface albedo a MODIS The calculation method is: the kernel coefficient f of each pixel in the red, green and blue bands of the MCD43A1 albedo product iso 、f vol 、f geo Bring it into the MODIS BRDF model for integration to obtain the black sky albedo and white sky albedo of these three bands, and then substitute the black sky albedo and white sky albedo of these three bands into the narrow band to wide band conversion relationship a of the MCD43A1 surface albedo product itself. MODIS =0.331a R +0.424a B +0.246a G In the MODIS Surface albedo; where a R 、a B 、a G It is the narrow-band albedo of the red, blue, and green light bands of the MCD43A1 albedo product.

Images