Models and Fitting (astropy.modeling)

Introduction

astropy.modeling provides a framework for representing models and performing model evaluation and fitting. It currently supports 1-D and 2-D models and fitting with parameter constraints.

It is designed to be easily extensible and flexible. Models do not reference fitting algorithms explicitly and new fitting algorithms may be added without changing the existing models (though not all models can be used with all fitting algorithms due to constraints such as model linearity).

The goal is to eventually provide a rich toolset of models and fitters such that most users will not need to define new model classes, nor special purpose fitting routines (while making it reasonably easy to do when necessary).

Note

astropy.modeling is currently a work-in-progress, and thus it is likely there will still be API changes in later versions of Astropy. Backwards compatibility support between versions will still be maintained as much as possible, but new features and enhancements are coming in future versions. If you have specific ideas for how it might be improved, feel free to let us know on the astropy-dev mailing list or at http://feedback.astropy.org

Getting started

The examples here use the predefined models and assume the following modules have been imported:

>>> import numpy as np
>>> from astropy.modeling import models, fitting

Using Models

The astropy.modeling package defines a number of models that are collected under a single namespace as astropy.modeling.models. Models behave like parametrized functions:

>>> from astropy.modeling import models
>>> g = models.Gaussian1D(amplitude=1.2, mean=0.9, stddev=0.5)
>>> print(g)
Model: Gaussian1D
Inputs: ('x',)
Outputs: ('y',)
Model set size: 1
Parameters:
    amplitude mean stddev
    --------- ---- ------
          1.2  0.9    0.5

Model parameters can be accessed as attributes:

>>> g.amplitude
Parameter('amplitude', value=1.2)
>>> g.mean
Parameter('mean', value=0.9)
>>> g.stddev  
Parameter('stddev', value=0.5, bounds=(1.1754943508222875e-38, None))

and can also be updated via those attributes:

>>> g.amplitude = 0.8
>>> g.amplitude
Parameter('amplitude', value=0.8)

Models can be evaluated by calling them as functions:

>>> g(0.1)
0.22242984036255528
>>> g(np.linspace(0.5, 1.5, 7))  
array([0.58091923, 0.71746405, 0.7929204 , 0.78415894, 0.69394278,
       0.54952605, 0.3894018 ])

As the above example demonstrates, in general most models evaluate array-like inputs according to the standard Numpy broadcasting rules for arrays.

Models can therefore already be useful to evaluate common functions, independently of the fitting features of the package.

Simple 1-D model fitting

In this section, we look at a simple example of fitting a Gaussian to a simulated dataset. We use the Gaussian1D and Trapezoid1D models and the LevMarLSQFitter fitter to fit the data:

import numpy as np
import matplotlib.pyplot as plt
from astropy.modeling import models, fitting

# Generate fake data
np.random.seed(0)
x = np.linspace(-5., 5., 200)
y = 3 * np.exp(-0.5 * (x - 1.3)**2 / 0.8**2)
y += np.random.normal(0., 0.2, x.shape)

# Fit the data using a box model
t_init = models.Trapezoid1D(amplitude=1., x_0=0., width=1., slope=0.5)
fit_t = fitting.LevMarLSQFitter()
t = fit_t(t_init, x, y)

# Fit the data using a Gaussian
g_init = models.Gaussian1D(amplitude=1., mean=0, stddev=1.)
fit_g = fitting.LevMarLSQFitter()
g = fit_g(g_init, x, y)

# Plot the data with the best-fit model
plt.figure(figsize=(8,5))
plt.plot(x, y, 'ko')
plt.plot(x, t(x), label='Trapezoid')
plt.plot(x, g(x), label='Gaussian')
plt.xlabel('Position')
plt.ylabel('Flux')
plt.legend(loc=2)

(png, svg, pdf)

../_images/index-11.png

As shown above, once instantiated, the fitter class can be used as a function that takes the initial model (t_init or g_init) and the data values (x and y), and returns a fitted model (t or g).

Simple 2-D model fitting

Similarly to the 1-D example, we can create a simulated 2-D data dataset, and fit a polynomial model to it. This could be used for example to fit the background in an image.

import warnings
import numpy as np
import matplotlib.pyplot as plt
from astropy.modeling import models, fitting

# Generate fake data
np.random.seed(0)
y, x = np.mgrid[:128, :128]
z = 2. * x ** 2 - 0.5 * x ** 2 + 1.5 * x * y - 1.
z += np.random.normal(0., 0.1, z.shape) * 50000.

# Fit the data using astropy.modeling
p_init = models.Polynomial2D(degree=2)
fit_p = fitting.LevMarLSQFitter()

with warnings.catch_warnings():
    # Ignore model linearity warning from the fitter
    warnings.simplefilter('ignore')
    p = fit_p(p_init, x, y, z)

# Plot the data with the best-fit model
plt.figure(figsize=(8, 2.5))
plt.subplot(1, 3, 1)
plt.imshow(z, origin='lower', interpolation='nearest', vmin=-1e4, vmax=5e4)
plt.title("Data")
plt.subplot(1, 3, 2)
plt.imshow(p(x, y), origin='lower', interpolation='nearest', vmin=-1e4,
           vmax=5e4)
plt.title("Model")
plt.subplot(1, 3, 3)
plt.imshow(z - p(x, y), origin='lower', interpolation='nearest', vmin=-1e4,
           vmax=5e4)
plt.title("Residual")

(png, svg, pdf)

../_images/index-21.png

A list of models is provided in the Reference/API section. The fitting framework includes many useful features that are not demonstrated here, such as weighting of datapoints, fixing or linking parameters, and placing lower or upper limits on parameters. For more information on these, take a look at the Fitting Models to Data documentation.

