InSAR coupled ERT frozen soil internal hydrothermal inversion method and system

Abstract

The invention discloses an InSAR (Interferometric Synthetic Aperture Radar) coupled ERT (Enhanced Reverse Transcription) frozen soil internal hydrothermal inversion method and system, and relates to the technical field of frozen soil internal hydrothermal characteristic inversion. Comprising the following steps: separating freeze-thaw deformation data according to InSAR earth surface deformation data; according to the multi-stage ERT data, determining a dynamic freeze-thaw interface depth, and inverting a first water content; according to the dynamic freeze-thaw interface depth, the thickness of an active layer in the frost heaving theoretical model is corrected, and a corrected frost heaving theoretical model is obtained; according to the freeze-thaw deformation data, water content inversion is carried out based on the corrected frost heaving theoretical model, and a second water content is obtained; optimizing the hydrothermal coupling numerical model by using the first water content; and inputting the InSAR earth surface deformation data and the surface layer temperature into the optimized dynamic hydrothermal coupling numerical model, and outputting the temperatures and the third water content of different depths in the frozen soil in the research area in different periods. And finally, high-precision dynamic inversion of the internal hydrothermal characteristics of the regional scale frozen soil is realized.

Description

Technical Field

[0001] The present invention relates to the technical field of inversion of hydrothermal characteristics inside frozen soil, and in particular to an InSAR-coupled ERT hydrothermal inversion method and system for frozen soil. Background Art

[0002] The internal hydrothermal characteristics of frozen soil are the key factors to reveal the process of frozen soil degradation under climate change and predict the evolution of frozen soil disasters. High-precision inversion of the internal hydrothermal characteristics of frozen soil is of great significance. At present, the main inversion methods include field measurement method, model simulation method and satellite remote sensing inversion method: the field measurement method obtains single-point vertical hydrothermal profile data through borehole temperature measurement, time domain reflectometer (TDR) and other equipment. Although the single-point detection accuracy can reach centimeter-level resolution, it is necessary to intensively deploy and consume a lot of resources to obtain regional data; the model simulation method is mainly based on the law of conservation of matter and energy, mainly based on the land process model, simulating the energy cycle process of atmosphere-vegetation-soil. Although it can obtain hydrothermal product data above the kilometer level, the spatial resolution is low and the data update cycle is long, which cannot meet the application requirements; satellite remote sensing inversion methods, such as thermal infrared remote sensing and passive microwave remote sensing, can realize large-scale surface hydrothermal monitoring on a spatial scale, but its detection depth is limited to the surface (<10cm), and the hydrothermal inversion of the active layer of frozen soil (1-3m inside) is completely ineffective. Therefore, the method of obtaining internal hydrothermal data of permafrost at the regional scale has not yet been effectively solved.

[0003] In view of this, this application is hereby filed. Summary of the invention

[0004] The technical problem to be solved by the present invention is that the existing hydrothermal inversion methods for the interior of frozen soil have problems such as low inversion accuracy and poor universality. The purpose of the present invention is to provide a hydrothermal inversion method and system for the interior of frozen soil by InSAR coupling ERT, which obtains water content by monitoring the surface deformation of InSAR and combining the frost heave theory, obtains the vertical distribution function of resistivity and the dynamic value of the depth of the freeze-thaw interface by ERT resistivity tomography, and constructs a hydrothermal coupling numerical model considering time-varying characteristics, and finally realizes high-precision dynamic inversion of the hydrothermal characteristics of the interior of frozen soil at a regional scale, providing key data and solutions for frozen soil engineering safety, frozen soil disaster warning, etc.

[0005] Among them, InSAR (Interferometry Synthetic Aperture Radar) refers to synthetic aperture radar interferometry, and ERT (Electrical Resistance Tomography) refers to electrical resistance tomography technology.

[0006] The present invention is achieved through the following technical solutions: In a first aspect, the present invention provides a method for inversion of internal hydrothermal data of frozen soil by InSAR coupled with ERT, the method comprising: Obtain InSAR surface deformation data, ERT data and surface layer temperature in the same study area; According to the InSAR surface deformation data, freeze-thaw deformation data are separated based on the improved InSAR frozen soil deformation separation model; Based on multiple ERT data, the dynamic freeze-thaw interface depth is determined and the first water content is inverted; According to the depth of the dynamic freeze-thaw interface, the thickness of the active layer in the frost heave theoretical model is corrected to obtain the corrected frost heave theoretical model; According to the freeze-thaw deformation data, the water content is inverted based on the modified frost heave theory model to obtain the second water content; The second water content and the surface temperature are used as initial conditions, and the dynamic freeze-thaw interface depth is used as the boundary condition to input the hydrothermal coupling numerical model; and the hydrothermal coupling numerical model is optimized using the first water content to obtain an optimized dynamic hydrothermal coupling numerical model; The InSAR surface deformation data and surface temperature are input into the optimized dynamic hydrothermal coupling numerical model, and the temperature and third water content at different depths and periods inside the frozen soil in the study area are output.

[0007] Furthermore, InSAR surface deformation data is obtained by processing SAR image data and DEM data based on InSAR time-series interferometry technology; ERT data is obtained by deploying resistivity tomography systems at typical points in the study area to obtain multi-period resistivity profile data ρ(z,t) covering the entire year; The surface temperature is the time series data Ts(t) of the surface temperature obtained through remote sensing product data, meteorological data or ground monitoring stations.

