Reference

Model

Model(x, *, weights=None, name='dirbe', gauss_quad_degree=50, extrapolate=False, ephemeris='builtin')

Main interface to ZodiPy.

Initialize a zodiacal light model.

PARAMETER	DESCRIPTION
x	Wavelength or frequency. If x is a sequence it is assumed to be a the points corresponding to an instrument bandpass and the corresponding weights argument must be provided. TYPE: Quantity[micron | GHz]
weights	Bandpass weights corresponding the the frequencies/wavelengths in x. The weights are assumed to represent a normalized instrument response in units of spectral radiance [Jy/sr]. TYPE: ArrayLike | None DEFAULT: None
name	Zodiacal light model to use. See the docs for list of available models. Defaults to 'dirbe'. TYPE: str DEFAULT: 'dirbe'
gauss_quad_degree	Order of the Gaussian-legendre quadrature representing the number of discrete points along each line-of-sight. Default is 50 points. TYPE: int DEFAULT: 50
extrapolate	If True all spectral quantities in the selected model are extrapolated to the requested frequencies/wavelengths. Else, an exception is raised on values of x outside of the valid model range. Default is False. TYPE: bool DEFAULT: False
ephemeris	Ephemeris used in Astropy's solar_system_ephemeris to compute the positions of Earth and optionally the observer. See the Astropy documentation for all available ephemerides. Defaults to 'builtin'. TYPE: str DEFAULT: 'builtin'

Source code in zodipy/model.py

def __init__(
    self,
    x: units.Quantity[units.micron | units.GHz],
    *,
    weights: npt.ArrayLike | None = None,
    name: str = "dirbe",
    gauss_quad_degree: int = 50,
    extrapolate: bool = False,
    ephemeris: str = "builtin",
) -> None:
    """Initialize a zodiacal light model.

    Args:
        x: Wavelength or frequency. If `x` is a sequence it is assumed to be a the points
            corresponding to an instrument bandpass and the corresponding `weights` argument
            must be provided.
        weights: Bandpass weights corresponding the the frequencies/wavelengths in `x`. The
            weights are assumed to represent a normalized instrument response in units of
            spectral radiance [Jy/sr].
        name: Zodiacal light model to use. See the
            [docs](https://cosmoglobe.github.io/zodipy/introduction/) for list of available
            models. Defaults to 'dirbe'.
        gauss_quad_degree: Order of the Gaussian-legendre quadrature representing the number of
            discrete points along each line-of-sight. Default is 50 points.
        extrapolate: If `True` all spectral quantities in the selected model are extrapolated to
            the requested frequencies/wavelengths. Else, an exception is raised on values of `x`
            outside of the valid model range. Default is `False`.
        ephemeris: Ephemeris used in Astropy's `solar_system_ephemeris` to compute the positions
            of Earth and optionally the observer. See the
            [Astropy documentation](https://docs.astropy.org/en/stable/coordinates/solarsystem.html)
            for all available ephemerides. Defaults to 'builtin'.

    """
    try:
        if not x.isscalar and weights is None:
            msg = "Bandpass weights must be provided for non-scalar `x`."
            raise ValueError(msg)
    except AttributeError as error:
        msg = "The input 'x' must be an astropy Quantity."
        raise TypeError(msg) from error
    if x.isscalar and weights is not None:
        msg = "Bandpass weights should not be provided for scalar `x`."
        raise ValueError(msg)

    self._ipd_model = model_registry.get_model(name)

    if not extrapolate and not self._ipd_model.is_valid_at(x):
        msg = (
            "The requested frequencies are outside the valid range of the model. "
            "If this was intended, set the extrapolate argument to True."
        )
        raise ValueError(msg)

    # Bandpass is provided rather than a delta wavelength or frequency.
    if weights is not None:
        weights = np.asarray(weights)
        if x.size != weights.size:
            msg = "Number of wavelengths and weights must be the same in the bandpass."
            raise ValueError(msg)
        normalized_weights = weights / integrate.trapezoid(weights, x)
    else:
        normalized_weights = None

    self._x = x
    self._normalized_weights = normalized_weights
    self._b_nu_table = tabulate_blackbody_emission(self._x, self._normalized_weights)

    quad_points, quad_weights = np.polynomial.legendre.leggauss(gauss_quad_degree)
    self._integrate_leggauss = functools.partial(
        integrate_leggauss,
        points=quad_points,
        weights=quad_weights,
    )

    self._ephemeris = ephemeris

    # Make mypy happy by declaring types of to-be initialized attributes.
    self._number_density_partials: dict[ComponentLabel, functools.partial]
    self._interped_comp_params: dict[ComponentLabel, dict]
    self._interped_shared_params: dict

    self._init_ipd_model_partials()

evaluate

evaluate(skycoord, *, obspos='earth', return_comps=False, nprocesses=1)

Return the simulated zodiacal light.

The zodiacal light is simulated for all sky coordinates present in the skycoord argument. If an obstime and obspos value is not provided for each coordinate value, all coordinates are assumed to be observed at an instant from the same position.