Model sets

In some cases it is necessary to describe many models of the same type but with different sets of parameter values. This could be done simply by instantiating as many instances of a Model as are needed. But that can be inefficient for a large number of models. To that end, all model classes in astropy.modeling can also be used to represent a model set which is a collection of models of the same type, but with different values for their parameters.

To instantiate a model set, use argument n_models=N where N is the number of models in the set when constructing the model. The value of each parameter must be a list or array of length N, such that each item in the array corresponds to one model in the set:

>>> g = models.Gaussian1D(amplitude=[1, 2], mean=[0, 0],
...                       stddev=[0.1, 0.2], n_models=2)
>>> print(g)
Model: Gaussian1D
Inputs: ('x',)
Outputs: ('y',)
Model set size: 2
Parameters:
    amplitude mean stddev
    --------- ---- ------
          1.0  0.0    0.1
          2.0  0.0    0.2

This is equivalent to two Gaussians with the parameters amplitude=1, mean=0, stddev=0.1 and amplitude=2, mean=0, stddev=0.2 respectively. When printing the model the parameter values are displayed as a table, with each row corresponding to a single model in the set.

The number of models in a model set can be determined using the len builtin:

>>> len(g)
2

Single models have a length of 1, and are not considered a model set as such.

When evaluating a model set, by default the input must be the same length as the number of models, with one input per model:

>>> g([0, 0.1])  
array([1.        , 1.76499381])

The result is an array with one result per model in the set. It is also possible to broadcast a single value to all models in the set:

>>> g(0)  
array([1., 2.])

Model sets are used primarily for fitting, allowing a large number of models of the same type to be fitted simultaneously (and independently from each other) to some large set of inputs. For example, fitting a polynomial to the time response of each pixel in a data cube. This can greatly speed up the fitting process, especially for linear models.

Compound models

New in version 1.0: This feature is experimental and expected to see significant further development, but the basic usage is stable and expected to see wide use.

While the Astropy modeling package makes it very easy to define new models either from existing functions, or by writing a Model subclass, an additional way to create new models is by combining them using arithmetic expressions. This works with models built into Astropy, and most user-defined models as well. For example, it is possible to create a superposition of two Gaussians like so:

>>> from astropy.modeling import models
>>> g1 = models.Gaussian1D(1, 0, 0.2)
>>> g2 = models.Gaussian1D(2.5, 0.5, 0.1)
>>> g1_plus_2 = g1 + g2

The resulting object g1_plus_2 is itself a new model. Evaluating, say, g1_plus_2(0.25) is the same as evaluating g1(0.25) + g2(0.25):

>>> g1_plus_2(0.25)  
0.5676756958301329
>>> g1_plus_2(0.25) == g1(0.25) + g2(0.25)
True

This model can be further combined with other models in new expressions. It is also possible to define entire new model classes using arithmetic expressions of other model classes. This allows general compound models to be created without specifying any parameter values up front. This more advanced usage is explained in more detail in the compound model documentation.

These new compound models can also be fitted to data, like most other models (though this currently requires one of the non-linear fitters):

import numpy as np
import matplotlib.pyplot as plt
from astropy.modeling import models, fitting

# Generate fake data
np.random.seed(42)
g1 = models.Gaussian1D(1, 0, 0.2)
g2 = models.Gaussian1D(2.5, 0.5, 0.1)
x = np.linspace(-1, 1, 200)
y = g1(x) + g2(x) + np.random.normal(0., 0.2, x.shape)

# Now to fit the data create a new superposition with initial
# guesses for the parameters:
gg_init = models.Gaussian1D(1, 0, 0.1) + models.Gaussian1D(2, 0.5, 0.1)
fitter = fitting.SLSQPLSQFitter()
gg_fit = fitter(gg_init, x, y)

# Plot the data with the best-fit model
plt.figure(figsize=(8,5))
plt.plot(x, y, 'ko')
plt.plot(x, gg_fit(x))
plt.xlabel('Position')
plt.ylabel('Flux')

(png, svg, pdf)

../_images/index-3.png

This works for 1-D models, 2-D models, and combinations thereof, though there are some complexities involved in correctly matching up the inputs and outputs of all models used to build a compound model. You can learn more details in the Compound Models documentation.

Astropy models also support convolution through the function convolve_models, which returns a compound model.

For instance, the convolution of two Gaussian functions is also a Gaussian function in which the resulting mean (variance) is the sum of the means (variances) of each Gaussian.

import numpy as np
import matplotlib.pyplot as plt
from astropy.modeling import models
from astropy.convolution import convolve_models

g1 = models.Gaussian1D(1, -1, 1)
g2 = models.Gaussian1D(1, 1, 1)
g3 = convolve_models(g1, g2)

x = np.linspace(-3, -3, 50)
plt.plot(x, g1(x), 'k-')
plt.plot(x, g2(x), 'k-')
plt.plot(x, g3(x), 'k-')

(png, svg, pdf)

../_images/index-4.png

Reference/API

astropy.modeling Package

This subpackage provides a framework for representing models and performing model evaluation and fitting. It supports 1D and 2D models and fitting with parameter constraints. It has some predefined models and fitting routines.

Functions

custom_model(*args, **kwargs) Create a model from a user defined function.

Classes

Fittable1DModel Base class for one-dimensional fittable models.
Fittable2DModel Base class for two-dimensional fittable models.
FittableModel Base class for models that can be fitted using the built-in fitting algorithms.
InputParameterError Used for incorrect input parameter values and definitions.
Model Base class for all models.
ModelDefinitionError Used for incorrect models definitions
Parameter([name, description, default, …]) Wraps individual parameters.
ParameterError Generic exception class for all exceptions pertaining to Parameters.