[0008] Furthermore, based on the improved InSAR frozen soil deformation separation model, freeze-thaw deformation data are separated, including: According to the multi-period SAR image data and DEM data, after image registration and terrain phase removal, the preset time baseline and space baseline threshold are selected to form an interference pair to generate a differential interferometric image; According to the differential interference image, selecting high coherence points as a target point set; Using a minimum cost flow algorithm, phase unwrapping is performed on the target point set to obtain an unwrapped phase; According to the unwrapped phase, the traditional frozen soil deformation separation model is improved by combining the internal hydrothermal factor to simulate the seasonal phase, and the improved InSAR frozen soil deformation separation model is obtained. Based on the improved InSAR frozen soil deformation separation model, the seasonal deformation of frozen soil and the deformation of the slope itself are separated; and the seasonal deformation of frozen soil is converted into freeze-thaw deformation data.

[0009] Furthermore, the freeze-thaw deformation data is the surface freeze-thaw deformation d(t) that changes with time in the study area, where d(t)=[d(t1),d(t2),…,d(t N )],t1~t N is a time series, d(t N ) is the tth N Surface freeze-thaw deformation over time.

[0010] Furthermore, based on the multi-period ERT data, the dynamic freeze-thaw interface depth is determined and the first water content is inverted, including: Based on multi-period ERT data, the depth of the dynamic freeze-thaw interface at different times is determined by identifying the extreme values ​​of resistivity gradient; Combined with Archie's law, the first water content θ at different depths at different times is obtained by inversion: ERT (z,t).

[0011] Furthermore, the second moisture content is an average soil moisture content corresponding to the freezing depth or the thawing depth.

[0012] Furthermore, the water content inversion is performed based on the modified frost heave theoretical model to obtain the second water content, including: The inversion formula of water content during freezing period is: ; The inversion formula of water content during the melting period is: ; In the formula, For time t Changes in freezing depth or thawing depth The average soil moisture content, i.e. the second moisture content; It is the residual water content of soil after freezing; is the temperature correction factor; is the soil freezing temperature; is the surface temperature; is the freeze-thaw deformation change, equal to d ( t )-d( t -1), d ( t ) is the time t The surface freeze-thaw deformation under d ( t -1) is time t The time before t-1 freeze-thaw deformation; is the ice resistivity; is the resistivity of water; is the change in soil freezing depth or melting depth, equal to z f ( t )-z f ( t -1), z f ( t ) is the time t Depth of freeze-thaw interface below, z f ( t -1) is time t The time before t The freeze-thaw interface depth of -1 is provided by ERT data.

[0013] Furthermore, the optimization steps of the hydrothermal coupling numerical model are: Based on the prior parameters, generate N Parameter Set ; is the parameter vector to be optimized, i is a positive integer, is the soil effective thermal conductivity parameter, is the volume heat capacity of soil, is the water permeability, is the residual unfrozen water content, is the expansion coefficient; [ ] T is the transpose of a vector; For each parameter set , the depth is obtained by the hydrothermal coupling numerical model z time t Predicted moisture content ; The first moisture content of the inversion is calculated by time t Generates the observation vector , θ ERT ( z m , t ) is the ERT reverse performance time t Lower depth z m The first moisture content at According to the water content prediction value and the observation vector, the objective function is constructed. Calculate the global residual sum of squares, N t is the observation time, N z is the maximum freezing or thawing thickness; θ ERT ( z , t) is the ERT reverse performance time t The first water content at the lower depth z; According to the global residual square sum, the minimum residual is obtained through gradient calculation to update the hydrothermal characteristic parameters and obtain an optimized dynamic hydrothermal coupling numerical model.

[0014] The above technical solution corrects the water content prediction value based on the observation vector and optimizes the water-heat coupling numerical model.

[0015] In a second aspect, the present invention further provides an InSAR coupled ERT frozen soil internal hydrothermal inversion system, the system comprising: The acquisition unit is used to obtain InSAR surface deformation data, ERT data and surface layer temperature in the same study area; A separation unit is used to separate freeze-thaw deformation data according to the InSAR surface deformation data based on an improved InSAR frozen soil deformation separation model; A determination unit is used to determine the depth of the dynamic freeze-thaw interface and invert the first water content based on multiple ERT data; A frost heave model correction unit is used to correct the thickness of the active layer in the frost heave theoretical model according to the depth of the dynamic freeze-thaw interface to obtain a corrected frost heave theoretical model; A water content inversion unit is used to perform water content inversion based on freeze-thaw deformation data and a modified frost heave theory model to obtain a second water content; The hydrothermal model optimization unit is used to input the hydrothermal coupling numerical model with the second water content and the surface temperature as initial conditions and the dynamic freeze-thaw interface depth as boundary conditions; and optimize the hydrothermal coupling numerical model with the first water content to obtain an optimized dynamic hydrothermal coupling numerical model; The hydrothermal inversion unit inside the frozen soil is used to input the InSAR surface deformation data and the surface layer temperature into the optimized dynamic hydrothermal coupling numerical model, and output the temperature and the third water content at different depths in the frozen soil of the study area at different times.

[0016] Further, the determination unit includes: Determine subunits for determining the depth of the dynamic freeze-thaw interface at different times by identifying the extreme values ​​of resistivity gradient based on multi-period ERT data; The inversion subunit is used to combine Archie's law to invert the first water content at different depths at different times.

