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)[source]

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
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)[source]

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)[source]

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)[source]

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, unit='adu')[source]

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.CCDData` 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.

  • 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) – a Measurement 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) – 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)[source]

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

Write

Auto-identify

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

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