Last modified: 25 Jan 2022

URL: https://cxc.cfa.harvard.edu/ciao/threads/match_binning/

Match the Binning of an Image

CIAO 4.16 Science Threads


Overview

Synopsis:

There are times that you would like to create an image so that it matches an already-existing image. One example would be creating a full-resolution exposure map of a single detector chip.

Purpose:

To use the get_sky_limits script to find the binning specification of an image for input into other CIAO tools, e.g. dmcopy or mkexpmap.

Related Links:

Last Update: 25 Jan 2022 - Review for CIAO 4.14. Updated for Repro-5 and CALDB 4.9.6.


Contents


Get Started

Download the sample data: 1838 (ACIS-S, G21.5-09)

unix% download_chandra_obsid 1838 evt2

Create an image

As an example, we select a region of the observation that we are interested in:

unix% dmcopy \
      "acisf01838N004_evt2.fits[sky=box(4182,4392,1172,1184)][bin x=::4,y=::4]" img.fits

What region does the image cover?

The region represented by the image can be found in one of several ways:

A. dmlist opt=cols

The cols option of dmlist produces a report on the various coordinate systems stored for the image. In this case we find:

unix% dmlist img.fits cols
 
--------------------------------------------------------------------------------
Columns for Image Block EVENTS_IMAGE
--------------------------------------------------------------------------------
 
ColNo  Name                 Unit        Type             Range
   1   EVENTS_IMAGE[294,297]              Int2(294x297)  -                    
 
--------------------------------------------------------------------------------
Physical Axis Transforms for Image Block EVENTS_IMAGE
--------------------------------------------------------------------------------
 
