A landslide hidden danger identification method based on SBAS-InSAR and optical remote sensing

Abstract

The present application relates to the technical field of geological disaster monitoring and identification, and particularly relates to a landslide hidden danger identification method based on SBAS-InSAR and optical remote sensing. The method comprises the following steps: S1, obtaining SAR images within the scope of the research area, and using SBAS-InSAR technology to process the SAR images, and then using the standard deviation method to divide the unstable area to obtain the initial unstable area set of the research area; S2, establishing landslide element feature interpretation marks, obtaining optical remote sensing images of the initial unstable area set, and obtaining a landslide hidden danger area set after eliminating the human activity area, and then identifying the historical deformation area landslide hidden danger and the deformation area landslide hidden danger to obtain a landslide hidden danger list within the scope of the research area. The optical remote sensing technology is introduced to perform secondary screening on the unstable area circled by SBAS-InSAR, the human activity area is eliminated, and the historical deformation area landslide hidden danger and the deformation area landslide hidden danger are accurately identified, which significantly improves the accuracy of landslide hidden danger identification.

Description

Technical Field

[0001] This invention relates to the field of geological disaster monitoring and identification technology, specifically to a landslide hazard identification method based on SBAS-InSAR and optical remote sensing, which is applicable to the accurate detection of landslide hazards in complex terrain areas such as high mountains and canyons. Background Technology

[0002] Landslides, as a common geological hazard, pose a serious threat to the ecological environment, infrastructure, and the safety of people's lives and property. Especially in areas with crisscrossing high mountains and deep valleys, such as the Lijiaxia to Gongboxia area in the upper reaches of the Yellow River, landslide hazards are characterized by being high-altitude, hidden, and complex, which brings great challenges to the investigation of landslide hazards.

[0003] Traditional methods for landslide hazard investigation, such as manual on-site surveys and the deployment of ground monitoring equipment, are limited by complex terrain. They are not only difficult and costly to carry out, but also difficult to achieve comprehensive detection over large areas and in all weather conditions. They also suffer from many blind spots and low accuracy, and cannot meet the needs for efficient identification of landslide hazards in such areas.

[0004] The development of remote sensing technology has provided new technical pathways for landslide hazard identification. InSAR technology can acquire surface deformation information, while SBAS-InSAR, as an optimized version of InSAR, can reduce the impact of unwrapping errors, atmospheric delay effects, and DEM errors through multi-temporal SAR image differential processing, enabling large-area surface deformation monitoring across all times and weather. Optical remote sensing images can visually present surface morphology, vegetation cover, and texture features, facilitating detailed observation and analysis of surface areas. However, current technologies often use either SBAS-InSAR or optical remote sensing alone for landslide hazard identification. Using SBAS-InSAR alone can easily misclassify surface deformation caused by human activities (such as farmland and village construction) as landslide hazards, while using optical remote sensing alone struggles to capture early landslide deformation information that is highly concealed and lacks obvious morphological features. Both suffer from insufficient accuracy and cannot fully meet the needs for precise identification of landslide hazards in complex terrain areas. Summary of the Invention

[0005] The purpose of this invention is to overcome the shortcomings of insufficient accuracy in landslide hazard identification when using SBAS-InSAR or optical remote sensing alone in existing technologies. It provides a landslide hazard identification method based on SBAS-InSAR and optical remote sensing, combining the advantages of SBAS-InSAR in monitoring surface deformation with the advantages of optical remote sensing in identifying surface features to achieve accurate identification of landslide hazards. To achieve the above objective, this invention adopts the following technical solution: A landslide hazard identification method based on SBAS-InSAR and optical remote sensing includes the following steps: S1. Define the initial set of unstable regions: Use synthetic aperture radar satellites to acquire SAR images within the study area, and use SBAS-InSAR technology to process the SAR images. Then, use the standard deviation method to divide the unstable regions and obtain the initial set of unstable regions in the study area. S2. Drafting a list of landslide hazards: Establish interpretation markers for landslide element characteristics, acquire optical remote sensing images of the initial unstable area set, and remove human activity areas from the images to obtain a set of landslide hazard areas. Then, identify landslide hazards in historical deformation areas and landslide hazards in deforming areas within the set of landslide hazard areas to obtain a list of landslide hazards within the study area.

