"PDRTPY"("1") manual page

Name

pdrtpy - pdrtpy 2.0.7 \$1 \n[an-margin] level \n[rst2man-indent-level] level margin: \n[rst2man-indent\n[rst2man-indent-level]] - \n[rst2man-indent0] \n[rst2man-indent1] \n[rst2man-indent2]

Reliable astrophysics at everyday low, low prices! �®

----

Astrophysics Source Code Library 1102.022Project Status: Active - The project has reached a stable, usable state and is being actively developed.Python versionGNU GPL v3 License
pdrtpy is the new and improved version of the classic PhotoDissociation Region Toolbox, rewritten in Python with new capabilities and giving more flexibility to end users. (The Perl/CGI version of PDRT is deprecated and no longer supported).
The new PDR Toolbox will cover many more 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 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 u13d C chemistry along with the emitted line intensities for [u13d C II] and u13d CO.

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, 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

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:

PDRT_Example_Measurements.ipynb
PDRT_Example_ModelSets.ipynb
PDRT_Example_Model_Plotting.ipynb
PDRT_Example_Find_n_G0_Single_Pixel.ipynb

PDRT_Example_Make_n_G0_maps.ipynb

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. Your FITS images can be in intensity units, equivalent to {\rm erg~s^{-1}~cm^{-2}~sr^{-1}} or can be in K km/s. PDRT will do appropriate conversion as necessary when it uses your images (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)

Bases: astropy.nddata.ccddata.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

data (numpy.ndarray-like) – The actual data contained in this Measurement object. Note that the data will always be saved by reference, so you should make a copy of the data before passing it in if that’s the desired behavior.

uncertainty (astropy.nddata.StdDevUncertainty, astropy.nddata.VarianceUncertainty, astropy.nddata.InverseVariance or numpy.ndarray) – Uncertainties on the data. If the uncertainty is a numpy.ndarray, it assumed to be, and stored as, a astropy.nddata.StdDevUncertainty. Required.

unit (astropy.units.Unit or str) – The units of the data. Required.

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

title (str) – formatted string (e.g. LaTeX) describing this observation that can be used for plotting, e.g., r’$^{13}$CO(3-2)’

bmaj (astropy.units.Quantity) – [optional] beam major axis diameter. This will be converted to degrees for storage in FITS header

bmin (astropy.units.Quantity) – [optional] beam minor axis diameter. This will be converted to degrees for storage in FITS header

bpa (astropy.units.Quantity) – [optional] beam position angle.

This will be converted to degrees for storage in FITS header

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 axis 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 (flux/error)
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
flux: Return the underlying flux data array
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
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.
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.

identifier(id)

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

make_measurement(fluxfile, error, outfile[, …])

Create a FITS files with 2 HDUS, the first being the flux and the 2nd being the flux 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 flux in 1st HDU and error in 2nd HDU.

property SN: Return the signal to noise ratio (flux/error)
Return type

numpy.ndarray

add(other): Add this Measurement to another, propagating errors, units, and updating identifiers. Masks are logically or’d.
Parameters: other (Measurement) –

a Measurement to add

divide(other): Divide this Measurement by another, propagating errors, units, and updating identifiers. Masks are logically or’d.
Parameters: other (Measurement) –

a Measurement to divide

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

property flux: Return the underlying flux data array
Return type

numpy.ndarray

property id: Return the string ID of this measurement, e.g., CO_10
Return type

str

identifier(id): Set the string ID of this measurement, e.g., CO_10
Parameters: id (str) –

the identifier

property levels

static make_measurement(fluxfile, error, outfile, rms=None, masknan=True, overwrite=False)

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

Parameters

fluxfile (str) – The FITS file containing the flux data as a function of spatial coordinates

error (str) –

The errors on the flux data Possible values for error are:

a filename with the same shape as fluxfile 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 fluxfile

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 flux.
masknan (bool) – Whether to mask any pixel where the flux or the error is NaN. Default:true
overwrite (bool) – If True, overwrite the output file if it exists. Default: False.

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=aq10%aq,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=aqrmsaq,outfile="my_outfile.fits",units="K
km/s",overwrite=True)

multiply(other): Multiply this Measurement by another, propagating errors, units, and updating identifiers. Masks are logically or’d.
Parameters: other (Measurement) –

