template<typename _realT, typename aosysT>
struct mx::AO::analysis::fourierTemporalPSD< _realT, aosysT >

Class to manage the calculation of temporal PSDs of the Fourier modes in atmospheric turbulence.

Works with both basic (sines/cosines) and modified Fourier modes.

Template Parameters

realT	is a real floating point type for calculations. Currently must be double due to gsl_integration.
aosysT	is an AO system type, usually of type ao_system.

Todo:

Split off the integration parameters in a separate structure.

once integration parameters are in a separate structure, make this a class with protected members.

GSL error handling

Definition at line 114 of file fourierTemporalPSD.hpp.

#include <ao/analysis/fourierTemporalPSD.hpp>

Public Types
typedef _realT	realT
	The type for arithmetic.

typedef std::complex< realT >	complexT
	The complex type for arithmetic.

Public Member Functions
	fourierTemporalPSD ()
	Default c'tor.

	fourierTemporalPSD (bool allocate)
	Constructor with workspace allocation.

	~fourierTemporalPSD ()
	Destructor.

realT	fastestPeak (int m, int n)
	Determine the frequency of the highest V-dot-k peak.

int	singleLayerPSD (std::vector< realT > &PSD, std::vector< realT > &freq, realT m, realT n, int layer_i, int p, realT fmax=0)
	Calculate the temporal PSD for a Fourier mode for a single layer.

template<bool parallel = true>
int	multiLayerPSD (std::vector< realT > &PSD, std::vector< realT > &freq, realT m, realT n, int p, realT fmax=0)
	Calculate the temporal PSD for a Fourier mode in a multi-layer model.

void	makePSDGrid (const std::string &dir, int mnMax, realT dFreq, realT maxFreq, realT fmax=0)
	Calculate PSDs over a grid of spatial frequencies.

int	analyzePSDGrid (const std::string &subDir, const std::string &psdDir, int mnMax, int mnCon, realT gfixed, int lpNc, realT lpRegPrecision, std::vector< realT > &mags, int lifetimeTrials=0, bool uncontrolledLifetimes=false, bool writePSDs=false, bool writeXfer=false)
	Analyze a PSD grid under closed-loop control.

int	intensityPSD (const std::string &subDir, const std::string &psdDir, const std::string &CvdPath, int mnMax, int mnCon, const std::string &si_or_lp, std::vector< realT > &mags, int lifetimeTrials, bool writePSDs)

