 ## “Non-Causal“ Digital Chebyshev Filter with Linear Phase Response. C Implementation and Octave Script

### 2. Introduction

##### Often in Digital Signal Processing applications it is necessary to use filters with a linear phase characteristic. In the case is possible to avoid phase distortion of the signal. It is impossible to construct a stable IIR filter with an ideal linear phase characteristic, since the following condition must be met: ##### The (1) means that all the filter poles must be inverted: ##### Figure 2-1: Non-Causal IIR Filter with Null Phase/Delay Response

##### (7) matches the filter shown in Figure 2-2. Figure 2-2: Equivalent filter for Figure 2-1

### 3. Time Reverse (TR) Implementation from 

##### See the Figure 2-1. First the input signal x(n) is passed through the first filter H(z) after which the g(n) must be reversed in time and then g(-n) is passed through the second filter H(z) again. We divide g(n) into blocks of length N and process each block separately. See the Figure 3-1. Figure 3-1: Block Definitions

##### Blocks in reverse time you can see in the Figure 3-2. Figure 3-2: Reverse Time for g(-n) Blocks

##### At the output second H(z) the y (- n) is got . See the Figure 3-3. Figure 3-3: Reverse Time for y(-n) Blocks

##### The subsequent time reversion for y(- n) blocks gets the output signal y(n). See the Figure 3-4. Figure 3-4: Output y(n) Blocks

• ##### Processing of the input signal x(n) is possible if the (N + Noverlap) samples is accumulated and then the algorithm shown in Figure 2-1 can be started. Thus the filter has the delay: ##### The formula (10) is valid only for the first / starting block and then only the next N samples should be accumulated, since the next block already has Noverlap samples from previous step and the filter time delay is reduced: #####  is defined another time delay of the filter: ##### It means that h(n) decreases with time and one can find a Noverlap where the values of h(n) can be neglected without a significant error: ### 4. Estimate the Transient Process of the Analog Filter Prototype

##### To begin with consider the simplest RC circuit: Figure 4-1: RC Analog Low Pass Filter

##### For example: the RC circuit (low pass filter) is described with the first order differential equation: ##### or in the Laplace form: ##### Now let’s feed the Heaviside function 1(t) to the RC circuit input ##### Laplace transform of the 1(t): ##### Then using (16) and (18), we get the filter response at the output: ##### This corresponds to the signal: ##### The transition process takes approximately 5τ: #####  is defined a generalized time constant as follows: ##### Let’s check the definition (22) for the RC filter: ##### We will calculate the generalized time constant (22) for the circuit: ##### Let’s use two properties of Laplace Transform: ##### Now the step function Heaviside 1(t) will be fed to the input (24). Then the signal at the filter output is: ##### Using (25) and (26) can be written: ##### Note: (28) is true for the case n>m. If n=m then H (+∞)=bn/an ##### Now we’ll calculate: ##### Suppose t → +∞, which corresponds to s → 0 in (30), then (22) can be written using (26): ##### Calculating (31) for (24): ##### Using the 5τ rule for (32): ##### Let’s check the formula (33) for (16): ##### Then: ### 5. Chebyshev Analog Filter Prototype

##### Then ##### Analog Chebyshev Filter (40) has a magnitude response in Figure 5.1, a phase response in Figure 5-2 and a step response in Figure 5-3. Figure 5-1: Magnitude Response of the Analog Low Pass Chebyshev Filter, Fband=1kHz, 0.5dB Ripple, 4th order

##### The Magnitude response meets the filter requirements: Fband=1kHz, Ripple = 0.5dB Figure 5-2: Phase Response of the Analog Low Pass Chebyshev Filter, Fband=1kHz, 0.5dB Ripple, 4th order

##### The phase response is non-linear. Figure 5-3: Step Response of the Analog Low Pass Chebyshev Filter, Fband=1kHz, 0.5dB Ripple, 4th order

##### Now let’s represent (40) as (24) and use the formula (33): ##### Figure 5-3 shows that after 2.4 ms there are two more fluctuations with levels of 0.97 (3%) and 1.016 (1.6%), so it is more reasonable to choose the transition time: ### 6. Digital Chebyshev Filter

