Tutorial 1: Single File Data Reduction

Welcome to the WIRC+POL data reduction pipeline (DRP) tutorial. The notebook will give you a basic run through of how to reduce WIRC+POL data and give you some details about the setup of the pipeline.

Pipeline Description

The pipeline has been designed so that once read-in, each data file will be held in an object called a wirc_data object. Each wircpol_data object will contain, the raw (or calibrated) 2D data itself, relevant information about the data and a number of functions that will calibrate and process the data. It will also contain a list of wircpol_source objects, where each wircpol_source will correspond to a source in the field. The source object will contain thumbnails of the spectral traces, as well as the extracted spectrum and polarized spectrum for each source.

Note: The current implementation has been designed with WIRC+POL data in mind, but we hope to extent this to WIRC+Spec data soon.

Assumptions

We assume you have installed the DRP following the instructions in the github README (including the installation of all the dependencies): https://github.com/WIRC-Pol/wirc_drp

This Tutorial

This tutorial will give an example of reducing a raw WIRC+POL data image from start to finish, and will give you an idea of the implementation of the pipeline along the way. We’ll note that many of the data reduction steps are works in progress. The purpose of this tutorial is more to demonstrate the features of the pipeline rather than to teach you how to fully calibrate your data.

Step 0 - Make your calibration files

The very first step will be to create your own master darks and master flats. For the sake of this tutorial we have already done this for you, but we provide some example code commented-out below so you can do it yourself. Already-made master darks and flats have been provided in the tutorial directory.

#Import the important things
%matplotlib inline
import wirc_drp
import wirc_drp.wirc_object as wo
from wirc_drp.utils import calibration
import matplotlib.pyplot as plt
import numpy as np

import warnings; warnings.simplefilter('ignore')

Create a master dark

First you will need a list of files that will go into the master dark. It can be created like this:

$ ls dark*.fits > dark_list.dat

Then you read it in and create your mater dark by taking the median of all the files in the list. A hotpixel map will also be generated by default.

### Get the file list
# import astropy.io.ascii as asci
# import os
#os.chdir('/hcig1-nas/wircpol/data/20170611')
#dark_list_fn = "dark_list.dat" #The file name for your list of dark files
#dark_list1s = (asci.read(dark_list_fn, format = 'fast_no_header'))['col1'] #Read in the list

### Create the master dark and a bad pixel map.
### This function creates a new fits file based on the last filename in your dark_list and appends "_master_dark.fits"
### The hot pixel map will be the same, except with a "_bp_map.fits" suffix.
#dark_name1s, bp_name = calibration.masterDark(dark_list1s) # The output of this function is a the filename
                                                         # of the master dark and bad pixel maps

Create a master flat

First you will need a list of files that will go into the master flat. It can be created like this:

$ ls flat*.fits > flat_list.dat

Then you read it in and create your mater flat by subtracting the master dark from each image, taking the median of all the files in the list and then normalize by either the median (default) or the mode, set by the normalize keyword. A bad pixel map will also be generated by default, which you can supplement with the hot pixel map output by creating the master dark.

### Get the file list
# flat_list_fname = "flat_list.dat" #The file name for your list of flat files
# flat_list = (asci.read(flat_list_fname, format = 'fast_no_header'))['col1']

### Create the master dark and a bad pixel map. You will need the filename
### This function creates a new fits file based on the last filename in your flat list and appends "_master_flat.fits"

# flat_name, bp_name = calibration.masterFlat(flat_list, dark_name1, hotp_map_fname = None)

Step 1 - Read in and calibrate your data

Assuming that you now have a master dark and master flat file, we can now read in a wirc_data file and perform the calibration.

#First we'll set up all the directories and filenames:
wircpol_dir = wirc_drp.__path__[0]+"/../"
tutorial_dir = wircpol_dir + "Tutorial/sample_data/Single_File_Tutorial/"

raw_fn =tutorial_dir+"wirc1586.fits"
flat_fn = tutorial_dir+"wirc2012_master_flat.fits"
dark_fn = tutorial_dir+"wirc0141_master_dark.fits"
bp_fn = tutorial_dir+"wirc2012_bp_map.fits"

#Now we'll create the wirc_data object, passing in the filenames for the master dark, flat and bad pixel maps
raw_data = wo.wirc_data(raw_filename=raw_fn, flat_fn = flat_fn, dark_fn = dark_fn, bp_fn = bp_fn,verbose=True)