[0006] Step S1 includes the following steps: S11. Acquire SAR images and auxiliary data within the study area: Collect Sentinel-1A ascending and descending SAR images within the study area, and acquire the corresponding POD precise orbit data and external DEM data. S12 and SBAS-InSAR data processing: Data processing is performed on the rising-orbit SAR image and the falling-orbit SAR image to obtain the rising-orbit surface deformation rate map and the falling-orbit surface deformation rate map within the study area. S13. Comparative analysis of the deformation results of the ascending and descending orbits: The range of surface deformation rate of the ascending orbit and the surface deformation rate of the descending orbit, as well as the number of effective deformation information points, are statistically analyzed. The consistency of the distribution of deformation areas of the ascending and descending orbits is analyzed, and the positive and negative characteristics of the deformation data caused by the difference in the flight azimuth of Sentinel-1A in the same area are clarified. S14. Delineation of Unstable Regions: The standard deviation method is used to delineate unstable regions, resulting in an initial set of unstable regions in the study area.

[0007] In step S12, the SARScape module in ENVI software is used to process the ascending-orbit SAR image and the descending-orbit SAR image to obtain the ascending-orbit surface deformation rate map and the descending-orbit surface deformation rate map within the study area.

[0008] Step S12 includes the following steps: S121. Data import and cropping: Import the Sentinel-1A ascending SAR image and descending SAR image into ENVI software to generate intensity maps and strip index files. Then import POD precise orbit data to correct orbit information. Finally, use the SampleSelection SAR Geometry Data tool to perform vector cropping of the ascending SAR image and descending SAR image according to the study area. S122. Baseline estimation: Perform interferometric relative registration on the cropped ascending-orbit SAR image and descending-orbit SAR image respectively, automatically select the super master image, and register the remaining images to the super master image. Then set the temporal baseline and spatial baseline to generate a spatiotemporal baseline connection map. S123. Interference workflow: The registered interferometric image pairs are processed for coherence generation, flattening, filtering and phase unwrapping. The Goldstein adaptive filtering method is used to remove noise from the interferogram, and the least squares method is used for phase unwrapping to remove phase pairs with low coherence. S124, Track Refining and Re-flattening: Select control points to remove phase residuals caused by track errors, and further remove residual phase shifts and trends in the interferogram; S125, Two inversions: The first inversion estimates the deformation rate and residual topography, and performs a second unwrapping to optimize the quality of the interferogram. The second inversion is based on the deformation rate obtained from the first inversion, and estimates and removes atmospheric phase delay through atmospheric high-pass filtering and atmospheric low-pass filtering. S126. Geocoding: The deformation results after removing atmospheric phase delay are transformed from the SAR coordinate system to the WGS-84 mapping reference system to obtain the surface deformation rate maps of the ascending orbit and descending orbit within the study area.

[0009] In step S122, the time baseline is set to 365 days, and the spatial baseline is 40% of the critical baseline.

[0010] The standard deviation calculation formula used in step S14 is as follows: , In the formula, The deformation of each cell, The mean of the deformation of all pixels. n Total number of pixels; The method for delineating unstable regions using the standard deviation method is as follows: Calculate the mean and standard deviation of the surface deformation rate for both ascending and descending orbits, and then determine the deformation rate within the range specified by the standard deviation. The regions where pixels are located outside the range are identified as unstable regions. The unstable regions defined by the ascending orbit and the unstable regions defined by the descending orbit are merged to obtain the initial set of unstable regions in the study area.

[0011] Step S2 includes the following steps: S21. Establish interpretation markers for landslide elements: take landslide walls, landslide steps, and landslide tongues as core interpretation markers, and clarify that landslide hazards include landslide hazards in historical deformation zones and landslide hazards in deforming zones. Among them, landslide hazards in historical deformation zones have complete landslide walls, landslide steps, and landslide tongues, while landslide hazards in deforming zones have arc-shaped tension cracks. S22. Acquire optical remote sensing image data: Collect GF-2 optical remote sensing images of the initial unstable area set, combine them with Tianditu online land use data, remove human activity areas from the initial unstable area set to obtain a landslide hazard area set, then identify historical deformation area landslide hazards and currently deforming area landslide hazards in the landslide hazard area set to obtain a list of landslide hazards within the study area.

[0012] The landslide wall appears as an arc-shaped or semi-circular steep slope in optical remote sensing images, the landslide steps appear as flat terrain distributed in a stepped manner in optical remote sensing images, and the slope tongue appears as a tongue-shaped accumulation with a protruding leading edge in optical remote sensing images.

