Quickstart with XSNAP

Prerequisites

You’ll need to know a bit of Python, NumPy, Matplotlib, Astropy, and XSPEC (or PyXspec). Please refer to each of their own documentation for details.

Additionally, a little understanding on basic astronomical photometry, spectroscopy, and stellar physics would help.

Installation

Details are available on Installation but a quick command below can be used.

pip install xsnap

Important

You must remember to first install HEASoft + XSPEC and PyXspec. Details of the required dependencies are here.

Before doing anything with XSNAP, make sure to initialize HEASoft:

If you have Bourne shell variants (bash/sh/zsh):

>>> export HEADAS=/path/to/heasoft/platform/
>>> source $HEADAS/headas-init.sh

Or if you have C-shell variants (tcsh/csh):

>>> setenv HEADAS /path/to/heasoft/platform/
>>> source $HEADAS/headas-init.csh

Hint

It is best to use a new virtual environment, particularly conda environment, for XSNAP.

In addition, it may be a good idea to make an alias for HEASoft in your ~/.bashrc or ~/.zshrc by doing:

1. Open the ~/.bashrc or ~/.zshrc file:

   >>> nano ~/.bashrc       # or ~/.zshrc
   

2. Add these two lines:

   >>> export HEADAS=/path/to/heasoft/platform/
   >>> alias heainit='. $HEADAS/headas-init.sh'
   

3. Reload the shell config:

   >>> source ~/.bashrc     # or ~/.zshrc
   

4. Now, you can initialize HEASoft by:

   >>> heainit
   

Please adjust accordingly if you have C-shell variants (tcsh/csh).

Quick CLI Usage

Once installed, XSNAP command-line interface (CLI) scripts are easy and ready-to-use.

Initialize HEASoft in the terminal first:

>>> export HEADAS=/path/to/heasoft/software
>>> source $HEADAS/headas-init.sh

A quick hint is also available above.

Then, you can immediately use in terminal, e.g.,

>>> extract-xmm ./0882480901/odf/ sou_physical.reg bkg_physical.reg

If you want to run it in an .ipynb file:

>>> ! extract-xmm ./0882480901/odf/ sou_physical.reg bkg_physical.reg

You can also run it in vanilla Python:

>>> import subprocess
>>> cmd = ["make-region",
...        "/path/to/event/file",
...        "169.59", "-32.83"]
>>> subprocess.run(cmd)

Note

Make sure to have install the dependencies of each script. Details of each script is also available here. Some examples are also made here.

Quick API Usage

Importing XSNAP

Fundamentally, you can import XSNAP by:

import xsnap

However, it can be tricky as it heavily relies on PyXspec, which is not available in PyPI.

In this page, there will be two walkthroughs on importing XSNAP, i.e. through Jupyter and Visual Studio Code (VS Code)

Importing through Jupyter

1. Activate the virtual environment that you installed XSNAP in:

   >>> conda activate xsnap-venv
   

2. Initialize HEASoft in the environment:

   >>> export HEADAS=/path/to/heasoft/software
   >>> source $HEADAS/headas-init.sh
   

3. Run Jupyter lab or Jupyter Notebook

   >>> jupyter lab
   

   or

   >>> jupyter notebook
   

4. Import XSNAP through the .ipynb file:

   >>> import xsnap
   

Importing through VS Code

1. Make sure to add VS Code to the PATH environment variable, such that you can open VS Code using code in terminal. Details can be found in the VS Code setup page for Linux, macOS, and Windows.

2. Activate the virtual environment that you installed XSNAP in:

   >>> conda activate xsnap-venv
   

3. Initialize HEASoft in the environment:

   >>> export HEADAS=/path/to/heasoft/software
   >>> source $HEADAS/headas-init.sh
   

4. Open VS Code through the exact same terminal:

   >>> code .
   

5. Import XSNAP through the .ipynb file:

   >>> import xsnap
   

Spectrum Fitting

Below is a minimal example to fit a spectrum.

Hint

The SpectrumFit assumes you have a grouped spectrum file myspectrum.pi or at least have the same name as the background and response file

