SignalUtils.DcFilter Method

Overload List

#	Signature	Description
1	function DcFilter(NewValue: Double; var State: TCplx; alpha: Double): Double;	State parameter holds the filter state.
2	procedure DcFilter(alpha: Double; Num: TVec; Den: TVec);	Design a DC-blocking (high-pass) IIR filter from alpha.
3	procedure DcFilter(TransitionBandwidth: Double; FS: Double; Num: TVec; Den: TVec);	Design a DC filter.

Overload 1: function DcFilter(NewValue: Double; var State: TCplx; alpha: Double): Double;

State parameter holds the filter state.

#	Name	Type	Description
1	NewValue	Double	scalar
2	State	TCplx
3	alpha	Double	scalar

Returns: Double

Remarks:

NewValue is the next sample and alpha is typically between 0.99 and 0.9999 and must be < 1. Big alpha will cause longer filter delay and more ringing. State should be initialized to zero before the routine is called for the first time.

Overload 2: procedure DcFilter(alpha: Double; Num: TVec; Den: TVec);

Design a DC-blocking (high-pass) IIR filter from alpha.

#	Name	Type	Description
1	alpha	Double	scalar
2	Num	TVec	source TVec
3	Den	TVec	source TVec

Result: stored in self (calling object)

Remarks:

Builds the transfer function of a first-order DC blocker and stores it in num / den for use with IirInit / IirFilter; the filtered signal is centred around zero (its mean is removed). With unity-DC-gain correction g = alpha/big(alpha + 1/2(1-alpha)big) the filter is

H(z) = g (1 - z^(-1))/(1 - alpha z^(-1))

i.e. a differentiator/integrator pair with a zero at z=1 (DC) and a pole at z=alpha. Domain: 0 < alpha < 1 (typically 0.99-0.9999); alpha nearer 1 pushes the cut-off lower and lengthens the transient. A finite alpha in range yields finite coefficients.

Examples

uses MtxExpr, Math387, MtxVec, SignalUtils, MtxVecTee, MtxVecEdit;

procedure TForm1.Button1Click(Sender: TObject);
var b,c,Response,num,den,X: Vector;
n,i: integer;
IirState: TIirState;
DCState: TCplx;
begin
    Tone(b,300,5/300,0,1); //generate sine with 5 periods in 300 samples

    //   Alternative:
    //          b.RandGauss;

    b := b + 2;
    c.Copy(b);
    n := 10;
    DcFilter(0.95,num,den);
    IirInit(num,den,IirState);
    for i := 0 to b.Length div n-1 do   //only to test the streaming
begin
    b.SetSubRange(i*n,n);
    c.SetSubRange(i*n,n);
    IirFilter(b,c,IirState);
end;
b.SetFullRange;
c.SetFullRange;
DrawIt([b,c],['Unfiltered','Filtered'],'DC filter');

FrequencyResponse(num,den,Response,64);
DrawIt(Response,'Frequency response');

DcFilter(0.05,2,num,den);
IirInit(num,den,IirState);
for i := 0 to b.Length div n-1 do  // the loop is only to test the streaming
begin
    b.SetSubRange(i*n,n);
    c.SetSubRange(i*n,n);
    IirFilter(b,c,IirState);
end;
b.SetFullRange;
c.SetFullRange;
DrawIt([b,c],['Unfiltered','Filtered']);

FrequencyResponse(num,den,Response,64);
DrawIt(Response,'Frequency response');

DCState := C_ZERO;
for i := 0 to b.Length -1 do
c.Values[i] := DcFilter(b.Values[i],DcState,0.95);

DrawIt([b,c],['Unfiltered','Filtered']);
end;

Overload 3: procedure DcFilter(TransitionBandwidth: Double; FS: Double; Num: TVec; Den: TVec);

Design a DC filter.

#	Name	Type	Description
1	TransitionBandwidth	Double	scalar
2	FS	Double	scalar
3	Num	TVec	source TVec
4	Den	TVec	source TVec

Result: stored in self (calling object)

Remarks:

Design a DC filter with TransitionBandwidth and place the transfer function in Num (numerator) and Den (denominator). You can then use this transfer function to initialize an IIR filter with a call to IirInit.

A DC filtered signal will be centered around zero. This DC filter is a simple differentiator/integrator pair. Transition bandwidth is the width of the frequency band where the amplitude is not yet completely attenuated. With DC filters, the transition band starts at 0 Hz. Narrow transition band (TransitionBandwidth/FS ratio is small) will result in filters with longer delays. FS is the sampling frequency. The filter implements the following difference equation:

y[i] = x[i] - x[i-1] + alpha y[i-1]

x.. input signal
y.. output signal
alpha.. parameter

Alpha paremeter can control the 3dB frequency of the transition bandwidth:

alpha := 1-(TransitionBandwidth/FS)*Pi;

FS.. sampling frequency
TransitionBandwidth.. frequency up to which will the filter have more
                      then 3dB attenuation. Must be less then FS/2.

Examples

alpha := 1-(TransitionBandwidth/FS)*Pi;

FS.. sampling frequency
TransitionBandwidth.. frequency up to which will the filter have more
then 3dB attenuation. Must be less then FS/2.