FFT inside gen?


    Jul 18 2014 | 5:07 pm
    Hi! I'm recently getting into using gen to provide exported code to C programmers. For me this is a very unknown terrain. Now, I know, if they ask me to do a compressor etc, how to do it, but what if something si supposed to happen in the frequency domain? Can I program a fast fourier transform+ inverse *inside* gen? Has anybiody done that? Or should I tell these programmers to send a freq. domain signal to my program, since they probably use some existing code to do the transforms? Any hints tips thoughts very welcome, thanks a lot!

    • Jul 21 2014 | 5:00 pm
      it's fair to assume that gen~ is not designed to operate a fourrier transformation. It works each sample after the other. It's more designed to be a subpatcher inside a pfft~ if you want to do that kind of things, so you may want to ask these programmers to send you a frequency domain signal indeed. this might be tricky of course... but if you don't do that you'll have to create a buffer that fills at every vector, then as you process calculation on what is inside that buffer ; you'll have to fill the next buffer...at the same time... and knowing that the calculation can't go faster that sample rate, it is indeed not very practical....
    • Jul 21 2014 | 5:13 pm
      hi vichug! Thanks for your thoughts!
      and knowing that the calculation can’t go faster that sample rate
      This is not quite true right? I mean I have done for loops with multiple iterations per sample. Maybe I misunderstand you. I'm not sure about all this, but to me it seems at the moment, that not being able to export something that can work with frequency domain signals s a big liability that I'd like to get rid of. Right now it seems to me that there is quite a market for gen exports :) so I'd also invest some time to get some homebrew FFT going.. if it's possible(in a more or less efficient way) that is. Which I'm not a hundred percent sure about.
    • Aug 08 2014 | 7:23 am
      100% possible!!!
      If you use a window of N samples, which is base 2 power (2^ß=N), then you will need ONLY (N*ß)/2 complex multiplications and N*ß complex additions. in "normal" DFT you will need N^2 complex mul's! 2N/ß ratio! You can do 2 signals simultaneously if you only use real values input (like most audio sig). Complexity is O(ε log N) for the classic Cooley-Tukey algorithm.
      For example, take N=2^12=4096 samples window, at 44100 sr (lowest freq bin is ~11Hz): you will need ca. 50K mul's in a time duration of ca. 100ms… thats 12 mul's to calculate between samples…
      so you can even add some overlapping windows and run HD Jitter patch at that time.
      @VICHUG, making a buffer fft code, is exactly what the fft msp objects do in a normal patch. classic FFT algorithm (Cooley-Tukey) are sample by sample DSP, to minimise memory requirements.
      Im just working on making a simple FFT in gen~ for getting high resolution spectra of speech modulation frequencies (less then 12Hz).
      If you wanna work on it together it could be interesting, and we could post it here for gen~erations to come.
      Just pick-up any beginners books on FFT with some basic easy ways to implement it in code (I hope you understand FORTRAN :))
      Keep us updated! Cheers!
    • Dec 05 2015 | 3:24 am
      Did you three ever get anywhere with this? Was looking into it but it looks as if the FFT implementations I've seen require passing arrays as arguments which is impossible in gen.
    • Dec 05 2015 | 1:51 pm
      I'm interested in this project, too. Please let us know of your progress!
    • Dec 05 2015 | 5:34 pm
      I was thinking of using a buffer that was 1024*1024 samples and then the function would just write into it at varying places and return the offset and the length. I'll try it later.
      -Matt
    • Dec 06 2015 | 4:10 pm
      Buffer is a kind of Array. Try it out w stereo buffer to function as complex array.
      I will go back to this. Completely forgot about it, and it can be very nice Object to have + great DSP tutorial...
    • Dec 06 2015 | 6:13 pm
      So... here's what I think it should be but there's a really un-useful error. I have no clue how to fix it. But basically this is an exact replica of http://introcs.cs.princeton.edu/java/97data/FFT.java.html
      The only difference is instead of creating and passing arrays it uses it gets passed the buffer location of the beginning of the array, the length of the array and the next "free" piece of real estate in the buffer.
      I have no clue if it should work in gen but it seems pretty spot on in terms of math.
    • Apr 26 2020 | 8:30 pm
      I'm pretty sure you're getting that error because gen cannot perform recursion.
    • May 07 2020 | 7:16 pm
      While function recursion in that exact form is not feasible in gen~, most recursive algorithms can be restated as iterative loops. An example that is structurally similar to the FFT is the gen~ Haar wavelet analysis/resynthesis example included in Max (search for "gen~.haar.wavelet.maxpat"). This is an essentially recursive algorithm (progressively dividing a problem in halves until the lowest level is reached) but the recursion is expressed there through a while() loop. I think it should be possible to express the FFT in a similar way; after all, Haar is to square as Fourier is to sine. If I have chance I may try and take a shot at this. There's no way that this will have respectable performance relative to a proper C++ FFT library, but it might be useful for educational purposes.
    • May 07 2020 | 7:23 pm
      Graham. I was trying to figure out the For() form of it and... my brain broke. You’ve got wavelets happening in max? I’ll have to check that out.
      can you use arbitrary kernels? I‘m really trying to put together a spike encoder/decoder for max. Would be super useful.
    • May 07 2020 | 7:42 pm
      No no no, it's just a super-basic Haar wavelet demo -- arbitrary kernels would be quite a bit harder. It's also more useful as an educational tool; but *can* do some dirty crunchy multi-band things to your audio if you want it. It's part of Max already so open it up and take a look :-)
      Going to give myself an hour to try the FFT thing.
    • May 07 2020 | 9:00 pm
      most recursive algorithms can be restated as iterative loops
      I would go even further: all recursive functions can be expressed in an iterative way.
    • May 07 2020 | 9:05 pm
      I would go even further: all recursive functions can be expressed in an iterative way.
      The trick is doing it with out banjaxing your brain. I just found this if it helps. http://blog.moertel.com/posts/2013-05-11-recursive-to-iterative.html
      Going to give myself an hour to try the FFT thing.
      Morgan Freeman: Some say he worked on it for twenty five weeks...
    • May 07 2020 | 9:16 pm
      One hour flying blind got me partway through before my brain started to melt -- got the even/odd dissection working but got a bit lost with the butterfly part. A quick look at wikipedia however confirms that this should be possible -- see the "iterative-fft" pseudocode here: https://en.wikipedia.org/wiki/Cooley%E2%80%93Tukey_FFT_algorithm#Data_reordering,_bit_reversal,_and_in-place_algorithms
      Leaving that link here for me to come back and try again when I have more coffee and free time, unless someone else wants to try it first!
    • May 09 2020 | 12:08 am
      Despite pitching at windmills I feel like this is getting somewhere near. The result kind of looks like a spectrum, with more energy in the lower bins, and the expected double-repeat of the spectrum. But... somehow not quite right. I mean, as far as I remember, a sine that is an integer harmonic of the fft frame (as generated in a buffer~ via 'fill sin N') really ought to give a spectrum with a very nice single spike in the corresponding harmonic bin and zero elsewhere, shouldn't it?
      The bit reversal stage seems to be correct, and of course the real, imag => mag, phase conversion is no doubt solid, so I guess something must be off with the inner loops of the iterative FFT. Reaching out to the mighty grey matter among you to see the error, as I'm fried now.
    • May 11 2020 | 12:23 am
      Got it working -- both forward and inverse!
      Next steps: windowing & overlap...
    • May 11 2020 | 5:47 am
      That's HUGE! 🙌 Praise to Graham 👏👏 (I don't understand any of it, but going to keep poking through it til I do) THANK YOU SO MUCH!
    • May 11 2020 | 10:52 am
      Magic.
      soooooo coool!!!
    • May 11 2020 | 1:26 pm
      Praise to Graham, indeed! Beautiful!
    • May 13 2020 | 8:33 pm
      That is sweet! I thought about trying this a few years ago but gave up rather quickly. Looking forward to looking at it a bit more later!
    • May 13 2020 | 9:07 pm
      @cycling: include this in the examples!
    • May 13 2020 | 10:14 pm
      ☝️+1 👍 (also, add +1 buttons to the forums-interface near the bottom-right 'share' button on posts here 😁)
    • May 14 2020 | 6:33 pm
      Well thanks for the support :-)
      Made more progress, now have a windowed 2x overlap version working (2nd overlap uses additional channels of the respective buffers). Demo below includes a very basic spectral gate lowshelf (zeroes out magnitude of bins above a certain threshold).
      A couple of things I'm not yet happy with:
      • I can hear some weird aliasing noise when the cutoff bin is very low, not sure where this is coming from.
      • I think I'm doing 2x the processing that I really need to; the polar spectrum is symmetric (with the phase inverted in the 2nd half), which suggests that I could just compute half of it and resynthesize from that.
      • Some hard-coded numbers, a bunch of things that could probably be cached, etc., but it is getting there.
      • The structure of the code could probably also be improved so it is more obvious/easy to start inserting spectral processing.
      Any ideas/suggestions very welcome!!
      Graham
      PS. I really really never thought I was going to implement the FFT in my life. Weird.
    • May 15 2020 | 7:51 pm
      Amazing stuff Graham! thanks