simple LPF and HPF (Low/High Pass Filter) without FFT

[Froggerprogger]

Here's an example how to make a simple LPF/HPF (and a Gain) that is not using the CPU-Power-intensive way via FFT/DSP/IFFT

princip:
The actual audiobufferdata is sent to a procedure by a FMOD-DSP-callback.
The procedure changes the data.

High value-differences from one samplevalue to another means a 'loud' frequency of Nyquist-FQ (22050 Hz at 44100 Hz) When we lower this difference logarithmic by a constant factor it should work similar to a LPF with center-frequency [nyquist-frequency] and -6dB damp/oktave.
So if this LPF is set to -18dB, then e.g. 11050Hz is still damped by 12dB, 5525 Hz by 6dB and <= 2762,5 by 0dB

The example needs the actuel fmod-import.
(at the moment you'll get it here: http://www.fmod.de/files/FMOD371_DLL_WRAPPER_PB.zip )

You can find source and an EXE here, too: http://fmod.2mal2mal.de

Code: Select all

; If this program crashes, try setting the debugger off!
;
; /////////////////// HighPassFilter and LowPassFilter (without FFT) and Gain version 1.0
; /////////////////// realized as DSP-callbacks using fmod's DSP-chain
; /////////////////// This program needs the PB-FMOD-Import for the direct function-calls
; /////////////////// and the fmod.dll in version 3.71
; /////////////////// 
; /////////////////// (v1.2 by Froggerprogger, 20.01.04)

; These algorithms aren't 'real' High/Low Pass Filters. Their influence on the spectrum is
; different than the classic filters'.
; Here only the delta-values of two continguous samplevalues are 'strechted' or 'pressed'.
; Higher frequencies have more often rapid valuechanges, so they are more affected by this
; method than lower frequencies. But - in opposite to a 'real' filter - all frequencies are
; >always< affected.
; Anyway, this filters manipulate the sound in a similar way and can be seen as an own kind of
; HPF/LPF-filter - and they are VERY CPU-friendly. 


;- PROCEDURE declaration
Declare.l myFMOD_HPF(*originalbuffer.l, *outputbuffer.l, dspBufferlengthSamples.l, *hpfRatio.l)
Declare.l myFMOD_LPF(*originalbuffer.l, *outputbuffer.l, dspBufferlengthSamples.l, *lpfRatio.l)
Declare.l myFMOD_GAIN(*originalbuffer.l, *outputbuffer.l, dspBufferlengthSamples.l, *Gain_dB.l)
Declare.f dB2Dec(dB.f)
Declare.f Dec2dB(dec.f)

;- CONST declaration
;diverse
#Window_Title.s = "HPF + LPF (without FFT) and Gain v1.2 by Froggerprogger"

; maths consts
#e.f = 2.7182818

; PB-Memory-IDs
#MemID_lastHpfSampleBuffer = 1
#MemID_lastLpfSampleBuffer = 2

; Gadget-IDs
#Gadget_HpfTrackbar = 1
#Gadget_HpfText = 2
#Gadget_HpfOnOff = 3
#Gadget_LpfTrackbar = 4
#Gadget_LpfOnOff = 5
#Gadget_LpfText = 6
#Gadget_GainTrackbar = 7
#Gadget_GainText = 8
#Gadget_DistortionOnText = 9

; FMOD-consts
#myFSOUND_StreamInit = #FSOUND_LOOP_NORMAL | #FSOUND_16BITS | #FSOUND_STEREO | #FSOUND_SIGNED

;- GLOABL VAR declaration
Global soundfile.s

Global dspNumChannels.l     : dspNumChannels = 2        ; fmod's dsp-chain is always in stereo, so this one is nearly redundant
Global dspBitrate.l         : dspBitrate = 16           ; fmod's dsp-chain has always 16 significant bits per sample, so same as above
Global dspBitalign.l        : dspBitalign = Int(Round(dspBitrate/8,1)) ; (independant of dspNumChannels!)

Global masterSamplerate.l   : masterSamplerate = 44100  ; the samplerate for FMOD_Init()
Global masterNumChannels.l  : masterNumChannels = 2     ; the number of channels for FMOD_Init()

Global *dsp_hpf.l
Global setHpfRatio.f        : setHpfRatio = -1.0
Global *dsp_lpf.l
Global setLpfRatio.f        : setLpfRatio = -1.0
Global *dsp_Gain.l
Global setGainRatio.f       : setGainRatio = -1.0

Global mainresume.l         : mainresume = #True
;-
;-
;-                                         MAIN PROGRAM
;-
;-

;_____________ check fmod.dll's version

If FSOUND_GetVersion() <> 3.71
  If MessageRequester("","This program is for fmod.dll in version 3.71"+Chr(13)+"You have the fmod.dll in version "+StrF(FSOUND_GetVersion(), 2)+Chr(13)+"Continue anyway ?", #MB_ICONWARNING | #MB_YESNO) = #IDNO
    End
  EndIf
EndIf


;_____________ open a window and put some gadgets onto it
OpenWindow(1,0,0,670,210,#PB_Window_SystemMenu | #PB_Window_ScreenCentered, #Window_Title)
CreateGadgetList(WindowID(1))
TrackBarGadget(#Gadget_HpfTrackbar,40,10,580,25,0,960,0)
    SetGadgetState(#Gadget_HpfTrackbar,Int(setHpfRatio-1.0))
    DisableGadget(#Gadget_HpfTrackbar,1)
CheckBoxGadget(#Gadget_HpfOnOff, 10, 10, 25, 25, "")
TextGadget(#Gadget_HpfText,40, 40, 160, 25, "HPF: -0.00 dB")
    DisableGadget(#Gadget_HpfText,1)
TrackBarGadget(#Gadget_LpfTrackbar,40,80,580,25,0,960,0)
    SetGadgetState(#Gadget_LpfTrackbar,Int(setLpfRatio-1.0))
    DisableGadget(#Gadget_LpfTrackbar,1)
CheckBoxGadget(#Gadget_LpfOnOff, 10, 80, 25, 25, "")
TextGadget(#Gadget_LpfText,40,110,160,25,"LPF: -0.00 dB")
    DisableGadget(#Gadget_LpfText,1)
TrackBarGadget(#Gadget_GainTrackbar,40,150,580,25,0,1000,0)
    SetGadgetState(#Gadget_GainTrackbar,500)
TextGadget(#Gadget_GainText,40,180,160,25,"Gain: +0.00 dB")    
TextGadget(#Gadget_DistortionOnText,520,180,100,25,"",#PB_Text_Right)    

;_____________ load an audio-file
soundfile = OpenFileRequester("","","Audiofiles (*.mp3, *.wav, *.mp2, *.ogg, *.raw)|*.mp3;*.wav;*.mp2;*.ogg;*.raw|All files|*.*",0,0) 

;_____________ initialize FMOD
FSOUND_SetBufferSize(125)
If FSOUND_Init(masterSamplerate, masterNumChannels, 0)
  
  ;_____________ create and activate the DSP-Procedures
  *dsp_lpf = FSOUND_DSP_Create(@myFMOD_LPF(),880, @setLpfRatio)
  *dsp_hpf = FSOUND_DSP_Create(@myFMOD_HPF(),881, @setHpfRatio)
  *dsp_Gain = FSOUND_DSP_Create(@myFMOD_GAIN(),882, @setGainRatio)
  FSOUND_DSP_SetActive(*dsp_lpf, #True)
  FSOUND_DSP_SetActive(*dsp_hpf, #True)
  FSOUND_DSP_SetActive(*dsp_Gain, #True)

  ;_____________ start playing the soundfile
  *hStream=FSOUND_Stream_Open(soundfile, #myFSOUND_StreamInit, 0, 0)
  If *hStream = #False : MessageRequester("","Couldn't load audio-file "+soundfile,0) : End : EndIf
  FSOUND_Stream_Play(1, *hStream)
Else
  MessageRequester("", "Couldn't initialize fmod.dll", #MB_ICONERROR )
  End
EndIf

;_____________ MAIN LOOP
While mainresume
    event = WaitWindowEvent()
    Select event
        Case #PB_Event_Gadget
          Select EventGadgetID()
            Case #Gadget_HpfTrackbar
                setHpfRatio = GetGadgetState(#Gadget_HpfTrackbar)/10
                SetGadgetText(#Gadget_HpfText, "HPF: -"+StrF(setHpfRatio,2)+" dB")
            
            Case #Gadget_LpfTrackbar
              setLpfRatio = GetGadgetState(#Gadget_LpfTrackbar)/10
              SetGadgetText(#Gadget_LpfText, "LPF: -"+StrF(setLpfRatio,2)+" dB")
              
            Case #Gadget_GainTrackbar
              setGainRatio = (GetGadgetState(#Gadget_GainTrackbar)-500)/10
              If setGainRatio > 0
                SetGadgetText(#Gadget_GainText, "Gain: +"+StrF(setGainRatio,2)+" dB")
              Else
                SetGadgetText(#Gadget_GainText, "Gain: "+StrF(setGainRatio,2)+" dB")
              EndIf
              
            Case #Gadget_HpfOnOff
              DisableGadget(#Gadget_HpfTrackbar, 1!GetGadgetState(#Gadget_HpfOnOff))
              DisableGadget(#Gadget_HpfText, 1!GetGadgetState(#Gadget_HpfOnOff))
              If GetGadgetState(#Gadget_HpfOnOff)
                setHpfRatio = GetGadgetState(#Gadget_HpfTrackbar)/10
                SetGadgetText(#Gadget_HpfText, "HPF: -"+StrF(setHpfRatio,2)+" dB")
              Else
                setHpfRatio = -1.0
              EndIf
            
            Case #Gadget_LpfOnOff
              DisableGadget(#Gadget_LpfTrackbar, 1!GetGadgetState(#Gadget_LpfOnOff))
              DisableGadget(#Gadget_LpfText, 1!GetGadgetState(#Gadget_LpfOnOff))
              If GetGadgetState(#Gadget_LpfOnOff)
                setLpfRatio = GetGadgetState(#Gadget_LpfTrackbar)/10
                SetGadgetText(#Gadget_LpfText, "LPF: -"+StrF(setLpfRatio,2)+" dB")
              Else
                setLpfRatio = -1.0
              EndIf
          EndSelect
        
        Case #PB_Event_CloseWindow
          ;_____________ shut down FMOD
          FSOUND_DSP_SetActive(*dsp_lpf, #False)
          FSOUND_DSP_SetActive(*dsp_hpf, #False)
          FSOUND_DSP_SetActive(*dsp_Gain, #False)
          FSOUND_DSP_Free(*dsp_lpf)
          FSOUND_DSP_Free(*dsp_hpf)
          FSOUND_DSP_Free(*dsp_Gain)
          FSOUND_Stream_Stop(*hStream)
          FSOUND_Close()
          mainresume = 0
          
    EndSelect

Wend
End



;-
;-
;-                                          PROCEDURES
;-
;-



;-
;- ___ myFMOD_HPF
;-
Procedure.l myFMOD_HPF(*originalbuffer.l, *outputbuffer.l, dspBufferlengthSamples.l, *hpfRatio.l)
  Shared *lasthpfSampleBuffer.l    : ; stores the previous samplevalues
  
  Protected actCHSample.l          : actCHSample = 0 ; CH means Channel, so it's value alters between all mutlichannel-sample-values e.g. L | R for stereo
  Protected CHSwitch.l             : CHSwitch = 0    ; a variable for channelIDs, [left (0), right (1) for stereo]
  Protected dspBufferlengthBytes.l : dspBufferlengthBytes = dspNumChannels * dspBitalign * dspBufferlengthSamples ; bufferlength in bytes, not samples
  Protected hpfRatio.f             
  Protected hpfFactor.f            
  Protected diff.f ; contains the difference of two contiguous samplevalues
  Protected tempL.l
  Protected maxval.l               : maxval = Pow(2, dspBitalign * 8 - 1) - 1
  Protected minval.l               : minval = -1 * maxval - 1
  
  If PeekF(*hpfRatio) <> -1
    hpfRatio = dspBitalign * 48 - PeekF(*hpfRatio)  ; 2 * 8(=16 Bits) * 6 (dB)
  Else
    hpfRatio = -1
  EndIf
  hpfFactor = 1/dB2Dec(hpfRatio)
  
  If *lasthpfSampleBuffer = 0 ; allocate new memory just the first time
    *lasthpfSampleBuffer = AllocateMemory(#MemID_lastHpfSampleBuffer, dspNumChannels * 4)
  EndIf
  
  If dspBitalign = 2 And hpfRatio <> -1 ; 16 Bit signal processing only here, because of PeekW in the following  
    
    While actCHSample < dspBufferlengthBytes-1 ; process the whole buffer
      hpfCHSample = PeekW(*outputbuffer+actCHSample) ; read in the actual 16Bit-sample-value
      tempL = hpfCHSample 
      hpfLastCHSample = PeekL(*lasthpfSampleBuffer+4*CHSwitch) ; read in the actual previous 16Bit-sample-value
      
      diff = hpfCHSample - hpfLastCHSample ; compare it with the previous one
      
      diff = diff * hpfFactor ; lower the difference between them by hpfFactor
      hpfCHSample = hpfLastCHSample + diff ; calculate the hpf-corrected sample-value
      
      PokeL(*lasthpfSampleBuffer+4*CHSwitch, hpfCHSample) ; the actual one will be the previous one for left channel
      
      hpfCHSample = tempL - hpfLastCHSample
      
      If hpfCHSample > maxval ; actually this could happen here
        hpfCHSample = maxval
      ElseIf hpfCHSample < minval
        hpfCHSample = minval
      EndIf
      
      PokeW(*outputbuffer + actCHSample, hpfCHSample) ; write it into the output-buffer
      
      CHSwitch + 1 : If CHSwitch = dspNumChannels : CHSwitch = 0 : EndIf ; circulate CHSwitch between all Channels
      actCHSample + dspBitalign ; increases the byteposition by bitalign (2 at 16 bit-audio)
    Wend
  EndIf  
  ProcedureReturn *outputbuffer ; overgive a pointer to our outputbuffer to the next dsp-station (will be *outputbuffer there)
EndProcedure
; needs during FSOUND_DSP_Create() for the *hpfRatio-parameter a pointer to a float-variable,
; whose value is in the range [0.0; 6 * Fmod-Bitrate (so normally 6*16 = 96)] and can be
; live updated during playback.
; It represents the factor for reducing the low freqencies in dB.
; Setting the float-variable to -1 sets the DSP to bypass.




;-
;- ___ myFMOD_LPF
;-
Procedure.l myFMOD_LPF(*originalbuffer.l, *outputbuffer.l, dspBufferlengthSamples.l, *lpfRatio.l)
  Shared *lastlpfSampleBuffer.l    : ; stores the previous samplevalues
  
  Protected actCHSample.l          : actCHSample = 0 ; CH means Channel, so it's value alters between all mutlichannel-sample-values e.g. L | R for stereo
  Protected CHSwitch.l             : CHSwitch = 0    ; a variable for channelIDs, [left (0), right (1) for stereo]
  Protected dspBufferlengthBytes.l : dspBufferlengthBytes = dspNumChannels * dspBitalign * dspBufferlengthSamples ; bufferlength in bytes, not samples
  Protected lpfRatio.f             : lpfRatio = PeekF(*lpfRatio)
  Protected lpfFactor.f            : lpfFactor = 1/dB2Dec(lpfRatio)
  Protected diff.f ; contains the difference of two contiguous samplevalues
  
  If *lastlpfSampleBuffer = 0 ; allocate new memory just the first time
    *lastlpfSampleBuffer = AllocateMemory(#MemID_lastLpfSampleBuffer, dspNumChannels * 4)
  EndIf
  
  If dspBitalign = 2 And lpfRatio <> -1 ; 16 Bit signal processing only here, because of PeekW in the following  
    
    While actCHSample < dspBufferlengthBytes-1 ; process the whole buffer
      lpfCHSample = PeekW(*outputbuffer+actCHSample) ; read in the actual 16Bit-sample-value
      lpfLastCHSample = PeekL(*lastlpfSampleBuffer+4*CHSwitch) ; read in the actual previous 16Bit-sample-value
      
      diff = lpfCHSample - lpfLastCHSample ; compare it with the previous one
      
      diff = diff * lpfFactor ; lower the difference between them by lpfFactor
      lpfCHSample = lpfLastCHSample + diff ; calculate the lpf-corrected sample-value
      
      PokeL(*lastlpfSampleBuffer+4*CHSwitch, lpfCHSample) ; the actual one will be the previous one for left channel
      PokeW(*outputbuffer + actCHSample, lpfCHSample) ; write it into the output-buffer
      
      CHSwitch + 1 : If CHSwitch = dspNumChannels : CHSwitch = 0 : EndIf ; circulate CHSwitch between all Channels
      actCHSample + dspBitalign ; increases the byteposition by bitalign (2 at 16 bit-audio)
    Wend
  EndIf  
  ProcedureReturn *outputbuffer ; overgive a pointer to our outputbuffer to the next dsp-station (will be *outputbuffer there)
EndProcedure
; needs during FSOUND_DSP_Create() for the *lpfRatio-parameter a pointer to a float-variable,
; whose value is in the range [0.0; +oo] and can be live updated during playback.
; It represents the factor for reducing the high freqencies in dB.
; Setting the float-variable to -1 sets the DSP to bypass.



;-
;- ___ myFMOD_GAIN
;-
Procedure.l myFMOD_GAIN(*originalbuffer.l, *outputbuffer.l, dspBufferlengthSamples.l, *Gain_dB.l)
  Protected actCHSample.l          : actCHSample = 0 ; CH means Channel, so it's value alters between all mutlichannel-sample-values e.g. L | R for stereo
  Protected dspBufferlengthBytes.l : dspBufferlengthBytes = dspNumChannels * dspBitalign * dspBufferlengthSamples ; bufferlength in bytes, not samples
  Protected Gain_Dec.f             : Gain_Dec = dB2Dec(PeekF(*Gain_dB))
  Protected actSampleVal.l
  Protected distortionOn.l

  If dspBitalign = 2 And PeekF(*Gain_dB) <> 0 And PeekF(*Gain_dB) <= 100 And PeekF(*Gain_dB) >= -100; 16 Bit signal processing only here, because of PeekW in the following  
      While actCHSample < dspBufferlengthBytes-1 ; process the whole buffer
        actSampleVal = PeekW(*outputbuffer+actCHSample)
        actSampleVal * Gain_Dec
        
        If actSampleVal > 32767 : actSampleVal = 32767 : distortionOn = #True

        ElseIf actSampleVal < -32768 : actSampleVal = -32768 : distortionOn = #True

        EndIf

        PokeW(*outputbuffer + actCHSample, actSampleVal) ; write it into the output-buffer
        
        actCHSample + dspBitalign ; increases the byteposition by bitalign (2 at 16 bit-audio)
      Wend
  EndIf  

  If distortionOn = #True
    SetGadgetText(#Gadget_DistortionOnText, ">> distortion << ")
  Else
    SetGadgetText(#Gadget_DistortionOnText, "")
  EndIf
  
  ProcedureReturn *outputbuffer ; overgive a pointer to our outputbuffer to the next dsp-station (will be *outputbuffer there)
EndProcedure
; needs during FSOUND_DSP_Create() for the *lpfRatio-parameter a pointer to a float-variable,
; whose value is in the range [-100.0; +100.0] and can be live updated during playback.
; This one gives the gain in dB.
; Setting it to -1 sets the DSP to bypass.
; Here the Textgadget #Gadget_DistortionOnText is needed to show whether the actual samplevalue overrides or not.



;-
;- ___ dB2Dec
;-
Procedure.f dB2Dec(dB.f)
  ProcedureReturn Pow(10.0, dB / 20.0)
EndProcedure



;-
;- ___ Dec2dB
;-
Procedure.f Dec2dB(dec.f)
  ProcedureReturn 20.0 * Log10(dec)
EndProcedure

[blueznl]

what would be faster in your opinion: multiple low pass filters as you described, 'subtracting' the results to split up an audio signal in bands, or use real fft to turn audio signals into something resembling an old fashioned 'light show'?

[Froggerprogger]

This 'Simple LPF' will not be useful for that for that, because it would not seperate the bands strong enough (only -6dB damp/oktave). I'm not sure if this is changeable, I'll have to take a closer look at the algo again (it's from July last year)
Further I think that it is much slower to do 512 post-LPF-subtractions with an array of about 1024 samplevalues per channel than doing one FFT onto these samples. (Fmod's FFT is very fast & good for this purpose - but of course just for 'reading' data, not for modifiing it)

For multiband-EQ the concept: sampledata -> FFT -> change frequency-amplitudes -> IFFT -> sampledata would be the best I think.

But I hadn't the time for understanding & programming FFT / IFFT until now.

[Froggerprogger]

I added HPF and Gain to the above code.

Replacing the line:

Code: Select all

diff = diff * lpfFactor

inside the LPF by

Code: Select all

diff = Pow(diff, lpfFactor)

and playing with the LPF-control-bar results in a cool sound-accident

And the moment I think it's not possible to increase the damp/oktave.
All linear changes made to the diff result in -6dB. I think this depends on the probability distribution of stepper curves in higher frequencies.

All nonlinear changes do affect the signal as if we would change the lpf-gain each sample an sound something like the above one.

[Road Runner]

Blueznl,
FFT is always faster when you want the amplitude of lots of frequencies within the signal.


Froggerprogger,
there are better ways to low/band/high pass filter data.

FFT is useful for complex filters or for times when you want to know how much of each frequency is present in a signal.

But, there is a whole class of time domain filters which operate directly on the data without the need for FFT.
These filters allow all sorts of characteristics such as a direct mimic of Butterworth or Chebychev filters, any number of poles/zeroes in the design allowing cut off's of well in excess of 6dB/octave.

Roughly, the processing required is 2N multiplies and 2N additions for a filter of order N.

The theory behind it is a bit involved but you can get software to design the filter for you and then just plug the resulting equation into your program. There's a good website at:
http://www-users.cs.york.ac.uk/~fisher/ ... /trad.html


The final page produced by that site gives you a "Recurrence relation". That's the bit you need.

y[n] is the nth output sample i.e. the nth filtered value.
x[n] is the nth input sample (i.e. nth sample in the original data)

That formula is very easy to plug into software to give "real" filters with well defined characteristics.

I hope that's of interest.

[Froggerprogger]

Thank you for the amazing link!
I mean especially the great stuff beyond this tool, the pages of Tony Fischer!

[dell_jockey]

FroggerProgger,

amazing link indeed, thanks!

Since you're obviously into signal processing I have a question for you: I want to remove the higher frequency components of a signal and want the remaining trend signal to be as close to the original signal, i.e. with delay reduced as much as possible. The problem is, I can't take numerous samples, I only have the last 30 to 50 samples to work with.
To this end I worked out a scheme where I take an exponential moving average of the original signal. After that I use this smoothed signal and take an exponential moving average of this signal as well. To reduce delay, I take the delta between these two ema's, and add it back to the first ema. That modified first ema is the smoothed signal I use for further processing.

Do you know better ways of getting a low-lagged smoothed version of a signal?

thanks!

[Road Runner]

dell_jockey,
use the link above. Plug in your values and a simple equation pops out which will filter your data for you.

The number of samples you have to work with is not that important. It just sets how closely you can control the cut-off frequency

[dell_jockey]

RoadRunner,

there is a relation between the number of samples available and what a filter can do. Nyquist's theorem tells us that if we have only the last, say, 50 samples, there's an upper limit of frequencies that can be detected. Surely this must be reflected in filter design.

Well, all of this wasn't PB related, but interesting nonetheless. Thanks again for the link!

[Road Runner]

dell_jockey,
that's not quite right. Nyquist relates sample rate to maximum frequency.
The number of samples is different variable.

Roughly, if you have 50 samples at sample rate of 1000Hz then you'll be able to specify filter cut-off frequencies to approximately 1000Hz/(50/2) = 40Hz. (rule of thumb only)

The maximum frequency you can handle is 1000Hz/2 = 500Hz. (that's Nyquist)


Now, take 100 samples (twice as many) then you can get twice the accuracy in specifying cut-off frequency (i.e about 20Hz) but the maximum frequency you can handle remains at 1000/2=500Hz.