iOS加速低通FFT濾波器鏡像結果

我正嘗試使用Accelerate vDSP框架將現有的基於FFT的低通濾波器移植到iOS。iOS加速低通FFT濾波器鏡像結果

似乎FFT的預期效果約爲樣本的前1/4。但之後的結果似乎是錯誤的，甚至更奇怪的是鏡像（與信號的後半部分鏡像上半年的大部分）。

您可以在下面看到測試應用程序的結果。首先是繪製原始採樣數據，然後是預期濾波結果（濾除高於15Hz的信號）的示例，然後是我當前FFT代碼的結果（請注意，所需結果和示例FFT結果的比例不同於原始數據）：

FFT Results

爲我的低通濾波器的實際代碼如下：

double *lowpassFilterVector(double *accell, uint32_t sampleCount, double lowPassFreq, double sampleRate) 
{ 
    double stride = 1; 

    int ln = log2f(sampleCount); 
    int n = 1 << ln; 

    // So that we get an FFT of the whole data set, we pad out the array to the next highest power of 2. 
    int fullPadN = n * 2; 
    double *padAccell = malloc(sizeof(double) * fullPadN); 
    memset(padAccell, 0, sizeof(double) * fullPadN); 
    memcpy(padAccell, accell, sizeof(double) * sampleCount); 

    ln = log2f(fullPadN); 
    n = 1 << ln; 

    int nOver2 = n/2; 

    DSPDoubleSplitComplex A; 
    A.realp = (double *)malloc(sizeof(double) * nOver2); 
    A.imagp = (double *)malloc(sizeof(double) * nOver2); 

    // This can be reused, just including it here for simplicity. 
    FFTSetupD setupReal = vDSP_create_fftsetupD(ln, FFT_RADIX2); 

    vDSP_ctozD((DSPDoubleComplex*)padAccell,2,&A,1,nOver2); 

    // Use the FFT to get frequency counts 
    vDSP_fft_zripD(setupReal, &A, stride, ln, FFT_FORWARD); 


    const double factor = 0.5f; 
    vDSP_vsmulD(A.realp, 1, &factor, A.realp, 1, nOver2); 
    vDSP_vsmulD(A.imagp, 1, &factor, A.imagp, 1, nOver2); 
    A.realp[nOver2] = A.imagp[0]; 
    A.imagp[0] = 0.0f; 
    A.imagp[nOver2] = 0.0f; 

    // Set frequencies above target to 0. 

    // This tells us which bin the frequencies over the minimum desired correspond to 
    NSInteger binLocation = (lowPassFreq * n)/sampleRate; 

    // We add 2 because bin 0 holds special FFT meta data, so bins really start at "1" - and we want to filter out anything OVER the target frequency 
    for (NSInteger i = binLocation+2; i < nOver2; i++) 
    { 
     A.realp[i] = 0; 
    } 

    // Clear out all imaginary parts 
    bzero(A.imagp, (nOver2) * sizeof(double)); 
    //A.imagp[0] = A.realp[nOver2]; 


    // Now shift back all of the values 
    vDSP_fft_zripD(setupReal, &A, stride, ln, FFT_INVERSE); 

    double *filteredAccell = (double *)malloc(sizeof(double) * fullPadN); 

    // Converts complex vector back into 2D array 
    vDSP_ztocD(&A, stride, (DSPDoubleComplex*)filteredAccell, 2, nOver2); 

    // Have to scale results to account for Apple's FFT library algorithm, see: 
    // http://developer.apple.com/library/ios/#documentation/Performance/Conceptual/vDSP_Programming_Guide/UsingFourierTransforms/UsingFourierTransforms.html#//apple_ref/doc/uid/TP40005147-CH202-15952 
    double scale = (float)1.0f/fullPadN;//(2.0f * (float)n); 
    vDSP_vsmulD(filteredAccell, 1, &scale, filteredAccell, 1, fullPadN); 

    // Tracks results of conversion 
    printf("\nInput & output:\n"); 
    for (int k = 0; k < sampleCount; k++) 
    { 
     printf("%3d\t%6.2f\t%6.2f\t%6.2f\n", k, accell[k], padAccell[k], filteredAccell[k]); 
    } 


    // Acceleration data will be replaced in-place. 
    return filteredAccell; 
}

在原始代碼庫被處理的非乘方的兩種尺寸輸入數據;在我的加速代碼中，我將輸入填充到最接近的兩個冪。在下面的樣本測試中，原始樣本數據是1000個樣本，因此填充到1024.我不認爲這會影響結果，但爲了可能的差異，我將其包括在內。

如果你想一個解決實驗，你可以下載，在這裏生成的圖表（在FFTTest文件夾）示例項目：

FFT Example Project code

感謝任何見解，我不工作與FFT的之前，所以我覺得我失去了一些關鍵的東西。

來源

2013-05-27 Kendall Helmstetter Gelner

一大問題，你的方法（而不是你可能有什麼問題，實現）正當權力是試圖在頻率申請磚牆過濾器像這樣的域會在時域中產生巨大的文物。您需要使用窗口方法來避免這種情況。 –

您是否找到了解決方案？ – NTNT