Surface reflectance uncertainty analysis method based on atmospheric terrain coupling model

Abstract

In order to solve the problem of hyperspectral surface reflectance product quality evaluation under the non-uniform and undulating surface condition, the invention provides a surface reflectance uncertainty analysis method based on an atmospheric topography coupling model, and on the basis of an atmospheric topography coupling radiation transmission model, two uncertainty quantification methods of a formula transfer method and a Monte Carlo method are combined, so that the quality of a hyperspectral surface reflectance product is evaluated. And the influence of the radiance product precision, the imaging geometric parameter precision, the terrain parameter precision and the atmospheric parameter precision on the surface reflectance is evaluated. The uncertainty of the parameters and the probability distribution of the parameters serve as input of uncertainty analysis, and the uncertainty is transmitted through a formula transmission method for the process with analytical formula transmission on a surface reflectance inversion processing chain; under the condition that no available analytical formula exists, uncertainty propagation is simulated by applying a Monte Carlo method, the uncertainty of the surface reflectance under the condition of non-uniform and undulating surface is obtained, and uncertainty contribution quantitative analysis sorting and main uncertainty traceability are carried out.

Description

Technical Field

[0001] The present invention relates to a surface reflectance uncertainty analysis method based on an atmosphere-terrain coupling model, which can be used in the atmosphere-terrain joint correction process of satellite-borne or airborne hyperspectral remote sensing data, thereby serving the quality assessment of reflectance data obtained by remote sensors, and is a supporting technology in the field of earth observation. Background Art

[0002] In hyperspectral remote sensing, surface reflectance is the basic data for subsequent applications. Due to the high data dimension, weak signal energy and long information chain, the accuracy of surface reflectance has become a key factor affecting quantitative applications. With the development of satellite technology, the remote sensing field has attached more and more importance to metrological research. In 2010, the Committee on Earth Observation Satellites established the "Quality Assurance Framework for Earth Observation", introduced the main methods of the "Guide to the Expression of Uncertainty in Measurement" (GUM) into the field of earth observation, and applied GUM to the uncertainty analysis of earth observation for the first time. The uncertainty assessment of hyperspectral remote sensing reflectance products has attracted widespread attention from researchers at home and abroad.

[0003] The existing accuracy evaluation method in the field of quantitative remote sensing generally uses the measured reflectance of a flat and uniform ground as a reference value, and evaluates the accuracy level of the remote sensing inversion reflectance data by calculating the average error or root mean square error of multiple verification points. However, in the case of non-uniform and undulating surfaces, the lack of quantification of the accuracy of surface reflectance under the influence of atmosphere-topography coupling makes it difficult to trace the data errors and improve the quality. At present, only Ma Zhaoning from the Beijing University of Aeronautics and Astronautics has published a paper entitled "Uncertainty Analysis for Topographic Correction of Hyperspectral Remote Sensing Images". This document discloses an uncertainty analysis method for terrain correction, but the analyzed and processed data remain at radiance data, and the uncertainty of the reflectance product needs further analysis and discussion. The Institute of Space Information Innovation of the Chinese Academy of Sciences has applied for a patent for "Method and Device for Realizing Quality Traceability of Land Surface Remote Sensing Surface Reflectance Products" (CN 115824971 A). This application quantitatively measures the uncertainty of hyperspectral remote sensing reflectance products, but does not take into account factors such as terrain undulation and non-uniformity of land objects. Therefore, the uncertainty measurement of this application is inaccurate and is not adaptively adjusted based on actual land objects and actual terrain conditions. Therefore, in order to evaluate the quality and stability of hyperspectral remote sensing reflectance products under non-uniform and undulating terrain conditions, trace the main sources of uncertainty, support algorithm improvement and data accuracy improvement, it is urgent to develop uncertainty quantification and analysis methods for surface reflectance under non-uniform and undulating surface conditions. Summary of the invention

[0004] In order to solve the problem of quality assessment of hyperspectral surface reflectance products under non-uniform and undulating surfaces, the present invention proposes a surface reflectance uncertainty analysis method based on an atmosphere-terrain coupling model. Based on the radiation transmission model of atmosphere-terrain coupling, the formula transfer method and the Monte Carlo method are combined to evaluate the influence of radiance product accuracy, imaging geometry parameter accuracy, terrain parameter accuracy and atmospheric parameter accuracy on surface reflectance. The uncertainty of the above parameters and their related probability distribution are used as inputs for uncertainty analysis. In the processing chain of surface reflectance inversion, for the process with analytical formula transmission, the uncertainty is transmitted by the formula transfer method. In the absence of an available analytical formula, the Monte Carlo method is used to simulate uncertainty propagation to obtain the uncertainty of hyperspectral surface reflectance under non-uniform and undulating surfaces, and the uncertainty contribution quantitative analysis and ranking and main uncertainty traceability are performed.

