PageRenderTime 12ms CodeModel.GetById 8ms app.highlight 1ms RepoModel.GetById 1ms app.codeStats 0ms

/WWDirectComputeTest/SincConvolution.hlsl

http://bitspersampleconv2.googlecode.com/
text | 329 lines | 267 code | 62 blank | 0 comment | 0 complexity | be779ff4a4f5ed8927baf75a86213629 MD5 | raw file
  1// 日本語UTF-8
  2
  3/*
  4
  5OutputBuffer[t+convolutionN] = Σ[sample[t+x] * sinc(πx + XBuffer[t])]
  6CONV_START <= x < CONV_END
  7を計算する
  8
  9convolutionN = 256
 10sampleN = 100
 11の場合
 12
 13CONV_START = -256
 14CONV_END   =  256
 15CONV_COUNT =  512
 16SAMPLE_N   =  100
 17GROUP_THREAD_COUNT 2の乗数
 18を#defineしてCS5.0 DirectCompute シェーダーとしてコンパイルする。
 19
 20// シェーダー定数を渡す
 21shaderParams.c_convOffs = 0
 22shaderParams.c_dispatchCount = convolutionN*2/GROUP_THREAD_COUNT;
 23ComputeShaderのrun(shaderParams, sampleN, 1, 1);
 24
 25する。
 26
 27用意するデータ
 28
 29①SampleDataBuffer…前後を水増しされたサンプルデータsample[t]
 30SampleDataBuffer[0]~SampleDataBuffer[convolutionN-1]…0を詰める
 31SampleDataBuffer[convolutionN]~SampleDataBuffer[convolutionN + sampleN-1]…サンプルデータsample[t]
 32SampleDataBuffer[convolutionN+SampleN]~SampleDataBuffer[convolutionN*2 + sampleN-1]…0を詰める
 33
 34②SinxBuffer リサンプル地点のsin(x) 適当に作る
 35SinxBuffer[0]~SinxBuffer[sampleN-1] sin(x)の値
 36
 37③XBuffer リサンプル地点x
 38XBuffer[0]~XBuffer[sampleN-1] xの値
 39
 40④出力バッファー
 41OutputBuffer[0]~OutputBuffer[sampleN-1]
 42OutputBuffer[]はsampleN個用意する
 43
 44*/
 45#ifdef HIGH_PRECISION
 46// 主にdouble精度
 47
 48StructuredBuffer<float>   g_SampleDataBuffer : register(t0);
 49StructuredBuffer<double>  g_SinxBuffer       : register(t1);
 50StructuredBuffer<float>   g_XBuffer          : register(t2);
 51RWStructuredBuffer<float> g_OutputBuffer     : register(u0);
 52
 53/// 定数。16バイトの倍数のサイズの構造体。
 54cbuffer consts {
 55    /// 畳み込み要素オフセット値。n * GROUP_THREAD_COUNTの飛び飛びの値が渡る。
 56    uint c_convOffs;
 57    /// Dispatch繰り返し回数。
 58    uint c_dispatchCount;
 59    uint c_reserved1;
 60    uint c_reserved2;
 61};
 62
 63inline double
 64SincF(double sinx, float x)
 65{
 66    if (-0.000000001f < x && x < 0.000000001f) {
 67        return 1.0;
 68    } else {
 69        // 割り算ができないので、ここで精度落ちる。残念。
 70        return sinx * rcp(x);
 71    }
 72}
 73
 74#define PI_F 3.141592653589793238462643f
 75
 76// TGSM
 77groupshared double s_scratch[GROUP_THREAD_COUNT];
 78groupshared double s_sinX;
 79groupshared float  s_xOffs;
 80
 81/// 畳み込み計算要素1回実行。
 82/// sample[t+x] * sinc(πx + XBuffer[t])
 83inline double
 84ConvolutionElemValue(uint pos, uint convOffs)
 85{
 86    const int offs = c_convOffs + convOffs;
 87    const float x = mad(PI_F, offs + CONV_START, s_xOffs);
 88    return ((double)g_SampleDataBuffer[offs + pos]) * SincF(s_sinX, x);
 89}
 90
 91// スレッドグループとTGSMを使用して、GPUメモリからの読み出し回数を減らす最適化。
 92
 93// groupIdXYZはDispatch()のパラメータXYZ=(nx,1,1)の場合(0,0,0)~(nx-1, 0, 0)。
 94// スレッドグループが作られ、tid==0~groupDim_x-1までのtidを持ったスレッドが同時に走る。
 95[numthreads(GROUP_THREAD_COUNT, 1, 1)]
 96void
 97CSMain(
 98        uint  tid:        SV_GroupIndex,
 99        uint3 groupIdXYZ: SV_GroupID)
