Topographic relief area tree height estimation method and system based on SAR scattering intensity correction

Abstract

The invention discloses a topographic relief region tree height estimation method and system based on SAR scattering intensity correction. The method comprises the steps that the irradiation area of SAR image pixels in a geographic coordinate system is converted into an SAR pixel coordinate system through inverse distance weighting; calculating the SAR backscattering intensity after the irradiation area correction by combining the image SAR resolution and the irradiation area of the SAR image pixel in the SAR pixel coordinate system; the actual measurement forest height is used as priori knowledge, forest height estimation is carried out based on a random volume-earth surface model, an optimal fitting coefficient of the random volume-earth surface model is obtained by using a nonlinear regression analysis method, and then the forest height of a research area is obtained. According to the method, the problem of estimation deviation caused by insufficient compensation of the terrain effect in the rugged terrains such as hills and mountainous regions in the existing forest height estimation technology based on single SAR backscattering is solved, and large-range and high-robustness forest height mapping can be realized.

Description

Technical Field

[0001] The present invention belongs to the technical field of radar remote sensing and forest parameter inversion, and particularly relates to a method and system for estimating tree height in a terrain undulation area based on SAR scattering intensity correction. Background Art

[0002] Forest height is one of the most basic structural parameters of the forest ecosystem, and its accurate estimation is of great significance for comprehensively understanding the forest structure and function and realizing the sustainable management of natural resources. As an active microwave remote sensing system, Synthetic Aperture Radar (SAR) can penetrate the forest canopy and obtain forest vertical structure information. Currently, the widely studied Polarimetric Interferometric SAR (PolInSAR) and SAR Tomography (TomoSAR) technologies based on the interference characteristics of SAR respectively introduce multi-polarization and multi-orbit observation information to extract more refined forest vertical structure distributions, and then achieve high-precision forest estimation. However, they rely on repeat orbit interferometric measurements and are susceptible to temporal decoherence in forest-covered areas, and even suffer from severe coherence loss, resulting in estimation failure.

[0003] Different from the above two technologies, considering the significant correlation between the SAR backscattering intensity and biophysical parameters such as forest height and biomass, the technology based on single-pass SAR backscattering intensity can establish statistical empirical models and semi-empirical models for different forest scenarios only through single-pass observations, and then can carry out forest height estimation work. Moreover, it has low computational complexity, is suitable for large-scale spatial processing, and has strong multi-source data fusion ability. However, due to the side-looking imaging characteristics of the SAR sensor, the measured signal is easily affected by terrain undulations. The same ground object will show different signal characteristics on the SAR image due to different local terrains. Therefore, when using SAR for ground object interpretation or quantitative parameter inversion, the influence of terrain is an inevitable problem. Especially for mountain forests with complex terrains, the accuracy of the statistical empirical model based on single-pass SAR scattering intensity significantly decreases in complex terrain areas, specifically manifested as: 1) the change of irradiation area, and it is difficult for traditional homomorphic correction methods (such as radiation correction based on local incident angle) to accurately reflect the abnormal relationship between geography and SAR slant range geometry; 2) the angle change effect, and terrain changes will change the forest scattering mechanism and microwave penetration path. However, the existing methods have insufficient compensation effects on terrain influence and cannot take into account both the change of irradiation area and the angle variation effect at the same time. To sum up, it is difficult for the existing technologies to achieve high-precision and robust forest height estimation in hilly areas, and there is an urgent need for a new method that combines terrain slope correction and scattering mechanism optimization. Summary of the Invention

