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];
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 bands | Amplitude response | Weighting function | Comments |
---|---|---|---|
Matrix | Union{Vector,Matrix} | Union{Vector,Matrix} | Vectors are interpreted as constant values over the frequency band. |
Vector | Vector | Vector | Vectors are interpreted as frequency knotpoints and values at those knotpoints. |
Matrix | Union{Vector,Matrix} | N/A | Vectors are interpreted as constant values over the frequency band. |
Vector | Vector | N/A | Vectors are interpreted as frequency knotpoints and values at those knotpoints. |
Functions
FIRLSFilterDesign._to_impulse_response
— MethodFor type III FIR filters:
\[ h = \bigg[ a[M+1] \quad a[M] \quad \cdots \quad a[2] \quad 0 \quad -a[2] \quad \cdots \quad -a[M+1] \bigg]^T\]
FIRLSFilterDesign._to_impulse_response
— MethodFor type II FIR filters:
\[ h = \bigg[ a[M+1] \quad a[M] \quad \cdots \quad a[2] \quad a[1] \quad a[1] \quad a[2] \quad \cdots \quad a[M+1] \bigg]^T\]
FIRLSFilterDesign._to_impulse_response
— MethodFor type IV FIR filters:
\[ h = \bigg[ a[M+1] \quad a[M] \quad \cdots \quad a[2] \quad 0 \quad 0 \quad -a[2] \quad \cdots \quad -a[M+1] \bigg]^T\]
FIRLSFilterDesign._to_impulse_response
— Method_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\]
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}}
...
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}
...
FIRLSFilterDesign.bn!
— Methodbn!(_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.
...
FIRLSFilterDesign.bn!
— Methodbn!(_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 to2/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.
...
FIRLSFilterDesign.bn_n0!
— Methodbn_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
.
FIRLSFilterDesign.constants_b
— Methodconstants_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 inf
. 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 inf
. The first and second columns indicate the weighting at the lower and upper bound of the frequency bands, interpolated linearly in between.
...
FIRLSFilterDesign.firls_design
— Methodfirls_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 inbands_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 inbands_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 returnsa
.
Outputs
h
: a vector of linear-phase FIR filter coefficients.
FIRLSFilterDesign.firls_design
— Methodfirls_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 inbands_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 returnsa
.
Outputs
h
: a vector of linear-phase FIR filter coefficients.
FIRLSFilterDesign.firls_design
— Methodfirls_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 inbands_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 inbands_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 returnsa
.
Outputs
h
: a vector of linear-phase FIR filter coefficients.
FIRLSFilterDesign.firls_design
— Methodfirls_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 inknotpoints_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 returnsa
.
Outputs
h
: a vector of linear-phase FIR filter coefficients.
FIRLSFilterDesign.firls_design
— Methodfirls_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 inknotpoints_DW
.W::Vector
: a vector of size(N,)
which contains weighting function values for the frequency knotpoints inknotpoints_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 returnsa
.
Outputs
h
: a vector of linear-phase FIR filter coefficients.
FIRLSFilterDesign.get_Q
— Methodget_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 inf
. 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$.
FIRLSFilterDesign.get_Q
— Methodget_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$.
FIRLSFilterDesign.get_Q
— Methodget_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$.
FIRLSFilterDesign.get_b
— Methodget_b(M, f, D, W, fir_type)
Finds the vector $b$ used in the equation $Qa = b$. ...
Arguments
M::Integer
: Size of theb
vector isM+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 inf
. 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 inf
. 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$.
...
FIRLSFilterDesign.get_flength_M
— Methodget_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 coefficientsM::Integer
: the number of unique amplitude coefficients needed to form the filter is equal to $M+1$
FIRLSFilterDesign.get_q
— Methodget_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 inf
. 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.
...
FIRLSFilterDesign.idx2n_b
— Methodidx2n_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
...
FIRLSFilterDesign.idx2n_b
— Methodidx2n_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
...
FIRLSFilterDesign.infer_fir_type
— Methodinfer_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:
FIR_I
, when filter length is odd and the filter is not antisymmetricFIR_II
, when filter length is even and the filter is not antisymmetricFIR_III
, when filter length is odd and the filter is antisymmetricFIR_IV
, when filter length is even and the filter is antisymmetric
Arguments
is_odd::Bool
is_antisymmetric::Bool
Outputs
fir_type::FIR