Leaf area index retrieval method based on terrain correction

Abstract

The present application relates to the technical field of leaf area index inversion, solves the problem that the prior art is easy to cause serious redundancy of the characteristic vector, and affects the recognition degree of the model to the real difference of the complex mountainous area and reduces the inversion precision, and particularly relates to a leaf area index inversion method based on terrain correction, obtains multi-band remote sensing images of a target area and carries out pretreatment, generates a standard image sequence, synchronously obtains a digital elevation model and an observation geometric parameter, extracts terrain topological features and illumination state features according to pixels, and differentiates and schedules the pixels in the standard image sequence according to the terrain topological features and the illumination state features.The present application reduces the spectral abnormal fluctuation caused by shadow shielding and terrain undulation, reduces the uncertainty of the remote sensing signal caused by the multiple scattering effect of the complex terrain, reduces the radiation distortion error, improves the recognition degree of the model to the real difference of the complex mountainous area and the precision of the leaf area index inversion.

Description

Technical Field

[0001] This invention relates to the field of leaf area index inversion technology, and in particular to a leaf area index inversion method based on topographic correction. Background Technology

[0002] Leaf area index (LAI) is a key biophysical parameter characterizing forest canopy structure and ecosystem function. In mountainous regions, topographic effects often significantly reduce the accuracy of high-resolution LAI retrieval. Currently, large-scale LAI retrieval using multi-band remote sensing imagery has become the mainstream engineering method in this field. However, when the observation field of remote sensing technology is turned to vast and dramatically undulating mountainous forest areas, the complex topography causes severe unevenness in solar illumination angles and drastic deviations in the bidirectional reflectance distribution function signal. This greatly affects the accuracy of the remote sensing radiation signal, thus severely reducing the accuracy of high-resolution LAI retrieval.

[0003] Existing inversion systems typically incorporate terrain correction algorithms to overcome terrain effects. However, the spectral reflectance of existing technologies has usually already had its radiation distortion removed through terrain correction when inputting into the model. If terrain topology features are then used as input variables, it can easily lead to severe redundancy in the feature vectors. Consequently, existing ensemble learning networks lack crucial terrain ecological boundary explanatory power during the feature fusion stage, affecting the model's ability to identify the real differences in complex mountainous terrain and reducing the inversion accuracy of leaf area index. Summary of the Invention

[0004] To address the shortcomings of existing technologies, this invention provides a leaf area index inversion method based on terrain correction, which solves the technical problem that existing technologies easily lead to serious redundancy in feature vectors, thus affecting the model's ability to identify the real differences in complex mountainous terrain and reducing inversion accuracy.

[0005] To solve the above-mentioned technical problems, the present invention provides the following technical solution: a leaf area index inversion method based on terrain correction, the method comprising the following steps: S1. Acquire multi-band remote sensing images of the target area and preprocess them to generate a standard image sequence. Simultaneously acquire digital elevation models and observation geometric parameters, and extract topographic features and illumination features by pixel. S2. Based on the topographic features and illumination characteristics, perform differentiated scheduling of pixels in the standard image sequence; If a pixel meets the first illumination state and the terrain topology features do not exceed the preset slope threshold, it is assigned to the through branch, and the original reflectance is retained and output. If a pixel meets the second illumination state or the terrain topology features exceed the preset slope threshold, it is assigned to the compensation branch. If a pixel is in a transitional region of illumination, it is assigned to the dynamic scheduling branch; S3. In the compensation branch, the total incident radiation intercepted by the pixel is analyzed for radiative transmission attenuation based on the dual-view factor to separate the multi-source radiation components and reconstruct the topographically modulated observation reflectance to generate the corrected reflectance. S4. Synchronously aggregate the reflectance outputs of the direct branch, compensation branch, and dynamic scheduling branch, and extract the corresponding vegetation feature parameters. Perform feature-level fusion of reflectance, vegetation feature parameters, topographic features, and meteorological variables to construct a multi-dimensional joint feature set. S5. Input the multidimensional joint feature set into the pre-trained ensemble learning network, and use the ensemble learning network to perform a nonlinear mapping mechanism to coordinate the radiation distortion compensation features and terrain topology features to generate the leaf area index of the target region.

