Intel® oneAPI Math Kernel Library
Ask questions and share information with other developers who use Intel® Math Kernel Library.

How to use mkl for FFT?

薛_清_
Beginner
1,973 Views

I am a new student in mkl and c language.

now when i replace FFT function with DftiComputeForward, i can not get the correct out put,why?

Realft(info.x,info.fftN);  // FFT info.x is  float in-out x[1] is real ,x[2] is complex.

 

  MKL_Complex8* x_cmplx = 0;
 x_cmplx = (MKL_Complex8*)mkl_malloc((info.fftN/2 + 1) * sizeof(MKL_Complex8), 64);
 DftiComputeForward(hand, info.x,x_cmplx);

0 Kudos
6 Replies
Jing_Xu
Employee
1,973 Views

Please check https://software.intel.com/en-us/node/521971 to view how FFT results are stored.

0 Kudos
薛_清_
Beginner
1,973 Views

Jing X. (Intel) wrote:

Please check https://software.intel.com/en-us/node/521971 to view how FFT results are stored.

Hi ,Jing X

     Thanks for your help. I have read the instructions you show me,but still have some questions

     status = DftiSetValue(hand, DFTI_INPUT_STRIDES, <real_strides>);
     status = DftiSetValue(hand, DFTI_OUTPUT_STRIDES, <complex_strides>);

     How can i fill the  <real_strides> and <complex_strides>?

0 Kudos
Jing_Xu
Employee
1,973 Views

Please check https://software.intel.com/en-us/node/521967 for information about strides.

0 Kudos
薛_清_
Beginner
1,973 Views

Jing X. (Intel) wrote:

Please check https://software.intel.com/en-us/node/521967 for information about strides.

thanks for your help,i have learnt a lot.

now I want to konw that

how to get the power conveniently 

,like ek=Σ real*real+imag*imag

0 Kudos
Jing_Xu
Employee
1,973 Views

I don't think there is a direct MKL API function that works on this task.

0 Kudos
薛_清_
Beginner
1,973 Views

Jing X. (Intel) wrote:

I don't think there is a direct MKL API function that works on this task.

thank you ,and if I want to use fft function for many times ,  how can it take less time?

what must be repeat for every time when sue fft function  ?

       status = DftiCreateDescriptor(&hand, DFTI_SINGLE, DFTI_REAL,1, (long int)info.fftN);
       status = DftiSetValue(hand, DFTI_PLACEMENT, DFTI_NOT_INPLACE);
       status = DftiSetValue(hand, DFTI_CONJUGATE_EVEN_STORAGE, DFTI_COMPLEX_COMPLEX);
       status = DftiCommitDescriptor(hand);

To be honest,when i use mkl fft function ,my  running time of the program gets longer. I really want to know why   thank you.

void FFT(float* s, int frSize, int invert)
{
   int ii,jj,n,nn,limit,m,j,inc,i;
  float wx,wr,wpr,wpi,wi,theta;
  float xre,xri,x;
   
   //n=VectorSize(s);
   n=frSize;
   nn=n / 2; j = 1;
   for (ii=1;ii<=nn;ii++) {
      i = 2 * ii - 1;
      if (j>i) {
         xre = s; xri = s[j + 1];
         s = s;  s[j + 1] = s[i + 1];
         s = xre; s[i + 1] = xri;
      }
      m = n / 2;
      while (m >= 2  && j > m) {
         j -= m; m /= 2;
      }
      j += m;
   };
   limit = 2;
   while (limit < n) {
      inc = 2 * limit; theta = TPI / limit;
      if (invert) theta = -theta;
      x = sin(0.5 * theta);
      wpr = -2.0 * x * x; wpi = sin(theta); 
      wr = 1.0; wi = 0.0;
      for (ii=1; ii<=limit/2; ii++) {
         m = 2 * ii - 1;
         for (jj = 0; jj<=(n - m) / inc;jj++) {
            i = m + jj * inc;
            j = i + limit;
            xre = wr * s - wi * s[j + 1];
            xri = wr * s[j + 1] + wi * s;
            s = s - xre; s[j + 1] = s[i + 1] - xri;
            s = s + xre; s[i + 1] = s[i + 1] + xri;
         }
         wx = wr;
         wr = wr * wpr - wi * wpi + wr;
         wi = wi * wpr + wx * wpi + wi;
      }
      limit = inc;
   }
   if (invert)
      for (i = 1;i<=n;i++) 
         s = s / nn;   
}

/* EXPORT-> Realft: apply fft to real s */

void Realft (float* s,int fftSize)
{

   int n, n2, i, i1, i2, i3, i4;
   float xr1, xi1, xr2, xi2, wrs, wis;
   float yr, yi, yr2, yi2, yr0, theta, x;

   //n=VectorSize(s) / 2; 
   n=fftSize / 2; 
   n2 = n/2;
   theta = PI / n;


   FFT(s, fftSize, FALSE);//    3.27.22.35 XQW-MKL


   x = sin(0.5 * theta);
   yr2 = -2.0 * x * x;
   yi2 = sin(theta); yr = 1.0 + yr2; yi = yi2;
   for (i=2; i<=n2; i++) {
      i1 = i + i - 1;      i2 = i1 + 1;
      i3 = n + n + 3 - i2; i4 = i3 + 1;
      wrs = yr; wis = yi;
      xr1 = (s[i1] + s[i3])/2.0; xi1 = (s[i2] - s[i4])/2.0;
      xr2 = (s[i2] + s[i4])/2.0; xi2 = (s[i3] - s[i1])/2.0;
      s[i1] = xr1 + wrs * xr2 - wis * xi2;
      s[i2] = xi1 + wrs * xi2 + wis * xr2;
      s[i3] = xr1 - wrs * xr2 + wis * xi2;
      s[i4] = -xi1 + wrs * xi2 + wis * xr2;
      yr0 = yr;
      yr = yr * yr2 - yi  * yi2 + yr;
      yi = yi * yr2 + yr0 * yi2 + yi;
   }
   xr1 = s[1];
   s[1] = xr1 + s[2];
   s[2] = 0.0;
}

 

0 Kudos
Reply