CN110109100A - A kind of more baseline least square phase unwrapping methods based on Quality Map weighting - Google Patents
A kind of more baseline least square phase unwrapping methods based on Quality Map weighting Download PDFInfo
- Publication number
- CN110109100A CN110109100A CN201910271422.5A CN201910271422A CN110109100A CN 110109100 A CN110109100 A CN 110109100A CN 201910271422 A CN201910271422 A CN 201910271422A CN 110109100 A CN110109100 A CN 110109100A
- Authority
- CN
- China
- Prior art keywords
- phase
- winding
- orientation
- gradient
- distance
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
- G01S13/00—Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
- G01S13/88—Radar or analogous systems specially adapted for specific applications
- G01S13/89—Radar or analogous systems specially adapted for specific applications for mapping or imaging
- G01S13/90—Radar or analogous systems specially adapted for specific applications for mapping or imaging using synthetic aperture techniques, e.g. synthetic aperture radar [SAR] techniques
- G01S13/9021—SAR image post-processing techniques
- G01S13/9023—SAR image post-processing techniques combined with interferometric techniques
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
- G01S7/00—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
- G01S7/02—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00
- G01S7/41—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00 using analysis of echo signal for target characterisation; Target signature; Target cross-section
Landscapes
- Engineering & Computer Science (AREA)
- Remote Sensing (AREA)
- Radar, Positioning & Navigation (AREA)
- Physics & Mathematics (AREA)
- Computer Networks & Wireless Communication (AREA)
- General Physics & Mathematics (AREA)
- Electromagnetism (AREA)
- Magnetic Resonance Imaging Apparatus (AREA)
- Image Analysis (AREA)
Abstract
The invention discloses a kind of more baseline least square phase unwrapping methods based on Quality Map weighting, this method calculates the single order phase gradient main value of more baseline InSAR winding phases first, and by its mirror-symmetric extension, calculate second order phase gradient, then Quality Map judgement winding phase masses are utilized, phase compensation is carried out to the residue points in winding phase as prior information, and weighting is optimized based on second order phase gradient of the Quality Map to winding phase, finally carry out more baseline least square phase unwrappings, it obtains solution and twines phase, this method enriches the interferometric phase information of target scene by more base line measurements, avoid the residue points phase error problem of transmission in traditional single baseline least square phase unwrapping algorithm, improve phase unwrapping precision.
Description
Technical field
The invention belongs to Radar Technology field, its in particular to interference synthetic aperture radar (Interferometric
Synthetic Aperture Radar, InSAR) phase unwrapping technical field in hypsographic survey.
Background technique
Interference synthetic aperture radar (InSAR) technology is developed in the early 1970s, it is by radio interferometry skill
Art and SAR imaging technique combine, and can carry out round-the-clock, round-the-clock interferometry to large-scale terrain, obtain two dimension
While SAR image, the three-dimensional information of earth's surface can be extracted using the phase information of SAR image.The working principle of InSAR is logical
It crosses two slave antennas while observation or single antenna difference was navigated observation, obtain the SAR image pair of ground Same Scene, extract SAR figure
As interferometric phase information, digital elevation model (Digital Elevation is finally inversed by conjunction with the geometrical relationship of image scene
Model, DEM), to realize high-precision three-dimensional mapping, be widely used in Natural calamity monitoring, mapping and provide naturally
The fields such as source investigation.More baseline InSAR (multi-baseline InSAR) technologies are the extensions of single baseline InSAR technology,
The contradiction that can overcome single baseline InSAR system altimetry precision and phase unwrapping difficulty that can not take into account to a certain extent.Using more
Baseline InSAR measuring technique fusion multi-frame interferometry image can identify the biggish region of actual landform mesorelief well, thus grind
It is increasingly important in high-precision three-dimensional acquisition of information field to study carefully more baseline InSAR Systems Theorys and technology.It is detailed in document " airborne Interference
Data of synthetic aperture radar processing technique ", red violent wind of fourth etc. are write, Science Press.
It is corresponding with landform altitude true that phase unwrapping (Phase Unwrapping, PU), which refers to winding phase recovery,
The process of reality position is step most important in InSAR data process flow.Two width SAR image conjugate multiplications are interfered
Phase be wrapped in [- π, π) in section, the integral multiple of 2 π is differed with its true phase, phase unwrapping is to obtain digital elevation letter
The main source of error of breath is directly related to the precision of elevation information extraction.Minimum two based on Fast Fourier Transform (FFT) (FFT)
Multiplying phase unwrapping algorithm is a kind of algorithm globally optimal, can obtain the optimal solution under least square meaning, is current practice
In common phase unwrapping method, basic thought is the phase of the phase gradient of the neighbor pixel of the true phase of original scene
Position main value is identical as the phase main value of phase gradient of the winding phase neighbor pixel obtained from echo data.It is detailed in document
“Zhang Z,Guo Z.Study on the unweighted least-squares phase unwrapping
algorithm[J].Proceedings of SPIE-The International Society for Optical
Engineering,2010,7848(1):30-36.”。
Fold existing for single baseline InSAR system cover, phase unwrapping essence caused by the effects such as shade and rolling topography
Low problem is spent, needs to improve traditional phase unwrapping algorithm, more baseline InSAR phase unwrapping technologies are single baselines
More baseline InSAR phase unwrapping algorithms are studied in the extension of InSAR phase unwrapping technology, in conjunction with the different multiple interference phases in visual angle
Position information carries out Combined Treatment, and can effectively expand more baseline InSAR systems does not obscure elevation range, improves rolling topography scene
The precision of lower phase unwrapping, obtaining high-precision three-dimensional information to more baseline InSAR has significant application value.
Summary of the invention
The invention discloses a kind of more baseline least square phase unwrapping methods based on Quality Map weighting, and this method is first
The single order phase gradient main value of more baseline InSAR winding phases is calculated, and by its mirror-symmetric extension, calculates second order phase ladder
Degree carries out phase to the residue points in winding phase as prior information then using Quality Map judgement winding phase masses
Position compensation, and weighting is optimized based on second order phase gradient of the Quality Map to winding phase, finally carry out more baseline minimums two
Multiply phase unwrapping, obtains solution and twine phase, this method is enriched the interferometric phase information of target scene by more base line measurements, avoided
Residue points phase error problem of transmission in traditional list baseline least square phase unwrapping algorithm, improves phase unwrapping precision.
In order to facilitate the description contents of the present invention, make following term definition first:
Define 1, interference synthetic aperture radar (InSAR)
Interference synthetic aperture radar, which refers to, utilizes two width or two width obtained in same observation scene difference observation angle
The above SAR image carries out interference imaging processing, extracts interferometric phase information, several then in conjunction with radar system parameters, radar platform
The synthetic aperture radar technique of what location parameter and observation terrain information inverting Terrain Elevation and elevation change information.It is detailed in document
" synthetic aperture radar image-forming principle ", skin, which also rings etc., to be write, publishing house, University of Electronic Science and Technology.
Define 2, baseline
Baseline refers to the distance between two antennas in InSAR system, and InSAR system baseline length is denoted as B in the present invention.
It is detailed in document " synthetic aperture radar image-forming principle ", skin, which also rings, etc. writes, publishing house, University of Electronic Science and Technology.
Define 3, winding phase:
Winding phase refer to practical interference processing in, interferometric phase be wrapped in [- π, π) section, with its true phase
Differ the integral multiple of 2 π.It is detailed in document " spaceborne InSAR ", king is superfine to write, Science Press.
Define 4, phase unwrapping
The process that phase recovery is true phase corresponding with landform altitude will be wound and be known as phase unwrapping.It is detailed in document
" airborne Interference synthetic aperture radar data processing technique ", red violent wind of fourth etc. are write, Science Press.
It defines 5, solution and twines phase:
Solution twines phase and refers to the true phase restored after phase unwrapping is handled in interference processing from winding phase.It is detailed in
Document " spaceborne InSAR ", king is superfine to write, Science Press.
Define 6, Fast Fourier Transform (FFT) and inverse fast Fourier transform
The basic thought of Fast Fourier Transform (FFT) (Fast Fourier Transform, FFT) is original N point sequence
Column, resolve into a series of short sequence, make full use of the symmetry and periodicity of exponential factor in discrete Fourier transform, pass through
These appropriately combined short corresponding discrete Fourier transforms of sequence, reach deletion and compute repeatedly, and reduce multiplying and simplify knot
The purpose of structure.The inverse process of Fast Fourier Transform (FFT) is denoted as inverse fast Fourier transform (Inverse Fast Fourier
Transform, IFFT).It is detailed in document " Digital Signal Processing study course ", Cheng Peiqing writes, publishing house, Tsinghua University.
Define 7, Quality Map
Quality Map is the evaluation model of phase data quality, and each pixel corresponding phase quality is good in characterization phase diagram
It is bad.Quality Map is generally divided into three kinds: coherence factor figure, pseudo- coherence factor figure and phase derivative variation diagram.Coherence factor figure is derivative
In radar plural number image, show the coherence of different zones；Pseudo- coherence factor figure refers to that answer image true that can not obtain two width
In the case where information, it is assumed that amplitude and the variation of its phase information are consistent；Phase derivative variation diagram characterizes phase data quality
The fast-changing region of phase is defined as low quality region by the degree of " bad ", when hypsography is larger, phase derivative variation
Figure can be used as the first choice of Quality Map.It is detailed in document " solution of optimization fiber end face surface shape measurement twines algorithm research ", Ru Feng, Zhejiang
Jiang great Xue master thesis.
Define 8, residue points
Residue points, also referred to as singular value point refer to that winding phase gradient carries out ring product in 2 × 2 forms of current pixel point
Point, the point that integral result is 2 π defines the residue points that are positive, and the point that integral result is -2 π defines the residue points that are negative, and phase noise is dirty
The factors such as dye, violent hypsography can all cause the presence of residue points.Document is detailed in " at airborne Interference synthetic aperture radar data
Reason technology ", red violent wind of fourth etc. are write, Science Press.
Define 9, main value function
For trigonometric function, inverse function is multivalued function, in order to express determining monotropic function relationship, to multivalue
Antitrigonometric function is limited, take specific codomain with meet function define in monambiguity requirement, this restriction relation be known as it is more
It is worth the main value function of antitrigonometric function, is denoted as W [].It is detailed in document " mathematics guide: practical mathematics handbook ", Ai Bohadecai
De Le etc. writes, Science Press.
Define 10, mirror-symmetric extension method
Mirror-symmetric extension method refers to the method that image is realized to periodic extension in such a way that mirror symmetry turns down.In detail
See document " Digital Image Processing ", Paul Gonzales etc. are write, Electronic Industry Press.
Define 11, normalization
Normalization, also referred to as deviation standardize, and are the linear transformations to initial data, result is made to be mapped to [0,1] section, become
Exchange the letters number is expressed asWherein xmaxFor the maximum value of sample data, xminFor the minimum value of sample data.It is detailed in
Document " mathematics guide: practical mathematics handbook ", Ai Bohadecaidele etc. writes, Science Press.
The present invention provides a kind of more baseline least square phase unwrapping methods based on Quality Map weighting, it includes following
Several steps:
Ginseng needed for more baseline InSAR least square phase unwrapping methods that step 1, initialization are weighted based on Quality Map
Number:
Parameter needed for initializing the more baseline InSAR least square phase unwrapping methods weighted based on Quality Map, comprising:
More baseline InSAR wind phase data group number, are denoted as L, wind the distance of phase to points, be denoted as M, wind the orientation of phase
Points, are denoted as N, baseline length is denoted as Bl, more baseline InSAR l groups winding phases are denoted as φl(a, b), wherein l=1,2 ...,
L, a=1,2 ..., M, b=1,2 ..., N, a indicate distance to a-th point, b indicates b-th point of orientation；Wind phase
Make FFT transform corresponding pixel coordinate into frequency domain to be denoted as (j, k), wherein j=1,2 ..., M, k=1,2 ..., N, j indicate away from
J-th point of descriscent frequency domain, k indicate that k-th point of orientation frequency domain, FFT are the fast Fourier change for defining traditional standard in 6
It changes；
Step 2, the single order phase gradient main value for calculating winding phase distance descriscent and orientation:
Using formula Δl,x(a, b)=W [φl(a+1,b)-φl(a, b)], l group winding phase is calculatedl(a,
B) distance to single order phase gradient main value, be denoted as Δl,x(a,b)；Using formula Δl,y(a, b)=W [φl(a,b+1)-φl
(a, b)], l group winding phase is calculatedlThe single order phase gradient main value of (a, b) orientation, is denoted as Δl,y(a, b), l
=1,2 ..., L, a=1,2 ..., M, b=1,2 ..., N, wherein x indicates the distance of winding phase to y indicates winding phase
Orientation, L indicate that the winding phase data group number initialized in step 1, M indicate the winding phase distance initialized in step 1
To points, N indicates the winding phase orientation initialized in step 1 points, and W [] is the main value letter for defining traditional standard in 9
Number, φl(a, b) be step 1 defined in l group wind phase, a indicate distance to a-th point, b expression orientation b
A point；
Step 3, mirror-symmetric extension single order phase gradient main value matrix calculate second order phase gradient:
Using formulaL group is twined
Around phasel(a, b) distance is to single order phase gradient main value Δl,x(a, b) carry out mirror-symmetric extension, after obtaining continuation away from
Descriscent single order phase gradient main value, is denoted as ρl,x(a,b)；Using formulaPhase is wound to l grouplThe orientation (a, b)
To single order phase gradient main value Δl,y(a, b) carries out mirror-symmetric extension, the orientation single order phase gradient master after obtaining continuation
Value, is denoted as ρl,y(a, b), l=1,2 ..., L, a=1,2 ..., 2M, b=1,2 ..., 2N, wherein L indicates to initialize in step 1
Winding phase data group number, M indicate step 1 in initialize winding phase distance descriscent points, N indicate step 1 in initialize
Winding phase orientation points, Δl,x(a, b) and Δl,y(a, b) is respectively the l group winding phase calculated in step 2l
For (a, b) distance to the single order phase gradient main value with orientation, mirror-symmetric extension is the mirror image pair for defining traditional standard in 10
Claim continuation method；
Using formula ρl(a, b)=[ρl,x(a+1,b)-ρl,x(a,b)]+[ρl,y(a,b+1)-ρl,y(a, b)], it is calculated
L group winds phaselThe second order phase gradient of (a, b), is denoted as ρl(a, b), l=1,2 ..., L, a=1,2 ..., 2M, b=1,
2 ..., 2N, wherein L indicates that the winding phase data group number initialized in step 1, M indicate the winding phase initialized in step 1
For distance to points, N indicates the winding phase orientation initialized in step 1 points, ρl,x(a,b)、ρl,y(a, b), which is respectively indicated, to be prolonged
The single order phase gradient main value of winding phase distance descriscent, orientation after opening up；
Residue points in step 4, identification winding phase:
Using formulaA-th of distance is calculated in winding phase to b-th of side
Position to pixel (a, b) First-order Gradient main value in 2 × 2 neighborhoods ring integral, be denoted as R (a, b), wherein ρi(i=1,2,
It 3,4) is winding phase First-order Gradient main value, ρ1=ρl,x(a,b)、ρ3=-ρl,x(a, b+1) be in step 3 after continuation distance to
Single order phase gradient main value, ρ2=ρl,y(a+1,b)、ρ4=-ρl,y(a, b) be in step 3 continuation back side to single order phase ladder
Main value is spent, l=1,2 ..., L, a=1,2 ..., 2M, b=1,2 ..., 2N, L indicate the winding phase data initialized in step 1
Group number, M indicate the winding phase distance descriscent initialized in step 1 points, and N indicates the winding phase orientation initialized in step 1
To points；
The methods of residue points is in identification winding phase: when winding phase pixel point (a, b) First-order Gradient main value 2 ×
Ring integral result in 2 neighborhoods | R (a, b) | when=2 π, (a, b) is residue points, and as R (a, b)=0, (a, b) is continuity point；
Step 5 carries out phase compensation to residue points based on Quality Map:
Using formulaIt calculates
To the Quality Map of l group winding phase, a-th of distance is to the phase masses of the pixel (a, b) of b-th of orientation are denoted as zl
(a, b), l=1,2 ..., L, a=1,2 ..., 2M, b=1,2 ..., 2N, wherein L indicates the winding phase initialized in step 1
Data group number, M indicate the winding phase distance descriscent initialized in step 1 points, and N indicates the winding phase initialized in step 1
Orientation points, ρl,x(a,b)、ρl,y(a, b) respectively indicates single order phase ladder of the distance in step 3 after continuation to, orientation
Main value is spent,Respectively indicate distance of the current pixel point in m × m (2≤m≤M) window after continuation to,
The mean value of the single order phase gradient main value of orientation；
Using formulaPhase compensating factor is calculated, is denoted as ni(i=1,2,3,4), wherein
zi(i=1,2,3,4) is phase masses, z1Indicate the corresponding phase masses of residue points (a, b) in winding phase, z2Indicate point (a+
1, b) corresponding phase masses, z3Indicate point (a+1, b+1) corresponding phase masses, z4Indicate point (a, b+1) corresponding phase matter
Amount, R (a, b) are ring of the First-order Gradient main value for the winding phase pixel point (a, b) being calculated in step 4 in 2 × 2 neighborhoods
Integral；Using formulaCarry out residue points phase compensation, a=1,2 ..., 2M, b=1,2 ..., 2N, wherein ρi(i
=1,2,3,4) it is single order phase gradient main value in step 4,It is the single order phase gradient after phase compensation
Main value, M indicate the winding phase distance descriscent initialized in step 1 points, and N indicates the winding phase orientation initialized in step 1
To points；
Step 6 is based on Quality Map to second order phase gradient weighted optimization:
By the winding phase masses z after compensation residue points phase in step 5l(a, b), using traditional standard in definition 11
Deviation standardized method normalizes in [0,1] section, using formula Zl,x(a, b)=min (zl(a-1,b),zl(a, b)), meter
Calculation obtains distance to gradient weight, is denoted as Zl,x(a,b)；Using formula Zl,y(a, b)=min (zl(a,b-1),zl(a, b)), meter
Calculation obtains orientation weight gradient, is denoted as Zl,y(a, b), l=1,2 ..., L, a=1,2 ..., 2M, b=1,2 ..., 2N, wherein
zl(a, b) is that the l group being calculated in step 5 winds a-th of distance of phase to the phase matter of b-th of orientation pixel
Amount, L indicate that the winding phase data group number initialized in step 1, M indicate the winding phase distance descriscent point initialized in step 1
Number, N indicate the winding phase orientation initialized in step 1 points, and function min (parameter 1, parameter 2) indicates calculating parameter column
Minimum value in table；
Using formulaIt calculates
The second order phase gradient compensation factor of m group winding phase is obtained, is denoted asL=1,2 ..., L, m=1,2 ..., L,
A=1,2 ..., 2M, b=1,2 ..., 2N, wherein L indicates that the winding phase data group number initialized in step 1, M indicate step 1
The winding phase distance descriscent of middle initialization is counted, and N indicates the winding phase orientation initialized in step 1 points, Zl,x(a,b)
And Zl,y(a, b) is respectively gradient weight of the distance to gradient weight and orientation, BlIndicate the l articles baseline defined in step 1
Length, BmIndicate the length of the m articles baseline defined in step 1；
Using formula
M group after weighted optimization is calculated winds phase second order phase gradient, is denoted asM=1,2 ..., L, a=1,
2 ..., 2M, b=1,2 ..., 2N, wherein L indicates that the winding phase data group number initialized in step 1, M indicate in step 1 just
The winding phase distance descriscent of beginningization is counted, and N indicates the winding phase orientation initialized in step 1 points,For m
The second order phase gradient compensation factor of group winding phase, ρm(a, b) is the second order phase for the winding phase being calculated in step 3
Gradient；
Step 7 carries out phase unwrapping using second order phase gradient after weighted optimization:
Phase second order phase gradient is wound to the m group after weighted optimization in step 6It is traditional in 6 using defining
Winding phase second order phase gradient is calculated in the Fast Fourier Transform (FFT) method of standardTwo-dimensional FFT as a result, note
For P (j, k), using formulaSolution after continuation is calculated twines phaseTwo-dimensional FFT as a result, be denoted as Φ (j, k), j=1,2 ..., 2M, k=1,2 ..., 2N, j indicates the jth apart from frequency domain
A, k indicates that k-th point of orientation frequency domain, M indicate the winding phase distance descriscent initialized in step 1 points, and N indicates step 1
The winding phase orientation of middle initialization is counted；To Φ (j, k) using the inverse fast Fourier transform side for defining traditional standard in 6
Method, the solution after continuation is calculated twine phase, are denoted asA=1 at this time, 2 ..., 2M, b=1,2 ..., 2N takeMiddle a=1,2 ..., M, b=1, the value of 2 ..., N obtain final solution and twine phase, be denoted as ψ (a, b).
The innovation of the invention consists in that proposing a kind of more baseline least squares phase unwrapping side based on Quality Map weighting
Method, the present invention compensate the residue points in winding phase using Quality Map as prior information, avoid traditional single baseline minimum two
Multiply residue points phase error problem of transmission in phase unwrapping algorithm, winding phase second order phase gradient optimizing is added based on Quality Map
Power twines error to reduce solution, improves phase unwrapping precision.
The advantage of the invention is that realizing that simple, precision is high, applicability is good, pass through more baseline InSAR technological incorporation multiple groups
Phase information is wound, can preferably identify the region that elevation rises and falls in real terrain；Residue points in compensation winding phase, avoid
Error propagation problem caused by phase residual error point in traditional algorithm；Weighted optimization winds phase second order phase gradient, improves phase
Position solution twines precision.
Detailed description of the invention
Fig. 1 is the schematic process flow diagram of method provided by the present invention；
Fig. 2 is more baseline InSAR phase unwrapping simulation parameters of method provided by the present invention.
Specific embodiment
The present invention can be verified using the method for emulation experiment, and all steps, conclusion are all soft in MATLAB R2017b
It is verified on part correct.Specific implementation step is as follows:
Ginseng needed for more baseline InSAR least square phase unwrapping methods that step 1, initialization are weighted based on Quality Map
Number:
Parameter needed for initializing the more baseline InSAR least square phase unwrapping methods weighted based on Quality Map, comprising:
More baseline InSAR wind phase data group number, are denoted as L=5, wind the distance of phase to points, are denoted as M=500, wind phase
Orientation points, be denoted as N=500, baseline length is respectively B1=3m, B2=7m, B3=13m, B4=19m, B5=31m is more
Baseline InSAR l group winds phase, is denoted as φl(a, b), wherein l=1,2 ..., 5, a=1,2 ..., 500, b=1,2 ...,
500, a indicate distances to a-th point, b indicates b-th point of orientation；It is corresponding into frequency domain that winding phase makees FFT transform
Pixel coordinate is denoted as (j, k), wherein j=1, and 2 ..., 500, k=1,2 ..., 500, j indicate at j-th point apart from frequency domain, k table
K-th point for showing orientation frequency domain；
Step 2, the single order phase gradient main value for calculating winding phase distance descriscent and orientation:
Using formula Δl,x(a, b)=W [φl(a+1,b)-φl(a, b)], l group winding phase is calculatedl(a,
B) distance to single order phase gradient main value, be denoted as Δl,x(a,b)；Using formula Δl,y(a, b)=W [φl(a,b+1)-φl
(a, b)], l group winding phase is calculatedlThe single order phase gradient main value of (a, b) orientation, is denoted as Δl,y(a, b), l
=1,2 ..., L, a=1,2 ..., M, b=1,2 ..., N, L=5, M=500, N=500, wherein L is to initialize in step 1
Phase data group number is wound, M is the winding phase distance descriscent points initialized in step 1, and N is the winding initialized in step 1
Phase orientation points, x indicate the distance of winding phase to y indicates the orientation of winding phase, and W [] is to define tradition in 9
The main value function of standard, φl(a, b) be step 1 defined in l group wind phase, a expression distance to a-th point, b table
B-th point for showing orientation,；
Step 3, mirror-symmetric extension single order phase gradient main value matrix calculate second order phase gradient:
Using formulaL group is twined
Around phasel(a, b) distance is to single order phase gradient main value Δl,x(a, b) carry out mirror-symmetric extension, after obtaining continuation away from
Descriscent single order phase gradient main value, is denoted as ρl,x(a,b)；Using formulaPhase is wound to l grouplThe orientation (a, b)
To single order phase gradient main value Δl,y(a, b) carries out mirror-symmetric extension, the orientation single order phase gradient master after obtaining continuation
Value, is denoted as ρl,y(a, b), l=1,2 ..., L, a=1,2 ..., 2M, b=1,2 ..., 2N, L=5, M=500, N=500, wherein
L is the winding phase data group number initialized in step 1, and M is the winding phase distance descriscent points initialized in step 1, and N is step
The winding phase orientation points initialized in rapid 1, Δl,x(a, b) and Δl,y(a, b) is respectively the l group calculated in step 2
Wind phaselFor (a, b) distance to the single order phase gradient main value with orientation, mirror-symmetric extension is to define tradition mark in 10
Quasi- mirror-symmetric extension method；
Using formula ρl(a, b)=[ρl,x(a+1,b)-ρl,x(a,b)]+[ρl,y(a,b+1)-ρl,y(a, b)], it is calculated
L group winds phaselThe second order phase gradient of (a, b), is denoted as ρl(a, b), l=1,2 ..., L, a=1,2 ..., 2M, b=1,
2 ..., 2N, L=5, M=500, N=500, wherein L indicates that the winding phase data group number initialized in step 1, M indicate step
The winding phase distance descriscent points initialized in 1, N indicate the winding phase orientation initialized in step 1 points, ρl,x(a,
b)、ρl,y(a, b) respectively indicates the single order phase gradient main value of the winding phase distance descriscent after continuation, orientation；
Residue points in step 4, identification winding phase:
Using formulaA-th of distance is calculated in winding phase to b-th of side
Position to pixel (a, b) First-order Gradient main value in 2 × 2 neighborhoods ring integral, be denoted as R (a, b), wherein ρi(i=1,2,
It 3,4) is winding phase First-order Gradient main value, ρ1=ρl,x(a,b)、ρ3=-ρl,x(a, b+1) be in step 3 after continuation distance to
Single order phase gradient main value, ρ2=ρl,y(a+1,b)、ρ4=-ρl,y(a, b) be in step 3 continuation back side to single order phase ladder
Degree main value, l=1,2 ..., L, a=1,2 ..., 2M, b=1,2 ..., 2N, L=5, M=500, N=500, L are indicated in step 1
The winding phase data group number of initialization, M indicate the winding phase distance descriscent initialized in step 1 points, and N is indicated in step 1
The winding phase orientation of initialization is counted；
The method of residue points is in identification winding phase: when the First-order Gradient main value of pixel (a, b) in winding phase is 2
Ring integral result in × 2 neighborhoods | R (a, b) | when=2 π, (a, b) is residue points, and as R (a, b)=0, (a, b) is continuous
Point；
Step 5 carries out phase compensation to residue points based on Quality Map:
Using formulaIt calculates
To the Quality Map of l group winding phase, a-th of distance is to the phase masses of the pixel (a, b) of b-th of orientation are denoted as zl
(a, b), l=1,2 ..., L, a=1,2 ..., 2M, b=1,2 ..., 2N, L=5, M=500, N=500, wherein L indicates step 1
The winding phase data group number of middle initialization, M indicate the winding phase distance descriscent initialized in step 1 points, and N indicates step 1
The winding phase orientation of middle initialization is counted, ρl,x(a,b)、ρl,y(a, b) respectively indicate the distance in step 3 after continuation to,
The single order phase gradient main value of orientation,Current pixel point is respectively indicated in 2 × 2 windows after continuation
Single order phase gradient main value from distance to, orientation mean value；
Using formulaPhase compensating factor is calculated, is denoted as ni(i=1,2,3,4), wherein
zi(i=1,2,3,4) is phase masses, z1Indicate the corresponding phase masses of residue points (a, b) in winding phase, z2Indicate point (a+
1, b) corresponding phase masses, z3Indicate point (a+1, b+1) corresponding phase masses, z4Indicate point (a, b+1) corresponding phase matter
Amount, R (a, b) are ring of the First-order Gradient main value for the winding phase pixel point (a, b) being calculated in step 4 in 2 × 2 neighborhoods
Integral；Using formulaCarry out residue points phase compensation, a=1,2 ..., 2M, b=1,2 ..., 2N, M=500, N
=500, wherein M indicates the winding phase distance descriscent initialized in step 1 points, and N indicates the winding phase initialized in step 1
Orientation points, ρi(i=1,2,3,4) is the single order phase gradient main value in step 4,After being phase compensation
Single order phase gradient main value；
Step 6 is based on Quality Map to more baseline second order phase gradient weighted optimizations:
By the winding phase masses z after compensation residue points phase in step 5l(a, b), using deviation traditional in definition 11
Standardized method normalizes in [0,1] section, using formula Zl,x(a, b)=min (zl(a-1,b),zl(a, b)), it calculates
To distance to gradient weight, it is denoted as Zl,x(a,b)；Using formula Zl,y(a, b)=min (zl(a,b-1),zl(a, b)), it calculates
To orientation weight gradient, it is denoted as Zl,y(a, b), l=1,2 ..., L, a=1,2 ..., 2M, b=1,2 ..., 2N, L=5, M=
500, N=500, wherein L indicates that the winding phase data group number initialized in step 1, M indicate the winding initialized in step 1
Phase distance descriscent points, N indicate the winding phase orientation initialized in step 1 points, zl(a, b) is to calculate in step 5
The l group arrived winds a-th of distance of phase to the phase matter of b-th of orientation pixel, min (parameter 1, parameter 2) indicates meter
Calculate the minimum value in parameter list；
Using formulaIt calculates
The second order phase gradient compensation factor of m group winding phase is obtained, is denoted asL=1,2 ..., L, m=1,2 ..., L,
A=1,2 ..., 2M, b=1,2 ..., 2N, L=5, M=500, N=500, wherein L indicates the winding phase initialized in step 1
Data group number, M indicate the winding phase distance descriscent initialized in step 1 points, and N indicates the winding phase initialized in step 1
Orientation points, Zl,x(a, b) and Zl,y(a, b) is respectively gradient weight of the distance to gradient weight and orientation, baseline length
Respectively B1=3m, B2=7m, B3=13m, B4=19m, B5=31m；
Using formulaMeter
The m group winding phase second order phase gradient after obtaining weighted optimization is calculated, is denoted asM=1,2 ..., L, a=1,
2 ..., 2M, b=1,2 ..., 2N, L=5, M=500, N=500, wherein L indicates the winding phase data initialized in step 1
Group number, M indicate the winding phase distance descriscent initialized in step 1 points, and N indicates the winding phase orientation initialized in step 1
To points,The second order phase gradient compensation factor of phase, ρ are wound for m groupm(a, b) is to be calculated in step 3
Winding phase second order phase gradient；
Step 7 carries out phase unwrapping using second order phase gradient after weighted optimization:
Phase second order phase gradient is wound to the m group after weighted optimization in step 6It is traditional in 6 using defining
Winding phase second order phase gradient is calculated in the Fast Fourier Transform (FFT) method of standardTwo-dimensional FFT as a result, note
For P (j, k), using formulaSolution after continuation is calculated twines phaseTwo-dimensional FFT as a result, be denoted as Φ (j, k), j=1,2 ..., 2M, k=1,2 ..., 2N, M=500, N=500, j is indicated
Apart from j-th point of frequency domain, k indicates that k-th point of orientation frequency domain, M indicate the winding phase distance descriscent point initialized in step 1
Number, N indicate the winding phase orientation initialized in step 1 points；To Φ (j, k) using define 6 in traditional standard it is quick
Inverse Fourier transform method, the solution after continuation is calculated twine phase, are denoted asA=1 at this time, 2 ..., 2M, b=1,
2 ..., 2N takeMiddle a=1,2 ..., M, b=1, the value of 2 ..., N obtain final solution and twine phase, be denoted as ψ (a, b).
Claims (1)
1. it is a kind of based on Quality Map weighting more baseline least square phase unwrapping methods, it is characterized in that it the following steps are included:
Parameter needed for more baseline InSAR least square phase unwrapping methods that step 1, initialization are weighted based on Quality Map:
Parameter needed for initializing the more baseline InSAR least square phase unwrapping methods weighted based on Quality Map, comprising: more bases
Line InSAR winds phase data group number, is denoted as L, winds the distance of phase to points, be denoted as M, winds the orientation point of phase
Number, is denoted as N, baseline length is denoted as Bl, more baseline InSAR l groups winding phases are denoted as φl(a, b), wherein l=1,2 ..., L, a
=1,2 ..., M, b=1,2 ..., N, a indicate distance to a-th point, b indicates b-th point of orientation；Phase is wound to make
FFT transform corresponding pixel coordinate into frequency domain is denoted as (j, k), wherein j=1,2 ..., M, k=1, and 2 ..., N, j indicate distance
To j-th point of frequency domain, k indicates that k-th point of orientation frequency domain, FFT are the Fast Fourier Transform (FFT) of traditional standard；
Step 2, the single order phase gradient main value for calculating winding phase distance descriscent and orientation:
Using formula Δl,x(a, b)=W [φl(a+1,b)-φl(a, b)], l group winding phase is calculatedl(a, b) away from
The single order phase gradient main value of descriscent, is denoted as Δl,x(a,b)；Using formula Δl,y(a, b)=W [φl(a,b+1)-φl(a,
B)], l group winding phase is calculatedlThe single order phase gradient main value of (a, b) orientation, is denoted as Δl,y(a, b), l=1,
2 ..., L, a=1,2 ..., M, b=1,2 ..., N, wherein x indicates the distance of winding phase to y indicates the orientation of winding phase
To L indicates that the winding phase data group number initialized in step 1, M indicate the winding phase distance descriscent point initialized in step 1
Number, N indicate the winding phase orientation initialized in step 1 points, and W [] is the main value function of traditional standard, φl(a,b)
Wind phase for the l group that initializes in step 1, a indicate distance to a-th point, b indicates b-th point of orientation；
Step 3, mirror-symmetric extension single order phase gradient main value matrix calculate second order phase gradient:
Using formulaPhase is wound to l group
Position φl(a, b) distance is to single order phase gradient main value Δl,x(a, b) carries out mirror-symmetric extension, distance after obtaining continuation to
Single order phase gradient main value, is denoted as ρl,x(a,b)；Using formulaPhase is wound to l grouplThe orientation (a, b)
To single order phase gradient main value Δl,y(a, b) carries out mirror-symmetric extension, the orientation single order phase gradient master after obtaining continuation
Value, is denoted as ρl,y(a, b), l=1,2 ..., L, a=1,2 ..., 2M, b=1,2 ..., 2N, wherein L indicates to initialize in step 1
Winding phase data group number, M indicate step 1 in initialize winding phase distance descriscent points, N indicate step 1 in initialize
Winding phase orientation points, Δl,x(a, b) and Δl,y(a, b) is respectively the l group winding phase calculated in step 2l
For (a, b) distance to the single order phase gradient main value with orientation, mirror-symmetric extension is the mirror-symmetric extension side of traditional standard
Method；
Using formula ρl(a, b)=[ρl,x(a+1,b)-ρl,x(a,b)]+[ρl,y(a,b+1)-ρl,y(a, b)], l group is calculated
Wind phaselThe second order phase gradient of (a, b), is denoted as ρl(a, b), l=1,2 ..., L, a=1,2 ..., 2M, b=1,2 ...,
2N, wherein L indicates that the winding phase data group number initialized in step 1, M indicate the winding phase distance initialized in step 1
To points, N indicates the winding phase orientation initialized in step 1 points, ρl,x(a,b)、ρl,yAfter (a, b) respectively indicates continuation
Winding phase distance descriscent, orientation single order phase gradient main value；
Residue points in step 4, identification winding phase:
Using formulaA-th of distance is calculated in winding phase to b-th of orientation
Pixel (a, b) First-order Gradient main value in 2 × 2 neighborhoods ring integral, be denoted as R (a, b), wherein ρi(i=1,2,3,4)
To wind phase First-order Gradient main value, ρ1=ρl,x(a,b)、ρ3=-ρl,x(a, b+1) be in step 3 after continuation distance to single order
Phase gradient main value, ρ2=ρl,y(a+1,b)、ρ4=-ρl,y(a, b) be in step 3 continuation back side to single order phase gradient master
Value, l=1,2 ..., L, a=1,2 ..., 2M, b=1,2 ..., 2N, L indicate the winding phase data group initialized in step 1
Number, M indicate the winding phase distance descriscent initialized in step 1 points, and N indicates the winding phase orientation initialized in step 1
Points；
The method of residue points is in identification winding phase: when the First-order Gradient main value of winding phase pixel point (a, b) is in 2 × 2 neighbours
Ring integral result in domain | R (a, b) | when=2 π, (a, b) is residue points, and as R (a, b)=0, (a, b) is continuity point；
Step 5 carries out phase compensation to residue points based on Quality Map:
Using formulaL is calculated
The Quality Map of group winding phase, a-th of distance is to the phase masses of the pixel (a, b) of b-th of orientation are denoted as zl(a, b),
L=1,2 ..., L, a=1,2 ..., 2M, b=1,2 ..., 2N, wherein L indicates the winding phase data group initialized in step 1
Number, M indicate the winding phase distance descriscent initialized in step 1 points, and N indicates the winding phase orientation initialized in step 1
Points, ρl,x(a,b)、ρl,y(a, b) respectively indicates single order phase gradient main value of the distance in step 3 after continuation to, orientation,Distance of the current pixel point in m × m (2≤m≤M) window after continuation is respectively indicated to, orientation
Single order phase gradient main value mean value；
Using formulaPhase compensating factor is calculated, is denoted as ni(i=1,2,3,4), wherein zi(i=
It 1,2,3,4) is phase masses, z1Indicate the corresponding phase masses of residue points (a, b) in winding phase, z2Indicate that point (a+1, b) is right
The phase masses answered, z3Indicate point (a+1, b+1) corresponding phase masses, z4Indicate point (a, b+1) corresponding phase masses, R
(a, b) is ring product of the First-order Gradient main value for the winding phase pixel point (a, b) being calculated in step 4 in 2 × 2 neighborhoods
Point；Using formulaCarry out residue points phase compensation, a=1,2 ..., 2M, b=1,2 ..., 2N, wherein ρi(i=
1,2,3,4) it is single order phase gradient main value in step 4,It is the single order phase gradient master after phase compensation
Value, M indicate the winding phase distance descriscent initialized in step 1 points, and N indicates the winding phase orientation initialized in step 1
Points；
Step 6 is based on Quality Map to second order phase gradient weighted optimization:
By the winding phase masses z after compensation residue points phase in step 5l(a, b), using the deviation standardization side of traditional standard
Method normalizes in [0,1] section, using formula Zl,x(a, b)=min (zl(a-1,b),zl(a, b)), be calculated distance to
Gradient weight is denoted as Zl,x(a,b)；Using formula Zl,y(a, b)=min (zl(a,b-1),zl(a, b)), orientation is calculated
Weight gradient, is denoted as Zl,y(a, b), l=1,2 ..., L, a=1,2 ..., 2M, b=1,2 ..., 2N, wherein zl(a, b) is step
The l group being calculated in 5 winds a-th of distance of phase to the phase masses of b-th of orientation pixel, L indicates step 1
The winding phase data group number of middle initialization, M indicate the winding phase distance descriscent initialized in step 1 points, and N indicates step 1
The winding phase orientation of middle initialization is counted, and function min (parameter 1, parameter 2) indicates the minimum value in calculating parameter list；
Using formulaIt is calculated
M group winds the second order phase gradient compensation factor of phase, is denoted asL=1,2 ..., L, m=1,2 ..., L, a=1,
2 ..., 2M, b=1,2 ..., 2N, wherein L indicates that the winding phase data group number initialized in step 1, M indicate in step 1 just
The winding phase distance descriscent of beginningization is counted, and N indicates the winding phase orientation initialized in step 1 points, Zl,x(a, b) and Zl,y
(a, b) is respectively gradient weight of the distance to gradient weight and orientation, BlIndicate the l articles baseline initialized in step 1
Length, BmIndicate the length of the m articles baseline initialized in step 1；
Using formulaMeter
The m group winding phase second order phase gradient after obtaining weighted optimization is calculated, is denoted asM=1,2 ..., L, a=1,
2 ..., 2M, b=1,2 ..., 2N, wherein L indicates that the winding phase data group number initialized in step 1, M indicate in step 1 just
The winding phase distance descriscent of beginningization is counted, and N indicates the winding phase orientation initialized in step 1 points,For m
The second order phase gradient compensation factor of group winding phase, ρm(a, b) is the second order phase for the winding phase being calculated in step 3
Gradient；
Step 7 carries out phase unwrapping using second order phase gradient after weighted optimization:
Phase second order phase gradient is wound to the m group after weighted optimization in step 6Using quick Fu of traditional standard
In leaf transformation method, be calculated winding phase second order phase gradientTwo-dimensional FFT as a result, be denoted as P (j, k), use
FormulaSolution after continuation is calculated twines phaseTwo dimension
FFT result is denoted as Φ (j, k), j=1,2 ..., 2M, k=1,2 ..., 2N, and j indicates at j-th point apart from frequency domain, the expression side k
K-th point of position frequency domain, M indicate the winding phase distance descriscent initialized in step 1 points, and N indicates that is initialized in step 1 twines
It counts around phase orientation；The inverse fast Fourier transform method that traditional standard is used to Φ (j, k), after continuation is calculated
Solution twines phase, is denoted asA=1 at this time, 2 ..., 2M, b=1,2 ..., 2N takeMiddle a=1,2 ..., M, b=1,
The value of 2 ..., N obtain final solution and twine phase, is denoted as ψ (a, b).
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910271422.5A CN110109100B (en) | 2019-04-04 | 2019-04-04 | Multi-baseline least square phase unwrapping method based on quality map weighting |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910271422.5A CN110109100B (en) | 2019-04-04 | 2019-04-04 | Multi-baseline least square phase unwrapping method based on quality map weighting |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN110109100A true CN110109100A (en) | 2019-08-09 |
| CN110109100B CN110109100B (en) | 2022-05-03 |
Family
ID=67485256
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201910271422.5A Active CN110109100B (en) | 2019-04-04 | 2019-04-04 | Multi-baseline least square phase unwrapping method based on quality map weighting |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN110109100B (en) |
Cited By (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN110793463A (en) * | 2019-09-25 | 2020-02-14 | 西安交通大学 | Unwrapped phase error detection and correction method based on phase distribution |
| CN110907932A (en) * | 2019-11-26 | 2020-03-24 | 上海卫星工程研究所 | Distributed InSAR satellite height measurement precision influence factor analysis method and system |
| CN113238227A (en) * | 2021-05-10 | 2021-08-10 | 电子科技大学 | Improved least square phase unwrapping method and system combined with deep learning |
| CN114265062A (en) * | 2021-11-11 | 2022-04-01 | 电子科技大学 | InSAR phase unwrapping method based on phase gradient estimation network |
Citations (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20020041263A1 (en) * | 2000-08-28 | 2002-04-11 | Seiko Epson Corporation | System and method for providing an image processing circuit that improves image quality |
| US20040062453A1 (en) * | 1999-02-25 | 2004-04-01 | Ludwig Lester F. | Relative optical path phase reconstruction in the correction of misfocused images using fractional powers of the fourier transform |
| US20050273257A1 (en) * | 2004-06-03 | 2005-12-08 | Hager James R | Methods and systems for enhancing accuracy of terrain aided navigation systems |
| EP2413158A1 (en) * | 2010-07-26 | 2012-02-01 | Consorci Institut de Geomatica | A method for monitoring terrain and man-made feature displacements using ground-based synthetic aperture radar (GBSAR) data |
| CN107730491A (en) * | 2017-10-13 | 2018-02-23 | 北京工业大学 | A kind of phase unwrapping package method based on Quality Map |
| CN108008383A (en) * | 2017-11-09 | 2018-05-08 | 电子科技大学 | A kind of four FFT phase unwrapping methods of more baseline high accuracy of iteration |
-
2019
- 2019-04-04 CN CN201910271422.5A patent/CN110109100B/en active Active
Patent Citations (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20040062453A1 (en) * | 1999-02-25 | 2004-04-01 | Ludwig Lester F. | Relative optical path phase reconstruction in the correction of misfocused images using fractional powers of the fourier transform |
| US20020041263A1 (en) * | 2000-08-28 | 2002-04-11 | Seiko Epson Corporation | System and method for providing an image processing circuit that improves image quality |
| US20050273257A1 (en) * | 2004-06-03 | 2005-12-08 | Hager James R | Methods and systems for enhancing accuracy of terrain aided navigation systems |
| EP2413158A1 (en) * | 2010-07-26 | 2012-02-01 | Consorci Institut de Geomatica | A method for monitoring terrain and man-made feature displacements using ground-based synthetic aperture radar (GBSAR) data |
| CN107730491A (en) * | 2017-10-13 | 2018-02-23 | 北京工业大学 | A kind of phase unwrapping package method based on Quality Map |
| CN108008383A (en) * | 2017-11-09 | 2018-05-08 | 电子科技大学 | A kind of four FFT phase unwrapping methods of more baseline high accuracy of iteration |
Non-Patent Citations (3)
| Title |
|---|
| M.D. PRITT 等: ""Least-squares two-dimensional phase unwrapping using FFT"s"", 《IEEE TRANSACTIONS ON GEOSCIENCE AND REMOTE SENSING》 * |
| YU YONG 等: ""A phase unwrapping method based on minimum cost flows method in irregular network"", 《IEEE INTERNATIONAL GEOSCIENCE AND REMOTE SENSING SYMPOSIUM 》 * |
| 蒋留兵 等: ""基于误差迭代补偿的干涉合成孔雷达四向加权最小二乘相位解缠法"", 《科学技术与工程 》 * |
Cited By (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN110793463A (en) * | 2019-09-25 | 2020-02-14 | 西安交通大学 | Unwrapped phase error detection and correction method based on phase distribution |
| CN110907932A (en) * | 2019-11-26 | 2020-03-24 | 上海卫星工程研究所 | Distributed InSAR satellite height measurement precision influence factor analysis method and system |
| CN113238227A (en) * | 2021-05-10 | 2021-08-10 | 电子科技大学 | Improved least square phase unwrapping method and system combined with deep learning |
| CN114265062A (en) * | 2021-11-11 | 2022-04-01 | 电子科技大学 | InSAR phase unwrapping method based on phase gradient estimation network |
| CN114265062B (en) * | 2021-11-11 | 2023-11-10 | 电子科技大学 | InSAR phase unwrapping method based on phase gradient estimation network |
Also Published As
| Publication number | Publication date |
|---|---|
| CN110109100B (en) | 2022-05-03 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN110109100A (en) | A kind of more baseline least square phase unwrapping methods based on Quality Map weighting | |
| CN104123464B (en) | Method for inversion of ground feature high elevation and number of land subsidence through high resolution InSAR timing sequence analysis | |
| CN102681033B (en) | Sea surface wind measurement method based on X-band marine radar | |
| CN104459633B (en) | Wavelet field InSAR interferometric phase filtering method in conjunction with local frequency estimation | |
| Zhou et al. | Deep learning-based branch-cut method for InSAR two-dimensional phase unwrapping | |
| CN103454636B (en) | Differential interferometric phase estimation method based on multi-pixel covariance matrixes | |
| CN103235301A (en) | Polarimetric synthetic aperture radar interferometry (POLInSAR) vegetation height inversion method based on complex field adjustment theory | |
| CN102621549A (en) | Multi-baseline/multi-frequency-band interference phase unwrapping frequency domain quick algorithm | |
| CN104698462B (en) | Synthetic aperture radar Ocean Wind-field fusion method based on variation | |
| CN104730519A (en) | High-precision phase unwrapping method adopting error iteration compensation | |
| CN104808203A (en) | Multi-baseline InSAR phase unwrapping method by iterating maximum likelihood estimation | |
| CN107064933B (en) | SAR chromatography building height method based on cyclic spectrum estimation | |
| CN103809180B (en) | For InSAR topographic Pre-Filter processing method | |
| CN103630903B (en) | The method of flow field, sea radial velocity is measured based on straight rail interference SAR | |
| CN112685819A (en) | Data post-processing method and system for monitoring dam and landslide deformation GB-SAR | |
| CN103630898B (en) | To the method that multi-baseline interference SAR phase bias is estimated | |
| CN108132468A (en) | A kind of more baseline polarimetric SAR interferometry depth of building extracting methods | |
| Budd et al. | The scaling and skewness of optimally transported meshes on the sphere | |
| Wu et al. | A deep learning based method for local subsidence detection and InSAR phase unwrapping: Application to mining deformation monitoring | |
| CN106097404B (en) | Utilize the method for non-linear vector Surface Construction InSAR phase image model | |
| Li et al. | An adaptive phase optimization algorithm for distributed scatterer phase history retrieval | |
| CN111025294A (en) | InSAR phase unwrapping method based on mean square volume Kalman filter | |
| CN103940834B (en) | Synthetic Aperture Radar Technique is adopted to measure the method for soil moisture | |
| CN111291316B (en) | Multi-scale resistivity inversion method and system based on wavelet transformation | |
| CN104484880A (en) | SAR image segmentation method based on coherence map migration and clustering |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |