Creating a Mosaic World Coordinate System

Table of Contents

• Introduction
• Input Data
• Creating a WCS From Scratch
• Updating the WCS
• Applying a WCS Database to Other Exposures
• CCMAP Verses MSCTPEAK and Logical Verses Physical Coordinates
• Verifying a WCS

Introduction

This document describe how to create or update a world coordinate system (WCS) for mosaic data using IRAF and the MSCRED package. A WCS is the header description that relates the image pixels to equatorial celestial coordinates. For a mosaic each extension requires its own WCS to account for the relative orientations of the CCDs as well as for the optical distortions from the focal plane to the sky. The unifying aspect of a mosaic WCS is that in all extensions the WCS uses the same reference point on the sky. This is normally the optical axis of the telescope and also the approximate center of the mosaic detector. Figure 1 illustrates a reference point in a mosaic of CCDs. The WCS for different exposures are generally the same (excluding binning and subraster readouts) except for the celestial coordinate tied to the reference point.

While it is possible to derive the WCS without having any header information it helps considerably to be able to use MSCDISPLAY and CCDPROC. In this discussion we assume that the headers contain the information needed to display and process the mosaic exposures. A description of how to set this information is given in Using MSCRED with Generic Mosaic Data.

In principle one would want to derive WCS for each filter. This has been done to some extent for the KPNO Mosaic Imager. The data acquisition system then adds the WCS for the appropriate filter and telescope (this camera is used at two telescopes). However, there is a stage in the mosaic data reduction where the WCS is adjusted to account for zeropoint offsets, small rotations, and refraction effects and scale changes. This is done by registering the WCS to a set of coordinates, either relative to one exposure or from a catalog of star coordinates. This registration provides the small adjustments to match filters to the same coordinate system.

Input Data

To begin you need a calibration exposure of a field with many astrometric stars per CCD. Process the data to subtract bias and trim to just the CCD data, that is trim away any prescan and overscan. The example used here to illustrate deriving and updating a WCS is an exposure from the NOAO CTIO Mosaic . In fact this is the exposure taken early in the commissioning of this mosaic detector used to derive the first WCS.

Astrometric stars are those for which you can obtain accurate right ascension and declination in J2000 coordinates and epoch of observation (i.e. proper motion corrected coordinates). Of course you can use less accurate coordinates and then the WCS will be less accurate. In the example the calibration exposure is of USNO standard field K. The coordinates from the USNO-A2 catalog were obtain from http://www.noao.edu/gateway/catalogs/usno.html . Another source of USNO coordinates is http://cadcwww.dao.nrc.ca/cadcbin/getusno2 . For the KPNO Mosaic star clusters with published astrometric coordinates were used.

The required coordinate list is a simple text file of RA in hours and Dec in degrees. The IRAF routines allow values to be input in sexigesimal notation as well as in decimal format. Only the first two columns are required but the remaining columns, in the example the red and blue magnitudes, are useful for creating subsets of stars such as by magnitude and for identification.

Figure 1 shows the example calibration mosaic exposure displayed with MSCDISPLAY. A subset of USNO stars brighter than a red magnitude of 14th is marked. This is done with MSCTVMARK as discussed in the section Verifying a WCS . This figure illustrates the goal of the WCS derivation. Initially you will have the star coordinate list but you may not have a WCS to allow marking the stars. If you have an approximate or prior WCS you could attempt to mark the stars to check your starting coordinate system.

Figure 1: Example Calibration Field

Creating a WCS From Scratch

This section describes creating an initial, rough, WCS for data which does not have a starting one. It is much easier to update a prior WCS than it is to begin without one. If your data has even a crude WCS or one that needs to be updated, either because the mosaic detector was worked on in a way that might have moved the CCDs or because you want to redo and possibly improve the current WCS, then you can skip this section and go to Updating the WCS.

