CN109212522B

CN109212522B - High-precision DEM inversion method and device based on double-base satellite-borne SAR

Info

Publication number: CN109212522B
Application number: CN201810525764.0A
Authority: CN
Inventors: 秦小芳; 张衡; 王宇; 邓云凯; 张华春; 刘华有
Original assignee: Institute of Electronics of CAS
Current assignee: Institute of Electronics of CAS
Priority date: 2018-05-28
Filing date: 2018-05-28
Publication date: 2020-01-21
Anticipated expiration: 2038-05-28
Also published as: CN109212522A

Abstract

The invention discloses a high-precision DEM inversion method, a device, a storage medium and an information processing device based on a double-base satellite-borne SAR, which are used for acquiring a first main image and a first auxiliary image which are obtained by respectively imaging an imaging target by a first satellite-borne Synthetic Aperture Radar (SAR) and a second satellite-borne SAR in a first track, and acquiring a second main image and a second auxiliary image which are obtained by respectively imaging the imaging target by the first satellite-borne SAR and the second satellite-borne SAR in a second track; generating a first interference pattern corresponding to the first track and a second interference pattern corresponding to the second track; then, carrying out filtering processing, and determining a coherence coefficient corresponding to each pixel on the two interferograms; respectively converting the filtered first interferogram and the filtered second interferogram into a first single-baseline Digital Elevation Model (DEM) and a second single-baseline DEM; and fusing the first single-baseline DEM and the second single-baseline DEM according to the corresponding coherence coefficient of each pixel to obtain a fused DEM of the imaging target.

Description

High-precision DEM inversion method and device based on double-base satellite-borne SAR

Technical Field

The invention relates to a Synthetic Aperture Radar (SAR) technology, in particular to a high-precision DEM inversion method and device based on a double-base spaceborne SAR.

Background

A Digital Elevation Model (DEM) is a three-dimensional Digital Model that describes the shape of the earth's surface and consists of a series of datasets containing geographical plane coordinates and elevations; DEM has important application value in scientific research, economic construction and military fields. In specific application occasions such as seismic deformation extraction, terrain monitoring and the like, the high-resolution high-precision DEM is particularly important; however, DEM extraction is very complex. Interferometric Synthetic Aperture Radar (InSAR) is one of the mature methods for accurately acquiring digital maps due to its all-time and all-weather working characteristics and as an active sensor. However, in a complex terrain area with large surface relief, such as a mountainous area or an urban area with high building height and dense distribution, in the process of acquiring an echo image pair by radar side-view imaging, an area with a large number of data holes, such as overlapping masks and shadows, is easily formed, and the quality of the generated DEM is seriously influenced.

A first satellite-borne double-satellite distributed system in the world is formed by a Terras SAR-X satellite transmitted in 2007 and a TanDEM-X satellite transmitted in 2010, namely TSX/TDX (Terras SAR-X/TanDEM-X); the primary purpose of the TSX/TDX is to obtain a global high precision DEM. Compared with the traditional repeated orbit satellite, the TSX/TDX has the main characteristic of acquiring the interference pair of the '0' time base line and the high-precision orbit, and can well overcome the atmospheric delay and the orbit error. In the generated DEM precision analysis, most of the generated DEMs are aimed at mountainous areas with steep terrain or residential areas with sparsely distributed buildings; quantitative experiments and researches are lacked in the dense urban area of the high building stand, and the generated DEM is easy to have abnormal imaging areas such as perspective shrinkage, overlapping and shadow areas. Meanwhile, the existing data fusion method is mostly based on data in a TSX/TDX strip mode, the spatial resolution is 12 meters, and the interpretation of the topographic information of urban areas is far insufficient.

Therefore, how to improve the imaging quality of the DEM in the dense area of the high-fall object is a problem to be solved urgently.

Disclosure of Invention

In view of this, embodiments of the present invention are expected to provide a high-precision DEM inversion method and apparatus based on double-base satellite-borne SAR, which can improve the resolution and precision of a DEM in a high-fall object dense area, and improve imaging quality.

In order to achieve the purpose, the technical scheme of the invention is realized as follows:

the embodiment of the invention provides a high-precision DEM inversion method based on a double-base spaceborne SAR, which comprises the following steps:

acquiring a first main image and a first auxiliary image which are obtained by respectively imaging an imaging target by a first satellite-borne SAR and a second satellite-borne SAR in a first track, and acquiring a second main image and a second auxiliary image which are obtained by respectively imaging the imaging target by the first satellite-borne SAR and the second satellite-borne SAR in a second track;

processing the first main image and the first auxiliary image by adopting a preset image processing rule to generate a first interference image corresponding to the first track, and processing the second main image and the second auxiliary image to generate a second interference image corresponding to the second track;

according to a preset filtering rule, respectively carrying out filtering processing on the first interference image and the second interference image, and determining a coherence coefficient corresponding to each pixel on the first interference image and the second interference image;

according to a preset conversion rule, respectively converting the filtered first interferogram and the filtered second interferogram into a first single-baseline digital elevation model DEM and a second single-baseline DEM;

and fusing the first single-baseline DEM and the second single-baseline DEM according to the corresponding coherence coefficient of each pixel on the first interference image and the second interference image by adopting a preset fusion rule to obtain a fused DEM of the imaging target.

In the foregoing solution, the processing the first main image and the first auxiliary image according to the preset image processing rule to generate the first interference pattern corresponding to the first track, and processing the second main image and the second auxiliary image to generate the second interference pattern corresponding to the second track includes:

carrying out image registration on the first main image and the first auxiliary image, and carrying out conjugate multiplication on the first main image and the first auxiliary image after registration processing to obtain a first interference image;

and carrying out image registration on the second main image and the second auxiliary image, and carrying out conjugate multiplication on the second main image and the second auxiliary image after registration processing to obtain a second interference image.

In the foregoing solution, the performing filtering processing on the first interferogram and the second interferogram respectively according to a preset filtering rule, and determining a coherence coefficient corresponding to each pixel on the first interferogram and the second interferogram includes:

respectively carrying out iterative estimation on each pixel in the first interference image and the second interference image by adopting a non-local filtering method;

according to the similarity weight of each pixel of the first interference image in iterative estimation, a Goodman model of the first main image and the first auxiliary image after registration processing is adopted to obtain a corresponding coherence coefficient of each pixel of the first interference image;

and obtaining a corresponding coherence coefficient of each pixel of the second interference image by adopting a Goodman model of the second main image and the second auxiliary image after registration processing according to the similarity weight of each pixel of the second interference image in iterative estimation.

In the above scheme, the first filtered interferogram and the second filtered interferogram are converted into a first single-baseline DEM and a second single-baseline DEM respectively according to a preset conversion rule; the method comprises the following steps:

differentiating the preset DEM simulation interference phase and the interference phase of the first interference pattern to obtain a first differential phase;

masking the first differential phase of a region with the coherence coefficient lower than a first preset coherence threshold value in the first interferogram, and unwrapping the masked first differential phase;

