Atmospheric turbulence module (arte.atmo)

Introduction

arte.atmo provides several functions to perform computation related to the atmospheric turbulence and its optical effects.

The main class to represent the turbulence is cn2_profile().

Random phase screen can be generated with is phase_screen_generator()

Kolmogorov and Von Karman spectra are available von_karman_psd()

Covariances and cross-power-spectral-densities of Von Karman turbulence are computed in von_karman_covariance_calculator()

Submodules

cn2_profile module

This module manages atmospheric Cn2 profiles.

The class Cn2Profile is the base class: it allows to describe a profile as a set of layers each one with specified J, outer scale \(L_0\), wind speed and wind direction and located at a specified altitude h above the ground

Cn2Profile allows to specify the zenith angle z and the wavelength and

The following relations are used, where the index i identifies the i-th layer and the airmass is \(X=\sec z\).

\begin{eqnarray} k & = & \frac{2 \pi}{\lambda} \\ J_i & = & {C_n^2}_i \, dh_i \\ {r_0}_i & = & (0.423 \, X \, k^2 J_i)^{-3/5} \\ r_0 & = & \big( \sum{{{r_0}_i}^{-5/3}} \big)^{-3/5} \\ & = & \big( \sum{0.423 X k^2 J_i} \big)^{-3/5} \\ \theta_0 & = & \big( 2.914 k^2 X^{8/3} \sum{J_i h_i^{5/3}} \big)^{-3/5} \end{eqnarray}
class arte.atmo.cn2_profile.Cn2Profile(layersJs, layersL0, layersAltitude, layersWindSpeed, layersWindDirection)

Bases: object

Parameters:
  • layersJs (ndarray) – array of layers Js in meters**(1/3)
  • layersL0 (ndarray) – array of layers outer-scale L0 in meters
  • layersAltitude (ndarray) – array of layers altitude at Zenith in meters
  • layersWindSpeed (ndarray) – array of layers wind speed in meters per second
  • layersWindDirection (ndarray) – array of layers wind direction in degrees, clockwise from North

Every array must be 1D and have the same size, defining the number of layers of the profile

All parameters must be defined at zenith

DEFAULT_AIRMASS = 1.0
DEFAULT_LAMBDA = 5e-07
airmass()
Returns:airmass – airmass at specified zenith angle
Return type:float
fractional_j()
classmethod from_fractional_j(r0AtZenith, layersFractionalJ, layersL0, layersAltitude, layersWindSpeed, layersWindDirection)

Cn2 profile constructor from total r0 at zenith and fractional J of each layer

Parameters:
  • r0AtZenith (float) – overall r0 at zenith [m]
  • layersFractionalJ (ndarray) – array of J values for each layer. Array must sum up to 1
  • layersL0 (ndarray) – array of layers outer-scale L0 in meters
  • layersAltitude (ndarray) – array of layers altitude at Zenith in meters
  • layersWindSpeed (ndarray) – array of layers wind speed in meters per second
  • layersWindDirection (ndarray) – array of layers wind direction in degrees, clockwise from North

Every array must be 1D and have the same size, defining the number of layers of the profile

All parameters must be defined at zenith

classmethod from_r0s(layersR0, layersL0, layersAltitude, layersWindSpeed, layersWindDirection)

Cn2 profile constructor from r0 values of each layer

Parameters:
  • layersR0 (ndarray) – array of layers r0 in meters at 500nm
  • layersL0 (ndarray) – array of layers outer-scale L0 in meters
  • layersAltitude (ndarray) – array of layers altitude at Zenith in meters
  • layersWindSpeed (ndarray) – array of layers wind speed in meters per second
  • layersWindDirection (ndarray) – array of layers wind direction in degrees, clockwise from North

Every array must be 1D and have the same size, defining the number of layers of the profile

All parameters must be defined at zenith