[0006] Furthermore, step S2 specifically includes the following steps: S21. Based on terrain topological features and observed geometric parameters, analyze the three-dimensional physical illumination state of pixels, and calculate the cosine of the solar incidence angle as an illumination state feature. The calculation formula is as follows: In the above formula, Indicates the characteristics of illumination conditions; These represent the slope parameter and aspect parameter in the topographic features, respectively. These represent the solar zenith angle and solar azimuth angle in the observation geometric parameters, respectively. S22. Configure an adaptive routing mechanism to perform differentiated scheduling of pixels in the standard image sequence based on illumination state characteristics and slope parameters; If the illumination state feature is greater than zero and the slope parameter is less than the preset slope threshold, then the pixel is determined to meet the first illumination state and has not exceeded the slope threshold. It is then assigned to the through branch, and the original reflectance is retained and output. If the illumination state feature is less than or equal to zero, or the slope parameter is greater than or equal to the preset slope threshold, then the pixel is determined to meet the second illumination state or exceed the slope threshold, and is assigned to the compensation branch. S23. If a pixel is in the transition region of mixed illumination, it is assigned to the dynamic scheduling branch for dynamic correction.

[0007] Furthermore, the dual-view factor includes a sky visibility factor and a terrain visibility factor; the multi-source radiation component includes a direct radiation component, a sky diffuse radiation component, and an adjacent terrain reflected radiation component.

[0008] Step S3 specifically includes the following steps: S31. In the compensation branch, based on the observed geometric parameters and atmospheric optical parameters, the direct radiation attenuation process of the pixels is analyzed, and the horizontal reference direct radiation component and the slope modulation direct radiation component are quantified respectively. The calculation formula is as follows: In the above formula, Indicates the direct radiation component of the horizontal reference; Indicates the direct radiation component modulated by the slope; Indicates solar irradiance at the top of the atmosphere; Indicates atmospheric optical thickness; Represents a binary function for shadow occlusion; These represent the slope parameter, aspect parameter, solar zenith angle, and solar azimuth angle, respectively. S32. Based on terrain topology features, quantify the sky visibility factor. Combining the sky visibility factor with the physical characteristics of anisotropic sky scattering, quantify the horizontal reference sky diffuse component and the slope-modulated sky diffuse component respectively. The calculation formula is as follows: In the above formula, Indicates the horizontal reference sky diffuse component; This indicates the slope modulated sky diffuse component. Represents the anisotropy factor. Indicates average sky radiance. Indicates the visibility factor of the sky; S33. Quantify the adjacent terrain reflectance components of the target pixel, and perform total radiation integration and inverse reconstruction on the various components obtained by analysis to generate the corrected reflectance.

[0009] Furthermore, the formula for calculating the corrected reflectivity is: In the above formula, Indicates corrected reflectivity; Indicates the original reflectivity; This represents the reflected radiation component from adjacent terrain under flat terrain conditions; This represents the reflected radiation component of adjacent terrain under slope conditions.

[0010] Furthermore, step S4 specifically includes the following steps: S41. Synchronously converge the reflectance outputs of the direct branch, compensation branch, and dynamic scheduling branch to generate a global reflectance matrix; S42. Based on the multispectral bands in the global reflectivity matrix, perform band operations to extract vegetation feature parameters; S43. Extract the topographic features and meteorological variables that are registered with the global reflectance matrix on a spatiotemporal scale. Then, splice the global reflectance matrix, vegetation feature parameters, topographic features and meteorological variables along the feature dimension to perform feature-level heterogeneous fusion and construct a multidimensional joint feature set.

[0011] Furthermore, the vegetation characteristic parameters include the atmospheric resistance vegetation index and the simple ratio index, calculated using the following formula: In the above formula, A simple ratio index for characterizing canopy structure, The atmospheric vegetation resistance index represents the ability to eliminate aerosol interference. This represents the reflectance of the corresponding pixel in the near-infrared band within the global reflectance matrix; This represents the reflectance of the corresponding pixel in the red band of the global reflectance matrix; This represents the reflectance of the blue light band of the corresponding pixel in the global reflectance matrix.

