CN108663678B - Multi-baseline InSAR phase unwrapping algorithm based on mixed integer optimization model - Google Patents

Abstract

The invention discloses a multi-baseline InSAR phase unwrapping algorithm based on a mixed integer optimization model. The implementation steps of the invention are: (1) inputting primary and auxiliary SAR complex image data; (2) generating a complex interference phase map through interference processing; (3) Obtain the local optimal window; (4) Construct a mixed integer optimization model according to the InSAR system parameters and the local optimal window parameters; (5) Calculate the fuzzy integer of the complex interferometric phase diagram; (6) Calculate the absolute value of each interferometric phase diagram Phase; (7) Output the DEM of the entire scene. In the case of multiple baselines, the invention can complete high-quality phase unwrapping of the interferometric InSAR image after registration, meet the actual engineering performance requirements of high-quality interferometric synthetic aperture radar InSAR processing, and obtain high-quality surveying and mapping products.

Description

Multibaseline InSAR Phase Unwrapping Algorithm Based on Mixed Integer Optimization Model

technical field

The invention belongs to the technical field of communication, and further relates to data processing for generating a digital elevation model of the earth surface by an Interferometric Synthetic Aperture Radar (InSAR) in the technical field of radar detection. The present invention can be used in multi-baseline InSAR, because SAR images with different incident angles generate two or more interferograms in the same area, and then process the generated interferometric phase images to extract the terrain height of each scattering unit in the observation scene value.

Background technique

Interferometric Synthetic Aperture Radar combines two or more coherent SAR images from different viewing angles in the same scene on the basis of the traditional synthetic aperture radar (Synthetic Aperture Radar, SAR) obtaining two-dimensional information in the range dimension and azimuth dimension of the actual three-dimensional scene, through Interference processing technology to obtain 3D terrain information of the target scene.

Conventional single-baseline interferometry systems with two antenna phase centers have been shown to be capable of providing high-precision DEMs for general terrain during interferometric synthetic aperture radar InSAR data acquisition. However, for complex terrain, including steep slopes or discontinuous surfaces (such as man-made buildings, canyons, and steep mountains, etc.), the measurement performance drops due to severe phase undersampling (i.e., the adjacent phase difference of the interference may be greater than π). Seriously, these regions will be "blind zones" for single-baseline interferometry, making phase unwrapping impossible or reducing the accuracy of DEM generation. In order to solve this problem, a multi-baseline interferometric synthetic aperture radar phase unwrapping technique based on mixed integer optimization model is proposed. By combining the center pixel with its adjacent position pixels, it is another strategy to improve the performance of CRT estimation. The phase unwrapping accuracy and the acquisition of more accurate absolute phase values can improve the measurement performance of the InSAR system.

Xu, in the paper "Phase-unwrapping of SAR interferograms with multi-frequencyor multi-baseline" (IEEE Int.Geosci.and Remote Sens.Symp., Pasadena, CA, 1994, pp.730-732) first proposed the Chinese residual Theorem method (GRT), projection method and linear combination method are three basic methods to improve the accuracy of phase unwrapping and improve the noise suppression effect.

Liu Huitao, in the literature "ANovel Mixed-norm Multibaseline Phase Unwrapping Algorithm Based on Linear Programming" (IEEE Geoscience and Remote Sensing Letters, 2015, 12(5), pp.1086-1090) proposed an effective baseline interference method, which overcomes the The disadvantage of the phase continuation assumption, and the absolute phase of complex terrain can be generated.

Yu Hanwen, in the document "Multi-baseline InSAR Phase Unwrapping Technique Using L~1 Norm" (Journal of Xidian University, 2013, 40(04): 37-41), proposed a method based on L1 norm. Multibaseline phase unwrapping algorithm. The algorithm integrates the relationship between multiple interferometric phase images obtained by multi-baseline interferometric synthetic aperture radar into the traditional single-baseline L1-norm phase unwrapping optimization model by using the L1-norm optimization model, which can be applied to Mapping of complex terrain.

Among the above algorithms, the single-pixel-based phase unwrapping method is inevitably seriously affected by interference phase noise, and the algorithm is very sensitive to noise, so it cannot achieve a good application in practice. Although the multi-baseline phase unwrapping method of the proposed method can overcome the influence of noise to a certain extent, the influence of the local terrain on the model accuracy is not considered in the model construction. With the increase of the phase noise level, its performance also drops sharply.

SUMMARY OF THE INVENTION

The present invention is aimed at the engineering application of the interferometric synthetic aperture radar InSAR data processing technology to terrain measurement in the prior art. The present invention can be used for multi-baseline InSAR phase unwrapping to further improve the accuracy of DEM establishment on the earth's surface. The existing interferometric phase diagram When the entanglement method cannot meet the required performance requirements, a multi-baseline InSAR phase unwrapping algorithm based on a mixed integer optimization model is proposed. By obtaining the optimal window, constructing the optimization model, calculating the fuzzy integer of the interferometric phase map, and finally using the unwrapped phase to output the DEM of the entire scene.

For achieving the above object, the main steps of the present invention are as follows:

