const dctSize = 512; // Default DCT size, read from header // --- Utility Functions for Audio Processing --- // JavaScript equivalent of C++ idct_512 function javascript_idct_512(input) { const output = new Float32Array(dctSize); const PI = Math.PI; const N = dctSize; for (let n = 0; n < N; ++n) { let sum = input[0] / 2.0; for (let k = 1; k < N; ++k) { sum += input[k] * Math.cos((PI / N) * k * (n + 0.5)); } output[n] = sum * (2.0 / N); } return output; } // Hanning window for smooth audio transitions (JavaScript equivalent) function hanningWindow(size) { const window = new Float32Array(size); const PI = Math.PI; for (let i = 0; i < size; i++) { window[i] = 0.5 * (1 - Math.cos((2 * PI * i) / (size - 1))); } return window; } const hanningWindowArray = hanningWindow(dctSize); // Pre-calculate window // ============================================================================ // FFT-based DCT/IDCT Implementation // ============================================================================ // Fast Fourier Transform using Radix-2 Cooley-Tukey algorithm // This implementation MUST match the C++ version in src/audio/fft.cc exactly // Bit-reversal permutation (in-place) // Reorders array elements by reversing their binary indices function bitReversePermute(real, imag, N) { // Calculate number of bits needed let temp_bits = N; let num_bits = 0; while (temp_bits > 1) { temp_bits >>= 1; num_bits++; } for (let i = 0; i < N; i++) { // Compute bit-reversed index let j = 0; let temp = i; for (let b = 0; b < num_bits; b++) { j = (j << 1) | (temp & 1); temp >>= 1; } // Swap if j > i (to avoid swapping twice) if (j > i) { const tmp_real = real[i]; const tmp_imag = imag[i]; real[i] = real[j]; imag[i] = imag[j]; real[j] = tmp_real; imag[j] = tmp_imag; } } } // In-place radix-2 FFT (after bit-reversal) // direction: +1 for forward FFT, -1 for inverse FFT function fftRadix2(real, imag, N, direction) { const PI = Math.PI; // Butterfly operations for (let stage_size = 2; stage_size <= N; stage_size *= 2) { const half_stage = stage_size / 2; const angle = direction * 2.0 * PI / stage_size; // Precompute twiddle factors for this stage let wr = 1.0; let wi = 0.0; const wr_delta = Math.cos(angle); const wi_delta = Math.sin(angle); for (let k = 0; k < half_stage; k++) { // Apply butterfly to all groups at this stage for (let group_start = k; group_start < N; group_start += stage_size) { const i = group_start; const j = group_start + half_stage; // Complex multiplication: (real[j] + i*imag[j]) * (wr + i*wi) const temp_real = real[j] * wr - imag[j] * wi; const temp_imag = real[j] * wi + imag[j] * wr; // Butterfly operation real[j] = real[i] - temp_real; imag[j] = imag[i] - temp_imag; real[i] = real[i] + temp_real; imag[i] = imag[i] + temp_imag; } // Update twiddle factor for next k (rotation) const wr_old = wr; wr = wr_old * wr_delta - wi * wi_delta; wi = wr_old * wi_delta + wi * wr_delta; } } } // Forward FFT: Time domain → Frequency domain function fftForward(real, imag, N) { bitReversePermute(real, imag, N); fftRadix2(real, imag, N, +1); } // Inverse FFT: Frequency domain → Time domain function fftInverse(real, imag, N) { bitReversePermute(real, imag, N); fftRadix2(real, imag, N, -1); // Scale by 1/N const scale = 1.0 / N; for (let i = 0; i < N; i++) { real[i] *= scale; imag[i] *= scale; } } // DCT-II via FFT using double-and-mirror method (matches C++ dct_fft) // This is a more robust algorithm that avoids reordering issues function javascript_dct_fft(input, N) { const PI = Math.PI; // Allocate arrays for 2N-point FFT const M = 2 * N; const real = new Float32Array(M); const imag = new Float32Array(M); // Pack input: [x[0], x[1], ..., x[N-1], x[N-1], x[N-2], ..., x[1]] // This creates even symmetry for real-valued DCT for (let i = 0; i < N; i++) { real[i] = input[i]; } for (let i = 0; i < N; i++) { real[N + i] = input[N - 1 - i]; } // imag is already zeros (Float32Array default) // Apply 2N-point FFT fftForward(real, imag, M); // Extract DCT coefficients // DCT[k] = Re{FFT[k] * exp(-j*pi*k/(2*N))} * normalization // Note: Need to divide by 2 because we doubled the signal length const output = new Float32Array(N); for (let k = 0; k < N; k++) { const angle = -PI * k / (2.0 * N); const wr = Math.cos(angle); const wi = Math.sin(angle); // Complex multiplication: (real + j*imag) * (wr + j*wi) // Real part: real*wr - imag*wi const dct_value = real[k] * wr - imag[k] * wi; // Apply DCT-II normalization (divide by 2 for double-length FFT) if (k === 0) { output[k] = dct_value * Math.sqrt(1.0 / N) / 2.0; } else { output[k] = dct_value * Math.sqrt(2.0 / N) / 2.0; } } return output; } // IDCT (Inverse DCT-II) via FFT using double-and-mirror method (matches C++ idct_fft) function javascript_idct_fft(input, N) { const PI = Math.PI; // Allocate arrays for 2N-point FFT const M = 2 * N; const real = new Float32Array(M); const imag = new Float32Array(M); // Prepare FFT input from DCT coefficients // IDCT = Re{IFFT[DCT * exp(j*pi*k/(2*N))]} * 2 for (let k = 0; k < N; k++) { const angle = PI * k / (2.0 * N); // Positive for inverse const wr = Math.cos(angle); const wi = Math.sin(angle); // Apply inverse normalization let scaled_input; if (k === 0) { scaled_input = input[k] * Math.sqrt(N) * 2.0; } else { scaled_input = input[k] * Math.sqrt(N / 2.0) * 2.0; } // Complex multiplication: DCT[k] * exp(j*theta) real[k] = scaled_input * wr; imag[k] = scaled_input * wi; } // Fill second half with conjugate symmetry (for real output) for (let k = 1; k < N; k++) { real[M - k] = real[k]; imag[M - k] = -imag[k]; } // Apply inverse FFT fftInverse(real, imag, M); // Extract first N samples (real part only, imag should be ~0) const output = new Float32Array(N); for (let i = 0; i < N; i++) { output[i] = real[i]; } return output; } // Convenience wrappers for dctSize = 512 (backward compatible) function javascript_dct_512_fft(input) { return javascript_dct_fft(input, dctSize); } function javascript_idct_512_fft(input) { return javascript_idct_fft(input, dctSize); }