Creating a new wirc_data object from file /Users/maxwellmb/Dropbox (Personal)/Library/Python/wirc_drp/wirc_drp/../Tutorial/sample_data/Single_File_Tutorial/wirc1586.fits
Found a J-band filter in the header of file /Users/maxwellmb/Dropbox (Personal)/Library/Python/wirc_drp/wirc_drp/../Tutorial/sample_data/Single_File_Tutorial/wirc1586.fits

The wirc_data object holds the 2D data image in the property full_image. We can take a look at it now:

plt.figure(figsize=(8,8))
plt.imshow(raw_data.full_image, vmin=0, vmax=5000)

<matplotlib.image.AxesImage at 0x1c2fe1ce48>

We can see (even in this zoomed out image) that there are a bunch of bad pixels. Let’s run the calibration step. It will subtract the dark, divide by the flat, and interpolate over the bad pixels.

raw_data.calibrate()

/Users/maxwellmb/Dropbox (Personal)/Library/Python/wirc_drp/wirc_drp/wirc_object.py:229: RuntimeWarning: divide by zero encountered in true_divide
  self.full_image = self.full_image/master_flat
/Users/maxwellmb/Dropbox (Personal)/Library/Python/wirc_drp/wirc_drp/wirc_object.py:229: RuntimeWarning: invalid value encountered in true_divide
  self.full_image = self.full_image/master_flat
/Users/maxwellmb/Dropbox (Personal)/Library/Python/wirc_drp/wirc_drp/utils/calibration.py:736: RuntimeWarning: invalid value encountered in less_equal
  bad_pixel_map = np.logical_or(bad_pixel_map, redux_science <= 0)

plt.figure(figsize=(8,8))
plt.imshow(raw_data.full_image, vmin=0, vmax=5000)

<matplotlib.image.AxesImage at 0x1c2ffafb00>

There are a few remaining artifacts in the image, but overall, not bad!

We’re now done with calibration and so we want to save our newly calibrated file. Up until now these steps can apply to both spec mode and pol mode data. However, the rest of the tutorial will just demonstrate how to reduce POL data.

raw_data.save_wirc_object(tutorial_dir+"calibrated.fits",verbose=True)

Saving a wirc_object to /Users/maxwellmb/Dropbox (Personal)/Library/Python/wirc_drp/wirc_drp/../Tutorial/sample_data/Single_File_Tutorial/calibrated.fits

Step 2 - Extract the thumbnails of the spectral traces

In this step we will create s wircpol_source object and append. Each source object will contain cutout thumbnails of each spectral trace, and eventually spectra and polarized spectra. We will then append the wircpol_source object to the wirc object’s “source_list”.

#First we'll create a new data object, mostly just to demonstrate how to read in an existing wirc_data object.
calibrated_data = wo.wirc_data(wirc_object_filename=tutorial_dir+"calibrated.fits")

Loading a wirc_data object from file /Users/maxwellmb/Dropbox (Personal)/Library/Python/wirc_drp/wirc_drp/../Tutorial/sample_data/Single_File_Tutorial/calibrated.fits
Didn't find any background files

calibrated_data.generate_bkg()

Create a wircpol_source object

#This specific object will be for a souce in the slit (with coordinates [1063,1027]).
wp_source = wo.wircpol_source([1063,1027],'1',0) #The second argument indicates that this source is in the middle slit

We’ll now get the trace cutouts for this source

wp_source.get_cutouts(calibrated_data.full_image, calibrated_data.DQ_image, 'J',calibrated_data.bkg_image)
wp_source.plot_cutouts(figsize=(10,6))

Now let’s add this source to the calibrated_data source list.

calibrated_data.source_list.append(wp_source)
calibrated_data.n_sources += 1

We’ll now manually add another source in the field, that isn’t in the slit

calibrated_data.source_list.append(wo.wircpol_source([1035,640],'slitless',calibrated_data.n_sources+1))
calibrated_data.n_sources += 1
calibrated_data.source_list[1].get_cutouts(calibrated_data.full_image, calibrated_data.DQ_image,'J',
                                           calibrated_data.bkg_image)