[0012] Furthermore, step S5 specifically includes the following steps: S51. The multi-dimensional joint feature set is used as an input vector and input into a pre-trained ensemble learning network. The ensemble learning network contains multiple parallel decision tree evaluators. Each decision tree evaluator independently performs nonlinear feature mapping based on the input vector and infers and generates the corresponding single-leaf area index prediction value. The calculation formula is as follows: In the above formula, This represents the predicted leaf area index value of a single leaf output by the i-th decision tree evaluator; Let represent the nonlinear mapping function constructed by the i-th decision tree evaluator; This represents the input multidimensional joint feature set, and , This represents the total number of decision tree evaluators in an ensemble learning network; S52. Using an ensemble aggregation mechanism, perform a global mean calculation on the predicted leaf area index values ​​of individual leaves generated by all decision tree evaluators. Then, by coordinating radiation distortion compensation features and terrain topology features at the decision layer, infer the leaf area index of the output target region. The calculation formula is as follows: In the above formula, This represents the leaf area index of the target region.

[0013] By employing the above technical solution, the present invention provides a leaf area index inversion method based on terrain correction, which has at least the following beneficial effects: 1. This invention utilizes an adaptive routing mechanism that includes a compensation branch and introduces a reverse reconstruction algorithm based on dual-viewpoint factors. This reduces the severe unevenness of solar illumination angles caused by complex terrain, reduces spectral anomalies caused by shadow occlusion and terrain undulations, reduces the uncertainty of remote sensing signals caused by multiple scattering effects in complex terrain, reduces radiation distortion errors, and improves the accuracy of leaf area index inversion in complex mountainous environments.

[0014] 2. When constructing a multidimensional joint feature set, this invention performs feature-level heterogeneous fusion of topographic features and distortion-free spectral reflectance, reducing multicollinearity interference caused by variable overlap, compensating for the lack of key topographic ecological boundary explanatory power in the nonlinear inference stage of the ensemble learning network, and improving the model's ability to identify the real differences in complex mountainous terrain.

[0015] 3. Based on the illumination characteristics and slope parameters of pixels, this invention performs differentiated scheduling of multi-band standard image sequences, avoiding the blind use of high-computational-complexity terrain correction models in flat, sunny slope areas, reducing ineffective computing power loss, lowering the computing power consumption and time cost in processing large amounts of high-resolution remote sensing data, and improving the overall efficiency of continuous leaf area index inversion. Attached Figure Description

[0016] The accompanying drawings, which are included to provide a further understanding of this application and form part of this application, illustrate exemplary embodiments of this application and are used to explain this application, but do not constitute an undue limitation of this application. In the drawings: Figure 1 This is a flowchart of the leaf area index inversion method of the present invention. Detailed Implementation

[0017] To make the above-mentioned objects, features, and advantages of the present invention more apparent and understandable, the present invention will be further described in detail below with reference to the accompanying drawings and specific embodiments. This will allow for a full understanding of how the present application uses technical means to solve technical problems and achieve technical effects, and to facilitate its implementation.

[0018] Current techniques for leaf area index (LAI) inversion often result in significant feature vector redundancy, affecting the model's ability to identify real-world differences in complex mountainous terrain and reducing inversion accuracy. To improve the model's ability to identify real-world differences in complex mountainous environments and the accuracy of LAI inversion, this invention proposes a terrain-corrected LAI inversion method, such as... Figure 1 As shown, the method includes the following steps: S1. By acquiring multi-band remote sensing images and digital elevation models of the target area, and simultaneously acquiring observation geometric parameters, this embodiment employs a time-series-based reflectivity adjustment (TRA) algorithm to coordinate the radiometric response differences of image sequences from different sensors (such as the Landsat series). Quality assessment (QA) bands are used in conjunction with time-series filters to mask and remove clouds, cloud shadows, and residual noise, eliminating radiometric response differences from heterogeneous sensors and atmospheric interference. This generates a standard image sequence with unified reference radiometric accuracy. Based on the digital elevation model and observation geometric parameters, the spatial normal vector and ray projection vector of the three-dimensional physical surface are analyzed at the pixel level to extract topographic features and illumination characteristics. Through pre-processing systematic deviation coordination and temporal noise removal, the radiometric reference drift caused by multi-source sensor heterogeneity and atmospheric cloud pollution is effectively reduced, lowering the uncertainty of the input data and improving the data stability of the inversion algorithm.

