tessl/pypi-astropy
Astronomy and astrophysics core library providing comprehensive tools for astronomical computations and data handling
—
Pending
Quality
Pending
Does it follow best practices?
Impact
Pending
No eval scenarios have been run
modeling.mddocs/
Modeling and Fitting
Model fitting framework with astronomical models (Gaussian, Sérsic, etc.) and fitting algorithms for data analysis and parameter estimation.
Capabilities
Base Model Classes
Foundation classes for creating and manipulating models with parameters, constraints, and evaluation methods.
class Model:
"""
Base class for all models.
Parameters:
- *args: positional arguments for model parameters
- **kwargs: keyword arguments for model parameters and metadata
"""
def __init__(self, *args, **kwargs): ...
def __call__(self, *args, **kwargs):
"""
Evaluate model at given inputs.
Parameters:
- *args: input coordinate arrays
- **kwargs: additional evaluation options
Returns:
ndarray: model values
"""
def evaluate(self, *args, **kwargs):
"""Evaluate model (implementation method)."""
def copy(self):
"""Create a copy of the model."""
def rename(self, name):
"""Rename the model."""
@property
def param_names(self):
"""List of parameter names."""
@property
def parameters(self):
"""Array of parameter values."""
@parameters.setter
def parameters(self, values):
"""Set parameter values."""
def __add__(self, other):
"""Add models together."""
def __sub__(self, other):
"""Subtract models."""
def __mul__(self, other):
"""Multiply models."""
def __truediv__(self, other):
"""Divide models."""
def __pow__(self, other):
"""Raise model to power."""
class Fittable1DModel(Model):
"""Base class for fittable 1D models."""
def __init__(self, *args, **kwargs): ...
class Fittable2DModel(Model):
"""Base class for fittable 2D models."""
def __init__(self, *args, **kwargs): ...
class CompoundModel(Model):
"""Model composed of multiple sub-models."""
def __init__(self, op, left, right, name=None): ...Model Parameters
Parameter class for managing model parameters with bounds, constraints, and units.
class Parameter:
"""
Model parameter with value, bounds, and constraints.
Parameters:
- name: parameter name
- default: default value
- description: parameter description
- unit: parameter unit
- min: minimum allowed value
- max: maximum allowed value
- bounds: (min, max) bounds tuple
- tied: tie parameter to another parameter or function
- fixed: whether parameter is fixed during fitting
"""
def __init__(self, name='', default=0, description='', unit=None,
min=None, max=None, bounds=None, tied=False, fixed=False): ...
@property
def value(self):
"""Parameter value."""
@value.setter
def value(self, val):
"""Set parameter value."""
@property
def bounds(self):
"""Parameter bounds (min, max)."""
@property
def tied(self):
"""Parameter tied constraint."""
@property
def fixed(self):
"""Whether parameter is fixed."""
def copy(self):
"""Create copy of parameter."""1D Astronomical Models
One-dimensional models commonly used in astronomical data analysis.
class Gaussian1D(Fittable1DModel):
"""
1D Gaussian model.
Parameters:
- amplitude: Gaussian amplitude
- mean: Gaussian center
- stddev: Gaussian standard deviation
"""
def __init__(self, amplitude=1, mean=0, stddev=1, **kwargs): ...
class Lorentz1D(Fittable1DModel):
"""
1D Lorentzian model.
Parameters:
- amplitude: Lorentzian amplitude
- x_0: Lorentzian center
- fwhm: full width at half maximum
"""
def __init__(self, amplitude=1, x_0=0, fwhm=1, **kwargs): ...
class Voigt1D(Fittable1DModel):
"""
1D Voigt profile (convolution of Gaussian and Lorentzian).
Parameters:
- x_0: center position
- amplitude_L: Lorentzian amplitude
- fwhm_L: Lorentzian FWHM
- fwhm_G: Gaussian FWHM
"""
def __init__(self, x_0=0, amplitude_L=1, fwhm_L=1, fwhm_G=1, **kwargs): ...
class Moffat1D(Fittable1DModel):
"""
1D Moffat model for stellar profiles.
Parameters:
- amplitude: profile amplitude
- x_0: profile center
- gamma: profile width parameter
- alpha: profile power index
"""
def __init__(self, amplitude=1, x_0=0, gamma=1, alpha=1, **kwargs): ...
class Sersic1D(Fittable1DModel):
"""
1D Sérsic profile for galaxy surface brightness.
Parameters:
- amplitude: profile amplitude
- r_eff: effective radius
- n: Sérsic index
"""
def __init__(self, amplitude=1, r_eff=1, n=4, **kwargs): ...
class PowerLaw1D(Fittable1DModel):
"""
1D power law model.
Parameters:
- amplitude: amplitude at reference point
- x_0: reference point
- alpha: power law index
"""
def __init__(self, amplitude=1, x_0=1, alpha=1, **kwargs): ...
class BlackBody(Fittable1DModel):
"""
Blackbody spectrum model.
Parameters:
- temperature: blackbody temperature
- bolometric_flux: total bolometric flux
"""
def __init__(self, temperature=5000*u.K, bolometric_flux=1*u.erg/u.s/u.cm**2, **kwargs): ...
class Polynomial1D(Fittable1DModel):
"""
1D polynomial model.
Parameters:
- degree: polynomial degree
- c0, c1, c2, ...: polynomial coefficients
"""
def __init__(self, degree, **kwargs): ...2D Astronomical Models
Two-dimensional models for image analysis and surface brightness modeling.
class Gaussian2D(Fittable2DModel):
"""
2D Gaussian model.
Parameters:
- amplitude: Gaussian amplitude
- x_mean: x-coordinate of center
- y_mean: y-coordinate of center
- x_stddev: standard deviation in x
- y_stddev: standard deviation in y
- theta: rotation angle
"""
def __init__(self, amplitude=1, x_mean=0, y_mean=0, x_stddev=1, y_stddev=1, theta=0, **kwargs): ...
class Moffat2D(Fittable2DModel):
"""
2D Moffat model for stellar PSF.
Parameters:
- amplitude: profile amplitude
- x_0: x-coordinate of center
- y_0: y-coordinate of center
- gamma: profile width parameter
- alpha: power index
"""
def __init__(self, amplitude=1, x_0=0, y_0=0, gamma=1, alpha=1, **kwargs): ...
class Sersic2D(Fittable2DModel):
"""
2D Sérsic profile for galaxy modeling.
Parameters:
- amplitude: surface brightness at effective radius
- r_eff: effective radius
- n: Sérsic index
- x_0: x-coordinate of center
- y_0: y-coordinate of center
- ellip: ellipticity
- theta: position angle
"""
def __init__(self, amplitude=1, r_eff=1, n=4, x_0=0, y_0=0, ellip=0, theta=0, **kwargs): ...
class AiryDisk2D(Fittable2DModel):
"""
2D Airy disk for diffraction-limited PSF.
Parameters:
- amplitude: peak amplitude
- x_0: x-coordinate of center
- y_0: y-coordinate of center
- radius: Airy disk radius
"""
def __init__(self, amplitude=1, x_0=0, y_0=0, radius=1, **kwargs): ...
class Disk2D(Fittable2DModel):
"""
2D disk (top-hat) model.
Parameters:
- amplitude: disk amplitude
- x_0: x-coordinate of center
- y_0: y-coordinate of center
- R_0: disk radius
"""
def __init__(self, amplitude=1, x_0=0, y_0=0, R_0=1, **kwargs): ...
class Polynomial2D(Fittable2DModel):
"""
2D polynomial model.
Parameters:
- degree: polynomial degree
- c0_0, c1_0, c0_1, ...: polynomial coefficients
"""
def __init__(self, degree, **kwargs): ...Fitting Algorithms
Various fitting algorithms for parameter estimation and model optimization.
class Fitter:
"""Base class for fitting algorithms."""
def __call__(self, model, x, y, weights=None, **kwargs):
"""
Fit model to data.
Parameters:
- model: model to fit
- x: independent variable data
- y: dependent variable data
- weights: data weights
- **kwargs: fitter-specific options
Returns:
Model: fitted model instance
"""
class LinearLSQFitter(Fitter):
"""Linear least squares fitter for linear models."""
def __init__(self, calc_uncertainties=False): ...
class LevMarLSQFitter(Fitter):
"""Levenberg-Marquardt least squares fitter."""
def __init__(self, calc_uncertainties=False): ...
class SLSQPLSQFitter(Fitter):
"""Sequential Least Squares Programming fitter."""
def __init__(self, calc_uncertainties=False): ...
class SimplexLSQFitter(Fitter):
"""Simplex algorithm fitter."""
def __init__(self, calc_uncertainties=False): ...
class DifferentialEvolutionLSQFitter(Fitter):
"""Differential evolution global optimizer."""
def __init__(self, calc_uncertainties=False): ...
class LMLSQFitter(Fitter):
"""Levenberg-Marquardt fitter using scipy.optimize.leastsq."""
def __init__(self, calc_uncertainties=False): ...Model Utilities
Utility functions for model evaluation, parameter handling, and model combinations.
def custom_model(func, fit_deriv=None):
"""
Create custom model from function.
Parameters:
- func: function that defines the model
- fit_deriv: function that computes derivatives
Returns:
Model: custom model class
"""
def bind_bounding_box(model, bounding_box, order='C'):
"""Bind a bounding box to a model."""
def fix_inputs(model, fixed_inputs):
"""Fix certain inputs of a model."""Usage Examples
Basic Model Fitting
from astropy.modeling import models, fitting
import numpy as np
import matplotlib.pyplot as plt
# Generate synthetic data
np.random.seed(42)
x = np.linspace(-3, 3, 100)
y_true = 2 * np.exp(-0.5 * ((x - 0.5) / 0.8)**2) + 0.1
noise = np.random.normal(0, 0.1, len(x))
y = y_true + noise
# Create and fit Gaussian model
gaussian_model = models.Gaussian1D(amplitude=2, mean=0.5, stddev=0.8)
fitter = fitting.LevMarLSQFitter()
fitted_model = fitter(gaussian_model, x, y)
print(f"Fitted parameters:")
print(f" Amplitude: {fitted_model.amplitude.value:.3f}")
print(f" Mean: {fitted_model.mean.value:.3f}")
print(f" Stddev: {fitted_model.stddev.value:.3f}")
# Evaluate fitted model
y_fit = fitted_model(x)2D Surface Brightness Modeling
from astropy.modeling import models, fitting
import numpy as np
# Create 2D coordinate grids
y, x = np.mgrid[0:100, 0:100]
# Generate synthetic galaxy with Sérsic profile
true_model = models.Sersic2D(amplitude=1000, r_eff=15, n=2,
x_0=50, y_0=45, ellip=0.3, theta=0.5)
galaxy_data = true_model(x, y)
# Add noise
noise = np.random.poisson(galaxy_data + 100) - 100 # Poisson noise
galaxy_data_noisy = galaxy_data + noise
# Fit Sérsic model
initial_model = models.Sersic2D(amplitude=800, r_eff=10, n=1,
x_0=48, y_0=47, ellip=0.2, theta=0)
fitter = fitting.LevMarLSQFitter()
fitted_galaxy = fitter(initial_model, x, y, galaxy_data_noisy)
print("Galaxy fitting results:")
print(f" Effective radius: {fitted_galaxy.r_eff.value:.2f} pixels")
print(f" Sérsic index: {fitted_galaxy.n.value:.2f}")
print(f" Ellipticity: {fitted_galaxy.ellip.value:.3f}")
print(f" Position angle: {fitted_galaxy.theta.value:.3f} radians")Spectral Line Fitting
# Fit multiple spectral lines
wavelength = np.linspace(6540, 6580, 200)
# Create composite model with multiple Gaussian lines
h_alpha = models.Gaussian1D(amplitude=100, mean=6562.8, stddev=0.5, name='H_alpha')
nii_6548 = models.Gaussian1D(amplitude=30, mean=6548.1, stddev=0.4, name='NII_6548')
nii_6583 = models.Gaussian1D(amplitude=90, mean=6583.5, stddev=0.4, name='NII_6583')
continuum = models.Polynomial1D(degree=1, c0=10, c1=0, name='continuum')
# Combine models
line_model = continuum + h_alpha + nii_6548 + nii_6583
# Generate synthetic spectrum
spectrum = line_model(wavelength)
spectrum += np.random.normal(0, 2, len(wavelength)) # Add noise
# Fit the composite model
fitter = fitting.LevMarLSQFitter()
fitted_lines = fitter(line_model, wavelength, spectrum)
# Extract individual line fluxes
h_alpha_flux = fitted_lines['H_alpha'].amplitude * fitted_lines['H_alpha'].stddev * np.sqrt(2*np.pi)
nii_ratio = fitted_lines['NII_6583'].amplitude / fitted_lines['H_alpha'].amplitude
print(f"H-alpha flux: {h_alpha_flux:.1f}")
print(f"[NII]/H-alpha ratio: {nii_ratio:.3f}")Parameter Constraints and Bounds
# Model with parameter constraints
psf_model = models.Moffat2D(amplitude=1000, x_0=50, y_0=50, gamma=2, alpha=2)
# Set parameter bounds
psf_model.amplitude.bounds = (100, 10000)
psf_model.x_0.bounds = (48, 52)
psf_model.y_0.bounds = (48, 52)
psf_model.gamma.bounds = (0.5, 5)
psf_model.alpha.bounds = (1, 10)
# Fix certain parameters
psf_model.x_0.fixed = True
psf_model.y_0.fixed = True
# Tie parameters together (force circular PSF)
def tie_gamma_to_alpha(model):
return model.alpha * 0.5
psf_model.gamma.tied = tie_gamma_to_alpha
# Fit with constraints
fitter = fitting.LevMarLSQFitter()
# fitted_psf = fitter(psf_model, x, y, psf_data)Custom Model Creation
from astropy.modeling import Fittable1DModel, Parameter
import numpy as np
class ExponentialCutoff(Fittable1DModel):
\"\"\"Custom exponentially cutoff power law model.\"\"\"
amplitude = Parameter(default=1, bounds=(0, None))
x_0 = Parameter(default=1, bounds=(0, None))
alpha = Parameter(default=1)
cutoff = Parameter(default=1, bounds=(0, None))
@staticmethod
def evaluate(x, amplitude, x_0, alpha, cutoff):
\"\"\"Model evaluation function.\"\"\"
return amplitude * (x / x_0)**(-alpha) * np.exp(-x / cutoff)
# Use custom model
x_data = np.logspace(0, 3, 50)
custom_model = ExponentialCutoff(amplitude=100, x_0=10, alpha=2, cutoff=200)
y_model = custom_model(x_data)
print(f"Custom model at x=50: {custom_model(50):.2f}")Model Arithmetic and Combinations
# Create individual models
gaussian = models.Gaussian1D(amplitude=10, mean=0, stddev=1)
linear = models.Linear1D(slope=0.5, intercept=2)
constant = models.Const1D(amplitude=5)
# Combine models using arithmetic
combined_additive = gaussian + linear + constant
combined_multiplicative = gaussian * linear
combined_complex = (gaussian + constant) * linear
# Evaluate combined models
x_test = np.linspace(-5, 5, 100)
y_combined = combined_additive(x_test)
print(f"Number of parameters in combined model: {len(combined_additive.parameters)}")
print(f"Parameter names: {combined_additive.param_names}")Advanced Fitting with Uncertainties
# Fitting with parameter uncertainties
fitter_with_errs = fitting.LevMarLSQFitter(calc_uncertainties=True)
# Generate data with known uncertainties
x_data = np.linspace(0, 10, 50)
true_model = models.Polynomial1D(degree=2, c0=1, c1=2, c2=-0.1)
y_true = true_model(x_data)
y_errors = 0.1 + 0.05 * np.abs(y_true) # Heteroscedastic errors
y_data = y_true + np.random.normal(0, y_errors)
# Fit with error weights
weights = 1.0 / y_errors**2
poly_model = models.Polynomial1D(degree=2)
fitted_poly = fitter_with_errs(poly_model, x_data, y_data, weights=weights)
# Extract parameter uncertainties (if supported by fitter)
if hasattr(fitter_with_errs, 'fit_info') and fitter_with_errs.fit_info.get('param_cov') is not None:
param_errors = np.sqrt(np.diag(fitter_with_errs.fit_info['param_cov']))
print("Parameter uncertainties:")
for name, value, error in zip(fitted_poly.param_names, fitted_poly.parameters, param_errors):
print(f" {name}: {value:.4f} ± {error:.4f}")