OE Letters

Unfolding wrapped phase

[+] Author Affiliations
Carlos Gerardo Treviño-Palacios

Instituto Nacional de Astrofísica, Óptica y Electrónica, Luis Enrique Erro 1, Sta Ma Tonantzintla, Puebla 72840, Mexico

Opt. Eng. 54(11), 110503 (Nov 20, 2015). doi:10.1117/1.OE.54.11.110503
History: Received August 11, 2015; Accepted October 16, 2015
Text Size: A A A

Open Access Open Access

Abstract.  Phase unwrapping is the final step in phase extraction methods, which consists of recovering the correct phase from the wrapped phase by removing 2π discontinuities. The difference between the correct phase and the wrapped phase is the phase wrapping map. A new method for phase unwrapping is presented by identifying the phase wrapping map as a sequence of binary valued intermediate wrapping maps and iteratively removing them producing the correct phase by phase-wrapped unfolding. A path-following algorithm is presented to exemplify the phase wrapped unfolding method.

Figures in this Article

Phase extraction methods are used in different coherent signal processing areas, such as optical and microwave interferometry,1 magnetic resonance imaging,2 synthetic aperture radar,3 and Fourier transform optical coherent tomography,4 among others. In these areas, the signal modulates a harmonic function and phase extraction algorithms are used to recover the signal. The recovered signal within the principal branch of the harmonic function exhibits a series of discontinuities, which need to be removed to obtain the real signal. This process, known as phase unwrapping, is the final phase extraction step. We collect these discontinuities in a wrapping map, normally overlooked based on the principle that the proposed methods aim to recover the real phase map knowing its expected continuity properties. Since phase unwrapping was first showed in optical interferometry phase extraction,5,6 there have been extensive studies to recover the correct phase map. Because of the presence of noise, signal discontinuities, and undersampling, phase unwrapping is considered to be one of the most difficult problems in mathematics and engineering.7 To date, there are two main approaches to the phase unwrapping problem, which can be classified in path following methods of minimum-norm methods.810 Independent of the approach, a definite solution for phase unwrapping remains an open problem.

On this note, an alternative solution to the phase unwrapping problem focusing on the properties of the phase wrapping map instead of the real phase signal properties is presented.

Phase extraction methods retrieve a wrapped phase from a harmonic function confined in the main branch of the function, which is normally a 2π band in either [0,2π) or [π,π). Because phase is a relative quantity, we can focus only on the [0,2π) band without loss of generality. In order to obtain an insight on the phase unwrapping problem instead of looking at the problem as the final step in a phase extraction methods as is customary, let us look at the problem from the solution point of view. In other words, starting from the correct phase map U [Fig. 1(a)], the wrapped phase Uw is obtained [Fig. 1(b)] by adding multiples of 2π [Fig. 1(c)] to contain Uw in the [0,2π) band. The form of the integer-valued correcting field [Fig. 1(c)], or phase wrapping map M, is dependent on the properties of U. Display Formula

Uw=U2πM.(1)

Graphic Jump Location
Fig. 1
F1 :

For a 15π amplitude Gaussian curve (a) genuine unwrapped phase U, (b) corresponding wrapped phase Uw within [0,2π) values, and (c) integer-valued 2π multiple correcting field, or phase wrapping map 2πM.

The aim of the phase unwrapping methods aforementioned is to add a correct multiple of 2π to the wrapped phase such that a continuous unwrapped phase map is reconstructed. These methods focus on the properties of the expected phase map U and the challenge is to restore a function free of the 2π jumps.

As mentioned, instead of focusing in the properties of U, here we investigate the properties of the phase wrapping map M. To be more precise, focus on an alternative construction of M. To do such and bearing in mind that if the phase map U is sampled with scaled harmonic functions with a multiple of 2π periodicity, in principle, we can obtain the correct solution after an adequate phase unwrapping procedure. Considering using powers of two multiples of 2π (2nπ) is of particular interest. By using this approach, an insight on the phase wrapping map is obtained and it is found that it can be properly decomposed by simpler phase wrapping maps.

The first task is to obtain the minimum power of two multiples of 2π, which contains U. This is, look for N, which is the minimum integer such that Display Formula

|U|2Nπ.(2)

