使用Matlab计算大数据（16GB）的FFT

Question

我正在尝试计算从一个文本文件导入的大块数据的快速傅里叶变换，该文件大小约为16 GB。我试图想出一种在Matlab中计算其FFT的方法，但由于我的计算机内存（8 GB），它给出了内存不足的错误。我尝试使用memmap、textscan，但无法应用于获取组合数据的FFT。

请问有人能够指导我如何获得傅里叶变换吗？我还尝试在远程服务器上使用C++代码获取傅里叶变换（使用定义），但执行时间很长。请问有人能够给我适当的见解，告诉我如何处理这个大数据量？

- Nehal P

---

要么尝试将数据分成较小的块，以便您的计算机可以处理它，要么获取更多内存。 - Some programmer dude

你是在实现自己的FFT还是使用内置的fft函数fft(x)？ - shieldfoss

我真的很想知道你打算用那个FFT做什么？ - Roger Rowland

你真的想要对数十亿个点进行FFT吗？因为通常在更长的信号上，你会计算一个频谱图（它只是在短时间间隔内重复计算FFT）。 - MSalters

@MSalters 只是一个应用示例，请看这个链接。该应用程序是用来搜索脉冲星，而在那篇论文中谈到的数据范围从2^29到2^33不等。我不知道 OP 的应用程序，但有时需要对大数据集进行 FFT 处理。 - triple_r

Answer 1

最好使用自己的代码实现FFT。

FFT算法有一个“蝴蝶”操作。因此，您可以将整个步骤分成较小的块。

文件大小对于典型的个人电脑来说太大了。但是FFT并不需要一次性处理所有数据。它始终可以从2点（也许8点更好）FFT开始，并且您可以通过级联阶段来逐步构建。这意味着您可以一次只读取少量点，进行一些计算，然后保存数据到磁盘。下一次迭代时，您可以从磁盘中读取保存的数据。

根据您如何构建数据结构，您可以将所有数据存储在一个单独的文件中，并使用指针（在Matlab中仅是一个数字）进行读/写；或者您可以将每个单独的点存储在一个单独的文件中，生成数十亿个文件，并通过文件名加以区分。

这个想法是你可以把计算结果存入磁盘而不是内存中。当然，这需要更为可行的磁盘空间来完成。

---

我可以为您翻译一段伪代码。根据原始数据的数据结构（那个16GB的文本文件），实现方式可能会有所不同，但您可以轻松地操作它。我将从2点FFT开始，并在this wikipedia picture中使用8点采样进行操作。

1. 对x进行2点FFT，生成y，即从左侧开始的白色圆圈的第三列。

read x[0], x[4] from file 'origin'
y[0] = x[0] + x[4]*W(N,0);
y[1] = x[0] - x[4]*W(N,0);
save y[0], y[1] to file 'temp'
remove x[0], x[4], y[0], y[1] from memory
read x[2], x[6] from file 'origin'
y[2] = x[2] + x[6]*W(N,0);
y[3] = x[2] - x[6]*W(N,0);
save y[2], y[3] to file 'temp'
remove x[2], x[6], y[2], y[3] from memory
....

2. 对 y 进行2点FFT，生成白色圆圈的第5列 z。

3. 对 z 进行2点FFT，生成最终结果 X。

基本上 Cooley–Tukey FFT algorithm 的设计是为了让你可以将数据分割并逐个计算，因此可以处理大量数据。我知道这不是一种常规方式，但如果您能查看中文版的维基百科页面，您可能会发现许多图片，这些图片可能有助于您理解它如何拆分点。

- Yvon

---

谢谢Yvon，我仍然不清楚如何分离数据，应用FFT并重新组合它们。您能否为此举个例子？ - Nehal P

我向您展示了从x到y的第一次迭代。您可以类似地执行以下迭代。 - Yvon

Answer 2

这取决于您需要的FFT分辨率。如果您只需要一个1024点的FFT，那么您可以将数据重新整形为或按顺序读取N x 1024块。一旦您以这种格式获得了数据，就可以将每个FFT结果的输出添加到一个1024点复累加器中。

如果您需要相同的FFT分辨率，则需要更多的内存或一个特殊的fft程序。Matlab中未包含此程序（但我不确定通过缓冲小块进行完整分辨率的部分FFT是否在数学上可行）。

- learnvst

---

谢谢，您能详细解释一下吗？我不太明白您的意思。我对Matlab还很陌生，谢谢。 - Nehal P

Answer 3

我遇到了同样的问题。最终在一篇论文中找到了解决方案：扩展有效卷积算法的尺寸。它基本上涉及加载较短的块，乘以相位因子并进行FFT，然后在系列中加载下一个块。这给出了完整信号的总FFT的采样。然后使用不同的相位因子多次重复此过程以填充其余点。我将在此处尝试总结（从论文中的表II改编）：

1. 对于长度为 N 的总信号 f(j)，决定一个数字 m 或更短的块，每个块的长度为 N/m，可以将它们存储在内存中（如果需要，用零填充信号，使得 N 是 m 的倍数）。

2. 对于 beta = 0, 1, 2, ... ,m - 1，执行以下操作：

3. 将新序列分成 m 个长度为 N/m 的子区间。

4. 对于每个子区间，将每个第 j 个元素乘以 exp(i*2*pi*j*beta/N)。这里，j 根据点相对于整个数据流中的第一个点的位置进行索引。

5. 将每个子区间的第一个元素相加，得到一个单一的数字，将第二个元素相加，依此类推。这可以在从文件读取点时完成，因此无需将全部的 N 个点存储在内存中。

6. 对结果序列进行傅里叶变换，其中包含 N/m 个点。

7. 这将给出 k = ml + beta 时的 F(k)，其中 l = 0, ..., N/m-1。将这些值保存到磁盘。

8. 转到步骤2，继续下一个 beta 的值。

- Matthew Reichert