FiveLayerAtmosphereDensityProfile#
- class ctapipe.atmosphere.FiveLayerAtmosphereDensityProfile(table: Table)[source]#
Bases:
AtmosphereDensityProfile
CORSIKA 5-layer atmosphere model
Layers 1-4 are modeled with:
\[T(h) = a_i + b_i \exp{-h/c_i}\]Layer 5 is modeled with:
\[T(h) = a_5 - b_5 \frac{h}{c_5}\]References
- [corsika-user] D. Heck and T. Pierog, “Extensive Air Shower Simulation with CORSIKA:
A User’s Guide”, 2021, Appendix F
Methods Summary
__call__
(height)from_array
(array[, metadata])construct from a 5x5 array as provided by eventio
from_table
(table)return a subclass of AtmosphereDensityProfile from a serialized table
height_from_overburden
(overburden)Get the height a.s.l. from the mass overburden in the atmosphere.
height_from_slant_depth
(slant_depth[, ...])Calculates height a.s.l. in the atmosphere from traversed slant depth
integral
(height)Integral of the profile along the height axis, i.e. the atmospheric depth \(X\).
peek
()Draw quick plot of profile
slant_depth_from_height
(height[, ...])Line-of-sight integral from height to infinity, along the direction specified by the zenith angle.
Methods Documentation
- classmethod from_array(array: ndarray, metadata: dict | None = None)[source]#
construct from a 5x5 array as provided by eventio
- classmethod from_table(table: Table)#
return a subclass of AtmosphereDensityProfile from a serialized table
- height_from_overburden(overburden) Quantity [source]#
- Get the height a.s.l. from the mass overburden in the atmosphere.
Inverse of the integral function
- Returns:
- u.Quantity[“m”]:
Height a.s.l. for given overburden
- height_from_slant_depth(slant_depth: ~astropy.units.quantity.Quantity, zenith_angle=<Quantity 0. deg>, output_units=Unit("m"))#
- Calculates height a.s.l. in the atmosphere from traversed slant depth
taking into account the shower zenith angle.
- Parameters:
- slant_depth: u.Quantity[“grammage”]
line-of-site distance from observer to point
- zenith_angle: u.Quantity[“angle”]
zenith angle of observation
- output_units: u.Unit
unit to output (must be convertible to m)
- integral(height) Quantity [source]#
Integral of the profile along the height axis, i.e. the atmospheric depth \(X\).
\[X(h) = \int_{h}^{\infty} \rho(h') dh'\]- Returns:
- u.Quantity[“g/cm2”]:
Integral of the density from height h to infinity
- peek()#
Draw quick plot of profile
- slant_depth_from_height(height: ~astropy.units.quantity.Quantity, zenith_angle=<Quantity 0. deg>, output_units=Unit("g / cm2"))#
Line-of-sight integral from height to infinity, along the direction specified by the zenith angle. This is sometimes called the slant depth. The atmosphere here is assumed to be Cartesian, the curvature of the Earth is not taken into account. This approximation breaks down for large zenith angles (>70 deg), in which case this function does not give correct results. Inverse of height_from_slant_depth function.
\[X(h, \Psi) = \int_{h}^{\infty} \rho(h') dh' / \cos{\Psi}\]- Parameters:
- height: u.Quantity[“length”]
height a.s.l. at which to start integral
- zenith_angle: u.Quantity[“angle”]
zenith angle of observation
- output_units: u.Unit
unit to output (must be convertible to g/cm2)