Skip to content

Commit

Permalink
Merge pull request #200 from mkelley/dust-and-dataclasses
Browse files Browse the repository at this point in the history
Dust and dataclasses
  • Loading branch information
mkelley committed Jul 22, 2019
2 parents 6e28595 + e47b4e7 commit 9255705
Show file tree
Hide file tree
Showing 3 changed files with 23 additions and 94 deletions.
100 changes: 17 additions & 83 deletions sbpy/activity/dust.py
Original file line number Diff line number Diff line change
Expand Up @@ -16,16 +16,14 @@
-------
Afrho - Coma dust quantity for scattered light.
Efrho - Coma dust quantity for thermal emission.
Syndynes - Dust dynamical model for zero-ejection velocities.
"""

__all__ = [
'phase_HalleyMarcus',
'Afrho',
'Efrho',
'Syndynes'
'Efrho'
]


Expand All @@ -35,15 +33,20 @@
import numpy as np
import astropy.units as u

from astropy.utils.exceptions import AstropyWarning
from .. import bib
from ..calib import Sun
from ..spectroscopy import BlackbodySource
from ..data import Ephem
from .. import data as sbd
from .. import exceptions as sbe
from ..spectroscopy.sources import SinglePointSpectrumError
from .core import Aperture, rho_as_length


@bib.cite({
'Halley-Marcus phase function': '2011AJ....141..177S',
'Halley phase function': '1998Icar..132..397S',
'Marcus phase function': '2007ICQ....29...39M'
})
def phase_HalleyMarcus(phase):
"""Halley-Marcus composite dust phase function.
Expand Down Expand Up @@ -105,15 +108,6 @@ def phase_HalleyMarcus(phase):
except ImportError:
scipy = None

bib.register(
'activity.dust.phase_HalleyMarcus',
{
'Halley-Marcus phase function': '2011AJ....141..177S',
'Halley phase function': '1998Icar..132..397S',
'Marcus phase function': '2007ICQ....29...39M'
}
)