[0019] S2. Based on terrain topology and illumination characteristics, perform differentiated scheduling of pixels in the standard image sequence: If a pixel meets the first illumination condition and the terrain topology features do not exceed the preset slope threshold, it is assigned to the direct branch, and the original reflectance is retained and output. If a pixel meets the second illumination condition or the terrain topology features exceed the preset slope threshold, it is assigned to the compensation branch. If a pixel is in the transition area of ​​the illumination condition, it is assigned to the dynamic scheduling branch.

[0020] In this embodiment, as a preferred implementation, the specific steps are as follows: based on terrain topological features and observation geometric parameters, the three-dimensional physical illumination state of the pixels is analyzed, and the cosine value of the solar incidence angle is calculated as the illumination state feature. The calculation formula is: In the above formula, Indicates the characteristics of illumination conditions; These represent the slope parameter and aspect parameter in the topographic features, respectively. These represent the solar zenith angle and solar azimuth angle in the observation geometry parameters, respectively.

[0021] Based on illumination characteristics and slope parameters, pixels in the standard image sequence are differentiated and scheduled. During this differentiated scheduling process: like If the illumination state feature is greater than zero, that is, the slope parameter is less than the preset slope threshold, and the slope threshold is set to 15 degrees, then the pixel is determined to meet the first illumination state and has not exceeded the slope threshold. It is then assigned to the through branch, and the original reflectance is retained and output.

[0022] like If the illumination state feature is less than or equal to zero, or the slope parameter is greater than or equal to the preset slope threshold, then the pixel is determined to meet the second illumination state or exceed the slope threshold. The slope threshold is set to 15 degrees, and the pixel is assigned to the compensation branch.

[0023] If a pixel is in a transitional region under mixed lighting conditions, it is assigned to a dynamic scheduling branch for dynamic correction. Dynamic scheduling reduces data fragmentation and spatial inconsistency noise at the boundaries of different processing branches, thereby improving the overall data processing efficiency of the system.

[0024] Based on the physical characteristic that sky radiance decreases with increasing wavelength according to Rayleigh scattering theory, an average sky radiance weighting coefficient with exponentially decreasing values ​​is preset for different multispectral bands. Combined with the maximum sky angle calculated by using the elevation difference and horizontal distance of surrounding mountain peaks, and then integrating the dual-view factor, a nonlinear correction is applied to the multi-source radiation components received by the slope in different wavelength ranges, thereby eliminating the vegetation index deviation caused by the natural difference in the reflection compensation mechanism between the red band and the near-infrared band in adjacent terrain.

[0025] S3. In the compensation branch, based on the dual-view factor representing spatial occlusion (including the sky visibility factor and the terrain visibility factor), radiative transmission attenuation is calculated on the total incident radiation intercepted by the pixel to separate multi-source radiation components, including direct radiation components, sky diffuse radiation components, and adjacent terrain reflected radiation components. The observed reflectance modulated by terrain is then reconstructed to generate a distortion-free corrected reflectance. This avoids severe spectral distortion caused by uneven solar illumination and deep shadow effects in complex mountainous terrain, reduces optical noise caused by multiple reflections from adjacent terrain, and improves the fidelity and spatial continuity of canopy surface reflectance under complex undulating terrain.

[0026] In this embodiment, as a preferred implementation method, the implementation steps are as follows: In the compensation branch, based on the observed geometric parameters and atmospheric optical parameters, the direct radiation attenuation process of the pixels is analyzed, and the horizontal reference direct radiation component and the slope modulation direct radiation component are quantified separately. The calculation formula is as follows: In the above formula, Indicates the direct radiation component of the horizontal reference; Indicates the direct radiation component modulated by the slope; Indicates solar irradiance at the top of the atmosphere; Indicates atmospheric optical thickness; Represents a binary function for shadow occlusion; These represent the slope parameter, aspect parameter, solar zenith angle, and solar azimuth angle, respectively.

[0027] The sky visibility factor is quantified based on terrain topology features. Combining the sky visibility factor with the physical characteristics of anisotropic sky scattering, the horizontal reference sky diffuse component and the slope-modulated sky diffuse component are quantified separately. The calculation formula is as follows: In the above formula, Indicates the horizontal reference sky diffuse component; This indicates the slope modulated sky diffuse component. Represents the anisotropy factor. Indicates average sky radiance. This represents the visibility factor of the sky.

