在Python中舍入很小/接近零的复数值

Question

考虑以下序列/信号：

import numpy as np
from scipy.fft import fft

# %% Discrete data

x = np.array([0.7966875, 0.7966875, 0.7966875, 0.7966875, 0.7966875, 0.7966875,
       0.7966875, 0.7966875, 0.7966875, 0.7966875, 0.7966875, 0.7966875,
       0.7966875, 0.7966875, 0.7966875, 0.7966875, 0.7966875, 0.7966875,
       0.7966875, 0.7966875, 0.7966875, 0.7966875, 0.7966875, 0.7966875,
       0.7966875, 0.7966875, 0.7966875, 0.7966875, 0.7966875, 0.7966875,
       0.7966875, 0.7966875, 0.7966875, 0.7966875, 0.7966875, 0.7966875,
       0.7966875, 0.7966875, 0.7966875, 0.7966875, 0.7966875, 0.7966875,
       0.7966875, 0.7966875, 0.7966875, 0.7966875, 0.7966875, 0.7966875])

我定义了以下函数来计算离散傅立叶变换：

def DFT(x):
    N = len(x)
    n = np.arange(0, N)
    k = n.reshape((N, 1))
    omega = 2*np.pi*k/N
    e = np.exp(-1j*omega*n)
    X = np.dot(x, e)
    X_round = np.round(X, 6)
    return X, X_round, N, n

X, X_round, N, n = DFT(x)

我们的

X

对于非零频率具有非常小的值（正如预期的那样，信号

x

不会改变）。然后，我尝试使用

np.round(X, 6)

将其四舍五入，返回：

array([38.241+0.j, -0.   +0.j, -0.   -0.j, -0.   -0.j,  0.   -0.j,
        0.   +0.j, -0.   -0.j,  0.   -0.j, -0.   +0.j, -0.   -0.j,
       -0.   -0.j, -0.   -0.j, -0.   +0.j,  0.   -0.j, -0.   -0.j,
       -0.   -0.j, -0.   -0.j,  0.   -0.j,  0.   -0.j,  0.   -0.j,
        0.   +0.j,  0.   +0.j, -0.   +0.j, -0.   +0.j,  0.   +0.j,
       -0.   -0.j,  0.   +0.j, -0.   +0.j, -0.   -0.j, -0.   -0.j,
       -0.   -0.j,  0.   -0.j,  0.   -0.j,  0.   -0.j, -0.   -0.j,
       -0.   -0.j,  0.   -0.j,  0.   -0.j,  0.   -0.j,  0.   +0.j,
       -0.   +0.j, -0.   -0.j,  0.   -0.j, -0.   +0.j, -0.   -0.j,
        0.   +0.j,  0.   -0.j,  0.   -0.j])

但是，如果我尝试使用

X_round

获取

X_round_angle = np.angle(X_round)

的角度，它不会像我预期的那样为零（您可以尝试使用 scipy 中的

fft

函数）。这是

X_round_angle

的输出：

array([ 0.        ,  3.14159265, -3.14159265, -3.14159265, -0.        ,
        0.        , -3.14159265, -0.        ,  3.14159265, -3.14159265,
       -3.14159265, -3.14159265,  3.14159265, -0.        , -3.14159265,
       -3.14159265, -3.14159265, -0.        , -0.        , -0.        ,
        0.        ,  0.        ,  3.14159265,  3.14159265,  0.        ,
       -3.14159265,  0.        ,  3.14159265, -3.14159265, -3.14159265,
       -3.14159265, -0.        , -0.        , -0.        , -3.14159265,
       -3.14159265, -0.        , -0.        , -0.        ,  0.        ,
        3.14159265, -3.14159265, -0.        ,  3.14159265, -3.14159265,
        0.        , -0.        , -0.        ])

使用

np.angle(fft(x))

返回所有频率的正确相位值为零。我怎样才能纠正我的

DFT

功能？

Answer 1

您观察到的行为并非完全不合理，因为无法定义复数 0 的相位角。

np.angle()

为您提供 0、pi 和 -pi 的混合，所有这些都同样正确。
np.angle()的

文档

这样说：

虽然复数0的角度没有定义， numpy.angle(0) 返回值 0。

如果我们完全按照说明进行测试，我们会得到预期结果 0。但是，当我们将参数

z

从 0（整数）调整到 0.0（浮点数）时，您所看到的行为就会发生。从技术上讲，这是一个未记录的边缘情况，因此我们不能指望那里有任何一致性。我们可以测试正负 0.0、有或没有 0.0j 的不同组合，结果是好坏参半。
从数学角度来说，这都是合理的，但是您可能希望相位对于所有类型的零都保持一致（例如，用于绘图），这肯定是有原因的。在这种情况下，我只需使用推导式将有问题的条目替换为 0。

phase = np.angle(X_round)
phase = np.array([0.0 if np.abs(a)<0.000001 else theta for a, theta in zip(X_round, phase)])