[0017] Compared with the prior art, the present invention has the following advantages and beneficial effects: 1. The InSAR-coupled ERT method and system for inversion of water and heat inside frozen soil of the present invention obtains water content by monitoring surface deformation with InSAR and combining frost heave theory, obtains the vertical distribution function of resistivity and the dynamic value of freeze-thaw interface depth by ERT resistivity tomography, and constructs a water-heat coupling numerical model considering time-varying characteristics, ultimately achieving high-precision dynamic inversion of water and heat characteristics inside frozen soil at a regional scale, with good universality, and providing key data and solutions for frozen soil engineering safety, frozen soil disaster warning, etc.; 2. The present invention uses high-resolution SAR images such as Sentinel-1, LT-1, and GF-3, and uses time-series InSAR time-series interferometry technology to obtain seasonal surface deformation, and then inverts high-temporal and spatial resolution water content (≤12 days, ≤30m), replacing traditional meteorological station interpolation and traditional remote sensing data; 3. Traditional remote sensing methods can only invert the surface water content (<10cm). This invention combines InSAR with ERT data to build a water-heat coupling model that takes into account time-varying characteristics, and realizes the depth z of the freeze-thaw interface. f (t) dynamic calibration, real-time tracking of the freezing and thawing process of frozen soil, and dynamic inversion of temperature and water content parameters at different depths inside the frozen soil; 4. Soil hydrothermal parameters input by traditional hydrothermal model (parameter soil effective thermal conductivity , soil volume heat capacity The present invention can obtain high-resolution (<0.1m) vertical water content by combining in-situ observation with ERT, optimize the parameters of the traditional hydrothermal coupling model, and further improve the accuracy of the hydrothermal model. 5. The present invention integrates multi-source data such as InSAR, ERT, and in-situ observations, combines time series analysis and hydrothermal coupling numerical models, and can provide technical solutions and effective data for inverting the internal hydrothermal field at the regional level. BRIEF DESCRIPTION OF THE DRAWINGS

[0018] The drawings described herein are used to provide a further understanding of the embodiments of the present invention, constitute a part of this application, and do not constitute a limitation of the embodiments of the present invention. In the drawings: Figure 1 This is a flow chart of the hydrothermal inversion method for frozen soil internal using InSAR coupled with ERT of the present invention; Figure 2 This is a detailed flow chart of the hydrothermal inversion method for frozen soil internal using InSAR coupled with ERT of the present invention; Figure 3 This is a graph showing the monitoring results of the surface deformation rate in the region of the present invention; Figure 4 This is the ERT inversion water content result diagram of the present invention; Figure 5 It is the regional deformation-water content inversion diagram of the present invention; FIG6 (a) is a temperature distribution diagram of the present invention at a depth of 0.3 m; FIG6 (b) is a temperature distribution diagram of the present invention at a depth of 1 m; Figure 6 (c) is a temperature distribution diagram of the present invention at a depth of 2 m; Figure 6 (d) is a temperature distribution diagram of the present invention at a depth of 2.5 m; Figure 7 (a) is a water content distribution diagram of the present invention at a depth of 0.3 m; Figure 7 (b) is a water content distribution diagram of the present invention at a depth of 1 m; Figure 7 (c) is a water content distribution diagram of the present invention at a depth of 2 m; Figure 7 (d) is a water content distribution diagram of the present invention at a depth of 2.5 m; Figure 8 This is a structural block diagram of the InSAR-coupled ERT frozen soil internal hydrothermal inversion system of the present invention. DETAILED DESCRIPTION

[0019] In order to make the objectives, technical solutions and advantages of the present invention more clearly understood, the present invention is further described in detail below in conjunction with embodiments and drawings. The exemplary embodiments of the present invention and their description are only used to explain the present invention and are not intended to limit the present invention.

[0020] The innovation of this invention is that InSAR is coupled with ERT to achieve high-precision dynamic inversion of the internal hydrothermal characteristics of frozen soil at a regional scale, providing key data and solutions for frozen soil engineering safety, frozen soil disaster warning, etc.

[0021] First, by utilizing the high temporal and spatial resolution observation characteristics of InSAR time-series interferometry technology, combined with the frost heave theory and the high vertical resolution of ERT, we can break through the blind spots of existing remote sensing technology for internal water and heat detection of frozen soil, and achieve high-precision inversion of internal water and heat fields in multiple dimensions of time (≤12 days), space (≤30m), and vertical (≤3m). Specifically, the combined use of InSAR and ERT: InSAR can provide surface deformation data with high temporal and spatial resolution, while ERT can provide relatively accurate resistivity data over a long period of time, which is used to determine the depth of the dynamic freeze-thaw interface and the internal soil moisture content. The combination of the two through the frost heave theory can well compensate for the vertical distribution characteristics of soil moisture content inverted by remote sensing, and provide multi-dimensional internal moisture content information of frozen soil.

[0022] Second, by dynamically updating the freeze-thaw interface depth through ERT, the low regional hydrothermal inversion accuracy and poor universality caused by factors such as unclear internal structure, inaccurate freeze-thaw process parameters, and uncertain freeze-thaw interface in the traditional hydrothermal coupling numerical model can be solved, thereby solving the problem of low traditional accuracy. Specifically, the combination of ERT and hydrothermal coupling numerical model: introducing ERT to obtain the time-varying function of the freeze-thaw interface depth, correcting the traditional hydrothermal coupling theory, solving the uncertainty of hydrothermal inversion parameters, and improving the inversion accuracy of internal water content.

[0023] The present invention obtains InSAR surface deformation data by processing multiple SAR images, separates freeze-thaw deformation data (i.e., surface freeze-thaw deformation) d(t), and inverts the water content corresponding to seasonal deformation in combination with frost heave theory as the initial condition of water-heat coupling inside frozen soil; uses time-series resistivity tomography ERT to obtain the resistivity vertical distribution function ρ(z, t), and obtains the dynamic value z of the freeze-thaw interface depth of frozen soil. f (t) is used as the internal hydrothermal inversion boundary condition to construct a hydrothermal coupling numerical model considering the time-varying characteristics of frozen soil. The first water content θ ERT (z, t) is used as the verification data of the parameters of the hydrothermal coupling numerical model. Through the optimized model, the surface temperature and InSAR surface deformation data are input to invert the temperature T(z, t) and the third water content θ(z, t) of the hydrothermal distribution inside the regional permafrost.

[0024] Example 1 like Figure 1 and Figure 2 As shown, Figure 1 This is a flow chart of the hydrothermal inversion method for frozen soil internal using InSAR coupled with ERT of the present invention; Figure 2 The following is a detailed flow chart of the method for hydrothermal inversion of frozen soil by InSAR coupled with ERT according to the present invention. The method for hydrothermal inversion of frozen soil by InSAR coupled with ERT according to the present invention comprises: S1, obtain InSAR surface deformation data, ERT data and surface layer temperature in the same study area; Specifically, InSAR surface deformation data is obtained by processing SAR image data and DEM data based on InSAR time-series interferometry technology; ERT data is obtained by deploying resistivity tomography systems at typical points in the study area to obtain multi-period resistivity profile data ρ(z,t) covering the entire year; The surface temperature is the time series data Ts(t) of the surface temperature obtained through remote sensing product data, meteorological data or ground monitoring stations.

[0025] In this embodiment, Sentinel-1SAR image data and SRTM30 m resolution DEM data covering the study area from 2019 to 2020 are obtained. The time span needs to cover the entire freeze-thaw cycle, and image registration and terrain phase removal are performed to eliminate the errors caused by uneven terrain. The weighted least squares method is used to complete the precise registration of the primary and secondary images. The standard deviation of the offset discreteness of the same-name point offset is used to evaluate the registration accuracy, which is generally less than 0.001 pixel.

[0026] The resistivity tomography system is deployed on the same section in the study area to obtain resistivity profile data covering multiple periods throughout the year (it is necessary to ensure that the measurement data covers the key periods in the freeze-thaw cycle, the initial freezing period, the deepest period of the freezing interface, the initial melting period, and the maximum melting period), namely, ERT data ρ(z,t); Soil temperature and humidity sensors deployed in the study area obtain the surface soil temperature Ts(t).

[0027] S2, based on the InSAR surface deformation data and the improved InSAR frozen soil deformation separation model, freeze-thaw deformation data are separated; Specifically, step S2 includes: S21, based on the multi-period SAR image data and DEM data, after image registration and terrain phase removal, a preset time baseline and space baseline threshold are selected to form an interference pair to generate a differential interferometric image; S22, selecting high coherence points as target point sets according to the differential interference image; S23, using the minimum cost flow algorithm, performing phase unwrapping on the target point set to obtain the unwrapped phase; S24, based on the unwrapped phase, the traditional frozen soil deformation separation model is improved by combining the internal hydrothermal factor to simulate the seasonal phase, and the improved InSAR frozen soil deformation separation model is obtained; S25, based on the improved InSAR frozen soil deformation separation model, separates the seasonal deformation of frozen soil and the deformation of the slope itself; and converts the seasonal deformation of frozen soil into freeze-thaw deformation data. The freeze-thaw deformation data is the surface freeze-thaw deformation d(t) that changes with time in the study area, that is, the long-term surface freeze-thaw deformation dataset of the study area; where d(t)=[d(t1),d(t2),…,d(t N )],t1~t N is a time series, d(t N ) is the tth N Surface freeze-thaw deformation over time.

[0028] In this embodiment, in step S21, appropriate time baseline (≤24 days) and space baseline (≤200m) thresholds are selected to form an interference pair to generate a differential interference pattern; step S22 sets the optimal threshold by analyzing the characteristics of amplitude stability, phase stability and timing correlation, and selects stable high coherence points as the target point set for subsequent phase extraction. Step S23 uses the minimum cost flow algorithm to untangle the target point set and restore the true phase. The formula is: ; Among them, Δ φ is the unwrapped deformation phase, φ top is the phase related to elevation, φ atm is the phase related to atmospheric disturbance, φ seasonal is the seasonal deformation of frozen soil, φ def is its own deformation driven by gravity, φ noise For other noise.

[0029] The frozen soil deformation separation in step S24 includes: ①Seasonal deformation separation and noise removal; Considering that frozen soil experiences seasonal subsidence and uplift, the frozen soil deformation model can be expressed as: ; Where λ is the radar wavelength, Δh is the height difference, φ atm-h and φ atm-tur : The atmospheric component of altitude and turbulence path delay caused by the complex mountain effect, which is related to the change of the propagation path of the radar signal when it passes through the atmosphere, is removed by spatial and temporal filtering; φ top (Δh,t) is the elevation phase associated with Δh at time t; φ atm-h (Δh,t) is the atmospheric phase associated with Δh at time t; Δφ atm-tur (t) is the delayed atmospheric component of the turbulent path at time t; φ noise By removing spatial filtering, φ def (t) is the slope deformation component at a specific time t, A is the phase constant, N is the number of internal environmental factors; βn is the coefficient related to the internal environmental factor (n=1,2,3,...), Bn(t) is the phase within the time interval t in different interference patterns; through measured data, it is determined that the parameters affecting frozen soil deformation are soil temperature and humidity, so N=2. After removing noise, the phase can be expressed as: ; in, β 1 is the correlation coefficient of internal environmental factor 1, B1(t) is the phase of internal environmental factor 1 at time t; β2 is the correlation coefficient of internal environmental factor 2, B2(t) is the phase of internal environmental factor 2 at time t; ②Phase conversion and solution; The unknown parameter matrix Substitute into the regression equation and solve it using the following matrix equations: ; ; Where B is a matrix that defines the small baseline interferogram combination. The singular value decomposition (SVD) method can be used to obtain a unique solution. The inversion accuracy is evaluated by the phase residual to obtain φ seasonal Seasonal deformation of frozen soil and φ def is its own deformation driven by gravity, and finally φ seasonal The phase information is converted into freeze-thaw deformation, and the surface freeze-thaw deformation d(t) that changes with time in the target area can be obtained. Figure 3 Shown is the monitoring results of regional surface deformation rate from February to April 2020. The unit of deformation rate is mm / d (millimeter / day), and “-” represents melt settlement.

[0030] S3, based on multi-period ERT data, determine the depth of the dynamic freeze-thaw interface and invert the first water content; Specifically, based on multiple ERT data, the dynamic freeze-thaw interface depth is determined and the first water content is inverted, including: Based on multi-period ERT data, the depth of the dynamic freeze-thaw interface at different times is determined by identifying the extreme values ​​of resistivity gradient; Combined with Archie's law, the first water content θ at different depths at different times is obtained by inversion: ERT (z,t).

[0031] In this embodiment, step S3 is specifically as follows: (1) Dynamic calibration of dynamic freeze-thaw interface depth: For each period of ERT data ρ(z,t k ), the depth z is calculated according to the following, and the time t k Downward resistivity gradient , z is the depth, and Δz is the depth step.

[0032] ; Among them, ρ(z,t k ) is time t k The resistivity at depth z, ρ(z+Δz,t k ) is time t k The resistivity at the depth z+Δz; ρ(z-Δz,t k ) is time tk The resistivity at the lower depth z-Δz; time t k Depth of lower freeze-thaw interface is the maximum point of resistivity gradient, and the dynamic value of freeze-thaw interface depth changing with time is obtained according to the following formula.

[0033] ; when hour( T f is the soil freezing temperature), then the corrected freeze-thaw interface depth is: ; in, Z f is the depth of freeze-thaw interface; is time t k Lower depth z f The temperature at is the temperature gradient; T is temperature; z is depth; is the change of temperature with depth; is the corrected freeze-thaw interface depth.

[0034] (2) Water content inversion based on Archie model Based on Archie's law and combined with the dynamic freeze-thaw interface depth, the internal resistivity-water content relationship of frozen soil is constructed: In unfrozen areas, the soil is mainly liquid water, has good conductivity, and the resistivity-water content relationship is: ; in, ϕ(z) is the soil porosity at depth z, z is the depth, ρ w is the resistivity of water, To reflect the empirical coefficient of soil particle characteristics, m is the porosity index, s is the soil saturation, is the resistivity at depth z at time t.

[0035] In the frozen area, the relationship between resistivity and unfrozen water content is: ; in, ρ i is the ice resistivity, is the residual water content of soil after freezing, η is the influence coefficient of unfrozen water, is the resistivity at depth z at time t.

[0036] (3) Parameter optimization Layered soil temperature and humidity sensors are deployed in the typical monitoring area profile to obtain the measured data of soil temperature and humidity. sensor (z,t) and θ sensor (z, t), and determine the resistivity-water content inversion parameters by the least square method, etc. , m, s, and generate resistivity inversion water content distribution map of typical areas.

[0037] like Figure 4 The figure shows the ERT inversion water content result of a certain period during the thawing of frozen soil. In this embodiment, , m=1, ϕ=0.5~0.6, s=3. The current freeze-thaw interface position z is obtained from the ERT detection results. f like Figure 4 Indicated by the red dotted line.

[0038] S4, correcting the thickness of the active layer in the frost heave theoretical model according to the depth of the dynamic freeze-thaw interface to obtain a corrected frost heave theoretical model; In this embodiment, considering that the InSAR seasonal deformation is affected by the internal structure of frozen soil and the depth of the freeze-thaw interface, the dynamic freeze-thaw interface depth observed by ERT is introduced to correct the active layer thickness in the frost heave theoretical model.

[0039] S5, according to the freeze-thaw deformation data, the water content is inverted based on the modified frost heave theory model to obtain the second water content; Specifically, the second moisture content is an average soil moisture content corresponding to the freezing depth or the thawing depth.

[0040] In this example, the surface freeze-thaw deformation d(t) is mainly caused by the freezing expansion or melting settlement of water in the soil, and the deformation is related to the freezing or melting thickness and water content. Considering the inconsistency of heat conduction during the freezing and melting processes of frozen soil, the surface temperature conducts unidirectionally from top to bottom during the freezing period and conducts bidirectionally from top and bottom to the middle during the melting period. The freeze-thaw interface depth z is introduced by ERT in-situ observation. f (t), it can be concluded that during the freezing / thawing period, the period t Thickness change Δz under freezing or thawing f Average soil moisture content (t) , the relationship is expressed as: The inversion formula of water content during freezing period is: ; The inversion formula of water content during the melting period is: ; In the formula, For time t Changes in freezing depth or thawing depth The average soil moisture content, i.e. the second moisture content; It is the residual water content of soil after freezing; is the temperature correction factor; is the soil freezing temperature; is the surface temperature; is the freeze-thaw deformation change, equal to d ( t )-d( t -1), d ( t ) is the time t The surface freeze-thaw deformation under d ( t -1) is time t The time before t -1 freeze-thaw deformation; is the ice resistivity; is the resistivity of water; is the change in soil freezing depth or melting depth, equal to z f ( t )-z f ( t -1), z f ( t ) is the time t Depth of freeze-thaw interface below, z f ( t -1) is time t The time before t -1 freeze-thaw interface depth, provided by ERT data; is the water-ice density ratio, which takes a value of 1.09.

[0041] Based on the regional InSAR surface deformation data obtained above, the regional deformation variable - water content is obtained. Figure 5 As shown, the water content here is the second water content.

[0042] S6, using the second water content and the surface temperature as initial conditions and the dynamic freeze-thaw interface depth as boundary conditions, inputting the hydrothermal coupling numerical model; and optimizing the hydrothermal coupling numerical model using the first water content to obtain an optimized dynamic hydrothermal coupling numerical model; the initial water content is the water content observed in the InSAR surface deformation data; Specifically, the temperature gradient change of frozen soil follows heat conduction, and the water migration follows Darcy's law. The dynamic freeze-thaw interface obtained by ERT is combined to build a hydrothermal coupling numerical model considering time-varying characteristics. The initial water content and surface temperature obtained by InSAR deformation observation are used as initial conditions, and the depth of the dynamic freeze-thaw interface observed by ERT is used as the boundary condition. The hydrothermal coupling numerical model is input and the first water content θ inverted by ERT is used. ERT(z, t) to optimize the model parameters (thermal conductivity k, water diffusion coefficient D, expansion coefficient α, etc.), and iteratively calculate to obtain the optimal hydrothermal numerical solution. Specifically: 1) Establish a hydrothermal coupling numerical model considering time-varying characteristics The water content in frozen soil is not statically distributed. During the freezing / thawing process, the water will undergo a phase change (liquid water ↔ ice), accompanied by latent heat release / absorption. In addition, during the freezing / thawing process, the depth of the freeze-thaw interface z f It is not fixed. Based on the above, an improved hydrothermal coupling numerical model is established as follows: Heat conduction equation: ; Moisture migration equation: ; The latent heat flux of phase change is: ; In the formula, is the effective thermal conductivity of soil, is the volume heat capacity of soil, is the latent heat of water-ice phase change, θ i is the ice volume content, K w is the water infiltration rate, which is related to soil type and temperature; and are the thermal conductivity of the frozen zone and the non-frozen zone, is the non-freezing zone temperature, is the water source item, including factors such as precipitation and evaporation; To solve the partial differential, Solve for the gradient; is the temperature at depth z at time t; is the volume content of ice at depth z at time t; is the water content at depth z at time t; z f (t) is the depth of the freeze-thaw interface at time t.

[0043] 2) Determine the boundary conditions of the equation The upper boundary is the surface layer (z=0), satisfying ; is the surface temperature at time t; The surface temperature obtained by monitoring sensors / weather stations, etc. time t Lower melting depth Δ z f ( t ) corresponds to the average water content , obtained by step S5; Depth of freeze-thaw interface zf Department: ; is the depth position z of the freeze-thaw interface at time t f The temperature at The soil freezing temperature.

[0044] 3) Optimization of hydrothermal coupling process parameters Specifically, the optimization steps of the hydrothermal coupling numerical model are: Based on the prior parameters (soil effective thermal conductivity , soil volume heat capacity , water permeability , residual unfrozen water content , expansion coefficient ),generate N Parameter Set ; is the parameter vector to be optimized, i is a positive integer, [ ] T is the transpose of a vector; For each parameter set , the depth is obtained by the hydrothermal coupling numerical model z time t Predicted moisture content ; The first moisture content of the inversion is calculated by time t Generates the observation vector , θ ERT ( z m , t ) is the ERT reverse performance time t Lower depth z m The first moisture content at According to the water content prediction value and the observation vector, the objective function is constructed. Calculate the global residual sum of squares, N t is the observation time, N z is the maximum freezing or thawing thickness; θ ERT ( z , t ) is the ERT reverse performance time t The first water content at the lower depth z; According to the global residual sum of squares, the minimum residual is obtained through gradient calculation to update the hydrothermal characteristic parameters and obtain the optimized dynamic hydrothermal coupling numerical model.

[0045] The above technical solution corrects the water content prediction value based on the observation vector and optimizes the water-heat coupling numerical model.

[0046] S7, input the InSAR surface deformation data and surface layer temperature into the optimized dynamic hydrothermal coupling numerical model, and output the temperature T(z, t) and the third water content θ(z, t) at different depths and periods inside the frozen soil in the study area.

[0047] In this embodiment, the optimized dynamic water-heat coupling numerical model is used, based on the premise that the internal structure of the regional frozen soil is consistent, and the InSAR surface deformation data and the surface layer temperature are input to obtain the temperature T(z, t) and the third water content θ(z, t) at different depths in the regional frozen soil at different periods. As shown in Figures 6 (a)-6 (d) and 7 (a)-7 (d), Figures 6 (a)-6 (d) and 7 (a)-7 (d) are the regional internal water-heat content maps obtained according to the above steps, reflecting that the method can track the freezing and melting process of frozen soil in real time and dynamically invert the temperature and water content parameters at different depths of the frozen soil in the region, providing an effective solution for inverting the internal water-heat field at the regional level.

[0048] Example 2 like Figure 8 As shown, the difference between this embodiment and embodiment 1 is that this embodiment provides an InSAR coupled ERT frozen soil internal hydrothermal inversion system, which corresponds one-to-one to the InSAR coupled ERT frozen soil internal hydrothermal inversion method of embodiment 1; the system includes: The acquisition unit is used to obtain InSAR surface deformation data, ERT data and surface layer temperature in the same study area; A separation unit is used to separate freeze-thaw deformation data according to the InSAR surface deformation data based on an improved InSAR frozen soil deformation separation model; A determination unit is used to determine the depth of the dynamic freeze-thaw interface and invert the first water content based on multiple ERT data; A frost heave model correction unit is used to correct the thickness of the active layer in the frost heave theoretical model according to the depth of the dynamic freeze-thaw interface to obtain a corrected frost heave theoretical model; A water content inversion unit is used to perform water content inversion based on freeze-thaw deformation data and a modified frost heave theory model to obtain a second water content; The hydrothermal model optimization unit is used to input the hydrothermal coupling numerical model with the second water content and the surface temperature as initial conditions and the dynamic freeze-thaw interface depth as boundary conditions; and optimize the hydrothermal coupling numerical model with the first water content to obtain an optimized dynamic hydrothermal coupling numerical model; The hydrothermal inversion unit inside the frozen soil is used to input the InSAR surface deformation data and the surface layer temperature into the optimized dynamic hydrothermal coupling numerical model, and output the temperature T(z, t) and the third water content θ(z, t) at different depths and periods inside the frozen soil in the study area.

[0049] In this embodiment, the determining unit includes: Determine subunits for determining the depth of the dynamic freeze-thaw interface at different times by identifying the extreme values ​​of resistivity gradient based on multi-period ERT data; The inversion subunit is used to combine Archie's law to invert the first water content θ at different times and depths. ERT (z,t).

[0050] The execution process of each unit can be performed according to the process steps of the InSAR coupled ERT frozen soil internal hydrothermal inversion method in Example 1, and will not be repeated in this embodiment.

[0051] Those skilled in the art will appreciate that the embodiments of the present application may be provided as methods, systems, or computer program products. Therefore, the present application may adopt the form of a complete hardware embodiment, a complete software embodiment, or an embodiment combining software and hardware. Moreover, the present application may adopt 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 codes.

[0052] The present application is described with reference to the flowcharts and / or block diagrams of the methods, devices (systems), and computer program products according to the embodiments of the present application. It should be understood that each process and / or box in the flowchart and / or block diagram, as well as the combination of the processes and / or boxes in the flowchart and / or block diagram, can be implemented by computer program instructions. These computer program instructions can be provided to a processor of a general-purpose computer, a special-purpose computer, an embedded processor, or other programmable data processing device to generate a machine, so that the instructions executed by the processor of the computer or other programmable data processing device generate instructions for implementing the processes in the flowchart and / or block diagram. Figure 1 A process or multiple processes and / or boxes Figure 1 A device that provides the functions specified in a block or multiple blocks.