in the unwrapped first differential phase, interpolating the mask region by adopting a preset interpolation rule, and adding the interpolated unwrapped first differential phase to a preset DEM simulation interference phase to obtain an unwrapped phase of a first interference pattern;

determining the first single-baseline DEM corresponding to the first track according to the unwrapping phase of the first interferogram by adopting an interferometric synthetic aperture radar InSAR height measurement model;

differentiating the preset DEM simulation interference phase and the interference phase of the second interference pattern to obtain a second differential phase;

masking the second differential phase of the region with the coherence coefficient lower than a first preset coherence threshold value in the second interference pattern, and unwrapping the masked second differential phase;

in the unwrapped second differential phase, interpolating the mask region by adopting a preset interpolation rule, and adding the interpolated unwrapped second differential phase to a preset DEM simulation interference phase to obtain an unwrapped phase of a second interference pattern;

and determining the second single-baseline DEM corresponding to the second track according to the unwrapping phase of the second interferogram by adopting the InSAR height measurement model.

In the above scheme, the first single-baseline DEM and the second single-baseline DEM are fused by adopting a preset fusion rule according to the coherence coefficient corresponding to each pixel on the first interferogram and the second interferogram to obtain a fused DEM of the imaging target; the method comprises the following steps:

respectively masking an imaging abnormal region and a region with a corresponding coherence coefficient lower than a second preset coherence threshold value in the first single-baseline DEM and the second single-baseline DEM;

carrying out weighted fusion on the first single-baseline DEM non-mask area and the second single-baseline DEM non-mask area, supplementing the first single-baseline DEM normal imaging area with the second single-baseline DEM mask area at the corresponding position, and supplementing the second single-baseline DEM normal imaging area with the first single-baseline DEM mask area at the corresponding position to obtain a fused DEM of an imaging target;

determining the proportion of the corresponding coherence coefficient of the first single-baseline DEM non-mask area to the sum of the corresponding coherence coefficients of the first single-baseline DEM non-mask area and the second single-baseline DEM non-mask area as a weight for performing weighted fusion on the first single-baseline DEM non-mask area;

and determining the weight of the second single-baseline DEM non-mask area for carrying out weighted fusion according to the ratio of the corresponding coherence coefficient of the second single-baseline DEM non-mask area to the sum of the corresponding coherence coefficients of the first single-baseline DEM non-mask area and the second single-baseline DEM non-mask area.

In the above scheme, the method further comprises:

determining the coherence coefficient of each area in the first interference pattern by adopting a coherence coefficient determination rule according to the coherence coefficient of each pixel of the first interference pattern;

each region in the first interferogram includes more than one pixel;

determining the coherence coefficient of each area in the second interference pattern by adopting a coherence coefficient determination rule according to the coherence coefficient of each pixel of the second interference pattern;

each region in the second interferogram includes more than one pixel.

In the foregoing solution, before performing the filtering process on the first and second interferograms, the method further includes:

and respectively carrying out land leveling effect removal treatment on the first interference pattern and the second interference pattern.

In the above scheme, the radar irradiation direction of the first satellite-borne SAR and the second satellite-borne SAR when the first track images the target forms a preset included angle with the radar irradiation direction of the second track when the second track images the target.

The embodiment of the invention also provides a high-precision DEM inversion device based on the double-base satellite-borne SAR, which comprises: the device comprises an acquisition module, an image processing module, a filtering module, a conversion module and a fusion module; wherein the content of the first and second substances,

the acquisition module is used for acquiring a first main image and a first auxiliary image which are respectively obtained by the first satellite-borne SAR and the second satellite-borne SAR in the first track for imaging the imaging target, and acquiring a second main image and a second auxiliary image which are respectively obtained by the first satellite-borne SAR and the second satellite-borne SAR in the second track for imaging the imaging target;

the image processing module is configured to process the first main image and the first auxiliary image by using a preset image processing rule, generate a first interference pattern corresponding to the first track, process the second main image and the second auxiliary image, and generate a second interference pattern corresponding to the second track;

the filtering module is used for respectively carrying out filtering processing on the first interference pattern and the second interference pattern according to a preset filtering rule and determining a coherent coefficient corresponding to each pixel on the first interference pattern and the second interference pattern;

the conversion module is used for respectively converting the filtered first interferogram and the filtered second interferogram into a first single-base-line DEM and a second single-base-line DEM according to a preset conversion rule;

and the fusion module is used for fusing the first single-baseline DEM and the second single-baseline DEM according to the coherence coefficient corresponding to each pixel on the first interference image and the second interference image by adopting a preset fusion rule to obtain a fused DEM of the imaging target.

In the foregoing solution, the image processing module is specifically configured to:

In the foregoing solution, the filtering module is specifically configured to:

In the foregoing scheme, the conversion module is specifically configured to:

determining the first single-baseline DEM corresponding to the first track according to the unwrapping phase of the first interferogram by adopting an InSAR height measurement model;

and determining the second single-baseline DEM corresponding to the second track according to the unwrapping phase of the second interferogram by adopting an InSAR height measurement model.

In the foregoing solution, the fusion module is specifically configured to:

and determining the proportion of the corresponding coherence coefficient of the second single-baseline DEM non-mask area to the sum of the corresponding coherence coefficients of the first single-baseline DEM and the second single-baseline DEM non-mask area as the weight of the second single-baseline DEM non-mask area for weighting fusion.

In the foregoing solution, the filtering module is further configured to:

each region in the first interferogram includes more than one pixel;

each region in the second interferogram includes more than one pixel.

In the foregoing solution, the image processing module is further configured to:

The embodiment of the invention also provides a storage medium, wherein an executable program is stored on the storage medium, and when the executable program is executed by a processor, the steps of any one of the methods for high-precision DEM inversion based on the double-base spaceborne SAR are realized.

The embodiment of the invention also provides a high-precision DEM inversion device based on the double-base spaceborne SAR, which comprises a processor, a memory and an executable program which is stored on the memory and can be run by the processor, wherein when the executable program is run by the processor, the step of any one of the methods is executed.

The high-precision DEM inversion method and device based on the double-base satellite-borne SAR provided by the embodiment of the invention are used for acquiring a first main image and a first auxiliary image which are obtained by respectively imaging an imaging target by a first satellite-borne SAR and a second satellite-borne SAR in a first track, and acquiring a second main image and a second auxiliary image which are obtained by respectively imaging the imaging target by the first satellite-borne SAR and the second satellite-borne SAR in a second track; processing the first main image and the first auxiliary image by adopting a preset image processing rule to generate a first interference image corresponding to the first track, and processing the second main image and the second auxiliary image to generate a second interference image corresponding to the second track; according to a preset filtering rule, respectively carrying out filtering processing on the first interference image and the second interference image, and determining a coherence coefficient corresponding to each pixel on the first interference image and the second interference image; according to a preset conversion rule, the filtered first interferogram and the filtered second interferogram are converted into a first single-base-line DEM and a second single-base-line DEM respectively; and fusing the first single-baseline DEM and the second single-baseline DEM according to the corresponding coherence coefficient of each pixel on the first interference image and the second interference image by adopting a preset fusion rule to obtain a fused DEM of the imaging target. Therefore, the single-baseline DEM of the two tracks is fused according to the coherence coefficients of the main image and the auxiliary image on the two tracks, the blank area in single-track imaging is filled, and the imaging quality of the DEM in the high-fall object dense area can be improved.

