Merging Data from Multiple Imaging Observations

[CIAO 3.4 Science Threads]

Overview

Last Update: 16 Apr 2007 - merge_all v3.6: "[subspace -expno]" filter was added to the dmmerge command in merge_all as a workaround for the extra GTIs in merged event files problem; sign of "delta x" changed to negative and reproject_events call corrected (removed match file) in Correcting Aspect Prior to Merging section

Synopsis:

How to merge event lists from observations with similar, but not necessarily identical, target coordinates into a single data file.

Purpose:

To combine one multiple datasets, the event sky coordinates (X,Y in the event list) must be recomputed; this is often referred to as "regridding" or "reprojecting events to a different tangent point." It is necessary because the sky coordinates are computed relative to the nominal R.A. and Dec. for each individual observation; these coordinates will likely differ from observation to observation of the same target.

Read this thread if:

you want to combine ACIS or HRC imaging data from multiple observations, e.g. to search for faint sources or low-surface-brightness extended features. The merged event list should not be used for spectral analysis, since it does not contain sufficient information to generate correct response files. The recommended technique for the spectral analysis case is to generate separate PHA, RMF, and ARF files for each observation and to analyze them simultaneously.

Related Links:

Proceed to the HTML or hardcopy (PDF: A4 | letter) version of the thread.

Getting Started

Sample ObsIDs used: 1842 (ACIS-I, G21.5-09); 1843 (ACIS-I, G21.5-09)

File types needed: evt2; asol1

This thread uses the merge_all script; for information about the script, consult the help file ("ahelp merge_all"). The most recent version of merge_all is v3.6 (11 April 2007):

unix% egrep '(Version|Date)' `which merge_all`
# Version: 3.6
# Date:    11 April 2007

Please check that you are using the most recent version before continuing. If you do not have the script installed or need to update to a newer version, please refer to the Scripts page.

In order to reduce the event lists to manageable size for the examples, they were filtered to include data from chips I0-I3 only:

unix% dmcopy "acisf01842N001_evt2.fits[ccd_id=0:3]" 1842_evt2.fits
unix% dmcopy "acisf01843N001_evt2.fits[ccd_id=0:3]" 1843_evt2.fits

Caveat: extra GTIs in merged event files

When event files that have different ranges of exposure numbers (expno column) are merged, the output file will have more than one GTI block for each chip. This causes incorrect exposure values later in the analysis, as the CIAO tools cannot use multiple GTI blocks for a single chip. Further details are available on the dmmerge bug page.

A "[subspace -expno]" filter was added to the dmmerge command in merge_all v3.6 as a workaround for this problem.

Note that the addition of the subspace filter means that any user who intends to create lightcurves binned on exposure number from the merged output cannot use merge_all, since that information is eliminated from the subspace. (In general, lightcurves are binned on time.)

Merging about a Common Tangent Point

Each event list contains all the World Coordinate System (WCS) information required to transform between sky coordinates and celestial coordinates in its header. If the WCSs are accurate, one may merge event lists simply by converting sky coordinates in one event list to celestial, using that list's WCS, and then converting those celestial coordinates to sky coordinates in the second event list using the second WCS (see the ds9 thread coordinate definitions for more information).

The level=2 event files may be reprojected to a common tangent point in one of two ways:

Using WCS Information to Merge: all the files are reprojected to match one file.
Merging about a Specific R.A. and Dec.: the files are reprojected to a new tangent point that is different than any of the input observations.

Using WCS Information to Merge

The sky coordinates of 1842_evt2.fits are converted to the system of 1843_evt2.fits, then the two lists are merged; all four CCDs (0-3) are included in the final file.