layers_distance()
mean_wind_speed()
number_of_layers()
outer_scale()
r0()
Returns:r0 – Fried parameter at defined wavelength and zenith angle
Return type:Quantity equivalent to meters
r0s()
Returns:r0 – Fried parameter of each layer at defined wavelength and zenith angle
Return type:Quantity equivalent to meters of ndarray
seeing()
Returns:seeing – seeing value at specified lambda and zenith angle defined as 0.98 * lambda / r0
Return type:Quantity equivalent to arcsec
set_wavelength(wavelengthInMeters)
set_wind_direction(windDirectionInDeg)
set_wind_speed(windSpeed)
set_zenith_angle(zenithAngleInDeg)
tau0()
Returns:tau0 at specified lambda and zenith angle [sec]
Return type:tau0 (float)
theta0()
Returns:
isoplanatic angle at specified lambda and zenith
angle [arcsec]
Return type:theta0 (float)
wavelength()
wind_direction()
Returns:wind – wind direction of each layer [deg, clockwise from N]
Return type:Quantity containing an array.
wind_speed()
Returns:windspeed of each layer in m/s
Return type:(array)
zenith_angle()
static zenith_angle_to_airmass(zenithAngleInDeg)
class arte.atmo.cn2_profile.EsoEltProfiles

Bases: object

L0 = 25
classmethod Median(*args, **kwargs)
classmethod Q1()
classmethod Q2()
classmethod Q3()
classmethod Q4()
class arte.atmo.cn2_profile.MaoryStereoScidarProfiles2021

Bases: object

L0 = 25
classmethod P10()
classmethod P25()
classmethod P50()
classmethod P75()
classmethod P90()
class arte.atmo.cn2_profile.MaorySteroScidarProfiles

Bases: object

classmethod P10()
classmethod P25()
classmethod P50()
classmethod P75()
classmethod P90()
class arte.atmo.cn2_profile.MiscellaneusProfiles

Bases: object

classmethod ERIS()

From Doc. No.: VLT-SPE-ESO-11250-4110

classmethod LBT()

From G. Agapito, C. Arcidiacono, F. Quiros-Pacheco, S. Esposito, “Adaptive optics at short wavelengths - Expected performance and sky coverage of the FLAO system going toward visible wavelengths”, doi: 10.1007/s10686-014-9380-7

classmethod MaunaKea()

From Brent L. Ellerbroek, Francois J. Rigaut, “Scaling multiconjugate adaptive optics performance estimates to extremely large telescopes,” Proc. SPIE 4007, Adaptive Optical Systems Technology, (7 July 2000); doi: 10.1117/12.390314

phase_screen_generator module

class arte.atmo.phase_screen_generator.PhaseScreenGenerator(screenSizeInPixels, screenSizeInMeters, outerScaleInMeters, seed=None)

Bases: object

generate_normalized_phase_screens(numberOfScreens)
get_in_meters()
get_in_radians_at(wavelengthInMeters)
rescale_to(r0At500nm)

utils module

class arte.atmo.utils.Seeing(seeingAt500nm)

Bases: object

to_r0(wavelength=<Quantity 500. nm>)

von_karman_psd module

@author: giuliacarla

class arte.atmo.von_karman_psd.VonKarmanPsd(fried_param, outer_scale)

Bases: object

This class computes the spatial Power Spectral Density (PSD) of turbulent phase assuming the Von Karman spectrum. The PSD is obtained from the following expression:

\[ \begin{align}\begin{aligned}PSD(f,h) = c * r_0(h)^{-5/3} * ( f^2+ \frac{1}{L_0^2})^{-11/6}\\c = ( \frac{24}{5} \Gamma(6/5) )^{5/6} \frac{\Gamma(11/6)^2}{2 \pi^{11/3} } = 0.0228955\end{aligned}\end{align} \]
Parameters:
  • fried_param (float) – Fried parameter characterizing the atmospheric turbulence [m].
  • outer_scale (float) – Outer scale of the atmospheric turbulence [m].

