CN107977939A - A kind of weighted least-squares phase unwrapping computational methods based on reliability - Google Patents
A kind of weighted least-squares phase unwrapping computational methods based on reliability Download PDFInfo
- Publication number
- CN107977939A CN107977939A CN201711226954.4A CN201711226954A CN107977939A CN 107977939 A CN107977939 A CN 107977939A CN 201711226954 A CN201711226954 A CN 201711226954A CN 107977939 A CN107977939 A CN 107977939A
- Authority
- CN
- China
- Prior art keywords
- mrow
- phase
- msup
- msub
- reliability
- 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
- 238000000205 computational method Methods 0.000 title claims abstract description 20
- 238000004364 calculation method Methods 0.000 claims description 16
- 238000010587 phase diagram Methods 0.000 claims description 8
- 238000000034 method Methods 0.000 claims description 7
- 239000000203 mixture Substances 0.000 claims description 4
- 230000009466 transformation Effects 0.000 claims description 4
- 238000001228 spectrum Methods 0.000 claims description 3
- 206010061619 Deformity Diseases 0.000 claims 1
- 230000001815 facial effect Effects 0.000 claims 1
- 238000009499 grossing Methods 0.000 abstract description 9
- 230000000694 effects Effects 0.000 abstract description 5
- 230000005540 biological transmission Effects 0.000 abstract description 3
- 230000008030 elimination Effects 0.000 abstract 1
- 238000003379 elimination reaction Methods 0.000 abstract 1
- 238000005516 engineering process Methods 0.000 description 5
- 230000006870 function Effects 0.000 description 2
- 238000005305 interferometry Methods 0.000 description 2
- 238000005259 measurement Methods 0.000 description 2
- 239000004411 aluminium Substances 0.000 description 1
- XAGFODPZIPBFFR-UHFFFAOYSA-N aluminium Chemical compound [Al] XAGFODPZIPBFFR-UHFFFAOYSA-N 0.000 description 1
- 229910052782 aluminium Inorganic materials 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000008859 change Effects 0.000 description 1
- 239000002131 composite material Substances 0.000 description 1
- 238000009434 installation Methods 0.000 description 1
- 238000012804 iterative process Methods 0.000 description 1
- 238000000691 measurement method Methods 0.000 description 1
- 230000003287 optical effect Effects 0.000 description 1
- 230000008569 process Effects 0.000 description 1
- 230000035945 sensitivity Effects 0.000 description 1
Classifications
-
- G06T5/70—
Abstract
The invention discloses a kind of weighted least-squares phase unwrapping computational methods based on reliability.By the speckle interference figure of industrial camera equipment shooting, collecting determinand, obtain the two-dimensional phase comprising determinand three-dimensional information by image procossing and wrap up figure；The reliability for obtaining every bit in phase parcel figure is calculated using the second differnce of wrapped phase value；Calculate and determine reliability threshold value, calculate and obtain the binaryzation mask factor as weighted least-squares phase unwrapping weight, be iterated calculating, obtain final true phase.The present invention not only have the advantages that weighted least-squares computational methods provide smoothing solution, and have the advantages that calculating speed is fast, precision is high, it is effective suppress noise transmission and elimination smoothing effect.
Description
Technical field
The present invention relates to digital speckle interference technical field, and in particular to a kind of weighted least-squares phase based on reliability
Position unfolding calculation method.
Background technology
Digital speckle interference technology (Digital Speckle Pattern Interferometry, DSPI) is a kind of sharp
The optical interferometry technology that light technology is combined with Digital Image Processing, compared with conventional measurement techniques, has high-precision, non-
The advantages that contacting harmless, high sensitivity and measurement of full field, is widely used to space flight and aviation, biomedicine, shipbuilding etc.
Field.Phase unwrapping is step the most key in the application of DSPI technologies, directly affects measurement accuracy.In recent years, learn both at home and abroad
Many phase unwrapping computational methods have been proposed in person.Wherein, least square phase unwrapping computational methods utilize wrapped phase
The criterion of the difference minimum of discrete partial differential and the discrete partial differential of true phase, which carrys out unpacking, can eliminate bracing wire phenomenon, be put down
Sliding solution, is that a kind of simple robust, arithmetic speed be fast and the most commonly used computational methods small to request memory.But work as wrapped phase
In there are when noise spot, null point, low-key system point, least square phase unwrapping computational methods can not be solved correctly.Add
Power least square phase unwrapping computational methods are overcome the influence that noise etc. is brought by weight, but the weight introduced is fitted
When whether can also influence the correctness of unpacking, while also reduce calculating speed, and cannot suppress to make an uproar and sonic propagation and disappear
The error brought except smoothing effect.
The content of the invention
In order to solve the above technical problem, the present invention provides a kind of weighted least-squares phase unwrapping based on reliability
Computational methods, calculating speed is fast, can effectively suppress noise transmission, and eliminate smoothing effect.
The present invention is achieved through the following technical solutions：
Step 1：By the speckle interference figure of industrial camera equipment shooting, collecting determinand, wrapped by image procossing
The two-dimensional phase parcel figure that the size of the three-dimensional information containing determinand is M × N；
Step 2：The reliability for obtaining every bit in phase parcel figure is calculated using the second differnce of wrapped phase value；
Step 3：Determine reliability threshold value, calculate and obtain the binaryzation mask factor, and be used as weighted least-squares phase exhibition
Open weight；
Step 4：Calculating is iterated according to weighted least-squares phase unwrapping weight, obtains final true phase.
Step 4 includes：Initialize iterations k and absolute phase φ0, calculate the wrapped phase each put's
Weight discrete partial differential；According to iterative formula, carry out kth time least square phase unwrapping and obtain φk+1；Judge φk+1It is whether full
The sufficient condition of convergence, if satisfied, then obtaining true phase φk+1；If not satisfied, then iterations k=k+1, and return and continue to change
Generation.
The present invention be directed to plate class determinand, collection be plate class determinand areal deformation speckle interference figure.
Plate class determinand passes through external force stress thin plate class determinand table using plate made of aluminium, composite material etc.
The deformation of indent or evagination occurs for face, then is detected using the method for the present invention.
The step 1, is specially：By the speckle interference figure before and after industrial camera equipment shooting, collecting composition deformation to be measured,
Then the two-dimensional phase parcel figure that the size comprising determinand 3 D deformation information is M × N is obtained by image procossing.
The step 1 is handled with spatial carrier method and obtains wrapped phase figure：Determinand is shot with industrial camera apparatus
Two width speckle interference figures are carried out Fourier transformation and obtained by the two width speckle interference figures with carrier wave amount before and after areal deformation respectively
Two amplitude-frequency spectrograms are obtained, select the positive level-one frequency spectrum in two amplitude-frequency spectrograms to obtain the positive level-one spectrogram of two width respectively, then to two width
Positive level-one spectrogram, which does inverse Fourier transform and negates tangent, obtains two amplitude phase diagrams, finally subtracts each other to obtain band by two amplitude phase diagrams
Have the phase parcel figure of determinand areal deformation information, phase parcel figure is filtered and denoising after to obtain size be M × N
Two-dimensional phase parcel figure
The step 2, is specially：
The second differnce of every bit wrapped phase value in two-dimensional phase parcel figure is first calculated using the following formula：
Wherein, V (i, j) be horizontal direction on second differnce, H (i, j) be vertical direction on second differnce, D1(i,
And D j)2(i, j) is the second differnce in two oblique line directions, and two oblique line directions are respectively in 45 degree with image level direction
Both direction positioned at both sides, wrap { } they are parcel computing, operation result is limited to (- π, π] in the range of；DIF2(i,j)
Represent the second differnce of coordinate position point (i, j) in two-dimensional phase parcel figure,Expression two-dimensional phase parcel figure midpoint (i,
J) wrapped phase value, point (i, j) represents abscissa i-th, j-th point of ordinate, 0≤i≤M-1,0≤j≤N-1；
Then, the reliability R (i, j) of the wrapped phase value of each point is calculated using the following formula：
Wherein, R (i, j) represents the reliability of the wrapped phase value of two-dimensional phase parcel figure midpoint (i, j).
The step 3, is specially：
The reliability of all the points in two-dimentional wrapped phase figure is sorted to obtain one-dimensional sequence q from small to large, is taken in sequence q
The reliability value of M × N × 5% point is as reliability threshold θ, the downward rounding if M × N × 5% is not integer, further according to
The following formula calculates the binaryzation mask factor W (i, j) for obtaining and each putting：
Wherein, W (i, j) represents the binaryzation mask factor of two-dimensional phase parcel figure midpoint (i, j).
The step 4, is specially：
4.1) Initialize installation iterations k and absolute phase φ0, the parcel phase each put first is calculated using the following formula
Place valueThe corresponding discrete partial differential c (i, j) of weighting：
Wherein,Represent the first-order difference after the parcel computing of wrapped phase vertical direction,Represent bag
Wrap up in the first-order difference after the parcel computing in phase levels direction；Two first-order differencesWithIt is calculated as respectively：
4.2) least square phase unwrapping is carried out in the following ways, obtains the true phase iterated to calculate every time：
By the true phase φ of kth time iterative calculationkSubstitute into the following formula and obtain the true phase of+1 iterative calculation of kth
Second order local derviation ρk+1：
ρk+1=c-F (φk)
Wherein, F (φk) represent true phase φkThe discrete partial differential function of weighting, c represent two-dimensional phase parcel figure in institute
The vector of discrete partial differential c (i, the j) composition of weighting a little；
True phase φkThe discrete partial differential function F (φ of weightingk) such as the following formula：
Wherein,Represent the first-order difference of true phase vertical direction,Represent true phase level side
To first-order difference, two first-order differencesWithIt is calculated as respectively：
During initial calculation, iterations k=0, true phase φ are set in specific implementation0For absolute phase r, φ0=r.
4.3) by the second order local derviation ρ of the true phase of+1 iterative calculation of kthk+1The Poisson side represented as the following formula
The input of journey, the Poisson's equation that the following formula expression is then solved with discrete cosine transform (DCT) obtain+1 iterative calculation of kth
True phase φk+1：
4.4) the true phase φ of+1 iterative calculation of kth is judgedk+1Whether the condition of convergence is met：
If satisfied, then with true phase φk+1As final true phase；
If not satisfied, then repeating the above steps 4.2)~4.3) carry out the true phase that processing obtains next iteration calculating
Position, iterations k=k+1.So as to be iterated processing by above-mentioned steps, until whether true phase meets the condition of convergence.
The condition of convergence such as the following formula represents, meets that the following formula then meets the condition of convergence：
Wherein, ε represents convergence threshold, ε=10-2；φk+1Represent the true phase of+1 iterative calculation of kth, φkRepresent the
The true phase of k iterative calculation.
Compared with prior art, the beneficial effects of the invention are as follows：
The present invention obtains binaryzation mask as the weight of weighted least-squares method by the use of RELIABILITY INDEX, only carries out once
The calculating (weight in i.e. all iterative process is consistent) of weight, does not further relate to the meter again of weight in successive iterations
Calculate, and by iterative solution true phase, can effectively suppress to make an uproar sonic propagation and eliminate smoothing effect, and improve calculating
Speed.
The present invention not only have the advantages that weighted least-squares computational methods provide smoothing solution, and with calculating speed soon,
The advantages that effectively suppressing noise transmission and eliminating smoothing effect.
Brief description of the drawings
Fig. 1 is phase unwrapping computational methods flow chart of the present invention；
Fig. 2 is the filtered two-dimensional phase parcel figure of embodiment；
Fig. 3 is the composition figure for the weighted factor that embodiment is obtained based on reliability；
Phase result figure is unfolded for embodiment in Fig. 4；
Fig. 5 is that the rewind of phase is unfolded around figure in embodiment.
Embodiment
The invention will be further described with reference to the accompanying drawings and examples.
The embodiment of the present invention as shown in the flowchart of fig.1, comprises the following steps that：
Step 1：Obtain by the front and rear speckle interference figure of industrial camera equipment shooting deformation and by respective image processing
The two-dimensional phase parcel figure that size comprising determinand 3 D deformation information is M × N
The present embodiment specifically uses spatial carrier method, and detailed process is：Carried before and after being deformed with the shooting of industrial camera apparatus
Two width speckle interference figures of carrier wave amount, carry out two width speckle interference figures Fourier transformation and obtain two amplitude-frequency spectrograms, respectively respectively
Select the positive level-one frequency spectrum in two amplitude-frequency spectrograms to obtain the positive level-one spectrogram of two width, inverse Fu is then to the positive level-one spectrogram of two width
In leaf transformation and negate tangent and obtain two amplitude phase diagrams, two amplitude phase diagrams finally are subtracted each other to obtain the phase bag with deformation information
Wrap up in figure, phase parcel figure is filtered and denoising after obtain the two-dimensional phase parcel figure for treating unpacking as shown in Figure 2, obtain
Obtained the wrapped phase value of each pointWherein 0≤i≤M-1,0≤j≤N-1.
Step 2：The reliability R (i, j) of the every bit of phase parcel figure is calculated using the second differnce of phase value.
The second differnce of every bit phase value is calculated by formula (1) in phase diagram：
WillSubstitute into formula (1) and the second differnce of each phase value is calculated, then each phase is calculated by formula (2)
The reliability R (i, j) of place value：
Step 3：Determine reliability threshold value, obtain the binaryzation mask factor and calculated as weighted least-squares phase unwrapping
The weight of method.
The reliability of two-dimentional every, wrapped phase figure is sorted to obtain one-dimensional sequence q from small to large, takes the M in sequence q
The reliability value of × N × 5% point is as reliability threshold θ, if M × N × 5% is not integer, to its downward rounding, further according to
Formula (3), which calculates, obtains binaryzation mask factor W (i, j), as shown in Figure 3.
Step 4：Calculate wrapped phaseThe discrete partial differential of weighting.
WillSubstitute into the discrete partial differential c (i, j) of weighting that formula (4) calculates each wrapped phase：
Step 5：Initialize iterations k and absolute phase φ0, make iterations k=0, absolute phase φ0=r, foundation
Iterative formula carries out kth time least square phase unwrapping and obtains φk+1。
By φkSubstitute into iterative formula ρk+1=c-F (φk) obtain ρk+1, wherein F is the operation as shown in formula (5)：
Again by ρk+1As the input of Poisson's equation shown in formula (6), then the party is solved with discrete cosine transform (DCT)
Journey obtains true phase φk+1。
Step 6：The condition of convergence is judged, if satisfied, then obtaining true phase φk+1；If not satisfied, then iterations k=k
+ 1, and return to step five, continue iteration.
The condition of convergence is by formula (7) Suo Shi, and convergence threshold ε=10-2.If meeting the condition of convergence, final true phase is
φk+1；If being unsatisfactory for the condition of convergence, iterations k=k+1, and return to step five continues iteration.
The results are shown in Figure 4 for the present embodiment phase unwrapping, it can be seen that present invention expansion smoothing pseudorange, is illustrated in figure 5
Phase rewind is unfolded around figure in it, it can be seen that it is consistent with original parcel phase diagram striped around bar graph to be unfolded the rewind of phase, guarantor
Effectiveness of the invention is demonstrate,proved.
It is the problem of present invention is difficult to determine for weight in weighted least-square solution parcel computational methods, true using reliability
Binaryzation mask is determined as weight, and true phase is obtained by iteration, improves the reliability of weighting.The present invention has minimum
Two multiply the advantages of computational methods provide smoothing solution, eliminate bracing wire phenomenon, while determine that weight solves using RELIABILITY INDEX
The error that smooth belt is come, improves precision, and inhibit sonic propagation of making an uproar.
Claims (6)
1. a kind of weighted least-squares phase unwrapping computational methods based on reliability, it is characterised in that include the following steps：
Step 1：By the speckle interference figure of industrial camera equipment shooting, collecting determinand, obtain to include by image procossing and treat
The two-dimensional phase parcel figure that the size for surveying thing three-dimensional information is M × N；
Step 2：The reliability for obtaining every bit in phase parcel figure is calculated using the second differnce of wrapped phase value；
Step 3：Determine reliability threshold value, calculate and obtain the binaryzation mask factor, and as weighted least-squares phase unwrapping power
Weight；
Step 4：Calculating is iterated according to weighted least-squares phase unwrapping weight, obtains final true phase.
2. a kind of weighted least-squares phase unwrapping computational methods based on reliability according to claim 1, its feature
It is：The step 1 is handled with spatial carrier method and obtains wrapped phase figure：Determinand table is shot with industrial camera apparatus
Two width speckle interference figures are carried out Fourier transformation acquisition by the two width speckle interference figures with carrier wave amount before and after facial disfigurement respectively
Two amplitude-frequency spectrograms, select the positive level-one frequency spectrum in two amplitude-frequency spectrograms to obtain the positive level-one spectrogram of two width, then to two width just respectively
Level-one spectrogram, which does inverse Fourier transform and negates tangent, obtains two amplitude phase diagrams, finally subtracts each other two amplitude phase diagrams and is carried
The phase parcel figure of determinand areal deformation information, phase parcel figure is filtered and denoising after to obtain size be M × N
Two-dimensional phase parcel figure
3. a kind of weighted least-squares phase unwrapping computational methods based on reliability according to claim 1, its feature
It is：The step 2, is specially：
The second differnce of every bit wrapped phase value in two-dimensional phase parcel figure is first calculated using the following formula：
<mrow>
<msub>
<mi>DIF</mi>
<mn>2</mn>
</msub>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<msqrt>
<mrow>
<mi>V</mi>
<msup>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>)</mo>
</mrow>
<mn>2</mn>
</msup>
<mo>+</mo>
<mi>H</mi>
<msup>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>)</mo>
</mrow>
<mn>2</mn>
</msup>
<mo>+</mo>
<msub>
<mi>D</mi>
<mn>1</mn>
</msub>
<msup>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>)</mo>
</mrow>
<mn>2</mn>
</msup>
<mo>+</mo>
<msub>
<mi>D</mi>
<mn>2</mn>
</msub>
<msup>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>)</mo>
</mrow>
<mn>2</mn>
</msup>
</mrow>
</msqrt>
</mrow>
Wherein, V (i, j) be horizontal direction on second differnce, H (i, j) be vertical direction on second differnce, D1(i, j) and D2
(i, j) is the second differnce in two oblique line directions, and wrap { } is to wrap up computing, DIF2(i, j) represents two-dimensional phase parcel
The second differnce of coordinate position point (i, j) in figure,The wrapped phase value of expression two-dimensional phase parcel figure midpoint (i, j), 0
≤ i≤M-1,0≤j≤N-1；
Then, the reliability R (i, j) of the wrapped phase value of each point is calculated using the following formula：
<mrow>
<mi>R</mi>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mfrac>
<mn>1</mn>
<mrow>
<msub>
<mi>DIF</mi>
<mn>2</mn>
</msub>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>)</mo>
</mrow>
</mrow>
</mfrac>
<mo>.</mo>
</mrow>
Wherein, R (i, j) represents the reliability of the wrapped phase value of two-dimensional phase parcel figure midpoint (i, j).
4. a kind of weighted least-squares phase unwrapping computational methods based on reliability according to claim 1, its feature
It is：The step 3, is specially：
The reliability of all the points in two-dimentional wrapped phase figure is sorted to obtain one-dimensional sequence q from small to large, takes the M in sequence q
The reliability value of × N × 5% point is as reliability threshold θ, the downward rounding if M × N × 5% is not integer, further according to
Lower formula calculates the binaryzation mask factor W (i, j) for obtaining and each putting：
<mrow>
<mi>W</mi>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mfenced open = "{" close = "">
<mtable>
<mtr>
<mtd>
<mn>1</mn>
<mo>,</mo>
<mi>R</mi>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>)</mo>
<mo>&GreaterEqual;</mo>
<mi>&theta;</mi>
</mtd>
</mtr>
<mtr>
<mtd>
<mn>0</mn>
<mo>,</mo>
<mi>R</mi>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>)</mo>
<mo><</mo>
<mi>&theta;</mi>
</mtd>
</mtr>
</mtable>
</mfenced>
</mrow>
Wherein, W (i, j) represents the binaryzation mask factor of two-dimensional phase parcel figure midpoint (i, j).
5. a kind of weighted least-squares phase unwrapping computational methods based on reliability according to claim 1, its feature
It is：The step 4, is specially：
4.1) the wrapped phase value each put first is calculated using the following formulaThe corresponding discrete partial differential c (i, j) of weighting：
Wherein,Represent the first-order difference after the parcel computing of wrapped phase vertical direction,Represent wrapped phase
First-order difference after the parcel computing of horizontal direction；Two first-order differencesWithIt is calculated as respectively：
4.2) least square phase unwrapping is carried out in the following ways, obtains the true phase iterated to calculate every time：
By the true phase φ of kth time iterative calculationkSubstitute into the following formula obtains the true phase of+1 iterative calculation of kth two
Rank local derviation ρk+1：
ρk+1=c-F (φk)
Wherein, F (φk) represent true phase φkThe discrete partial differential function of weighting, c represent two-dimensional phase parcel figure in all the points
Weighting discrete partial differential c (i, j) composition vector；
True phase φkThe discrete partial differential function F (φ of weightingk) such as the following formula：
<mfenced open = "" close = "">
<mtable>
<mtr>
<mtd>
<mrow>
<mi>F</mi>
<mrow>
<mo>(</mo>
<msub>
<mi>&phi;</mi>
<mi>k</mi>
</msub>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mrow>
<mo>(</mo>
<mi>min</mi>
<mo>(</mo>
<mrow>
<msup>
<mrow>
<mo>(</mo>
<mrow>
<mi>W</mi>
<mrow>
<mo>(</mo>
<mrow>
<mi>i</mi>
<mo>+</mo>
<mn>1</mn>
<mo>,</mo>
<mi>j</mi>
</mrow>
<mo>)</mo>
</mrow>
</mrow>
<mo>)</mo>
</mrow>
<mn>2</mn>
</msup>
<mo>,</mo>
<msup>
<mrow>
<mo>(</mo>
<mrow>
<mi>W</mi>
<mrow>
<mo>(</mo>
<mrow>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
</mrow>
<mo>)</mo>
</mrow>
</mrow>
<mo>)</mo>
</mrow>
<mn>2</mn>
</msup>
</mrow>
<mo>)</mo>
<mo>-</mo>
<mn>1</mn>
<mo>)</mo>
</mrow>
<msubsup>
<mi>&Delta;</mi>
<mi>&phi;</mi>
<mi>x</mi>
</msubsup>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>)</mo>
</mrow>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mi>min</mi>
<mo>(</mo>
<mrow>
<msup>
<mrow>
<mo>(</mo>
<mrow>
<mi>W</mi>
<mrow>
<mo>(</mo>
<mrow>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
</mrow>
<mo>)</mo>
</mrow>
</mrow>
<mo>)</mo>
</mrow>
<mn>2</mn>
</msup>
<mo>,</mo>
<msup>
<mrow>
<mo>(</mo>
<mrow>
<mi>W</mi>
<mrow>
<mo>(</mo>
<mrow>
<mi>i</mi>
<mo>-</mo>
<mn>1</mn>
<mo>,</mo>
<mi>j</mi>
</mrow>
<mo>)</mo>
</mrow>
</mrow>
<mo>)</mo>
</mrow>
<mn>2</mn>
</msup>
</mrow>
<mo>)</mo>
<mo>-</mo>
<mn>1</mn>
<mo>)</mo>
</mrow>
<msubsup>
<mi>&Delta;</mi>
<mi>&phi;</mi>
<mi>x</mi>
</msubsup>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>-</mo>
<mn>1</mn>
<mo>,</mo>
<mi>j</mi>
<mo>)</mo>
</mrow>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<mo>+</mo>
<mrow>
<mo>(</mo>
<mi>min</mi>
<mo>(</mo>
<mrow>
<msup>
<mrow>
<mo>(</mo>
<mrow>
<mi>W</mi>
<mrow>
<mo>(</mo>
<mrow>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>+</mo>
<mn>1</mn>
</mrow>
<mo>)</mo>
</mrow>
</mrow>
<mo>)</mo>
</mrow>
<mn>2</mn>
</msup>
<mo>,</mo>
<msup>
<mrow>
<mo>(</mo>
<mrow>
<mi>W</mi>
<mrow>
<mo>(</mo>
<mrow>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
</mrow>
<mo>)</mo>
</mrow>
</mrow>
<mo>)</mo>
</mrow>
<mn>2</mn>
</msup>
</mrow>
<mo>)</mo>
<mo>-</mo>
<mn>1</mn>
<mo>)</mo>
</mrow>
<msubsup>
<mi>&Delta;</mi>
<mi>&phi;</mi>
<mi>y</mi>
</msubsup>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>)</mo>
</mrow>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mi>min</mi>
<mo>(</mo>
<mrow>
<msup>
<mrow>
<mo>(</mo>
<mrow>
<mi>W</mi>
<mrow>
<mo>(</mo>
<mrow>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
</mrow>
<mo>)</mo>
</mrow>
</mrow>
<mo>)</mo>
</mrow>
<mn>2</mn>
</msup>
<mo>,</mo>
<msup>
<mrow>
<mo>(</mo>
<mrow>
<mi>W</mi>
<mrow>
<mo>(</mo>
<mrow>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
<mo>)</mo>
</mrow>
</mrow>
<mo>)</mo>
</mrow>
<mn>2</mn>
</msup>
</mrow>
<mo>)</mo>
<mo>-</mo>
<mn>1</mn>
<mo>)</mo>
</mrow>
<msubsup>
<mi>&Delta;</mi>
<mi>&phi;</mi>
<mi>y</mi>
</msubsup>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>-</mo>
<mn>1</mn>
<mo>)</mo>
</mrow>
</mrow>
</mtd>
</mtr>
</mtable>
</mfenced>
Wherein,Represent the first-order difference of true phase vertical direction,Represent true phase horizontal direction
First-order difference, two first-order differencesWithIt is calculated as respectively：
<mrow>
<msubsup>
<mi>&Delta;</mi>
<mi>&phi;</mi>
<mi>x</mi>
</msubsup>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mfenced open = "{" close = "">
<mtable>
<mtr>
<mtd>
<mi>&phi;</mi>
<mo>(</mo>
<mi>i</mi>
<mo>+</mo>
<mn>1</mn>
<mo>,</mo>
<mi>j</mi>
<mo>)</mo>
<mo>-</mo>
<mi>&phi;</mi>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>)</mo>
<mo>,</mo>
<mo>(</mo>
<mi>i</mi>
<mo>=</mo>
<mn>0</mn>
<mo>...</mo>
<mi>M</mi>
<mo>-</mo>
<mn>2</mn>
<mo>,</mo>
<mi>j</mi>
<mo>=</mo>
<mn>0</mn>
<mo>...</mo>
<mi>N</mi>
<mo>-</mo>
<mn>1</mn>
<mo>)</mo>
</mtd>
</mtr>
<mtr>
<mtd>
<mn>0</mn>
<mo>,</mo>
<mi>o</mi>
<mi>t</mi>
<mi>h</mi>
<mi>e</mi>
<mi>r</mi>
<mi>w</mi>
<mi>i</mi>
<mi>s</mi>
<mi>e</mi>
</mtd>
</mtr>
</mtable>
</mfenced>
</mrow>
<mrow>
<msubsup>
<mi>&Delta;</mi>
<mi>&phi;</mi>
<mi>y</mi>
</msubsup>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mfenced open = "{" close = "">
<mtable>
<mtr>
<mtd>
<mi>&phi;</mi>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>+</mo>
<mn>1</mn>
<mo>)</mo>
<mo>-</mo>
<mi>&phi;</mi>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>)</mo>
<mo>,</mo>
<mo>(</mo>
<mi>i</mi>
<mo>=</mo>
<mn>0</mn>
<mo>...</mo>
<mi>M</mi>
<mo>-</mo>
<mn>1</mn>
<mo>,</mo>
<mi>j</mi>
<mo>=</mo>
<mn>0</mn>
<mo>...</mo>
<mi>N</mi>
<mo>-</mo>
<mn>2</mn>
<mo>)</mo>
</mtd>
</mtr>
<mtr>
<mtd>
<mn>0</mn>
<mo>,</mo>
<mi>o</mi>
<mi>t</mi>
<mi>h</mi>
<mi>e</mi>
<mi>r</mi>
<mi>w</mi>
<mi>i</mi>
<mi>s</mi>
<mi>e</mi>
</mtd>
</mtr>
</mtable>
</mfenced>
</mrow>
4.3) by the second order local derviation ρ of the true phase of+1 iterative calculation of kthk+1The Poisson's equation represented as the following formula
Input, the Poisson's equation that the following formula expression is then solved with discrete cosine transform (DCT) obtain the true of+1 iterative calculation of kth
Real phasek+1：
<mrow>
<mfrac>
<msup>
<mo>&part;</mo>
<mn>2</mn>
</msup>
<mrow>
<mo>&part;</mo>
<msup>
<mi>x</mi>
<mn>2</mn>
</msup>
</mrow>
</mfrac>
<msub>
<mi>&phi;</mi>
<mrow>
<mi>k</mi>
<mo>+</mo>
<mn>1</mn>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>)</mo>
</mrow>
<mo>+</mo>
<mfrac>
<msup>
<mo>&part;</mo>
<mn>2</mn>
</msup>
<mrow>
<mo>&part;</mo>
<msup>
<mi>y</mi>
<mn>2</mn>
</msup>
</mrow>
</mfrac>
<msub>
<mi>&phi;</mi>
<mrow>
<mi>k</mi>
<mo>+</mo>
<mn>1</mn>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<msub>
<mi>&rho;</mi>
<mrow>
<mi>k</mi>
<mo>+</mo>
<mn>1</mn>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>)</mo>
</mrow>
</mrow>
4.4) the true phase φ of+1 iterative calculation of kth is judgedk+1Whether the condition of convergence is met：
If satisfied, then with true phase φk+1As final true phase；
If not satisfied, then repeating the above steps 4.2)~4.3) carry out the true phase that processing obtains next iteration calculating.
6. a kind of weighted least-squares phase unwrapping computational methods based on reliability according to claim 5, its feature
It is：The condition of convergence such as the following formula represents, meets that the following formula then meets the condition of convergence：
<mrow>
<msqrt>
<mrow>
<mfrac>
<mn>1</mn>
<mrow>
<mi>M</mi>
<mo>&times;</mo>
<mi>N</mi>
</mrow>
</mfrac>
<munderover>
<mo>&Sigma;</mo>
<mrow>
<mi>i</mi>
<mo>=</mo>
<mn>0</mn>
</mrow>
<mrow>
<mi>M</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</munderover>
<munderover>
<mo>&Sigma;</mo>
<mrow>
<mi>j</mi>
<mo>=</mo>
<mn>0</mn>
</mrow>
<mrow>
<mi>N</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</munderover>
<msup>
<mrow>
<mo>&lsqb;</mo>
<msub>
<mi>&phi;</mi>
<mrow>
<mi>k</mi>
<mo>+</mo>
<mn>1</mn>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>)</mo>
</mrow>
<mo>-</mo>
<msub>
<mi>&phi;</mi>
<mi>k</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>)</mo>
</mrow>
<mo>&rsqb;</mo>
</mrow>
<mn>2</mn>
</msup>
</mrow>
</msqrt>
<mo>&le;</mo>
<mi>&epsiv;</mi>
</mrow>
Wherein, ε represents convergence threshold, ε=10-2；φk+1Represent the true phase of+1 iterative calculation of kth, φkRepresent kth time
The true phase of iterative calculation.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201711226954.4A CN107977939B (en) | 2017-11-29 | 2017-11-29 | Reliability-based weighted least square phase unwrapping calculation method |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201711226954.4A CN107977939B (en) | 2017-11-29 | 2017-11-29 | Reliability-based weighted least square phase unwrapping calculation method |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN107977939A true CN107977939A (en) | 2018-05-01 |
| CN107977939B CN107977939B (en) | 2021-11-16 |
Family
ID=62008573
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201711226954.4A Active CN107977939B (en) | 2017-11-29 | 2017-11-29 | Reliability-based weighted least square phase unwrapping calculation method |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN107977939B (en) |
Cited By (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN109472834A (en) * | 2018-10-23 | 2019-03-15 | 桂林电子科技大学 | A kind of Kalman filtering phase developing method based on wavelet transformation |
| CN112665529A (en) * | 2021-01-19 | 2021-04-16 | 浙江理工大学 | Object three-dimensional shape measuring method based on stripe density area segmentation and correction |
| CN112797917A (en) * | 2021-01-19 | 2021-05-14 | 浙江理工大学 | High-precision digital speckle interference phase quantitative measurement method |
Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN102721374A (en) * | 2012-06-04 | 2012-10-10 | 中国科学院光电技术研究所 | Planar sub-aperture splicing method based on weighted least square method |
| CN103063155A (en) * | 2012-12-12 | 2013-04-24 | 浙江师范大学 | Rapid package removal method of digital microscopic holographic phase diagram |
| CN104407332A (en) * | 2014-11-25 | 2015-03-11 | 沈阳建筑大学 | Correction method for updating DEM (digital elevation model) by aid of ground based SAR (synthetic aperture radar) |
| CN107014313A (en) * | 2017-05-16 | 2017-08-04 | 深圳大学 | The method and system of weighted least-squares phase unwrapping based on S-transformation ridge value |
-
2017
- 2017-11-29 CN CN201711226954.4A patent/CN107977939B/en active Active
Patent Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN102721374A (en) * | 2012-06-04 | 2012-10-10 | 中国科学院光电技术研究所 | Planar sub-aperture splicing method based on weighted least square method |
| CN103063155A (en) * | 2012-12-12 | 2013-04-24 | 浙江师范大学 | Rapid package removal method of digital microscopic holographic phase diagram |
| CN104407332A (en) * | 2014-11-25 | 2015-03-11 | 沈阳建筑大学 | Correction method for updating DEM (digital elevation model) by aid of ground based SAR (synthetic aperture radar) |
| CN107014313A (en) * | 2017-05-16 | 2017-08-04 | 深圳大学 | The method and system of weighted least-squares phase unwrapping based on S-transformation ridge value |
Non-Patent Citations (6)
| Title |
|---|
| DENNIS C. GHIGLIA等: ""Robust two-dimensional weighted and unweighted phase unwrapping that uses fast transforms and iterative methods"", 《J.OPT. SOC.AM.A》 * |
| MIGUEL AREVALLILO HERRA´EZ ET.AL.: ""Fast two-dimensional phase-unwrapping algorithm based on sorting by reliability following a noncontinuous path"", 《APPLIED OPTICS》 * |
| XIAN WANG, SUPING FANG, XINDONG ZHU: ""Weighted least-squares phase unwrapping algorithm based on a non-interfering image of an object"", 《APPLIED OPTICS》 * |
| 刘景峰等: ""含噪声包裹相位图的加权最小二乘相位展开算法研究"", 《光学技术》 * |
| 殷爱菡等: ""一种改进的加权最小二乘相位展开方法研究"", 《计算机工程与应用》 * |
| 黄文宇等: ""基于调制度分析的加权最小二乘位相展开方法"", 《北京理工大学学报》 * |
Cited By (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN109472834A (en) * | 2018-10-23 | 2019-03-15 | 桂林电子科技大学 | A kind of Kalman filtering phase developing method based on wavelet transformation |
| CN109472834B (en) * | 2018-10-23 | 2023-04-14 | 桂林电子科技大学 | Kalman filtering phase expansion method based on wavelet transformation |
| CN112665529A (en) * | 2021-01-19 | 2021-04-16 | 浙江理工大学 | Object three-dimensional shape measuring method based on stripe density area segmentation and correction |
| CN112797917A (en) * | 2021-01-19 | 2021-05-14 | 浙江理工大学 | High-precision digital speckle interference phase quantitative measurement method |
| CN112665529B (en) * | 2021-01-19 | 2022-06-24 | 浙江理工大学 | Object three-dimensional shape measuring method based on stripe density area segmentation and correction |
Also Published As
| Publication number | Publication date |
|---|---|
| CN107977939B (en) | 2021-11-16 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN104091064B (en) | PS-DInSAR ground surface deformation measurement parameter estimation method based on optimal solution space search method | |
| CN102854533B (en) | A kind of denoising method improving seismic data signal to noise ratio (S/N ratio) based on wave field separation principle | |
| CN104123464B (en) | Method for inversion of ground feature high elevation and number of land subsidence through high resolution InSAR timing sequence analysis | |
| CN104050681B (en) | A kind of road vanishing Point Detection Method method based on video image | |
| CN107977939A (en) | A kind of weighted least-squares phase unwrapping computational methods based on reliability | |
| CN102621549A (en) | Multi-baseline/multi-frequency-band interference phase unwrapping frequency domain quick algorithm | |
| CN103454636B (en) | Differential interferometric phase estimation method based on multi-pixel covariance matrixes | |
| CN106093939A (en) | A kind of InSAR image phase unwrapping method based on phase contrast statistical model | |
| CN104933673A (en) | Interference SAR (Synthetic Aperture Radar) image precise registration method based on resolution search sub-pixel offset | |
| CN108680901A (en) | A kind of novel sound bearing localization method | |
| CN102831588A (en) | De-noising processing method for three-dimensional seismic images | |
| CN109556797A (en) | The pipeline leakage detection and location method with convolutional neural networks is decomposed based on spline local mean value | |
| CN105469398A (en) | Deformation speckle generation method based on reverse mapping method | |
| CN103308914B (en) | One-station fixed bistatic interference synthetic aperture radar (SAR) processing method | |
| CN111580163A (en) | Full waveform inversion method and system based on non-monotonic search technology | |
| CN105809649A (en) | Variation multi-scale decomposing based SAR image and visible light image integration method | |
| CN111950438A (en) | Depth learning-based effective wave height inversion method for Tiangong No. two imaging altimeter | |
| CN106289051A (en) | The direction of big change density of electronic speckle interference fringe pattern and density processing method | |
| CN104730521A (en) | SBAS-DInSAR method based on nonlinear optimization strategy | |
| CN104076324A (en) | Method for estimating high-accuracy arrival direction without knowing information source number | |
| CN104331857A (en) | Phase position difference iteration compensation method in light intensity transmission equation phase retrieval | |
| CN107909606A (en) | A kind of SAR image registration communication center elimination of rough difference method | |
| CN104240230A (en) | Method for improving matching accuracy of phase correlation algorithm | |
| CN102903084B (en) | Wavelet field image noise variance method of estimation under a kind of α stable model | |
| CN104765063A (en) | Oil-gas detection method and device for calculating absorption attenuation attribute based on frequency spectrum |
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 |