Applying Customized Background Regions to LETG/HRC-S Observations

Introduction to tgextract2

tgextract2 is an alternative to tgextract, the CIAO tool which filters and bins spectra. tgextract2 provides support for variable-shaped background regions which do not need to have a constant width ratio with the source region. tgextract2 computes a vector "backscale"—the scaling factor which relates the background count rate at any wavelegth to the source extraction region width at that wavelength. tgextract2 also sums the background extractions from each side of the spectrum, whereas tgextract writes two counts arrays, one from the "_UP" and one from the "_DOWN" sides of the spectrum.

Since tgextract2 is not backward-compatible with tgextract, it is not used in standard processing. It is instead recommended for special cases where some detailed and interactive analysis is required, involving some scientific judgement by the analyst.

Getting Started

Download the sample data: 29 (LETG/HRC-S, Alpha Cen)

unix% download_chandra_obsid 29 evt2

The standard HRC grating region file, found in CalDB, is also used:

ln -s $CALDB/data/chandra/hrc/tgmask2/letgD1999-07-22regN0002.fits .

Open the event file in DS9:

unix% ds9 hrcf00029N006_evt2.fits &

Next, we set the height and width to something convenient for grating spectra, and adjust the colorscale and colormap in the DS9 menu.

• Frame → Frame Parameters → Set Display Size → X=1024, Y=400
• Scale → Log
• Color → BB

Figure 1 shows the resulting image.

Figure 1: Source Image

[Version: full-size]

Figure 1: Source Image

Source image with BB colormap and log-scale selected

Next, we set the binning parameters in the DS9 menu Bin → Binning Parameters (Figure 2):

Figure 2: Grating Coordinates Binning Parameters

Figure 2: Grating Coordinates Binning Parameters

Binning parameter settings applied in DS9 by going to the menu and choosing Bin → Binning Parameters.

The binning was iteratively chosen such that the field fits the screen and it has enough counts/bin to give good statistics. Fiddling with parameters, we start off with knowing that the tg_r coordinate range is -1 to 1 degrees in about 10000 bins, or 0.0002 deg/bin. Binning by a factor of 10, we have a tg_d range of about -0.002 to 0.002 deg/bin.

After applying the binning parameters, the field image in grating dispersion coordinates shows two spectra (Figure 3):

Figure 3: Grating Spectrum of the Two Field Sources

[Version: full-size]

Figure 3: Grating Spectrum of the Two Field Sources

This is the image in grating dispersion coordinates with tg_r (in-dispersion) on the horizontal axis, and tg_d (cross-dispersion) on the vertical axis (both in degrees). We can clearly see the two sources zero-orders near physical coordinates (0,0) and their dispersed and line-rich spectra, with negative-orders to the left and positive-orders to the right. The vertical dark regions are gaps between the HRC-S microchannel plates. Overplotted is the standard LETG extraction region, where, in this case, the lower background region includes the spectrum of the second source and therefore should not be used.

Customizing the Extraction Regions

The standard HRC/LETG extraction region file must be edited for the current field. This is most easily done in DS9 as long as the information about which region represents source and up and down background is preserved. The tgmask2reg and reg2tgmask scripts allow the use of DS9 to define custom extraction regions, and use the DS9 region "tagging" feature to encode the grating information along with the extraction region coordinates.

The first step in the editing process requires running tgmask2reg to convert the CALDB tgmask file into a DS9 format. The additional (non-coordinate) columns in the CALDB files are added as "tags" so that later tools know each regions' purpose.

unix% tgmask2reg infile=CALDB outfile=tg.reg

