Skip to content

Commit

Permalink
add tutorial template. Update the notebooks to follow the template
Browse files Browse the repository at this point in the history
  • Loading branch information
@zoghbi-a
zoghbi-a committed Feb 28, 2024
1 parent f6f254c commit 2b1ac77
Show file tree
Hide file tree
Showing 18 changed files with 1,857 additions and 1,331 deletions.
87 changes: 87 additions & 0 deletions _files/template.md
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,87 @@
---
jupyter:
jupytext:
text_representation:
extension: .md
format_name: markdown
format_version: '1.3'
jupytext_version: 1.16.0
kernelspec:
display_name: Python 3 (ipykernel)
language: python
name: python3
---

# Descriptive Title of the Content of the Notebook
---
Header Section: include the following information.

- **Description:** A template on how to write a notebook tutorial on sciserver.
- **Level:** Beginner | Intermediate | Advanced.
- **Data:** Descirbe what data, if any will be use. If None, write: NA
- **Requirements:** Describe what is needed to run the notebooks. For example: "Run in the (heasof) conda environment on Sciserver". Or "python packages: [`heasoftpy`, `astropy`, `numpy`]".
- **Credit:** Who wrote the notebebook and when.
- **Support:** How to get help.
- **Last verified to run:** (00/00/000) When was this last tested.

---


## 1. Introduction
Describe the content. It can contain 0plain text, bullets, and/or images as needed.
Use `Markdown` when writing.

The following are suggested subsections. Not all are needed:
- Motivation / Science background.
- Learning goals.
- Details about the requirements, and on running the notebook outside Sciserver.
- Type of outcome or end product.

You may want to include the following section on how to run the notebook outise sciserver.
<div style='color: #333; background: #ffffdf; padding:20px; border: 4px solid #fadbac'>
<b>Running On Sciserver:</b><br>
When running this notebook inside Sciserver, make sure the HEASARC data drive is mounted when initializing the Sciserver compute container. <a href='https://heasarc.gsfc.nasa.gov/docs/sciserver/'>See details here</a>.
<br><br>
<b>Running Outside Sciserver:</b><br>
This notebook runs in the heasoftpy conda environment on Sciserver.
If running outside Sciserver, some changes will be needed, including:<br>
&bull; Make sure heasoftpy and heasoft are correctly installed (<a herf='https://heasarc.gsfc.nasa.gov/docs/software/lheasoft/'>Download and Install heasoft</a>).<br>
&bull; Unlike on Sciserver, where the data is available locally, you will need to download the data to your machine.<br>
</div>

The following gives an example sections.


## 2. Imports

```python
# add imports here
import heasoftpy as hsp
```

## 3.0 Define Input if needed
This section will include things like:
- obsIDs
- Plot settings
- Work directory
- Detector settings
- etc

```python
obsid = '000000000'
```

## 4. Data Access
How is the data used here can be found, and accessed (e.g. copied, downloaded etc.)


## 5. Data Processing or Extraction
Include this section if needed (choose a relevant title)


## 6. Data Products
If the notebook produces some data products, describe them here, and how the user will be using them.

```python

```
115 changes: 70 additions & 45 deletions ixpe_intro.md → analysis-ixpe-example.md
View file
Original file line number Diff line number Diff line change
Expand Up @@ -5,33 +5,70 @@ jupyter:
extension: .md
format_name: markdown
format_version: '1.3'
jupytext_version: 1.15.2
jupytext_version: 1.16.0
kernelspec:
display_name: (heasoft)
language: python
name: heasoft
---

# Getting Started with IXPE data
# Getting Started with IXPE Data
<hr style="border: 2px solid #fadbac" />