Group# Axis# 
   1   1,2    sky(x) = (+3595.50) +(+4.0)* ((#1)-(+0.50))
                 (y)   (+3799.50)  (+4.0)  ((#2) (+0.50))
 
--------------------------------------------------------------------------------
World Coordinate Axis Transforms for Image Block EVENTS_IMAGE
--------------------------------------------------------------------------------
 
Group# Axis# 
   1   1,2    EQPOS(RA ) = (+278.3850)[deg] +TAN[(-0.000136667)* (sky(x)-(+4096.50))]
                   (DEC)   (-10.5884 )           (+0.000136667)  (   (y) (+4096.50)) 

B. dmlist opt=subspace

The subspace option of dmlist describes the filters that have been applied to the data:

unix% dmlist img.fits subspace | less
 
 --- Component 1 --- 
 --- Component 1 --- 
   1 sky                  Real8               Box(4182,4392,1172,1184)
   1 sky                  Real8               Field area = 1.39709e+06 Region area = 1.38765e+06

   1 sky                  [ 1] x                   3596.0:     4768.0 
   1 sky                  [ 2] y                   3800.0:     4984.0 
   2 time                 Real8               TABLE GTI7
                                              
                                               84245787.5549994707: 84253739.5550068617
   3 ccd_id               Int2                7:7 
   4 node_id              Int2                0:3 
   5 expno                Int4                0:2147483647 
   6 chipx                Int2                1:1024 
   7 chipy                Int2                1:1024 
   8 tdetx                Int2                1:8192 
   9 tdety                Int2                1:8192 
  10 detx                 Real4                   0.50:     8192.50 
  11 dety                 Real4                   0.50:     8192.50 
  12 pha                  Int4                0:36855 
  13 pha_ro               Int4                0:36855 
  14 energy               Real4                                  0:  1000000.0 
  15 pi                   Int4                1:1024 
  16 fltgrade             Int2                0:255 
  17 grade                Int2                0:0,2:2,3:3,4:4,6:6 
  18 phas                 Int2                -4096:4095 
... (other components listed) ...

C. dmcoords

The dmcoords tool can be used to calculate the edges of the image in a variety of coordinate systems. Here we use it in the non-interactive mode to find the bottom-left and top-right corners:

unix% punlearn dmcoords

unix% dmcoords img.fits opt=logical logicalx=0.5 logicaly=0.5 mode=h
unix% echo "x = `pget dmcoords x`  y = `pget dmcoords y`"
x = 3595.5  y = 3799.5

unix% punlearn dmcoords
unix% dmcoords img.fits opt=logical logicalx=294.5 logicaly=297.5 mode=h
unix% echo "x = `pget dmcoords x`  y = `pget dmcoords y`"
x = 4771.5  y = 4987.5

The parameter files after running these two commands can be found below: bottom left and top right.


Run the script (get_sky_limits)

The get_sky_limits script returns the text you would use in both a DM-binning and mkexpmap xygrid specification:

unix% get_sky_limits img.fits verbose=1
Running: get_sky_limits
  version: 25 March 2011
Checking binning of image: img.fits
  Image has 294 x 297 pixels
  Pixel size is 4 by 4
  Lower left (0.5,0.5) corner is x,y= 3595.5, 3799.5
  Upper right (294.5,297.5) corner is x,y= 4771.5, 4987.5
  DM filter is:
    x=3595.5:4771.5:#294,y=3799.5:4987.5:#297
  mkexpmap xygrid value is:
    3595.5:4771.5:#294,3799.5:4987.5:#297

The binning specifications are also stored in the parameter file in the dmfilter and xygrid parameters:

unix% pget get_sky_limits dmfilter
x=3595.5:4771.5:#294,y=3799.5:4987.5:#297

unix% pget get_sky_limits xygrid
3595.5:4771.5:#294,3799.5:4987.5:#297

As discussed in the Caveats section, you may get better results if you use the X,Y range from the data subspace, rather than using dmcoords.


Caveats

As discussed in the dmbinning ahelp page, care must be taken when binning real-values columns (such as the SKY column of an event file), otherwise the edges of the files may not match up exactly. This can be seen in the above example it you use the dmfilter expression to create a copy of the image and subtract it from the original. In this case, five pixels at the bottom edge of the frame are different.

unix% dmcopy \
      "acisf01838N004_evt2.fits[bin x=3595.5:4771.5:#294,y=3799.5:4987.5:#297]" img2.fits
unix% dmimgcalc img.fits img2.fits diff.fits sub
unix% dmstat diff.fits centroid=no
diff.fits
    min:        -1            @:        ( 3757.5 3801.5 )
    max:        0             @:        ( 3597.5 3801.5 )
   mean:        -5.7261962024e-05
  sigma:        0.0075669467483
    sum:        -5
   good:        87318
   null:        0

The following example shows a larger difference:

unix% dmcopy "acisf01838N004_evt2.fits[bin x=3570:4874:16,y=3650:4980:16]" zoom16_1.fits
unix% dmcopy "acisf01838N004_evt2.fits[x=3570:4874,y=3650:4980][bin x=::16,y=::16]" zoom16_2.fits
unix% dmimgcalc zoom16_1.fits zoom16_2.fits zoom_diff.fits sub 
unix% dmstat zoom_diff.fits centroid=no
zoom_diff.fits
    min:	-104 	      @:	( 4058 4250 )
    max:	164 	      @:	( 4058 4234 )
   mean:	0 
  sigma:	2.8512788255 
    sum:	0 
   good:	6888 
   null:	0 

The resulting image (zoom_diff.fits) is shown in Figure 1, and the subspace option of dmlist reports the X,Y ranges of the files to be:

zoom16_1.fits: x = 3570.0:4882.0  y = 3650.0:4994.0 
zoom16_2.fits: x = 3570.0:4874.0  y = 3650.0:4980.0 

Figure 1: Offset due to different binning filters

[Thumbnail image: The subtracted image highlights the pixels which differed in the input files.]

[Version: full-size]

[Print media version: The subtracted image highlights the pixels which differed in the input files.]

Figure 1: Offset due to different binning filters

Care must be taken when binning real-values columns (such as the SKY column of an event file), otherwise the edges of the files may not match up exactly.



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


        infile = img.fits         Input dataset/block specification
      asolfile = none             Input aspect solution file
#
# Position of photon in different coord systems
#
       chip_id = 6                Chip ID number
         chipx = 995.6866592050034 Chip X [pixel]
         chipy = -107.7026386592746 Chip Y [pixel]
         tdetx = 0                TDETX [pixel]
         tdety = 0                TDETY [pixel]
          detx = 3827.861836131362 FPC X [pixel]
          dety = 4613.262553706417 FPC Y [pixel]
             x = 3595.5           Sky X [pixel]
             y = 3799.5           Sky Y [pixel]
      logicalx = 0.5              X coordinate in binned image [pixel]
      logicaly = 0.5              Y coordinate in binned image [pixel]
            ra = 18:33:49.497     RA [deg or hh:mm:ss]
           dec = -10:37:47.46     Dec [deg or dd:mm:ss]
         theta = 4.775819411656162 Off axis angle [arcmin]
           phi = 117.4675919434541 Azimuthal angle [deg]
         order = 0                Grating order
        energy = 1                Energy [keV]
    wavelength = 0                Wavelength [A]
         ra_zo = 18:33:49.497     RA of zero order
        dec_zo = -10:37:47.46     Dec of zero order
     (asolfile = none)		  Input aspect solution file
       (option = )                Conversion option
#
# Override setup for observation
# All parameters here are strings so that they can
# be set blank, in which case the data file value is used
#
       (celfmt = hms)             RA and Dec format [deg or hms] (xx.xx or xx:xx:xx.x)
     (detector = )                Detector (ACIS or HRC-I or HRC-S)
      (grating = )                Grating
        (fpsys = )                FP convention
          (sim = )                SIM position (eg 0.0 0.0 -190.6)
     (displace = )                STF displacement (X,Y,Z,AX,AY,AZ)
       (ra_nom = )                Nominal pointing RA [deg or hh:mm:ss]
      (dec_nom = )                Nominal dec [deg or dd:mm:ss]
     (roll_nom = )                Nominal roll [deg]
       (ra_asp = )                Instantaneous pointing RA [deg]
      (dec_asp = )                Instantaneous pointing Dec [deg]
     (roll_asp = )                Instantaneous Aspect roll [deg]
#
      (geompar = geom)            Parameter file for Pixlib Geometry files
      (verbose = 0)               Debug Level
         (mode = ql)              
    


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


        infile = img.fits         Input dataset/block specification
      asolfile = none             Input aspect solution file
#
# Position of photon in different coord systems
#
       chip_id = 8                Chip ID number
         chipx = 424.3637197084861 Chip X [pixel]
         chipy = 1155.542356460029 Chip Y [pixel]
         tdetx = 0                TDETX [pixel]
         tdety = 0                TDETY [pixel]
          detx = 5339.917834531454 FPC X [pixel]
          dety = 3351.096693883851 FPC Y [pixel]
             x = 4771.5           Sky X [pixel]
             y = 5379.5           Sky Y [pixel]
      logicalx = 294.5            X coordinate in binned image [pixel]
      logicaly = 395.5            Y coordinate in binned image [pixel]
            ra = 18:33:10.266     RA [deg or hh:mm:ss]
           dec = -10:24:50.08     Dec [deg or dd:mm:ss]
         theta = 11.88772992568937 Off axis angle [arcmin]
           phi = 329.0581834129313 Azimuthal angle [deg]
         order = 0                Grating order
        energy = 1                Energy [keV]
    wavelength = 0                Wavelength [A]
         ra_zo = 18:33:10.266     RA of zero order
        dec_zo = -10:24:50.08     Dec of zero order
     (asolfile = )		  Input aspect solution file
       (option = )                Conversion option
#
# Override setup for observation
# All parameters here are strings so that they can
# be set blank, in which case the data file value is used
#
       (celfmt = hms)             RA and Dec format [deg or hms] (xx.xx or xx:xx:xx.x)
     (detector = )                Detector (ACIS or HRC-I or HRC-S)
      (grating = )                Grating
        (fpsys = )                FP convention
          (sim = )                SIM position (eg 0.0 0.0 -190.6)
     (displace = )                STF displacement (X,Y,Z,AX,AY,AZ)
       (ra_nom = )                Nominal pointing RA [deg or hh:mm:ss]
      (dec_nom = )                Nominal dec [deg or dd:mm:ss]
     (roll_nom = )                Nominal roll [deg]
       (ra_asp = )                Instantaneous pointing RA [deg]
      (dec_asp = )                Instantaneous pointing Dec [deg]
     (roll_asp = )                Instantaneous Aspect roll [deg]
#
      (geompar = geom)            Parameter file for Pixlib Geometry files
      (verbose = 0)               Debug Level
         (mode = ql)              
    

History

04 Jan 2005 updated for CIAO 3.2: minor change to dmcoords parameter file
16 Dec 2005 reviewed for CIAO 3.3: no changes
01 Dec 2006 reviewed for CIAO 3.4: no changes
16 Jan 2008 updated for CIAO 4.0: get_sky_limits v1.9 (internal code updates needed for software changes; the case of the axis names, as reported on-screen and saved to the dmfilter parameter, may change); filename and screen output updated for reprocessed data (version N004 event file)
24 Oct 2008 get_sky_limits v1.13 (fixes a rare segmentation fault and adds the pixel size in sky coordinates to the screen output.)
04 Feb 2009 updated for CIAO 4.1: images are inline
20 Apr 2009 updated for CIAO 4.1.2: the asolfile parameter in dmcoords has changed from hidden to automatic (updated parameter file listing)
06 May 2009 check the version of the CIAO scripts package instead of the individual script
05 Feb 2010 reviewed for CIAO 4.2: no changes
13 Jan 2011 reviewed for CIAO 4.3: no changes
04 Apr 2011 updated for 04 Apr scripts package release: get_sky_limits script prints the version at verbose > 0.
20 Jul 2011 required software updates are listed in Synopsis
09 Jan 2012 reviewed for CIAO 4.4: no changes
03 Dec 2012 Review for CIAO 4.5; file version update, fix time DSS (extra new-line)
03 Dec 2013 Review for CIAO 4.6; no changes.
22 Dec 2013 Review for CIAO 4.7; added link to apply_fov_limits help file.
25 Jan 2022 Review for CIAO 4.14. Updated for Repro-5 and CALDB 4.9.6.