[0005] The uncertainty analysis method of surface reflectivity based on the atmosphere-terrain coupling model is characterized by comprising the following steps: Step 1: Establish an atmosphere-terrain joint correction model to complete the value transfer of hyperspectral remote sensing radiance data to surface reflectance data; Step 2: Based on the value transfer process, the uncertainty is traced back to the uncertainty input of the model, which is the accuracy of the radiance product, the accuracy of the imaging geometric parameters, the accuracy of the terrain parameters, and the accuracy of the atmospheric parameters; Step 3: For each input variable obtained by tracing back in step 2, estimate the uncertainty by using the type A or type B assessment method and determine the probability distribution of these uncertain inputs; Step 4: Use the uncertainty assessment results of the input variables obtained in step 3 to transfer and synthesize the uncertainty, and transfer the uncertainty of the formula transfer method to the links with analytical formulas and transfer the uncertainty of the Monte Carlo method to the links without analytical formulas, and calculate the uncertainty assessment results of the surface reflectance synthesis standard; Step 5: Analyze the uncertainty components in the synthetic surface reflectance, and rank and compare the uncertainty contributions of each link band by band, and further trace the main uncertainties to determine the main source factors affecting the uncertainty of the surface reflectance, guide parameter optimization and algorithm improvement, and achieve improved data accuracy.

[0006] The core formula for inverting the reflectivity of the undulating and non-uniform surface in step 1 is:

[0007] In this formula, is the direct solar radiation, is the atmospheric scattered radiation, is the ground-atmosphere coupled radiation, Radiation reflected by surrounding objects, is the radiance caused by the atmospheric proximity effect, is the atmospheric radiation, is the radiance data, is the atmospheric upward transmittance, is the reflectivity of the target object, Represents the number of rows, columns and elevation respectively. is the wavelength, and the uncertainty introduced by the radiation transfer model error and the surface Lambertian assumption is represented by “+0”.

[0008] Among them, the sources of uncertainty traced back in step 2 specifically include: the accuracy of the radiation calibration and spectral calibration of the radiance data product; the accuracy of imaging geometric parameters, including the central longitude and latitude, sensor altitude, observation zenith angle, azimuth, imaging time, and image spatial resolution; the accuracy of atmospheric parameters, including water vapor content and aerosol optical thickness; the accuracy of terrain parameters, including slope , Slope 、Cosine of solar incident angle , Sky Visibility Factor and the shadow occlusion factor .

[0009] Among them, the uncertainty of the radiance data product and imaging geometric parameters in step 3 is determined by the Class B evaluation method based on the relevant literature on hyperspectral payloads. The uncertainty of atmospheric parameters including the uncertainty of aerosol optical thickness and water vapor content is obtained by Monte Carlo simulation. The terrain parameters including slope , Slope 、Cosine of solar incident angle , Sky Visibility Factor Calculated by formula transfer method, shadow shading factor The propagation is simulated by Monte Carlo method, and the terrain parameters are calculated based on the input digital surface model. , is the elevation change rate of the terrain surface in the east-west and north-south directions, ( ) is the area around the target pixel The grid, is the grid pixel resolution, the solar zenith angle and azimuth angle are , , is the direction in hemispherical space, The pixel location is The maximum sky angle in the direction, The value is between 0 and The specific calculation and uncertainty transfer formula are as follows:

[0010]

[0011]

[0012]

[0013]

[0014]

[0015]

[0016]

[0017]

[0018]

[0019] The shadow shielding factor is calculated using a depth buffer algorithm by inputting the sun's illumination angle and digital surface model data. Since there is no analytical formula, the Monte Carlo method is used to simulate uncertainty propagation, and the uncertainty of all terrain parameters is obtained.

[0020] Among them, step 4 uses the uncertainty source of step 3 to propagate the uncertainty through the Monte Carlo method, and obtains the uncertainty of each radiation component through Monte Carlo calculation. , , , , , , And uncertainty related items, and then use the formula transfer method to calculate the formula for the composite standard uncertainty of reflectivity as follows:

[0021] In the formula For the The components and The correlation coefficient between the components.

