Mercurial > hg > Members > yuuhi > OpenCL
view fft_Example/fft_base_kernels.cl @ 2:ccea4e6a1945
add OpenCL example
author | Yuhi TOMARI <yuhi@cr.ie.u-ryukyu.ac.jp> |
---|---|
date | Tue, 22 Jan 2013 23:19:41 +0900 |
parents | |
children |
line wrap: on
line source
#ifndef M_PI #define M_PI 0x1.921fb54442d18p+1 #endif #define complexMul(a,b) ((float2)(mad(-(a).y, (b).y, (a).x * (b).x), mad((a).y, (b).x, (a).x * (b).y))) #define conj(a) ((float2)((a).x, -(a).y)) #define conjTransp(a) ((float2)(-(a).y, (a).x)) #define fftKernel2(a,dir) { float2 c = (a)[0]; (a)[0] = c + (a)[1]; (a)[1] = c - (a)[1]; } define fftKernel2S(d1,d2,dir) { float2 c = (d1); (d1) = c + (d2); (d2) = c - (d2); } #define fftKernel4(a,dir) { fftKernel2S((a)[0], (a)[2], dir); fftKernel2S((a)[1], (a)[3], dir); fftKernel2S((a)[0], (a)[1], dir); (a)[3] = (float2)(dir)*(conjTransp((a)[3])); fftKernel2S((a)[2], (a)[3], dir); float2 c = (a)[1]; (a)[1] = (a)[2]; (a)[2] = c; } #define fftKernel4s(a0,a1,a2,a3,dir) { fftKernel2S((a0), (a2), dir); fftKernel2S((a1), (a3), dir); fftKernel2S((a0), (a1), dir); (a3) = (float2)(dir)*(conjTransp((a3))); fftKernel2S((a2), (a3), dir); float2 c = (a1); (a1) = (a2); (a2) = c; } #define bitreverse8(a) { float2 c; c = (a)[1]; (a)[1] = (a)[4]; (a)[4] = c; c = (a)[3]; (a)[3] = (a)[6]; (a)[6] = c; } #define fftKernel8(a,dir) { const float2 w1 = (float2)(0x1.6a09e6p-1f, dir*0x1.6a09e6p-1f); const float2 w3 = (float2)(-0x1.6a09e6p-1f, dir*0x1.6a09e6p-1f); float2 c; fftKernel2S((a)[0], (a)[4], dir); fftKernel2S((a)[1], (a)[5], dir); fftKernel2S((a)[2], (a)[6], dir); fftKernel2S((a)[3], (a)[7], dir); (a)[5] = complexMul(w1, (a)[5]); (a)[6] = (float2)(dir)*(conjTransp((a)[6])); (a)[7] = complexMul(w3, (a)[7]); fftKernel2S((a)[0], (a)[2], dir); fftKernel2S((a)[1], (a)[3], dir); fftKernel2S((a)[4], (a)[6], dir); fftKernel2S((a)[5], (a)[7], dir); (a)[3] = (float2)(dir)*(conjTransp((a)[3])); (a)[7] = (float2)(dir)*(conjTransp((a)[7])); fftKernel2S((a)[0], (a)[1], dir); fftKernel2S((a)[2], (a)[3], dir); fftKernel2S((a)[4], (a)[5], dir); fftKernel2S((a)[6], (a)[7], dir); bitreverse8((a)); } #define bitreverse4x4(a) { float2 c; c = (a)[1]; (a)[1] = (a)[4]; (a)[4] = c; c = (a)[2]; (a)[2] = (a)[8]; (a)[8] = c; c = (a)[3]; (a)[3] = (a)[12]; (a)[12] = c; c = (a)[6]; (a)[6] = (a)[9]; (a)[9] = c; c = (a)[7]; (a)[7] = (a)[13]; (a)[13] = c; c = (a)[11]; (a)[11] = (a)[14]; (a)[14] = c; } #define fftKernel16(a,dir) { const float w0 = 0x1.d906bcp-1f; const float w1 = 0x1.87de2ap-2f; const float w2 = 0x1.6a09e6p-1f; fftKernel4s((a)[0], (a)[4], (a)[8], (a)[12], dir); fftKernel4s((a)[1], (a)[5], (a)[9], (a)[13], dir); fftKernel4s((a)[2], (a)[6], (a)[10], (a)[14], dir); fftKernel4s((a)[3], (a)[7], (a)[11], (a)[15], dir); (a)[5] = complexMul((a)[5], (float2)(w0, dir*w1)); (a)[6] = complexMul((a)[6], (float2)(w2, dir*w2)); (a)[7] = complexMul((a)[7], (float2)(w1, dir*w0)); (a)[9] = complexMul((a)[9], (float2)(w2, dir*w2)); (a)[10] = (float2)(dir)*(conjTransp((a)[10])); (a)[11] = complexMul((a)[11], (float2)(-w2, dir*w2)); (a)[13] = complexMul((a)[13], (float2)(w1, dir*w0)); (a)[14] = complexMul((a)[14], (float2)(-w2, dir*w2)); (a)[15] = complexMul((a)[15], (float2)(-w0, dir*-w1)); fftKernel4((a), dir); fftKernel4((a) + 4, dir); fftKernel4((a) + 8, dir); fftKernel4((a) + 12, dir); bitreverse4x4((a)); } #define bitreverse32(a) { float2 c1, c2; c1 = (a)[2]; (a)[2] = (a)[1]; c2 = (a)[4]; (a)[4] = c1; c1 = (a)[8]; (a)[8] = c2; c2 = (a)[16]; (a)[16] = c1; (a)[1] = c2; c1 = (a)[6]; (a)[6] = (a)[3]; c2 = (a)[12]; (a)[12] = c1; c1 = (a)[24]; (a)[24] = c2; c2 = (a)[17]; (a)[17] = c1; (a)[3] = c2; c1 = (a)[10]; (a)[10] = (a)[5]; c2 = (a)[20]; (a)[20] = c1; c1 = (a)[9]; (a)[9] = c2; c2 = (a)[18]; (a)[18] = c1; (a)[5] = c2; c1 = (a)[14]; (a)[14] = (a)[7]; c2 = (a)[28]; (a)[28] = c1; c1 = (a)[25]; (a)[25] = c2; c2 = (a)[19]; (a)[19] = c1; (a)[7] = c2; c1 = (a)[22]; (a)[22] = (a)[11]; c2 = (a)[13]; (a)[13] = c1; c1 = (a)[26]; (a)[26] = c2; c2 = (a)[21]; (a)[21] = c1; (a)[11] = c2; c1 = (a)[30]; (a)[30] = (a)[15]; c2 = (a)[29]; (a)[29] = c1; c1 = (a)[27]; (a)[27] = c2; c2 = (a)[23]; (a)[23] = c1; (a)[15] = c2; } #define fftKernel32(a,dir) { fftKernel2S((a)[0], (a)[16], dir); fftKernel2S((a)[1], (a)[17], dir); fftKernel2S((a)[2], (a)[18], dir); fftKernel2S((a)[3], (a)[19], dir); fftKernel2S((a)[4], (a)[20], dir); fftKernel2S((a)[5], (a)[21], dir); fftKernel2S((a)[6], (a)[22], dir); fftKernel2S((a)[7], (a)[23], dir); fftKernel2S((a)[8], (a)[24], dir); fftKernel2S((a)[9], (a)[25], dir); fftKernel2S((a)[10], (a)[26], dir); fftKernel2S((a)[11], (a)[27], dir); fftKernel2S((a)[12], (a)[28], dir); fftKernel2S((a)[13], (a)[29], dir); fftKernel2S((a)[14], (a)[30], dir); fftKernel2S((a)[15], (a)[31], dir); (a)[17] = complexMul((a)[17], (float2)(0x1.f6297cp-1f, dir*0x1.8f8b84p-3f)); (a)[18] = complexMul((a)[18], (float2)(0x1.d906bcp-1f, dir*0x1.87de2ap-2f)); (a)[19] = complexMul((a)[19], (float2)(0x1.a9b662p-1f, dir*0x1.1c73b4p-1f)); (a)[20] = complexMul((a)[20], (float2)(0x1.6a09e6p-1f, dir*0x1.6a09e6p-1f)); (a)[21] = complexMul((a)[21], (float2)(0x1.1c73b4p-1f, dir*0x1.a9b662p-1f)); (a)[22] = complexMul((a)[22], (float2)(0x1.87de2ap-2f, dir*0x1.d906bcp-1f)); (a)[23] = complexMul((a)[23], (float2)(0x1.8f8b84p-3f, dir*0x1.f6297cp-1f)); (a)[24] = complexMul((a)[24], (float2)(0x0p+0f, dir*0x1p+0f)); (a)[25] = complexMul((a)[25], (float2)(-0x1.8f8b84p-3f, dir*0x1.f6297cp-1f)); (a)[26] = complexMul((a)[26], (float2)(-0x1.87de2ap-2f, dir*0x1.d906bcp-1f)); (a)[27] = complexMul((a)[27], (float2)(-0x1.1c73b4p-1f, dir*0x1.a9b662p-1f)); (a)[28] = complexMul((a)[28], (float2)(-0x1.6a09e6p-1f, dir*0x1.6a09e6p-1f)); (a)[29] = complexMul((a)[29], (float2)(-0x1.a9b662p-1f, dir*0x1.1c73b4p-1f)); (a)[30] = complexMul((a)[30], (float2)(-0x1.d906bcp-1f, dir*0x1.87de2ap-2f)); (a)[31] = complexMul((a)[31], (float2)(-0x1.f6297cp-1f, dir*0x1.8f8b84p-3f)); fftKernel16((a), dir); fftKernel16((a) + 16, dir); bitreverse32((a)); } __kernel void clFFT_1DTwistInterleaved(__global float2 *in, unsigned int startRow, unsigned int numCols, unsigned int N, unsigned int numRowsToProcess, int dir) { float2 a, w; float ang; unsigned int j; unsigned int i = get_global_id(0); unsigned int startIndex = i; if(i < numCols) { for(j = 0; j < numRowsToProcess; j++) { a = in[startIndex]; ang = 2.0f * M_PI * dir * i * (startRow + j) / N; w = (float2)(native_cos(ang), native_sin(ang)); a = complexMul(a, w); in[startIndex] = a; startIndex += numCols; } } } __kernel void clFFT_1DTwistSplit(__global float *in_real, __global float *in_imag , unsigned int startRow, unsigned int numCols, unsigned int N, unsigned int numRowsToProcess, int dir) { float2 a, w; float ang; unsigned int j; unsigned int i = get_global_id(0); unsigned int startIndex = i; if(i < numCols) { for(j = 0; j < numRowsToProcess; j++) { a = (float2)(in_real[startIndex], in_imag[startIndex]); ang = 2.0f * M_PI * dir * i * (startRow + j) / N; w = (float2)(native_cos(ang), native_sin(ang)); a = complexMul(a, w); in_real[startIndex] = a.x; in_imag[startIndex] = a.y; startIndex += numCols; } } }