[0013] Step S22 includes the following steps: S221. Removal of human activity areas: By observing the surface features of each initial unstable area through optical remote sensing images and combining them with the online land use data of Tianditu, areas located in cultivated land areas, village flatlands and without slope structures are identified as deformation areas caused by human activities and removed from the set of initial unstable areas to obtain a set of landslide hazard areas. S222. Historical Deformation Zone Identification: Optical remote sensing visual interpretation is performed on each landslide hazard area in the landslide hazard area set to identify areas with complete landslide walls, landslide steps, and landslide tongues, which are identified as historical deformation zone landslide hazards, and their relevant information is recorded. S223. Identification of Deformation Zones: For areas without complete landslide features but with deformation, the boundaries of complex terrain areas are redefined according to the direction of slope sliding, based on the topographic features presented by optical remote sensing images, to obtain potential landslide hazards in deformation zones. S224. Draft a list of potential landslide hazards: The list of potential landslide hazards includes landslide hazards in historical deformation zones and landslide hazards in areas currently undergoing deformation.

[0014] It also includes step S3, verification of the accuracy of landslide hazard identification results: spatially overlaying the set of landslide hazard areas with historical landslide records within the study area to verify the accuracy of the landslide hazard identification results.

[0015] The beneficial effects of this invention are as follows: 1. This invention achieves all-day, all-weather monitoring of large-area surface deformation in the study area using SBAS-InSAR technology. It objectively delineates unstable areas using the standard deviation method, avoiding the subjectivity of manual delineation. Compared with traditional manual survey methods, it significantly improves the efficiency and coverage of unstable area detection, and is especially suitable for complex terrain areas such as high mountains and canyons.

[0016] 2. This invention introduces optical remote sensing technology to perform secondary screening of unstable areas delineated by SBAS-InSAR. By establishing clear landslide interpretation markers, it effectively eliminates non-landslide deformation areas caused by human activities. At the same time, it accurately identifies landslide hazards in historical deformation areas and landslide hazards in deforming areas, solving the problem of easy misjudgment when using SBAS-InSAR alone, and significantly improving the accuracy of landslide hazard identification.

[0017] 3. This invention eliminates the influence of satellite flight azimuth differences on deformation data by processing and comparing SAR image data from ascending and descending orbits. Combined with the analysis of terrain features from optical remote sensing images, it achieves refined identification of landslide hazards in complex terrain areas, providing accurate basic data support for subsequent landslide susceptibility prediction, risk assessment, and disaster prevention and mitigation work. Attached Figure Description

[0018] Figure 1 This is a flowchart of the SBAS-InSAR data processing according to Embodiment 2 of the present invention; Figure 2 The above are surface deformation rate diagrams of the lifting rail in Embodiment 2 of the present invention, wherein (a) is the surface deformation rate diagram of the lifting rail and (b) is the surface deformation rate diagram of the lowering rail. Figure 3 This is a delineation diagram of the unstable region of the lifting rail in Embodiment 2 of the present invention, wherein (a) is a delineation diagram of the unstable region of the lifting rail and (b) is a delineation diagram of the unstable region of the lowering rail. Figure 4 This is a morphological feature diagram of a mature landslide according to Embodiment 2 of the present invention; Figure 5 The images show the initial unstable area map and land use map of Embodiment 2 of the present invention, wherein (a) is the initial unstable area map and (b) is the land use map; Figure 6 This is a landslide hazard diagram of the historical deformation zone in Hedong Village, according to Embodiment 2 of the present invention. Figure 7 This is a deformation delineation and redistribution area map of Embodiment 2 of the present invention, wherein (a) is an initial unstable area map spanning the ridge, and (b) is a landslide hazard map of the deforming zone after the boundary has been redefined; Figure 8 This is a spatial overlay diagram of the landslide hazard area set and historical landslide records in Embodiment 2 of the present invention. Detailed Implementation

[0019] The technical solution of the present invention is further described below, but the scope of protection is not limited to what is described.

[0020] Example 1: The present invention discloses a landslide hazard identification method based on SBAS-InSAR and optical remote sensing, wherein SBAS-InSAR refers to small baseline aggregator aperture radar interferometry, and includes the following steps: S1. Define the initial set of unstable regions: Use synthetic aperture radar satellites to acquire SAR images within the study area, and use SBAS-InSAR technology to process the SAR images. Then, use the standard deviation method to divide the unstable regions and obtain the initial set of unstable regions in the study area.

[0021] S11. Acquire SAR images and auxiliary data within the study area: Collect Sentinel-1A ascending and descending SAR images within the study area, and simultaneously acquire the corresponding POD precise orbit data and external DEM data.