[0022] Among them, in step 5, the eight uncertainty components of the above synthetic reflectivity are divided by the uncertainty of the total reflectivity, and a band-by-band comparison is performed to achieve quantitative sorting of uncertainty contributions and traceability of major uncertainties. The specific calculation is as follows: ; ; ; ; ; ; ; ; It is After obtaining the main uncertainty contribution components, the contribution ratio of each uncertainty component can be further traced back to the input uncertainty parameters, which can then guide parameter optimization and algorithm improvement, and ultimately improve the accuracy of surface reflectance data products.

[0023] Beneficial effects achieved by the present invention: (1) An uncertainty transfer model for joint correction of hyperspectral remote sensing atmosphere and terrain was established, and an uncertainty transfer analysis method for surface reflectance under non-uniform and undulating surfaces was proposed, thus achieving quality assessment of surface reflectance under non-uniform and undulating surfaces; (2) Obtain the quantitative uncertainty assessment results of the surface reflectance for each pixel and band, conduct the main uncertainty source analysis, realize the quantitative ranking of uncertainty contributions, locate the main sources of uncertainty, and guide the improvement of algorithms and data accuracy; (3) This application targets undulating and non-uniform surfaces and designs terrain parameter factors to reflect the characteristics of the undulating surface, including slope , Slope , and then introduce it into the quantitative analysis of uncertainty, realizing the reflectivity uncertainty analysis of non-flat ground; at the same time, the reflected radiation of surrounding objects is introduced to reflect the characteristics of non-uniform surface, and then the terrain parameters and non-uniform object parameters are considered at the same time. The two are coupled and act together on the uncertainty quantitative evaluation calculation of surface reflectivity, avoiding the idealization and simplification of the surface calculation in the existing technology, greatly improving the calculation accuracy of the reflectivity uncertainty of the real surface environment, and providing a reliable and solid theoretical basis for the subsequent improvement of the calculation of surface reflectivity; (4) This application takes into account a variety of factors that affect the surface reflectivity, and the above factors are not isolated and have no influence on each other. In order to calculate the uncertainty of the surface reflectivity more realistically, comprehensively and accurately, this application deeply explores the interaction between various factors. For example, terrain factors will affect the sky visibility factor and will also affect the shadow shielding factor. BRIEF DESCRIPTION OF THE DRAWINGS

[0024] Figure 1 Flow chart for uncertainty assessment of surface reflectance under complex surface; Figure 2 is the uncertainty tree diagram of surface reflectivity; Figure 3 The quantitative ranking results of the uncertainty contribution of vegetation pixels in surface reflectance of GF-5 hyperspectral satellite data are presented. DETAILED DESCRIPTION

[0025] The uncertainty analysis method of surface reflectivity based on the atmosphere-terrain coupling model is characterized by comprising the following steps: Step 1: Establish an atmosphere-terrain joint correction model to complete the value transfer of hyperspectral remote sensing radiance data to surface reflectance data; Specifically, the core formula for complex subsurface reflectivity inversion is used: ;

[0026] In this formula, is the direct solar radiation, is the atmospheric scattered radiation, is the ground-atmosphere coupled radiation, The parameter reflects the radiation reflected by the surrounding objects. It reflects the influence of the non-uniform surface. is the radiance caused by the atmospheric proximity effect, is the atmospheric radiation, is the radiance data, is the atmospheric upward transmittance, is the reflectivity of the target object, Represents the number of rows, columns and elevation respectively. is the wavelength, and the uncertainty introduced by the radiation transfer model error and the surface Lambertian assumption is represented by “+ 0”.

[0027] Step 2: Based on the value transfer process, the uncertainty is traced back to the uncertainty input of the model, which is the accuracy of the radiance product, the accuracy of the imaging geometric parameters, the accuracy of the terrain parameters, and the accuracy of the atmospheric parameters; Specifically, the uncertainty is traced and the uncertainty transfer link is drawn. Starting from the core formula of the joint correction of atmosphere and terrain, the uncertainty of each variable in the model is traced back to the input variables, such as Figure 2 Shown by: Figure 2 The uncertainty transfer chain shows that the uncertainty traceability parameters of reflectivity include: the radiometric calibration and spectral calibration accuracy of radiance data products; the accuracy of imaging geometric parameters, including central longitude and latitude, sensor altitude, observation zenith angle, azimuth, imaging time, and image spatial resolution; the accuracy of atmospheric parameters, including water vapor content and aerosol optical thickness; the accuracy of terrain parameters, including slope , Slope 、Cosine of solar incident angle , Sky Visibility Factor and the shadow occlusion factor .