th = np.arange(181)
ph = np.array(
[1.0000e+00, 9.5960e-01, 9.2170e-01, 8.8590e-01,
Expand Down Expand Up @@ -168,7 +162,7 @@ def phase_HalleyMarcus(phase):
if scipy:
Phi = splev(_phase, splrep(th, ph))
else:
warn(AstropyWarning(
warn(sbe.OptionalPackageUnavailable(
'scipy is not present, using linear interpolation.'))
Phi = np.interp(_phase, th, ph)

Expand Down Expand Up @@ -199,6 +193,7 @@ def __new__(cls, value, unit=None, dtype=None, copy=None):
def from_fluxd(cls, wfb, fluxd, aper, eph, **kwargs):
"""Initialize from spectral flux density.
Parameters
----------
wfb : `~astropy.units.Quantity`, `~synphot.SpectralElement`, list
Expand All @@ -214,8 +209,8 @@ def from_fluxd(cls, wfb, fluxd, aper, eph, **kwargs):
or angular units), or as an `~sbpy.activity.Aperture`.
eph: dictionary-like, `~sbpy.data.Ephem`
Ephemerides of the comet, requires heliocentric distance
and observer-comet distance.
Ephemerides of the comet. Required fields: 'rh', 'delta'.
Optional: 'phase'.
**kwargs
Keyword arguments for `~to_fluxd`.
Expand All @@ -232,6 +227,7 @@ def from_fluxd(cls, wfb, fluxd, aper, eph, **kwargs):

return coma

@sbd.dataclass_input(eph=sbd.Ephem)
def to_fluxd(self, wfb, aper, eph, unit=None, **kwargs):
"""Express as spectral flux density in an observation.
Expand All @@ -250,8 +246,8 @@ def to_fluxd(self, wfb, aper, eph, unit=None, **kwargs):
or angular units), or as an sbpy `~sbpy.activity.Aperture`.
eph: dictionary-like, `~sbpy.data.Ephem`
Ephemerides of the comet, requires heliocentric distance
and observer-comet distance.
Ephemerides of the comet. Required fields: 'rh', 'delta'.
Optional: 'phase'.
unit : `~astropy.units.Unit`, string, optional
The flux density unit for the output.
Expand All @@ -261,10 +257,6 @@ def to_fluxd(self, wfb, aper, eph, unit=None, **kwargs):
# This method handles the geometric quantities. Sub-classes
# will handle the photometric quantities in `_source_fluxd`.

# validate eph
if not isinstance(eph, Ephem):
eph = Ephem.from_dict(eph)

# rho = effective circular aperture radius at the distance of
# the comet. Keep track of array dimensionality as Ephem
# objects can needlessly increase the number of dimensions.
Expand Down Expand Up @@ -443,11 +435,9 @@ def to_fluxd(self, wfb, aper, eph, unit=None, phasecor=False,
"""

@bib.cite({'model': '1984AJ.....89..579A'})
def _source_fluxd(self, wfb, eph, unit=None, phasecor=False,
Phi=None, **kwargs):
bib.register('activity.dust.Afrho', {
'model': '1984AJ.....89..579A'})

# get solar flux density
sun = Sun.from_default()
try:
Expand Down Expand Up @@ -617,10 +607,8 @@ def to_fluxd(self, wfb, aper, eph, unit=None, Tscale=1.1,
"""

@bib.cite({'model': '2013Icar..225..475K'})
def _source_fluxd(self, wfb, eph, unit=None, Tscale=1.1, T=None, B=None):
bib.register('activity.dust.Efrho',
{'model': '2013Icar..225..475K'})

if T is None:
T = Tscale * 278 / np.sqrt(eph['rh'].to('au').value)

Expand All @@ -639,57 +627,3 @@ def _source_fluxd(self, wfb, eph, unit=None, Tscale=1.1, T=None, B=None):
'flux density, e.g., W/m2/μm or W/m2/Hz')

return B


class Syndynes:
"""Syndyne and synchrone generator.
Parameters
----------
orbit : `~sbpy.data.Orbit`
Comet orbital parameters.
date : `~astropy.time.Time`
Date of observation.
location : various, optional
IAU observatory code as a string; a location on the Earth as
an `~astropy.coordinates.EarthLocation`; or, a location in the
Solar System via :mod"`~astropy.coordinates` objects, e.g.,
`~astropy.coordinates.HeliocentricTrueEcliptic`.
"""

def __init__(self, orb, date, location='500'):
self.orb = orb
self.date = date
raise NotImplementedError

def plot_syndynes(self, beta, age):
"""Plot syndynes for an observer.
Parameters
----------
beta: array-like
Test particle β values.
age: `~astropy.units.Quantity`
Test particle ages.
Returns
-------
ax: `~matplotlib.pyplot.Axes`
Examples
--------
TBD
not yet implemented
"""

pass
15 changes: 5 additions & 10 deletions sbpy/activity/gas/core.py
Original file line number Diff line number Diff line change
Expand Up @@ -15,7 +15,7 @@
Classes
-------
Activity - Abstract base class for gas coma models.
GasComa - Abstract base class for gas coma models.
Haser - Haser coma model for gas (Haser 1957).
Vectorial - Vectorial coma model for gas (Festou 1981).
Expand Down Expand Up @@ -48,7 +48,7 @@
from astropy.utils.exceptions import AstropyWarning
from ... import bib
from .. core import (Aperture, RectangularAperture, GaussianAperture,
AnnularAperture, CircularAperture)
AnnularAperture, CircularAperture)
from .. core import rho_as_length


Expand Down Expand Up @@ -600,11 +600,10 @@ class Haser(GasComa):
"""

@bib.cite({'model': '1957BSRSL..43..740H'})
def __init__(self, Q, v, parent, daughter=None):
super().__init__(Q, v)

bib.register('activity.gas.Haser', {'model': '1957BSRSL..43..740H'})

if not parent.unit.is_equivalent(u.m):
raise ValueError('parent must have units of length')
self.parent = parent
Expand Down Expand Up @@ -643,10 +642,8 @@ def _K1(self, x):
raise AstropyWarning('scipy is not present, cannot continue.')
return special.k1(x.decompose().value)

@bib.cite({'model': '1978Icar...35..360N'})
def column_density(self, rho, eph=None):
bib.register('activity.gas.Haser.column_density',
{'model': '1978Icar...35..360N'})

r = rho_as_length(rho, eph=eph)
x = 0 if self.parent is None else (r / self.parent).decompose()
y = 0 if self.daughter is None else (r / self.daughter).decompose()
Expand All @@ -662,10 +659,8 @@ def column_density(self, rho, eph=None):
return sigma.decompose()
column_density.__doc__ = GasComa.column_density.__doc__

@bib.cite({'model': '1978Icar...35..360N'})
def total_number(self, aper, eph=None):
bib.register('activity.gas.Haser.total_number',
{'model': '1978Icar...35..360N'})

# Inspect aper and handle as appropriate
if isinstance(aper, Aperture):
aper = aper.as_length(eph)
Expand Down
2 changes: 1 addition & 1 deletion sbpy/data/__init__.py
Original file line number Diff line number Diff line change
Expand Up @@ -94,7 +94,7 @@ class Conf():
'provenance': ['ephem', 'obs'],
'dimension': 'angle'},
{'description': 'Solar Phase Angle',
'fieldnames': ['alpha', 'phaseangle', 'Phase'],
'fieldnames': ['alpha', 'phaseangle', 'Phase', 'phase'],
'provenance': ['ephem', 'obs'],
'dimension': 'angle'},
{'description': 'Solar Elongation Angle',
Expand Down

0 comments on commit 9255705

Please sign in to comment.