#include <filter.h>

Public Member Functions | |
| FIR_filt () | |
| ~FIR_filt () | |
| void | clear (void) |
| void | set_coeffs (const arma::Mat< T2 > &_b) |
| void | set_coeffs (const arma::Col< T2 > &_b_col) |
| arma::Col< T2 > | get_coeffs () |
| void | update_coeffs (const arma::Mat< T2 > &_b) |
| T3 | operator() (const T1 &in) |
| arma::Mat< T3 > | filter (const arma::Mat< T1 > &in) |
| arma::Col< T3 > | filter (const arma::Col< T1 > &in) |
| void | setup_lms (const arma::uword _N, const double _mu, const arma::uword _L=1) |
| void | lms_adapt (const T3 _err) |
| void | setup_nlms (const arma::uword _N, const double _mu, const T2 _c, const arma::uword _L=1) |
| void | nlms_adapt (const T3 _err) |
| void | setup_newt (const arma::uword _N, const double _mu, const T2 _c, const arma::uword _L=1) |
| void | newt_adapt (const T3 _err) |
| void | setup_rls (const arma::uword _N, const double _lmd, const double _P0) |
| void | rls_adapt (const T3 _err) |
| void | setup_kalman (const arma::uword _N, const double _P0, const double _Q0, const double _R0) |
| void | kalman_adapt (const T3 _err) |
| double | get_step_size (void) |
| arma::Mat< T1 > | get_P (void) |
| arma::Mat< T1 > | get_K (void) |
| void | set_step_size (const double _mu) |
| void | adapt_enable (void) |
| void | adapt_disble (void) |
Private Attributes | |
| arma::uword | M |
| arma::uword | cur_p |
| arma::Mat< T1 > | buf |
| arma::Mat< T2 > | b |
| double | mu |
| arma::uword | L |
| arma::uword | blk_ctr |
| T2 | c |
| arma::Mat< T1 > | P |
| arma::Mat< T1 > | Q |
| arma::Mat< T1 > | R |
| arma::Mat< T1 > | K |
| double | lmd |
| arma::Mat< T1 > | X_tpz |
| arma::uword | do_adapt |
Implements FIR/MA filter functions as
\[ y(n) = \sum_{k=0}^{M-1}{b_kx(n-k)}=b_0x(n)+b_1x(n-1)+...+b_{M-1}x(n-(M-1))\]
where M is the number of taps in the FIR filter. The filter order is M-1. Adaptive update of filter is possible with LMS or NLMS algorithms
|
inline |
|
inline |
|
inline |
Stop adapt.
Definition at line 489 of file filter.h.
References sp::FIR_filt< T1, T2, T3 >::do_adapt.
|
inline |
Start adapt.
Definition at line 481 of file filter.h.
References sp::FIR_filt< T1, T2, T3 >::do_adapt.
|
inline |
Clears the internal states and pointer.
Definition at line 59 of file filter.h.
References sp::FIR_filt< T1, T2, T3 >::buf, and sp::FIR_filt< T1, T2, T3 >::cur_p.
Referenced by sp::FIR_filt< T1, T2, T3 >::set_coeffs().
|
inline |
Definition at line 147 of file filter.h.
References sp::FIR_filt< T1, T2, T3 >::filter().