Drawings

FIG. 1 is a schematic flow chart of a high-precision DEM inversion method based on a bistatic spaceborne SAR in an embodiment of the invention;

FIG. 2 is a schematic diagram of a data processing framework of a high-precision DEM inversion method based on a bistatic spaceborne SAR in the embodiment of the invention;

FIG. 3 is a schematic diagram of weights calculated based on block search according to an embodiment of the present invention;

FIG. 4 is an InSAR altimetric model according to an embodiment of the present invention;

FIG. 5 is a schematic diagram of interference after rail-lifting registration and land-leveling according to an embodiment of the present invention;

FIG. 6 is a schematic diagram of interference after non-local filtering for rail lifting according to an embodiment of the present invention;

FIG. 7 is a schematic diagram of an elevated DEM after masking a low coherence region according to an embodiment of the present invention;

FIG. 8 is a diagram illustrating a reduced-track DEM after a low coherence region mask, in accordance with an embodiment of the present invention;

FIG. 9 is a schematic diagram of a merged DEM according to an embodiment of the present invention;

fig. 10 is a schematic structural diagram of a high-precision DEM inversion device based on bistatic spaceborne SAR in the embodiment of the present invention.

Detailed Description

In the embodiment of the invention, a first main image and a first auxiliary image obtained by respectively imaging an imaging target by a first satellite-borne SAR and a second satellite-borne SAR in a first track are obtained, and a second main image and a second auxiliary image obtained by respectively imaging the imaging target by the first satellite-borne SAR and the second satellite-borne SAR in a second track are obtained; processing the first main image and the first auxiliary image by adopting a preset image processing rule to generate a first interference image corresponding to the first track, and processing the second main image and the second auxiliary image to generate a second interference image corresponding to the second track; according to a preset filtering rule, respectively carrying out filtering processing on the first interference image and the second interference image, and determining a coherence coefficient corresponding to each pixel on the first interference image and the second interference image; according to a preset conversion rule, the filtered first interferogram and the filtered second interferogram are converted into a first single-base-line DEM and a second single-base-line DEM respectively; and fusing the first single-baseline DEM and the second single-baseline DEM according to the corresponding coherence coefficient of each pixel on the first interference image and the second interference image by adopting a preset fusion rule to obtain a fused DEM of the imaging target.

The high-precision DEM inversion method based on the double-base spaceborne SAR provided by the embodiment of the invention is shown in figure 1 and comprises the following steps:

step 110: acquiring a first main image and a first auxiliary image which are obtained by respectively imaging an imaging target by a first satellite-borne SAR and a second satellite-borne SAR in a first track, and acquiring a second main image and a second auxiliary image which are obtained by respectively imaging the imaging target by the first satellite-borne SAR and the second satellite-borne SAR in a second track;

here, the first and second on-board SAR may be on-board SAR of two satellites in a two-satellite distributed system, such as on-board SAR of TerraSAR-X satellite and on-board SAR of TanDEM-X satellite; the first satellite and the second satellite may be a primary satellite and a secondary satellite in a two-satellite distributed system; the dual-satellite system composed of the main satellite and the auxiliary satellite flies in a close formation mode to image the ground, for example, the distance between a Terras SAR-X satellite and a TanDEM-X satellite is 200-300; the imaging target may be a city, etc.;

imaging of the spaceborne SAR is Single view Complex (SLC); the SAR with the first main image and the second main image as main stars are respectively imaged on a first track and a second track to obtain an SLC; the first auxiliary image and the second auxiliary image are SLCs obtained by imaging SAR of an auxiliary satellite on a first track and a second track respectively;

in the DEM data processing framework shown in fig. 2, the first track is located relatively high and the second track is located relatively low, the image at the first track is called an ascending track SLC and the image at the second track is called a descending track SLC. The main satellite and the auxiliary satellite respectively image the target on the first track to obtain an ascending rail SLC1 and an ascending rail SLC2, and respectively image the target on the second track to obtain a descending rail SLC1 and a descending rail SLC 2; after the SLC is obtained, all images can be transmitted to a main satellite for processing, and can also be transmitted to an auxiliary satellite or a ground station for processing.

Furthermore, a preset included angle is formed between the radar irradiation direction of the first satellite-borne SAR and the radar irradiation direction of the second satellite-borne SAR when the first orbit images the target and the radar irradiation direction when the second orbit images the target;

the first orbit and the second orbit refer to different orbit heights and positions of the double-star system running at different times; due to the fact that the earth rotates and the satellite flies around the earth, when the double-satellite system flies over the imaging target at different time points, the orbit heights are different, and the directions of the imaging target irradiated by the satellite radar are also different; when the imaging target is imaged in the first orbit, the radar irradiation direction of the first spaceborne SAR and the second spaceborne SAR when the target is imaged can be from east to west; when the imaging target is imaged by the second orbit, the radar irradiation direction when the target is imaged by the first and second spaceborne SAR may be from west to east. Therefore, the imaging target can be imaged in two directions, the situation that imaging can only be carried out on the slope-facing surface such as a mountain peak and a tall building due to the included angle between the track and the imaging target is avoided, and the probability of imaging on the back slope surface is increased.

Step 120: processing the first main image and the first auxiliary image by adopting a preset image processing rule to generate a first interference image corresponding to the first track, and processing the second main image and the second auxiliary image to generate a second interference image corresponding to the second track;

here, the preset image processing rule may be set according to imaging conditions of the first main image and the first auxiliary image, and the second main image and the second auxiliary image, such as a positional relationship between a satellite and an imaging target, and an imaging size; the interference pattern can be generated after the image registration of the main image and the auxiliary image. The registration can be carried out by adopting methods such as a maximum interference spectrum method, a phase difference image average fluctuation function method and the like; the first main image and the first auxiliary image after the registration may be conjugate multiplied to obtain a first interference image corresponding to the first track, and the second main image and the second auxiliary image after the registration may be conjugate multiplied to obtain a second interference image corresponding to the second track.

Further, carrying out image registration on the first main image and the first auxiliary image, and carrying out conjugate multiplication on the first main image and the first auxiliary image after registration processing to obtain a first interference image; carrying out image registration on the second main image and the second auxiliary image, and carrying out conjugate multiplication on the second main image and the second auxiliary image after registration processing to obtain a second interference image; here, the first main image and the first auxiliary image may be registered and the second main image and the second auxiliary image may be registered by using a correlation coefficient method, and then the first interference image and the second interference image may be generated;

specifically, the correlation coefficient method adopts a window-based registration method, and a matching window and a search window are respectively established on the main image and the auxiliary image for registration. And selecting some control points from the main image and the auxiliary image, and establishing a matching window and a searching window by taking the control points as centers. Sliding the matching window to obtain a complex correlation function of the window of the matching window and the window of the searching window so as to obtain an offset; the complex correlation function may be represented by expression (1):