PARAMETER	DESCRIPTION
skycoord	astropy.coordinates.SkyCoord object representing the coordinates for which to simulate the zodiacal light. The obstime attribute must be specified, and correspond to a either a single, or a sequence of observational times, one for each coordinate in skycoord. The coordinate frame, provided through the frame keyword in the the astropy.coordinates.SkyCoord object (defaults to astropy.coordinates.ICRS), must be convertible to the astropy.coordinates.BarycentricMeanEcliptic frame. TYPE: SkyCoord
obspos	The heliocentric ecliptic position of the observer, or a string representing an observer supported by the astropy.coordinates.solar_system_ephemeris. If an explicit position is given, it must either be a single, or a sequence of positions, one for each coordinate. Defaults to 'earth'. TYPE: Quantity[AU] | str DEFAULT: 'earth'
return_comps	If True, the emission is returned component-wise. Defaults to False. TYPE: bool DEFAULT: False
nprocesses	Number of cores to use. If nprocesses >= 1, the line-of-sight integrals are distributed and computed in parallel using the multiprocessing module. Defaults to 1. TYPE: int DEFAULT: 1

RETURNS	DESCRIPTION
emission	Simulated zodiacal light [MJy/sr]. TYPE: Quantity[MJy / sr]

Source code in zodipy/model.py

def evaluate(
    self,
    skycoord: coords.SkyCoord,
    *,
    obspos: units.Quantity[units.AU] | str = "earth",
    return_comps: bool = False,
    nprocesses: int = 1,
) -> units.Quantity[units.MJy / units.sr]:
    """Return the simulated zodiacal light.

    The zodiacal light is simulated for all sky coordinates present in the `skycoord` argument.
    If an obstime and obspos value is not provided for each coordinate value, all coordinates
    are assumed to be observed at an instant from the same position.

    Args:
        skycoord: `astropy.coordinates.SkyCoord` object representing the coordinates for which
            to simulate the zodiacal light. The `obstime` attribute must be specified, and
            correspond to a either a single, or a sequence of observational times, one for each
            coordinate in `skycoord`. The coordinate frame, provided through the `frame` keyword
            in the the `astropy.coordinates.SkyCoord` object (defaults to
            `astropy.coordinates.ICRS`), must be convertible to the
            `astropy.coordinates.BarycentricMeanEcliptic` frame.
        obspos: The heliocentric ecliptic position of the observer, or a string representing an
            observer supported by the `astropy.coordinates.solar_system_ephemeris`. If an
            explicit position is given, it must either be a single, or a sequence of positions,
            one for each coordinate. Defaults to 'earth'.
        return_comps: If `True`, the emission is returned component-wise. Defaults to `False`.
        nprocesses: Number of cores to use. If `nprocesses >= 1`, the line-of-sight integrals
            are distributed and computed in parallel using the `multiprocessing` module.
            Defaults to 1.

    Returns:
        emission: Simulated zodiacal light [MJy/sr].

    """
    try:
        if skycoord.obstime is None:
            msg = "The `obstime` attribute of the `SkyCoord` object is not set."
            raise ValueError(msg)
    except AttributeError as error:
        msg = "The input coordinates must be an astropy SkyCoord object."
        raise TypeError(msg) from error

    try:
        if not (obspos_isstr := isinstance(obspos, str)) and (
            (obspos.ndim > 1 and skycoord.obstime.size != obspos.shape[-1])
            or (obspos.ndim == 1 and skycoord.obstime.size != 1)
        ):
            msg = "The number of obstime (ncoords) and obspos (3, ncoords) does not match."
            raise ValueError(msg)
    except AttributeError as error:
        msg = "The observer position is not a string or an astropy Quantity."
        raise TypeError(msg) from error

    if skycoord.obstime.size > skycoord.size:
        msg = "The size of obstime must be either 1 or ncoords."
        raise ValueError(msg)

    if skycoord.obstime.size == 1:
        interp_obstimes = None
    else:
        interp_obstimes = arrange_obstimes(skycoord.obstime[0].mjd, skycoord.obstime[-1].mjd)

    dist_coords_to_cores = skycoord.size > nprocesses and nprocesses > 1
    if dist_coords_to_cores:
        skycoord_splits = np.array_split(skycoord, nprocesses)
        obspos_splits = (
            itertools.repeat(obspos, nprocesses)
            if obspos_isstr
            else np.array_split(obspos, nprocesses, axis=-1)
        )
        interp_obstime_splits = itertools.repeat(interp_obstimes, nprocesses)
        with multiprocessing.get_context(PLATFORM_METHOD).Pool(nprocesses) as pool:
            emission_splits = [
                pool.apply_async(self._evaluate, args=(skycoord, obspos, obstime_lims))
                for skycoord, obspos, obstime_lims in zip(
                    skycoord_splits, obspos_splits, interp_obstime_splits
                )
            ]
            emission = np.concatenate([split.get() for split in emission_splits], axis=-1)
    else:
        emission = self._evaluate(skycoord, obspos, interp_obstimes)

    emission <<= units.MJy / units.sr
    return emission if return_comps else emission.sum(axis=0)