[0022] S12, SBAS-InSAR Data Processing: The SARScape module in ENVI software was used to process the ascending-orbit SAR images and descending-orbit SAR images to obtain the ascending-orbit surface deformation rate map and the descending-orbit surface deformation rate map within the study area.

[0023] S121. Data Import and Cropping: Import Sentinel-1A ascending and descending SAR images into ENVI software to generate intensity maps and strip index files. Then, import POD precise orbit data to correct orbit information. Finally, use the SampleSelection SAR Geometry Data tool to perform vector cropping of the ascending and descending SAR images according to the study area.

[0024] S122. Baseline estimation: Perform interferometric relative registration on the cropped ascending-orbit SAR image and descending-orbit SAR image respectively, automatically select the super master image, and register the remaining images to the super master image. Then set the temporal baseline and spatial baseline to generate a spatiotemporal baseline connection map.

[0025] The time baseline is set at 365 days, and the spatial baseline is 40% of the critical baseline.

[0026] S123. Interference Workflow: The registered interferometric image pairs are processed for coherence generation, flattening, filtering, and phase unwrapping. The Goldstein adaptive filtering method is used to remove noise from the interferogram, and the least squares method is used for phase unwrapping to remove phase pairs with low coherence.

[0027] S124. Orbit Refinement and Re-flattening: Select control points to remove phase residuals caused by orbital errors, and further remove residual phase shifts and trends in the interferogram. Specifically, select stable control points (such as GCPs) to remove phase residuals caused by orbital errors.

[0028] S125. Two inversions: The first inversion estimates the deformation rate and residual topography, and performs a second unwrapping to optimize the quality of the interferogram. The second inversion is based on the deformation rate obtained from the first inversion, and estimates and removes atmospheric phase delay through atmospheric high-pass filtering and atmospheric low-pass filtering.

[0029] S126. Geocoding: The deformation results after removing atmospheric phase delay are transformed from the SAR coordinate system to the WGS-84 mapping reference system to obtain the surface deformation rate maps of the ascending orbit and descending orbit within the study area.

[0030] S13. Comparative Analysis of Deformation Results of Ascending and Descending Orbits: The range and number of effective deformation information points of the surface deformation rate of ascending and descending orbits were statistically analyzed. The consistency of the distribution of deformation areas of ascending and descending orbits was analyzed, and the positive and negative characteristics of deformation data caused by the difference in the flight azimuth of Sentinel-1A in the same area were clarified.

[0031] S14. Delineation of Unstable Regions: The standard deviation method is used to delineate unstable regions, resulting in an initial set of unstable regions in the study area.

[0032] The formula for calculating the standard deviation used in the standard deviation method is: , In the formula, The deformation of each cell, The mean of the deformation of all pixels. n Total number of pixels; The method for delineating unstable regions using the standard deviation method is as follows: Calculate the mean and standard deviation of the surface deformation rate for both ascending and descending orbits, and then determine the deformation rate within the range specified by the standard deviation. The regions where pixels are located outside the range are identified as unstable regions. The unstable regions defined by the ascending orbit and the unstable regions defined by the descending orbit are merged to obtain the initial set of unstable regions in the study area.

[0033] S2. Drafting a list of landslide hazards: Establish interpretation markers for landslide element characteristics, acquire optical remote sensing images of the initial unstable area set, and remove human activity areas from the images to obtain a set of landslide hazard areas. Then, identify landslide hazards in historical deformation areas and landslide hazards in deforming areas within the set of landslide hazard areas to obtain a list of landslide hazards within the study area.

[0034] S21. Establish interpretation markers for landslide elements: take landslide walls, landslide steps, and landslide tongues as core interpretation markers, and clarify that landslide hazards include landslide hazards in historical deformation zones and landslide hazards in deformation zones. Among them, landslide hazards in historical deformation zones have complete landslide walls, landslide steps, and landslide tongues, while landslide hazards in deformation zones have a small number of arc-shaped tension cracks.

[0035] The landslide wall appears as an arc-shaped or semi-circular steep slope in optical remote sensing images, the landslide steps appear as flat terrain distributed in a stepped manner in optical remote sensing images, and the slope tongue appears as a tongue-shaped accumulation with a protruding leading edge in optical remote sensing images.

[0036] S22. Acquire optical remote sensing image data: Collect GF-2 optical remote sensing images of the initial unstable area set, combine them with Tianditu online land use data, remove human activity areas from the initial unstable area set to obtain a landslide hazard area set, then identify historical deformation area landslide hazards and currently deforming area landslide hazards in the landslide hazard area set to obtain a list of landslide hazards within the study area.

