InSAR Interferometric Fringe Matching Method Based on Wavelet Transform

Abstract

This invention discloses an InSAR interferometric fringe matching method based on wavelet transform, which is used for InSAR matching navigation. It matches the InSAR fringes acquired in real time by the aircraft with the InSAR fringes obtained by inversion from the onboard digital elevation model to obtain relevant state parameters such as range and azimuth positioning offsets. Interferometric synthetic aperture radar image fringe matching is used to correct the systematic errors of the InSAR system, achieving autonomous, long-term, and accurate positioning navigation.

Description

Technical Field

[0001] This invention relates to the field of navigation technology, and specifically to an InSAR interferometric fringe matching method based on wavelet transform. Background Technology

[0002] Inertial Navigation Systems (INS), as autonomous navigation systems independent of the external environment, are widely used in military and civilian fields. However, INS drift errors accumulate over time. To eliminate INS positioning errors and achieve long-term autonomous navigation, navigation systems often use other navigation methods to assist INS. Currently, the main autonomous integrated navigation methods in practical applications are terrain matching navigation and Synthetic Aperture Radar (SAR) scene matching navigation.

[0003] Terrain-matching navigation is a process where, during flight, an aircraft uses a radar altimeter to measure the elevation of terrain profiles along its flight path. The real-time elevation data is then compared with pre-stored baseline terrain profile data, and the aircraft's geographic location is determined based on the optimal match. This precise geographic location information is then fused with information from the airborne navigation system (INS) to eliminate INS errors and achieve accurate navigation.

[0004] SAR scene matching navigation leverages the all-weather, all-day observation capabilities and high-resolution two-dimensional capabilities of SAR, providing high-resolution images similar to those from an optical camera even in low-visibility environments. Real-time SAR image information is matched with map data from the corresponding mapping zone in an onboard digital map database to obtain precise geographic location information at the current moment. This allows for the calculation of the INS positioning error, which is then used as the observation measurement and passed through a Kalman filter to estimate the INS error. This error is then used to correct the INS error, resulting in long-term accurate geographic location information. Simultaneously, this information can be used for motion compensation and view-of-sight positioning parameter calculations on the SAR, thereby determining the aircraft's state information and improving navigation accuracy.

[0005] Existing terrain matching navigation mainly involves searching for terrain profiles measured by radar altimeters on a reference terrain map stored in the aircraft. This method is prone to mismatches in flat areas, so it is necessary to select areas with significant terrain slope features for matching. Moreover, terrain matching navigation is a one-dimensional matching method without two-dimensional resolution, so the matching accuracy is not high.

[0006] SAR image matching, with its all-weather, all-day operation and high-resolution two-dimensional capabilities, has been successfully applied in the military field. However, it also has some drawbacks. For example, SAR image matching is susceptible to seasonal changes in ground features, requiring the selection of typical terrain features to ensure matching stability. Furthermore, single-channel SAR cannot detect changes in aircraft attitude, thus hindering effective motion compensation, affecting imaging accuracy, and consequently reducing matching accuracy. Summary of the Invention

[0007] To address the aforementioned technical problems, this invention provides an InSAR interferometric fringe matching method based on wavelet transform. This method uses interferometric synthetic aperture radar (InSAR) image fringe matching to correct the systematic errors of INS, thereby achieving autonomous, long-term, and accurate positioning and navigation.

[0008] To achieve the above objectives, the technical solution adopted by the present invention is as follows:

[0009] An InSAR interferometric fringe matching method based on wavelet transform includes the following steps:

[0010] Step 1, obtaining interference fringes, specifically includes:

[0011] (1) First, the registration of single-view complex image s1 and single-view complex image s2 is performed based on the actual InSAR echo data, including the calculation of the range and azimuth offset values ​​between single-view complex image s1 and single-view complex image s2 and the resampling of the sub-image to the main image.

[0012] (2) Interference processing is performed on the two registered images, that is, one image is multiplied by the conjugate of the other image, i.e., s1.s * 2. Then, extract the phase information;

[0013] (3) Remove the flat phase to obtain the interference fringes to be matched;

[0014] Step 2, separation of interference fringe features, specifically includes:

