339 #define fft2d_alloc_error_check(p) { \
341 fprintf(stderr, "fft2d memory allocation error\n"); \
347 #ifdef USE_FFT2D_PTHREADS
348 #define USE_FFT2D_THREADS
349 #ifndef FFT2D_MAX_THREADS
350 #define FFT2D_MAX_THREADS 4
352 #ifndef FFT2D_THREADS_BEGIN_N
353 #define FFT2D_THREADS_BEGIN_N 65536
356 #define fft2d_thread_t pthread_t
357 #define fft2d_thread_create(thp,func,argp) { \
358 if (pthread_create(thp, NULL, func, (void *) (argp)) != 0) { \
359 fprintf(stderr, "fft2d thread error\n"); \
363 #define fft2d_thread_wait(th) { \
364 if (pthread_join(th, NULL) != 0) { \
365 fprintf(stderr, "fft2d thread error\n"); \
372 #ifdef USE_FFT2D_WINTHREADS
373 #define USE_FFT2D_THREADS
374 #ifndef FFT2D_MAX_THREADS
375 #define FFT2D_MAX_THREADS 4
377 #ifndef FFT2D_THREADS_BEGIN_N
378 #define FFT2D_THREADS_BEGIN_N 131072
381 #define fft2d_thread_t HANDLE
382 #define fft2d_thread_create(thp,func,argp) { \
384 *(thp) = CreateThread(NULL, 0, (LPTHREAD_START_ROUTINE) (func), (LPVOID) (argp), 0, &thid); \
386 fprintf(stderr, "fft2d thread error\n"); \
390 #define fft2d_thread_wait(th) { \
391 WaitForSingleObject(th, INFINITE); \
397 void cdft2d(
int n1,
int n2,
int isgn,
double **a,
double *t,
400 void makewt(
int nw,
int *ip,
double *w);
401 void cdft(
int n,
int isgn,
double *a,
int *ip,
double *w);
402 void cdft2d_sub(
int n1,
int n2,
int isgn,
double **a,
double *t,
404 #ifdef USE_FFT2D_THREADS
405 void xdft2d0_subth(
int n1,
int n2,
int icr,
int isgn,
double **a,
407 void cdft2d_subth(
int n1,
int n2,
int isgn,
double **a,
double *t,
410 int n, itnull, nthread, nt, i;
416 if (n > (ip[0] << 2)) {
417 makewt(n >> 2, ip, w);
423 #ifdef USE_FFT2D_THREADS
424 nthread = FFT2D_MAX_THREADS;
426 nt = 8 * nthread * n1;
427 if (n2 == 4 * nthread) {
429 }
else if (n2 < 4 * nthread) {
432 t = (
double *) malloc(
sizeof(
double) * nt);
433 fft2d_alloc_error_check(t);
435 #ifdef USE_FFT2D_THREADS
436 if ((
double) n1 * n2 >= (
double) FFT2D_THREADS_BEGIN_N) {
437 xdft2d0_subth(n1, n2, 0, isgn, a, ip, w);
438 cdft2d_subth(n1, n2, isgn, a, t, ip, w);
442 for (i = 0; i < n1; i++) {
443 cdft(n2, isgn, a[i], ip, w);
445 cdft2d_sub(n1, n2, isgn, a, t, ip, w);
453 void rdft2d(
int n1,
int n2,
int isgn,
double **a,
double *t,
456 void makewt(
int nw,
int *ip,
double *w);
457 void makect(
int nc,
int *ip,
double *c);
458 void rdft(
int n,
int isgn,
double *a,
int *ip,
double *w);
459 void cdft2d_sub(
int n1,
int n2,
int isgn,
double **a,
double *t,
461 void rdft2d_sub(
int n1,
int n2,
int isgn,
double **a);
462 #ifdef USE_FFT2D_THREADS
463 void xdft2d0_subth(
int n1,
int n2,
int icr,
int isgn,
double **a,
465 void cdft2d_subth(
int n1,
int n2,
int isgn,
double **a,
double *t,
468 int n, nw, nc, itnull, nthread, nt, i;
480 if (n2 > (nc << 2)) {
482 makect(nc, ip, w + nw);
488 #ifdef USE_FFT2D_THREADS
489 nthread = FFT2D_MAX_THREADS;
491 nt = 8 * nthread * n1;
492 if (n2 == 4 * nthread) {
494 }
else if (n2 < 4 * nthread) {
497 t = (
double *) malloc(
sizeof(
double) * nt);
498 fft2d_alloc_error_check(t);
500 #ifdef USE_FFT2D_THREADS
501 if ((
double) n1 * n2 >= (
double) FFT2D_THREADS_BEGIN_N) {
503 rdft2d_sub(n1, n2, isgn, a);
504 cdft2d_subth(n1, n2, isgn, a, t, ip, w);
506 xdft2d0_subth(n1, n2, 1, isgn, a, ip, w);
508 cdft2d_subth(n1, n2, isgn, a, t, ip, w);
509 rdft2d_sub(n1, n2, isgn, a);
515 rdft2d_sub(n1, n2, isgn, a);
516 cdft2d_sub(n1, n2, isgn, a, t, ip, w);
518 for (i = 0; i < n1; i++) {
519 rdft(n2, isgn, a[i], ip, w);
522 cdft2d_sub(n1, n2, isgn, a, t, ip, w);
523 rdft2d_sub(n1, n2, isgn, a);
532 void rdft2dsort(
int n1,
int n2,
int isgn,
double **a)
539 for (i = n1h + 1; i < n1; i++) {
540 a[i][0] = a[i][n2 + 1];
544 a[n1h][1] = a[n1h][n2];
546 for (i = n1h + 1; i < n1; i++) {
552 a[n1 - i][n2 + 1] = -y;
553 a[i][0] = a[n1 - i][0];
554 a[i][1] = -a[n1 - i][1];
559 a[n1h][n2] = a[n1h][1];
566 void ddct2d(
int n1,
int n2,
int isgn,
double **a,
double *t,
569 void makewt(
int nw,
int *ip,
double *w);
570 void makect(
int nc,
int *ip,
double *c);
571 void ddct(
int n,
int isgn,
double *a,
int *ip,
double *w);
572 void ddxt2d_sub(
int n1,
int n2,
int ics,
int isgn,
double **a,
573 double *t,
int *ip,
double *w);
574 #ifdef USE_FFT2D_THREADS
575 void ddxt2d0_subth(
int n1,
int n2,
int ics,
int isgn,
double **a,
577 void ddxt2d_subth(
int n1,
int n2,
int ics,
int isgn,
double **a,
578 double *t,
int *ip,
double *w);
580 int n, nw, nc, itnull, nthread, nt, i;
594 makect(nc, ip, w + nw);
600 #ifdef USE_FFT2D_THREADS
601 nthread = FFT2D_MAX_THREADS;
603 nt = 4 * nthread * n1;
604 if (n2 == 2 * nthread) {
606 }
else if (n2 < 2 * nthread) {
609 t = (
double *) malloc(
sizeof(
double) * nt);
610 fft2d_alloc_error_check(t);
612 #ifdef USE_FFT2D_THREADS
613 if ((
double) n1 * n2 >= (
double) FFT2D_THREADS_BEGIN_N) {
614 ddxt2d0_subth(n1, n2, 0, isgn, a, ip, w);
615 ddxt2d_subth(n1, n2, 0, isgn, a, t, ip, w);
619 for (i = 0; i < n1; i++) {
620 ddct(n2, isgn, a[i], ip, w);
622 ddxt2d_sub(n1, n2, 0, isgn, a, t, ip, w);
630 void ddst2d(
int n1,
int n2,
int isgn,
double **a,
double *t,
633 void makewt(
int nw,
int *ip,
double *w);
634 void makect(
int nc,
int *ip,
double *c);
635 void ddst(
int n,
int isgn,
double *a,
int *ip,
double *w);
636 void ddxt2d_sub(
int n1,
int n2,
int ics,
int isgn,
double **a,
637 double *t,
int *ip,
double *w);
638 #ifdef USE_FFT2D_THREADS
639 void ddxt2d0_subth(
int n1,
int n2,
int ics,
int isgn,
double **a,
641 void ddxt2d_subth(
int n1,
int n2,
int ics,
int isgn,
double **a,
642 double *t,
int *ip,
double *w);
644 int n, nw, nc, itnull, nthread, nt, i;
658 makect(nc, ip, w + nw);
664 #ifdef USE_FFT2D_THREADS
665 nthread = FFT2D_MAX_THREADS;
667 nt = 4 * nthread * n1;
668 if (n2 == 2 * nthread) {
670 }
else if (n2 < 2 * nthread) {
673 t = (
double *) malloc(
sizeof(
double) * nt);
674 fft2d_alloc_error_check(t);
676 #ifdef USE_FFT2D_THREADS
677 if ((
double) n1 * n2 >= (
double) FFT2D_THREADS_BEGIN_N) {
678 ddxt2d0_subth(n1, n2, 1, isgn, a, ip, w);
679 ddxt2d_subth(n1, n2, 1, isgn, a, t, ip, w);
683 for (i = 0; i < n1; i++) {
684 ddst(n2, isgn, a[i], ip, w);
686 ddxt2d_sub(n1, n2, 1, isgn, a, t, ip, w);
697 void cdft2d_sub(
int n1,
int n2,
int isgn,
double **a,
double *t,
700 void cdft(
int n,
int isgn,
double *a,
int *ip,
double *w);
704 for (j = 0; j < n2; j += 8) {
705 for (i = 0; i < n1; i++) {
707 t[2 * i + 1] = a[i][j + 1];
708 t[2 * n1 + 2 * i] = a[i][j + 2];
709 t[2 * n1 + 2 * i + 1] = a[i][j + 3];
710 t[4 * n1 + 2 * i] = a[i][j + 4];
711 t[4 * n1 + 2 * i + 1] = a[i][j + 5];
712 t[6 * n1 + 2 * i] = a[i][j + 6];
713 t[6 * n1 + 2 * i + 1] = a[i][j + 7];
715 cdft(2 * n1, isgn, t, ip, w);
716 cdft(2 * n1, isgn, &t[2 * n1], ip, w);
717 cdft(2 * n1, isgn, &t[4 * n1], ip, w);
718 cdft(2 * n1, isgn, &t[6 * n1], ip, w);
719 for (i = 0; i < n1; i++) {
721 a[i][j + 1] = t[2 * i + 1];
722 a[i][j + 2] = t[2 * n1 + 2 * i];
723 a[i][j + 3] = t[2 * n1 + 2 * i + 1];
724 a[i][j + 4] = t[4 * n1 + 2 * i];
725 a[i][j + 5] = t[4 * n1 + 2 * i + 1];
726 a[i][j + 6] = t[6 * n1 + 2 * i];
727 a[i][j + 7] = t[6 * n1 + 2 * i + 1];
730 }
else if (n2 == 4) {
731 for (i = 0; i < n1; i++) {
733 t[2 * i + 1] = a[i][1];
734 t[2 * n1 + 2 * i] = a[i][2];
735 t[2 * n1 + 2 * i + 1] = a[i][3];
737 cdft(2 * n1, isgn, t, ip, w);
738 cdft(2 * n1, isgn, &t[2 * n1], ip, w);
739 for (i = 0; i < n1; i++) {
741 a[i][1] = t[2 * i + 1];
742 a[i][2] = t[2 * n1 + 2 * i];
743 a[i][3] = t[2 * n1 + 2 * i + 1];
745 }
else if (n2 == 2) {
746 for (i = 0; i < n1; i++) {
748 t[2 * i + 1] = a[i][1];
750 cdft(2 * n1, isgn, t, ip, w);
751 for (i = 0; i < n1; i++) {
753 a[i][1] = t[2 * i + 1];
759 void rdft2d_sub(
int n1,
int n2,
int isgn,
double **a)
766 for (i = 1; i < n1h; i++) {
768 xi = a[i][0] - a[j][0];
771 xi = a[j][1] - a[i][1];
776 for (i = 1; i < n1h; i++) {
778 a[j][0] = 0.5 * (a[i][0] - a[j][0]);
780 a[j][1] = 0.5 * (a[i][1] + a[j][1]);
787 void ddxt2d_sub(
int n1,
int n2,
int ics,
int isgn,
double **a,
788 double *t,
int *ip,
double *w)
790 void ddct(
int n,
int isgn,
double *a,
int *ip,
double *w);
791 void ddst(
int n,
int isgn,
double *a,
int *ip,
double *w);
795 for (j = 0; j < n2; j += 4) {
796 for (i = 0; i < n1; i++) {
798 t[n1 + i] = a[i][j + 1];
799 t[2 * n1 + i] = a[i][j + 2];
800 t[3 * n1 + i] = a[i][j + 3];
803 ddct(n1, isgn, t, ip, w);
804 ddct(n1, isgn, &t[n1], ip, w);
805 ddct(n1, isgn, &t[2 * n1], ip, w);
806 ddct(n1, isgn, &t[3 * n1], ip, w);
808 ddst(n1, isgn, t, ip, w);
809 ddst(n1, isgn, &t[n1], ip, w);
810 ddst(n1, isgn, &t[2 * n1], ip, w);
811 ddst(n1, isgn, &t[3 * n1], ip, w);
813 for (i = 0; i < n1; i++) {
815 a[i][j + 1] = t[n1 + i];
816 a[i][j + 2] = t[2 * n1 + i];
817 a[i][j + 3] = t[3 * n1 + i];
820 }
else if (n2 == 2) {
821 for (i = 0; i < n1; i++) {
826 ddct(n1, isgn, t, ip, w);
827 ddct(n1, isgn, &t[n1], ip, w);
829 ddst(n1, isgn, t, ip, w);
830 ddst(n1, isgn, &t[n1], ip, w);
832 for (i = 0; i < n1; i++) {
840 #ifdef USE_FFT2D_THREADS
841 struct fft2d_arg_st {
853 typedef struct fft2d_arg_st fft2d_arg_t;
856 void xdft2d0_subth(
int n1,
int n2,
int icr,
int isgn,
double **a,
859 void *xdft2d0_th(
void *p);
860 fft2d_thread_t th[FFT2D_MAX_THREADS];
861 fft2d_arg_t ag[FFT2D_MAX_THREADS];
864 nthread = FFT2D_MAX_THREADS;
868 for (i = 0; i < nthread; i++) {
869 ag[i].nthread = nthread;
878 fft2d_thread_create(&th[i], xdft2d0_th, &ag[i]);
880 for (i = 0; i < nthread; i++) {
881 fft2d_thread_wait(th[i]);
886 void cdft2d_subth(
int n1,
int n2,
int isgn,
double **a,
double *t,
889 void *cdft2d_th(
void *p);
890 fft2d_thread_t th[FFT2D_MAX_THREADS];
891 fft2d_arg_t ag[FFT2D_MAX_THREADS];
894 nthread = FFT2D_MAX_THREADS;
896 if (n2 == 4 * FFT2D_MAX_THREADS) {
898 }
else if (n2 < 4 * FFT2D_MAX_THREADS) {
902 for (i = 0; i < nthread; i++) {
903 ag[i].nthread = nthread;
909 ag[i].t = &t[nt * i];
912 fft2d_thread_create(&th[i], cdft2d_th, &ag[i]);
914 for (i = 0; i < nthread; i++) {
915 fft2d_thread_wait(th[i]);
920 void ddxt2d0_subth(
int n1,
int n2,
int ics,
int isgn,
double **a,
923 void *ddxt2d0_th(
void *p);
924 fft2d_thread_t th[FFT2D_MAX_THREADS];
925 fft2d_arg_t ag[FFT2D_MAX_THREADS];
928 nthread = FFT2D_MAX_THREADS;
932 for (i = 0; i < nthread; i++) {
933 ag[i].nthread = nthread;
942 fft2d_thread_create(&th[i], ddxt2d0_th, &ag[i]);
944 for (i = 0; i < nthread; i++) {
945 fft2d_thread_wait(th[i]);
950 void ddxt2d_subth(
int n1,
int n2,
int ics,
int isgn,
double **a,
951 double *t,
int *ip,
double *w)
953 void *ddxt2d_th(
void *p);
954 fft2d_thread_t th[FFT2D_MAX_THREADS];
955 fft2d_arg_t ag[FFT2D_MAX_THREADS];
958 nthread = FFT2D_MAX_THREADS;
960 if (n2 == 2 * FFT2D_MAX_THREADS) {
962 }
else if (n2 < 2 * FFT2D_MAX_THREADS) {
966 for (i = 0; i < nthread; i++) {
967 ag[i].nthread = nthread;
974 ag[i].t = &t[nt * i];
977 fft2d_thread_create(&th[i], ddxt2d_th, &ag[i]);
979 for (i = 0; i < nthread; i++) {
980 fft2d_thread_wait(th[i]);
985 void *xdft2d0_th(
void *p)
987 void cdft(
int n,
int isgn,
double *a,
int *ip,
double *w);
988 void rdft(
int n,
int isgn,
double *a,
int *ip,
double *w);
989 int nthread, n0, n1, n2, icr, isgn, *ip, i;
992 nthread = ((fft2d_arg_t *) p)->nthread;
993 n0 = ((fft2d_arg_t *) p)->n0;
994 n1 = ((fft2d_arg_t *) p)->n1;
995 n2 = ((fft2d_arg_t *) p)->n2;
996 icr = ((fft2d_arg_t *) p)->ic;
997 isgn = ((fft2d_arg_t *) p)->isgn;
998 a = ((fft2d_arg_t *) p)->a;
999 ip = ((fft2d_arg_t *) p)->ip;
1000 w = ((fft2d_arg_t *) p)->w;
1002 for (i = n0; i < n1; i += nthread) {
1003 cdft(n2, isgn, a[i], ip, w);
1006 for (i = n0; i < n1; i += nthread) {
1007 rdft(n2, isgn, a[i], ip, w);
1014 void *cdft2d_th(
void *p)
1016 void cdft(
int n,
int isgn,
double *a,
int *ip,
double *w);
1017 int nthread, n0, n1, n2, isgn, *ip, i, j;
1020 nthread = ((fft2d_arg_t *) p)->nthread;
1021 n0 = ((fft2d_arg_t *) p)->n0;
1022 n1 = ((fft2d_arg_t *) p)->n1;
1023 n2 = ((fft2d_arg_t *) p)->n2;
1024 isgn = ((fft2d_arg_t *) p)->isgn;
1025 a = ((fft2d_arg_t *) p)->a;
1026 t = ((fft2d_arg_t *) p)->t;
1027 ip = ((fft2d_arg_t *) p)->ip;
1028 w = ((fft2d_arg_t *) p)->w;
1029 if (n2 > 4 * nthread) {
1030 for (j = 8 * n0; j < n2; j += 8 * nthread) {
1031 for (i = 0; i < n1; i++) {
1033 t[2 * i + 1] = a[i][j + 1];
1034 t[2 * n1 + 2 * i] = a[i][j + 2];
1035 t[2 * n1 + 2 * i + 1] = a[i][j + 3];
1036 t[4 * n1 + 2 * i] = a[i][j + 4];
1037 t[4 * n1 + 2 * i + 1] = a[i][j + 5];
1038 t[6 * n1 + 2 * i] = a[i][j + 6];
1039 t[6 * n1 + 2 * i + 1] = a[i][j + 7];
1041 cdft(2 * n1, isgn, t, ip, w);
1042 cdft(2 * n1, isgn, &t[2 * n1], ip, w);
1043 cdft(2 * n1, isgn, &t[4 * n1], ip, w);
1044 cdft(2 * n1, isgn, &t[6 * n1], ip, w);
1045 for (i = 0; i < n1; i++) {
1047 a[i][j + 1] = t[2 * i + 1];
1048 a[i][j + 2] = t[2 * n1 + 2 * i];
1049 a[i][j + 3] = t[2 * n1 + 2 * i + 1];
1050 a[i][j + 4] = t[4 * n1 + 2 * i];
1051 a[i][j + 5] = t[4 * n1 + 2 * i + 1];
1052 a[i][j + 6] = t[6 * n1 + 2 * i];
1053 a[i][j + 7] = t[6 * n1 + 2 * i + 1];
1056 }
else if (n2 == 4 * nthread) {
1057 for (i = 0; i < n1; i++) {
1058 t[2 * i] = a[i][4 * n0];
1059 t[2 * i + 1] = a[i][4 * n0 + 1];
1060 t[2 * n1 + 2 * i] = a[i][4 * n0 + 2];
1061 t[2 * n1 + 2 * i + 1] = a[i][4 * n0 + 3];
1063 cdft(2 * n1, isgn, t, ip, w);
1064 cdft(2 * n1, isgn, &t[2 * n1], ip, w);
1065 for (i = 0; i < n1; i++) {
1066 a[i][4 * n0] = t[2 * i];
1067 a[i][4 * n0 + 1] = t[2 * i + 1];
1068 a[i][4 * n0 + 2] = t[2 * n1 + 2 * i];
1069 a[i][4 * n0 + 3] = t[2 * n1 + 2 * i + 1];
1071 }
else if (n2 == 2 * nthread) {
1072 for (i = 0; i < n1; i++) {
1073 t[2 * i] = a[i][2 * n0];
1074 t[2 * i + 1] = a[i][2 * n0 + 1];
1076 cdft(2 * n1, isgn, t, ip, w);
1077 for (i = 0; i < n1; i++) {
1078 a[i][2 * n0] = t[2 * i];
1079 a[i][2 * n0 + 1] = t[2 * i + 1];
1086 void *ddxt2d0_th(
void *p)
1088 void ddct(
int n,
int isgn,
double *a,
int *ip,
double *w);
1089 void ddst(
int n,
int isgn,
double *a,
int *ip,
double *w);
1090 int nthread, n0, n1, n2, ics, isgn, *ip, i;
1093 nthread = ((fft2d_arg_t *) p)->nthread;
1094 n0 = ((fft2d_arg_t *) p)->n0;
1095 n1 = ((fft2d_arg_t *) p)->n1;
1096 n2 = ((fft2d_arg_t *) p)->n2;
1097 ics = ((fft2d_arg_t *) p)->ic;
1098 isgn = ((fft2d_arg_t *) p)->isgn;
1099 a = ((fft2d_arg_t *) p)->a;
1100 ip = ((fft2d_arg_t *) p)->ip;
1101 w = ((fft2d_arg_t *) p)->w;
1103 for (i = n0; i < n1; i += nthread) {
1104 ddct(n2, isgn, a[i], ip, w);
1107 for (i = n0; i < n1; i += nthread) {
1108 ddst(n2, isgn, a[i], ip, w);
1115 void *ddxt2d_th(
void *p)
1117 void ddct(
int n,
int isgn,
double *a,
int *ip,
double *w);
1118 void ddst(
int n,
int isgn,
double *a,
int *ip,
double *w);
1119 int nthread, n0, n1, n2, ics, isgn, *ip, i, j;
1122 nthread = ((fft2d_arg_t *) p)->nthread;
1123 n0 = ((fft2d_arg_t *) p)->n0;
1124 n1 = ((fft2d_arg_t *) p)->n1;
1125 n2 = ((fft2d_arg_t *) p)->n2;
1126 ics = ((fft2d_arg_t *) p)->ic;
1127 isgn = ((fft2d_arg_t *) p)->isgn;
1128 a = ((fft2d_arg_t *) p)->a;
1129 t = ((fft2d_arg_t *) p)->t;
1130 ip = ((fft2d_arg_t *) p)->ip;
1131 w = ((fft2d_arg_t *) p)->w;
1132 if (n2 > 2 * nthread) {
1133 for (j = 4 * n0; j < n2; j += 4 * nthread) {
1134 for (i = 0; i < n1; i++) {
1136 t[n1 + i] = a[i][j + 1];
1137 t[2 * n1 + i] = a[i][j + 2];
1138 t[3 * n1 + i] = a[i][j + 3];
1141 ddct(n1, isgn, t, ip, w);
1142 ddct(n1, isgn, &t[n1], ip, w);
1143 ddct(n1, isgn, &t[2 * n1], ip, w);
1144 ddct(n1, isgn, &t[3 * n1], ip, w);
1146 ddst(n1, isgn, t, ip, w);
1147 ddst(n1, isgn, &t[n1], ip, w);
1148 ddst(n1, isgn, &t[2 * n1], ip, w);
1149 ddst(n1, isgn, &t[3 * n1], ip, w);
1151 for (i = 0; i < n1; i++) {
1153 a[i][j + 1] = t[n1 + i];
1154 a[i][j + 2] = t[2 * n1 + i];
1155 a[i][j + 3] = t[3 * n1 + i];
1158 }
else if (n2 == 2 * nthread) {
1159 for (i = 0; i < n1; i++) {
1160 t[i] = a[i][2 * n0];
1161 t[n1 + i] = a[i][2 * n0 + 1];
1164 ddct(n1, isgn, t, ip, w);
1165 ddct(n1, isgn, &t[n1], ip, w);
1167 ddst(n1, isgn, t, ip, w);
1168 ddst(n1, isgn, &t[n1], ip, w);
1170 for (i = 0; i < n1; i++) {
1171 a[i][2 * n0] = t[i];
1172 a[i][2 * n0 + 1] = t[n1 + i];
1174 }
else if (n2 == nthread) {
1175 for (i = 0; i < n1; i++) {
1179 ddct(n1, isgn, t, ip, w);
1181 ddst(n1, isgn, t, ip, w);
1183 for (i = 0; i < n1; i++) {