0
$\begingroup$

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!

(left) input wavefield (right) Hilbert transformed wavefield (left) instantaneous phase phi (right) derived propagation angle using the gradient of the instantaneous phase

$\endgroup$

1 Answer 1

1
$\begingroup$

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.

$\endgroup$
5
  • $\begingroup$ 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? $\endgroup$ Sep 23, 2014 at 13:04
  • $\begingroup$ does the time derivative get computated accordingly? $\endgroup$ Sep 23, 2014 at 15:09
  • $\begingroup$ @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$. $\endgroup$
    – Matt L.
    Sep 23, 2014 at 15:30
  • $\begingroup$ Thanks a lot - it works! Do you have a reference for this alternative way to compute $\phi'$? $\endgroup$ Sep 24, 2014 at 16:12
  • $\begingroup$ @MichaelScott: You can find more information in this paper. $\endgroup$
    – Matt L.
    Sep 24, 2014 at 17:32

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.