GSL Integration Tolerances
For good results it seems that absolute tolerance (absTol) needs to be 1e-10. Lower tolerances cause some frequencies to drop out, etc. Relative tolerance (relTol) seems to be less sensitive, and 1e-4 works on cases tested as of 1 Jan, 2017. See the documentation for the GSL Library integrators at (https://www.gnu.org/software/gsl/manual/htmlm_node/QAGI-adaptive-integration-on-infinite-intervals.html)
void	absTol (realT at)
	Set absolute tolerance.

realT	absTol ()
	Get the current absolute tolerance.

void	relTol (realT rt)
	Set relative tolerance.

realT	relTol ()
	Get the current relative tolerance.

Disk Storage
These methods handle writing to and reading from disk. The calculated PSDs are store in the mx::BinVector binary format. A grid of PSDs is specified by its directory name. The directory contains one frequency file (freq.binv), and a set of PSD files, named according to psd_<m>_<n>_.binv.
int	getGridFreq (std::vector< realT > &freq, const std::string &dir)
	Get the frequency scale for a PSD grid.

int	getGridPSD (std::vector< realT > &psd, const std::string &dir, int m, int n)
	Get a single PSD from a PSD grid.

int	getGridPSD (std::vector< realT > &freq, std::vector< realT > &psd, const std::string &dir, int m, int n)
	Get both the frequency scale and a single PSD from a PSD grid.

Public Attributes
aosysT *	m_aosys { nullptr }
	Pointer to an AO system structure.

realT	m_f { 0 }
	the current temporal frequency

realT	m_m { 0 }
	the spatial frequency m index

realT	m_n { 0 }
	the spatial frequency n index

realT	m_cq { 0 }
	The cosine of the wind direction.

realT	m_sq { 0 }
	The sine of the wind direction.

realT	m_spatialFilter { false }
	Flag indicating if a spatial filter is applied.

realT	m_f0 { 0 }
	the Berdja boiling parameter

int	m_p { 1 }

int	_layer_i
	The index of the current layer.

int	_useBasis

gsl_integration_workspace *	_w
	Workspace for the gsl integrators, allocated to WSZ if constructed as worker (with allocate == true).

realT	_absTol
	The absolute tolerance to use in the GSL integrator.

realT	_relTol
	The relative tolerance to use in the GSL integrator.

int	m_mode_i
	Projected basis mode index.

Eigen::Array< realT, -1, -1 >	m_modeCoeffs
	Coeeficients of the projection onto the Fourier modes.

Protected Member Functions
void	initialize ()
	Initialize parameters to default values.

Member Typedef Documentation

◆ complexT

template<typename _realT , typename aosysT >

typedef std::complex<realT> mx::AO::analysis::fourierTemporalPSD< _realT, aosysT >::complexT

The complex type for arithmetic.

Definition at line 120 of file fourierTemporalPSD.hpp.

◆ realT

template<typename _realT , typename aosysT >

typedef _realT mx::AO::analysis::fourierTemporalPSD< _realT, aosysT >::realT

The type for arithmetic.

Definition at line 117 of file fourierTemporalPSD.hpp.

Constructor & Destructor Documentation

◆ fourierTemporalPSD() [1/2]

template<typename realT , typename aosysT >

mx::AO::analysis::fourierTemporalPSD< realT, aosysT >::fourierTemporalPSD ( )

Default c'tor.

Definition at line 421 of file fourierTemporalPSD.hpp.

◆ fourierTemporalPSD() [2/2]

template<typename realT , typename aosysT >

mx::AO::analysis::fourierTemporalPSD< realT, aosysT >::fourierTemporalPSD ( bool allocate )

explicit

Constructor with workspace allocation.

Parameters

allocate if true, then the workspace for GSL integrators is allocated.

Definition at line 428 of file fourierTemporalPSD.hpp.

◆ ~fourierTemporalPSD()

template<typename realT , typename aosysT >

mx::AO::analysis::fourierTemporalPSD< realT, aosysT >::~fourierTemporalPSD ( )

Destructor.

Frees GSL workspace if it was allocated.

Definition at line 440 of file fourierTemporalPSD.hpp.

Member Function Documentation

◆ absTol() [1/2]

template<typename realT , typename aosysT >

realT mx::AO::analysis::fourierTemporalPSD< realT, aosysT >::absTol ( )

Get the current absolute tolerance.

Returns: _absTol

Definition at line 465 of file fourierTemporalPSD.hpp.

◆ absTol() [2/2]

template<typename realT , typename aosysT >

void mx::AO::analysis::fourierTemporalPSD< realT, aosysT >::absTol ( realT at )

Set absolute tolerance.

Parameters

at	is the new absolute tolerance.

Definition at line 459 of file fourierTemporalPSD.hpp.

◆ analyzePSDGrid()

template<typename realT , typename aosysT >

int mx::AO::analysis::fourierTemporalPSD< realT, aosysT >::analyzePSDGrid	(	const std::string &	subDir,
		const std::string &	psdDir,
		int	mnMax,
		int	mnCon,
		realT	gfixed,
		int	lpNc,
		realT	lpRegPrecision,
		std::vector< realT > &	mags,
		int	lifetimeTrials = `0`,
		bool	uncontrolledLifetimes = `false`,
		bool	writePSDs = `false`,
		bool	writeXfer = `false`
	)

Analyze a PSD grid under closed-loop control.

This always analyzes the simple integrator, and can also analyze the linear predictor controller. Outputs maps of optimum gains, predictor coefficients, variances, and contrasts for a range of guide star magnitudes. Optionally calculates speckle lifetimes. Optionally writes the closed-loop PSDs and transfer functions.

Todo:: need upstream NCP and CP NCP and NCP NCP

Parameters

[out]	subDir	the sub-directory of psdDir where to write the results. Is created.
[in]	psdDir	the directory containing the grid of PSDs.
[in]	mnMax	the maximum value of m and n in the grid.
[in]	mnCon	the maximum value of m and n which can be controlled.
[in]	gfixed	if > 0 then this fixed gain is used in the SI.
[in]	lpNc	the number of linear predictor coefficients to analyze. If 0 then LP is not analyzed.
[in]	lpRegPrecision	the initial precision for the LP regularization algorithm. Normal value is 2. Higher is faster. Decrease if getting stuck in local minima.
[in]	mags	the guide star magnitudes to analyze for.
[in]	lifetimeTrials	[optional] number of trials used for calculating speckle lifetimes. If 0, lifetimes are not calculated.
[in]	uncontrolledLifetimes	[optional] flag controlling whether lifetimes are calculated for uncontrolled modes.
[in]	writePSDs	[optional] flag controlling if resultant PSDs are saved
[in]	writeXfer	[optional] flag controlling if resultant Xfer functions are saved

Definition at line 784 of file fourierTemporalPSD.hpp.

References mx::ioutils::convertToString(), mx::ioutils::createDirectories(), mx::ipc::ompLoopWatcher< _outputT, _printPretty, _printLoops, _printPercent, _printNLine, _time >::incrementAndOutputStatus(), mx::sigproc::psdVar(), mx::AO::analysis::clAOLinearPredictor< _realT >::regularizeCoefficients(), mx::AO::analysis::speckleAmpPSD(), mx::fits::fitsFile< dataT >::write(), and mx::ioutils::writeBinVector().

◆ fastestPeak()

template<typename realT , typename aosysT >

realT mx::AO::analysis::fourierTemporalPSD< realT, aosysT >::fastestPeak	(	int	m,
		int	n
	)

Determine the frequency of the highest V-dot-k peak.

Parameters

m	the spatial frequency u index
n	the spatial frequency v index

Returns: the frequency of the fastest peak

Definition at line 483 of file fourierTemporalPSD.hpp.

◆ getGridFreq()

template<typename realT , typename aosysT >

int mx::AO::analysis::fourierTemporalPSD< realT, aosysT >::getGridFreq	(	std::vector< realT > &	freq,
		const std::string &	dir
	)

Get the frequency scale for a PSD grid.

Parameters

[out]	freq	the vector to populate with the frequency scale.
[in]	dir	specifies the directory containing the grid.

Definition at line 2200 of file fourierTemporalPSD.hpp.

References mx::ioutils::readBinVector().

◆ getGridPSD() [1/2]

template<typename realT , typename aosysT >

int mx::AO::analysis::fourierTemporalPSD< realT, aosysT >::getGridPSD	(	std::vector< realT > &	freq,
		std::vector< realT > &	psd,
		const std::string &	dir,
		int	m,
		int	n
	)

Get both the frequency scale and a single PSD from a PSD grid.

Parameters

[out]	freq	the vector to populate with the frequency scale.
[out]	psd	the vector to populate with the PSD.
[in]	dir	specifies the directory containing the grid.
[in]	m	specifies the u component of spatial frequency.
[in]	n	specifies the v component of spatial frequency.

Definition at line 2216 of file fourierTemporalPSD.hpp.

◆ getGridPSD() [2/2]

template<typename realT , typename aosysT >

int mx::AO::analysis::fourierTemporalPSD< realT, aosysT >::getGridPSD	(	std::vector< realT > &	psd,
		const std::string &	dir,
		int	m,
		int	n
	)

Get a single PSD from a PSD grid.

Parameters

[out]	psd	the vector to populate with the PSD.
[in]	dir	specifies the directory containing the grid.
[in]	m	specifies the u component of spatial frequency.
[in]	n	specifies the v component of spatial frequency.

Definition at line 2208 of file fourierTemporalPSD.hpp.

References mx::ioutils::convertToString(), and mx::ioutils::readBinVector().

◆ initialize()

template<typename realT , typename aosysT >

void mx::AO::analysis::fourierTemporalPSD< realT, aosysT >::initialize ( )

protected

Initialize parameters to default values.

Definition at line 449 of file fourierTemporalPSD.hpp.

◆ intensityPSD()

template<typename realT , typename aosysT >

int mx::AO::analysis::fourierTemporalPSD< realT, aosysT >::intensityPSD	(	const std::string &	subDir,
		const std::string &	psdDir,
		const std::string &	CvdPath,
		int	mnMax,
		int	mnCon,
		const std::string &	si_or_lp,
		std::vector< realT > &	mags,
		int	lifetimeTrials,
		bool	writePSDs
	)

Todo:: should be able to scale the covar by r0, and possibly D

Definition at line 1514 of file fourierTemporalPSD.hpp.

References mx::sigproc::augment1SidedPSD(), mx::sigproc::augment1SidedPSDFreq(), backward, mx::ioutils::convertToString(), mx::sigproc::window::hann(), mx::improc::eigenCube< dataT >::image(), mx::ipc::ompLoopWatcher< _outputT, _printPretty, _printLoops, _printPercent, _printNLine, _time >::incrementAndOutputStatus(), mx::sigproc::normPSD(), mx::sigproc::psdVar(), mx::fits::fitsFile< dataT >::read(), mx::ioutils::readBinVector(), mx::sigproc::averagePeriodogram< realT >::resize(), mx::math::randomT< typeT, _ranengT, _randistT >::seed(), mx::math::six_fifths(), mx::sigproc::averagePeriodogram< realT >::size(), mx::math::vectorMeanSub(), mx::math::vectorVariance(), mx::sigproc::averagePeriodogram< realT >::win(), and mx::fits::fitsFile< dataT >::write().

◆ makePSDGrid()

template<typename realT , typename aosysT >

void mx::AO::analysis::fourierTemporalPSD< realT, aosysT >::makePSDGrid	(	const std::string &	dir,
		int	mnMax,
		realT	dFreq,
		realT	maxFreq,
		realT	fmax = `0`
	)

Calculate PSDs over a grid of spatial frequencies.

The grid of spatial frequencies is square, set by the maximum value of m and n.

The PSDs are written as mx::binVector binary files to a directory. We do not use FITS since this adds overhead and cfitisio handles parallelization poorly due to the limitation on number of created file pointers.

Parameters

[in]	dir	the directory for output of the PSDs
[in]	mnMax	the maximum value of m and n in the grid.
[in]	dFreq	the temporal frequency spacing.
[in]	maxFreq	the maximum temporal frequency to calculate
[in]	fmax	[optional] set the maximum temporal frequency for the calculation. The PSD is filled in with a -17/3 power law past this frequency. If 0, then it is taken to be 150 Hz + 2*fastestPeak(m,n).

Definition at line 696 of file fourierTemporalPSD.hpp.

References mx::ioutils::convertToString(), mx::ipc::ompLoopWatcher< _outputT, _printPretty, _printLoops, _printPercent, _printNLine, _time >::incrementAndOutputStatus(), mx::math::vectorScale(), and mx::ioutils::writeBinVector().

◆ multiLayerPSD()

template<typename realT , typename aosysT >

template<bool parallel>

int mx::AO::analysis::fourierTemporalPSD< realT, aosysT >::multiLayerPSD	(	std::vector< realT > &	PSD,
		std::vector< realT > &	freq,
		realT	m,
		realT	n,
		int	p,
		realT	fmax = `0`
	)

Calculate the temporal PSD for a Fourier mode in a multi-layer model.

Template Parameters

parallel controls whether layers are calculated in parallel. Default is true. Set to false if this is called inside a parallelized loop, as in makePSDGrid.

Todo:

implement error checking

handle reports of convergence failures form singleLayerPSD when implemented.

Parameters

[out]	PSD	the calculated PSD
[in]	freq	the populated temporal frequency grid defining at which frequencies the PSD is calculated
[in]	m	the first index of the spatial frequency
[in]	n	the second index of the spatial frequency
[in]	p	sets which mode is calculated (if basic modes, p = -1 for sine, p = +1 for cosine)
[in]	fmax	[optional] set the maximum temporal frequency for the calculation. The PSD is filled in with a -17/3 power law past this frequency. If 0, then it is taken to be 150 Hz + 2*fastestPeak(m,n).

Definition at line 680 of file fourierTemporalPSD.hpp.

◆ relTol() [1/2]

template<typename realT , typename aosysT >

realT mx::AO::analysis::fourierTemporalPSD< realT, aosysT >::relTol ( )

Get the current relative tolerance.

Returns: _relTol

Definition at line 477 of file fourierTemporalPSD.hpp.

◆ relTol() [2/2]

template<typename realT , typename aosysT >

void mx::AO::analysis::fourierTemporalPSD< realT, aosysT >::relTol ( realT rt )

Set relative tolerance.

Parameters

rt	is the new relative tolerance.

Definition at line 471 of file fourierTemporalPSD.hpp.

◆ singleLayerPSD()

template<typename realT , typename aosysT >

int mx::AO::analysis::fourierTemporalPSD< realT, aosysT >::singleLayerPSD	(	std::vector< realT > &	PSD,
		std::vector< realT > &	freq,
		realT	m,
		realT	n,
		int	layer_i,
		int	p,
		realT	fmax = `0`
	)

Calculate the temporal PSD for a Fourier mode for a single layer.

Todo:

implement error checking.

need a way to track convergence failures in integral without throwing an error.

need better handling of averaging for the -17/3 extension.

Parameters

[out]	PSD	the calculated PSD
[in]	freq	the populated temporal frequency grid defining the frequencies at which the PSD is calculated
[in]	m	the first index of the spatial frequency
[in]	n	the second index of the spatial frequency
[in]	layer_i	the index of the layer, for accessing the atmosphere parameters
[in]	p	sets which mode is calculated (if basic modes, p = -1 for sine, p = +1 for cosine)
[in]	fmax	[optional] set the maximum temporal frequency for the calculation. The PSD is filled in with a -17/3 power law past this frequency.