Skip to content

📝 Add docstrings to improve-numerical-accuracy #52

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Closed
coderabbitai wants to merge 1 commit into from
Closed

📝 Add docstrings to improve-numerical-accuracy #52

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
39 changes: 39 additions & 0 deletions include/constfilt/analog_filter.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -79,13 +79,28 @@ class AnalogFilter : public Filter<T, N + 1u, N + 1u>
"method-tag constructor with constfilt::prewarp(wc_hz)");
}

/**
* @brief Discretizes a continuous-time transfer function.
*
* @param continuous_tf Continuous-time transfer function to convert.
* @param sample_rate_hz Sampling rate in hertz.
* @param method_tag Discretization method tag.
*/
constexpr AnalogFilter(TransferFunction<T, N + 1u, N + 1u> continuous_tf,
T sample_rate_hz, BoundMethod method_tag)
: AnalogFilter(checked_discretize(continuous_tf.b, continuous_tf.a,
sample_rate_hz, method_tag))
{
}

/**
* @brief Discretizes an analog transfer function using factored pole and zero data.
*
* @param continuous_tf Continuous-time transfer function to discretize.
* @param factored_tf Factored pole and zero representation used during discretization.
* @param sample_rate_hz Sampling rate in hertz.
* @param method_tag Discretization method tag.
*/
constexpr AnalogFilter(TransferFunction<T, N + 1u, N + 1u> continuous_tf,
const FactoredTF<T, N> &factored_tf,
T sample_rate_hz, BoundMethod method_tag)
Expand All @@ -96,12 +111,26 @@ class AnalogFilter : public Filter<T, N + 1u, N + 1u>
}

private:
/**
* @brief Initializes an analog filter from a digital transfer function.
*
* @param digital_tf Digital transfer function to store in the base filter.
*/
constexpr explicit AnalogFilter(
TransferFunction<T, N + 1u, N + 1u> digital_tf)
: Filter<T, N + 1u, N + 1u>(digital_tf.b, digital_tf.a)
{
}

/**
* @brief Discretizes an analog transfer function.
*
* @param b_c Continuous-time numerator coefficients.
* @param a_c Continuous-time denominator coefficients.
* @param sample_rate_hz Sampling rate in hertz.
* @param method_tag Discretization method tag.
* @return TransferFunction<T, N + 1u, N + 1u> Discrete transfer function.
*/
static constexpr TransferFunction<T, N + 1u, N + 1u> checked_discretize(
const T (&b_c)[N + 1u], const T (&a_c)[N + 1u], T sample_rate_hz,
BoundMethod method_tag)
Expand Down Expand Up @@ -137,6 +166,16 @@ class AnalogFilter : public Filter<T, N + 1u, N + 1u>
b_c, a_c, static_cast<T>(1) / sample_rate_hz, method_tag);
}