unix% punlearn merge_all
unix% pset merge_all evtfile="1842_evt2.fits,1843_evt2.fits"
unix% pset merge_all refcoord=1843_evt2.fits
unix% pset merge_all chip=0:3
unix% pset merge_all merged=merge.fits
unix% merge_all
Input event file(s) (1842_evt2.fits,1843_evt2.fits): 
Input asol file(s); time-sorted (): 
Input Live Time Correction list files for HRC-I (): 
ACIS CCD ID(s), comma-separated or range; or 'HRC-I' or 'HRC-S' (0:3): 
Reference coordinates or evt2 file for reproject_events (1843_evt2.fits): 
XY grid (e.g., 0.5:8192.5:#1024,0.5:8192.5:#1024) (): 
Energy or spectrum file for mkinstmap (): 
Output merged event file (merge.fits): 
Output merged exposure map (): 
Output exposure-corrected image (): 

Reprojecting 1842_evt2.fits to 1843_evt2.fits

Merging data...

All done!

The contents of the parameter file may be checked using plist merge_all.

Display the merged image:

unix% ds9 merge.fits &

as shown in Figure 1 .

The source photon distribution in the merged event list appears double-peaked, indicating that there is in fact a small (2-3 pixels) systematic offset in the two WCSs. This may be minimized by first correcting the aspect, as described in the Correcting Aspect Prior to Merging section.

Merging about a Specific R.A. and Dec.

It is also possible to reproject the events to a new tangent point that is different than any of the input observations. In this example, both event lists are reprojected to the arbitrarily chosen point (R.A.,Dec) = (278.086,-10.58625) prior to merging. Also, only data from CCDs 1 and 3 are merged.

unix% punlearn merge_all
unix% pset merge_all evtfile="1842_evt2.fits,1843_evt2.fits"
unix% pset merge_all chip="0,1"
unix% pset merge_all refcoord="278.086,-10.58625"
unix% pset merge_all merged="merge_match.fits"
unix% merge_all
Input event file(s) (1842_evt2.fits,1843_evt2.fits): 
Input asol file(s); time-sorted (): 
Input Live Time Correction list files for HRC-I (): 
ACIS CCD ID(s), comma-separated or range; or 'HRC-I' or 'HRC-S' (0,1): 
Reference coordinates or evt2 file for reproject_events (278.086,-10.58625): 
XY grid (e.g., 0.5:8192.5:#1024,0.5:8192.5:#1024) (): 
Energy or spectrum file for mkinstmap (): 
Output merged event file (merge_match.fits): 
Output merged exposure map (): 
Output exposure-corrected image (): 

Reprojecting 1842_evt2.fits to 278.086,-10.58625

Reprojecting 1843_evt2.fits to 278.086,-10.58625

Merging data...

All done!

The contents of the parameter file may be checked using plist merge_all.

Display the merged image:

unix% ds9 merge_match.fits &

Figure 2 shows merge_match.fits on the left and the WCS-merged merge.fits on the right. Again, the double-peaked source photon distribution is visible in the merged event list.

Incorporating the Aspect Solution

The individual aspect solutions may also be used in the reprojection. When aspect solutions are provided, (DETX,DETY) are recalculated from (CHIPX,CHIPY), then (X,Y) are calculated from (DETX,DETY) using the new tangent point.

The following two examples both include the aspect solution in the analysis:

Applying Aspect Solution: use the aspect solution obtained from the Archive in reprojecting the event files.
Correcting Aspect Prior to Merging: shift one aspect solution to make it consistent with the other before reprojecting the files.

Applying Aspect Solution

In this example, the aspect solution from the Archive is used to recalculate the coordinates for each observation before merging them.

Run reproject_events on each data set:

unix% punlearn reproject_events
unix% pset reproject_events infile=1842_evt2.fits 
unix% pset reproject_events outfile=new_1842.fits
unix% pset reproject_events match="278.086 -10.58625"
unix% pset reproject_events aspect=pcadf084280882N001_asol1.fits 
unix% reproject_events
Input dataset/block specification (1842_evt2.fits):
Output dataset/block specification (new_1842.fits): 
Match file (278.086 -10.58625): 

unix% pset reproject_events infile=1843_evt2.fits 
unix% pset reproject_events outfile=new_1843.fits
unix% pset reproject_events match="278.086 -10.58625"
unix% pset reproject_events aspect=pcadf084271087N001_asol1.fits
unix% reproject_events
Input dataset/block specification (1843_evt2.fits): 
Output dataset/block specification (new_1843.fits): 
Match file (278.086 -10.58625):

The contents of the parameter file may be checked using plist reproject_events.

Now merge the reprojected event lists:

unix% punlearn dmmerge
unix% pset dmmerge infile="new_1842.fits,new_1843.fits" 
unix% pset dmmerge outfile=merge_new.fits
unix% dmmerge
Input dataset/block specification (new_1842.fits,new_1843.fits): 
Output dataset name (merge_new.fits): 
warning: DS_IDENT has different value...Merged...
warning: OBS_ID has different value...Merged...
omit - ROLL_NOM values different more than 1.000000
warning: SEQ_NUM has different value...Merged...
omit - SIM_Z values different more than 0.100000

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 contents of the parameter file may be checked using plist dmmerge.

Display the merged image:

unix% ds9 merge_new.fits &

and you will see something like Figure 3 . Since we did not alter the aspect solution from the pipeline, the double-peak source distribution is still evident.

Correcting Aspect Prior to Merging

In order to improve on the results of the previous examples, we must address the problem of the systematic shift in the coordinate systems of the two event lists. To do this, we will shift the aspect solution of one to make it consistent with the other.

If you choose to use this method, the shifted aspect solution should then be used throughout the rest of the analysis in order to obtain correct results.

First we need to determine what shift to apply to the data. This is done by finding the centroid of a bright source that is in both files with dmstat, then find the error in the centroid using the (RA,DEC) coordinates. Figure 4 shows the source region defined for each of these files. The regions are saved as source_1842.reg and source_1843.reg:

unix% cat source_1842.reg 
# Region file format: CIAO version 1.0
circle(4425.125,4098.125,3.937505)

unix% cat source_1843.reg
# Region file format: CIAO version 1.0
circle(3815.2273,3986.1253,4)

The centroid of the image is simply the mean of the event coordinates (EQPOS column):

unix% dmstat "1842_evt2.fits[sky=region(source_1842.reg)][cols eqpos]"
EQPOS(RA, DEC)
    min:        ( 278.08635883 -10.586861679 )        @:        ( 53 67 )
    max:        ( 278.08720587 -10.586258083 )        @:        ( 16 12 )
   mean:        ( 278.08678919 -10.58657894 )
  sigma:        ( 0.0002011502262 0.00013137615646 )
    sum:        ( 23359.290292 -889.27263093 )
   good:        ( 84 84 )
   null:        ( 0 0 )


unix% dmstat "1843_evt2.fits[sky=region(source_1843.reg)][cols eqpos]"
EQPOS(RA, DEC)
    min:        ( 278.08600209 -10.586755728 )        @:        ( 45 43 )
    max:        ( 278.08694789 -10.585734185 )        @:        ( 27 44 )
   mean:        ( 278.08648462 -10.586224615 )
  sigma:        ( 0.00014211719079 0.00019158941942 )
    sum:        ( 33370.378154 -1270.3469538 )
   good:        ( 120 120 )
   null:        ( 0 0 )

Since 1843_evt2.fits is the smaller of the two lists, it will be shifted to match the 1842 observation. To find the shift:

calculate the difference between the two centroids:

delta RA  = 278.08678919-278.08648462   =  0.00030457 deg

delta Dec = -10.58657894-(-10.586224615)= -0.00035433 deg

correct the delta RA by "cos(Dec)" to get the true spherical distance:

-10.586 deg * 2pi rad / 360 deg = -0.184761 rad

0.00030457 deg * cos(-0.184761) =  0.000299 deg = delta RA

convert the values to sky pixels. The delta x conversion is done with -1.3667E-4 deg/pixel because the x value is opposite in sign from RA. Dec and y have the same sign, so 1.3667E-4 deg/pixel is used for that calculation.
```
1 pixel = 0.492" = 1.3667E-4 deg 

delta x =  0.000299 deg / -1.3667E-4 deg/pixel  = -2.19 pixels

delta y = -0.00035433 deg / 1.3667E-4 deg/pixel = -2.59 pixels
```

The tool wcs_update is used to apply the offsets. The shifts are input to the deltax and deltay parameters in units of sky pixels. Additionally, the WCS information from the event file for the match observation (new_1842.fits, created in the previous section) is supplied in the wcsfile parameter. The new aspect solution file will be named 1843_aspect.fits.

unix% punlearn wcs_update
unix% pset wcs_update infile=pcadf084271087N001_asol1.fits
unix% pset wcs_update outfile=1843_aspect.fits
unix% pset wcs_update wcsfile=new_1842.fits
unix% pset wcs_update deltax= -2.19
unix% pset wcs_update deltay= -2.59
unix% wcs_update
Input coordinate transform file (): 
Either input asol file, or file with WCS to be updated (pcadf084271087N001_asol1.fits): 
Output asol file (1843_aspect.fits):

The contents of the parameter file may be checked with plist wcs_update.

The new aspect file is then used to adjust the reprojected event list (new_1843.fits):

unix% punlearn reproject_events
unix% pset reproject_events infile=new_1843.fits
unix% pset reproject_events outfile=newasp_1843.fits
unix% pset reproject_events match=none
unix% pset reproject_events aspect=1843_aspect.fits
unix% reproject_events
Input dataset/block specification (new_1843.fits): 
Output dataset/block specification (newasp_1843.fits): 
Match file (none):

The contents of the parameter file may be checked using plist reproject_events.

Finally, merge newasp_1843.fits and new_1842.fits:

unix% punlearn dmmerge
unix% pset dmmerge infile="new_1842.fits,newasp_1843.fits"
unix% pset dmmerge outfile=merge_asol_match.fits
unix% dmmerge 
Input dataset/block specification (new_1842.fits,newasp_1843.fits): 
Output dataset name (merge_asol_match.fits): 
warning: DS_IDENT has different value...Merged...
warning: OBS_ID has different value...Merged...
omit - ROLL_NOM values different more than 1.000000
warning: SEQ_NUM has different value...Merged...
omit - SIM_Z values different more than 0.100000

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 contents of the parameter file may be checked using plist dmmerge.

Display the merged image:

unix% ds9 merge_asol_match.fits &

As shown in Figure 5 , there is no longer a noticable shift between the events of the two observations since we corrected the aspect solution before combining them. The left frame is the result of correcting the aspect solution before merging the event files (merge_asol_match.fits). The right frame is the merged file from the "Applying the Aspect Solution" example, merge_new.fits.


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


       evtfile = 1842_evt2.fits,1843_evt2.fits Input event file(s)
          asol =                  Input asol file(s); time-sorted
       dtffile =                  Input Live Time Correction list files for HRC-I
          chip = 0:3              ACIS CCD ID(s), comma-separated or range; or 'HRC-I' or 'HRC-S'
      refcoord = 1843_evt2.fits Reference coordinates or evt2 file for reproject_events
        xygrid =                  XY grid (e.g., 0.5:8192.5:#1024,0.5:8192.5:#1024)
        energy =                  Energy or spectrum file for mkinstmap
        merged = merge.fits       Output merged event file
        expmap =                  Output merged exposure map
       expcorr =                  Output exposure-corrected image
       (intdir = /tmp)            Directory for intermediate files
      (clobber = no)              OK to overwrite existing output file(s)?
         (mode = ql)


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


       evtfile = 1842_evt2.fits,1843_evt2.fits Input event file(s)
          asol =                  Input asol file(s); time-sorted
       dtffile =                  Input Live Time Correction list files for HRC-I
          chip = 0,1              ACIS CCD ID(s), comma-separated or range; or 'HRC-I' or 'HRC-S'
      refcoord = 278.086,-10.58625 Reference coordinates or evt2 file for reproject_events
        xygrid =                  XY grid (e.g., 0.5:8192.5:#1024,0.5:8192.5:#1024)
        energy =                  Energy or spectral weights file for mkinstmap
        merged = merge_match.fits Output merged event file
        expmap =                  Output merged exposure map
       expcorr =                  Output exposure-corrected image
       (intdir = /tmp)            Directory for intermediate files
      (clobber = no)              OK to overwrite existing output file(s)?
         (mode = ql)


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


        infile = 1843_evt2.fits Input dataset/block specification
       outfile = new_1843.fits    Output dataset/block specification
         match = 278.086 -10.58625 Match file
       (aspect = pcadf084271087N001_asol1.fits) Aspect file
       (random = -1)              random seed (0 use time, -1 no randomize)
      (geompar = geom)            Parameter file for Pixlib Geometry files
      (verbose = 0)               Debug Level(0-5)
      (clobber = no)              Clobber existing file
         (mode = ql)


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


        infile = new_1842.fits,new_1843.fits Input dataset/block specification
       outfile = merge_new.fits   Output dataset name
     (outBlock = )                Output block name(1st blkname to be duplicated)
    (lookupTab = ${ASCDS_CALIB}/dmmerge_header_lookup.txt -> /soft/ciao/data/dmmerge_header_lookup.txt) lookup table
   (columnList = )                Column list
      (clobber = no)              Clobber existing file[y/n]
      (verbose = 0)               Tool chatter level
         (mode = ql)


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


        infile = pcadf084271087N001_asol1.fits Either input asol file, or file with WCS to be updated
       outfile = 1843_aspect.fits Output asol file
 transformfile =                  Input coordinate transform file
      (wcsfile = new_1842.fits)  Input reference WCS file
      (logfile = STDOUT)          debug log file ( STDOUT | stdout | <filename>)
       (deltax = -2.19)            transform delta_x value (sky pixels)
       (deltay = -2.59)           transform delta_y value (sky pixels)
       (rotang = 0)               transform rotation angle (degrees)
     (scalefac = 1)               transform scale factor
      (clobber = no)              Overwrite existing output dataset with same name?
      (verbose = 0)               debug level (0-5)
         (mode = ql)


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


        infile = new_1843.fits    Input dataset/block specification
       outfile = newasp_1843.fits Output dataset/block specification
         match = none             Match file
       (aspect = 1843_aspect.fits) Aspect file
       (random = -1)              random seed (0 use time, -1 no randomize)
      (geompar = geom)            Parameter file for Pixlib Geometry files
      (verbose = 0)               Debug Level(0-5)
      (clobber = no)              Clobber existing file
         (mode = ql)


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


        infile = new_1842.fits,newasp_1843.fits Input dataset/block specification
       outfile = merge_asol_match.fits Output dataset name
     (outBlock = )                Output block name(1st blkname to be duplicated)
    (lookupTab = ${ASCDS_CALIB}/dmmerge_header_lookup.txt -> /soft/ciao/data/dmmerge_header_lookup.txt) lookup table
   (columnList = )                Column list
      (clobber = no)              Clobber existing file[y/n]
      (verbose = 0)               Tool chatter level
         (mode = ql)

History

04 Jan 2005	reviewed for CIAO 3.2: no changes
20 Dec 2005	updated for CIAO 3.3: the `merge_all` script has been updated to version 3.3 for the new `asphist` tool syntax; minor changes to `dmmerge` parameter file (`outBlock`, `lookupTab`, and `columnList` are hidden parameters; added `verbose` parameter)
06 Jun 2006	Correcting Aspect Prior to Merging section updated to use `wcs_update` instead of `dmtcalc`
01 Dec 2006	updated for CIAO 3.4: version 3.4 of `merge_all` (removed "punlearn mkinstmap" to allow new CIAO 3.4 `mkinstmap` parameters (pbkfile, dafile) to be set before running the script, if desired; parameters supplied on the command line will not be prompted)
26 Feb 2007	merge_all v3.5: added a test to determine whether "energy" value is monochromatic or a spectrum file; temporary filenames have a timestamp value as the root, rather than hardwired names.
03 Mar 2007	added Caveat: extra GTIs in merged event files
16 Apr 2007	merge_all v3.6: "[subspace -expno]" filter was added to the dmmerge command in `merge_all` as a workaround for the extra GTIs in merged event files problem; sign of "delta x" changed to negative and `reproject_events` call corrected (removed match file) in Correcting Aspect Prior to Merging section