Manual

Installation

To install FIRLSFilterDesign.jl, open up a Julia REPL and do:

pkg> update
pkg> add FIRLSFilterDesign

Usage

Unweighted least-squares FIR design

Let's say we want to design a bandpass filter with a passband from 1 Hz to 2 Hz, with a sampling frequency of 6 Hz, a filter order of 10, and symmetric filter coefficients.

We can do this as follows:

julia> using FIRLSFilterDesign;

julia> fs = 6;

julia> filter_order = 10;

julia> antisymmetric = false;

Set the frequency-band matrix:

julia> freq_bands = [0 1; 1 2; 2 3];
Note

Frequency bands must not overlap and must completely cover the range $[0, f_s/2]$.

Define the amplitude response:

julia> D = [0. 0.; 1. 1.; 0. 0.];

Now we can design the filter:

h = firls_design(filter_order, freq_bands, D, antisymmetric; fs = fs)
11-element Vector{Float64}:
  7.424814296795033e-17
  0.13783222385544808
 -3.898171832519375e-17
 -0.27566444771089604
  3.222575649650622e-17
  0.3333333333333333
  3.222575649650622e-17
 -0.27566444771089604
 -3.898171832519375e-17
  0.13783222385544808
  7.424814296795033e-17

Weighted least-squares FIR design

Now let's say we want to design the same filter, but this time we want to place more weight on errors in the passband. We can do this by defining the following weighting coefficient matrix:

julia> W = [1 2; 2 2; 2 1];

Now we can design the filter with the weight function:

h = firls_design(filter_order, bands_D, D, W, antisymmetric; fs = fs)
11-element Vector{Float64}:
  5.851941012279437e-17
  0.1403757865290026
 -7.540357228407273e-18
 -0.27528418302327057
  2.972899058715575e-17
  0.3342918305637282
  2.972899058715575e-17
 -0.27528418302327057
 -7.540357228407273e-18
  0.1403757865290026
  5.851941012279437e-17

Different input forms

You can use the firls_design function with several different input shapes for the frequency, amplitude response, and weighting function values. They are listed in the table below:

Frequency bandsAmplitude responseWeighting functionComments
MatrixUnion{Vector,Matrix}Union{Vector,Matrix}Vectors are interpreted as constant values over the frequency band.
VectorVectorVectorVectors are interpreted as frequency knotpoints and values at those knotpoints.
MatrixUnion{Vector,Matrix}N/AVectors are interpreted as constant values over the frequency band.
VectorVectorN/AVectors are interpreted as frequency knotpoints and values at those knotpoints.

Functions

FIRLSFilterDesign._to_impulse_responseMethod
_to_impulse_response(a, fir_type)

Creates a linear phase FIR filter based on fir_type and the coefficients in vector $a$, which was obtained by solving the linear equation $Qa = b$.

...

Arguments

-a::Vector : a vector of size (M+1,) with coefficients.

  • fir_type::FIR : indicates the type of FIR filter.

Outputs

-h::Vector : a vector of size (filter_order+1,) with the filter coefficients. ...

For type I FIR filters:

\[ h = \bigg[ a[M+1] \quad a[M] \quad \cdots \quad a[2] \quad a[1] \quad a[2] \quad \cdots \quad a[M+1] \bigg]^T\]

source
FIRLSFilterDesign._update_trig_arg_b!Method
_update_trig_arg_b!(_αn, n, fir_type::Union{FIR_I,FIR_II})

Updates the argument of the trigonometric functions in bn! by multiplying with the current n.

...

Arguments

  • _αn::Vector : a vector of size (J,)
  • n::Real
  • fir_type{Union{FIR_I,FIR_II}}

...

source
FIRLSFilterDesign._update_trig_arg_b!Method
_update_trig_arg_b!(_αn, n, fir_type::Union{FIR_III,FIR_IV})

Updates the argument of the trigonometric functions in bn! by multiplying with $n$ and subtracting $\pi$. The subtraction of $\pi$ is necessary because when the filter is antisymmetric (type III and IV FIR filters), the filter response is a sum of sines instead of cosines and $\sin(x) = \cos(x - \pi/2)$ (see page 12 and 13 of this).

