/normxcorr/trunk

To get this branch, use:
bzr branch http://suren.me/webbzr/normxcorr/trunk

« back to all changes in this revision

Viewing changes to cuda/normxcorr_hw.cu

  • Committer: Suren A. Chilingaryan
  • Date: 2009-12-06 01:52:56 UTC
  • Revision ID: csa@dside.dyndns.org-20091206015256-evn0sne8d18ovm8o
A little more computations are moved to CUDA

Show diffs side-by-side

added added

removed removed

Lines of Context:
cuda/normxcorr_hw.cu
10
10
#include "normxcorr_hw_msg.h"
11
11
#include "normxcorr_hw_kernel.cu"
12
12
 
 
13
#define max2(a,b) ((a>b)?a:b)
 
14
#define min2(a,b) ((a<b)?a:b)
 
15
 
13
16
 
14
17
static TProcessingState *pstate = NULL;
15
18
 
16
19
static void fftFree(TProcessingState *ps) {
17
20
    if (ps->cuda_base_buffer) {
18
 
 
19
21
        cudaFree(ps->cuda_lsum_temp);
20
22
        
21
23
        cudaFree(ps->cuda_lsum_buffer);
27
29
        cudaFree(ps->cuda_data_buffer);
28
30
        cudaFree(ps->cuda_base_buffer);
29
31
        cudaFree(ps->cuda_input_buffer);
 
32
        cudaFreeHost(ps->input_buffer);
 
33
        
 
34
        cudaFree(ps->cuda_cp);
 
35
 
 
36
        cudaFreeHost(ps->data_x);
 
37
        cudaFreeHost(ps->data_y);
30
38
        
31
39
        ps->cuda_base_buffer = NULL;
32
40
    }
55
63
    cufftResult cufft_err;
56
64
    cudaError cuda_err;
57
65
 
 
66
    int size;
58
67
    int lsum_alloc_size2 = ps->lsum_alloc_size * ps->lsum_alloc_size;
 
68
    int side_alloc_size2 = ps->side_alloc_size * ps->side_alloc_size;
59
69
    
60
70
    fftFree(ps);
61
71
 
98
108
        return ERROR_CUDA_MALLOC;
99
109
    }
100
110
 
101
 
    cuda_err = cudaMalloc((void**)&ps->cuda_data_buffer, ps->fft_alloc_size * sizeof(cufftComplex));
 
111
 
 
112
    size = max2(
 
113
        ps->fft_alloc_size * sizeof(cufftComplex),              /* FFT multiplication */
 
114
        2 * CP_BLOCK * ps->side_alloc_size * sizeof(int)        /* Sum, Std computations */
 
115
    );
 
116
 
 
117
    cuda_err = cudaMalloc((void**)&ps->cuda_data_buffer, size);
102
118
    if (cuda_err) {
103
119
        reportError("Device memory allocation of %u*cufftComplex bytes for cuda_data_buffer is failed", ps->fft_alloc_size);
104
120
        fftFree(ps);
105
121
        return ERROR_CUDA_MALLOC;
106
122
    }
107
123
 
108
 
    cuda_err = cudaMalloc((void**)&ps->cuda_input_buffer, ps->ncp * ps->fft_alloc_size * sizeof(uint8_t));
