Correcting Absolute Astrometry with reproject_aspect
![[CXC Logo]](../../imgs/cxc-logo.gif){kind=logo}
CIAO 4.4 Science Threads
Overview
Synopsis:
reproject_aspect applies corrections for translation, scale, and rotation to the WCS of a file by comparing two sets of source lists from the same sky region. If three or more sources are found to be a close match in position, the tool calculates the transformation that relates the two files and updates the WCS mapping from SKY(X,Y) to (RA,Dec) either by modifying the aspect solution or by revising the WCS keywords, depending on how the parameters are set.
reproject_aspect is actually a script which runs two tools: wcs_match and wcs_update. These tools may be run individually for slightly more flexibility; see the help files for details.
Purpose:
-
To eliminate the errors in the absolute astrometry between two files before merging.
-
To improve the absolute astrometry of a Chandra dataset using external information (e.g. a position from optical data).
Related Links:
- Overview: Reprojecting Files
-
Correcting Aspect Prior to Merging: a section of the Merging Data from Multiple Imaging Observations thread which shows how to correct the WCS by hand in the case where this thread cannot be used (i.e. "Cannot find at least 3 source matches").
Last Update: 11 Jan 2012 - reviewed for CIAO 4.4: minor change to output from wavdetect (34 common sources found by reproject_aspect instead of 33)
Contents
- Get Started
- Creating Source Lists (wavdetect)
- Running reproject_aspect: Which Method to Use
- Method 1: Create a New Aspect Solution
- Method 2: Update the Input File
- "ERROR: Cannot find at least 3 source matches"
- Parameter files:
- History
- Images
Get Started
Download the sample data: 441 (ACIS-I, Chandra Deep Field South); 581 (ACIS-I, Chandra Deep Field South)
unix% download_chandra_obsid 441,581 evt2,asol
In this thread, we assume that all relevant files are in the same working directory.
Creating Source Lists (wavdetect)
The reproject_aspect tool needs an input source list which is in FITS format and contains RA and Dec columns. The source list can be generated by CIAO, created with another data analysis program, or be obtained from a catalog. The source list from the archive (src2.fits file) is not accurate enough for use with this tool.
This thread creates a new source list by running one of the detect tools; wavdetect is recommended because of its ability to separate closely-spaced point sources.
wavdetect operates on an image, so we begin by binning the event files. When selecting the image region and binning, keep in mind that you need good position determinations from a set of sources well distributed across the field. However, you also need good source position measurements to see the differences, which implies wavdetect run on full-resolution data. Due to the heavy computational load, it is recommended to run wavdetect images no larger than 1024x1024.
To accomodate all these requirements, we create an image of the full field of view, binned by a factor of 2. While not full resolution, this is a small enough binning to yield good source positions:
unix% dmcopy "acisf00441N003_evt2.fits[bin x=2425:6095:2,y=1620:5770:2]" \ 441_img.fits unix% dmcopy "acisf00581N004_evt2.fits[bin x=2425:6095:2,y=1620:5770:2]" \ 581_img.fits
Next, run wavdetect on each image to create a source list.
unix% punlearn wavdetect unix% pset wavdetect infile=441_img.fits unix% pset wavdetect outfile=441_src.fits unix% pset wavdetect scellfile=441_scell.fits unix% pset wavdetect imagefile=441_imgfile.fits unix% pset wavdetect defnbkgfile=441_nbgd.fits unix% pset wavdetect regfile=441_src.reg unix% wavdetect Input file name (441_img.fits): Output source list file name (441_src.fits): Output source cell image file name (441_scell.fits): Output reconstructed image file name (441_imgfile.fits): Output normalized background file name (441_nbgd.fits): unix% punlearn wavdetect unix% pset wavdetect infile=581_img.fits unix% pset wavdetect outfile=581_src.fits unix% pset wavdetect scellfile=581_scell.fits unix% pset wavdetect imagefile=581_imgfile.fits unix% pset wavdetect defnbkgfile=581_nbgd.fits unix% pset wavdetect regfile=581_src.reg unix% wavdetect Input file name (581_img.fits): Output source list file name (581_src.fits): Output source cell image file name (581_scell.fits): Output reconstructed image file name (581_imgfile.fits): Output normalized background file name (581_nbgd.fits):
Figure 1 shows the detections for each observation. The contents of the parameter file may be checked with plist wavdetect for ObsID 441 and ObsID 581.
![[Thumbnail image: The two event files are displayed side-by-side in ds9 with the source lists overlaid in green.]](detect.png)
[Version: full-size]
Figure 1: Source detections from wavdetect
ObsID 441 is in the left frame and ObsID 581 is in the right frame. The two frames are registered so the WCS coordinates are matched.
reproject_aspect requires that the input source lists have RA and Dec. columns, so the FITS format files - 441_src.fits and 581_src.fits - must be used.
Running reproject_aspect: Which Method to Use
There are two output methods for this tool, depending on what type of file is provided in the updfile parameter:
-
Create a new aspect solution: If updfile is an aspect solution file, a new aspect solution file is created with changes applied to the RA, Dec, and roll columns.
-
Update the input file: If updfile contains an image or event file with a WCS associated with it, the WCS elements within that file are updated. No new output file is created.
If you are planning on reprocessing the data, i.e. running acis_process_events or hrc_process_events to apply new calibration, then you should run reproject_aspect on the aspect solution file(s). Then use the new aspect solution when reprocessing the data; see the "Apply the new aspect solution" section for further information. The modified aspect files should then be used for any further analysis that requires the files, such as running asphist when creating ARFs or exposure maps.
Do not mix and match the original and reprojected files. If you create new aspect solution files, you must reprocess the event file with them. You will get incorrect results if you create new aspect solution files and pair them with the archived event data.
In the case where you have created images already, it is easier to update the WCS in those files. Make sure to update all the images with which you are working, though. If you have a separate counts image, background image, and exposure image, apply the same correction to all three files.
Method 1: Create a New Aspect Solution
Run reproject_aspect
The source list for ObsID 441 is used as the reference file, and the WCS information is taken from the image of that observation. The aspect solution file for ObsID 581 is given in the updfile parameter, which means that a new aspect solution file will be created.
If there is more than one aspect solution file for the observation, it is possible to reproject all of them with a single run of reproject_aspect. There must be a one-to-one match between the number of files listed in the updfile and outfile parameters; see ahelp stack for more information on using stack inputs in CIAO. Here we used:
unix% cat pcad_asol1.lis pcadf056322870N004_asol1.fits unix% cat outfiles.lis 056322870N003_aspect.fits unix% punlearn reproject_aspect unix% pset reproject_aspect infile=581_src.fits unix% pset reproject_aspect refsrcfile=441_src.fits unix% pset reproject_aspect updfile=@pcad_asol1.lis unix% pset reproject_aspect outfile=@outfiles.lis unix% pset reproject_aspect wcsfile=441_img.fits
Running the tool with verbose=2 prints information about the transform to the screen:
unix% reproject_aspect verbose=2 Input file with duplicate srcs (581_src.fits): Input file with reference srcs (441_src.fits): Either input asol file, or file with WCS to be updated (@pcad_asol1.lis): Output asol file (outfiles.lis): ... 34 common sources found between: 441_src.fits 581_src.fits After deleting poor matches, 33 sources remain ... num_asolInFiles: 1 num_asolOutFiles: 1 Transform scale_factor: 0.999976 Transform rotation angle (deg): 0.008637 Transform x translation (pixels): 0.441327 Transform y translation (pixels): -0.254642
The tool goes through a few iterations of calculating the transform, then prints the final scale, rotation, and translation elements to the screen. The new aspect file is 056322870N003_aspect.fits. The contents of the parameter file may be checked with plist reproject_aspect.
Compare the original and updated files
The dmdiff tool can be used to compare the new and old aspect solution files:
unix% dmdiff 056322870N003_aspect.fits pcadf056322870N004_asol1.fits keys=no | less Infile 1: 056322870N003_aspect.fits Infile 2: pcadf056322870N004_asol1.fits -------------------------- SUBSPACE VALUE DIFFERENCES -------------------------- Message: Column: Row Value(s): Diff: -------- ------- --- --------- ----- TABLE NAME: aspsol ----------------------- TABLE VALUE DIFFERENCES ----------------------- Message: Row:Cell: Column: Value(s): Diff: -------- --------- ------- --------- ----- Values are not equal 1 ra 53.1193922377255 53.1194604232829 +6.81856e-05 (+0.000128 %) Values are not equal 1 dec -27.805515733692 -27.8054809290395 +3.48047e-05 (-0.000125 %) Values are not equal 1 roll 47.2904982221841 47.2818611562033 -0.00863707 (-0.0183 %) Values are not equal 2 ra 53.1195099069574 53.1195780926213 +6.81857e-05 (+0.000128 %) Values are not equal 2 dec -27.8057456325813 -27.8057108278711 +3.48047e-05 (-0.000125 %) Values are not equal 2 roll 47.2827244138623 47.2740873478816 -0.00863707 (-0.0183 %) Values are not equal 3 ra 53.1195311108065 53.119599296466 +6.81857e-05 (+0.000128 %) Values are not equal 3 dec -27.8057493070577 -27.805714502337 +3.48047e-05 (-0.000125 %) Values are not equal 3 roll 47.2827554404384 47.2741183744577 -0.00863707 (-0.0183 %) Values are not equal 4 ra 53.1195433545475 53.1196115402043 +6.81857e-05 (+0.000128 %) Values are not equal 4 dec -27.8057512279805 -27.8057164232538 +3.48047e-05 (-0.000125 %) Values are not equal 4 roll 47.2827822848197 47.274145218839 -0.00863707 (-0.0183 %) Values are not equal 5 ra 53.1195645584019 53.1196327440542 +6.81857e-05 (+0.000128 %) Values are not equal 5 dec -27.8057549024379 -27.8057200977008 +3.48047e-05 (-0.000125 %) Values are not equal 5 roll 47.2828133114 47.2741762454193 -0.00863707 (-0.0183 %) ...
Apply the new aspect solution
The event file for ObsID 581 now needs to be reprocessed with the new aspect solution in order to apply the updated information. For instructions on how to do so, follow the Create a New Level=2 Event File thread.
Here we compare two versions of the level=2 event file reprocessed with acis_process_events. 581_newaspect_evt2.fits was reprocessed with the new aspect solution file created by reproject_aspect, while the aspect solution from the Archive was applied to 581_oldaspect_evt2.fits.
unix% dmdiff 581_newaspect_evt2.fits 581_oldaspect_evt2.fits keys=no | less Infile 1: 581_newaspect_evt2.fits Infile 2: 581_oldaspect_evt2.fits ... ----------------------- TABLE VALUE DIFFERENCES ----------------------- Message: Row:Cell: Column: Value(s): Diff: -------- --------- ------- --------- ----- Values are not equal 1 x 3334.86254882812 3332.91088867188 -1.95166 (-0.0585 %) Values are not equal 1 y 4630.58349609375 4629.17919921875 -1.4043 (-0.0303 %) Values are not equal 2 x 2872.9736328125 2871.14331054688 -1.83032 (-0.0637 %) Values are not equal 2 y 4010.57666015625 4009.08178710938 -1.49487 (-0.0373 %) Values are not equal 3 x 3480.75708007812 3478.85571289062 -1.90137 (-0.0546 %) Values are not equal 3 y 4373.736328125 4372.36083984375 -1.37549 (-0.0314 %) Values are not equal 4 x 3488.25708007812 3486.392578125 -1.8645 (-0.0535 %) Values are not equal 4 y 4185.38232421875 4184.00830078125 -1.37402 (-0.0328 %)
Remember that the modified aspect files should be used for any further analysis that requires the files, such as running asphist when creating ARFs or exposure maps.
A note on combining the data
reproject_aspect makes its correction relative to the WCS supplied in the infile. The WCS in the infile and refsrcfile may be completely different, e.g. one may be TAN and one could be some other projection.
If the WCSs of the two files are different, and you want to use them as input to dmmerge (or dmimgcalc) they still need to be reprojected to have the same projection (e.g. be on the same tangent plane. To do this, run reproject_events, reproject_image, or reproject_image_grid after reprocessing the event files.
Method 2: Update the Input File
Run reproject_aspect
Again, the source list for ObsID 441 is used as the reference file, and the WCS information is taken from the image of that observation. A copy of the Archive event file for ObsID 581 is given in the updfile parameter, so the WCS information will be directly updated in that file. (If you choose this method, make sure the file has write permission so it can be updated by the tool.)
unix% cp acisf00581N004_evt2.fits 581_wcs_evt2.fits unix% chmod +w 581_wcs_evt2.fits unix% punlearn reproject_aspect unix% pset reproject_aspect infile=581_src.fits unix% pset reproject_aspect refsrcfile=441_src.fits unix% pset reproject_aspect updfile=581_wcs_evt2.fits unix% pset reproject_aspect wcsfile=441_img.fits unix% reproject_aspect Input file with duplicate srcs (581_src.fits): Input file with reference srcs (441_src.fits): Either input asol file, or file with WCS to be updated (581_wcs_evt2.fits): Output asol file ():
The contents of the parameter file may be checked with plist reproject_aspect.
Compare the original and updated files
The cols option of dmlist shows the WCS information for the event file:
unix% dmlist acisf00581N004_evt2.fits cols
...
----------------------------------------------------------------------------
World Coord Transforms for Columns in Table Block EVENTS
----------------------------------------------------------------------------
ColNo Name
...
8: EQPOS(RA ) = (+53.1226)[deg] +TAN[(-0.000136667)* (sky(x)-(+4096.50))]
(DEC) (-27.8061) (+0.000136667) ( (y) (+4096.50))
unix% dmlist 581_wcs_evt2.fits cols
...
----------------------------------------------------------------------------
World Coord Transforms for Columns in Table Block EVENTS
----------------------------------------------------------------------------
ColNo Name
...
8: EQPOS(RA ) = (+53.1227)[deg] +TAN[(-0.000136497)* (sky(x)-(+4096.50))]
(DEC) (-27.8060) (+0.000136497) ( (y) (+4096.50))
Some rounding is introduced when these values are displayed. For greater precision, look at the RA---TAN and DEC--TAN keywords in the raw header of the event files:
unix% dmlist acisf00581N004_evt2.fits header,raw ... Key 534: C *MTYPE7 = EQPOS / DM Keyword: Descriptor name. Key 535: C *MFORM7 = RA,DEC / [deg] Key 536: C *TCTYP11 = RA---TAN / Key 537: R *TCRVL11 = 53.12262703 / Key 538: R *TCRPX11 = 4096.50000000 / Key 539: R *TCDLT11 = -0.000136667 / Key 540: C *TCUNI11 = deg / Key 541: C *TCTYP12 = DEC--TAN / Key 542: R *TCRVL12 = -27.80606296 / Key 543: R *TCRPX12 = 4096.50000000 / Key 544: R *TCDLT12 = 0.000136667 / Key 545: C *TCUNI12 = deg / ... unix% dmlist 581_wcs_evt2.fits header,raw ... Key 535: C *MTYPE7 = EQPOS / DM Keyword: Descriptor name. Key 536: C *MFORM7 = RA,DEC / [deg] Key 537: C *TCTYP11 = RA---TAN / Key 537: R *TCRVL11 = 53.12266791 / Key 538: R *TCRPX11 = 4096.50000000 / Key 539: R *TCDLT11 = -0.000136497 / Key 541: C *TCUNI11 = deg / Key 542: C *TCTYP12 = DEC--TAN / Key 542: R *TCRVL12 = -27.80603845 / Key 543: R *TCRPX12 = 4096.50000000 / Key 544: R *TCDLT12 = 0.000136497 / Key 546: C *TCUNI12 = deg / ...
The keyword number will change depending on if you have ACIS or HRC data, and whether a grating was used.
Using the updated file in analysis
No further processing related to the aspect correction is required. The updated event or image file may now be used in your data analysis.
"ERROR: Cannot find at least 3 source matches"
reproject_aspect needs at least three sources with a close match in position to be present in the infile and refsrcfile. If the input files don't meet this requirement, the tool exits with errors from wcs_match and wcs_update:
# wcs_match (CIAO 4.4): ERROR: Cannot find at least 3 source matches between reference and duplicate source files. # wcs_update (CIAO 4.4): ERROR: Cannot read transform data from file. # wcs_update (CIAO 4.4): ERROR: Could not close wcsUpdBlock. # wcs_update (CIAO 4.4): ERROR: Could not close wcsUpdDataset.
Users that encounter this problem should correct the WCS by hand as shown in the Correcting Aspect Prior to Merging section of the Merging Data from Multiple Imaging Observations thread.
Parameters for /home/username/cxcds_param/wavdetect.par
#
# parameter file for wavdetect
#
#
# input
#
infile = 441_img.fits Input file name
#
# output
#
outfile = 441_src.fits Output source list file name
scellfile = 441_scell.fits Output source cell image file name
imagefile = 441_imgfile.fits Output reconstructed image file name
defnbkgfile = 441_nbgd.fits Output normalized background file name
#
# scales
#
scales = 2.0 4.0 wavelet scales (pixels)
#
# end of wtransform parameters
#
########################################################################
########################################################################
#
# wrecon parameters
#
#
# PSF size parameters
#
psffile = Image of the size of the PSF
(regfile = 441_src.reg) ASCII regions output file
#
# output options
#
(clobber = no) Overwrite existing outputs?
(ellsigma = 3.0) Size of output source ellipses (in sigmas)
(interdir = ${ASCDS_WORK_PATH} -> /tmp) Directory for intermediate outputs
#
#########################################################################
#
# wtransform parameters
#
#
# optional input
#
(bkginput = ) Input background file name
(bkgerrinput = no) Use bkginput[2] for background error
#
# output info
#
(outputinfix = ) Output filename infix
#
# output content options
#
(sigthresh = 1e-06) Threshold significance for output source pixel list
(bkgsigthresh = 0.001) Threshold significance when estimating bkgd only
(falsesrc = -1.0) Allowed number of false sources per image
(sigcalfile = ${ASCDS_CALIB}/wtsimresult.fits -> /soft/ciao-4.4/data/wtsimresult.fits) Significance calibration file
#
# exposure info
#
(exptime = 0) Exposure time (if zero, estimate from map itself
(expfile = ) Exposure map file name (blank=none)
(expthresh = 0.1) Minimum relative exposure needed in pixel to analyze it
#
# background
#
(bkgtime = 0) Exposure time for input background file
#
# iteration info
#
(maxiter = 2) Maximum number of source-cleansing iterations
(iterstop = 0.0001) Min frac of pix that must be cleansed to continue
#
# end of wrecon parameters
#
########################################################################
#
# run log verbosity and content
#
(log = no) Make a log file?
(verbose = 0) Log verbosity
#
# mode
#
(mode = ql)
Parameters for /home/username/cxcds_param/wavdetect.par
#
# parameter file for wavdetect
#
#
# input
#
infile = 581_img.fits Input file name
#
# output
#
outfile = 581_src.fits Output source list file name
scellfile = 581_scell.fits Output source cell image file name
imagefile = 581_imgfile.fits Output reconstructed image file name
defnbkgfile = 581_nbgd.fits Output normalized background file name
#
# scales
#
scales = 2.0 4.0 wavelet scales (pixels)
#
# end of wtransform parameters
#
########################################################################
########################################################################
#
# wrecon parameters
#
#
# PSF size parameters
#
psffile = Image of the size of the PSF
(regfile = 581_src.reg) ASCII regions output file
#
# output options
#
(clobber = no) Overwrite existing outputs?
(ellsigma = 3.0) Size of output source ellipses (in sigmas)
(interdir = ${ASCDS_WORK_PATH} -> /tmp) Directory for intermediate outputs
#
#########################################################################
#
# wtransform parameters
#
#
# optional input
#
(bkginput = ) Input background file name
(bkgerrinput = no) Use bkginput[2] for background error
#
# output info
#
(outputinfix = ) Output filename infix
#
# output content options
#
(sigthresh = 1e-06) Threshold significance for output source pixel list
(bkgsigthresh = 0.001) Threshold significance when estimating bkgd only
(falsesrc = -1.0) Allowed number of false sources per image
(sigcalfile = ${ASCDS_CALIB}/wtsimresult.fits -> /soft/ciao-4.4/data/wtsimresult.fits) Significance calibration file
#
# exposure info
#
(exptime = 0) Exposure time (if zero, estimate from map itself
(expfile = ) Exposure map file name (blank=none)
(expthresh = 0.1) Minimum relative exposure needed in pixel to analyze it
#
# background
#
(bkgtime = 0) Exposure time for input background file
#
# iteration info
#
(maxiter = 2) Maximum number of source-cleansing iterations
(iterstop = 0.0001) Min frac of pix that must be cleansed to continue
#
# end of wrecon parameters
#
########################################################################
#
# run log verbosity and content
#
(log = no) Make a log file?
(verbose = 0) Log verbosity
#
# mode
#
(mode = ql)
Parameters for /home/username/cxcds_param/reproject_aspect.par
infile = 581_src.fits Input file with duplicate srcs
refsrcfile = 441_src.fits Input file with reference srcs
updfile = @pcad_asol1.lis Either input asol file, or file with WCS to be updated
outfile = @outfiles.lis Output asol file
(wcsfile = 441_img.fits) Input file with WCS used in transform
(radius = 12) radius used to match sources (arcsec)
(residlim = 2) src pairs with residuals > residlim are dropped (arcsec)
(residfac = 0) src pairs with residuals > residfac * position error are dropped
(residtype = 0) residfac applies to: (0) each residual, (1) avg residuals
(logfile = STDOUT) debug log file ( STDOUT | stdout | <filename>)
(clobber = no) Overwrite existing output dataset with same name?
(verbose = 0) debug level (0-5)
(mode = ql)
Parameters for /home/username/cxcds_param/reproject_aspect.par
infile = 581_src.fits Input file with duplicate srcs
refsrcfile = 441_src.fits Input file with reference srcs
updfile = 581_wcs_evt2.fits Either input asol file, or file with WCS to be updated
outfile = Output asol file
(wcsfile = 441_img.fits) Input file with WCS used in transform
(radius = 12) radius used to match sources (arcsec)
(residlim = 2) src pairs with residuals > residlim are dropped (arcsec)
(residfac = 0) src pairs with residuals > residfac * position error are dropped
(residtype = 0) residfac applies to: (0) each residual, (1) avg residuals
(logfile = STDOUT) debug log file ( STDOUT | stdout | <filename>)
(clobber = no) Overwrite existing output dataset with same name?
(verbose = 0) debug level (0-5)
(mode = ql)
History
| 24 Apr 2006 | new for CIAO 3.3: original version |
| 01 Dec 2006 | updated for CIAO 3.4: CIAO version in error messages |
| 23 Jan 2008 | updated for CIAO 4.0: kernel parameter removed from wavdetect; filename and screen output updated for reprocessed data (version N003 event file for 441) |
| 23 Jun 2008 | updated image display to place figures inline with text |
| 13 Feb 2009 | updated for CIAO 4.1: minor change to screen output due to CIAO 4.1 change in reproject_aspect |
| 22 May 2009 | added A note on combining the data |
| 22 Jul 2009 | clarified that the source list does not have to be made by wavdetect; a list from another data analysis tool or catalog may be used |
| 08 Feb 2010 | reviewed for CIAO 4.2: no changes |
| 13 Jan 2011 | reviewed for CIAO 4.3: no changes |
| 11 Jan 2012 | reviewed for CIAO 4.4: minor change to output from wavdetect (34 common sources found by reproject_aspect instead of 33) |
