Chandra X-Ray Observatory (CXC)
Skip to the navigation links
Search the CIAO website or contact the CXC HelpDesk
Last modified: 3 Dec 2012

URL: http://cxc.harvard.edu/ciao/threads/extended/thread.html
Return to Threads Page: Top | All | Imag Spec

Extract Spectrum and Response Files for an Extended Source

CIAO 4.5 Science Threads

Overview

Synopsis:

Using a combination of CIAO tools, we extract source and background spectra for an extended source. The background spectrum is grouped, if desired. The appropriate Response Matrix Files (RMFs) and Ancillary Response Files (ARFs) are created for both source and background.

The specextract script automates these steps for extended and pointlike sources observed with the ACIS detector.

Purpose:

To generate source and background spectra of an extended ACIS source and build the proper RMFs and ARFs.

Software & Calibration Updates:

This thread requires the following updates to the standard CIAO 4.5 installation

Related Links:

Last Update: 3 Dec 2012 - Review for CIAO 4.5; minor formatting

Contents

Get Started

Download the sample data: 869 (ACIS-S, ARP 220)

unix% download_chandra_obsid 869 evt2,asol,bpix,pbk,msk

How specextract Works

specextract runs the following tools for the extended source case:

  • dmextract: to extract source and (optionally) background spectra. This tool also creates the WMAP used as input to mkacisrmf.
  • sky2tdet: to create the WMAP input for mkwarf.
  • mkwarf: to create weighted ARF(s).
  • mkrmf or mkacisrmf: to build the RMF(s), depending on which is appropriate for the data and the calibration; see the Creating ACIS RMFs why topic for details.
  • dmgroup: to group the source spectrum and/or background spectrum.
  • dmhedit: to update the BACKFILE, RESPFILE and ANCRFILE keys in the source and background spectrum files.

Build Source and Background Regions

We need to define two regions, one for the source and another for the background. To do this, first display the file: 

unix% ds9 acisf00869N004_evt2.fits &

Refer to the Using CIAO Regions thread for information on creating region files. The CIAO-format files for this example are: 

unix% cat simple.reg 
# Region file format: CIAO version 1.0
ellipse(4026,4138.9,50,40,0)

unix% cat simple_bkg.reg
# Region file format: CIAO version 1.0
circle(3975,4233,20)

The regions are shown displayed on the event file in Figure 1; the source is a dotted white line and the background is a solid green line.

Figure 1: Extraction regions on the event file

[Thumbnail image: The source region and the background region are each a single ellipse.]

[Version: full-size]

[Print media version: The source region and the background region are each a single ellipse.]

Figure 1: Extraction regions on the event file

The background was chosen from a source-free area of the same chip for this example, but it may also be chosen from a different chip or different event file.

Choosing a background file

The examples in this thread use a background region that was chosen from a source-free area of the same chip as the source region. The source or field of sources may cover the chip, however, making it impossible to select a local background. There are several options in this case:

  • Choose the background from another same-type chip that was on in the observation. If your source is on ACIS-S3 and ACIS-S1 (the other back-illuminated chip) was also turned on, define the background on ACIS-S1.

  • If using another chip from the observation is not an option, use one of the ACIS "blank-Sky" background files which are distributed in the CALDB.

    If you plan on using one of the ACIS "blank-sky" background files from the CALDB and intend to subtract the background spectrum (i.e. not fit it), then you do not need to create a background RMF and ARF. Set the bkgresp parameter to "no" and the script will only extract a background spectrum.

  • Omit the background completely (i.e. leave the bkgfile blank).

Using specextract

1. Set the input parameters

Input the event file with the appropriate region file for the source and background, as well as the aspect solution files. In many cases, there will be more than one aspect solution file (pcad_asol1.fits) for an observation. All the files must be input to the asp parameter, either as a list or as a stack. Here we use: 

unix% cat pcad_asol1.lis 
pcadf078246822N003_asol1.fits

unix% punlearn specextract
unix% pset specextract infile="acisf00869N004_evt2.fits[sky=region(simple.reg)]"
unix% pset specextract outroot=simple
unix% pset specextract bkgfile="acisf00869N004_evt2.fits[sky=region(simple_bkg.reg)]"
unix% pset specextract asp=@pcad_asol1.lis
unix% pset specextract pbkfile=acisf078247287N003_pbk0.fits
unix% pset specextract mskfile=acisf00869_000N003_msk1.fits
unix% pset specextract badpixfile=acisf00869_000N003_bpix1.fits
unix% pset specextract weight=yes correct=no
unix% pset specextract grouptype=NUM_CTS binspec=15

Grouping the output spectra

We choose to use the default grouping values: the source spectrum will be grouped to a minimum number of 15 counts per new channel and the background spectrum will be ungrouped. The grouptype and binspec parameters are used to specify the source spectrum grouping, and the bkg_grouptype and bkg_binspec parameters specify the background spectrum grouping.

Setting the energy range for the WMAP

The energy_wmap parameter determines the energy band used to create the weight map when weight=yes.

If you use the soft band (e.g. the default value of 300:2000 keV), then you are weighting by the emission from the source (unless you have a very hard spectrum). To weight by the effective area of the telescope (again, unless you have a hard source), set the range to, for example, 2000:7000 keV.

Unless the spectrum of your source is unusual, adjusting the WMAP energy range will not affect the ARF - or, therefore, the spectral fit - noticeably.

2. Run the script

Running the tool with verbose=2 shows what it is doing: 

Source event file(s) (acisf00869N003_evt2.fits[sky=region(simple.reg)]): 
Output directory path + root name for output files (simple): 
Should response files be weighted? (yes): 
Apply point source aperture correction to ARF? (no): 
Combine ungrouped output spectra and responses? (no): 
Background event file(s) (acisf00869N003_evt2.fits[sky=region(simple_bkg.reg)]):
Create background ARF and RMF? (yes): 
Source aspect solution or histogram file(s) (@pcad_asol1.lis): 
pbkfile input to mkwarf (acisf078247287N003_pbk0.fits): 
mskfile input to mkwarf (acisf00869_000N003_msk1.fits): 
  
Running: specextract
  Version:    14 February 2012

Checking parameter block file(s) for readability...

Checking dead area file(s) for readability...

Checking mask file(s) for readability...

Checking bad pixel file(s) for readability...

Checking aspect file(s) for readability...

Checking source file(s) for readability...

Checking background file(s) for readability...
Running tool dmlist using:
>>> dmlist "acisf00869N003_evt2.fits[sky=region(simple.reg)]" counts

Setting bad pixel file for item 1 of 1 in input list

Running tool acis_set_ardlib using:
>>> acis_set_ardlib acisf00869_000N003_bpix1.fits verbose=2

Extracting src spectra for item 1 of 1 in input list

Running tool dmextract using:
>>> dmextract "acisf00869N003_evt2.fits[sky=region(simple.reg)][bin PI]" simple.pi wmap="[energy=300:2000][bin tdet=8]" verbose=2

Creating src ARF for item 1 of 1 in input list

Running tool dmlist using:
>>> dmlist "acisf00869N003_evt2.fits[sky=region(simple.reg)]" counts
Running tool dmlist using:
>>> dmlist "acisf00869N003_evt2.fits[sky=region(simple.reg)]" counts
Running tool asphist using:
>>> asphist @pcad_asol1.lis simple_asphist7.fits "acisf00869N003_evt2.fits[ccd_id=7]" dtffile= verbose=2
Running tool sky2tdet using:
>>> sky2tdet "acisf00869N003_evt2.fits[sky=region(simple.reg)][energy=300:2000][bin sky]" simple_asphist7.fits "simple_tdet.fits[wmap]" verbose=2 clobber=yes
Running tool mkwarf using:
>>> mkwarf "simple_tdet.fits[wmap]" simple.warf simple.wfef spectrumfile= egridspec=0.3:11.0:0.01 pbkfile=acisf078247287N003_pbk0.fits mskfile=acisf00869_000N003_msk1.fits verbose=2

Creating src RMF for item 1 of 1 in input list

Searching for P2_RESP calibration file...
Running tool calquiz using:
>>> calquiz "simple.pi[WMAP]" telescope= instrument= product=SC_MATRIX calfile="CALDB(ccd_id=7)" outfile=y echo=yes verbose=2

Using mkacisrmf...

Running tool mkacisrmf using:
>>> mkacisrmf CALDB simple.wrmf "simple.pi[WMAP]" 0.3:11.0:0.01 1:1024:1 PI 0 INDEF INDEF CALDB verbose=2
Running tool dmgroup using:
>>> dmgroup "simple.pi[SPECTRUM]" simple_grp.pi NUM_CTS 15.0 binspec= xcolumn=channel ycolumn=counts verbose=2
Updating header of simple.pi with RESPFILE and ANCRFILE keywords.
Running tool dmhedit using:
>>> dmhedit simple.pi none add RESPFILE simple.wrmf verbose=2
Running tool dmhedit using:
>>> dmhedit simple.pi none add ANCRFILE simple.warf verbose=2
Updating header of simple_grp.pi with RESPFILE and ANCRFILE keywords.
Running tool dmhedit using:
>>> dmhedit simple_grp.pi none add RESPFILE simple.wrmf verbose=2
Running tool dmhedit using:
>>> dmhedit simple_grp.pi none add ANCRFILE simple.warf verbose=2
Running tool dmlist using:
>>> dmlist "acisf00869N003_evt2.fits[sky=region(simple_bkg.reg)]" counts

Setting bad pixel file for item 1 of 1 in input list

Running tool acis_set_ardlib using:
>>> acis_set_ardlib acisf00869_000N003_bpix1.fits verbose=2

Extracting bkg spectra for item 1 of 1 in input list

Running tool dmextract using:
>>> dmextract "acisf00869N003_evt2.fits[sky=region(simple_bkg.reg)][bin PI]" simple_bkg.pi wmap="[energy=300:2000][bin tdet=8]" verbose=2

Creating bkg ARF for item 1 of 1 in input list

Running tool dmlist using:
>>> dmlist "acisf00869N003_evt2.fits[sky=region(simple_bkg.reg)]" counts
Running tool dmlist using:
>>> dmlist "acisf00869N003_evt2.fits[sky=region(simple_bkg.reg)]" counts
Running tool asphist using:
>>> asphist @pcad_asol1.lis simple_bkg_asphist7.fits "acisf00869N003_evt2.fits[ccd_id=7]" dtffile= verbose=2
Running tool sky2tdet using:
>>> sky2tdet "acisf00869N003_evt2.fits[sky=region(simple_bkg.reg)][energy=300:2000][bin sky]" simple_bkg_asphist7.fits "simple_bkg_tdet.fits[wmap]" verbose=2 clobber=yes
Running tool mkwarf using:
>>> mkwarf "simple_bkg_tdet.fits[wmap]" simple_bkg.warf simple_bkg.wfef spectrumfile= egridspec=0.3:11.0:0.01 pbkfile=acisf078247287N003_pbk0.fits mskfile=acisf00869_000N003_msk1.fits verbose=2

Creating bkg RMF for item 1 of 1 in input list

Searching for P2_RESP calibration file...
Running tool calquiz using:
>>> calquiz "simple_bkg.pi[WMAP]" telescope= instrument= product=SC_MATRIX calfile="CALDB(ccd_id=7)" outfile=y echo=yes verbose=2

Using mkacisrmf...

Running tool mkacisrmf using:
>>> mkacisrmf CALDB simple_bkg.wrmf "simple_bkg.pi[WMAP]" 0.3:11.0:0.01 1:1024:1 PI 0 INDEF INDEF CALDB verbose=2
Updating header of simple_bkg.pi with RESPFILE and ANCRFILE keywords.
Running tool dmhedit using:
>>> dmhedit simple_bkg.pi none add RESPFILE simple_bkg.wrmf verbose=2
Running tool dmhedit using:
>>> dmhedit simple_bkg.pi none add ANCRFILE simple_bkg.warf verbose=2
Updating header of simple.pi with BACKFILE keyword.
Running tool dmhedit using:
>>> dmhedit simple.pi none add BACKFILE simple_bkg.pi verbose=2
Updating header of simple_grp.pi with BACKFILE keyword.
Running tool dmhedit using:
>>> dmhedit simple_grp.pi none add BACKFILE simple_bkg.pi verbose=2

The contents of the parameter file may be checked with plist specextract.

3. Examining the Output Files

The number of files created depends on if a background event file was provided and if source and/or background grouping was specified. In this case, the output files are: 

Source:
simple.pi         ungrouped spectrum
simple_tdet.fits  WMAP from sky2tdet
simple.warf       weighted ARF
simple.wfef       FEF weight file
simple.wrmf       weighted RMF
simple_grp.pi     grouped spectrum

Background:
simple_bkg.pi         ungrouped spectrum
simple_bkg_tdet.fits  WMAP from sky2tdet
simple_bkg.warf       weighted ARF
simple_bkg.wfef       FEF weight file
simple_bkg.wrmf       weighted RMF

The FEF weight files (.wfef) are required as input to mkrmf; they are created by mkwarf before the scripts determines whether to use mkrmf or mkacisrmf. After the RMFs are created, these files are no longer needed.

Step-by-Step Guide

This section explains how to manually run the analysis steps that are used by specextract.

Make sure that you have set up ardlib to use the bad pixel file for your observation before following this section of the thread.

1. Extract Source and Background Spectra

First, we extract spectra from the events enclosed by the source and background regions.

To create a WMAP block in a PHA file, dmextract must be given a binning specification for its wmap option. This WMAP will be used as input to mkacisrmf. It must be in detector coordinates, and the suggested binning factor is 8; an optional energy filter is included for the WMAP. 

unix% punlearn dmextract
unix% pset dmextract infile="acisf00869N004_evt2.fits[sky=region(simple.reg)][bin pi]"
unix% pset dmextract outfile=simple_steps.pi
unix% pset dmextract wmap="[energy=300:2000][bin tdet=8]"
unix% dmextract
Input event file  (acisf00869N004_evt2.fits[sky=region(simple.reg)][bin pi]): 
Enter output file name (simple_steps.pi):

And for the background spectrum: 

unix% pset dmextract infile="acisf00869N004_evt2.fits[sky=region(simple_bkg.reg)][bin pi]"
unix% pset dmextract outfile=simple_steps_bkg.pi
unix% dmextract
Input event file  (acisf00869N004_evt2.fits[sky=region(simple_bkg.reg)][bin pi]): 
Enter output file name (simple_steps_bkg.pi):

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

2. Calculate the ARFs

Compute the Aspect Histogram (asphist)

We now need to create the aspect histogram, which is a binned representation of aspect motion during the observation.

In many cases, there will be more than one aspect solution file (pcad_asol1.fits) for an observation. All the files must be input to the infile parameter, either as a list or as a stack. Here we use: 

unix% cat pcad_asol1.lis 
pcadf078246822N003_asol1.fits

We also need to know what chip the source and background regions are on: 

unix% dmstat "acisf00869N004_evt2.fits[sky=region(simple.reg)][cols ccd_id]"
ccd_id
    min:        7             @:        1 
    max:        7             @:        1 
   mean:        7 
  sigma:        0 
    sum:        17752 
   good:        2536 
   null:        0 

unix% dmstat "acisf00869N004_evt2.fits[sky=region(simple_bkg.reg)][cols ccd_id]"
    min:        7             @:        1 
    max:        7             @:        1 
   mean:        7 
  sigma:        0 
    sum:        1715 
   good:        245 
   null:        0

Each region is on ccd_id=7. This means we only need to make one aspect histogram; if the regions were on different chips, we would make one for each chip. Now run the tool: 

unix% punlearn asphist
unix% pset asphist infile=@pcad_asol1.lis
unix% pset asphist outfile=simple_steps.asphist
unix% pset asphist evtfile="acisf00869N004_evt2.fits[ccd_id=7]"
unix% asphist
Aspect Solution List Files (@pcad_asol1.lis): 
Aspect Histogram Output File (simple_steps.asphist): 
Event List Files (acisf00869N004_evt2.fits[ccd_id=7]):

Create a WMAP for the ARF (sky2tdet)

The WMAP created by sky2tdet properly weights the ARF based on how much of the source flux fell onto the bad pixels, columns, or a node boundary, and which bad pixels are actually exposed. Without accounting for these effects, the ARF is significantly over-estimated. 

unix% punlearn sky2tdet
unix% pset sky2tdet infile="acisf00869N004_evt2.fits[sky=region(simple.reg)][energy=300:2000][bin sky=1]"
unix% pset sky2tdet asphistfile="simple_steps.asphist"
unix% pset sky2tdet outfile="simple_steps_tdet.fits[wmap]"
unix% sky2tdet
Input image in sky (x,y) coordinates (acisf00869N004_evt2.fits[sky=region(simple.reg)][energy=300:2000][bin sky=1]): 
Input aspect histogram file (simple_steps.asphist): 
Output TDET  WMAP file (simple_steps_tdet.fits[wmap]):

If you made a separate aspect histogram for the background, set it in the asphistfile parameter before running the tool for the background spectrum. 

unix% pset sky2tdet infile="acisf00869N004_evt2.fits[sky=region(simple_bkg.reg)][energy=300:2000][bin sky=1]"
unix% pset sky2tdet outfile="simple_steps_bkg_tdet.fits[wmap]"
unix% sky2tdet
Input image in sky (x,y) coordinates (acisf00869N004_evt2.fits[sky=region(simple_bkg.reg)][energy=300:2000][bin sky=1]): 
Input aspect histogram file (simple_steps.asphist): 
Output TDET  WMAP file (simple_steps_bkg_tdet.fits[wmap]):

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

Compute the ARFs (mkwarf)

mkwarf creates a weighted ARF file based on the WMAP that was created with sky2tdet. It also produces a set of weights that are used by mkrmf to create the associated weighted RMF.

Run mkwarf for the source: 

unix% punlearn mkwarf
unix% pset mkwarf infile="simple_steps_tdet.fits[wmap]"
unix% pset mkwarf outfile=simple_steps.warf
unix% pset mkwarf weightfile=simple_steps.wfef
unix% pset mkwarf egridspec=0.3:11.0:0.01 
unix% pset mkwarf pbkfile=acisf078247287N003_pbk0.fits 
unix% pset mkwarf mskfile=acisf00869_000N003_msk1.fits

unix% mkwarf
Input detector WMAP (simple_steps_tdet.fits[wmap]): 
Output weighted ARF file (simple_steps.warf): 
Output FEF weights (simple_steps.wfef): 
Input Spectral weighting file (<filename>|NONE) (): 
Output energy grid [kev] (0.3:11.0:0.01): 
Parameter block file (acisf078247287N003_pbk0.fits):

and for the background: 

unix% pset mkwarf infile="simple_steps_bkg_tdet.fits[wmap]"
unix% pset mkwarf outfile=simple_steps_bkg.warf
unix% pset mkwarf weightfile=simple_steps_bkg.wfef
unix% mkwarf 
Input detector WMAP (simple_steps_bkg_tdet.fits[wmap]): 
Output weighted ARF file (simple_steps_bkg.warf): 
Output FEF weights (simple_steps_bkg.wfef): 
Input Spectral weighting file (<filename>|NONE) (): 
Output energy grid [kev] (0.3:11.0:0.01): 
Parameter block file (acisf078247287N003_pbk0.fits):

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

3. Calculate the RMFs

The observation used in this thread (ObsID 869) was taken at the -120 C focal plane temperature. Therefore, mkacisrmf should be used to create the RMF file.

The syntax for both mkacisrmf and mkrmf are given in this section. Users should choose the appropriate tool for the data and calibration.

Using mkacisrmf (mkacisrmf)

The mkacisrmf why topic has more details on using the mkacisrmf tool.

The RMF doesn't change quickly enough to require the detail of the sky2tdet WMAP, and its precision will cause mkacisrmf to run very slowly. Nearly identical results are obtained with the dmextract WMAP with a more reasonable runtime.

For the source: 

unix% punlearn mkacisrmf
unix% pset mkacisrmf infile=CALDB
unix% pset mkacisrmf outfile=simple_steps_mkacisrmf.wrmf
unix% pset mkacisrmf energy=0.3:11.0:0.01
unix% pset mkacisrmf channel=1:1024:1
unix% pset mkacisrmf chantype=PI
unix% pset mkacisrmf wmap="simple_steps.pi[WMAP]"

unix% mkacisrmf
scatter/rsp matrix file (CALDB): 
RMF output file (simple_steps_mkacisrmf.wrmf): 
WMAP file (simple_steps.pi[WMAP]): 
energy grid in keV (lo:hi:bin) (0.3:11.0:0.01): 
channel grids in pixel (min:max:bin) (1:1024:1): 
channel type (PI|PHA) (PI): 
gain file (CALDB): 

Total 18 regions to be processed:
    1> reg# 1263  processed
    2> reg# 1264  processed
    3> reg# 1265  processed
    4> reg# 1266  processed
    5> reg# 1294  processed
    6> reg# 1295  processed
    7> reg# 1296  processed
    8> reg# 1297  processed
    9> reg# 1298  processed
   10> reg# 1326  processed
   11> reg# 1327  processed
   12> reg# 1328  processed
   13> reg# 1329  processed
   14> reg# 1330  processed
   15> reg# 1359  processed
   16> reg# 1360  processed
   17> reg# 1361  processed
   18> reg# 1391  processed

For the background, repeat the command with the background spectrum: 

unix% pset mkacisrmf outfile=simple_steps_bkg_mkacisrmf.wrmf
unix% pset mkacisrmf wmap="simple_steps_bkg.pi[WMAP]"

unix% mkacisrmf
scatter/rsp matrix file (CALDB): 
RMF output file (simple_steps_mkacisrmf.wrmf): 
WMAP file (simple_steps_bkg.pi[WMAP]): 
energy grid in keV (lo:hi:bin) (0.3:11.0:0.01): 
channel grids in pixel (min:max:bin) (1:1024:1): 
channel type (PI|PHA) (PI): 
gain file (CALDB): 

Total 6 regions to be processed:
    1> reg# 1393  processed
    2> reg# 1394  processed
    3> reg# 1424  processed
    4> reg# 1425  processed
    5> reg# 1426  processed
    6> reg# 1457  processed

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

You can now skip to the Update the Spectrum Files section.

Using mkrmf (mkrmf)

The mkrmf tool uses the weightfile from mkwarf to create the weighted RMF. The tool queries the CALDB for the correct FEF file and the energy axis grid is ignored since the grid is read from the weights file. (Although the axis1 parameter is ignored, a syntactically-correct value must be entered; this example uses axis1="energy=0:1".)

For the source: 

unix% punlearn mkrmf
unix% pset mkrmf infile=CALDB
unix% pset mkrmf outfile=simple_steps_mkrmf.wrmf
unix% pset mkrmf axis1="energy=0:1"
unix% pset mkrmf axis2="pi=1:1024:1"
unix% pset mkrmf weights=simple_steps.wfef
unix% mkrmf
name of FEF input file (CALDB): 
name of RMF output file (simple_steps_mkrmf.wrmf): 
axis-1(name=lo:hi:btype) (energy=0:1): 
axis-2(name=lo:hi:btype) (pi=1:1024:1):

and for the background: 

unix% pset mkrmf outfile=simple_steps_bkg_mkrmf.wrmf
unix% pset mkrmf weights=simple_steps_bkg.wfef

unix% mkrmf
name of FEF input file (CALDB): 
name of RMF output file (simple_steps_bkg_mkrmf.wrmf): 
axis-1(name=lo:hi:btype) (energy=0:1): 
axis-2(name=lo:hi:btype) (pi=1:1024:1):

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

4. Update the Spectrum Files

Group the Source Spectrum (dmgroup)

The source spectrum is grouped to have a minimum number of 15 counts per new channel: 

unix% punlearn dmgroup
unix% pset dmgroup infile=simple_steps.pi
unix% pset dmgroup outfile=simple_steps_grp.pi
unix% pset dmgroup grouptype=NUM_CTS 
unix% pset dmgroup grouptypeval=15
unix% pset dmgroup xcolumn=channel
unix% pset dmgroup ycolumn=counts
unix% dmgroup
Input dataset name (simple_steps.pi): 
Output dataset name (simple_steps_grp.pi): 
Grouping type (NONE|BIN|SNR|NUM_BINS|NUM_CTS|ADAPTIVE|ADAPTIVE_SNR|BIN_WIDTH|MIN_SLOPE|MAX_SLOPE|BIN_FILE) (NUM_CTS): 
Grouping type value (15): 
Binning specification (): 
Name of x-axis (channel): 
Name of y-axis (counts):

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

Update File Headers (dmhedit)

Finally, update the appropriate header keywords in the both the ungrouped and grouped source files.

Note that the ungrouped background file name is used as the BACKFILE header keyword value. The source grouping is applied to the background grouping when fitting in Sherpa. For fitting only background data, or simultaneous fitting of source and background data, Sherpa can group background dynamically with the group functions; see the Fitting section of this thread for more information. 

unix% punlearn dmhedit
unix% dmhedit infile=simple_steps.pi filelist="" operation=add key=BACKFILE value=simple_steps_bkg.pi
unix% dmhedit infile=simple_steps.pi filelist="" operation=add key=RESPFILE value=simple_steps_mkacisrmf.wrmf
unix% dmhedit infile=simple_steps.pi filelist="" operation=add key=ANCRFILE value=simple_steps_corr.arf

unix% dmhedit infile=simple_steps_grp.pi filelist="" operation=add key=BACKFILE value=simple_steps_bkg.pi
unix% dmhedit infile=simple_steps_grp.pi filelist="" operation=add key=RESPFILE value=simple_steps_mkacisrmf.wrmf
unix% dmhedit infile=simple_steps_grp.pi filelist="" operation=add key=ANCRFILE value=simple_steps.arf

in the ungrouped background file: 

unix% dmhedit infile=simple_steps_bkg.pi filelist="" operation=add key=RESPFILE value=simple_steps_bkg_mkacisrmf.wrmf
unix% dmhedit infile=simple_steps_bkg.pi filelist="" operation=add key=ANCRFILE value=simple_steps_bkg.arf

Extracting Multiple Spectra

In this example, we show how specextract can create multiple output spectra from a single run of the script.

Build Source Regions and Background Regions

This example uses two observations of G21.5-09. Both observations will be processed by specextract at the same time, producing two sets of output files; this is explained further in the Run specextract section.

We define have defined a source and background region for each observation: 

unix% cat 1842_src.reg
# Region file format: CIAO version 1.0
circle(2249.5,4221.5,102.0092)

unix% cat 1842_bg.reg
# Region file format: CIAO version 1.0
circle(2565.5,4129.5,40)

unix% cat 1843_src.reg
# Region file format: CIAO version 1.0
circle(1635.5,4113.5,135.11408)

unix% cat 1843_bg.reg
# Region file format: CIAO version 1.0
circle(2129.5,4007.5,40)

The regions are shown displayed on the event files in Figure 2; ObsID 1842 is in the top frame ObsID 1843 is in the bottom frame. The source region is a dotted white ellipse and the background region is a solid green ellipse in each frame.

Figure 2: Extraction regions on the event files

[Thumbnail image: The two event files used for spectral extraction are displayed in ds9 with the source and background regions overlaid.]

[Version: full-size]

[Print media version: The two event files used for spectral extraction are displayed in ds9 with the source and background regions overlaid.]

Figure 2: Extraction regions on the event files

ObsID 1842 is displayed in the top frame and ObsID 1843 is displayed in the bottom frame.

Run specextract

The event files are input to the script as a stack; this syntax means that a separate spectrum will be created for each of the file/region pairs: 

unix% cat multi_src.lis
acisf01842N002_evt2.fits[sky=region(1842_src.reg)]
acisf01843N002_evt2.fits[sky=region(1843_src.reg)]

When working with stack inputs to specextract, the source and background stacks must contain the same number of items, unless you are not extracting background spectra. Make sure that the background file definitions are in the same order as the source files: 

unix% cat multi_bg.lis
acisf01842N002_evt2.fits[sky=region(1842_bg.reg)]
acisf01843N002_evt2.fits[sky=region(1843_bg.reg)]

When applying the ACIS dead area correction, the same number of parameter block files is also required: 

unix% cat multi_pbk.lis
acisf084281294N002_pbk0.fits
acisf084272477N002_pbk0.fits

We also need to create a stack of output root names, maskfiles, aspect solution, and bad pixel files. These could be supplied to the parameters as a comma-separate list instead: 

unix% cat multi_out.lis
1842
1843

unix% cat multi_mask.lis
acisf01842_000N002_msk1.fits
acisf01843_000N002_msk1.fits

unix% cat multi_asp.lis
pcadf084280882N002_asol1.fits
pcadf084271087N002_asol1.fits

unix% cat multi_bpix.lis
acisf01842_000N002_bpix1.fits
acisf01843_000N002_bpix1.fits

If you prefer, you may just give a string as the outroot and specextract will create output files designated as "src1", "src2", "bkg1", "bkg2", etc. 

unix% pset specextract outroot=multi

Set the parameters: 

unix% punlearn specextract
unix% pset specextract infile=@multi_src.lis
unix% pset specextract outroot=@multi_out.lis
unix% pset specextract bkgfile=@multi_bg.lis
unix% pset specextract asp=@multi_asp.lis
unix% pset specextract pbkfile=@multi_pbk.lis
unix% pset specextract mskfile=@multi_mask.lis
unix% pset specextract badpixfile=@multi_bpix.lis
unix% pset specextract grouptype=BIN binspec=10
unix% pset specextract weight=yes correct=no

The source spectra will be grouped into bins of 10 channels each.

Note that this method allows you to input as many event file/region file pairs as you like, but the same grouping will be applied to all of them. The tool is run with verbose=0 for no screen output: 

unix% specextract 
Source event file(s) (@multi_src.lis): 
Output directory path + root name for output files (@multi_out_test.lis): 
Should response files be weighted? (yes): 
Apply point source aperture correction to ARF? (no): 
Combine ungrouped output spectra and responses? (no): 
Background event file(s) (@multi_bg.lis): 
Create background ARF and RMF? (yes): 
Source aspect solution or histogram file(s) (@multi_asp.lis): 
pbkfile input to mkwarf (@multi_pbk.lis): 
mskfile input to mkwarf (@multi_mask.lis):

The contents of the parameter file may be checked with plist specextract.

Examining the Output Files

The number of files created depends on if a background event file was provided and if source and/or background grouping was specified. In this case, the output files are: 

Source 1 (1842_src.reg):
1842.pi         ungrouped spectrum
1842_tdet.fits  WMAP from sky2tdet
1842.warf       weighted ARF
1842.wfef       FEF weight file
1842.wrmf       weighted RMF
1842_grp.pi     grouped spectrum

Source 2 (1843_src.reg):
1843.pi
1843_tdet.fits
1843.warf
1843.wfef
1843.wrmf
1843_grp.pi

Background 1 (1842_bg.reg):
1842_bkg.pi
1842_bkg_tdet.fits
1842_bkg.warf
1842_bkg.wfef
1842_bkg.wrmf

Background 2 (1843_bg.reg):
1843_bkg.pi
1843_bkg_tdet.fits
1843_bkg.warf
1843_bkg.wfef
1843_bkg.wrmf

If one string is provided for the outroot parameter, the script simply numbers the output to match the order in which the files were input: src1 = 1842_src.reg and src2 = 1843_src.reg. Similarly, the background files are bkg1 and bkg2.

Fitting

If you would like to fit the background-subtracted source spectrum using a common RMF and ARF for source and background, simply read the source spectrum FITS file into Sherpa, subtract the background, and fit. See the Introduction to Fitting PHA Spectra thread for details.

To fit source and background spectra simultaneously with distinct RMFs and ARFs, follow the Independent Background Responses thread.

The RESPFILE and ANCRFILE header keywords in the source and background spectra have been updated, as well as the BACKFILE in the source spectra. This means that when the spectra are read into Sherpa, all the supporting files will automatically be read as well; the background (if available) will be defined, as will the source and background response files.

Analysis Caveats

Analysis at the edges of ACIS CCDs

Users should be cautious about analyzing the data for sources near the edges of the ACIS CCDs.

  1. For X-rays passing through the mirrors, the very bottom of each CCD is obscured by the frame store. As a result, some of the events in rows with CHIPY <= 8 are not detected. (The set of rows affected varies from CCD to CCD.) Since the CIAO tools do not compensate for this effect, the ARFs and exposure maps for sources in these regions may be inaccurate.

  2. For sources within about thirty-two pixels of any edge of a CCD, the source may be dithered off the CCD during part of an observation. The aspect histogram, which is used to create ARFs and exposure maps, is designed to compensate for this effect.

  3. An ARF calculated at the edge of a chip will not be accurate since the response tools for spectral extraction (specifically the ARF) assume that 100% of the PSF is enclosed - i.e. on the chip - all the time, which may not be the case. The amount of error introduced depends on how close the source is to the edge, the morphology of the source, and the characteristics of the PSF, which depends on the source spectrum.

  4. A contaminant has accumulated on the optical-blocking filters of the ACIS detectors, as described in the ACIS QE Contamination why topic. Since there is a gradient in the temperature across the filters (the edges are colder), there is a gradient in the amount of material on the filters. (The contaminant is thicker at the edges.) Within about 100 pixels of the outer edges of the ACIS-I and ACIS-S arrays, the gradient is relatively steep. Therefore, the effective low-energy (' 1 keV) detection efficiency may vary within the dither pattern in this region. The ARF and instrument map tools are designed to read a calibration file which describes this spatial dependence. 


Parameters for /home/username/cxcds_param/specextract.par


        infile = acisf00869N004_evt2.fits[sky=region(simple.reg)] Source event file(s)
       outroot = simple           Output directory path + root name for output files
        weight = yes              Should response files be weighted?
       correct = no               Apply point source aperture correction to ARF?
       bkgfile = acisf00869N004_evt2.fits[sky=region(simple_bkg.reg)] Background event file(s)
       bkgresp = yes              Create background ARF and RMF?
           asp = @pcad_asol1.lis  Source aspect solution or histogram file(s)
       combine = no               Combine ungrouped output spectra and responses?
       pbkfile = acisf078247287N003_pbk0.fits pbkfile input to mkwarf
       mskfile = acisf00869_000N003_msk1.fits mskfile input to mkwarf
      (rmffile = CALDB)           rmffile input for CALDB
        (ptype = PI)              PI or PHA
    (grouptype = NUM_CTS)         Spectrum grouping type (same as grouptype in dmgroup)
      (binspec = 15)              Spectrum grouping specification (NONE,1:1024:10,etc)
(bkg_grouptype = NONE)            Background spectrum grouping type (NONE, BIN, SNR, NUM_BINS, NUM_CTS, or ADAPTIVE)
  (bkg_binspec = )                Background spectrum grouping specification (NONE,10,etc)
       (energy = 0.3:11.0:0.01)   Energy grid
      (channel = 1:1024:1)        RMF binning attributes
  (energy_wmap = 300:2000)        Energy range for (dmextract) WMAP input to mkacisrmf
      (binwmap = tdet=8)          Binning factor for (dmextract) WMAP input to mkacisrmf
   (binarfwmap = 1)               Binning factor for (sky2tdet) WMAP input to mkwarf
       (dafile = CALDB)           dafile input to mkwarf
   (badpixfile = acisf00869_000N003_bpix1.fits) Bad pixel file for the observation
      (clobber = no)              OK to overwrite existing output file?
      (verbose = 1)               Debug Level(0-5)
         (mode = ql)              
    

Parameters for /home/username/cxcds_param/dmextract.par


        infile = acisf00869N004_evt2.fits[sky=region(simple_bkg.reg)][bin pi] Input event file 
       outfile = simple_steps_bkg.pi Enter output file name
          (bkg = )                Background region file or fixed background (counts/pixel/s) subtraction
        (error = gaussian)        Method for error determination(gaussian|gehrels|<variance file>)
     (bkgerror = gaussian)        Method for background error determination(gaussian|gehrels|<variance file>)
      (bkgnorm = 1.0)             Background normalization
          (exp = )                Exposure map image file
       (bkgexp = )                Background exposure map image file
      (sys_err = 0)               Fixed systematic error value for SYS_ERR keyword
          (opt = pha1)            Output file type 
     (defaults = ${ASCDS_CALIB}/cxo.mdb -> /soft/ciao/data/cxo.mdb) Instrument defaults file
         (wmap = [energy=300:2000][bin tdet=8]) WMAP filter/binning (e.g. det=8 or default)
      (clobber = no)              OK to overwrite existing output file(s)?
      (verbose = 0)               Verbosity level
         (mode = ql)              
    

Parameters for /home/username/cxcds_param/sky2tdet.par


        infile = acisf00869N004_evt2.fits[sky=region(simple_bkg.reg)][energy=300:2000][bin sky=1] Input image in sky (x,y) coordinates
   asphistfile = simple_steps.asphist Input aspect histogram file
       outfile = simple_steps_bkg_tdet.fits[wmap] Output TDET  WMAP file
          (bin = 1)               Binning factor
      (geompar = geom)            Pixlib geometry file
      (verbose = 0)               Verbosity
      (clobber = no)              Remove existing files?
         (mode = ql)              
    

Parameters for /home/username/cxcds_param/mkwarf.par


        infile = simple_steps_bkg_tdet.fits[wmap] Input detector WMAP
       outfile = simple_steps_bkg.warf Output weighted ARF file
    weightfile = simple_steps_bkg.wfef Output FEF weights
  spectrumfile =                  Input Spectral weighting file (<filename>|NONE)
     egridspec = 0.3:11.0:0.01    Output energy grid [kev]
       pbkfile = acisf078247287N003_pbk0.fits Parameter block file
    (threshold = 0)               Percent threshold cut for FEF regions
      (feffile = CALDB)           FEF file
      (mskfile = acisf00869_000N003_msk1.fits) Mask file
     (asolfile = )                Stack of aspect solution files
       (mirror = HRMA)            ARDLIB Mirror specification
 (detsubsysmod = )                Detector sybsystem modifier
       (dafile = CALDB)           Dead area file
    (ardlibpar = ardlib)          Parameter file for ARDLIB files
      (geompar = geom)            Parameter file for Pixlib Geometry files
      (clobber = no)              Clobber existing outputs
      (verbose = 0)               Tool chatter level
         (mode = ql)              
    

Parameters for /home/username/cxcds_param/mkacisrmf.par


        infile = CALDB            scatter/rsp matrix file
       outfile = simple_steps_mkacisrmf.rmf RMF output file
          wmap = simple_steps_bkg.pi[WMAP] WMAP file
        energy = 0.3:11.0:0.01    energy grid in keV (lo:hi:bin)
       channel = 1:1024:1         channel grids in pixel (min:max:bin)
      chantype = PI               channel type
        ccd_id = 0                filter CCD-ID
         chipx = 0                filter chipx in pixel
         chipy = 0                filter chipy in pixel
          gain = CALDB            gain file
     (asolfile = )                aspect solution file or a stack of asol files
      (obsfile = )wmap -> simple_steps_bkg.pi[WMAP]) obs file
      (logfile = )                log file
      (contlvl = 100)             # contour level
      (geompar = geom)            pixlib geometry parameter file
       (thresh = 1e-06)           low threshold of energy cut-off probability
      (clobber = no)              overwrite existing output file (yes|no)?
      (verbose = 0)               verbosity level (0 = no display)
         (mode = ql)              
    

Parameters for /home/username/cxcds_param/mkrmf.par


        infile = CALDB            name of FEF input file
       outfile = simple_steps_bkg_mkrmf.wrmf name of RMF output file
         axis1 = energy=0:1       axis-1(name=lo:hi:btype)
         axis2 = pi=1:1024:1      axis-2(name=lo:hi:btype)
      (logfile = STDOUT)          name of log file
      (weights = simple_steps_bkg.wfef) name of weight file
       (thresh = 1e-5)            low threshold of energy cut-off probability
       (outfmt = legacy)          RMF output format (legacy|cxc)
      (clobber = no)              overwrite existing output file (yes|no)?
      (verbose = 0)               verbosity level (0 = no display)
        (axis3 = none)            axis-3(name=lo:hi:btype)
        (axis4 = none)            axis-4(name=lo:hi:btype)
        (axis5 = none)            axis-5(name=lo:hi:btype)
         (mode = ql)              
    

Parameters for /home/username/cxcds_param/dmgroup.par


        infile = simple_steps.pi  Input dataset name
       outfile = simple_steps_grp.pi Output dataset name
     grouptype = NUM_CTS          Grouping type
  grouptypeval = 15               Grouping type value
       binspec =                  Binning specification
       xcolumn = channel          Name of x-axis
       ycolumn = counts           Name of y-axis
      (tabspec = )                Tab specification
    (tabcolumn = )                Name of tab column
     (stopspec = )                Stop specification
   (stopcolumn = )                Name of stop column
    (errcolumn = )                Name of error column
      (clobber = no)              Clobber existing output file?
      (verbose = 0)               Verbosity level
    (maxlength = 0)               Maximum size of groups (in channels)
         (mode = ql)              
    

Parameters for /home/username/cxcds_param/specextract.par


        infile = @multi_src.lis   Source event file(s)
       outroot = @multi_out.lis   Output directory path + root name for output files
        weight = yes              Should response files be weighted?
       correct = no               Apply point source aperture correction to ARF?
       bkgfile = @multi_bg.lis    Background event file(s)
       bkgresp = yes              Create background ARF and RMF?
           asp = @multi_asp.lis   Source aspect solution or histogram file(s)
       combine = no               Combine ungrouped output spectra and responses?
       pbkfile = @multi_pbk.lis   pbkfile input to mkwarf
       mskfile = @multi_mask.lis  mskfile input to mkwarf
      (rmffile = CALDB)           rmffile input for CALDB
        (ptype = PI)              PI or PHA
    (grouptype = BIN)             Spectrum grouping type (same as grouptype in dmgroup)
      (binspec = 10)              Spectrum grouping specification (NONE,1:1024:10,etc)
(bkg_grouptype = NONE)            Background spectrum grouping type (NONE, BIN, SNR, NUM_BINS, NUM_CTS, or ADAPTIVE)
  (bkg_binspec = )                Background spectrum grouping specification (NONE,10,etc)
       (energy = 0.3:11.0:0.01)   Energy grid
      (channel = 1:1024:1)        RMF binning attributes
  (energy_wmap = 300:2000)        Energy range for (dmextract) WMAP input to mkacisrmf
      (binwmap = tdet=8)          Binning factor for (dmextract) WMAP input to mkacisrmf
   (binarfwmap = 1)               Binning factor for (sky2tdet) WMAP input to mkwarf
       (dafile = CALDB)           dafile input to mkwarf
   (badpixfile = @multi_bpix.lis) Bad pixel file for the observation
      (clobber = no)              OK to overwrite existing output file?
      (verbose = 1)               Debug Level(0-5)
         (mode = ql)

History

01 Feb 2006 new for CIAO 3.3
15 Feb 2006 created Running mkacisrmf Independently section
31 Mar 2006 specextract use update added to Overview
05 Apr 2006 In light of the specextract usage change, the thread has been rewritten to use extended sources in the examples
14 Apr 2006 added Analysis at the edges of ACIS CCDs caveat
24 May 2006 added new information to Using the ACIS "Blank-Sky" Background Files caveat
14 Jun 2006 corrected link in "Calibration Updates"; clarified information on GRADED mode data
18 Dec 2006 updated for CIAO 3.4: new calibration files in CALDB 3.3.0; Extracting Multiple Spectra section uses a stack of output file roots (new feature in CIAO 3.4); output files in one-output case no longer have "src1" or "bkg1" included in the filename; mkrmf no longer prints messages at verbose=0; CIAO version in warnings
06 Mar 2007 added ACIS dead area correction section
22 Jan 2008 updated for CIAO 4.0: ACIS dead area correction parameters added to the specextract.par file: pbkfile and dafile (dead area correction is turned on by default); new ACIS blank-sky background file in CALDB 3.4.0 eliminate the header keyword issue; available links point to the Sherpa Beta website; removed outdated calibration updates
31 Mar 2008 updated for CALDB 3.4.3: use mkacisrmf for -110 BI chips if TGAIN calibration has been applied
26 Jun 2008 updated Analysis at the edges of ACIS CCDs caveats (aspect information not taken into account by specextract)
04 Feb 2009 updated for CIAO 4.1: images are inline; Sherpa link updated to 4.1 website; screen output changed with CALDB 4; input data must have a CTI_APP keyword
17 Feb 2009 added "for Extended Sources" to the title
12 Jan 2010 updated for CIAO 4.2: specextract uses a CALDB query to decide which RMF tool should be used; calibration update - the ACIS QE contamination model has been upgraded to vN0005.
05 Mar 2010 added additional information to the Choosing a background file section
09 Mar 2010 The ACIS detector is calibrated over the range 0.224004 - 12 keV; choosing values outside this range will result in errors from specextract.
22 Mar 2010 added "Special case: -110 C data on a front-illuminated chip" to the Creating RMFs: mkrmf vs mkacisrmf section
15 Dec 2010 updated for CIAO 4.3: new ACIS contamination calibration file; specextract is part of the scripts package; the script has been updated to run the sky2tdet tool when creating weighted responses, so an aspect histogram must be provided; the mask file is needed as input for mkwarf
22 Dec 2010 specextract bug fixes released in 22 Dec 2010 scripts package update
25 Feb 2011 updated for 25 Feb scripts package release: the asp parameter can take aspect solution files directly (don't have to run asphist make an aspect histogram first); if the bad pixel file is supplied in the badpixfile parameter, the ardlib parameters are set by the script
01 Mar 2011 CALDB 4.4.2 release: fix to the header of the ACIS QE contamination file. Prior to this release, CIAO would fail when trying to look up the contamination model correction for chips ACIS-8 (S4) and ACIS-9 (S5).
04 Apr 2011 updated for 04 Apr scripts package release: acis_fef_lookup script prints the version at verbose > 0. specextract: user is prompted for the bkgfile parameter (previously was hidden); new parameter, bkgresp, determines whether a background ARF and RMF should be created; additional changes are outlined in detail in the script release notes for 04 Apr 2011.
12 Apr 2011 updated for 12 Apr scripts package release: specextract bug fix - regions may be specified on the command line or in a file, i.e. "sky=circle(344,435,10)" vs. "sky=region(src.reg)".
26 Apr 2011 install version 2 of the tools package for CIAO 4.3 to fix the mkrmf bug; updated for 25 Apr scripts package release: specextract runs the remove_extra_time_keywords script to work around a sky2tdet bug
21 Jun 2011 added a "step-by-step" section
06 Jul 2011 specextract bug fixes were released in the 07 July 2011 scripts package; required software updates are listed in Synopsis
06 Sep 2011 a specextract bug fix was released in the 06 Sep 2011 scripts package
06 Oct 2011 a specextract bug fix and a workaround for a known CIAO bug were released in the 06 Oct 2011 scripts package
15 Dec 2011 reviewed for CIAO 4.4: the remove_extra_time_keywords script is no longer needed, since the sky2tdet bug is fixed in CIAO 4.4
16 Feb 2012 a specextract bug fix was released in the 16 Feb 2012 scripts package: the script will not error out if you set outroot to be a stack and have correct=yes and weight=no.
03 Dec 2012 Review for CIAO 4.5; minor formatting
Return to Threads Page: Top | All | Imag Spec

Last modified: 3 Dec 2012
CXC logo

The Chandra X-Ray Center (CXC) is operated for NASA by the Smithsonian Astrophysical Observatory. 60 Garden Street, Cambridge, MA 02138 USA.   Email: cxcweb@head.cfa.harvard.edu Smithsonian Institution, Copyright © 1998-2012. All rights reserved.