mxlib
c++ tools for analyzing astronomical data and other tasks by Jared R. Males. [git repo]
|
A class to manage optimizing closed-loop gains.
_realT | the real floating point type in which to do all arithmetic. Is used to define the complex type as well. |
Definition at line 63 of file clGainOpt.hpp.
#include <ao/analysis/clGainOpt.hpp>
Public Types | |
typedef _realT | realT |
The real data type. | |
typedef std::complex< _realT > | complexT |
The complex data type. | |
Public Member Functions | |
clGainOpt () | |
Default c'tor. | |
clGainOpt (realT Ti, realT tau) | |
C'tor setting the loop timings. | |
void | init () |
Initialize this instance. | |
int | N () |
Get the number of integrations in the (optional) moving average. | |
void | N (int newN) |
Set the number of integrations in the moving average. | |
realT | Ti () |
Get the loop sampling interval. | |
void | Ti (realT newTi) |
Set the loop sampling interval. | |
realT | tau () |
Get the loop delay. | |
void | tau (realT newTau) |
Set the loop delay. | |
void | b (const std::vector< realT > &newB) |
Set the vector of FIR coefficients. | |
void | b (const Eigen::Array< realT, -1, -1 > &newB) |
Set the vector of FIR coefficients. | |
realT | b (size_t i) |
Get a single FIR coefficient. | |
const std::vector< realT > & | b () |
Get the vector of FIR coefficients. | |
void | a (const std::vector< realT > &newA) |
Set the vector of IIR coefficients. | |
void | a (const Eigen::Array< realT, -1, -1 > &newA) |
Set the vector of IIR coefficients. | |
realT | a (size_t i) |
Get a single IIR coefficient. | |
const std::vector< realT > & | a () |
Get the vector of IIR coefficients. | |
void | remember (const realT &rem) |
Set the remember factor for a leaky integrator. | |
realT | remember () |
Get the remember factor. | |
void | setLeakyIntegrator (realT remember) |
Set the FIR and IIR coefficients so that the control law is a leaky integrator. | |
void | f (realT *newF, size_t nF) |
Set the vector of frequencies. | |
void | f (std::vector< realT > &newF) |
Set the vector of frequencies. | |
size_t | f_size () |
Get the size of the frequency vector. | |
realT | f (size_t i) |
Get the i-th value of frequency. | |
complexT | olXfer (int fi, complexT &H_dm, complexT &H_del, complexT &H_con) |
Calculate the open-loop transfer function. | |
complexT | olXfer (int fi) |
Calculate the open-loop transfer function. | |
complexT | clETF (int fi, realT g) |
Return the closed loop error transfer function (ETF) at frequency f for gain g. | |
realT | clETFPhase (int fi, realT g) |
Return the closed loop error transfer function (ETF) phase at frequency f for gain g. | |
realT | clETF2 (int fi, realT g) |
Return the norm of the closed loop error transfer function (ETF) at frequency f for gain g. | |
complexT | clNTF (int fi, realT g) |
Return the closed loop noise transfer function (NTF) at frequency f for gain g. | |
realT | clNTF2 (int fi, realT g) |
Return the norm of the closed loop noise transfer function (NTF) at frequency f for gain g. | |
void | clTF2 (realT &ETF, realT &NTF, int fi, realT g) |
Return the norm of the closed loop transfer functions at frequency f for gain g. | |
realT | clVariance (realT &varErr, realT &varNoise, const std::vector< realT > &PSDerr, const std::vector< realT > &PSDnoise, realT g) |
Calculate the closed loop variance given open-loop PSDs and gain. | |
realT | clVariance (const std::vector< realT > &PSDerr, const std::vector< realT > &PSDnoise, realT g) |
Calculate the closed loop variance given open-loop PSDs and gain. | |
realT | maxStableGain (realT &ll, realT &ul) |
Find the maximum stable gain for the loop parameters. | |
realT | maxStableGain (const realT &ll, const realT &ul) |
Find the maximum stable gain for the loop parameters. | |
realT | maxStableGain () |
Find the maximum stable gain for the loop parameters. | |
realT | optGainOpenLoop (realT &var, const std::vector< realT > &PSDerr, const std::vector< realT > &PSDnoise, bool gridSearch) |
Return the optimum closed loop gain given an open loop PSD. | |
realT | optGainOpenLoop (realT &var, const std::vector< realT > &PSDerr, const std::vector< realT > &PSDnoise, realT &gmax, bool gridSearch) |
Return the optimum closed loop gain given an open loop PSD. | |
int | pseudoOpenLoop (std::vector< realT > &PSD, realT g) |
Calculate the pseudo open-loop PSD given a closed loop PSD. | |
Protected Attributes | |
int | m_N |
Number of integrations in the (optional) moving average. Default is 1. | |
realT | m_Ti |
The loop sampling interval. | |
realT | m_tau |
The loop delay. | |
realT | m_remember { 1.0 } |
The leaky integrator forget factor. | |
std::vector< realT > | m_b |
Vector of FIR coefficients. | |
std::vector< realT > | m_a |
Vector of IIR coefficients. | |
std::vector< realT > | m_f |
Vector of frequencies. | |
bool | m_fChanged { true } |
True if frequency or max size of m_a and m_b changes. | |
bool | m_changed { true } |
True if any of the members which make up the basic transfer functions are changed. | |
typedef std::complex<_realT> mx::AO::analysis::clGainOpt< _realT >::complexT |
The complex data type.
Definition at line 66 of file clGainOpt.hpp.
typedef _realT mx::AO::analysis::clGainOpt< _realT >::realT |
The real data type.
Definition at line 65 of file clGainOpt.hpp.
mx::AO::analysis::clGainOpt< realT >::clGainOpt | ( | ) |
Default c'tor.
Definition at line 440 of file clGainOpt.hpp.
mx::AO::analysis::clGainOpt< realT >::clGainOpt | ( | realT | Ti, |
realT | tau | ||
) |
C'tor setting the loop timings.
[in] | Ti | the desired loop sampling interval. |
[in] | tau | the desired loop delay. |
Definition at line 446 of file clGainOpt.hpp.
|
inline |
Get the vector of IIR coefficients.
Definition at line 212 of file clGainOpt.hpp.
References mx::AO::analysis::clGainOpt< _realT >::m_a.
void mx::AO::analysis::clGainOpt< realT >::a | ( | const Eigen::Array< realT, -1, -1 > & | newA | ) |
Set the vector of IIR coefficients.
[in] | newA | a column-vector Eigen::Array of coefficients, which is copied to m_a. |
Definition at line 583 of file clGainOpt.hpp.
void mx::AO::analysis::clGainOpt< realT >::a | ( | const std::vector< realT > & | newA | ) |
Set the vector of IIR coefficients.
[in] | newA | a vector of coefficients, which is copied to m_a. |
Definition at line 571 of file clGainOpt.hpp.
Referenced by mx::AO::analysis::clAOLinearPredictor< _realT >::_regularizeCoefficients(), mx::AO::analysis::fourierTemporalPSD< _realT, aosysT >::analyzePSDGrid(), and mx::AO::analysis::clAOLinearPredictor< _realT >::regularizeCoefficients().
|
inline |
Get a single IIR coefficient.
Definition at line 204 of file clGainOpt.hpp.
References mx::AO::analysis::clGainOpt< _realT >::m_a.
|
inline |
Get the vector of FIR coefficients.
Definition at line 183 of file clGainOpt.hpp.
References mx::AO::analysis::clGainOpt< _realT >::m_b.
void mx::AO::analysis::clGainOpt< realT >::b | ( | const Eigen::Array< realT, -1, -1 > & | newB | ) |
Set the vector of FIR coefficients.
[in] | newB | a column-vector Eigen::Array of coefficients, which is copied to m_b. |
Definition at line 542 of file clGainOpt.hpp.
void mx::AO::analysis::clGainOpt< realT >::b | ( | const std::vector< realT > & | newB | ) |
Set the vector of FIR coefficients.
[in] | newB | a vector of coefficients, which is copied to m_b. |
Definition at line 530 of file clGainOpt.hpp.
Referenced by mx::AO::analysis::clAOLinearPredictor< _realT >::_regularizeCoefficients(), mx::AO::analysis::fourierTemporalPSD< _realT, aosysT >::analyzePSDGrid(), and mx::AO::analysis::clAOLinearPredictor< _realT >::regularizeCoefficients().
|
inline |
Get a single FIR coefficient.
[in] | i | the index of the FIR coefficient |
Definition at line 175 of file clGainOpt.hpp.
References mx::AO::analysis::clGainOpt< _realT >::m_b.
std::complex< realT > mx::AO::analysis::clGainOpt< realT >::clETF | ( | int | fi, |
realT | g | ||
) |
Return the closed loop error transfer function (ETF) at frequency f for gain g.
[in] | fi | the index of the frequency at which to calculate the ETF |
[in] | g | the loop gain. |
Definition at line 844 of file clGainOpt.hpp.
Referenced by mx::AO::analysis::fourierTemporalPSD< _realT, aosysT >::analyzePSDGrid().
realT mx::AO::analysis::clGainOpt< realT >::clETF2 | ( | int | fi, |
realT | g | ||
) |
Return the norm of the closed loop error transfer function (ETF) at frequency f for gain g.
[in] | fi | the index of the frequency at which to calculate the ETF. |
[in] | g | the loop gain. |
Definition at line 872 of file clGainOpt.hpp.
realT mx::AO::analysis::clGainOpt< realT >::clETFPhase | ( | int | fi, |
realT | g | ||
) |
Return the closed loop error transfer function (ETF) phase at frequency f for gain g.
[in] | fi | the index of the frequency at which to calculate the ETF |
[in] | g | the loop gain. |
Definition at line 858 of file clGainOpt.hpp.
std::complex< realT > mx::AO::analysis::clGainOpt< realT >::clNTF | ( | int | fi, |
realT | g | ||
) |
Return the closed loop noise transfer function (NTF) at frequency f for gain g.
[in] | fi | the index of the frequency at which to calculate the NTF |
[in] | g | the loop gain. |
Definition at line 886 of file clGainOpt.hpp.
Referenced by mx::AO::analysis::fourierTemporalPSD< _realT, aosysT >::analyzePSDGrid().
realT mx::AO::analysis::clGainOpt< realT >::clNTF2 | ( | int | fi, |
realT | g | ||
) |
Return the norm of the closed loop noise transfer function (NTF) at frequency f for gain g.
[in] | fi | the index of the frequency at which to calculate the NTF |
[in] | g | the loop gain. |
Definition at line 904 of file clGainOpt.hpp.
void mx::AO::analysis::clGainOpt< realT >::clTF2 | ( | realT & | ETF, |
realT & | NTF, | ||
int | fi, | ||
realT | g | ||
) |
Return the norm of the closed loop transfer functions at frequency f for gain g.
Calculates both the error transfer function (ETF) and the noise transfer function (NTF). This minimizes the various complex number operations, compared to calling both clETF2 and clNTF2.
[out] | ETF | is set to the ETF at f and g |
[out] | NTF | is set to the NTF at f and g |
[in] | fi | is the index of the frequency at which to calculate the ETF |
[in] | g | is the loop gain. |
Definition at line 924 of file clGainOpt.hpp.
Referenced by mx::AO::analysis::fourierTemporalPSD< _realT, aosysT >::analyzePSDGrid().
realT mx::AO::analysis::clGainOpt< realT >::clVariance | ( | const std::vector< realT > & | PSDerr, |
const std::vector< realT > & | PSDnoise, | ||
realT | g | ||
) |
Calculate the closed loop variance given open-loop PSDs and gain.
Overload of clVariance without the varErr and varNoise output parameters.
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts.
[in] | PSDerr | the open-loop process error PSD. |
[in] | PSDnoise | the open-loop measurement noise PSD. |
[in] | g | the gain. |
Definition at line 990 of file clGainOpt.hpp.
realT mx::AO::analysis::clGainOpt< realT >::clVariance | ( | realT & | varErr, |
realT & | varNoise, | ||
const std::vector< realT > & | PSDerr, | ||
const std::vector< realT > & | PSDnoise, | ||
realT | g | ||
) |
Calculate the closed loop variance given open-loop PSDs and gain.
Calculates the following quantities.
\[ \sigma_{err}^2 = \sum_i \left| ETF(f_i) \right|^2 PSD_{err}(fi) \Delta f\\ \sigma_{noise}^2 = \sum_i \left| NTF(f_i) \right|^2 PSD_{noise}(fi) \Delta f\\ \sigma^2 = \sigma_{err}^2 + \sigma_{noise}^2 \]
\( \sigma^2 \) is returned, and \( \sigma_{err}^2 \) and \( \sigma_{noise}^2 \) are available as the optional arguments varErr and varNoise.
[out] | varErr | the variance in the residual process error. |
[out] | varNoise | the variance in the residual measurement noise. |
[in] | PSDerr | the open-loop process error PSD. |
[in] | PSDnoise | the open-loop measurement noise PSD. |
[in] | g | the gain. |
Definition at line 955 of file clGainOpt.hpp.
Referenced by mx::AO::analysis::fourierTemporalPSD< _realT, aosysT >::analyzePSDGrid().
void mx::AO::analysis::clGainOpt< realT >::f | ( | realT * | newF, |
size_t | nF | ||
) |
Set the vector of frequencies.
[in] | newF | a pointer to an array containing the new frequencies |
[in] | nF | the number of elements of size(realT) in newF. |
Definition at line 658 of file clGainOpt.hpp.
Referenced by mx::AO::analysis::fourierTemporalPSD< _realT, aosysT >::analyzePSDGrid().
realT mx::AO::analysis::clGainOpt< realT >::f | ( | size_t | i | ) |
Get the i-th value of frequency.
No range checks are conducted.
[in] | i | the index of the frequency to return |
Definition at line 680 of file clGainOpt.hpp.
void mx::AO::analysis::clGainOpt< realT >::f | ( | std::vector< realT > & | newF | ) |
Set the vector of frequencies.
[in] | newF | a vector containing the new frequencies |
Definition at line 671 of file clGainOpt.hpp.
|
inline |
Get the size of the frequency vector.
Definition at line 248 of file clGainOpt.hpp.
References mx::AO::analysis::clGainOpt< _realT >::m_f.
void mx::AO::analysis::clGainOpt< realT >::init | ( | ) |
Initialize this instance.
Definition at line 455 of file clGainOpt.hpp.
realT mx::AO::analysis::clGainOpt< realT >::maxStableGain | ( | ) |
Find the maximum stable gain for the loop parameters.
Conducts a search along the Nyquist contour of the open-loop transfer function to find the most-negative crossing of the real axis.
This version uses m_maxFindMin for the lower limit and no upper limit.
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts.
Definition at line 1038 of file clGainOpt.hpp.
realT mx::AO::analysis::clGainOpt< _realT >::maxStableGain | ( | const realT & | ll, |
const realT & | ul | ||
) |
Find the maximum stable gain for the loop parameters.
Conducts a search along the Nyquist contour of the open-loop transfer function to find the most-negative crossing of the real axis.
This version allows constant arguments. This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts.
[in] | ll | the lower limit used for the search |
[in] | ul | the upper limit used for hte search |
realT mx::AO::analysis::clGainOpt< realT >::maxStableGain | ( | realT & | ll, |
realT & | ul | ||
) |
Find the maximum stable gain for the loop parameters.
ll | [in.out] the lower limit used for the search |
ul | [in.out] the upper limit used for hte search |
Definition at line 999 of file clGainOpt.hpp.
Referenced by mx::AO::analysis::clAOLinearPredictor< _realT >::_regularizeCoefficients(), and mx::AO::analysis::clAOLinearPredictor< _realT >::regularizeCoefficients().
int mx::AO::analysis::clGainOpt< realT >::N | ( | ) |
Get the number of integrations in the (optional) moving average.
Definition at line 476 of file clGainOpt.hpp.
void mx::AO::analysis::clGainOpt< realT >::N | ( | int | newN | ) |
Set the number of integrations in the moving average.
[in] | newN | the value of m_N. |
Definition at line 482 of file clGainOpt.hpp.
std::complex< realT > mx::AO::analysis::clGainOpt< realT >::olXfer | ( | int | fi | ) |
Calculate the open-loop transfer function.
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts.
[in] | fi | the index of the frequency |
Definition at line 687 of file clGainOpt.hpp.
std::complex< realT > mx::AO::analysis::clGainOpt< realT >::olXfer | ( | int | fi, |
complexT & | H_dm, | ||
complexT & | H_del, | ||
complexT & | H_con | ||
) |
Calculate the open-loop transfer function.
[in] | fi | the index of the frequency |
[out] | H_dm | the transfer function of the DM |
[out] | H_del | the delay transfer function |
[out] | H_con | the controller transfer function. |
Definition at line 701 of file clGainOpt.hpp.
References mx::math::six_fifths().
realT mx::AO::analysis::clGainOpt< realT >::optGainOpenLoop | ( | realT & | var, |
const std::vector< realT > & | PSDerr, | ||
const std::vector< realT > & | PSDnoise, | ||
bool | gridSearch | ||
) |
Return the optimum closed loop gain given an open loop PSD.
Uses _gmax for the upper limit.
[out] | var | the variance at the optimum gain |
[in] | PSDerr | open loop error PSD |
[in] | PSDnoise | open loop measurement noise PSD |
[in] | gridSearch | flag controlling whether an intial grid search is done to find the global minimum |
Definition at line 1136 of file clGainOpt.hpp.
References gmax.
Referenced by mx::AO::analysis::clAOLinearPredictor< _realT >::_regularizeCoefficients(), mx::AO::analysis::fourierTemporalPSD< _realT, aosysT >::analyzePSDGrid(), and mx::AO::analysis::clAOLinearPredictor< _realT >::regularizeCoefficients().
realT mx::AO::analysis::clGainOpt< realT >::optGainOpenLoop | ( | realT & | var, |
const std::vector< realT > & | PSDerr, | ||
const std::vector< realT > & | PSDnoise, | ||
realT & | gmax, | ||
bool | gridSearch | ||
) |
Return the optimum closed loop gain given an open loop PSD.
[out] | var | the variance at the optimum gain |
[in] | PSDerr | open loop error PSD |
[in] | PSDnoise | open loop measurement noise PSD |
[in] | gmax | maximum gain to consider. If 0, then _gmax is used. |
[in] | gridSearch | flag controlling whether an intial grid search is done to find the global minimum |
Definition at line 1146 of file clGainOpt.hpp.
References gmax.
int mx::AO::analysis::clGainOpt< realT >::pseudoOpenLoop | ( | std::vector< realT > & | PSD, |
realT | g | ||
) |
Calculate the pseudo open-loop PSD given a closed loop PSD.
PSD | [in.out] input closed loop PSD, on output contains the pseudo open loop error PSD | |
[in] | g | the loop gain when PSD was measured. |
Definition at line 1204 of file clGainOpt.hpp.
realT mx::AO::analysis::clGainOpt< realT >::remember | ( | ) |
Get the remember factor.
Definition at line 625 of file clGainOpt.hpp.
void mx::AO::analysis::clGainOpt< realT >::remember | ( | const realT & | rem | ) |
Set the remember factor for a leaky integrator.
Definition at line 612 of file clGainOpt.hpp.
void mx::AO::analysis::clGainOpt< realT >::setLeakyIntegrator | ( | realT | remember | ) |
Set the FIR and IIR coefficients so that the control law is a leaky integrator.
Set remember to 1.0 for a pure integrator control law.
[in] | remember | a number usually close to 1 setting the amount "remembered" from previous iterations. |
Definition at line 631 of file clGainOpt.hpp.
realT mx::AO::analysis::clGainOpt< realT >::tau | ( | ) |
Get the loop delay.
Definition at line 512 of file clGainOpt.hpp.
void mx::AO::analysis::clGainOpt< realT >::tau | ( | realT | newTau | ) |
Set the loop delay.
[in] | newTau | the new value of m_tau. |
Definition at line 518 of file clGainOpt.hpp.
realT mx::AO::analysis::clGainOpt< realT >::Ti | ( | ) |
Get the loop sampling interval.
Definition at line 494 of file clGainOpt.hpp.
void mx::AO::analysis::clGainOpt< realT >::Ti | ( | realT | newTi | ) |
Set the loop sampling interval.
[in] | newTi | the new value of m_Ti. |
Definition at line 500 of file clGainOpt.hpp.
|
protected |
Vector of IIR coefficients.
Definition at line 75 of file clGainOpt.hpp.
Referenced by mx::AO::analysis::clGainOpt< _realT >::a(), and mx::AO::analysis::clGainOpt< _realT >::a().
|
protected |
Vector of FIR coefficients.
Definition at line 74 of file clGainOpt.hpp.
Referenced by mx::AO::analysis::clGainOpt< _realT >::b(), and mx::AO::analysis::clGainOpt< _realT >::b().
|
protected |
True if any of the members which make up the basic transfer functions are changed.
Definition at line 81 of file clGainOpt.hpp.
|
protected |
Vector of frequencies.
Definition at line 77 of file clGainOpt.hpp.
Referenced by mx::AO::analysis::clGainOpt< _realT >::f_size().
|
protected |
True if frequency or max size of m_a and m_b changes.
Definition at line 79 of file clGainOpt.hpp.
realT mx::AO::analysis::clGainOpt< _realT >::m_maxFindMin |
The Minimum value for the maximum stable gain finding algorithm.
Parameters for stability analysis
Definition at line 97 of file clGainOpt.hpp.
int mx::AO::analysis::clGainOpt< _realT >::m_minFindBits |
The bits of precision to use for minimum finding. Defaults to std::numeric_limits<realT>::digits.
Definition at line 107 of file clGainOpt.hpp.
realT mx::AO::analysis::clGainOpt< _realT >::m_minFindMaxFact |
The maximum value, as a multiplicative factor of maximum gain.
Definition at line 106 of file clGainOpt.hpp.
uintmax_t mx::AO::analysis::clGainOpt< _realT >::m_minFindMaxIter |
The maximum iterations allowed for minimization.
Definition at line 109 of file clGainOpt.hpp.
realT mx::AO::analysis::clGainOpt< _realT >::m_minFindMin |
The Minimum value for the minimum finding algorithm.
Parameters for minimization finding
Definition at line 105 of file clGainOpt.hpp.
|
protected |
Number of integrations in the (optional) moving average. Default is 1.
Definition at line 69 of file clGainOpt.hpp.
|
protected |
The leaky integrator forget factor.
Definition at line 73 of file clGainOpt.hpp.
|
protected |
The loop delay.
Definition at line 71 of file clGainOpt.hpp.
|
protected |
The loop sampling interval.
Definition at line 70 of file clGainOpt.hpp.