Style fixes for FFT code.

This commit is contained in:
Oleksii Trekhleb 2018-08-15 13:27:50 +03:00
parent a88c45a329
commit ac9920a0f3
2 changed files with 65 additions and 51 deletions

View File

@ -2,50 +2,65 @@ import fastFourierTransform from '../fastFourierTransform';
import ComplexNumber from '../../complex-number/ComplexNumber'; import ComplexNumber from '../../complex-number/ComplexNumber';
/** /**
* @param {ComplexNumber[]} [seq1] * @param {ComplexNumber[]} sequence1
* @param {ComplexNumber[]} [seq2] * @param {ComplexNumber[]} sequence2
* @param {Number} [eps] * @param {Number} delta
* @return {boolean} * @return {boolean}
*/ */
function approximatelyEqual(seq1, seq2, eps) { function sequencesApproximatelyEqual(sequence1, sequence2, delta) {
if (seq1.length !== seq2.length) { return false; } if (sequence1.length !== sequence2.length) {
return false;
}
for (let i = 0; i < seq1.length; i += 1) { for (let numberIndex = 0; numberIndex < sequence1.length; numberIndex += 1) {
if (Math.abs(seq1[i].real - seq2[i].real) > eps) { return false; } if (Math.abs(sequence1[numberIndex].re - sequence2[numberIndex].re) > delta) {
if (Math.abs(seq1[i].complex - seq2[i].complex) > eps) { return false; } return false;
}
if (Math.abs(sequence1[numberIndex].im - sequence2[numberIndex].im) > delta) {
return false;
}
} }
return true; return true;
} }
const delta = 1e-6;
describe('fastFourierTransform', () => { describe('fastFourierTransform', () => {
it('should calculate the radix-2 discrete fourier transform after zero padding', () => { it('should calculate the radix-2 discrete fourier transform after zero padding', () => {
const eps = 1e-6; const input = [new ComplexNumber({ re: 0, im: 0 })];
const in1 = [new ComplexNumber({ re: 0, im: 0 })]; const expectedOutput = [new ComplexNumber({ re: 0, im: 0 })];
const expOut1 = [new ComplexNumber({ re: 0, im: 0 })]; const output = fastFourierTransform(input);
const out1 = fastFourierTransform(in1); const invertedOutput = fastFourierTransform(output, true);
const invOut1 = fastFourierTransform(out1, true);
expect(approximatelyEqual(expOut1, out1, eps)).toBe(true);
expect(approximatelyEqual(in1, invOut1, eps)).toBe(true);
const in2 = [ expect(sequencesApproximatelyEqual(expectedOutput, output, delta)).toBe(true);
expect(sequencesApproximatelyEqual(input, invertedOutput, delta)).toBe(true);
});
it('should calculate the radix-2 discrete fourier transform after zero padding', () => {
const input = [
new ComplexNumber({ re: 1, im: 2 }), new ComplexNumber({ re: 1, im: 2 }),
new ComplexNumber({ re: 2, im: 3 }), new ComplexNumber({ re: 2, im: 3 }),
new ComplexNumber({ re: 8, im: 4 }), new ComplexNumber({ re: 8, im: 4 }),
]; ];
const expOut2 = [ const expectedOutput = [
new ComplexNumber({ re: 11, im: 9 }), new ComplexNumber({ re: 11, im: 9 }),
new ComplexNumber({ re: -10, im: 0 }), new ComplexNumber({ re: -10, im: 0 }),
new ComplexNumber({ re: 7, im: 3 }), new ComplexNumber({ re: 7, im: 3 }),
new ComplexNumber({ re: -4, im: -4 }), new ComplexNumber({ re: -4, im: -4 }),
]; ];
const out2 = fastFourierTransform(in2);
const invOut2 = fastFourierTransform(out2, true);
expect(approximatelyEqual(expOut2, out2, eps)).toBe(true);
expect(approximatelyEqual(in2, invOut2, eps)).toBe(true);
const in3 = [ const output = fastFourierTransform(input);
const invertedOut = fastFourierTransform(output, true);
expect(sequencesApproximatelyEqual(expectedOutput, output, delta)).toBe(true);
expect(sequencesApproximatelyEqual(input, invertedOut, delta)).toBe(true);
});
it('should calculate the radix-2 discrete fourier transform after zero padding', () => {
const input = [
new ComplexNumber({ re: -83656.9359385182, im: 98724.08038374918 }), new ComplexNumber({ re: -83656.9359385182, im: 98724.08038374918 }),
new ComplexNumber({ re: -47537.415125808424, im: 88441.58381765135 }), new ComplexNumber({ re: -47537.415125808424, im: 88441.58381765135 }),
new ComplexNumber({ re: -24849.657029355192, im: -72621.79007878687 }), new ComplexNumber({ re: -24849.657029355192, im: -72621.79007878687 }),
@ -58,7 +73,7 @@ describe('fastFourierTransform', () => {
new ComplexNumber({ re: -39327.43830818355, im: 30611.949874562706 }), new ComplexNumber({ re: -39327.43830818355, im: 30611.949874562706 }),
]; ];
const expOut3 = [ const expectedOutput = [
new ComplexNumber({ re: -203215.3322151, im: -100242.4827503 }), new ComplexNumber({ re: -203215.3322151, im: -100242.4827503 }),
new ComplexNumber({ re: 99217.0805705, im: 270646.9331932 }), new ComplexNumber({ re: 99217.0805705, im: 270646.9331932 }),
new ComplexNumber({ re: -305990.9040412, im: 68224.8435751 }), new ComplexNumber({ re: -305990.9040412, im: 68224.8435751 }),
@ -77,9 +92,10 @@ describe('fastFourierTransform', () => {
new ComplexNumber({ re: -179002.5662573, im: 239821.0124341 }), new ComplexNumber({ re: -179002.5662573, im: 239821.0124341 }),
]; ];
const out3 = fastFourierTransform(in3); const output = fastFourierTransform(input);
const invOut3 = fastFourierTransform(out3, true); const invertedOutput = fastFourierTransform(output, true);
expect(approximatelyEqual(expOut3, out3, eps)).toBe(true);
expect(approximatelyEqual(in3, invOut3, eps)).toBe(true); expect(sequencesApproximatelyEqual(expectedOutput, output, delta)).toBe(true);
expect(sequencesApproximatelyEqual(input, invertedOutput, delta)).toBe(true);
}); });
}); });

View File

@ -10,10 +10,15 @@ import bitLength from '../bits/bitLength';
*/ */
function reverseBits(input, bitsCount) { function reverseBits(input, bitsCount) {
let reversedBits = 0; let reversedBits = 0;
for (let i = 0; i < bitsCount; i += 1) { for (let i = 0; i < bitsCount; i += 1) {
reversedBits *= 2; reversedBits *= 2;
if (Math.floor(input / (1 << i)) % 2 === 1) { reversedBits += 1; }
if (Math.floor(input / (1 << i)) % 2 === 1) {
reversedBits += 1;
}
} }
return reversedBits; return reversedBits;
} }
@ -21,8 +26,8 @@ function reverseBits(input, bitsCount) {
* Returns the radix-2 fast fourier transform of the given array. * Returns the radix-2 fast fourier transform of the given array.
* Optionally computes the radix-2 inverse fast fourier transform. * Optionally computes the radix-2 inverse fast fourier transform.
* *
* @param {ComplexNumber[]} [inputData] * @param {ComplexNumber[]} inputData
* @param {Boolean} [inverse] * @param {boolean} [inverse]
* @return {ComplexNumber[]} * @return {ComplexNumber[]}
*/ */
export default function fastFourierTransform(inputData, inverse = false) { export default function fastFourierTransform(inputData, inverse = false) {
@ -30,48 +35,41 @@ export default function fastFourierTransform(inputData, inverse = false) {
const N = 1 << bitsCount; const N = 1 << bitsCount;
while (inputData.length < N) { while (inputData.length < N) {
inputData.push(new ComplexNumber({ inputData.push(new ComplexNumber());
real: 0,
imaginary: 0,
}));
} }
const output = []; const output = [];
for (let i = 0; i < N; i += 1) { output[i] = inputData[reverseBits(i, bitsCount)]; } for (let i = 0; i < N; i += 1) {
output[i] = inputData[reverseBits(i, bitsCount)];
}
for (let blockLength = 2; blockLength <= N; blockLength *= 2) { for (let blockLength = 2; blockLength <= N; blockLength *= 2) {
let phaseStep; const imaginarySign = inverse ? -1 : 1;
if (inverse) { const phaseStep = new ComplexNumber({
phaseStep = new ComplexNumber({ re: Math.cos(2 * Math.PI / blockLength),
real: Math.cos(2 * Math.PI / blockLength), im: imaginarySign * Math.sin(2 * Math.PI / blockLength),
imaginary: -1 * Math.sin(2 * Math.PI / blockLength), });
});
} else {
phaseStep = new ComplexNumber({
real: Math.cos(2 * Math.PI / blockLength),
imaginary: Math.sin(2 * Math.PI / blockLength),
});
}
for (let blockStart = 0; blockStart < N; blockStart += blockLength) { for (let blockStart = 0; blockStart < N; blockStart += blockLength) {
let phase = new ComplexNumber({ let phase = new ComplexNumber({ re: 1, im: 0 });
real: 1,
imaginary: 0,
});
for (let idx = blockStart; idx < blockStart + blockLength / 2; idx += 1) { for (let idx = blockStart; idx < blockStart + blockLength / 2; idx += 1) {
const upd1 = output[idx].add(output[idx + blockLength / 2].multiply(phase)); const upd1 = output[idx].add(output[idx + blockLength / 2].multiply(phase));
const upd2 = output[idx].subtract(output[idx + blockLength / 2].multiply(phase)); const upd2 = output[idx].subtract(output[idx + blockLength / 2].multiply(phase));
output[idx] = upd1; output[idx] = upd1;
output[idx + blockLength / 2] = upd2; output[idx + blockLength / 2] = upd2;
phase = phase.multiply(phaseStep); phase = phase.multiply(phaseStep);
} }
} }
} }
if (inverse) { if (inverse) {
for (let idx = 0; idx < N; idx += 1) { for (let idx = 0; idx < N; idx += 1) {
output[idx] /= N; output[idx] /= N;
} }
} }
return output; return output;
} }