[0004] To solve the problem that existing forest height estimation techniques based on single - frame SAR backscattering result in estimation biases in hilly, mountainous, and other undulating terrains due to insufficient compensation for terrain effects, and to achieve large - scale and highly robust forest height mapping, in view of the deficiencies of the existing technologies, the present invention proposes a method for estimating tree height in undulating terrain based on SAR scattering intensity correction, including the following steps: Step 1: Obtain the HV - polarized SAR single - look complex image, digital elevation model, and forest / non - forest classification data covering the study area; Step 2: Calculate the irradiation area of SAR image pixels in the geographic coordinate system using the SAR system orbit parameters, SAR image resolution, and DEM data; Step 3: According to the conversion relationship between the geographic coordinate system and the SAR pixel coordinate system, obtain the irradiation area of SAR image pixels in the SAR pixel coordinate system using inverse distance weighting; Step 4: Combine the SAR resolution of the image and the irradiation area of SAR image pixels in the SAR pixel coordinate system to calculate the SAR backscattering intensity after irradiation area correction; Step 5: Using the measured forest height as prior knowledge and combining the SAR backscattering intensity after irradiation area correction, obtain the best - fitting coefficients of the random volume - surface model using the nonlinear regression analysis method, and then obtain the forest height of the study area.

[0005] Further, in the geographic coordinate system in Step 2, the calculation formula for the irradiation area of each SAR image pixel is as follows: (1) In the formula, is the irradiation area of the SAR image pixel in the geographic coordinate system, is the azimuth resolution of the SAR image, is the slant - range resolution of the SAR image, is the cosine of the projection angle, obtained by the dot - product of the terrain surface normal vector and the radar slant - range plane normal vector: (2) (3) (4) In the formula, is the terrain surface normal vector, is the change rate of elevation in the azimuth direction of the SAR image, is the change rate of elevation in the ground - range direction of the SAR image, is the radar slant - range plane normal vector, is the velocity vector of the sensor, is the slant - range vector.

[0006] Further, in step 3, the irradiation area of the SAR image pixels in the geographic coordinate system is transformed to the SAR pixel coordinate system by using the SAR geocoding lookup table, and the irradiation area in the radar pixel coordinates is assigned to the SAR image pixels and will be deformed. The irradiation area in the radar pixel coordinates is assigned to the 4 neighboring pixels around the SAR image pixels ; Suppose any pixel 、 、 、 in the SAR image will be assigned N irradiation areas in the geographic coordinate system. The corresponding 、 、 、 ; Suppose any pixel in the SAR image will be assigned N irradiation areas in the geographic coordinate system. The corresponding irradiation areas in the geographic coordinate system are weighted and summed to obtain the irradiation area of pixel in the radar pixel coordinates : (5) (6) (7) (8) In the formula, is the k-th irradiation area in the geographic coordinate system, is the weight of the k-th irradiation area, is the SRA image pixel point corresponding to the k-th irradiation area, is the distance between and the surrounding integer-indexed pixels. floor represents rounding down, ceil represents rounding up, represents the azimuth integer index obtained by rounding down

[0007]

[0007] Further, the SAR backscattering intensity after irradiation area correction in step 4 is calculated by the following formula: (9) In the formula, is the SAR backscattering intensity after irradiation area correction, is the azimuth resolution of the SAR image, is the resolution in the slant range direction of the SAR image, is the SLC intensity value after radiometric calibration, is the SAR image pixel is the irradiation area at the SAR pixel coordinates.

[0008] Furthermore, the random volume-surface model in step 5 is expressed as: (10) In the formula, represents the SAR backscattering intensity, is the height value in the direction perpendicular to the inclined surface, is the surface elevation in the direction perpendicular to the inclined surface, is the actual forest height, , represents the "pseudo" forest height caused by the terrain range slope , is the scattering intensity of the volume scattering mechanism, is the scattering intensity of the double-bounce scattering mechanism, is the extinction coefficient, is the local incidence angle.

[0009] Substitute the SAR backscattering intensity after correcting the irradiation area obtained in step 4 into formula (10) to establish its relationship with the forest height, that is: (11) Where: (12) In the formula, is the predictor variable, is the response variable, is the fitting coefficient.

[0010] Using the geocoding lookup table, convert the SAR backscattering intensity after correcting the irradiation area to the geographic coordinate system. Use the forest / non-forest classification data to mask the non-forest areas in the backscattering intensity image in the geographic coordinate system. Calculate the "pseudo" forest height value caused by the terrain range slope based on some measured forest height values, and use it as a prior value. Combine the backscattering intensity image of the forest area, and use the nonlinear regression analysis technique to obtain the best fitting coefficient of formula (11). Substitute the best fitting coefficient and the SAR backscattering intensity after correcting the irradiation area of the study area into formula (11), and then obtain the forest height of the study area.