To get an initial WCS requires you to create a text file for each CCD having lines with pixel coordinates (X and Y) and matching celestial coordinates (RA and Dec). There must be at least three coordinate pairs, usually the pixel position and celestial coordinates of stars, in the image. A finder chart is clearly useful for this. Making the text file can be accomplished in various ways including just using a text file editor.

The tools you need are an image display, such as XIMTOOL, and some way to get the image pixel coordinates using the cursor. This might be as simple as looking at the coordinate readout box in the display, or use tasks such as RIMCURSOR or IMEXAMINE. While centroided positions can be used it is not really necessary for the determining an initial WCS. The initial WCS will be refined later using centroiding and more coordinates. The pixel coordinates can be recorded directly to a file or edited by hand. Similarly the celestial coordinates can be added after the pixel coordinates are determined or as the pixel coordinates are measured.

Here we illustrate using RIMCURSOR . This task is simply a cursor reading loop that writes the keystrokes typed in the display to the standard output, which can be redirected to a file. First display a particular CCD (or single amplifier readout) with DISPLAY (not MSCDISPLAY which is for displaying all the CCDs) using the "fill" option. Remember to specify the image extension. Then run RIMCURSOR with the output directed to a file. When you have matched an object in the image with an RA and Dec in your finding chart and list, position the cursor on the object and type colon followed by the RA and Dec.

Figure 2 shows the first extension in the example calibration exposure displayed and the three coordinates identified with RIMCURSOR. The postions marked are shown for illustration though no marks are usually generated by RIMCURSOR. Note that to exit the cursor reading loop use your EOF keystroke (typically either ^D or ^Z).

Figure 2: Using RIMCURSOR to collect coordinates.

ms> display example[im1] 1 fill+
ms> rimcur > example1.dat
:15:23:32 -0:19:07
:15:23:41 -0:16:13
:15:23:54 -0:15:22
[EOF]

The resulting file looks as shown in Figure 3. One could edit the file to remove the "101 :" but in the next step we have the option of specifying which columns are used.

Figure 3: Result of identifying stars with RIMCURSOR.

ms> type example1.dat
300.412 1156.412 101 : 15:23:32 -0:19:07
980.436 1668.428 101 : 15:23:41 -0:16:13
1188.444 2412.452 101 : 15:23:54 -0:15:22

The task that computes a WCS is CCMAP. This task is basically a fitting task that fits a transformation between pixel coordinates and celestial coordinates and records the result in both a database file and in the image header as a WCS. This task could be used directly to make a final WCS if the the input is generated with many centroided stars and celestial coordinates. This is discussed further in CCMAP Verses MSCTPEAK and Logical Verses Physical Coordinates .

Below shows the result of using this task with the example of three stars. Note that we need to change the "lngcol" and "latcol" parameters. We use an output database of "dev$null" since for now we only want to update the header. In the interactive fitting you just need to type 'q'. Later you will do more with this interactive fitting mode. You could also turn off the interactive fitting if desired. Another possiblity for doing multiple fitts in one step is to use a list of images (again given with the extension syntax) and a matching list of files with the coordinates.

Figure 4: Using CCMAP to make the initial WCS with three stars.

ms> ccmap example1.dat dev$null image=example[im1] lngcol=5 latcol=6 update+

Coords File: example1.dat  Image: example[im1]
    Database: dev$null  Solution: example[im1]
Refsystem: j2000  Coordinates: equatorial FK5
    Equinox: J2000.000 Epoch: J2000.00000000 MJD: 51544.50000
Insystem: j2000  Coordinates: equatorial FK5
    Equinox: J2000.000 Epoch: J2000.00000000 MJD: 51544.50000
Coordinate mapping status
    Ra/Dec or Long/Lat fit rms: 1.06E-12  6.32E-13   (arcsec  arcsec)
Coordinate mapping parameters
    Sky projection geometry: tan
    Reference point: 15:23:42.333  -0:16:54.00  (hours  degrees)
    Reference point: 823.098  1745.763  (pixels  pixels)
    X and Y scale: 0.259  0.262  (arcsec/pixel  arcsec/pixel)
    X and Y axis rotation: 270.332  90.829  (degrees  degrees)