[0028] By quantizing the adjacent terrain reflectance components of the target pixel, and performing total radiation integration and inverse reconstruction on the various components obtained from the analysis, nonlinear terrain distortion is eliminated, and corrected reflectance is generated. The calculation formula is as follows: In the above formula, Indicates corrected reflectivity; Indicates the original reflectivity; This represents the reflected radiation component from adjacent terrain under flat terrain conditions; This represents the reflected radiation component of adjacent terrain under slope conditions.

[0029] S4. To overcome the impact of high canopy closure in summer, the extracted parameters not only include the simple ratio index (SR) but also specifically incorporate the Atmospheric Anti-Vegetation Index (ARVI) constructed using the blue light band. Reflectance outputs from the direct branch, compensation branch, and dynamic scheduling branch are simultaneously aggregated, and corresponding vegetation feature parameters are extracted. Reflectance, vegetation feature parameters, topographic features, and meteorological variables are fused at the feature level to construct a multi-dimensional joint feature set. This reduces multiple scattering noise and atmospheric aerosol interference under high canopy closure, minimizes spectral saturation in leaf area index inversion, and, by heterogeneously fusing distorted reflectance with topographic parameters, overcomes a situation generally considered by the industry to lead to severe redundancy.

[0030] In this embodiment, as a preferred implementation method, the implementation steps are as follows: The reflectance outputs from the direct branch, compensation branch, and dynamic scheduling branch are simultaneously aggregated to generate a global reflectance matrix. Based on the multispectral bands in the global reflectance matrix, band operations are performed to extract vegetation feature parameters. These vegetation feature parameters include the atmospheric resistance vegetation index and the simple ratio index, calculated using the following formula: In the above formula, A simple ratio index for characterizing canopy structure, The atmospheric vegetation resistance index represents the ability to eliminate aerosol interference. This represents the reflectance of the corresponding pixel in the near-infrared band within the global reflectance matrix; This represents the reflectance of the corresponding pixel in the red band of the global reflectance matrix; This represents the reflectance of the blue light band of the corresponding pixel in the global reflectance matrix.

[0031] Topographic features and meteorological variables that are registered with the global reflectance matrix at the spatiotemporal scale are extracted. The global reflectance matrix, vegetation feature parameters, topographic features and meteorological variables are spliced ​​along the feature dimension to perform feature-level heterogeneous fusion and construct a multi-dimensional joint feature set.

[0032] S5. The multi-dimensional joint feature set is input into a pre-trained ensemble learning network (such as a random forest architecture). Through a nonlinear mapping mechanism, the ensemble learning network coordinates the radiation distortion compensation features and the terrain topology features to generate the leaf area index (LAI) of the target region. This reduces the underestimation of high LAI values ​​in extreme seasons (such as summer) and steep, shady slope environments, and improves the accuracy of the high-resolution leaf area index inversion method in wide-area complex terrain applications.

[0033] In this embodiment, as a preferred implementation method, the specific implementation steps are as follows: The multidimensional joint feature set is used as an input vector and fed into a pre-trained ensemble learning network. The ensemble learning network contains multiple parallel decision tree evaluators, each of which independently performs nonlinear feature mapping based on the input vector, and respectively infers and generates the corresponding single-leaf area index prediction value. The calculation formula is as follows: In the above formula, This represents the predicted leaf area index value of a single leaf output by the i-th decision tree evaluator; Let represent the nonlinear mapping function constructed by the i-th decision tree evaluator; This represents the input multidimensional joint feature set, and , This represents the total number of decision tree evaluators in an ensemble learning network; Using an ensemble aggregation mechanism, the global mean of the leaf area index predictions of individual leaves generated by all decision tree evaluators is calculated. By coordinating radiation distortion compensation features and terrain topology features at the decision layer, the leaf area index of the output target region is inferred. The calculation formula is as follows: In the above formula, This represents the leaf area index of the target region.

