Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
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
4 changes: 2 additions & 2 deletions README.md
Comment thread
gabrielfrasantos marked this conversation as resolved.
Comment thread
gabrielfrasantos marked this conversation as resolved.
Original file line number Diff line number Diff line change
Expand Up @@ -16,7 +16,7 @@ Refer to the documentation to quickly integrate and utilize the library's signal

| Category | Description |
|--------------------------------------------------------------------|----------------------------------------------------------------------|
| [Analysis](doc/analysis/README.md) | FFT, Real-Input FFT (RFFT), Power Spectral Density, DCT, Window Functions, Signal Detectors, Convolution & Correlation, Goertzel Algorithm |
| [Analysis](doc/analysis/README.md) | FFT, Real-Input FFT (RFFT), Power Spectral Density, DCT, Window Functions, Signal Detectors, Convolution & Correlation, Goertzel Algorithm, Decibels |
| [Control Analysis](doc/control_analysis/README.md) | Frequency Response, Root Locus, Controllability/Observability Matrices & Gramians |
| [Controllers](doc/controllers/README.md) | Bang-Bang/Hysteresis, PID, LQR, LQI (Integral/Servo State Feedback), MPC, Saturation, Rate Limiter, Slew-Limited Saturation, Feedforward/2-DOF, Gain-Scheduled Controller, Lead-Lag Compensator, Luenberger Observer |
| [Dynamics](doc/dynamics/README.md) | Euler-Lagrange, Newton-Euler, Recursive Newton-Euler, ABA |
Expand All @@ -26,7 +26,7 @@ Refer to the documentation to quickly integrate and utilize the library's signal
| [Neural Network](doc/neural_network/README.md) | Layers, activations, losses, model |
| [Optimization](doc/optimization/README.md) | Gradient Descent |
| [Regularization](doc/regularization/README.md) | L1 (Lasso), L2 (Ridge) |
| [Math](doc/math/README.md) | CORDIC, Quaternion |
| [Math](doc/math/README.md) | CORDIC, Quaternion, Step Response Metrics |
| [Solvers](doc/solvers/README.md) | Gaussian Elimination, Levinson-Durbin, Durand-Kerner, Cholesky, DARE, Runge-Kutta ODE Integrators (RK4 + Dormand-Prince) |
| [Performance Optimization](doc/performance-optimization/README.md) | Compiler optimizations, SIMD |

