hayati ayguen | b2d2936 | 2020-01-02 00:06:09 +0100 | [diff] [blame] | 1 | /* |
| 2 | Copyright (c) 2013 Julien Pommier. |
| 3 | |
| 4 | Small test & bench for PFFFT, comparing its performance with the scalar FFTPACK, FFTW, and Apple vDSP |
| 5 | |
| 6 | How to build: |
| 7 | |
| 8 | on linux, with fftw3: |
| 9 | gcc -o test_pffft -DHAVE_FFTW -msse -mfpmath=sse -O3 -Wall -W pffft.c test_pffft.c fftpack.c -L/usr/local/lib -I/usr/local/include/ -lfftw3f -lm |
| 10 | |
| 11 | on macos, without fftw3: |
| 12 | clang -o test_pffft -DHAVE_VECLIB -O3 -Wall -W pffft.c test_pffft.c fftpack.c -L/usr/local/lib -I/usr/local/include/ -framework Accelerate |
| 13 | |
| 14 | on macos, with fftw3: |
| 15 | clang -o test_pffft -DHAVE_FFTW -DHAVE_VECLIB -O3 -Wall -W pffft.c test_pffft.c fftpack.c -L/usr/local/lib -I/usr/local/include/ -lfftw3f -framework Accelerate |
| 16 | |
| 17 | as alternative: replace clang by gcc. |
| 18 | |
| 19 | on windows, with visual c++: |
| 20 | cl /Ox -D_USE_MATH_DEFINES /arch:SSE test_pffft.c pffft.c fftpack.c |
| 21 | |
| 22 | build without SIMD instructions: |
| 23 | gcc -o test_pffft -DPFFFT_SIMD_DISABLE -O3 -Wall -W pffft.c test_pffft.c fftpack.c -lm |
| 24 | |
| 25 | */ |
| 26 | |
hayati ayguen | ca11241 | 2020-04-13 00:19:40 +0200 | [diff] [blame] | 27 | #ifdef PFFFT_ENABLE_FLOAT |
hayati ayguen | b2d2936 | 2020-01-02 00:06:09 +0100 | [diff] [blame] | 28 | #include "pffft.h" |
| 29 | |
hayati ayguen | ca11241 | 2020-04-13 00:19:40 +0200 | [diff] [blame] | 30 | typedef float pffft_scalar; |
| 31 | #else |
| 32 | /* |
| 33 | Note: adapted for double precision dynamic range version. |
| 34 | */ |
| 35 | #include "pffft_double.h" |
| 36 | |
| 37 | typedef double pffft_scalar; |
| 38 | #endif |
| 39 | |
hayati ayguen | b2d2936 | 2020-01-02 00:06:09 +0100 | [diff] [blame] | 40 | #include <math.h> |
| 41 | #include <stdio.h> |
| 42 | #include <stdlib.h> |
| 43 | #include <time.h> |
| 44 | #include <assert.h> |
| 45 | #include <string.h> |
| 46 | |
hayati ayguen | c974c1d | 2020-03-29 03:39:30 +0200 | [diff] [blame] | 47 | /* define own constants required to turn off g++ extensions .. */ |
| 48 | #ifndef M_PI |
| 49 | #define M_PI 3.14159265358979323846 /* pi */ |
| 50 | #endif |
hayati ayguen | b2d2936 | 2020-01-02 00:06:09 +0100 | [diff] [blame] | 51 | |
| 52 | /* EXPECTED_DYN_RANGE in dB: |
| 53 | * single precision float has 24 bits mantissa |
| 54 | * => 24 Bits * 6 dB = 144 dB |
| 55 | * allow a few dB tolerance (even 144 dB looks good on my PC) |
| 56 | */ |
hayati ayguen | ca11241 | 2020-04-13 00:19:40 +0200 | [diff] [blame] | 57 | #ifdef PFFFT_ENABLE_FLOAT |
hayati ayguen | b2d2936 | 2020-01-02 00:06:09 +0100 | [diff] [blame] | 58 | #define EXPECTED_DYN_RANGE 140.0 |
hayati ayguen | ca11241 | 2020-04-13 00:19:40 +0200 | [diff] [blame] | 59 | #else |
| 60 | #define EXPECTED_DYN_RANGE 215.0 |
| 61 | #endif |
hayati ayguen | b2d2936 | 2020-01-02 00:06:09 +0100 | [diff] [blame] | 62 | |
| 63 | /* maximum allowed phase error in degree */ |
| 64 | #define DEG_ERR_LIMIT 1E-4 |
| 65 | |
| 66 | /* maximum allowed magnitude error in amplitude (of 1.0 or 1.1) */ |
| 67 | #define MAG_ERR_LIMIT 1E-6 |
| 68 | |
| 69 | |
| 70 | #define PRINT_SPEC 0 |
| 71 | |
| 72 | #define PWR2LOG(PWR) ( (PWR) < 1E-30 ? 10.0*log10(1E-30) : 10.0*log10(PWR) ) |
| 73 | |
| 74 | |
hayati ayguen | b2d2936 | 2020-01-02 00:06:09 +0100 | [diff] [blame] | 75 | |
| 76 | int test(int N, int cplx, int useOrdered) { |
| 77 | int Nfloat = (cplx ? N*2 : N); |
hayati ayguen | ca11241 | 2020-04-13 00:19:40 +0200 | [diff] [blame] | 78 | #ifdef PFFFT_ENABLE_FLOAT |
| 79 | pffft_scalar *X = pffft_aligned_malloc((unsigned)Nfloat * sizeof(pffft_scalar)); |
| 80 | pffft_scalar *Y = pffft_aligned_malloc((unsigned)Nfloat * sizeof(pffft_scalar)); |
| 81 | pffft_scalar *R = pffft_aligned_malloc((unsigned)Nfloat * sizeof(pffft_scalar)); |
| 82 | pffft_scalar *Z = pffft_aligned_malloc((unsigned)Nfloat * sizeof(pffft_scalar)); |
| 83 | pffft_scalar *W = pffft_aligned_malloc((unsigned)Nfloat * sizeof(pffft_scalar)); |
| 84 | #else |
| 85 | pffft_scalar *X = pffftd_aligned_malloc((unsigned)Nfloat * sizeof(pffft_scalar)); |
| 86 | pffft_scalar *Y = pffftd_aligned_malloc((unsigned)Nfloat * sizeof(pffft_scalar)); |
| 87 | pffft_scalar *R = pffftd_aligned_malloc((unsigned)Nfloat * sizeof(pffft_scalar)); |
| 88 | pffft_scalar *Z = pffftd_aligned_malloc((unsigned)Nfloat * sizeof(pffft_scalar)); |
| 89 | pffft_scalar *W = pffftd_aligned_malloc((unsigned)Nfloat * sizeof(pffft_scalar)); |
| 90 | #endif |
| 91 | pffft_scalar amp = (pffft_scalar)1.0; |
hayati ayguen | b2d2936 | 2020-01-02 00:06:09 +0100 | [diff] [blame] | 92 | double freq, dPhi, phi, phi0; |
| 93 | double pwr, pwrCar, pwrOther, err, errSum, mag, expextedMag; |
hayati ayguen | ca11241 | 2020-04-13 00:19:40 +0200 | [diff] [blame] | 94 | int k, j, m, iter, kmaxOther, retError = 0; |
hayati ayguen | b2d2936 | 2020-01-02 00:06:09 +0100 | [diff] [blame] | 95 | |
hayati ayguen | ca11241 | 2020-04-13 00:19:40 +0200 | [diff] [blame] | 96 | #ifdef PFFFT_ENABLE_FLOAT |
hayati ayguen | e6cffc9 | 2020-02-28 19:57:30 +0100 | [diff] [blame] | 97 | assert( pffft_is_power_of_two(N) ); |
hayati ayguen | b2d2936 | 2020-01-02 00:06:09 +0100 | [diff] [blame] | 98 | PFFFT_Setup *s = pffft_new_setup(N, cplx ? PFFFT_COMPLEX : PFFFT_REAL); |
hayati ayguen | ca11241 | 2020-04-13 00:19:40 +0200 | [diff] [blame] | 99 | #else |
| 100 | assert( pffftd_is_power_of_two(N) ); |
| 101 | PFFFTD_Setup *s = pffftd_new_setup(N, cplx ? PFFFT_COMPLEX : PFFFT_REAL); |
| 102 | #endif |
hayati ayguen | b2d2936 | 2020-01-02 00:06:09 +0100 | [diff] [blame] | 103 | assert(s); |
| 104 | if (!s) { |
| 105 | printf("Error setting up PFFFT!\n"); |
| 106 | return 1; |
| 107 | } |
| 108 | |
| 109 | for ( k = m = 0; k < (cplx? N : (1 + N/2) ); k += N/16, ++m ) |
| 110 | { |
hayati ayguen | ca11241 | 2020-04-13 00:19:40 +0200 | [diff] [blame] | 111 | amp = (pffft_scalar)( ( (m % 3) == 0 ) ? 1.0 : 1.1 ); |
hayati ayguen | b2d2936 | 2020-01-02 00:06:09 +0100 | [diff] [blame] | 112 | freq = (k < N/2) ? ((double)k / N) : ((double)(k-N) / N); |
| 113 | dPhi = 2.0 * M_PI * freq; |
| 114 | if ( dPhi < 0.0 ) |
| 115 | dPhi += 2.0 * M_PI; |
| 116 | |
| 117 | iter = -1; |
| 118 | while (1) |
| 119 | { |
| 120 | ++iter; |
| 121 | |
| 122 | if (iter) |
| 123 | printf("bin %d: dphi = %f for freq %f\n", k, dPhi, freq); |
| 124 | |
| 125 | /* generate cosine carrier as time signal - start at defined phase phi0 */ |
| 126 | phi = phi0 = (m % 4) * 0.125 * M_PI; /* have phi0 < 90 deg to be normalized */ |
| 127 | for ( j = 0; j < N; ++j ) |
| 128 | { |
| 129 | if (cplx) { |
hayati ayguen | ca11241 | 2020-04-13 00:19:40 +0200 | [diff] [blame] | 130 | X[2*j] = amp * (pffft_scalar)cos(phi); /* real part */ |
| 131 | X[2*j+1] = amp * (pffft_scalar)sin(phi); /* imag part */ |
hayati ayguen | b2d2936 | 2020-01-02 00:06:09 +0100 | [diff] [blame] | 132 | } |
| 133 | else |
hayati ayguen | ca11241 | 2020-04-13 00:19:40 +0200 | [diff] [blame] | 134 | X[j] = amp * (pffft_scalar)cos(phi); /* only real part */ |
hayati ayguen | b2d2936 | 2020-01-02 00:06:09 +0100 | [diff] [blame] | 135 | |
| 136 | /* phase increment .. stay normalized - cos()/sin() might degrade! */ |
| 137 | phi += dPhi; |
| 138 | if ( phi >= M_PI ) |
| 139 | phi -= 2.0 * M_PI; |
| 140 | } |
| 141 | |
| 142 | /* forward transform from X --> Y .. using work buffer W */ |
hayati ayguen | ca11241 | 2020-04-13 00:19:40 +0200 | [diff] [blame] | 143 | #ifdef PFFFT_ENABLE_FLOAT |
hayati ayguen | b2d2936 | 2020-01-02 00:06:09 +0100 | [diff] [blame] | 144 | if ( useOrdered ) |
| 145 | pffft_transform_ordered(s, X, Y, W, PFFFT_FORWARD ); |
| 146 | else |
| 147 | { |
hayati ayguen | 61ec6da | 2020-03-29 08:48:01 +0200 | [diff] [blame] | 148 | pffft_transform(s, X, R, W, PFFFT_FORWARD ); /* use R for reordering */ |
| 149 | pffft_zreorder(s, R, Y, PFFFT_FORWARD ); /* reorder into Y[] for power calculations */ |
hayati ayguen | b2d2936 | 2020-01-02 00:06:09 +0100 | [diff] [blame] | 150 | } |
hayati ayguen | ca11241 | 2020-04-13 00:19:40 +0200 | [diff] [blame] | 151 | #else |
| 152 | if ( useOrdered ) |
| 153 | pffftd_transform_ordered(s, X, Y, W, PFFFT_FORWARD ); |
| 154 | else |
| 155 | { |
| 156 | pffftd_transform(s, X, R, W, PFFFT_FORWARD ); /* use R for reordering */ |
| 157 | pffftd_zreorder(s, R, Y, PFFFT_FORWARD ); /* reorder into Y[] for power calculations */ |
| 158 | } |
| 159 | #endif |
hayati ayguen | b2d2936 | 2020-01-02 00:06:09 +0100 | [diff] [blame] | 160 | |
| 161 | pwrOther = -1.0; |
| 162 | pwrCar = 0; |
| 163 | |
| 164 | |
| 165 | /* for positive frequencies: 0 to 0.5 * samplerate */ |
| 166 | /* and also for negative frequencies: -0.5 * samplerate to 0 */ |
| 167 | for ( j = 0; j < ( cplx ? N : (1 + N/2) ); ++j ) |
| 168 | { |
| 169 | if (!cplx && !j) /* special treatment for DC for real input */ |
| 170 | pwr = Y[j]*Y[j]; |
| 171 | else if (!cplx && j == N/2) /* treat 0.5 * samplerate */ |
| 172 | pwr = Y[1] * Y[1]; /* despite j (for freq calculation) we have index 1 */ |
| 173 | else |
| 174 | pwr = Y[2*j] * Y[2*j] + Y[2*j+1] * Y[2*j+1]; |
| 175 | if (iter || PRINT_SPEC) |
| 176 | printf("%s fft %d: pwr[j = %d] = %g == %f dB\n", (cplx ? "cplx":"real"), N, j, pwr, PWR2LOG(pwr) ); |
| 177 | if (k == j) |
| 178 | pwrCar = pwr; |
| 179 | else if ( pwr > pwrOther ) { |
| 180 | pwrOther = pwr; |
| 181 | kmaxOther = j; |
| 182 | } |
| 183 | } |
| 184 | |
| 185 | if ( PWR2LOG(pwrCar) - PWR2LOG(pwrOther) < EXPECTED_DYN_RANGE ) { |
| 186 | printf("%s fft %d amp %f iter %d:\n", (cplx ? "cplx":"real"), N, amp, iter); |
| 187 | printf(" carrier power at bin %d: %g == %f dB\n", k, pwrCar, PWR2LOG(pwrCar) ); |
| 188 | printf(" carrier mag || at bin %d: %g\n", k, sqrt(pwrCar) ); |
| 189 | printf(" max other pwr at bin %d: %g == %f dB\n", kmaxOther, pwrOther, PWR2LOG(pwrOther) ); |
| 190 | printf(" dynamic range: %f dB\n\n", PWR2LOG(pwrCar) - PWR2LOG(pwrOther) ); |
| 191 | retError = 1; |
| 192 | if ( iter == 0 ) |
| 193 | continue; |
| 194 | } |
| 195 | |
| 196 | if ( k > 0 && k != N/2 ) |
| 197 | { |
| 198 | phi = atan2( Y[2*k+1], Y[2*k] ); |
| 199 | if ( fabs( phi - phi0) > DEG_ERR_LIMIT * M_PI / 180.0 ) |
| 200 | { |
| 201 | retError = 1; |
| 202 | printf("%s fft %d bin %d amp %f : phase mismatch! phase = %f deg expected = %f deg\n", |
| 203 | (cplx ? "cplx":"real"), N, k, amp, phi * 180.0 / M_PI, phi0 * 180.0 / M_PI ); |
| 204 | } |
| 205 | } |
| 206 | |
| 207 | expextedMag = cplx ? amp : ( (k == 0 || k == N/2) ? amp : (amp/2) ); |
| 208 | mag = sqrt(pwrCar) / N; |
| 209 | if ( fabs(mag - expextedMag) > MAG_ERR_LIMIT ) |
| 210 | { |
| 211 | retError = 1; |
| 212 | printf("%s fft %d bin %d amp %f : mag = %g expected = %g\n", (cplx ? "cplx":"real"), N, k, amp, mag, expextedMag ); |
| 213 | } |
| 214 | |
| 215 | |
| 216 | /* now convert spectrum back */ |
hayati ayguen | ca11241 | 2020-04-13 00:19:40 +0200 | [diff] [blame] | 217 | #ifdef PFFFT_ENABLE_FLOAT |
hayati ayguen | 61ec6da | 2020-03-29 08:48:01 +0200 | [diff] [blame] | 218 | if (useOrdered) |
| 219 | pffft_transform_ordered(s, Y, Z, W, PFFFT_BACKWARD); |
| 220 | else |
| 221 | pffft_transform(s, R, Z, W, PFFFT_BACKWARD); |
hayati ayguen | ca11241 | 2020-04-13 00:19:40 +0200 | [diff] [blame] | 222 | #else |
| 223 | if (useOrdered) |
| 224 | pffftd_transform_ordered(s, Y, Z, W, PFFFT_BACKWARD); |
| 225 | else |
| 226 | pffftd_transform(s, R, Z, W, PFFFT_BACKWARD); |
| 227 | #endif |
hayati ayguen | b2d2936 | 2020-01-02 00:06:09 +0100 | [diff] [blame] | 228 | |
| 229 | errSum = 0.0; |
| 230 | for ( j = 0; j < (cplx ? (2*N) : N); ++j ) |
| 231 | { |
| 232 | /* scale back */ |
| 233 | Z[j] /= N; |
| 234 | /* square sum errors over real (and imag parts) */ |
| 235 | err = (X[j]-Z[j]) * (X[j]-Z[j]); |
| 236 | errSum += err; |
| 237 | } |
| 238 | |
| 239 | if ( errSum > N * 1E-7 ) |
| 240 | { |
| 241 | retError = 1; |
| 242 | printf("%s fft %d bin %d : inverse FFT doesn't match original signal! errSum = %g ; mean err = %g\n", (cplx ? "cplx":"real"), N, k, errSum, errSum / N); |
| 243 | } |
| 244 | |
| 245 | break; |
| 246 | } |
| 247 | |
| 248 | } |
hayati ayguen | ca11241 | 2020-04-13 00:19:40 +0200 | [diff] [blame] | 249 | #ifdef PFFFT_ENABLE_FLOAT |
hayati ayguen | b2d2936 | 2020-01-02 00:06:09 +0100 | [diff] [blame] | 250 | pffft_destroy_setup(s); |
| 251 | pffft_aligned_free(X); |
| 252 | pffft_aligned_free(Y); |
| 253 | pffft_aligned_free(Z); |
hayati ayguen | 61ec6da | 2020-03-29 08:48:01 +0200 | [diff] [blame] | 254 | pffft_aligned_free(R); |
hayati ayguen | b2d2936 | 2020-01-02 00:06:09 +0100 | [diff] [blame] | 255 | pffft_aligned_free(W); |
hayati ayguen | ca11241 | 2020-04-13 00:19:40 +0200 | [diff] [blame] | 256 | #else |
| 257 | pffftd_destroy_setup(s); |
| 258 | pffftd_aligned_free(X); |
| 259 | pffftd_aligned_free(Y); |
| 260 | pffftd_aligned_free(Z); |
| 261 | pffftd_aligned_free(R); |
| 262 | pffftd_aligned_free(W); |
| 263 | #endif |
hayati ayguen | b2d2936 | 2020-01-02 00:06:09 +0100 | [diff] [blame] | 264 | |
| 265 | return retError; |
| 266 | } |
| 267 | |
hayati ayguen | ca11241 | 2020-04-13 00:19:40 +0200 | [diff] [blame] | 268 | /* small functions inside pffft.c that will detect (compiler) bugs with respect to simd instructions */ |
| 269 | void validate_pffft_simd(); |
| 270 | int validate_pffft_simd_ex(FILE * DbgOut); |
| 271 | void validate_pffftd_simd(); |
| 272 | int validate_pffftd_simd_ex(FILE * DbgOut); |
| 273 | |
hayati ayguen | b2d2936 | 2020-01-02 00:06:09 +0100 | [diff] [blame] | 274 | |
| 275 | |
| 276 | int main(int argc, char **argv) |
| 277 | { |
hayati ayguen | ca11241 | 2020-04-13 00:19:40 +0200 | [diff] [blame] | 278 | int N, result, resN, resAll, i, k, resNextPw2, resIsPw2, resFFT; |
hayati ayguen | b2d2936 | 2020-01-02 00:06:09 +0100 | [diff] [blame] | 279 | |
hayati ayguen | e6cffc9 | 2020-02-28 19:57:30 +0100 | [diff] [blame] | 280 | int inp_power_of_two[] = { 1, 2, 3, 4, 5, 6, 7, 8, 9, 511, 512, 513 }; |
| 281 | int ref_power_of_two[] = { 1, 2, 4, 4, 8, 8, 8, 8, 16, 512, 512, 1024 }; |
| 282 | |
hayati ayguen | ca11241 | 2020-04-13 00:19:40 +0200 | [diff] [blame] | 283 | for ( i = 1; i < argc; ++i ) { |
| 284 | |
| 285 | if (!strcmp(argv[i], "--test-simd")) { |
| 286 | #ifdef PFFFT_ENABLE_FLOAT |
| 287 | int numErrs = validate_pffft_simd_ex(stdout); |
| 288 | #else |
| 289 | int numErrs = validate_pffftd_simd_ex(stdout); |
| 290 | #endif |
| 291 | fprintf( ( numErrs != 0 ? stderr : stdout ), "validate_pffft_simd_ex() returned %d errors!\n", numErrs); |
| 292 | return ( numErrs > 0 ? 1 : 0 ); |
| 293 | } |
| 294 | } |
| 295 | |
hayati ayguen | e6cffc9 | 2020-02-28 19:57:30 +0100 | [diff] [blame] | 296 | resNextPw2 = 0; |
| 297 | resIsPw2 = 0; |
| 298 | for ( k = 0; k < (sizeof(inp_power_of_two)/sizeof(inp_power_of_two[0])); ++k) { |
hayati ayguen | ca11241 | 2020-04-13 00:19:40 +0200 | [diff] [blame] | 299 | #ifdef PFFFT_ENABLE_FLOAT |
hayati ayguen | e6cffc9 | 2020-02-28 19:57:30 +0100 | [diff] [blame] | 300 | N = pffft_next_power_of_two(inp_power_of_two[k]); |
hayati ayguen | ca11241 | 2020-04-13 00:19:40 +0200 | [diff] [blame] | 301 | #else |
| 302 | N = pffftd_next_power_of_two(inp_power_of_two[k]); |
| 303 | #endif |
hayati ayguen | e6cffc9 | 2020-02-28 19:57:30 +0100 | [diff] [blame] | 304 | if (N != ref_power_of_two[k]) { |
| 305 | resNextPw2 = 1; |
| 306 | printf("pffft_next_power_of_two(%d) does deliver %d, which is not reference result %d!\n", |
| 307 | inp_power_of_two[k], N, ref_power_of_two[k] ); |
| 308 | } |
| 309 | |
hayati ayguen | ca11241 | 2020-04-13 00:19:40 +0200 | [diff] [blame] | 310 | #ifdef PFFFT_ENABLE_FLOAT |
hayati ayguen | e6cffc9 | 2020-02-28 19:57:30 +0100 | [diff] [blame] | 311 | result = pffft_is_power_of_two(inp_power_of_two[k]); |
hayati ayguen | ca11241 | 2020-04-13 00:19:40 +0200 | [diff] [blame] | 312 | #else |
| 313 | result = pffftd_is_power_of_two(inp_power_of_two[k]); |
| 314 | #endif |
hayati ayguen | e6cffc9 | 2020-02-28 19:57:30 +0100 | [diff] [blame] | 315 | if (inp_power_of_two[k] == ref_power_of_two[k]) { |
| 316 | if (!result) { |
| 317 | resIsPw2 = 1; |
| 318 | printf("pffft_is_power_of_two(%d) delivers false; expected true!\n", inp_power_of_two[k]); |
| 319 | } |
| 320 | } else { |
| 321 | if (result) { |
| 322 | resIsPw2 = 1; |
| 323 | printf("pffft_is_power_of_two(%d) delivers true; expected false!\n", inp_power_of_two[k]); |
| 324 | } |
| 325 | } |
| 326 | } |
| 327 | if (!resNextPw2) |
| 328 | printf("tests for pffft_next_power_of_two() succeeded successfully.\n"); |
| 329 | if (!resIsPw2) |
| 330 | printf("tests for pffft_is_power_of_two() succeeded successfully.\n"); |
| 331 | |
| 332 | resFFT = 0; |
hayati ayguen | b2d2936 | 2020-01-02 00:06:09 +0100 | [diff] [blame] | 333 | for ( N = 32; N <= 65536; N *= 2 ) |
| 334 | { |
| 335 | result = test(N, 1 /* cplx fft */, 1 /* useOrdered */); |
| 336 | resN = result; |
hayati ayguen | e6cffc9 | 2020-02-28 19:57:30 +0100 | [diff] [blame] | 337 | resFFT |= result; |
hayati ayguen | b2d2936 | 2020-01-02 00:06:09 +0100 | [diff] [blame] | 338 | |
| 339 | result = test(N, 0 /* cplx fft */, 1 /* useOrdered */); |
| 340 | resN |= result; |
hayati ayguen | e6cffc9 | 2020-02-28 19:57:30 +0100 | [diff] [blame] | 341 | resFFT |= result; |
hayati ayguen | b2d2936 | 2020-01-02 00:06:09 +0100 | [diff] [blame] | 342 | |
| 343 | result = test(N, 1 /* cplx fft */, 0 /* useOrdered */); |
| 344 | resN |= result; |
hayati ayguen | e6cffc9 | 2020-02-28 19:57:30 +0100 | [diff] [blame] | 345 | resFFT |= result; |
hayati ayguen | b2d2936 | 2020-01-02 00:06:09 +0100 | [diff] [blame] | 346 | |
| 347 | result = test(N, 0 /* cplx fft */, 0 /* useOrdered */); |
| 348 | resN |= result; |
hayati ayguen | e6cffc9 | 2020-02-28 19:57:30 +0100 | [diff] [blame] | 349 | resFFT |= result; |
hayati ayguen | b2d2936 | 2020-01-02 00:06:09 +0100 | [diff] [blame] | 350 | |
| 351 | if (!resN) |
| 352 | printf("tests for size %d succeeded successfully.\n", N); |
| 353 | } |
| 354 | |
hayati ayguen | ca11241 | 2020-04-13 00:19:40 +0200 | [diff] [blame] | 355 | if (!resFFT) { |
| 356 | #ifdef PFFFT_ENABLE_FLOAT |
| 357 | printf("all pffft transform tests (FORWARD/BACKWARD, REAL/COMPLEX, float) succeeded successfully.\n"); |
| 358 | #else |
| 359 | printf("all pffft transform tests (FORWARD/BACKWARD, REAL/COMPLEX, double) succeeded successfully.\n"); |
| 360 | #endif |
| 361 | } |
hayati ayguen | e6cffc9 | 2020-02-28 19:57:30 +0100 | [diff] [blame] | 362 | |
| 363 | resAll = resNextPw2 | resIsPw2 | resFFT; |
hayati ayguen | b2d2936 | 2020-01-02 00:06:09 +0100 | [diff] [blame] | 364 | if (!resAll) |
| 365 | printf("all tests succeeded successfully.\n"); |
hayati ayguen | e6cffc9 | 2020-02-28 19:57:30 +0100 | [diff] [blame] | 366 | else |
| 367 | printf("there are failed tests!\n"); |
hayati ayguen | b2d2936 | 2020-01-02 00:06:09 +0100 | [diff] [blame] | 368 | |
| 369 | return resAll; |
| 370 | } |
| 371 | |