a Measurement to multiply

subtract(other): Subtract another Measurement from this one, propagating errors, units, and updating identifiers. Masks are logically or’d.
Parameters: other (Measurement) –

a Measurement to subtract

property title: A formatted title (e.g., LaTeX) that can be in plotting.
Return type

str or None

write(filename, **kwd)

Write this Measurement to a FITS file with flux 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

T} T{

Write

T} T{

Auto-identify

T} T{

pdrtpy.measurement.fits_measurement_reader(filename, hdu=0, unit=None, hdu_uncertainty=aqUNCERTaq, hdu_mask=aqMASKaq, hdu_flags=None, key_uncertainty_type=aqUTYPEaq, **kwd)

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 aqUNCERTaq.

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 aqMASKaq.

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. The current default are the Wolfire/Kaufman 2006 models, which have both constant density and constant thermal pressure versions, for metallicities z=1 and z=3 (limited spectral lines). Models are stored as ratios of intensities as a function of radiation field G_0 and hydrogen nucleus volume density n . We expect to update these soon with new physics and a wider range of lines and metallicities. Any PDR models can be used if they are stored in the correct FITS format. We are currently working with Marcus Rollig to import the Kosma-\tau models.

For example how to use ModelSets, see the notebook

PDRT_Example_ModelSets.ipynb

----

Manage pre-computed PDR models

class pdrtpy.modelset.ModelSet(name, z)

Bases: object

Class for computed PDR Model Sets. ModelSet will interface with a directory containing the model FITS files.

Parameters

name (str) – identifier.

z (float) – metallicity in solar units.

Raises: ValueError – If identifier or z not recognized/found.
Attributes
description: The description of this model
metallicity: The metallicity of this ModelSet
name: The name of this model
supported_intensities: Table of lines and continuum that are included in ratios models of this ModelSet.
supported_lines: Table of lines 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.
version: The version of this model
z

The metallicity of this ModelSet

Methods

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()

List the names and descriptions of available models (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.

property description: The description of this model
Return type

str

find_files(m, ext=aqfitsaq)

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): 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=aqfitsaq): 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=aqratioaq, ext=aqfitsaq)

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

static list(): List the names and descriptions of available models (not just this one)

property metallicity: The metallicity of this ModelSet
Return type

float

model_intensities(m): 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): 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 model
Return type

str

ratiocount(m): 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 included in ratios models of this ModelSet.
Return type

astropy.table.Table

property supported_lines: Table of lines 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 version: The version of this model
Return type

str

property z: The metallicity of this ModelSet
Return type

float

Utilities: Various constants and methods used by the Toolbox

pdrtpy.pdrutils.addkey(key, value, image)

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

Parameters

key (str) – The keyword to add to the header

value (any native type) – the value for the keyword

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

pdrtpy.pdrutils.check_units(input_unit, compare_to)

Check if the input unit is equivalent to another.

Parameters

input_unit (astropy.units.Unit, astropy.units.Quantity or str) – the unit to check.

compare_unit (astropy.units.Unit, astropy.units.Quantity or str) –

the unit to check against

Returns: True if the input unit is equivalent to compare unit,

False otherwise

pdrtpy.pdrutils.comment(value, image)

Add a comment to an image header

Parameters

value (str) – the value for the comment

image (astropy.io.fits.ImageHDU, astropy.nddata.CCDData, or Measurement.) –

The image which to add the comment to

pdrtpy.pdrutils.convert_if_necessary(image): 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)

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

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 header BUNIT keyword. It’s units must be K km/s

wavelength (astropy.units.Quantity) – the wavelength of the observation.

The default is to determine wavelength from the image header RESTFREQ keyword

Returns

an image with converted values and units

pdrtpy.pdrutils.dataminmax(image): 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): 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): Return the “first” key in a dictionary
Parameters: d (dict) –

the dictionary

pdrtpy.pdrutils.fliplabel(label): 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.get_rad(key): 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=aqipacaq, path=None)

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”

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
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 (FUV radiation field flux) 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): 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

savefig(fname, **kwargs)

Save the current figure to a file.

usetex(use)

Control whether plots use LaTeX formatting in axis labels and other text components.