[0034] This inversion method acquires and preprocesses multi-band remote sensing images of the target area to generate a standard image sequence. Simultaneously, it acquires digital elevation models and observational geometric parameters. Topographic and illumination features are extracted at the pixel level. Based on these features, pixels in the standard image sequence are differentiated and assigned to direct-access, compensation, and dynamic scheduling branches. In the compensation branch, radiative transmission attenuation analysis is performed on the total incident radiation intercepted by the pixels using a dual-viewpoint factor to separate multi-source radiation components. The topographically modulated observed reflectance is reconstructed to generate corrected reflectance. The reflectance outputs from the direct-access, compensation, and dynamic scheduling branches are simultaneously converged, and corresponding vegetation feature parameters are extracted. Reflectance, vegetation feature parameters, topographic features, and meteorological variables are fused at the feature level to construct a multi-dimensional joint feature set, which is then input into a pre-trained ensemble learning network. Through a nonlinear mapping mechanism, the ensemble learning network coordinates radiative distortion compensation features with topographic features to generate the leaf area index of the target area.

[0035] Those skilled in the art will understand that all or part of the steps in the methods of the above embodiments can be implemented by a program instructing related hardware. Therefore, this application can take the form of a completely hardware embodiment, a completely software embodiment, or an embodiment combining software and hardware aspects. Moreover, this application can take the form of a computer program product implemented on one or more computer-usable storage media (including but not limited to disk storage, CD-ROM, optical storage, etc.) containing computer-usable program code.

[0036] The various embodiments in this specification are described in a progressive manner, with each embodiment focusing on its differences from other embodiments. Similar or identical parts between embodiments can be referred to interchangeably. Since the above embodiments are substantially similar to the method embodiments, their descriptions are relatively simple; relevant parts can be referred to the descriptions of the method embodiments.

[0037] The above embodiments provide a detailed description of the present invention. Specific examples have been used to illustrate the principles and implementation methods of the present invention. The descriptions of the above embodiments are only for the purpose of helping to understand the method and core ideas of the present invention. At the same time, for those skilled in the art, there will be changes in the specific implementation methods and application scope based on the ideas of the present invention. Therefore, the content of this specification should not be construed as a limitation of the present invention.

Claims

1. A method for retrieving leaf area index based on terrain correction, characterized in that, The method includes the following steps: S1. Acquire multi-band remote sensing images of the target area and preprocess them to generate a standard image sequence. Simultaneously acquire digital elevation models and observation geometric parameters, and extract topographic features and illumination features by pixel. S2. Based on the topographic features and illumination characteristics, perform differentiated scheduling of pixels in the standard image sequence; If a pixel meets the first illumination state and the terrain topology features do not exceed the preset slope threshold, it is assigned to the through branch, and the original reflectance is retained and output. If a pixel meets the second illumination state or the terrain topology features exceed the preset slope threshold, it is assigned to the compensation branch. If a pixel is in a transitional region of illumination, it is assigned to the dynamic scheduling branch; S3. In the compensation branch, the total incident radiation intercepted by the pixel is analyzed for radiative transmission attenuation based on the dual-view factor to separate the multi-source radiation components and reconstruct the topographically modulated observation reflectance to generate the corrected reflectance. S4. Synchronously aggregate the reflectance outputs of the direct branch, compensation branch, and dynamic scheduling branch, and extract the corresponding vegetation feature parameters. Perform feature-level fusion of reflectance, vegetation feature parameters, topographic features, and meteorological variables to construct a multi-dimensional joint feature set. S5. Input the multidimensional joint feature set into the pre-trained ensemble learning network, and use the ensemble learning network to perform a nonlinear mapping mechanism to coordinate the radiation distortion compensation features and terrain topology features to generate the leaf area index of the target region.

2. The inversion method according to claim 1, characterized in that, Step S2 specifically includes the following steps: S21. Based on terrain topological features and observed geometric parameters, analyze the three-dimensional physical illumination state of pixels, and calculate the cosine of the solar incidence angle as an illumination state feature. The calculation formula is as follows: In the above formula, Indicates the characteristics of illumination conditions; These represent the slope parameter and aspect parameter in the topographic features, respectively. These represent the solar zenith angle and solar azimuth angle in the observation geometric parameters, respectively. S22. Configure an adaptive routing mechanism to perform differentiated scheduling of pixels in the standard image sequence based on illumination state characteristics and slope parameters; If the illumination state feature is greater than zero and the slope parameter is less than the preset slope threshold, then the pixel is determined to meet the first illumination state and has not exceeded the slope threshold. It is then assigned to the through branch, and the original reflectance is retained and output. If the illumination state feature is less than or equal to zero, or the slope parameter is greater than or equal to the preset slope threshold, then the pixel is determined to meet the second illumination state or exceed the slope threshold, and is assigned to the compensation branch. S23. If a pixel is in the transition region of mixed illumination, it is assigned to the dynamic scheduling branch for dynamic correction.