109
 
    if (cuda_err) {
110
 
        reportError("Device memory allocation of %u*%u*uint8 bytes for cuda_input_buffer is failed", ps->ncp, ps->fft_alloc_size);
 
124
    cuda_err = cudaHostAlloc((void**)&ps->data_x, ps->ncp * sizeof(float), 0);
 
125
    if (!cuda_err) cuda_err = cudaHostAlloc((void**)&ps->data_y, ps->ncp * sizeof(float), 0);
 
126
    if (cuda_err) {
 
127
        reportError("Host memory allocation of 2*%u*float bytes for control points is failed", ps->ncp);
 
128
        fftFree(ps);
 
129
        return ERROR_CUDA_MALLOC;
 
130
    }
 
131
 
 
132
    cuda_err = cudaMalloc((void**)&ps->cuda_cp, 2 * CP_BLOCK * sizeof(float));
 
133
    if (cuda_err) {
 
134
        reportError("Device memory allocation of %u*%u*uint8 bytes for cuda_input_buffer is failed", ps->ncp, side_alloc_size2);
 
135
        fftFree(ps);
 
136
        return ERROR_CUDA_MALLOC;
 
137
    }
 
138
 
 
139
    cuda_err = cudaMalloc((void**)&ps->cuda_input_buffer, ps->ncp * side_alloc_size2 * sizeof(uint8_t));
 
140
    if (cuda_err) {
 
141
        reportError("Device memory allocation of %u*%u*uint8 bytes for cuda_input_buffer is failed", ps->ncp, side_alloc_size2);
 
142
        fftFree(ps);
 
143
        return ERROR_CUDA_MALLOC;
 
144
    }
 
145
 
 
146
    cuda_err = cudaHostAlloc((void**)&ps->input_buffer, ps->ncp * ps->fft_alloc_size * sizeof(uint8_t), cudaHostAllocWriteCombined);
 
147
    if (cuda_err) {
 
148
        reportError("Host memory allocation of %u*%u*uint8 bytes for input_buffer is failed", ps->ncp, ps->fft_alloc_size);
111
149
        fftFree(ps);
112
150
        return ERROR_CUDA_MALLOC;
113
151
    }
120
158
    }
121
159
    cudaMemset((void*)ps->cuda_final_buffer, 0, ps->ncp * ps->fft_alloc_size * sizeof(float));
122
160
 
123
 
    cuda_err = cudaMalloc((void**)&ps->cuda_result_buffer, ps->fft_alloc_size * sizeof(cufftReal));
 
161
    cuda_err = cudaMalloc((void**)&ps->cuda_result_buffer, ps->ncp*ps->fft_alloc_size * sizeof(cufftReal));
124
162
    if (cuda_err) {
125
163
        reportError("Device memory allocation of %u*cufftReal bytes for cuda_result_buffer is failed", ps->fft_alloc_size);
126
164
        fftFree(ps);
204
242
    
205
243
    int size = ps->fft_size;
206
244
    int alloc_size = ps->fft_alloc_size;
207
 
    
 
245
 
208
246
    int N = mxGetM(data);
209
247
    int N2 = N * N;
210
248
 
 
249
    int side_alloc = ps->side_alloc_size;
 
250
    int side_alloc2 = side_alloc * side_alloc;
 
251
 
211
252
    dim3 input_block_dim(N, 1, 1);
212
253
    dim3 input_grid_dim(N, 1, 1);
213
254
 
214
 
    uint8_t *cudaInputPtr = ps->cuda_input_buffer + icp * alloc_size;
 
255
    uint8_t *cudaInputPtr = ps->cuda_input_buffer + icp * side_alloc2;
215
256
    cufftComplex *cudaPtr = ps->cuda_base_buffer + icp * alloc_size;
216
257
    cufftReal *cudaRealPtr = ps->cuda_temp_buffer;
217
258
 
241
282
    return cudaPtr;
242
283
}
243
284
 
244
 
