Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

Python numpy unwrap function

I am hoping to convert a array of radians into range [0, 2*pi) and numpy unwrap function is exactly what I need

However, when I run the following code to input a = [pi, 2*pi, 3*pi]:

import numpy as np

a = np.array([np.pi, 2*np.pi, 3*np.pi])
np.unwrap(a)

I expect the results to be close to [pi, 0, pi]. However, the output is still:

array([ 3.14159265,  6.28318531,  9.42477796])

It is not unwrapped. However, if I instead run the following without using the numpy.pi

a = np.array([3.14159265,  6.28318531,  9.42477796])
np.unwrap(a)

The output is correct:

array([  3.14159265e+00,   2.82041412e-09,   3.14159265e+00])

what is going on?

like image 985
SunnyIsaLearner Avatar asked Jun 20 '17 22:06

SunnyIsaLearner


2 Answers

Although the accepted answer gives you the result you want, I don't think it gets to the heart of the problem which, if I'm interpreting your question correctly, is that you actually want wrap your phase, not unwrap it.

The reason np.unwrap works in this case, with small changes to your data, is actually a consequence of the naive way that np.unwrap computes its result; it simply looks for local discontinuities in your data and adjusts accordingly. Getting the result you're looking for in this way is a result of sampling errors. In other words, if you improve your sampling by interpolating to get a = np.array([np.pi, 3*np.pi/2, 2*np.pi, 5*np.pi/2, 3*np.pi]), adjusting the data won't work any more.

A more sophisticated method of phase unwrapping, such as a Fourier transform method, will leave your data unwrapped, even if the sampling is poor.

If you really want to constrain your data to [0, 2*pi), np.unwrap is the inverse of what you want. The simplest way I can think of to wrap your phase is with the modulo operator:

import numpy as np

a = np.array([np.pi, 2 * np.pi, 3 * np.pi])
a_wrapped = a % (2 * np.pi)
print (a_wrapped)

Of course, because of the sampling errors, np.unwrap(a_wrapped) does not return your original a, so it may not be clear that this is the inverse. However, if you improve your sampling, it does indeed return the original a:

import numpy as np

a = np.arange(0, 4 * np.pi, np.pi/10)
print (a)
a_wrapped = a % (2 * np.pi)
print (a_wrapped)
a = np.unwrap(a_wrapped)
print (a)
like image 75
Zac-K Avatar answered Sep 28 '22 10:09

Zac-K


Taken from np.unwrap docs:

Unwrap radian phase p by changing absolute jumps greater than discont to their 2*pi complement along the given axis.

Where discont = np.pi (by default). When

a = np.array([np.pi, 2*np.pi, 3*np.pi])

The jumps a[1] - a[0] = np.pi and a[2] - a[1] = np.pi are not greater than np.pi and therefore not 'unwraped'. However, if

a = np.array([3.14159265,  6.28318531,  9.42477796])

you have a[1] - a[0] = 3.1415926600000001 greater than np.pi, thus the function unwraps the values.

like image 38
rodgdor Avatar answered Sep 28 '22 09:09

rodgdor