# Instantaneous phase calculation problems

I am currently struggling with a problem when calculating the instantaneous phase of a wavefield using a 1D Hilbert transformation.

The scheme is as follows.

1D Hilbert transformation of wave field $U(x,z)$ $$q_z(x,z) = H_z[U(x,z)] = \int_R U(x-\xi,z)\frac{d\xi}{\xi}$$

Calculation of the instantaneous phase $\phi(x,z)$

$$\phi(x,z) = arctan\left(\frac{q_z(x,z)}{U(x,z}\right),$$ where $x$ and $z$ are the coordinates and $H_z$ is the 1D Hilbert transformation over the vertical direction.

My problem is, that $\phi(x,z)$ which is implemented using the intrinsic function $atan2()$ contains phase shifts [$\pi \rightarrow -\pi$].

Using the instantaneous phase for further calculations in which its gradient is derived hence leads to imaging artifacts as shown in the attached pictures:

TOP LEFT: input wavefield $U(x,z)$

TOP RIGHT: Hilbert transformed wavefield $q_z(x,z) = H_z[U(x,z)]$

BOTTOM LEFT: instantaneous phase $\phi(x,z)$

BOTTOM RIGHT: derived propagation angle with artifacts

.

My question: I am able to overcome this problem using Gauss filtering or Median filtering. Yet this does not wotrk properly for more complicated wavefields and higher frequencies. Is there another way to overcome those sharp transitions in the area of my wavefield? (not the noise around)

Thanks alot!

One way to overcome this problem is phase unwrapping, which is also implemented by the Matlab function unwrap. But I believe that the best option for obtaining the derivative of the phase is to avoid a direct computation, but instead start out from the complex function itself. To simplify notation, let's use $U$ (without arguments):

$$U=re^{j\phi}=a+ib$$ where $r$ is the magnitude and $\phi$ is the phase. The derivative of $U$ is given by

$$U'=r'e^{j\phi}+rj\phi'e^{j\phi}\\ \frac{U'}{U}=\frac{r'}{r}+j\phi'$$

Since $r'/r$ is real-valued, we get

$$\phi'=\text{Im}\left\{\frac{U'}{U}\right\}=\frac{ab'-a'b}{a^2+b^2}\tag{1}$$

where $a$ and $b$ are the real and imaginary parts of $U$, respectively. Using (1) instead of a direct derivative of the phase should solve the problem.

• Thanks for your reply @Matt! Just to get this right: Since I deal with a 2D problem $\phi '$ has to be derived for both $x$ and $z$ right? And $b$ would be my Hilbert transform of a? – MichaelScott Sep 23 '14 at 13:04
• does the time derivative get computated accordingly? – MichaelScott Sep 23 '14 at 15:09
• @MichaelScott: Yes, you can replace the derivative sign with all derivatives that you need. And $b$ is indeed the imaginary part of the analytic signal, so it's the Hilbert transform of $a$. – Matt L. Sep 23 '14 at 15:30
• Thanks a lot - it works! Do you have a reference for this alternative way to compute $\phi'$? – MichaelScott Sep 24 '14 at 16:12
• @MichaelScott: You can find more information in this paper. – Matt L. Sep 24 '14 at 17:32