Efrho

class sbpy.activity.Efrho(value, unit=None, dtype=None, copy=True)[source]

Bases: DustComaQuantity

Coma dust quantity for thermal emission.

Efrho behave like Quantity objects with units of length.

Parameters:
valuenumber, Quantity

The value(s).

unitstr, Unit, optional

The unit of the input value. Strings must be parseable by units package.

dtypedtype, optional

See Quantity.

copybool, optional

See Quantity.

Notes

εfρ is the product of dust emissivity, dust filling factor, and circular aperture radius. It is nominally a constant for a steady-state coma in free expansion, and is the thermal emission equivalent for the Afρ quantity. See Kelley et al. (2013) for details.

References

A’Hearn et al. 1984, AJ 89, 579-591. Kelley et al. 2013, Icarus 225, 475-494.

Examples

>>> from sbpy.activity import Efrho
>>> print(Efrho(1000.0, 'cm'))
1000.0 cm

Methods Summary

from_fluxd(wfb, fluxd, aper, eph, **kwargs)

Initialize from spectral flux density.

to_fluxd(wfb, aper, eph[, unit, Tscale, T, B])

Express as spectral flux density in an observation.

Methods Documentation

classmethod from_fluxd(wfb, fluxd, aper, eph, **kwargs)[source]

Initialize from spectral flux density.

Parameters:
wfbQuantity, SpectralElement, list

Wavelengths, frequencies, bandpass, or list of bandpasses of the observation. Bandpasses require synphot.

fluxdQuantity

Flux density per unit wavelength or frequency.

aperQuantity or Aperture

Aperture of the observation as a circular radius (length or angular units), or as an Aperture.

eph: dictionary-like, `~sbpy.data.Ephem`

Ephemerides of the comet. Required fields: ‘rh’, ‘delta’. Optional: ‘phase’.

**kwargs

Keyword arguments for to_fluxd.

to_fluxd(wfb, aper, eph, unit=None, Tscale=1.1, T=None, B=None)[source]

Express as spectral flux density in an observation.

Assumes the small angle approximation.

Parameters:
wfbQuantity, SpectralElement, list

Wavelengths, frequencies, bandpass, or list of bandpasses of the observation. Bandpasses require synphot. Ignored if S is provided.

aper: `~astropy.units.Quantity`, `~sbpy.activity.Aperture`

Aperture of the observation as a circular radius (length or angular units), or as an sbpy Aperture.

eph: dictionary-like, `~sbpy.data.Ephem`

Ephemerides of the comet. Required fields: ‘rh’, ‘delta’. Optional: ‘phase’.

unitUnit, string, optional

The flux density unit for the output.

Tscalefloat, optional

Scale factor for blackbody in LTE with sunlight. Ignored if T or B is provided.

TQuantity, optional

Blackbody temperature. Ignored if B is provided.

BQuantity, optional

Observed spectral flux density from a blackbody sphere, i.e., pi * Planck function. Overrides T and Tscale.

Returns:
fluxdQuantity

Spectral flux density.

Examples

>>> from sbpy.activity import Efrho
>>> import astropy.units as u
>>> efrho = Efrho(120.0, 'cm')
>>> freq = 15.8 * u.um
>>> aper = 11.1 * u.arcsec
>>> eph = {'rh': 4.42 * u.au, 'delta': 4.01 * u.au}
>>> fluxd = efrho.to_fluxd(freq, aper, eph=eph, unit='Jy')
>>> print(fluxd)                           
0.006519545281786034 Jy