3. The inversion method according to claim 1, characterized in that, The dual-view factor includes the sky visibility factor and the terrain visibility factor; The multi-source radiation components include direct radiation components, diffuse sky radiation components, and radiation components reflected from adjacent terrain.

4. The inversion method according to claim 3, characterized in that, Step S3 specifically includes the following steps: S31. In the compensation branch, based on the observed geometric parameters and atmospheric optical parameters, the direct radiation attenuation process of the pixels is analyzed, and the horizontal reference direct radiation component and the slope modulation direct radiation component are quantified respectively. The calculation formula is as follows: In the above formula, Indicates the direct radiation component of the horizontal reference; Indicates the direct radiation component modulated by the slope; Indicates solar irradiance at the top of the atmosphere; Indicates atmospheric optical thickness; Represents a binary function for shadow occlusion; These represent the slope parameter, aspect parameter, solar zenith angle, and solar azimuth angle, respectively. S32. Based on terrain topology features, quantify the sky visibility factor. Combining the sky visibility factor with the physical characteristics of anisotropic sky scattering, quantify the horizontal reference sky diffuse component and the slope-modulated sky diffuse component respectively. The calculation formula is as follows: In the above formula, Indicates the horizontal reference sky diffuse component; This indicates the slope modulated sky diffuse component. Represents the anisotropy factor. Indicates average sky radiance. Indicates the visibility factor of the sky; S33. Quantify the adjacent terrain reflectance components of the target pixel, and perform total radiation integration and inverse reconstruction on the various components obtained by analysis to generate the corrected reflectance.

5. The inversion method according to claim 4, characterized in that, The formula for calculating the corrected reflectivity is: In the above formula, Indicates corrected reflectivity; Indicates the original reflectivity; This represents the reflected radiation component from adjacent terrain under flat terrain conditions; This represents the reflected radiation component of adjacent terrain under slope conditions.

6. The inversion method according to claim 1, characterized in that, Step S4 specifically includes the following steps: S41. Synchronously converge the reflectance outputs of the direct branch, compensation branch, and dynamic scheduling branch to generate a global reflectance matrix; S42. Based on the multispectral bands in the global reflectivity matrix, perform band operations to extract vegetation feature parameters; S43. Extract the topographic features and meteorological variables that are registered with the global reflectance matrix on a spatiotemporal scale. Then, splice the global reflectance matrix, vegetation feature parameters, topographic features and meteorological variables along the feature dimension to perform feature-level heterogeneous fusion and construct a multidimensional joint feature set.

7. The inversion method according to claim 6, characterized in that, The vegetation characteristic parameters include the atmospheric resistance vegetation index and the simple ratio index, calculated using the following formula: In the above formula, It is a simple ratio index. Indicates the atmospheric resistance vegetation index; This represents the reflectance of the corresponding pixel in the near-infrared band within the global reflectance matrix; This represents the reflectance of the corresponding pixel in the red band of the global reflectance matrix; This represents the reflectance of the blue light band of the corresponding pixel in the global reflectance matrix.

8. The inversion method according to claim 1, characterized in that, Step S5 specifically includes the following steps: S51. The multi-dimensional joint feature set is used as an input vector and input into a pre-trained ensemble learning network. The ensemble learning network contains multiple parallel decision tree evaluators. Each decision tree evaluator independently performs nonlinear feature mapping based on the input vector and infers and generates the corresponding single-leaf area index prediction value. The calculation formula is as follows: In the above formula, This represents the predicted leaf area index value of a single leaf output by the i-th decision tree evaluator; Let represent the nonlinear mapping function constructed by the i-th decision tree evaluator; This represents the input multidimensional joint feature set, and , This represents the total number of decision tree evaluators in an ensemble learning network; S52. Using an ensemble aggregation mechanism, perform a global mean calculation on the predicted leaf area index values ​​of individual leaves generated by all decision tree evaluators. Then, by coordinating radiation distortion compensation features and terrain topology features at the decision layer, infer the leaf area index of the output target region. The calculation formula is as follows: In the above formula, This represents the leaf area index of the target region.

Images