(1) Input SAR complex image data;

(1a) Input the main image data obtained by the InSAR main antenna of the interferometric synthetic aperture radar;

(1b) Input the auxiliary image data obtained by the interferometric synthetic aperture radar InSAR auxiliary antenna that has been fully registered with the main image;

(1c) Input the processing parameters and system parameters of the InSAR imaging processing of the interferometric synthetic aperture radar;

(2) Complex interference phase diagram:

(2a) Interferometric processing is performed on the input primary and secondary SAR images to obtain several interferometric phase diagrams;

(3) Obtain the optimal local window;

(4) Build a mixed integer optimization model;

(5) Calculate the fuzzy integer of the complex interference phase diagram;

(6) Calculate the unwrapped phase of each interferogram;

(7) Output the DEM of the entire scene.

Compared with the prior art, the present invention has the following advantages:

First, the present invention proposes a multi-baseline InSAR phase unwrapping algorithm based on a mixed integer optimization model. By combining the center and its adjacent pixels, under the premise that the pixels in the local window can be approximated to the small slope terrain, A mixed-integer optimization model is jointly constructed, which improves the robustness of images compared with the state-of-the-art. The method overcomes the influence of many phase noises in the interference phase image due to the method of using a single pixel in the prior art, reduces the difficulty of phase unwrapping of the filtered interference phase image, and can obtain a more real absolute phase value.

Second, the present invention utilizes an optimization strategy, which further reduces the influence of noise on the solution estimation and improves its application accuracy compared with the interference phase unwrapping method in the prior art interferometric synthetic aperture radar InSAR data processing.

Description of drawings (tables)

Fig. 1 is the flow chart of the present invention;

Figure 2 simulated terrain image;

Fig. 3 simulated interference image;

Fig. 4 Reconstructed topographic map of different methods;

Figure 5. The relationship between the standard deviation and the height root mean square error of different scenes.

Detailed ways

The present invention will be further described below with reference to the accompanying drawings.

With reference to accompanying drawing 1, the specific implementation steps of the present invention are as follows:

Step 1, input SAR image data and auxiliary parameters.

The main image data obtained by the main antenna of the interferometric synthetic aperture radar InSAR and the auxiliary image data obtained by the auxiliary InSAR antenna of the interferometric synthetic aperture radar that have been fully registered with the main image, and the assistants related to system parameters and imaging processing used in the processing process The parameters are input into the system, and the input primary and secondary SAR images must meet the quality requirements of interferometric processing in terms of coherence and imaging quality.

Step 2, complex interference phase map.

Perform interferometric processing on the input primary and secondary SAR images to obtain a complex interferometric phase map;

Step 3, obtain the optimal local window.

Taking the InSAR operation of two baselines as an example, an InSAR system with three-antenna phase centers (that is, two independent baselines) can provide two independent interferences with different baseline lengths, assuming that both interferences have undergone fine registration . The size of the local window (ie, the number of pixels in the azimuth and range directions) is a key parameter of the proposed mixed-integer optimization model. A large window size will result in an increased error between the actual terrain height and the assumed sloped terrain. On the other hand, if the local window is too small, fewer pixels will be included in the optimized model. This will affect the performance of resolving ambiguities. Therefore, there should be an optimized model with an optimal window size. The invention estimates the optimal window size according to the deviation of the interference packet phase and the assumed linear terrain. First estimate the size of the local optimum window of the complex interferometric phase map:

Among them, P represents the number of pixels in the local window azimuth, S ₁ and S ₂ represent the registered complex SAR image pair, auxiliary image data, (m, n) represents the coordinates of the center pixel, ( ) ^* represents Take the conjugate operation.

The number of pixels in the distance direction of the local window can be obtained by a similar method.

Step 4, build a mixed integer optimization model.

First, a mixed integer optimization model is constructed using the InSAR system parameters and the locally optimal window parameters obtained in step (3). In order to improve the robustness of noise, the local pixel method is adopted under the premise that the pixels of the local window satisfy the facet condition, that is, the terrain of the local window can be approximated as a slope. In fact, the DEM of the Earth's surface can be thought of as the entire surface consisting of many slopes, even in areas with complex topography. The optimal model is constructed by:

和

表示各个干涉图中通过精准SAR图像的共轭相乘法直接测量到的干涉相位，

和

表示有效的跨轨道的基线长度，

和

表示每个相位图的模糊数。where |ε ₁ | and |ε ₂ | represent the maximum distance between the true height and the inclined plane, |ε ₃ | represents the height difference between the individual interferograms, and a, b and c are the intrinsic parameters for constructing an inclined plane equation , X _i,j and R _i,j are the azimuth and distance directions relative to the center of the pixel, respectively, α ₁ and α ₂ are the scale factors of the absolute phase and height of each interferogram,

and

represents the interferometric phase directly measured by the conjugate multiplication of the precise SAR image in each interferogram,

and

represents the effective cross-track baseline length,

and

Represents the ambiguity number for each phase map.

Step 5: Calculate the fuzzy integer of the complex interference phase map.