Wcs mapping status
    Ra/Dec or Long/Lat wcs rms: 1.06E-12  6.32E-13   (arcsec  arcsec)
Updating image header wcs

Updating the WCS

Starting with an approximately correct WCS what is done is to register a list of celestial coordinates to the objects in the image, fit a WCS of an appropriate order to describe the optical distortions, and update the header and write a WCS database file to be used to apply the WCS to other images. This is done separately for each CCD though the WCS are constrained to have the same celestial reference point as noted earlier. This requires that there are RA and Dec keywords in the header (in hours and degrees respectively) for each extension which are the same and which point to the optical center of the exposure. The keywords may be added if not present with HEDIT .

The task we will use is MSCTPEAK in the MSCFINDER subpackage of MSCRED. This task was derived from an older program, TFINDER , for making plate solutions using the Guide Star Catalog. The MSCFINDER version instead accepts any simple text file of RA and Dec coordinates, understands the multiextension format, and produces a database appropriate for application to other mosaic exposures. There is currently no specific help for this task but it is fairly self-guiding and the help for TFINDER can be used as a reference.

The single list of coordinates includes all objects within the field of the entire mosaic and may contain objects outside the field. The task will use the starting WCS to select those stars in or near each element of the mosaic for display and for you to register to the image. The example list, part of which is shown below, is taken from the USNO catalog. As before only the first two columns are used and other columns may be present for other purposes. In this case the list has been previously sorted by red magnitude though the order does not matter for MSCTPEAK.

Figure 5: Input coordinate list of USNO objects for the example.

ms> head usnoextract
15:23:41.251 -0:16:17.79 8.6 11.6
15:26:00.732 0:17:46.26 8.6 9.4
15:22:58.946 0:05:55.95 10.4 13.9
15:25:54.554 -0:17:02.79 10.4 11.8
15:23:19.939 -0:08:41.78 10.6 12.2
15:25:34.787 -0:01:02.79 10.6 12.3
15:25:40.378 -0:01:11.50 10.6 11.5
15:24:57.694 -0:06:16.80 10.9 12.1
15:23:05.066 0:04:17.01 11.0 12.1
15:24:09.960 0:03:04.16 11.0 11.9
15:24:23.282 0:03:03.20 11.1 12.7
15:22:57.130 -0:20:26.01 11.3 12.8

Figure 6 shows the parameters for MSCTPEAK that will update the WCS for our example. The "extname" parameter can be used to work only on one or a subset of extensions. Otherwise if you abort the task or the task crashes you would need to go through all the extensions again. The "projection" in the example is "tnx" which is a tangent plane projection with distortion polynomials that is understood by IRAF tasks. The orders are for terms up to a power of three. This is needed for the NOAO 4meter telescopes. You should use "tan" if that is sufficient for your telescope distortions otherwise use the lowest representation you can. In future other functions might be used and will incorporate the developing FITS WCS standards. You will notice that there are some parameters you saw with CCMAP. That is because the fitting engine in this task is also CCMAP.

Figure 6: Parameters for the MSCTPEAK task.

    ms> epar msctpeak
			      I R A F  
	       Image Reduction and Analysis Facility
    PACKAGE = mscfinder
       TASK = msctpeak

    images  =     example  List of WCS calibrated Mosaic images
    coordina= usnoextract  List of ra(hr), dec(deg), optional id
    database=  example.db  Database for astrometric fit
    (extname=            ) Extensions
    (epoch  =       2000.) Coordinate epoch
    (update =         yes) Update image header WCS following fit?
    (autocen=          no) Center catalog coords when entering task?
    (boxsize=           9) Centering box fullwidth

    (project=         tnx) Sky projection geometry
    (fitgeom=     general) Fitting geometry
    (functio=  polynomial) Surface type
    (xxorder=           4) Order of xi fit in x
    (xyorder=           4) Order of xi fit in y
    (xxterms=        half) Xi fit cross terms type

    (yxorder=           4) Order of eta fit in x
    (yyorder=           4) Order of eta fit in y
    (yxterms=        half) Eta fit cross terms type?

    (reject =       INDEF) Rejection limit in sigma units

    (interac=         yes) Enter interactive image cursor loop?
    (frame  =           1) Display frame number
    (marker =      circle) Marker type
    (omarker=        plus) Overlay marker type
    (goodcol=        blue) Color of good marker
    (badcolo=         red) Color of bad marker
    (fdimage=            )
    (fdcoord=            )
    (mode   =           q)