[0015] (1) Decomposition and reconstruction of interference fringes;

[0016] (2) Perform interference fringe texture filtering;

[0017] (3) Extract the stripe profile line features;

[0018] (4) Rapidly extract the interference fringe features in the shaded area;

[0019] (5) Rapidly extract interference fringe features from overlapping regions;

[0020] Step 3, multi-scale correlation matching, specifically includes:

[0021] Let I(x,y) be a reference interference fringe image or feature at a certain scale. The actual acquired interference fringes and extracted features are used as a template image J(x,y). The template image J(x,y) is overlaid on the reference interference fringe image or feature I(x,y) at a certain scale for a translational sliding window search. A normalized cross-correlation matching operator is used to measure the similarity between the reference interference fringes and the actual interference fringes. The normalized correlation coefficient is:

[0022]

[0023] Where, μ I and μ J These are the means of the two images, σ and σ', respectively. I and σ J These are the standard deviations of the two images, and N represents the number of pixels in the image. The normalized cross-correlation range is [0,1].

[0024] Further, step (1) of step 2 includes:

[0025] The interference fringes are decomposed at the j-1 scale using the two-dimensional Mallat wavelet decomposition algorithm:

[0026] C j =H az H rg C j-1

[0027]

[0028]

[0029]

[0030] Among them, H az For azimuth low-pass decomposition filters, H rg For a range-directed low-pass decomposition filter, G az For azimuth high-pass decomposition filters, G rg For range-oriented high-pass decomposition filters; Corresponding to stripe image C j-1 The low-frequency components, high-frequency components in the vertical direction, high-frequency components in the horizontal direction, and high-frequency components in the diagonal direction;

[0031] The interference fringes at the original resolution were reconstructed using the Mallat wavelet reconstruction algorithm.

[0032]

[0033] in, These are the weighting coefficients of the four components in the j-th layer. When all these coefficients are equal to 1, the original interference fringe image is completely recovered; H az * H rg * G az * and G rg * H respectively az H rg G az and G rg . conjugate.

[0034] Furthermore, step (2) of step 2 specifically includes:

[0035] The coherence distance between the center pixel and its neighboring pixels is defined as:

[0036] D=|γ i -γ j |

[0037] Where, γ i and γ j These are the coherence coefficients of the center pixel i and its neighboring pixel j, respectively.

[0038] Based on the relationship between the correlation coefficient distance and the neighborhood energy value, the improved neighborhood energy form is as follows:

[0039]

[0040] Where, x i Indicates pixel segmentation label, x j It is x i The nearest neighbor segmentation label; α is used to control the curve shape, and β is a parameter greater than zero, used to adjust the weights of likelihood energy and neighborhood energy.

[0041] Furthermore, step (3) of step 2 specifically includes:

[0042] Mallat wavelet decomposition and reconstruction preserve the high-frequency component D in the vertical direction of the azimuth distance. j 1 and the high-frequency component D in the horizontal direction j 2 The line edge structure can be completely reconstructed using two orthogonal components, horizontal and vertical, while also preserving low-frequency components for the extraction of equiphase lines and maximum gradient profiles. The line feature wavelet reconstruction equation is expressed as:

[0043]

[0044] After wavelet decomposition and reconstruction, high-frequency noise information and high-frequency change edges in the interference phase map are removed to obtain a smoothed interference phase map. Then, based on the definition of equiphase lines and maximum gradient profiles, equiphase lines and maximum gradient profiles are extracted from both the actual interference fringes and the reference simulated interference fringes, and feature matching is performed to obtain the matching results.

[0045] Furthermore, step (4) of step 2 specifically includes:

[0046] Within the framework of wavelet multi-scale analysis, Mallat wavelet decomposition and reconstruction should, as far as possible, preserve high-frequency components in the diagonal direction. The shadow region is quickly detected using the diagonal direction component. The wavelet reconstruction equation for the shadow region is:

[0047]

