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.
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)[source]¶
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 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 acold
\(T_{ex}\).Safely replace an existing intensity Measurement.
run
()Runs the tool.
excitation_diagram
plot_column_densities
plot_intensities
x_lin
- 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}}\)
- average_column_density(norm, x, y, xsize, ysize, line)[source]¶
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
- colden(intensity)[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
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
- column_densities(norm=False)[source]¶
The computed upper state column densities of stored intensities
- fit_excitation(**kwargs)[source]¶
Fit the \(log N_u-E\) diagram with two excitation temperatures, a
warm
\(T_{ex}\) and acold
\(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
- property intensities¶
The stored intensities. See
add_measurement()
- Return type
list of
Measurement
- 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}}\)
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)[source]¶
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 keysmodelset
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 a Measurement to internal dictionary used to compute ratios.
chisq
([min])The computed chisquare value(s).
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).
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)[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
- compute_density_radiation_field()[source]¶
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
- 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
- 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 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
- property ratiocount¶
The number of ratios that match models available in the current
ModelSet
given the current set of measurements- Return type
- 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
- remove_measurement(id)[source]¶
Delete a measurement from the internal dictionary used to compute ratios.
- run(mask=['mad', 1.0])[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
[‘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