[0011] The present invention also provides a system for estimating tree height in a terrain undulating area based on SAR scattering intensity correction, which is used to implement the method for estimating tree height in a terrain undulating area based on SAR scattering intensity correction as described above.

[0012] Furthermore, the method comprises a processor and a memory, wherein the memory is used to store program instructions, and the processor is used to call the stored instructions in the memory to execute the above-mentioned method for estimating tree height in a terrain undulating area based on SAR scattering intensity correction.

[0013] Alternatively, it comprises a readable storage medium having a computer program stored thereon, and when the computer program is executed, the method for estimating tree height in a terrain undulating area based on SAR scattering intensity correction as described above is implemented.

[0014] Compared with the prior art, the present invention has the following advantages: 1) The present invention takes into account the difference in radar wave irradiation area caused by changes in terrain slope, as well as the angular dependence of SAR forest backscatter correction in the vertical direction. On the basis of the amplitude correction related to the ground irradiation area in the first stage, the backscattering semi-empirical method with slope correction of the ground scattering mechanism is further used in the second stage to extract the forest height, thereby weakening the SAR signal distortion in the forest area caused by terrain undulations, and suppressing the interference of terrain slope changes on the SAR forest height estimation performance, thereby improving the applicability in complex terrain areas; 2) The present invention uses the spatially sparsely distributed LiDAR high-precision tree height data as the key "ground truth" anchor point, implements the collaborative fusion strategy of "SAR continuous surface + LiDAR discrete point", and can generate enhanced forest height products with both high LiDAR point accuracy and wide SAR surface coverage; 3) The present invention only relies on single-view SAR backscatter intensity images, avoiding the complex interferometric measurement analysis process, and is expected to provide key algorithm support for the production of commercial, large-scale, and high-precision forest height products; 4) The present invention breaks the limitation of traditional methods in areas with significant terrain undulations, and can effectively improve the forest height mapping effect and spatial continuity in mountainous and hilly areas. BRIEF DESCRIPTION OF THE DRAWINGS

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

[0016] Figure 1This is the flowchart of the tree height estimation method in the terrain undulation area based on SAR scattering intensity correction according to the embodiments of the present invention.

[0017] Figure 2 This is the comparison diagram of the airborne SAR backscattering intensity before and after correction according to the embodiments of the present invention, where (a) is the airborne SAR backscattering intensity map before correction, and (b) is the airborne SAR backscattering intensity map after correction.

[0018] Figure 3 This is the statistical distribution histogram of the airborne SAR backscattering intensity before and after correction according to the embodiments of the present invention.

[0019] Figure 4 This is the relationship diagram of the response variable and the predictor variable in the regression model according to the embodiments of the present invention, where (a) is the relationship diagram of the predictor variable and the response variable before correction, and (b) is the relationship diagram of the predictor variable and the response variable after correction.

[0020] Figure 5 This is the comparison diagram of the forest height product and the LiDAR forest height reference data before and after correction according to the embodiments of the present invention, where (a) is the forest height product before correction, (b) is the forest height product after correction, (c) is the LiDAR forest height reference data, and (d) is the RMSE distribution of the forest height before and after correction and the LiDAR forest height reference data. Detailed implementation manners

[0021] To make the objectives, technical solutions and advantages of the present invention clearer, the technical solutions of the present invention will be further described below with reference to the accompanying drawings and embodiments. Obviously, the described embodiments are some but not all of the embodiments of the present invention. All other embodiments obtained by those of ordinary skill in the art based on the embodiments of the present invention without creative efforts shall fall within the protection scope of the present invention.

[0022] Embodiment 1 As Figure 1 shown, the embodiments of the present invention provide a tree height estimation method in the terrain undulation area based on SAR scattering intensity correction, including the following steps: Step 1, obtain the HV polarization SAR single-look complex (SLC) image, digital elevation model (DEM) and forest / non-forest classification data covering the study area.

