synthesizer.grid

The Grid class containing tabulated spectral and line data.

This object underpins all of Synthesizer function which generates synthetic spectra and line luminosities from models or simulation outputs. The Grid object contains attributes and methods for interfacing with spectral and line grids.

The grids themselves use a standardised HDF5 format which can be generated using the synthesizer-grids sister package.

Example usage:

from synthesizer import Grid

# Load a grid grid = Grid(“bc03_salpeter”)

# Get the axes of the grid print(grid.axes)

# Get the spectra grid print(grid.spectra)

Classes

class synthesizer.grid.Grid(grid_name, grid_dir=None, spectra_to_read=None, read_lines=True, new_lam=None, lam_lims=())[source]

The Grid class containing tabulated spectral and line data.

This object contains attributes and methods for reading and manipulating the spectral grids which underpin all spectra/line generation in synthesizer.

grid_dir

The directory containing the grid HDF5 file.

Type:

str

grid_name

The name of the grid (as defined by the file name) with no extension.

Type:

str

grid_ext

The grid extension. Either “.hdf5” or “.h5”. If the passed grid_name has no extension then “.hdf5” is assumed.

Type:

str

grid_filename

The full path to the grid file.

Type:

str

available_lines

A list of lines on the Grid.

Type:

bool/list

available_spectra

A list of spectra on the Grid.

Type:

bool/list

lam

The wavelengths at which the spectra are defined.

Type:

Quantity, float

spectra

The spectra array from the grid. This is an N-dimensional grid where N is the number of axes of the SPS grid. The final dimension is always wavelength.

Type:

dict, array-like, float

line_lams

A dictionary of line wavelengths.

Type:

dict, dist, float

line_lums

A dictionary of line luminosities.

Type:

dict, dict, float

line_conts

A dictionary of line continuum luminosities.

Type:

dict, dict, float

parameters

A dictionary containing the grid’s parameters used in its generation.

Type:

dict

axes

A list of the names of the spectral grid axes.

Type:

list, str

naxes

The number of axes the spectral grid has.

logQ10

A dictionary of ionisation Q parameters. (DEPRECATED)

Type:

dict

log10_specific_ionising_luminosity

A dictionary of log10 specific ionising luminosities.

Type:

dict

<grid_axis>

A Grid will always contain 1D arrays corresponding to the axes of the spectral grid. These are read dynamically from the HDF5 file so can be anything but usually contain at least stellar ages and stellar metallicity.

Type:

array-like, float

animate_grid(show=False, save_path=None, fps=30, spectra_type='incident')[source]

Create an animation of the grid stepping through wavelength.

Each frame of the animation is a wavelength bin.

Parameters:
  • show (bool) – Should the animation be shown?

  • save_path (str, optional) – Path to save the animation. If not specified, the animation is not saved.

  • fps (int, optional) – the number of frames per second in the output animation. Default is 30 frames per second.

Returns:

The animation object.

Return type:

matplotlib.animation.FuncAnimation

get_delta_lambda(spectra_id='incident')[source]

Calculate the delta lambda for the given spectra.

Parameters:

spectra_id (str) – Identifier for the spectra (default is “incident”).

Returns:

tuple

A tuple containing the list of wavelengths and delta lambda.

get_grid_point(**kwargs)[source]

Identify the nearest grid point for a tuple of values.

Parameters:

**kwargs (dict) – Pairs of axis names and values for the desired grid point, e.g. log10ages=9.3, log10metallicities=-2.1.

Returns:

tuple

A tuple of integers specifying the closest grid point.

get_line(grid_point, line_id, spectra_type='nebular')[source]

Create a Line object for a given line_id and grid_point.

Parameters:
  • grid_point (tuple) – A tuple of integers specifying the closest grid point.

  • line_id (str) – The id of the line.

  • spectra_type (str) – The spectra type to extract the line from. Default is “nebular”, all other spectra will have line luminosities of 0 by definition.

Returns:

line (synthesizer.line.Line)

A synthesizer Line object.

get_lines(grid_point, line_ids=None)[source]

Create a LineCollection for multiple lines.

Parameters:
  • grid_point (tuple) – A tuple of the grid point indices.

  • line_ids (list) – A list of lines, if None use all available lines.

