//+------------------------------------------------------------------+ //| FFT.mq5 | //| Copyright 2000-2025, MetaQuotes Ltd. | //| https://www.mql5.com | //+------------------------------------------------------------------+ #property copyright "Copyright 2000-2025, MetaQuotes Ltd." #property link "https://www.mql5.com" #property version "1.00" #include #include #resource "Kernels/fft.cl" as string cl_program #define kernel_init "fft_init" #define kernel_stage "fft_stage" #define kernel_scale "fft_scale" #define NUM_POINTS 16384 #define FFT_DIRECTION 1 //+------------------------------------------------------------------+ //| Fast Fourier transform and its inverse (both recursively) | //| Copyright (C) 2004, Jerome R. Breitenbach. All rights reserved. | //| Reference: | //| Matthew Scarpino, "OpenCL in Action: How to accelerate graphics | //| and computations", Manning, 2012, Chapter 14. | //+------------------------------------------------------------------+ //| Recursive direct FFT transform | //+------------------------------------------------------------------+ void fft(const int N,float &x_real[],float &x_imag[],float &X_real[],float &X_imag[]) { //--- prepare temporary arrays float XX_real[],XX_imag[]; ArrayResize(XX_real,N); ArrayResize(XX_imag,N); //--- calculate FFT by a recursion fft_rec(N,0,1,x_real,x_imag,X_real,X_imag,XX_real,XX_imag); } //+------------------------------------------------------------------+ //| Recursive inverse FFT transform | //+------------------------------------------------------------------+ void ifft(const int N,float &x_real[],float &x_imag[],float &X_real[],float &X_imag[]) { int N2=N/2; // half the number of points in IFFT //--- calculate IFFT via reciprocity property of DFT fft(N,X_real,X_imag,x_real,x_imag); x_real[0]=x_real[0]/N; x_imag[0]=x_imag[0]/N; x_real[N2]=x_real[N2]/N; x_imag[N2]=x_imag[N2]/N; for(int i=1; i0) fft(N,data_real,data_imag,XX_real,XX_imag); else ifft(N,XX_real,XX_imag,data_real,data_imag); //--- CPU calculation finished time_cpu=ulong((GetMicrosecondCount()-time_cpu)); //--- copy calculated data ArrayCopy(data_real,XX_real,0,0,WHOLE_ARRAY); ArrayCopy(data_imag,XX_imag,0,0,WHOLE_ARRAY); //--- return(true); } //+------------------------------------------------------------------+ //| FFT_GPU | //+------------------------------------------------------------------+ bool FFT_GPU(int direction,int power,float &data_real[],float &data_imag[],ulong &time_gpu) { //--- calculate the number of points int num_points=1; for(int i=0;inum_points) points_per_group=num_points; //--- set kernel arguments OpenCL.SetArgumentBuffer(0,0,0); OpenCL.SetArgumentBuffer(0,1,1); OpenCL.SetArgumentLocalMemory(0,2,local_mem_size); OpenCL.SetArgument(0,3,points_per_group); OpenCL.SetArgument(0,4,num_points); OpenCL.SetArgument(0,5,direction); //--- OpenCL execute settings int task_dimension=1; uint global_size=(uint)((num_points/points_per_group)*local_size); uint global_work_offset[1]={0}; uint global_work_size[1]; global_work_size[0]=global_size; uint local_work_size[1]; local_work_size[0]=local_size; //--- GPU calculation start time_gpu=GetMicrosecondCount(); //-- execute kernel fft_init if(!OpenCL.Execute(0,task_dimension,global_work_offset,global_work_size,local_work_size)) { PrintFormat("fft_init: Error in CLExecute. Error code=%d",GetLastError()); return(false); } //-- further stages of the FFT if(num_points>points_per_group) { //--- set arguments for kernel 1 OpenCL.SetArgumentBuffer(1,0,1); OpenCL.SetArgument(1,2,points_per_group); OpenCL.SetArgument(1,3,direction); for(int stage=2; stage<=num_points/points_per_group; stage<<=1) { OpenCL.SetArgument(1,1,stage); //-- execute kernel fft_stage if(!OpenCL.Execute(1,task_dimension,global_work_offset,global_work_size,local_work_size)) { PrintFormat("fft_stage: Error in CLExecute. Error code=%d",GetLastError()); return(false); } } } //--- scale values if performing the inverse FFT if(direction<0) { OpenCL.SetArgumentBuffer(2,0,1); OpenCL.SetArgument(2,1,points_per_group); OpenCL.SetArgument(2,2,num_points); //-- execute kernel fft_scale if(!OpenCL.Execute(2,task_dimension,global_work_offset,global_work_size,local_work_size)) { PrintFormat("fft_scale: Error in CLExecute. Error code=%d",GetLastError()); return(false); } } //--- read the results from GPU memory if(!OpenCL.BufferRead(1,data,0,0,2*num_points)) { PrintFormat("Error in BufferRead for data_buffer2. Error code=%d",GetLastError()); return(false); } //--- GPU calculation finished time_gpu=ulong((GetMicrosecondCount()-time_gpu)); //--- copy calculated data and release OpenCL handles for(int i=0; i