|
| 1 | +use num_complex::{Complex64, c64}; |
| 2 | + |
| 3 | +use crate::{ |
| 4 | + systems::Tf, |
| 5 | + utils::traits::{Mag2Db, Rad2Deg, Time}, |
| 6 | +}; |
| 7 | +/// Generates a linearly spaced iterator between `start` and `end`, inclusive. |
| 8 | +/// |
| 9 | +/// # Arguments |
| 10 | +/// - `start`: The starting value of the sequence. |
| 11 | +/// - `end`: The ending value of the sequence. |
| 12 | +/// - `n`: The number of points to generate (must be greater than 1). |
| 13 | +/// |
| 14 | +/// # Returns |
| 15 | +/// - An iterator producing `n` evenly spaced `f64` values from `start` to |
| 16 | +/// `end`. |
| 17 | +/// |
| 18 | +/// # Panics |
| 19 | +/// - Panics if `n` is less than or equal to 1. |
| 20 | +pub fn lin_space( |
| 21 | + start: f64, |
| 22 | + end: f64, |
| 23 | + n: usize, |
| 24 | +) -> impl ExactSizeIterator<Item = f64> { |
| 25 | + assert!(n >= 1, "n must be greater than or equal to one"); |
| 26 | + let step = (end - start) / (n as f64 - 1.0); |
| 27 | + (0..n).map(move |i| start + step * i as f64) |
| 28 | +} |
| 29 | + |
| 30 | +/// Generates a logarithmically spaced iterator between `start` and `end`, using |
| 31 | +/// the specified logarithmic base. |
| 32 | +/// |
| 33 | +/// # Arguments |
| 34 | +/// - `start`: The starting value of the sequence (must be greater than 0). |
| 35 | +/// - `end`: The ending value of the sequence (must be greater than 0). |
| 36 | +/// - `n`: The number of points to generate. |
| 37 | +/// - `base`: The logarithmic base to use for spacing. |
| 38 | +/// |
| 39 | +/// # Returns |
| 40 | +/// - An iterator producing `n` logarithmically spaced `f64` values from `start` |
| 41 | +/// to `end`. |
| 42 | +/// |
| 43 | +/// # Panics |
| 44 | +/// - Panics if `start` or `end` is less than or equal to 0. |
| 45 | +pub fn log_space( |
| 46 | + start: f64, |
| 47 | + end: f64, |
| 48 | + n: usize, |
| 49 | + base: usize, |
| 50 | +) -> impl ExactSizeIterator<Item = f64> { |
| 51 | + assert!( |
| 52 | + start > 0., |
| 53 | + "logarithm of negative numbers are not implemented" |
| 54 | + ); |
| 55 | + assert!( |
| 56 | + end > 0., |
| 57 | + "logarithm of negative numbers are not implemented" |
| 58 | + ); |
| 59 | + assert!(base > 1, "log_base() must be well defined"); |
| 60 | + let start_log = start.log(base as f64); |
| 61 | + let end_log = end.log(base as f64); |
| 62 | + |
| 63 | + let nums = lin_space(start_log, end_log, n); |
| 64 | + |
| 65 | + nums.map(move |x| (base as f64).powf(x)) |
| 66 | +} |
| 67 | + |
| 68 | +impl<U: Time> Tf<f64, U> { |
| 69 | + /// Computes the Bode plot (magnitude and phase) for a transfer function |
| 70 | + /// over a frequency range. |
| 71 | + /// |
| 72 | + /// # Arguments |
| 73 | + /// - `sys`: The transfer function of the system to evaluate. |
| 74 | + /// - `min_freq`: The minimum frequency for the plot. |
| 75 | + /// - `max_freq`: The maximum frequency for the plot. |
| 76 | + /// |
| 77 | + /// # Returns |
| 78 | + /// - A vector of `[magnitude (dB), phase (degrees), frequency]` tuples for |
| 79 | + /// each evaluated frequency. |
| 80 | + pub fn bode(&self, min_freq: f64, max_freq: f64) -> Vec<[f64; 3]> { |
| 81 | + let freqs = log_space(min_freq, max_freq, 1000, 10); |
| 82 | + self.bode_freqs(freqs) |
| 83 | + } |
| 84 | + |
| 85 | + /// Computes the Bode plot (magnitude and phase) for a transfer function |
| 86 | + /// over a given set of frequencies. |
| 87 | + /// |
| 88 | + /// # Arguments |
| 89 | + /// - `sys`: The transfer function of the system to evaluate. |
| 90 | + /// - `freqs`: An iterator of frequencies to evaluate the system at. |
| 91 | + /// |
| 92 | + /// # Returns |
| 93 | + /// - A vector of `[magnitude (dB), phase (degrees), frequency]` tuples for |
| 94 | + /// each evaluated frequency. |
| 95 | + pub fn bode_freqs( |
| 96 | + &self, |
| 97 | + freqs: impl Iterator<Item = f64>, |
| 98 | + ) -> Vec<[f64; 3]> { |
| 99 | + let mut mag_phase_freq_vec = Vec::with_capacity(freqs.size_hint().0); |
| 100 | + |
| 101 | + for omega in freqs { |
| 102 | + let c = c64(0., omega); |
| 103 | + let sys_val = self.eval(&c); |
| 104 | + mag_phase_freq_vec.push([ |
| 105 | + sys_val.norm().mag2db(), |
| 106 | + sys_val.arg().rad2deg(), |
| 107 | + omega, |
| 108 | + ]); |
| 109 | + } |
| 110 | + mag_phase_freq_vec |
| 111 | + } |
| 112 | + |
| 113 | + /// Computes the Nyquist plot for a transfer function over a frequency |
| 114 | + /// range. |
| 115 | + /// |
| 116 | + /// # Arguments |
| 117 | + /// - `sys`: The transfer function of the system to evaluate. |
| 118 | + /// - `min_freq`: The minimum frequency for the plot. |
| 119 | + /// - `max_freq`: The maximum frequency for the plot. |
| 120 | + /// |
| 121 | + /// # Returns |
| 122 | + /// - A vector of complex numbers representing the Nyquist plot. |
| 123 | + pub fn nyquist(&self, min_freq: f64, max_freq: f64) -> Vec<Complex64> { |
| 124 | + let freqs = log_space(min_freq, max_freq, 1000, 10); |
| 125 | + self.nyquist_freqs(freqs) |
| 126 | + } |
| 127 | + |
| 128 | + /// Computes the Nyquist plot for a transfer function over a given set of |
| 129 | + /// frequencies. |
| 130 | + /// |
| 131 | + /// # Arguments |
| 132 | + /// - `sys`: The transfer function of the system to evaluate. |
| 133 | + /// - `freqs`: An iterator of frequencies to evaluate the system at. |
| 134 | + /// |
| 135 | + /// # Returns |
| 136 | + /// - A vector of complex numbers representing the Nyquist plot. |
| 137 | + pub fn nyquist_freqs( |
| 138 | + &self, |
| 139 | + freqs: impl Iterator<Item = f64>, |
| 140 | + ) -> Vec<Complex64> { |
| 141 | + let mut pos_vals = Vec::with_capacity(freqs.size_hint().0); |
| 142 | + let mut neg_vals = Vec::with_capacity(freqs.size_hint().0); |
| 143 | + |
| 144 | + for freq in freqs { |
| 145 | + pos_vals.push(self.eval(&c64(0., freq))); |
| 146 | + neg_vals.push(self.eval(&c64(0., -freq))); |
| 147 | + } |
| 148 | + |
| 149 | + pos_vals.extend(neg_vals.iter().rev()); |
| 150 | + pos_vals |
| 151 | + } |
| 152 | +} |
0 commit comments