|
inline |
Filter function.
| in | Input vector |
Definition at line 138 of file filter.h.
Referenced by sp::FIR_filt< T1, T2, T3 >::filter().
|
inline |
Get coefficients from FIR filter.
Definition at line 94 of file filter.h.
References sp::FIR_filt< T1, T2, T3 >::b.
Referenced by main().
|
inline |
|
inline |
|
inline |
Get step size.
Definition at line 446 of file filter.h.
References sp::FIR_filt< T1, T2, T3 >::mu.
|
inline |
Kalman Filter update function.
The Kalman filter is updated as
\( \mathbf{b(n)} = \mathbf{b(n-1)}+\mathbf{Kx(n)}\)
where
\( \mathbf{K} =\frac{\mathbf{Px}}{R+\mathbf{x^TPx}} \)
and
\( \mathbf{P^+} =\mathbf{P^-+xP^-x^T }+Q \)
| _err | Feedback error |
Definition at line 417 of file filter.h.
References sp::FIR_filt< T1, T2, T3 >::b, sp::FIR_filt< T1, T2, T3 >::buf, sp::FIR_filt< T1, T2, T3 >::cur_p, sp::FIR_filt< T1, T2, T3 >::do_adapt, sp::FIR_filt< T1, T2, T3 >::K, sp::FIR_filt< T1, T2, T3 >::M, sp::FIR_filt< T1, T2, T3 >::P, sp::FIR_filt< T1, T2, T3 >::Q, and sp::FIR_filt< T1, T2, T3 >::R.
|
inline |
LMS Filter update function.
The LMS filter is updated as
\( \mathbf{b(n)} = \mathbf{b(n-1)}+2\mu\mathbf{x(n)}err(n) \)
where
\( err(n) = d(n)-\mathbf{b(n-1)^Tx(n)} \)
| _err | Feedback error |
Definition at line 183 of file filter.h.
References sp::FIR_filt< T1, T2, T3 >::b, sp::FIR_filt< T1, T2, T3 >::blk_ctr, sp::FIR_filt< T1, T2, T3 >::buf, sp::FIR_filt< T1, T2, T3 >::cur_p, sp::FIR_filt< T1, T2, T3 >::do_adapt, sp::FIR_filt< T1, T2, T3 >::K, sp::FIR_filt< T1, T2, T3 >::L, sp::FIR_filt< T1, T2, T3 >::M, and sp::FIR_filt< T1, T2, T3 >::mu.
|
inline |
LMS-Newton Filter update function.
The LMS-Newton filter is updated as
\( \mathbf{b(n)} =
\mathbf{b(n-1)}+2\mu\frac{\mathbf{x(n)}err(n)}{c+\mathbf{R_{xx}}} \)
where
\( err(n) = d(n)-\mathbf{b(n-1)^Tx(n)} \)
and \(
\mathbf{R_{xx}} \) is the correlation matrix
| _err | Feedback error |
Definition at line 299 of file filter.h.
References sp::FIR_filt< T1, T2, T3 >::b, sp::FIR_filt< T1, T2, T3 >::blk_ctr, sp::FIR_filt< T1, T2, T3 >::buf, sp::FIR_filt< T1, T2, T3 >::c, sp::FIR_filt< T1, T2, T3 >::cur_p, sp::FIR_filt< T1, T2, T3 >::do_adapt, sp::FIR_filt< T1, T2, T3 >::K, sp::FIR_filt< T1, T2, T3 >::L, sp::FIR_filt< T1, T2, T3 >::M, sp::FIR_filt< T1, T2, T3 >::mu, and sp::FIR_filt< T1, T2, T3 >::X_tpz.
|
inline |
NLMS Filter update function.
The NLMS filter is updated as
\( \mathbf{b(n)} =
\mathbf{b(n-1)}+2\mu\frac{\mathbf{x(n)}err(n)}{c+\mathbf{x(n)^Tx(n)}} \)
where
\( err(n) = d(n)-\mathbf{b(n-1)^Tx(n)} \)
| _err | Feedback error |
Definition at line 239 of file filter.h.
References sp::FIR_filt< T1, T2, T3 >::b, sp::FIR_filt< T1, T2, T3 >::blk_ctr, sp::FIR_filt< T1, T2, T3 >::buf, sp::FIR_filt< T1, T2, T3 >::c, sp::FIR_filt< T1, T2, T3 >::cur_p, sp::FIR_filt< T1, T2, T3 >::do_adapt, sp::FIR_filt< T1, T2, T3 >::K, sp::FIR_filt< T1, T2, T3 >::L, sp::FIR_filt< T1, T2, T3 >::M, and sp::FIR_filt< T1, T2, T3 >::mu.
Referenced by main().
|
inline |
Filter operator.
| in | Input sample |
Definition at line 114 of file filter.h.
References sp::FIR_filt< T1, T2, T3 >::b, sp::FIR_filt< T1, T2, T3 >::buf, sp::FIR_filt< T1, T2, T3 >::cur_p, and sp::FIR_filt< T1, T2, T3 >::M.
|
inline |
RLS Filter update function.
The RLS filter is updated as
\( \mathbf{b(n)} = \mathbf{b(n-1)}+\mathbf{Kx(n)}\)
where
\( \mathbf{K} =\frac{\mathbf{Px}}{\lambda+\mathbf{x^TPx}} \)
and
\( \mathbf{P^+} =\frac{\mathbf{P^-+xP^-x^T }}{\lambda} \)
| _err | Feedback error |
Definition at line 360 of file filter.h.
References sp::FIR_filt< T1, T2, T3 >::b, sp::FIR_filt< T1, T2, T3 >::buf, sp::FIR_filt< T1, T2, T3 >::cur_p, sp::FIR_filt< T1, T2, T3 >::do_adapt, sp::FIR_filt< T1, T2, T3 >::K, sp::FIR_filt< T1, T2, T3 >::lmd, sp::FIR_filt< T1, T2, T3 >::M, and sp::FIR_filt< T1, T2, T3 >::P.
|
inline |
Sets coefficients in FIR filter (col format) The internal state and pointers are cleared.
| _b_col | Filter coefficients \( [b_0 ..b_{M-1}]^T \) |
Definition at line 84 of file filter.h.
References sp::FIR_filt< T1, T2, T3 >::set_coeffs().