Example

Compute the total variance of Kolmogorov over a 10m telescope and compare variance [rad2] with Noll(‘76): \(\Delta_1= 1.029 (D/r_0)^{5/3}\)

>>> R = 5
>>> r0 = 0.1
>>> L0 = np.inf
>>> psd = von_karman_psd.VonKarmanPsd(r0, L0)
>>> freqs = np.logspace(-8, 4, 1000)
>>> bess = scipy.special.jv(1, 2*np.pi*R*freqs)
>>> psdPistonRem = psd.spatial_psd(freqs) * (1 - (bess/(np.pi*R*freqs))**2)
>>> varInRad2 = np.trapz(psdPistonRem*2*np.pi*freqs, freqs)
>>> varInRad2Noll = 1.029*(2*R/r0)**(5./3)
>>> print("%g %g" % (varInRad2, varInRad2Noll))
2214.36 2216.91

or use shortcut function von_karman_psd.rms()

NUM_CONST = 0.02289558710855519
plot_von_karman_psd_vs_frequency(freqs, idx=None)
spatial_psd(freqs)

Spatial Power Spectral Density of Von Karman turbulence

Parameters:freqs (ndarray) – Spatial frequencies vector[m^-1].
Returns:psd – power spectral density computed at the specified frequencies
Return type:ndarray
arte.atmo.von_karman_psd.rms(diameter: Unit("m"), wavelength: Unit("nm"), fried_param: Unit("m"), outer_scale: Unit("m"), freqs=None)

Von Karman wavefront rms value over a circular aperture

Parameters:
  • diameter (Quantity equivalent to meter) – Aperture diameter
  • wavelength (Quantity equivalent to nanometer) – wavelength
  • fried_param (Quantity equivalent to meter) – Fried parameter r0 defined at the specified wavelength
  • outer_scale (Quantity equivalent to meter) – Outer scale L0. Use np.inf for Kolmogorov spectrum
Other Parameters:
 

freqs (array of Quantity equivalent to 1/meter) – spatial frequencies array. Default logspace(-8, 4, 1000) m^-1

Returns:

rms – wavefront rms for the specified von Karman turbulence

Return type:

Quantity equivalent to nm

von_karman_covariance_calculator module

@author: giuliacarla

class arte.atmo.von_karman_covariance_calculator.VonKarmanSpatioTemporalCovariance(source1, source2, aperture1, aperture2, cn2_profile, spat_freqs, logger=<Logger VK_COVARIANCE (WARNING)>)

Bases: object

Covariance of Von Karman atmospheric turbulence

This class computes the covariance and its Fourier transform, the Cross Power Spectral Density (CPSD), between Zernike coefficients describing the turbulence induced phase aberrations of two sources seen by two different circular apertures. The CPSD of the phase is also computed.

References

Plantet et al. (2020) - “Spatio-temporal statistics of the turbulent Zernike coefficients and piston-removed phase from two distinct beams”.

Plantet et al. (2018) - “LO WFS of MAORY: performance and sky coverage assessment.”

Whiteley et al. (1998) - “Temporal properties of the Zernike expansion coefficients of turbulence-induced phase aberrations for aperture and source motion”.