Now run the task. You must have an image display tool running. The task will display, with the "fill" option, each image extension in turn. We will only look at the first extension in this example. It will draw circles where the current WCS predicts the objects to be. If you see no circles when you expect to then there is a problem with the WCS. This might be due to not having the CRVALn keywords correct for the field or the initial WCS is not correct. Figure 7 shows the first part of running the example. Figure 8 shows the cursor help and figure 9 shows a zoomed portion of the display.

Figure 7: Running MSCTPEAK.

ms> msctpeak
List of WCS calibrated Mosaic images (example): 
List of ra(hr), dec(deg), optional id (usnoextract): 
Database for astrometric fit (example.db): 

Marking uncentered catalog sources in red...
  all 219 sources will be marked...done

Marking centered catalog sources in blue...
  no sources to mark

Marking objects in green...
  no objects to mark

You will be prompted with the image cursor for input of various commands. The '?' command will page a list of the commands, see figure 8, in the terminal window.

Figure 8:

			     TPEAK Commands
			     --------------

			Cursor Keystroke Summary

		  "good" = centered, "bad" = uncentered

 a All source toggle			 j Center from the current coords
 b Redisplay only bad sources		 k Center with one keystroke
 c Recenter good sources		 l Center with two keystrokes
 d Delete source(s)			 o Overlay the raw coordinates
 f Fit good sources, reposition bad	 r Redisplay good and bad sources
 g Redisplay only good sources		 u Undelete source(s)
 i Initialize to the raw coordinates

 q Exit the task			 ? Get this help

The important keys are 'a', 'f', 'k', 'l', and 'r'. We could attempt to center all objects using one object with the 'a'+'k' command. But because of the rotation and rough nature of the initial WCS only stars in the vicinity of the reference object would be correctly centered. Instead we move around the field typing 'k' on the image of some of the stars from which we will do a first WCS fit. Note that the 'k' key will center the nearest red circle. If there is a possiblity of confusion the 'l' key is used by typing 'l' on the red circle and then 'l' again on the object to be centered.

When moving around the field feel free to zoom the display. Figure 9 below shows two smaller zoomed regions with a few of the stars marked. Notice how the red circles are off in different directions in the two regions because the initial WCS was only approximately correct. You should center a few stars in the corners and the edge, roughly 12 stars. The centered stars are shown in blue. Because the circles can't be erased except by reloading the image the red circles for the centered objects are still visible.

Figure 9: MSCTPEAK display with objects marked.

Once some stars have been centered type 'f' to fit a new WCS. This will bring up the interactive graphical fitting of CCMAP. You can type '?' to get a reminder of the keys that can be used. The initial figure is something like an x,y plot of the positions to show the distribution of the objects to be fit. Most of the work is done by using the 'x', 'y', 'r', and 's' to display residuals as a function of x and y. In these plots you might use 'd' to delete misidentified, miscentered, or stars with large and uncorrected proper motions. After deleting stars type 'f' to redo the fit. For the first pass don't worry too much about the details. We will come back to the fit after adding in the rest of the objects. Type 'q' to go back to the image. Figure 10 shows three graphs from the first pass at fitting a WCS. These are the initial x/y plot and the residual plots produced by the 'x' and 'y' keys.

Figure 10: First pass at fitting the WCS with fewer stars.

Type 'r' to redisplay the image with the circles now at the position of the new WCS. Figure 11 shows this for one of the regions shown previously.

