c# - 在C#中使用FFTW计算HilbertTransform

我想在C＃中实现希尔伯特变换。从this文章中，我看到最快的FFT开源实现似乎是FFTW，因此我下载了该示例并将其用于学习如何将fftw包装器用于C＃。

我有200.000点的电流信号用于测试。通过fft获得Hilbert转换相对简单：


计算fft。
除直流分量和奈奎斯特分量外，将所有正频率乘以2（如果样本大小为偶数，则为0和n / 2 +1）。
将所有负频率（[n / 2 + 1，n]）乘以0。
计算逆ftf。


到目前为止，我已经完成了所有步骤。唯一的问题是逆fft。使用fftw无法获得与Matlab的ifft相同的结果。

我的密码

    RealArray _input;
    ComplexArray _fft;

    void ComputeFFT()
    {
        _fft = new ComplexArray(_length / 2 + 1);
        _input.Set(Data);
        _plan = Plan.Create1(_length, _input, _fft, Options.Estimate);
        _plan.Execute();
    }

到目前为止，我只使用了正频率。因此，我不需要将负频率乘以零：它们甚至不存在。使用以下代码，我可以得到原始信号：

    double[] ComputeIFFT(ComplexArray input)
    {
        double[] temp = new double[_length];
        RealArray output = new RealArray(_length);

        _plan = Plan.Create1(_length, input, output, Options.Estimate);
        _plan.Execute();
        temp = output.ToArray();

        for (int i = 0; i < _length; ++i)
        {
            temp[i] /= _length;
        }

        return temp;
    }

当我尝试从信号中获得一个复杂的逆时，问题就来了。

    void ComputeHilbert()
    {
        double[] fft = FFT.ToArray();
        double[] h = new double[_length / 2 + 1];
        double[] temp = new double[_length * 2];

        bool fftLengthIsOdd = (_length | 1) == 1;

        h[0] = 1;

        for (int i = 1; i < _length / 2; i++) h[i] = 2;

        if (!fftLengthIsOdd) h[(_length / 2)] = 1;

        for (int i = 0; i <= _length / 2; i++)
        {
            temp[2 * i] = fft[2*i] * h[i];
            temp[2 * i + 1] = fft[2*i + 1] * h[i];
        }

        ComplexArray _tempHilbert = new ComplexArray(_length);

        _tempHilbert.Set(temp);

        _hilbert = ComputeIFFT(_tempHilbert);
        _hilbertComputed = true;
    }

重要的是要注意，当我将ToArray()方法应用于ComplexArray对象时，得到的double[]长度是原始数组的两倍，并且实部和虚部连续。就这样，对于包含“ 3 + 1i”的ComplexArray对象，我将得到一个带有[3，1]的双精度矢量。

因此，目前，我所拥有的是：

[DC Frequency, 2*positive frequencies, Nyquist Frequency, zeros]
如果将此数据导出到Matlab并计算IFFT，则得到的结果与其hilbert（signal）相同。

但是，如果我尝试应用fftw提供的IFFT，则会从奈奎斯特频率到末尾得到怪异的值（也就是说，零与fftw混为一谈）。

这是我正在使用的ifft：

double[] ComputeIFFT(ComplexArray input)
{
    double[] temp;
    ComplexArray output = new ComplexArray(_length);

    _plan = Plan.Create1(_length, input, output, Direction.Backward, Options.Estimate);
    _plan.Execute();
    temp = output.ToArray();

    for (int i = 0; i < _length; ++i)
    {
        temp[i] /= _length;
    }

    return temp;
}

因此，总而言之，我的问题是我用来计算ifft的方式。它似乎不能很好地与零。也许Matlab有能力理解它必须应用一些不同的方法，而我应该手动进行，但是我不知道该怎么做。

非常感谢您的提前帮助，非常感谢！

最佳答案

因此问题出在ComputeIFFT函数上。在for循环中，我正在执行i
这就是为什么我只有一半正确的价值观。

正确的代码是：

double[] ComputeIFFT(ComplexArray input)
{
    double[] temp;
    ComplexArray output = new ComplexArray(_length);

    _plan = Plan.Create1(_length, input, output, Direction.Backward, Options.Estimate);
    _plan.Execute();
    temp = output.ToArray();

    for (int i = 0; i < temp.Length; ++i)
    {
        temp[i] /= _length;
    }

    return temp;
}

我希望这对尝试通过C＃通过FFTW实现希尔伯特变换的人有用。

关于c# - 在C#中使用FFTW计算HilbertTransform，我们在Stack Overflow上找到一个类似的问题：https://stackoverflow.com/questions/37988344/