Search code examples
c++fftintel-ipp

Intel integrated performance primitives Fourier Transform magnitudes


When I am using Intel IPP's ippsFFTFwd_RToCCS_64f and then ippsMagnitude_64fc I get a massive peak at zero index in magnitudes array.

My sine wave is long and main component I am interested is somewhere between 0.15 Hz and 0.25 Hz. I take the sample with 500Hz sampling frequency. If I reduce mean from the signal before FFT I get really small zero component not that peak anymore. Below is a pic of magnitudes array head:

enter image description here

Also the magnitude scaling seems to be 10 times the magnitude I see in the time series of the signal e.g. if amplitude is 29 in magnitudes it is 290.

I Am not sure why this is so and my question is 1. Do I really need to address the zero index peak with mean reduction and 2. Where does this scale of 10 come?

void CalculateForwardTransform(array<double> ^signal, array<double> ^transformedSignal, array<double> ^magnitudes)
{ 
    // source signal
    pin_ptr<double> pinnedSignal = &signal[0];
    double *pSignal = pinnedSignal;
    int order = (int)Math::Round(Math::Log(signal->Length, 2));

    // get sizes
    int sizeSpec = 0, sizeInit = 0, sizeBuf = 0;
    int status = ippsFFTGetSize_R_64f(order, IPP_FFT_DIV_INV_BY_N, ippAlgHintNone, &sizeSpec, &sizeInit, &sizeBuf);

    // memory allocation
    IppsFFTSpec_R_64f* pSpec;
    Ipp8u *pSpecMem = (Ipp8u*)ippMalloc(sizeSpec);
    Ipp8u *pMemInit = (Ipp8u*)ippMalloc(sizeInit);
    
    //  FFT specification structure initialized
    status = ippsFFTInit_R_64f(&pSpec, order, IPP_FFT_DIV_INV_BY_N, ippAlgHintNone, pSpecMem, pMemInit);

    // transform
    pin_ptr<double> pinnedTransformedSignal = &transformedSignal[0];
    double *pDst = pinnedTransformedSignal;
    Ipp8u *pBuffer = (Ipp8u*)ippMalloc(sizeBuf);
    status = ippsFFTFwd_RToCCS_64f(pSignal, pDst, pSpec, pBuffer);

    // get magnitudes
    pin_ptr<double> pinnedMagnitudes = &magnitudes[0];
    double *pMagn = pinnedMagnitudes;
    status = ippsMagnitude_64fc((Ipp64fc*)pDst, pMagn, magnitudes->Length); // magnitudes is half of signal len

    // free memory
    ippFree(pSpecMem);
    ippFree(pMemInit);
    ippFree(pBuffer);
}

Solution

  • Do I really need to address the zero index peak with mean reduction?

    For low frequency signal analysis a small bias can really interfere (especially due to spectral leakage). For sake of illustration, consider the following low-frequency signal tone and another one with a constant bias tone_with_bias:

    fs = 1;
    f0 = 0.15;
    for (int i = 0; i < N; i++)
    {
      tone[i] = 0.001*cos(2*M_PI*i*f0/fs);
      tone_with_bias[i] = 1 + tone[i];
    }
    

    If we plot the frequency spectrum of a 100 sample window of these signals, you should notice that the spectrum of tone_with_bias completely drowns out the spectrum of tone: Spectrum of tone & tone with bias

    So yes it's better if you can remove that bias. However, it should be emphasized that this is possible provided that you know the nature of the bias. If you know that the bias is indeed a constant, removing it will reveal the low-frequency component. Otherwise, removing the mean from the signal may not achieve the desired effect if the bias is just a very low-frequency component.

    Where does this scale of 10 come?

    Scaling of the magnitude by the FFT should be expected, as described in this answer of approximately 0.5*N (where N is the FFT size). If you were processing a small chunk of 20 samples, then you should get such a factor of 10 scaling. If you scale the output of the FFT by 2/N (or equivalently scale by 2 while also using the IPP_FFT_DIV_FWD_BY_N flag) you should get results that have similar magnitudes as the time-domain signal.