where ρ is_cRepresenting a complex correlation function, s₁And s₂Representing the main and secondary images, (u, v) the coordinates of the center point of the sliding matching window, (m, n) the coordinates of the center point of the pixel, i.e. the coordinates of the control point, M, N the size of the image, and x represents taking the complex conjugate. | □ | represents taking magnitude values, "□" in | □ | may represent any expression. The size of the matching window may be determined according to processing speed and accuracy, such as 20x 20;

complex correlation function rho_cThe maximum value of (2) is the best matching point, the maximum value of the complex correlation function can be interpolated, for example, 0.5 pixel is used as a unit interpolation to obtain more accurate offset, the registration offset of other pixels is obtained according to the registration offset corresponding to the control point, and the registration offset can be obtained by using a second-order polynomial (2):

wherein (x, y) represents the pixel coordinates in the main image, (Δ x, Δ y) represents the offset of the auxiliary image relative to the main image, and a and b represent fitting coefficients, which can be obtained by fitting the offset obtained according to equation 1, wherein a₀And b₀Representing the coefficient of residency, the subscripts representing the dimension and order; a is₁₀Represents a 1-dimensional 0-order coefficient, a₂₁2-dimensional 1-order coefficients are shown, and so on, and will not be described herein.

Acquiring a first main image and a first auxiliary image after registration, and a second main image and a second auxiliary image, and performing conjugate multiplication on the first main image and the first auxiliary image after registration to obtain a first interference image corresponding to a first track; and carrying out conjugate multiplication on the registered second main image and the second auxiliary image to obtain a second interference image corresponding to the second track.

Further, before the filtering process is performed on the first and second interferograms, the method further comprises: respectively removing the land leveling effect on the first interference pattern and the second interference pattern;

here, the ground leveling effect can be removed by an orbital method; taking TSX/TDX as an example, TSX/TDX has accurate track parameters and other information, can calculate the flat ground phase according to the track parameters and the scene ground incident angle, and removes the flat ground effect by conjugate multiplication with an interferogram.

Step 130: according to a preset filtering rule, respectively carrying out filtering processing on the first interference image and the second interference image, and determining a coherence coefficient corresponding to each pixel on the first interference image and the second interference image;

due to the inherent characteristics of a coherent imaging system, system noise of a satellite-borne SAR and the like, different noises such as speckle noise, system noise and the like exist in the first interferogram and the second interferogram; therefore, the interference pattern needs to be filtered; the first interferogram and the second interferogram can be filtered by mean filtering, Goldstein filtering and the like;

the coherence coefficient corresponding to each pixel on the first interferogram is a coherence coefficient reflecting the coherence relationship between each pixel of the registered first main image and the first auxiliary image corresponding to each pixel on the first interferogram; the coherence coefficient corresponding to each pixel on the second interferogram is a coherence coefficient reflecting the coherence relationship between each pixel of the registered second main image and the second auxiliary image corresponding to each pixel on the second interferogram; the coherence coefficient corresponding to each pixel on the first interference image can be obtained by adopting methods such as multi-view processing and the like on adjacent pixels at corresponding positions of the first main image and the first auxiliary image after registration; the coherence coefficient corresponding to each pixel on the second interference pattern can be obtained by using methods such as multi-view processing and the like on adjacent pixels at corresponding positions of the second main image and the second auxiliary image after registration.

Further, respectively carrying out iterative estimation on each pixel in the first interference image and the second interference image by adopting a non-local filtering method; determining a coherence coefficient corresponding to each pixel of the first interference image by adopting a Goodman model of the first main image and the first auxiliary image after registration processing according to the similarity weight of each pixel of the first interference image in iterative estimation; determining a coherence coefficient corresponding to each pixel of the second interference image by adopting a Goodman model of the second main image and the second auxiliary image after registration processing according to the similarity weight of each pixel of the second interference image in iterative estimation;

specifically, the pixel model of each point of the two registered SLC complex images can be described by three parameters: reflection coefficient, true interference phase and coherence coefficient, assuming that z and z 'are corresponding pixel points on the two registration images, and according to the Goodman model, z and z' obey expression (3):

where Σ represents a 2 × 2 coherence coefficient matrix, which can be expressed by expression (4):

wherein, R and R' represent reflection coefficients, beta represents a real phase, D represents a coherence coefficient, and E represents an expected value;

the three parameters R, β, D can be jointly estimated by non-local filtering, assuming that the scattering properties of the two images are the same, i.e. R ═ R ', this assumption holds in regions with high coherence, while let a ═ z |, a ' ═ z ' |; phi ═ arg (zz'^*) For the phase after flattening, there is an expression (5):

Non-Local Interferometric synthetic aperture Radar (NL-InSAR, Non Local Interferometric synthetic Aperture Radar) estimates weights using Weight Maximum Likelihood Estimation (WMLE), in blocks s, the parameter estimation can be represented by expression (6):

let R ═ R', extend the formula into weighted maximum likelihood estimation, let Θ_s＝(R_s,β_s,D_s) And

the estimated R, β, D can be represented by expression (7):

wherein the content of the first and second substances,

w (s, t) is a similarity weight based on block search calculation in the filtering of the interferogram, namely, by adopting an iterative estimation method; the calculation method for implementing the filtering of the interferogram by using the iterative estimation method can be as shown in fig. 3:

suppose that a pixel point s to be estimated on an image I can be represented by expression (8):

v＝{v_s|s∈I} (8)

judging whether a pixel point s to be estimated is similar to a similar pixel point t, selecting a window with a certain size by taking s and t as centers, and judging whether the window with s and t as centers is similar; taking a pixel window with a pixel point s as a center as a target window 31, taking a pixel window with a pixel point t as a center as a similar window 32, searching by a preset route 34 in a selected search window 33, determining the similarity between each similar window 32 and the target window 31, wherein an estimated value of the pixel point s can be represented by an expression (9):