[0028] Step 3: For each input variable obtained by tracing back in step 2, estimate the uncertainty by using the type A or type B assessment method and determine the probability distribution of these uncertain inputs; Specifically, the uncertainty of the radiance data product and imaging geometric parameters in step 3 is determined by the Class B assessment method based on the relevant literature on hyperspectral payloads and is considered to be normally distributed. Since atmospheric parameters and terrain parameters cannot be obtained through repeated measurements, both atmospheric parameters and terrain parameters are assessed through Class B uncertainty assessment and are considered to be normally distributed. The uncertainty of atmospheric parameters, including the uncertainty of aerosol optical thickness and water vapor content, is propagated through Monte Carlo simulation, and the uncertainty of terrain parameters, including slope , Slope 、Cosine of solar incident angle , Sky Visibility Factor Calculated by formula transfer method, shadow shading factor The propagation is simulated by Monte Carlo method, and the terrain parameters are calculated based on the input digital surface model. , is the elevation change rate of the terrain surface in the east-west and north-south directions, ( ) is the area around the target pixel The grid, is the grid pixel resolution, the solar zenith angle and azimuth angle are , , is the direction in hemispherical space, The pixel location is The maximum sky angle in the direction, The value is between 0 and The specific calculation and uncertainty transfer formula are as follows:

[0029]

[0030]

[0031]

[0032]

[0033]

[0034]

[0035]

[0036]

[0037]

[0038] in That is, the height uncertainty is evaluated by consulting the literature for Class B uncertainty, and is considered to be normally distributed. That is, the uncertainty of spatial resolution is evaluated by Class B uncertainty and is considered to be uniformly distributed. It can be assumed that the error limit is ,Right now The shadow shielding factor is calculated by using the depth buffer algorithm. It is calculated by inputting the solar illumination angle and the digital surface model data. There is no analytical formula. Therefore, the uncertainty of the shadow shielding factor is obtained by simulating the propagation through the Monte Carlo method. So far, the uncertainty of all terrain parameters has been obtained. The terrain parameter appears in the calculation formula , here reflects the quantification method of the impact of terrain factors on sky visibility factor, which has not been considered in the prior art. At the same time, terrain factors will also affect the calculation of shadow shielding factors, and this application takes them into comprehensive consideration.

[0039] Step 4: Use the uncertainty assessment results of the input variables obtained in step 3 to transfer and synthesize the uncertainty, and transfer the uncertainty of the formula transfer method to the links with analytical formulas and transfer the uncertainty of the Monte Carlo method to the links without analytical formulas, and calculate the uncertainty assessment results of the surface reflectance synthesis standard; Specifically, in step 4, under the uncertainty source of step 3, the uncertainty is propagated by Monte Carlo method, and the uncertainty of each radiation component is obtained by Monte Carlo calculation. , , , , , , And uncertainty related items, and then use the formula transfer method to calculate the formula for the composite standard uncertainty of reflectivity as follows:

[0040] In the formula For the The components and The correlation coefficient between the components.

[0041] Step 5: Analyze the uncertainty components in the synthetic surface reflectance, and perform a comparative analysis of the uncertainty contribution of each link on a band-by-band basis.

[0042] Specifically, the eight uncertainty components of the above synthetic reflectivity are divided by the uncertainty of the total reflectivity, and a band-by-band comparison is performed to achieve quantitative sorting of uncertainty contributions and traceability of major uncertainties. The specific calculation is as follows: ; ; ; ; ; ; ; ;

[0043] It is After obtaining the main uncertainty contribution components, the contribution ratio of each uncertainty component can be further traced back to the input uncertainty parameters, which can then guide parameter optimization and algorithm improvement, and ultimately improve the accuracy of surface reflectance data products.

[0044] In order to illustrate the proposed uncertainty analysis method of surface reflectance based on the atmosphere-terrain coupling model, the Guanting Reservoir area with undulating and non-uniform terrain was selected as the test area to carry out uncertainty transfer assessment. Figure 3 The quantitative ranking results of the uncertainty contribution of vegetation pixels in the surface reflectance of GF-5 hyperspectral satellite data are given. It can be seen that the uncertainty of the direct solar component and radiance data has a greater impact on the reflectance. In the future, the traceability parameters can be optimized and the parameter calculation algorithm can be improved to improve the accuracy of remote sensing image surface reflectance data products in a targeted and justified manner.

Claims