[0053] These computer program instructions may also be stored in a computer-readable memory capable of directing a computer or other programmable data processing device to operate in a specific manner, so that the instructions stored in the computer-readable memory produce an article of manufacture comprising an instruction device, which implements the process Figure 1 A process or multiple processes and / or boxes Figure 1 A function specified in one or more boxes.

[0054] These computer program instructions can also be loaded onto a computer or other programmable data processing device so that a series of operating steps are executed on the computer or other programmable device to produce a computer-implemented process, thereby providing instructions for implementing the process. Figure 1 A process or multiple processes and / or boxes Figure 1 The steps for the functions specified in one or more boxes.

[0055] The specific implementation methods described above further illustrate the objectives, technical solutions and beneficial effects of the present invention in detail. It should be understood that the above description is only a specific implementation method of the present invention and is not intended to limit the scope of protection of the present invention. Any modifications, equivalent substitutions, improvements, etc. made within the spirit and principles of the present invention should be included in the scope of protection of the present invention.

Claims

1. The method of hydrothermal inversion of frozen soil internal by InSAR coupled with ERT is characterized by: The method includes: Obtain InSAR surface deformation data, ERT data and surface layer temperature in the same study area; According to the InSAR surface deformation data, freeze-thaw deformation data are separated based on the improved InSAR frozen soil deformation separation model; Based on multiple ERT data, the dynamic freeze-thaw interface depth is determined and the first water content is inverted; According to the depth of the dynamic freeze-thaw interface, the thickness of the active layer in the frost heave theoretical model is corrected to obtain the corrected frost heave theoretical model; According to the freeze-thaw deformation data, the water content is inverted based on the modified frost heave theory model to obtain the second water content; The second water content and the surface temperature are used as initial conditions, and the dynamic freeze-thaw interface depth is used as the boundary condition to input the hydrothermal coupling numerical model; and the hydrothermal coupling numerical model is optimized using the first water content to obtain an optimized dynamic hydrothermal coupling numerical model; The InSAR surface deformation data and surface temperature are input into the optimized dynamic hydrothermal coupling numerical model, and the temperature and third water content at different depths and periods inside the frozen soil in the study area are output.

2. The method for hydrothermal inversion of frozen soil interior by InSAR coupled with ERT according to claim 1 is characterized in that: The InSAR surface deformation data is surface deformation data obtained by processing SAR image data and DEM data based on InSAR time series interferometry technology; The ERT data is obtained by deploying a resistivity tomography system in the study area to obtain multi-period resistivity profile data covering the entire year; The surface layer temperature is time series data of the surface temperature obtained through remote sensing product data, meteorological data or ground monitoring stations.

