Using 2D DFT Transformation in Intel IPP 8.2 (Image DFT)

[Royi]

Hello,

I'm comparing the results of the magnitude of DFT Transform of Intel DFTFwd to MATLAB's fft.

The input is Real Image in 32fC1 format.
The code snippet is:

int sizeSpec;
int sizeInit;
int sizeBuffer;
/// get sizes for required buffers
	IppStatus status = ippiDFTGetSize_R_32f(mI.size_, IPP_FFT_DIV_FWD_BY_N, ippAlgHintNone,
                   &sizeSpec, &sizeInit, &sizeBuffer);

    /// allocate memory for required buffers
    IppiDFTSpec_R_32f* pMemSpec = (IppiDFTSpec_R_32f*)ippMalloc(sizeSpec);

Ipp8u* pMemInit   = NULL;
Ipp8u* pMemBuffer = NULL;
    if (sizeInit > 0)
         pMemInit = (Ipp8u*)ippMalloc(sizeInit);

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

    /// initialize DFT specification structure
    status = ippiDFTInit_R_32f(mI.size_, IPP_FFT_DIV_FWD_BY_N, ippAlgHintNone, pMemSpec, pMemInit);

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

    /// perform forward DFT to put source data to frequency domain
status = ippiDFTFwd_RToPack_32f_C1IR(mI, mI.bytes_step_[0], pMemSpec, pMemBuffer);
status = ippiMagnitudePack_32f_C1R(mI, mI.bytes_step_[0], mI, mI.bytes_step_[0], mI.size_);
//status = ippiPhasePack_32f_C1R(mI, mI.bytes_step_[0], mI, mI.bytes_step_[0], mI.size_);

/// ...

    /// free buffers
    if (sizeBuffer > 0)
        ippFree(pMemBuffer);

    ippFree(pMemSpec);
}

The result is different from MATLAB's (In the bottom right quadrant).
I also tried the phase and got it is different as well (Again, all perfect but certain quadrant).

Any idea what's happening?
Have I missed something?

Thank You.

[Igor_A_Intel]

Hi Royi,

do you use Pack format for Matlab code too? Do you know the difference between Pack, Perm, CCS and complex FFT/DFT result representation? It is described in details in the IPP manual.

PS Matlab uses IPP FFT/DFT implementation internally.

Regards, Igor

[Royi]

Hi Igor,

I don't compare pack format.
I compare Magnitude after using ippiMagnitudePack_32f_C1R().
So the output is only a floating point array of numbers.

At the end I want array at the size of the image with the magnitude of the DFT.

How could one achieve that?

Thank You.

[Igor_A_Intel]

Hi Royi,

take a look at the IPP manual, please, - it states (for MagnitudePack): "This function computes magnitude of elements of the source image pSrc given in Real - Complex Packed (RCPack2D) Format and stores results in the destination image pDst."

Pack format is based on DFT/FFT symmetry for real input and stores only half of result as the second half is just the mirror reflection of first half. If you need the same result as Matlab provides - you should use "RToC" DFT/FFT and magnitude for complex "image", or should unpack packed format to its complex representation with PackToCplxExtend
function.

regards, Igor

[Royi]

Igor,

Either the function doesn't make sense or you're assuming I didn't read.

Let's say our image is 100 x 100.
So its DFT is also 100 x 100 array of complex numbers.

In order to calculate the magnitude you could run on each element in the array (Complex Element) and calculate it.
Namely the Magnitude of the DFT of an 100 x 100 image is a 100 x 100 array of real numbers.

Now, the packed format says if the image is real there is a redundancy in its DFT (Conjugate Symmetry).
Hence one could store the DFT of real image using real numbers in the same size of array (100 x 100).

But the magnitude calculation must yield the same result as above.
This is why you created 'ippiMagnitudePack_32f_C1R()`.
According to documentation it should create an array of 100 x 100 (Same as the original image) of the image DFT Magnitude as one would do the "Long Way". Otherwise, what's the point?

Igor,
Please feel free to share code with us which just give us what we need using IPP.
We followed the example in the site. The result doesn't make sense.
Pay attention we are not comparing a packed array of Intel IPP to MATLAB.
We compare the Magnitude Array of both (Which should be identical up to numerical issues) and Intel's has only 1 quadrant which is weird.

We read about the pack format, but it is just a format, if functions are adapted to it, you expect the result to be independent of the format.

Thank You.

[Royi]

Hello,

Anyone could assist with this?

It seems the function Magnitude, though is documented to be aware of the packed format does something it shouldn't and mess with the quadrant.

Any assistance?

[Igor_A_Intel]

Hi Royi,

according to your pictures it seems to be a real bug. We will investigate this issue and I'll provide you an update (in ~1-3 weeks).

Regards, Igor.

[Royi]

Hi Igor,

Yea, we also think it is some kind of a bug.

It's important for us so please let us know as soon as possible.

Thank You.

[Royi]

Hi,

Any update on this?

[Igor_A_Intel]

Hi Royi,

Safe and correct in-place operation in IPP is guaranteed only for functions with suffix "_I". If function doesn't have such suffix - we don't guarantee correct behavior for pSrc == pDst. In your case you've faced with such issue. If you use this function in correct mode - separated non-overlapping source and destination - the result will be correct.

regards, Igor.

[Royi]

Hi Igor,

We tried that as well and the result is the same.

[Igor_A_Intel]

Hi Royi,

please prove your last statement (some reproducer appreciated): we've performed double check (comparison with Matlab and comparison with CToC flavor) - MagnitudePack doesn't show any issues.

regards, Igor

[Royi]

We moved to Intel IPP 2018 Update 1 and things seems to work (The error is 1e-4, compared to MATLAB).
We're not sure about the previous version as we don't have access to it anymore.

There is a different issue bothering us.
Once we have the data we decompose it into Magnitude and Phase using:

ippiDFTFwd_RToPack_32f_C1IR(mI, mI.bytesStep, pMemSpec, pMemBuffer);
 ippiMagnitudePack_32f_C1R  (mI, mI.bytesStep, mO.Plane(0), mO.bytesStep, mO.size_, pMagBuffer);
 ippiPhasePack_32f_C1R      (mI, mI.bytesStep, mO.Plane(1), mO.bytesStep, mI.size_, pPhaseBuffer);

The problem is once we have data in form of Phase and Magnitude we have no way going back.
Assume we keep the properties of the data such that the inverse is real, we still don't have a built in function to go back.

What's the recommended solution?

Thank You.

[Royi]

Igor,

Any assistance with the above?

[Royi]

Igor / Anyone form Intel,

Any assistance on this?

Thank You.