diff options
author | 2006-02-09 14:55:48 +0000 | |
---|---|---|
committer | 2006-02-09 14:55:48 +0000 | |
commit | 2ae43be7b9d4118335c9d2cef6e098f9b9f807fe (patch) | |
tree | bbb5e49ccbf7e3614966571acc317f8d318fecad /test/cminor/fft.cm |
Initial import of compcert
git-svn-id: https://yquem.inria.fr/compcert/svn/compcert/trunk@1 fca1b0fc-160b-0410-b1d3-a4f43f01ea2e
Diffstat (limited to 'test/cminor/fft.cm')
-rw-r--r-- | test/cminor/fft.cm | 149 |
1 files changed, 149 insertions, 0 deletions
diff --git a/test/cminor/fft.cm b/test/cminor/fft.cm new file mode 100644 index 0000000..b4e6a40 --- /dev/null +++ b/test/cminor/fft.cm @@ -0,0 +1,149 @@ +/********************************************************/ +/* A Duhamel-Hollman split-radix dif fft */ +/* Ref: Electronics Letters, Jan. 5, 1984 */ +/* Complex input and output data in arrays x and y */ +/* Length is n. */ +/********************************************************/ + +"dfft"(x, y, np): int /*float ptr*/ -> int /*float ptr*/ -> int -> int +{ + var px, py, /*float ptr*/ + i, j, k, m, n, i0, i1, i2, i3, is, id, n1, n2, n4, a, e, a3, /*int*/ + cc1, ss1, cc3, ss3, r1, r2, s1, s2, s3, xt, tpi /*float*/; + + px = x - 8; + py = y - 8; + i = 2; + m = 1; + + {{ loop { + if (! (i < np)) exit; + i = i+i; + m = m+1; + } }} + + n = i; + + if (n != np) + { + i = np + 1; + {{ loop { + if (! (i <= n)) exit; + float64[px + i*8] = 0.0; + float64[py + i*8] = 0.0; + i = i + 1; + } }} + } + + n2 = n+n; + tpi = 2.0 *f 3.14159265358979323846; + k = 1; + {{ loop { + if (! (k <= m - 1)) exit; + n2 = n2 / 2; + n4 = n2 / 4; + e = tpi /f floatofint(n2); + a = 0.0; + + j = 1; + {{ loop { + if (! (j <= n4)) exit; + cc1 = "cos_static"(a) : float -> float; + ss1 = "sin_static"(a) : float -> float; + a3 = 3.0 *f a; + cc3 = "cos_static"(a3) : float -> float; + ss3 = "sin_static"(a3) : float -> float; + a = e *f floatofint(j); + is = j; + id = 2 * n2; + + {{ loop { + if (! ( is < n )) exit; + i0 = is; + {{ loop { + /* DEBUG "trace"(); */ + if (! (i0 <= n - 1)) exit; + i1 = i0 + n4; /*DEBUG "print_int"(i1); */ + i2 = i1 + n4; /*DEBUG "print_int"(i2); */ + i3 = i2 + n4; /*DEBUG "print_int"(i3); */ + r1 = float64[px+i0*8] -f float64[px+i2*8]; + /* DEBUG "print_float"(r1); */ + float64[px+i0*8] = float64[px+i0*8] +f float64[px+i2*8]; + r2 = float64[px+i1*8] -f float64[px+i3*8]; + /* DEBUG "print_float"(r2); */ + float64[px+i1*8] = float64[px+i1*8] +f float64[px+i3*8]; + s1 = float64[py+i0*8] -f float64[py+i2*8]; + /* DEBUG "print_float"(s1); */ + float64[py+i0*8] = float64[py+i0*8] +f float64[py+i2*8]; + s2 = float64[py+i1*8] -f float64[py+i3*8]; + /* DEBUG "print_float"(s2); */ + float64[py+i1*8] = float64[py+i1*8] +f float64[py+i3*8]; + s3 = r1 -f s2; r1 = r1 +f s2; + s2 = r2 -f s1; r2 = r2 +f s1; + float64[px+i2*8] = (r1 *f cc1) -f (s2 *f ss1); + float64[py+i2*8] = (-f s2 *f cc1) -f (r1 *f ss1); + float64[px+i3*8] = (s3 *f cc3) +f (r2 *f ss3); + float64[py+i3*8] = (r2 *f cc3) -f (s3 *f ss3); + i0 = i0 + id; + } }} + is = 2 * id - n2 + j; + id = 4 * id; + } }} + j = j + 1; + } }} + k = k + 1; + } }} + +/************************************/ +/* Last stage, length=2 butterfly */ +/************************************/ + is = 1; + id = 4; + + {{ loop { + if (! ( is < n)) exit; + i0 = is; + {{ loop { + if (! (i0 <= n)) exit; + i1 = i0 + 1; + r1 = float64[px+i0*8]; + float64[px+i0*8] = r1 +f float64[px+i1*8]; + float64[px+i1*8] = r1 -f float64[px+i1*8]; + r1 = float64[py+i0*8]; + float64[py+i0*8] = r1 +f float64[py+i1*8]; + float64[py+i1*8] = r1 -f float64[py+i1*8]; + i0 = i0 + id; + } }} + is = 2*id - 1; + id = 4 * id; + } }} + +/*************************/ +/* Bit reverse counter */ +/*************************/ + j = 1; + n1 = n - 1; + + i = 1; + {{ loop { + if (! (i <= n1)) exit; + if (i < j) { + xt = float64[px+j*8]; + float64[px+j*8] = float64[px+i*8]; + float64[px+i*8] = xt; + xt = float64[py+j*8]; + float64[py+j*8] = float64[py+i*8]; + float64[py+i*8] = xt; + } + k = n / 2; + {{ loop { + if (! (k < j)) exit; + j = j - k; + k = k / 2; + } }} + j = j + k; + i = i + 1; + } }} + + return(n); +} |