After determining this number, consider wrapping planes in subpowers of two multiples of 2π:2N1π, 2N2π,,4π, 2π (Fig. 2) in such a way that iterative wrapping maps M1,M2,,MN1, MN [Figs. 3(b), 3(d), 3(f), and 3(h)] are obtained by wrapping the immediate previous wrapped function (U+M1·2N1π+M2·2N2π++Mn·2Nnπ) by the next subpower of two multiples of 2π. By doing so the n’th wrapping map has binary values: zero or 2Nn+1π. Notice that by construction, the position, but not necessarily the sign, of the phase jumps for each wrapping map is contained in the next one. Thus, the following intermediate wrapped functions are constructed [Figs. 3(a), 3(c), 3(e), and 3(g)]: Display Formula

{Uw,0=U;|Uw,0|2NπUw,1=Uw,0M1·2N1π;|Uw,1|=2N1πUw,N2=UMN2·4π=Uw,N3MN2·4π;|Uw,N2|=4πUw=UMN1·2π=Uw,N2MN1·2π;|Uw|=2π.(3)

Graphic Jump Location
Fig. 2
F2 :

Alternative wrapping planes separated in a 2Nnπ sequence for the function on Fig. 1. In this case, N=4 and planes for M1,M2, and M3 are considered.

Graphic Jump Location
Fig. 3
F3 :

Intermediate wrapped functions and intermediate wrapped maps for the example in Fig. 1: (a) Uw,0=U, (b) 8πM1 used to fold Uw,0, (c) Uw,1=U08πM1, (d) 4πM2 used to fold Uw,1, (e) Uw,2=Uw,14πM2, (f) 2πM3 used to fold Uw,2, (g) Uw,3=Uw,22πM3; this is the final wrapped phase Uw within [0,2π) values, and (h) πM3 used as constructing wrapped function in the phase unfolding algorithm.

The wrapping map down to the main branch (2π) is therefore the sum of the binary intermediate wrapping maps. Display Formula

M=i=1N1Mi·2Nnπ.(4)

Using this approach instead of finding the correct multiple of 2π for phase unwrapping, the solution is reduced to reconstruct each step backward. In other words, instead of unwrapping the phase, we iteratively unfold the phase. Display Formula

Uw,n=Uw,n1ϕnMn,(5)
where ϕn=2Nnπ. Unfolding the phase consists on finding and removing iteratively Mn. The benefit is that the partial wrapping map Mn consists of only zeroes and ones.

In spite of the apparent advantage using this binary approach, the unfolding problem is as complex as the unwrapping approach. This complexity is made evident by manipulating an intermediate wrapped function Uw,n placed on top of itself, thus mimicking one unfolding step. Using this construction, it is found that each point has two possible solutions (Fig. 4): Display Formula

Uw,n1I=Uw,n1,(6)
Display Formula
Uw,n1II=Uw,n12Mnϕn+ϕn.(7)
Equation (6) is the correct solution [Fig. 4(b) solid line], whereas Eq. (7) is a false solution, which increases the amount of wrapping [Fig. 4(b) dotted line]. Nevertheless, the correct solution is found among only two possible solutions and phase unwrapping has been reduced to find the correct one among these two.

Graphic Jump Location
Fig. 4
F4 :

Graphical principle of the phase unfold method. (a) Original Uw,n intermediate wrapped phase and (b) Uw,n duplicated and placed on top of itself showing the correct solution Uw,n1 (solid line) and false solution in Eq. (7) (dotted line).

Different approaches for unfolding the phase have been explored. Here, a straightforward path following method is presented. In order to implement the basic path following phase unfolding process, notice that if the intermediate phase wrapped function Un is shifted by ϕn/2, the following folding step produces the same wrapped function Un+1. This property is used to simplify the unfolding process. Thus, the path following phase unfolding process can be performed by splitting the task down into the following steps:

  1. Start with an intermediate wrapped function Uw,n=Uw,n1+Mnϕn and roll the phase within its [0,ϕn] domain by ϕn/2. This is done by adding ϕn/2 to values lower than half the step amplitude (ϕn) and subtracting ϕn/2 to the upper half values [Figs. 5(a) and 5(b)]. The resulting intermediate wrapped phase is Uw,n1+Onϕn, where On is the correspondent phase wrapping map for the ϕn/2 shifted Uw,n function.
    • for the first iteration, use ϕN1=2π and Uw,n=Uw.
  2. Subtract the previous two functions [Fig. 5(c)]: Display Formula
    (Uw,n1Onϕn)(Uw,n1Mnϕn)=(MnOn)ϕn=Mn+1ϕnϕn/2(8)
    This is a binary function, which by construction contains the position of the unknown phase jumps of Mn. It also contains intercalated the phase jumps of On
  3. Initialize an unfolding variable Sw to zero.
  4. Starting from the first sample of the left, calculate the difference in the function resulting from Eq. (8) between the current sample (j) and its direct adjacent right-hand neighbor (j+1). Only pay attention to nonzero differences.
    • For the first nonzero value, which can be either positive or negative, assign 1 or 1 to an unfolding sign variable sgnw, correspondingly.
  5. For the nonzero difference points identified in step 4, if the difference between adjacent values of Un added to Sw is larger than ϕn/2{Un(j)Un(j+1)+Sw>ϕn/2?}:
    • if Sw is zero, assign sgnwφn to it.
    • if Sw is nonzero, assign zero to it.
    This step selects between the two possible paths Mn or On.
  6. In parallel to steps 4 and 5, starting from the first sample on the left, add Sw to Uw,n (Uw,n1=Uw,n+Sw).
  7. Continue steps 4, 5, and 6 to the right-hand samples until all the samples have been processed. This concludes an unfolding step and Uw,n1 has been obtained.
  8. Double the value of ϕn (ϕn1=2ϕn) and repeat steps 1 through 7 using the obtained Uw,n1.
  9. Continue until the difference in Eq. (8) is constant and/or smaller than a given threshold (ϕ/4 for instance) and exit in step 2.

Graphic Jump Location
Fig. 5
F5 :

(a) Intermediate wrapped function Uw,n=Uw,n1Mnϕn, the arrows show the required displacements to roll the phase by ϕn/2, as described in the text, (b) complementary intermediate wrapped function after a ϕn/2 roll, Uw,nOnϕn and (c) difference between Mnϕn and Onϕn signifying Mn+1ϕn.

Using this procedure after N1 iterations, in which the final value of ϕ is 2Nπ, the correct phase map U is recovered. The algorithm performs the unfolding procedure by locating phase jumps on a binary function Mn+1 governed by Eq. (8). Therefore, in contrast to regular phase unwrapping techniques instead of searching for discontinuities on a stepwise continuous function Un at each point to execute the unwrapping decision like ice walking the unfolding decision is made only on abrupt changes in Mn+1 in a cliff hopping manner.

Phase unfolding produces N1 times more operations than phase unwrapping methods but has the advantage that the number of discontinuity correcting decisions is reduced to the phase jumps clearly identified on Mn+1 [Eq. (8)], which decrease on each iteration. Phase unfolding also has the same restrictions than conventional phase unwrapping: is limited by signal noise, under sampling beyond Shannon–Nyquist limit and discontinuities in the wrapped signal. Furthermore, the first iteration in the unfolding algorithm requires an additional downfolding the wrapped phase to half of the main branch [Fig. 3(h)], thus increasing the sampling limitations to fulfill the Shannon–Nyquist limit. On the other hand, errors due to false wrapping do not propagate as it does in the path following phase unwrapping because of the local phase jump detection in Eq. (7).

As an example of the phase unfolding algorithm, a unitary circle function closely related to third order spherical aberration with 15π amplitude is analyzed (Fig. 6). Figure 6(a) shows the expected result of a phase extraction method confined in the [0,2π) band. In this case, N is 4 (24π>15π>23π), therefore unfolding phase unwrapping is achieved in three steps as shown. On each consecutive step [Figs. 6(b)6(d)], the phase unfolds as a folding cup opening until the correct phase map U is recovered in Fig. 6(d).

Graphic Jump Location
Fig. 6
F6 :

Phase unfolding algorithm operating on the unitary circle function 10π {6(x2+y2)26(x2+y2)}. Only half of the function is shown to better appreciate the effect of the phase unfolding algorithm. (a) Original wrapped function Uw=Uw,3 in the [0,2π) domain, (b) first unfolding iteration Uw,2 in the [0,4π) domain, (c) second iteration Uw,1 in the [0,8π) domain, and (d) third and final iteration recovering the correct phase map Uw,0=U with a 15π amplitude.

For the simple path following the method demonstrated, phase unfolding presents some advantages with respect to phase unwrapping, namely the correct identification of phase jumps, a decrease number of discontinuity correcting decisions and splitting the problem to binary maps [Eq. (5)]. It also presents disadvantages such as tighter limitation in the sampling period and the number of computational steps. Both approaches are similarly sensitive to signal noise or discontinuities in the signal. Comparable computational time has been observed using either approach. Further research is required to incorporate simultaneously the properties of both the real phase map U and the phase wrapping map M in phase unwrapping or unfolding algorithms. Additional studies also need to be conducted to fully determine the effect of the wrapping map decomposition [Eqs. (4) and (5)] on either path following or minimum norm unwrapping methods.

In conclusion, an alternative method for phase unwrapping, named phase unfolding, has been presented by exploring the properties of the phase wrapping map decomposed in a series of binary intermediate wrapping maps. An effective iterative algorithm for path following phase unfolding was also presented, identifying some advantages and disadvantages with respect to standard phase unwrapping.

The author would like to thank Dr. Blas M. Rodriguez-Lara for fruitful discussions. This work was partially funded by grant CONACyT CB-2011-01-168558.

Schofield  M. A., and Zhu  Y., “Fast phase unwrapping algorithm for interferometric applications,” Opt. Lett.. 28, (14 ), 1194 –1196 (2003). 0146-9592 CrossRef
Ying  L., “Phase unwrapping,” in Wiley Encyclopedia of Biomedical Engineering. , and Akay  M., Ed.,  John Wiley & Sons ,  Hoboken, New Jersey  (2006).
Osmanoglu  B.  et al., “On the importance of path for phase unwrapping in synthetic aperture radar interferometry,” Appl. Opt.. 50, (19 ), 3205 –3220 (2011). 0003-6935 CrossRef
Nassif  N. A.  et al., “In vivo high-resolution video-rate spectral-domain optical coherence tomography of the human retina and optic nerve,” Opt. Express. 12, (3 ), 367 –376 (2004). 1094-4087 CrossRef
Takeda  M., , Ina  H., and Kobayashi  S., “Fourier-transform method of fringe-pattern analysis for computer-based topography and interferometry,” J. Opt. Soc. Am.. 72, (1 ), 156 –160 (1982). 0030-3941 CrossRef
Itoh  K., “Analysis of the phase unwrapping algorithm,” Appl. Opt.. 21, (14 ), 2470 –2470 (1982). 0003-6935 CrossRef
Ghiglia  D., and Pritt  M., Two-Dimensional Phase Unwrapping Theory, Algorithms and Software. ,  John Wiley & Sons ,  Hoboken, New Jersey  (1998).
Volkov  V. V., and Zhu  Y., “Deterministic phase unwrapping in the presence of noise,” Opt. Lett.. 28, (22 ), 2156 –2158 (2003). 0146-9592 CrossRef
Servin  M.  et al., “Phase unwrapping with a regularized phase-tracking system,” Appl. Opt.. 37, (10 ), 1917 –1923 (1998). 0003-6935 CrossRef
Huang  H. Y. H.  et.al., “Path-independent phase unwrapping using phase gradient and total-variation (TV) denoising,” Opt. Express. 20, (13 ), 14075 –14089 (2012). 1094-4087 CrossRef
© The Authors. Published by SPIE under a Creative Commons Attribution 3.0 Unported License. Distribution or reproduction of this work in whole or in part requires full attribution of the original publication, including its DOI.

Citation

Carlos Gerardo Treviño-Palacios
"Unfolding wrapped phase", Opt. Eng. 54(11), 110503 (Nov 20, 2015). ; http://dx.doi.org/10.1117/1.OE.54.11.110503


Figures

Graphic Jump Location
Fig. 1
F1 :

For a 15π amplitude Gaussian curve (a) genuine unwrapped phase U, (b) corresponding wrapped phase Uw within [0,2π) values, and (c) integer-valued 2π multiple correcting field, or phase wrapping map 2πM.

Graphic Jump Location
Fig. 2
F2 :

Alternative wrapping planes separated in a 2Nnπ sequence for the function on Fig. 1. In this case, N=4 and planes for M1,M2, and M3 are considered.

Graphic Jump Location
Fig. 3
F3 :

Intermediate wrapped functions and intermediate wrapped maps for the example in Fig. 1: (a) Uw,0=U, (b) 8πM1 used to fold Uw,0, (c) Uw,1=U08πM1, (d) 4πM2 used to fold Uw,1, (e) Uw,2=Uw,14πM2, (f) 2πM3 used to fold Uw,2, (g) Uw,3=Uw,22πM3; this is the final wrapped phase Uw within [0,2π) values, and (h) πM3 used as constructing wrapped function in the phase unfolding algorithm.

Graphic Jump Location
Fig. 4
F4 :

Graphical principle of the phase unfold method. (a) Original Uw,n intermediate wrapped phase and (b) Uw,n duplicated and placed on top of itself showing the correct solution Uw,n1 (solid line) and false solution in Eq. (7) (dotted line).

Graphic Jump Location
Fig. 5
F5 :

(a) Intermediate wrapped function Uw,n=Uw,n1Mnϕn, the arrows show the required displacements to roll the phase by ϕn/2, as described in the text, (b) complementary intermediate wrapped function after a ϕn/2 roll, Uw,nOnϕn and (c) difference between Mnϕn and Onϕn signifying Mn+1ϕn.

Graphic Jump Location
Fig. 6
F6 :

Phase unfolding algorithm operating on the unitary circle function 10π {6(x2+y2)26(x2+y2)}. Only half of the function is shown to better appreciate the effect of the phase unfolding algorithm. (a) Original wrapped function Uw=Uw,3 in the [0,2π) domain, (b) first unfolding iteration Uw,2 in the [0,4π) domain, (c) second iteration Uw,1 in the [0,8π) domain, and (d) third and final iteration recovering the correct phase map Uw,0=U with a 15π amplitude.

Tables

References

Schofield  M. A., and Zhu  Y., “Fast phase unwrapping algorithm for interferometric applications,” Opt. Lett.. 28, (14 ), 1194 –1196 (2003). 0146-9592 CrossRef
Ying  L., “Phase unwrapping,” in Wiley Encyclopedia of Biomedical Engineering. , and Akay  M., Ed.,  John Wiley & Sons ,  Hoboken, New Jersey  (2006).
Osmanoglu  B.  et al., “On the importance of path for phase unwrapping in synthetic aperture radar interferometry,” Appl. Opt.. 50, (19 ), 3205 –3220 (2011). 0003-6935 CrossRef
Nassif  N. A.  et al., “In vivo high-resolution video-rate spectral-domain optical coherence tomography of the human retina and optic nerve,” Opt. Express. 12, (3 ), 367 –376 (2004). 1094-4087 CrossRef
Takeda  M., , Ina  H., and Kobayashi  S., “Fourier-transform method of fringe-pattern analysis for computer-based topography and interferometry,” J. Opt. Soc. Am.. 72, (1 ), 156 –160 (1982). 0030-3941 CrossRef
Itoh  K., “Analysis of the phase unwrapping algorithm,” Appl. Opt.. 21, (14 ), 2470 –2470 (1982). 0003-6935 CrossRef
Ghiglia  D., and Pritt  M., Two-Dimensional Phase Unwrapping Theory, Algorithms and Software. ,  John Wiley & Sons ,  Hoboken, New Jersey  (1998).
Volkov  V. V., and Zhu  Y., “Deterministic phase unwrapping in the presence of noise,” Opt. Lett.. 28, (22 ), 2156 –2158 (2003). 0146-9592 CrossRef
Servin  M.  et al., “Phase unwrapping with a regularized phase-tracking system,” Appl. Opt.. 37, (10 ), 1917 –1923 (1998). 0003-6935 CrossRef
Huang  H. Y. H.  et.al., “Path-independent phase unwrapping using phase gradient and total-variation (TV) denoising,” Opt. Express. 20, (13 ), 14075 –14089 (2012). 1094-4087 CrossRef

Some tools below are only available to our subscribers or users with an online account.

Related Content

Customize your page view by dragging & repositioning the boxes below.

Related Book Chapters

Topic Collections

PubMed Articles
Advertisement
  • Don't have an account?
  • Subscribe to the SPIE Digital Library
  • Create a FREE account to sign up for Digital Library content alerts and gain access to institutional subscriptions remotely.
Access This Article
Sign in or Create a personal account to Buy this article ($20 for members, $25 for non-members).
Access This Proceeding
Sign in or Create a personal account to Buy this article ($15 for members, $18 for non-members).
Access This Chapter

Access to SPIE eBooks is limited to subscribing institutions and is not available as part of a personal subscription. Print or electronic versions of individual SPIE books may be purchased via SPIE.org.