Change signature

Author: dy

benchmark/run.js

@@ -62,6 +62,11 @@ const { rfft: kissRfft } = await import('kissfft-wasm')
 const kissInput = new Float64Array(input)
 const kissRun = () => kissRfft(kissInput)
 
+// als-fft
+const AlsFFT = require('als-fft')
+const alsInput = Array.from(input)
+const alsRun = () => AlsFFT.fft(alsInput)
+
 // fft-js (naive Cooley-Tukey, educational)
 const fftjs2 = require('fft-js')
 const fftjs2Input = Array.from(input)
@@ -77,6 +82,7 @@ const libs = [
 	['ndarray-fft', ndfftRun],
 	['ooura', oouraRun],
 	['kissfft-wasm', kissRun],
+	['als-fft', alsRun],
 	['fft-js', fftjs2Run],
 ]
 

changelog.md

@@ -11,7 +11,7 @@
 
 ## Added
 
-- `fft()` named export — returns complex DFT as `{ re, im }`, N/2+1 bins, unnormalized.
+- `fft()` named export — returns complex DFT as `[re, im]`, N/2+1 bins, unnormalized.
 - Optional output buffer parameter for both `rfft()` and `fft()`.
 
 ## Performance

index.js

@@ -17,7 +17,7 @@ function init(N) {
 	const spectrum = new Float64Array(half)
 	const im = new Float64Array(half + 1)
 	const re = x.subarray(0, half + 1) // zero-copy view into x
-	const complex = { re, im }
+	const complex = [re, im]
 	const bSi = 2 / N
 
 	// Precompute bit-reversal permutation table
@@ -221,8 +221,8 @@ export default function rfft(input, output) {
 /**
  * Compute complex spectrum of real-valued input (unnormalized DFT).
  * @param {ArrayLike<number>} input - length must be power of 2 (>= 2).
- * @param {{re: Float64Array, im: Float64Array}} [output] - Optional buffers (length N/2+1 each). If omitted, returns internal view.
- * @returns {{re: Float64Array, im: Float64Array}} Complex spectrum, N/2+1 bins (DC through Nyquist).
+ * @param {[Float64Array, Float64Array]} [output] - Optional [re, im] buffers (length N/2+1 each). If omitted, returns internal view.
+ * @returns {[Float64Array, Float64Array]} Complex spectrum [re, im], N/2+1 bins (DC through Nyquist).
  */
 export function fft(input, output) {
 	const entry = transform(input)
@@ -231,15 +231,15 @@ export function fft(input, output) {
 	const { x, complex } = entry
 
 	if (output) {
-		const re = output.re, im = output.im
+		const re = output[0], im = output[1]
 		for (let k = 0; k <= half; k++) re[k] = x[k]
 		im[0] = 0; im[half] = 0
 		for (let k = 1; k < half; k++) im[k] = x[N - k]
 		return output
 	}
 
 	// re is already a zero-copy view of x[0..half] — no copy needed
-	const im = complex.im
+	const im = complex[1]
 	im[0] = 0; im[half] = 0
 	for (let k = 1; k < half; k++) im[k] = x[N - k]
 

readme.md

@@ -15,7 +15,7 @@ const spectrum = rfft(waveform)
 import { fft } from 'fourier-transform'
 
 // complex DFT output (N/2+1 bins, unnormalized)
-const { re, im } = fft(waveform)
+const [re, im] = fft(waveform)
 ```
 
 ## API
@@ -32,9 +32,9 @@ Normalization: a unit-amplitude cosine at frequency bin *k* produces `spectrum[k
 
 ### `fft(input, output?)` — named export
 
-Returns complex DFT as `{ re, im }`, each `Float64Array` of length N/2+1 (DC through Nyquist).
+Returns complex DFT as `[re, im]`, each `Float64Array` of length N/2+1 (DC through Nyquist).
 
-- `output` — optional `{ re: Float64Array(N/2+1), im: Float64Array(N/2+1) }`.
+- `output` — optional `[Float64Array(N/2+1), Float64Array(N/2+1)]`.
 - Unnormalized: `X[k] = sum( x[n] * e^(-j*2*pi*k*n/N) )`.
 - DC and Nyquist bins always have `im = 0` (real input).
 
@@ -52,17 +52,18 @@ rfft(signal, out) // safe to keep
 N=4096 real-valued FFT, complex output, 20k iterations (lower is better):
 
 ```
-fft.js (indutny)          16.2µs  ×1.0  — radix-4, interleaved output
-fourier-transform         17.3µs  ×1.1  — split-radix, separate re/im
-ooura                     23.1µs  ×1.4  — Ooura C port
-ml-fft                    36.0µs  ×2.2
-dsp.js                    47.1µs  ×2.9  — our split-radix ancestor
-kissfft-wasm              49.4µs  ×3.1  — WASM KissFFT
-ndarray-fft               62.6µs  ×3.9
-fft-js                  2297.4µs  ×142   — naive recursive
+fft.js (indutny)          16.5µs  ×1.0  — radix-4, interleaved output
+fourier-transform         17.8µs  ×1.1  — split-radix, separate re/im
+ooura                     23.6µs  ×1.4  — Ooura C port
+ml-fft                    37.0µs  ×2.2
+dsp.js                    48.1µs  ×2.9  — our split-radix ancestor
+kissfft-wasm              49.4µs  ×3.0  — WASM KissFFT
+ndarray-fft               63.1µs  ×3.8
+als-fft                 2311.4µs  ×140
+fft-js                  2329.2µs  ×141  — naive recursive
 ```
 
-Raw transform speed is identical to fft.js. The ×1.1 gap is entirely the cost of returning separate `re`/`im` arrays vs interleaved output.
+Raw transform speed is identical to fft.js. The gap is the cost of returning separate `re`/`im` arrays vs interleaved output.
 
 `npm run benchmark` to reproduce.
 

test.js

@@ -174,14 +174,14 @@ test('internal buffer overwritten on repeated calls (view semantics)', () => {
 
 test('fft: returns N/2+1 complex bins', () => {
 	const N = 64
-	const { re, im } = fft(new Float32Array(N))
+	const [re, im] = fft(new Float32Array(N))
 	assert.equal(re.length, N / 2 + 1)
 	assert.equal(im.length, N / 2 + 1)
 })
 
 test('fft: DC signal → real-only X[0]=N', () => {
 	const N = 64
-	const { re, im } = fft(new Float32Array(N).fill(1))
+	const [re, im] = fft(new Float32Array(N).fill(1))
 	assert(Math.abs(re[0] - N) < EPSILON, `DC re: ${re[0]}`)
 	assert(Math.abs(im[0]) < EPSILON, `DC im: ${im[0]}`)
 	for (let k = 1; k <= N / 2; k++) {
@@ -195,7 +195,7 @@ test('fft: cosine → real component X[k] = N/2', () => {
 	for (const k of [1, 5, 32]) {
 		const input = new Float32Array(N)
 		for (let n = 0; n < N; n++) input[n] = Math.cos(2 * Math.PI * k * n / N)
-		const { re, im } = fft(input)
+		const [re, im] = fft(input)
 		assert(Math.abs(re[k] - N / 2) < EPSILON, `cos k=${k}: re=${re[k]}, expected ${N / 2}`)
 		assert(Math.abs(im[k]) < EPSILON, `cos k=${k}: im=${im[k]}, expected 0`)
 	}
@@ -206,7 +206,7 @@ test('fft: sine → imaginary component X[k] = -jN/2', () => {
 	for (const k of [1, 5, 63]) {
 		const input = new Float32Array(N)
 		for (let n = 0; n < N; n++) input[n] = Math.sin(2 * Math.PI * k * n / N)
-		const { re, im } = fft(input)
+		const [re, im] = fft(input)
 		assert(Math.abs(re[k]) < EPSILON, `sin k=${k}: re=${re[k]}, expected 0`)
 		assert(Math.abs(im[k] - (-N / 2)) < EPSILON, `sin k=${k}: im=${im[k]}, expected ${-N / 2}`)
 	}
@@ -216,7 +216,7 @@ test('fft: DC and Nyquist always have zero imaginary', () => {
 	for (const N of [4, 64, 512]) {
 		const input = new Float32Array(N)
 		for (let i = 0; i < N; i++) input[i] = Math.sin(i * 1.7) + Math.cos(i * 0.3)
-		const { im } = fft(input)
+		const [, im] = fft(input)
 		assert.equal(im[0], 0, `N=${N}: DC im must be 0`)
 		assert.equal(im[N / 2], 0, `N=${N}: Nyquist im must be 0`)
 	}
@@ -228,7 +228,7 @@ test('fft: magnitude matches rfft', () => {
 	for (let i = 0; i < N; i++) input[i] = Math.sin(i * 0.7) + Math.cos(i * 1.3)
 
 	const mag = new Float64Array(rfft(input))
-	const { re, im } = fft(input)
+	const [re, im] = fft(input)
 
 	const bSi = 2 / N
 	assert(Math.abs(mag[0] - Math.abs(bSi * re[0])) < EPSILON, `DC mismatch`)
@@ -240,10 +240,10 @@ test('fft: magnitude matches rfft', () => {
 
 test('fft: output buffer parameter', () => {
 	const N = 64
-	const out = { re: new Float64Array(N / 2 + 1), im: new Float64Array(N / 2 + 1) }
+	const out = [new Float64Array(N / 2 + 1), new Float64Array(N / 2 + 1)]
 	const ret = fft(new Float32Array(N).fill(1), out)
 	assert.equal(ret, out)
-	assert(Math.abs(out.re[0] - N) < EPSILON)
+	assert(Math.abs(out[0][0] - N) < EPSILON)
 })
 
 test('fft: view overwritten on repeated calls', () => {
@@ -254,10 +254,10 @@ test('fft: view overwritten on repeated calls', () => {
 		b[i] = Math.sin(2 * Math.PI * 10 * i / N)
 	}
 	const ra = fft(a)
-	assert(Math.abs(ra.re[5] - N / 2) < EPSILON)
+	assert(Math.abs(ra[0][5] - N / 2) < EPSILON)
 
 	fft(b) // overwrites ra since same N
 
-	assert(Math.abs(ra.re[5]) < EPSILON, 'ra.re[5] should be overwritten')
-	assert(Math.abs(ra.im[10] - (-N / 2)) < EPSILON, 'ra now reflects b')
+	assert(Math.abs(ra[0][5]) < EPSILON, 'ra[0][5] should be overwritten')
+	assert(Math.abs(ra[1][10] - (-N / 2)) < EPSILON, 'ra now reflects b')
 })