[0037] S221. Removal of human activity areas: By observing the surface features of each initial unstable area through optical remote sensing images and combining them with the online land use data from Tianditu, areas located in cultivated land areas (such as terraced fields) and flat land near villages without sloping structures are identified as deformation areas caused by human activities and removed from the set of initial unstable areas to obtain a set of landslide hazard areas.

[0038] S222. Identification of Historical Deformation Zones: Optical remote sensing visual interpretation is used to identify landslide hazard areas in the landslide hazard area set. Areas with intact landslide walls, landslide steps, and landslide tongues are identified as historical deformation zone landslide hazards, and their relevant information is recorded. The relevant information for historical deformation zone landslide hazards includes their boundaries, area, and main sliding direction.

[0039] S223. Identification of Deformation Zones: For areas without complete landslide features but with deformation, the boundaries of complex terrain areas are redefined according to the direction of slope sliding, based on the topographic features presented by optical remote sensing images, to obtain potential landslide hazards in deformation zones. S224. Draft a list of potential landslide hazards: The list of potential landslide hazards includes landslide hazards in historical deformation zones and landslide hazards in areas currently undergoing deformation.

[0040] It also includes step S3, verification of the accuracy of landslide hazard identification results: spatially overlaying the set of landslide hazard areas with historical landslide records within the study area to verify the accuracy of the landslide hazard identification results.

[0041] Example 2: like Figures 1 to 8As shown, this embodiment takes the area from Lijiaxia to Gongboxia in the upper reaches of the Yellow River as the study area, and uses the method described in this invention to identify landslide hazards. The specific steps are as follows: S1. Define the initial set of unstable regions: Use synthetic aperture radar satellites to acquire SAR images within the study area, and use SBAS-InSAR technology to process the SAR images. Then, use the standard deviation method to divide the unstable regions and obtain the initial set of unstable regions in the study area.

[0042] S11. Acquire SAR images and auxiliary data within the study area: Collect 33 Sentinel-1A ascending orbit SAR images and 54 descending orbit SAR images within the study area from August 2022 to August 2024, and simultaneously acquire the corresponding POD precise orbit data and 30m resolution external DEM data.

[0043] S12 and SBAS-InSAR Data Processing: The SARScape 5.4 module in ENVI 5.6 software was used to process the ascending and descending SAR images, obtaining ascending and descending surface deformation rate maps for the study area. The data processing workflow is as follows: Figure 1 As shown, the specific steps include: S121. Data Import and Cropping: Import Sentinel-1A ascending-orbit SAR images and descending-orbit SAR images into ENVI software to generate intensity maps and strip index files. Then, import POD precise orbit data to correct orbit information. Finally, use the SampleSelection SAR Geometry Data tool to load the vector boundary of the study area and crop the ascending-orbit SAR images and descending-orbit SAR images to remove redundant data outside the study area.

[0044] S122. Baseline estimation: Perform interferometric relative registration on the cropped ascending SAR images, automatically select the image from May 10, 2023 as the super master image, and register the remaining 32 ascending SAR images to the super master image. Set the time baseline to 365 days and the spatial baseline to 40% of the critical baseline to generate an ascending orbit spatiotemporal baseline connection map. Similarly, the de-orbiting SAR images were processed. The image from June 2, 2023 was selected as the super master image, and the remaining 53 de-orbiting SAR images were registered to this super master image. The time baseline was set to 365 days and the spatial baseline was 40% of the critical baseline to generate a de-orbiting spatiotemporal baseline connection map.

[0045] S123. Interferometric Workflow: Coherence generation is performed on ascending-orbit interferometric image pairs to obtain a coherence coefficient map; the coherence coefficient map is flattened based on external DEM data to generate a differential interferogram; the differential interferogram is filtered using the Goldstein adaptive filtering method to remove speckle noise; phase unwrapping is performed using the least squares method to obtain a phase unwrapped map; phase pairs with low coherence are removed with a coherence coefficient threshold of 0.3. Descending-orbit interferometric image pairs are processed using the same workflow.

[0046] S124. Orbit Refinement and Re-leveling: Select 10 stable control points (such as exposed bedrock areas) in the study area, remove the phase residual caused by orbit error through orbit refinement, and then perform re-leveling to eliminate residual phase shift and trend.