[0048] After wavelet domain decomposition and reconstruction, threshold segmentation is further performed by combining the intensity information on the image grayscale. Interference phase noise is removed by using graphics erosion and dilation methods, and the shadow area detected by small area shadow detection is obtained to obtain the detection result of larger area shadow area. The edge detection algorithm is used to extract the edge of the shadow area, and the boundary information of this shadow area is matched with the shadow edge vector stored in the reference image to obtain the matching result of the shadow area.

[0049] Furthermore, step (5) of step 2 specifically includes:

[0050] Thresholding segmentation is performed by combining the intensity information on the grayscale of the image. The overlapping false alarm regions are removed by using the methods of erosion and dilation in graphics. The interference fringes of the overlapping regions have certain texture features. The interference fringes of the overlapping regions are used as surface features to be directly matched with the overlapping regions of the baseline simulated interference fringes to obtain the matching results of the overlapping regions. At the same time, the boundary of the overlapping regions is mapped to the geographical location of the pixels in the baseline library. Therefore, the boundary of the overlapping regions is used as the matching result of the overlapping regions.

[0051] Furthermore, step 3 specifically includes:

[0052] Step (1): Perform j-level wavelet decomposition on the actual interference fringes and the reference simulated interference fringes respectively to obtain the interference phase diagram after each level of decomposition;

[0053] Step (2): Based on the typical features of texture features, line features, shadow interferograms, and overlay interferograms of the interferometric phase map, under the wavelet multi-scale decomposition architecture, the interferometric fringe map is processed according to the wavelet decomposition and reconstruction method in step (1) of step 2. According to the definition of fringe features, the interferometric phase surface features, equiphase line features, maximum gradient profile line features, shadow lines and surface features, and overlay surface features in the interferometric fringes are extracted.

[0054] Step (3): On layer j, the cross-correlation similarity measurement template matching method is used to perform coarse matching of the features extracted from the two images at the lowest level of the scale space to obtain the best matching region at that scale.

[0055] Step (4): At layer j-1, in the best matching region of the previous step, perform high-resolution cross-correlation matching calculation at the j-1 scale to filter out mismatched features.

[0056] Step (5) Repeat steps (3) and (4) until the recursion reaches the top level of the scale space to obtain the final matching result.

[0057] Beneficial effects:

[0058] (1) InSAR fringe matching and traditional SAR scene matching have significant differences in image features and observable motion parameters. SAR scene matching is easily affected by changes in ground features, and single-channel SAR cannot sense the attitude information of the aircraft. In contrast, the interferometric fringes after InSAR is removed from flat ground are sensitive to long-term stable terrain slope features, and InSAR has a stronger position and attitude inversion capability.

[0059] (2) Existing scene matching navigation only uses the position of the matching point to correct inertial navigation errors. InSAR fringe matching, compared with SAR image matching, not only obtains position information through interferometric fringe matching, but also can invert platform attitude parameters with high accuracy. InSAR interferometric phase contains roll angle information, Doppler frequency contains yaw and pitch angle attitude information, and point target phase history contains relative motion information. These motion parameters include absolute three-dimensional position, radar slant range error, time-varying baseline error, and Doppler frequency error information. By fusing this information with INS, the platform's attitude information can be inverted with high accuracy, improving positioning and navigation accuracy. Attached Figure Description

[0060] Figure 1 This is a flowchart of an InSAR interferometric fringe matching method based on wavelet transform according to the present invention.

[0061] Figure 2 This is a schematic diagram of the interference fringe texture information;

[0062] Figure 3 These are the phase lines of the interference fringes;

[0063] Figure 4 The maximum gradient profile of the interference fringes;

[0064] Figure 5 This is a schematic diagram illustrating the influence of the shaded area and the characteristics of the interference fringes.

[0065] Figure 6 This is a schematic diagram of the overlapping area image and interference fringe characteristics. Detailed Implementation

[0066] To make the objectives, technical solutions, and advantages of this invention clearer, the invention will be further described in detail below with reference to the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are merely illustrative and not intended to limit the invention. Furthermore, the technical features involved in the various embodiments of this invention described below can be combined with each other as long as they do not conflict with each other.

[0067] This invention discloses an InSAR interferometric fringe matching method based on wavelet transform, which is an interferometric fringe matching method used in InSAR matching navigation. It matches the InSAR fringes acquired by the aircraft in real time with the InSAR fringes obtained by inversion from the digital elevation model stored in the aircraft to obtain relevant state parameters such as range and azimuth positioning offset.