[0023] The SAR data covering the study area was collected by the UAVSAR (Uninhabited Aerial Vehicle Synthetic Aperture Radar), an airborne synthetic aperture radar system developed by the Jet Propulsion Laboratory (JPL) of the National Aeronautics and Space Administration (NASA). It operates in the L-band (wavelength approximately 24 cm), has full-polarization observation capabilities, can penetrate the vegetation canopy, and is sensitive to surface and subsurface structures. The digital elevation model uses the Copernicus Digital Elevation Model without forests and buildings (FABDEM_V1-2), and the forest / non-forest classification data uses the National Land Cover Dataset (NLCD).

[0024] Step 2: Calculate the irradiation area of each SAR image pixel in the geographic coordinate system using the SAR system orbital parameters, SAR image resolution, and DEM data.

[0025] In the geographic coordinate system, the formula for calculating the irradiation area of each SAR image pixel is as follows: (1) where is the irradiation area of the SAR image pixel in the geographic coordinate system, is the azimuth resolution of the SAR image, is the range resolution of the SAR image, is the cosine of the projection angle, obtained by the dot product of the terrain surface normal vector and the radar range plane normal vector: (2) (3) (4) where is the terrain surface normal vector, is the rate of change of elevation in the azimuth direction of the SAR image, is the rate of change of elevation in the range direction of the SAR image, is the radar range plane normal vector, is the velocity vector of the sensor, is the range vector.

[0026] Step 3: According to the conversion relationship between the geographic coordinate system and the SAR pixel coordinate system, use inverse distance weighting to obtain the irradiation area of the SAR image pixel in the SAR pixel coordinate system.

[0027] The geographic coordinate system regards the Earth as an approximate sphere, with the geocenter as the origin, the equatorial plane as the reference plane, and the plane where the prime meridian lies as the starting meridian plane. The coordinates of each point on the Earth are determined by the intersection line of these two planes and the Earth's axis of rotation. Among them, longitude represents the angle between a certain point and the plane where the prime meridian lies, and latitude represents the angle between a certain point and the equatorial plane.

[0028] The SAR pixel coordinate system takes the upper left corner of the SAR image as the origin, the horizontal direction as the x-axis, with the right direction as the positive direction, and the vertical direction as the y-axis, with the downward direction as the positive direction. Each pixel point has a unique coordinate value (slant range coordinate, azimuth coordinate) in this coordinate system, which is used to represent its position in the radar image.

[0029] The conversion relationship between SAR image pixels in the geographic coordinate system and the SAR pixel coordinate system can be obtained through the SAR geocoding lookup table. SAR image pixels The irradiated area of SAR image pixels in the geographic coordinate system will be distorted after being converted to the SAR pixel coordinate system. Therefore, the irradiated area in the radar pixel coordinates is assigned to the SAR image pixels The 4 neighboring pixels around 、 、 、 , among which, 、 、 、 . Assume that any pixel on the SAR image will be assigned N irradiated areas in the geographic coordinate system. The corresponding irradiated areas in the geographic coordinate system are weighted and summed to obtain the irradiated area of pixel in the radar pixel coordinates : (5) (6) (7) (8) In the formula, is the k-th irradiated area in the geographic coordinate system, is the weight of the k-th irradiated area, is the SRA image pixel corresponding to the k-th irradiated area, is the distance between and the surrounding integer-indexed pixels. floor represents rounding down, ceil represents rounding up, represents the azimuth integer index obtained by rounding down respectively represent The azimuth integer index obtained by rounding up, indicating the azimuth integer index obtained by rounding down, indicating the azimuth integer index obtained by rounding up.

[0030] Step 4: Combine the SAR resolution of the image and the irradiation area of the SAR image pixels in the SAR pixel coordinate system to calculate the SAR backscattering intensity after irradiation area correction.

[0031] The SAR backscattering intensity after irradiation area correction is calculated by the following formula: (9) In the formula, is the SAR backscattering intensity after irradiation area correction, is the resolution in the azimuth direction of the SAR image, is the resolution in the slant range direction of the SAR image, is the SLC intensity value after radiometric calibration, is the SAR image pixel irradiation area in the SAR pixel coordinates.

