PhotoDissociation Region Toolbox — Python

Reliable astrophysics at everyday low, low prices! ®


Astrophysics Source Code Library 1102.022 Project Status: Active - The project has reached a stable, usable state and is being actively developed. Python version GNU GPL v3 License Documentation status Contributor Covenant Code of Conduct Integration test statusCode coverageCodacy quality grade

pdrtpy is the new and improved version of the formerly web-based PhotoDissociation Region Toolbox, rewritten in Python with new capabilities and giving more flexibility to end users. (The web-based /CGI version of PDRT is deprecated and no longer supported).

The PDR Toolbox is a science-enabling tool for the community, designed to help astronomers determine the physical parameters of photodissociation regions from observations. Typical observations of both Galactic and extragalactic PDRs come from ground- and space-based millimeter, submillimeter, and far-infrared telescopes such as ALMA, SOFIA, JWST, Spitzer, and Herschel. Given a set of observations of spectral line or continuum intensities, PDR Toolbox can compute best-fit FUV incident intensity and cloud density based on our models of PDR emission.

The PDR Toolbox will cover a wide range of spectral lines and metallicities and allows map-based analysis so users can quickly compute spatial images of density and radiation field from map data. We provide Jupyter Example Notebooks for data analysis. It also can support models from other PDR codes enabling comparison of derived properties between codes.

The underlying PDR model code has improved physics and chemistry. Critical updates include those discussed in Neufeld & Wolfire 2016, plus photo rates from Heays et al. 2017, oxygen chemistry rates from Kovalenko et al. 2018 and Tran et al. 2018, and carbon chemistry rates from Dagdigian 2019. We have also implemented new collisional excitation rates for [O I] from Lique et al. 2018 (and Lique private communication) and have included 13C chemistry along with the emitted line intensities for [13C II] and 13CO.

We also support fitting of temperatures and column densities to H2 excitation diagrams.

Up to date documentation can be found at pdrtpy.readthedocs.io.

What is a PDR?

Photodissociation regions (PDRs) include all of the neutral gas in the ISM where far-ultraviolet (FUV) photons dominate the chemistry and/or heating. In regions of massive star formation, PDRS are created at the boundaries between the HII regions and neutral molecular cloud, as photons with energies 6 eV < h nu < 13.6 eV. photodissociate molecules and photoionize other elements. The gas is heated from photo-electrons and cools mostly through far-infrared fine structure lines like [O I] and [C II].

For a full review of PDR physics and chemistry, see Hollenbach & Tielens 1997.

Getting Started

Installation

Requirements

pdrtpy requires Python 3 and recent versions of astropy, numpy, scipy, lmfit, and matplotlib. If you want to run the Example Notebooks, you also need jupyter.

First make sure you are using Python 3:

python --version

should show e.g., 3.7.6.

Install the package

With pip

Python has numerous ways to install packages; the easiest is with pip. The code is hosted at the Python Packaging Index, so you can type:

pip install pdrtpy

If you do not have permission to install into your Python system package area, you will need to do a user-install, which will install the package locally.

pip install --user pdrtpy

Then go ahead and install the Example Notebooks.

Example Notebooks

We have prepared Jupyter iPython notebooks with examples of how to use pdrtpy. You can download these as follows.

git clone https://github.com/mpound/pdrtpy-nb.git

If you don’t have git, you can download a zip file of the repository.

To familiarize yourself with the capabilities of pdrtpy, we suggest you do the notebooks in this order:

Getting Help & Giving Feedback

If you have a question or wish to give feedback about using PDR Toolbox or about the example notebooks, head on over to our PDR Toolbox online forum. There you can post your question and engage in discussion with the developers and other users. Feature requests from the community are welcome.

Reporting Issues

If you find a bug or something you think is in error, please report it on the github issue tracker. (You must have a Github account to submit an issue). If you aren’t sure if something is a bug or not, or if you don’t wish to create a Github account, you can post to the PDR Toolbox forum.

Contribute Code or Documentation

We welcome contributions and ideas to improve the PDR Toolbox! All contributors agree to follow our Code of Conduct . Please look at our Roadmap of Functionality to see the main new features we want to build. You can help out with those or suggest new features.

For Developers

If you plan to tinker with the code, you should fork the repo and work on your own fork. Point your browser to https://github.com/mpound/pdrtpy and click on fork in the upper right corner. After you have made your changes, create a pull request to merge them into the master branch.

You may want to use a virtual environment to protect from polluting your daily working environment (especially if you have a stable version of pdrtpy installed).

sudo apt-get install python3-venv
python -m venv ~/pdrtpy_venv
source ~/pdrtpy_venv/bin/activate[.csh]
cd pdrtpy
pip install -r requirements.txt
pip install -e .

Module Descriptions and APIs

Measurements: How you put observations to the Toolbox

To use PDR Toolbox, you need to create Measurements from your observations. A Measurement consists of a value and an error. These can be single-valued or an array of values. In the typical case of an image, the Measurement is a representation of a FITS file with two HDUs, the first HDU is the spatial map of intensity and the 2nd HDU is the spatial map of the errors. It is based on astropy’s CCDData if you are familiar with that. Typical sub-millimeter maps we get from telescopes don’t have the error plane, but PDRT makes it easy for you to create one if you know the magnitude of the error. Typical FITS images will be in intensity units, equivalent to \({\rm erg~s^{-1}~cm^{-2}~sr^{-1}}\), or in \({\rm K~km~s^{-1}}\). For the latter, PDRT will do appropriate conversion as necessary when it uses your images (the original Measurement remains unchanged).

For example how to use Measurements, see the notebook PDRT_Example_Measurements.ipynb.


Manage spectral line or continuum observations

class pdrtpy.measurement.Measurement(*args, **kwargs)[source]

Bases: CCDData

Measurement represents one or more observations of a given spectral line or continuum. It is made up of a value array, an uncertainty array, units, and a string identifier. It is based on astropy.nddata.CCDData. It can represent a single pixel observation or an image. Mathematical operations using Measurements will correctly propagate errors.

Typically, Measurements will be instantiated from a FITS file by using the the read() or make_measurement() methods. For a list of recognized spectral line identifiers, see supported_lines().

Parameters:
Raises:

TypeError – if beam parameters are not Quantities

Measurements can also be instantiated by the read(\*args, \**kwargs), to create an Measurement instance based on a FITS file. This method uses fits_measurement_reader() with the provided parameters. Example usage:

from pdrtpy.measurement import Measurement

my_obs = Measurement.read("file.fits",identifier="CII_158")
my_other_obs = Measurement.read("file2.fits",identifier="CO2_1",
                                 unit="K km/s",
                                 bmaj=9.3*u.arcsec,
                                 bmin=14.1*u.arcsec,
                                 bpa=23.2*u.degrees)

By default image axes with only a single dimension are removed on read. If you do not want this behavior, used read(squeeze=False). See also: astropy.nddata.CCDData.

Attributes:
SN

Return the signal to noise ratio (value/error)

beam

Return the beam parameters as astropy Quantities or None if beam is not set

data

~numpy.ndarray-like : The stored dataset.

dtype

numpy.dtype of this object’s data.

error

Return the underlying error array

filename

The FITS file that created this measurement, or None if it didn’t originate from a file

flags
header
id

Return the string ID of this measurement, e.g., CO_10

levels
mask

any type : Mask for the dataset, if any.

meta

dict-like : Additional meta information about the dataset.

ndim

integer dimensions of this object’s data

psf

Image representation of the PSF for the dataset.

shape

shape tuple of this object’s data.

size

integer size of this object’s data.

title

A formatted title (e.g., LaTeX) that can be in plotting.

uncertainty

any type : Uncertainty in the dataset, if any.

unit

~astropy.units.Unit : Unit for the dataset, if any.

value

Return the underlying data array

wcs

any type : A world coordinate system (WCS) for the dataset, if any.

Methods

add(other)

Add this Measurement to another, propagating errors, units, and updating identifiers.

convert_unit_to(unit[, equivalencies])

Returns a new NDData object whose values have been converted to a new unit.

copy()

Return a copy of the CCDData object.

divide(other)

Divide this Measurement by another, propagating errors, units, and updating identifiers.

from_table(filename[, format, array])

Table file reader for Measurement class.

get(world_x, world_y[, log])

Get the value(s) at the give world coordinates

get_pixel(world_x, world_y)

Return the nearest pixel coordinates to the input world coordinates

identifier(id)

Set the string ID of this measurement, e.g., CO_10

is_ratio()

Indicate if this Measurement is a ratio.

is_single_pixel()

Is this Measurement a single value? :returns: True if a single value (pixel) :rtype: bool

make_measurement(datafile, error, outfile[, ...])

Create a FITS files with 2 HDUS, the first being the datavalue and the 2nd being the data uncertainty.

multiply(other)

Multiply this Measurement by another, propagating errors, units, and updating identifiers.

read

subtract(other)

Subtract another Measurement from this one, propagating errors, units, and updating identifiers.

to_hdu([hdu_mask, hdu_uncertainty, ...])

Creates an HDUList object from a CCDData object.

write(filename, **kwd)

Write this Measurement to a FITS file with value in 1st HDU and error in 2nd HDU.

property SN

Return the signal to noise ratio (value/error)

Return type:

numpy.ndarray

add(other)[source]

Add this Measurement to another, propagating errors, units, and updating identifiers. Masks are logically or’d.

Parameters:

other (Measurement or number) – a Measurement or number to add

property beam

Return the beam parameters as astropy Quantities or None if beam is not set

divide(other)[source]

Divide this Measurement by another, propagating errors, units, and updating identifiers. Masks are logically or’d.

Parameters:

other (Measurement or number) – a Measurement or number to divide by

property error

Return the underlying error array

Return type:

numpy.ndarray

property filename

The FITS file that created this measurement, or None if it didn’t originate from a file