##### Using the analog prototype of the Chebyshev filter (40) can be calculated the digital filter. The bilinear transformation shall be used for it: ##### Choose the sampling frequency: ##### Digital Chebyshev Filter (45) has a magnitude response in Figure 6.1, a phase response in Figure 6-2 and a step response in Figure 6-3. Figure 6-1: Magnitude Response of the Digital Low Pass Chebyshev Filter, Fband=1kHz, 0.5dB Ripple, 4th order, fs=12kHz (Ts=83.3μs) Figure 6-2: Phase Response of the Digital Low Pass Chebyshev Filter, Fband=1kHz, 0.5dB Ripple, 4th order, fs=12kHz (Ts=83.3μs) Figure 6-3: Step Response of the Digital Low Pass Chebyshev Filter, Fband=1kHz, 0.5dB Ripple, 4th order

### 7. “Non-Causal“ Digital Chebyshev Filter with the Linear Phase

##### Now let’s analyze and implement the filter shown in Figure 2-1. Using (9), we get the Magnitude Response of the filter. See the Figure 7-1. Figure 7-1: Magnitude Response of the non-causal Digital Low Pass Chebyshev Filter, Fband=1kHz, 0.5dB Ripple, 4th order, fs=12kHz (Ts=83.3μs). See the Figure 2-1

##### Figure 7-2: Step Response of the Non-Causal Low Pass Chebyshev Filter 4th Order. Blue – input, red – output

##### Now we will feed the sum of two sinusoids: 500Hz (inside the bandwidth) and 5kHz (outside the bandwidth) to the filter input. See the Figure 7-3. Only the 500Hz signal remains at the filter output. See the Figure 7-4 Figure 7-3: Sum of the two sines 500Hz and 5kHz. Filter Input Signal. Figure 7-4: Response on the Sum of the two sines 500Hz and 5kHz. Non-Causal Low Pass Chebyshev Filter 4th Order. Filter Output Signal.

##### Comparing Figure 7-3 and Figure 7-4 we can see that the phase of the input and output signals of 500Hz does not change, which indicates a “zero” signal delay or a linear phase response. In reality the filter has a delay but it is related to the fact that the output signal can be calculated only when N + Noverlap input samples are accumulated. See the formula (10), (11). If we neglect the Tdsp time then for our example: ### 8. Octave GNU Script DigitalChebyshevFilter.m

##### Basic Functions:

% Plotting the Magnitude Response for the 4th order analog filter: 1/(a11*s^2+a12*s+a13)(a21*s^2+a22*s+a23)
% Input parameters:
% SA1 – denominator of the second order section
% SA2 – denominator of the second order section
% FhighInHz – Frequency space in Hz
function AnalogMagnitudeResponse(SA1, SA2, FhighInHz)

% Plotting the Phase Response for the 4th order analog filter: 1/(a11*s^2+a12*s+a13)(a21*s^2+a22*s+a23)
% Input parameters:
% SA1 – denominator 1 of the second order section
% SA2 – denominator 2 of the second order section
% FhighInHz – Frequency space in Hz
function AnalogPhaseResponse(SA1, SA2, FhighInHz)

% Plotting the Step Response for the 4th order analog filter: 1/(a11*s^2+a12*s+a13)(a21*s^2+a22*s+a23)
% Input parameters:
% SA1 – denominator 1 of the second order section
% SA2 – denominator 2 of the second order section
% TimeVector – Time space
function plot_StepResponseAnalogFilter(SA1, SA2, TimeVector)

% Plotting the Magnitude Response for the 4th order digital filter: (b11+b12*z^-1+b13*z^-2)(b21+b22*z^-1+b23*z^-2)/(a11+a12*z^-1+a13*z^-2)(a21+a22*z^-1+a23*z^-2)
% Input parameters:
% ZB1 – numerator 1 of the second order section
% ZA1 – denominator 1 of the second order section
% ZB2 – numerator 2 of the second order section
% ZA2 – denominator 2 of the second order section
% Tsample – sampling time in seconds
function plot_DigitalMagnitudeResponse(ZB1, ZA1, ZB2, ZA2, Tsample)

