mxlib
c++ tools for analyzing astronomical data and other tasks by Jared R. Males. [git repo]
mx::AO::sim::turbAtmosphere< _aoSystemT > Struct Template Reference

template<typename _aoSystemT>
struct mx::AO::sim::turbAtmosphere< _aoSystemT >

A turbulent atmosphere simulator.

Generate multi-layer phase screens using PSD filtering with the FFT. Manages the shifting and combination of the phase screens (without propagation).

The turbSubHarmonic class provides low-frequency sub-harmonics if desired

The turbLayer class manages the actual interpolation.

Todo:

add facility for switching interpolators

need to include wrap detection when using subharmonics (e.g. frames needs to be right)

Definition at line 82 of file turbAtmosphere.hpp.

#include <ao/sim/turbAtmosphere.hpp>

+ Inheritance diagram for mx::AO::sim::turbAtmosphere< _aoSystemT >:

Public Member Functions

int wfPS (realT ps)
 dummy function for simulatedAOSystem. Does nothing. More...
 
Construction and Configuration
 turbAtmosphere ()
 Default c'tor. More...
 
void setup (uint32_t wfSz, uint32_t buffSz, aoSystemT *aosys, uint32_t shLevel)
 Setup the overall atmosphere. More...
 
void wfSz (uint32_t ws)
 Set the wavefront size. More...
 
uint32_t wfSz ()
 Get the wavefront size. More...
 
void buffSz (uint32_t bs)
 Set the interpolation buffer size. More...
 
uint32_t buffSz ()
 Get the interpolation buffer size. More...
 
void aosys (aoSystemT *aos)
 Set the pointer to an AO system object. More...
 
aoSystemT * aosys ()
 Get the pointer to an AO system object. More...
 
void shLevel (uint32_t shl)
 Set the subharmonic level. More...
 
uint32_t shLevel ()
 Get the subharmonic level. More...
 
void shPreCalc (bool shp)
 Set whether subharmonic modes are pre-calculated. More...
 
bool shPreCalc ()
 Get whether subharmonic modes are pre-calculated. More...
 
void retain (bool rtn)
 Set whether memory for screen generation is retained. More...
 
bool retain ()
 Get whether memory for screen generation is retained. More...
 
void forceGen (bool fg)
 Set whether new screen generation is forced. More...
 
bool forceGen ()
 Get whether new screen generation is forced. More...
 
void dataDir (const std::string &dd)
 Set the data directory for saving phase screens. More...
 
std::string dataDir ()
 Get the data directory for saving phase screens. More...
 
void setLayers (const std::vector< size_t > &scrnSz)
 Setup the layers and prepare them for phase screen generation. More...
 
void setLayers (const size_t scrnSz)
 Setup the layers and prepare them for phase screen generation. More...
 
uint32_t nLayers ()
 Get the number of layers. More...
 
turbLayer< aoSystemT > & layer (uint32_t n)
 Get a layer. More...
 
math::normDistT< realT > & normVar ()
 Get the random number generator. More...
 
Screen Generation
void genLayers ()
 Generate all phase screens. More...
 
- Public Member Functions inherited from mx::base::changeable
void changed ()
 Increment the counter. More...
 
changeT change ()
 Get the value of the counter. More...
 

Protected Attributes

Configuration Data
uint32_t m_wfSz {0}
 Size of the wavefront in pixels. More...
 
uint32_t m_buffSz {0}
 Buffer to apply around wavefront for interpolation. More...
 
aoSystemT * m_aosys {nullptr}
 The AO system object containing all system parameters. See aoSystem. More...
 
uint32_t m_shLevel {0}
 
bool m_shPreCalc {true}
 
bool m_retain {false}
 
bool m_forceGen {true}
 
std::string m_dataDir
 
imageT * m_pupil {0}
 A pointer to the pupil mask. More...
 
realT m_timeStep {0}
 Length of each iteration, in seconds. More...
 
Internal State
std::vector< turbLayer< aoSystemT > > m_layers
 Vector of turbulent layers. More...
 
std::vector< turbSubHarmonic< turbAtmosphere > > m_subHarms
 Sub-harmonic layers. More...
 
size_t m_frames {0}
 Length of the turbulence sequence. More...
 
int m_nWf {0}
 Number of iterations which have occurred. More...
 
math::normDistT< realT > m_normVar
 Normal random deviate generator. This seeded in the constructor. More...
 
improc::eigenImage< realT > m_phaseWork
 