Class Inheritance Diagram

Inheritance diagram of astropy.modeling.core.Fittable1DModel, astropy.modeling.core.Fittable2DModel, astropy.modeling.core.FittableModel, astropy.modeling.parameters.InputParameterError, astropy.modeling.core.Model, astropy.modeling.core.ModelDefinitionError, astropy.modeling.parameters.Parameter, astropy.modeling.parameters.ParameterError

astropy.modeling.functional_models Module

Mathematical models.

Classes

AiryDisk2D Two dimensional Airy disk model.
Moffat1D One dimensional Moffat model.
Moffat2D Two dimensional Moffat model.
Box1D One dimensional Box model.
Box2D Two dimensional Box model.
Const1D One dimensional Constant model.
Const2D Two dimensional Constant model.
Ellipse2D A 2D Ellipse model.
Disk2D Two dimensional radial symmetric Disk model.
BaseGaussian1D Base class for 1D Gaussian models.
Gaussian1D One dimensional Gaussian model.
GaussianAbsorption1D .
Gaussian2D Two dimensional Gaussian model.
Linear1D One dimensional Line model.
Lorentz1D One dimensional Lorentzian model.
MexicanHat1D One dimensional Mexican Hat model.
MexicanHat2D Two dimensional symmetric Mexican Hat model.
RedshiftScaleFactor One dimensional redshift scale factor model.
Scale Multiply a model by a factor.
Sersic1D One dimensional Sersic surface brightness profile.
Sersic2D Two dimensional Sersic surface brightness profile.
Shift Shift a coordinate.
Sine1D One dimensional Sine model.
Trapezoid1D One dimensional Trapezoid model.
TrapezoidDisk2D Two dimensional circular Trapezoid model.
Ring2D Two dimensional radial symmetric Ring model.
Voigt1D One dimensional model for the Voigt profile.

Class Inheritance Diagram

Inheritance diagram of astropy.modeling.functional_models.AiryDisk2D, astropy.modeling.functional_models.Moffat1D, astropy.modeling.functional_models.Moffat2D, astropy.modeling.functional_models.Box1D, astropy.modeling.functional_models.Box2D, astropy.modeling.functional_models.Const1D, astropy.modeling.functional_models.Const2D, astropy.modeling.functional_models.Ellipse2D, astropy.modeling.functional_models.Disk2D, astropy.modeling.functional_models.BaseGaussian1D, astropy.modeling.functional_models.Gaussian1D, astropy.modeling.functional_models.GaussianAbsorption1D, astropy.modeling.functional_models.Gaussian2D, astropy.modeling.functional_models.Linear1D, astropy.modeling.functional_models.Lorentz1D, astropy.modeling.functional_models.MexicanHat1D, astropy.modeling.functional_models.MexicanHat2D, astropy.modeling.functional_models.RedshiftScaleFactor, astropy.modeling.functional_models.Scale, astropy.modeling.functional_models.Sersic1D, astropy.modeling.functional_models.Sersic2D, astropy.modeling.functional_models.Shift, astropy.modeling.functional_models.Sine1D, astropy.modeling.functional_models.Trapezoid1D, astropy.modeling.functional_models.TrapezoidDisk2D, astropy.modeling.functional_models.Ring2D, astropy.modeling.functional_models.Voigt1D

astropy.modeling.powerlaws Module

Power law model variants

Classes

PowerLaw1D One dimensional power law model.
BrokenPowerLaw1D One dimensional power law model with a break.
SmoothlyBrokenPowerLaw1D One dimensional smoothly broken power law model.
ExponentialCutoffPowerLaw1D One dimensional power law model with an exponential cutoff.
LogParabola1D One dimensional log parabola model (sometimes called curved power law).

Class Inheritance Diagram

Inheritance diagram of astropy.modeling.powerlaws.PowerLaw1D, astropy.modeling.powerlaws.BrokenPowerLaw1D, astropy.modeling.powerlaws.SmoothlyBrokenPowerLaw1D, astropy.modeling.powerlaws.ExponentialCutoffPowerLaw1D, astropy.modeling.powerlaws.LogParabola1D

astropy.modeling.blackbody Module

Model and functions related to blackbody radiation.

Blackbody Radiation

Blackbody flux is calculated with Planck law (Rybicki & Lightman 1979):

\[ \begin{align}\begin{aligned}B_{\lambda}(T) = \frac{2 h c^{2} / \lambda^{5}}{exp(h c / \lambda k T) - 1}\\B_{\nu}(T) = \frac{2 h \nu^{3} / c^{2}}{exp(h \nu / k T) - 1}\end{aligned}\end{align} \]

where the unit of \(B_{\lambda}(T)\) is \(erg \; s^{-1} cm^{-2} \mathring{A}^{-1} sr^{-1}\), and \(B_{\nu}(T)\) is \(erg \; s^{-1} cm^{-2} Hz^{-1} sr^{-1}\). blackbody_lambda() and blackbody_nu() calculate the blackbody flux for \(B_{\lambda}(T)\) and \(B_{\nu}(T)\), respectively.

For blackbody representation as a model, see BlackBody1D.

Examples
>>> import numpy as np
>>> from astropy import units as u
>>> from astropy.modeling.blackbody import blackbody_lambda, blackbody_nu

Calculate blackbody flux for 5000 K at 100 and 10000 Angstrom while suppressing any Numpy warnings:

>>> wavelengths = [100, 10000] * u.AA
>>> temperature = 5000 * u.K
>>> with np.errstate(all='ignore'):
...     flux_lam = blackbody_lambda(wavelengths, temperature)
...     flux_nu = blackbody_nu(wavelengths, temperature)
>>> flux_lam  
<Quantity [  1.27452545e-108,  7.10190526e+005] erg / (Angstrom cm2 s sr)>
>>> flux_nu  
<Quantity [  4.25135927e-123,  2.36894060e-005] erg / (cm2 Hz s sr)>