Parameters:
  • source1 (GuideSource) – Geometry of the first source. We consider rho as the angle in arcsec wrt the z-axis and theta as the angle in degrees wrt the x-axis. (e.g. source1 = arte.types.guide_source.GuideSource((1,90), 9e3)
  • source2 (GuideSource) – Geometry of the second source. Same conventions as source1.
  • aperture1 (CircularOpticalAperture) – Geometry of the first optical aperture. (e.g. aperture1 = arte.types.aperture.CircularOpticalAperture( 10, (0, 0, 0)))
  • aperture2 (CircularOpticalAperture) – Geometry of the second optical aperture.
  • cn2_profile (cn2_profile) – Cn2 profile. (e.g. cn2_eso = arte.atmo.cn2_profile.EsoEltProfiles.Q1() e.g. cn2_invented = arte.atmo.cn2_profile.Cn2Profile.from_r0s( [0.16], [25], [10e3], [0.1], [0]))
  • spat_freqs (numpy.ndarray) – Range of spatial frequencies that are used in Zernike covariance, Zernike CPSD and phase CPSD computation.
aperture1()
aperture2()
cn2_profile()
getGeneralPhaseCPSD(temp_freqs)
getGeneralZernikeCPSD(j, k, temp_freqs)

Return the generalized expression of Zernike CPSD that we get from ‘getZernikeCPSD’ function. This expression is needed when we have to integrate the Zernike CPSD in the temporal frequency range from -infinity to +infinity. Instead of this computation, we can obtain the same result performing the integral of the generalized Zernike CPSD in the temporal frequency range from 0 to infinity. This is what we need, for example, when we want to compare the Zernike covariance with the Zernike CPSD integrated in the temporal frequencies’ domain.

Parameters:
  • j (int or list) – Index of Zernike coefficients related to source1 on aperture1.
  • k (int or list) – Index of Zernike coefficients related to source2 on aperture2.
  • temp_freqs (numpy.ndarray) – Temporal frequencies array.
Returns:

cpsdTotal – General Zernike CPSD in [rad**2/Hz].

Return type:

Quantity

getPhaseCPSD(temp_freqs)

Return the Cross Power Spectral Density (CPSD) of the turbulent phase seen by aperture1 and aperture2 observing, respectively, source1 and source2. The CPSD is a function of temporal frequency.

Parameters:temp_freqs (numpy.ndarray) – Temporal frequencies array.
Returns:phaseCPSD – Phase CPSD in [rad**2/Hz].
Return type:Quantity
getPhaseCovariance()

Return the covariance between the phase seen from source1 on aperture1 and the phase seen from source2 on aperture2.

Returns:phaseCovariance – Covariance between phase1 and phase2 in [rad**2].
Return type:Quantity
getZernikeCPSD(j, k, temp_freqs)

Return the Cross Power Spectral Density (CPSD) of the Zernike coefficients with index j and k describing the phase seen on aperture1 and aperture2 observing, respectively, source1 and source2. The CPSD is a function of the temporal frequency.

Parameters:
  • j (int) – Index of the Zernike coefficient (related to source1 on aperture1).
  • k (int) – Index of the Zernike coefficient (related to source2 on aperture2).
  • temp_freqs (numpy.ndarray) – Temporal frequencies array in Hz.
Returns:

zernikeCPSD – Zernike CPSD or matrix of Zernike CPSDs in [rad**2/Hz].

Return type:

Quantity

getZernikeCovariance(j, k)

Return the covariance between two Zernike coefficients with index j and k describing the phase seen, respectively, on aperture1 from source1 and on aperture2 from source2.

Parameters:
  • j (int or list) – Index of Zernike coefficients related to source1 on aperture1.
  • k (int or list) – Index of Zernike coefficients related to source2 on aperture2.
Returns:

zernikeCovariance – Zernike covariance or covariance matrix (matrix of shape nxm if n and m are, respectively, the dimension of j and k) in [rad**2].

Return type:

Quantity

integrandOfGeneralPhaseCPSD(nLayer, temp_freqs)
integrandOfGeneralZernikeCPSD(j, k, nLayer, temp_freqs)
integrandOfPhaseCPSD(nLayer, temp_freqs)
integrandOfZernikeCPSD(j, k, nLayer, temp_freqs)
plotCPSD(cpsd, func_part, scale, legend)
setAperture1(ap1)
setAperture2(ap2)
setCn2Profile(cn2_profile)
setSource1(s1)
setSource2(s2)
setSpatialFrequencies(freqs)
source1()
source2()
spatial_frequencies()
temporal_frequencies()
useGPU()