diff --git a/Makefile b/Makefile index c42a863..6046917 100644 --- a/Makefile +++ b/Makefile @@ -16,7 +16,7 @@ EXTRALIBS:=-lusb-1.0 colorchord : os_generic.o main.o dft.o decompose.o filter.o color.o sort.o notefinder.o util.o outdrivers.o $(RAWDRAW) $(SOUND) $(OUTS) parameters.o chash.o gcc -o $@ $^ $(CFLAGS) $(LDLIBS) $(EXTRALIBS) $(RAWDRAWLIBS) -embeddedcc : os_generic.c embeddedcc.c dft.c +embeddedcc : os_generic.c embeddedcc.c dft.c embeddednf.c gcc -o $@ $^ $(CFLAGS) -DCCEMBEDDED $(LDFLAGS) $(EXTRALIBS) $(RAWDRAWLIBS) runembedded : embeddedcc diff --git a/dft.c b/dft.c index bfa3d02..75aab1c 100644 --- a/dft.c +++ b/dft.c @@ -559,13 +559,13 @@ void HandleProgressiveIntSkippy( int8_t sample1 ) //Add TS and TC to the datspace stuff. (24 instructions) tmp1 = (*ds); //Read out, sin component. 4 Accurate. // tmp1 -= tmp1>>4; //Subtract from the MSB (with carry) 2 -> 6 AS/IS: 7+7 = 14 - tmp1 += ts>>4; //Add MSBs with carry 2 -> 6 AS/IS: 6 + tmp1 += ts>>SHIFT_ADD_DETAIL; //Add MSBs with carry 2 -> 6 AS/IS: 6 *(ds++) = tmp1; //Store values back 4 tmp1 = *ds; //Read out, sin component. 4 // tmp1 -= tmp1>>4; //Subtract from the MSB (with carry) 2 -> 6 AS/IS: 7+7 = 14 - tmp1 += tc>>4; //Add MSBs with carry 2 -> 6 AS/IS: 6 + tmp1 += tc>>SHIFT_ADD_DETAIL; //Add MSBs with carry 2 -> 6 AS/IS: 6 *ds++ = tmp1; //Store values back 4 diff --git a/dft.h b/dft.h index 1589c22..531d0cc 100644 --- a/dft.h +++ b/dft.h @@ -56,6 +56,13 @@ void Push8BitIntegerSkippy( int8_t dat ); //Call this to push on new frames of s #define FIXBINS (FIXBPERO*OCTAVES) #define BINCYCLE (1< - +#include "embeddednf.h" #include "dft.h" -#define DFREQ 8000 -#define BASE_FREQ 55.0 - -const float bf_table[24] = { - 1.000000, 1.029302, 1.059463, 1.090508, 1.122462, 1.155353, - 1.189207, 1.224054, 1.259921, 1.296840, 1.334840, 1.373954, - 1.414214, 1.455653, 1.498307, 1.542211, 1.587401, 1.633915, - 1.681793, 1.731073, 1.781797, 1.834008, 1.887749, 1.943064 }; - -#define BUILD_BUG_ON(condition) ((void)sizeof(char[1 - 2*!!(condition)])) - -/* The above table was generated using the following code: - -#include -#include - -int main() -{ - int i; - #define FIXBPERO 24 - printf( "const float bf_table[%d] = {", FIXBPERO ); - for( i = 0; i < FIXBPERO; i++ ) - { - if( ( i % 6 ) == 0 ) - printf( "\n\t" ); - printf( "%f, ", pow( 2, (float)i / (float)FIXBPERO ) ); - } - printf( "};\n" ); - return 0; -} -*/ - -void UpdateFreqs() -{ - uint16_t fbins[FIXBPERO]; - int i; - - BUILD_BUG_ON( sizeof(bf_table) != FIXBPERO*4 ); - - for( i = 0; i < FIXBPERO; i++ ) - { - float frq = ( bf_table[i] * BASE_FREQ ); - fbins[i] = ( 65536.0 ) / ( DFREQ ) * frq * 16; - } - - UpdateBinsForProgressiveIntegerSkippyInt( fbins ); -} - -void Init() -{ - //Step 1: Initialize the Integer DFT. - SetupDFTProgressiveIntegerSkippy(); - - //Step 2: Set up the frequency list. - UpdateFreqs(); -} - -void HandleFrameInfo() -{ - uint16_t folded_bins[FIXBPERO]; - - int i, j, k = 0; - - for( i = 0; i < FIXBPERO; i++ ) - folded_bins[i] = 0; - - for( j = 0; j < OCTAVES; j++ ) - { - for( i = 0; i < FIXBPERO; i++ ) - { - folded_bins[i] += embeddedbins[k++]; - } - } - - - - //XXX TODO Taper the first and last octaves. -// for( i = 0; i < freqbins; i++ ) -// { -// nf->outbins[i] *= (i+1.0)/nf->freqbins; -// } -// for( i = 0; i < freqbins; i++ ) -// { -// nf->outbins[freqs-i-1] *= (i+1.0)/nf->freqbins; -// } - - //We now have the system folded into one - - for( i = 0; i < FIXBPERO; i++ ) - { - printf( "%5d ", folded_bins[i] ); - } - printf( "\n" ); -} int main() { diff --git a/embeddednf.c b/embeddednf.c new file mode 100644 index 0000000..07e7930 --- /dev/null +++ b/embeddednf.c @@ -0,0 +1,342 @@ +#include "embeddednf.h" +#include + +uint16_t folded_bins[FIXBPERO]; +uint16_t fuzzed_bins[FIXBINS]; +uint8_t note_peak_freqs[MAXNOTES]; +uint16_t note_peak_amps[MAXNOTES]; +uint16_t note_peak_amps2[MAXNOTES]; + + +static const float bf_table[24] = { + 1.000000, 1.029302, 1.059463, 1.090508, 1.122462, 1.155353, + 1.189207, 1.224054, 1.259921, 1.296840, 1.334840, 1.373954, + 1.414214, 1.455653, 1.498307, 1.542211, 1.587401, 1.633915, + 1.681793, 1.731073, 1.781797, 1.834008, 1.887749, 1.943064 }; + +/* The above table was generated using the following code: + +#include +#include + +int main() +{ + int i; + #define FIXBPERO 24 + printf( "const float bf_table[%d] = {", FIXBPERO ); + for( i = 0; i < FIXBPERO; i++ ) + { + if( ( i % 6 ) == 0 ) + printf( "\n\t" ); + printf( "%f, ", pow( 2, (float)i / (float)FIXBPERO ) ); + } + printf( "};\n" ); + return 0; +} +*/ + + +#define BUILD_BUG_ON(condition) ((void)sizeof(char[1 - 2*!!(condition)])) + + +void UpdateFreqs() +{ + uint16_t fbins[FIXBPERO]; + int i; + + BUILD_BUG_ON( sizeof(bf_table) != FIXBPERO*4 ); + + //Warning: This does floating point. Avoid doing this frequently. If you + //absolutely cannot have floating point on your system, you may precompute + //this and store it as a table. It does preclude you from changing + //BASE_FREQ in runtime. + + for( i = 0; i < FIXBPERO; i++ ) + { + float frq = ( bf_table[i] * BASE_FREQ ); + fbins[i] = ( 65536.0 ) / ( DFREQ ) * frq * 16; + } + + UpdateBinsForProgressiveIntegerSkippyInt( fbins ); +} + +void Init() +{ + int i; + //Set up and initialize arrays. + for( i = 0; i < MAXNOTES; i++ ) + { + note_peak_freqs[i] = 255; + note_peak_amps[i] = 0; + note_peak_amps2[i] = 0; + } + + for( i = 0; i < FIXBPERO; i++ ) + { + folded_bins[i] = 0; + } + + for( i = 0; i < FIXBINS; i++ ) + { + fuzzed_bins[i] = 0; + } + + //Step 1: Initialize the Integer DFT. + SetupDFTProgressiveIntegerSkippy(); + + //Step 2: Set up the frequency list. You could do this multiple times + //if you want to change the loadout of the frequencies. + UpdateFreqs(); +} + +void HandleFrameInfo() +{ + int i, j, k; + + //Copy out the bins from the DFT to our fuzzed bins. + for( i = 0; i < FIXBINS; i++ ) + { + fuzzed_bins[i] = (fuzzed_bins[i] + (embeddedbins[i]>>FUZZ_IIR_BITS) - + (fuzzed_bins[i]>>FUZZ_IIR_BITS)); + } + + //Taper first octave + for( i = 0; i < FIXBPERO; i++ ) + { + uint32_t taperamt = (65536 / FIXBPERO) * i; + fuzzed_bins[i] = (taperamt * fuzzed_bins[i]) >> 16; + } + + //Taper last octave + for( i = 0; i < FIXBPERO; i++ ) + { + int newi = FIXBINS - i - 1; + uint32_t taperamt = (65536 / FIXBPERO) * i; + fuzzed_bins[newi] = (taperamt * fuzzed_bins[newi]) >> 16; + } + + //Fold the bins from fuzzedbins into one octave. + for( i = 0; i < FIXBPERO; i++ ) + folded_bins[i] = 0; + + k = 0; + for( j = 0; j < OCTAVES; j++ ) + { + for( i = 0; i < FIXBPERO; i++ ) + { + folded_bins[i] += fuzzed_bins[k++]; + } + } + + //Now, we must blur the folded bins to get a good result. + //Sometimes you may notice every other bin being out-of + //line, and this fixes that. We may consider running this + //more than once, but in my experience, once is enough. + { + //Extra scoping because this is a large on-stack buffer. + uint16_t folded_out[FIXBPERO]; + uint8_t adjLeft = FIXBPERO-1; + uint8_t adjRight = 1; + for( i = 0; i < FIXBPERO; i++ ) + { + uint16_t lbin = folded_bins[adjLeft]>>2; + uint16_t rbin = folded_bins[adjRight]>>2; + uint16_t tbin = folded_bins[i]>>1; + folded_out[i] = lbin + rbin + tbin; + + //We do this funny dance to avoid a modulus operation. On some + //processors, a modulus operation is slow. This is cheap. + adjLeft++; if( adjLeft == FIXBPERO ) adjLeft = 0; + adjRight++; if( adjRight == FIXBPERO ) adjRight = 0; + } + + for( i = 0; i < FIXBPERO; i++ ) + { + folded_bins[i] = folded_out[i]; + } + } + + //Next, we have to find the peaks, this is what "decompose" does in our + //normal tool. As a warning, it expects that the values in foolded_bins + //do NOT exceed 32767. + { + uint8_t adjLeft = FIXBPERO-1; + uint8_t adjRight = 1; + for( i = 0; i < FIXBPERO; i++ ) + { + int16_t prev = folded_bins[adjLeft]; + int16_t next = folded_bins[adjRight]; + int16_t this = folded_bins[i]; + uint8_t thisfreq = i< this || next > this ) continue; + if( prev == this && next == this ) continue; + + //i is at a peak... + int32_t totaldiff = (( this - prev ) + ( this - next )); + int32_t porpdiffP = ((this-prev)<<16) / totaldiff; //close to 0 = + //closer to this side, 32768 = in the middle, 65535 away. + int32_t porpdiffN = ((this-next)<<16) / totaldiff; + + if( porpdiffP < porpdiffN ) + { + //Closer to prev. + offset = -(32768 - porpdiffP); + } + else + { + //Closer to next + offset = (32768 - porpdiffN); + } + + //Need to round. That's what that extra +(15.. is in the center. + thisfreq += (offset+(1<<(15-SEMIBITSPERBIN)))>>(16-SEMIBITSPERBIN); + + //In the event we went 'below zero' need to wrap to the top. + if( thisfreq > 255-(1< ((1<<(SEMIBITSPERBIN-1))*FIXBPERO) ) + { + distance = ((1<<(SEMIBITSPERBIN))*FIXBPERO) - distance; + } + + if( distance < closest_note_distance ) + { + closest_note_id = j; + closest_note_distance = distance; + } + } + + int8_t marked_note = -1; + + if( closest_note_distance <= MAX_JUMP_DISTANCE ) + { + //We found the note we need to augment. + //XXX: TODO: Should we be IIRing this? + + note_peak_freqs[closest_note_id] = thisfreq; + marked_note = closest_note_id; + } + //The note was not found. + else if( lowest_found_free_note != -1 ) + { + note_peak_freqs[lowest_found_free_note] = thisfreq; + marked_note = lowest_found_free_note; + } + + if( marked_note != -1 ) + { + if( note_peak_amps[marked_note] <= this ) + note_peak_amps[marked_note] = this; + if( note_peak_amps2[marked_note] <= this ) + note_peak_amps2[marked_note] = this; + } + } + } + + //Now we need to handle combining notes. + for( i = 0; i < MAXNOTES; i++ ) + for( j = 0; j < i; j++ ) + { + //We'd be combining nf2 (j) into nf1 (i) if they're close enough. + uint8_t nf1 = note_peak_freqs[i]; + uint8_t nf2 = note_peak_freqs[j]; + int16_t distance = nf1 - nf2; + + if( nf1 == 255 || nf2 == 255 ) continue; + + if( distance < 0 ) distance = -distance; + + if( distance > ((1<<(SEMIBITSPERBIN-1))*FIXBPERO) ) + { + distance = ((1<<(SEMIBITSPERBIN))*FIXBPERO) - distance; + } + + if( distance > MAX_JUMP_DISTANCE * 2 ) + { + continue; + } + + //We need to combine the notes. We need to move the new note freq + //towards the stronger of the two notes. + int16_t amp1 = note_peak_amps[i]; + int16_t amp2 = note_peak_amps[j]; + + //0 to 32768 porportional to how much of amp1 we want. + uint32_t porp = (amp1<<15) / (amp1+amp2); + uint16_t newnote = (nf1 * porp + nf2 * (32768-porp))>>15; + + note_peak_amps[i] = amp1 + amp2; + note_peak_amps[j] = 0; + note_peak_freqs[i] = newnote; + note_peak_freqs[j] = 255; + note_peak_amps2[i] += note_peak_amps2[j]; + } + + + for( i = 0; i < MAXNOTES; i++ ) + { + if( note_peak_freqs[i] == 255 ) continue; + + note_peak_amps[i] -= note_peak_amps[i]>>AMP_1_NERFING_BITS; + note_peak_amps2[i] -= note_peak_amps2[i]>>AMP_2_NERFING_BITS; + + if( note_peak_amps[i] < MINIMUM_AMP_FOR_NOTE_TO_DISAPPEAR ) + { + note_peak_freqs[i] = 255; + note_peak_amps[i] = 0; + note_peak_amps2[i] = 0; + } + } + + //We now have notes!!! +/* + for( i = 0; i < MAXNOTES; i++ ) + { + if( note_peak_freqs[i] == 255 ) continue; + printf( "(%3d %4d %4d) ", note_peak_freqs[i], note_peak_amps[i], note_peak_amps2[i] ); + } + printf( "\n") ; +*/ + +/* + for( i = 0; i < FIXBPERO; i++ ) + { + printf( "%5d ", folded_bins[i] ); + } + printf( "\n" );*/ +} + + diff --git a/embeddednf.h b/embeddednf.h new file mode 100644 index 0000000..a28e930 --- /dev/null +++ b/embeddednf.h @@ -0,0 +1,53 @@ +#ifndef _EMBEDDEDNF_H +#define _EMBEDDEDNF_H + +#include "dft.h" + +#define DFREQ 8000 +#define BASE_FREQ 55.0 // You may make this a float. + +//The higher the number the slackier your FFT will be come. +#define FUZZ_IIR_BITS 1 + +//Notes are the individually identifiable notes we receive from the sound. +//We track up to this many at one time. Just because a note may appear to +//vaporize in one frame doesn't mean it is annihilated immediately. +#define MAXNOTES 10 + +//Determines bit shifts for where notes lie. We represent notes with an uint8_t +//We have to define all of the possible locations on the note line in this. +//note_frequency = 0..((1<