import xsnap

# 1) Create the fit object
spec = xsnap.SpectrumFit(abund="aspl")

# 2) Load your PHA spectrum file
spec.load_data("myspectrum.pi")

# 3) Define a model (e.g. absorbed power law)
spec.set_model("tbabs*pow",
                TBabs_nH="0.05 -1",
                powerlaw_PhoIndex=2)
# The model has 3 parameters here:
# TBabs.nH, powerlaw.PhoIndex, and powerlaw.norm

# 4) Fit!
spec.fit(nIterations=500)

# 5) Plot (if you want to see the plot in svg for example)
spec.set_plot("ldata", device="/svg")

# 6) Getting best-fit parameters for powerlaw PhoIndex and norm.
# This means that we're fitting for 1σ for parameter 2 and 3 (PhoIndex and norm)
spec.get_params("1.0 2 3")

# 7) Getting observation times and count rates
spec.get_counts()
spec.get_time()

# 8) Getting flux and luminosity
df_flux = spec.get_fluxes()
spec.get_lumin(df_flux['unabsorbed'], redshift=0.1)

Handling a Collection of Spectrum

Below is a minimal example to plot light curves and parameter evolutions from multiple spectra.

import xsnap

# 1) Create manager object
manager = xsnap.SpectrumManager()

# 2) Analyze a spectrum
spec1 = xsnap.SpectrumFit(abund="aspl")
spec1.load_data("spectrum1.pha")
spec1.set_model("tbabs*pow",
                TBabs_nH="0.05 -1",
                powerlaw_PhoIndex=2)
spec1.fit()
spec1.get_params("1.0 2 3")
spec1.get_counts()
spec1.get_time()
df_flux = spec1.get_fluxes()
spec1.get_lumin(df_flux['unabsorbed'], redshift=0.1)

# 3) Load the spectrum to manager
# (Optionally) you can also specify
# the instrument by manager.load(spec1, instrument="XMM")
# Otherwise, the manager will try to find the instrument
# from the spectrum file
manager.load(spec1)

# 4) Repeat analysis of other spectrums
# Due to the nature of PyXspec, it is
# best to load spectrum each time after analysis

spec2 = xsnap.SpectrumFit(abund="aspl")
# clear existing data before loading a new one by doing:
spec2.load_data("spectrum2.pha", clear=True)
# you can also clear by importing xspec and run:
# xspec.AllData.clear()
spec2.set_model("tbabs*pow",
                TBabs_nH="0.05 -1",
                powerlaw_PhoIndex=2)
spec2.fit()
spec2.get_params("1.0 2 3")
spec2.get_counts()
spec2.get_time()
df_flux = spec2.get_fluxes()
spec2.get_lumin(df_flux['unabsorbed'], redshift=0.1)

# 5) Load the spectrum again to manager
manager.load(spec2)

# 6) Plot flux and luminosity light curves
manager.plot_flux()
manager.plot_lumin()

# 7) Plot parameter evolution
manager.plot_params()

Detecting Sources

Below is a minimal example to detect a source given an event file and the source coordinates.

import glob
from xsnap import SourceDetection

# 1) Define the paths to the event files and exposure image files
# (exposure image are optional but helpful and recommended to use)
# Here, I use glob to list all the .evt and .img file inside of a directory

evt_paths = sorted(glob.glob("./parent/dir/of/event/*.evt"))
img_paths = sorted(glob.glob("./parent/dir/of/image/*.img"))

# Keep in mind that the amount of event path and image must be the same
# i.e. len(evt_paths) == len(img_paths)
# This is because they will be paired up when used for detecting, i.e.
# evt_paths[0] with img_paths[0]

# 2) Create the detection object and load the files
detect = SourceDetection()
detect.load(evt_paths, img_paths)

# 3) Define the source RA and dec
# Make sure it's in decimal degrees
RA = 169.59
dec = -32.83

# 4) Detect the source in all event files!
detect.detect_all(RA, dec)

More Examples

More examples are available in the Examples page (more notebooks will be made soon!).