[0047] S125. Two inversions: The first inversion is based on the processed interferometric image pairs to estimate the surface deformation rate and residual topography of the ascending and descending orbits, and at the same time, a second unwrapping is performed to optimize the quality of the interferogram; the second inversion is based on the deformation rate obtained from the first inversion, and atmospheric high-pass filtering and atmospheric low-pass filtering are used to estimate and remove atmospheric phase delay.

[0048] S126. Geocoding: Using the WGS-84 coordinate system and combining external DEM data, the ascent and descent deformation results are transformed from the SAR coordinate system to the geographic coordinate system to obtain an ascent surface deformation rate map (e.g., Figure 2 (as shown in (a)), and the surface deformation rate diagram of the descending orbit (as shown in (a)). Figure 2 (as shown in (b)).

[0049] S13. Comparative analysis of the deformation results of the ascending and descending rails: The deformation rate of the ascending rail ranged from -174 mm / a to 67 mm / a, and 2,362,936 effective deformation information points (30m×30m pixels) were detected, accounting for 73.81% of the total area of ​​the study area; the deformation rate of the descending rail ranged from -149 mm / a to 148 mm / a, and 2,325,509 effective deformation information points were detected, accounting for 72.64% of the total area of ​​the study area; the deformation areas of the ascending and descending rails were distributed in the same way, but the deformation data of the same area were opposite in sign.

[0050] S14. Delineation of Unstable Regions: Calculate the mean deformation rate of the lifting orbit. =-28.5mm / a, standard deviation is =35.2mm / a; the average deformation rate of the lower rail is =-28.5mm / a, standard deviation is =35.2 mm / a; regions in the ascending orbit with deformation rates less than -98.9 mm / a and greater than 41.9 mm / a, and regions in the descending orbit with deformation rates less than -84 mm / a and greater than = 86.4 mm / a, were identified as unstable regions. Ultimately, 66 unstable regions were delineated in the ascending orbit (e.g., Figure 3 As shown in (a), the descent orbit delineates 128 unstable regions (as shown in...). Figure 3 As shown in (b), after merging, 154 initial unstable regions were obtained (as shown in [the diagram]). Figure 5 As shown in (a), this is the initial set of unstable regions in the study area.

[0051] S2. Drafting a list of landslide hazards: Establish interpretation markers for landslide element characteristics, acquire optical remote sensing images of the initial unstable area set, and remove human activity areas from the images to obtain a set of landslide hazard areas. Then, identify landslide hazards in historical deformation areas and landslide hazards in deforming areas within the set of landslide hazard areas to obtain a list of landslide hazards within the study area.

[0052] S21. Establish interpretation markers for landslide feature characteristics, such as Figure 4 As shown: The core interpretation markers are the landslide wall (arc-shaped / semi-circular steep slope), landslide step (stepped flat terrain), and landslide tongue (tongue-shaped deposits at the leading edge); landslide hazards are divided into landslide hazards in the deformation zone (a small number of arc-shaped tension cracks) and landslide hazards in the historical deformation zone (possessing the three elements of a complete landslide and having ancient landslide deposits).

[0053] S22. Acquire optical remote sensing image data: Collect GF-2 optical remote sensing images of the study area, combined with Tianditu online land use data (land use details such as...). Figure 5 As shown in (b), human activity areas are removed from the initial unstable area set to obtain the landslide hazard area set. Then, the landslide hazards in the historical deformation area and the landslide hazards in the deformation area are identified in the landslide hazard area set to obtain the landslide hazard list within the study area.

[0054] S221. Removal of Human Activity Areas: Observation of optical images revealed that the initial unstable areas in Shalianbao Township and Dehenglong Township were located in cultivated land areas, specifically terraced fields with a slope of <15°, indicating deformation caused by cultivated land. In addition, some initial unstable areas in Cuo'en Township were located on flat land near villages, with a deformation rate of approximately -27 mm / a, but without sloping structures, indicating deformation caused by village construction. A total of 32 human activity areas were removed, resulting in a set of landslide hazard areas. S222. Historical Deformation Zone Identification: Visual interpretation of the remaining 122 unstable areas revealed an area in Hedong Village with an arc-shaped steep slope (landslide wall), three levels of stepped terrain (landslide steps), and a tongue-shaped deposit at the leading edge (landslide tongue), covering an area of ​​23,300 m². 2(180m long, 140m wide), identified as a potential landslide hazard in a historical deformation zone (e.g.) Figure 6 (as shown) S223. Identification of Deformation Zones: Targeting initially unstable areas spanning ridges (e.g., Figure 7 As shown in (a), combined with the slope (30°~45°) and slope aspect (sunny slope) characteristics of the optical image, the boundaries were redefined according to the sliding direction of the slopes on both sides of the ridge (sliding southwest on the west side and southeast on the east side), resulting in 3 potential landslide sites in the deformation zone (e.g. Figure 7 (as shown in (b)); a total of 113 deformation zones were identified. S224. Drafting a list of potential landslide hazards: As can be seen from step S222, 122 potential landslide hazards have been identified in the area from Lijiaxia to Gongboxia in the upper reaches of the Yellow River, which constitutes the list of potential landslide hazards.

[0055] S3. Verification of the accuracy of landslide hazard identification results: The accuracy of the landslide hazard identification results was verified by spatially overlaying the 122 landslide hazard areas (i.e., the set of landslide hazard areas) with historical landslide records within the study area. For example... Figure 8 As shown in the figure, most of the historical landslide points fall near the spatial locations of the landslide points identified by InSAR and optical remote sensing. This indicates that these two technologies have a certain degree of scientific validity and practicality in landslide hazard identification. They can capture key information related to landslides, such as topography and surface deformation, and the extracted landslide information has a high degree of consistency with the actual situation.

[0056] The above embodiments demonstrate that the present invention, by combining SBAS-InSAR with optical remote sensing, successfully identified 122 potential landslide sites in the Lijiaxia to Gongboxia area of ​​the upper Yellow River, significantly improving the efficiency and accuracy of landslide hazard identification in complex terrain areas such as high mountains and canyons, and providing reliable data support for disaster prevention and mitigation work.

Claims

1. A method for identifying landslide hazards based on SBAS-InSAR and optical remote sensing, characterized in that: Includes the following steps: S1. Define the initial set of unstable regions: Use synthetic aperture radar satellites to acquire SAR images within the study area, and use SBAS-InSAR technology to process the SAR images. Then, use the standard deviation method to divide the unstable regions and obtain the initial set of unstable regions in the study area. S2. Drafting a list of landslide hazards: Establish interpretation markers for landslide element characteristics, acquire optical remote sensing images of the initial unstable area set, and remove human activity areas from the images to obtain a set of landslide hazard areas. Then, identify landslide hazards in historical deformation areas and landslide hazards in deforming areas within the set of landslide hazard areas to obtain a list of landslide hazards within the study area.

2. The landslide hazard identification method based on SBAS-InSAR and optical remote sensing as described in claim 1, characterized in that: Step S1 includes the following steps: S11. Acquire SAR images and auxiliary data within the study area: Collect Sentinel-1A ascending and descending SAR images within the study area, and acquire the corresponding POD precise orbit data and external DEM data. S12 and SBAS-InSAR data processing: Data processing is performed on the rising-orbit SAR image and the falling-orbit SAR image to obtain the rising-orbit surface deformation rate map and the falling-orbit surface deformation rate map within the study area. S13. Comparative analysis of the deformation results of the ascending and descending orbits: The range of surface deformation rate of the ascending orbit and the surface deformation rate of the descending orbit, as well as the number of effective deformation information points, are statistically analyzed. The consistency of the distribution of deformation areas of the ascending and descending orbits is analyzed, and the positive and negative characteristics of the deformation data caused by the difference in the flight azimuth of Sentinel-1A in the same area are clarified. S14. Delineation of Unstable Regions: The standard deviation method is used to delineate unstable regions, resulting in an initial set of unstable regions in the study area.

3. The landslide hazard identification method based on SBAS-InSAR and optical remote sensing as described in claim 2, characterized in that: In step S12, the SARScape module in ENVI software is used to process the ascending-orbit SAR image and the descending-orbit SAR image to obtain the ascending-orbit surface deformation rate map and the descending-orbit surface deformation rate map within the study area.

4. The landslide hazard identification method based on SBAS-InSAR and optical remote sensing as described in claim 3, characterized in that: Step S12 includes the following steps: S121. Data import and cropping: Import the Sentinel-1A ascending SAR image and descending SAR image into ENVI software to generate intensity maps and strip index files. Then import POD precise orbit data to correct orbit information. Finally, use the SampleSelection SAR Geometry Data tool to perform vector cropping of the ascending SAR image and descending SAR image according to the study area. S122. Baseline estimation: Perform interferometric relative registration on the cropped ascending-orbit SAR image and descending-orbit SAR image respectively, automatically select the super master image, and register the remaining images to the super master image. Then set the temporal baseline and spatial baseline to generate a spatiotemporal baseline connection map. S123. Interference workflow: The registered interferometric image pairs are processed for coherence generation, flattening, filtering and phase unwrapping. The Goldstein adaptive filtering method is used to remove noise from the interferogram, and the least squares method is used for phase unwrapping to remove phase pairs with low coherence. S124, Track Refining and Re-flattening: Select control points to remove phase residuals caused by track errors, and further remove residual phase shifts and trends in the interferogram; S125, Two inversions: The first inversion estimates the deformation rate and residual topography, and performs a second unwrapping to optimize the quality of the interferogram. The second inversion is based on the deformation rate obtained from the first inversion, and estimates and removes atmospheric phase delay through atmospheric high-pass filtering and atmospheric low-pass filtering. S126. Geocoding: The deformation results after removing atmospheric phase delay are transformed from the SAR coordinate system to the WGS-84 mapping reference system to obtain the surface deformation rate maps of the ascending orbit and descending orbit within the study area.

5. The landslide hazard identification method based on SBAS-InSAR and optical remote sensing as described in claim 4, characterized in that: In step S122, the time baseline is set to 365 days, and the spatial baseline is 40% of the critical baseline.

6. The landslide hazard identification method based on SBAS-InSAR and optical remote sensing as described in claim 2, characterized in that: The standard deviation calculation formula used in step S14 is as follows: ， In the formula, The deformation of each cell, The mean of the deformation of all pixels. n Total number of pixels; The method for delineating unstable regions using the standard deviation method is as follows: Calculate the mean and standard deviation of the surface deformation rate for both ascending and descending orbits, and then determine the deformation rate within the range specified by the standard deviation. The regions where pixels are located outside the range are identified as unstable regions. The unstable regions defined by the ascending orbit and the unstable regions defined by the descending orbit are merged to obtain the initial set of unstable regions in the study area.

7. The landslide hazard identification method based on SBAS-InSAR and optical remote sensing as described in claim 1, characterized in that: Step S2 includes the following steps: S21. Establish interpretation markers for landslide elements: take landslide walls, landslide steps, and landslide tongues as core interpretation markers, and clarify that landslide hazards include landslide hazards in historical deformation zones and landslide hazards in deforming zones. Among them, landslide hazards in historical deformation zones have complete landslide walls, landslide steps, and landslide tongues, while landslide hazards in deforming zones have arc-shaped tension cracks. S22. Acquire optical remote sensing image data: Collect GF-2 optical remote sensing images of the initial unstable area set, combine them with Tianditu online land use data, remove human activity areas from the initial unstable area set to obtain a landslide hazard area set, then identify historical deformation area landslide hazards and currently deforming area landslide hazards in the landslide hazard area set to obtain a list of landslide hazards within the study area.

8. The landslide hazard identification method based on SBAS-InSAR and optical remote sensing as described in claim 7, characterized in that: The landslide wall appears as an arc-shaped or semi-circular steep slope in optical remote sensing images, the landslide steps appear as flat terrain distributed in a stepped manner in optical remote sensing images, and the slope tongue appears as a tongue-shaped accumulation with a protruding leading edge in optical remote sensing images.

9. The landslide hazard identification method based on SBAS-InSAR and optical remote sensing as described in claim 7, characterized in that: Step S22 includes the following steps: S221. Removal of human activity areas: By observing the surface features of each initial unstable area through optical remote sensing images and combining them with the online land use data of Tianditu, areas located in cultivated land areas, village flatlands and without slope structures are identified as deformation areas caused by human activities and removed from the set of initial unstable areas to obtain a set of landslide hazard areas. S222. Historical Deformation Zone Identification: Optical remote sensing visual interpretation is performed on each landslide hazard area in the landslide hazard area set to identify areas with complete landslide walls, landslide steps, and landslide tongues, which are identified as historical deformation zone landslide hazards, and their relevant information is recorded. S223. Identification of Deformation Zones: For areas without complete landslide features but with deformation, the boundaries of complex terrain areas are redefined according to the direction of slope sliding, based on the topographic features presented by optical remote sensing images, to obtain potential landslide hazards in deformation zones. S224. Draft a list of potential landslide hazards: The list of potential landslide hazards includes landslide hazards in historical deformation zones and landslide hazards in areas currently undergoing deformation.

10. The landslide hazard identification method based on SBAS-InSAR and optical remote sensing as described in claim 1, characterized in that: It also includes step S3, verification of the accuracy of landslide hazard identification results: spatially overlaying the set of landslide hazard areas with historical landslide records within the study area to verify the accuracy of the landslide hazard identification results.

Images