Plot a blackbody spectrum for 5000 K:

(png, svg, pdf)

../_images/index-5.png

Note that an array of temperatures can also be given instead of a single temperature. In this case, the Numpy broadcasting rules apply: for instance, if the frequency and temperature have the same shape, the output will have this shape too, while if the frequency is a 2-d array with shape (n, m) and the temperature is an array with shape (m,), the output will have a shape (n, m).

See Also

Rybicki, G. B., & Lightman, A. P. 1979, Radiative Processes in Astrophysics (New York, NY: Wiley)

Functions

blackbody_nu(in_x, temperature) Calculate blackbody flux per steradian, \(B_{\nu}(T)\).
blackbody_lambda(in_x, temperature) Like blackbody_nu() but for \(B_{\lambda}(T)\).

Classes

BlackBody1D One dimensional blackbody model.

Class Inheritance Diagram

Inheritance diagram of astropy.modeling.blackbody.BlackBody1D

astropy.modeling.polynomial Module

This module contains models representing polynomials and polynomial series.

Classes

Chebyshev1D Univariate Chebyshev series.
Chebyshev2D Bivariate Chebyshev series.
Hermite1D Univariate Hermite series.
Hermite2D Bivariate Hermite series.
InverseSIP Inverse Simple Imaging Polynomial
Legendre1D Univariate Legendre series.
Legendre2D Bivariate Legendre series.
Polynomial1D 1D Polynomial model.
Polynomial2D 2D Polynomial model.
SIP Simple Imaging Polynomial (SIP) model.
OrthoPolynomialBase This is a base class for the 2D Chebyshev and Legendre models.
PolynomialModel Base class for polynomial models.

Class Inheritance Diagram

Inheritance diagram of astropy.modeling.polynomial.Chebyshev1D, astropy.modeling.polynomial.Chebyshev2D, astropy.modeling.polynomial.Hermite1D, astropy.modeling.polynomial.Hermite2D, astropy.modeling.polynomial.InverseSIP, astropy.modeling.polynomial.Legendre1D, astropy.modeling.polynomial.Legendre2D, astropy.modeling.polynomial.Polynomial1D, astropy.modeling.polynomial.Polynomial2D, astropy.modeling.polynomial.SIP, astropy.modeling.polynomial.OrthoPolynomialBase, astropy.modeling.polynomial.PolynomialModel

astropy.modeling.projections Module

Implements projections–particularly sky projections defined in WCS Paper II [R115].

All angles are set and and displayed in degrees but internally computations are performed in radians. All functions expect inputs and outputs degrees.

References

[R115]Calabretta, M.R., Greisen, E.W., 2002, A&A, 395, 1077 (Paper II)

Classes

