Biquad Filter

Share your advanced PureBasic knowledge/code with the community.
AndyMK
Enthusiast
Enthusiast
Posts: 540
Joined: Wed Jul 12, 2006 4:38 pm
Location: UK

Biquad Filter

Post by AndyMK »

Code: Select all

;************************************************************
;                                                           *
;  Date Created 14/08/2022                                  *
;  Author AndyCY                                            *
;  Original Source https://github.com/hosackm/BiquadFilter  *
;                                                           *
;************************************************************

Enumeration
	#LOWPASS
	#HIGHPASS
	#BANDPASS
	#NOTCH
	#PEAK
	#LOWSHELF
	#HIGHSHELF
EndEnumeration

Structure biquad
	a0.f
	a1.f
	a2.f
	b0.f
	b1.f
	b2.f
	prev_input_1.f
	prev_input_2.f
	prev_output_1.f
	prev_output_2.f
	type.s
EndStructure

Procedure bq_load_coefficients(*bq.biquad, filter_type.l,	A.f, omega.f, sn.f, cs.f, alpha.f, beta.f)
  
  Select filter_type
    Case #LOWPASS
      *bq\b0 = (1.0 - cs) /2.0
      *bq\b1 = 1.0 - cs
      *bq\b2 = (1.0 - cs) /2.0
      *bq\a0 = 1.0 + alpha
      *bq\a1 = -2.0 * cs
      *bq\a2 = 1.0 - alpha
      *bq\type = "LOWPASS"
      
    Case #HIGHPASS
      *bq\b0 = (1 + cs) /2.0
      *bq\b1 = -(1 + cs)
      *bq\b2 = (1 + cs) /2.0
      *bq\a0 = 1 + alpha
      *bq\a1 = -2 * cs
      *bq\a2 = 1 - alpha
      *bq\type = "HIGHPASS"
      
    Case #BANDPASS
      *bq\b0 = alpha
      *bq\b1 = 0
      *bq\b2 = -alpha
      *bq\a0 = 1 + alpha
      *bq\a1 = -2 * cs
      *bq\a2 = 1 - alpha
      *bq\type = "BANDPASS"
      
    Case #NOTCH
      *bq\b0 = 1
      *bq\b1 = -2 * cs
      *bq\b2 = 1
      *bq\a0 = 1 + alpha
      *bq\a1 = -2 * cs
      *bq\a2 = 1 - alpha
      *bq\type = "NOTCH"
      
    Case #PEAK
      *bq\b0 = 1 + (alpha * A)
      *bq\b1 = -2 * cs
      *bq\b2 = 1 - (alpha * A)
      *bq\a0 = 1 + (alpha /A)
      *bq\a1 = -2 * cs
      *bq\a2 = 1 - (alpha /A)
      *bq\type = "PEAK"
      
    Case #LOWSHELF
      *bq\b0 = A * ((A + 1) - (A - 1) * cs + beta * sn)
      *bq\b1 = 2 * A * ((A - 1) - (A + 1) * cs)
      *bq\b2 = A * ((A + 1) - (A - 1) * cs - beta * sn)
      *bq\a0 = (A + 1) + (A - 1) * cs + beta * sn
      *bq\a1 = -2 * ((A - 1) + (A + 1) * cs)
      *bq\a2 = (A + 1) + (A - 1) * cs - beta * sn
      *bq\type = "LOWSHELF"
      
    Case #HIGHSHELF
      *bq\b0 = A * ((A + 1) + (A - 1) * cs + beta * sn)
      *bq\b1 = -2 * A * ((A - 1) + (A + 1) * cs)
      *bq\b2 = A * ((A + 1) + (A - 1) * cs - beta * sn)
      *bq\a0 = (A + 1) - (A - 1) * cs + beta * sn
      *bq\a1 = 2 * ((A - 1) - (A + 1) * cs)
      *bq\a2 = (A + 1) - (A - 1) * cs - beta * sn
      *bq\type = "HIGHSHELF"
      
    Default
      
  EndSelect
  
EndProcedure

Procedure.l bq_new(filter_type.l , frequency.f , q.f, dbGain.f,	sample_rate.l)
  
 	*tmp.biquad = AllocateStructure(biquad)
  	A.f = Pow(10, dbGain / 20); convert to db
  	omega.f = (2 * #PI * frequency) / sample_rate
  	sn.f = Sin(omega)
  	cs.f = Cos(omega)
  	alpha.f = sn / (2*Q)
  	beta.f = Sqr(A + A)

	bq_load_coefficients(*tmp, filter_type, A, omega, sn, cs, alpha, beta);
	
	*tmp\a1 / *tmp\a0
	*tmp\a2 / *tmp\a0
	*tmp\b0 / *tmp\a0
	*tmp\b1 / *tmp\a0
	*tmp\b2 / *tmp\a0
	
	*tmp\prev_input_1 = 0.0
	*tmp\prev_input_2 = 0.0
	*tmp\prev_output_1 = 0.0
	*tmp\prev_output_2 = 0.0

	ProcedureReturn *tmp
		
EndProcedure


Procedure.f bq_process(*bq.biquad, input.f)
  
 	output.f = (*bq\b0 * input) + (*bq\b1 * *bq\prev_input_1) + (*bq\b2 * *bq\prev_input_2) - (*bq\a1 * *bq\prev_output_1) -	(*bq\a2 * *bq\prev_output_2)
	*bq\prev_input_2 = *bq\prev_input_1
	*bq\prev_input_1 = input
	*bq\prev_output_2 = *bq\prev_output_1
	*bq\prev_output_1 = output
	
	ProcedureReturn output
	
EndProcedure

Procedure bq_destroy(*bq)
  
 	FreeMemory(*bq)
  
EndProcedure
User avatar
idle
Always Here
Always Here
Posts: 5040
Joined: Fri Sep 21, 2007 5:52 am
Location: New Zealand

Re: Biquad Filter

Post by idle »

thanks for sharing
Post Reply