Return type:

str or None

static from_table(filename, format='ipac', array=False)[source]

Table file reader for Measurement class. Create one or more Measurements from a table. The input table header must contain the columns:

data - the data value

uncertainty - the error on the data, can be absolute error or percent. If percent, the header unit row entry for this column must be “%”

identifier - the identifier of this Measurement which should match a model in the ModelSet you are using, e.g., “CII_158” for [C II] 158 $\mu$m

The following columns are optional:

bmaj - beam major axis size

bmin - beam minor axis size

bpa - beam position angle

The table must specify the units of each column, e.g. a unit row in the header for IPAC format. Leave column entry blank if unitless. Units of value and error should be the same or conformable. Units must be transformable to a valid astropy.unit.Unit.

Parameters:
  • filename (str) – Name of table file.

  • formatAstropy Table format format. e.g., ascii, ipac, votable. Default is IPAC format

  • array (bool) – Controls whether a list of Measurements or a single Measurement is returned. If array is True, one Measurement instance will be created for each row in the table and a Python list of Measurements will be returned. If array is False, one Measurement containing all the points in the data member will be returned. If array is False, the identifier and beam parameters of the first row will be used. If feeding the return value to a plot method such as phasespace(), choose array=False. Default:False.

Return type:

Measurement or list of Measurement

get(world_x, world_y, log=False)[source]

Get the value(s) at the give world coordinates

Parameters:
  • world_x (float or array-like) – the x value in world units of naxis1

  • world_y (float or array-lke) – the y value in world units of naxis2

  • log (bool) – True if the input coords are logarithmic Default:False

Returns:

The value(s) of the Measurement at input coordinates

Return type:

float

get_pixel(world_x, world_y)[source]

Return the nearest pixel coordinates to the input world coordinates

Parameters:
  • world_x (float) – The horizontal world coordinate

  • world_y (float) – The vertical world coordinate

property id

Return the string ID of this measurement, e.g., CO_10

Return type:

str

identifier(id)[source]

Set the string ID of this measurement, e.g., CO_10

Parameters:

id (str) – the identifier

is_ratio()[source]

Indicate if this Measurement is a ratio.. This method looks for the ‘/’ past the first character of the` Measurement` identifier, such as “CII_158/CO_32” See also pdrutils.is_ratio(string)

Returns:

True if the Measurement is a ratio, False otherwise

Return type:

bool

is_single_pixel()[source]

Is this Measurement a single value? :returns: True if a single value (pixel) :rtype: bool

property levels
static make_measurement(datafile, error, outfile, rms=None, masknan=True, overwrite=False, unit='adu')[source]

Create a FITS files with 2 HDUS, the first being the datavalue and the 2nd being the data uncertainty. This format makes allows the resulting file to be read into the underlying :class:’~astropy.nddata.CCDData` class.

Parameters:
  • datafile (str) – The FITS file containing the data as a function of spatial coordinates

  • error (str) –

    The errors on the data Possible values for error are:

    • a filename with the same shape as datafile containing the error values per pixel

    • a percentage value ‘XX%’ must have the “%” symbol in it

    • ’rms’ meaning use the rms parameter if given, otherwise look for the RMS keyword in the FITS header of the datafile

  • outfile (str) – The output file to write the result in (FITS format)

  • rms (float or astropy.units.Unit) – If error == ‘rms’, this value may give the rms in same units as data (e.g ‘erg s-1 cm-2 sr-1’).

  • masknan (bool) – Whether to mask any pixel where the data or the error is NaN. Default:true

  • overwrite (bool) – If True, overwrite the output file if it exists. Default: False.

  • unit (astropy.units.Unit or str) – Intensity unit to use for the data, this will override BUNIT in header if present.

Raises:
  • Exception – on various FITS header issues

  • OSError – if overwrite is False and the output file exists.

Example usage:

# example with percentage error
Measurement.make_measurement("my_infile.fits",error='10%',outfile="my_outfile.fits")

# example with measurement in units of K km/s and error
# indicated by RMS keyword in input file.
Measurement.make_measurement("my_infile.fits",error='rms',outfile="my_outfile.fits",unit="K km/s",overwrite=True)
multiply(other)[source]

Multiply this Measurement by another, propagating errors, units, and updating identifiers. Masks are logically or’d.

Parameters:

other (Measurement or number) – a Measurement or number to multiply

subtract(other)[source]

Subtract another Measurement from this one, propagating errors, units, and updating identifiers. Masks are logically or’d.

Parameters:

other (Measurement or number) – a Measurement or number to subtract

property title

A formatted title (e.g., LaTeX) that can be in plotting.

Return type:

str or None

property value

Return the underlying data array

Return type:

numpy.ndarray

write(filename, **kwd)[source]

Write this Measurement to a FITS file with value in 1st HDU and error in 2nd HDU. See astropy.nddata.CCDData.write().

Parameters:
  • filename (str) – Name of file.

  • kwd – All additional keywords are passed to astropy.io.fits

The available built-in formats are:

Format

Read

Write

Auto-identify

pdrtpy.measurement.fits_measurement_reader(filename, hdu=0, unit=None, hdu_mask='MASK', hdu_flags=None, key_uncertainty_type='UTYPE', **kwd)[source]

FITS file reader for Measurement class, which will be called by Measurement.read().

Parameters:
  • filename (str) – Name of FITS file.

  • identifier (str) – string indicating what this is an observation of, e.g., “CO_10” for CO(1-0)

  • squeeze (bool) – If True, remove single dimension axes from the input image. Default: True

  • hdu (int, optional) – FITS extension from which Measurement should be initialized. If zero and and no data in the primary extension, it will search for the first extension with data. The header will be added to the primary header. Default is 0.

  • unit (astropy.units.Unit, optional) – Units of the image data. If this argument is provided and there is a unit for the image in the FITS header (the keyword BUNIT is used as the unit, if present), this argument is used for the unit. Default is None.

  • hdu_uncertainty (str or None, optional) – FITS extension from which the uncertainty should be initialized. If the extension does not exist the uncertainty of the Measurement is None. Default is 'UNCERT'.

  • hdu_mask (str or None, optional) – FITS extension from which the mask should be initialized. If the extension does not exist the mask of the Measurement is None. Default is 'MASK'.

  • hdu_flags (str or None, optional) – Currently not implemented. Default is None.

  • kwd – Any additional keyword parameters are passed through to the FITS reader in astropy.io.fits

Raises:

TypeError – If the conversion from CCDData to Measurement fails

ModelSets: The interface to models in the Toolbox

PDRT supports a variety of PDR models to be used to fit your data. These are represented in the Python class ModelSet. Broadly three classes are available:

  1. Wolfire/Kaufman 2020 models for constant density media (metallicities Z=0.5,1)

  2. Wolfire/Kaufman 2006 models for constant density media (Z=0.1,1,3)

  3. Kosma-\(\tau\) 2013 models for clumpy and non-clumpy media (Z=1)

Models are stored in FITS format as ratios of intensities as a function of radiation field and hydrogen nucleus volume density.

For example how to use ModelSets, see the notebook PDRT_Example_ModelSets.ipynb


Manage pre-computed PDR models

class pdrtpy.modelset.ModelSet(name, z, medium='constant density', mass=None, modelsetinfo=None, format='ipac')[source]

Bases: object

Class for computed PDR Model Sets. ModelSet provides interface to a directory containing the model FITS files and the ability to query details about.

Parameters:
  • name (str) – identifying name, e.g., ‘wk2006’

  • z (float) – metallicity in solar units.

  • medium (str) – medium type, e.g. ‘constant density’, ‘clumpy’, ‘non-clumpy’

  • mass – maximum clump mass (for KosmaTau models). Default:None (appropriate for Wolfire/Kaufman models)

Params modelsetinfo:

For adding user specified ModelSet, this parameter specifies the file with information about the ModelSet. It can be a pathname to an externalt tabular file in an astropy-recognized or an astropy Table object. If an external tabular file, its format should be give by the format keyword. The columns are:

[‘PDR code’,’name’, ‘version’,’path’,’filename’,’medium’,’z’,’mass’,’description’]

z and mass should be floats, the rest should be strings. The ‘name`, version, medium, z, and mass columns contain the values available in the input ModelSet as described above. PDR code is the originator of the code e.g., “KOSMA-tau”, version is the code version, path is the full-qualified pathname to the FITS files, filename is the tabular file which contains descriptions of individual FITS files. This must also be in the format specified by the format keyword. description is a description of the input ModelSet. :type modelsetinfo: str or :class:~astropy.Table :format: If modelsetinfo is a file, give the format of file. Must be an astropy recognized Table format. e.g., ascii, ipac, votable. Default is IPAC format :type format: str :raises ValueError: If model set not recognized/found.

Attributes:
code

The PDR code that created this ModelSet, e.g.

description

The description of this model

identifiers

Table of lines and continuum that are included in ratio models of this ModelSet.

is_wk2006

method to indicate this is a wk2006 model, to deal with quirks

medium

The medium type of this model, e.g.

metallicity

The metallicity of this ModelSet

name

The name of this ModelSet

supported_intensities

Table of lines and continuum that are covered by this ModelSet and have models separate from the any ratio model they might be in.

supported_ratios

The emission ratios that are covered by this ModelSet

table

The table containing details of the models in this ModelSet.

user_added_models

Show which models have been added to this ModelSet by the user

version

The version of this ModelSet

z

The metallicity of this ModelSet

Methods

add_model(identifier, model, title[, overwrite])

Add your own model to this ModelSet.

all_sets()

Return a table of the names and descriptions of available ModelSets (not just this one)

find_files(m[, ext])

Find the valid model ratios files in this ModelSet for a given list of measurement IDs.

find_pairs(m)

Find the valid model ratios labels in this ModelSet for a given list of measurement IDs