[0068] like Figure 1 As shown, the InSAR interferometric fringe matching method based on wavelet transform of the present invention specifically includes the following steps:

[0069] Step 1, obtaining interference fringes, specifically includes:

[0070] (1) First, the registration of single-view complex image s1 and single-view complex image s2 is performed based on the actual InSAR echo data, including the calculation of the range and azimuth offset values ​​between single-view complex image s1 and single-view complex image s2 and the resampling of the sub-image to the main image.

[0071] (2) Interference processing is performed on the two registered images, that is, one image is multiplied by the conjugate of the other image, i.e., s1.s * 2. Then, extract the phase information;

[0072] (3) Remove the flat phase to obtain the interference fringes to be matched.

[0073] Step 2, interference fringe feature separation, specifically includes:

[0074] (1) Decomposition and reconstruction of interference fringes:

[0075] The interference fringes are decomposed at the j-1 scale using the two-dimensional Mallat wavelet decomposition algorithm:

[0076] C j =H az H rg C j-1

[0077]

[0078]

[0079]

[0080] Among them, H az For azimuth low-pass decomposition filters, H rg For a range-directed low-pass decomposition filter, G az For azimuth high-pass decomposition filters, G rg For range-oriented high-pass decomposition filters; Corresponding to stripe image C j-1 The interference fringes are composed of low-frequency components, high-frequency components in the vertical direction, high-frequency components in the horizontal direction, and high-frequency components in the diagonal direction. The original resolution interference fringes can be reconstructed using the Mallat wavelet reconstruction algorithm.

[0081]

[0082] in, These are the weighting coefficients of the four components in the j-th layer. When all these coefficients are equal to 1, the original interference fringe image is completely recovered; H az * H rg * G az * and G rg * H respectively az H rg G az and G rg . conjugate.

[0083] Wavelet multiresolution analysis can decompose signals at different scales and break down mixed signals composed of various frequencies into sub-signals of different frequency bands. The specific steps for denoising using wavelet decomposition and reconstruction are: decomposing the noisy signal into different frequency bands at a certain scale, then setting the frequency band containing the noise to zero, and performing wavelet reconstruction to achieve the separation of different features.

[0084] The slowly varying terrain information in the interferometric phase will be preserved in the low-frequency component C. j In the mid-mountain and alpine striped winding edge information, the main focus is on the horizontal and vertical high-frequency components. In the shadowed region, the interference fringes exhibit high-frequency noise characteristics, primarily consisting of high-frequency components in the diagonal direction. In this study, by controlling this weighting coefficient, low-frequency interference phase, linear edge features, high-frequency random noise in the shadowed area, and edge features of the overlapping region in the original interference fringes can be separated.

[0085] (2) Perform interference fringe texture filtering, the specific method is as follows:

[0086] In mountainous areas with large topographic relief, such as high mountains and the Hengduan Mountains, even after the flattening process, the interference fringes still exhibit a 2π-wound phase, and the fringes contain a large number of fringe texture features, such as... Figure 2 As shown, these texture features are actually periodic interleaving features of the interference phase, with the edge interference phase values ​​between the stripes being 0 or 2π.

[0087] In mountainous areas with significant topographic relief, interference fringe entanglement is a frequent phenomenon, with the entanglement characteristics exhibiting a [0 2π] periodic variation from the mountaintop to the valley floor. Furthermore, compared to the baseline simulated interference fringes, the actual interference fringe map contains a large amount of ground feature texture information, which needs to be removed. This invention, within the framework of wavelet domain decomposition, first denoises the interference phase, proposing a Bayesian method for denoising using a wavelet domain coherent-Markov model. Utilizing the strong correlation between the true values ​​of the interference phase in the interference fringes, and the similarly strong correlation between wavelet coefficients after wavelet transform, noise filtering of the interference fringes and preservation of texture information are achieved by processing the local interaction relationships between wavelet coefficients.