...

Arguments

  • _αn::Vector : a vector of size (J,)
  • n::Real
  • fir_type::Union{FIR_III,FIR_IV}

...

source
FIRLSFilterDesign.bn!Method
bn!(_bn, k, f, a, b, c, d)

Special case for when $n = 0$, since then the integral is simplified:

\[g_j(f,0) = (c_j f+d_j) (a_j f+b_j) \cos(\pi k f 0) = (c_j f+d_j) (a_j f+b_j)\]

And the antiderivative becomes:

\[G_j(f,0) = a_j c_j \frac{f^3}{3} + (a_j d_j + b_j c_j)\frac{f^2}{2} + b_j d_j f\]

...

Arguments

  • _bn::Vector : a vector of size (J,) that is used to store the intermdediate values.
  • k::Real : equal to $2/f_s$.
  • f::Matrix : a matrix of size (J,2) which contains rows of sequential frequency bands, spanning [0, fs/2].
  • a::Vector : a vector of size (J,) with the $a_j$ values in the equation $a_j f + b_j$ that equates to the linear function that describes the desired frequency response in the $j^{th}$ frequency band.
  • b::Vector : a vector of size (J,) with the $b_j$ values in the equation $a_j f + b_j$ that equates to the linear function that describes the desired frequency response in the $j^{th}$ frequency band.
  • c::Vector : a vector of size (J,) with the $c_j$ values in the equation $c_j f + d_j$ that equates to the linear function that describes the error weighting function in the $j^{th}$ frequency band.
  • d::Vector : a vector of size (J,) with the $d_j$ values in the equation $c_j f + d_j$ that equates to the linear function that describes the error weighting function in the $j^{th}$ frequency band.

...

source
FIRLSFilterDesign.bn!Method
bn!(_bn, n, k, _αn, _βn², γ, _δn, fir_type)

Calculates the elements of the b-vector, which are equal to:

\[b[i] = \frac{2}{f_s} \int_0^{f_s/2} W(f) D(f) cos(\pi \frac{2}{f_s} n f) df, \quad i = 1, 2, \cdots, M+1\]

Both $D(f)$ and $W(f)$ are piecewise linear functions, and $n$ is calculated by idx2n_b. Using $k = 2/f_s$, this integral becomes:

\[b[i] = k \sum_{j=1}^J \int_{F_{j,1}}^{F_{j,2}} g_j(f,n) df = k \sum_{j=1}^J \int_{F_{j,1}}^{F_{j,2}} \big(c_j f+d_j\big) \big(a_j f+b_j\big) \cos(\pi k n f) df, \quad i = 1, 2, \cdots, M+1\]

Where:

  • $F_{j,1}$ and $F_{j,2}$ are the lower and upper bound of the $j^{th}$ frequency band.
  • $a_j$ and $b_j$ are the parameters of the linear function that describes the desired frequency response in the $j^{th}$ frequency band.
  • $c_j$ and $d_j$ are the parameters of the linear function that describes the error weighting function in the $j^{th}$ frequency band.

The antiderivative of $g(f,n)$ is equal to:

\[G_j(f,n) = \frac{1}{\pi^3 k^3 n^3} \bigg(\sin\big(\alpha(f) n\big)\Big(\beta_j(f) n^2 + \gamma_j(f)\Big) + \delta_j(f) n \cos\big(\alpha(f) n\big)\bigg) + constants\]

Where:

  • $\alpha(f) = \pi kf$,
  • $\beta(f) = \pi^2 k^2 \big(acf^2 + (ad+bc)f + bd\big)$,
  • $\gamma(f) = -2ac$,
  • $\delta(f) = \pi k\big(2acf + ad + bc\big)$.

Note that the subscripts have been dropped here for clarity.

...

Arguments

  • _bn::Vector : a vector of size (J,) that is used to store the intermediate values.
  • n::Integer : integer denoting the current cosine mode.
  • k::Real : equal to 2/fs.
  • _αn::Vector : a vector of size (J,) that holds the values of $\alpha n$.
  • _βn²::Vector : a vector of size (J,) that holds the values of $\beta n^2$.
  • γ::Vector : a vector of size (J,) that holds the values of $\gamma$.
  • _δn::Vector : a vector of size (J,) that holds the values of $\delta n$.
  • fir_type::FIR : indicates the type of FIR filter.