Projection Base class for all sky projections.
Pix2SkyProjection Base class for all Pix2Sky projections.
Sky2PixProjection Base class for all Sky2Pix projections.
Zenithal Base class for all Zenithal projections.
Cylindrical Base class for Cylindrical projections.
PseudoCylindrical Base class for pseudocylindrical projections.
Conic Base class for conic projections.
PseudoConic Base class for pseudoconic projections.
QuadCube Base class for quad cube projections.
HEALPix Base class for HEALPix projections.
AffineTransformation2D Perform an affine transformation in 2 dimensions.
Pix2Sky_ZenithalPerspective Zenithal perspective projection - pixel to sky.
Sky2Pix_ZenithalPerspective Zenithal perspective projection - sky to pixel.
Pix2Sky_SlantZenithalPerspective Slant zenithal perspective projection - pixel to sky.
Sky2Pix_SlantZenithalPerspective Zenithal perspective projection - sky to pixel.
Pix2Sky_Gnomonic Gnomonic projection - pixel to sky.
Sky2Pix_Gnomonic Gnomonic Projection - sky to pixel.
Pix2Sky_Stereographic Stereographic Projection - pixel to sky.
Sky2Pix_Stereographic Stereographic Projection - sky to pixel.
Pix2Sky_SlantOrthographic Slant orthographic projection - pixel to sky.
Sky2Pix_SlantOrthographic Slant orthographic projection - sky to pixel.
Pix2Sky_ZenithalEquidistant Zenithal equidistant projection - pixel to sky.
Sky2Pix_ZenithalEquidistant Zenithal equidistant projection - sky to pixel.
Pix2Sky_ZenithalEqualArea Zenithal equidistant projection - pixel to sky.
Sky2Pix_ZenithalEqualArea Zenithal equidistant projection - sky to pixel.
Pix2Sky_Airy Airy projection - pixel to sky.
Sky2Pix_Airy Airy - sky to pixel.
Pix2Sky_CylindricalPerspective Cylindrical perspective - pixel to sky.
Sky2Pix_CylindricalPerspective Cylindrical Perspective - sky to pixel.
Pix2Sky_CylindricalEqualArea Cylindrical equal area projection - pixel to sky.
Sky2Pix_CylindricalEqualArea Cylindrical equal area projection - sky to pixel.
Pix2Sky_PlateCarree Plate carrée projection - pixel to sky.
Sky2Pix_PlateCarree Plate carrée projection - sky to pixel.
Pix2Sky_Mercator Mercator - pixel to sky.
Sky2Pix_Mercator Mercator - sky to pixel.
Pix2Sky_SansonFlamsteed Sanson-Flamsteed projection - pixel to sky.
Sky2Pix_SansonFlamsteed Sanson-Flamsteed projection - sky to pixel.
Pix2Sky_Parabolic Parabolic projection - pixel to sky.
Sky2Pix_Parabolic Parabolic projection - sky to pixel.
Pix2Sky_Molleweide Molleweide’s projection - pixel to sky.
Sky2Pix_Molleweide Molleweide’s projection - sky to pixel.
Pix2Sky_HammerAitoff Hammer-Aitoff projection - pixel to sky.
Sky2Pix_HammerAitoff Hammer-Aitoff projection - sky to pixel.
Pix2Sky_ConicPerspective Colles’ conic perspective projection - pixel to sky.
Sky2Pix_ConicPerspective Colles’ conic perspective projection - sky to pixel.
Pix2Sky_ConicEqualArea Alber’s conic equal area projection - pixel to sky.
Sky2Pix_ConicEqualArea Alber’s conic equal area projection - sky to pixel.
Pix2Sky_ConicEquidistant Conic equidistant projection - pixel to sky.
Sky2Pix_ConicEquidistant Conic equidistant projection - sky to pixel.
Pix2Sky_ConicOrthomorphic Conic orthomorphic projection - pixel to sky.
Sky2Pix_ConicOrthomorphic Conic orthomorphic projection - sky to pixel.
Pix2Sky_BonneEqualArea Bonne’s equal area pseudoconic projection - pixel to sky.
Sky2Pix_BonneEqualArea Bonne’s equal area pseudoconic projection - sky to pixel.
Pix2Sky_Polyconic Polyconic projection - pixel to sky.
Sky2Pix_Polyconic Polyconic projection - sky to pixel.
Pix2Sky_TangentialSphericalCube Tangential spherical cube projection - pixel to sky.
Sky2Pix_TangentialSphericalCube Tangential spherical cube projection - sky to pixel.
Pix2Sky_COBEQuadSphericalCube COBE quadrilateralized spherical cube projection - pixel to sky.
Sky2Pix_COBEQuadSphericalCube COBE quadrilateralized spherical cube projection - sky to pixel.
Pix2Sky_QuadSphericalCube Quadrilateralized spherical cube projection - pixel to sky.
Sky2Pix_QuadSphericalCube Quadrilateralized spherical cube projection - sky to pixel.
Pix2Sky_HEALPix HEALPix - pixel to sky.
Sky2Pix_HEALPix HEALPix projection - sky to pixel.
Pix2Sky_HEALPixPolar HEALPix polar, aka “butterfly” projection - pixel to sky.
Sky2Pix_HEALPixPolar HEALPix polar, aka “butterfly” projection - pixel to sky.
Pix2Sky_AZP alias of astropy.modeling.projections.Pix2Sky_ZenithalPerspective.
Sky2Pix_AZP alias of astropy.modeling.projections.Sky2Pix_ZenithalPerspective.
Pix2Sky_SZP alias of astropy.modeling.projections.Pix2Sky_SlantZenithalPerspective.
Sky2Pix_SZP alias of astropy.modeling.projections.Sky2Pix_SlantZenithalPerspective.
Pix2Sky_TAN alias of astropy.modeling.projections.Pix2Sky_Gnomonic.
Sky2Pix_TAN alias of astropy.modeling.projections.Sky2Pix_Gnomonic.
Pix2Sky_STG alias of astropy.modeling.projections.Pix2Sky_Stereographic.
Sky2Pix_STG alias of astropy.modeling.projections.Sky2Pix_Stereographic.
Pix2Sky_SIN alias of astropy.modeling.projections.Pix2Sky_SlantOrthographic.
Sky2Pix_SIN alias of astropy.modeling.projections.Sky2Pix_SlantOrthographic.
Pix2Sky_ARC alias of astropy.modeling.projections.Pix2Sky_ZenithalEquidistant.
Sky2Pix_ARC alias of astropy.modeling.projections.Sky2Pix_ZenithalEquidistant.
Pix2Sky_ZEA alias of astropy.modeling.projections.Pix2Sky_ZenithalEqualArea.
Sky2Pix_ZEA alias of astropy.modeling.projections.Sky2Pix_ZenithalEqualArea.
Pix2Sky_AIR alias of astropy.modeling.projections.Pix2Sky_Airy.
Sky2Pix_AIR alias of astropy.modeling.projections.Sky2Pix_Airy.
Pix2Sky_CYP alias of astropy.modeling.projections.Pix2Sky_CylindricalPerspective.
Sky2Pix_CYP alias of astropy.modeling.projections.Sky2Pix_CylindricalPerspective.
Pix2Sky_CEA alias of astropy.modeling.projections.Pix2Sky_CylindricalEqualArea.
Sky2Pix_CEA alias of astropy.modeling.projections.Sky2Pix_CylindricalEqualArea.
Pix2Sky_CAR alias of astropy.modeling.projections.Pix2Sky_PlateCarree.
Sky2Pix_CAR alias of astropy.modeling.projections.Sky2Pix_PlateCarree.
Pix2Sky_MER alias of astropy.modeling.projections.Pix2Sky_Mercator.
Sky2Pix_MER alias of astropy.modeling.projections.Sky2Pix_Mercator.
Pix2Sky_SFL alias of astropy.modeling.projections.Pix2Sky_SansonFlamsteed.
Sky2Pix_SFL alias of astropy.modeling.projections.Sky2Pix_SansonFlamsteed.
Pix2Sky_PAR alias of astropy.modeling.projections.Pix2Sky_Parabolic.
Sky2Pix_PAR alias of astropy.modeling.projections.Sky2Pix_Parabolic.
Pix2Sky_MOL alias of astropy.modeling.projections.Pix2Sky_Molleweide.
Sky2Pix_MOL alias of astropy.modeling.projections.Sky2Pix_Molleweide.
Pix2Sky_AIT alias of astropy.modeling.projections.Pix2Sky_HammerAitoff.
Sky2Pix_AIT alias of astropy.modeling.projections.Sky2Pix_HammerAitoff.
Pix2Sky_COP alias of astropy.modeling.projections.Pix2Sky_ConicPerspective.
Sky2Pix_COP alias of astropy.modeling.projections.Sky2Pix_ConicPerspective.
Pix2Sky_COE alias of astropy.modeling.projections.Pix2Sky_ConicEqualArea.
Sky2Pix_COE alias of astropy.modeling.projections.Sky2Pix_ConicEqualArea.
Pix2Sky_COD alias of astropy.modeling.projections.Pix2Sky_ConicEquidistant.
Sky2Pix_COD alias of astropy.modeling.projections.Sky2Pix_ConicEquidistant.
Pix2Sky_COO alias of astropy.modeling.projections.Pix2Sky_ConicOrthomorphic.
Sky2Pix_COO alias of astropy.modeling.projections.Sky2Pix_ConicOrthomorphic.
Pix2Sky_BON alias of astropy.modeling.projections.Pix2Sky_BonneEqualArea.
Sky2Pix_BON alias of astropy.modeling.projections.Sky2Pix_BonneEqualArea.
Pix2Sky_PCO alias of astropy.modeling.projections.Pix2Sky_Polyconic.
Sky2Pix_PCO alias of astropy.modeling.projections.Sky2Pix_Polyconic.
Pix2Sky_TSC alias of astropy.modeling.projections.Pix2Sky_TangentialSphericalCube.
Sky2Pix_TSC alias of astropy.modeling.projections.Sky2Pix_TangentialSphericalCube.
Pix2Sky_CSC alias of astropy.modeling.projections.Pix2Sky_COBEQuadSphericalCube.
Sky2Pix_CSC alias of astropy.modeling.projections.Sky2Pix_COBEQuadSphericalCube.
Pix2Sky_QSC alias of astropy.modeling.projections.Pix2Sky_QuadSphericalCube.
Sky2Pix_QSC alias of astropy.modeling.projections.Sky2Pix_QuadSphericalCube.
Pix2Sky_HPX alias of astropy.modeling.projections.Pix2Sky_HEALPix.
Sky2Pix_HPX alias of astropy.modeling.projections.Sky2Pix_HEALPix.
Pix2Sky_XPH alias of astropy.modeling.projections.Pix2Sky_HEALPixPolar.
Sky2Pix_XPH alias of astropy.modeling.projections.Sky2Pix_HEALPixPolar.

