Making an exposure-corrected mosaic
CIAO 4.16 Science Threads
Overview
Synopsis:
The merge_obs script replaces the merge_all script and lets you easily create an exposure-corrected image of a group of observations. These steps can also be taken manually: the reproject_image tool maps an image in one WCS reference frame to the WCS in another image. If you do not have a reference image available, reproject_image_grid can reproject the image to a user-defined WCS grid instead.
Purpose:
To create an exposure-corrected mosaic image from multiple observations of SN1006.
Related Links:
Last Update: 7 Feb 2022 - Review for CIAO 4.14. No changes.
Contents
- Get Started
- Using merge_obs to automate the process
- Manual reduction
- Improving the Combined Image
- Parameter files:
- History
-
Images
- Figure 1: Broad-band exposure-corrected image of SN1006 created by merge_obs
- Figure 2: Broad-band image and exposure map of SN1006 created by merge_obs
- Figure 3: Three color, exposure-corrected image of SN1006 created by merge_obs
- Figure 4: Images of each observation
- Figure 5: Input images reprojected to the same WCS
- Figure 6: Reprojected exposure maps
- Figure 7: Exposure-corrected image
- Figure 8: Comparing exposure-corrected images with and without the S3 chip
Get Started
Download the sample data: 3838; 4385; 4386; 4387; 4388; 4389; 4390; 4391; 4392; 4393; 4394
unix% download_chandra_obsid 3838,4385,4386,4387,4388,4389,4390,4391,4392,4393,4394
All of these data are ACIS-I observations of SN1006, a supernova remnant.
In this thread we assume that the observations have been reprocessed by chandra_repro and can be found in obsid/repro/; e.g.
unix% chandra_repro 3838,4385,4386,4387,4388,4389,4390,4391,4392,4393,4394 outdir= verbose=1 ...
The verbose parameter is set to 1 to provide some feedback on the reprocessing (the default value of 0 provides no on-screen updates).
Using merge_obs to automate the process
The merge_obs script automates the construction of an exposure-corrected image from a set of observations. It replaces the merge_all script.
To create exposure-corrected images in the broad, soft, medium, and hard bands using a binning of 4, you just need to say the following (and wait a while):
unix% punlearn merge_obs unix% merge_obs "*/repro/*evt*[ccd_id=0:3]" sn1006/ bin=4 bands=broad,csc Running merge_obs Version: 05 November 2021 Verifying 11 observations. Using CSC ACIS broad science energy band. Using CSC ACIS soft science energy band. Using CSC ACIS medium science energy band. Using CSC ACIS hard science energy band. Calculating new tangent point. New tangent point: RA=15h 2m 56.505s Dec=-41d 56' 15.66" Observations to be reprojected: Obsid Obs Date Exp DETNAM SIM_Z FP Sepn PA (ks) (mm) (K) (') (deg) ---------------------------------------------------------------- 1 3838 2003-04-08 20.1 ACIS-01237 -233.587 153.4 14.2 +80 2 4385 2003-04-08 19.8 ACIS-01237 -233.587 153.4 13.6 +125 3 4386 2003-04-08 19.8 ACIS-01237 -233.587 153.4 13.9 +172 4 4387 2003-04-09 19.8 ACIS-01237 -233.587 153.4 14.3 -142 5 4388 2003-04-09 19.8 ACIS-01237 -233.587 153.4 15.1 -101 6 4389 2003-04-09 19.8 ACIS-01237 -233.587 153.6 13.7 -56 7 4390 2003-04-09 19.8 ACIS-01237 -233.587 153.9 13.8 -5 8 4391 2003-04-10 19.9 ACIS-01237 -233.587 153.4 15.8 +36 9 4392 2003-04-10 19.6 ACIS-01237 -233.587 153.4 4.6 +54 10 4393 2003-04-11 19.6 ACIS-01237 -233.587 153.4 4.1 -80 11 4394 2003-04-11 19.6 ACIS-01237 -233.587 153.4 4.5 +177 Running tasks in parallel with 4 processors. Reprojecting 11 event files to a common tangent point. Merging reprojected events files to sn1006/merged_evt.fits Calculating the output grid The merged images will have 1619 by 1579 pixels, a pixel size of 1.968 arcsec, and cover x=844.5:7320.5:4, y=928.5:7244.5:4. Creating the fluxed images for 11 observations in parallel. Creating 4 aspect histograms for obsid 3838 Creating 4 aspect histograms for obsid 4385 Creating 4 aspect histograms for obsid 4386 Creating 4 aspect histograms for obsid 4387 Creating 4 aspect histograms for obsid 4388 Creating 4 aspect histograms for obsid 4389 Creating 4 aspect histograms for obsid 4390 Creating 4 aspect histograms for obsid 4391 Creating 4 aspect histograms for obsid 4392 Creating 4 aspect histograms for obsid 4393 Creating 4 aspect histograms for obsid 4394 Creating 16 instrument maps for obsid 3838 Creating 16 instrument maps for obsid 4385 Creating 16 instrument maps for obsid 4386 Creating 16 instrument maps for obsid 4387 Creating 16 instrument maps for obsid 4388 Creating 16 instrument maps for obsid 4389 Creating 16 instrument maps for obsid 4390 Creating 16 instrument maps for obsid 4391 Creating 16 instrument maps for obsid 4392 Creating 16 instrument maps for obsid 4393 Creating 16 instrument maps for obsid 4394 Creating 16 exposure maps for obsid 3838 Creating 16 exposure maps for obsid 4385 Creating 16 exposure maps for obsid 4386 Creating 16 exposure maps for obsid 4387 Creating 16 exposure maps for obsid 4388 Creating 16 exposure maps for obsid 4389 Creating 16 exposure maps for obsid 4390 Creating 16 exposure maps for obsid 4391 Creating 16 exposure maps for obsid 4392 Creating 16 exposure maps for obsid 4393 Creating 16 exposure maps for obsid 4394 Combining 4 exposure maps for 4 bands (obsid 3838) Combining 4 exposure maps for 4 bands (obsid 4385) Combining 4 exposure maps for 4 bands (obsid 4386) Combining 4 exposure maps for 4 bands (obsid 4387) Combining 4 exposure maps for 4 bands (obsid 4388) Combining 4 exposure maps for 4 bands (obsid 4389) Combining 4 exposure maps for 4 bands (obsid 4390) Combining 4 exposure maps for 4 bands (obsid 4391) Combining 4 exposure maps for 4 bands (obsid 4392) Combining 4 exposure maps for 4 bands (obsid 4393) Combining 4 exposure maps for 4 bands (obsid 4394) Thresholding data for obsid 3838 Thresholding data for obsid 4385 Thresholding data for obsid 4386 Thresholding data for obsid 4387 Thresholding data for obsid 4388 Thresholding data for obsid 4389 Thresholding data for obsid 4390 Thresholding data for obsid 4391 Exposure-correcting 4 images for obsid 3838 Exposure-correcting 4 images for obsid 4385 Exposure-correcting 4 images for obsid 4386 Exposure-correcting 4 images for obsid 4387 Exposure-correcting 4 images for obsid 4388 Exposure-correcting 4 images for obsid 4389 Exposure-correcting 4 images for obsid 4390 Thresholding data for obsid 4392 Thresholding data for obsid 4393 Exposure-correcting 4 images for obsid 4391 Exposure-correcting 4 images for obsid 4392 Thresholding data for obsid 4394 Exposure-correcting 4 images for obsid 4393 Exposure-correcting 4 images for obsid 4394 Combining 11 observations. The following files were created: The co-added clipped counts images are: sn1006/broad_thresh.img sn1006/soft_thresh.img sn1006/medium_thresh.img sn1006/hard_thresh.img The co-added clipped exposure maps are: sn1006/broad_thresh.expmap sn1006/soft_thresh.expmap sn1006/medium_thresh.expmap sn1006/hard_thresh.expmap The combined FOVs are: sn1006/merged.fov The co-added exposure-corrected images are: sn1006/broad_flux.img sn1006/soft_flux.img sn1006/medium_flux.img sn1006/hard_flux.img Warning: the merged event file sn1006/merged_evt.fits should not be used to create ARF/RMF/exposure maps because the RA_NOM keyword varies by 0.6461 (limit is 0.0003) the DEC_NOM keyword varies by 0.4584 (limit is 0.0003)
The contents of the parameter file may be checked with plist merge_obs or using dmhistory:
unix% dmhistory sn1006/broad_flux.img merge_obs expand=no merge_obs infiles="*/repro/*evt*[ccd_id=0:3]" outroot="sn1006/" bands="broad,csc" xygrid="" maxsize="INDEF" binsize="4" asolfiles="" badpixfiles="" maskfiles="" dtffiles="None" refcoord="" units="default" expmapthresh="1.5%" background="default" parallel="yes" nproc="INDEF" tmpdir="/tmp" cleanup="yes" clobber="no" verbose="1"
The script performs all the tasks needed to create the exposure corrected image (Figure 1), including:
- validating each observation (e.g. to check that they are all compatible);
- calculating the position for the reprojected data (this can be over-ridden by setting the refcoord parameter);
- time sort the observations and display a summary;
- automatically find the supporting files for each observation, such as aspect solution, bad-pixel, DTF file for HRC, mask file for ACIS (these files can be explicitly included or ignored if required using the asolfiles, badpixfiles, maskfiles, and dtffiles parameters);
- create the output directory if required;
- reproject each observation to the new reference position;
- merge the reprojected event files (keeping only the common columns);
- run fluximage on each reprojected observation;
- combine the per-observation images and exposure maps, and use these to create the exposure-corrected image (or images, if multiple bands are given);
- and summarize any observational differences that make using the merged event file with CIAO tools and scripts like specextract or mkwarf problematic.
[Version: full-size]
Figure 1: Broad-band exposure-corrected image of SN1006 created by merge_obs
[Version: full-size]
Figure 2: Broad-band image and exposure map of SN1006 created by merge_obs
Since we set bands=broad,csc, we can also create a three-color image of the remnant (see Figure 3).
[Version: full-size]
Figure 3: Three color, exposure-corrected image of SN1006 created by merge_obs
Manual reduction
In this section we run fluximage on each observation and then reproject the images to combine them. This is slightly different to the approach used by merge_obs above, which does the reprojection first, on the event files.
Users should be careful with the images created with reproject_image_grid (below). While the tool does preserve total flux, it does not preserve integer pixel values. Partial pixels with integer counts in the input image are mapped into non-integer, real valued pixels in the output image. This can lead to some aliasing effects especially when the input and output bin sizes are nearly the same.
The approach taken by merge_obs, to dmmerge the individual event files is the best way to combine images.
Create an Image and Exposure Maps for Each Observation (fluximage)
To run this thread, we want a source image and a congruent exposure map for each of the observations. We choose to create both by running the fluximage script. The image and exposure maps may also be created by following the Multiple-Chip ACIS Exposure Map thread.
In this example, we filter out the ACIS-S3 chip (ccd_id=7) from the images. Since the back-illuminated chips (ACIS-S1 and S3) have different background than the front-illuminated chips, eliminating S3 will reduce the chip-to-chip variations in the final image. If you choose to filter out a chip from your image, make sure that doing so does not remove important features from the data.
The image is binned by a factor of 4, and the preset "broad" energy band filter (500-7000 eV) is applied:
unix% punlearn fluximage unix% pset binsize=4 unix% fluximage "3838/repro/*evt*[ccd_id=0:3]" manual/3838 Running fluximage Version: 04 November 2021 Using event file 3838/repro/acisf03838_repro_evt2.fits[ccd_id=0:3] Using CSC ACIS broad science energy band. Aspect solution 3838/repro/pcadf03838_000N001_asol1.fits found. Bad-pixel file 3838/repro/acisf03838_repro_bpix1.fits found. Mask file 3838/repro/acisf03838_000N004_msk1.fits found. The output images will have 740 by 741 pixels, pixel size of 1.968 arcsec, and cover x=2552.5:5512.5:4,y=2604.5:5568.5:4. Running tasks in parallel with 4 processors. Creating 4 aspect histograms for obsid 3838 Creating 4 instrument maps for obsid 3838 Creating 4 exposure maps for obsid 3838 Combining 4 exposure maps for obsid 3838 Thresholding data for obsid 3838 Exposure-correcting image for obsid 3838 The following files were created: The clipped counts image is: manual/3838_broad_thresh.img The observation FOV is: manual/3838.fov The clipped exposure map is: manual/3838_broad_thresh.expmap The exposure-corrected image is: manual/3838_broad_flux.img
Repeat the fluximage command for each of the ObsIDs, changing the infile and outroot:
unix% foreach obsid (4385 4386 4387 4388 4389 4390 4391 4392 4393 4394) foreach? fluximage "${obsid}/repro/*evt*[ccd_id=0:3]" manual/${obsid} foreach? end ...
A number of output files are created. The images that will be combined are:
unix% ls -1 manual/*_broad_thresh.img manual/3838_broad_thresh.img manual/4385_broad_thresh.img manual/4386_broad_thresh.img manual/4387_broad_thresh.img manual/4388_broad_thresh.img manual/4389_broad_thresh.img manual/4390_broad_thresh.img manual/4391_broad_thresh.img manual/4392_broad_thresh.img manual/4393_broad_thresh.img manual/4394_broad_thresh.img
Figure 4 shows the eleven individual images.
[Version: full-size]
Figure 4: Images of each observation
Reproject the Images to a Common Tangent Point (reproject_image_grid)
reproject_image_grid maps an image - or many images - from the current WCS to a specified grid. This tool is used to reproject all the images of SN1006 to a common frame of reference, creating a mosaic of the object from the individual observations.
The filenames are input to reproject_image_grid as a stack file, images.lis:
unix% cd manual unix% ls *_broad_thresh.img > images.lis
The parameters are set to create a 1500 square pixel image, centered at (225.7, -41.9). Each sky pixel will correspond to 2"; read the caveat about defining a pixel size in arcseconds before running reproject_image_grid.
unix% punlearn reproject_image_grid unix% reproject_image_grid @images.lis sn1006.image 1500 1500 225.7 -41.9 0 Pixel size (): 2" BTIMDRFT values are different...FAIL... BTIMNULL values are different...FAIL... BTIMRATE values are different...FAIL... omit - DEC_NOM values different more than 0.000300 omit - DEC_PNT values different more than 0.000300 warning: DS_IDENT has different value...Merged... warning: OBS_ID has different value...Merged... omit - RA_NOM values different more than 0.000300 omit - RA_PNT values different more than 0.000300
The messages are related to how the tool merges the header information in the input files. The merging_rules ahelp file explains the rules and how they affect the output file header.
Figure 5 shows the output image. The contents of the parameter file may be checked with plist reproject_image_grid or with dmhistory (which does not suffer the problems with converting 2" to 2'):
unix% dmhistory sn1006.image reproject_image_grid expand=no reproject_image_grid infile="@images.lis" outfile="sn1006.image" xsize="1500" ysize="1500" xcenter="225.7" ycenter="-41.9" theta="0" pixelsize="2"" projection="tan" resolution="1" method="sum" coord_sys="world" lookupTab="/soft/ciao/data/dmmerge_header_lookup.txt" clobber="no" verbose="0"
[Version: full-size]
Figure 5: Input images reprojected to the same WCS
Combining the Exposure Maps (reproject_image)
In order to exposure-correct the images, we need to combine the exposure map images into a single file that has the same resolution and WCS as the mosaic image. reproject_image maps an image - or many images - in one WCS reference frame to the WCS in a match image. The tool is used here to reproject the exposure maps to match the mosaic image.
unix% ls -1 *broad_thresh.expmap > expmaps.lis unix% head -5 expmaps.lis 3838_broad_thresh.expmap 4385_broad_thresh.expmap 4386_broad_thresh.expmap 4387_broad_thresh.expmap 4388_broad_thresh.expmap
It is important to note that the method parameter is set to average; this is because we want the average exposure value for each pixel in the combined image, not the sum of all of them.
unix% punlearn reproject_image unix% reproject_image @expmaps.lis sn1006.image sn1006.expmap method=average BTIMDRFT values are different...FAIL... BTIMNULL values are different...FAIL... BTIMRATE values are different...FAIL... omit - DEC_NOM values different more than 0.000300 omit - DEC_PNT values different more than 0.000300 warning: DS_IDENT has different value...Merged... warning: OBS_ID has different value...Merged... omit - RA_NOM values different more than 0.000300 omit - RA_PNT values different more than 0.000300
The messages are related to how the tool merges the header information in the input files. The merging_rules ahelp file explains the rules and how they affect the output file header.
The reprojected exposure maps are shown in Figure 6. The contents of the parameter file may be checked with dmhistory:
unix% dmhistory sn1006.expmap reproject_image expand- reproject_image infile="@expmaps.lis" matchfile="sn1006.image" outfile="sn1006.expmap" resolution="1" method="average" coord_sys="world" lookupTab="/soft/ciao/data/dmmerge_header_lookup.txt" clobber="no" verbose="0"
[Version: full-size]
Figure 6: Reprojected exposure maps
Exposure-correct the Image Mosaic (dmimgthresh & dmimgcalc)
The strongly variable exposure near the edge of a dithered field may produce "hot" pixels when divided into an image. While technically proper, these hot pixels can be an eyesore, drawing attention to a noisy, uninteresting portion of the image. The images and exposure maps used above (those containing the label _thresh) have already been "thresholded" to remove these pixels. If you want to do this yourself you can use the dmimgthresh; for example
unix% punlearn dmimgthresh unix% dmimgthresh sn1006.image sn1006_thresh.image expfile=sn1006.expmap cut=1.5% unix% dmimgthresh sn1006.expmap sn1006_thresh.expmap cut=1.5%
which sets those pixels in sn1006.image whose corresponding exposure value (from sn1006.expmap) is less than 1.5% the peak value to 0.0 in sn1006_thresh.image. You may want to adjust these values for your own observation.
We use dmimgcalc to divide the image by the exposure map, creating the exposure-corrected image:
unix% punlearn dmimgcalc unix% dmimgcalc sn1006.image sn1006.expmap sn1006.flux div warning: CONTENT has 1 different values.
The exposure-corrected image is shown in Figure 7.
[Version: full-size]
Figure 7: Exposure-corrected image
Improving the Combined Image
There are several additional methods that may be used to improve the quality of the final image.
Background Subtraction
As discussed in the Create an Image of Each Observation section, we filtered out the ACIS-S3 chip to reduce the chip-to-chip exposure variations. Figure 8 compares this fluxed image to the results if the S3 chip is not removed. This is not an option for all observations, e.g. if there is significant data on the back-illuminated chips. If you must include both front- and back-illuminated chips, try using background subtraction to compensate for the different chip exposures.
[Version: full-size]
Figure 8: Comparing exposure-corrected images with and without the S3 chip
Energy Filtering
This thread uses monochromatic exposure maps created for the full energy range of the data files. To more accurately represent the spectrum of the data you can use narrower energy bands - such as the soft, medium, and hard bands used to create Figure 3 - or use a spectral weights file to weight the exposure map by a source spectrum.
Smoothing
Adaptively smoothing the data with aconvolve will create a nice image for publication; see the ahelp file for examples. Note that care must be taken in the scientific interpretation of the smoothed image as it is highly processed. The An Image of Diffuse Emission thread describes one way of removing point sources before smoothing data.
Parameters for /home/username/cxcds_param/merge_obs.par infiles = */repro/*evt*[ccd_id=0:3] Input events files outroot = sn1006/ Root of output files (bands = default) Energy bands, comma-separated list, min:max:center in keV or ultrasoft, soft, medium, hard, broad, wide, CSC (xygrid = ) xygrid for output or filename (maxsize = INDEF) Maximum image width or height in pixels (binsize = INDEF) Image binning factor (asolfiles = ) Input aspect solutions (badpixfiles = ) Input bad pixel files (maskfiles = ) Input mask files (dtffiles = ) Input dtf files for HRC observations (refcoord = ) Reference coordinates or evt2 file (units = default) Units for the exposure map (expmapthresh = 1.5%) Remove low-exposure regions? '2%' excludes pixels where exposure is < 2% of the maximum (background = default) Method for background removal (HRC-I) (parallel = yes) Run processes in parallel? (nproc = INDEF) Number of processors to use (tmpdir = ${ASCDS_WORK_PATH} -> /tmp) Directory for temporary files (cleanup = yes) Delete intermediary files? (clobber = no) OK to overwrite existing output file? (verbose = 1) Verbosity level (mode = ql)
Parameters for /home/username/cxcds_param/reproject_image_grid.par infile = @images.lis Input image file name outfile = sn1006.image Output file name xsize = 1500 X-Size of output image in image pixels ysize = 1500 Y-Size of output image in image pixels xcenter = 225.7 X-Center of image in coord_sys ycenter = -41.9 Y-Center of image in coord_sys theta = 0 Angle between world coord north and output y axis pixelsize = 2' Pixel size (projection = tan) Projection from world to physical (resolution = 1) Number of point per side to evaluate (method = sum) Average value (coord_sys = world) Coordinate system to match images in (lookupTab = ${ASCDS_CALIB}/dmmerge_header_lookup.txt -> /soft/ciao/data/dmmerge_header_lookup.txt) lookup table (clobber = no) Clobber existing files (verbose = 0) Tool verbosity (mode = ql)
Parameters for /home/username/cxcds_param/dmimgcalc.par infile = sn1006_image_clean.fits Input file #1 infile2 = sn1006_expmap.fits Input file #2 outfile = sn1006_fluxed.fits output file operation = div arithmetic operation (weight = 1) weight for first image (weight2 = 1) weight for second image (lookupTab = ${ASCDS_CALIB}/dmmerge_header_lookup.txt -> /soft/ciao/data/dmmerge_header_lookup.txt) lookup table (clobber = no) delete old output (verbose = 0) output verbosity (mode = ql)
History
24 Jan 2012 | reviewed for CIAO 4.4: use fluximage to create the images and the exposure maps; renamed "Create an Image of Each Observation" section to Create an Image and Exposure Maps for Each Observation |
15 Oct 2012 | The thread has been updated to use the new merge_obs script and to account for updates to fluximage in the 15 October 2012 scripts package release. |
03 Dec 2012 | Review for CIAO 4.5; mkexpmap chatter removed and dmhistory bug now fixed. |
03 Dec 2013 | Review for CIAO 4.6: the pbkfiles parameter has been removed from the merge_obs script. |
22 Dec 2014 | Reviewed for CIAO 4.7. Added a caveat about using reproject_image_grid with integer counts images. |
07 Feb 2022 | Review for CIAO 4.14. No changes. |