Additional Inherited Members

- Public Types inherited from mx::base::changeable
typedef uint64_t changeT
 The integer type of the counter. More...
 

Constructor & Destructor Documentation

◆ turbAtmosphere()

template<typename aoSystemT >
mx::AO::sim::turbAtmosphere< aoSystemT >::turbAtmosphere

Default c'tor.

Definition at line 322 of file turbAtmosphere.hpp.

Member Function Documentation

◆ aosys() [1/2]

template<typename aoSystemT >
aoSystemT * mx::AO::sim::turbAtmosphere< aoSystemT >::aosys

Get the pointer to an AO system object.

Returns
the current value of m_aosys

Definition at line 386 of file turbAtmosphere.hpp.

◆ aosys() [2/2]

template<typename aoSystemT >
void mx::AO::sim::turbAtmosphere< aoSystemT >::aosys ( aoSystemT *  aos)

Set the pointer to an AO system object.

see m_aosys

Parameters
[in]aosthe new pointer to an AO system object

Definition at line 376 of file turbAtmosphere.hpp.

◆ buffSz() [1/2]

template<typename aoSystemT >
uint32_t mx::AO::sim::turbAtmosphere< aoSystemT >::buffSz

Get the interpolation buffer size.

Returns
the current value of m_buffSz

Definition at line 370 of file turbAtmosphere.hpp.

◆ buffSz() [2/2]

template<typename aoSystemT >
void mx::AO::sim::turbAtmosphere< aoSystemT >::buffSz ( uint32_t  bs)

Set the interpolation buffer size.

see m_buffSz

Parameters
[in]bsthe new buffer size

Definition at line 359 of file turbAtmosphere.hpp.

◆ dataDir() [1/2]

template<typename aoSystemT >
std::string mx::AO::sim::turbAtmosphere< aoSystemT >::dataDir

Get the data directory for saving phase screens.

Returns
the current value of m_dataDir

Definition at line 468 of file turbAtmosphere.hpp.

◆ dataDir() [2/2]

template<typename aoSystemT >
void mx::AO::sim::turbAtmosphere< aoSystemT >::dataDir ( const std::string &  dd)

Set the data directory for saving phase screens.

see m_dataDir

Parameters
[in]ddthe new data directory for phase screen saving

Definition at line 457 of file turbAtmosphere.hpp.

◆ forceGen() [1/2]

template<typename aoSystemT >
bool mx::AO::sim::turbAtmosphere< aoSystemT >::forceGen

Get whether new screen generation is forced.

Returns
the current value of m_forceGen

Definition at line 451 of file turbAtmosphere.hpp.

◆ forceGen() [2/2]

template<typename aoSystemT >
void mx::AO::sim::turbAtmosphere< aoSystemT >::forceGen ( bool  fg)

Set whether new screen generation is forced.

see m_forceGen

Parameters
[in]fgthe new value of the flag controlling whether screen generation is forced

Definition at line 440 of file turbAtmosphere.hpp.

◆ genLayers()

template<typename aoSystemT >
void mx::AO::sim::turbAtmosphere< aoSystemT >::genLayers

Generate all phase screens.

Loads them from disk if possible and m_forceGen == false.

Deallocates working memory when done, unless m_retain == true.

Definition at line 556 of file turbAtmosphere.hpp.

References mx::fits::fitsFile< dataT >::read(), and mx::fits::fitsFile< dataT >::write().

◆ layer()

template<typename aoSystemT >
turbLayer< aoSystemT > & mx::AO::sim::turbAtmosphere< aoSystemT >::layer ( uint32_t  n)

Get a layer.

Returns
a reference to on entry in m_layers

Definition at line 539 of file turbAtmosphere.hpp.

◆ nLayers()

template<typename aoSystemT >
uint32_t mx::AO::sim::turbAtmosphere< aoSystemT >::nLayers

Get the number of layers.

Returns
the size of the layers vector m_layers

Definition at line 533 of file turbAtmosphere.hpp.

◆ normVar()

template<typename aoSystemT >
math::normDistT< typename aoSystemT::realT > & mx::AO::sim::turbAtmosphere< aoSystemT >::normVar

Get the random number generator.

Returns
a reference to m_normVar

Definition at line 550 of file turbAtmosphere.hpp.

◆ retain() [1/2]

template<typename aoSystemT >
bool mx::AO::sim::turbAtmosphere< aoSystemT >::retain

