From e8a3b901c235ebeae910d9e4b73f52cf4aed4738 Mon Sep 17 00:00:00 2001 From: Cameron Taylor Date: Mon, 27 Sep 2021 22:30:38 -0400 Subject: [PATCH] FFT IN PROGRESS LOL --- CHANGELOG.md | 1 + source/ChartingState.hx | 4 + source/PlayState.hx | 28 ++++--- source/SpectogramSprite.hx | 152 +++++++++++++++++++++++++++++++++++- source/dsp/Complex.hx | 80 +++++++++++++++++++ source/dsp/FFT.hx | 154 +++++++++++++++++++++++++++++++++++++ source/dsp/OffsetArray.hx | 78 +++++++++++++++++++ source/dsp/Signal.hx | 110 ++++++++++++++++++++++++++ 8 files changed, 593 insertions(+), 14 deletions(-) create mode 100644 source/dsp/Complex.hx create mode 100644 source/dsp/FFT.hx create mode 100644 source/dsp/OffsetArray.hx create mode 100644 source/dsp/Signal.hx diff --git a/CHANGELOG.md b/CHANGELOG.md index 893637f59..f33bd2b12 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -11,6 +11,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 - 3 AWESOME PICO VS. DARNELL SONGS!! - Character offset editor / spritesheet viewer ## Changed +- Lerp'd the healthbar - Resetting from game over and "restart song" should be faster - Health gain is different depending on how accurate you hit notes! - slight less health gained on sustain notes diff --git a/source/ChartingState.hx b/source/ChartingState.hx index 13dd88cbc..5874d769c 100644 --- a/source/ChartingState.hx +++ b/source/ChartingState.hx @@ -3,6 +3,7 @@ package; import Conductor.BPMChangeEvent; import Section.SwagSection; import Song.SwagSong; +import dsp.FFT; import flixel.FlxSprite; import flixel.FlxStrip; import flixel.addons.display.FlxGridOverlay; @@ -36,6 +37,7 @@ import openfl.media.Sound; import openfl.net.FileReference; import openfl.utils.ByteArray; +using Lambda; using StringTools; using flixel.util.FlxSpriteUtil; @@ -399,6 +401,7 @@ class ChartingState extends MusicBeatState musSpec.x += 70; musSpec.daHeight = FlxG.height / 2; musSpec.scrollFactor.set(); + musSpec.visType = FREQUENCIES; add(musSpec); // trace(audioBuf.data.length); @@ -418,6 +421,7 @@ class ChartingState extends MusicBeatState { var vocalSpec:SpectogramSprite = new SpectogramSprite(voc, FlxG.random.color(0xFFAAAAAA, FlxColor.WHITE, 100)); vocalSpec.x = 70 - (50 * index); + vocalSpec.visType = FREQUENCIES; vocalSpec.daHeight = musSpec.daHeight; vocalSpec.y = vocalSpec.daHeight; vocalSpec.scrollFactor.set(); diff --git a/source/PlayState.hx b/source/PlayState.hx index eed3fb04a..67a43866b 100644 --- a/source/PlayState.hx +++ b/source/PlayState.hx @@ -1824,11 +1824,6 @@ class PlayState extends MusicBeatState { healthDisplay = FlxMath.lerp(healthDisplay, health, 0.15); - #if !debug - perfectMode = false; - #else - if (FlxG.keys.justPressed.H) - camHUD.visible = !camHUD.visible; if (needsReset) { resetCamFollow(); @@ -1843,13 +1838,8 @@ class PlayState extends MusicBeatState vocals.pause(); FlxG.sound.music.time = 0; - regenNoteData(); + regenNoteData(); // loads the note data from start health = 1; - - // resyncVocals(); - - // FlxG.sound.music.play(); - restartCountdownTimer(); needsReset = false; @@ -1865,6 +1855,22 @@ class PlayState extends MusicBeatState */ // sys.io.File.saveContent('./swag.png', png.readUTFBytes(png.length)); } + + #if !debug + perfectMode = false; + #else + if (FlxG.keys.justPressed.H) + camHUD.visible = !camHUD.visible; + if (FlxG.keys.justPressed.K) + { + @:privateAccess + var funnyData:Array = cast FlxG.sound.music._channel.__source.buffer.data; + + funnyData.reverse(); + + @:privateAccess + FlxG.sound.music._channel.__source.buffer.data = cast funnyData; + } #end // do this BEFORE super.update() so songPosition is accurate diff --git a/source/SpectogramSprite.hx b/source/SpectogramSprite.hx index 54848b920..03415ba3f 100644 --- a/source/SpectogramSprite.hx +++ b/source/SpectogramSprite.hx @@ -1,5 +1,6 @@ package; +import dsp.FFT; import flixel.FlxSprite; import flixel.group.FlxGroup; import flixel.group.FlxSpriteGroup.FlxTypedSpriteGroup; @@ -10,6 +11,7 @@ import flixel.system.FlxSound; import flixel.util.FlxColor; import lime.utils.Int16Array; +using Lambda; using flixel.util.FlxSpriteUtil; class SpectogramSprite extends FlxTypedSpriteGroup @@ -48,14 +50,21 @@ class SpectogramSprite extends FlxTypedSpriteGroup } var setBuffer:Bool = false; - var audioData:Int16Array; + + public var audioData:Int16Array; + var numSamples:Int = 0; override function update(elapsed:Float) { - if (visType == UPDATED) + switch (visType) { - updateVisulizer(); + case UPDATED: + updateVisulizer(); + + case FREQUENCIES: + updateFFT(); + default: } // if visType is static, call updateVisulizer() manually whenever you want to update it! @@ -125,6 +134,78 @@ class SpectogramSprite extends FlxTypedSpriteGroup } } + public function updateFFT() + { + if (daSound != null) + { + var remappedShit:Int = 0; + + checkAndSetBuffer(); + + if (setBuffer) + { + if (daSound.playing) + remappedShit = Std.int(FlxMath.remapToRange(daSound.time, 0, daSound.length, 0, numSamples)); + else + remappedShit = Std.int(FlxMath.remapToRange(Conductor.songPosition, 0, daSound.length, 0, numSamples)); + + var i = remappedShit; + var prevLine:FlxPoint = new FlxPoint(); + + var swagheight:Int = 200; + + var fftSamples:Array = []; + + for (sample in remappedShit...remappedShit + lengthOfShit) + { + var left = audioData[i] / 32767; + var right = audioData[i + 1] / 32767; + + var balanced = (left + right) / 2; + + i += 2; + + // var remappedSample:Float = FlxMath.remapToRange(sample, remappedShit, remappedShit + lengthOfShit, 0, lengthOfShit - 1); + fftSamples.push(balanced); + } + + var freqShit = funnyFFT(fftSamples); + + for (i in 0...group.members.length) + { + // var sampleApprox:Int = Std.int(FlxMath.remapToRange(i, 0, group.members.length, startingSample, startingSample + samplesToGen)); + var remappedFreq:Int = Std.int(FlxMath.remapToRange(i, 0, group.members.length, 0, freqShit.length - 1)); + + group.members[i].x = prevLine.x; + group.members[i].y = prevLine.y; + + var freqIDK:Float = FlxMath.remapToRange(freqShit[remappedFreq], 0, 0.002, 0, 20); + + prevLine.x = (freqIDK * swagheight / 2 + swagheight / 2) + x; + prevLine.y = (i / group.members.length * daHeight) + y; + + // var line = FlxVector.get(prevLine.x - group.members[i].x, prevLine.y - group.members[i].y); + + // group.members[i].setGraphicSize(Std.int(Math.max(line.length, 1)), Std.int(1)); + // group.members[i].angle = line.degrees; + } + /* + for (freq in 0...freqShit.length) + { + var remappedFreq:Float = FlxMath.remapToRange(freq, 0, freqShit.length, 0, lengthOfShit - 1); + + group.members[Std.int(remappedFreq)].x = prevLine.x; + group.members[Std.int(remappedFreq)].y = prevLine.y; + + var freqShit:Float = FlxMath.remapToRange(freqShit[freq], 0, 0.002, 0, 20); + + prevLine.x = (freqShit * swagheight / 2 + swagheight / 2) + x; + prevLine.y = (Math.ceil(remappedFreq) / lengthOfShit * daHeight) + y; + }*/ + } + } + } + public function updateVisulizer():Void { if (daSound != null) @@ -172,10 +253,75 @@ class SpectogramSprite extends FlxTypedSpriteGroup } } } + + function funnyFFT(samples:Array):Array + { + var fs:Float = 44100; // sample rate shit? + + final fftN = 2048; + final halfN = Std.int(fftN / 2); + final overlap = 0.5; + final hop = Std.int(fftN * (1 - overlap)); + + // window function to compensate for overlapping + final a0 = 0.50; // => Hann(ing) window + final window = (n:Int) -> a0 - (1 - a0) * Math.cos(2 * Math.PI * n / fftN); + + // helpers, note that spectrum indexes suppose non-negative frequencies + final binSize = fs / fftN; + final indexToFreq = (k:Int) -> 1.0 * k * binSize; // we need the `1.0` to avoid overflows + + // "melodic" band-pass filter + final minFreq = 32.70; + final maxFreq = 4186.01; + final melodicBandPass = function(k:Int, s:Float) + { + final freq = indexToFreq(k); + final filter = freq > minFreq - binSize && freq < maxFreq + binSize ? 1 : 0; + return s * filter; + }; + + var freqOutput:Array = []; + + var c = 0; // index where each chunk begins + while (c < samples.length) + { + // take a chunk (zero-padded if needed) and apply the window + final chunk = [ + for (n in 0...fftN) + (c + n < samples.length ? samples[c + n] : 0.0) * window(n) + ]; + + // compute positive spectrum with sampling correction and BP filter + final freqs = FFT.rfft(chunk).map(z -> z.scale(1 / fs).magnitude).mapi(melodicBandPass); + + // find spectral peaks and their instantaneous frequencies + for (k => s in freqs) + { + final time = c / fs; + final freq = indexToFreq(k); + final power = s; + if (FlxG.keys.justPressed.N) + { + haxe.Log.trace('${time};${freq};${power}', null); + } + if (freq < 4200) + freqOutput.push(power); + // + } + // haxe.Log.trace("", null); + + // move to next (overlapping) chunk + c += hop; + } + + return freqOutput; + } } enum VISTYPE { STATIC; UPDATED; + FREQUENCIES; } diff --git a/source/dsp/Complex.hx b/source/dsp/Complex.hx new file mode 100644 index 000000000..d99b4490a --- /dev/null +++ b/source/dsp/Complex.hx @@ -0,0 +1,80 @@ +package dsp; + +/** + Complex number representation. +**/ +@:forward(real, imag) @:notNull @:pure +abstract Complex({ + final real:Float; + final imag:Float; +}) +{ + public inline function new(real:Float, imag:Float) + this = {real: real, imag: imag}; + + /** + Makes a Complex number with the given Float as its real part and a zero imag part. + **/ + @:from + public static inline function fromReal(r:Float) + return new Complex(r, 0); + + /** + Complex argument, in radians. + **/ + public var angle(get, never):Float; + + inline function get_angle() + return Math.atan2(this.imag, this.real); + + /** + Complex module. + **/ + public var magnitude(get, never):Float; + + inline function get_magnitude() + return Math.sqrt(this.real * this.real + this.imag * this.imag); + + @:op(A + B) + public inline function add(rhs:Complex):Complex + return new Complex(this.real + rhs.real, this.imag + rhs.imag); + + @:op(A - B) + public inline function sub(rhs:Complex):Complex + return new Complex(this.real - rhs.real, this.imag - rhs.imag); + + @:op(A * B) + public inline function mult(rhs:Complex):Complex + return new Complex(this.real * rhs.real - this.imag * rhs.imag, this.real * rhs.imag + this.imag * rhs.real); + + /** + Returns the complex conjugate, does not modify this object. + **/ + public inline function conj():Complex + return new Complex(this.real, -this.imag); + + /** + Multiplication by a real factor, does not modify this object. + **/ + public inline function scale(k:Float):Complex + return new Complex(this.real * k, this.imag * k); + + public inline function copy():Complex + return new Complex(this.real, this.imag); + + /** + The imaginary unit. + **/ + public static final im = new Complex(0, 1); + + /** + The complex zero. + **/ + public static final zero = new Complex(0, 0); + + /** + Computes the complex exponential `e^(iw)`. + **/ + public static inline function exp(w:Float) + return new Complex(Math.cos(w), Math.sin(w)); +} diff --git a/source/dsp/FFT.hx b/source/dsp/FFT.hx new file mode 100644 index 000000000..b07b9def4 --- /dev/null +++ b/source/dsp/FFT.hx @@ -0,0 +1,154 @@ +package dsp; + +import dsp.Complex; + +// these are only used for testing, down in FFT.main() +using dsp.OffsetArray; +using dsp.Signal; + + +/** + Fast/Finite Fourier Transforms. +**/ +class FFT { + /** + Computes the Discrete Fourier Transform (DFT) of a `Complex` sequence. + + If the input has N data points (N should be a power of 2 or padding will be added) + from a signal sampled at intervals of 1/Fs, the result will be a sequence of N + samples from the Discrete-Time Fourier Transform (DTFT) - which is Fs-periodic - + with a spacing of Fs/N Hz between them and a scaling factor of Fs. + **/ + public static function fft(input:Array) : Array + return do_fft(input, false); + + /** + Like `fft`, but for a real (Float) sequence input. + + Since the input time signal is real, its frequency representation is + Hermitian-symmetric so we only return the positive frequencies. + **/ + public static function rfft(input:Array) : Array { + final s = fft(input.map(Complex.fromReal)); + return s.slice(0, Std.int(s.length / 2) + 1); + } + + /** + Computes the Inverse DFT of a periodic input sequence. + + If the input contains N (a power of 2) DTFT samples, each spaced Fs/N Hz + from each other, the result will consist of N data points as sampled + from a time signal at intervals of 1/Fs with a scaling factor of 1/Fs. + **/ + public static function ifft(input:Array) : Array + return do_fft(input, true); + + // Handles padding and scaling for forwards and inverse FFTs. + private static function do_fft(input:Array, inverse:Bool) : Array { + final n = nextPow2(input.length); + var ts = [for (i in 0...n) if (i < input.length) input[i] else Complex.zero]; + var fs = [for (_ in 0...n) Complex.zero]; + ditfft2(ts, 0, fs, 0, n, 1, inverse); + return inverse ? fs.map(z -> z.scale(1 / n)) : fs; + return fs; + } + + // Radix-2 Decimation-In-Time variant of Cooley–Tukey's FFT, recursive. + private static function ditfft2( + time:Array, t:Int, + freq:Array, f:Int, + n:Int, step:Int, inverse: Bool + ) : Void { + if (n == 1) { + freq[f] = time[t].copy(); + } else { + final halfLen = Std.int(n / 2); + ditfft2(time, t, freq, f, halfLen, step * 2, inverse); + ditfft2(time, t + step, freq, f + halfLen, halfLen, step * 2, inverse); + for (k in 0...halfLen) { + final twiddle = Complex.exp((inverse ? 1 : -1) * 2 * Math.PI * k / n); + final even = freq[f + k].copy(); + final odd = freq[f + k + halfLen].copy(); + freq[f + k] = even + twiddle * odd; + freq[f + k + halfLen] = even - twiddle * odd; + } + } + } + + // Naive O(n^2) DFT, used for testing purposes. + private static function dft(ts:Array, ?inverse:Bool) : Array { + if (inverse == null) inverse = false; + final n = ts.length; + var fs = new Array(); + fs.resize(n); + for (f in 0...n) { + var sum = Complex.zero; + for (t in 0...n) { + sum += ts[t] * Complex.exp((inverse ? 1 : -1) * 2 * Math.PI * f * t / n); + } + fs[f] = inverse ? sum.scale(1 / n) : sum; + } + return fs; + } + + /** + Finds the power of 2 that is equal to or greater than the given natural. + **/ + static function nextPow2(x:Int) : Int { + if (x < 2) return 1; + else if ((x & (x-1)) == 0) return x; + var pow = 2; + x--; + while ((x >>= 1) != 0) pow <<= 1; + return pow; + } + + // testing, but also acts like an example + static function main() { + // sampling and buffer parameters + final Fs = 44100.0; + final N = 512; + final halfN = Std.int(N / 2); + + // build a time signal as a sum of sinusoids + final freqs = [5919.911]; + final ts = [for (n in 0...N) freqs.map(f -> Math.sin(2 * Math.PI * f * n / Fs)).sum()]; + + // get positive spectrum and use its symmetry to reconstruct negative domain + final fs_pos = rfft(ts); + final fs_fft = new OffsetArray( + [for (k in -(halfN - 1) ... 0) fs_pos[-k].conj()].concat(fs_pos), + -(halfN - 1) + ); + + // double-check with naive DFT + final fs_dft = new OffsetArray( + dft(ts.map(Complex.fromReal)).circShift(halfN - 1), + -(halfN - 1) + ); + final fs_err = [for (k in -(halfN - 1) ... halfN) fs_fft[k] - fs_dft[k]]; + final max_fs_err = fs_err.map(z -> z.magnitude).max(); + if (max_fs_err > 1e-6) haxe.Log.trace('FT Error: ${max_fs_err}', null); + // else for (k => s in fs_fft) haxe.Log.trace('${k * Fs / N};${s.scale(1 / Fs).magnitude}', null); + + // find spectral peaks to detect signal frequencies + final freqis = fs_fft.array.map(z -> z.magnitude) + .findPeaks() + .map(k -> (k - (halfN - 1)) * Fs / N) + .filter(f -> f >= 0); + if (freqis.length != freqs.length) { + trace('Found frequencies: ${freqis}'); + } else { + final freqs_err = [for (i in 0...freqs.length) freqis[i] - freqs[i]]; + final max_freqs_err = freqs_err.map(Math.abs).max(); + if (max_freqs_err > Fs / N) trace('Frequency Errors: ${freqs_err}'); + } + + // recover time signal from the frequency domain + final ts_ifft = ifft(fs_fft.array.circShift(-(halfN - 1)).map(z -> z.scale(1 / Fs))); + final ts_err = [for (n in 0...N) ts_ifft[n].scale(Fs).real - ts[n]]; + final max_ts_err = ts_err.map(Math.abs).max(); + if (max_ts_err > 1e-6) haxe.Log.trace('IFT Error: ${max_ts_err}', null); + // else for (n in 0...ts_ifft.length) haxe.Log.trace('${n / Fs};${ts_ifft[n].scale(Fs).real}', null); + } +} diff --git a/source/dsp/OffsetArray.hx b/source/dsp/OffsetArray.hx new file mode 100644 index 000000000..ab5dd2b22 --- /dev/null +++ b/source/dsp/OffsetArray.hx @@ -0,0 +1,78 @@ +package dsp; + +/** + A view into an Array with an indexing offset. + + Usages include 1-indexed sequences or zero-centered buffers with negative indexing. +**/ +@:forward(array, offset) +abstract OffsetArray({ + final array : Array; + final offset : Int; +}) { + public inline function new(array:Array, offset:Int) + this = { array: array, offset: offset }; + + public var length(get,never) : Int; + inline function get_length() + return this.array.length; + + @:arrayAccess + public inline function get(index:Int) : T + return this.array[index - this.offset]; + + @:arrayAccess + public inline function set(index:Int, value:T) : Void + this.array[index - this.offset] = value; + + /** + Iterates through items in their original order while providing the altered indexes as keys. + **/ + public inline function keyValueIterator() : KeyValueIterator + return new OffsetArrayIterator(this.array, this.offset); + + @:from + static inline function fromArray(array:Array) + return new OffsetArray(array, 0); + + @:to + inline function toArray() + return this.array; + + /** + Makes a shifted version of the given `array`, where elements are in the + same order but shifted by `n` positions (to the right if positive and to + the left if negative) in **circular** fashion (no elements discarded). + **/ + public static function circShift(array:Array, n:Int) : Array { + if (n < 0) return circShift(array, array.length + n); + + var shifted = new Array(); + + n = n % array.length; + for (i in array.length - n ... array.length) shifted.push(array[i]); + for (i in 0 ... array.length - n) shifted.push(array[i]); + + return shifted; + } +} + +private class OffsetArrayIterator { + private final array : Array; + private final offset : Int; + private var enumeration : Int; + + public inline function new(array:Array, offset:Int) { + this.array = array; + this.offset = offset; + this.enumeration = 0; + } + + public inline function next() : {key:Int, value:T} { + final i = this.enumeration++; + return { key: i + this.offset, value: this.array[i] }; + } + + public inline function hasNext() : Bool + return this.enumeration < this.array.length; +} diff --git a/source/dsp/Signal.hx b/source/dsp/Signal.hx new file mode 100644 index 000000000..dcc752f16 --- /dev/null +++ b/source/dsp/Signal.hx @@ -0,0 +1,110 @@ +package dsp; + +using Lambda; + + +/** + Signal processing miscellaneous utilities. +**/ +class Signal { + /** + Returns a smoothed version of the input array using a moving average. + **/ + public static function smooth(y:Array, n:Int) : Array { + if (n <= 0) { + return null; + } else if (n == 1) { + return y.copy(); + } else { + var smoothed = new Array(); + smoothed.resize(y.length); + for (i in 0...y.length) { + var m = i + 1 < n ? i : n - 1; + smoothed[i] = sum(y.slice(i - m, i + 1)); + } + return smoothed; + } + } + + /** + Finds indexes of peaks in the order they appear in the input sequence. + + @param threshold Minimal peak height wrt. its neighbours, defaults to 0. + @param minHeight Minimal peak height wrt. the whole input, defaults to global minimum. + **/ + public static function findPeaks( + y:Array, + ?threshold:Float, + ?minHeight:Float + ) : Array { + threshold = threshold == null ? 0.0 : Math.abs(threshold); + minHeight = minHeight == null ? Signal.min(y) : minHeight; + + var peaks = new Array(); + + final dy = [for (i in 1...y.length) y[i] - y[i-1]]; + for (i in 1...dy.length) { + // peak: function growth positive to its left and negative to its right + if ( + dy[i-1] > threshold && dy[i] < -threshold && + y[i] > minHeight + ) { + peaks.push(i); + } + } + + return peaks; + } + + /** + Returns the sum of all the elements of a given array. + + This function tries to minimize floating-point precision errors. + **/ + public static function sum(array:Array) : Float { + // Neumaier's "improved Kahan-Babuska algorithm": + + var sum = 0.0; + var c = 0.0; // running compensation for lost precision + + for (v in array) { + var t = sum + v; + c += Math.abs(sum) >= Math.abs(v) + ? (sum - t) + v // sum is bigger => low-order digits of v are lost + : (v - t) + sum; // v is bigger => low-order digits of sum are lost + sum = t; + } + + return sum + c; // correction only applied at the very end + } + + /** + Returns the average value of an array. + **/ + public static function mean(y:Array) : Float + return sum(y) / y.length; + + /** + Returns the global maximum. + **/ + public static function max(y:Array) : Float + return y.fold(Math.max, y[0]); + + /** + Returns the global maximum's index. + **/ + public static function maxi(y:Array) : Int + return y.foldi((yi, m, i) -> yi > y[m] ? i : m, 0); + + /** + Returns the global minimum. + **/ + public static function min(y:Array) : Float + return y.fold(Math.min, y[0]); + + /** + Returns the global minimum's index. + **/ + public static function mini(y:Array) : Int + return y.foldi((yi, m, i) -> yi < y[m] ? i : m, 0); +}