[0032] Figure 2 is the comparison chart of the airborne SAR backscattering intensity before and after irradiation area correction, where (a) is the airborne SAR backscattering intensity map before irradiation area correction, and (b) is the airborne SAR backscattering intensity map after irradiation area correction. From Figure 2 (a), it can be seen that the area facing the radar incident direction appears brighter than the area facing away from the incident direction. This is because the irradiation area is calculated based on the geoid model, ignoring the non-uniform correspondence between the SAR slant range image and the geographic coordinate space, resulting in insufficient terrain correction. In addition, the viewing angle range of the airborne SAR is very wide, about 20° - 70°. The perspective shortening effect of the slope facing the radar at close range is more obvious and appears brighter, while the slope facing away from the radar at a long distance is relatively less affected. Therefore, there are more shadow areas on the rear slope. From Figure 2 (b), it can be seen that after irradiation area correction, the difference in radar scattering intensity between the upslope and the downslope becomes smaller, and the texture features are more uniform, especially significantly improved at close range. Figure 3 shows the statistical distribution of the airborne SAR backscattering intensity before and after correction. The variances of the backscattering intensity before and after correction are 8.45 and 4.35 respectively, and the variance is reduced by 48.5%. This indicates that after the correction operation, the brightness difference of the image pixel points becomes smaller.

[0033] Step 5: Using the measured forest height as prior knowledge, combining with the SAR backscattering intensity after irradiance area correction, and using the non-linear regression analysis method to obtain the best fitting coefficients of the random volume-ground model, and then obtaining the forest height of the study area.

[0034] In the area with terrain undulation, the change of slope will change the penetration path of microwave signal in the forest scene, thus changing the scattering mechanism. At this time, the scattering process of the forest scene can be described by the Slope-based Random Volume over Ground (Slope-RVoG) model, which can be expressed as: (10) In the formula, represents the SAR backscattering intensity, is the height value along the direction perpendicular to the inclined ground surface, is the ground elevation along the direction perpendicular to the inclined ground surface, is the actual forest height, represents the "pseudo" forest height caused by the terrain range slope , is the scattering intensity of the volume scattering mechanism, is the scattering intensity of the double-bounce scattering mechanism, is the extinction coefficient, which describes the one-way power loss of microwave propagation in the forest canopy, is the local incident angle.

[0035] Substitute the SAR backscattering intensity after irradiance area correction obtained in Step 4 into formula (10) to establish its relationship with the forest height, that is: (11) Among them: (12) In the formula, is the predictor variable, is the response variable, is the fitting coefficient.

[0036] Using a geocoding lookup table, the SAR backscattering intensity after irradiance area correction is converted to the geographic coordinate system. Using GAMMA software, based on some measured forest height values, the "pseudo" forest height value caused by the topographic range slope is calculated. After masking the non-forest areas in the backscattering intensity image using forest / non-forest classification data, combining some "pseudo" forest height values as prior knowledge, the best fitting coefficients of formula (11) are obtained using non-linear regression analysis technology. Substituting the best fitting coefficients and the SAR backscattering intensity after irradiance area correction in the study area into formula (11), the forest height of the study area is then obtained.

[0037] In this embodiment, the goodness of fit of the regression model is evaluated by describing the degree of dispersion of the observed data points around the regression line - the Standard Error of Regression (SER): (13) In the formula, is the actual response variable value, is the model predicted value, is the number of samples, is the number of predictor variables in the model.

[0038] The smaller the SER, the higher the goodness of fit of the regression equation to the data and the stronger the prediction ability. After determining the parameters of the non-linear regression model, the model transforms from a generalized framework to a model with specific prediction ability for the current study area and data characteristics. By applying it to the entire study area, a spatial distribution map of forest height covering all forest pixels can be generated. In this embodiment, the model fitting coefficient corresponding to the minimum SER is taken as the optimal fitting coefficient.

[0039] Figure 4 shows the relationship between the response variable and the predictor variable in the regression model and the model fitting results, where the predictor variable in (a) is the SAR backscattering intensity without irradiance area correction, and the predictor variable in (b) is the SAR backscattering intensity after irradiance area correction. From Figure 4 it can be seen that the dynamic range of the corrected predictor variable and the response variable is smaller and the sensitivity is higher. The SER index decreases from 4.31 to 3.40, indicating that the corrected regression equation has a higher goodness of fit to the data and better prediction ability.

