X ray data analysis

Notes on X-ray data analysis tasks relating to spectroscopy and timing for different missions.

Specific observatories

XMM-Newton

• One starting point could be the Streamlined ABC Guide provided by HEASARC.
• Also see ESA's XMM-SAS user guide.

NuSTAR

• See the HEASARC NuSTAR page which includes a link to the analysis guide page which has a quick start guide.

Obtaining data from the archive

• You can search for data and download it using the "HEASARC Browse" interface

Data preparation

You can group the data to a minimum number of counts per bin, for example, using grppha (see the help file for more information about FTOOLS in general). Here's an example of grouping data to a minimum of 20 counts per bin and setting some of the header keywords.

grppha
input_file.fits
output_file_grp.fits
chkey BACKFILE background_pi.fits
chkey RESPFILE rmf.fits
chkey ANCRFILE arf.fits
group min 20
exit

When you read the output_grp.fits file into whatever fitting program you are using it should read in the grouping flags, as well as the background, RMF, and ARF files in the header.

Fitting algorithms

A theoretical spectral model can be compared to X-ray data by "folding" the model through the instrument response. This needs to be done because the observed spectram cannot be uniquely inverted (for example, different incoming photon energies might result in the same detector response because it is probabalistic). The folded model can be comapred to the data to generate a statistic such as a $\chi^2$ value where a smaller value of the statistic represents a better fit (although the data can be "overfit" if the model is implausible good given the noise in the data). The job of a fitting algorithm is to adjust the model's knobs (e.g., photon index and normalisation for a powerlaw) to find a range of parameters that provide a statistically acceptable description of the spectrum. There are many ways of doing this. The simplest models are gradient descent models that go "downhill" in the statistic space unitl they arrive at a minimum, but this might not be a global minimum. More sophisticated algorithms are better at not getting trapped in loacl minima and more rigorously exploring a multi-dimensional parameters space. To make things even more difficult, model evaluations can be quite slow, so efficiency is important. Bayesian and Markov-Chain Monte-Carlo methods are often exmployed.

MCMC

Here are some Markov-Chain Monte-Carlo resources.

• A paper by Hogg & Foreman-Mackey (2017), "Data analysis recipes: Using Markov Chain Monte Carlo"
• An emcee tutorial on Fitting a model to data
• An nDspec example fit
• A SpectralFitting.jl MCMC walkthrough

XSPEC

• Many astronomers use XSPEC for fitting spectra. Here's an example walkthrough for loading some XMM data and fitting it with a simple absorbed power-law model.
• To run XSPEC on the group's computers see our local documentation (UOB access only) for information on how to set up the paths and initialise HEASOFT.

Here's a walk through assuming your data set is called pn_spec_grp.fits and its header already points to the associated background, RMF, and ARF. The specific data set here is for the AGN MCG--6-30-15. The ignore commands exclude fitting the very low and high energy ranges where the instrument calibration is not good. cps /xs sets the plot device, and setp e sets the x-scale to be energy. mo tbabs * pow sets up the absorbed powerlaw model. You can use, e.g., fit 100 for it to try 100 iterations before asking you to continue. pl ld de plots the log of the data with the $\Delta \chi$ residuals.

xspec
data pn_spec_grp.fits
ignore 0.-0.5
ignore 10.-**
cpd /xs
setp e
mo tbabs * pow

Hit return a few times to choose default values.

fit

When XSPEC says "Number of trials exceeded: continue fitting?" hit y.

pl ld de

You should get a plot like this

Not a great fit, but this shows that your data extraction and XSPEC are working.

SpectralFitting.jl

Fergus Baker has been leading the development of SpectralFitting.jl that is a Julia-based spectral fitting package. See the documentation.

nDspec

Another options is nDspec which is a Python based fitting package.