Figure 11: Objects marked after the first pass fit.

You can now center all the objects quickly with 'a'+'k' on a star which has a blue circle, this centers without any initial offset. After this you will find most of the red cicles are blue and correctly centered on an object. In some cases the centering may fail and the circle remains red. Also if you have wrong centering you can delete a source with 'd'. Once most of the stars have been recentered you are ready to get the final WCS solution. Type 'f' again. Figure 12 shows the results with all the USNO stars on this field.

Figure 12: Fitting of the final WCS with CCMAP.

After quitting the final fit and returning to the image displays type 'q' to quit the current extension and go on to the next one. When this is done the image header is updated and a database record is written.

The final database record will look something like what is shown in figure 13. The database record name will be the same as the extension name.

Figure 13: Example CCMAP WCS database.

ms> page example.db
# Mon 12:25:01 02-Aug-99
begin   im1
        xrefmean        1044.292343750001
        yrefmean        2156.800455357142
        lngmean         15.39719017981151
        latmean         -0.2653082589285715
        pixsystem       physical
        coosystem       j2000
        projection      tnx
        lngref          15.407730555556
        latref          -0.03330555555555601
        lngunits        hours
        latunits        degrees
        xpixref         4244.325766925275
        ypixref         4304.248081888623
        geometry        general
        function        polynomial
        xishift         -1141.614713161611
        etashift        -1101.716071777617
        xmag            0.2639284056373001
        ymag            0.2655069305191946
        xrotation       269.9466255882454
        yrotation       90.92663219116762
        wcsxirms        0.1472868023213945
        wcsetarms       0.160997703816063
        xirms           0.1472868023213879
        etarms          0.1609977038160636
        surface1        11
                        3.      3.
                        2.      2.
                        2.      2.
                        0.      0.
                        -23.    -23.
                        2112.   2112.
                        1.      1.
                        4096.   4096.
                        -1141.614713161611      -1101.716071777617
                        -2.458649045090564E-4   0.2639282911184323
                        0.2654722084256873      -0.004293798704277087
        surface2        18
                        3.      3.
                        4.      4.
                        4.      4.
                        2.      2.
                        -23.    -23.
                        2112.   2112.
                        1.      1.
                        4096.   4096.
                        6.977981629782491       4.83341022221393
                        -0.004048419139473676   -0.005838836166964449
                        9.986220347252658E-7    2.009699573532316E-6
                        -9.720500785419318E-11  -1.103085078407953E-10
                        -0.00780283211407008    -0.003276602257239252
                        1.364528047712044E-6    1.287934476820951E-6
                        -1.462311976812420E-10  3.907533256511147E-11
                        2.333215504725742E-6    5.414511921261557E-7
                        9.518409839069242E-12   -1.571748447126591E-10
                        -2.004754169679203E-10  2.479352519503381E-11

Applying a WCS Database to Other Exposures

The end result of creating and updating a WCS is that the calibration image will have the right WCS in its header and a database file will be produced. To apply the WCS database to other images requires two things. First the WCS keywords must be added to the image extensions based on the extension names. Second the position of the reference point in the WCS, the CRVALn keywords, must be set to the right point on the sky. The task MSCSETWCS, whose parameters are shown in figure 14, does this.

Figure 14: Parameters for MSCSETWCS.

                           I R A F  
            Image Reduction and Analysis Facility
PACKAGE = mscred
   TASK = mscsetwcs

images  =     newimage  Mosaic images
database=   example.db  WCS database
(ra     =           ra) Right ascension keyword (hours)
(dec    =          dec) Declination keyword (degrees)
(equinox=      equinox) Epoch keyword (years)
(ra_offs=           0.) RA offset (arcsec)
(dec_off=           0.) Dec offset (arcsec)
(extlist=             )