This notebook is a tutorial on accessing IXPE data on Sciserver and getting started with analysing the data. You will learn to download the data, extract the source and background regions and perform spectro-polarimetric fits.
- **Description:** An example on analysing IXPE data on Sciserver using heasoftpy.
- **Level:** Intermediate.
- **Data:** IXPE observation of blazar **Mrk 501** (ObsID 01004701).
- **Requirements:** `heasoftpy`, `xspec`, `matplotlib`
- **Credit:** Kavitha Arur (Jun 2023).
- **Support:** Contact the [IXPE Guest Observer Facility (GOF)](https://heasarc.gsfc.nasa.gov/cgi-bin/Feedback?selected=ixpe) or the [HEASARC helpdesk](https://heasarc.gsfc.nasa.gov/cgi-bin/Feedback).
- **Last verified to run:** 01/26/2024.

<hr style="border: 2px solid #fadbac" />


## 1. Introduction

This notebook is a tutorial on accessing IXPE data on Sciserver and getting started with analysing them. You will learn to download the data, extract the source and background regions and perform spectro-polarimetric fits.

It also highly recommended that new users read the IXPE Quick Start Guide ([linked here](https://heasarc.gsfc.nasa.gov/docs/ixpe/analysis/IXPE_quickstart.pdf)) and the recommended practices for statistical treatment of IXPE results [here](https://heasarcdev.gsfc.nasa.gov/docs/ixpe/analysis/IXPE_Stats-Advice.pdf).

<div style='color: #333; background: #ffffdf; padding:20px; border: 4px solid #fadbac'>
<b>Running On Sciserver:</b><br>
When running this notebook inside Sciserver, make sure the HEASARC data drive is mounted when initializing the Sciserver compute container. <a href='https://heasarc.gsfc.nasa.gov/docs/sciserver/'>See details here</a>.
<br>
Also, this notebook requires <code>heasoftpy</code>, which is available in the (heasoft) conda environment. You should see (heasoft) at the top right of the notebook. If not, click there and select it.

<b>Running Outside Sciserver:</b><br>
If running outside Sciserver, some changes will be needed, including:<br>
&bull; Make sure heasoftpy and heasoft are installed (<a herf='https://heasarc.gsfc.nasa.gov/docs/software/lheasoft/'>Download and Install heasoft</a>).<br>
&bull; Unlike on Sciserver, where the data is available locally, you will need to download the data to your machine.<br>
</div>


## 2. Module Imports
We need the following python modules:

<div style='color: #333; background: #ffffdf; padding:20px; border: 4px solid #fadbac'>
In this example, reprocessing the data is not required. Instead the level 2 data products are sufficient. If you need to reprocess the data, the IXPE tools are available with <code>from heasoftpy import ixpe</code>.
</div>

```python
import glob
import matplotlib.pyplot as plt
import heasoftpy as hsp
import xspec
```

## Finding and exploring data
## 3. Finding the Data

All the heasarc data is mounted into the compute under /FTP/, so once we have the path to the data, we can directly access it without the need to download it.
On Sciserver, all the HEASARC data ire mounted locally under `/FTP/`, so once we have the path to the data, we can directly access it without the need to download it.

For our exploratory data analysis, we will use an observation of the blazar Mrk 501 (ObsID 01004701).
For our exploratory data analysis, we will use an observation of the blazar **Mrk 501** (ObsID 01004701).

(For more information on how to locate datasets of interest, see the [data access notebook](data_access.ipynb).)
You can also see the [Getting Started](getting-started.md), [Data Access](data-access.md) and [Finding and Downloading Data](data-find-download.md) tutorials for examples on how to find data.

```python
paths_txt = """
/FTP/ixpe/data/obs/01/01004701
"""
paths = paths_txt.split('\n')[1:-1]
data_path = "/FTP/ixpe/data/obs/01/01004701"
```

Check the contents of this folder
Expand All @@ -44,26 +81,22 @@ It should contain the standard IXPE data files, which include:
For a complete description of data formats of the level 1, level 2 and calibration data products, see the support documentation on the [IXPE Website](https://heasarc.gsfc.nasa.gov/docs/ixpe/analysis/#supportdoc)

```python
import glob
glob.glob(f'{paths[0]}/*')
glob.glob(f'{data_path}/*')
```

---
## Analyzing The Data
## 4. Exploring The Data
To Analyze the data within the notebook, we use `heasoftpy`.

In the folder for each observation, check for a README file. This file is included with a description of known issues (if any) with the processing for that observation.
In the folder for each observation, check for a `README` file. This file is included with a description of known issues (if any) with the processing for that observation.

In this *IXPE* example, it is not necessary to reprocess the data. Instead the level 2 data products can be analysed directly.

```python
import heasoftpy as hsp

# set some input
indir = paths[0]
indir = data_path
obsid = indir.split('/')[-1]

filelist = glob.glob(f'{paths[0]}/event_l2/*')
filelist = glob.glob(f'{indir}/event_l2/*')
filelist
```

Expand All @@ -74,16 +107,15 @@ det1_fits = filelist[0]
det2_fits = filelist[1]
det3_fits = filelist[2]

#print the file structure for event 1 dectector file
#print the file structure for event 1 detector file
out = hsp.fstruct(infile=det1_fits).stdout
print(out)

```

---
## Extracting the spectro polarimetric data
## 5. Extracting the spectro polarimetric data

### Defining the source and background regions
### 5.1 Defining the Source and Background Regions

To obtain the source and background spectra from the Level 2 files, we need to define a source region and background region for the extraction. This can also be done using `ds9`.

Expand All @@ -101,7 +133,7 @@ f.write('annulus(16:53:51.766,+39:45:44.41,132.000",252.000")')
f.close()
```

### Running the extractor tools
### 5.2 Running the extractor tools

The `extractor` tool from FTOOLS, can now be used to extract I,Q and U spectra from IXPE Level 2
event lists as shown below.
Expand Down Expand Up @@ -165,21 +197,20 @@ if out.returncode != 0:
raise Exception('extractor for det3 failed!')
```

---
### Obtaining the response files
### 5.3 Obtaining the Response Files

For the I spectra, you will need to include the RMF (Response Matrix File), and
the ARF (Ancillary Response File).

For the Q and U spectra, you will need to include the RMF and MRF (Modulation Response File). The MRF is defined as defined by the product of the energy-dependent modulation factor, $\mu$(E) and the ARF.
For the Q and U spectra, you will need to include the RMF and MRF (Modulation Response File). The MRF is defined by the product of the energy-dependent modulation factor, $\mu$(E) and the ARF.

The location of the calibration files can be obtained through the `hsp.quzcif` tool. Type in `hsp.quzcif?` to get more information on this function.

Note that the output of the `hsp.quzcif` gives the path to more than one file. This is because there are 3 sets of response files, corresponding to the different weighting schemes.

For the 'NEFF' weighting, use 'alpha07_02'.
For the 'SIMPLE' weighting, use 'alpha075simple_02'.
For the 'UNWEIGHTED' version, use '20170101_02'.
- For the 'NEFF' weighting, use 'alpha07_02'.
- For the 'SIMPLE' weighting, use 'alpha075simple_02'.
- For the 'UNWEIGHTED' version, use '20170101_02'.

```python
# hsp.quzcif?
Expand Down Expand Up @@ -233,11 +264,10 @@ res = hsp.quzcif(mission='ixpe', instrument='gpd',detector='DU3',
mrf3 = [x.split()[0] for x in res.output if 'alpha075_02' in x][0]
```

---
### Load data into PyXSPEC and start fitting

### 5.4 Load data into PyXSPEC and start fitting

```python
import xspec

rmf_list = [rmf1,rmf2,rmf3]
mrf_list = [mrf1,mrf2,mrf3]
Expand All @@ -251,6 +281,7 @@ for (du, rmf_file, mrf_file, arf_file) in zip(du_list, rmf_list, mrf_list, arf_l

#Load the I data
xspec.AllData("%i:%i ixpe_det%i_src_I.pha"%(du, du+x, du))
xspec.AllData(f"{du}:{du+x} ixpe_det{du}_src_I.pha")
s = xspec.AllData(du+x)

# #Load response and background files
Expand Down Expand Up @@ -318,16 +349,11 @@ xspec.AllModels.show()
xspec.Fit.perform()
```

---
### Plotting the results
### 5.5 Plotting the results

This is done through matplotlib.
This is done through `matplotlib`.

```python
import matplotlib.pyplot as plt


%matplotlib inline
xspec.Plot.area=True
xspec.Plot.xAxis='keV'
xspec.Plot('lda')
Expand Down Expand Up @@ -366,8 +392,7 @@ ax.set_xlabel('Energy (keV)')
ax.set_ylabel(r'Polangle')
```

---
## Interpreting the results from XSPEC
## 6. Interpreting the results from XSPEC

There are two parameters of interest in our example. These given by the polarization fraction, A,
and polarization angle, $\psi$. The XSPEC error (or uncertainty) command can be used
Expand Down Expand Up @@ -414,7 +439,7 @@ plt.xlabel('A')
plt.errorbar(m1.polconst.A.values[0],m1.polconst.psi.values[0],fmt='+')
```

### Determining the flux and calculating MDP
### 6.1 Determining the flux and calculating MDP


Note that the detection is deemed "highly probable" (confidence C > 99.9%) as
Expand All @@ -439,7 +464,7 @@ mean of exposure time of 97243 s gives an MDP99 of 5.70% meaning that, for an un
This is consistent with the highly probable detection deduced here of a polarization fraction of 7.45$\pm$1.8%.


## Additional Resources
## 7. Additional Resources

Visit the IXPE [GOF Website](https://heasarcdev.gsfc.nasa.gov/docs/ixpe/analysis/) and the IXPE [Project Website at MSFC](https://ixpe.msfc.nasa.gov/for_scientists/index.html) for more resources.

Expand Down
Loading

0 comments on commit 2b1ac77

Please sign in to comment.