Afrho¶
- class sbpy.activity.Afrho(value, unit=None, dtype=None, copy=True)[source]¶
Bases:
sbpy.activity.dust.DustComaQuantity
Coma dust quantity for scattered light.
Afrho
objects behave likeQuantity
objects with units of length.- Parameters
Notes
Afρ is the product of dust albedo, dust filling factor, and circular aperture radius. It is nominally a constant for a steady-state coma in free expansion. See A’Hearn et al. (1984) for details.
References
A’Hearn et al. 1984, AJ 89, 579-591.
Examples
>>> from sbpy.activity import Afrho >>> print(Afrho(1000, 'cm')) 1000.0 cm
Methods Summary
from_fluxd
(wfb, fluxd, aper, eph, **kwargs)Initialize from spectral flux density.
to_fluxd
(wfb, aper, eph[, unit, phasecor, Phi])Express as spectral flux density in an observation.
to_phase
(to_phase, from_phase[, Phi])Scale to another phase angle.
Methods Documentation
- classmethod from_fluxd(wfb, fluxd, aper, eph, **kwargs)[source]¶
Initialize from spectral flux density.
- Parameters
- wfb
Quantity
,SpectralElement
, list Wavelengths, frequencies, bandpass, or list of bandpasses of the observation. Bandpasses require
synphot
.- fluxd
Quantity
Flux density per unit wavelength or frequency.
- aper
Quantity
orAperture
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
.
- wfb
- to_fluxd(wfb, aper, eph, unit=None, phasecor=False, Phi=None)[source]¶
Express as spectral flux density in an observation.
Assumes the small angle approximation.
- Parameters
- wfb
Quantity
,SpectralElement
, list Wavelengths, frequencies, bandpass, or list of bandpasses of the observation. Bandpasses require
synphot
. Ignored ifS
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’.
- unit
Unit
, string, optional The flux density unit for the output.
- phasecor: bool, optional
Scale the result by the phase function
Phi
, assumingAfrho
is quoted for 0° phase. Requires phase angle ineph
.- Phicallable, optional
Phase function, see
to_phase()
.- **kwargs
Keyword arguments for
observe
.
- wfb
- Returns
- fluxd
Quantity
Spectral flux density.
- fluxd
Examples
>>> from sbpy.activity import Afrho >>> import astropy.units as u >>> afrho = Afrho(1000 * u.cm) >>> wave = 0.55 * u.um >>> aper = 1 * u.arcsec >>> eph = dict(rh=1.5 * u.au, delta=1.0 * u.au) >>> fluxd = afrho.to_fluxd(wave, aper, eph) >>> print(fluxd) 6.730018324465526e-14 W / (m2 um)
With a phase correction: >>> eph[‘phase’] = 30 * u.deg >>> fluxd = afrho.to_fluxd(wave, aper, eph, phasecor=True) >>> print(fluxd) # doctest: +FLOAT_CMP 2.8017202649540757e-14 W / (m2 um)
In magnitudes through the Johnson V filter: >>> import sbpy.units as sbu >>> from sbpy.photometry import bandpass >>> bp = bandpass(‘Johnson V’) >>> fluxd = afrho.to_fluxd(bp, aper, eph, unit=sbu.JMmag, … phasecor=True) >>> print(fluxd) # doctest: +FLOAT_CMP 15.321242371548918 mag(JM)
- to_phase(to_phase, from_phase, Phi=None)[source]¶
Scale to another phase angle.
- Parameters
- to_phase
Quantity
New target phase angle.
- from_phase
Quantity
Current target phase angle.
- Phicallable, optional
Phase function, a callable object that takes a single parameter, phase angle as a
Quantity
, and returns a scale factor. Default isphase_HalleyMarcus
.
- to_phase
- Returns
- afrho
Afrho
The scaled Afρ quantity.
- afrho
Examples
>>> from sbpy.activity import Afrho >>> afrho = Afrho(10 * u.cm).to_phase(15 * u.deg, 0 * u.deg) >>> print(afrho) 5.87201 cm