[0088] In traditional Markov random field models, the neighborhood energy is a fixed value when the center pixel and its neighboring pixels have the same label, and zero when the labels are different. This results in adjacent pixels having the same influence on the center pixel, thus failing to fully utilize the contextual information. Based on the interference coherence coefficient, a new neighborhood energy model for Markov is defined, and a detection method based on the coherence coefficient-Markov random field model is proposed to fully utilize the interference information of the context.

[0089] The coherence distance between the center pixel and its neighboring pixels is defined as:

[0090] D=|γ i -γ j |

[0091] Where, γ i and γ j are the coherence coefficients of the center pixel i and the neighboring pixel j, respectively.

[0092] When the center pixel and its neighboring pixels have the same segmentation label, the closer their coherence coefficients are, the more similar they are. Therefore, the corresponding neighborhood energy value should be lower, and the probability of them having the same label is greater. When the center pixel and its neighboring pixels have different segmentation labels, the closer their coherence coefficients are, the more similar they are. Therefore, the corresponding neighborhood energy value should be higher, and the probability of them having different labels is lower. Based on the above analysis, the improved neighborhood energy form is as follows:

[0093]

[0094] Where, x i Indicates pixel segmentation label, x j It is x i The nearest neighbor segmentation label.

[0095] Here, α controls the curve shape, and β is a parameter greater than zero used to adjust the weights of the likelihood energy and the neighborhood energy. In the wavelet domain, a coherent Markov prior model is constructed to reflect the gradually varying edge characteristics of the interference fringes as a criterion for wavelet coefficient reconstruction, thereby reducing interference phase noise while preserving the fringe texture at large scales.

[0096] (3) Extract the stripe profile line features. The specific method is as follows:

[0097] Interference fringe texture feature extraction is a surface matching method, but it is computationally intensive. Therefore, it is necessary to use simpler line feature representations to reduce the amount of image data to be matched and improve matching efficiency. Based on the variation characteristics of interference fringes, this invention proposes an interference fringe matching method based on two line features: equiphase lines and maximum gradient profiles.

[0098] Within the wavelet multi-scale analysis framework, line features are a type of edge feature. Mallat wavelet decomposition and reconstruction are used to preserve the high-frequency component D in the vertical direction of the azimuth distance. j 1 and the high-frequency component D in the horizontal direction j 2 The line edge structure can be completely reconstructed using two orthogonal components, horizontal and vertical, while also preserving low-frequency components for the extraction of equiphase lines and maximum gradient profiles. The line feature wavelet reconstruction equation can be expressed as:

[0099]

[0100] Equiphase lines are defined as the curves connecting adjacent points with equal interference phase values ​​in an interferometric fringe pattern. Within a 2π interference phase period, equiphase lines represent contour lines, typically appearing as closed curves in the interferometric fringe pattern, as shown in Figure 3. Equiphase lines also represent terrain slope; denser equiphase lines indicate a steeper slope, while a greater horizontal distance and sparser arrangement indicate a gentler slope. In particular, equiphase lines change direction when passing through ridges or valleys, revealing significant variations in terrain slope.

[0101] The maximum gradient profile line is defined as the curve connecting points in the azimuth and distance directions of the interference fringe gradient pattern, such as... Figure 4As shown, the maximum gradient profile represents the region with the greatest change in terrain slope between two orthogonal directions. It is the steepest part of the fringe pattern and further simplifies the characteristics of the isophase lines compared to the isophase lines, thus further refining the most significant features of the interference fringe pattern.

[0102] After wavelet decomposition and reconstruction, high-frequency noise and high-frequency variation edges in the interferometric phase map are removed, resulting in a smoothed interferometric phase map. Then, based on the definitions of equiphase lines and maximum gradient profiles, equiphase lines and maximum gradient profiles are extracted from both the actual interference fringes and the reference simulated interference fringes, and feature matching is performed to obtain the matching results.

[0103] (4) Perform rapid extraction of interference fringe features in the shaded area. The specific method is as follows:

[0104] Because shadows obstruct the radar's line of sight, the back echo cannot reach the receiver. Therefore, the interference phase within the shadowed area exhibits entirely as random noise, such as... Figure 5 As shown. Due to the geometry of InSAR side-view imaging, the edge features of shadows are closely related to the radar incident angle; the same mountain will exhibit different shadow features when observed from different angles. Interference fringe matching in the shadow area is achieved by combining random phase features and shadow edge features, increasing the number of matching points. Conventional shadow region detection methods typically use image brightness and statistical coherence coefficients to determine shadows. In fact, in the interferometric phase map, the shadow region appears as chaotic high-frequency white noise. Within the wavelet multi-scale analysis framework, Mallat wavelet decomposition and reconstruction should, as far as possible, preserve the high-frequency components in the diagonal direction. The diagonal component can be used to quickly detect shadow regions. The wavelet reconstruction equation for the shadow region is as follows:

[0105]

[0106] After wavelet domain decomposition and reconstruction, threshold segmentation is further performed by combining the intensity information in the image grayscale. Interference phase noise and shadow regions detected in smaller shadow areas are removed using graphical erosion and dilation methods, resulting in larger shadow region detection results. Shadow regions exhibit random noise, and point-to-point matching is not possible within these regions. However, the boundaries of the shadow regions can be mapped to the geographical locations of pixels in the benchmark database. Therefore, an edge detection algorithm is used to extract the edges of the shadow regions. This shadow boundary information is then matched with the shadow edge vectors stored in the benchmark image to obtain the matching results for the shadow regions. Figure 5 As shown.

[0107] (5) Perform rapid extraction of interference fringe features in the overlapping region. The specific method is as follows:

[0108] The mechanism of superposition is that radar echoes from the mountaintop and the base of the mountain arrive at the radar receiver simultaneously, superimposing on the same pixel unit in the image. This causes a phase jump and reverse bias in the interference, which appears as a bright area in the SAR image, such as... Figure 6 As shown. Since the coherence of the overlapping region is not low, it is impossible to use coherence to determine the overlapping region completely. Therefore, the commonly used method is still the image brightness threshold segmentation method combined with a semi-automatic detection method with manual intervention.

[0109] The phase of the overlapping region in the interferometric fringe pattern is reversed. Within the overall framework of wavelet multi-scale analysis, the overlapping region is automatically detected using the extraction results of the aforementioned isophase lines and maximum gradient profile. The overlapping region is the area with the steepest terrain slope and will intersect with the maximum gradient profile; at the same time, the isophase lines in the overlapping region will undergo jumps and truncation.

[0110] In practical feature extraction, threshold segmentation is performed by combining intensity information in the image grayscale. Graphical erosion and dilation methods are used to remove overlapping false alarm regions. The interference fringes in the overlapping regions have certain texture features. These interference fringes are used as surface features for direct matching with the overlapping regions of the baseline simulated interference fringes, yielding the matching results for the overlapping regions. Simultaneously, the boundaries of the overlapping regions can be mapped to the geographical locations of pixels in the baseline database; therefore, the boundaries of the overlapping regions are used as the matching results for the overlapping regions.

[0111] Step 3, multi-scale correlation matching, specifically includes:

[0112] Considering the real-time requirements of navigation applications and to improve matching efficiency, this invention proposes a wavelet multi-scale interferometric fringe correlation matching method. After wavelet decomposition, a target feature vector in scale space is established for the image. To reduce the number of matching searches, matching starts from low resolution, where the low-resolution features represent the overall texture features. Based on this matching, high-resolution feature matching is further performed. The search strategy utilizes low-resolution features for layer-by-layer filtering, thereby improving the efficiency of fringe matching search.

[0113] Let I(x,y) be a reference interference fringe image or feature at a certain scale. The actual acquired interference fringes and extracted features are used as a template image J(x,y). The template image J(x,y) is overlaid on the reference interference fringe image or feature I(x,y) at the certain scale for a translational sliding window search. Since the interference fringe features include line, surface, and interference phase values, a normalized cross-correlation matching operator is used to measure the similarity between the reference interference fringes and the actual interference fringes. The normalized correlation coefficient is as follows:

[0114]

[0115] Where, μ I and μJ These are the means of the two images, σ and σ', respectively. I and σ J These are the standard deviations of the two images, and N represents the number of pixels in the image. The normalized cross-correlation range is [0,1].