Outputs

  • bn : $n^{th}$ element in the $b$ vector.

...

source
FIRLSFilterDesign.bn_n0!Method
bn_n0!(_bn, n, k, f, a, b, c, d, _αn, _βn², γ, _δn, fir_type)

Dispatches to the correct function to calculate the first element of the b vector, based on the type of FIR filter. Needed because for type I FIR filters the value of n at the first iteration is 0.

source
FIRLSFilterDesign.constants_bMethod
constants_b(f, D, W)

Calculates data that is reused at every evaluation of bn!.

...

Arguments

  • f::Matrix : a matrix of size (N,2) which contains rows of sequential frequency bands, spanning the interval [0, fs/2].
  • D::Matrix : a matrix of size (N,2) which contains rows of desired frequency response values for the frequency bands in f. The first and second columns indicate the desired response at the lower and upper bound of the frequency bands, interpolated linearly in between.
  • W::Matrix : a matrix of size (N,2) which contains rows of weighting coefficients for the frequency bands in f. The first and second columns indicate the weighting at the lower and upper bound of the frequency bands, interpolated linearly in between.

...

source
FIRLSFilterDesign.firls_designMethod
firls_design(filter_order::Integer, bands_DW::Matrix, D::Matrix, W::Matrix, antisymmetric::Bool; fs::Real = 1, solver::Function = \)

Designs a linear-phase FIR filter.

Arguments

  • filter_order::Integer : the order of the FIR filter.
  • bands_DW::Matrix : a matrix of size (N,2) which contains rows of sequential frequency bands, spanning [0, fs/2].
  • D::Matrix : a matrix of size (N,2) which contains rows of amplitude responses for the frequency bands in bands_DW. The first and second columns indicate the amplitude response at the lower and upper bound of the frequency bands, interpolated linearly in between.
  • W::Matrix : a matrix of size (N,2) which contains rows of weighting coefficients for the frequency bands in bands_DW. The first and second columns indicate the weighting at the lower and upper bound of the frequency bands, interpolated linearly in between.
  • antisymmetric::Bool : a Boolean that signifies whether the filter coefficients will be anti-symmetric, as used in type III and IV FIR filters.
  • fs::Real : the sampling frequency.
  • solver::Function : the function that is called to solve the equation $Qa = b$, with the function call: solver(Q,b) which returns a.

Outputs

  • h : a vector of linear-phase FIR filter coefficients.
source
FIRLSFilterDesign.firls_designMethod
firls_design(filter_order::Integer, bands_DW::Matrix, D::Union{Vector,Matrix}, antisymmetric::Bool; fs::Real = 1, solver::Function = \)

Arguments

  • filter_order::Integer : the order of the FIR filter.
  • bands_DW::Matrix : a matrix of size (N,2) which contains rows of sequential frequency bands, spanning [0, fs/2].
  • D::Union{Vector,Matrix} : a matrix of size (N,2), or a vector of size (N,), which amplitude responses for the frequency bands in bands_DW. In the case of a matrix, the first and second columns indicate the amplitude response at the lower and upper bound of the frequency bands, interpolated linearly in between. In the case of a vector the elements the amplitude response is piecewise constant.
  • antisymmetric::Bool : a Boolean that signifies whether the filter coefficients will be anti-symmetric, as used in type III and IV FIR filters.
  • fs::Real : the sampling frequency.
  • solver::Function : the function that is called to solve the equation $Qa = b$, with the function call: solver(Q,b) which returns a.

Outputs

  • h : a vector of linear-phase FIR filter coefficients.
source
FIRLSFilterDesign.firls_designMethod
firls_design(filter_order::Integer, bands_DW::Matrix, D::Union{Vector,Matrix}, W::Union{Vector,Matrix}, antisymmetric::Bool; fs::Real = 1, solver::Function = \)