static inline mxArray *fftCompute(TProcessingState *ps, int icp, const mxArray *data, float sum, float denom) {
245
 
    uint8_t *dataPtr;
246
 
    double *ar;
247
 
    mxArray *res;
 
285
 
 
286
/*
 
287
static int fftComputeFragment(TProcessingState *ps, int icp, const mxArray *data, float sum, float denom) {
 
288
//    uint8_t *dataPtr;
 
289
//    double *ar;
 
290
//    mxArray *res;
248
291
 
249
292
    int size = ps->fft_size;
250
293
    int size2 = size * size;
251
294
    int alloc_size = ps->fft_alloc_size;
252
295
 
253
 
    int N = mxGetM(data);
 
296
    int half_size = ps->corr_size;
 
297
    int N = 2 * half_size + 1;
 
298
 
 
299
//    int N = mxGetM(data);
254
300
    int N2 = N * N;
255
301
 
 
302
 
256
303
    dim3 input_block_dim(N, 1, 1);
257
304
    dim3 input_grid_dim(N, 1, 1);
258
305
 
262
309
    dim3 output_block_dim(size, 1, 1);
263
310
    dim3 output_grid_dim(size, 1, 1);
264
311
 
265
 
    uint8_t *cudaInputPtr = ps->cuda_input_buffer + icp * alloc_size;
266
 
    cufftComplex *cudaPtr = ps->cuda_base_buffer + icp * alloc_size;
 
312
//    uint8_t *cudaInputPtr = ps->cuda_input_buffer + icp * alloc_size;
 
313
//    cufftComplex *cudaPtr = ps->cuda_base_buffer + icp * alloc_size;
267
314
    cufftReal *cudaRealPtr = ps->cuda_temp_buffer;
268
 
    cufftComplex *cudaDataPtr = ps->cuda_data_buffer;
269
 
    float *cudaResultPtr = ps->cuda_final_buffer;
 
315
//    cufftComplex *cudaDataPtr = ps->cuda_data_buffer;
 
316
    float *cudaResultPtr = ps->cuda_final_buffer + icp * alloc_size;
 
317
 
270
318
 
271
319
    dataPtr = (uint8_t*)mxGetData(data);
272
320
    cudaMemcpy(cudaInputPtr, dataPtr, N2*sizeof(uint8_t), cudaMemcpyHostToDevice);
273
 
    vecPack<<<input_grid_dim, input_block_dim>>>(cudaRealPtr, size, cudaInputPtr, N);
 
321
 
 
322
   vecPack<<<input_grid_dim, input_block_dim>>>(cudaRealPtr, size, cudaInputPtr, N);
274
323
 
275
324
    cufftExecR2C(ps->cufft_r2c_plan, cudaRealPtr, cudaDataPtr);
276
325
 
279
328
    cudaRealPtr = ps->cuda_result_buffer;
280
329
    cufftExecC2R(ps->cufft_c2r_plan, cudaDataPtr, cudaRealPtr);
281
330
 
282
 
    float *cudaDenom = ps->cuda_denom_buffer + icp*alloc_size;
283
 
    float *cudaLSum = ps->cuda_lsum_buffer + icp*alloc_size;
284
 
 
285
 
    vecCompute<<<output_grid_dim, output_block_dim>>>(
286
 
            cudaResultPtr,
287
 
            cudaRealPtr, 1./(size2 * (N2 - 1)),
288
 
            cudaLSum, sum / (N2 * (N2 - 1)),
289
 
            cudaDenom, denom,
290
 
            size
291
 
    );
292
 
    
293
 
    res = mxCreateNumericMatrix(size, size, mxSINGLE_CLASS, mxREAL);
294
 
    ar = mxGetPr(res);
 
331
 
 
332
    cudaRealPtr = ps->cuda_result_buffer + icp*alloc_size;
 
333
 
 
334
    float *cudaDenom = ps->cuda_denom_buffer + icp*alloc_size;
 
335
    float *cudaLSum = ps->cuda_lsum_buffer + icp*alloc_size;
 
336
 
 
337
    vecCompute<<<output_grid_dim, output_block_dim>>>(
 
338
            cudaResultPtr,
 
339
            cudaRealPtr, 1./(size2 * (N2 - 1)),
 
340
            cudaLSum, sum / (N2 * (N2 - 1)),
 
341
            cudaDenom, denom,
 
342
            size
 
343
    );
 
344
*/
 
345
/*
 
346
    int size = ps->fft_size;
 
347
    int size2 = size * size;
 
348
    int alloc_size = ps->fft_alloc_size;
 
349
 
 
350
    int half_size = ps->corr_size;
 
351
    int N = 2 * half_size + 1;
 
352
    int N2 = N * N;
 
353
 
 
354
    dim3 output_block_dim(size, 1, 1);
 
355
    dim3 output_grid_dim(size, 1, 1);
 
356
 
 
357
    cufftReal *cudaRealPtr = ps->cuda_result_buffer + icp*alloc_size;
 
358
    float *cudaResultPtr = ps->cuda_final_buffer + icp * alloc_size;
 
359
    float *cudaDenom = ps->cuda_denom_buffer + icp*alloc_size;
 
360
    float *cudaLSum = ps->cuda_lsum_buffer + icp*alloc_size;
 
361
 
 
362
    vecCompute<<<output_grid_dim, output_block_dim>>>(
 
363
            cudaResultPtr,
 
364
            cudaRealPtr, 1./(size2 * (N2 - 1)),
 
365
            cudaLSum, sum / (N2 * (N2 - 1)),
 
366
            cudaDenom, denom,
 
367
            size
 
368
    );
 
369
    return 0;
 
370
}
 
371
*/
 
372
 
 
373
static inline mxArray *fftCompute(TProcessingState *ps, int icp) {
 
374
    int size = ps->fft_size;
 
375
    int size2 = size * size;
 
376
    int alloc_size = ps->fft_alloc_size;
 
377
    float *cudaResultPtr = ps->cuda_final_buffer + icp * alloc_size;
 
378
 
 
379
//    fftComputeFragment(ps, icp, data, sum, denom);
 
380
 
 
381
    mxArray *res = mxCreateNumericMatrix(size, size, mxSINGLE_CLASS, mxREAL);
 
382
    float *ar = (float*)mxGetPr(res);
295
383
 
296
384
    cudaMemcpy(ar, cudaResultPtr, size2*sizeof(cufftReal), cudaMemcpyDeviceToHost);
297
385
 
298
386
    return res;
299
387
}
300
388
 
 
389
int fftLoadFragment(TProcessingState *ps, int icp, int ncp, const mxArray *image) {
 
390
    int width = mxGetN(image);
 
391
    int height = mxGetM(image);
 
392
 
 
393
    int half_size = ps->corr_size;
 
394
    int size = 2 * half_size + 1;
 
395
    int size2 = size * size;
 
396
 
 
397
    int fft_size = ps->fft_size;
 
398
    int fft_size2 = fft_size * fft_size;
 
399
    int alloc_size = ps->fft_alloc_size;
 
400
    int side_alloc = ps->side_alloc_size;
 
401
    int side_alloc2 = side_alloc * side_alloc;
 
402
 
 
403
    uint8_t *fullimg = ((uint8_t*)mxGetData(image));
 
404
    uint8_t *img = ps->input_buffer;
 
405
 
 
406
    cufftComplex *cudaDataPtr = (cufftComplex*)ps->cuda_data_buffer;
 
407
    cufftReal *cudaRealPtr = ps->cuda_temp_buffer;
 
408
 
 
409
    dim3 input_block_dim(size, 1, 1);
 
410
    dim3 input_grid_dim(size, 1, 1);
 
411
    dim3 block_dim(fft_size / 2 + 1, 1, 1);
 
412
    dim3 grid_dim(fft_size, 1, 1);
 
413
 
 
414
    for (int i = 0;i < ncp;i++) {
 
415
        float x = ps->data_x[i+icp] - 1;
 
416
        float y = ps->data_y[i+icp] - 1;
 
417
    
 
418
        int xstart = roundf(x) - half_size;
 
419
        int ystart = roundf(y) - half_size;
 
420
    
 
421
        int xend = xstart + size;
 
422
        int yend = xstart + size;
 
423
 
 
424
        if ((xstart < 0)||(ystart < 0)||(xend >= width)||(yend >= height)) {
 
425
                // Somehow mark we have skipped it
 
426
            continue;
 
427
        }
 
428
 
 
429
        cudaMemcpy2D(
 
430
            img,
 
431
            size * sizeof(uint8_t),
 
432
            fullimg + (xstart * height + ystart),
 
433
            height * sizeof(uint8_t),
 
434
            size * sizeof(uint8_t),
 
435
            size,
 
436
            cudaMemcpyHostToHost
 
437
        );
 
438
 
 
439
 
 
440
        cufftComplex *cudaBasePtr = ps->cuda_base_buffer + (i+icp) * alloc_size;
 
441
        cufftReal *cudaResultPtr = ps->cuda_result_buffer + (i+icp) * alloc_size;
 
442
        
 
443
        uint8_t *cudaInputPtr = ps->cuda_input_buffer + i*side_alloc2;
 
444
 
 
445
        cudaMemcpy2D(
 
446
            cudaInputPtr, side_alloc * sizeof(uint8_t),
 
447
            img, size * sizeof(uint8_t),
 
448
            size * sizeof(uint8_t), size, cudaMemcpyHostToDevice
 
449
        );
 
450
 
 
451
        vecPack<<<input_grid_dim, input_block_dim>>>(cudaRealPtr, fft_size, cudaInputPtr, side_alloc, size);
 
452
 
 
453
        cufftExecR2C(ps->cufft_r2c_plan, cudaRealPtr, cudaDataPtr);
 
454
 
 
455
        vecMul<<<grid_dim,block_dim>>>(cudaDataPtr, cudaBasePtr, fft_size/2+1);
 
456
 
 
457
        cufftExecC2R(ps->cufft_c2r_plan, cudaDataPtr, cudaResultPtr);
 
458
    }
 
459
 
 
460
 
 
461
    int cp_blocks, side_blocks;
 
462
    if (ncp%CP_BLOCK_SIZE) cp_blocks = (ncp / CP_BLOCK_SIZE) + 1;
 
463
    else cp_blocks = ncp / CP_BLOCK_SIZE;
 
464
 
 
465
    if (size%SIDE_BLOCK_SIZE) side_blocks = (size / SIDE_BLOCK_SIZE) + 1;
 
466
    else side_blocks = size / SIDE_BLOCK_SIZE;
 
467
 
 
468
 
 
469
    int *stat_buf = (int*)ps->cuda_data_buffer;
 
470
 
 
471
    float *sumbuf = ps->cuda_cp;// + 2*ps->ncp + icp;
 
472
    float *stdbuf = ps->cuda_cp + CP_BLOCK;
 
473
 
 
474
    dim3 joint_block_dim(SIDE_BLOCK_SIZE, CP_BLOCK_SIZE, 1);
 
475
    dim3 joint_grid_dim(side_blocks, cp_blocks, 1);
 
476
    
 
477
    stat1<<<joint_grid_dim, joint_block_dim>>>(stat_buf, stat_buf + side_alloc * CP_BLOCK, ps->cuda_input_buffer, side_alloc2, side_alloc, size);
 
478
 
 
479
    if (ncp%BLOCK_SIZE_1D) cp_blocks = (ncp / BLOCK_SIZE_1D) + 1;
 
480
    else cp_blocks = ncp / BLOCK_SIZE_1D;
 
481
 
 
482
    stat2<<<cp_blocks, BLOCK_SIZE_1D>>>(sumbuf, stdbuf, stat_buf, stat_buf + side_alloc * CP_BLOCK, size);
 
483
 
 
484
/*
 
485
    float *cp = (float*)malloc(2*CP_BLOCK*sizeof(float));
 
486
    cudaMemcpy(cp, ps->cuda_cp, 2*CP_BLOCK*sizeof(float), cudaMemcpyDeviceToHost);
 
487
 
 
488
    float *sumbuf1 = cp;
 
489
    float *stdbuf1 = cp + CP_BLOCK;
 
490
    free(cp);
 
491
*/
 
492
 
 
493
    dim3 output_block_dim(fft_size, 1, 1);
 
494
    dim3 output_grid_dim(fft_size, 1, 1);
 
495
 
 
496
 
 
497
    for (int i = 0;i < ncp;i++) {
 
498
        float *cudaDenom = ps->cuda_denom_buffer + (i+icp)*alloc_size;
 
499
        float *cudaLSum = ps->cuda_lsum_buffer + (i+icp)*alloc_size;
 
500
        cufftReal *cudaRealPtr = ps->cuda_result_buffer + (i+icp)*alloc_size;
 
501
        float *cudaResultPtr = ps->cuda_final_buffer + (i+icp)*alloc_size;
 
502
 
 
503
        vecCompute<<<output_grid_dim, output_block_dim>>>(
 
504
            cudaResultPtr,
 
505
            cudaRealPtr, 1./(fft_size2 * (size2 - 1)),
 
506
            cudaLSum, sumbuf+i, 1. / (size2 * (size2 - 1)),
 
507
            cudaDenom, stdbuf+i,
 
508
            fft_size
 
509
        );
 
510
    }
 
511
 
 
512
 
 
513
    return 0;
 
514
}
 
515
 
301
516
 
302
517
/*
303
518
static inline double fftDownloadData(TProcessingState *ps, mxArray *data) {
341
556
    const mxArray *nonzero;
342
557
#endif /* VALIDATE_LSUM */
343
558
 
 
559
    const mxArray *x, *y;
 
560
 
344
561
    double input_sum;
345
562
    double input_denom;
346
563
 
460
677
        }
