Mercurial > hg > js-dsp-test

annotate fft/nayuki/fft-test.html @ 40:223f770b5341 kissfft-double tip

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Try a double-precision kissfft
author Chris Cannam
date Wed, 07 Sep 2016 10:40:32 +0100
parents d7c216b6a84f
children
rev   line source
Chris@0 1 <!--
Chris@0 2 - FFT and convolution test (JavaScript)
Chris@0 3 -
Chris@0 4 - Copyright (c) 2014 Project Nayuki
Chris@0 5 - http://www.nayuki.io/page/free-small-fft-in-multiple-languages
Chris@0 6 -
Chris@0 7 - (MIT License)
Chris@0 8 - Permission is hereby granted, free of charge, to any person obtaining a copy of
Chris@0 9 - this software and associated documentation files (the "Software"), to deal in
Chris@0 10 - the Software without restriction, including without limitation the rights to
Chris@0 11 - use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of
Chris@0 12 - the Software, and to permit persons to whom the Software is furnished to do so,
Chris@0 13 - subject to the following conditions:
Chris@0 14 - * The above copyright notice and this permission notice shall be included in
Chris@0 15 - all copies or substantial portions of the Software.
Chris@0 16 - * The Software is provided "as is", without warranty of any kind, express or
Chris@0 17 - implied, including but not limited to the warranties of merchantability,
Chris@0 18 - fitness for a particular purpose and noninfringement. In no event shall the
Chris@0 19 - authors or copyright holders be liable for any claim, damages or other
Chris@0 20 - liability, whether in an action of contract, tort or otherwise, arising from,
Chris@0 21 - out of or in connection with the Software or the use or other dealings in the
Chris@0 22 - Software.
Chris@0 23 -->
Chris@0 24 <!DOCTYPE html>
Chris@0 25 <html>
Chris@0 26 <head>
Chris@0 27 <meta charset="UTF-8">
Chris@0 28 <title>FFT and convolution test (JavaScript)</title>
Chris@0 29 </head>
Chris@0 30 <body>
Chris@0 31 <script src="fft.js" type="application/javascript"></script>
Chris@0 32 <pre>
Chris@0 33 <script>
Chris@0 34 /* Main and test functions */
Chris@0 35
Chris@0 36 function main() {
Chris@0 37 // Test power-of-2 size FFTs
Chris@0 38 for (var i = 0; i <= 12; i++)
Chris@0 39 testFft(1 << i);
Chris@0 40
Chris@0 41 // Test small size FFTs
Chris@0 42 for (var i = 0; i < 30; i++)
Chris@0 43 testFft(i);
Chris@0 44
Chris@0 45 // Test diverse size FFTs
Chris@0 46 var prev = 0;
Chris@0 47 for (var i = 0; i <= 100; i++) {
Chris@0 48 var n = Math.round(Math.pow(1500, i / 100.0));
Chris@0 49 if (n > prev) {
Chris@0 50 testFft(n);
Chris@0 51 prev = n;
Chris@0 52 }
Chris@0 53 }
Chris@0 54
Chris@0 55 // Test power-of-2 size convolutions
Chris@0 56 for (var i = 0; i <= 12; i++)
Chris@0 57 testConvolution(1 << i);
Chris@0 58
Chris@0 59 // Test diverse size convolutions
Chris@0 60 prev = 0;
Chris@0 61 for (var i = 0; i <= 100; i++) {
Chris@0 62 var n = Math.round(Math.pow(1500, i / 100.0));
Chris@0 63 if (n > prev) {
Chris@0 64 testConvolution(n);
Chris@0 65 prev = n;
Chris@0 66 }
Chris@0 67 }
Chris@0 68
Chris@0 69 document.write("\nMax log err = " + maxLogError.toFixed(1));
Chris@0 70 document.write("\nTest " + (maxLogError < -10 ? "passed" : "failed"));
Chris@0 71 }
Chris@0 72
Chris@0 73
Chris@0 74 function testFft(size) {
Chris@0 75 var inputreal = randomReals(size);
Chris@0 76 var inputimag = randomReals(size);
Chris@0 77
Chris@0 78 var refoutreal = new Array(size);
Chris@0 79 var refoutimag = new Array(size);
Chris@0 80 naiveDft(inputreal, inputimag, refoutreal, refoutimag, false);
Chris@0 81
Chris@0 82 var actualoutreal = inputreal.slice(0);
Chris@0 83 var actualoutimag = inputimag.slice(0);
Chris@0 84 transform(actualoutreal, actualoutimag);
Chris@0 85
Chris@0 86 document.write("fftsize=" + size + " logerr=" + log10RmsErr(refoutreal, refoutimag, actualoutreal, actualoutimag).toFixed(1) + "\n");
Chris@0 87 }
Chris@0 88
Chris@0 89
Chris@0 90 function testConvolution(size) {
Chris@0 91 var input0real = randomReals(size);
Chris@0 92 var input0imag = randomReals(size);
Chris@0 93
Chris@0 94 var input1real = randomReals(size);
Chris@0 95 var input1imag = randomReals(size);
Chris@0 96
Chris@0 97 var refoutreal = new Array(size);
Chris@0 98 var refoutimag = new Array(size);
Chris@0 99 naiveConvolve(input0real, input0imag, input1real, input1imag, refoutreal, refoutimag);
Chris@0 100
Chris@0 101 var actualoutreal = new Array(size);
Chris@0 102 var actualoutimag = new Array(size);
Chris@0 103 convolveComplex(input0real, input0imag, input1real, input1imag, actualoutreal, actualoutimag);
Chris@0 104
Chris@0 105 document.write("convsize=" + size + " logerr=" + log10RmsErr(refoutreal, refoutimag, actualoutreal, actualoutimag).toFixed(1) + "\n");
Chris@0 106 }
Chris@0 107
Chris@0 108
Chris@0 109 /* Naive reference computation functions */
Chris@0 110
Chris@0 111 function naiveDft(inreal, inimag, outreal, outimag, inverse) {
Chris@0 112 if (inreal.length != inimag.length || inreal.length != outreal.length || outreal.length != outimag.length)
Chris@0 113 throw "Mismatched lengths";
Chris@0 114
Chris@0 115 var n = inreal.length;
Chris@0 116 var coef = (inverse ? 2 : -2) * Math.PI;
Chris@0 117 for (var k = 0; k < n; k++) { // For each output element
Chris@0 118 var sumreal = 0;
Chris@0 119 var sumimag = 0;
Chris@0 120 for (var t = 0; t < n; t++) { // For each input element
Chris@0 121 var angle = coef * (t * k % n) / n; // This is more accurate than t * k
Chris@0 122 sumreal += inreal[t]*Math.cos(angle) - inimag[t]*Math.sin(angle);
Chris@0 123 sumimag += inreal[t]*Math.sin(angle) + inimag[t]*Math.cos(angle);
Chris@0 124 }
Chris@0 125 outreal[k] = sumreal;
Chris@0 126 outimag[k] = sumimag;
Chris@0 127 }
Chris@0 128 }
Chris@0 129
Chris@0 130
Chris@0 131 function naiveConvolve(xreal, ximag, yreal, yimag, outreal, outimag) {
Chris@0 132 if (xreal.length != ximag.length || xreal.length != yreal.length || yreal.length != yimag.length || xreal.length != outreal.length || outreal.length != outimag.length)
Chris@0 133 throw "Mismatched lengths";
Chris@0 134
Chris@0 135 var n = xreal.length;
Chris@0 136 for (var i = 0; i < n; i++) {
Chris@0 137 var sumreal = 0;
Chris@0 138 var sumimag = 0;
Chris@0 139 for (var j = 0; j < n; j++) {
Chris@0 140 var k = (i - j + n) % n;
Chris@0 141 sumreal += xreal[k] * yreal[j] - ximag[k] * yimag[j];
Chris@0 142 sumimag += xreal[k] * yimag[j] + ximag[k] * yreal[j];
Chris@0 143 }
Chris@0 144 outreal[i] = sumreal;
Chris@0 145 outimag[i] = sumimag;
Chris@0 146 }
Chris@0 147 }
Chris@0 148
Chris@0 149
Chris@0 150 /* Utility functions */
Chris@0 151
Chris@0 152 var maxLogError = Number.NEGATIVE_INFINITY;
Chris@0 153
Chris@0 154 function log10RmsErr(xreal, ximag, yreal, yimag) {
Chris@0 155 if (xreal.length != ximag.length || xreal.length != yreal.length || yreal.length != yimag.length)
Chris@0 156 throw "Mismatched lengths";
Chris@0 157
Chris@0 158 var err = 0;
Chris@0 159 for (var i = 0; i < xreal.length; i++)
Chris@0 160 err += (xreal[i] - yreal[i]) * (xreal[i] - yreal[i]) + (ximag[i] - yimag[i]) * (ximag[i] - yimag[i]);
Chris@0 161 err = Math.sqrt(err / Math.max(xreal.length, 1)); // Now this is a root mean square (RMS) error
Chris@0 162 err = err > 0 ? Math.log(err) / Math.log(10) : -99;
Chris@0 163 maxLogError = Math.max(err, maxLogError);
Chris@0 164 return err;
Chris@0 165 }
Chris@0 166
Chris@0 167
Chris@0 168 function randomReals(size) {
Chris@0 169 var result = new Array(size);
Chris@0 170 for (var i = 0; i < result.length; i++)
Chris@0 171 result[i] = Math.random() * 2 - 1;
Chris@0 172 return result;
Chris@0 173 }
Chris@0 174
Chris@0 175
Chris@0 176 main();
Chris@0 177 </script>
Chris@0 178 </pre>
Chris@0 179 </body>
Chris@0 180 </html>