Arguments

  • filter_order::Integer : the order of the FIR filter.
  • bands_DW::Matrix : a matrix of size (N,2) which contains rows of sequential frequency bands, spanning [0, fs/2].
  • D::Union{Vector,Matrix} : a matrix of size (N,2), or a vector of size (N,), which amplitude responses for the frequency bands in bands_DW. In the case of a matrix, the first and second columns indicate the amplitude response at the lower and upper bound of the frequency bands, interpolated linearly in between. In the case of a vector the elements the amplitude response is piecewise constant.
  • W::Union{Vector,Matrix} : a matrix of size (N,2), or a vector of size (N,), which contains weighting function values for the frequency bands in bands_DW. In the case of a matrix, the first and second columns indicate the weighting function values at the lower and upper bound of the frequency bands, interpolated linearly in between. In the case of a vector the elements the weighting function is piecewise constant.
  • antisymmetric::Bool : a Boolean that signifies whether the filter coefficients will be anti-symmetric, as used in type III and IV FIR filters.
  • fs::Real : the sampling frequency.
  • solver::Function : the function that is called to solve the equation $Qa = b$, with the function call: solver(Q,b) which returns a.

Outputs

  • h : a vector of linear-phase FIR filter coefficients.
source
FIRLSFilterDesign.firls_designMethod
firls_design(filter_order::Integer, knotpoints_D::Vector, D::Vector, antisymmetric::Bool; fs::Real = 1, solver::Function = \)

Arguments

  • filter_order::Integer : the order of the FIR filter.
  • knotpoints_D::Vector : a vector of size (N,) which contains frequency knotpoints, spanning [0, fs/2].
  • D::Vector : a vector of size (N,) which contains amplitude response values for the frequency knotpoints in knotpoints_D::Vector.
  • antisymmetric::Bool : a Boolean that signifies whether the filter coefficients will be anti-symmetric, as used in type III and IV FIR filters.
  • fs::Real : the sampling frequency.
  • solver::Function : the function that is called to solve the equation $Qa = b$, with the function call: solver(Q,b) which returns a.

Outputs

  • h : a vector of linear-phase FIR filter coefficients.
source
FIRLSFilterDesign.firls_designMethod
firls_design(filter_order::Integer, knotpoints_DW::Vector, D::Vector, W::Vector, antisymmetric::Bool; fs::Real = 1, solver::Function = \)

Arguments

  • filter_order::Integer : the order of the FIR filter.
  • knotpoints_DW::Vector : a vector of size (N,) which contains frequency knotpoints, spanning [0, fs/2].
  • D::Vector : a vector of size (N,) which contains amplitude response values for the frequency knotpoints in knotpoints_DW.
  • W::Vector : a vector of size (N,) which contains weighting function values for the frequency knotpoints in knotpoints_DW.
  • antisymmetric::Bool : a Boolean that signifies whether the filter coefficients will be anti-symmetric, as used in type III and IV FIR filters.
  • fs::Real : the sampling frequency.
  • solver::Function : the function that is called to solve the equation $Qa = b$, with the function call: solver(Q,b) which returns a.

Outputs

  • h : a vector of linear-phase FIR filter coefficients.
source
FIRLSFilterDesign.get_QMethod
get_Q(M, f, W, fir_type)

Constructs the matrix $Q$ used in the equation $Qa = b$, based on a set of weights.

Arguments

  • M::Integer : indicator of the amount of elements needed.
  • f::Matrix : a matrix of size (N,2) which contains rows of sequential frequency bands, spanning [0, fs/2].
  • W::Matrix : a matrix of size (N,2) which contains rows of weighting coefficients for the frequency bands in f. The first and second columns indicate the weighting at the lower and upper bound of the frequency bands, interpolated linearly in between.
  • fir_type::FIR : indicates the type of FIR filter.

Outputs

  • Q::Matrix : the matrix $Q$ used in the equation $Qa = b$.
source
FIRLSFilterDesign.get_QMethod
get_Q(M, fir_type)

Constructs the matrix $Q$ used in the equation $Qa = b$, when there are no weights. Which results in $Q$ being the identity matrix.

Arguments

  • M::Integer : indicator of the amount of elements needed.
  • fir_type::FIR : indicates the type of FIR filter.

Outputs

  • Q::Matrix : the matrix $Q$ used in the equation $Qa = b$.