[0040] Figure 5 In (a) and (b) respectively represent the forest height results estimated using the SAR backscattering intensity without irradiance area correction and the SAR backscattering intensity after irradiance area correction. Figure 5 In (c) is the LiDAR forest height reference data. From Figure 5As can be seen from (a), (b), and (c) in the figure, the forest height estimation results without considering terrain undulations show a trend of underestimation at close distances and overestimation at far distances. However, when using the forest height estimation strategy considering terrain undulations proposed in the present invention, the obtained results are consistent with the spatial distribution of LiDAR forest height. To further carry out quantitative evaluation, the LiDAR forest height values are binned and statistically analyzed in units of 1 m, and the corresponding root mean square error (RMSE) is calculated as the forest height evaluation index. As can be seen from Figure 5 (d) in the figure, the error trends estimated before and after correction are similar at different forest height levels. Compared with the overall average root mean square error of 6.20 m before correction, the overall average root mean square error after correction is 4.49 m, and the estimation accuracy of forest height is improved by 27.6%.

[0041] Embodiment 2 Based on the same inventive concept, the present invention also provides a tree height estimation system for terrain undulation areas based on SAR scattering intensity correction, including a processor and a memory. The memory is used to store program instructions, and the processor is used to call the program instructions in the memory to execute the above-mentioned tree height estimation method for terrain undulation areas based on SAR scattering intensity correction.

[0042] Embodiment 3 Based on the same inventive concept, the present invention also provides a tree height estimation system for terrain undulation areas based on SAR scattering intensity correction, including a readable storage medium. A computer program is stored on the readable storage medium, and when the computer program is executed, it implements the above-mentioned tree height estimation method for terrain undulation areas based on SAR scattering intensity correction.

[0043] Specifically in implementation, the method proposed by the technical solution of the present invention can be automatically run by those skilled in the art using computer software technology. The system device for implementing the method, such as a computer-readable storage medium storing the corresponding computer program of the technical solution of the present invention and a computer device including running the corresponding computer program, should also be within the protection scope of the present invention.

[0044] The specific embodiments described herein are merely illustrative of the spirit of the present invention. Those skilled in the art of the present invention can make various modifications or supplements or use similar methods to replace the described specific embodiments, but will not deviate from the spirit of the present invention or exceed the scope defined by the appended claims.

Claims

1. A method for estimating tree height in terrain undulation areas based on SAR scattering intensity correction, characterized in that, It includes the following steps: Step 1: Obtain the HV polarized SAR single-look complex image, digital elevation model, and forest / non-forest classification data covering the study area; Step 2: Calculate the irradiation area of SAR image pixels in the geographic coordinate system using the SAR system orbit parameters, SAR image resolution, and DEM data; Step 3: According to the conversion relationship between the geographic coordinate system and the SAR pixel coordinate system, use inverse distance weighting to obtain the irradiation area of SAR image pixels in the SAR pixel coordinate system; Step 4: Combine the SAR resolution of the image and the irradiation area of SAR image pixels in the SAR pixel coordinate system to calculate the SAR backscattering intensity after irradiation area correction; Step 5: Using the measured forest height as prior knowledge, combine the SAR backscattering intensity after irradiation area correction, and use the non-linear regression analysis method to obtain the best fitting coefficients of the random volume-surface model, and then obtain the forest height of the study area.

2. The method for estimating tree height in a terrain undulation area based on SAR scattering intensity correction according to claim 1, wherein: In the geographic coordinate system in Step 2, the calculation formula for the irradiation area of each SAR image pixel is as follows: （1） Wherein, is the irradiation area of the SAR image pixel in the geographic coordinate system, is the resolution in the azimuth direction of the SAR image, is the resolution in the slant range direction of the SAR image, is the cosine of the projection angle, obtained by the dot product of the terrain surface normal vector and the radar slant range plane normal vector: （2） （3） （4） In the formula, is the terrain surface normal vector, is the change rate of elevation in the azimuth direction of the SAR image, is the change rate of elevation in the range direction of the SAR image, is the radar slant range plane normal vector, is the velocity vector of the sensor, is the slant range vector.