1. The uncertainty analysis method of surface reflectivity based on the atmosphere-terrain coupling model is characterized by: The following steps are involved: Step 1: Establish an atmosphere-terrain joint correction model to complete the value transfer of hyperspectral remote sensing radiance data to surface reflectance data; Step 2: Based on the value transfer process, the uncertainty is traced back to the uncertainty input of the model, which is the accuracy of the radiance product, the accuracy of the imaging geometric parameters, the accuracy of the terrain parameters, and the accuracy of the atmospheric parameters; Step 3: For each input variable obtained by tracing back in step 2, estimate the uncertainty by using the type A or type B assessment method and determine the probability distribution of these uncertain inputs; Step 4: Use the uncertainty assessment results of the input variables obtained in step 3 to transfer and synthesize the uncertainty, and transfer the uncertainty of the formula transfer method to the links with analytical formulas and transfer the uncertainty of the Monte Carlo method to the links without analytical formulas, and calculate the uncertainty assessment results of the surface reflectance synthesis standard; Step 5: Analyze the uncertainty components in the synthetic surface reflectance, and rank and compare the uncertainty contributions of each link band by band, and further trace the main uncertainties to determine the main source factors affecting the uncertainty of the surface reflectance, guide parameter optimization and algorithm improvement, and achieve improved data accuracy.

2. The method according to claim 1, characterized in that The core formula for inverting the reflectivity of the undulating and non-uniform surface in step 1 is: ； In this formula, is the direct solar radiation, is the atmospheric scattered radiation, is the ground-atmosphere coupled radiation, Radiation reflected by surrounding objects, is the radiance caused by the atmospheric proximity effect, is the atmospheric radiation, is the radiance data, is the atmospheric upward transmittance, is the reflectivity of the target object, Represents the number of rows, columns and elevation respectively. is the wavelength, and the uncertainty introduced by the radiation transfer model error and the surface Lambertian assumption is represented by "+0".

3. The method according to claim 1, characterized in that The sources of uncertainty traced in step 2 include: the accuracy of radiometric and spectral calibration of radiance data products; the accuracy of imaging geometric parameters, including central longitude and latitude, sensor altitude, observation zenith angle, azimuth, imaging time, and image spatial resolution; the accuracy of atmospheric parameters, including water vapor content and aerosol optical depth; the accuracy of terrain parameters, including slope , Slope 、Cosine of solar incident angle , Sky Visibility Factor and the shadow occlusion factor .

4. The method according to claim 1, characterized in that The uncertainty of the radiance data product and imaging geometric parameters in step 3 is determined by the Class B evaluation method based on the literature related to hyperspectral payloads. The uncertainty of atmospheric parameters including aerosol optical thickness and water vapor content is simulated by the Monte Carlo method. The terrain parameters including slope , Slope 、Cosine of solar incident angle , Sky Visibility Factor Calculated by formula transfer method, shadow shading factor The propagation is simulated by Monte Carlo method, and the terrain parameters are calculated based on the input digital surface model. , is the elevation change rate of the terrain surface in the east-west and north-south directions, ( ) is the area around the target pixel The grid, is the grid pixel resolution, the solar zenith angle and azimuth angle are , , is the direction in hemispherical space, The pixel location is The maximum sky angle in the direction, The value is between 0 and The specific calculation and uncertainty transfer formula are as follows: ； ； ； ； ； ； ； ； ； ； The shadow shielding factor is calculated using a depth buffer algorithm by inputting the sun's illumination angle and digital surface model data. Since there is no analytical formula, the Monte Carlo method is used to simulate uncertainty propagation, and the uncertainty of all terrain parameters is obtained.

5. The method according to claim 1, characterized in that Step 4: Based on the uncertainty source of step 3, the uncertainty is propagated by Monte Carlo method. The uncertainty of each radiation component is obtained by Monte Carlo calculation. , , , , , , And uncertainty related items, and then use the formula transfer method to calculate the formula for the composite standard uncertainty of reflectivity as follows: ； In the formula For the The components and The correlation coefficient between the components.

6. The method according to claim 1, characterized in that In step 5, the eight uncertainty components of the above synthetic reflectivity are divided by the uncertainty of the total reflectivity, and a band-by-band comparison is performed to achieve quantitative sorting of uncertainty contributions and traceability of major uncertainties. The specific calculation is as follows: ； ； ； ； ； ； ； ； It is After obtaining the main uncertainty contribution components, the contribution ratio of each uncertainty component can be further traced back to the input uncertainty parameters, which can then guide parameter optimization and algorithm improvement, and ultimately improve the accuracy of surface reflectance data products.

Images