get_model(identifier[, unit, ext])

Get a specific model by its identifier

get_models(identifiers[, model_type, ext])

get the models from thie ModelSet that match the input list of identifiers

list()

Print the names and descriptions of available ModelSets (not just this one)

model_intensities(m)

Return the model intensities in this ModelSet that match the input Measurement ID list.

model_ratios(m)

Return the model ratios that match the input Measurement ID list.

ratiocount(m)

The number of valid ratios in this ModelSet, given a list of observation (Measurement) identifiers.

add_model(identifier, model, title, overwrite=False)[source]

Add your own model to this ModelSet.

param identifier:

a Measurement ID. It can be an intensity or a ratio, e.g., “CII_158”,”CI_609/FIR”.

type identifier:

str

param model:

the model to add. If a string, this must be the fully-qualified path of a FITS file. If a Measurement it must have the same CTYPEs and CUNITs as the models in the ModelSet(?).

type model:

str or Measurement

param title:

A formatted string (e.g., LaTeX) describing this observation that can be used for plotting. Python r-strings are accepted, e.g., r’$^{13}$CO(3-2)’ would give \(^{13}{\rm CO(3-2)}\).

static all_sets()[source]

Return a table of the names and descriptions of available ModelSets (not just this one)

Return type:

Table

property code

The PDR code that created this ModelSet, e.g. ‘KOSMA-tau’

Return type:

str

property description

The description of this model

Return type:

str

find_files(m, ext='fits')[source]

Find the valid model ratios files in this ModelSet for a given list of measurement IDs. See id()

Parameters:
  • m (list) – list of string Measurement IDs, e.g. [“CII_158”,”OI_145”,”FIR”]

  • ext (str) – file extension. Default: “fits”

Returns:

An iterator of model ratio files for the given list of Measurement IDs

Return type:

iterator

find_pairs(m)[source]

Find the valid model ratios labels in this ModelSet for a given list of measurement IDs

Parameters:

m (list) – list of string Measurement IDs, e.g. [“CII_158”,”OI_145”,”FIR”]

Returns:

An iterator of model ratios labels for the given list of measurement IDs

Return type:

iterator

get_model(identifier, unit=None, ext='fits')[source]

Get a specific model by its identifier

Parameters:

identifier (str) – a Measurement ID. It can be an intensity or a ratio, e.g., “CII_158”,”CI_609/FIR”.

Returns:

The model matching the identifier

Return type:

Measurement

Raises:

KeyError if identifier not found in this ModelSet

get_models(identifiers, model_type='ratio', ext='fits')[source]

get the models from thie ModelSet that match the input list of identifiers

Parameters:
  • identifiers (list) – list of string Measurement IDs, e.g., [“CII_158”,”OI_145”,”CS_21”]

  • model_type (str) – indicates which type of model is requested one of ‘ratio’ or ‘intensity’

Returns:

The matching models as a list of Measurement.

Return type:

list

Raises:

KeyError if identifiers not found in this ModelSet

property identifiers

Table of lines and continuum that are included in ratio models of this ModelSet. Only lines and continuum that are part of ratios are included in this list. For a separate list of line and continuum intensity models see supported_intensities().

Return type:

astropy.table.Table

property is_wk2006

method to indicate this is a wk2006 model, to deal with quirks of that modelset

Returns:

True if it is.

static list()[source]

Print the names and descriptions of available ModelSets (not just this one)

property medium

The medium type of this model, e.g. ‘constant density’, ‘clumpy’, ‘non-clumpy’

Return type:

str

property metallicity

The metallicity of this ModelSet

Return type:

float

model_intensities(m)[source]

Return the model intensities in this ModelSet that match the input Measurement ID list. This method will return the intersection of the input list and the list of supported lines.

Parameters:

m (list) – list of string Measurement IDs, e.g., [“CII_158”,”OI_145”,”CS_21”]

Returns:

list of string identifiers of ratios IDs, e.g., [‘CII_158’,’OI_145’]

Return type:

list

model_ratios(m)[source]

Return the model ratios that match the input Measurement ID list. You must provide at least 2 Measurements IDs

Parameters:

m (list) – list of string Measurement IDs, e.g., [“CII_158”,”OI_145”,”FIR”]

Returns:

list of string identifiers of ratios IDs, e.g., [‘OI_145/CII_158’, ‘OI_145+CII_158/FIR’]

Return type:

list

property name

The name of this ModelSet

Return type:

str

ratiocount(m)[source]

The number of valid ratios in this ModelSet, given a list of observation (Measurement) identifiers.

Parameters:

m (list) – list of string Measurement IDs, e.g. [“CII_158”,”OI_145”,”FIR”]

Returns:

The number of model ratios found for the given list of measurement IDs

Return type:

int

property supported_intensities

Table of lines and continuum that are covered by this ModelSet and have models separate from the any ratio model they might be in.

Return type:

astropy.table.Table

property supported_ratios

The emission ratios that are covered by this ModelSet

Return type:

astropy.table.Table

property table

The table containing details of the models in this ModelSet.

Return type:

astropy.table.Table

property user_added_models

Show which models have been added to this ModelSet by the user

Return type:

list

property version

The version of this ModelSet

Return type:

str

property z

The metallicity of this ModelSet

Return type:

float

Plotting Tools: Display models and data

The plot module provides mechanisms for plotting models, observations, and model fits.

The ModelPlot class can be used plotting models and observations without any \(\chi^2\) fitting. An example notebook for using ModelPlot is PDRT_Example_ModelPlotting.ipynb .

Some classes are paired with analysis tools in the tool module. LineRatioPlot which is used to plot the results of LineRatioFit, and H2ExcitationPlot that is used in H2Excitation. All plot classes are derived from PlotBase.

Plot Keywords

To manage the plots, the methods in Plot classes take keywords (**kwargs) that turn on or off various options, specify plot units, or map to matplotlib’s plot(), imshow(), contour() keywords. The methods have reasonable defaults, so try them with no keywords to see what they do before modifying keywords.

  • units (str or astropy.units.Unit) image data units to use in the plot. This can be either a string such as, ‘cm^-3’ or ‘Habing’, or it can be an astropy.units.Unit. Data will be converted to the desired unit. Note these are not the axis units, but the image data units. Modifying axis units is implemented via the xaxis_unit and yaxis_unit keywords.

  • image (bool) whether or not to display the image map (imshow).

  • show (str) which quantity to display in the Measurement, one of ‘data’, ‘error’, ‘mask’. For example, this can be used to plot the errors in observed ratios. Default: ‘data’

  • cmap (str) colormap name, Default: ‘plasma’

  • colorbar (str) whether or not to display colorbar

  • colors (str) color of the contours. Default: ‘whitecolor of the contours. Default: ‘white’

  • contours (bool), whether or not to plot contours

  • label (bool), whether or not to label contours

  • linewidths (float or sequence of float), the line width in points, Default: 1.0

  • legend (bool) Draw a legend on the plot. If False, a title is drawn above the plot with the value of the title keyword

  • bbox_to_anchor (tuple) The matplotlib legend keyword for controlling the placement of the legend. See the matplotlib Legend Guide

  • loc (str) The matplotlib legend keyword for controlling the location of the legend. See legend().

  • levels (int or array-like) Determines the number and positions of the contour lines / regions. If an int n, use n data intervals; i.e. draw n+1 contour lines. The level heights are automatically chosen. If array-like, draw contour lines at the specified levels. The values must be in increasing order.

  • measurements (array-like) A list of single pixel Measurements that can be contoured over a model ratio or intensity map.

  • meas_color (array of str) A list of colors to use when overlaying Measurement contours. There should be one color for each Measurement in the measurement keyword. The Default of None will use a color-blind friendly color cycle.

  • norm (str or astropy.visualization normalization object) The normalization to use in the image. The string ‘simple’ will normalize with simple_norm() and ‘zscale’ will normalize with IRAF’s zscale algorithm. See ZScaleInterval.

  • stretch (str) {‘linear’, ‘sqrt’, ‘power’, ‘log’, ‘asinh’}. The stretch function to apply to the image for simple_norm. The Default is ‘linear’.

  • aspect (str) aspect ratio, ‘equal’ or ‘auto’ are typical defaults.

  • origin (str) Origin of the image. Default: ‘lower’

  • title (str) A title for the plot. LaTeX allowed.

  • vmin (float) Minimum value for colormap normalization

  • vmax (float) Maximum value for colormap normalization

  • xaxis_unit (str or astropy.units.Unit) X axis (density) units to use when plotting models, such as in overlay_all_ratios() or modelratio(). If None, the native model axis units are used.

  • yaxis_unit (str or astropy.units.Unit) Y axis (density) units to use when plotting models, such as in overlay_all_ratios() or modelratio(). If None, the native model axis units are used.

The following keywords are available, but you probably won’t touch.

  • nrows (int) Number of rows in the subplot

  • ncols (int) Number of columns in the subplot

  • index (int) Index of the subplot

  • reset (bool) Whether or not to reset the figure.

Providing keywords other than these has undefined results, but may just work!


PlotBase

class pdrtpy.plot.plotbase.PlotBase(tool)[source]

Bases: object

Base class for plotting.

Parameters:

tool (Any class derived from ToolBase) – Reference to a tool object or None. This is used for classes that inherit from PlotBase and are coupled to a specific tool, e.g. LineRatioPlot and LineRatioFit.

Attributes:
axis

The last axis that was drawn.

figure

The last figure that was drawn.

Methods

colorcycle(colorcycle)

Set the plot color cycle for multi-trace plots.

reset_colorcycle()

Reset the color cycle to the default color-blind friendly one

savefig(fname, **kwargs)

Save the current figure to a file.

text(x, y, s[, fontdict])

Add text to the Axes.

usetex(use)

Control whether plots delegate rendering of fancy text components in axis labels and elsewhere to the system version of LaTeX or use matplotlib's rendering.