[0116] Specifically, the steps of the wavelet multi-scale interferometric fringe matching algorithm are as follows:

[0117] Step (1) Perform j-level wavelet decomposition on the actual interference fringes and the reference simulated interference fringes respectively to obtain the interference phase diagram after each level of decomposition.

[0118] Step (2) For the typical features of the interferometric phase map, such as texture features, line features, shadow interferogram, and overlapping interferogram, the interferometric fringe map is processed under the wavelet scale decomposition framework and the wavelet decomposition reconstruction method proposed above. According to the definition of fringe features, the interferometric phase surface features, equiphase line features, maximum gradient profile line features, shadow lines and surface features, and overlapping surface features are extracted from the interferometric fringes.

[0119] In step (3), on layer j, the cross-correlation similarity measurement template matching method is used to perform coarse matching of the features extracted from the two images at the lowest level of the scale space to obtain the best matching region at that scale.

[0120] Step (4) At layer j-1, in the best matching region of the previous step, perform high-resolution cross-correlation matching calculation at the j-1 scale to filter out mismatched features.

[0121] Step (5) Repeat steps (3) and (4) until the recursion reaches the top level of the scale space to obtain the final matching result.

[0122] Those skilled in the art will readily understand that the above description is merely a preferred embodiment of the present invention and is not intended to limit the present invention. Any modifications, equivalent substitutions, and improvements made within the spirit and principles of the present invention should be included within the scope of protection of the present invention.

Claims

1. A wavelet transform based InSAR interferogram fringe matching method, characterized in that, Specifically comprising the following steps: Step 1, the acquisition of interference fringes, specifically comprising: (1) first, according to the actual InSAR echo data, the registration of single view complex image s1 and single view complex image s2, including the calculation of the distance and azimuth offset values between single view complex image s1 and single view complex image s2 and the resampling of the sub-image to the main image; (2) The two registered images are interfered, i.e. one image is multiplied by the conjugate of the other, i.e. s1.s * 2, after which the phase information is taken out; (3) remove the flat phase, get the interference fringes to be matched; Step 2, the separation of interference fringe features, specifically comprising: (1) decomposition and reconstruction of interference fringes; (2) interference fringe texture filtering is carried out; (3) the feature extraction of fringe profile line is carried out; (4) the fast extraction of interference fringe features in shadow area is carried out; (5) the fast extraction of interference fringe features in overlap area is carried out; Step 3, multi-scale correlation matching, specifically comprising: Let the reference interference fringe image or feature I(x,y) under a certain scale, the actual acquired interference fringe and the extracted feature as the template graph J(x,y), and the template graph J(x,y) is overlaid on the reference interference fringe image or feature I(x,y) under a certain scale to carry out translation sliding window search; the normalized cross correlation matching operator is used to realize the similarity measurement between the reference interference fringe and the actual interference fringe; the normalized correlation coefficient is: where μ I and μ J are the mean values of the two images, σ I and σ J are the standard deviations of the two images, and N represents the number of pixels in the image. The normalized cross-correlation has a value range of [0, 1].

2. The wavelet transform based InSAR interferogram fringe matching method according to claim 1, characterized in that: The step 2 (1) comprises: The interference fringe is decomposed into: C j =H az H rg C j-1 where H az is an azimuth low-pass decomposition filter, H rg is a range low-pass decomposition filter, G az is an azimuth high-pass decomposition filter, G rg is a range high-pass decomposition filter; C j , correspond to the low frequency component, the high frequency component in the vertical direction, the high frequency component in the horizontal direction and the high frequency component in the diagonal direction of the fringe image C j-1 , respectively. The original resolution interference fringe is reconstructed by the Mallat wavelet reconstruction algorithm: wherein, are the weighting coefficients of the four components of the jth layer, respectively, and when these coefficients are all equal to 1, the original interference fringe pattern is completely recovered; H az * ,H rg * ,G az * and G rg * are the conjugates of H az ,H rg ,G az and G rg , respectively.