/**
* @brief Discretizes an analog transfer function using factored pole-zero data.
*
* @param b_c Continuous-time numerator coefficients.
* @param a_c Continuous-time denominator coefficients.
* @param factored_tf Factored pole-zero representation of the continuous transfer function.
* @param sample_rate_hz Sampling rate in hertz.
* @param method_tag Discretization method tag.
* @return TransferFunction<T, N + 1u, N + 1u> Discrete transfer function coefficients.
*/
static constexpr TransferFunction<T, N + 1u, N + 1u>
checked_discretize_factored(const T (&b_c)[N + 1u], const T (&a_c)[N + 1u],
const FactoredTF<T, N> &factored_tf,
Expand Down
31 changes: 27 additions & 4 deletions include/constfilt/butterworth.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -31,7 +31,12 @@ class Butterworth
using BoundMethod = typename bind_method<T, Method>::type;

public:
// Construct from filter specification; all math is constexpr.
/**
* @brief Constructs a Butterworth filter from its cutoff and sample rate.
*
* @param cutoff_hz Cutoff frequency in hertz.
* @param sample_rate_hz Sampling frequency in hertz.
*/
constexpr Butterworth(T cutoff_hz, T sample_rate_hz)
: AnalogFilter<T, N, BoundMethod>(
compute_continuous_tf(cutoff_hz),
Expand Down Expand Up @@ -78,7 +83,13 @@ class Butterworth
//
// Denominator: a[k] = p[N-k] * wc^k. The normalized Butterworth
// coefficients appear in reversed order, each scaled by wc^k (p[0] = 1
// so a[0] = 1, keeping the denominator monic).
/**
* @brief Computes the continuous-time transfer function for a high-pass Butterworth filter.
*
* @param wc Cutoff angular frequency in radians per second.
* @param b Numerator coefficients in ascending power order.
* @param a Denominator coefficients in ascending power order.
*/
static constexpr void continuous_tf(T wc, T (&b)[N + 1u], T (&a)[N + 1u],
HighPass)
{
Expand All @@ -92,7 +103,13 @@ class Butterworth
}
}

// LP: poles at wc*exp(j*theta_k), no finite zeros, gain = wc^N.
/**
* @brief Computes the factored transfer function for a low-pass Butterworth filter.
*
* @param cutoff_hz Cutoff frequency in hertz.
* @param LowPass Low-pass filter selector.
* @return FactoredTF<T, N> Factored transfer function with no finite zeros, gain set to the cutoff frequency raised to the filter order, and poles on the Butterworth circle.
*/
static constexpr FactoredTF<T, N> compute_factored_tf(T cutoff_hz, LowPass)
{
const T wc = static_cast<T>(2) * static_cast<T>(GCEM_PI) * cutoff_hz;
Expand All @@ -110,7 +127,13 @@ class Butterworth
return factored_tf;
}

// HP: poles at wc*exp(-j*theta_k) (magnitude wc), N zeros at s=0, gain=1.
/**
* @brief Computes the factored transfer function for a high-pass Butterworth filter.
*
* @param cutoff_hz Cutoff frequency in hertz.
* @param HighPass High-pass filter tag.
* @return FactoredTF<T, N> Factored transfer function with N zeros at the origin, unit gain, and poles on the Butterworth pole circle.
*/
static constexpr FactoredTF<T, N> compute_factored_tf(T cutoff_hz, HighPass)
{
using Complex = consteig::Complex<T>;
Expand Down
93 changes: 93 additions & 0 deletions include/constfilt/discretize.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -123,6 +123,13 @@ template <typename T, consteig::Size N> struct FactoredTF
// exponential factors, compute V * diag(exp_factors) * V^{-1} and return
// the real part.
template <typename T, consteig::Size N>
/**
* @brief Forms a real matrix exponential from an eigenvector basis and per-mode exponential factors.
*
* @param V Complex eigenvector matrix.
* @param exp_factors Exponential value for each eigenmode.
* @return Real part of V * diag(exp_factors) * V^{-1}.
*/
constexpr consteig::Matrix<T, N, N> spectral_matrix_exp(
const consteig::Matrix<consteig::Complex<T>, N, N> &V,
const consteig::Complex<T> (&exp_factors)[N])
Expand Down Expand Up @@ -173,6 +180,11 @@ constexpr consteig::Matrix<T, N, N> spectral_matrix_exp(
// where V = eigenvectors, lam_i = eigenvalues (complex).
// Real part extracted at the end (imaginary parts cancel for real A).
template <typename T, consteig::Size N>
/**
* @brief Computes the matrix exponential of a square matrix.
*
* @return consteig::Matrix<T, N, N> The matrix exponential of @p A.
*/
constexpr consteig::Matrix<T, N, N> matrix_exp(
const consteig::Matrix<T, N, N> &A)
{
Expand Down Expand Up @@ -205,6 +217,13 @@ constexpr consteig::Matrix<T, N, N> matrix_exp(
// V[r][i] = p_i^r (Vandermonde, built from UNSCALED poles)
// exp factor = exp(Ts * p_i) (spectral scaling)
template <typename T, consteig::Size N>
/**
* @brief Computes the state transition matrix for a controllable canonical form system.
*
* @param factored_tf Analog poles used to form the canonical eigenvector basis.
* @param Ts Sampling period.
* @return Matrix exponential of the continuous-time controllable canonical form state matrix.
*/
constexpr consteig::Matrix<T, N, N> matrix_exp_ccf(
const FactoredTF<T, N> &factored_tf, T Ts)
{
Expand Down Expand Up @@ -241,6 +260,13 @@ constexpr consteig::Matrix<T, N, N> matrix_exp_ccf(
// discrete matrices for callers that need them (e.g. state estimation, observer
// design).
template <typename T, consteig::Size N>
/**
* @brief Discretizes a continuous-time state-space model using zero-order hold.
*
* @param sys_c Continuous-time system.
* @param Ts Sampling period.
* @return State-space model of the zero-order-hold discretization.
*/
constexpr StateSpace<T, N> zoh_discretize(const StateSpace<T, N> &sys_c, T Ts,
ZOH /*tag*/)
{
Expand Down Expand Up @@ -272,6 +298,14 @@ constexpr StateSpace<T, N> zoh_discretize(const StateSpace<T, N> &sys_c, T Ts,
// Bypasses both the QR eigenvalue search and eigenvector inverse iteration
// by exploiting the Vandermonde structure of the CCF eigenvectors.
template <typename T, consteig::Size N>
/**
* @brief Discretizes a continuous-time state-space system with ZOH using factored poles.
*
* @param sys_c Continuous-time system.
* @param Ts Sampling period.
* @param factored_tf Analog pole data used to compute the state transition matrix.
* @return State-space system with discrete-time A and B matrices and unchanged C and D.
*/
constexpr StateSpace<T, N> zoh_discretize_with_evals(
const StateSpace<T, N> &sys_c, T Ts, const FactoredTF<T, N> &factored_tf)
{
Expand Down Expand Up @@ -456,6 +490,16 @@ constexpr StateSpace<T, N> tf_to_ss(const T (&b)[N + 1u], const T (&a)[N + 1u])
// and sample period, maps to z-domain, pads missing zeros at z=-1, matches
// gain, and assembles the discrete TF.
template <typename T, consteig::Size N>
/**
* @brief Assembles a matched-Z discrete transfer function from analog poles, zeros, and gain.
*
* @param poles Analog poles to map to the z-plane.
* @param zeros Analog finite zeros to map to the z-plane.
* @param nz Number of valid entries in @p zeros.
* @param k_c Analog gain used for gain matching.
* @param Ts Sampling period.
* @return TransferFunction<T, N + 1u, N + 1u> Discrete transfer function with matched poles, zeros, and gain.
*/
constexpr TransferFunction<T, N + 1u, N + 1u> matched_z_assemble(
const consteig::Complex<T> (&poles)[N],
const consteig::Complex<T> (&zeros)[N], consteig::Size nz, T k_c, T Ts)
Expand Down Expand Up @@ -604,6 +648,15 @@ constexpr TransferFunction<T, N + 1u, N + 1u> matched_z_assemble(
// Steps 1-3: extract poles and zeros from polynomial coefficients via companion
// matrix eigenvalues, then delegate to matched_z_assemble.
template <typename T, consteig::Size N>
/**
* @brief Discretizes an analog transfer function using the matched-Z method.
*
* @param b_c Continuous-time numerator coefficients.
* @param a_c Continuous-time denominator coefficients.
* @param Ts Sampling period.
* @param tag Matched-Z method tag.
* @return TransferFunction<T, N + 1u, N + 1u> Discrete-time transfer function.
*/
constexpr TransferFunction<T, N + 1u, N + 1u> matched_z_discretize_tf(
const T (&b_c)[N + 1u], const T (&a_c)[N + 1u], T Ts, MatchedZ /*tag*/)
{
Expand Down Expand Up @@ -695,6 +748,13 @@ constexpr TransferFunction<T, N + 1u, N + 1u> matched_z_discretize_tf(
// Matched-Z discretization from a FactoredTF: poles and zeros are known
// analytically, so companion-matrix eigendecompositions are not needed.
template <typename T, consteig::Size N>
/**
* @brief Discretizes a factored transfer function using matched-Z mapping.
*
* @param factored_tf Analog poles, zeros, gain, and zero count.
* @param Ts Sampling period.
* @return TransferFunction<T, N + 1u, N + 1u> Discrete transfer function coefficients.
*/
constexpr TransferFunction<T, N + 1u, N + 1u> matched_z_discretize_factored(
const FactoredTF<T, N> &factored_tf, T Ts)
{
Expand Down Expand Up @@ -819,6 +879,12 @@ constexpr TransferFunction<T, N + 1u, N + 1u> analog_to_digital(
}

template <typename T, consteig::Size N>
/**
* @brief Discretizes a continuous-time state-space system using prewarped Tustin and converts it to a transfer function.
*
* @param tag Prewarp data containing the warp frequency.
* @return Discrete-time transfer function obtained from the prewarped bilinear transform.
*/
constexpr TransferFunction<T, N + 1u, N + 1u> analog_to_digital(
const StateSpace<T, N> &sys_c, T Ts, TustinPWData<T> tag)
{
Expand All @@ -830,6 +896,20 @@ constexpr TransferFunction<T, N + 1u, N + 1u> analog_to_digital(
// path.

template <typename T, consteig::Size N>
/**
* @brief Discretizes a transfer function with ZOH using factored poles and zeros.
*
* Converts the analog transfer function to controllable canonical form, applies
* zero-order-hold discretization with the stored pole factors, and converts the
* result back to transfer-function form.
*
* @param b_c Continuous-time numerator coefficients.
* @param a_c Continuous-time denominator coefficients.
* @param factored_tf Analog pole/zero/gain data used to accelerate discretization.
* @param Ts Sampling interval.
*
* @return Discrete-time transfer function.
*/
constexpr TransferFunction<T, N + 1u, N + 1u> discretize_with_factored(
const T (&b_c)[N + 1u], const T (&a_c)[N + 1u],
const FactoredTF<T, N> &factored_tf, T Ts, ZOH)
Expand All @@ -839,6 +919,13 @@ constexpr TransferFunction<T, N + 1u, N + 1u> discretize_with_factored(
}

template <typename T, consteig::Size N>
/**
* @brief Discretizes a factored analog transfer function with matched-z mapping.
*
* @param factored_tf Analog poles, zeros, and gain.
* @param Ts Sampling period.
* @return TransferFunction<T, N + 1u, N + 1u> Discrete transfer function.
*/
constexpr TransferFunction<T, N + 1u, N + 1u> discretize_with_factored(
const T (& /*b_c*/)[N + 1u], const T (& /*a_c*/)[N + 1u],
const FactoredTF<T, N> &factored_tf, T Ts, MatchedZ)
Expand All @@ -847,6 +934,12 @@ constexpr TransferFunction<T, N + 1u, N + 1u> discretize_with_factored(
}

template <typename T, consteig::Size N, typename M>
/**
* @brief Discretizes a continuous transfer function using the selected method.
*
* @param method_tag Discretization method tag.
* @return TransferFunction<T, N + 1u, N + 1u> Discrete transfer function coefficients.
*/
constexpr TransferFunction<T, N + 1u, N + 1u> discretize_with_factored(
const T (&b_c)[N + 1u], const T (&a_c)[N + 1u],
const FactoredTF<T, N> & /*factored_tf*/, T Ts, M method_tag)
Expand Down
Loading
Loading