property axis

The last axis that was drawn.

Return type:

matplotlib.axes._subplots.AxesSubplot

colorcycle(colorcycle)[source]

Set the plot color cycle for multi-trace plots. The default color cycle is optimized for color-blind users.

Parameters:

colorcycle (list) – List of colors to use, typically a list of hex color strings. This list will be passed to matplotlib.pyplot.rc() as the axes prop_cycle parameter using matplotlib.cycler.

property figure

The last figure that was drawn.

Return type:

matplotlib.figure.Figure

reset_colorcycle()[source]

Reset the color cycle to the default color-blind friendly one

savefig(fname, **kwargs)[source]

Save the current figure to a file.

Parameters:

fname (str) – filename to save in

Keyword Arguments:

Additional arguments (**kwargs) are passed to matplotlib.pyplot.savefig(). e.g. bbox_inches=’tight’ for a tight layout.

text(x, y, s, fontdict=None, **kwargs)[source]

Add text to the Axes. Add the text s to the Axes at location x, y in data coordinates. This calls through to matplotlib.pyplot.text().

Parameters:
  • x (float) – the horizontal coordinate for the text

  • y (float) – the vertical coordinate for the text

  • s (str) – the text

  • fontdict (dict) – A dictionary to override the default text properties. If fontdict is None, the defaults are determined by rcParams.

  • **kwargs – Other miscellaneous Text parameters.

usetex(use)[source]

Control whether plots delegate rendering of fancy text components in axis labels and elsewhere to the system version of LaTeX or use matplotlib’s rendering. This method sets matplotlib parameter rcParams[“text.usetex”] in the local pyplot instance. Note: You must have LaTeX installed on your system if setting this to True or an exception will be raised when you try to plot.

Parameters:

use (bool) – whether to use LaTeX or not

ExcitationPlot

class pdrtpy.plot.excitationplot.ExcitationPlot(tool, label)[source]

Bases: PlotBase

ExcitationPlot creates excitation diagrams using the results of H2ExcitationFit. It can plot the observed excitation diagram with or without fit results, and allows averaging over user-given spatial areas.

Attributes:
axis

The last axis that was drawn.

figure

The last figure that was drawn.

Methods

colorcycle(colorcycle)

Set the plot color cycle for multi-trace plots.

column_density(component[, log])

Plot the column density of hot or cold gas component, or total column density.

ex_diagram([position, size, norm, show_fit])

Plot the excitation diagram.

explore([data, interaction_type])

Explore the fitted parameters of a map.

opr(**kwargs)

Plot the ortho-to-para ratio.

reset_colorcycle()

Reset the color cycle to the default color-blind friendly one

savefig(fname, **kwargs)

Save the current figure to a file.

temperature(component, **kwargs)

Plot the temperature of hot or cold gas component.

text(x, y, s[, fontdict])

Add text to the Axes.

usetex(use)

Control whether plots delegate rendering of fancy text components in axis labels and elsewhere to the system version of LaTeX or use matplotlib's rendering.

column_density(component, log=True, **kwargs)[source]

Plot the column density of hot or cold gas component, or total column density.

Parameters:
  • component (str) – ‘hot’, ‘cold’, or ‘total

  • log – take the log10 of the column density before plotting

ex_diagram(position=None, size=None, norm=True, show_fit=False, **kwargs)[source]

Plot the excitation diagram. For maps of excitation parameters, a position and optional size are required. To examine the entire map, use explore().

Parameters:
  • position (tuple) – The spatial position of excitation diagram. For spatial averaging this is the cutout array’s center with respect to the data array. The position is specified as a (x, y) tuple of pixel coordinates.

  • size (int, array_like, or astropy.units.Quantity) – The size of the cutout array along each axis. If size is a scalar number or a scalar Quantity, then a square cutout of size will be created. If size has two elements, they should be in (ny, nx) order. Scalar numbers in size are assumed to be in units of pixels. size can also be a Quantity object or contain Quantity objects. Such Quantity objects must be in pixel or angular units. For all cases, size will be converted to an integer number of pixels, rounding the the nearest integer. See Cutout2D

  • norm (bool) – if True, normalize the column densities by the statistical weight of the upper state, \(g_u\).

  • show_fit (bool) – Show the most recent fit done the the associated H2ExcitationFit tool.

explore(data=None, interaction_type='click', **kwargs)[source]

Explore the fitted parameters of a map. A user-requested map is displayed in the left panel and in the right panel is the fitted excitation diagram for a point selected by the user. The user clicks on a point in the left panel and the right panel will update with the excitation diagram for that point.

Parameters:
  • data (Measurement) – A reference image to use for the left panel, e.g. the total column density, the cold temperature, etc. This should be a reference results in the H2Excitation tool used for this ExcitationPlot (e.g., htool.temperature[‘cold’])

  • interaction_type (str) – whether to use mouse click or mouse move to update the right hand panel. Valid values are ‘click’ or ‘move’.

  • **kwargs

    Other parameters passed to _plot(), ex_diagram(), or matplotlib methods.

    • units,image, contours, label, title, norm, figsize – See the general Plot Keywords documentation

    • show_fit - show the fit in the excitation diagram, Default: True

    • log - plot the log10 of the image, can be useful for column density, Default: False

    • markersize - size of the marker displayed where clicked, in points, Default: 20

    • fmt - matplotlib format for the marker, Default:. ‘r+’

opr(**kwargs)[source]

Plot the ortho-to-para ratio. This will be a map if the input data are a map, otherwise a float value is returned.

temperature(component, **kwargs)[source]

Plot the temperature of hot or cold gas component.

Parameters:

component (str) – ‘hot’ or ‘cold’

LineRatioPlot

class pdrtpy.plot.lineratioplot.LineRatioPlot(tool)[source]

Bases: PlotBase

LineRatioPlot plots the results of LineRatioFit. It can plot maps of fit results, observations with errors on top of models, chi-square and confidence intervals and more.

Keyword Arguments:

The methods of this class can take a variety of optional keywords. See the general Plot Keywords documentation

Attributes:
axis

The last axis that was drawn.

figure

The last figure that was drawn.

Methods

chisq(**kwargs)

Plot the \(\chi^2\) map that was computed by the LineRatioFit tool.

colorcycle(colorcycle)

Set the plot color cycle for multi-trace plots.

confidence_intervals(**kwargs)

Plot the confidence intervals from the \(\chi^2\) map computed by the LineRatioFit tool.

density(**kwargs)

Plot the hydrogen nucleus volume density map that was computed by LineRatioFit tool.

modelintensity(id, **kwargs)

Plot one of the model intensities

modelratio(id, **kwargs)

Plot one of the model ratios

observedratio(id, **kwargs)

Plot one of the observed ratios

overlay_all_ratios(**kwargs)

Overlay all the measured ratios and their errors on the \((n,F_{FUV})\) space.

radiation_field(**kwargs)

Plot the radiation field map that was computed by LineRatioFit tool.

ratios_on_models(**kwargs)

Overlay all the measured ratios and their errors on the individual models for those ratios.

reduced_chisq(**kwargs)

Plot the reduced \(\chi^2\) map that was computed by the LineRatioFit tool.

reset_colorcycle()

Reset the color cycle to the default color-blind friendly one

savefig(fname, **kwargs)

Save the current figure to a file.

show_both([units])

Plot both radiation field and volume density maps computed by the LineRatioFit tool in a 1x2 panel subplot.

text(x, y, s[, fontdict])

Add text to the Axes.

usetex(use)

Control whether plots delegate rendering of fancy text components in axis labels and elsewhere to the system version of LaTeX or use matplotlib's rendering.

chisq(**kwargs)[source]

Plot the \(\chi^2\) map that was computed by the LineRatioFit tool.

confidence_intervals(**kwargs)[source]

Plot the confidence intervals from the \(\chi^2\) map computed by the LineRatioFit tool. Default levels: [50., 68., 80., 95., 99.]

Currently only works for single-pixel Measurements

density(**kwargs)[source]

Plot the hydrogen nucleus volume density map that was computed by LineRatioFit tool. Default units: cm \(^{-3}\)

modelintensity(id, **kwargs)[source]

Plot one of the model intensities