The first thing the task does is apply the WCS database specified by the "database" parameter. If no database is specified then this task can be used only for the second step of updating the WCS reference coordinate based on the coordinate keywords in the header. This step consists of transforming the keywords specified by the "ra" and "dec" parameters, presumably set by the data acquisition system, to the WCS keywords "crval1" and "crval2". The keywords give the RA and Dec of the observation in hours and degrees. Ideally these fairly obvious keywords will be a pointing put into the data by the data acquisition system based on the telescope pointing. The "equinox" keyword, if specified and found in the header, will precess the header coordinates to the equinox of the J2000 equinox of the WCS. The offset parameters may be used if there is always a known offset between the RA and Dec values in the header and the WCS reference point.

Currently this task and the method of updating the WCS from a database assumes that there is no rotation or only a small rotation between the WCS calibration and later observations. The small rotations, as well as scale and atmospheric refraction effects are removed later as a global adjustment to the WCS using the task MSCCMATCH.

CCMAP Verses MSCTPEAK and Logical Verses Physical Coordinates

MSCTPEAK is a convenient routine. However you can also do everything with other tools which produce the input to CCMAP, the list of accurate pixel and celestial coordinates, and then use CCDMAP directly. The important points to note in the context of a mosaic and MSCRED is that the "pixsystem" should be specified as "physical" and the database WCS names should be just the image extension names. This is needed to allow MSCSETWCS, see Applying a WCS Database to Other Exposures , to update the headers of mosaic exposures which have a different pointing and possibly are untrimmed (have overscan and prescan regions), trimmed to a smaller region, binned, or the result of a subregion readout.

The idea is that the WCS database describes a transformation between "physical" coordinates and celestial coordinates. Physical coordinates are coordinates tied to some fundamental system. For mosaic data this should be the unbinned CCD coordinates. The "logical" coordinate system is the pixel coordinates of the actual image raster. By determining the WCS and specifying the database WCS in physical coordinates, then when the WCS is applied to other exposures first the logical pixels of the exposure are transformed to the same physical system (unbinned CCD pixels) regardless of overscan/prescan regions, binning, and subraster readouts. Another way to think of this is that the WCS database in physical coordinates is a WCS tied to the CCD pixels. The transformation between logical and physical coordinates is specified by the LTV and LTM keywords found in the image headers.

Earlier we recommended that that the calibration exposure be trimmed to the full unbinned CCD data. This means that the headers should all have LTM1_1/LTM2_2 = 1 and LTV1/LTV2 missing or zero. CCDPROC will update or create these keywords based on the CCDSEC and NSUM keywords. However, the database WCS is best applied as early as possible in the raw data. If the data have been binned then LTM1_1/LTM2_2 will be different than one and if there is an overscan on the left of the image, say 64 pixels of overscan, then the raw data should have LTV1=-64. MSCSETWCS will then apply the WCS correctly in this case.

Verifying a WCS

During any operation involving the WCS it is useful to be able to see if the coordinate system is correct. A picture is worth a thosand words and can help you diagnose problems. The tasks MSCZERO and MSCTVMARK are useful for this. You need to have a list of coordinates that you expect to overlay on your field. This might be easily obtained from some catalog such as the USNO-A. Then you display the image either with MSCDISPLAY or MSCZERO. The latter does MSCDISPLAY if the image is not already loaded and enters an interactive cursor loop. Then you would use MSCTVMARK or the 'm' key in MSCZERO which excutes this task. This takes the celestial coordinates, determines in which, if any, of the extensions the star falls based on the WCS and draws marks at the appropriate logical pixel coordinate.

Figure 15 provides an example. This is the same example field and, in this case, there is a large list of USNO stars. However any size list is useful. The two images show the full field and then zoomed on the region we have looked at previously. This tells us the WCS is correct.

Figure 15: Using MSCZERO and MSCTVMARK to mark coordinates.

ms> msczero example
[Type 'm'...]
List of coordinates (usnoextract14): usnoextract
Coordinate type (logical|physical|world) (logical|physical|world) (world): 
Radii of concentric circles (50): 40
Gray level of marks to be drawn (205):