Class Inheritance Diagram

Inheritance diagram of astropy.modeling.projections.Projection, astropy.modeling.projections.Pix2SkyProjection, astropy.modeling.projections.Sky2PixProjection, astropy.modeling.projections.Zenithal, astropy.modeling.projections.Cylindrical, astropy.modeling.projections.PseudoCylindrical, astropy.modeling.projections.Conic, astropy.modeling.projections.PseudoConic, astropy.modeling.projections.QuadCube, astropy.modeling.projections.HEALPix, astropy.modeling.projections.AffineTransformation2D, astropy.modeling.projections.Pix2Sky_ZenithalPerspective, astropy.modeling.projections.Sky2Pix_ZenithalPerspective, astropy.modeling.projections.Pix2Sky_SlantZenithalPerspective, astropy.modeling.projections.Sky2Pix_SlantZenithalPerspective, astropy.modeling.projections.Pix2Sky_Gnomonic, astropy.modeling.projections.Sky2Pix_Gnomonic, astropy.modeling.projections.Pix2Sky_Stereographic, astropy.modeling.projections.Sky2Pix_Stereographic, astropy.modeling.projections.Pix2Sky_SlantOrthographic, astropy.modeling.projections.Sky2Pix_SlantOrthographic, astropy.modeling.projections.Pix2Sky_ZenithalEquidistant, astropy.modeling.projections.Sky2Pix_ZenithalEquidistant, astropy.modeling.projections.Pix2Sky_ZenithalEqualArea, astropy.modeling.projections.Sky2Pix_ZenithalEqualArea, astropy.modeling.projections.Pix2Sky_Airy, astropy.modeling.projections.Sky2Pix_Airy, astropy.modeling.projections.Pix2Sky_CylindricalPerspective, astropy.modeling.projections.Sky2Pix_CylindricalPerspective, astropy.modeling.projections.Pix2Sky_CylindricalEqualArea, astropy.modeling.projections.Sky2Pix_CylindricalEqualArea, astropy.modeling.projections.Pix2Sky_PlateCarree, astropy.modeling.projections.Sky2Pix_PlateCarree, astropy.modeling.projections.Pix2Sky_Mercator, astropy.modeling.projections.Sky2Pix_Mercator, astropy.modeling.projections.Pix2Sky_SansonFlamsteed, astropy.modeling.projections.Sky2Pix_SansonFlamsteed, astropy.modeling.projections.Pix2Sky_Parabolic, astropy.modeling.projections.Sky2Pix_Parabolic, astropy.modeling.projections.Pix2Sky_Molleweide, astropy.modeling.projections.Sky2Pix_Molleweide, astropy.modeling.projections.Pix2Sky_HammerAitoff, astropy.modeling.projections.Sky2Pix_HammerAitoff, astropy.modeling.projections.Pix2Sky_ConicPerspective, astropy.modeling.projections.Sky2Pix_ConicPerspective, astropy.modeling.projections.Pix2Sky_ConicEqualArea, astropy.modeling.projections.Sky2Pix_ConicEqualArea, astropy.modeling.projections.Pix2Sky_ConicEquidistant, astropy.modeling.projections.Sky2Pix_ConicEquidistant, astropy.modeling.projections.Pix2Sky_ConicOrthomorphic, astropy.modeling.projections.Sky2Pix_ConicOrthomorphic, astropy.modeling.projections.Pix2Sky_BonneEqualArea, astropy.modeling.projections.Sky2Pix_BonneEqualArea, astropy.modeling.projections.Pix2Sky_Polyconic, astropy.modeling.projections.Sky2Pix_Polyconic, astropy.modeling.projections.Pix2Sky_TangentialSphericalCube, astropy.modeling.projections.Sky2Pix_TangentialSphericalCube, astropy.modeling.projections.Pix2Sky_COBEQuadSphericalCube, astropy.modeling.projections.Sky2Pix_COBEQuadSphericalCube, astropy.modeling.projections.Pix2Sky_QuadSphericalCube, astropy.modeling.projections.Sky2Pix_QuadSphericalCube, astropy.modeling.projections.Pix2Sky_HEALPix, astropy.modeling.projections.Sky2Pix_HEALPix, astropy.modeling.projections.Pix2Sky_HEALPixPolar, astropy.modeling.projections.Sky2Pix_HEALPixPolar, astropy.modeling.projections.Pix2Sky_ZenithalPerspective, astropy.modeling.projections.Sky2Pix_ZenithalPerspective, astropy.modeling.projections.Pix2Sky_SlantZenithalPerspective, astropy.modeling.projections.Sky2Pix_SlantZenithalPerspective, astropy.modeling.projections.Pix2Sky_Gnomonic, astropy.modeling.projections.Sky2Pix_Gnomonic, astropy.modeling.projections.Pix2Sky_Stereographic, astropy.modeling.projections.Sky2Pix_Stereographic, astropy.modeling.projections.Pix2Sky_SlantOrthographic, astropy.modeling.projections.Sky2Pix_SlantOrthographic, astropy.modeling.projections.Pix2Sky_ZenithalEquidistant, astropy.modeling.projections.Sky2Pix_ZenithalEquidistant, astropy.modeling.projections.Pix2Sky_ZenithalEqualArea, astropy.modeling.projections.Sky2Pix_ZenithalEqualArea, astropy.modeling.projections.Pix2Sky_Airy, astropy.modeling.projections.Sky2Pix_Airy, astropy.modeling.projections.Pix2Sky_CylindricalPerspective, astropy.modeling.projections.Sky2Pix_CylindricalPerspective, astropy.modeling.projections.Pix2Sky_CylindricalEqualArea, astropy.modeling.projections.Sky2Pix_CylindricalEqualArea, astropy.modeling.projections.Pix2Sky_PlateCarree, astropy.modeling.projections.Sky2Pix_PlateCarree, astropy.modeling.projections.Pix2Sky_Mercator, astropy.modeling.projections.Sky2Pix_Mercator, astropy.modeling.projections.Pix2Sky_SansonFlamsteed, astropy.modeling.projections.Sky2Pix_SansonFlamsteed, astropy.modeling.projections.Pix2Sky_Parabolic, astropy.modeling.projections.Sky2Pix_Parabolic, astropy.modeling.projections.Pix2Sky_Molleweide, astropy.modeling.projections.Sky2Pix_Molleweide, astropy.modeling.projections.Pix2Sky_HammerAitoff, astropy.modeling.projections.Sky2Pix_HammerAitoff, astropy.modeling.projections.Pix2Sky_ConicPerspective, astropy.modeling.projections.Sky2Pix_ConicPerspective, astropy.modeling.projections.Pix2Sky_ConicEqualArea, astropy.modeling.projections.Sky2Pix_ConicEqualArea, astropy.modeling.projections.Pix2Sky_ConicEquidistant, astropy.modeling.projections.Sky2Pix_ConicEquidistant, astropy.modeling.projections.Pix2Sky_ConicOrthomorphic, astropy.modeling.projections.Sky2Pix_ConicOrthomorphic, astropy.modeling.projections.Pix2Sky_BonneEqualArea, astropy.modeling.projections.Sky2Pix_BonneEqualArea, astropy.modeling.projections.Pix2Sky_Polyconic, astropy.modeling.projections.Sky2Pix_Polyconic, astropy.modeling.projections.Pix2Sky_TangentialSphericalCube, astropy.modeling.projections.Sky2Pix_TangentialSphericalCube, astropy.modeling.projections.Pix2Sky_COBEQuadSphericalCube, astropy.modeling.projections.Sky2Pix_COBEQuadSphericalCube, astropy.modeling.projections.Pix2Sky_QuadSphericalCube, astropy.modeling.projections.Sky2Pix_QuadSphericalCube, astropy.modeling.projections.Pix2Sky_HEALPix, astropy.modeling.projections.Sky2Pix_HEALPix, astropy.modeling.projections.Pix2Sky_HEALPixPolar, astropy.modeling.projections.Sky2Pix_HEALPixPolar