Parameters:
  • id (str) – the intensity identifier, such as CO_32`.

  • **kwargs – see class documentation above

Raises:

KeyError – if is id not in existing model intensities

modelratio(id, **kwargs)[source]

Plot one of the model ratios

Parameters:
  • id (str) – the ratio identifier, such as CII_158/CO_32.

  • **kwargs – see class documentation above

Raises:

KeyError – if is id not in existing model ratios

observedratio(id, **kwargs)[source]

Plot one of the observed ratios

Parameters:

id – the ratio identifier, such as CII_158/CO_32.

Raises:

KeyError – if id is not in existing observed ratios

overlay_all_ratios(**kwargs)[source]

Overlay all the measured ratios and their errors on the \((n,F_{FUV})\) space.

This only works for single-valued Measurements; an overlay for multi-pixel doesn’t make sense.

radiation_field(**kwargs)[source]

Plot the radiation field map that was computed by LineRatioFit tool. Default units: Habing.

ratios_on_models(**kwargs)[source]

Overlay all the measured ratios and their errors on the individual models for those ratios. Plots are displayed in multi-column format, controlled the ncols keyword. Default: ncols=2

Currently only works for single-pixel Measurements

reduced_chisq(**kwargs)[source]

Plot the reduced \(\chi^2\) map that was computed by the LineRatioFit tool.

show_both(units=['Habing', 'cm^-3'], **kwargs)[source]

Plot both radiation field and volume density maps computed by the LineRatioFit tool in a 1x2 panel subplot. Default units: [‘Habing’,’cm^-3’]

ModelPlot

class pdrtpy.plot.modelplot.ModelPlot(modelset, figure=None, axis=None)[source]

Bases: PlotBase

ModelPlot is a tool for exploring sets of models. It can plot individual intensity or ratio models, phase-space diagrams, and optionally overlay observations. Units are seamlessly transformed, so you can plot in Habing units, Draine units, or any conformable quantity. ModelPlot does not require model fitting with LineRatioFit first.

Keyword Arguments:

The methods of this class can take a variety of optional keywords. See the general Plot Keywords documentation.

Attributes:
axis

The last axis that was drawn.

figure

The last figure that was drawn.

Methods

colorcycle(colorcycle)

Set the plot color cycle for multi-trace plots.

intensity(identifier, **kwargs)

Plot a model ratio

isoplot(identifier, plotnaxis[, nax_clip])

Plot lines of constant model parameter as a function of the other model parameter and a model intensity or ratio

overlay(measurements, **kwargs)

Overlay one or more single-pixel measurements in the model space \((n,F_{FUV})\).

phasespace(identifiers[, nax1_clip, ...])

Plot lines of constant density and radiation field on a ratio-ratio, ratio-intensity, or intensity-intensity map

plot(identifier, **kwargs)

Plot a model intensity or ratio

ratio(identifier, **kwargs)

Plot a model ratio

reset_colorcycle()

Reset the color cycle to the default color-blind friendly one

savefig(fname, **kwargs)

Save the current figure to a file.

text(x, y, s[, fontdict])

Add text to the Axes.

usetex(use)

Control whether plots delegate rendering of fancy text components in axis labels and elsewhere to the system version of LaTeX or use matplotlib's rendering.

intensity(identifier, **kwargs)[source]

Plot a model ratio

Parameters:

identifier (str) – Identifier tag for the model to plot, e.g., “OI_63”, “CII_158”, “CO_10”]

See also

supported_intensities() for a list of available identifer tags

isoplot(identifier, plotnaxis, nax_clip=<Quantity [  1000., 100000.] K>, **kwargs)[source]

Plot lines of constant model parameter as a function of the other model parameter and a model intensity or ratio

Parameters:
  • identifier (str) – identifier tag for the model to plot, e.g., “OI_63/CO_21” or “CII_158”

  • plotnaxis (int) – Which NAXIS to use to compute lines of constant value. Since models have two axes, this must be either 1 or 2

  • nax_clip (array-like, must contain Quantity) – The range of model parameters on NAXIS{plotnaxis} to show in the plot. Must be given as a range of astropy quanitities, e.g. [10,1E7]*Unit(“cm-3”). Default is None which means plot full range of the parameter.

  • step (int) – Allows skipping of lines of constant value, e.g. plot every step-th value. Useful when parameter space is crowded, and a cleaner looking plot is desired. Default: 1 – plot every value

overlay(measurements, **kwargs)[source]

Overlay one or more single-pixel measurements in the model space \((n,F_{FUV})\).

Parameters:
  • measurements (list) – a list of one or more Measurement to overlay.

  • shading (float) – Controls how measurements and errors are drawn. If shading is zero, Measurements will be drawn in solid contour for the value and dashed for the +/- errors. If shading is between 0 and 1, Measurements are drawn with as filled contours representing the size of the errors (see matplotlib.pyplot.contourf()) with alpha set to the shading value. Default value: 0.4

phasespace(identifiers, nax1_clip=<Quantity [1.e+01, 1.e+07] 1 / cm3>, nax2_clip=<Quantity [1.e+01, 1.e+06] Habing>, reciprocal=[False, False], **kwargs)[source]

Plot lines of constant density and radiation field on a ratio-ratio, ratio-intensity, or intensity-intensity map

Parameters:
  • identifiers (list of str) – list of two identifier tags for the model to plot, e.g., [“OI_63/CO_21”, “CII_158”]

  • nax1_clip (array-like, must contain Quantity) – The range of model densities on NAXIS1 to show in the plot. For most model NAXIS1 is hydrogen number density $n_H$ in cm$^{-3}$. For ionized gas models, it is electron temperature $T_e$ in K. Must be given as a range of astropy quanitities. Default: [10,1E7]*Unit(“cm-3”)

  • nax2_clip (array-like, must contain Quantity) – The range of model parameters on NAXIS2 to show in the plot. For most models NAXIS2 is radiation field intensities in Habing or cgs units. For ionized gas models, it is electron volume density $n_e$. Must be given as a range of astropy quantities. Default: nax1_clip=[10,1E6]*utils.habing_unit.

  • reciprocal (array-like bool) – Whether or not the plot the reciprocal of the model on each axis. Given as a pair of booleans. e.g. [False,True] means don’t flip the quantity X axis, but flip quantity the Y axis. i.e. if the model is “CII/OI”, and reciprocal=True then the axis will be “OI/CII”. Default: [False, False]

The following keywords are supported as **kwargs:

Parameters:
  • measurements (array-like of Measurement.) – A list of two Measurement, one for each identifier, that will be used to plot a data point on the grid. At least two Measurements, one for x and one for y, must be given. Subsequent Measurements must also be paired since they represent x and y, e.g [m1x, m1y, m2x, m2y,…]. Measurement data and uncertainty members may be arrays. Default: None

  • errorbar (bool) – Plot error bars when given measurements. Default: True

  • fmt (array of str) – The format to use when plotting Measurement data. There should be one for each pair of Measurements. See matplotlib.axes.Axes.plot() for examples. Default is ‘sk’ for all points.

  • label (array of str) – The label(s) to use the Measurement data legend. There should be one for each pair of Measurements. Default is ‘data’ for all points.

  • legend (bool) – Draw a legend on the plot. Default: True

  • title (str) – Title to draw on the plot. Default: None

  • linewidth (float) – line width

  • grid (bool) – show grid or not, Default: True

  • figsize (2-tuple of floats) – Figure dimensions (width, height) in inches. Default: (8,5)

  • capsize (float) – end cap length of errorbars if shown, in points. Default: 3.

  • markersize (float) – size of data point marker in points. Default: 8

plot(identifier, **kwargs)[source]

Plot a model intensity or ratio

Parameters:

identifier (str) – Identifier tag for the model to plot, e.g., “CII_158”,”OI_145”,”CO_43/CO_21’]

See also

supported_lines() for a list of available identifer tags

ratio(identifier, **kwargs)[source]

Plot a model ratio

Parameters:

identifier (str) – Identifier tag for the model to plot, e.g., “OI_63+CII_158/FIR”, “CO_43/CO_21’]

See also

supported_ratios() for a list of available identifer tags

Analysis Tools: Fit models to data

The tool module contains the analysis tools in the PDR Toolbox. All tools are derived from ToolBase.

For examples how to use LineRatioFit, see the notebooks PDRT_Example_Find_n_G0_Single_Pixel.ipynb and PDRT_Example_Make_n_G0_maps.ipynb.

For an example how to use H2ExcitationFit and ExcitationPlot see the notebook PDRT_Example_H2_Excitation.ipynb.


ToolBase

The base class of all tools. Tools have a run() method both of which subclasses must define.

class pdrtpy.tool.toolbase.ToolBase[source]

Bases: object

Base class object for PDR Toolbox tools. This class implements a simple interface with a run method. Tools will generally do some set up such as reading in observational data before run() can be invoked.

Attributes:
has_maps

Are the Measurements used map-based?.

has_scalar

Are the Measurements used scalars.

has_vectors

Are the Measurements used a Nx1 vector, e.g.

Methods

run()

Runs the tool.

property has_maps

Are the Measurements used map-based?. (i.e., have 2 spatial axes)

Returns:

True, if the observational inputs are spatial maps, False otherwise

Return type:

bool

property has_scalar

Are the Measurements used scalars.

Returns:

True, if the observational inputs are scalars, False otherwise

Return type:

bool

property has_vectors

Are the Measurements used a Nx1 vector, e.g. read in from a table with from_table().

Returns:

True, if the observational inputs are a vector, False otherwise

Return type:

bool

run()[source]

Runs the tool. Each subclass Tool must implement its own run() method.

Excitation Diagram Fitting

H2ExcitationFit is a tool for fitting temperature, column density, and ortho-to-para ratio in \(H_2\) excitation diagrams. A two temperature model is assumed, and the fit will find \(T_{hot}, T_{cold}, N_{hot}(H_2), N_{cold}(H_2),\) and optionally OPR. The base class ExcitationFit can be used to create a tool to fit a different molecule.

class pdrtpy.tool.h2excitation.ExcitationFit(measurements=None, constantsfile=None)[source]

Bases: ToolBase

Base class for creating excitation fitting tools for various species.

Parameters:

measurements (list of Measurement.) – Input measurements to be fit.

Attributes:
has_maps

Are the Measurements used map-based?.

has_scalar

Are the Measurements used scalars.

has_vectors

Are the Measurements used a Nx1 vector, e.g.

Methods

add_measurement(m)

Add an intensity Measurement to internal dictionary used to compute the excitation diagram.

remove_measurement(identifier)

Delete a measurement from the internal dictionary used to compute column densities.

replace_measurement(m)

Safely replace an existing intensity Measurement.

run()

Runs the tool.

add_measurement(m)[source]

Add an intensity Measurement to internal dictionary used to compute the excitation diagram. This method can also be used to safely replace an existing intensity Measurement.

Parameters:

m – A Measurement instance containing intensity in units equivalent to \({\rm erg~cm^{-2}~s^{-1}~sr^{-1}}\)

remove_measurement(identifier)[source]

Delete a measurement from the internal dictionary used to compute column densities. Any associated column density will also be removed.

Parameters:

identifier (str) – the measurement identifier

Raises:

KeyError – if identifier not in existing Measurements

replace_measurement(m)[source]

Safely replace an existing intensity Measurement. Do not change a Measurement in place, use this method. Otherwise, the column densities will be inconsistent.

Parameters:

m – A Measurement instance containing intensity in units equivalent to \({\rm erg~cm^{-2}~s^{-1}~sr^{-1}}\)

class pdrtpy.tool.h2excitation.H2ExcitationFit(measurements=None, constantsfile='RoueffEtAl.tab')[source]

Bases: ExcitationFit

Tool for fitting temperatures, column densities, and ortho-to-para ratio(OPR) from an \(H_2\) excitation diagram. It takes as input a set of \(H_2\) rovibrational line observations with errors represented as Measurement.

Often, excitation diagrams show evidence of both “hot” and “cold” gas components, where the cold gas dominates the intensity in the low J transitions and the hot gas dominates in the high J transitions. Given data over several transitions, one can fit for \(T_{cold}, T_{hot}, N_{total} = N_{cold}+ N_{hot}\), and optionally OPR. One needs at least 5 points to fit the temperatures and column densities (slope and intercept \(\times 2\)), though one could compute (not fit) them with only 4 points. To additionally fit OPR, one should have 6 points (5 degrees of freedom).

Once the fit is done, ExcitationPlot can be used to view the results.

Parameters:

measurements (list of Measurement.) – Input \(H_2\) measurements to be fit.

Attributes:
cold_colden

The fitted cold gas total column density

fit_result

The result of the fitting procedure which includes fit statistics, variable values and uncertainties, and correlations between variables.

has_maps

Are the Measurements used map-based?.

has_scalar

Are the Measurements used scalars.

has_vectors

Are the Measurements used a Nx1 vector, e.g.

hot_colden

The fitted hot gas total column density

intensities

The stored intensities.

numcomponents

Number of temperature components in the fit

opr

The ortho-to-para ratio (OPR)

opr_fitted

Was the ortho-to-para ratio fitted?

tcold

The fitted cold gas excitation temperature

temperature

The fitted gas temperatures, returned in a dictionary with keys ‘hot’ and ‘cold’.

thot

The fitted hot gas excitation temperature

total_colden

The fitted total column density

Methods

add_measurement(m)

Add an intensity Measurement to internal dictionary used to compute the excitation diagram.

average_column_density([position, size, ...])

Compute the average column density over a spatial box.

colden(component)

The column density of hot or cold gas component, or total column density.

column_densities([norm, unit, line])

The computed upper state column densities of stored intensities

energies([line])

Upper state energies of stored intensities, in K.

gu(id, opr)

Get the upper state statistical weight $g_u$ for the given transition identifer, and, if the transition is odd-$J$, scale the result by the given ortho-to-para ratio.

intensity(colden)

Given an upper state column density \(N_u\), compute the intensity \(I\).

remove_measurement(identifier)

Delete a measurement from the internal dictionary used to compute column densities.

replace_measurement(m)

Safely replace an existing intensity Measurement.

run([position, size, fit_opr])

Fit the \(log N_u-E\) diagram with two excitation temperatures, a hot \(T_{ex}\) and a cold \(T_{ex}\).

upper_colden(intensity, unit)

Compute the column density in upper state \(N_u\), given an intensity \(I\) and assuming optically thin emission.

average_column_density(position=None, size=None, norm=True, unit=Unit("1 / cm2"), line=True, clip=<Quantity -1.e+40 1 / cm2>)[source]

Compute the average column density over a spatial box. The box is created using astropy.nddata.utils.Cutout2D.

Parameters:
  • position (tuple) – The position of the cutout array’s center with respect to the data array. The position can be specified either as a (x, y) tuple of pixel coordinates.

  • size (int, array_like`) – The size of the cutout array along each axis. If size is a scalar number or a scalar Quantity, then a square cutout of size will be created. If size has two elements, they should be in (nx,ny) order [this is the opposite of Cutout2D signature]. Scalar numbers in size are assumed to be in units of pixels. Default value of None means use all pixels (position is ignored)

  • norm (bool) – if True, normalize the column densities by the statistical weight of the upper state, \(g_u\). For ortho-$H_2$ $g_u = OPR times (2J+1)$, for para-$H_2$ $g_u=2J+1$. In LTE, $OPR = 3$.

  • unit (str or astropy.units.Unit) – The units in which to return the column density. Default: \({\rm cm}^{-2}\)

  • line (bool) – if True, the returned dictionary index is the Line name, otherwise it is the upper state \(J\) number.

  • clip (astropy.units.Quantity) – Column density value at which to clip pixels. Pixels with column densities below this value will not be used in the average. Default: a large negative number, which translates to no clipping.