Returns:

lines (lines.LineCollection)

static get_nearest_index(value, array)[source]

Calculate the closest index in an array for a given value.

TODO: What is this doing here!?

Parameters:
  • value (float/unyt_quantity) – The target value.

  • array (np.ndarray/unyt_array) – The array to search.

Returns:

int

The index of the closet point in the grid (array)

get_sed(spectra_type)[source]

Get the spectra grid as an Sed object.

This enables grid wide use of Sed methods for flux, photometry, indices, ionising photons, etc.

Parameters:

spectra_type (string) – The key of the spectra grid to extract as an Sed object.

Returns:

Sed

The spectra grid as an Sed object.

get_spectra(grid_point, spectra_id='incident')[source]

Create an Sed object for a specific grid point.

Parameters:
  • grid_point (tuple) – A tuple of integers specifying the closest grid point.

  • spectra_id (str) – The name of the spectra (in the grid) that is desired.

Returns:

synthesizer.sed.Sed

A synthesizer Sed object

property has_lines

Return whether the Grid has lines.

property has_spectra

Return whether the Grid has spectra.

interp_spectra(new_lam, loop_grid=False)[source]

Interpolates the spectra grid onto the provided wavelength grid.

NOTE: this will overwrite self.lam and self.spectra, overwriting the attributes loaded from the grid file. To get these back a new grid will need to instantiated with no lam argument passed.

Parameters:
  • new_lam (unyt_array/array-like, float) – The new wavelength array to interpolate the spectra onto.

  • loop_grid (bool) – flag for whether to do the interpolation over the whole grid, or loop over the first axes. The latter is less memory intensive, but slower. Defaults to False.

property lines_available

Flag for whether line emission exists.

This will only access the file the first time this property is accessed.

Returns:

True if lines are available, False otherwise.

Return type:

bool

property log10ages

Return the log10 age axis.

This is an alias to provide a pluralised version of the log10age attribute.

property log10metallicities

Return the log10 metallicity axis.

plot_specific_ionising_lum(ion='HI', hsize=3.5, vsize=None, cmap=<matplotlib.colors.ListedColormap object>, vmin=None, vmax=None, max_log10age=None)[source]

Make a simple plot of the specific ionising photon luminosity.

The resulting figure will show the (log) specific ionsing photon luminosity as a function of (log) age and metallicity for a given grid and ion.

Parameters:

ion (str) –

The desired ion. Most grids only have HI and HII calculated by

default.

hsize (float)

The horizontal size of the figure

vsize (float)

The vertical size of the figure

cmap (object/str)

A colourmap object or string defining the matplotlib colormap.

vmin (float)

The minimum specific ionising luminosity used in the colourmap

vmax (float)

The maximum specific ionising luminosity used in the colourmap

max_log10age (float)

The maximum log10(age) to plot

Returns:

matplotlib.Figure

The created figure containing the axes.

matplotlib.Axis

The axis on which to plot.

property reprocessed

Flag for whether grid has been reprocessed through cloudy.

This will only access the file the first time this property is accessed.

Returns:

True if reprocessed, False otherwise.

property shape

Return the shape of the grid.

truncate_grid_lam(min_lam, max_lam)[source]

Truncate the grid to a specific wavelength range.

If out of range wavlengths are requested, the grid will be truncated to the nearest wavelength within the grid.

Parameters:
  • min_lam (unyt_quantity) – The minimum wavelength to truncate the grid to.

  • max_lam (unyt_quantity) – The maximum wavelength to truncate the grid to.

class synthesizer.grid.Template(lam, lnu, fesc=0.0, unify_with_grid=None, **kwargs)[source]

A simplified grid contain only a single template spectra.

The model is different to all other emission models in that it scales a template by bolometric luminosity.

sed

The template spectra for the AGN.

Type:

Sed

normalisation

The normalisation for the spectra. In reality this is the bolometric luminosity.

Type:

unyt_quantity

get_spectra(bolometric_luminosity)[source]

Calculate the blackhole spectra by scaling the template.

Parameters:

bolometric_luminosity (float) – The bolometric luminosity of the blackhole(s) for scaling.