461
678
*/
462
679
        
463
 
        icp = (unsigned int)mxGetScalar(prhs[2]);
 
680
        icp = (unsigned int)mxGetScalar(prhs[2]) - 1;
464
681
/*      if (icp >= ps->ncp) {
465
682
            reportError("The control point (%i) is out of range (0-%u)", icp, ps->ncp - 1);
466
683
            return;
467
684
        }
468
 
*/
 
685
 
469
686
        input = prhs[3];
470
 
/*    
 
687
    
471
688
        if (mxGetNumberOfDimensions(input) != 2) {
472
689
            reportError("Invalid dimensionality of input matrix, 2D matrix is expected");
473
690
            return;
478
695
            return;
479
696
        }
480
697
 
481
 
*/
482
698
        input_sum = mxGetScalar(prhs[4]);
483
699
        input_denom = mxGetScalar(prhs[5]);
484
 
 
485
 
        plhs[0] = fftCompute(ps, icp, input, input_sum, input_denom);
 
700
*/
 
701
 
 
702
        plhs[0] = fftCompute(ps, icp);
 
703
     break;
 
704
     case ACTION_COMPUTE:
 
705
        if (nrhs != 3) {
 
706
            reportError("This action expects 1 argument, but %i is passed", nrhs - 2);
 
707
            return;
 
708
        }
 