Get whether memory for screen generation is retained.

Returns
the current value of m_retain

Definition at line 434 of file turbAtmosphere.hpp.

◆ retain() [2/2]

template<typename aoSystemT >
void mx::AO::sim::turbAtmosphere< aoSystemT >::retain ( bool  rtn)

Set whether memory for screen generation is retained.

see m_retain

Parameters
[in]rtnthe new value of the flag controlling whether memory is retained

Definition at line 424 of file turbAtmosphere.hpp.

◆ setLayers() [1/2]

template<typename aoSystemT >
void mx::AO::sim::turbAtmosphere< aoSystemT >::setLayers ( const size_t  scrnSz)

Setup the layers and prepare them for phase screen generation.

This sets all layer phase screens to be the same size.

Parameters
[in]scrnSzthe screen sizes for all layers.

Definition at line 520 of file turbAtmosphere.hpp.

◆ setLayers() [2/2]

template<typename aoSystemT >
void mx::AO::sim::turbAtmosphere< aoSystemT >::setLayers ( const std::vector< size_t > &  scrnSz)

Setup the layers and prepare them for phase screen generation.

The number of entries in scrnSz must mach the number of layers in the atmosphere of m_aosys.

Parameters
[in]scrnSzthe screen sizes for each layer.

Definition at line 474 of file turbAtmosphere.hpp.

◆ setup()

template<typename aoSystemT >
void mx::AO::sim::turbAtmosphere< aoSystemT >::setup ( uint32_t  wfSz,
uint32_t  buffSz,
aoSystemT *  aosys,
uint32_t  shLevel 
)

Setup the overall atmosphere.

Parameters
[in]wfSzThe size of the wavefront in pixels.
[in]buffSzThe size of the interpolation buffer to use.
[in]aosysPointer to an AO System. See m_aosys.
[in]shLevelnumber of subharmonic levels to use. 0 turns off subharmonics. See m_shLevel.

Definition at line 328 of file turbAtmosphere.hpp.

◆ shLevel() [1/2]

template<typename aoSystemT >
uint32_t mx::AO::sim::turbAtmosphere< aoSystemT >::shLevel

Get the subharmonic level.

Returns
the current value of m_shLevel

Definition at line 402 of file turbAtmosphere.hpp.

◆ shLevel() [2/2]

template<typename aoSystemT >
void mx::AO::sim::turbAtmosphere< aoSystemT >::shLevel ( uint32_t  shl)

Set the subharmonic level.

see m_shLevel

Parameters
[in]shlthe new subharmonic level

Definition at line 392 of file turbAtmosphere.hpp.

◆ shPreCalc() [1/2]

template<typename aoSystemT >
bool mx::AO::sim::turbAtmosphere< aoSystemT >::shPreCalc

Get whether subharmonic modes are pre-calculated.

Returns
the current value of m_shPreCalc

Definition at line 418 of file turbAtmosphere.hpp.

◆ shPreCalc() [2/2]

template<typename aoSystemT >
void mx::AO::sim::turbAtmosphere< aoSystemT >::shPreCalc ( bool  shp)

Set whether subharmonic modes are pre-calculated.

see m_shPreCalc

Parameters
[in]shpthe new value of flag controlling subharmonic mode precalculation

Definition at line 408 of file turbAtmosphere.hpp.

◆ wfPS()

template<typename aoSystemT >
int mx::AO::sim::turbAtmosphere< aoSystemT >::wfPS ( realT  ps)

dummy function for simulatedAOSystem. Does nothing.

Definition at line 687 of file turbAtmosphere.hpp.

◆ wfSz() [1/2]

template<typename aoSystemT >
uint32_t mx::AO::sim::turbAtmosphere< aoSystemT >::wfSz

Get the wavefront size.

Returns
the current value of m_wfSz

Definition at line 353 of file turbAtmosphere.hpp.

◆ wfSz() [2/2]

template<typename aoSystemT >
void mx::AO::sim::turbAtmosphere< aoSystemT >::wfSz ( uint32_t  ws)

Set the wavefront size.

see m_wfSz

Parameters
[in]wsthe new wavefront size

Definition at line 341 of file turbAtmosphere.hpp.

Member Data Documentation

◆ m_aosys

template<typename _aoSystemT >
aoSystemT* mx::AO::sim::turbAtmosphere< _aoSystemT >::m_aosys {nullptr}
protected

The AO system object containing all system parameters. See aoSystem.