3. The wavelet transform based InSAR interferogram matching method according to claim 2, characterized in that: The step 2 (2) specifically comprises: The coherent coefficient distance of the center pixel and the neighborhood pixel is defined as: D = | γ i - γ j | where γ i and γ j are the coherence factors of the central pixel i and the neighboring pixel j, respectively. Based on the relationship between the correlation coefficient distance and the domain energy value, the improved neighborhood energy form is as follows: where x i represents the pixel segmentation label, x j is the neighboring segmentation label of x i ; a is used to control the curve shape, and β is a parameter greater than zero used to adjust the weight of the likelihood energy and the neighborhood energy.

4. The wavelet transform based InSAR interferogram matching method according to claim 3, characterized in that: The step 2 (3) specifically comprises: The high frequency components D in the vertical direction and the high frequency components D in the horizontal direction are reserved by using Mallat wavelet decomposition and reconstruction j 1 j 2 The line edge structure can be completely reconstructed by using the horizontal and vertical orthogonal components, and the low frequency components are also reserved for the extraction of the equal phase line and the maximum gradient profile. The line feature wavelet reconstruction equation is represented as:​ After wavelet decomposition and reconstruction, the high-frequency noise information and the high-frequency change edge in the interference phase graph are removed to obtain the smoothed interference phase graph; then, according to the equal phase line and the maximum gradient profile definition, the equal phase line and the maximum gradient profile line are further extracted in the actual interference fringe and the reference simulation interference fringe, and the feature matching is carried out to obtain the matching result.

5. The wavelet transform based InSAR interferogram matching method according to claim 4, characterized in that: The step 2 (4) specifically comprises: In the framework of wavelet multi-scale analysis, Mallat wavelet decomposition and reconstruction should be used to retain the high frequency components in the diagonal direction as much as possible The shadow region is quickly detected by using the diagonal direction component, and the wavelet reconstruction equation of the shadow region is: After wavelet domain decomposition and reconstruction, the threshold segmentation processing is further carried out combined with the intensity information on the image gray, the shadow area is removed by using the methods of graph erosion and expansion, the interference phase noise and the small area shadow detection in the shadow area are removed, and the detection result of the large area shadow area is obtained; the edge detection algorithm is used to extract the edge of the shadow area, the boundary information of the shadow area is matched with the stored shadow edge vector in the reference graph, and the matching result of the shadow area is obtained.

6. The wavelet transform based InSAR interferogram matching method according to claim 5, characterized in that: The step 2 (5) specifically comprises: The intensity information on the image gray scale is used to make threshold segmentation processing, the false alarm area of the overlay is removed by using the method of erosion and expansion of graph theory, the interference fringes of the overlay area have certain texture characteristics, the interference fringes of the overlay area are directly matched with the reference simulation interference fringes of the overlay area as the surface features, and the matching result of the overlay area is obtained; meanwhile, the boundary of the overlay area is mapped with the geographical position of the pixel points in the reference library, so the boundary of the overlay area is taken as the matching result of the overlay area.

7. The wavelet transform based InSAR interferogram matching method according to claim 6, characterized in that: The step 3 specifically comprises: Step (1), j-layer wavelet decomposition is performed on the actual interference fringes and the reference simulation interference fringes respectively, and the interference phase patterns after each level of decomposition are obtained; Step (2), according to the texture characteristics, line characteristics, shadow interference pattern and typical characteristics of the overlay interference pattern of the interference phase pattern, under the wavelet multi-scale decomposition architecture, the interference fringe pattern is processed according to the wavelet decomposition and reconstruction method in the first step in the step 2, and according to the definition of the fringe characteristics, the interference phase surface features, the equal phase line features, the maximum gradient profile line features, the shadow line and surface features and the overlay surface features in the interference fringes are extracted; Step (3), on the j-layer, the cross-correlation similarity measurement template matching method is used to perform coarse matching on the features extracted from the two images at the bottom layer of the scale space, and the best matching area on the scale is obtained; Step (4), on the j-1 layer, the j-1 scale high-resolution cross-correlation matching calculation is further performed on the best matching area in the last step, and the false matching features are screened; Step (5), steps (3) and (4) are repeated until the top layer of the scale space is reached, and the final matching result is obtained.

Images