3. The method for hydrothermal inversion of frozen soil interior by InSAR coupled with ERT according to claim 1 is characterized in that: Based on the improved InSAR frozen soil deformation separation model, freeze-thaw deformation data are separated, including: According to the multi-period SAR image data and DEM data, after image registration and terrain phase removal, the preset time baseline and space baseline threshold are selected to form an interference pair to generate a differential interferometric image; According to the differential interference image, selecting high coherence points as a target point set; Using a minimum cost flow algorithm, phase unwrapping is performed on the target point set to obtain an unwrapped phase; According to the unwrapped phase, the traditional frozen soil deformation separation model is improved by combining the internal hydrothermal factor to simulate the seasonal phase, and the improved InSAR frozen soil deformation separation model is obtained. Based on the improved InSAR frozen soil deformation separation model, the seasonal deformation of frozen soil and the deformation of the slope itself are separated; and the seasonal deformation of frozen soil is converted into freeze-thaw deformation data.

4. The method for hydrothermal inversion of frozen soil interior by InSAR coupled with ERT according to claim 3 is characterized in that: The freeze-thaw deformation data is the surface freeze-thaw deformation d(t) that changes with time in the study area, where d(t)=[d(t1),d(t2),…,d(t N )],t1~t N is a time series, d(t N ) is the tth N Surface freeze-thaw deformation over time.

5. The method for hydrothermal inversion of frozen soil interior by InSAR coupled with ERT according to claim 1 is characterized in that: Based on multiple ERT data, the dynamic freeze-thaw interface depth is determined and the first water content is inverted, including: Based on multi-period ERT data, the depth of the dynamic freeze-thaw interface at different times is determined by identifying the extreme values ​​of resistivity gradient; Combined with Archie's law, the first water content at different depths at different times is inverted.

6. The method for hydrothermal inversion of frozen soil interior by InSAR coupled with ERT according to claim 1 is characterized in that: The second moisture content is the average soil moisture content corresponding to the freezing depth or the thawing depth.

7. The method for hydrothermal inversion of frozen soil interior by InSAR coupled with ERT according to claim 6 is characterized in that: Based on the modified frost heave theory model, the water content is inverted to obtain the second water content, including: The inversion formula of water content during freezing period is: ； The inversion formula of water content during the melting period is: ； In the formula, For time t Changes in freezing depth or thawing depth The average soil moisture content, i.e. the second moisture content; It is the residual water content of soil after freezing; is the temperature correction factor; is the soil freezing temperature; is the surface temperature; is the freeze-thaw deformation change, equal to d ( t )-d( t -1), d ( t ) is the time t The surface freeze-thaw deformation under d ( t -1) is time t The time before t -1 freeze-thaw deformation; is the ice resistivity; is the resistivity of water; is the change in soil freezing depth or melting depth, equal to z f ( t )-z f ( t -1), z f ( t ) is the time t Depth of freeze-thaw interface below, z f ( t -1) is time t The time before t The freeze-thaw interface depth of -1 is provided by ERT data.

8. The method for hydrothermal inversion of frozen soil interior by InSAR coupled with ERT according to claim 1 is characterized in that: The optimization steps of the hydrothermal coupling numerical model are: Based on the prior parameters, generate N Parameter Set ; is the parameter vector to be optimized, i is a positive integer, is the soil effective thermal conductivity parameter, is the volume heat capacity of soil, is the water permeability, is the residual unfrozen water content, is the expansion coefficient; [ ] T is the transpose of a vector; For each parameter set , the depth is obtained by the hydrothermal coupling numerical model z time t Predicted moisture content ; The first moisture content of the inversion is calculated by time t Generates the observation vector , θ ERT ( z m , t ) is the ERT reverse performance time t Lower depth z m The first moisture content at According to the water content prediction value and the observation vector, the objective function is constructed. Calculate the global residual sum of squares, N t is the observation time, N z is the maximum freezing or thawing thickness; θ ERT ( z , t ) is the ERT reverse performance time t The first water content at the lower depth z; According to the global residual square sum, the minimum residual is obtained through gradient calculation to update the hydrothermal characteristic parameters and obtain an optimized dynamic hydrothermal coupling numerical model.

9. The InSAR coupled ERT frozen soil internal water and heat inversion system is characterized by: The system includes: The acquisition unit is used to obtain InSAR surface deformation data, ERT data and surface layer temperature in the same study area; A separation unit is used to separate freeze-thaw deformation data based on the InSAR surface deformation data and the improved InSAR frozen soil deformation separation model; A determination unit is used to determine the depth of the dynamic freeze-thaw interface and invert the first water content based on multiple ERT data; A frost heave model correction unit is used to correct the thickness of the active layer in the frost heave theoretical model according to the depth of the dynamic freeze-thaw interface to obtain a corrected frost heave theoretical model; A water content inversion unit is used to perform water content inversion based on freeze-thaw deformation data and a modified frost heave theory model to obtain a second water content; The hydrothermal model optimization unit is used to input the hydrothermal coupling numerical model with the second water content and the surface temperature as initial conditions and the dynamic freeze-thaw interface depth as boundary conditions; and optimize the hydrothermal coupling numerical model with the first water content to obtain an optimized dynamic hydrothermal coupling numerical model; The hydrothermal inversion unit inside the frozen soil is used to input the InSAR surface deformation data and the surface layer temperature into the optimized dynamic hydrothermal coupling numerical model, and output the temperature and the third water content at different depths in the frozen soil of the study area at different times.

10. The InSAR-coupled ERT frozen soil internal hydrothermal inversion system according to claim 9, characterized in that: The determining unit comprises: Determine subunits for determining the depth of the dynamic freeze-thaw interface at different times by identifying the extreme values ​​of resistivity gradient based on multi-period ERT data; The inversion subunit is used to combine Archie's law to invert the first water content at different depths at different times.

Images