我是一名从事电信项目的计算机程序员
在我们的项目中，我必须将一系列复数更改为它们的傅立叶变换。因此，我需要一个用于C89标准的有效FFT代码
我正在使用以下代码，它运行得很好：
short FFT(short int dir,long m,double *x,double *y)
{
long n,i,i1,j,k,i2,l,l1,l2;
double c1,c2,tx,ty,t1,t2,u1,u2,z;
/* Calculate the number of points */
n = 1;
for (i=0;i<m;i++)
n *= 2;
/* Do the bit reversal */
i2 = n >> 1;
j = 0;
for (i=0;i<n-1;i++) {
if (i < j) {
tx = x[i];
ty = y[i];
x[i] = x[j];
y[i] = y[j];
x[j] = tx;
y[j] = ty;
}
k = i2;
while (k <= j) {
j -= k;
k >>= 1;
}
j += k;
}
/* Compute the FFT */
c1 = -1.0;
c2 = 0.0;
l2 = 1;
for (l=0;l<m;l++) {
l1 = l2;
l2 <<= 1;
u1 = 1.0;
u2 = 0.0;
for (j=0;j<l1;j++) {
for (i=j;i<n;i+=l2) {
i1 = i + l1;
t1 = u1 * x[i1] - u2 * y[i1];
t2 = u1 * y[i1] + u2 * x[i1];
x[i1] = x[i] - t1;
y[i1] = y[i] - t2;
x[i] += t1;
y[i] += t2;
}
z = u1 * c1 - u2 * c2;
u2 = u1 * c2 + u2 * c1;
u1 = z;
}
c2 = sqrt((1.0 - c1) / 2.0);
if (dir == 1)
c2 = -c2;
c1 = sqrt((1.0 + c1) / 2.0);
}
/* Scaling for forward transform */
if (dir == 1) {
for (i=0;i<n;i++) {
x[i] /= n;
y[i] /= n;
}
}
return(true);
}
但这段代码，只支持2^m大小的数组。就像CLRS书中的代码一样
我们应该转换的数组不在这个范围内，加零将是昂贵的，所以我正在寻找另一种解决方案，它可以帮助我获得任何大小的输入
类似于IT++和matlab所做的事情。但正如我们在纯C中所希望的那样，使用它们是不可能的。此外，当我检查
IT++代码被阻止如果您在任何大众市场计算平台上工作（Intel with Windows或OS X、iOS等），则供应商或制造商提供了高性能FFT实现。
否则，应评估FFTW。
为二次幂以外的大小编写高性能FFT是一项复杂的任务。
如果你要使用自己的实现，那么，关于两种大小的幂：
您展示的实现在FFT期间计算sqrt。大多数高性能FFT实现提前计算常数并将其存储在表中。
缩放包含x[i] /= n和y[i] /= n中的除法运算。编译器可能会将这些指令作为除法指令来实现。除法在普通处理器上通常是一条缓慢的指令。最好计算一次CCD_ 11并乘以CCD_ 12，而不是除以CCD_ 13。
更好的做法是完全省略比额表。变换通常在没有缩放的情况下很有用，或者可以从单个变换中省略缩放，并将其作为聚合缩放应用一次。（例如，不进行两次缩放操作，一次在正向FFT中，另一次在反向FFT中，将缩放操作排除在FFT例程之外，在正向FFT和反向FFT之后只进行一次。）
如果频域数据的比特反转顺序是可以接受的，那么您可以省略比特反转排列。
如果保留比特反转排列，则可以对其进行优化。实现这一点的技术依赖于平台。一些平台具有用于反转整数中的位的指令（例如，ARM具有rbit）。如果您的平台没有，您可能希望将位反转索引保存在表中，或者研究比当前代码更快地计算它们的方法。
如果同时保留位反转排列和缩放，则应考虑同时执行它们。排列使用了大量的内存运动，缩放使用了处理器的算术单元。大多数现代处理器可以同时完成这两项操作，因此您可以从重叠操作中获得一些好处。
您当前的代码使用基数-2蝴蝶。Radix-4通常更好，因为它可以通过改变使用的数据段和将一些加法改为减法来实现与i的乘法，反之亦然，从而获得一些优势。
如果您的数组长度接近处理器上一级内存缓存的大小，则FFT实现的部分将对缓存进行猛烈抨击并显著降低速度，除非您设计适当的代码来处理此问题（通常是通过将数组的部分临时复制到缓冲区中）。
如果你的目标处理器有SIMD功能，你绝对应该在FFT中使用这些功能；它们极大地提高了FFT性能。
以上内容应该告诉您，编写高效的FFT是一项复杂的任务。除非您想在这方面花费大量精力，否则最好使用FFTW或其他现有实现。
在您的实现中，我关心的是以下代码：
z = u1 * c1 - u2 * c2;
u2 = u1 * c2 + u2 * c1;
u1 = z;
当l1较大时，（u1，u2）会经历大量的累积舍入误差。您可能会得到不准确的转换。
我赞同FFTW的建议，但我相信它是高度针对平台的。（大多数FFT库都是。）[EDIT：不，它实际上是直C89。这正是你所说的你想要的。]
维基百科FFT页面列出了一组适用于大小怪异的输入数组的算法。我不知道它们是如何工作的，但我相信一般的想法是，你可以使用Rader算法或Bluestein算法来处理素数大小的输入，并使用Cooley Tukey将复合大小的变换简化为一组素数大小的变换。
对于FFTW的替代方案，请查看我的混合基数FFT，它也处理非2^N的FFT长度。C源可从http://www.corix.dk/Mix-FFT/mix-fft.html.