% Plotting the Phase Response for the 4th order digital filter: (b11+b12*z^-1+b13*z^-2)(b21+b22*z^-1+b23*z^-2)/(a11+a12*z^-1+a13*z^-2)(a21+a22*z^-1+a23*z^-2)
% Input parameters:
% ZB1 – numerator 1 of the second order section
% ZA1 – denominator 1 of the second order section
% ZB2 – numerator 2 of the second order section
% ZA2 – denominator 2 of the second order section
% Tsample – sampling time in seconds
function plot_DigitalPhaseResponse(ZB1, ZA1, ZB2, ZA2, Tsample)

% Plotting the Step Response for the 4th order digital filter: (b11+b12*z^-1+b13*z^-2)(b21+b22*z^-1+b23*z^-2)/(a11+a12*z^-1+a13*z^-2)(a21+a22*z^-1+a23*z^-2)
% Input parameters:
% ZB1 – numerator 1 of the second order section
% ZA1 – denominator 1 of the second order section
% ZB2 – numerator 2 of the second order section
% ZA2 – denominator 2 of the second order section
% SampleNumber – Sample number of the step function
function plot_StepResponseDigitalFilter(ZB1, ZA1, ZB2, ZA2, SampleNumber)

% Plotting the non-causal Magnitude Response for the 4th order digital filter: (b11+b12*z^-1+b13*z^-2)(b21+b22*z^-1+b23*z^-2)/(a11+a12*z^-1+a13*z^-2)(a21+a22*z^-1+a23*z^-2)
% Input parameters:
% ZB1 – numerator 1 of the second order section
% ZA1 – denominator 1 of the second order section
% ZB2 – numerator 2 of the second order section
% ZA2 – denominator 2 of the second order section
% Tsample – sampling time in seconds
function plot_non_causal_DigitalMagnitudeResponse(ZB1, ZA1, ZB2, ZA2, Tsample)

% Plotting the non-causal Magnitude Response for the 4th order digital filter: (b11+b12*z^-1+b13*z^-2)(b21+b22*z^-1+b23*z^-2)/(a11+a12*z^-1+a13*z^-2)(a21+a22*z^-1+a23*z^-2)
% Input parameters:
% ZB – filter numerator
% ZA – filter denominator
% Xinput – input samples
% Nblock – block size
% Nov – overlap size
function Youtput = non_causal_digitalFilter(ZB, ZA, Xinput, Nblock, Nov)

% Integer Test Signal to File
% Support C-Code
% Input parameters:
% SignalAmpl – amplitude of the test signal
% SignalFreqInHz – frequency of the test signal in Hz
% SignalLength – table length
% Tsample – sampling time in s
% FileNameString – output file name
% Output parameters:
% TestSignal – test signal vector
function TestSignal = IntegerTestSignals2file(SignalAmpl, SignalFreqInHz, SignalLength, Tsample, FileNameString)

% Integer Table to File
% Support C-Code
% Input parameters:
% ParameterVector – data array
% FileNameString – output file name
function IntegerParamVector2file(ParameterVector, FileNameString)

% Float Table to File
% Support C-Code
% Input parameters:
% ParameterVector – data array
% FileNameString – output file name
function FloatParamVector2file(ParameterVector, FileNameString)

% IIR Filter Coefficients Quantization
% Input parameters:
% ZA, ZB – IIR coefficients
% ScaleFactor – integer scale
% Output parameters:
% ZA_quantization – integer format
% ZB_quantization – integer format
function [ZA_quantization ZB_quantization] = IirFilterQuantization(ZA, ZB, ScaleFactor)

### 9. C Language non_causal_filter.c/h

##### non_causal_filter_integer.c/h – fixed point arithmeticnon_causal_filter_float.c/h – float arithmeticBase functions:

/*
IIR Non-Causal Filter Initialization
Clear the both histories:
iir_history_first and iir_history_second
Input
void
Return
void
*/
void iir_non_causal_init(void);

/*
IIR Non-Causal Filter
Input
sint16 buffer[N_BLOCK_SIZE+N_OVERLAP_SIZE] – pointer on the input sample buffer
buffer[0…N_OVERLAP_SIZE-1] – from previous processing step: first H(z) was already done
Return
buffer[N_BLOCK_SIZE] – filter output buffer
*/
void iir__non_causal_filter(sint16* buffer)

/* Test */
int main(void); – Check the filter with the step signal and sum of the two sines 500Hz and 5kHz