Returns:

dictionary of column density Measurements, with keys as \(J\) number or Line name

Return type:

dict

property cold_colden

The fitted cold gas total column density

Return type:

Measurement

colden(component)[source]

The column density of hot or cold gas component, or total column density.

Parameters:

component (str) – ‘hot’, ‘cold’, or ‘total

Return type:

Measurement

column_densities(norm=False, unit=Unit('1 / cm2'), line=True)[source]

The computed upper state column densities of stored intensities

Parameters:
  • norm (bool) – if True, normalize the column densities by the statistical weight of the upper state, \(g_u\). Default: False

  • unit (str or astropy.units.Unit) – The units in which to return the column density. Default: \({\\rm }cm^{-2}\)

  • line (bool) – if True, the dictionary index is the Line name, otherwise it is the upper state \(J\) number. Default: False

Returns:

dictionary of column densities indexed by upper state \(J\) number or Line name. Default: False means return indexed by \(J\).

Return type:

dict

energies(line=True)[source]

Upper state energies of stored intensities, in K.

Parameters:

line (bool) – if True, the dictionary index is the Line name, otherwise it is the upper state \(J\) number. Default: False

Returns:

dictionary indexed by upper state \(J\) number or Line name. Default: False means return indexed by \(J\).

Return type:

dict

property fit_result

The result of the fitting procedure which includes fit statistics, variable values and uncertainties, and correlations between variables.

Return type:

lmfit.model.ModelResult

gu(id, opr)[source]

Get the upper state statistical weight $g_u$ for the given transition identifer, and, if the transition is odd-$J$, scale the result by the given ortho-to-para ratio. If the transition is even-$J$, the LTE value is returned.

Parameters:
  • id (str) – the measurement identifier

  • opr (float) –

Raises:

KeyError – if id not in existing Measurements

Return type:

float

property hot_colden

The fitted hot gas total column density

Return type:

Measurement

property intensities

The stored intensities. See add_measurement()

Return type:

list of Measurement

intensity(colden)[source]

Given an upper state column density \(N_u\), compute the intensity \(I\).

\[I = {A \Delta E~N_u \over 4\pi}\]

where \(A\) is the Einstein A coefficient and \(\Delta E\) is the energy of the transition.

Parameters:

colden (Measurement) – upper state column density

Returns:

optically thin intensity

Return type:

Measurement

property numcomponents

Number of temperature components in the fit

Return type:

int

property opr

The ortho-to-para ratio (OPR)

Returns:

The fitted OPR is it was determined in the fit, otherwise the canonical LTE OPR

Return type:

Measurement

property opr_fitted

Was the ortho-to-para ratio fitted?

Returns:

True if OPR was fitted, False if canonical LTE value was used

Return type:

bool

run(position=None, size=None, fit_opr=False, **kwargs)[source]

Fit the \(log N_u-E\) diagram with two excitation temperatures, a hot \(T_{ex}\) and a cold \(T_{ex}\).

If position and size are given, the data will be averaged over a spatial box before fitting. The box is created using astropy.nddata.utils.Cutout2D. If position or size is None, the data are averaged over all pixels. If the Measurements are single values, these arguments are ignored.