unix% more tg.reg 
# Region file format: DS9 version 4.1
# Filename:
global color=green dashlist=8 3 width=1 font="helvetica 10 normal roman" select=1 highlite=1 dash=0 fixed=0 edit=1 move=1 delete=0 include=1 source=1
physical
polygon(-0.0057804,0.000531,-0.29867,0.000531,-1.1549,0.002114,-1.1549,-0.002114,-0.29867,-0.000531,-0.0057804,-0.000531,-0.0057804,0.000531) # tag={row_idx=1} tag={rowid=SOURCE} tag={tg_part=3} tag={tg_srcid=1} tag={tg_m=-1} tag={backscal=1.0}
polygon(-0.0057804,0.00731,-0.29867,0.00731,-1.1549,0.02314,-1.1549,0.002,-0.29867,0.002,-0.0057804,0.002,-0.0057804,0.00731) # tag={row_idx=2} tag={rowid=BACKGROUND_UP} tag={tg_part=3} tag={tg_srcid=1} tag={tg_m=-1} tag={backscal=5.0}
polygon(-0.0057804,-0.00731,-0.29867,-0.00731,-1.1549,-0.02314,-1.1549,-0.002,-0.29867,-0.002,-0.0057804,-0.002,-0.0057804,-0.00731) # tag={row_idx=3} tag={rowid=BACKGROUND_DOWN} tag={tg_part=3} tag={tg_srcid=1} tag={tg_m=-1} tag={backscal=5.0}
polygon(0.0057804,0.000531,0.29867,0.000531,1.1549,0.002114,1.1549,-0.002114,0.29867,-0.000531,0.0057804,-0.000531,0.0057804,0.000531) # tag={row_idx=1} tag={rowid=SOURCE} tag={tg_part=3} tag={tg_srcid=1} tag={tg_m=1} tag={backscal=1.0}
polygon(0.0057804,0.00731,0.29867,0.00731,1.1549,0.02314,1.1549,0.002,0.29867,0.002,0.0057804,0.002,0.0057804,0.00731) # tag={row_idx=2} tag={rowid=BACKGROUND_UP} tag={tg_part=3} tag={tg_srcid=1} tag={tg_m=1} tag={backscal=5.0}
polygon(0.0057804,-0.00731,0.29867,-0.00731,1.1549,-0.02314,1.1549,-0.002,0.29867,-0.002,0.0057804,-0.002,0.0057804,-0.00731) # tag={row_idx=3} tag={rowid=BACKGROUND_DOWN} tag={tg_part=3} tag={tg_srcid=1} tag={tg_m=1} tag={backscal=5.0}

Once the region is loaded into DS9 (Region → Load Regions → tg.reg), it can be customized interactivelty as needed either by moving points or adding new ones both for the source and the background areas. Note that polygon regions cannot added or deleted. They can only be modified. This ensures that the needed meta-data is preserved.

Note that in DS9 regions can be inspected to see the extra information associated with each region by going to Region → Groups . This will provide a list of the different tags in use. Selecting one of the tags will show and select the corresponding region(s) in the main ds9 window.

The final customized set of regions must be saved in DS9 format in physical coordinates (Region → Save Regions → ds9.reg) (see Figure 4).

Figure 4: Customized Extraction Regions

[Version: full-size]

Figure 4: Customized Extraction Regions

This is the modified extraction region for source and background: the lower background regions have been moved well away from the second source spectrum while the source region has been left unaltered.

Finally the reg2tgmask tool converts the modified region ("ds9.reg") back into a FITS format ("tg_mask.fits") compatible with tgextract2.

unix% reg2tgmask infile=ds9.reg srcfile=CALDB out=tg_mask.fits 

Running tgextract2

With the background regions appropriately defined, we are now ready to run tgextract2. The opt parameter is set to "pha1" so that the tool creates separate files for the +1 and -1 spectra (Type 1 PHA file), rather than both spectra be written as rows in the same file (Type 2 PHA file); see the grating spectra data products guide for more information these formats. The files will also contain the background extractions and the "areas" (BACKSCAL) of the source and background regions for each spectral bin.

The outfile parameter value (a_) is used as the root name of the output files. The files will be named as "root<srcid><part><order>_pha2.fits", one for each source ID, part, and order combination.

unix% punlearn tgextract2
unix% pset tgextract2 infile=hrcf00029N006_evt2.fits
unix% pset tgextract2 outfile=a_ 
unix% pset tgextract2 region_file=tg_mask.fits

unix% pset tgextract2 opt=pha1

unix% tgextract2 
Input event file (output event file from L1.5 processing) (hrcf00029N006_evt2.fits):
Output pha file (for pha2, enter full file name or '.'; for pha1, enter rootname.) (a_):
Source ID's to process: 'all', comma list, @file (1):
Grating parts to process: HETG, HEG, MEG, LETG, header_value (HETG|HEG|MEG|LETG|header_value) (header_value):
Grating diffraction orders to process: 'default', comma list, range list, @file (default):

This produces two new files, a_01Lm1_pha2.fits and a_01Lp1_pha2.fits, which are -1 order and +1 order LEG spectra, respectively.

You can check the parameter file that was used with plist tgextract2.

After we run tgextract2, we will copy the background counts to a separate file, one for each order. This is because various analysis softwares typically need to have the background spectrum formatted the same way as the source (i.e. having the same column names).

unix% dmcopy 'a_01Lm1_pha2.fits[cols CHANNEL,BIN_LO,BIN_HI,COUNTS=bg_counts,BACKSCAL=bg_area]' a_01Lm1_pha2_bg.fits
unix% dmcopy 'a_01Lp1_pha2.fits[cols CHANNEL,BIN_LO,BIN_HI,COUNTS=bg_counts,BACKSCAL=bg_area]' a_01Lp1_pha2_bg.fits
unix% ls a_*
a_01Lm1_pha2.fits     a_01Lm1_pha2_bg.fits  a_01Lp1_pha2.fits     a_01Lp1_pha2_bg.fits

Inspecting the Counts and Background Spectra with Sherpa

We will now inspect the counts and background spectra. Using Sherpa, we will create a strip-chart displaying the plus- and minus-order counts spectra, with background over-plotted. We begin by loading the data and background with the minus-order spectrum associated with dataset ID=1 and positive-order spectrum associated with dataset ID=2.

unix% sherpa

sherpa In []: load_pha(1,"a_01Lm1_pha2.fits")
sherpa In []: load_pha(2,"a_01Lp1_pha2.fits")
sherpa In []: load_bkg(1,"a_01Lm1_pha2_bg.fits")
sherpa In []: load_bkg(2,"a_01Lp1_pha2_bg.fits")

Next we set Sherpa plotting preferences. The key parameter is set_analysis, which plots the spectrum as counts vs. wavelength.

sherpa In []: set_analysis("wave",type="counts")

Loading a two plot strip chart, we first add the minus-order counts spectrum to the upper plot:

sherpa In []: %matplotlib
sherpa In []: import matplotlib.pylab as plt
sherpa In []: plt.subplot(2,1,1)
sherpa In []: plt.subplots_adjust(hspace=0.0 )
sherpa In []: plt.plot(get_data_plot(1).x,get_data_plot(1).y,marker="None", color="black")

Next, we over-plot the minus-order count spectrum with the minus-order background in red:

sherpa In []: bs1=get_bkg_scale(1)
sherpa In []: ct1=get_counts(1,bkg_id=1)

sherpa In []: plt.plot(get_data_plot(1).x,(ct1*bs1),marker="None", color="red")

and add some labels to plot

sherpa In []: plt.ylabel("Counts") 
sherpa In []: plt.figtext(0.93,0.75,"Negative-Order",rotation=270)
sherpa In []: plt.title('Counts and Background (red) for Negative- and Positive-Orders')
sherpa In []: plt.xlim(0, 180)
sherpa In []: plt.gca().set_xticks([]) # hide x-ticks

Now we populate the lower plot with the plus-order counts and background spectra in the same manner

sherpa In []: plt.subplot(2,1,2)
sherpa In []: plt.plot(get_data_plot(2).x,get_data_plot(2).y, marker="None",color="black")
sherpa In []: 
sherpa In []: bs2=get_bkg_scale(2)
sherpa In []: ct2=get_counts(2,bkg_id=1) # 1st (only) background associated with dataset 2
sherpa In []: plt.plot(get_data_plot(2).x,(ct2*bs2),marker="None", color="red")
sherpa In []: 
sherpa In []: plt.xlabel(r"$\lambda (\AA)$")
sherpa In []: plt.ylabel("Counts")
sherpa In []: plt.figtext(0.93,0.35,"Positive-Order",rotation=270)
sherpa In []: plt.xlim(0, 180)
sherpa In []: plt.ylim(0,35)

The results of this plot can be seen in Figure 5.

Figure 5: Spectra of Field Sources with Background Spectra over-plotted in red

[Version: full-size]

Figure 5: Spectra of Field Sources with Background Spectra over-plotted in red

Strip chart of the spectra of counts (white) and background (red) with negative orders in the upper frame and postive orders in the lower frame.

Before exiting Sherpa, you may want to write the plot to a file. In this example, we use the interactive print dialog.

sherpa In []: plt.savefig("letgs_spec_bg.png")

Events in Wavelength Coordinates

To confirm that the edits made to the region were well-chosen, we inspect the events in wavelength coordinates. With the event file displayed in DS9:

unix% ds9 hrcf00029N006_evt2.fits &

we change the parameters to bin on wavelength coordinates (tg_lam,tg_d) (Bin → Binning Parameters). Again, we have to fiddle with our binning parameters, but we have some starting points. We can keep the tg_d binning of 0.00015 since the cross-dispersion ranges between -1 to 1 degree over about 10000 bins. The nominal wave-band coverage for the LETGS when using HRC-S is ≅150-175 Å. In wavelength, LETGS has 0.015 Å/bin; binning up by a factor of 10, we can get a good plot.

Figure 6: Wavelength Coordinates Binning Parameters

Figure 6: Wavelength Coordinates Binning Parameters

Binning parameter settings applied in DS9 by going to the menu and choosing Bin → Binning Parameters.

Then load the region which is stored in the PHA file, a_01Lp1_pha2.fits[REGION]. Figure 7 shows the events and extraction regions in wavelength coordinates.

Figure 7: Grating Spectrum of the Two Field Sources in Wavelength Coordinates

[Version: full-size]

Figure 7: Grating Spectrum of the Two Field Sources in Wavelength Coordinates

This is the field image in wavelength coordinates with tg_lam (wavelength) on the horizontal axis, and tg_d (cross-dispersion) on the vertical axis. The region from a_01Lp1_pha2.fits is overlaid, showing the source and background.

Scripting It

The file, tgextract2.py, performs the Sherpa commands used above.