使用Intel IPP FFT处理非2的幂次数组时零填充结果异常问题

我在项目中使用Intel IPP的FFT（快速傅里叶变换）函数，但输入数据大小可变。IPP的FFT函数仅支持长度为2的幂次的数组，因此我尝试通过零填充将数组补为2的幂次长度，但得到的结果与MATLAB等程序的结果不符。

注：当输入数组长度N为2的幂次时结果正常，且我清楚零填充后结果会有差异，已采用正确方式查看结果。

请问是我的零填充方式有误？还是Intel IPP的FFT即使零填充也不支持原长度非2的幂次的数组？是否必须改用Intel MKL等其他库？

相关代码如下：

#include <ipp.h>

// N is input array size

void IPP_fft(int N){

Ipp32fc* src;  //FFT input array
Ipp32fc* dst;  //FFT output array

int powerOf2Sz = pow(2,ceil(log2(N))); //Finding input's upper power of 2 value
                                      //if it's equal to power of 2 it's equal to N

int orderOfFFT = (int)(log((double)powerOf2Sz) / log(2.0));

src = ippsMalloc_32fc(powerOf2Sz);
dst = ippsMalloc_32fc(powerOf2Sz);

//set input with variable

for(int i = 0; i<N;i++){
   src[i].re = i+1;
   src[i].im = 0;
}

//do zero padding for rest of size

for(int i = N; i<powerOf2Sz;i++){
   src[i].re = 0;
   src[i].im = 0;
}

IppsFFTSpec_C_32fc* pSpec=0;

Ipp8u* pMemSpec   = 0;
Ipp8u* pMemInit   = 0;
Ipp8u* pMemBuffer = 0;

int sizeSpec   = 0;
int sizeInit   = 0;
int sizeBuffer = 0;

int flag = IPP_FFT_NODIV_BY_ANY;

//get sizes for required buffer
ippsFFTGetSize_C_32fc(FFTOrder, flag, ippAlgHintNone, &sizeSpec, &sizeInit, &sizeBuffer);

//allocate memory for buffer
pMemSpec = (Ipp8u*)ippMalloc(sizSpec);

if(sizeInit > 0){
pMemInit = (Ipp8u*)ippMalloc(sizInit);
}

if(sizeBuffer > 0){
pMemBuffer = (Ipp8u*)ippMalloc(sizBuffer);
}

//Initialize FFT Specification Structures
ippsFFTInit_C_32fc(&pSpec, FFTOrder, flag, ippAlgHintNone, pMemSpec, pMemInit);

//free initialization buffer
if(sizeInit > 0){
ippFree(pMemInit);
}

//perform Forward FFT
ippsFFTFwd_CToC_32fc(src,dst,pSpec,pMemBuffer);

//you can use dst output array

if(sizeBuffer > 0){
ippFree(pMemBuffer);
}

ippFree(pMemSpec);

}

你的问题既不是零填充方式错误，也不是IPP不支持零填充后的非2幂次原数组，核心问题是代码里的变量拼写错误导致FFT初始化逻辑异常，最终结果偏差。具体问题和修正方案如下：

1. 核心错误点 #

• 未定义变量FFTOrder：代码中实际定义的FFT阶数变量是orderOfFFT，但后续调用ippsFFTGetSize_C_32fc和ippsFFTInit_C_32fc时误用了大小写错误的FFTOrder，触发未定义行为。
• 变量名拼写错误：sizSpec和sizInit应为sizeSpec和sizeInit（少了字母e），会导致内存分配时使用错误的大小值，引发内存越界或初始化失败。
• 浮点计算精度隐患：用pow(2,ceil(log2(N)))计算2的幂次可能因浮点精度问题出错，比如极端情况下log2(N)的结果会因精度偏差被ceil多算一位。

2. 修正后的代码 #

#include <ipp.h>

// 整数方法计算大于等于N的最小2的幂次，避免浮点误差
int nextPowerOf2(int N) {
    if (N <= 1) return 1;
    int sz = 1;
    while (sz < N) sz <<= 1;
    return sz;
}

void IPP_fft(int N){
    Ipp32fc* src;  // FFT输入数组
    Ipp32fc* dst;  // FFT输出数组

    int powerOf2Sz = nextPowerOf2(N);
    int orderOfFFT = 0;
    int temp = powerOf2Sz;
    // 整数方法计算FFT阶数
    while (temp >>= 1) orderOfFFT++;

    src = ippsMalloc_32fc(powerOf2Sz);
    dst = ippsMalloc_32fc(powerOf2Sz);

    // 填充输入数据
    for(int i = 0; i < N; i++){
        src[i].re = i + 1;
        src[i].im = 0;
    }

    // 零填充剩余部分
    for(int i = N; i < powerOf2Sz; i++){
        src[i].re = 0;
        src[i].im = 0;
    }

    IppsFFTSpec_C_32fc* pSpec = 0;
    Ipp8u* pMemSpec = 0;
    Ipp8u* pMemInit = 0;
    Ipp8u* pMemBuffer = 0;

    int sizeSpec = 0;
    int sizeInit = 0;
    int sizeBuffer = 0;

    int flag = IPP_FFT_NODIV_BY_ANY;

    // 获取所需缓冲区大小
    ippsFFTGetSize_C_32fc(orderOfFFT, flag, ippAlgHintNone, &sizeSpec, &sizeInit, &sizeBuffer);

    // 分配内存
    pMemSpec = (Ipp8u*)ippMalloc(sizeSpec);
    if(sizeInit > 0){
        pMemInit = (Ipp8u*)ippMalloc(sizeInit);
    }
    if(sizeBuffer > 0){
        pMemBuffer = (Ipp8u*)ippMalloc(sizeBuffer);
    }

    // 初始化FFT规格结构
    ippsFFTInit_C_32fc(&pSpec, orderOfFFT, flag, ippAlgHintNone, pMemSpec, pMemInit);

    // 释放初始化缓冲区
    if(sizeInit > 0){
        ippFree(pMemInit);
    }

    // 执行正向FFT
    ippsFFTFwd_CToC_32fc(src, dst, pSpec, pMemBuffer);

    // 使用dst输出数组...

    // 释放所有资源
    if(sizeBuffer > 0){
        ippFree(pMemBuffer);
    }
    ippFree(pMemSpec);
    ippsFree(src);
    ippsFree(dst);
}

3. 额外注意事项 #

• 资源释放：原代码未释放src和dst，会导致内存泄漏，修正后的代码补充了这部分逻辑。
• 缩放因子：你使用的IPP_FFT_NODIV_BY_ANY标志和MATLAB的fft行为一致（正向FFT不缩放），修正后结果会和MATLAB零填充后的FFT结果对齐。
• 无需改用MKL，IPP完全支持零填充到2的幂次后的FFT计算，只要代码逻辑正确即可。

内容的提问来源于stack exchange，提问作者Hansov