Source 1's traces will hit the vertical bar of doom, compensating by subtracting the median of the edges of each row
Source 1's traces will hit the vertical bar of doom, compensating by subtracting the median of the edges of each row

We can also look at what the background subtracted data look like:

calibrated_data.source_list[1].plot_cutouts(figsize=(10,6),plot_bkg_sub=True,vmin=-200,vmax=400)

Step 3 - Extract the spectra from each thumbnail

In this step we will extract the specta from one of the source objects we found in Step 2.

#Here you can enable "plot" to see where the traces were "found"
calibrated_data.source_list[0].extract_spectra(verbose=False)

/Users/maxwellmb/Dropbox (Personal)/Library/Python/wirc_drp/wirc_drp/utils/spec_utils.py:553: RuntimeWarning: invalid value encountered in true_divide
  flux_opt_final = np.nan_to_num(sum1/sum3)
/Users/maxwellmb/Dropbox (Personal)/Library/Python/wirc_drp/wirc_drp/utils/spec_utils.py:554: RuntimeWarning: invalid value encountered in true_divide
  variance_opt_final = np.nan_to_num(sumP/sum3)
/Users/maxwellmb/Dropbox (Personal)/Library/Python/wirc_drp/wirc_drp/utils/spec_utils.py:457: RuntimeWarning: divide by zero encountered in true_divide
  P_0 = (data - background)/flux_0 #this is dividing each column (x) in data-background by the sum in that column
/Users/maxwellmb/Dropbox (Personal)/Library/Python/wirc_drp/wirc_drp/utils/spec_utils.py:457: RuntimeWarning: invalid value encountered in true_divide
  P_0 = (data - background)/flux_0 #this is dividing each column (x) in data-background by the sum in that column
/Users/maxwellmb/Dropbox (Personal)/Library/Python/wirc_drp/wirc_drp/utils/spec_utils.py:462: RuntimeWarning: invalid value encountered in less
  P_0[P_0 < 0] = 0
/Users/maxwellmb/Dropbox (Personal)/Library/Python/wirc_drp/wirc_drp/utils/spec_utils.py:533: RuntimeWarning: invalid value encountered in less
  Mask = (data - background - flux_0*P_0)**2 < sig_clip**2*variance_opt

#And we'll plot the data!
calibrated_data.source_list[0].plot_trace_spectra(figsize=(12,8))

Now we’ll do a rough wavelength calibration. Currently it uses the edge of the filters to shift and scale the spectra. This will eventually need to be refined as we’re not super happy with it.

## Now we'll do a rough wavelength calibration and plot the data again
calibrated_data.source_list[0].rough_lambda_calibration(method=2)
#Note method=1 is broken.

#And we'll plot the data again.
calibrated_data.source_list[0].plot_trace_spectra(figsize=(10,6),xlow=1.16,xhigh=1.34)

Step 4 - Calculate the polarized spectra

We now combine the 4 spectra to calculate Stokes Q and U.

We’ll note that at this point it is a “single-differenced” spectrum, and is subject to instrumental biases. These can be easily removed in data taken after Feb 2019 when we installed a modulator into WIRC+Pol. Tutorial 3 covers how to process data taken with the modulator (though don’t skip Tutorial 2!).

These instrumental biases can also be calibrated, but to a worse precision, in data before Feb 2019 date using an unpolarized standard taken at a similar telescope pointing.

calibrated_data.source_list[0].compute_polarization(cutmin=20, cutmax=150)

/Users/maxwellmb/Dropbox (Personal)/Library/Python/wirc_drp/wirc_drp/utils/spec_utils.py:1751: RuntimeWarning: invalid value encountered in true_divide
  u = (Up-Um)/(Up+Um)

calibrated_data.source_list[0].plot_Q_and_U(figsize=(8,3),xlow=1.17,xhigh=1.32,ylow=-0.1,yhigh=0.1)

Step 5 - Save the wirc object

We now save the new information and tables to a new file Note how it initiates the table columns even if the data hasn’t been computed (for example, we haven’t computed source_list[1] yet)

calibrated_data.save_wirc_object(tutorial_dir+"calibrated.fits")

Saving a wirc_object to /Users/maxwellmb/Dropbox (Personal)/Library/Python/wirc_drp/wirc_drp/../Tutorial/sample_data/Single_File_Tutorial/calibrated.fits