Expand Down
19 changes: 0 additions & 19 deletions ROADMAP.md
Original file line number Diff line number Diff line change
Expand Up @@ -48,8 +48,6 @@ Difficulty legend:
| 45 | IIR filter design (Butterworth/Chebyshev + bilinear) | `filters/passive` | ★★★★★ |
| 46 | H∞ state-feedback control | `robust_control` (new) | ★★★★★ |
| 47 | Model Reference Adaptive Control (MRAC) | `nonlinear_control` (new) | ★★★★★ |
| 48 | Decibel & magnitude-response helpers | `analysis` | ★☆☆☆☆ |
| 49 | Step / transient-response metrics | `math` | ★★☆☆☆ |
| 50 | Matrix norms & condition number | `math` | ★★★☆☆ |
| 51 | Spectral radius / discrete stability margin | `math` | ★★★☆☆ |
| 52 | Estimator consistency metrics (NEES / NIS) | `estimators` | ★★★☆☆ |
Expand Down Expand Up @@ -561,23 +559,6 @@ and [`control_analysis/FrequencyResponse`](numerical/control_analysis/FrequencyR
magnitude/phase. The items below are the missing pieces. All are **float-only**, no-heap, and operate
on bounded `math::Vector`/`math::Matrix` inputs; tests are `TEST_F` on `float`.

### 48. Decibel & magnitude-response helpers ★☆☆☆☆ — `analysis`
- **What:** `ToDecibels(ratio)` / `FromDecibels`, plus convenience magnitude(-in-dB) and attenuation
helpers layered over `control_analysis::FrequencyResponse`.
- **Metric value:** M2 (frequency response) — pass-band ripple / stop-band attenuation in dB, the
natural unit for filter and controller tests.
- **Algorithm:** `20·log10(·)`; guard the zero/`-inf` case with a floor.
- **Reuses:** `math::TrigonometricFunctions`/`std::log10`, existing `FrequencyResponse`.

### 49. Step / transient-response metrics ★★☆☆☆ — `math`
- **What:** From a bounded step-response `Vector`: `RiseTime` (10–90 %), `SettlingTime` (±band),
`PercentOvershoot`, `PeakTime`, `SteadyStateError`.
- **Metric value:** M3 (time-response) — the core acceptance criteria for every controller
(`controllers/`), `LinearTimeInvariant`, and IIR filter.
- **Algorithm:** single forward pass over the sampled response against the reference/steady value;
standard control-systems definitions.
- **Reuses:** `math::Vector`, `math::Statistics` for the steady-state estimate.

### 50. Matrix norms & condition number ★★★☆☆ — `math`
- **What:** `FrobeniusNorm`, `OneNorm`, `InfinityNorm` on `Matrix`; `Vector` `Norm`/`Normalize`;
`ConditionNumber` estimate.
Expand Down
98 changes: 98 additions & 0 deletions doc/analysis/Decibels.md
Original file line number Diff line number Diff line change
@@ -0,0 +1,98 @@
# Decibels & Magnitude-Response Helpers

## Overview & Motivation

Audio, filter, and control-system specifications express signal levels and filter performance in decibels (dB) because the human auditory system and most engineering metrics scale logarithmically with amplitude ratio. Specifying a stop-band attenuation of 60 dB or a pass-band ripple of 0.1 dB is natural and compact; the equivalent linear ratios (1 000 : 1 and 1.01161 : 1) are not. A small set of conversion primitives — `ToDecibels`, `FromDecibels`, and two derived helpers for attenuation and ripple — centralises this conversion and eliminates scattered, error-prone inline `20·log10` expressions throughout the rest of the library.

## Mathematical Theory

### Magnitude Decibel Conversion

For a positive amplitude ratio $r > 0$, the equivalent level in decibels is:

$$L_{\mathrm{dB}} = 20 \log_{10}(r)$$

The factor 20 (rather than 10) reflects the voltage/pressure convention: power is proportional to the square of amplitude, so a doubling of amplitude ($r = 2$) gives a 6 dB increase, matching the $10 \log_{10}(4) = 6.02$ dB power equivalent.

### Inverse Conversion

$$r = 10^{L_{\mathrm{dB}}/20}$$

This inverse is exact for all finite $L_{\mathrm{dB}}$; no guard is needed on the output side.

### Zero and Negative Input Guard

$\log_{10}(0) = -\infty$; negative ratios are physically meaningless. Both cases are mapped to a finite floor value $L_{\min}$ chosen well below any engineering specification of interest:

$$L_{\mathrm{dB}} = \max\!\left(20\log_{10}(r),\; L_{\min}\right), \quad r > 0$$
$$L_{\mathrm{dB}} = L_{\min}, \quad r \leq 0$$

A floor of $-160\,\mathrm{dB}$ corresponds to an amplitude ratio below $10^{-8}$, safely beyond the dynamic range of any practical floating-point computation in 32-bit single precision.

### Derived Helpers

**Stop-band attenuation** between a pass-band ratio $r_p$ and a stop-band ratio $r_s$:

$$A = L_{\mathrm{dB}}(r_p) - L_{\mathrm{dB}}(r_s)$$

**Pass-band ripple** between the maximum and minimum in-band ratios $r_{\max}$ and $r_{\min}$:

$$\Delta = L_{\mathrm{dB}}(r_{\max}) - L_{\mathrm{dB}}(r_{\min})$$

Both are simple differences in decibel space, exploiting the logarithm identity $\log(a/b) = \log a - \log b$.

## Complexity Analysis

| Operation | Time | Space | Notes |
|-----------------|------|-------|-------------------------------------|
| `ToDecibels` | O(1) | O(1) | One `log10` + one `max` + one `mul` |
| `FromDecibels` | O(1) | O(1) | One `pow` |
| `AttenuationDb` | O(1) | O(1) | Two `ToDecibels` + one subtraction |
| `RippleDb` | O(1) | O(1) | Two `ToDecibels` + one subtraction |

No state, no buffers. All operations are pure functions.

## Step-by-Step Walkthrough

Converting a ratio of 10 to decibels:

1. Input $r = 10$; guard passes ($r > 0$).
2. Compute $20 \cdot \log_{10}(10) = 20 \cdot 1 = 20$.
3. Apply floor: $\max(20, -160) = 20$.
4. Output: $20\,\mathrm{dB}$.

Round-trip for $r = 0.5$:

1. `ToDecibels(0.5)` = $20 \cdot \log_{10}(0.5) \approx -6.0206\,\mathrm{dB}$.
2. `FromDecibels(-6.0206)` = $10^{-6.0206/20} \approx 0.5$.

## Pitfalls & Edge Cases

- Passing $r = 0$ produces $-\infty$ from `log10`; the floor guard prevents propagation into downstream computations.
- Negative ratios indicate a programming error (signed sample values must not be passed directly as ratios without taking absolute value first); they are silently floored rather than raising an exception, consistent with the no-exception policy.
- With `fast-math` enabled, the compiler may fuse or reorder floating-point operations. The `log10` result is still monotone and the floor remains correct because it uses `std::max`, which is not reordered away.
- `FromDecibels` has no floor; at very large positive dB values the result overflows to `+inf` in float — this is expected behaviour for out-of-range inputs.

## Variants & Generalizations

- Power decibels use $10 \log_{10}(\cdot)$ (factor 10 rather than 20). The amplitude convention used here ($\times 20$) is correct for voltage, pressure, and filter transfer-function magnitude.
- Field-quantity vs. power-quantity disambiguation: IEEE 60268 / IEC 61672 mandate $20 \log_{10}$ for sound pressure level; the same convention applies to filter magnitude response.
- The floor can be parameterised if a stricter or looser sentinel is required; the default of $-160\,\mathrm{dB}$ is conservative for 32-bit float.

## Applications

- Filter specification: stop-band attenuation and pass-band ripple in dB are the primary acceptance criteria for IIR/FIR designs.
- Frequency response plots: `FrequencyResponse::Calculate()` already returns magnitude in dB using $20 \log_{10}$; these helpers provide the same conversion for ad-hoc analysis.
- Controller gain margin is expressed in dB; converting from a linear ratio with `ToDecibels` avoids duplication.
- Audio dynamic processing (compressor thresholds, limiter ceilings) and acoustic measurement both use dB natively.

## Connections to Other Algorithms

- `control_analysis::FrequencyResponse` internally applies $20 \log_{10}(\|H\|)$ on its magnitude output vector; these helpers are the scalar equivalent exposed for library consumers.
- Pass-band ripple computed by `RippleDb` feeds directly into filter-design acceptance testing alongside the step/transient-response metrics in the evaluation primitives family.

## References & Further Reading

- Proakis, J. & Manolakis, D., "Digital Signal Processing", 4th ed., Prentice Hall, 2007 — Appendix A (decibel notation).
- Zolzer, U., "DAFX: Digital Audio Effects", 2nd ed., Wiley, 2011 — Chapter 2 (level and gain in dB).
- IEC 61672-1:2013, "Electroacoustics — Sound level meters — Part 1: Specifications."
3 changes: 2 additions & 1 deletion doc/analysis/README.md
Comment thread
gabrielfrasantos marked this conversation as resolved.
Original file line number Diff line number Diff line change
Expand Up @@ -11,7 +11,8 @@ Signal analysis algorithms for frequency-domain decomposition and spectral estim
| [Power Spectral Density](PowerDensitySpectrum.md) | Estimation of signal power distribution across frequencies using Welch's method |
| [Discrete Cosine Transform](DiscreteCosineTransform.md) | Real-valued frequency decomposition via cosine basis functions, computed through FFT |
| [Signal Detectors](SignalDetectors.md) | Peak hold, zero-crossing counter, and RMS envelope detectors for real-time signal monitoring |
| [Goertzel Algorithm](GoertzelAlgorithm.md) | Single-bin DFT via a second-order recurrence for O(N) tone detection with O(1) memory |
| [Decibels](Decibels.md) | `ToDecibels` / `FromDecibels` conversion helpers with zero-floor guard, plus attenuation and ripple utilities |
| [Goertzel Algorithm](GoertzelAlgorithm.md) | Single-bin DFT via a second-order recurrence for O(N) tone detection with O(1) memory |

## Sub-domains

Expand Down
5 changes: 3 additions & 2 deletions doc/math/README.md
Original file line number Diff line number Diff line change
Expand Up @@ -6,5 +6,6 @@ Core mathematical primitives for numerical computation.

| Algorithm | Description |
|-----------------------------|-----------------------------------------------------------------------------------------------|
| [CORDIC](Cordic.md) | Iterative shift-add engine for sin/cos, atan2, magnitude, and vector rotation — no multiplier |
| [Quaternion](Quaternion.md) | Unit-quaternion rotation type: Hamilton product, SLERP, rotation-matrix and Euler conversions |
| [CORDIC](Cordic.md) | Iterative shift-add engine for sin/cos, atan2, magnitude, and vector rotation — no multiplier |
| [Quaternion](Quaternion.md) | Unit-quaternion rotation type: Hamilton product, SLERP, rotation-matrix and Euler conversions |
| [Step Response Metrics](StepResponseMetrics.md) | Rise time, settling time, percent overshoot, peak time, and steady-state error from a bounded step-response vector |
Comment thread
gabrielfrasantos marked this conversation as resolved.
107 changes: 107 additions & 0 deletions doc/math/StepResponseMetrics.md
Original file line number Diff line number Diff line change
@@ -0,0 +1,107 @@
# Step / Transient-Response Metrics

## Overview & Motivation

When a control system or filter receives a step input, its output traces a transient trajectory before settling at the final value. Quantifying that trajectory with standardised scalar metrics — rise time, settling time, percent overshoot, peak time, and steady-state error — is the primary acceptance test for any closed-loop design. These metrics translate the raw sample sequence into the language of control specifications, allowing automated pass/fail decisions without manual inspection of time-domain plots.

## Mathematical Theory

### Definitions

Let $y[k]$, $k = 0, \ldots, N-1$ be the sampled step response and $y_{ss}$ the steady-state value. The sample period is $\Delta t$.

**Rise Time** $T_r$

The elapsed time for the response to travel from 10 % to 90 % of steady state:

$$T_r = (k_{90} - k_{10})\,\Delta t$$

where $k_{10} = \min\{k : y[k] \ge 0.1\,y_{ss}\}$ and $k_{90} = \min\{k \ge k_{10} : y[k] \ge 0.9\,y_{ss}\}$.

**Settling Time** $T_s$

The first time after which the response remains permanently inside the band $[(1-\delta)y_{ss},\,(1+\delta)y_{ss}]$ (typically $\delta = 0.02$):

$$T_s = (k^* + 1)\,\Delta t, \quad k^* = \max\{k : |y[k] - y_{ss}| > \delta\,|y_{ss}|\}$$

**Percent Overshoot** $\%OS$

$$\%OS = 100\,\frac{y_{\max} - y_{ss}}{y_{ss}}, \quad y_{\max} = \max_k y[k]$$

For an underdamped second-order system with damping ratio $\zeta$:

$$\%OS = 100\,\exp\!\left(-\frac{\pi\zeta}{\sqrt{1-\zeta^2}}\right)$$

**Peak Time** $T_p$

$$T_p = k_p\,\Delta t, \quad k_p = \arg\max_k y[k]$$

For a continuous underdamped second-order system with natural frequency $\omega_n$:

$$T_p = \frac{\pi}{\omega_n\sqrt{1-\zeta^2}}$$

**Steady-State Error** $e_{ss}$

$$e_{ss} = r - \bar{y}_{\text{tail}}$$

where $r$ is the reference (command) value and $\bar{y}_{\text{tail}}$ is the mean of the final quarter of the response buffer, providing a robust estimate of the achieved steady state.

## Complexity Analysis

| Case | Time | Space | Notes |
|---------|----------|--------|--------------------------------------------|
| All | $O(N)$ | $O(1)$ | Single forward pass; no auxiliary storage |
Comment thread
gabrielfrasantos marked this conversation as resolved.

Each metric requires at most one traversal of the $N$-element vector. The tail-mean for steady-state error adds a constant-fraction second scan of the same data — still $O(N)$ total.

## Step-by-Step Walkthrough

Consider a 10-sample ramp to $y_{ss} = 1$ followed by a constant plateau (N = 20):

```
k: 0 1 2 3 4 5 6 7 8 9 10 11 …
y: 0 .1 .2 .3 .4 .5 .6 .7 .8 .9 1 1 …
```

- **Rise Time:** $k_{10} = 1$ (first sample $\ge 0.1$), $k_{90} = 9$ (first sample $\ge 0.9$). $T_r = 8\,\Delta t$.
- **Settling Time:** With $\delta = 0.02$, last sample outside the band is $k = 9$. $T_s = 10\,\Delta t$.
- **Percent Overshoot:** $y_{\max} = 1.0 = y_{ss}$, so $\%OS = 0$.
- **Peak Time:** $k_p = 10$ (first occurrence of max). $T_p = 10\,\Delta t$.
- **Steady-State Error:** Tail mean $= 1.0$, reference $= 1.0$. $e_{ss} = 0$.

## Pitfalls & Edge Cases

**Zero steady state.** Division by $y_{ss}$ in percent overshoot is guarded; the function returns zero when $y_{ss} = 0$ to avoid a NaN.

**Non-monotone ramp.** If the response crosses 90 % before 10 % (e.g., DC offset or wrong initial condition), $k_{10}$ may be found after the first 90 % crossing. The implementation returns the first pair that satisfies the threshold order.

**Oscillatory settling.** Settling time is defined as the last time the trajectory leaves the band, not the first time it enters it. Repeated crossings near the boundary extend the metric correctly.

**Finite buffer.** With a bounded vector of length $N$, if the response has not yet settled by the final sample, `SettlingTime` returns $N\,\Delta t$ and `RiseTime` returns $(N-1)\,\Delta t$ as conservative bounds.

**Tail-mean length.** Using the last $\lfloor N/4 \rfloor + 1$ samples for the steady-state estimate assumes the transient has decayed to within numerical noise by that point. Poorly chosen $N$ relative to the system time constant degrades the estimate.

## Variants & Generalizations

- **Delay Time** $T_d$: the time to reach 50 % of steady state — obtainable with the same threshold-scan pattern.
- **Band-relative rise time**: using a band other than 10–90 % (e.g., 20–80 %) is a trivial parameter change.
- **Multi-channel:** applying the scalar functions element-wise to each row of a response matrix generalises to MIMO systems without algorithmic change.

## Applications

- Automated controller tuning acceptance: verify that a PID or LQR design meets specification ($T_r < T_{r,\text{spec}}$, $\%OS < \%OS_{\text{spec}}$, etc.).
- Filter characterisation: measure the transient of a step fed through an IIR or FIR filter.
- Hardware-in-the-loop test harnesses: compute metrics directly from sampled actuator responses.

## Connections to Other Algorithms

- **Statistics** (this library): the tail-mean for steady-state error replicates the `Mean` function on a sub-range.
- **LinearTimeInvariant**: the primary source of step responses whose metrics are evaluated here.
- **Filters/active** (Kalman, EKF): step-excitation tests use these metrics to validate estimator transient behaviour.
- **Controllers**: PID and LQR tuning loops iterate until all five metrics satisfy design targets.

## References & Further Reading

- K. J. Åström and R. M. Murray, *Feedback Systems: An Introduction for Scientists and Engineers*, Princeton University Press, 2008. Chapter 10.
- G. F. Franklin, J. D. Powell, and A. Emami-Naeini, *Feedback Control of Dynamic Systems*, 8th ed., Pearson, 2019. Chapter 3.
- N. S. Nise, *Control Systems Engineering*, 8th ed., Wiley, 2019. Chapter 4.
1 change: 1 addition & 0 deletions numerical/analysis/CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -13,6 +13,7 @@ target_link_libraries(numerical.analysis ${NUMERICAL_VISIBILITY}

target_sources(numerical.analysis PRIVATE
ConvolutionCorrelation.hpp
Decibels.hpp
DiscreteCosineTransform.hpp
FastFourierTransform.hpp
FastFourierTransformRadix2Impl.hpp
Expand Down
Loading
Loading