|
inline |
Sets coefficients in FIR filter. The internal state and pointers are cleared.
| _b | Filter coefficients \( [b_0 ..b_{M-1}]^T \) |
Definition at line 70 of file filter.h.
References sp::FIR_filt< T1, T2, T3 >::b, sp::FIR_filt< T1, T2, T3 >::buf, sp::FIR_filt< T1, T2, T3 >::clear(), and sp::FIR_filt< T1, T2, T3 >::M.
Referenced by main(), sp::resampling< T1 >::resampling(), sp::resampling< T1 >::resampling(), and sp::FIR_filt< T1, T2, T3 >::set_coeffs().

|
inline |
Set step size.
| _mu | Step size mu |
Definition at line 473 of file filter.h.
References sp::FIR_filt< T1, T2, T3 >::mu.
|
inline |
Kalman Filter function setup.
| _N | Number of filter taps |
| _P0 | Inverse corr matrix initializer |
| _Q0 | Process noise matrix initializer |
| _R0 | Measurement noise matrix initializer |
Definition at line 388 of file filter.h.
References sp::FIR_filt< T1, T2, T3 >::b, sp::FIR_filt< T1, T2, T3 >::buf, sp::FIR_filt< T1, T2, T3 >::cur_p, sp::FIR_filt< T1, T2, T3 >::do_adapt, sp::FIR_filt< T1, T2, T3 >::K, sp::FIR_filt< T1, T2, T3 >::L, sp::FIR_filt< T1, T2, T3 >::M, sp::FIR_filt< T1, T2, T3 >::P, sp::FIR_filt< T1, T2, T3 >::Q, and sp::FIR_filt< T1, T2, T3 >::R.
|
inline |
LMS Filter function setup.
| _N | Number of filter taps |
| _mu | Step size |
| _L | Block length |
Definition at line 159 of file filter.h.
References sp::FIR_filt< T1, T2, T3 >::b, sp::FIR_filt< T1, T2, T3 >::blk_ctr, sp::FIR_filt< T1, T2, T3 >::buf, sp::FIR_filt< T1, T2, T3 >::cur_p, sp::FIR_filt< T1, T2, T3 >::do_adapt, sp::FIR_filt< T1, T2, T3 >::K, sp::FIR_filt< T1, T2, T3 >::L, sp::FIR_filt< T1, T2, T3 >::M, and sp::FIR_filt< T1, T2, T3 >::mu.
|
inline |
LMS-Newton Filter function setup. (Affine Projection Algorithm)
| _N | Number of filter taps |
| _mu | Step size |
| _c | Regularization factor |
| _L | Block length |
Definition at line 270 of file filter.h.
References sp::FIR_filt< T1, T2, T3 >::b, sp::FIR_filt< T1, T2, T3 >::blk_ctr, sp::FIR_filt< T1, T2, T3 >::buf, sp::FIR_filt< T1, T2, T3 >::c, sp::FIR_filt< T1, T2, T3 >::cur_p, sp::FIR_filt< T1, T2, T3 >::do_adapt, sp::FIR_filt< T1, T2, T3 >::K, sp::FIR_filt< T1, T2, T3 >::L, sp::FIR_filt< T1, T2, T3 >::M, sp::FIR_filt< T1, T2, T3 >::mu, and sp::FIR_filt< T1, T2, T3 >::X_tpz.
|
inline |
NLMS Filter function setup.
| _N | Number of filter taps |
| _mu | Step size |
| _c | Regularization factor |
| _L | Block length |
Definition at line 213 of file filter.h.
References sp::FIR_filt< T1, T2, T3 >::b, sp::FIR_filt< T1, T2, T3 >::blk_ctr, sp::FIR_filt< T1, T2, T3 >::buf, sp::FIR_filt< T1, T2, T3 >::c, sp::FIR_filt< T1, T2, T3 >::cur_p, sp::FIR_filt< T1, T2, T3 >::do_adapt, sp::FIR_filt< T1, T2, T3 >::K, sp::FIR_filt< T1, T2, T3 >::L, sp::FIR_filt< T1, T2, T3 >::M, and sp::FIR_filt< T1, T2, T3 >::mu.
Referenced by main().
|
inline |
RLS Filter function setup.
| _N | Number of filter taps |
| _lmd | Lambda |
| _P0 | Inverse corr matrix initializer |
Definition at line 333 of file filter.h.
References sp::FIR_filt< T1, T2, T3 >::b, sp::FIR_filt< T1, T2, T3 >::buf, sp::FIR_filt< T1, T2, T3 >::cur_p, sp::FIR_filt< T1, T2, T3 >::do_adapt, sp::FIR_filt< T1, T2, T3 >::K, sp::FIR_filt< T1, T2, T3 >::L, sp::FIR_filt< T1, T2, T3 >::lmd, sp::FIR_filt< T1, T2, T3 >::M, and sp::FIR_filt< T1, T2, T3 >::P.
|
inline |
Updates coefficients in FIR filter without clearing the internal states.
| _b | Filter coefficients \( [b_0 ..b_{M-1}] \) |
Definition at line 104 of file filter.h.
References sp::FIR_filt< T1, T2, T3 >::b.
Referenced by main().
|
private |
Filter coefficients.
Definition at line 28 of file filter.h.
Referenced by sp::FIR_filt< T1, T2, T3 >::get_coeffs(), sp::FIR_filt< T1, T2, T3 >::kalman_adapt(), sp::FIR_filt< T1, T2, T3 >::lms_adapt(), sp::FIR_filt< T1, T2, T3 >::newt_adapt(), sp::FIR_filt< T1, T2, T3 >::nlms_adapt(), sp::FIR_filt< T1, T2, T3 >::operator()(), sp::FIR_filt< T1, T2, T3 >::rls_adapt(), sp::FIR_filt< T1, T2, T3 >::set_coeffs(), sp::FIR_filt< T1, T2, T3 >::setup_kalman(), sp::FIR_filt< T1, T2, T3 >::setup_lms(), sp::FIR_filt< T1, T2, T3 >::setup_newt(), sp::FIR_filt< T1, T2, T3 >::setup_nlms(), sp::FIR_filt< T1, T2, T3 >::setup_rls(), and sp::FIR_filt< T1, T2, T3 >::update_coeffs().
|
private |
Adaptive filter block length counter.
Definition at line 32 of file filter.h.
Referenced by sp::FIR_filt< T1, T2, T3 >::lms_adapt(), sp::FIR_filt< T1, T2, T3 >::newt_adapt(), sp::FIR_filt< T1, T2, T3 >::nlms_adapt(), sp::FIR_filt< T1, T2, T3 >::setup_lms(), sp::FIR_filt< T1, T2, T3 >::setup_newt(), and sp::FIR_filt< T1, T2, T3 >::setup_nlms().
|
private |
Signal buffer.
Definition at line 27 of file filter.h.
Referenced by sp::FIR_filt< T1, T2, T3 >::clear(), sp::FIR_filt< T1, T2, T3 >::kalman_adapt(), sp::FIR_filt< T1, T2, T3 >::lms_adapt(), sp::FIR_filt< T1, T2, T3 >::newt_adapt(), sp::FIR_filt< T1, T2, T3 >::nlms_adapt(), sp::FIR_filt< T1, T2, T3 >::operator()(), sp::FIR_filt< T1, T2, T3 >::rls_adapt(), sp::FIR_filt< T1, T2, T3 >::set_coeffs(), sp::FIR_filt< T1, T2, T3 >::setup_kalman(), sp::FIR_filt< T1, T2, T3 >::setup_lms(), sp::FIR_filt< T1, T2, T3 >::setup_newt(), sp::FIR_filt< T1, T2, T3 >::setup_nlms(), and sp::FIR_filt< T1, T2, T3 >::setup_rls().
|
private |
Adaptive filter NLMS regulation const.
Definition at line 33 of file filter.h.
Referenced by sp::FIR_filt< T1, T2, T3 >::newt_adapt(), sp::FIR_filt< T1, T2, T3 >::nlms_adapt(), sp::FIR_filt< T1, T2, T3 >::setup_newt(), and sp::FIR_filt< T1, T2, T3 >::setup_nlms().
|
private |
Pointer to current sample in buffer.
Definition at line 26 of file filter.h.
Referenced by sp::FIR_filt< T1, T2, T3 >::clear(), sp::FIR_filt< T1, T2, T3 >::kalman_adapt(), sp::FIR_filt< T1, T2, T3 >::lms_adapt(), sp::FIR_filt< T1, T2, T3 >::newt_adapt(), sp::FIR_filt< T1, T2, T3 >::nlms_adapt(), sp::FIR_filt< T1, T2, T3 >::operator()(), sp::FIR_filt< T1, T2, T3 >::rls_adapt(), sp::FIR_filt< T1, T2, T3 >::setup_kalman(), sp::FIR_filt< T1, T2, T3 >::setup_lms(), sp::FIR_filt< T1, T2, T3 >::setup_newt(), sp::FIR_filt< T1, T2, T3 >::setup_nlms(), and sp::FIR_filt< T1, T2, T3 >::setup_rls().
|
private |
Adaptive filter enable flag.
Definition at line 40 of file filter.h.
Referenced by sp::FIR_filt< T1, T2, T3 >::adapt_disble(), sp::FIR_filt< T1, T2, T3 >::adapt_enable(), sp::FIR_filt< T1, T2, T3 >::kalman_adapt(), sp::FIR_filt< T1, T2, T3 >::lms_adapt(), sp::FIR_filt< T1, T2, T3 >::newt_adapt(), sp::FIR_filt< T1, T2, T3 >::nlms_adapt(), sp::FIR_filt< T1, T2, T3 >::rls_adapt(), sp::FIR_filt< T1, T2, T3 >::setup_kalman(), sp::FIR_filt< T1, T2, T3 >::setup_lms(), sp::FIR_filt< T1, T2, T3 >::setup_newt(), sp::FIR_filt< T1, T2, T3 >::setup_nlms(), and sp::FIR_filt< T1, T2, T3 >::setup_rls().
|
private |
Adaptive filter gain vector.
Definition at line 37 of file filter.h.
Referenced by sp::FIR_filt< T1, T2, T3 >::get_K(), sp::FIR_filt< T1, T2, T3 >::kalman_adapt(), sp::FIR_filt< T1, T2, T3 >::lms_adapt(), sp::FIR_filt< T1, T2, T3 >::newt_adapt(), sp::FIR_filt< T1, T2, T3 >::nlms_adapt(), sp::FIR_filt< T1, T2, T3 >::rls_adapt(), sp::FIR_filt< T1, T2, T3 >::setup_kalman(), sp::FIR_filt< T1, T2, T3 >::setup_lms(), sp::FIR_filt< T1, T2, T3 >::setup_newt(), sp::FIR_filt< T1, T2, T3 >::setup_nlms(), and sp::FIR_filt< T1, T2, T3 >::setup_rls().
|
private |
Adaptive filter block length.
Definition at line 31 of file filter.h.
Referenced by sp::FIR_filt< T1, T2, T3 >::lms_adapt(), sp::FIR_filt< T1, T2, T3 >::newt_adapt(), sp::FIR_filt< T1, T2, T3 >::nlms_adapt(), sp::FIR_filt< T1, T2, T3 >::setup_kalman(), sp::FIR_filt< T1, T2, T3 >::setup_lms(), sp::FIR_filt< T1, T2, T3 >::setup_newt(), sp::FIR_filt< T1, T2, T3 >::setup_nlms(), and sp::FIR_filt< T1, T2, T3 >::setup_rls().
|
private |
Adaptive filter RLS forgetting factor.
Definition at line 38 of file filter.h.
Referenced by sp::FIR_filt< T1, T2, T3 >::rls_adapt(), and sp::FIR_filt< T1, T2, T3 >::setup_rls().
|
private |
Nr of filter taps.
Definition at line 25 of file filter.h.
Referenced by sp::FIR_filt< T1, T2, T3 >::kalman_adapt(), sp::FIR_filt< T1, T2, T3 >::lms_adapt(), sp::FIR_filt< T1, T2, T3 >::newt_adapt(), sp::FIR_filt< T1, T2, T3 >::nlms_adapt(), sp::FIR_filt< T1, T2, T3 >::operator()(), sp::FIR_filt< T1, T2, T3 >::rls_adapt(), sp::FIR_filt< T1, T2, T3 >::set_coeffs(), sp::FIR_filt< T1, T2, T3 >::setup_kalman(), sp::FIR_filt< T1, T2, T3 >::setup_lms(), sp::FIR_filt< T1, T2, T3 >::setup_newt(), sp::FIR_filt< T1, T2, T3 >::setup_nlms(), and sp::FIR_filt< T1, T2, T3 >::setup_rls().
|
private |
Adaptive filter step size.
Definition at line 30 of file filter.h.
Referenced by sp::FIR_filt< T1, T2, T3 >::get_step_size(), sp::FIR_filt< T1, T2, T3 >::lms_adapt(), sp::FIR_filt< T1, T2, T3 >::newt_adapt(), sp::FIR_filt< T1, T2, T3 >::nlms_adapt(), sp::FIR_filt< T1, T2, T3 >::set_step_size(), sp::FIR_filt< T1, T2, T3 >::setup_lms(), sp::FIR_filt< T1, T2, T3 >::setup_newt(), and sp::FIR_filt< T1, T2, T3 >::setup_nlms().
|
private |
Adaptive filter Inverse corr matrix (estimated accuracy)
Definition at line 34 of file filter.h.
Referenced by sp::FIR_filt< T1, T2, T3 >::get_P(), sp::FIR_filt< T1, T2, T3 >::kalman_adapt(), sp::FIR_filt< T1, T2, T3 >::rls_adapt(), sp::FIR_filt< T1, T2, T3 >::setup_kalman(), and sp::FIR_filt< T1, T2, T3 >::setup_rls().
|
private |
Adaptive filter Process noise.
Definition at line 35 of file filter.h.
Referenced by sp::FIR_filt< T1, T2, T3 >::kalman_adapt(), and sp::FIR_filt< T1, T2, T3 >::setup_kalman().
|
private |
Adaptive filter Measurement noise.
Definition at line 36 of file filter.h.
Referenced by sp::FIR_filt< T1, T2, T3 >::kalman_adapt(), and sp::FIR_filt< T1, T2, T3 >::setup_kalman().
|
private |
Adaptive filter Toeplitz for Corr matrix calc.
Definition at line 39 of file filter.h.
Referenced by sp::FIR_filt< T1, T2, T3 >::newt_adapt(), and sp::FIR_filt< T1, T2, T3 >::setup_newt().