astropy.modeling.rotations Module

Implements rotations, including spherical rotations as defined in WCS Paper II [R116]

RotateNative2Celestial and RotateCelestial2Native follow the convention in WCS Paper II to rotate to/from a native sphere and the celestial sphere.

The implementation uses EulerAngleRotation. The model parameters are three angles: the longitude (lon) and latitude (lat) of the fiducial point in the celestial system (CRVAL keywords in FITS), and the longitude of the celestial pole in the native system (lon_pole). The Euler angles are lon+90, 90-lat and -(lon_pole-90).

References

[R116]Calabretta, M.R., Greisen, E.W., 2002, A&A, 395, 1077 (Paper II)

Classes

RotateCelestial2Native Transform from Celestial to Native Spherical Coordinates.
RotateNative2Celestial Transform from Native to Celestial Spherical Coordinates.
Rotation2D Perform a 2D rotation given an angle.
EulerAngleRotation Implements Euler angle intrinsic rotations.

Class Inheritance Diagram

Inheritance diagram of astropy.modeling.rotations.RotateCelestial2Native, astropy.modeling.rotations.RotateNative2Celestial, astropy.modeling.rotations.Rotation2D, astropy.modeling.rotations.EulerAngleRotation

astropy.modeling.tabular Module

Tabular models.