wherein v is_sRepresenting the pixel point to be estimated, v_tRepresenting a center pixel of the similarity window; the weights w (s, t) are selected based on the similarity of the pixel windows around the s, t pixels, such that 0 ≦ w (s, t ≦ 1, and

the size of the target window and similar windows may be determined based on processing speed and accuracy, such as 3X 3;

the de-dot estimation value of each pixel is obtained, namely, the filtering is completed, wherein, the iterative estimation of the current pixel is a filtering processing process.

The NL-InSAR is used for filtering phase noise, the noise reduction performance is good, and the subsequent phase unwrapping is facilitated.

Further, according to the coherence coefficient of each pixel of the first interference pattern, determining the coherence coefficient of each area in the first interference pattern by adopting a coherence coefficient determination rule; each region in the first interferogram includes more than one pixel; determining the coherence coefficient of each area in the second interference pattern by adopting a coherence coefficient determination rule according to the coherence coefficient of each pixel of the second interference pattern; each region in the second interferogram comprises more than one pixel;

in order to improve the subsequent calculation speed, the coherence coefficient of different areas in the interferogram can be aimed at, wherein the area can be one pixel in the interferogram or the size of the area can be preset; the size of the region can be set according to the accuracy requirement of the interference pattern; by adopting a coherence coefficient determination rule, the determination can be performed according to the pixel relation between the region and the region, if the region is a pixel, the regional coherence coefficient is changed to be the coherence coefficient of the pixel, and if the region is larger than a pixel, the coherence coefficient of the region can perform arithmetic operation such as arithmetic mean and the like on the coherence coefficients of the pixels in the region.

Step 140: according to a preset conversion rule, the filtered first interferogram and the filtered second interferogram are converted into a first single-base-line DEM and a second single-base-line DEM respectively;

here, the DEM may be generated according to the interferogram and attitude data of the spaceborne SAR, and a single baseline DEM corresponding to the orbit may be generated according to the first interferogram and the second interferogram, respectively; the preset conversion rule can be set according to the track positions of the first satellite-borne SAR and the second satellite-borne SAR, the DEM precision requirement, the operation time requirement and the like; the DEM may be generated in a conventional manner.

Further, the preset DEM simulation interference phase and the interference phase of the first interference image can be differentiated to obtain a first differential phase; masking the first differential phase of a region with the coherence coefficient lower than a first preset coherence threshold value in the first interferogram, and unwrapping the masked first differential phase; in the unwrapped first differential phase, interpolating the mask region by adopting a preset interpolation rule, and adding the interpolated unwrapped first differential phase to a preset DEM simulation interference phase to obtain an unwrapped phase of a first interference pattern; determining the first single-baseline DEM corresponding to the first track according to the unwrapping phase of the first interferogram by adopting an InSAR height measurement model; differentiating the preset DEM simulation interference phase and the interference phase of the second interference pattern to obtain a second differential phase; masking the second differential phase of the region with the coherence coefficient lower than a first preset coherence threshold value in the second interference pattern, and unwrapping the masked second differential phase; in the unwrapped second differential phase, interpolating the mask region by adopting a preset interpolation rule, and adding the interpolated unwrapped second differential phase to a preset DEM simulation interference phase to obtain an unwrapped phase of a second interference pattern; determining the second single-baseline DEM corresponding to the second track according to the unwrapping phase of the second interferogram by adopting the InSAR height measurement model;

here, the preset DEM may be an existing DEM of an imaging target, and interference phases to the preset DEM are simulated according to the SAR system parameters of the first and second spaceborne SAR and the orbit parameters of the first and second spaceborne SAR; carrying out differential interference processing by adopting a simulation interference phase of a preset DEM and an interference pattern to obtain a winding differential phase; in order to ensure the reliability of phase unwrapping, a low-coherence region in an interferogram can be masked, and then a differential phase is unwrapped; a first preset coherence threshold value can be preset, areas with coherence coefficients lower than the first preset coherence threshold value in each area of the interference pattern are masked, and interpolation is carried out on the masked areas; adding the unwrapped differential phase after the mask processing to a preset simulated interference phase of the DEM to obtain an unwrapped phase of the interferogram, wherein the processing can be performed on both the first interferogram and the second interferogram to obtain the unwrapped phases of the first interferogram and the second interferogram respectively;

after the unwrapped phase is obtained, the height can be inverted by using an InSAR geometric relationship model, and the spatial geometric relationship of InSAR digital elevation reconstruction is shown in FIG. 4, wherein: a. the₁And A₂Respectively a first satellite-borne SAR and a second satellite-borne SAR; h denotes the platform height, A₁And A₂Respectively representing the phase centers of two antennas, alpha is a base line angle, B is the length of an interference base line, P is a certain target point on the ground, r₁And r₂＝r₁+ delta r is the distance from the phase center of two antennas to the target point P in the zero Doppler plane, theta is the visual angle from the phase center of the main antenna to the target, the terrain elevation is represented by h, y_PThe distance of the point P is shown, and the coordinate of the distance of the point Y is shown. From this spatial geometry, the elevation h of the terrain can be derived and expressed by expression (10):

and respectively obtaining elevation data corresponding to the first interferogram and the second interferogram through the unwrapping phases of the first interferogram and the second interferogram, and then coding together according to the geography of an imaging target to obtain a first single-baseline DEM and a second single-baseline DEM corresponding to the first track and the second track respectively.

Here, the external preset DEM is used for differential processing, so that the difficulty of phase unwrapping can be reduced, and the reliability of phase unwrapping can be ensured.

Step 150: fusing the first single-baseline DEM and the second single-baseline DEM according to the corresponding coherence coefficient of each pixel on the first interference image and the second interference image by adopting a preset fusion rule to obtain a fused DEM of the imaging target;

the first single-baseline DEM and the second single-baseline DEM have more data holes; cannot reflect the true topographic features; here, the first single baseline DEM and the second single baseline DEM may be merged; the preset fusion rule can be set according to different areas in the baseline DEM, such as mutual supplement of imaging shielding areas and the like, weighted fusion of normal imaging areas and the like;

further, masking an imaging abnormal region and a region with a corresponding coherence coefficient lower than a second preset coherence threshold value in the first single-baseline DEM and the second single-baseline DEM respectively; carrying out weighted fusion on the first single-baseline DEM non-mask area and the second single-baseline DEM non-mask area, supplementing the first single-baseline DEM normal imaging area with the second single-baseline DEM mask area at the corresponding position, and supplementing the second single-baseline DEM normal imaging area with the first single-baseline DEM mask area at the corresponding position to obtain a fused DEM of an imaging target; determining the proportion of the corresponding coherence coefficient of the first single-baseline DEM non-mask area to the sum of the corresponding coherence coefficients of the first single-baseline DEM non-mask area and the second single-baseline DEM non-mask area as a weight for performing weighted fusion on the first single-baseline DEM non-mask area; determining the proportion of the corresponding coherence coefficient of the second single-baseline DEM non-mask area to the sum of the corresponding coherence coefficients of the first single-baseline DEM non-mask area and the second single-baseline DEM non-mask area as a weight for performing weighted fusion on the second single-baseline DEM non-mask area;

specifically, the imaging abnormal region can be a perspective shrinkage, overlapping and shadow region in a single base line DEM; the corresponding positions refer to positions with the same geographic coordinates. Masking perspective shrinkage, overlapping and shading of the first single-baseline DEM and the second single-baseline DEM, and meanwhile masking an area with lower coherence, namely a lower coherence, by adopting an unbiased estimated coherence coefficient obtained by NL filtering; then, supplementing null value areas by using the characteristic that the first track and the second track have opposite irradiation directions; and weighting and fusing the non-null value areas to obtain a fused DEM. The weighted value in the weighted fusion can be determined according to the coherence coefficient corresponding to each pixel or area of the first single-baseline DEM and the second single-baseline DEM, for example, the coherence coefficient corresponding to each pixel or area of the first single-baseline DEM is D_AThe second sheetThe coherent coefficient corresponding to each pixel or area of the base line DEM is D_BThen the first single baseline DEM weighting value in the weighted fusion is thenThe second single baseline DEM weighting value isIn order to simplify the operation, the first preset coherence threshold and the second preset coherence threshold may take the same value;

therefore, images obtained in different track directions of the first track and the second track are used for complementing 'blank' areas such as overlapping and covering areas, for example, the area of the upper back slope surface of the descending track DEM is used for filling the blank area of the upper slope surface of the ascending track DEM, high-fall objects such as building dense areas and high-resolution high-precision urban DEM with rich detailed features are reconstructed. In addition, because the embodiment of the invention directly uses one-time multi-view processing operation, repeated iteration processing is avoided, the processing is simple and quick, and the error transmission problem during multiple iterations is avoided.

The positive effects produced by the present invention will be described in further detail with reference to specific examples below;

in the example, TSX/TDX lifting rail image data are adopted to verify the high-precision DEM inversion method based on the double-base satellite-borne SAR. FIG. 5 is an interferogram after two ascending rails SLC are registered and land removed; fig. 6 shows an interferogram obtained by non-local mean filtering, which can find that the noise of the modified image is well removed and the edge of the building is clear. Fig. 7 and 8 show an up-track DEM and a down-track DEM with a resolution of 4 m obtained by NL-InSAR, respectively, and the null area is indicated by a null value. Fig. 9 shows the DEM after the merging of the ascending and descending tracks, and it can be seen that compared with the DEM after the ascending or descending tracks, the invalid area of the DEM after the merging is obviously reduced, and the statistical result shows that the invalid point ratio is reduced to 1.34% from 4.93% and 4.52% of the ascending and descending tracks, and meanwhile, compared with the SRTM data, the Root Mean Square Error (RMSE) after the merging is reduced to 5.06 meters from 6.37 meters of the ascending track and 5.61 meters of the descending track.

As shown in fig. 10, the high-precision DEM inversion apparatus based on bistatic spaceborne SAR provided in the embodiment of the present invention includes: the system comprises an acquisition module 101, an image processing module 102, a filtering module 103, a conversion module 104 and a fusion module 105; wherein the content of the first and second substances,

the acquiring module 101 is configured to acquire a first main image and a first auxiliary image, which are obtained by respectively imaging an imaging target by a first satellite-borne SAR and a second satellite-borne SAR in a first orbit, and acquire a second main image and a second auxiliary image, which are obtained by respectively imaging the imaging target by the first satellite-borne SAR and the second satellite-borne SAR in a second orbit;

imaging of the spaceborne SAR is SLC; the SAR with the first main image and the second main image as main stars are respectively imaged on a first track and a second track to obtain an SLC; the first auxiliary image and the second auxiliary image are SLCs obtained by imaging SAR of an auxiliary satellite on a first track and a second track respectively;

The image processing module 102 is configured to process the first main image and the first auxiliary image by using a preset image processing rule, generate a first interference pattern corresponding to the first track, process the second main image and the second auxiliary image, and generate a second interference pattern corresponding to the second track;

specifically, the correlation coefficient method adopts a window-based registration method, and a matching window and a search window are respectively established on the main image and the auxiliary image for registration. And selecting some control points from the main image and the auxiliary image, and establishing a matching window and a searching window by taking the control points as centers. Sliding the matching window to obtain a complex correlation function of the window of the matching window and the window of the searching window so as to obtain an offset; the complex correlation function may be represented by expression (1);

complex correlation function rho_cThe maximum value of the complex correlation function is the optimal matching point, interpolation can be carried out on the maximum value of the complex correlation function, for example, interpolation is carried out by adopting 0.5 pixel as a unit to obtain more accurate offset, registration offset of other pixels is obtained according to registration offset corresponding to the control point, and the registration offset can be obtained by utilizing a second-order polynomial (2); wherein (x, y) represents the pixel coordinates in the main image, (Δ x, Δ y) represents the offset of the auxiliary image relative to the main image, and a and b represent fitting coefficients, which can be obtained by fitting the offset obtained according to equation 1, wherein a₀And b₀Representing the coefficient of residency, the subscripts representing the dimension and order; a is₁₀Represents a 1-dimensional 0-order coefficient, a₂₁2-dimensional 1-order coefficients are shown, and so on, and will not be described herein.

The filtering module 103 is configured to perform filtering processing on the first interference pattern and the second interference pattern respectively according to a preset filtering rule, and determine a coherence coefficient corresponding to each pixel on the first interference pattern and the second interference pattern;

specifically, the pixel model of each point of the two registered SLC complex images can be described by three parameters: the reflection coefficient, the real interference phase and the coherence coefficient, and if z and z 'are corresponding pixel points on the two registration images, according to the Goodman model, the z and z' obey an expression (3); where, Σ represents a 2 × 2 coherence coefficient matrix, which can be represented by expression (4); wherein, R and R' represent reflection coefficients, beta represents a real phase, D represents a coherence coefficient, and E represents an expected value;

the three parameters R, β, D can be jointly estimated by non-local filtering, assuming that the scattering properties of the two images are the same, i.e. R ═ R ', this assumption holds in regions with high coherence, while let a ═ z |, a ' ═ z ' |; phi ═ arg (zz'^*) The phase after the flat ground is removed, an expression (5) is provided;

