如何有效地计算（基本上）等于零的复数的相位角？

Question

我正在编写一个 C++ 程序，它采用包含 double 值的真实输入信号的 FFT 和 return 包含 std::complex<double> 值的向量 X。一旦我有了结果向量，我就会尝试计算结果的幅度和相位。

当其中一个输出为 "zero" 时，我运行遇到了计算相位角的问题。零在引号中是因为当计算结果为 0 return 时，returned 值将非常接近零，但不完全为零。

例如，在索引 3 处，我的输出数组具有计算出的 "zero" 值：

X[3] = 3.0531133177191805e-16 - i*5.5511151231257827e-17

我正在尝试使用标准库 std::arg 函数，该函数应该 return 复数的相角。 std::arg(X[3])

虽然 X[3] 本质上是 0，但它不完全是 0，并且计算相位的方式会导致问题，因为计算使用的是虚部除以实部的比率，即远离 0！

进行实际计算得出的结果远非理想。

我怎样才能让 C++ 意识到结果真的是 0，这样我才能得到正确的相位角？

我正在寻找一种比使用任意硬编码 "epsilon" 值来比较 double 更优雅的解决方案，但到目前为止在网上搜索我还没有找到更好的解决方案.

Answer 1

如果您正在计算输入信号的浮点 FFT，则该信号将包含噪声，因此具有信噪比，包括传感器噪声，热噪声、量化噪声、定时抖动噪声等

因此，将 FFT 结果丢弃为低于本底噪声的阈值很可能不是计算数学问题，而是物理或电子数据采集分析的一部分。您必须插入该数字，并将相位设置为 0.0 或 NaN 或任何您的默认标记值是无用的（处于或低于本底噪声）FFT 结果。

Answer 2

从复数到幅度和相位的映射在 0 处不连续。

这是由于您选择的坐标不连续造成的。

解决方案将取决于您在可能出现不连续点附近的值的情况下选择这些坐标的原因。

它不是 "really" 零。如果你正确地考虑了误差线，你的答案真的会是一个小量级（希望如此）和一个不受约束的角度。

Answer 3

我注意到，当 FFT 的输入被缩放后，我原来的答案将不起作用。我相信我现在有一个实际有效的解决方案... 原始答案保留在下面，以便评论仍然有意义。

从对这个答案和其他答案的评论中，我了解到用语言计算精确的舍入误差在技术上是可行的，但绝对不切实际。最实用的解决方案似乎是允许用户提供他们自己的噪声阈值（以 dB 为单位）并忽略任何功率水平低于该阈值的数据点。不可能针对所有情况提出通用阈值，但用户可以根据被分析信号的 signal-to-noise 比率提供合理的阈值并将其传入。

下面显示了一个通用相位计算函数，它计算复杂数据点向量的相位角。

std::vector<double> Phase(std::vector<std::complex<double>> X, double threshold, double amplitude)
{
    size_t N = X.size();
    std::vector<double> X_phase(N);

    std::transform(X.begin(), X.end(), X_phase.begin(), [threshold, amplitude](const std::complex<double>& value) {
        double level = 10.0 * std::log10(std::norm(value) / std::pow(amplitude, 2.0));
        return level > threshold ? std::arg(value) : 0.0;
    });

    return X_phase;
}

这个函数有 3 个参数：

要计算其相位的复杂信号数据向量。
一个合理的阈值 -- 可以根据用于捕获信号的任何测量设备的 signal-to-noise 比率计算得出。如果您的信号除了语言本身的舍入误差外不包含任何噪声，您可以将其设置为任意非常低的值，例如 -120dB。
输入信号的最大可能振幅。如果你的信号是计算出来的，这应该简单地设置为你的信号的幅度。如果您的信号被测量，则应将其设置为测量设备能够测量的最大幅度（如果您的信号来自读取音频文件，通常其数据将在 -1.0 和 1.0 之间归一化。在这种情况下，您只需将振幅值设置为 1.0)。

这个新实现仍然为我提供了正确的结果，但更加稳健。通过将阈值计算留给用户，他们可以根据用于测量其输入信号的测量设备的特性自行设置最合理的值。

如果您发现任何错误或任何我可以进一步改进设计的方法，请告诉我！

原答案

我找到了一个看起来足够通用的解决方案。

在#include <limits>header中，std::numeric_limits<double>::digits10有一个常量值。

根据文档：

The value of std::numeric_limits<T>::digits10 is the number of base-10 digits that can be represented by the type T without change, that is, any number with this many significant decimal digits can be converted to a value of type T and back to decimal form, without change due to rounding or overflow.

使用它我可以过滤掉幅度低于此限制的任何输出值：

计算X[3]的相位：

int N = X.size();

auto tmp = std::abs(X[3])/N > std::pow(10, -std::numeric_limits<double>::digits10)
         ? value
         : 0.0

double phase = std::arg(tmp);

这有效地过滤掉了由于 C++ 语言本身的舍入错误而导致的任何不精确为零的值。但是它不会过滤掉由输入信号中的噪声引起的垃圾数据。

将此添加到我的相位计算后，我得到了预期的结果。

如何有效地计算（基本上）等于零的复数的相位角？

How can I effectively calculate the phase angle of a complex number that is (essentially) equal to zero?

c++

fft

precision

division

phase

原答案