100{
101    uint offs = tid;
102
103    if (tid == 0) {
104        s_xOffs = g_XBuffer[groupIdXYZ.x];
105        s_sinX  = g_SinxBuffer[groupIdXYZ.x];
106    }
107    s_scratch[tid] = 0;
108
109    GroupMemoryBarrierWithGroupSync();
110
111    do {
112        s_scratch[tid] +=
113            ConvolutionElemValue(groupIdXYZ.x, offs) +
114            ConvolutionElemValue(groupIdXYZ.x, offs + GROUP_THREAD_COUNT);
115        offs += GROUP_THREAD_COUNT * 2;
116    } while (offs < CONV_COUNT);
117
118    GroupMemoryBarrierWithGroupSync();
119
120#if 1024 <= GROUP_THREAD_COUNT
121    if (tid < 512) { s_scratch[tid] += s_scratch[tid + 512]; }
122    GroupMemoryBarrierWithGroupSync();
123#endif
124
125#if 512 <= GROUP_THREAD_COUNT
126    if (tid < 256) { s_scratch[tid] += s_scratch[tid + 256]; }
127    GroupMemoryBarrierWithGroupSync();
128#endif
129
130#if 256 <= GROUP_THREAD_COUNT
131    if (tid < 128) { s_scratch[tid] += s_scratch[tid + 128]; }
132    GroupMemoryBarrierWithGroupSync();
133#endif
134
135#if 128 <= GROUP_THREAD_COUNT
136    if (tid < 64) { s_scratch[tid] += s_scratch[tid + 64]; }
137    GroupMemoryBarrierWithGroupSync();
138#endif
139
140#if 64 <= GROUP_THREAD_COUNT
141    if (tid < 32) { s_scratch[tid] += s_scratch[tid + 32]; }
142    //GroupMemoryBarrierWithGroupSync(); // これ以降要らないらしい。2260_GTC2010.pdf参照。
143#endif
144
145#if 32 <= GROUP_THREAD_COUNT
146    if (tid < 16) { s_scratch[tid] += s_scratch[tid + 16]; }
147    //GroupMemoryBarrierWithGroupSync();
148#endif
149
150#if 16 <= GROUP_THREAD_COUNT
151    if (tid < 8) { s_scratch[tid] += s_scratch[tid + 8]; }
152    //GroupMemoryBarrierWithGroupSync();
153#endif
154
155#if 8 <= GROUP_THREAD_COUNT
156    if (tid < 4) { s_scratch[tid] += s_scratch[tid + 4]; }
157   // GroupMemoryBarrierWithGroupSync();
158#endif
159
160#if 4 <= GROUP_THREAD_COUNT
161    if (tid < 2) { s_scratch[tid] += s_scratch[tid + 2]; }
162    //GroupMemoryBarrierWithGroupSync();
163#endif
164
165    if (tid == 0) {
166        s_scratch[0] += s_scratch[1];
167        g_OutputBuffer[groupIdXYZ.x] = (float)s_scratch[0];
168    }
169}
170
171#else
172
173// 主にfloat精度
174
175StructuredBuffer<float>   g_SampleDataBuffer : register(t0);
176StructuredBuffer<float>   g_SinxBuffer       : register(t1);
177StructuredBuffer<float>   g_XBuffer          : register(t2);
178RWStructuredBuffer<float> g_OutputBuffer     : register(u0);
179
180/// 定数。16バイトの倍数のサイズの構造体。
181cbuffer consts {
182    /// 畳み込み要素オフセット値。n * GROUP_THREAD_COUNTの飛び飛びの値が渡る。
183    uint c_convOffs;
184    /// Dispatch繰り返し回数。
185    uint c_dispatchCount;
186    uint c_reserved1;
187    uint c_reserved2;
188};
189
190inline float
191SincF(float sinx, float x)
192{
193    if (-0.000000001f < x && x < 0.000000001f) {
194        return 1.0f;
195    } else {
196        // どちらでも同じだった。
197#if 1
198        return sinx * rcp(x);
199#else
200        return sinx / x;
201#endif
202    }
203}
204
205#define PI_F 3.141592653589793238462643f
206
207// TGSM
208groupshared float s_scratch[GROUP_THREAD_COUNT];
209groupshared float s_sinX;
210groupshared float s_xOffs;
211
212/// 畳み込み計算要素1回実行。
213/// sample[t+x] * sinc(πx + XBuffer[t])
214inline float
215ConvolutionElemValue(uint pos, uint convOffs)
216{
217    const int offs = c_convOffs + convOffs;
218    const float x = mad(PI_F, offs + CONV_START, s_xOffs);
219    return g_SampleDataBuffer[offs + pos] * SincF(s_sinX, x);
220}
221
222// スレッドグループとTGSMを使用して、GPUメモリからの読み出し回数を減らす最適化。
223
224// groupIdXYZはDispatch()のパラメータXYZ=(nx,1,1)の場合(0,0,0)~(nx-1, 0, 0)。
225// スレッドグループが作られ、tid==0~groupDim_x-1までのtidを持ったスレッドが同時に走る。
226[numthreads(GROUP_THREAD_COUNT, 1, 1)]
227void
228CSMain(
229        uint  tid:        SV_GroupIndex,
230        uint3 groupIdXYZ: SV_GroupID)
231{
232    uint offs = tid;
233
234    if (tid == 0) {
235        s_xOffs = g_XBuffer[groupIdXYZ.x];
236#if 1
237        // 計算精度良好。
238        s_sinX  = g_SinxBuffer[groupIdXYZ.x];
239#else
240        // こうすると精度が落ちる。GPUのsin()の精度に問題あり。
241        s_sinX  = sin(s_xOffs);
242#endif
243    }
244    s_scratch[tid] = 0;
245
246    GroupMemoryBarrierWithGroupSync();
247
248    do {
249        s_scratch[tid] +=
250            ConvolutionElemValue(groupIdXYZ.x, offs) +
251            ConvolutionElemValue(groupIdXYZ.x, offs + GROUP_THREAD_COUNT);
252        offs += GROUP_THREAD_COUNT * 2;
253    } while (offs < CONV_COUNT);
254
255    GroupMemoryBarrierWithGroupSync();
256
257#if 1024 <= GROUP_THREAD_COUNT
258    if (tid < 512) { s_scratch[tid] += s_scratch[tid + 512]; }
259    GroupMemoryBarrierWithGroupSync();
260#endif
261
262#if 512 <= GROUP_THREAD_COUNT
263    if (tid < 256) { s_scratch[tid] += s_scratch[tid + 256]; }
264    GroupMemoryBarrierWithGroupSync();
265#endif
266
267#if 256 <= GROUP_THREAD_COUNT
268    if (tid < 128) { s_scratch[tid] += s_scratch[tid + 128]; }
269    GroupMemoryBarrierWithGroupSync();
270#endif
271
272#if 128 <= GROUP_THREAD_COUNT
273    if (tid < 64) { s_scratch[tid] += s_scratch[tid + 64]; }
274    GroupMemoryBarrierWithGroupSync();
275#endif
276
277#if 64 <= GROUP_THREAD_COUNT
278    if (tid < 32) { s_scratch[tid] += s_scratch[tid + 32]; }
279    //GroupMemoryBarrierWithGroupSync(); // これ以降要らないらしい。2260_GTC2010.pdf参照。
280#endif
281
282#if 32 <= GROUP_THREAD_COUNT
283    if (tid < 16) { s_scratch[tid] += s_scratch[tid + 16]; }
284    //GroupMemoryBarrierWithGroupSync();
285#endif
286
287#if 16 <= GROUP_THREAD_COUNT
288    if (tid < 8) { s_scratch[tid] += s_scratch[tid + 8]; }
289    //GroupMemoryBarrierWithGroupSync();
290#endif
291
292#if 8 <= GROUP_THREAD_COUNT
293    if (tid < 4) { s_scratch[tid] += s_scratch[tid + 4]; }
294   // GroupMemoryBarrierWithGroupSync();
295#endif
296
297#if 4 <= GROUP_THREAD_COUNT
298    if (tid < 2) { s_scratch[tid] += s_scratch[tid + 2]; }
299    //GroupMemoryBarrierWithGroupSync();
300#endif
301
302    if (tid == 0) {
303        s_scratch[0] += s_scratch[1];
304        g_OutputBuffer[groupIdXYZ.x] = s_scratch[0];
305    }
306}
307
308#if 0
309// 最適化前
310[numthreads(1, 1, 1)]
311void
312CSMain(uint3 groupIdXYZ  : SV_GroupID,
313       uint threadIdx : SV_GroupIndex)
314{
315    int   i;
316    float sinx  = SinxBuffer[c_pos];
317    float xOffs = XBuffer[c_pos];
318    float r = 0.0f;
319
320    for (i=CONV_START; i<CONV_END; ++i) {
321        float x = mad(PI, i, xOffs);
322        r = mad(SampleDataBuffer[c_pos+i+CONV_N], SincF(sinx, x), r);
323    }
324
325    OutputBuffer[c_pos] = r;
326}
327#endif // before optimization
328
329#endif // HIGH_PRECISION