3. The method for estimating tree height in a terrain undulation area based on SAR scattering intensity correction according to claim 1, characterized in that: In step 3, the SAR geocoding lookup table is used to convert the irradiation area of the SAR image pixels in the geodetic coordinate system to the SAR pixel coordinate system, and the irradiation area in the radar pixel coordinates is assigned to the SAR image pixels and its 4 surrounding neighboring pixels , , , , where , , , ; assume that any pixel on the SAR image is assigned N irradiation areas in the geodetic coordinate system, and the corresponding irradiation areas in the geodetic coordinate system are weighted and summed to obtain the irradiation area of the pixel in the radar pixel coordinates: （5） wherein, is the k-th irradiation area in the geographic coordinate system, is the weight of the k-th irradiation area.

4. The method for estimating tree height in a terrain undulation area based on SAR scattering intensity correction according to claim 3, wherein: The weight in Step 3 is calculated as follows: （6） （7） （8） Wherein, is the SRA image pixel corresponding to the k-th irradiation area, is the distance between and the surrounding integer-indexed pixels, floor represents rounding down, ceil represents rounding up, represents the azimuth integer index obtained by rounding down , respectively represent the azimuth integer index obtained by rounding up , represents the azimuth integer index obtained by rounding down , represents the azimuth integer index obtained by rounding up .

5. The method for estimating tree height in a terrain undulation area based on SAR scattering intensity correction according to claim 1, wherein: The SAR backscattering intensity after irradiation area correction in Step 4 is calculated by the following formula: （9） In the formula, is the SAR backscattering intensity after irradiation area correction, is the resolution in the azimuth direction of the SAR image, is the resolution in the slant range direction of the SAR image, is the SLC intensity value after radiometric calibration, is the SAR image pixel and is the irradiation area under the SAR pixel coordinates.

6. The method for estimating tree height in a terrain undulation area based on SAR scattering intensity correction according to claim 1, wherein: The random volume-surface model in Step 5 is expressed as: （10） In the formula, represents the SAR backscattering intensity, is the height value along the direction perpendicular to the inclined ground surface, is the ground elevation along the direction perpendicular to the inclined ground surface, is the actual forest height, , representing the "pseudo" forest height caused by the terrain distance to the slope , is the scattering intensity of the volume scattering mechanism, is the scattering intensity of the even-order scattering mechanism, is the extinction coefficient, is the local incident angle.

7. The method for estimating tree height in a terrain undulation area based on SAR scattering intensity correction according to claim 6, characterized in that: In Step 5, substitute the SAR backscattering intensity after irradiation area correction obtained in Step 4 into formula (10) to establish its relationship with the forest height, that is: （11） Where: （12） Wherein, is a predictor variable, is a response variable, is a fitting coefficient.

8. The method for estimating tree height in a terrain undulation area based on SAR scattering intensity correction according to claim 7, characterized in that: In Step 5, use the geocoding lookup table to convert the SAR backscattering intensity after irradiation area correction to the geographic coordinate system, mask the non-forest areas in the backscattering intensity image in the geographic coordinate system using the forest / non-forest classification data, calculate the "pseudo" forest height value caused by the terrain range slope according to some measured forest height values, use it as a prior value, combine the backscattering intensity image of the forest area, and use the non-linear regression analysis technique to obtain the best fitting coefficients of formula (11). Substitute the best fitting coefficients and the SAR backscattering intensity after irradiation area correction of the study area into formula (11), and then obtain the forest height of the study area.

9. A tree height estimation system for terrain undulation areas based on SAR scattering intensity correction, characterized in that, It includes a processor and a memory. The memory is used to store program instructions, and the processor is used to call the program instructions in the memory to execute a method for estimating tree height in a terrain undulation area based on SAR scattering intensity correction as described in any one of claims 1-8.

10. A tree height estimation system for undulating terrain areas based on SAR scattering intensity correction, characterized in that, It includes a readable storage medium, and a computer program is stored on the readable storage medium. When the computer program is executed, it implements a method for estimating tree height in a terrain undulation area based on SAR scattering intensity correction as described in any one of claims 1-8.

Images