Definition at line 101 of file turbAtmosphere.hpp.

◆ m_buffSz

template<typename _aoSystemT >
uint32_t mx::AO::sim::turbAtmosphere< _aoSystemT >::m_buffSz {0}
protected

Buffer to apply around wavefront for interpolation.

Definition at line 99 of file turbAtmosphere.hpp.

◆ m_dataDir

template<typename _aoSystemT >
std::string mx::AO::sim::turbAtmosphere< _aoSystemT >::m_dataDir
protected

Specifies a path where phase screens are stored. Leaving this unset or "" is equivalent to m_forceGen== true.

Definition at line 124 of file turbAtmosphere.hpp.

◆ m_forceGen

template<typename _aoSystemT >
bool mx::AO::sim::turbAtmosphere< _aoSystemT >::m_forceGen {true}
protected

Force generation of new screens if true. Note that this class does not presently do any checks of saved screens to verify they match the current configuration. Set to true with caution.

Todo:
need check for changes in parameters to decide if saved frames are valid per layer

Definition at line 117 of file turbAtmosphere.hpp.

◆ m_frames

template<typename _aoSystemT >
size_t mx::AO::sim::turbAtmosphere< _aoSystemT >::m_frames {0}
protected

Length of the turbulence sequence.

Definition at line 141 of file turbAtmosphere.hpp.

◆ m_layers

template<typename _aoSystemT >
std::vector<turbLayer<aoSystemT> > mx::AO::sim::turbAtmosphere< _aoSystemT >::m_layers
protected

Vector of turbulent layers.

Definition at line 137 of file turbAtmosphere.hpp.

◆ m_normVar

template<typename _aoSystemT >
math::normDistT<realT> mx::AO::sim::turbAtmosphere< _aoSystemT >::m_normVar
protected

Normal random deviate generator. This seeded in the constructor.

Definition at line 145 of file turbAtmosphere.hpp.

◆ m_nWf

template<typename _aoSystemT >
int mx::AO::sim::turbAtmosphere< _aoSystemT >::m_nWf {0}
protected

Number of iterations which have occurred.

Definition at line 143 of file turbAtmosphere.hpp.

◆ m_pupil

template<typename _aoSystemT >
imageT* mx::AO::sim::turbAtmosphere< _aoSystemT >::m_pupil {0}
protected

A pointer to the pupil mask.

Definition at line 128 of file turbAtmosphere.hpp.

◆ m_retain

template<typename _aoSystemT >
bool mx::AO::sim::turbAtmosphere< _aoSystemT >::m_retain {false}
protected

Whether or not to retain working memory after screen generation. One would set to true if many screens will be calculated in monte carlo fashion. For a standard case of one set of screens being shifted by wind velocity, this should be false to minimize memory use.

Definition at line 111 of file turbAtmosphere.hpp.

◆ m_shLevel

template<typename _aoSystemT >
uint32_t mx::AO::sim::turbAtmosphere< _aoSystemT >::m_shLevel {0}
protected

Number of subharmonic levels. 0 means no subharmonics. Generally only 1 level is needed to obtain good Tip/Tilt results. See turbSubHarmonic.

Definition at line 103 of file turbAtmosphere.hpp.

◆ m_shPreCalc

template<typename _aoSystemT >
bool mx::AO::sim::turbAtmosphere< _aoSystemT >::m_shPreCalc {true}
protected

Whether or not to pre-calculate the subharmonic modes. This is a trade of speed vs. memory use. See turbSubHarmonic.

Definition at line 107 of file turbAtmosphere.hpp.

◆ m_subHarms

template<typename _aoSystemT >
std::vector<turbSubHarmonic<turbAtmosphere> > mx::AO::sim::turbAtmosphere< _aoSystemT >::m_subHarms
protected

Sub-harmonic layers.

Definition at line 139 of file turbAtmosphere.hpp.

◆ m_timeStep

template<typename _aoSystemT >
realT mx::AO::sim::turbAtmosphere< _aoSystemT >::m_timeStep {0}
protected

Length of each iteration, in seconds.

Definition at line 130 of file turbAtmosphere.hpp.

◆ m_wfSz

template<typename _aoSystemT >
uint32_t mx::AO::sim::turbAtmosphere< _aoSystemT >::m_wfSz {0}
protected

Size of the wavefront in pixels.

Definition at line 97 of file turbAtmosphere.hpp.


The documentation for this struct was generated from the following file: