12. Advanced wavdetect: Using wtransform & wrecon

12.1 General Discussion

wavdetect has two parts, each of which can be run as a standalone tool:

wtransform performs the convolutions between the data and the various scaled wavelet functions;
wrecon uses the wtransform output products to construct source lists, measure parameters for each detection, and create various maps.

There are a few advantages to running the tools separately.

In wavdetect , the flux scales chosen to define the source properties are computed automatically by selecting which of the scales supplied by the user best matches the size of the point spread function (PSF) at that location on the detector. With wrecon however, the user has control over which of the wtransform wavelet scales should be included. The ability to select flux scales means that you are not tied to optimizing for unresolved structures.
In many instances, one may want to perform the convolution for a set of wavelet scales and then use those results to construct source lists based on different selection criteria. Instead of running wavdetect multiple times, one can run wtransform once, and then perform repeated runs of wrecon with altered conditions.
The Mexican Hat functions in wavdetect are circularly symetric, so there is a single scales parameter. In wtransform , however, the xscales and yscales parameters allows one to specify elliptical wavelet functions.

12.2 Running wtransform

This example illustrates running wavdetect using a Chandra HRC-I dataset of 3C 273 as input; see Chapter 2, Section 2.3 for further information about this dataset.

Herein we run wtransform with 5 scales, using autonaming for most of the output files. We then make two runs with wrecon to demonstrate how the output source properties change.

12.2.1 wtransform Input

The input parameters are set as follows:

unix% punlearn wtransform
unix% pset wtransform infile = \
      "data/hrcf00461N003_evt2.fits[EVENTS][sky=circle(16400,16400,1000)]"
unix% pset wtransform srclist = "."
unix% pset wtransform correl = "."
unix% pset wtransform nbkg = "."
unix% pset wtransform log = yes
unix% pset wtransform verbose = 2
unix% pset wtransform thresh = "."
unix% pset wtransform xscales = "1.0 2.0 4.0 8.0 16.0"
unix% pset wtransform yscales = "1.0 2.0 4.0 8.0 16.0"
unix% pset wtransform clobber = yes

A dot (".'') indicates that autonaming will be used, as explained in Section 1.3.2.

The wtransform tool is then run:

unix% wtransform
input file (data/hrcf00461N003_evt2.fits[EVENTS][sky=circle(16400,16400,1000)]): 
Source list output stack (.): 
Correlation image output stack (.): 
Normalized background image output stack (.): 
Warning: Keyword TLMVER is not in the header.
Warning: Keyword HDUSPEC is not in the header.
Warning: Keyword HDUDOC is not in the header.
Warning: Keyword HDUVERS is not in the header.
Warning: Keyword BTIMNULL is not in the header.
Warning: Keyword BTIMRATE is not in the header.
Warning: Keyword BTIMDRFT is not in the header.
Warning: Keyword BTIMCORR is not in the header.
Warning: Keyword FLSHTIME is not in the header.
Warning: Keyword OBI_NUM is not in the header.
Warning: Keyword READMODE is not in the header.
Warning: Keyword RA_PNT is not in the header.
Warning: Keyword DEC_PNT is not in the header.
Warning: Keyword ROLL_PNT is not in the header.

12.2.2 wtransform Output

A number of output files are created:

unix% ls -1
hrcf00461N003_correl_stk.fits[EVENTS][sky=circle(16400,16400,1000)]
hrcf00461N003_evt2_010_010.fits
hrcf00461N003_evt2_010_010_correl.fits
hrcf00461N003_evt2_010_010_nbkg.fits
hrcf00461N003_evt2_010_010_thresh.fits
hrcf00461N003_evt2_020_020.fits
hrcf00461N003_evt2_020_020_correl.fits
hrcf00461N003_evt2_020_020_nbkg.fits
hrcf00461N003_evt2_020_020_thresh.fits
hrcf00461N003_evt2_040_040.fits
hrcf00461N003_evt2_040_040_correl.fits
hrcf00461N003_evt2_040_040_nbkg.fits
hrcf00461N003_evt2_040_040_thresh.fits
hrcf00461N003_evt2_080_080.fits
hrcf00461N003_evt2_080_080_correl.fits
hrcf00461N003_evt2_080_080_nbkg.fits
hrcf00461N003_evt2_080_080_thresh.fits
hrcf00461N003_evt2_160_160.fits
hrcf00461N003_evt2_160_160_correl.fits
hrcf00461N003_evt2_160_160_nbkg.fits
hrcf00461N003_evt2_160_160_thresh.fits
hrcf00461N003_nbkg_stk.fits[EVENTS][sky=circle(16400,16400,1000)]
hrcf00461N003_src_stk.fits[EVENTS][sky=circle(16400,16400,1000)]
hrcf00461N003_thresh_stk.fits[EVENTS][sky=circle(16400,16400,1000)]
wtransform.log

As can be seen, for each scale size there are a number of output files:

hrcf00461N003_evt2_xscale#_yscale#.fits: the locations of correlation maxima.
hrcf00461N003_evt2_xscale#_yscale#_correl.fits: maps of the convolution of the Mexican Hat wavelet with the data; see Figure 12.1.
hrcf00461N003_evt2_xscale#_yscale#_nbkg.fits: the normalized background maps.
hrcf00461N003_evt2_xscale#_yscale#_thresh.fits: optional output maps which give the S/N available at each location.

The four output "stack" files (ASCII lists of the relevent files) have adopted a filename which includes the DM filter (e.g. hrcf00461N003_correl_stk.fits[EVENTS][sky=circle(16400,16400,1000)]). This is caused by the use of autonaming and is generally undesirable. It can be avoided by not using autonaming, or rename the files after the fact:

unix% mv 'hrcf00461N003_correl_stk.fits[EVENTS][sky=circle(16400,16400,1000)]' \
      hrcf00461N003_correl_stk.fits
unix% mv 'hrcf00461N003_nbkg_stk.fits[EVENTS][sky=circle(16400,16400,1000)]' \
      hrcf00461N003_nbkg_stk.fits
unix% mv 'hrcf00461N003_src_stk.fits[EVENTS][sky=circle(16400,16400,1000)]' \
      hrcf00461N003_src_stk.fits
unix% mv 'hrcf00461N003_thresh_stk.fits[EVENTS][sky=circle(16400,16400,1000)]' \
      hrcf00461N003_thresh_stk.fits

unix% more hrcf00461N003_correl_stk.fits
hrcf00461N003_evt2_010_010_correl.fits
hrcf00461N003_evt2_020_020_correl.fits
hrcf00461N003_evt2_040_040_correl.fits
hrcf00461N003_evt2_080_080_correl.fits
hrcf00461N003_evt2_160_160_correl.fits

**Figure 12.1:** A correlation map produced by *wtransform* for `scale = 4` pixels. The 3C 273 jet stands out. Red areas are positive (or just slightly negative), while the other colors correspond to significantly negative values. The map peaks around 200,000 and has negative regions as low as -25,000.
$\begin{figure}\centering \includegraphics*{plots/wav_recon_runwtransform_output_fig1.ps}\end{figure}$

12.3 Running wrecon

12.3.1 wrecon Input

wrecon is set up to use some of the wtransform output files:

unix% punlearn wrecon
unix% pset wrecon infile = \
      "data/hrcf00461N003_evt2.fits[EVENTS][sky=circle(16400,16400,1000)]"
unix% pset wrecon sourcefile = src1.fits
unix% pset wrecon scellfile = scell1.fits
unix% pset wrecon imagefile = image1.fits
unix% pset wrecon regfile = flux1.reg
unix% pset wrecon clobber = yes
unix% pset wrecon srclist = hrcf00461N003_src_stk.fits
unix% pset wrecon correl = hrcf00461N003_correl_stk.fits
unix% pset wrecon nbkg = hrcf00461N003_nbkg_stk.fits
unix% pset wrecon defnbkgfile = defnbkgfile1.fits
unix% pset wrecon fluxfile = flux1.stk
unix% pset wrecon xscales = "1.0 2 4 8 16"
unix% pset wrecon yscales = "1.0 2 4 8 16"
unix% pset wrecon fluxscales = 1
unix% pset wrecon log = yes
unix% pset wrecon verbose = 2

The wrecon tool is then run:

unix% wrecon
Input file name (data/hrcf00461N003_evt2.fits[EVENTS][sky=circle(16400,16400,1000)]): 
Output source list file name (src1.fits): 
Warning: Keyword TLMVER is not in the header.
Warning: Keyword HDUSPEC is not in the header.
Warning: Keyword HDUDOC is not in the header.
Warning: Keyword HDUVERS is not in the header.
Warning: Keyword BTIMNULL is not in the header.
Warning: Keyword BTIMRATE is not in the header.
Warning: Keyword BTIMDRFT is not in the header.
Warning: Keyword BTIMCORR is not in the header.
Warning: Keyword FLSHTIME is not in the header.
Warning: Keyword OBI_NUM is not in the header.
Warning: Keyword READMODE is not in the header.
Warning: Keyword RA_PNT is not in the header.
Warning: Keyword DEC_PNT is not in the header.
Warning: Keyword ROLL_PNT is not in the header.
Output source cell image file name (scell1.fits): 
Output image file name (image1.fits): 
unix%

Next, we repeat the wrecon run, changing only the output file names and the fluxscale parameter (from 1 to "3 4"). The logfile from the first run (wrecon.log) is also renamed so that it does not get overwritten.

unix% mv wrecon.log wrecon1.log
unix% pset wrecon sourcefile = src2.fits
unix% pset wrecon scellfile = scell2.fits
unix% pset wrecon imagefile = image2.fits
unix% pset wrecon regfile = flux2.reg
unix% pset wrecon defnbkgfile = "defnbkgfile2.fits"
unix% pset wrecon fluxfile = flux2.stk
unix% pset wrecon fluxscales = "3 4"
unix% 
unix% wrecon
...
unix% mv wrecon.log wrecon2.log
unix%

12.3.2 wrecon Output

The first run of wrecon produced:

defnbkgfile1.fits
flux1.reg
flux1.stk
image1.fits
scell1.fits
src1.fits
wrecon1.log

and the second run produced:

defnbkgfile2.fits
flux2.reg
flux2.stk
image2.fits
scell2.fits
src2.fits
wrecon2.log

Descriptions of the output files of wrecon are given in the wavdetect discussions (Chapters 11 and 14).

12.3.3 Comparing Two Runs of wrecon

Both runs made the same detections since they were based on the results from the same run of wtransform . However, the source properties, including location and counts, depend on the size and location of the source cell, so the actual values between the two wrecon runs differ somewhat. We have selected a few columns of interest from the output source lists, and show them here to demonstrate the sort of differences that occur:

POS(X,Y): the pixel coordinates of the centroid of the distribution within the cell.
NPIXSOU: number of pixels in the cell.
NET_COUNTS: net counts within the cell.
R[2]: the major and minor major axes of the count distribution within the cell, multiplied by the ellsigma parameter value.

Source Properties from fluxscale = 1

unix% dmlist "src1.fits[cols POS,NPIXSOU,NET_COUNTS,SRC_SIGNIFICANCE,R]" opt=data
--------------------------------------------------------------------------------
Data for Table Block SRCLIST
--------------------------------------------------------------------------------
ROW  POS(X,Y)           NPIXSOU  NET_COUNTS  SRC_SIGNIFICANCE  R[2]                    

1   (16478.13,16291.77) 806      213089.18   2627.78           [10.09 9.60]
2   (16544.38,16219.29) 170         395.94     72.64           [ 8.87 8.08]
3   (16556.39,16206.0)   67         146.50     39.68           [ 7.34 5.01]
4   (15914.40,16323.60)  16           4.93      2.59           [ 3.10 1.37]
5   (16565.68,16196.76)  37          35.42     12.54           [ 5.04 3.76]
6   (16580.92,16177.14)  33          26.67     10.92           [ 5.38 2.95]
7   (16417.0 ,16691.0)    5           0.97      0.51           [ 0    0]

Source Properties from fluxscale = "3 4"

unix% dmlist "src2.fits[cols POS,NPIXSOU,NET_COUNTS,SRC_SIGNIFICANCE,R]" opt=data
--------------------------------------------------------------------------------
Data for Table Block SRCLIST
--------------------------------------------------------------------------------
ROW  POS(X,Y)           NPIXSOU  NET_COUNTS  SRC_SIGNIFICANCE  R[2]                    

1   (16478.13,16291.76) 1244     213218.03   2405.44           [10.51 10.03]
2   (16544.58,16219.00)  497        498.73     60.22           [12.05 11.10]
3   (16554.74,16206.77)  385        352.63     50.30           [12.63 10.06]
4   (15913.19,16323.80)  342         19.85      8.36           [10.05  8.08]
5   (16567.28,16194.35)  385        150.90     25.21           [14.30 12.05]
6   (16580.52,16178.00)  410        100.61     20.05           [12.16 11.50]
7   (16417.40,16690.40)  103          4.40      2.03           [ 6.31  4.29]

The regions are shown in Figure 12.2.

**Figure 12.2:** The region files produced by two runs of *wrecon* . For the bright core of 3C 273, the regions are essentially the same size, but they differ for the jet detections. The smaller regions (cyan/blue) are for `fluxscales = 1`, and the larger regions (magenta) are for `fluxscales = "3 4"`.
$\begin{figure}\centering \includegraphics*{plots/wav_recon_runwrecon_comparetwo_fig1.ps}\end{figure}$

The cell maps are shown in Figures 12.3 and 12.4. Note how the cell sizes are larger for the jet detections when using larger flux scales.

**Figure 12.3:** The source cells for 3C 273 with `fluxscales = 1`. The map values denote the detection number.
$\begin{figure}\centering \includegraphics*{plots/wav_recon_runwrecon_comparetwo_fig2.ps}\end{figure}$

**Figure 12.4:** The source cells for 3C 273 with `fluxscales = "3 4"`. Compare the cell sizes to those in Figure 12.3.
$\begin{figure}\centering \includegraphics*{plots/wav_recon_runwrecon_comparetwo_fig3.ps} \end{figure}$

cxchelp@head.cfa.harvard.edu