Swift中的DFT结果与MATLAB不同

Question

import Cocoa
import Accelerate

let filePath = Bundle.main.path(forResource: "sinusoid", ofType: "txt")
let contentData = FileManager.default.contents(atPath: filePath!)
var content = NSString(data: contentData!, encoding: String.Encoding.utf8.rawValue) as? String

var idx = content?.characters.index(of: "\n")
idx = content?.index(after: idx!)

repeat {
    //let fromIndex = index(from: )
    content = content?.substring(from: idx!)
    idx = content?.characters.index(of: "\n")
    idx = content?.index(after: idx!)
} while content!.characters.contains("%")

let regex = try? NSRegularExpression(pattern: "[ ]+", options:[])

let delimiter = ","
var modifiedString = regex?.stringByReplacingMatches(in: content!, options: [], range: NSRange(location: 0, length: (content! as NSString).length), withTemplate: delimiter)

let lines = modifiedString?.components(separatedBy: "\n")

var s = [Double]()

for var line in lines! {
    if !line.isEmpty {
        let data = line.components(separatedBy: ",")
        s.append(Double(data[1])!)
    }
}

let length = vDSP_Length(pow(2, floor(log2(Float(s.count)))))
let L = Int(length)

// zrop or zop? 
// zrop covers real to complex, and zop covers complex
// length must be a power of 2 or specific multiples of powers of 2 if size is at least 4
let setup = vDSP_DFT_zrop_CreateSetupD(nil, length, vDSP_DFT_Direction.FORWARD)

var inputReal = UnsafeMutablePointer<Double>.allocate(capacity: L)
var inputImaginary = UnsafeMutablePointer<Double>.allocate(capacity: L)
var outputReal = UnsafeMutablePointer<Double>.allocate(capacity: L)
var outputImaginary = UnsafeMutablePointer<Double>.allocate(capacity: L)

for i in 0..<L {
    inputReal[i] = s[i]
    inputImaginary[i] = 0.0
}

vDSP_DFT_ExecuteD(setup!, inputReal, inputImaginary, outputReal, outputImaginary)

for i in 0..<L {
    print("\(outputReal[i]) + \(outputImaginary[i])i")
}

输入文件"sinusoid.txt"在以下链接中：https://dpaste.de/M1VD 输入文件数据包含两个正弦波，频率分别为50和120。Matlab代码生成了正确的输出结果，具体请参考以下链接：https://dpaste.de/2mdK 当从Matlab中得到的结果进行缩放并取幅值时，可以正确地显示出50Hz频率的振幅为0.7，而120Hz频率的振幅为1。

clear all; close all; clc;
data = load('sinusoid.txt');
S = data(:,2);
Fs = 1000;
Y = fft(S);
L = length(S);
P2 = abs(Y/L);
P1 = P2(1:L/2+1);
P1(2:end-1) = 2*P1(2:end-1);
f = Fs*(0:(L/2))/L;
plot(f,P1)
title('Single-Sided Amplitude Spectrum of X(t)')
xlabel('f (Hz)')
ylabel('|P1(f)|')

Swift代码的输出与Matlab输出完全不同，无论应用何种缩放因子以及应用实际到复杂或复杂到复杂的转换，这都是不可识别的。有什么想法吗？请参考以下链接：https://dpaste.de/MUwB。

- SwiftMatt

---

你输入了1000个数字，但只有512个被传递给DSP函数，剩下的488个值被简单地忽略了。这怎么可能？ - Martin R

尝试使用1024进行测试 - SwiftMatt

我对MATLAB不熟悉，但看起来你的MATLAB代码中1000个元素的数据集没有被截断到512。因此，你不能期望得到相同的结果。 - Martin R

好的，如果我不截断数据，而是将0填充到1024个数据元素，马丁R仍然没有得到相同的结果。 - SwiftMatt

据我了解，来自加速框架的FFT函数仅能处理形如2的n次幂乘以1、3、5或15的样本大小，而MATLAB可以处理任意大小的样本。相关问题：http://stackoverflow.com/q/12134208/1187415，https://dev59.com/KWgv5IYBdhLWcg3wSe-2和http://math.stackexchange.com/q/77118/42969来自MSE。正如您从最后一个参考所看到的那样，“任意大小”情况转换为“2的n次幂”情况并不容易，您不能简单地用零填充采样。 - Martin R

Answer 1

实际上，FFT结果不匹配的问题是由于输入大小不匹配引起的。将输入限制为2的幂的特定倍数会极大地限制加速框架中FFT的使用。建议是用0填充输入，直到其长度合适。无论您是用0填充还是截断输入，使其大小成为2的幂的特定倍数，加速框架的结果都与诸如MATLAB之类的程序的结果不同。解决此问题的方法是执行chirp-z变换，如Martin R指定的链接所述。chirp-z变换本身产生与FFT相同的结果，并且可以在任意大小的输入上执行。

- SwiftMatt

Answer 2

你的两个FFT长度不同，因此结果肯定不匹配。你还向两个FFT传递了不同数量的数据。

打印出FFT长度和输入数据向量以调试代码。在比较结果之前确保输入匹配。

另外，苹果的加速/vDSP FFT可以使用除2的幂外的其他长度（具有3或5的因子的长度也是允许的）。

另外，请注意Matlab将数组索引从1开始，而不是像C和Swift函数中更典型的从0开始。

- hotpaw2