For each local window, use the optimal model established in step (4) to solve the fuzzy integers of the two interferometric phase maps.

Step 6, calculate the unwrapped phase of each interferogram.

Step 7, output the DEM of the entire scene.

2. Simulation data processing experiment:

The simulation data experiment takes the distributed satellite as the platform, and the simulation parameters of the InSAR system are shown in the following table:

Fig. 2 is the simulated terrain, Fig. 3(a)(b) are the two ideal interferences measured by the long and short baselines, respectively, and Fig. 3(c)(d) are the long and , Short baseline interference. Compared with Fig. 3(a), the density of the interference phase fringes in Fig. 3(b) is much larger, because the longer the baseline length of the cross-track is, the smaller the blur height is. Also, due to discontinuous terrain, Figure 3(a) and (b) both contain some pixels whose adjacent phase difference exceeds π.

In order to facilitate the comparison with the simulated DEM, the recovered absolute phase values are converted into terrain heights, as shown in Fig. 4, where Fig. 4(a) is the terrain height map obtained by CRT. Figure 4(b) is the terrain height map obtained by the mixed-norm multi-baseline phase unwrapping algorithm based on linear programming. Figure 4(c) is a terrain height map obtained by the method of the present invention. It can be seen from FIG. 4 that the reconstructed terrain DEM obtained by the method of the present invention is more accurate.

The relationship between the phase standard deviation and the height root mean square error of the proposed method for different scenarios (mountain, slope, discontinuous step terrain) is shown in Figure 5, in which Figure 5(a) is the standard deviation of different methods in the mountain Plot versus height rms error. Figure 5(b) shows the relationship between the standard deviation of different slope methods and the root mean square error of height. Figure 5(c) shows the relationship between the standard deviation and the height root mean square error of different methods for discontinuous step terrain. Figure 5(d) shows the relationship between the standard deviation and the height root mean square error of different methods for the whole scene. It can be seen from Fig. 5 that the method proposed in the present invention has the ability to reduce the influence of noise and improves the robustness of baseline phase unwrapping.

It can be seen from the processing results in Fig. 4 and Fig. 5 that the method of the present invention can better realize the acquisition of terrain DEM.

Claims (3)

1. Multi-baseline InSAR phase unwrapping algorithm based on mixed integer optimization model, including the following steps: (1) Input SAR complex image data; (1a) Input the main image data obtained by the InSAR main antenna of the interferometric synthetic aperture radar; (1b) Input the auxiliary image data obtained by the interferometric synthetic aperture radar InSAR auxiliary antenna that has been fully registered with the main image; (1c) Input the processing parameters and system parameters of the InSAR imaging processing of the interferometric synthetic aperture radar; (2) Complex interference phase diagram: (2a) Interferometric processing of the input primary and secondary SAR images to obtain multiple interferometric phase images; (3) Obtain the local optimal window; The size of the local window is a key parameter of the mixed integer optimization model. According to the deviation of the interference packet phase and the assumed linear terrain, the size of the local optimal window is calculated; the size of the local window is the number of pixels in the azimuth and distance directions; Step ( 3) The described acquisition of the local optimal window is realized by the following formula:

Among them, P represents the number of pixels in the local window azimuth, S1 and S2 represent the registered complex SAR image pair, (m, n) represents the coordinates of the central pixel, and ( )* represents the conjugation operation; The number of pixels in the distance direction of the local window can be obtained by a similar method; (4) Build a mixed integer optimization model; Using the InSAR system parameters and the local optimal window parameters obtained in step (3), the terrain of the local window is approximated as a slope, and a mixed integer optimization model is constructed; (5) Calculate the fuzzy integer of the complex interference phase diagram; For each local window, use the mixed integer optimization model established in step (4) to solve the fuzzy parameters of the two interferometric phase diagrams; (6) Calculate the absolute phase of each interferogram; (7) Output the DEM of the entire scene. 2. the multi-baseline InSAR phase unwrapping algorithm based on mixed integer optimization model according to claim 1, is characterized in that: described in step (4), build mixed integer optimization model and utilize following formula to realize: minimising(|ε ₁ |+|ε ₂ |+|ε ₃ |) is subordinate to

where |ε ₁ | and |ε ₂ | represent the maximum distance between the true height and the inclined plane, |ε ₃ | represents the height difference between the individual interferograms, and a, b and c are the intrinsic parameters for constructing an inclined plane equation , X _i,j and Ri _,j are the azimuth and distance directions relative to the pixel center, respectively, α ₁ and α ₂ are the absolute phase and height scaling factors of each interferogram,

and

represents the interferometric phase directly measured by the conjugate multiplication of the precise SAR image in each interferogram,

and

represents the effective cross-track baseline length,

and

Represents the ambiguity number for each phase map. 3. the multi-baseline InSAR phase unwrapping algorithm based on mixed integer optimization model according to claim 1, is characterized in that: the fuzzy integer of calculating complex interferometric phase diagram described in step (5): Calculated by using the optimal window determined in claim 1 and the mixed integer optimization model described in claim 2.

Images