709
 
 
710
        input = prhs[2];
 
711
        
 
712
        if (mxGetClassID(input) != mxUINT8_CLASS) {
 
713
            reportError("Invalid type of image data, should be 8bit integers");
 
714
            return;
 
715
        }
 
716
        
 
717
        for (icp = 0; icp < ps->ncp; icp+=CP_BLOCK) {
 
718
            err = fftLoadFragment(ps, icp, min2(CP_BLOCK, ps->ncp - icp), input);
 
719
        }
 
720
        
486
721
     break;
487
722
     case ACTION_COMPUTE_BASE:
488
723
        if ((nrhs != 4)
494
729
            return;
495
730
        }
496
731
 
497
 
        icp = (unsigned int)mxGetScalar(prhs[2]);
 
732
        icp = (unsigned int)mxGetScalar(prhs[2]) - 1;
498
733
        if (icp >= ps->ncp) {
499
734
            reportError("The control point (%i) is out of range (0-%u)", icp, ps->ncp - 1);
500
735
            return;
546
781
#endif /* VALIDATE_LSUM */
547
782
 
548
783
     break;
 
784
     case ACTION_SET_POINTS:
 
785
        if (nrhs != 4) {
 
786
            reportError("SET_POINTS action expects two arrays with 'x' and 'y' coordinates of control points");
 
787
            return;
 
788
        }
 
789
 
 
790
        x = prhs[2];
 
791
        y = prhs[3];
 
792
        
 
793
        if (    (mxGetClassID(x) != mxSINGLE_CLASS)||
 
794
                (mxGetClassID(y) != mxSINGLE_CLASS)||
 
795
                (mxGetN(x)*mxGetM(x) != ps->ncp)||
 
796
                (mxGetN(y)*mxGetM(y) != ps->ncp)
 
797
        ) {
 
798
            reportError("Invalid control points are specified");
 
799
            return;
 
800
        }
 
801
        
 
802
        memcpy(ps->data_x, mxGetData(x), ps->ncp * sizeof(float));
 
803
        memcpy(ps->data_y, mxGetData(y), ps->ncp * sizeof(float));
 
804
     break;
549
805
     case ACTION_SETUP:
550
806
        if (nrhs != 4) {
551
807
            reportError("SETUP action expects 'ncp' and 'corrsize' parameters");
552
808
            return;
553
809
        }
554
810
 
555
 
        iprop = (int)mxGetScalar(prhs[2]);
556
 
        ps->ncp = iprop + 1;
 
811
        ps->ncp = (int)mxGetScalar(prhs[2]);
557
812
 
558
813
        iprop = (int)mxGetScalar(prhs[3]);
559
814
        ps->corr_size = iprop;
561
816
        ps->fft_size2 = ps->fft_size * ps->fft_size;
562
817
        ps->fft_inner_size = ps->fft_size * (ps->fft_size / 2 + 1);
563
818
 
 
819
        if (ps->fft_size % SIDE_BLOCK_SIZE) ps->side_alloc_size = (ps->fft_size / SIDE_BLOCK_SIZE + 1) * SIDE_BLOCK_SIZE;
 
820
        else ps->side_alloc_size = ps->fft_size;
 
821
 
 
822
        
564
823
        if (ps->fft_size2 % BLOCK_SIZE_1D) {
565
824
            ps->fft_alloc_size = ((ps->fft_size2 / BLOCK_SIZE_1D) + 1) * BLOCK_SIZE_1D;
566
825
        } else {