get_parameters

get_parameters()

Return a dictionary containing the zodiacal light model parameters.

RETURNS	DESCRIPTION
parameters	Zodiacal light model parameter dict. TYPE: dict

Source code in zodipy/model.py

def get_parameters(self) -> dict:
    """Return a dictionary containing the zodiacal light model parameters.

    Returns:
        parameters: Zodiacal light model parameter dict.
    """
    return self._ipd_model.to_dict()

update_parameters

update_parameters(parameters)

Update the zodiacal light model parameters from a parameter dictionary.

The structure of the input dictionary must match that of the output of the get_parameters method.

PARAMETER	DESCRIPTION
parameters	Zodiacal light model parameter dict. TYPE: dict

Source code in zodipy/model.py

def update_parameters(self, parameters: dict) -> None:
    """Update the zodiacal light model parameters from a parameter dictionary.

    The structure of the input dictionary must match that of the output of the `get_parameters`
    method.

    Args:
        parameters: Zodiacal light model parameter dict.

    """
    _dict = parameters.copy()
    _dict["comps"] = {}
    for key, value in parameters.items():
        if key == "comps":
            for comp_key, comp_value in value.items():
                _dict["comps"][ComponentLabel(comp_key)] = type(
                    self._ipd_model.comps[ComponentLabel(comp_key)]
                )(**comp_value)
        elif isinstance(value, dict):
            _dict[key] = {ComponentLabel(k): v for k, v in value.items()}

    self._ipd_model = self._ipd_model.__class__(**_dict)
    self._init_ipd_model_partials()

grid_number_density

grid_number_density(x, y, z, obstime, model='dirbe', ephemeris='builtin')

Return the component-wise tabulated densities of the zodiacal components for a given grid.

PARAMETER	DESCRIPTION
x	x-coordinates of a cartesian mesh grid. TYPE: Quantity[AU]
y	y-coordinates of a cartesian mesh grid. TYPE: Quantity[AU]
z	z-coordinates of a cartesian mesh grid. TYPE: Quantity[AU]
obstime	Time of observation. Required to compute the Earth position. TYPE: Time
model	String representing a built-in model or an explicit instance of a ZodiacalLightModel. Default is 'dirbe'. TYPE: str | ZodiacalLightModel DEFAULT: 'dirbe'
ephemeris	Solar system ephemeris to use. Default is 'builtin'. TYPE: str DEFAULT: 'builtin'

RETURNS	DESCRIPTION
number_density_grid	The tabulated zodiacal component densities. TYPE: NDArray[float64]

Source code in zodipy/number_density.py

def grid_number_density(
    x: units.Quantity[units.AU],
    y: units.Quantity[units.AU],
    z: units.Quantity[units.AU],
    obstime: time.Time,
    model: str | ZodiacalLightModel = "dirbe",
    ephemeris: str = "builtin",
) -> npt.NDArray[np.float64]:
    """Return the component-wise tabulated densities of the zodiacal components for a given grid.

    Args:
        x: x-coordinates of a cartesian mesh grid.
        y: y-coordinates of a cartesian mesh grid.
        z: z-coordinates of a cartesian mesh grid.
        obstime: Time of observation. Required to compute the Earth position.
        model: String representing a built-in model or an explicit instance of a
            `ZodiacalLightModel`. Default is 'dirbe'.
        ephemeris: Solar system ephemeris to use. Default is 'builtin'.

    Returns:
        number_density_grid: The tabulated zodiacal component densities.

    """
    if isinstance(model, str):
        zodiacal_light_model = model_registry.get_model(model)
    elif isinstance(model, ZodiacalLightModel):
        zodiacal_light_model = model
    else:
        msg = "model type must be a `str` or a `ZodiacalLightModel`."
        raise TypeError(msg)

    grid = np.asarray(np.meshgrid(x, y, z))

    earthpos_xyz = get_earthpos_inst(obstime, ephemeris=ephemeris)

    # broadcasting reshapes
    for comp in zodiacal_light_model.comps.values():
        comp.X_0 = comp.X_0.reshape(3, 1, 1, 1)

    density_partials = get_partial_number_density_func(
        comps=zodiacal_light_model.comps,
    )
    density_partials = update_partial_earth_pos(
        density_partials, earthpos_xyz[:, np.newaxis, np.newaxis, np.newaxis]
    )

    number_density_grid = np.zeros((len(zodiacal_light_model.comps), *grid.shape[1:]))
    for idx, number_density_callable in enumerate(density_partials.values()):
        number_density_grid[idx] = number_density_callable(grid)

    # Revert broadcasting reshapes
    for comp in zodiacal_light_model.comps.values():
        comp.X_0 = comp.X_0.reshape(3, 1)

    return number_density_grid