Merging Data from Multiple Imaging Observations
[CIAO 3.4 Science Threads]
OverviewLast 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: |
Contents
- Getting Started
- Merging about a Common Tangent Point
- Incorporating the Aspect Solution
- Parameter files:
- History
- Images
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 |