source
FIRLSFilterDesign.get_QMethod
get_Q(M, fir_type::FIR_I)

Constructs the matrix $Q$ used in the equation $Qa = b$, when there are no weights and the FIR filter is of type I.

Arguments

  • M::Integer : indicator of the amount of elements needed.
  • fir_type::FIR_I : indicates the type of FIR filter is I.

Outputs

  • Q::Matrix : the matrix $Q$ used in the equation $Qa = b$.
source
FIRLSFilterDesign.get_bMethod
get_b(M, f, D, W, fir_type)

Finds the vector $b$ used in the equation $Qa = b$. ...

Arguments

  • M::Integer : Size of the b vector is M+1 .
  • f::Matrix : a matrix of size (N,2) which contains rows of sequential frequency bands, spanning the interval [0, fs/2].
  • D::Matrix : a matrix of size (N,2) which contains rows of desired frequency response values for the frequency bands in f. The first and second columns indicate the desired response at the lower and upper bound of the frequency bands, interpolated linearly in between.
  • W::Matrix : a matrix of size (N,2) which contains rows of weighting coefficients for the frequency bands in f. The first and second columns indicate the weighting at the lower and upper bound of the frequency bands, interpolated linearly in between.
  • fir_type::FIR : indicates the type of FIR filter.

Outputs

  • b_out::Vector : a vector of size (M+1,), the b-vector used in the equation $Qa = b$.

...

source
FIRLSFilterDesign.get_flength_MMethod
get_flength_M(filter_order)

Determines the length of the FIR filter and the number of amplitude coefficients needed, based on its order.

Arguments

  • filter_order::Integer : the order of the FIR filter.

Outputs

  • filter_length::Integer : the number of filter coefficients
  • M::Integer : the number of unique amplitude coefficients needed to form the filter is equal to $M+1$
source
FIRLSFilterDesign.get_qMethod
get_q(M, f, W, fir_type)

Finds the vector $q$ which is used to populate the matrix $Q$. ...

Arguments

  • M::Integer : indicator of the amount of elements needed.
  • f::Matrix : a matrix of size (N,2) which contains rows of sequential frequency bands, spanning [0, fs/2].
  • W::Matrix : a matrix of size (N,2) which contains rows of weighting coefficients for the frequency bands in f. The first and second columns indicate the weighting at the lower and upper bound of the frequency bands, interpolated linearly in between.
  • fir_type::FIR : indicates the type of FIR filter.

Outputs

  • q_out::Vector : a vector of q-values that are used to fill in the Q-matrix.

...

source
FIRLSFilterDesign.idx2n_bMethod
idx2n_b(idx, fir_type::Union{FIR_I,FIR_III})

Determine the value of $n$ based on the index in the b-vector. For type I and III FIR filters the following holds: $n = i - 1$. Based on page 10 and 12 of this.

...

Arguments

  • i::Integer : index in the b-vector.
  • fir_type::Union{FIR_I,FIR_III}

Outputs

  • n::Real

...

source
FIRLSFilterDesign.idx2n_bMethod
idx2n_b(idx, fir_type::Union{FIR_II,FIR_IV})

Determine the value of $n$ based on the index in the b-vector, Where $n$ is used in For type II and IV FIR filters the following holds: $n = i - 1/2$. Based on page 11 and 13 of this.

...

Arguments

  • i::Integer : index in the b-vector.
  • fir_type::Union{FIR_II,FIR_IV}

Outputs

  • n::Real

...

source
FIRLSFilterDesign.infer_fir_typeMethod
infer_fir_type(is_odd, is_antisymmetric)

Determines the type of FIR filter to be designed, based on:

  • Whether the number of filter coefficients is odd (is_odd)
  • Whether the filter should be antisymmetric (is_antisymmetric)

The result is a fir_type, which can be:

  1. FIR_I, when filter length is odd and the filter is not antisymmetric
  2. FIR_II, when filter length is even and the filter is not antisymmetric
  3. FIR_III, when filter length is odd and the filter is antisymmetric
  4. FIR_IV, when filter length is even and the filter is antisymmetric

Arguments

  • is_odd::Bool
  • is_antisymmetric::Bool

Outputs

  • fir_type::FIR
source