property axis: The last axis that was drawn.
Return type

matplotlib.axes._subplots.AxesSubplot

property figure: The last figure that was drawn.
Return type

matplotlib.figure.Figure

savefig(fname, **kwargs): 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.

usetex(use): Control whether plots use LaTeX formatting in axis labels and other text components. This method sets matplotlib parameter rcParams[“text.usetex”] in the local pyplot instance.
Parameters: use (bool) –

whether to use LaTeX or not

H2ExcitationPlot

class pdrtpy.plot.h2excitationplot.H2ExcitationPlot(tool, **kwargs): Bases: pdrtpy.plot.plotbase.PlotBase
Class to plot various results from H2 Excitation diagram fitting.
Attributes
axis: The last axis that was drawn.
figure: The last figure that was drawn.

Methods

plot_diagram(x, y, xsize, ysize[, norm])

Plot the excitation diagram

savefig(fname, **kwargs)

Save the current figure to a file.

usetex(use)

Control whether plots use LaTeX formatting in axis labels and other text components.

plot_diagram(x, y, xsize, ysize, norm=True)

Plot the excitation diagram

Parameters

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

x (int) – bottom left corner x

y (int) – bottom left corner y

xsize (int) – box width, pixels

ysize (int) – box height, pixels

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

LineRatioPlot

class pdrtpy.plot.lineratioplot.LineRatioPlot(tool): Bases: pdrtpy.plot.plotbase.PlotBase
Class to plot various results from PDR Toolbox model fitting.
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.

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,G_0) 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.

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.

usetex(use)

Control whether plots use LaTeX formatting in axis labels and other text components.

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

confidence_intervals(**kwargs): 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): Plot the hydrogen nucleus volume density map that was computed by LineRatioFit tool. Default units: cm ^{-3}

modelintensity(id, **kwargs)

Plot one of the model intensities

Parameters

id (str) – the intensity identifier, such as CO_32ga.

**kwargs –

see class documentation above

Raises: KeyError –

if is id not in existing model intensities

modelratio(id, **kwargs)

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): 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): Overlay all the measured ratios and their errors on the (n,G_0) space.

Currently only works for single-pixel Measurements

radiation_field(**kwargs): Plot the radiation field map that was computed by LineRatioFit tool. Default units: Habing.

ratios_on_models(**kwargs): 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): Plot the reduced \chi^2 map that was computed by the LineRatioFit tool.

show_both(units=[aqHabingaq, aqcm^-3aq], **kwargs): 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): Bases: pdrtpy.plot.plotbase.PlotBase
Class to plot models and optionally Measurements. It does not require 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

intensity(identifier, **kwargs)

Plot a model ratio

overlay(measurements, **kwargs)

Overlay one or more single-pixel measurements in the model space ($n,G_0).

phasespace(identifiers[, dens_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

savefig(fname, **kwargs)

Save the current figure to a file.

usetex(use)

Control whether plots use LaTeX formatting in axis labels and other text components.

intensity(identifier, **kwargs): 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

overlay(measurements, **kwargs)

Overlay one or more single-pixel measurements in the model space ($n,G_0).

Parameters

measurements (list) – a list of one or more pdrtpy.measurement.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 contourf()) with alpha set to the shading value. Default value:

0.4

phasespace(identifiers, dens_clip=<Quantity [1.e+01, 1.e+07] 1 / cm3>, rad_clip=<Quantity [1.e+01, 1.e+06] Habing>, reciprocal=[False, False])

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”]

dens_clip (array like, must contain Quantity) – The range of model densities to show in the plot. Must be given as a range of astropy quanitities. Default: [10,1E7]*Unit(“cm-3”)

rad_clip (array like, must contain Quantity) – The range of model radiation field intensities to show in the plot. Must be given as a range of astropy quantities. Default: rad_clip=[10,1E6]*utils.habing_unit

reciprocal (list of 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]

plot(identifier, **kwargs): 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): 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.

----

ToolBase

The base class of all tools. Tools have a built-in plotter and a run() method both of which subclasses must define.

class pdrtpy.tool.toolbase.ToolBase

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.

Methods

run()

Runs the tool.

run()

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

H2Excitation

Tool for fitting temperatures in H_2 excitation diagrams. A two temperature model is assumed, T_{warm} and T_{cold}.