NL-InSAR uses WMLE to estimate weights, in block s, the parameter estimates can be represented by expression (6);

the estimated R, β, D can be represented by expression (7); wherein the content of the first and second substances, w (s, t) is a similarity weight based on block search calculation in the filtering of the interferogram, namely, by adopting an iterative estimation method; the calculation mode for implementing the filtering of the interferogram by using the iterative estimation method can be as shown in fig. 3; the pixel point s to be estimated on the image I can be represented by an expression (8); judging whether a pixel point s to be estimated is similar to a similar pixel point t, selecting a window with a certain size by taking s and t as centers, and judging whether the window with s and t as centers is similar; taking a pixel window with a pixel point s as a center as a target window 31, taking a pixel window with a pixel point t as a center as a similar window 32, searching by a preset route 34 in a selected search window 33, determining the similarity between each similar window 32 and the target window 31, and expressing the estimation value of the pixel point s by an expression (9); wherein v is_sRepresenting the pixel point to be estimated, v_tRepresenting a center pixel of the similarity window; the weights w (s, t) are selected based on the similarity of the pixel windows around the s, t pixels, such that 0 ≦ w (s, t ≦ 1, and

The conversion module 104 is configured to convert the filtered first interferogram and the filtered second interferogram into a first single-baseline DEM and a second single-baseline DEM, respectively, according to a preset conversion rule;

after the unwrapped phase is obtained, the height can be inverted by using an InSAR geometric relationship model, and the spatial geometric relationship of InSAR digital elevation reconstruction is shown in FIG. 4, wherein: a. the₁And A₂Respectively a first satellite-borne SAR and a second satellite-borne SAR; h denotes the platform height, A₁And A₂Respectively representing the phase centers of two antennas, alpha is a base line angle, B is the length of an interference base line, P is a certain target point on the ground, r₁And r₂＝r₁+ delta r is the distance from the phase center of two antennas to the target point P in the zero Doppler plane, theta is the visual angle from the phase center of the main antenna to the target, the terrain elevation is represented by h, y_PThe distance of the point P is shown, and the coordinate of the distance of the point Y is shown. According to the space tableThe relationship can deduce that the elevation h of the terrain can be expressed by an expression (10);

The fusion module 105 is configured to: fusing the first single-baseline DEM and the second single-baseline DEM according to the corresponding coherence coefficient of each pixel on the first interference image and the second interference image by adopting a preset fusion rule to obtain a fused DEM of the imaging target;

specifically, the imaging abnormal region can be a perspective shrinkage, overlapping and shadow region in a single base line DEM; the corresponding positions refer to positions with the same geographic coordinates. Masking perspective shrinkage, overlapping and shading of the first single-baseline DEM and the second single-baseline DEM, and meanwhile masking an area with lower coherence, namely a lower coherence, by adopting an unbiased estimated coherence coefficient obtained by NL filtering; then, supplementing null value areas by using the characteristic that the first track and the second track have opposite irradiation directions; and weighting and fusing the non-null value areas to obtain a fused DEM. The weighted value in the weighted fusion can be determined according to the coherence coefficient corresponding to each pixel or area of the first single-baseline DEM and the second single-baseline DEM, for example, the coherence coefficient corresponding to each pixel or area of the first single-baseline DEM is D_AThe coherent coefficient corresponding to each pixel or area of the second single-baseline DEM is D_BThen the first single baseline DEM weighting value in the weighted fusion is then

The second single baseline DEM weighting value is

In order to simplify the operation, the first preset coherence threshold and the second preset coherence threshold may take the same value;

In practical applications, the obtaining module 101, the image processing module 102, the filtering module 103, the transforming module 104 and the fusing module 105 may be implemented by a Central Processing Unit (CPU), a microprocessor unit (MPU), a Digital Signal Processor (DSP), a Field Programmable Gate Array (FPGA), or the like in a satellite system or a ground system.

The storage medium provided by the embodiment of the invention stores an executable program, and the executable program is executed by a processor to realize the step of the high-precision DEM inversion method based on the double-base spaceborne SAR; as shown in fig. 1, the high-precision DEM inversion method based on the bistatic spaceborne SAR includes the steps of:

complex correlation function rho_cThe maximum value of the complex correlation function is the optimal matching point, interpolation can be carried out on the maximum value of the complex correlation function, for example, interpolation is carried out by adopting 0.5 pixel as a unit to obtain more accurate offset, registration offset of other pixels is obtained according to registration offset corresponding to the control point, and the registration offset can be obtained by utilizing a second-order polynomial (2); wherein (x, y) represents pixel coordinates in the main image, (Δ x, Δ y) represents an offset amount of the sub image with respect to the main image, and a and b represent fitting coefficients that can be fitted according to the offset amount obtained by equation 1To, wherein, a₀And b₀Representing the coefficient of residency, the subscripts representing the dimension and order; a is₁₀Represents a 1-dimensional 0-order coefficient, a₂₁2-dimensional 1-order coefficients are shown, and so on, and will not be described herein.

the estimated R, β, D can be represented by expression (7); wherein the content of the first and second substances,

w (s, t) is a similarity weight based on block search calculation in the filtering of the interferogram, namely, by adopting an iterative estimation method; the calculation mode for implementing the filtering of the interferogram by using the iterative estimation method can be as shown in fig. 3; the pixel point s to be estimated on the image I can be represented by an expression (8); judging whether a pixel point s to be estimated is similar to a similar pixel point t, selecting a window with a certain size by taking s and t as centers, and judging whether the window with s and t as centers is similar; taking a pixel window with a pixel point s as a center as a target window 31, taking a pixel window with a pixel point t as a center as a similar window 32, searching by a preset route 34 in a selected search window 33, determining the similarity between each similar window 32 and the target window 31, and expressing the estimation value of the pixel point s by an expression (9); wherein v is_sRepresenting the pixel point to be estimated, v_tRepresenting a center pixel of the similarity window; the weights w (s, t) are selected based on the similarity of the pixel windows around the s, t pixels, such that 0 ≦ w (s, t ≦ 1, and

after the unwrapped phase is obtained, the height can be inverted by using an InSAR geometric relationship model, and the spatial geometric relationship of InSAR digital elevation reconstruction is shown in FIG. 4, wherein: a. the₁And A₂Respectively a first satellite-borne SAR and a second satellite-borne SAR; h denotes the platform height, A₁And A₂Respectively representing the phase centers of two antennas, alpha is a base line angle, B is the length of an interference base line, and P is the length on the groundA certain target point of r₁And r₂＝r₁+ delta r is the distance from the phase center of two antennas to the target point P in the zero Doppler plane, theta is the visual angle from the phase center of the main antenna to the target, the terrain elevation is represented by h, y_PThe distance of the point P is shown, and the coordinate of the distance of the point Y is shown. The elevation h of the terrain can be deduced according to the space geometric relation and can be expressed by an expression (10);

The second single baseline DEM weighting value is

The high-precision DEM inversion device based on the double-base spaceborne SAR comprises a processor, a memory and an executable program which is stored on the memory and can be operated by the processor, wherein the processor executes the high-precision DEM inversion method based on the double-base spaceborne SAR when the executable program is operated; as shown in fig. 1, the high-precision DEM inversion method based on the bistatic spaceborne SAR includes the steps of:

complex correlation function rho_cThe maximum value of the complex correlation function is the optimal matching point, interpolation can be carried out on the maximum value of the complex correlation function, for example, interpolation is carried out by adopting 0.5 pixel as a unit to obtain more accurate offset, registration offset of other pixels is obtained according to registration offset corresponding to the control point, and the registration offset can be obtained by utilizing a second-order polynomial (2); wherein (x, y) represents the pixel coordinates in the main image, (Δ x, Δ y) represents the offset of the sub-image with respect to the main image,a and b represent fitting coefficients, which can be obtained by fitting the offset obtained according to equation 1, wherein a₀And b₀Representing the coefficient of residency, the subscripts representing the dimension and order; a is₁₀Represents a 1-dimensional 0-order coefficient, a₂₁2-dimensional 1-order coefficients are shown, and so on, and will not be described herein.

let R ═ R'The formula is extended to the weighted maximum likelihood estimation, let Θ_s＝(R_s,β_s,D_s) And

after the unwrapped phase is obtained, the height can be inverted by using an InSAR geometric relationship model, and the spatial geometric relationship of InSAR digital elevation reconstruction is shown in FIG. 4, wherein: a. the₁And A₂Respectively a first satellite-borne SAR and a second satellite-borne SAR; h denotes the platform height, A₁And A₂Representing the phases of two antennas respectivelyPosition center, alpha is the base line angle, B is the length of the interference base line, P is a certain target point on the ground, r₁And r₂＝r₁+ delta r is the distance from the phase center of two antennas to the target point P in the zero Doppler plane, theta is the visual angle from the phase center of the main antenna to the target, the terrain elevation is represented by h, y_PThe distance of the point P is shown, and the coordinate of the distance of the point Y is shown. The elevation h of the terrain can be deduced according to the space geometric relation and can be expressed by an expression (10);

The second single baseline DEM weighting value is

The above description is only for the purpose of illustrating the preferred embodiments of the present invention and is not to be construed as limiting the scope of the present invention, and any modifications, equivalents and improvements made within the spirit and principle of the present invention should be included in the scope of the present invention.

Claims

1. A high-precision DEM inversion method based on double-base spaceborne SAR is characterized by comprising the following steps:

acquiring a first main image and a first auxiliary image which are obtained by respectively imaging an imaging target by a first satellite-borne Synthetic Aperture Radar (SAR) and a second satellite-borne SAR in a first track, and acquiring a second main image and a second auxiliary image which are obtained by respectively imaging the imaging target by the first satellite-borne SAR and the second satellite-borne SAR in a second track;

fusing the first single-baseline DEM and the second single-baseline DEM according to the corresponding coherence coefficient of each pixel on the first interference image and the second interference image by adopting a preset fusion rule to obtain a fused DEM of the imaging target; fusing the first single-baseline DEM and the second single-baseline DEM according to the coherence coefficient corresponding to each pixel on the first interference image and the second interference image by adopting a preset fusion rule to obtain a fused DEM of the imaging target; the method comprises the following steps: respectively masking an imaging abnormal region and a region with a corresponding coherence coefficient lower than a first preset coherence threshold value in the first single-baseline DEM and the second single-baseline DEM; carrying out weighted fusion on the first single-baseline DEM non-mask area and the second single-baseline DEM non-mask area, supplementing the first single-baseline DEM normal imaging area with the second single-baseline DEM mask area at the corresponding position, and supplementing the second single-baseline DEM normal imaging area with the first single-baseline DEM mask area at the corresponding position to obtain a fused DEM of an imaging target; determining the proportion of the corresponding coherence coefficient of the first single-baseline DEM non-mask area to the sum of the corresponding coherence coefficients of the first single-baseline DEM non-mask area and the second single-baseline DEM non-mask area as a weight for performing weighted fusion on the first single-baseline DEM non-mask area; and determining the proportion of the corresponding coherence coefficient of the second single-baseline DEM non-mask area to the sum of the corresponding coherence coefficients of the first single-baseline DEM non-mask area and the second single-baseline DEM non-mask area as a weight for performing weighted fusion on the second single-baseline DEM non-mask area.

2. The method according to claim 1, wherein the processing the first primary image and the first secondary image according to the preset image processing rule to generate a first interference map corresponding to the first track, and processing the second primary image and the second secondary image to generate a second interference map corresponding to the second track comprises:

3. The method according to claim 2, wherein the filtering the first and second interferograms respectively according to a preset filtering rule and determining the coherence coefficient corresponding to each pixel on the first and second interferograms comprises:

4. The method according to claim 3, wherein the filtered first and second interferograms are converted into first and second single-baseline DEMs, respectively, according to a preset conversion rule; the method comprises the following steps:

masking the first differential phase of a region with the coherence coefficient lower than a second preset coherence threshold value in the first interferogram, and unwrapping the masked first differential phase;

masking the second differential phase of a region with the coherence coefficient lower than a second preset coherence threshold value in the second interference pattern, and unwrapping the masked second differential phase;

5. The method of claim 1 or 4, further comprising:

each region in the first interferogram includes more than one pixel;

each region in the second interferogram includes more than one pixel.

6. The method according to any one of claims 1 to 4, wherein prior to the filtering the first and second interferograms, the method further comprises:

7. The method as claimed in any one of claims 1 to 4, wherein the radar irradiation direction of the first and second spaceborne SAR when the target is imaged in the first orbit is at a preset angle with respect to the radar irradiation direction when the target is imaged in the second orbit.

8. A high-precision DEM inversion device based on double-base satellite-borne SAR is characterized by comprising: the device comprises an acquisition module, an image processing module, a filtering module, a conversion module and a fusion module; wherein the content of the first and second substances,

the acquisition module is used for acquiring a first main image and a first auxiliary image which are obtained by respectively imaging an imaging target by a first satellite-borne SAR and a second satellite-borne SAR in a first track, and acquiring a second main image and a second auxiliary image which are obtained by respectively imaging the imaging target by the first satellite-borne SAR and the second satellite-borne SAR in a second track;

the fusion module is used for fusing the first single-baseline DEM and the second single-baseline DEM according to the coherence coefficient corresponding to each pixel on the first interference image and the second interference image by adopting a preset fusion rule to obtain a fused DEM of the imaging target; wherein, the fusion module is specifically configured to: respectively masking an imaging abnormal region and a region with a corresponding coherence coefficient lower than a first preset coherence threshold value in the first single-baseline DEM and the second single-baseline DEM; carrying out weighted fusion on the first single-baseline DEM non-mask area and the second single-baseline DEM non-mask area, supplementing the first single-baseline DEM normal imaging area with the second single-baseline DEM mask area at the corresponding position, and supplementing the second single-baseline DEM normal imaging area with the first single-baseline DEM mask area at the corresponding position to obtain a fused DEM of an imaging target; determining the proportion of the corresponding coherence coefficient of the first single-baseline DEM non-mask area to the sum of the corresponding coherence coefficients of the first single-baseline DEM non-mask area and the second single-baseline DEM non-mask area as a weight for performing weighted fusion on the first single-baseline DEM non-mask area; and determining the proportion of the corresponding coherence coefficient of the second single-baseline DEM non-mask area to the sum of the corresponding coherence coefficients of the first single-baseline DEM and the second single-baseline DEM non-mask area as the weight of the second single-baseline DEM non-mask area for weighting fusion.

9. The apparatus of claim 8, wherein the image processing module is specifically configured to:

10. The apparatus of claim 9, wherein the filtering module is specifically configured to:

11. The apparatus according to claim 10, wherein the conversion module is specifically configured to:

12. The apparatus of claim 8 or 11, wherein the filtering module is further configured to:

each region in the first interferogram includes more than one pixel;

each region in the second interferogram includes more than one pixel.

13. The apparatus of any of claims 8 to 11, wherein the image processing module is further configured to:

14. The apparatus as claimed in any one of claims 8 to 11, wherein the radar irradiation direction of the first and second SAR at the first orbit when imaging the target is at a predetermined angle with respect to the radar irradiation direction at the second orbit when imaging the target.

15. A storage medium having stored thereon an executable program, which when executed by a processor implements the steps of the double-base SAR-based high-precision DEM inversion method according to any one of claims 1 to 7.

16. A double-base spaceborne SAR-based high-precision DEM inversion apparatus comprising a processor, a memory and an executable program stored on the memory and capable of being executed by the processor, wherein the processor executes the executable program to perform the steps of the double-base spaceborne SAR-based high-precision DEM inversion method according to any one of claims 1 to 7.