Parameters:
  • position (tuple) – The position of the cutout array’s center with respect to the data array. The position can be specified either as a (x, y) tuple of pixel coordinates.

  • size (int, array_like`) – The size of the cutout array along each axis in pixels. If size is a scalar number or a scalar Quantity, then a square cutout of size will be created. If size has two elements, they should be in (nx, ny) order [this is the opposite of Cutout2D signature]. Scalar numbers in size are assumed to be in units of pixels. Default value of None means use all pixels (position is ignored)

  • fit_opr (bool) – Whether to fit the ortho-to-para ratio or not. If True, the OPR will be varied to determine the best value. If False, the OPR is fixed at the canonical LTE value of 3.

property tcold

The fitted cold gas excitation temperature

Return type:

Measurement

property temperature

The fitted gas temperatures, returned in a dictionary with keys ‘hot’ and ‘cold’. :rtype: dict

property thot

The fitted hot gas excitation temperature

Return type:

Measurement

property total_colden

The fitted total column density

Return type:

Measurement

upper_colden(intensity, unit)[source]

Compute the column density in upper state \(N_u\), given an intensity \(I\) and assuming optically thin emission. Units of \(I\) need to be equivalent to \({\rm erg~cm^{-2}~s^{-1}~sr^{-1}}\).

\[ \begin{align}\begin{aligned}I &= {A \Delta E~N_u \over 4\pi}\\N_u &= 4\pi {I\over A\Delta E}\end{aligned}\end{align} \]

where \(A\) is the Einstein A coefficient and \(\Delta E\) is the energy of the transition.

Parameters:
  • intensity (Measurement) – A Measurement instance containing intensity in units equivalent to \({\rm erg~cm^{-2}~s^{-1}~sr^{-1}}\)

  • unit (str or astropy.units.Unit) – The units in which to return the column density. Default: \({\rm }cm^{-2}\)

Returns:

a Measurement of the column density.

Return type:

Measurement

LineRatioFit

LineRatioFit is a tool for determining photodissociation region external radiation field and hydrogen nucleus density (commonly given as \(G_0\) and \(n\)) from measured spectral line intensity ratios.

class pdrtpy.tool.lineratiofit.LineRatioFit(modelset=<pdrtpy.modelset.ModelSet object>, measurements=None)[source]

Bases: ToolBase

LineRatioFit is a tool to fit observations of intensity ratios to a set of PDR models. It takes as input a set of observations with errors represented as Measurement and ModelSet for the models to which the data will be fitted. The observations should be spectral line or continuum intensities. They can be spatial maps or single pixel values. They should have the same spatial resolution.

The models to be fit are stored as intensity ratios. The input observations will be use to create ratios that correspond to models. From there a minimization fit is done to determine the density and radiation field that best fit the data. At least 3 observations are needed in order to make at least 2 ratios. With fewer ratios, no fitting can be done. More ratios generally means better determined density and radiation field, assuming the data are consistent with each other.

Once the fit is done, LineRatioPlot can be used to view the results.

Parameters:
  • modelset (ModelSet) – The set of PDR models to use for fitting.

  • measurements (list or dict of Measurement. If dict, the keys should be the Measurement identifiers.) – Input measurements to be fit.

Attributes:
density

The computed hydrogen nucleus density value(s).

fit_result

The result of the fitting procedure which includes fit statistics, variable values and uncertainties, and correlations between variables.

has_maps

Are the Measurements used map-based?.

has_scalar

Are the Measurements used scalars.

has_vectors

Are the Measurements used a Nx1 vector, e.g.

measurementIDs

The stored measurement IDs, which are strings.

measurements

The stored measurements as a dictionary with Measurement IDs as keys

modelset

The underlying ModelSet

observed_ratios

The list of the observed line ratios that have been input so far.

radiation_field

The computed radiation field value(s).

ratiocount

The number of ratios that match models available in the current ModelSet given the current set of measurements

table

Construct the table of input Measurements, and if the fit has been run, the density, radiation field, and \(\chi^2\) values

Methods

add_measurement(m)

Add a Measurement to internal dictionary used to compute ratios.

chisq([min])

The computed chisquare value(s).

read_models([unit])

Given a list of measurement IDs, find and open the FITS files that have matching ratios and populate the _modelratios dictionary.

reduced_chisq([min])

The computed reduced chisquare value(s).

remove_measurement(id)

Delete a measurement from the internal dictionary used to compute ratios.

run(**kwargs)

Run the full computation using all the observations added.

write_chisq([chi, rchi, overwrite])

Write the chisq and reduced-chisq data to a file

add_measurement(m)[source]

Add a Measurement to internal dictionary used to compute ratios. This measurement may be intensity units (erg \({\rm s}^{-1}\) \({\rm cm}^{-2}\)) or integrated intensity (K km/s).

Parameters:

m (Measurement.) – a Measurement instance to be added to this tool

chisq(min=False)[source]

The computed chisquare value(s).

Parameters:

min (bool) – If True return the minimum reduced \(\chi^2\). In the case of map inputs this will be a spatial map of mininum \(\chi^2\). If False with map inputs the entire \(\chi^2\) hypercube is returned. If True with single pixel inputs, a single value is returned. If False with single pixel inputs, \(\chi^2\) as a function of density and radiation field is returned.

Return type:

Measurement

property density

The computed hydrogen nucleus density value(s).

Return type:

Measurement

property fit_result

The result of the fitting procedure which includes fit statistics, variable values and uncertainties, and correlations between variables. For each pixel, there will be one instance of lmfit.minimizer.MinimizerResult.

Return type:

FitMap

property measurementIDs

The stored measurement IDs, which are strings.

Return type:

dict_keys

property measurements

The stored measurements as a dictionary with Measurement IDs as keys

Return type:

dict of Measurement

property modelset

The underlying ModelSet

property observed_ratios

The list of the observed line ratios that have been input so far.

Return type:

list of str

property radiation_field

The computed radiation field value(s).

Return type:

Measurement

property ratiocount

The number of ratios that match models available in the current ModelSet given the current set of measurements

Return type:

int

read_models(unit=Unit(dimensionless))[source]

Given a list of measurement IDs, find and open the FITS files that have matching ratios and populate the _modelratios dictionary. Uses pdrtpy.measurement.Measurement as a storage mechanism.

param m:

list of measurement IDS (string)

type m:

list

param unit:

units of the data

type unit:

string or astropy.Unit

reduced_chisq(min=False)[source]

The computed reduced chisquare value(s).

Parameters:

min (bool) – If True return the minimum reduced \(\chi_\nu^2\). In the case of map inputs this will be a spatial map of mininum \(\chi_\nu^2\). If False with map inputs the entire \(\chi_\nu^2\) hypercube is returned. If True with single pixel inputs, a single value is returned. If False with single pixel inputs, \(\chi_\nu^2\) as a function of density and radiation field is returned.

Return type:

Measurement

remove_measurement(id)[source]

Delete a measurement from the internal dictionary used to compute ratios.

Parameters:

id (str) – the measurement identifier

Raises:

KeyError – if id not in existing Measurements

run(**kwargs)[source]

Run the full computation using all the observations added. This will check compatibility of input observations (e.g., beam parameters, coordinate types, axes lengths) and raise exceptions if the observations don’t match each other.

param mask:
Indicate how to mask image observations (Measurements) before computing the density

and radiation field. Possible values are:

[‘mad’, multiplier] - compute standard deviation using median absolute deviation (astropy.mad_std),

and mask out values between +/- multiplier*mad_std. Example: [‘mad’,1.0]

[‘data’, (low,high)] - mask based on data values, mask out data between low and high

[‘clip’, (low,high)] - mask based on data values, mask out data below low and above high

[‘error’, (low,high)] - mask based on uncertainty plane, mask out data where the corresponding error pixel value

is below low or above high

None - Don’t mask data. This is the default.

type mask:

list or None

param method:

the fitting method to be used. The default is ‘leastsq’, which is Levenberg-Marquardt least squares. For other options see https://lmfit-py.readthedocs.io/en/latest/fitting.html#fit-methods-table

type method:

str

param nan_policy:

Specifies action if fit returns NaN values. One of: * ’raise’ : a ValueError is raised [Default] * ’propagate’ : the values returned from userfcn are un-altered * ’omit’ : non-finite values are filtered

type nan_policy:

str

raises Exception:

if no models match the input observations, observations are not compatible, or on unrecognized parameters, or NaN encountered.

property table

Construct the table of input Measurements, and if the fit has been run, the density, radiation field, and \(\chi^2\) values

Return type:

astropy.table.Table

write_chisq(chi='chisq.fits', rchi='rchisq.fits', overwrite=True)[source]

Write the chisq and reduced-chisq data to a file

Parameters:
  • chi (str) – FITS file to write the chisq map to.

  • rchi (str) – FITS file to write the reduced chisq map to.

FitMap

When fitting either single pixels or spatial maps, the fit results are stored per pixel in an NDData object that contains ModelResult objects for H2ExcitationFit or MinimizerResult objects for LineRatioFit. The user can thus examine in detail the fit at any pixel.

class pdrtpy.tool.fitmap.FitMap(data, *args, **kwargs)[source]

Bases: NDData

Attributes:
data

~numpy.ndarray-like : The stored dataset.

mask

any type : Mask for the dataset, if any.

meta

dict-like : Additional meta information about the dataset.

name

The name of this FitMap

psf

Image representation of the PSF for the dataset.

uncertainty

any type : Uncertainty in the dataset, if any.

unit

~astropy.units.Unit : Unit for the dataset, if any.

wcs

any type : A world coordinate system (WCS) for the dataset, if any.

Methods

get_pixel(world_x, world_y)

Return the nearest pixel coordinates to the input world coordinates

get_pixel(world_x, world_y)[source]

Return the nearest pixel coordinates to the input world coordinates

Parameters:
  • world_x (float) – The horizontal world coordinate

  • world_y (float) – The vertical world coordinate

property name

The name of this FitMap :rtype: str

Utilities: Various constants and methods used by the Toolbox

pdrtpy.pdrutils.addkey(key, value, image)[source]

Add a (FITS) keyword,value pair to the image header

Parameters:
pdrtpy.pdrutils.check_nb()[source]
pdrtpy.pdrutils.check_units(input_unit, compare_to)[source]

Check if the input unit is equivalent to another.

Parameters:
Returns:

True if the input unit is equivalent to compare unit, False otherwise

pdrtpy.pdrutils.comment(value, image)[source]

Add a comment to an image header

Parameters:
pdrtpy.pdrutils.convert_if_necessary(image)[source]

Helper method to convert integrated intensity units in an image or Measurement from \({\rm K~km~s}^{-1}\) to \({\rm erg~s^{-1}~cm^{-2}~sr^{-1}}\). If a conversion is necessary, the convert_integrated_intensity() is called. If not, the image is returned unchanged.

Parameters:

image (astropy.io.fits.ImageHDU, astropy.nddata.CCDData, or Measurement.) – the image to convert. It must have a numpy.ndarray data member and astropy.units.Unit unit member or a header BUNIT keyword. It’s units must be \({\rm K~km~s}^{-1}\). It must also have a header RESTFREQ keyword.

Returns:

an image with converted values and units

pdrtpy.pdrutils.convert_integrated_intensity(image, wavelength=None)[source]

Convert integrated intensity from \({\rm K~km~s}^{-1}\) to \({\rm erg~s^{-1}~cm^{-2}~sr^{-1}}\), assuming \(B_\lambda d\lambda = 2kT/\lambda^3 dV\) where \(T dV\) is the integrated intensity in K km/s and \(\lambda\) is the wavelength. The derivation:

\[B_\lambda = 2 h c^2/\lambda^5 {1\over{exp[hc/\lambda k T] - 1}}\]

The integrated line is \(B_\lambda d\lambda\) and for \(hc/\lambda k T << 1\):

\[B_\lambda d\lambda = 2c^2/\lambda^5 \times (\lambda kT/hc)~d\lambda\]

The relationship between velocity and wavelength, \(dV = \lambda/c~d\lambda\), giving

\[B_\lambda d\lambda = 2\times10^5~kT/\lambda^3~dV,\]

with \(\lambda\) in cm, the factor \(10^5\) is to convert \(dV\) in \({\rm km~s}^{-1}\) to \({\rm cm~s}^{-1}\).

Parameters:
Returns:

an image with converted values and units

pdrtpy.pdrutils.dataminmax(image)[source]

Set the data maximum and minimum in image header

Parameters:

image (astropy.io.fits.ImageHDU, astropy.nddata.CCDData, or Measurement.) – The image which to add the key,val to.

pdrtpy.pdrutils.draine_unit = Unit("Draine")

The Draine radiation field unit

\({\rm 1~Draine = 2.72\times10^{-3}~erg~s^{-1}~cm^{-2}}\)

pdrtpy.pdrutils.dropaxis(w)[source]

Drop the first single dimension axis from a World Coordiante System. Returns the modified WCS if it had a single dimension axis or the original WCS if not.

Parameters:

w (astropy.wcs.WCS) – a WCS

Return type:

astropy.wcs.WCS

pdrtpy.pdrutils.firstkey(d)[source]

Return the “first” key in a dictionary

Parameters:

d (dict) – the dictionary

pdrtpy.pdrutils.fliplabel(label)[source]

Given a label that has a numerator and a denominator separated by a ‘/’, return the reciprocal label. For example, if the input label is ‘(x+y)/z’ return ‘z/(x+y)’. This method simply looks for the ‘/’ and swaps the substrings before and after it.

Parameters:

label (str) – the label to flip

Returns:

the reciprocal label

Return type:

str

Raises:

ValueError – if the input label has no ‘/’

pdrtpy.pdrutils.float_formatter(quantity, precision)[source]
pdrtpy.pdrutils.get_rad(key)[source]

Get radiation field symbol (LaTeX) given radiation field unit. If key is unrecognized, ‘FUV Radiation Field’ is returned.

Parameters:

key (str or astropy.units.Unit) – input field unit name, e.g. ‘Habing’, ‘Draine’ or an astropy.units.Unit

Returns:

LaTeX string for the radiation field symbol e.g., \(G_0\), \(\chi\)

Return type:

str

pdrtpy.pdrutils.get_table(filename, format='ipac', path=None)[source]

Return an astropy Table read from the input filename.

Parameters:
  • filename (str) – input filename, no path

  • format (str) – file format, Default: “ipac”

  • path (str) – path to filename relative to models directory. Default of None means look in “tables” directory

Return type:

astropy.table.Table

pdrtpy.pdrutils.get_testdata(filename)[source]

Get fully qualified pathname to FITS test data file.

Parameters:

filename (str) – input filename, no path

pdrtpy.pdrutils.get_xy_from_wcs(data, quantity=False, linear=False)[source]

Get the x,y axis vectors from the WCS of the input image.

Parameters:
Returns:

The axis values as arrays. Values are center of pixel.

Return type:

numpy.ndarray or astropy.units.Quantity

pdrtpy.pdrutils.habing_unit = Unit("Habing")

The Habing radiation field unit

\({\rm 1~Habing = 1.6\times 10^{-3}~erg~s^{-1}~cm^{-2}}\)

pdrtpy.pdrutils.has_single_axis(w)[source]

Check if the input WCS has any single dimension axes

Parameters:

w (astropy.wcs.WCS) – a WCS

Returns:

True if the input WCS has any single dimension axes, False otherwise

Return type:

bool

pdrtpy.pdrutils.history(value, image)[source]

Add a history record to an image header

Parameters:
pdrtpy.pdrutils.is_even(number)[source]

Check if number is even

Parameters:

number (float) – a number

Returns:

True if even, False otherwise

Return type:

bool

pdrtpy.pdrutils.is_image(image)[source]

Check if a Measurement is an image. The be an image it must have a header with axes keywords and a WCS to be considered an image. This is to distiguish Measurements that have a data array with more than one member from a true image. :param image: the image to check. It must have a numpy.ndarray data member and astropy.units.Unit unit member or a header BUNIT keyword. :type image: astropy.io.fits.ImageHDU, astropy.nddata.CCDData, or Measurement. :return: True if it is an image, False otherwise.

pdrtpy.pdrutils.is_odd(number)[source]

Check if number is odd

Parameters:

number (float) – a number

Returns:

True if odd, False otherwise

Return type:

bool

pdrtpy.pdrutils.is_rad(input_unit)[source]
pdrtpy.pdrutils.is_ratio(identifier)[source]

Is the identifier a ratio (as opposed to an intensity)

Return type:

bool

pdrtpy.pdrutils.mask_union(arrays)[source]

Return the union mask (logical OR) of the input masked arrays. This is useful when doing arithmetic on images that don’t have identical masks and you want the most restrictive mask.

Parameters:

arrays (numpy.ma.masked_array) – masked arrays to unionize

Return type:

mask

pdrtpy.pdrutils.mathis_unit = Unit("Mathis")

The Mathis radiation field unit

\({\rm 1~Mathis = 1.81\times10^{-3}~erg~s^{-1}~cm^{-2}}\)

pdrtpy.pdrutils.model_dir()[source]

Project model directory, including trailing slash

Return type:

str

pdrtpy.pdrutils.now()[source]
Returns:

a string representing the current date and time in ISO format

pdrtpy.pdrutils.rescale_axis_units(x, from_unit, from_ctype, to_unit, loglabel=True)[source]
pdrtpy.pdrutils.root_dir()[source]

Project root directory, including trailing slash

Return type:

str

pdrtpy.pdrutils.root_path()[source]

Project root directory as path

Return type:

Path

pdrtpy.pdrutils.setkey(key, value, image)[source]

Set the value of an existing keyword in the image header

Parameters:
pdrtpy.pdrutils.signature(image)[source]

Add AUTHOR and DATE keywords to the image header Author is ‘PDR Toolbox’, date as returned by now()

Parameters:

image (astropy.io.fits.ImageHDU, astropy.nddata.CCDData, or Measurement.) – The image which to add the key,val to.

pdrtpy.pdrutils.squeeze(image)[source]

Remove single-dimensional entries from image data and WCS.

Parameters:

image (astropy.nddata.CCDData, or Measurement.) – the image to convert. It must have a numpy.ndarray data member and astropy.units.Unit unit member.

Returns:

an image with single axes removed

Return type:

astropy.nddata.CCDData, or Measurement as input

pdrtpy.pdrutils.table_dir()[source]

Project ancillary tables directory, including trailing slash

Return type:

str

pdrtpy.pdrutils.testdata_dir()[source]

Project test data directory, including trailing slash

Return type:

str

pdrtpy.pdrutils.to(unit, image)[source]

Convert the image values to another unit. While generally this is intended for converting radiation field strength maps between Habing, Draine, cgs, etc, it will work for any image that has a unit member variable. So, e.g., it would work to convert density from \({\rm cm ^{-3}}\) to \({\rm m^{-3}}\). If the input image is a Measurement, its uncertainty will also be converted.

Parameters:
Returns:

an image with converted values and units

pdrtpy.pdrutils.toDraine(image)[source]

Convert a radiation field strength image to Draine units (chi).

\({\rm 1~Draine = 2.72\times10^{-3}~erg~s^{-1}~cm^{-2}}\)

between 6eV and 13.6eV (912-2066 \(\unicode{xC5}\)). See Weingartner and Draine 2001, ApJS, 134, 263, section 4.1

Parameters:

image (astropy.io.fits.ImageHDU, astropy.nddata.CCDData, or Measurement.) – the image to convert. It must have a numpy.ndarray data member and astropy.units.Unit unit member.

Returns:

an image with converted values and units

pdrtpy.pdrutils.toHabing(image)[source]

Convert a radiation field strength image to Habing units \((G_0)\).

\({\rm G_0 \equiv 1~Habing = 1.6\times10^{-3}~erg~s^{-1}~cm^{-2}}\)

between 6eV and 13.6eV (912-2066 \(\unicode{xC5}\)). See Weingartner and Draine 2001, ApJS, 134, 263, section 4.1

Parameters:

image (astropy.io.fits.ImageHDU, astropy.nddata.CCDData, or Measurement.) – the image to convert. It must have a numpy.ndarray data member and astropy.units.Unit unit member.

Returns:

an image with converted values and units

pdrtpy.pdrutils.toMathis(image)[source]

Convert a radiation field strength image to Mathis units

\({\rm 1~Mathis = 1.81\times10^{-3}~erg~s^{-1}~cm^{-2}}\)

between 6eV and 13.6eV (912-2066 \(\unicode{xC5}\)). See Weingartner and Draine 2001, ApJS, 134, 263, section 4.1

Parameters:

image (astropy.io.fits.ImageHDU, astropy.nddata.CCDData, or Measurement.) – the image to convert. It must have a numpy.ndarray data member and astropy.units.Unit unit member.

Returns:

an image with converted values and units

pdrtpy.pdrutils.tocgs(image)[source]

Convert a radiation field strength image to \({\rm erg~s^{-1}~cm^{-2}}\).

Parameters:

image (astropy.io.fits.ImageHDU, astropy.nddata.CCDData, or Measurement.) – the image to convert. It must have a numpy.ndarray data member and astropy.units.Unit unit member.

Returns:

an image with converted values and units

pdrtpy.pdrutils.warn(cls, msg)[source]

Issue a warning

Parameters:
  • cls (Class) – The calling Class

  • msg (str) – The warning message

Indices

Credits

pdrtpy is developed by Marc Pound and Mark Wolfire. This project is supported by NASA Astrophysics Data Analysis Program grant 80NSSC19K0573; from JWST-ERS program ID 1288 provided through grants from the STScI under NASA contract NAS5-03127; and from the SOFIA C+ Legacy Project through a grant from NASA through award #208378 issued by USRA.