class pdrtpy.tool.h2excitation.H2Excitation(measurements=None)

Bases: pdrtpy.tool.toolbase.ToolBase

Tool for fitting temperatures to H_2 Excitation Diagrams

This tool is still under development

Attributes

intensities

The stored intensities.

Methods

add_measurement(m)

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

average_column_density(norm, x, y, xsize, …)

Compute the average column density over a spatial box.

colden(intensity)

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

column_densities([norm])

The computed upper state column densities of stored intensities

energies([line])

Upper state energies of stored intensities, in K.

fit_excitation(**kwargs)

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

replace_measurement(m)

Safely replace an existing intensity Measurement.

run()

Runs the tool.

excitation_diagram

plot_column_densities

plot_intensities

x_lin

add_measurement(m): 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}}

average_column_density(norm, x, y, xsize, ysize, line)

Compute the average column density over a spatial box.

Parameters

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

x (int) – bottom left corner x

y (int) – bottom left corner y

xsize (int) – box width, pixels

ysize (int) – box height, pixels

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

Returns: dictionary of column densities
Return type

dict

colden(intensity)

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}}.

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

N_u &= 4\pi {I\over A\Delta E}

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

Parameters

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

Returns

a Measurement of the column density.

Return type

Measurement

column_densities(norm=False): 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
Returns: dictionary of column densities indexed by Line name
Return type

dict

energies(line=False): 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.
Return type

dict

excitation_diagram(**kwargs)

fit_excitation(**kwargs): Fit the log N_u-E diagram with two excitation temperatures, a warm T_{ex} and a cold T_{ex}. A first pass guess is initially made using data partitioning and two linear fits.
Returns: The fit parameters as in scipy.optimize.curve_fit
Return type

list

property intensities: The stored intensities. See add_measurement()
Return type

list of Measurement

plot_column_densities(**kwargs)

plot_intensities(**kwargs)

replace_measurement(m): 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}}

x_lin(m1, n1, m2, n2)

LineRatioFit

Tool for determining photodissociation region external radiation field and density (commonly known as G_0 and n) from measured spectral line intensity ratios.

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

Bases: pdrtpy.tool.toolbase.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 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.

Attributes

density

The computed hydrogen nucleus density value(s).

has_maps

Are the Measurements used map-based or pixel-based?.

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

Methods

add_measurement(m)

Add a Measurement to internal dictionary used to compute ratios.

chisq([min])

The computed chisquare value(s).

compute_density_radiation_field()

Compute the best-fit density and radiation field spatial maps by searching for the minimum chi-squared at each spatial pixel.

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([mask])

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): 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): 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

compute_density_radiation_field(): Compute the best-fit density and radiation field spatial maps by searching for the minimum chi-squared at each spatial pixel.

property density: The computed hydrogen nucleus density value(s).
Return type

Measurement

property has_maps: Are the Measurements used map-based or pixel-based?.
Returns: True, if the observational inputs are spatial maps, False if they were single-pixel Measurements
Return type

bool

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)): 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): 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): 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(mask=[aqmadaq, 1.0]): 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

[‘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

type mask: list or None
raises Exception: if no models match the input observations, observations are not compatible,

or on unrecognized parameters

write_chisq(chi=aqchisq.fitsaq, rchi=aqrchisq.fitsaq, overwrite=True)

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.

for Developers

If you plan to tinker with the code, you may want to try this installation method.

git clone https://github.com/mpound/pdrtpy
sudo apt-get install python3-venv
python -m venv ~/venv
source ~/venv/bin/activate
cd pdrtpy
pip install -r requirements.txt

Indices

genindex
modindex

Name
What is a Pdr?
Getting Started

Installation
Requirements
Install the package
Example Notebooks

Module Descriptions and Apis

Measurements: How you put observations to the Toolbox
ModelSets: The interface to models in the Toolbox
Utilities: Various constants and methods used by the Toolbox
Plotting Tools: Display models and data
Plot Keywords
PlotBase
H2ExcitationPlot
LineRatioPlot
ModelPlot
Analysis Tools: Fit models to data
ToolBase
H2Excitation
LineRatioFit

for Developers
Indices
Credits
Author
Copyright