Skip to the navigation links
Last modified: 8 Feb 2010
Where are the PDFs?

Correcting Absolute Astrometry with reproject_aspect

CIAO 4.2 Science Threads



Overview

Last Update: 8 Feb 2010 - reviewed for CIAO 4.2: no changes

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.

Read this thread if:

Related Links:




Contents



Get Started

Sample ObsIDs used: 441 (ACIS-I, Chandra Deep Field South); 581 (ACIS-I, Chandra Deep Field South)

File types needed: evt2; asol1

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 "acisf00581N003_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): 

# wrecon (CIAO 4.2): WARNING: PSF size values not valid beyond 20.008001 arcmin off axis

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): 

# wrecon (CIAO 4.2): WARNING: PSF size values not valid beyond 20.008001 arcmin off axis

The warning indicates that some of the source detections are far off-axis; that is not a problem for this thread. 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.]

[Version: full-size]

[Print media version: The two event files are displayed side-by-side in ds9 with the source lists overlaid in green.]

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:

  1. 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.

  2. 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
pcadf056322870N003_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): 
...
33 common sources found between: 
441_src.fits
581_src.fits
After deleting poor matches, 32 sources remain
...
num_asolInFiles: 1      num_asolOutFiles: 1
Transform scale_factor: 0.998759
Transform rotation angle (deg): 0.039587
Transform x translation (pixels): -0.264562
Transform y translation (pixels): 0.179357

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 pcadf056322870N003_asol1.fits keys=no | less
Infile 1:  056322870N003_aspect.fits
Infile 2:  pcadf056322870N003_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.11954523052  53.11950435554    -4.0875e-05 (-7.69e-05 %)
Values are not equal                   1             dec -27.8055280682518  -27.8055525825571 -2.45143e-05 (+8.82e-05 %)
Values are not equal                   1            roll 47.3214685559945  47.2818816761885   -0.0395869 (-0.0837 %)
Values are not equal                   2              ra 53.1196653515724  53.1196244765415  -4.0875e-05 (-7.69e-05 %)
Values are not equal                   2             dec -27.8057354273435  -27.8057599416841 -2.45143e-05 (+8.82e-05 %)
Values are not equal                   2            roll 47.3190445238645  47.2794576440585   -0.0395869 (-0.0837 %)
Values are not equal                   3              ra 53.1196866080189  53.1196457329912  -4.0875e-05 (-7.69e-05 %)
Values are not equal                   3             dec -27.8057397015102  -27.8057642158571 -2.45143e-05 (+8.82e-05 %)
Values are not equal                   3            roll 47.3190516000757  47.2794647202697   -0.0395869 (-0.0837 %)
Values are not equal                   4              ra 53.1196989041701  53.1196580291443  -4.0875e-05 (-7.69e-05 %)
Values are not equal                   4             dec -27.8057422228761  -27.8057667372266 -2.45144e-05 (+8.82e-05 %)
Values are not equal                   4            roll 47.3190544940053  47.2794676141993   -0.0395869 (-0.0837 %)
Values are not equal                   5              ra 53.1197201606188  53.1196792855961  -4.0875e-05 (-7.69e-05 %)
Values are not equal                   5             dec -27.8057464970396  -27.8057710113964 -2.45144e-05 (+8.82e-05 %)
Values are not equal                   5            roll 47.3190615702197  47.2794746904137   -0.0395869 (-0.0837 %)
...

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_asepct 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 acisf00581N003_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 acisf00581N003_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 acisf00581N003_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.2): ERROR: Cannot find at least 3 source matches between reference and duplicate source files.
# wcs_update (CIAO 4.2): ERROR: Cannot read transform data from file.
# wcs_update (CIAO 4.2): ERROR: Could not close wcsUpdBlock.
# wcs_update (CIAO 4.2): 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


        infile = 441_img.fits     Input file name
       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 = 2.0 4.0          wavelet scales (pixels)
      (regfile = 441_src.reg)     ASCII regions output file
      (clobber = no)              Overwrite existing outputs?
     (ellsigma = 3.0)             Size of output source ellipses (in sigmas)
     (interdir = .)               Directory for intermediate outputs
     (bkginput = )                Input background file name
  (bkgerrinput = no)              Use bkginput[2] for background error
  (outputinfix = )                Output filename infix
    (sigthresh = 1e-06)           Threshold significance for output source pixel list
 (bkgsigthresh = 0.001)           Threshold significance when estimating bkgd only
      (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
      (bkgtime = 0)               Exposure time for input background file
      (maxiter = 2)               Maximum number of source-cleansing iterations
     (iterstop = 0.0001)          Min frac of pix that must be cleansed to continue
      (xoffset = INDEF)           Offset of x axis from optical axis
      (yoffset = INDEF)           Offset of y axis from optical axis
        (eband = 1.4967)          Energy band
      (eenergy = 0.393)           Encircled energy of PSF
     (psftable = ${ASCDS_CALIB}/psfsize20010416.fits -< /soft/ciao/data/psfsize20010416.fits) Table of PSF size data
          (log = no)              Make a log file?
      (verbose = 0)               Log verbosity
         (mode = ql)               
    


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


        infile = 581_img.fits     Input file name
       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 = 2.0 4.0          wavelet scales (pixels)
      (regfile = 581_src.reg)     ASCII regions output file
      (clobber = no)              Overwrite existing outputs?
     (ellsigma = 3.0)             Size of output source ellipses (in sigmas)
     (interdir = .)               Directory for intermediate outputs
     (bkginput = )                Input background file name
  (bkgerrinput = no)              Use bkginput[2] for background error
  (outputinfix = )                Output filename infix
    (sigthresh = 1e-06)           Threshold significance for output source pixel list
 (bkgsigthresh = 0.001)           Threshold significance when estimating bkgd only
      (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
      (bkgtime = 0)               Exposure time for input background file
      (maxiter = 2)               Maximum number of source-cleansing iterations
     (iterstop = 0.0001)          Min frac of pix that must be cleansed to continue
      (xoffset = INDEF)           Offset of x axis from optical axis
      (yoffset = INDEF)           Offset of y axis from optical axis
        (eband = 1.4967)          Energy band
      (eenergy = 0.393)           Encircled energy of PSF
     (psftable = ${ASCDS_CALIB}/psfsize20010416.fits -< /soft/ciao/data/psfsize20010416.fits) Table of PSF size data
          (log = no)              Make a log file?
      (verbose = 0)               Log verbosity
          (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

Return to Threads Page: Top | All | Imag

Where are the PDFs?
Last modified: 8 Feb 2010