Commit
This commit does not belong to any branch on this repository, and may belong to a fork outside of the repository.
Add
firls
function (least-squares minimization based FIR filter des…
…ign) Adds a new `firls` procedure that can be used to design a linear phase FIR filter using least-squares error criterion. This function is avaiable in `scpy` (as `scipy.signal.firls`). This version is similar but has some differences in its interface and also as some additional functionality (it can design a broader set of filter types). `firls` eturns the coefficients of the linear phase FIR filter of the requested order and type which best fullfills the desired frequency response "constraints" described by the `bands` (i.e. frequency), `desired` (i.e. gain) and `weights` (e.i. constraint weights) input tensors. Depending on the order and the value of the `symmetric` argument, the output filter coefficients will correspond to a Type I, II, III or IV FIR filter. The constraints are specified as a pair of tensors describing "constrained frequency bands" indicating the start and end frequencies of each band `k`, as well as the gain of the filter at those edge frequencies. Inputs: - fir_order: The "order" of the filter (which is 1 less than the length of the output tensor). - bands: 2-column matrix of non-overlapping, normalized frequency band edges in ascending order. Each row in the matrix corresponds to the edges of a constrained frequency band (the column contains the start frequency of the band and the second column contains the end frequency). Frequencies between the specified bands are "unconstrained" (i.e. ignored during the error minimization process). Frequency values must be in the [0.0, fs/2] range (i.e. they cannot be negative or exceed the Nyquist frequency). - desired: 2-column matrix that specifies the gains of the desired frequency response at the band "edges". Thus the lengths of the `bands` and `desired` tensors must match. If they do not, an exception is raised. For each band `k` the desired filter frequency response is such that its gain linearly changes from the `desired[k, 0]` value at the start of the band (i.e. at the `bands[k, 0]` frequency) to the value `desired[k, 1]` at the end of the band (i.e. at `bands[k, 1]`). - weights: Optional rank-1 Tensor of weights. Controls which frequency response "contraints" are given more "weight" during the least-squares error minimization process. The default is that all constraints are given the same weight. If provided, its length must be half the length of `bands` (i.e. there must be one constraint per band). An exception is raised otherwise. - symmetric: When `true` (the default), the result will be a symmetric FIR filter (Type I when `fir_order` is even and Type II when `fir_order` is odd). When `false`, the result will be an anti-symmetric FIR filter (Type III when `fir_order` is even and Type IV when `fir_order` is odd). - fs: The sampling frequency of the signal (as a float). Each frequency in `bands` must be between 0.0 and `fs/2` (inclusive). Default is 2.0, which means that by default the band frequencies are expected to be on the 0.0 to 1.0 range. Result: - A Tensor containing the `fir_order + 1` taps of the FIR filter that best approximates the desired frequency response (i.e. the filter that minimizes the least-squares error vs the given constraints). Notes: - Contrary to numpy's firls, the first argument is the FIR filter order, not the FIR length. The filter length is `filter_order + 1`. - Frequencies between "constrained bands" are considered "don't care" regions for which the error is not minimized. - When designing a filter with a gain other than zero at the Nyquist frequency (i.e. when `bands[^1] != 0.0`), such as high-pass and band-stop filters, the filter order must be even. An exception is raised otherwise. - The `bands` and `desired` can also be flat rank-1 tensors, as long as their length is even (so that they can be reshaped into 2 column matrices. Examples: ```nim # Example of an order 4, Type I, low-pass FIR filter which targets being # a pass-through in the 0.0 to 0.3 times Nyquist frequency range, while # filtering frequencies in the 0.4 to 1.0 Nyquist frequency range # Note that the range 0.3 to 0.4 is left unconstrained # Also note how the result length is 6 and the filter is symmetric # around the middle value: echo firls(4, [[0.0, 0.3], [0.4, 1.0]].toTensor, [[1.0, 1.0], [0.0, 0.0]].toTensor) # Tensor[system.float] of shape "[5]" on backend "Cpu" # 0.126496 0.278552 0.345068 0.278552 0.126496 # Same filter as above, but using rank-1, even length tensors as inputs echo firls(4, [0.0, 0.3, 0.4, 1.0].toTensor, [1.0, 1.0, 0.0, 0.0].toTensor) # Tensor[system.float] of shape "[5]" on backend "Cpu" # 0.126496 0.278552 0.345068 0.278552 0.126496 # Same filter as above, but using unnormalized frequencies and a sampling # frequency of 10.0 Hz (note how all the frequencies are 5 times greater # because the default fs is 2.0): echo firls(4, [[0.0, 1.5], [2.0, 5.0]].toTensor, [[1.0, 1.0], [0.0, 0.0]].toTensor, fs = 10.0) # Tensor[system.float] of shape "[5]" on backend "Cpu" # 0.126496 0.278552 0.345068 0.278552 0.126496 # Same kind of filter, but give more weight to the pass-through constraint echo firls(4, [[0.0, 0.3], [0.4, 1.0]].toTensor, [[1.0, 1.0], [0.0, 0.0]].toTensor, weights = [1.0, 0.5].toTensor) # Tensor[system.float] of shape "[5]" on backend "Cpu" # 0.110353 0.284473 0.360868 0.284473 0.110353 # A more complex, order 5, Type II, low-pass filter echo firls(5, [[0.0, 0.3], [0.3, 0.6], [0.6, 1.0]].toTensor, [[1.0, 1.0], [1.0, 0.2], [0.0, 0.0]].toTensor) # Tensor[system.float] of shape "[6]" on backend "Cpu" # -0.0560333 0.146107 0.430716 0.430716 0.146107 -0.0560333 # Example of an order 6 Type IV high-pass FIR filter: echo firls(6, [[0.0, 0.4], [0.6, 1.0]].toTensor, [[0.0, 0.0], [0.9, 1.0]].toTensor, symmetric = false) # Tensor[system.float] of shape "[7]" on backend "Cpu" # -0.13945 0.285186 -0.258596 0 0.258596 -0.285186 0.13945 Trying to design a high-pass filter with odd order generates an exception: echo firls(5, [0.0, 0.5, 0.6, 1.0].toTensor, [0.0, 0.0, 1.0, 1.0].toTensor, symmetric = false) # Filter order (5) must be even when the last # frequency is 1.0 and its gain is not 0.0 [ValueError] ```
- Loading branch information