summaryrefslogtreecommitdiff
path: root/test/files/specialized/fft.scala
diff options
context:
space:
mode:
Diffstat (limited to 'test/files/specialized/fft.scala')
-rw-r--r--test/files/specialized/fft.scala149
1 files changed, 149 insertions, 0 deletions
diff --git a/test/files/specialized/fft.scala b/test/files/specialized/fft.scala
new file mode 100644
index 0000000000..76029832a9
--- /dev/null
+++ b/test/files/specialized/fft.scala
@@ -0,0 +1,149 @@
+
+/*
+ * http://local.wasp.uwa.edu.au/~pbourke/miscellaneous/dft/
+ Modification of Paul Bourkes FFT code by Peter Cusack
+ to utilise the Microsoft complex type.
+
+ This computes an in-place complex-to-complex FFT
+ x and y are the real and imaginary arrays of 2^m points.
+ dir = 1 gives forward transform
+ dir = -1 gives reverse transform
+*/
+
+import Math.{sqrt, pow}
+
+/** Test that specialization handles tuples. Perform FFT transformation
+ * using pairs to represent complex numbers.
+ */
+object Test {
+ type Complex = (Double, Double)
+
+ def swap(x: Array[Complex], i: Int, j: Int) {
+ val tmp = x(i)
+ x(i) = x(j)
+ x(j) = tmp
+ }
+
+ def times(x: Complex, y: Complex): Complex =
+ (x._1 * y._1 - x._2 * y._2, x._1 * y._2 + x._2 * y._1)
+
+ def div(x: Complex, y: Complex): Complex = {
+ val num = pow(y._1, 2) + pow(y._2, 2)
+ ((x._1 * y._1 + x._2 * y._2)/num,
+ (x._2 * y._1 - x._1 * y._2)/num)
+ }
+
+ def div(x: Complex, y: Long) =
+ (x._1 / y, x._2 / y)
+
+ def add(x: Complex, y: Complex) =
+ (x._1 + y._1, x._2 + y._2)
+
+ def minus(x: Complex, y: Complex) =
+ (x._1 - y._1, x._2 - y._2)
+
+ def FFT(dir: Int, m: Long, x: Array[(Double, Double)]) {
+ var i, i1, i2,j, k, l, l1, l2, n = 0l
+// complex <double> tx, t1, u, c;
+ var tx, t1, u, c = (0.0, 0.0)
+
+ /*Calculate the number of points */
+ n = 1
+ for (i <- 0l until m)
+ n <<= 1
+
+ /* Do the bit reversal */
+ i2 = n >> 1
+ j = 0
+
+ for (i <- 0l until (n - 1)) {
+ if (i < j)
+ swap(x, i.toInt, j.toInt);
+
+ k = i2;
+
+ while (k <= j) {
+ j -= k;
+ k >>= 1;
+ }
+
+ j += k;
+ }
+
+ /* Compute the FFT */
+ // c.real(-1.0);
+ // c.imag(0.0);
+ c = (-1.0, 0.0)
+ l2 = 1
+ for (l <- 0l until m) {
+ l1 = l2
+ l2 <<= 1;
+ // u.real(1.0);
+ // u.imag(0.0);
+ u = (1.0, 0.0)
+
+ for (j <- 0l until l1) {
+ for (i <- j.until(n, l2)) {
+ i1 = i + l1;
+ t1 = times(u, x(i1.toInt))
+ x(i1.toInt) = minus(x(i.toInt), t1)
+ x(i.toInt) = add(x(i.toInt), t1)
+ }
+
+ u = times(u, c)
+ }
+
+ // c.imag(sqrt((1.0 - c.real()) / 2.0));
+ c = (c._1, sqrt( (1.0 - c._1) / 2.0 ))
+ // if (dir == 1)
+ // c.imag(-c.imag());
+ if (dir == 1)
+ c = (c._1, -c._2)
+
+ // c.real(sqrt((1.0 + c.real()) / 2.0));
+ c = (sqrt( (1.0 + c._1) / 2.0), c._2)
+ }
+
+ /* Scaling for forward transform */
+ if (dir == 1) {
+ for (i <- 0l until n)
+ x(i.toInt) = div(x(i.toInt), n)
+ }
+ }
+
+ def run() {
+ FFT(1, 16, data)
+ }
+ var data: Array[Complex] = null
+
+ def inputFileName = {
+ val cwd = System.getProperty("partest.cwd")
+ if (cwd ne null) cwd + java.io.File.separator + "input2.txt"
+ else "input2.txt"
+ }
+
+ def setUp {
+// print("Loading from %s.. ".format(inputFileName))
+ val f = io.Source.fromFile(inputFileName)
+ val lines = f.getLines
+ val n = lines.next.toInt
+ data = new Array[Complex](n)
+ var i = 0
+ for (line <- lines if line != "") {
+ val pair = line.trim.split(" ")
+ data(i) = (pair(0).trim.toDouble, pair(1).trim.toDouble)
+ i += 1
+ }
+// println("[loaded]")
+ println("Processing " + n + " items")
+ }
+
+ def main(args: Array[String]) {
+ setUp
+ run()
+
+ println("Boxed doubles: " + runtime.BoxesRunTime.doubleBoxCount)
+ println("Boxed ints: " + runtime.BoxesRunTime.integerBoxCount)
+ println("Boxed longs: " + runtime.BoxesRunTime.longBoxCount)
+ }
+}