Tabular models of any dimension can be created using tabular_model. For convenience Tabular1D and Tabular2D are provided.

Examples

>>> table = np.array([[ 3.,  0.,  0.],
...                  [ 0.,  2.,  0.],
...                  [ 0.,  0.,  0.]])
>>> points = ([1, 2, 3], [1, 2, 3])
>>> t2 = Tabular2D(points, lookup_table=table, bounds_error=False, fill_value=None, method='nearest')

Functions

tabular_model(dim[, name]) Make a Tabular model where n_inputs is based on the dimension of the lookup_table.
class astropy.modeling.tabular.Tabular1D

Tabular model in 1D. Returns an interpolated lookup table value.

Parameters:

points : array-like of float of ndim=1.

The points defining the regular grid in n dimensions.

lookup_table : array-like, of ndim=1.

The data in one dimensions.

method : str, optional

The method of interpolation to perform. Supported are “linear” and “nearest”, and “splinef2d”. “splinef2d” is only supported for 2-dimensional data. Default is “linear”.

bounds_error : bool, optional

If True, when interpolated values are requested outside of the domain of the input data, a ValueError is raised. If False, then fill_value is used.

fill_value : float, optional

If provided, the value to use for points outside of the interpolation domain. If None, values outside the domain are extrapolated. Extrapolation is not supported by method “splinef2d”.

Returns:

value : ndarray

Interpolated values at input coordinates.

Raises:

ImportError

Scipy is not installed.

Notes

Uses scipy.interpolate.interpn.

class astropy.modeling.tabular.Tabular2D

Tabular model in 2D. Returns an interpolated lookup table value.

Parameters:

points : tuple of ndarray of float, with shapes (m1, m2), optional

The points defining the regular grid in n dimensions.

lookup_table : array-like, shape (m1, m2)

The data on a regular grid in 2 dimensions.

method : str, optional

The method of interpolation to perform. Supported are “linear” and “nearest”, and “splinef2d”. “splinef2d” is only supported for 2-dimensional data. Default is “linear”.

bounds_error : bool, optional

If True, when interpolated values are requested outside of the domain of the input data, a ValueError is raised. If False, then fill_value is used.

fill_value : float, optional

If provided, the value to use for points outside of the interpolation domain. If None, values outside the domain are extrapolated. Extrapolation is not supported by method “splinef2d”.

Returns:

value : ndarray

Interpolated values at input coordinates.

Raises:

ImportError

Scipy is not installed.

Notes

Uses scipy.interpolate.interpn.

astropy.modeling.mappings Module

Special models useful for complex compound models where control is needed over which outputs from a source model are mapped to which inputs of a target model.

Classes

Mapping Allows inputs to be reordered, duplicated or dropped.
Identity Returns inputs unchanged.

Class Inheritance Diagram

Inheritance diagram of astropy.modeling.mappings.Mapping, astropy.modeling.mappings.Identity

astropy.modeling.fitting Module

This module implements classes (called Fitters) which combine optimization algorithms (typically from scipy.optimize) with statistic functions to perform fitting. Fitters are implemented as callable classes. In addition to the data to fit, the __call__ method takes an instance of FittableModel as input, and returns a copy of the model with its parameters determined by the optimizer.

Optimization algorithms, called “optimizers” are implemented in optimizers and statistic functions are in statistic. The goal is to provide an easy to extend framework and allow users to easily create new fitters by combining statistics with optimizers.

There are two exceptions to the above scheme. LinearLSQFitter uses Numpy’s lstsq function. LevMarLSQFitter uses leastsq which combines optimization and statistic in one implementation.

Classes

LinearLSQFitter() A class performing a linear least square fitting.
LevMarLSQFitter() Levenberg-Marquardt algorithm and least squares statistic.
FittingWithOutlierRemoval(fitter, outlier_func) This class combines an outlier removal technique with a fitting procedure.
SLSQPLSQFitter() SLSQP optimization algorithm and least squares statistic.
SimplexLSQFitter() Simplex algorithm and least squares statistic.
JointFitter(models, jointparameters, initvals) Fit models which share a parameter.
Fitter(optimizer, statistic) Base class for all fitters.

Class Inheritance Diagram

Inheritance diagram of astropy.modeling.fitting.LinearLSQFitter, astropy.modeling.fitting.LevMarLSQFitter, astropy.modeling.fitting.FittingWithOutlierRemoval, astropy.modeling.fitting.SLSQPLSQFitter, astropy.modeling.fitting.SimplexLSQFitter, astropy.modeling.fitting.JointFitter, astropy.modeling.fitting.Fitter

astropy.modeling.optimizers Module

Optimization algorithms used in fitting.

Classes

Optimization(opt_method) Base class for optimizers.
SLSQP() Sequential Least Squares Programming optimization algorithm.
Simplex() Neald-Mead (downhill simplex) algorithm.

Class Inheritance Diagram

Inheritance diagram of astropy.modeling.optimizers.Optimization, astropy.modeling.optimizers.SLSQP, astropy.modeling.optimizers.Simplex

astropy.modeling.statistic Module

Statistic functions used in fitting.

Functions

leastsquare(measured_vals, updated_model, …) Least square statistic with optional weights.